liboctave Wiki

連立微分方程式を解きたい

Tips

解説 †

連立微分方程式 \frac{d\vec{x}}{dt}=\vec{f}(\vec{x}) を与えられた初期値のもとで解きます。

サンプルプログラム †

\dot{x}_0=x_0+x_1
\dot{x}_1=-4~x_0-3~x_1

を初期値x_0(t=0)=0,x_1(t=0)=1として解きます。

filelsode.cpp

#include <iostream>
#include <octave/config.h>
#include <octave/Matrix.h>
#include <octave/LSODE.h> 

using namespace std; 

ColumnVector xdot(const ColumnVector& x, double t)
{
  ColumnVector dx(2);
  dx(0) =    x(0) + 2*x(1);
  dx(1) = -4*x(0) - 3*x(1);
  return dx;
}

int main()
{
  ColumnVector xi(2), x(2);
  ODEFunc odef(xdot); // set xdot  
 
  xi(0)=0.0;  // initial condition
  xi(1)=1.0;  //
  LSODE ls(xi,0.0,odef); // solver
  
  for(double t = 0; t <= 3.0; t+=0.1) {
    x=ls.do_integrate(t);
    cout << t << " " << x(0) << " " << x(1) << endl;
  }  
}

出力 †

0 0 1
0.1 0.179763 0.707037
0.2 0.318829 0.435272
0.3 0.418297 0.193126
0.4 0.480858 -0.0138417
0.5 0.510378 -0.182668
0.6 0.511514 -0.312648
(・・・中略・・・)
2.7 -0.0519342 0.0945891
2.8 -0.0383875 0.0855498
2.9 -0.0255641 0.0742883

解はx_0(t)=e^{-t}{\rm{sin}}(2t), x_0(t)=e^{-t}({\rm~cos}(2t)-{\rm~sin}(2t))となります。

plot.gif

コメント †

微分方程式に複素数が入った場合 †

管理人? (2005-01-05 (水) 16:42:53)

単純にColumnVectorをComplexColumnVectorにしただけでは駄目らしい。
複素数も扱えるソルバーは何？



コンパイルについて †

管理人? (2005-01-05 (水) 16:39:47)

私の環境(Windows,Cygwinにoctave-2.1.64をインストール)では、コンパイル時に

/usr/local/include/octave-2.1.64/octave/LSODE.h:34:24:\
LSODE-opts.h: No such file or directory

と出てコンパイルできませんでした。コンパイルのオプションとして

-I/usr/local/include/octave-2.1.64/octave

と付け加えることでコンパイルできました。

  • mkoctfile のほうでは特に問題なくコンパイルできました。(ld関連と思われるInfo:は出ましたが) -- 管理人? 2005-01-12 (水) 22:32:11

&ref(): File not found: "wikilogo1.png" at page "連立微分方程式を解きたい";

[edit]