how can I dissect a ode in 3 parts ! lets say I am starting from [0 t1] then stop the ode add a value to a parameter again start the ode[ t1 t2 ]then again change it back ?

조회 수: 2 (최근 30일)
sd()
Warning: RelTol has been increased to 2.22045e-14.
Warning: RelTol has been increased to 2.22045e-14.
%this is the code, lets say i dont add anything to the code then i should get back the original cycles
function sd
a = 0.021;
b = 0.025;
Dc = 100; % in microns
sigma = 1e1; % in MPa
mu_0 = 0.6;
Gmod = 3e4; % in MPa
v_shear = 3e9; % in microns/sec
rad = 0.5*Gmod/v_shear; % Radiation damping term
v_dyn = (a*sigma)/rad;
Kc = sigma*(b-a)/Dc; % in MPa/microns
v_init = 1e-2; % in microns/s
theta_init = Dc/v_init; % in seconds
tspan1=[0 2.5e7];
solopt = odeset('RelTol',1e-16);
t_step = 1e5; % in seconds
Im_stiff_osc = sigma*(sqrt(b)-sqrt(a))^2/Dc;
rat_oscillatory_sol = Im_stiff_osc/Kc;
stiff = 0.7*Im_stiff_osc; % Stiffness
B_s=mu_0*sigma;
law = 'A';
par = [a,b,Dc,sigma,stiff,rad];
taudot_step=[0 9*10e-9*B_s 2.6*10e-9*B_s];
[t1,y11] = ode45(@(t,y)ode5(t,y,par,v_init,t_step,law),...
tspan1,[theta_init;v_init],solopt);
theta_init=y11(end,1);
v_init=y11(end,2);
tspan2=tspan1+2.51e7;
[t2,y12]=ode45(@(t,y)ode5(t,y,par,v_init,t_step,law),...
tspan2,[theta_init;v_init],solopt);
figure()
semilogy(t1,y11(:,2),'r')
hold on
semilogy(t2,y12(:,2),'g')
ylabel ('$velocity$', 'Interpreter','latex')
xlabel ('time[s]')
hold off
end
function dydt = ode5(t,y,par,v_lp,t_step,law)
aa = par(1); bb = par(2); dc = par(3); sigma = par(4); stiff = par(5);
rad = par(6);
dydt = zeros(2,1);
omega = y(1)*y(2)/dc;
dydt(1) = 1.d0 - omega;
if t <= t_step
v_loc = v_lp;
elseif t > t_step
v_loc=1.000001*v_lp;
else
v_loc = 1.00000*v_lp;
end
tau_dot = stiff*(v_loc- y(2));
if strcmp(law,'A')
dydt(1) = 1.d0 - omega;
elseif strcmp(law,'S')
dydt(1) = -omega*log(omega);
end
dydt(2) = ((tau_dot/sigma - bb*dydt(1)/y(1)))/(aa/y(2)+rad/sigma);
end
  댓글 수: 4
SOUVIK DARIPA
SOUVIK DARIPA 2023년 4월 14일
I will, but first I have to ensure that I am getting back the same cycle Without adding anything else.
Torsten
Torsten 2023년 4월 14일
Yes, of course you will get the same plot if you don't change anything in the equations. But there seem to be singularities in your solution function.

댓글을 달려면 로그인하십시오.

답변 (0개)

카테고리

Help CenterFile 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