3
$\begingroup$

Can Mathematica solve (for b[k] and c[r]) the following system of recurrence relations? b[k] is defined for only odd natural numbers, and c[r] is defined for only even non-negative integers.

For n greater than or equal to 3; 

b[2n+1] = c[2n] + c[2n-2] + 2*a[n],
c[2n] = b[2n-1] + c[2n-2],
a[n] = .0794367(.460811)^n + .255972(-.675131)^n + .664591(3.21432)^n;

c[0]=1, b[1]=4, c[2]=5, b[3]=10, c[4]=15, b[5]=34.

I seemed to have tried "everything" in Mathematica 11, but nothing seems to work.... Any help is greatly appreciated!

UPDATE (September 4, 2016)

After consolidating everything, here is the question.

Find a "closed-form" solution (formula) to the following piecewise defined recurrence relation:

For integers k greater than or equal to 3:

b[k] = b[k-1] + b[k-3] + 2*[.0794367(.460811)^((k-1)/2) + \
.255972(-.675131)^((k-1)/2) + .664591(3.21432)^((k-1)/2)]`

If k is odd:

b[k] = b[k-1] + b[k-2], if k is even; b[0] = 1, b[1] = 4, b[2] = 5

This piecewise defined recurrence relation is the one that I am interested in. I verified (using this piecewise defined recurrence relation, by hand calculations) the correct values: 

b[0]=1, b[1]=4, b[2]=5, b[3]=10, b[4]=15, b[5]=34, b[6]=49, 
b[7]=108, b[8]=157, b[9]=348, b[10]=505, ...

Sincerely, Richard M. Low

Improve this question
$\endgroup$
3
  • 2
    $\begingroup$ Have you tried RSolve? What syntax did you use? Many problems are caused by subtle errors in syntax that someone could spot easily. $\endgroup$
    – mikado
    Sep 3, 2016 at 7:24
  • 2
    $\begingroup$ please show your efforts, tell us what you've tried. :) $\endgroup$
    – Wjx
    Sep 3, 2016 at 8:50
  • 2
    $\begingroup$ Please cross reference your posts here with the corresponding ones on Wolfram Community: community.wolfram.com/groups/-/m/t/917888 Just edit your post and include a link (on both websites). This is to prevent duplication of effort (i.e. two people posting the same answer in two places and not being aware of it). $\endgroup$
    – Szabolcs
    Sep 6, 2016 at 10:00

2 Answers 2

Reset to default
7
$\begingroup$

Eliminating b from the first two expressions in the question and redefining n yields

c[n] == 2 c[n - 1] + c[n - 2] + 2 a[n - 1]

Although RSolve can handle this expression with a as given in the question, it seems better to Rationalize it.

794367 10^-7 (460811 10^-6)^n+255972 10^-6 (-675131 10^-6)^n+664591 10^-6 (321432 10^-5)^n

Then, 

Simplify[RSolve[{c[n] == 2 c[n - 1] + c[n - 2] + 2 (794367 10^-7 (460811 10^-6)^(n - 1) + 
    255972 10^-6 (-675131 10^-6)^(n - 1) + 664591 10^-6 (321432 10^-5)^(n - 1)), 
    c[0] == 1, c[1] == 5}, c[n], n], n ∈ Integers] // First

(* c[n] -> *)

The first several values are

Table[Simplify[c[n] /. %], {n, 0, 8}]

(* {1, 5, 74999988159517/5000000000000, 244999948333156036527/5000000000000000000, 
    784999809145104350538564757/5000000000000000000000000, 
    2524999383279806241504389904195767/5000000000000000000000000000000, (4^(-4 - 3 n) 5^(-7 - 6 n) (10^(6 n) (1172409346315586312019869378750000 - 
   1088292064380910194582788537056291 Sqrt[2]) (1 - Sqrt[2])^n + 10^(6 n) (1 + Sqrt[2])^
   n (1172409346315586312019869378750000 + 1088292064380910194582788537056291 Sqrt[2]) + 
   1500000 (130809072583670204544034435447 240^n 13393^(1 + n) - 
   110652256667034006265872216 460811^(1 + n) - 756091031029961199112699840 (-1)^n 675131^
   (1 + n))))/89285790105620659634534766727097 
    8114998166092225878868708271269265029997/5000000000000000000000000000000000000, 
    26084994986029714868240982411237137250977319807/        
    5000000000000000000000000000000000000000000, 
    83844987673045565789165365810827032895161630505820037/
    5000000000000000000000000000000000000000000000000} *)

or, numericized,

N[%, 15]

(* {1.00000000000000, 5.00000000000000, 14.9999976319034, 48.9999896666312, 
    156.999961829021, 504.999876655961, 1622.99963321845, 5216.99899720594, 
    16768.9975346091} *)

That these numbers are very close to being integers suggests that a more precise expression for a would yield integers directly.

Edit: Computation of b

Obtaining b from c is straightforward. For convenience, define

fc[n_] := (4^(-4 - 3 n) 5^(-7 - 6 n) (10^(6 n) (1172409346315586312019869378750000 - 
   1088292064380910194582788537056291 Sqrt[2]) (1 - Sqrt[2])^n + 10^(6 n) (1 + Sqrt[2])^
   n (1172409346315586312019869378750000 + 1088292064380910194582788537056291 Sqrt[2]) + 
   1500000 (130809072583670204544034435447 240^n 13393^(1 + n) - 
   110652256667034006265872216 460811^(1 + n) - 756091031029961199112699840 (-1)^n 675131^
   (1 + n))))/89285790105620659634534766727097

Then

fb[n_] := Simplify[fc[n] + fc[n - 1] + 2 (794367 10^-7 (460811 10^-6)^n + 
    255972 10^-6 (-675131 10^-6)^n + 664591 10^-6 (321432 10^-5)^n)]

The two can be interleaved, as requested, by

fh[n_] := If[OddQ[n], fb[(n - 1)/2], fc[n/2]]
Table[fh[n], {n, 0, 11}] // N

(* {1., 4., 5., 10., 15., 34., 49., 108., 157., 348., 505., 1118.} *)

ListLogPlot[Table[fh[n], {n, 0, 30}] // N, PlotRange -> All]

enter image description here

Addendum

Another approach is to leave a undefined at first.

Simplify[RSolve[{c[n] == 2 c[n - 1] + c[n - 2] + 2 a[n], c[0] == 1, c[1] == 5},
   c[n], n], n ∈ Integers] // First;
Collect[% /. {K[1] -> k, K[2] -> k}, {a[0], a[1], _Sum}, Simplify[#, n ∈ Integers] &];
ReleaseHold[Hold[Evaluate[%]] /. {k, i1_, i2_} -> {i, i1 + 2, i2 + 2} /. 
    k -> k - 2 /. i -> k];
Simplify[% /. z1_ Sum[z2_, z3_] :> z1 Sum[z2, z3] + (z1 z2 /. k -> 0) + (z1 z2 /. k -> 1)]
    /. {k, 0, n + 1} -> {k, 2, n + 1}

(* {c[n] -> ((1 - Sqrt[2])^n - 2*Sqrt[2]*(1 - Sqrt[2])^n + (1 + Sqrt[2])^n + 
   2*Sqrt[2]*(1 + Sqrt[2])^n + 2*(1 + Sqrt[2])^n*
   Sum[((-1)^(2 - k)*(1 - Sqrt[2])^(-2 + k)*(-1 + Sqrt[2])*a[k])/Sqrt[2], {k, 2, 1 + n}]
   + 2*(1 - Sqrt[2])^n*Sum[((-1)^(2 - k)*(1 + Sqrt[2])^(-1 + k)*a[k])/
   Sqrt[2], {k, 2, 1 + n}])/2} *)

Results can be tested at each step of the simplification, in all cases yielding

Table[Simplify[c[n] /. % /. n -> i], {i, 0, 5}]

(* {1, 5, 11 + 2 a[2], 27 + 4 a[2] + 2 a[3], 65 + 10 a[2] + 4 a[3] + 2 a[4], 
   157 + 24 a[2] + 10 a[3] + 4 a[4] + 2 a[5]} *)

It is natural to ask whether the coefficients of a represent a known series, and they do.

2 FindSequenceFunction[Cases[Last@%, z_ a[_] -> z/2] // Reverse][n]

(* 2 Fibonacci[n, 2] *)

Thus, c[n] can be represented as

f[n_] := Simplify[1/2 ((1 - Sqrt[2])^n - 2 Sqrt[2] (1 - Sqrt[2])^n + (1 + Sqrt[2])^n + 
    2 Sqrt[2] (1 + Sqrt[2])^n)] + Plus @@ (2 Fibonacci[n + 1 - #, 2] a[#] & /@ Range[2, n])

One might hope that a[n], rounded to the nearest integer, also represents a known series, but, FindSequenceFunction returns unevaluated.

Table[Round@a[n], {n, 0, 6}] // FindSequenceFunction

(* FindSequenceFunction[{1, 2, 7, 22, 71, 228, 733}] *)

Second Addendum: a, b, and c in closed form

Vaclav Kotesovec made the very useful observation in a comment, that a[n], rounded to integers, has the recurrence relation,

a[n] == 3 a[n - 1] + a[n - 2] - a[n - 3]

It turns out that Mathematica can produce this result without difficulty, once it is realized from Kotesovec's comment that it is feasible to do so. With a[n] as given above,

Table[FullSimplify[a[n]], {n, 0, 11}] // N // Round // FindLinearRecurrence
(* {3, 1, -1} *)

This can be solved to obtain a in closed form.

Clear[a]
a[n] /. First@RSolve[{a[n] == 3 a[n - 1] + a[n - 2] - a[n - 3], 
    a[0] == 1, a[1] == 2, a[2] == 7}, a[n], n]
(* Root[1 - #1 - 3 #1^2 + #1^3 &, 2]^n Root[-1 + 18 #1 - 74 #1^2 + 74 #1^3 &, 1] + 
   Root[1 - #1 - 3 #1^2 + #1^3 &, 1]^n Root[-1 + 18 #1 - 74 #1^2 + 74 #1^3 &, 2] + 
   Root[1 - #1 - 3 #1^2 + #1^3 &, 3]^n Root[-1 + 18 #1 - 74 #1^2 + 74 #1^3 &, 3] *)

Sample results are as expected,

Table[FullSimplify[%], {n, 0, 11}]
(* {1, 2, 7, 22, 71, 228, 733, 2356, 7573, 24342, 78243, 251498} *)

This exact result can be numericized to whatever precision is desired using N. For $MachinePrecision, 

%% // N
(* 0.255972 (-0.675131)^n + 0.0794367 0.460811^n + 0.664591 3.21432^n *)

which is the approximate expression for a shown in the question. Now, all this would be of little interest, except that the same can be done for b and c!

Table[Simplify[fc[n]], {n, 0, 11}] // N // Round // FindLinearRecurrence
(* {3, 1, -1} *)
c[n] /. First@RSolve[{c[n] == 3 c[n - 1] + c[n - 2] - c[n - 3], 
    c[0] == 1, c[1] == 5, c[2] == 15}, c[n], n]
(* Root[1 - #1 - 3 #1^2 + #1^3 &, 1]^n Root[-1 - 25 #1 - 37 #1^2 + 37 #1^3 &, 1] + 
   Root[1 - #1 - 3 #1^2 + #1^3 &, 2]^n Root[-1 - 25 #1 - 37 #1^2 + 37 #1^3 &, 2] + 
   Root[1 - #1 - 3 #1^2 + #1^3 &, 3]^n Root[-1 - 25 #1 - 37 #1^2 + 37 #1^3 &, 3] *)
Table[FullSimplify[%], {n, 0, 11}]
(* {1, 5, 15, 49, 157, 505, 1623, 5217, 16769, 53901, 173255, 556897} *)
%% // N
(* -0.428786 (-0.675131)^n - 0.0428314 0.460811^n + 1.47162 3.21432^n *)

Table[Simplify[fb[n]], {n, 0, 11}] // N // Round // FindLinearRecurrence
(* {3, 1, -1} *)
b[n] /. First@RSolve[{b[n] == 3 b[n - 1] + b[n - 2] - b[n - 3], 
    b[0] == 4, b[1] == 10, b[2] == 34}, b[n], n]
(* Root[1 - #1 - 3 #1^2 + #1^3 &, 2]^n Root[-2 + 90 #1 - 148 #1^2 + 37 #1^3 &, 1] + 
   Root[1 - #1 - 3 #1^2 + #1^3 &, 1]^n Root[-2 + 90 #1 - 148 #1^2 + 37 #1^3 &, 2] + 
   Root[1 - #1 - 3 #1^2 + #1^3 &, 3]^n Root[-2 + 90 #1 - 148 #1^2 + 37 #1^3 &, 3] *)
Table[FullSimplify[%], {n, 0, 11}]
(* {4, 10, 34, 108, 348, 1118, 3594, 11552, 37132, 119354, 383642, 1233148} *)
N // %%
(* 0.718273 (-0.675131)^n + 0.0230942 0.460811^n + 3.25863 3.21432^n *)

Thus, we have been able to provide in a straightforward manner both exact and MachinePrecision expressions for a, b, and c.

Improve this answer
$\endgroup$
1
  • $\begingroup$ bbgodfrey, Your answer is very helpful! Thank you very much. $\endgroup$
    – richmlow
    Sep 5, 2016 at 22:32
3
$\begingroup$

If numerical values are desired:

a[n_] := .0794367 (.460811)^n + .255972 (-.675131)^
     n + .664591 (3.21432)^n;
b[j_] := c[j - 1] + c[j - 3] + 2 a[(j - 1)/2]
c[j_] := b[j - 1] + c[j - 2]
c[0] := 1
b[1] := 4
c[2] := 5
b[3] := 10
c[4] := 15
b[5] := 34.

Presentation:

fun[j_] := If[EvenQ[j], c[j], b[j]];
bcf[n_?OddQ] := 
 With[{r = {0 -> Red, 1 -> Blue}}, 
  Framed[TableForm[{Style[fun@#, Mod[#, 2] /. r], a@#} & /@ 
     Range[0, n], 
    TableHeadings -> {MapIndexed[
       Subscript[#1, 
         Style[ToString[#2[[1]] - 1], Mod[#2[[1]] - 1, 2] /. r]] &, 
       Flatten[Table[{Style["c", Red], Style["b", Blue]}, (n + 1)/
          2]]], {"value", "\!\(\*SubscriptBox[\(a\), \(n\)]\)"}}]]]

enter image description here

This does not scale well. I have not tried FindSequenceFunction etc. Others may wish to.

Apologies for any errors.

Improve this answer
$\endgroup$
4
  • 2
    $\begingroup$ Sequence a(n) is oeis.org/A030186 $\endgroup$
    – Vaclav Kotesovec
    Sep 3, 2016 at 13:29
  • $\begingroup$ Everybody has been quite helpful thus far! Thank you. My original question can be formulated in the following manner: Can Mathematica solve (ie., find a "closed-form" formula for b[k]) for the following "piecewise" defined recurrence relation? For integers k greater than or equal to 3, b[k] = b[k-1] + b[k-3] + 2*[.0794367(.460811)^((k-1)/2) + .255972(-.675131)^((k-1)/2) + .664591(3.21432)^((k-1)/2)], if k is odd; b[k] = b[k-1] + b[k-2], if k is even; b[0] = 1, b[1] = 4, b[2] = 5. Sincerely, Richard M. Low $\endgroup$
    – richmlow
    Sep 5, 2016 at 8:25
  • $\begingroup$ ubpdqn, Thank you for your reply. It is helping me to get to the final goal. $\endgroup$
    – richmlow
    Sep 5, 2016 at 22:38
  • $\begingroup$ @VaclavKotesovec Thank you for identifying the sequence a[n]. Knowing its recurrence relation, I then was able to determine the recurrence relations with Mathematica for b[n] and c[n]. Very helpful. $\endgroup$
    – bbgodfrey
    Sep 7, 2016 at 13:42

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

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