Re: How can I solve equation of motion?

[Przemek Klosowski]

   Now,I research to solve the equation of motion by using octave.
   The following is the example of equations.

   (dx/dt)''=F - (dx/dt)'

   #"F" is the externalforce which is measured by experiment.This value is
   refered to data of experiment every time.

   I came up with the program to solve this equation as follows.
   But,this program does not seem to work correctly.

   F=rand(5,1); # Now, I give randam "F"
   function dx = ex(x,t)
           global n;
           dx(1) = F(n+1,1)-x(1);
           dx(2) = x(1);
   endfunction

   global n;
   v=0; # set initial value
   X=0; # set initial value
   for n=0:4
           x0 = [v;X]   # initial value
           t = linspace(n,n+1,2)
           y = lsode("ex",x0,t)
           v=(2,1)              # new initial value
           X=(2,2)              # new initial value
   end

Traditionally, the Newton equation of motion is written as F = m x''
Your equation has inconsistent dimensions (the units on the left are
m/s^3, while the right side subtracts m/s^2 from kg*m/s^2 (assuming
that you define your force in 

You would solve F=mx'' by converting it into a vector first order
equation: X=[x,x']; then X'=[x',x''] = [x', F/m] = [0 1; 0 0] * X + [0 F/m]. 

X0=[0 0];  # X0 = [x,v];
T=linspace(0,10,100);
function xp = ex(x,t)
xp=[0 1 ; 0 0]*x + [0 3]'; # F/m = 3, for instance
endfunction
r=lsode("ex",X0,T);

gives what I think is the expected result (x = 3/2 t^2, x' = 3t). Note
that you don't need to split up the interval and solve separately---you 
can just modify F in the function F, depending on the interval.



-------------------------------------------------------------
Octave is freely available under the terms of the GNU GPL.

Octave's home on the web:  http://www.octave.org
How to fund new projects:  http://www.octave.org/funding.html
Subscription information:  http://www.octave.org/archive.html
-------------------------------------------------------------