Shadow Academy トップ > 我流・物理数学 > ファン・デル・ポール方程式のシミュレーション


ファン・デル・ポール方程式のシミュレーション

繰り返しの回数:n (n×Δt 秒後に終了)
時間の刻み幅:Δt [s] (n×Δt 秒後に終了)
ε ex). ε=0.1, 0.3, 1, 10, …etc.
ε=0で単振動(LC回路)と同じ
ω ex). ω2=1
初期条件:x0 x0p0は、少なくとも一方は、0以外。
初期条件:p0 x0p0は、少なくとも一方は、0以外。
データ出力形式:

データ出力欄:



解説

ファン・デル・ポール方程式は、二階微分方程式:

d2x   dx  

ε(1-x2)
ω2x=0
dt2   dt  
で表されるが、
dx  

p
dt  
と置くと、
dp  

ε(1-x2)pω2x
dt  
となるので、両者の連立方程式の形で表すことが出来る。
従って、ファン・デル・ポール方程式を
数値計算によってシミュレーションする場合、
上述の式を4次のルンゲ・クッタ法(RK4)で解けばよい。




スクリプトプログラムの使い方

各テキストフィールドに数値を入力後、「計算する」ボタンを押すと、
新しいページテキストエリアにデータが生成されるので、そこで、「全て選択」して、
テキストエディタを起動して、ファイルにコピーペーストし、名前を付けて保存する。
ここでは、Cドライブ以下のtempフォルダに「vanderpol.dat」という名前で保存する。
これをグラフ作成ツール「gnuplot」を用いて、グラフ化するには、
gnuplotを起動して、以下のコマンドを入力する。

横軸をt、縦軸をxとして出力

set terminal png
set output 'c:\temp\vanderpolxtgraph.png'
plot 'c:\temp\vanderpol.dat' using 1:2 with lines
exit

横軸をt、縦軸をpとして出力

set terminal png
set output 'c:\temp\vanderpolptgraph.png'
plot 'c:\temp\vanderpol.dat' using 1:3 with lines
exit

横軸をx、縦軸をpとして出力(位相空間)

set terminal png
set output 'c:\temp\vanderpolphasespace.png'
plot 'c:\temp\vanderpol.dat' using 2:3 with lines
exit




ソースコード(新版)

vanderpol.html

<form><!-- getElementByIdを使用する為、name属性を指定しない -->
<table border="1" style="margin-left:auto;margin-right:auto;text-align:center;">
  <tr>
    <td style="text-align:right;">繰り返しの回数:<i>n</i>=</td>
    <td style="text-align:left;"><input type="text" id="n" value="1000"></td>
    <td style="text-align:left;">(<i>n</i>&times;<i>&Delta;t</i> 秒後に終了)</td>
  </tr>
  <tr>
    <td style="text-align:right;">時間の刻み幅:<i>&Delta;t</i>=</td>
    <td style="text-align:left;"><input type="text" id="dt" value="0.01">[s]</td>
    <td style="text-align:left;">(<i>n</i>&times;<i>&Delta;t</i> 秒後に終了)</td>
  </tr>
  <tr>
    <td style="text-align:right;"><i>&epsilon;</i>=</td>
    <td style="text-align:left;"><input type="text" id="epsilon" value="0.3"></td>
    <td style="text-align:left;">ex). <i>&epsilon;</i>=0.1, 0.3, 1, 10, …etc.<br><i>&epsilon;</i>=0で単振動(LC回路)と同じ</td>
  </tr>
  <tr>
    <td style="text-align:right;"><i>&omega;</i>=</td>
    <td style="text-align:left;"><input type="text" id="omega" value="1"></td>
    <td style="text-align:left;">ex). <i>&omega;</i><sup>2</sup>=1</td>
  </tr>
  <tr>
    <td style="text-align:right;">初期条件:<i>x</i><sub>0</sub>=</td>
    <td style="text-align:left;"><input type="text" id="x0" value="1"></td>
    <td style="text-align:left;"><i>x</i><sub>0</sub>と<i>p</i><sub>0</sub>は、少なくとも一方は、0以外。</td>
  </tr>
  <tr>
    <td style="text-align:right;">初期条件:<i>p</i><sub>0</sub>=</td>
    <td style="text-align:left;"><input type="text" id="p0" value="0"></td>
    <td style="text-align:left;"><i>x</i><sub>0</sub>と<i>p</i><sub>0</sub>は、少なくとも一方は、0以外。</td>
  </tr>
  <tr>
    <td style="text-align:right;">データ出力形式:</td>
    <td style="text-align:left;" colspan="2">
      <label for="csv"><input type="radio" name="separator" id="csv">CSV(半角コンマ区切り)</label><br>
      <label for="ssv"><input type="radio" name="separator" id="ssv" checked>SSV(半角スペース区切り)</label><br>
      <input type="button" value="計算する" onClick="vdpcalc()" style="width:200px;height:40px;"></td>
  </tr>
  <tr>
    <td style="text-align:right;">データ出力欄:</td>
    <td style="text-align:left;" colspan="2"><textarea id="data" cols="50" rows="20"></textarea></td>
  </tr>
