Sum(j, 1, n, Sum(i, 1, j-1, p^(j-i)))

[mikldk]

I am trying to use Yacas to calculate the sum

Sum(j, 1, n, Sum(i, 1, j-1, p^(j-i)))

According to https://www.wolframalpha.com/input/?i=x:%3Dsum(sum(p%5E(j-i),i%3D1..j-1),j%3D1..n) this should be

-(p*(-p^n + n*(p-1) + 1))/((p-1)^2)

but Yacas cannot reduce the sum.

Is there anything I can do to help Yacas calculating this?

[grzegorzmazur]

Hi,

Thanks for the report. While I'm not at all surprised that yacas can't produce so nicely simplified result, the fact that it's not simplified at all is no good news at all. At the first glance this stems from the fact that IIRC there are no rules which could reduce the inner sum, Sum(i, 1, j-1, p^(j-i)). Whether the result of the inner sum could be later on simplified is another question, but definitely dealing with the inner one would be a start.

The support for sums in yacas is on one hand only rudimentary, but on the other a bit convoluted. The code is in share/yacas/sums.rep/code.ys and the rules are implicitly defined by creating instances of SumFunc(), eg

SumFunc(_k,0,_n,(r_IsFreeOf(k))^(_k), (1-r^(n+1))/(1-r) );

You're are more than welcome to try and come up with some more sum reduction rules, pull requests are always welcome :)

If you don't feel like this, could you please suggest some sums which are not implemented but would be good to have? And I'll do my best to implement the simplification rules for them.

Thanks a lot, Grzesiek

[mikldk]

I have tried making some additional rules together with simple tests (see below).

I apologise that it is not a pull request, but I don't think it is ready for that yet. It needs additional proofing and the rules probably needs more guards/predicates. And all must be checked again.

Rules:

// Geometric series starting from 1:
// [r1]
SumFunc(_k,1,_n,(r_IsFreeOf(k))^(_k), r*(1 - r^n)/(1-r) );

////////////////////////////////////////////////////////////////////////////
// Makes e.g. Sum(k, 0, n, k^2) work (from 0 instead as from 1), 
// including taking care of Sum(k, 0, n, k^0) = 0^0 + 1^0 + ... = Undefined
// [r2]
SumFunc(_k,0,_n,_k, n*(n+1)/2 );
SumFunc(_k,0,_n,_k^0, Undefined );

//Same body as when starting from 1, but using "Sum(k,1,n,k^p)" below doesn't work
// [r3]
SumFunc(_k,0,_n,_k^_p,(Bernoulli(p+1,n+1) - Bernoulli(p+1))/(p+1), IsPositiveInteger(p) );


// [r4]
SumFunc(_k,1,_n,1/_k, HarmonicNumber(n) );
SumFunc(_k,1,_n,r_IsFreeOf(k)/_k, r*HarmonicNumber(n) );


// [r5]
// FIXME: IsPositive(p), but gives error in yacas
SumFunc(_k,1,_n,1/_k^_p, HarmonicNumber(n, p) );
SumFunc(_k,1,_n,_k^(-_p), HarmonicNumber(n, p) );


// [r6]
// FIXME: IsPositive(r), but gives error in yacas
// https://www.wolframalpha.com/input/?i=Sum(k,+0,+n,+p%5E(-k))+%3D+Sum(k,+0,+n,+(1%2Fp)%5Ek)
// Using Geometric series with r replaces by 1/r
SumFunc(_k,0,_n,(r_IsFreeOf(k))^(-_k), (r - r^(-n)) / (r - 1) );
SumFunc(_k,1,_n,(r_IsFreeOf(k))^(-_k), (1 - r^(-n))/(r - 1) );


// [r7]
// FIXME: Take care of 0/0?
// FIXME: IsPositive(r), but gives error in yacas
SumFunc(_k,0,_n,r_IsFreeOf(k)^(_k-m_IsFreeOf(k)), Simplify( Sum(k, 0, n, r^k)/r^m ) );
SumFunc(_k,1,_n,r_IsFreeOf(k)^(_k-m_IsFreeOf(k)), Simplify( Sum(k, 1, n, r^k)/r^m ) );


// [r8]
// FIXME: Take care of 0/0?
// FIXME: IsPositive(r), but gives error in yacas
SumFunc(_k,0,_n,r_IsFreeOf(k)^(m_IsFreeOf(k)-_k), Simplify( Sum(k, 0, n, r^(-k) )*r^m ) );
SumFunc(_k,1,_n,r_IsFreeOf(k)^(m_IsFreeOf(k)-_k), Simplify( Sum(k, 1, n, r^(-k) )*r^m ) );


// Additional

// FIXME: Take care of 0/0?
// FIXME: r: IsPositive()
SumFunc(_k,0,_n,r_IsFreeOf(k)^(_k+m_IsFreeOf(k)), Simplify( Sum(k, 0, n, r^k)*r^m ) );
SumFunc(_k,1,_n,r_IsFreeOf(k)^(_k+m_IsFreeOf(k)), Simplify( Sum(k, 1, n, r^k)*r^m ) );

