ポテンシャル関数f(x,y)なら、電位分布は等...

　　もっと高速な方法はあるでしょうか？島村1243 2007/12/13 14:21:55 (修正1回) 　　└続き島村1243 2007/12/13 14:23:05 (修正1回) 　　　├(1)について（処理手順の改良）山中和義 2007/12/14 08:55:07 　　　│├!もったいない、整理すれば、４倍速。 SECOND 2007/12/14 21:34:04 　　　││└確かに速いです。島村1243 2007/12/14 22:33:52 　　　││　└どっちでしょう？ SECOND 2007/12/14 22:50:11 (修正1回) 　　　│└２０秒速くなりました。島村1243 2007/12/14 22:06:46 　　　└(3)について各点での流線山中和義 2007/12/15 09:43:05 　　　　├!等ポテンシャルを色の濃淡で表現する山中和義 2007/12/15 11:03:55 　　　　└これを利用すると流れの方向がわかるので、... 山中和義 2007/12/15 20:11:13 　　　　　└とても数理的で良いですね。島村1243 2007/12/15 23:15:50 (修正1回) 　　　　　　└接地導体球の作画に問題が有る様島村1243 2007/12/16 07:37:08 　　　　　　　└対策案できました。島村1243 2007/12/16 10:44:04 (修正1回) 　　　　　　　　└ポテンシャル関数f(x,y)なら、電位分布は等... 山中和義 2007/12/16 13:56:21

ポテンシャル関数f(x,y)なら、電位分布は等高線アルゴリズムで描けます。

!関数f(x,y)の等高線を描く

LET cEps=0.001 !許容誤差
LET h=0.001 !増分（偏微分係数）

SUB newton(fc,x,y, fx,fy,grad2) !ニュートン法
DO
LET ff=f(x,y) !∇f
LET fx=(f(x+h,y)-ff)/h
LET fy=(f(x,y+h)-ff)/h
LET grad2=fx*fx+fy*fy

IF grad2<1e-10 THEN
LET x=1e30
EXIT SUB
END IF

LET t=(fc-ff)/grad2
LET x=x+t*fx
LET y=y+t*fy
LOOP WHILE t*t*grad2>cEps*cEps
END SUB

SUB contour(fc,x,y,d) !等高線を描く
LET i=0
DO
CALL newton(fc,x,y, fx,fy,grad2)
IF ABS(x)+ABS(y)>1e10 THEN EXIT SUB
IF i=0 THEN
PLOT LINES: x,y; !始点
LET x0=x
LET y0=y
ELSE
PLOT LINES: x,y; !折れ線でつなげる
END IF
IF i>2 AND (x-x0)^2+(y-y0)^2<d*d THEN EXIT DO !始点近傍なら、終了

LET t=d/SQR(grad2)
LET x=x+fy*t
LET y=y-fx*t

LET i=i+1
LOOP
PLOT LINES: x0,y0 !閉じる
END SUB
!---------- ここまでがサブルーチン

!SET WINDOW -6,6,-6,6 !表示領域
SET WINDOW -15,15,-15,15 !表示領域
DRAW grid !座標

LET r=1 !導体球の半径
LET a1=4 !点電荷のx座標
LET a2=r^2/a1 !影像電荷のx座標
LET q1=1 !点電荷の大きさ
LET q2=-r/a1*q1 !影像電荷の大きさ

DEF f(x,y)=q1/SQR((x-a1)^2+y^2)+q2/SQR((x-a2)^2+y^2) !合成電位関数

SET AREA COLOR "red" !点電荷の位置に赤丸印を表示
DRAW disk WITH SCALE(0.1)*SHIFT(a1,0)

SET AREA COLOR "green" !接地導体球の位置に緑丸印を表示
DRAW disk WITH SCALE(r)*SHIFT(0,0)

LET d=0.015 !増分
FOR fc=0 TO 0.5 STEP d !等高線の関数値
CALL contour(fc,1,1,d)
NEXT fc

FOR fc=0 TO 0.06 STEP d !※左端
CALL contour(fc,-10,1,d)
NEXT fc

END

Re: 対策案できました。	返事を書くノートメニュー
山中和義 <drdlxujciw> 2007/12/16 13:56:21
ポテンシャル関数f(x,y)なら、電位分布は等高線アルゴリズムで描けます。 !関数f(x,y)の等高線を描く LET cEps=0.001 !許容誤差 LET h=0.001 !増分（偏微分係数） SUB newton(fc,x,y, fx,fy,grad2) !ニュートン法 DO LET ff=f(x,y) !∇f LET fx=(f(x+h,y)-ff)/h LET fy=(f(x,y+h)-ff)/h LET grad2=fxfx+fyfy IF grad2<1e-10 THEN LET x=1e30 EXIT SUB END IF LET t=(fc-ff)/grad2 LET x=x+tfx LET y=y+tfy LOOP WHILE ttgrad2>cEpscEps END SUB SUB contour(fc,x,y,d) !等高線を描く LET i=0 DO CALL newton(fc,x,y, fx,fy,grad2) IF ABS(x)+ABS(y)>1e10 THEN EXIT SUB IF i=0 THEN PLOT LINES: x,y; !始点 LET x0=x LET y0=y ELSE PLOT LINES: x,y; !折れ線でつなげる END IF IF i>2 AND (x-x0)^2+(y-y0)^2<dd THEN EXIT DO !始点近傍なら、終了 LET t=d/SQR(grad2) LET x=x+fyt LET y=y-fxt LET i=i+1 LOOP PLOT LINES: x0,y0 !閉じる END SUB !---------- ここまでがサブルーチン !SET WINDOW -6,6,-6,6 !表示領域 SET WINDOW -15,15,-15,15 !表示領域 DRAW grid !座標 LET r=1 !導体球の半径 LET a1=4 !点電荷のx座標 LET a2=r^2/a1 !影像電荷のx座標 LET q1=1 !点電荷の大きさ LET q2=-r/a1q1 !影像電荷の大きさ DEF f(x,y)=q1/SQR((x-a1)^2+y^2)+q2/SQR((x-a2)^2+y^2) !合成電位関数 SET AREA COLOR "red" !点電荷の位置に赤丸印を表示 DRAW disk WITH SCALE(0.1)SHIFT(a1,0) SET AREA COLOR "green" !接地導体球の位置に緑丸印を表示 DRAW disk WITH SCALE(r)*SHIFT(0,0) LET d=0.015 !増分 FOR fc=0 TO 0.5 STEP d !等高線の関数値 CALL contour(fc,1,1,d) NEXT fc FOR fc=0 TO 0.06 STEP d !※左端 CALL contour(fc,-10,1,d) NEXT fc END