レールガンコイルガン

コイルガン関連の話で、磁気双極子が電磁コイルに引き付けられる力の計算を少しやってみました。
(ゆっくり動かして、渦電流などは無視しています。)
コイル内半径は a1、外半径は 2*a1 と決め付けました。コイル長は 2*b1 で、中心は原点 (0,0) です。

電流素片が座標位置 (b, a) にあって、それに引きつけられる磁気双極子の位置を
P(x, 0) とします。 P は常に x 軸上を動きます。(電流素片は x 軸を取り巻く円周上にあり、
それらを総計しなくてはなりません。この辺の説明は飛ばしています。)
力のベクトルの方向は (b-x, a)
その x 成分を１に正規化すると (b-x)/sqrt((b-x)^2 + a^2)

コイル電流を I、巻き数を N としたとき (k は適当な係数)、力の大きさは
k*I*N / d^2, d = sqrt((b-x)^2 + a^2)

よって磁気双極子に加わる加速力は (b-x)/sqrt((b-x)^2 + a^2)^3 に比例します。
このへんから数式処理プログラム (MAXIMA) の力を借りました。

f1(x,a,b) := (b-x)/sqrt((b-x)^2 + a^2)^3;
fortran(defint(f1(x,a,b),b,-b1,b1));
f2(x,a) := 1/SQRT(x**2+2*B1*x+B1**2+a**2)-1/SQRT(x**2-2*B1*x+B1**2+a**2);
fortran(defint(f2(x,a),a,a1,2*a1));
(b1 > 0 かつ a1 > 0 です)

f3(x,a,b) := asinh(2*a/sqrt(x**2+2*b*x+b**2))-asinh(a/sqrt(x**2+2*b*x+b**2))
-asinh(2*a/sqrt(x**2-2*b*x+b**2))+asinh(a/sqrt(x**2-2*b*x+b**2));
plot2d(f3(x,0.5,1),[x,-3,3]);

得られたグラフは、次のようになりました。
ttp://radio.s56.xrea.com/radio/src/radio0630.png
あまり自信が持てません。よろしく御批判のほど、お願いします。
(これに間違いがなければ、コイル径・巻き幅を変えても加速力のより正確な見積もりができます。)