3
$\begingroup$

How can I implement this pseudocode in Mathematica?

Algorithm:

(* input: f, x0, xmax, y0, n  *)

h = (xmax-x0)/n;
u[1] = x0;
v[1] = y0;
For i = 2,n+1:
    u[i] = u[i-1] + h;
    v[i] = v[i-1] + h * f( u[i-1] , v[i-1] );
End for.
Return {{u[1],v[1]},{u[2],v[2]},...,{u[n+1],v[n+1]}};

This is my take on it but it's not working properly.

MyEuler[f, x0, xmax, y0, n] = Module[{},
  h = (xmax - x0)/n;
  u[1] = x0;
  v[1] = y0;
  For[i = 2, i <= n + 2, i++,
    u[i] = u[i - 1] + h;
    v[i] = v[i - 1] + h*f (u[i - 1], v[i - 1]);]
   Table[{u[k], v[k]}, {k, 1, n + 1}]]
$\endgroup$
5
  • 3
    $\begingroup$ f (u[i - 1], v[i - 1]) should probably be f[u[i - 1], v[i - 1]] if f is a function. And why end your loop at at n+2 instead of n+1 is a difference, but perhaps an unimportant one. Other than that, read a good introductory text on how to program Mathematica. You aren't taking good advantage of its features, other than its quasi-ability to pretend to be C or Java. $\endgroup$
    – Michael E2
    Nov 3 '21 at 12:18
  • 3
    $\begingroup$ Possible duplicates: (22042), (22413*), (41783), (43623), (206106), especially (22413). $\endgroup$
    – Michael E2
    Nov 3 '21 at 12:29
  • 2
    $\begingroup$ Further, arguments of a function are defined by underscore like: `MyEuler[f_,x0_, ....] $\endgroup$ Nov 3 '21 at 13:06
  • 1
    $\begingroup$ Are you aware of the fact that Euler method is already implemented into the standard NDSolve function and can be fixed by the option Method? It is just to make sure that we understood you correctly. $\endgroup$ Nov 4 '21 at 15:40
  • $\begingroup$ Yes I am aware,but my teacher specifically told us to solve the problem this way. $\endgroup$
    – Angry
    Nov 4 '21 at 22:18
4
$\begingroup$

There are several syntactical mistakes in your code.

  1. Function arguments are missing an underscore:

MyEuler[f_, x0_, xmax_, y0_, n_]

  1. You should use SetDelayed instead of Set, that is: you are missing a colon:

MyEuler[f_, x0_, xmax_, y0_, n_] :=

  1. Calling a function requires square brackets, not round ones:

f[u[i - 1], v[i - 1]]

  1. Separate For from Table by the additional semicolon:

For[ ... ]; Table[ ... ]

Put all of this together, and the code should work:

MyEuler[f_, x0_, xmax_, y0_, n_] := Module[{},
  h = (xmax - x0)/n;
  u[1] = x0;
  v[1] = y0;
  For[i = 2, i <= n + 2, i++,
   u[i] = u[i - 1] + h;
   v[i] = v[i - 1] + h*f[u[i - 1], v[i - 1]]
   ];
  Table[{u[k], v[k]}, {k, 1, n + 1}]
  ]

f[x_, y_] := x + y;

MyEuler[f, 0, 1, 0, 3]

(* {{0, 0}, {1/3, 0}, {2/3, 1/9}, {1, 10/27}} *)
$\endgroup$
2
  • $\begingroup$ We could put h = (xmax - x0)/n, u, v} inside the currently empty 1st argument of Module. And given that u and v are given initial values named x0 and y0, it would be sensible to change indexing to begin with 0 instead of 1 and to then run i only up to n. $\endgroup$
    – murray
    Nov 3 '21 at 14:53
  • $\begingroup$ @murray, yes, certainly there can be several improvements of the code (as pointed out in other answers). But since it looks like the OP is a beginner in Mathematica, I tried to make only the most necessary changes to make the code work. $\endgroup$
    – Domen
    Nov 3 '21 at 15:10
6
$\begingroup$

In a functional style, I recommend NestList, the key trick of which is to properly define an iteration function to be used repeatedly via pure/anonymous Function (&). And more can be found out in @MichaelE2 's comments.

Clear[myEuler]
myEuler[f : (_Symbol | _Function), x0_Real, xmax_Real, y0_Real, n_Integer] :=
Module[{h = (xmax - x0)/n, iterate},
    iterate = # + h {1, f @@ #} &;
    NestList[iterate, {x0, y0}, n]
]

Note that, considering the numeric nature of the problem, here I set explicitly type checking for each argument of myEuler. The consequence is that unless all the type checking returns True, myEuler will refuse to run. In this way would it be disentangled with MMA's symbolic ability to some extend, which is actually not needed here.

As an example, I set $ f(x, y) = x^2 + y^2, x_0 = 0., x_\text{max} = 1., y_0 = .1 $, and $ n = 100 $, and get the result plotted as below:

enter image description here

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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