2
$\begingroup$

I am trying to recursively define the following function:

G[t_, x_] := 
 G[t, x] = 
  X0[x] + t*V0[x] - 
   1/2*Integrate[
     Integrate[Sign[G[s, w] - G[s, x]], {w, 0, 1}], {s, 0, t}]

(X0 and V0 are previously defined piecewise constant functions on [0, 1]). Subsequently, I also have tried to set the base case by doing:

G[0, x_] = X0[x]

However, whenever I try to run, for example G[0.5, 0.5] I always run into errors. Typically,I get 3 recursion limits (of depth exceeding 1024), followed by a general error saying further recursion errors will be suppressed, and 2 iteration limits of 4096. Obviously running this recursive integral as is seems to be overloading the memory, so how can I fix this so even if it's just doing it approximately? (When I replace Integrate with NIntegrate the code just runs indefinitely and never terminates).

Thank you.

Edit: Here are the definitions for X0 and V0:

X0[x_] := 
 Piecewise[{{-85/4, 0 <= x <= 1/4}, {15/4, 1/4 < x <= 1/2}, {35/4, 
    1/2 < x <= 1}}]

V0[x_] := 
 Piecewise[{{4, 0 <= x <= 1/4}, {0, 1/4 < x <= 1/2}, {-2, 
    1/2 < x <= 1}}]
Improve this question
$\endgroup$
14
  • 1
    $\begingroup$ How could this ever work? G[t,x] would have to be defined everywhere initially for the integral to be evaluable. Recursion is good for integer-dependent problems. $\endgroup$ Commented Jun 30 at 19:52
  • 1
    $\begingroup$ Could you add an index 'n' to your function, which labels the order of the approximation, and then recurse over 'n'? For n=0, G = X0 + t*V0. This would work well only if the integral is small in some sense. $\endgroup$ Commented Jun 30 at 20:45
  • 2
    $\begingroup$ Would the numerical solution of the integral equation be an alternative to recursion ? $\endgroup$ Commented Jul 1 at 8:08
  • 1
    $\begingroup$ I'd think of this problem as a partial differential equation, not as a recursion. By differentiating your equation with respect to $t$, we get $\frac{\partial G(t,x)}{\partial t}=V_0(x)-\frac12\int_0^1 \text{sign}[G(t,w)-G(t,x)]dw$. Together with $G(0,x)=X_0(x)$ this could be used in a numerical method for a partial differential equation. The presence of the sign function makes it very unlikely that there is a closed-form solution, and I'd go straight to numerical solutions. $\endgroup$ Commented Jul 2 at 6:51
  • 1
    $\begingroup$ A simple forward integration of this differential equation would be, for a small $\tau>0$ and for a large set of points $x\in[0,1]$, to compute $G(\tau,x)=X_0(x)+\tau V_0(x)$, and then $G(2\tau,x)=G(\tau,x)+\tau\left[V_0(x)-\frac12\int_0^1\text{sign}[G(\tau,w)-G(\tau,x)]dw\right]$, where the integral in the last formula is done numerically over the list of computed $x$-values. And so forth to $G(3\tau,x)$ etc. Difficult to say more without setting up a concrete calculation with concrete $X_0(t)$ and $V_0(t)$. Also, a higher-order integration (Runge–Kutta or similar) will stabilize the method. $\endgroup$ Commented Jul 2 at 7:10

3 Answers 3

Reset to default
3
$\begingroup$

Extending my comments into a solution:

First, convert the integral equation into a differential equation:

$$ G(t,x)=X_0(x)+t V_0(x)-\frac12\int_0^t ds \int_0^1 dw\,\text{sign}\left[G(s,w)-G(s,x)\right] $$

$$ \frac{\partial G(t,x)}{\partial t}=V_0(x)-\frac12\int_0^1 dw\,\text{sign}\left[G(t,w)-G(t,x)\right] $$

Then, introduce a set of values $\xi_0,\xi_1,\ldots,\xi_n$ that cover the range $[0,1]$, so that we can do the $dw$ integral numerically. For example, we could set the $\xi_i$ from Gauss–Legendre integration. Here I'll do something simpler and use the trapezoidal rule, setting $\xi_i=i/n$.

The integral over $dw$ is therefore

$$ \int_0^1 dw\,\text{sign}\left[G(t,w)-G(t,x)\right] \approx \frac{1}{2n} \text{sign}\left[G(t,\xi_0)-G(t,x)\right] + \frac{1}{n} \sum_{i=1}^{n-1}\text{sign}\left[G(t,\xi_i)-G(t,x)\right] +\frac{1}{2n} \text{sign}\left[G(t,\xi_n)-G(t,x)\right] $$

The differential equation at point $\xi_j$ is thus

$$ \frac{\partial G(t,\xi_j)}{\partial t}\approx V_0(\xi_j)-\frac{1}{2n}\left[ \frac12 \text{sign}\left[G(t,\xi_0)-G(t,\xi_j)\right] + \sum_{i=1}^{n-1}\text{sign}\left[G(t,\xi_i)-G(t,\xi_j)\right] +\frac12 \text{sign}\left[G(t,\xi_n)-G(t,\xi_j)\right] \right] $$

This is a system of $n+1$ coupled differential equations for the functions $G_j(t)=G(t,\xi_j)$:

$$ G_j(0)= X(\xi_j)\\ G_j'(t)=V_0(\xi_j)-\frac{1}{2n}\left[ \frac12 \text{sign}\left[G_0(t)-G_j(t)\right] + \sum_{i=1}^{n-1}\text{sign}\left[G_i(t)-G_j(t)\right] +\frac12 \text{sign}\left[G_n(t)-G_j(t)\right] \right] $$

