調和関数表示用RC48Mのソースの紹介画面

// 2007年2月2日に作成する
// Excelで学ぶ物理数学(ナツメ社) 245ページを参考にしました。
// 中田　亨様のホームページを参考にしました。
// http://staff.aist.go.jp/toru-nakata/index-j.html

//　四元数を用いて、三次元空間で、回転操作を行って球面調和関数を描く

#　　t=δ, a=α, b = β, c = γ, s=θ, p=ψ

i@;　// 単位変数であることを宣言する
j@;
k@;

//　回転軸(α, β, γ)
α=0;
β=1;
γ=0;

δ = 0.4; // 回転角度(ラジアン)

// 座標(x, y, z)

r=1; // 半径
θ0 = 0;
θ16;
dθ=PI()/16;

ψ0 = 0.0;
ψ40;
dψ = 2PI()/40;

IF(1) { // Y10 Y20 では必要、極点近くでは細かく刻む
　　ALL(θ) {
　　　　θ(θ) = (1-COS(θ(θ)))/2*PI();
　　};
};

X3;

QUATERNION() { // 四元数の演算を宣言する
　　Q(i) = QUARTER(COS(δ/2), αSIN(δ/2), βSIN(δ/2), γSIN(δ/2));
};

ALL(ψ,θ) {
　　r = ABS(3COS(θ(θ))^2-1); 　　　// Y20
//　　r = ABS(COS(θ(θ)));　　　　　// Y10
//　　r = ABS(SIN(θ(θ)));　　　　　// Y11
//　　r = ABS(SIN(θ(θ))COS(θ(θ))); // Y21
　　x=rSIN(θ(θ))COS(ψ(ψ));
　　y=rSIN(θ(θ))SIN(ψ(ψ));
　　z=rCOS(θ(θ));

// 回転角度(δ)反時計回り

　　QUATERNION() { // 四元数の演算を宣言する
　　　　P(i) = QUARTER(0,x,y,z); // 四元数を生成する
　　　　A(i) = ~Q(i)P(i)Q(i); // δ/2を２回だけ回転させる
　　};
　　ALL(X) {
　　　　Z(ψ,θ,X:0) = ELEMENT(A(i), 0); // 要素を取り出す
　　　　Z(ψ,θ,X:1) = ELEMENT(A(i), 1);
　　　　Z(ψ,θ,X:2) = ELEMENT(A(i), 2);
　　　　Z(ψ,θ,X:3) = ELEMENT(A(i), 3);
　　};

};
MPRT(Z(ψ,θ,X)); // csvファイルに出力する

CSVファイルのExcelでの図表現の方法