Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

I have a sequence $(u_{n})$ $$u_1= 1, \quad u_2 = 2, \quad u_3 = 3, \quad u_{n}= -u_{n-3} + 3u _{n-2} +2 u_{n-1}, \quad \forall n \geqslant 4.$$ I want to list the first $20$ terms of this sequence and find the sum of all terms and sum of all even numbers of $20$ terms. I tried

Clear[u];
u[1] := 1; u[2] := 2;
u[3] := 3;
u[y_] := -u[y - 3] + 3u[y - 2] + 2u[y - 1]
A := Table[ u[k], {k, 1, 20}]

Sum[ A[[i]], {i, 1, Length[A]}]

I can't find sum of even numbers. How do I tell Mathematica to do that ?

share|improve this question
    
Do you want to sum all u[n] terms, i.e. up to infinity ? –  Artes Dec 23 '12 at 3:58
    
You might consider cacheing values of u to speed things up. Just change the first part of your code to: u[1] = 1; u[2] = 2; u[3] = 3; u[y_] := u[y] = -u[y - 3] + 3*u[y - 2] + 2*u[y - 1]; –  m_goldberg Dec 23 '12 at 4:23
    
Also, I note you can get the sum of all elements with Total@A. –  m_goldberg Dec 23 '12 at 4:27

4 Answers 4

up vote 9 down vote accepted

You can use Select to pull the elements meeting some specified criterion from a list. The selected elements are returned in a list that you can feed to Total to sum them. Total is better than Sum for your problem because you don't have to supply the length of the list and you have no need for indexing.

Clear[u];
u[1] := 1; u[2] := 2;
u[3] := 3;
u[y_] := -u[y - 3] + 3*u[y - 2] + 2*u[y - 1]
A = Table[u[k], {k, 1, 20}]
Total[Select[A, EvenQ]]

(* 278878570 *)
share|improve this answer

Working with RSolve we can find much more than only a few first terms, here is a general term of your function u[n] e.g. :

u[n_] = u[n] /. Flatten[ RSolve[{ u[1] == 1, u[2] == 2, u[3] == 3, 
                                  u[n] == -u[n-3] + 3 u[n-2] + 2 u[n-1]}, u[n], n]]
