How can I use Euler's method to solve this ODE problem in MATLAB?

6 次查看（过去 30 天）
Chitra Ram
Chitra Ram 2018-3-9
回答： James Tursa 2018-3-9
The well-know Michealis-Menten reaction equation for enzymatic reactions is a 2 parameter kinetic model,
dC/dt = - Vm * C/(Km + C)
where C is the reactant concentration, and Vm and Km are reaction rate constants.
For Vm = 20 mmol/L-min, Km = 10 mmol/L and at t = 0 min, Co=100 mmol/L,
use the Euler and Heun methods to estimate a numerical solutions for 0<t<10, plotting the concentration as a function of time for the Michealis-Menten equation for the following step sizes:
a) step size = 1 min
b) step size = 0.5 min
c) step size = 0.1 min
What do you observe about the accuracy/precision of the models?
  2 个评论
James Tursa
James Tursa 2018-3-9
What have you done so far? What specific problems are you having with your code? For starters, this appears to be a 1st order ODE, not a 2nd order ODE.
Chitra Ram
Chitra Ram 2018-3-9
编辑：James Tursa 2018-3-9
This is the template of code I was provided with but I don't know how to use my given information to find solutions.
if nargin <5 && nargin==4;
%error('myApp:argChk','Wrong number of input arguments. Read help comments of EULER_forward.m');
display('Wrong number of input arguments. Read help comments of EULER_forward.m');
pause
disp('Hit any key to continue anyway and wish to program make a decision for you')
disp('Default value for number of iterations is taken to be 10 !!! ');
Niter=10; % That is a default value for number of iterations
elseif nargin==3 && nargin<4;
disp('Do you know how many iterations to carry out? ')
Niter=input('Enter the number of iterations, i.e. you need to perform Niter = ');
disp('Do you know end value for t i.e., tend ? ');
tend=input('Enter tend = ');
elseif nargin==2;
disp('Do you know how many iterations to carry out? ')
Niter=input('Enter the number of iterations, i.e. you need to perform Niter = ');
disp('Do you know the end value for t? ');
tend=input('Enter tend = : ');
disp('Do you know initial value for y i.e., y0 ')
y0=input('Enter y0 = ');
elseif nargin==1;
disp('Do you know how many iterations to carry out? ')
Niter=input('Enter the number of iterations, Niter = ');
disp('Do you know the end value for t? ');
tend=input('Enter tend = : ');
disp('Do you know initial value for y i.e., y0 ')
y0=input('Enter y0 = ');
display('Do you know inital value of t i.e., t0');
t0=input('Enter t0 = ');
display('The entry is to be a function called by an anonymous function');
else
disp('Do you know how many iterations to carry out? ')
Niter=input('Enter the number of iterations, Niter = ');
disp('Do you know the end value for t? ');
tend=input('Enter tend = : ');
disp('Do you know initial value for y i.e., y0 ')
y0=input('Enter y0 = ');
display('Do you know inital value of t i.e., t0 ');
t0=input('Enter t0 = ');
disp('Create anonymous function f=@(t,y)expression');
f=input('Enter anonymous function: f = ');
end
hstep=(tend-t0)/Niter;
t=[t0, zeros(1,Niter)]; % Memory allocation
y=[y0, zeros(1,Niter)]; % Memory allocation
for jj=1:Niter
t(jj+1)=t(jj)+hstep;
y_NEW=y(jj)+hstep*(f(t(jj), y(jj)));
y(jj+1)=y(jj)+(hstep/2)*(f(t(jj), y(jj))+f(t(jj+1), y_NEW));
end
% Visualization of computed data
plot(t,y, 'rp--');
grid on;
title('Numerical solution of ODE using EULER modified Method');
axis tight;
end
if true
% code
end

请先登录，再进行评论。

请先登录，再回答此问题。

回答（1 个）

James Tursa
James Tursa 2018-3-9
Looks like you have been given the complete integrator code, and all you need to do is define the derivative function f. Although you don't show it, I am assuming that f is actually an input argument to this integrator routine. Can you confirm this? Assuming this is true, all you need to do is define a function handle in the calling routine that calculates the derivative given the time and the current state ... i.e. your derivative function needs to have the signature dy = f(t,y). So here is what you do with the information you have been given:
Given this ODE:
dC/dt = - Vm * C/(Km + C)
and given these constants:
Vm = 20 mmol/L-min, Km = 10 mmol/L
and given this initial condition:
at t = 0 min, Co=100 mmol/L,
and given this desired time span in minutes:
0<t<10
Define your constants in m-code:
Vm = 20; % mmol/L-min
Km = 10; % mmol/L
Construct your derivative function handle in m-code:
f = @(t,C) = - Vm * C/(Km + C);
Define your initial condition in m-code:
Co = 100; % mmol/L
And define your time span in m-code:
tspan = [0 10]; % min
It also appears that your integrator code needs the number of iterations, Niter. To get that for a desired step size, simply solve this equation that is found in the above code for Niter, so that given a desired hstep you can calculate Niter. The tend and t0 are from the desired time span.
hstep=(tend-t0)/Niter
These are the pieces that you need to solve your problem. Just have these pieces in your caller routine and pass them to your integrator routine as input arguments in the appropriate spots. Give it a try ...

类别

在 Help Center 和 File Exchange 中查找有关 Ordinary Differential Equations 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by