By choosing $n$ sufficiently large, we can get a high-resolution approximation for the $x$-dependence.

X0[x_] = Piecewise[{{-85/4, 0 <= x <= 1/4}, {15/4, 1/4 < x <= 1/2}, {35/4, 1/2 < x <= 1}}];
V0[x_] = Piecewise[{{4, 0 <= x <= 1/4}, {0, 1/4 < x <= 1/2}, {-2, 1/2 < x <= 1}}];

n = 5;
ξ[i_] = i/n;
sol[t_] = NDSolveValue[
  Join[Table[G[j][0] == X0[ξ[j]], {j, 0, n}], 
       Table[G[j]'[t] == V0[ξ[j]] - 
         1/(2 n) (1/2 Sign[G[0][t] - G[j][t]] + 
                  Sum[Sign[G[i][t] - G[j][t]], {i, 1, n - 1}] + 
                  1/2 Sign[G[n][t] - G[j][t]]), {j, 0, n}]], 
         Table[G[j][t], {j, 0, n}],
  {t, 0, 10}, 
  Method -> {"DiscontinuityProcessing" -> False}]

Plot[Evaluate@sol[t], {t, 0, 10}, PlotTheme -> "Detailed",
  FrameLabel -> {t, G}, PlotLegends -> (ξ /@ Range[0, n])]

enter image description here

Same with $n=64$ and rendered in 2D:

ListDensityPlot[Table[sol[t], {t, 0, 10, 1/10}],
  ColorFunction -> "Rainbow", PlotLegends -> Automatic, 
  DataRange -> {{0, 1}, {0, 10}}, FrameLabel -> {x, t}]

enter image description here

Improve this answer
$\endgroup$
3
  • $\begingroup$ Nice and clever solution(+1). I tried galerkin method and gaussquadrature but failed with NDSolve because I didn't know Method -> {"DiscontinuityProcessing" -> False}]. Without this option ode including Sign[] can't be solved $\endgroup$ Commented Jul 3 at 12:18
  • 1
    $\begingroup$ FWIW we could also approximate Sign[] with a scaled smooth function like LogisticSigmoid or ArcTan $\endgroup$ Commented Jul 4 at 3:57
  • $\begingroup$ Thank you! Apologies for the delay, this is very interesting. Really appreciate the help. :-) $\endgroup$ Commented Jul 15 at 1:48
1
$\begingroup$

Recursion might be defined using NestList:

iter=1;
soli = NestList[
X0[x] + t*V0[x] - 
1/2*Integrate[ Sign[#[s, w] - #[s, x]], {w, 0, 1} , {s, 0, t}, 
  Assumptions -> {0 < t < 1, 0 < x < 1 }] &, 
Function[{t, x}, X0[x] + t*V0[x]]
, iter]

enter image description here

The first recursion agrees with TomDickens'answer.

Plot3D [soli[[-1]] , {t, 0, 1}, {x, 0, 1}]

![enter image description here

Unfortunately Mathematica isn't able to evaluate the next recursions iter>1.

Probably the problem should be solved numerically (integral equation)

Improve this answer
$\endgroup$
1
$\begingroup$

This should not be an answer, since the method here does not completely work, but it's too long for a comment. I trust that a more skillful Mathematica reader will have suggestions to do this completely correctly.

I implemented an iterative technique that is recursive over the index n of the iteration number, with n=0 your initial condition. It appears that the result converges fairly quickly. It is fairly slow, I think because the piecewise function definition gets used over and over again, resulting in lots of memory use and branching.

My problem is that I do not know how force it to store complete definitions of earlier iterations without doing something like Plot[ ]. So, note that in the pointwise evaluations I show at the end, the code does not calculate a value for the iteration $n=3,$ the next one after the ones I already plotted. Also, it might be better to use NIntegrate[ ] because the exact results would be very unwieldy.

Maybe this can help you a bit.

Regards, Tom

X0[x_] := 
  Piecewise[{{-85/4, 0 <= x <= 1/4}, {15/4, 1/4 < x <= 1/2}, {35/4, 1/2 < x <= 1}}]

V0[x_] := 
  Piecewise[{{4, 0 <= x <= 1/4}, {0, 1/4 < x <= 1/2}, {-2, 1/2 < x <= 1}}]

 G[t_, x_, n_] := G[t, x, n] = 
    X0[x] + t*V0[x] - 
    1/2*Module[{w, s}, 
     Integrate[
      Integrate[Sign[G[s, w, n - 1] - G[s, x, n - 1]], {w, 0, 1}],{s,0,t}]]

 G[t_, x_, 0] = X0[x] + t*V0[x];

 Plot3D[G[t, x, 0], {t, 0, 1}, {x, 0, 1}]

initial value

 Plot3D[G[t, x, 1], {t, 0, 1}, {x, 0, 1}]

After one iteration

 Plot3D[G[t, x, 2], {t, 0, 1}, {x, 0, 1}]

After two iterations

 (* Illustration of values after several iterations - fails for n>2 that I have not already run *)
 G[.1, .2, 0]
 (* -20.85 *)

 G[.1, .2, 1]
 (* -20.8875 *)

 G[.1, .2, 2]
 (* -20.9 *)

 G[.1, .2, 3]
 (* Undefined *)
Improve this answer
$\endgroup$

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.