[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
How can I solve equation of motion?
From: 
John W. Eaton 
Subject: 
How can I solve equation of motion? 
Date: 
Mon, 2 Dec 2002 09:46:55 0600 
On 2Dec2002, NIITSU Taichi <address@hidden> wrote:
 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

 I would be very happy if you have any idea on it and let me know
 Thanks in advance,
The T argument in lsode can be a vector, so you can avoid the loop.
If I understand what you want the random F to do, I think you are
trying to do something like this:
global F
F = rand(5,1);
function dx = ex(x,t)
global F T
## Find value in F corresponding to current time interval.
r = F (max (find (t >= T)));
dx = [rx(1); x(1)];
endfunction
x0 = [0; 0];
global T
T = (0:4)';
## Since there are discontinuities in the derivative function, you
## should tell lsode where they are for both accuracy and speed.
Tcrit = T(2:4);
y = lsode("ex",x0,T,Tcrit)
jwe

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
