Numerical integration of iterative form

Question

I have the following code:

       u[0, t_, h_] := 0; For[k = 1, k < 11,u[k, t_, h_] = 
       u[k - 1, t, h] +h Integrate[(D[u[k - 1, \[Tau], h], \[Tau]] + 
       u[k - 1, \[Tau], h]^2 - 1), {\[Tau], 0, t}]; 
       Print["k=", k, "____", 0 Expand[u[k, t, h]]]; k++]
           

I want to calculate the integral numerically as it is time consuming by increasing the loop.

Answer 1

Use Memoization, it become much faster.

u[0, t_, h_] := 0;
u[k_, t_, h_] := u[k, t, h] = u[k - 1, t, h] + 
      h* Integrate[(D[u[k - 1, τ, h], τ] + u[k - 1, τ, h]^2 - 1), {τ, 0, t}]

data = Table[u[k, t, h], {k, 1, 5}]

To speed it up more, you can use numerical values for t,h like this

u[0,t_,h_]:=0;
u[k_,t_,h_]:=u[k,t,h]=u[k-1,t,h]+h* Integrate[(D[u[k-1,τ,h],τ]+u[k-1,τ,h]^2-1),{τ,0,t}];
t=1;h=3;
data=Table[u[k,t,h],{k,1,11}]//N

Which took less than one minute. For k=9 it took little over one second.

u[0, t_, h_] := 0;
u[k_, t_, h_] := 
  u[k, t, h] = 
   u[k - 1, t, h] + 
    h* Integrate[(D[u[k - 1, \[Tau], h], \[Tau]] + 
        u[k - 1, \[Tau], h]^2 - 1), {\[Tau], 0, t}];
t = 1; h = 3;
Timing[ Table[N@u[k, t, h], {k, 1, 9}]]

Answer 2

If you split up your integral, and make some observations about Integrate[u[-1 + k, t, h]^2, t] always being zero at t=0 you can simplify your expression to:

u[0, t_, h_] := 0;
u[k_, t_, h_] := u[k, t, h] =
  (1 + h) u[-1 + k, t, h] + 
   h ((Integrate[u[-1 + k, τ, h]^2, τ] /. τ -> t) - t)

Which takes a little over a minute to calculate to kMax = 9: (k = 10 seem to be a lot harder, the LeafCount of the DownValues stored in u is starting to get huge):

kMax = 9;
(*preload values*)
u[kMax, t, h]; // AbsoluteTiming

uDat = u[#, t, h] & /@ Range[0, kMax];

{72.2118, Null}

Compare this to the original:

uOrig[0, t_, h_] := 0; 
For[k = 1, k <= kMax, 
  uOrig[k, t_, h_] = 
   uOrig[k - 1, t, h] + 
    h Integrate[(D[uOrig[k - 1, \[Tau], h], \[Tau]] + 
        uOrig[k - 1, \[Tau], h]^2 - 1), {\[Tau], 0, t}];
  Print["k=", k, "____", 0 Expand[uOrig[k, t, h]]]; 
  k++] // AbsoluteTiming

{145.233, Null}

And of course verify our simplification is correct (I only verify up to k=8, because the k=9 case was taking too long to simplify):

origDat = uOrig[#, t, h] & /@ Range[0, kMax];
FullSimplify[Most@origDat == Most@uDat]

True