Implementierung des Euler-Cauchy-Verfahrens

Clear[EulerCauchy]
EulerCauchy[g_, {x_, x0_}, {y_, y0_}, h_, n_] := 
      Module[{j, plot, Achsen}, 
X[j_] := 
          x0 + j h; 
Y[0] = 
          y0; 
plot[0] = 
          Graphics[{PointSize[0.02], 
              Point[{X[0], Y[0]}]}]; 
Do[
          Y[j] = Y[j - 1] + 
              h (g /. {x -> X[j - 1], 
                      y -> Y[j - 1]}); 
plot[j] = 
            Graphics[{PointSize[0.02], 
                Point[{X[j], Y[j]}]}]
, {j, 1, 
            n}]; 
Achsen = 
          Plot[{}, {x, x0, n h}, 
            DisplayFunction -> Identity]; 
Show[
          Prepend[Table[plot[j], {j, 0, n}], Achsen], 
          DisplayFunction -> $DisplayFunction]
]

euler1 = EulerCauchy[(y^2), {x, 0}, {y, 1}, 0.1, 10]

Transpose[Table[{j, X[j], Y[j]}, {j, 0, 10}]]

dsol = 
      DSolve[{y'[x] == y [x]^2, y[0] == 1}, y[x], 
        x]

plot = Plot[y[x] /. dsol, {x, 0, 1}]

Show[euler1, plot, PlotRange -> {0, 7}]

euler2 = EulerCauchy[(y^2), {x, 0}, {y, 1}, 0.05, 20]

Table[{j, X[j], Y[j]}, {j, 0, 20}]

Show[euler2, plot, PlotRange -> {0, 10}]

euler5 = EulerCauchy[(y^2), {x, 0}, {y, 1}, 0.02, 50]

Show[euler5, plot, PlotRange -> {0, 20}]

Y[50]

Mathematicas eingebaute Prozedur zur numerischen Lösung von Anfangswertproblemen

sol = 
      NDSolve[{y'[x] == y[x]^2, y[0] == 1}, 
        y[x], {x, 0, 1}]

ndsolve = Plot[y[x] /. sol, {x, 0, 1}]

Show[ndsolve, plot]

Herleitung des Runge-Kutta-Verfahrens

RK = 
      y[x + h] - (y[
              x] + (∑_(j = 0)^(m - 1))A[j] k[j])

RK2 = 
      RK //. {m -> 2, k[0] -> h g[x, y [x ]], 
          k[1] -> 
            h g[x + α[1] h, 
                y [x ] + β[1, 0] k[0 ]]}

ser = Series[RK2, {h, 0, 3}]

ser = 
      ser //. Derivative[k_][y][x] -> 
          D[g[x, y[x]], {x, k - 1}]

coeff = CoefficientList[ser, h ]

Heunverfahren:

Solve[{Coefficient[coeff,Derivative[0, 1][g](x, y(x))]==0,Coefficient[coeff,Derivative[1, 0][g](x, y(x))]==0}/.{A[0] -> (1/2), A[1] -> (1/2)},{α[1], β [1, 0]}]

Modifiziertes Euler-Cauchy-Verfahren:

Solve[{Coefficient[coeff,Derivative[0, 1][g](x, y(x))]==0,Coefficient[coeff,Derivative[1, 0][g](x, y(x))]==0}/.{A[0] -> 0, A[1] -> 1},{α[1], β [1, 0]}]

Klassisches Runge-Kutta-Verfahren:

RK4 = 
      RK //. {m -> 4, k[0] -> h g[x, y [x ]], 
          k[1] -> 
            h g[x + α[1] h, y [x ] + β[1, 0] k[0 ]], 
          k[2] -> 
            h g[x + α[2] h, 
                y [x ] + β[2, 0] k[0 ] + β[2, 1] k[
                      1 ]], 
k[3] -> 
            h g[x + α[3] h, 
                y [x ] + β[3, 0] k[0 ] + β[3, 1] k[
                      1 ] + β[3, 2] k[2 ]]}

ser = Series[RK4, {h, 0, 4}]

ser = 
      ser //. Derivative[k_][y][x] -> 
          D[g[x, y[x]], {x, k - 1}]

coeff = CoefficientList[ser, h ]

coeff[[1]]

coeff[[2]]

eq1 = {coeff[[2]] == 0}

eq2 = coeff[[3]]

eq2 = Collect[eq2, g[x, y[x]]]

eq2 = Collect[eq2, Derivative[n_, m_][g][x, y[x]]]

eq2 = Map[# == 0 &, Apply[List, eq2]]

eq3 = Collect[coeff[[4]], g[x, y[x]]]

eq3 = Collect[eq3, Derivative[n_, m_][g][x, y[x]]]

eq3 = Map[# == 0 &, Apply[List, eq3]]

eq = Join[eq1, eq2, eq3]

eq = eq /. Derivative[n_, m_][g][x, y[x]] -> 1

eq = eq /. g[x, y[x]] -> 1

Length[eq]

wahl1 = {A[0] -> (1/6), A[1] -> (1/3), A[2] -> (1/3), 
        A[3] -> (1/6)}

EQ = eq /. wahl1

wahl2 = {β[2, 0] -> 0, β[3, 0] -> 
          0, β[3, 1] -> 0}

sol = 
      Solve[EQ /. 
          wahl2, {α[1], α[2], α[3], β[1, 
            0], β[2, 0], β[2, 1], β[3, 0], β[3, 
            1], β[3, 2]}]

coeff /. sol[[1]] /. wahl1 /. wahl2 // Simplify

coeff /. sol[[3]] /. wahl1 /. wahl2 // Simplify

coeff /. sol[[5]] /. wahl1 /. wahl2 // Simplify

Runge-Kutta-Daten

Union[Join[wahl1, wahl2, sol[[3]]]]