gabrielmp

5 Reputation

2 Badges

12 years, 201 days

MaplePrimes Activity


These are answers submitted by gabrielmp

Thank you very much for your help Preben. I probably wouldn't have done it without you. I appreciate it !

Thanks for your help. Preben. As you said, i changed things a little bit in the program. Now it looks like this: 

 

restart;

A := Matrix();

dt := 0.1e-2;

diff(xo(t), t) := vo;

diff(vo, t) := a(xo(t), vo);

xo(t) := 1; vo := 0;

for i to 1000 do

xn := xo(t)+vo*dt; vn := vo+a*dt;

dx1 := vn*dt; dv1 := a*dt;

dx2 := (vn+(1/2)*dv1)*dt; dv2 := a(x+(1/2)*dx1, (1/2)*dv1)*dt;

dx3 := (vn+(1/2)*dv2)*dt; dv3 := a(xn+(1/2)*dx2, vn+(1/2)*dv2)*dt;

dx4 := (vn+dv3)*dt; dv4 := a(xn+dx3, vn+dv3)*dt;

xn := evalf((xo(t)+dx1+2*dx2+2*dx3+dx4)*(1/6)); vn := evalf((vo+dv1+2*dv2+2*dv3+dv4)*(1/6));

A(i, 1) := xn; A(i, 2) := vn; A(i, 3) := i*dt;

xo(t) := xn;

vo := vn

end do;

ExportMatrix("RK4.txt", A)

 

One thing I didn't understand is what do you mean with "a" not being concrete.

Thanks once again

Page 1 of 1