Root[1 - 3 #1 - 2 #1^2 + #1^3 &, 3]^n Root[-45 + 457 #1 - 1028 #1^2 + 257 #1^3 &, 1] + 
Root[1 - 3 #1 - 2 #1^2 + #1^3 &, 1]^n Root[-45 + 457 #1 - 1028 #1^2 + 257 #1^3 &, 2] + 
Root[1 - 3 #1 - 2 #1^2 + #1^3 &, 2]^n Root[-45 + 457 #1 - 1028 #1^2 + 257 #1^3 &, 3]

Related : How do I work with Root objects.

Let's write the first 20 terms :

t = FullSimplify @ Table[ u[k], {k, 20}]
{1, 2, 3, 11, 29, 88, 252, 739, 2146, 6257, 18213, 53051, 154484, 449908, 1310217,  
 3815674, 11112091, 32360987, 94242573, 274456016}
Total @ t
417982742    

Let's sum only even terms :

Total @ Cases[ t, _?EvenQ]
278878570

One might understand as well the task this way :

Total @ FullSimplify @ Table[ u[2 k], {k, 10}]
311142733

Moreover we can find a general symbolic formula for partial sums s[n] of the first n terms of u[i] for 1<= i <= n :

s[n_] := FullSimplify[ Sum[ u[i], {i, n}], n ∈ Integers && n > 0]

s[n]
1/3 ( 3 +
  Root[1 - 3 #1 - 2 #1^2 + #1^3 &, 2]^n Root[405 - 1233 #1 + 771 #1^2 + 257 #1^3 &, 1] +
  Root[1 - 3 #1 - 2 #1^2 + #1^3 &, 1]^n Root[405 - 1233 #1 + 771 #1^2 + 257 #1^3 &, 2] +
  Root[1 - 3 #1 - 2 #1^2 + #1^3 &, 3]^n Root[405 - 1233 #1 + 771 #1^2 + 257 #1^3 &, 3])

To write s[n] in terms of more familiar expressions we have to make a trick applying ToRadicals and ComplexExpand to every Root object separately and then retrieve the original expression :

1/3 (3 + #1^n #2 + #3^n #4 + #5^n #6)& @@ ComplexExpand @ ToRadicals @ {
Root[1 - 3 #1 - 2 #1^2 + #1^3 &, 2], Root[405 - 1233 #1 + 771 #1^2 + 257 #1^3 &, 1], 
Root[1 - 3 #1 - 2 #1^2 + #1^3 &, 1], Root[405 - 1233 #1 + 771 #1^2 + 257 #1^3 &, 2], 
Root[1 - 3 #1 - 2 #1^2 + #1^3 &, 3], Root[405 - 1233 #1 + 771 #1^2 + 257 #1^3 &, 3]} // 
Simplify // TraditionalForm

enter image description here

It is not especially surprising that neither the series u[n] nor s[n] is convergent, e.g.

Limit[ s[n], n -> Infinity]

$\infty$

you could still find out certain interesting properties e.g. using the Regularization option of Sum. Regularization techniques are ubiquitous in modern theoretical physics, e.g. in quantum field theory one encounters this formula $\sum_{n=1}^{\infty} n = -\frac{1}{12}\;$:

Sum[ n, {n, Infinity}, Regularization -> "Dirichlet"]
-(1/12)

See e.g. this post The Euler-Maclaurin formula, Bernoulli numbers, the zeta function, and real-variable analytic continuation by Terrence Tao. Now we can figure out ( only "Borel" works here) the regularized sum of u[n]:

FullSimplify @ Sum[ u[n], {n, Infinity}, Regularization -> "Borel"]
1

For comparison one can regularize a geometric series of a general term a^n :

Sum[ a^n, {n, Infinity}, Regularization -> "Borel"]
-(a/(-1 + a))

We get 1 for a == 1/2.

share|improve this answer
    
I receive your help very much. Thank you. –  minthao_2011 Dec 23 '12 at 6:26
    
@minthao_2011 You are welcome. I added more advanced and interesting topics. –  Artes Dec 23 '12 at 6:58

If you don't explicitly need the array A you can use Sum in place of Table.

Also, you should be using memoization in your definition of u.

u[x : 1|2|3] := x;

u[y_] := u[y] = -u[y - 3] + 3 u[y - 2] + 2 u[y - 1]

The function evens returns n for any even integer n and zero for any odd integer.

evens = Mod[# - 1, 2] # &;

Sum[evens @ u @ k, {k, 20}]
278878570

Thanks to memoization much larger sums are practical:

Sum[evens @ u @ k, {k, 100000}] // IntegerLength
46421

That's an integer with over 46,000 digits.

share|improve this answer

I think RecurrenceTable deserves justice as a very efficient way to tabulate solutions of the recurrence equations.

Compare procedural approach with no dynamic programming used:

In[214]:= Block[{u, A},
  u[1] := 1; u[2] := 2;
  u[3] := 3;
  u[y_] := -u[y - 3] + 3*u[y - 2] + 2*u[y - 1];
  A = Table[u[k], {k, 1, 25}];
  Total[Select[A, EvenQ]]] // AbsoluteTiming

Out[214]= {4.617608, 7856903616}

to the use of RecurrenceTable, where all that is needed is to define the equation:

In[215]:= 
Select[RecurrenceTable[{u[1] == 1, u[2] == 2, u[3] == 3, 
     u[n] == -u[n - 3] + 3 u[n - 2] + 2 u[n - 1]}, u[n], {n, 1, 25}], 
   EvenQ] // Total // AbsoluteTiming

Out[215]= {0., 7856903616}

Also, it is worth pointing out that any linear recurrence equation of rank $r$ evaluated $\bmod p$ must be periodic, due to the recurrence relation and the finiteness of the state space. Indeed, our system turns out to have period of 7:

In[217]:= Partition[
 Mod[RecurrenceTable[{u[1] == 1, u[2] == 2, u[3] == 3, 
    u[n] == -u[n - 3] + 3 u[n - 2] + 2 u[n - 1]}, u[n], {n, 1, 32}], 
  2], 7]

Out[217]= {{1, 0, 1, 1, 1, 0, 0}, 
           {1, 0, 1, 1, 1, 0, 0}, 
           {1, 0, 1, 1, 1, 0, 0}, 
           {1, 0, 1, 1, 1, 0, 0}}

Thus $u_{7m+2}$, $u_{7m+5}$ and $u_{7m+6}$ are even.

We can now use the symbolic power of Mathematica to work out the recurrence relation for sums $v^{(2)}_m = \sum_{k=0}^{m} u_{7k+2}$, $v^{(5)}_m = \sum_{k=0}^{m} u_{7k+5}$ and $v^{(6)}_m = \sum_{k=0}^{m} u_{7k+6}$:

For example:

In[223]:= Sum[
 DifferenceRoot[
   Function[{u, n}, {u[n + 3] == -u[n] + 3 u[n + 1] + 2 u[n + 2], 
     u[1] == 1, u[2] == 2, u[3] == 3}]][7 k + 2], {k, 0, m}]

Out[223]= 
DifferenceRoot[
  Function[{y, n}, {-y[n] + 6318 y[1 + n] - 
      4544 y[2 + n] - 1774 y[3 + n] + y[4 + n] == 0, 
      y[0] == 0, y[1] == 2, y[2] == 2148, y[3] == 3817822}]][1 + m]

Now, evaluation for any $m$ is just a matter of substitution:

In[224]:= % /. m -> 20

Out[224]= \
118640843309870276692961346051090397680399199522287919469464837082

Working out the solution for $v^{(5)}_m$ and $v^{(6)}_m$ is an easy exercise.

share|improve this answer

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.