</table>
</form>

vanderpol.js

var vdpcalc = function(){
  var n = parseInt(document.getElementById("n").value); // 繰り返しの回数(n*dt秒後に終了)
  var dt = parseFloat(document.getElementById("dt").value); // 時間の刻み幅(n*dt秒後に終了)
  var epsilon = parseFloat(document.getElementById("epsilon").value);
  var omega = parseFloat(document.getElementById("omega").value);
  var x0 = parseFloat(document.getElementById("x0").value); // xの初期条件
  var p0 = parseFloat(document.getElementById("p0").value); // pの初期条件

  var t,x,p,i; // 時刻t, x, p及び、ループ変数を宣言
  var j1,j2,j3,j4,k1,k2,k3,k4; // 4次のルンゲクッタ法で使用する変数を宣言

  var data = ""; // データを表す変数を空の初期値と共に宣言
  var separator = ""; // データの区切り文字

  // データの区切り文字の種類を選択して取得
  if(document.getElementById("csv").checked == true){ separator=","; } // 半角コンマ区切り
  if(document.getElementById("ssv").checked == true){ separator=" "; } // 半角スペース区切り

  t = 0; // 時刻を初期化
  x = x0; // xにその初期条件を代入
  p = p0; // xにその初期条件を代入

  data = t + separator + Math.round(1000*x)/1000 + separator + Math.round(1000*p)/1000;
  // 時刻t=0での、xとp(小数点以下第3位未満を四捨五入)
  // をデータを表す変数dataに代入して、データを上書きする。

  for(i=0;i<n;++i){
    j1=dt*p;
    k1=dt*(epsilon*(1-Math.pow(x,2))*p-Math.pow(omega,2)*x);
    j2=dt*(p+k1/2);
    k2=dt*(epsilon*(1-Math.pow((x+j1/2),2))*(p+k1/2)-Math.pow(omega,2)*(x+j1/2));
    j3=dt*(p+k2/2);
    k3=dt*(epsilon*(1-Math.pow((x+j2/2),2))*(p+k2/2)-Math.pow(omega,2)*(x+j2/2));
    j4=dt*(p+k3);
    k4=dt*(epsilon*(1-Math.pow((x+j3),2))*(p+k3)-Math.pow(omega,2)*(x+j3));

    t=Math.round(1000*(t+dt))/1000;
    x=x+(j1+2*j2+2*j3+j4)/6;
    p=p+(k1+2*k2+2*k3+k4)/6;
    data += "\n" + t + separator + Math.round(1000*x)/1000 + separator + Math.round(1000*p)/1000;
    // 時刻tでの、xとp(小数点以下第3位未満を四捨五入)
    // をデータを表す変数dataに代入して、データを上書きする。
  }

  document.getElementById("data").value = data; // テキストエリアへ出力する
};



ソースコード(旧版)

vanderpol.html

