# Error in both numerical methods in Mathematica

Computed in Mathematica to obtain which method is efficient using residuals: Please what am I doing wrong because it not giving any output?

y[t], {t, 0, 0.5}, Method -> Fehlberg45]

s9 = NDSolve[{y'[t] == -100 y[t] - 99 Exp[2 t], y[0] == 0},
y, {t, 0, 0.5}, Method -> CRK4]

sys = {D[Y[t], t] == -100 Y[t] - 99 Exp[2 t]};
residuals = sys /. Equal -> Subtract;
LogPlot[Join[Abs[residuals /. s2], Abs[residuals /. s9]] //
Evaluate, {t, 0, 0.5}, PlotStyle -> {Red, Directive[Dashed, Green]},
PlotLegends -> {"CRK4", "Fehlberg45"}, Frame -> True, PlotRange -> All]`
• Welcome to the Mathematica Stack Exchange. Please include definitions for s2 and s9.
– Syed
Commented Feb 24, 2022 at 17:40
• To second @Syed comment, what you have included here is only the definition of equations, and it's missing the actual solution where you used NDSolve I assume? Commented Feb 24, 2022 at 17:50
• @MathX. It is been edited and can be viewed. Commented Feb 24, 2022 at 18:04

Your question is incomplete. A lot of your code is missing, I didn't see your method's definitions so I found them from the documentations. Change them if you are not happy with these definitions:

Fehlbergamat = {{1/4}, {3/32, 9/32}, {1932/2197, -7200/2197,
7296/2197}, {439/216, -8, 3680/513, -845/4104}, {-8/27,
2, -3544/2565, 1859/4104, -11/40}};
Fehlbergbvec = {25/216, 0, 1408/2565, 2197/4104, -1/5, 0};
Fehlbergcvec = {1/4, 3/8, 12/13, 1, 1/2};
Fehlbergevec = {-1/360, 0, 128/4275, 2197/75240, -1/50, -2/55};
FehlbergCoefficients[4, p_] :=
N[{Fehlbergamat, Fehlbergbvec, Fehlbergcvec, Fehlbergevec}, p];

Fehlberg45 = {"ExplicitRungeKutta",
"Coefficients" -> FehlbergCoefficients, "DifferenceOrder" -> 4,
"EmbeddedDifferenceOrder" -> 5, "StiffnessTest" -> False};
sfehlberg45 =
NDSolve[{D[y[t], t] == -100 y[t] - 99 Exp[2 t], y[0] == 0},
y, {t, 0, 0.5}, Method -> Fehlberg45]
CRK4[]["Step"[rhs_, h_, t_, x_, xp_]] :=
Module[{k0, k1, k2, k3}, k0 = h xp;
k1 = h rhs[t + h/2, x + k0/2];
k2 = h rhs[t + h/2, x + k1/2];
k3 = h rhs[t + h, x + k2];
(k0 + 2 k1 + 2 k2 + k3)/6]
CRK4[___]["StepInput"] = {"Function"["Time", "DependentVariables"],
"TimeStep", "Time", "DependentVariables", "TemporalDerivatives"};
CRK4[___]["StepOutput"] = "DependentVariablesIncrement";
CRK4[___]["DifferenceOrder"] := 4;
CRK4[___]["StepMode"] := Fixed;
scrk4 = NDSolve[{D[y[t], t] == -100 y[t] - 99 Exp[2 t], y[0] == 0},
y, {t, 0, 0.5}, Method -> CRK4]

And defining the residuals similar to documentations:

residual[t_] = D[y[t], t] + 100 y[t] + 99 Exp[2 t];
Plot[Evaluate[
RealExponent[{residual[x] /. scrk4,
residual[x] /. sfehlberg45}]], {x, 0, 0.5},
PlotStyle -> {GrayLevel[0], RGBColor[1, 0, 0]}, AxesOrigin -> {0, 0},
PlotLegends -> {{"CRK4", "Fehlberg45"}}]

• Thanks for the solution. Actually, I used the methods defined in the documents, I did not know it was necessary to put them here. However, I am confused about the image, as it stands, which method would say is the best fit or can be concluded to give the best approximation. I asked this because both graphs and numerical values are the same and I need residuals to be able to conclude which numerical technique is the best fit. I anticipate your response. Commented Feb 24, 2022 at 19:19
• You can remove the RealExponent function and plot in log scale with 'LogPlot. However, since residuals can be negative, in that case it would be better to define them as absolute value: residual[t_] = Abs[D[y[t], t] + 100 y[t] + 99 Exp[2 t]] and then LogPlot[Evaluate[{residual[x] /. scrk4, residual[x] /. sfehlberg45}], {x, 0, 0.5}, PlotStyle -> {GrayLevel[0], RGBColor[1, 0, 0]}, PlotLegends -> {{"CRK4", "Fehlberg45"}}]. Commented Feb 24, 2022 at 19:49
• Based on these results CRK4 is better because it has residuals with magnitude in orders of less than ~10^-9 whereas Fehlberg45 has residuals with magnitude in orders of less than ~10^-4. Commented Feb 24, 2022 at 19:51
• Thanks for the clarity, you just confirmed my conclusions. If you don't mind, please can you send me a mail so we can collaborate or for possible collaborations? [email protected] Commented Feb 24, 2022 at 21:13