Octave の使い方その４
常微分方程式とグラフ

□ 関数の定義と作図

    octave: > function y = f(x)
              y = x.^2 .- 4.*x .+ 6;
              endfunction
    octave: > f(2)
    octave: > f([2, 3])

    octave: > t = linspace(1, 3, 100);
    octave: > grid "on"
    octave: > plot(t, f(t))

    octave: > t = logspace(-1, 2, 100);
    octave: > semilogx(t, f(t))
    octave: > semilogy(t, f(t))
    octave: > loglog(t, f(t))

    octave: > closeplot

□ 常微分方程式の解

   x'' + 4*x' + 13*x = 0, x'(0) = 1, x(0) = 1
   の解を求める．一階の常微分方程式に書き改めると次のようになる．

   x1' = x2
   x2' = -4*x2 - 13*x1

    octave: > function xdot = f(x, t)
              xdot(1) = x(2);
              xdot(2) = -4*x(2) - 13*x(1);
              endfunction
    octave: > x0 = [1; 1]
    octave: > t = linspace(0, 4, 1000);
    octave: > y = lsode("f", x0, t);

    octave: > plot(t, y)
    octave: > grid "on"
    octave: > plot(t, y(:, 1))
    octave: > hold on
    octave: > plot(t, y(:, 2))
    octave: > hold off
    octave: > plot(y(:, 1), y(:, 2))
    octave: > closeplot

   強制振動系
   x'' + 2*x' + (4*π)^2*x = sin(4πt), x'(0) = 0, x(0) = 0
   の解を求める．一階の常微分方程式に書き改めると次のようになる．

   x1' = x2
   x2' = -2*x2 - (4*π)^2*x1 + sin(4πt)

    octave: > function xdot = f(x, t)
              xdot(1) = x(2);
              xdot(2) = -2*x(2) - (4*pi)^2*x(1) + sin(4*pi*t);
              endfunction
    octave: > x0 = [0; 0]
    octave: > t = linspace(0, 4, 1000);
    octave: > y = lsode("f", x0, t);

    octave: > plot(t, y)
    octave: > grid "on"
    octave: > plot(t, y(:, 1))
    octave: > plot(y(:, 1), y(:, 2))
    octave: > closeplot

□ ３次元プロット

   0 ≦ x ≦ 4 and 0 ≦ y ≦ 4 での sin(x)*sin(y) の図

    octave: > x = linspace(0, 4, 40);
    octave: > y = x;
    octave: > [X, Y] = meshgrid(x, y);
    octave: > Z = sin(X).*cos(Y);
    octave: > mesh(x, y, Z)
    octave: > closeplot

   -π ≦ x ≦ π and -π/2 ≦ y ≦ π/2 での sin(x*y) の図

    octave: > x = linspace(-pi, pi, 40);
    octave: > y = linspace(-pi/2, pi/2, 40);
    octave: > [X, Y] = meshgrid(x, y);
    octave: > Z = sin(X.*Y);
    octave: > mesh(x, y, Z)
    octave: > closeplot