solving the equation using rk-4 method, but couldn't get the right output, last time I forget to add the program please check this one.

1 次查看（过去 30 天）

SAHIL SAHOO 2022-7-20

评论： Torsten 2022-7-21

i want to solve this, where σ = 0.1, L0 = 0.5, k= 0.01, tc = 70E-9

format long

a = 0;

b = 4;

h = 0.1;

k = 0.01; % critical coupling strength

L = 0.5;

s = 0.1;

tc = 70E-9; % photon life time in cavity

n = (b-a)/h ;

t = a + (0:n)*h;

func1 = @(t,y) [(-(s^2)*k/tc)*sin(y - pi/2) + L*(s^2)/(2*tc)*sin(y + pi/2)/sqrt(1 + cos(y + pi/2))];

for i = 1:n

k1 = func1(y(i));

k2 = feval(func1, y(i)+(k1/2)*h);

k3 = feval(func1, y(i)+(k2/2)*h);

k4 = feval(func1, y(i)+(k3*h));

y(i+1) = y(i) + (h/6)*(k1+(2*(k2+k3))+k4);

end

plot(t,y)

Walter Roberson 2022-7-21

You should preallocate y as length n+1. With your current code if you just happen to have a longer y in your workspace you would have a problem.

Torsten 2022-7-21

Works for me (see above).

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

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

Translated by