言語ゲーム

とあるエンジニアが嘘ばかり書く日記

Twitter: @propella
トップ > maxima で遊ぶ

maxima で遊ぶ

マックに Grapher というかっこいいグラフ作成ツールが付いているんだけど、どうもかっこいいだけで操作が煩雑なので無料で遊べる maxima を試してみました。マックなら sudo port install maxima で簡単にインストール出来ます。maxima は数式処理ソフトと言って、電卓のように何が何でも一つの数字で答えを出すだけじゃなくて、方程式を解いたり出来るプログラムです。式の終わりに ; を打つと計算です。分数も頑張ってそれらしく表示します。% を先頭に付けると過去の式を参照します。

(%i1) 3 + 4;
(%o1)                                  7
(%i2) % * 6;
(%o2)                                 42
(%i3) %o2 / 56;
                                       3
(%o3)                                  -
                                       4

しかし私のようなマイコン中年は数式よりもプログラムに愛着があるので、プログラムっぽく表示してくれる方法を探しました。

(%i4) display2d : false;
(%o4) false
(%i5) 3 / 4;
(%o5) 3/4

display2d というのは変数名で、true の時に標準の二次元数式モード、false にすると硬派な一次元モードになります。変数の代入にはコロン : を使います。なお、このような(maxima の世界では)当たり前すぎる事をウェブや info で調べるのは結構大変です。maxima 内蔵の help コマンドが大変役に立ちました。? で普通のヘルプ。?? でもっとおおまかなヘルプです。

(%i6) ?? variable;

 0: Functions and Variables for Affine
 1: Functions and Variables for Arrays
 2: Functions and Variables for atensor
... ずらずらメニューが並ぶ

Enter space-separated numbers, `all' or `none': 88

 -- Function: define_variable (<name>, <default_value>, <mode>)
     Introduces a global variable into the Maxima environment.
     `define_variable' is useful in user-written packages, which are
... ずらずらヘルプが表示される。この中のサンプルでコロンが怪しそうだっ
たので、コロンについてヘルプに尋ねる。

(%i7) ? :;

 -- Operator: :
     Assignment operator.

     When the left-hand side is a simple variable (not subscripted),
     `:' evaluates its right-hand side and associates that value with
... ずらずらヘルプが表示される。

グラフ

やっと本題のグラフです。最近読んでる電子回路の問題集をやります。ここに書いてある事は全然間違ってる可能性があるので、違ったら教えてください。1uF のコンデンサの 1KHz の時のリアクタンスを求める問題をやります。こういう公式に代入がある問題をやるときは、ev 関数を使うと分かりやすいと思います。ev(式, 代入文かスイッチ) のように書くと代入した値で式を実行します。

(%i8) Xc: 1 / (2 * %pi * f * C); /* コンデンサのリアクタンスの公式 */
                                       1
(%o8)                              ---------
                                   2 %pi f C

(%i9) ev(Xc,
  C: 1e-6, /* C に 1u を代入 */
  f: 1e3,  /* f に 1K を代入 */
  numer); /* 小数で計算 */

(%o9)                      Rc = 159.1549430918953

コロンで代入した変数は remvalue で削除出来ます。remvalue(all) で全部の変数を削除する事も出来ます。

今度は周波数を色々変えてどうなるかグラフに描きます。グラフを描くときは答えが要らないのでセミコロン ; の代わりにダラー $ を置きます。



(%i10) /* まず単純に描きます */

plot2d (ev(Xc, C: 1e-6), [f, 0, 1000], [y, 0, 1000]) $
(%i11) /* 軸のラベルを付けてみます */

plot2d (ev(Xc, C: 1e-6), [f, 0, 1000], [y, 0, 1000],

  [xlabel, "frequency Hz"],

  [ylabel, "reactance ohm"]) $
(%i12) /* ファイルに出力します */

plot2d (ev(Xc, C: 1e-6), [f, 0, 1000], [y, 0, 1000],

  [xlabel, "frequency Hz"],

  [ylabel, "reactance ohms"],

  [gnuplot_term, "png size 300,200"],

  [gnuplot_out_file, "graph.png"]) $
graph

では次に共振の問題をやります。問題集によると、コイルとコンデンサを直列に繋いだときのリアクタンスはコイルからコンデンサの分を引いた物らしいです。という事は共振周波数はこの合成リアクタンスがゼロの時の周波数を求めれば良いです。

(%i13) Xl: 2 * %pi * f * L;  /* コイルのリアクタンスの公式 */
(%i14) X: Xl - Xc; /* コンデンサとコイルの合成 */
(%i15) solve(X = 0, f);
                               1                    1
(%o15)          [f = - ---------------, f = ---------------]
                        2 %pi sqrt(C L)      2 %pi sqrt(C L)

(%i16) /* 答えは二つ返るので、後ろの答えを得るのに [2], 右辺だけ取り出すのに rhs を使う。*/
Fr: rhs(solve(X = 0, f)[2]);
                                       1
(%o16)                          ---------------
                                2 %pi sqrt(C L)

例えば L: 100mH C: 1uF の時の共振周波数と、ついでに近辺のグラフも見てみます。



(%i17) Fr, L: 100e-3, C:1e-6, numer;

(%o17)                        503.2921210448703
(%i18) plot2d(ev(X, L: 100e-3, C: 1e-6), [f, 0, 4000], [y, -2000, 2000],

  [xlabel, "frequency Hz"],

  [ylabel, "reactance ohms"]) $
graph2

うふふふ。これは面白い。ではさらに直列に 500Ωの抵抗をつなげて 5V の電圧をかけて、周波数に応じた抵抗の周りの電圧を求めます。



(%i19) Z: sqrt(R^2 + X^2);　/* インピーダンスの公式 */

(%i20) Vout: Vin * R / Z; /* 抵抗の周りの電圧 */

(%i21) Vout, Vin: 5, L: 100e-3, C: 1e-6, R: 500;
(%i22) plot2d (ev(Vout, Vin: 5, L: 100e-3, C: 1e-6, R: 500), [f, 0, 2000], [y, 0, 5],

  [xlabel, "frequency Hz"],

  [ylabel, "Vout"]) $
graph3

かっこいい表示

見た目 maxima は mathematica に比べてださいです。マックでは wxMaxima やXmaxima でもっとかっこいい表示に出来ます。でもやってみた結果、結局 emacs で info を読みながら shell mode で使うのが最も楽でした。

参考