<form name="vdpform">
<table border="5" style="margin-left:auto;margin-right:auto;text-align:center;">
<tr>
<td style="text-align:right;">繰り返しの回数:<i>n</i>=</td>
<td style="text-align:left;"><input type="text" name="n"></td>
<td style="text-align:left;">(<i>n</i>×<i>Δt</i> 秒後に終了)</td>
</tr>
<tr>
<td style="text-align:right;">時間の刻み幅:<i>Δt</i>=</td>
<td style="text-align:left;"><input type="text" name="dt">[s]</td>
<td style="text-align:left;">(<i>n</i>×<i>Δt</i> 秒後に終了)</td>
</tr>
<tr>
<td style="text-align:right;"><i>ε</i>=</td>
<td style="text-align:left;"><input type="text" name="epsilon"></td>
<td style="text-align:left;">ex). <i>ε</i>=0で単振動(LC回路)と同じ</td>
</tr>
<tr>
<td style="text-align:right;"><i>ω</i>=</td>
<td style="text-align:left;"><input type="text" name="omega"></td>
<td style="text-align:left;">ex). <i>ω</i><sup>2</sup>=1</td>
</tr>
<tr>
<td style="text-align:right;">初期条件:<i>x</i><sub>0</sub>=</td>
<td style="text-align:left;"><input type="text" name="x0"></td>
<td style="text-align:left;">
<i>x</i><sub>0</sub>と<i>p</i><sub>0</sub>は、少なくとも一方は、0以外。</td>
</tr>
<tr>
<td style="text-align:right;">初期条件:<i>p</i><sub>0</sub>=</td>
<td style="text-align:left;"><input type="text" name="p0"></td>
<td style="text-align:left;">
<i>x</i><sub>0</sub>と<i>p</i><sub>0</sub>は、少なくとも一方は、0以外。</td>
</tr>
</table>
<input type="button" value="計算する" onClick="vdpcalc()">
</form>

vanderpol.js

function vdpcalc(){
  var n=parseInt(document.vdpform.n.value);
  /* 繰り返しの回数(n*dt秒後に終了) */
  var dt=parseFloat(document.vdpform.dt.value);
  /* 時間の刻み幅(n*dt秒後に終了) */
  var epsilon=parseFloat(document.vdpform.epsilon.value);

  var omega=parseFloat(document.vdpform.omega.value);

  var x0=parseFloat(document.vdpform.x0.value);
  /* xの初期条件 */
  var p0=parseFloat(document.vdpform.p0.value);
  /* pの初期条件 */
  var t,x,p,i;
  /* 時刻t, x, p及び、ループ変数を宣言 */
  var j1,j2,j3,j4,k1,k2,k3,k4;
  /* 4次のルンゲクッタ法で使用する変数を宣言 */

  t=0; /* 時刻を初期化 */
  x=x0; /* xにその初期条件を代入 */
  p=p0; /* xにその初期条件を代入 */

  document.writeln("<html>");
  document.writeln("<head>");
  document.writeln("<title>ファン・デル・ポール方程式</title>");
  document.writeln("</head>");
  document.writeln("<body bgcolor='#ffd700' text='#000000'>");

  document.writeln("<pre>");
  document.writeln(t," ",Math.round(1000*x)/1000," ",Math.round(1000*p)/1000);
  /* 時刻t=0での、xとpを出力(小数点以下第3位未満を四捨五入) */
  for(i=0;i<n;++i){
    j1=dt*p;
    k1=dt*(epsilon*(1-Math.pow(x,2))*p-Math.pow(omega,2)*x);
    j2=dt*(p+k1/2);
    k2=dt*(epsilon*(1-Math.pow((x+j1/2),2))*(p+k1/2)-Math.pow(omega,2)*(x+j1/2));
    j3=dt*(p+k2/2);
    k3=dt*(epsilon*(1-Math.pow((x+j2/2),2))*(p+k2/2)-Math.pow(omega,2)*(x+j2/2));
    j4=dt*(p+k3);
    k4=dt*(epsilon*(1-Math.pow((x+j3),2))*(p+k3)-Math.pow(omega,2)*(x+j3));

    t=Math.round(1000*(t+dt))/1000;
    x=x+(j1+2*j2+2*j3+j4)/6;
    p=p+(k1+2*k2+2*k3+k4)/6;
    document.writeln(t," ",Math.round(1000*x)/1000," ",Math.round(1000*p)/1000);
    /* 時刻tでの、xとpを出力(小数点以下第3位未満を四捨五入) */
  }
  document.writeln("</pre>");

  document.writeln("</body>");
  document.write("</html>");
}



参考文献

  1. 「非線形物理学―カオス・ソリトン・パターン―」(裳華房、2010年)



Shadow Academy トップへ戻る

inserted by FC2 system