1
$\begingroup$

The initial value problem is : $$P'(t)=0.7P(t)(1-\frac{P(t)}{750})-20, P(0)=30$$

The time step is set to $Δt=7$ days

For the algorithm we have:

$f(t,P)=0.7P(1-\frac{P}{750})-20$

with step $h=7$

$t_0=0$.

Euler's method reads in this case $$ P_0 = 30, \quad P_{i+1} = P_i + h \left( 0.7 P_i \left(1- \frac{P_i}{750}\right) -20\right) $$.

My target is to create a table and a plot that will change every time I will the step $Δt$

I have used the following commands

p'[t] == 7/10 p[t]*(1 - p[t]/750) - 20;
p[0]==30
h==7
[![NDSolve\[{p'\[t\] == 7/10 p\[t\]*(1 - p\[t\]/750) - 20, p\[0\]==30}, p\[t\], {t, 0, 1},
 Method -> "ExplicitEuler", "StartingStepSize" -> 7\]][1]][1]

My target is to create a table and a plot that will change every time I will the step $Δt$ and shows $P(T)$ as in the picture below. Is it possible?

enter image description here

Edit:

My professor shared with me the following Python code. However, I am not used to using Python, and I would like to run it in Mathematica. Is it possible? Could anyone help me?

# Program      : Euler's method
# Author       : MOOC team Mathematical Modelling Basics
# Created      : April, 2017

import numpy as np
import matplotlib.pyplot as plt

print("Solution for dP/dt = 0.7*P") # in Python 2.7: use no brackets

# Initializations

Dt = 0.1                                # timestep Delta t
P_init = 10                             # initial population 
t_init = 0                              # initial time
t_end = 5                               # stopping time
n_steps = int(round((t_end-t_init)/Dt)) # total number of timesteps

t_arr = np.zeros(n_steps + 1)           # create an array of zeros for t
P_arr = np.zeros(n_steps + 1)           # create an array of zeros for P
t_arr[0] = t_init                       # add the initial P to the array
P_arr[0] = P_init                       # add the initial t to the array

# Euler's method

for i in range (1, n_steps + 1):
    P = P_arr[i-1]
    t = t_arr[i-1]
    dPdt = 0.7*P                        # calculate the derivative 
    P_arr[i] = P + Dt*dPdt              # calculate P on the next time step
    t_arr[i] = t + Dt                   # adding the new t-value to the list

# Plot the results

fig = plt.figure()                      # create figure
plt.plot(t_arr, P_arr, linewidth = 4)   # plot population vs. time

plt.title('dP/dt = 0.7P, P(0)=10', fontsize = 25)  
plt.xlabel('t (in days)', fontsize = 20)
plt.ylabel('P(t)', fontsize = 20)

plt.xticks(fontsize = 15)
plt.yticks(fontsize = 15)
plt.grid(True)                          # show grid 
plt.axis([0, 5, 0, 200])                # define the axes
plt.show()                              # show the plot
# save the figure as .jpgde 
```
$\endgroup$
3
  • $\begingroup$ Here is your cleaned up code that produces a result: solution[t_] = p[t] /. NDSolve[{p'[t] == 7/10 p[t]*(1 - p[t]/750) - 20, p[0] == 30}, p[t], {t, 0, 1}, Method -> "ExplicitEuler", "StartingStepSize" -> 7][[1]] $\endgroup$ Dec 26, 2022 at 17:58
  • $\begingroup$ @DanielHuber Thanks! I have edited a bit my question with some Python Code my professor has been shared with me $\endgroup$ Dec 26, 2022 at 18:02
  • $\begingroup$ If you want to actually implement the Euler's method and not blindly use NDSolve, then there are several questions on this StackExchange about the implementation, such as this example. $\endgroup$
    – Domen
    Dec 26, 2022 at 18:13

2 Answers 2

6
$\begingroup$
ode = p'[t] == 7/10 p[t]*(1 - p[t]/750) - 20;
ic = p[0] == 30
sol = p /. First@NDSolve[{ode, ic}, p, {t, 0, 1}, Method -> "ExplicitEuler", StartingStepSize" -> 7];
exact = p[t] /. First@DSolve[{ode, ic}, p[t], t]
sol["Methods"]

Mathematica graphics

The above gives you all properties of the numerical solution. You can access the data used and compare to exact solution as follows

y = sol["ValuesOnGrid"] // Chop
x = First@sol["Coordinates"]
data = Transpose[{x, y}];
p1 = ListPlot[data, PlotStyle -> Red, GridLines -> Automatic, GridLinesStyle -> LightGray];
p2 = Plot[exact, {t, 0, 1}, PlotStyle -> Blue];
Show[p1, p2, PlotLabel -> "Exact vs. Euler"]

Mathematica graphics

Euler method gets worst with time. You can see this if you change t to 10 seconds instead of 1

Mathematica graphics

$\endgroup$
6
$\begingroup$

MaxStepFraction is needed if you want fewer than ten steps.

Manipulate[
 With[{sol = Quiet[
     NDSolveValue[{p'[t] == 7/10 p[t]*(1 - p[t]/750) - 20, 
       p[0] == 30}, p, {t, 0, n*h}, 
      Method -> {"FixedStep", Method -> "ExplicitEuler"}, 
      StartingStepSize -> h, MaxStepFraction -> 1, MaxSteps -> n],
     NDSolveValue::mxst]},
  foo = sol;
  Grid[{{Style["Results Euler's Method", 
      "Subsubsection"]}, {ListLinePlot[sol, Mesh -> All, 
      PlotRange -> {{0, 60}, {-10, 1000}}],
     Column[{
       Row[{"\[CapitalDelta]t = ", h, " days"}],
       "",
       Grid[{{"n", "t", HoldForm[P[t]], HoldForm[P'[t]]},
         {n, n*h, sol[n*h], sol'[n*h]}},
        Dividers -> All]
       }]}
    }]
  ],
 {{h, 7, "\[CapitalDelta]t"}, 1., 15},
 {{n, 3, "steps"}, 1, 60, 1}
 ]

Mathematica graphics

$\endgroup$
2
  • $\begingroup$ Thank you! This is what I need it !! I appreciate it $\endgroup$ Dec 26, 2022 at 18:47
  • $\begingroup$ @AthanasiosParaskevopoulos You're welcome $\endgroup$
    – Michael E2
    Dec 27, 2022 at 0:59

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