Maxima を使って，関数のグラフを描くことができます。また，数値データもグラフにすることができます。
Maxima の2次元グラフ作成は draw2d() でもできますが，ここでは，plot2d() について説明します。
plot2d() でプロットできるオブジェクトは以下のとおりです。
- 陽関数 $y = f(x)$:
plot2d(f(x), [x, xmin, xmax])
- 陰関数 $f(x, y) = 0$:
plot2d(f(x, y) = 0, [x, xmin, xmax], [y, ymin, ymax])
- 媒介変数表示 $x(t), y(t)$:
plot2d([parametric, x(t), y(t), [t, t0, t1]])
- 点，$x$ 座標 $y$ 座標の数値データ:
plot2d([discrete, [[x1, y1], ..., [xn, ..., yn]]])またはplot2d([discrete, [x1, ..., xn], [y1, ..., yn]])
陽関数のグラフ
例として，$y = \sin x$ のグラフを $−2\pi \leq x \leq 2\pi$ の範囲で描きます。基本的な定数の一つである円周率 $\pi$ は Maxima では %pi と書きます。
plot2d(sin(x), [x, -2*%pi, 2*%pi])$
いくつかオプションを設定する例。
/* オプションの設定例 */
plot2d(sin(x), [x, -2*%pi, 2*%pi],
/* 表示範囲 */
[x, -10, 10], [y, -1.1, 1.1],
/* グリッド */
grid2d,
/* 凡例 */
[legend, "正弦関数"],
/* 線の太さと色 */
[style, [lines, 2, red]]
)$
陰関数のグラフ
例として，$x^2 + y^2 = 1$ のグラフを描きます。$y$ について解いて陽関数の形にしたり，以下で述べるような媒介変数表示にしなくても，陰関数のまま，グラフに描くことができます。
陰関数 $f(x, y) = 0$ を [x, xmin, xmax], [y, ymin, ymax] の範囲で描く書式は…
plot2d(f(x, y) = 0, [x, xmin, xmax], [y, ymin, ymax])$
plot2d(x**2 + y**2 - 1 = 0, [x, -1.2, 1.2], [y, -1.2, 1.2])$
縦横比・サイズの設定
縦横の比 yx_ratio を設定して円らしく見えるようにします。size も変更します。弘大 JupyterHub ではグラフを svg でインライン表示させているので， gnuplot_svg_term_command で size 600, 600 などと設定します。
/* オプションの設定例 */
plot2d(x**2 + y**2 = 1, [x, -1.2, 1.2], [y, -1.2, 1.2],
/* 縦横比 [yx_ratio, 1] または... */
[same_xy],
[gnuplot_svg_term_command, "set term svg font \",14\" size 600,600"]
)$
媒介変数表示のグラフ
半径 $1$ の円の方程式は $x^2 +y^2 = 1$ です。この円を以下のような媒介変数表示にして描きます。
$$ x=\cos t, \quad y=\sin t \quad (0 \leq t \leq 2\pi) $$
書式は…
plot2d([parametric, x(t), y(t), [t, tmin, tmax]])$
plot2d([parametric, cos(t), sin(t), [t, 0, 2*%pi]]
)$
/* オプションの設定例 */
plot2d([parametric, cos(t), sin(t), [t, 0, 2*%pi]],
/* 縦横比 [same_xy] は縦横の表示範囲によらず円は円に見える */
[same_xy],
/* 表示範囲 */
[x, -1.5, 1.5], [y, -1.1, 1.1],
/* 横軸・縦軸のラベル */
[xlabel, "x"], [ylabel, "y"]
)$
点・数値データのグラフ
以下の例では，最初（左側）の数値を $x$ 座標の値，次（右側）の数値を $y$ 座標の値として 6 個の点 からなるリスト data の各点をつないだ折れ線グラフを描きます。
data: makelist([i, i^2], i, 0, 5);
plot2d([discrete, data])$
以下のように，$x$ 座標のみの数値データ，$y$ 座標のみの数値データをプロットしてもよいです。
X: makelist(i, i, 0, 5);
Y: makelist(i**2, i, 0, 5);
オプションを設定して，各点を線でつながずにプロットする例。
ここでは，[style, points]を設定して線でつながずに点で描きます。プロットオプションの style には，points の他にも，lines，linespoints， dots が設定できます。
/* オプションの設定例 */
plot2d([discrete, X, Y],
/* 表示範囲 */
[x, 0, 5], [y, 0, 28],
/* 点と大きさ */
[style, [points, 2]],
/* 色 blue, red, green, magenta, black, cyan から */
[color, red]
)$
複数のグラフを重ねて表示
Maxima で複数のグラフを重ねて表示する例を示します。
$y = x^2 – 1$ と $y = 4 x – 5$ の 2 つのグラフを重ねて描きます。
$x$ の範囲は $-5 \leq x \leq 5$ です。
以下の例でわかるように，重ねて描きたい関数を [x^2 - 1, 4*x - 5] のように [ と ] で囲んだリスト形式にして plot2d コマンドを使います。
plot2d([x^2 - 1, 4*x - 5], [x, -5, 5])$
system("cp ~/.maxplot.svg ./maxplot09.svg")$
いくつかオプションを設定して描く例です。
plot2d(
[x^2 - 1, 4*x - 5], [x, -5, 5],
/* 表示範囲 */
[y, -5, 10],
/* 線の色の設定 */
[color, red, black],
/* 線の太さの設定。最初の線の太さを 3 に */
[style, [lines, 3], lines],
/* 凡例の位置。グラフの線とかぶらないように */
[gnuplot_preamble, "set key bottom"]
)$
数値データファイルを読み込む
Maxima では，あらかじめ作成された数値データファイルを読み込んでグラフを描くこともできます。
/* 読み込むデータファイルを作成しておく */
tmp: makelist([i, i^2], i, 0, 5)$
write_data(tmp, "data.txt")$
データファイルを読み込む Maxima の関数は，read_nested_list を使います。
dat: read_nested_list("data.txt");
リスト dat に読み込んだ数値データを plot2d() でグラフにします。
plot2d([discrete, dat],
[x, 0, 6], [y, 0, 28], /* x, y の表示範囲を設定 */
[style, [points, 1.5]], [color, blue], /* 大きさ，色 */
[gnuplot_preamble, "set key top left"], /* 凡例の位置 */
[legend, "データ"] /* 凡例表示の設定 */
)$
数値データと理論曲線を重ねて表示
前節の数値データを読み込んだリスト dat と理論曲線 $y = x^2$ の 2 つのグラフを重ねて表示してみます。
/* もう少しオプションを指定して描く */
plot2d([[discrete, dat], x**2], [x, 0, 5],
[style, [points, 1.5], lines],
[color, blue, red],
[gnuplot_preamble, "set key top left"],
[legend, "データ", "理論曲線 y=x^2"]
)$
海王星と冥王星の軌道
焦点を原点とした楕円のグラフ
太陽からの万有引力を受けて運動する惑星は，太陽（二体問題では質量中心）を焦点の一つとした楕円運動を描いて運動します。
焦点の一つを原点とし，長半径 $a$，離心率 $e$ の楕円の方程式は，極座標 $r$，$\varphi$ を使って以下のように表すことができます。
$$ r= \frac{a (1−e^2)}{ 1 + e\cos\varphi}$$
さて，かつては第9惑星，現在では準惑星の一つである冥王星も楕円軌道を描きます。冥王星の軌道長半径 $𝑎_P=39.445 \mbox{au}$，離心率 $𝑒_P=0.250$
を使って冥王星の軌道を描きます。
まず，極座標表示の楕円の式を関数として定義します。
r(a, e, phi):= a*(1-e^2)/(1+e*cos(phi));
x(a, e, phi):= r(a, e, phi)*cos(phi);
y(a, e, phi):= r(a, e, phi)*sin(phi);
媒介変数表示でプロット
冥王星の軌道長半径と離心率を入れて，媒介変数表示でプロットしてみます。
aP: 39.445$
eP: 0.250$
plot2d([parametric, x(aP, eP, phi), y(aP, eP, phi), [phi, 0, 2*%pi]],
[same_xy],
[xlabel, "x"], [ylabel, "y"], [legend, "冥王星"],
[x,-50, 50], [y, -50, 50]
)$
では，この冥王星の軌道と，その内側をまわる海王星の軌道を重ねて描いてみましょう。
海王星の軌道長半径は $a_N=30.181 \mbox{au}$，離心率は $e_N=0.0097$ と小さいので簡単のために $e_N=0$ として扱います。
実際の2つの天体の軌道は同じ平面上にありませんが，太陽からの距離のみをみるために，ここでは同一平面上に描きます。
aN: 30.181$
eN: 0$
plot2d([[parametric, x(aP, eP, phi), y(aP, eP, phi), [phi, 0, 2*%pi]],
[parametric, x(aN, eN, phi), y(aN, eN, phi), [phi, 0, 2*%pi]]],
[same_xy],
[xlabel, ""], [ylabel, ""],
[legend, "冥王星", "海王星"],
[x,-50, 50], [y, -50, 50]
)$
月別平年気温の関数フィット
弘前市の月別平年気温のデータを使ってリスト HiroT を作ります。
HiroT: [
[1, -1.5],
[2, -1.0],
[3, 2.3],
[4, 8.6],
[5, 14.3],
[6, 18.3],
[7, 22.3],
[8, 23.5],
[9, 19.4],
[10, 12.9],
[11, 6.5],
[12, 0.8]
]$
このデータをグラフに描いてみます。
plot2d([discrete, HiroT], [style, [points, 1]])$
グラフを見ると，月別平年気温は 12 ヶ月周期の三角関数でフィットできるようにみえます。
では，以下のように関数フィットをしてみましょう。
まず，最小二乗法のパッケージ lsquares をロードします.これは lsquares_estimates を 初めて使う際，最初に 1 回だけ行えばよいです。
load(lsquares)$
次に，次のような関数を定義して，最小二乗法でパラメータ a, b, c を求めます。
$$f(x) = a + b \sin\left(\frac{2\pi}{12}\right)+ c \cos\left(\frac{2\pi}{12}\right)$$
f(x):= a + b*sin(%pi*x/6)+ c*cos(%pi*x/6);
最小二乗法を行う関数 lsquares_estimates には数値データを行列として入れてやる必要があります。
このために apply('matrix, HiroT) でリスト HiroT を行列に変換しています。
/* 表示は 5 桁くらいで... */
fpprintprec: 5$
lsquares_estimates(apply('matrix, HiroT), [x, y], y=f(x), [a, b, c])$
float(%);
最小二乗法によって求めた a, b, c の値を f(x) に代入して fit とし，数値データと合わせて，この関数のグラフも描きます。
fit: ev(f(x), %);
plot2d([[discrete, HiroT], fit], [x, 1, 12],
[style, [points, 1], lines],
/* 表示範囲 */ [x, 0.8, 12.2],
[color, blue, red],
[legend, "月別平年気温", "関数フィット"],
/* x 軸の目盛は 1 ごとに。凡例は center に */
[xtics, 1], [gnuplot_preamble, "set key center"],
[xlabel, "月"], [ylabel, "°C"]
)$
