Recursive computation of symbolic series

Question

I don't have experience in Mathematica and I would like to use the ability it has to do symbolic computations to find the terms of a series of functions recursively. The terms I want to find are given by:

$$S_{n+1}(x,t)=i\int_0^{t}\left(\frac{\partial^2S_n(x,\tau)}{\partial{x^2}}-A_n(x,\tau)\right)d\tau$$ $$A_n(x,t)=\sum_{i=0}^nS(x,t)_{i}S(x,t)_{n-i}$$

I have written the following:

S[0, x_, t_] = Exp[-(x)^2]
A[0, x_, t_] = S[0, x, t]*S[0, x, t]


S[n_Integer, r_, t_] := 
  S[n, r, t] = 
   Refine[Block[{e}, 
      I*Integrate[
        1/(2) D[S[n - 1, r, e], {r, 2}] - A[n - 1, r, e], {e, 0, t}]],
     Assumptions -> 
     Element[r | t | \[Sigma] | l | \[Alpha] | \[Beta], Reals]];

A[n_Integer, y_, t_] := 
  A[n, y, t] = 
   Refine[Sum[S[i, y, t]*S[n - i, y, t], {i, 0, n}], 
    Assumptions -> 
     Element[y | t | \[Sigma] | l | \[Alpha] | \[Beta], Reals]];


SN[n_Integer, r_, t_] := Sum[S[i, r, t], {i, 0, n}]


h = SN[5, r, t]

However, when I use this type of construction for more complicated equations or when I want to get higher-order terms, the running time starts to be of the order of hours. Is there a way to make this code more optimized?

Thank you all in advance.

Answer 1

The trick is to use partial memoization: only memoize the index $n$ but not the function parameters $(x,t)$,

Clear[S, A, SN];
S[0] = Function[{x, t}, Exp[-x^2]];
A[n_Integer?NonNegative] := A[n] = Function[{x, t},
  Evaluate@Sum[S[i][x,t] * S[n-i][x,t], {i, 0, n}]]
S[n_Integer?Positive] := S[n] = Function[{x, t},
  Evaluate[Expand[I*Integrate[D[S[n-1][x,τ], {x,2}] - A[n-1][x,τ], {τ, 0, t}]]]]
SN[n_Integer?NonNegative] := SN[n] = Function[{x, t},
  Evaluate@Sum[S[i][x,t], {i, 0, n}]]

Now evaluation is pretty quick:

SN[5][x, t]
(*    E^-x^2 - I E^(-2 x^2) t - 2 I E^-x^2 t - E^(-3 x^2) t^2 - 4 E^(-2 x^2) t^2 -
      6 E^-x^2 t^2 + I E^(-4 x^2) t^3 + 6 I E^(-3 x^2) t^3 + 56/3 I E^(-2 x^2) t^3 +
      20 I E^-x^2 t^3 + E^(-5 x^2) t^4 + 8 E^(-4 x^2) t^4 + 110/3 E^(-3 x^2) t^4 +
      96 E^(-2 x^2) t^4 + 70 E^-x^2 t^4 - I E^(-6 x^2) t^5 - 10 I E^(-5 x^2) t^5 -
      892/15 I E^(-4 x^2) t^5 - 3628/15 I E^(-3 x^2) t^5 - 2656/5 I E^(-2 x^2) t^5 -
      252 I E^-x^2 t^5 + 4 I E^-x^2 t x^2 + 12 E^(-2 x^2) t^2 x^2 + 24 E^-x^2 t^2 x^2 -
      68/3 I E^(-3 x^2) t^3 x^2 - 368/3 I E^(-2 x^2) t^3 x^2 - 120 I E^-x^2 t^3 x^2 -
      106/3 E^(-4 x^2) t^4 x^2 - 952/3 E^(-3 x^2) t^4 x^2 - 1008 E^(-2 x^2) t^4 x^2 -
      560 E^-x^2 t^4 x^2 + 248/5 I E^(-5 x^2) t^5 x^2 + 3088/5 I E^(-4 x^2) t^5 x^2 +
      17304/5 I E^(-3 x^2) t^5 x^2 + 38848/5 I E^(-2 x^2) t^5 x^2 + 2520 I E^-x^2 t^5 x^2 -
      8 E^-x^2 t^2 x^4 + 224/3 I E^(-2 x^2) t^3 x^4 + 80 I E^-x^2 t^3 x^4 +
      808/3 E^(-3 x^2) t^4 x^4 + 3776/3 E^(-2 x^2) t^4 x^4 + 560 E^-x^2 t^4 x^4 -
      9872/15 I E^(-4 x^2) t^5 x^4 - 91696/15 I E^(-3 x^2) t^5 x^4 - 74112/5 I E^(-2 x^2) t^5 x^4 -
      3360 I E^-x^2 t^5 x^4 - 32/3 I E^-x^2 t^3 x^6 - 320 E^(-2 x^2) t^4 x^6 -
      448/3 E^-x^2 t^4 x^6 + 6688/3 I E^(-3 x^2) t^5 x^6 + 38144/5 I E^(-2 x^2) t^5 x^6 +
      1344 I E^-x^2 t^5 x^6 + 32/3 E^-x^2 t^4 x^8 - 15872/15 I E^(-2 x^2) t^5 x^8 -
      192 I E^-x^2 t^5 x^8 + 128/15 I E^-x^2 t^5 x^10    *)

If you prefer an interface of the form S[n, x, t] etc., you can define these:

S[n_, x_, t_] := S[n][x, t]
A[n_, x_, t_] := A[n][x, t]
SN[n_, x_, t_] := SN[n][x, t]