// FIXME: Take care of 0/0?
// FIXME: r: IsPositive()
SumFunc(_k,0,_n,r_IsFreeOf(k)^(m_IsFreeOf(k)*_k), Simplify( Sum(k, 0, n, (r^m)^k) ) );
SumFunc(_k,1,_n,r_IsFreeOf(k)^(m_IsFreeOf(k)*_k), Simplify( Sum(k, 1, n, (r^m)^k) ) );

Tests:

// [r1]
Print("[r1]");
// https://www.wolframalpha.com/input/?i=Sum(k,+1,+n,+p%5Ek)
Print(TestYacas(  Sum(k, 1, n, p^k), p*(p^n - 1)/(p - 1)  ));

// [r2]
Print("[r2]");
// https://www.wolframalpha.com/input/?i=Sum(k,+0,+n,+k)
Print(TestYacas(  Sum(k, 0, n, k^0), Undefined )); // 0^0 is Undefined
Print(TestYacas(  Sum(k, 0, n, k), n*(n+1)/2  ));

// [r3]
Print("[r3]");
// https://www.wolframalpha.com/input/?i=Sum(k,+0,+n,+k%5E2)
Print(TestYacas(  Sum(k, 0, n, k^2), n*(n+1)*(2*n+1)/6  ));

// [r4]
Print("[r4]");
Print(TestYacas(  Sum(k, 1, n, 1/k), HarmonicNumber(n) ));
Print(TestYacas(  Sum(k, 1, n, 3/k), 3*HarmonicNumber(n) ));

// [r5]
Print("[r5]");
Print(TestYacas(  Sum(k, 1, n, 1/k^3), HarmonicNumber(n,3) ));
Print(TestYacas(  Sum(k, 1, n, k^(-3)), HarmonicNumber(n,3) ));

// [r6]
Print("[r6]");
// https://www.wolframalpha.com/input/?i=Sum(k,+0,+n,+p%5E(-k))
// https://www.wolframalpha.com/input/?i=Sum(k,+0,+n,+p%5E(-k))+%3D+(p*(1-(1%2Fp)%5E(n%2B1)))%2F(p-1)
Print(TestYacas(  Sum(k, 0, n, p^(-k)), (p-p^(-n))/(p-1) ));
// https://www.wolframalpha.com/input/?i=Sum(k,+1,+n,+p%5E(-k))
// https://www.wolframalpha.com/input/?i=Sum(k,+1,+n,+p%5E(-k))%3D(1-(1%2Fp)%5En)%2F(p-1)
Print(TestYacas(  Sum(k, 1, n, p^(-k)), (1-p^(-n))/(p-1) ));

// [r7]
Print("[r7]");
// https://www.wolframalpha.com/input/?i=Sum(k,+0,+n,+p%5E(k-3))
Print(TestYacas(  Sum(k, 0, n, p^(k-3)), (p^(n+1) - 1)/((p-1)*p^3) ));
// https://www.wolframalpha.com/input/?i=Sum(k,+1,+n,+p%5E(k-3))
Print(TestYacas(  Sum(k, 1, n, p^(k-3)), (p^n - 1)/((p-1)*p^2) ));
// https://www.wolframalpha.com/input/?i=Sum(k,+0,+n,+p%5E(k-m))
Print(TestYacas(  Sum(k, 0, n, p^(k-m)), (p^(n+1) - 1)*p^(-m)/(p-1) ));
// https://www.wolframalpha.com/input/?i=Sum(k,+1,+n,+p%5E(k-m))
// https://www.wolframalpha.com/input/?i=Sum(k,+1,+n,+p%5E(k-m))%3D(p*(1-p%5En))%2F((1-p)*p%5Em)
Print(TestYacas(  Sum(k, 1, n, p^(k-m)), (p*(1-p^n))/((1-p)*p^m) ));

// [r8]
Print("[r8]");
// https://www.wolframalpha.com/input/?i=Sum(k,+0,+n,+p%5E(3-k))
Print(TestYacas(  Sum(k, 0, n, p^(3-k)), (p^3*(p-p^(-n)))/(p-1) ));
// https://www.wolframalpha.com/input/?i=Sum(k,+1,+n,+p%5E(3-k))
Print(TestYacas(  Sum(k, 1, n, p^(3-k)), (p^3*(1-p^(-n)))/(p-1) ));
// https://www.wolframalpha.com/input/?i=Sum(k,0,+n,+p%5E(m-k))
Print(TestYacas(  Sum(k, 0, n, p^(m-k)), (p^m*(p-p^(-n)))/(p-1) ));
// https://www.wolframalpha.com/input/?i=Sum(k,1,+n,+p%5E(m-k))
Print(TestYacas(  Sum(k, 1, n, p^(m-k)), (p^m*(1-p^(-n)))/(p-1) ));

[mikldk]

The above rules still does not give the double sum, but the inner can now be done.