Talk:Euler–Maruyama method

Range of the index of Y_n

In the line

"recursively define Y_n for 1 ≤ n ≤ N"

shouldn't the index be 0 ≤ n ≤ N-1?

Otherwise Y₁ will not be computed, since Y_{n + 1} = Y₂ for n = 1.

Moreover, the last value to be computed will be Y_N+1, which is outside the interval of time [0, T] partitioned by 𝞃₀ < 𝞃₁ < ... < 𝞃_N. — Preceding unsigned comment added by Rabelais (talk • contribs) 16:49, 29 February 2020 (UTC)

Wrong time partition in the example implementations

In both Python and Matlab codes, Δt is defined as

dt = (t_end - t_init) / N

which, referring to the theoretical definition, implies a partition of [t_init t_end] in N+1 points.

But the partition is defined as

Python: ts = np.arange(t_init, t_end, dt)

which is a vector of N elements (should be N+1) given by [t_init, t_init+dt, t_init+2*dt, ..., t_end-dt], hence it does not contain the final time;

Matlab: ts = linspace(tBounds(1), tBounds(2), N)

which is a vector of N elements (should be N+1) given by [t_init, t_end/(N-1), 2*t_end/(N-1), ..., t_end], so it is a totally different time partition (excluding the extremes), even though, by the numerical point of view, it is not important what the partition contains but the number of the elements contained.

Then, the internal 'for' cycle devoted to the computation of the approximations is defined as

Python: for i in range(1, ts.size)

hence 'i' goes from 1 to N-1 (should be from 1 to N) and the approximation corresponding to t_end is not computed;

Matlab: for i = 2:numel(ts)

hence 'i' goes from 2 to N (should be from 2 to N+1) and the approximation corresponding to t_end is not computed.

The theoretical definition states that

• the time interval [t_init, t_end] is partitioned in N interval, i.e. in N+1 nodes

• the length of each subinterval is t_end/N

• an approximation for each node has to be computed, hence a total of N+1 approximations have to be computed (the first approximation Y_0 is given by the initial condition)

So either Δt should be changed to

dt = (t_end - t_init) / (N - 1)

or ts should be changed to

Python: ts = np.arange(t_init, t_end+dt, dt)

and to

Matlab: ts = linspace(tBounds(1), tBounds(2), N+1) (ts = tBounds(1):dt:tBounds(2) is fine too )

In the python code, the relation 'assert len(TS) == TS.size == N + 1' holds true and a for-loop starts at 0. I made everything explicit in the code now, with an assert. I.e. the python code uses a grid of 1001 points and does 1000 computations (as opposed to having a grid of 1000). I didn't touch the Matlab code, but I think it's fairly clear.

91.141.56.40 (talk) 18:36, 11 January 2023 (UTC)