Python Lab 3 - Cart1233/PHYS250 GitHub Wiki
Goal
The goal of this lab is to use the Euler method to model a system that is changing with time. To describe how the system changes with time, then we can determine how the system changes by stepping it forward in time using tiny intervals of time. This is the basic idea behind the Euler method.
Overview
I created a series a functions to use the Euler method to model a system that changes with time, and plotted my results.
Background
Ordinary Differential Equation
An ordinary differential equation (ODE) is a mathematical equation that relates a function of one variable to its derivatives with respect to that variable. The "ordinary" in "ordinary differential equation" distinguishes it from partial differential equations.
Taylor Expansion
A Taylor expansion is a math technique used to approximate a function using a polynomial series around a specific point. The usefulness of Taylor expansions lies in their ability to approximate complex functions with simpler polynomial functions, near a specific point.
Euler Method
The Euler Method is a numerical technique used to approximate the solution of an ordinary differential equation (ODE) with an initial value problem. The general form for the Euler Method for solving a first-order ODE is given as: $y_[n+1]=y_n+h*f(x_n,y_n)$
This method is particularly useful for approximating solutions to differential equations because it provides a simple way to obtain numerical solutions to ODEs, especially when analytical solutions are difficult or impossible to find.
Lab
How the Euler method was applied
We used the Euler method to calculate how many radioactive nuclei survive as a function of time. We first started by solving the ODE by separating the variables and integrating each side. This end up giving an analytic solution, which means we do not have to solve it using the Euler method but, it provides a nice way to introduce us to the Euler method and how it affects the accuracy of our result.
How the Euler method affected the accuracy
Our calculated values for the Euler method were very similar to those of the analytic solution as shown is the graphs below:
We were able to create the plot of the Euler method using our functions we created, with the first being our initial to set up the variables and give them specific equations.
def initialize(N0,tau,tmax,dt):
"""This function takes in variables N0,tau,tmax,dt and returns time and the number of surviving nuclei in an element """
nsteps=int(tmax/dt)
t = np.zeros(nsteps)
n_nuclei = np.zeros(nsteps)
n_nuclei[0]=N0
return t,n_nuclei
This equation was followed by our calculate function which will be used to give our variables their assigned values.
def calculate(t,n_nuclei,dt,tau):
for i in range(1,len(t)):
n_nuclei[i] = n_nuclei[i-1] - n_nuclei[i-1]/tau * dt
t[i] = t[i-1] + dt
return t,n_nuclei
How the time steps affected our results
We calculated the number of nuclei versus time using a time steps of 0.05 sec, 0.2 sec, and 0.5 sec. Time steps is our defined variable dt. We saved the output for each to a different name so we could plot all three together with the following code:
# code to get number of nuclei and time for dt = 0.05 s
dt1=.05
initial1=initialize(N0,tau,tmax,dt1)
calc1=calculate(initial1[0],initial1[1],dt1,tau)
# code to get number of nuclei and time for dt = 0.2 s
dt2=.2
initial2=initialize(N0,tau,tmax,dt2)
calc2=calculate(initial2[0],initial2[1],dt2,tau)
# code to get number of nuclei and time for dt = 0.5 s
dt3=.5
initial3=initialize(N0,tau,tmax,dt3)
calc3=calculate(initial3[0],initial3[1],dt3,tau)
Then using this code we created our plot using this code:
# code to plot the 3 different results
plt.figure()
plt.plot(calc1[0],calc1[1],'m*',label='dt=.05',alpha=.5)
plt.plot(calc2[0],calc2[1],'ko',label='dt=.2',alpha=.5)
plt.plot(calc3[0],calc3[1],'bs',label='dt=.5',alpha=.5)
plt.title('Number of nuclei vs. Time')
plt.xlabel('Time')
plt.ylabel('Number of Nuclei')
# code to plot the analytic solution
plt.plot(tline,N_t,'g-',label="Analytic solution")
plt.title('Number of nuclei vs. Time with the analytical solution')
plt.xlabel('Time')
plt.ylabel('Number of Nuclei')
plt.legend()
This gave us a plot that looks like this:
As you can see, when the time steps increase the slope decreases, meaning that there is more nuclei at a specific time with a lower time step compared to a higher time step. This affects our results juristically since with more time steps it becomes less accurate since there aren’t as many plottable points.