Динамічні властивості нелінійних локалізованих мод у лінійних молекулярних ланцюжках
//Для молекул рахуємо методом Рунґе-Кутта
//Функції v, w ,y, z стосуються невідомих f, fi, p, u
$tm1 = $t-$dtau;
$mp1 = $m+1;
$mm1 = $m-1;
$v1[$m][$tm1] =
f1($j, $fi[$mp1][$tm1], $fi[$mm1][$tm1], $g0, $u[$m][$tm1], $fi[$m][$tm1]);
$w1[$m][$tm1] =
f2($j, $f[$mp1][$tm1], $f[$mm1][$tm1], $g0,$u[$m][$tm1],$f[$m][$tm1]);
$y1[$m][$tm1] =
f3($g0,$c0, $fi[$m][$tm1], $f[$m][$tm1], $u[$m][$tm1], $hamma, $p[$m][$tm1]);
$z1[$m][$tm1] = f4($p[$m][$tm1]);
$v2[$m][$tm1] = f1($j, $fi[$mp1][$tm1]+$w1[$mp1][$tm1]/2, $fi[$mm1][$tm1]+$w1[$mm1][$tm1]/2, $g0, $u[$m][$tm1]+$z1[$m][$tm1]/2, $fi[$m][$tm1]+$w1[$m][$tm1]/2);
$w2[$m][$tm1] = f2($j, $f[$mp1][$tm1]+$v1[$mp1][$tm1]/2, $f[$mm1][$tm1]+$v1[$mm1][$tm1]/2, $g0, $u[$m][$tm1]+$z1[$m][$tm1]/2, $f[$m][$tm1]+$v1[$m][$tm1]/2);
$y2[$m][$tm1] = f3($g0,$c0, $fi[$m][$tm1]+$w1[$m][$tm1]/2, $f[$m][$tm1]+$v1[$m][$tm1]/2, $u[$m][$tm1]+$z1[$m][$tm1]/2, $hamma, $p[$m][$tm1]+$y1[$m][$tm1]/2);
$z2[$m][$tm1] = f4($p[$m][$tm1]+$y1[$m][$tm1]/2);
$v3[$m][$tm1] = f1($j, $fi[$mp1][$tm1]+$w2[$mp1][$tm1]/2, $fi[$mm1][$tm1]+$w2[$mm1][$tm1]/2, $g0, $u[$m][$tm1]+$z2[$m][$tm1]/2, $fi[$m][$tm1]+$w2[$m][$tm1]/2);
$w3[$m][$tm1] = f2($j, $f[$mp1][$tm1]+$v2[$mp1][$tm1]/2, $f[$mm1][$tm1]+$v2[$mm1][$tm1]/2, $g0, $u[$m][$tm1]+$z2[$m][$tm1]/2, $f[$m][$tm1]+$v2[$m][$tm1]/2);
$y3[$m][$tm1] = f3($g0,$c0, $fi[$m][$tm1]+$w2[$m][$tm1]/2, $f[$m][$tm1]+$v2[$m][$tm1]/2, $u[$m][$tm1]+$z2[$m][$tm1]/2, $hamma, $p[$m][$tm1]+$y2[$m][$tm1]/2);
$z3[$m][$tm1] = f4($p[$m][$tm1]+$y2[$m][$tm1]/2);
$v4[$m][$tm1] = f1($j, $fi[$mp1][$tm1]+$w3[$mp1][$tm1], $fi[$mm1][$tm1]+$w3[$mm1][$tm1], $g0, $u[$m][$tm1]+$z3[$m][$tm1], $fi[$m][$tm1]+$w3[$m][$tm1]);
$w4[$m][$tm1] = f2($j, $f[$mp1][$tm1]+$v3[$mp1][$tm1], $f[$mm1][$tm1]+$v3[$mm1][$tm1], $g0, $u[$m][$tm1]+$z3[$m][$tm1], $f[$m][$tm1]+$v3[$m][$tm1]);
$y4[$m][$tm1] = f3($g0,$c0, $fi[$m][$tm1]+$w3[$m][$tm1], $f[$m][$tm1]+$v3[$m][$tm1], $u[$m][$tm1]+$z3[$m][$tm1], $hamma, $p[$m][$tm1]+$y3[$m][$tm1]);
$z4[$m][$tm1] = f4($p[$m][$tm1]+$y3[$m][$tm1]);
//Маємо результати$f[$m][$t] = $f[$m][$tm1] + ($dtau/6)*($v1[$m][$tm1] + 2*$v2[$m][$tm1] + 2*$v3[$m][$tm1] + $v4[$m][$tm1]);