Octave でトーラスを描いてみた

Octave でトーラスを描いてみた。

Julia に翻訳--219　トーラス，3次元プロット 参照

Julia の range(), repeat() に相応する関数がないので，自前の関数を書く。

function [rng] = range_(begin_, end_, n)
    step = (end_ - begin_) / (n - 1);
    rng = begin_:step:end_;
end

function [res] = repeat(x, n)
    res = x;
    for i = 2:n
        res = [res, x];
    end
end

function torus(n=100, a=2, b=9)
    u = reshape(repeat(range_(0, 2*pi, length=n), n), n, n);
    v = u';
    x = (b + a * cos(u)) * cos(v);
    y = (b + a * cos(u)) * sin(v);
    z = a * sin(u);
    colormap(pink)
    surf(x, y, z)
    pbaspect([4, 4, 1])
end

torus(50, 5, 10)

matgeom パッケージにトーラスを描く関数があった。

標準では緑色なのだが，赤にしてもやはり浮き輪にしか見えない。赤道の位置（?）に明るい線があるからかな。

pkg load matgeom

figure;
drawTorus([50 50 50 30 10 30 45], "facecolor", "red");
axis equal; view([95 10]); light;

「三次元リング」という説明もあったので，金色にしてみたが，安っぽい色だ。真鍮の指輪みたいだな。

pkg load matgeom

figure;
drawTorus([0 0 0 30 10 0 0], "facecolor", [0xff 0xd3 0x26]/255, "nPhi", 240, "nTheta", 120);
axis equal; view([60 30]); light;

Octave で陰関数のプロット

陰関数で表された曲線の描画
https://blog.goo.ne.jp/r-de-r/e/ff2399cb40043afca226208e8fefaec9

の例を Octave で描いてみた。

clf;
hold on;
for c = [-50, -30, -20, -15, -10, 0, 1, 2, 5, 10, 20, 30, 50]
    ezplot(@(x, y) (x.^2 + y.^2 - 2 .* x).^3- 4*(x .* y).^2 - c, [-1.5, 3.5, -2.5, 2.5], n=1000)
end

clf;
hold on;
for c = -10:10
    ezplot(@(x, y) x.^2 - x .* y + y.^2 - c, [-4, 4], n=1000)
end

Jupyterlab で GNU Octave を使う @macOS

Jupyterlab 側で用意することは唯一，octave kernel のインストール

$ pip install octave_kernel

macOS 用の GNU Octave は

GNU Octave にも手を出してみる
https://blog.goo.ne.jp/r-de-r/e/a59385acd2252786199f2c45663ffaeb

に書いたように，

https://github.com/octave-app/octave-app/releases

で Octave 6.2.0 がインストールされる。

しかし，アプリケーションアイコンは Octave 自身ではなく Octave の IDE 環境（?）のアプレット applet を指している。

Octave 自身は，

/Applications/Octave-6.2.0.app/Contents/Resources/usr/Cellar/octave-octave-app@6.2.0/6.2.0/bin/octave-6.2.0

にある。これに対して /usr/local/bin の中にシンボリックリンク octave を張る。

$ ln -s /Applications/Octave-6.2.0.app/Contents/Resources/usr/Cellar/octave-octave-app@6.2.0/6.2.0/bin/octave-6.2.0 /usr/local/bin/octave

/usr/local/bin にはパスが通っているはず（通っていなければ通す）なので，コマンドライン（プロンプト）から

$ octave

と入力すれば，octave の REPL ウインドウが開く。

GNU Octave, version 6.2.0
Copyright (C) 2021 The Octave Project Developers.
This is free software; see the source code for copying conditions.
There is ABSOLUTELY NO WARRANTY; not even for MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE.  For details, type 'warranty'.

Octave was configured for "x86_64-apple-darwin18.7.0".

Additional information about Octave is available at https://www.octave.org.

Please contribute if you find this software useful.
For more information, visit https://www.octave.org/get-involved.html

Read https://www.octave.org/bugs.html to learn how to submit bug reports.
For information about changes from previous versions, type 'news'.

octave:1> 1+sqrt(8)
ans = 3.8284
octave:2> quit

これが終わったら Jupyterlab のカーネル Octave を選ぶだけだ。

 

GNU Octave: 二次元配列の要素の総和

>> x = [1 2 3 4 5; 6 7 8 9 10]
x =

    1    2    3    4    5
    6    7    8    9   10

のような場合，総和を求めるのに sum(sum(x)) としたのだが sum(x(:)) とも書けるそうだ。x(:) は行優先で一次元配列（ベクトル）に展開してくれる。

>> x(:)
ans =

    1
    6
    2
    7
    3
    8
    4
    9
    5
   10

>> sum(x(:))
ans = 55

>> sum(sum(x))
ans = 55

 

GNU Octave にも手を出してみる

MATLAB は GNU Octave と部分的な互換性があるようだ。

そこで，MATLAB プログラムを動かすときには GNU Octave をインストールして使えばよいようだ。

https://octave.org/download.html

Windows では 7.1.0 がダウンロードできるようだ。
https://octave.org/download.html#ms-windows

macOS ではちょっと面倒なようだ。
https://github.com/octave-app/octave-app/releases
6.2.0 で我慢する。

====

GNU Octave の実行は他のアプリケーションと同じように，アイコンをクリックするだけ。

IDE 環境ウインドウが開く。コマンドウインドウになにか入力すれば，それが実行され，結果が表示される。

必要なパッケージの準備

例えば，以下のプログラムでは 「chi2cdf は statistics パッケージにあるよ！」 と言われる。

コマンドラインで statistics パッケージの最新バージョンを調べて（たとえば statistics-1.4.3.tar.gz のとき） wget でダウンロードする。

% wget http://downloads.sourceforge.net/project/octave/Octave%20Forge%20Packages/Individual%20Package%20Releases/statistics-1.4.3.tar.gz

パッケージを octave コマンドウインドウでインストールする（1回だけ）

>> pkg install -global statistics-1.4.3.tar.gz

パッケージを使用するとき（セッションごと）

>> pkg load statistics

使用するプログラム。

>> function [chisq, df, p] = ChiSquareTest(x)
    [nrow , ncol] = size(x);
    df = (nrow - 1) * (ncol - 1);
    colsums = sum(x, 1);
    rowsums = sum(x, 2);
    n = sum(rowsums);
    expectations = (rowsums * colsums) / n;
    Zscores = (x .- expectations) .^ 2 ./ expectations;
    chisq = sum(sum(Zscores));
    p = 1 - chi2cdf(chisq, df);
end

>> x = [55 22 16 7; 40 32 24 4];
[chisq, df, p] = ChiSquareTest(x)
chisq = 6.6385
df = 3
p = 0.084359

>> [chisq, df, p] = ChiSquareTest([10 5; 6 11])
chisq = 3.1373
df = 1
p = 0.076523

* プログラム中に日本語があると変なことが起きることがあるようだ（プログラムをコピペできない）。
* MATLAB のすべてが octave で利用できるというものではない。例えば，icdf('chi2', x, df) はないので，1 - chi2cdf(x, df) とするなど，エラー回避が必要になることもある。