!２重振子Chaos
LET g= 9.8 !m/s^2
LET m1=0.1 !kg
LET m2=0.1 !kg
LET L1= 5 !m
LET L2= 5 !m
!
LET dt=0.05 !sec. 演算ピッチ。
!
LET μ2=m2/(m1+m2)
LET L21=L2/L1
DEF ss1(w2,θ1,θ2)=-g/L1*SIN(θ1) -μ2*L21*w2^2*SIN(θ1-θ2)
DEF ss2(w1,θ1,θ2)=-g/L2*SIN(θ2) +w1^2*SIN(θ1-θ2)/L21
DEF D(θ1,θ2)=1-μ2*COS(θ1-θ2)^2
DEF α1(w1,w2,θ1,θ2)=( ss1(w2,θ1,θ2) -L21*μ2*COS(θ1-θ2)*ss2(w1,θ1,θ2) )/D(θ1,θ2)
DEF α2(w1,w2,θ1,θ2)=(-ss1(w2,θ1,θ2)*COS(θ1-θ2)/L21 +ss2(w1,θ1,θ2) )/D(θ1,θ2)
SUB RK(θ1,θ2,w1,w2)
LET w11=w1
LET w12=w2
LET α11=α1(w1,w2,θ1,θ2)
LET α12=α2(w1,w2,θ1,θ2)
!
LET w21=w1+α11*dt/2
LET w22=w2+α12*dt/2
LET α21=α1(w21,w22,θ1+w11*dt/2,θ2+w12*dt/2)
LET α22=α2(w21,w22,θ1+w11*dt/2,θ2+w12*dt/2)
!
LET w31=w1+α21*dt/2
LET w32=w2+α22*dt/2
LET α31=α1(w31,w32,θ1+w21*dt/2,θ2+w22*dt/2)
LET α32=α2(w31,w32,θ1+w21*dt/2,θ2+w22*dt/2)
!
LET w41=w1+α31*dt
LET w42=w2+α32*dt
LET α41=α1(w41,w42,θ1+w31*dt,θ2+w32*dt)
LET α42=α2(w41,w42,θ1+w31*dt,θ2+w32*dt)
!
LET θ1=θ1+(w11+2*w21+2*w31+w41)*dt/6
LET θ2=θ2+(w12+2*w22+2*w32+w42)*dt/6
LET w1=w1+(α11+2*α21+2*α31+α41)*dt/6
LET w2=w2+(α12+2*α22+2*α32+α42)*dt/6
END SUB
!----init.
LET a_1=PI*0.8 !初期角度１
LET a_2=PI*0.9 ! ～ ２
LET a_3=0 ! 初期角速度１
LET a_4=0 ! ～ ２
!
LET b_1=-a_1+0.001
LET b_2=-a_2
LET b_3=0
LET b_4=0
!
LET c_1=a_1
LET c_2=a_2+0.002
LET c_3=0
LET c_4=0
!
LET d_1=-a_1
LET d_2=-a_2+0.003
LET d_3=0
LET d_4=0
!
!----run
LET w=14
SET WINDOW -w,w,-w,w
SET LINE width 2 !4
SET LINE COLOR 2 !43
LET r1=SQR(m1)
LET r2=SQR(m2)
LET t0=TIME
DO
LET t=TIME
IF dt=< ABS(t-t0) THEN
SET DRAW mode hidden
CLEAR
PLOT TEXT,AT 0.25*w,0.9*w:"マウス 右ボタンで、終了。"
PLOT TEXT,AT -0.98*w,0.93*w,USING"演算ピッチ=#.### 秒":dt
PLOT TEXT,AT -0.98*w,0.87*w,USING"描画ピッチ=#.### 秒":t-t0
LET t0=t
SET AREA COLOR 15
DRAW disk WITH SCALE(3.86,4.67)
SET AREA COLOR 1
DRAW PDL1X2(a_1,a_2) WITH ROTATE(a_1)*SHIFT(-3,3)
DRAW PDL1X2(b_1,b_2) WITH ROTATE(b_1)*SHIFT(3,3)
DRAW PDL1X2(c_1,c_2) WITH ROTATE(c_1)*SHIFT(-3,-3)
DRAW PDL1X2(d_1,d_2) WITH ROTATE(d_1)*SHIFT(3,-3)
CALL RK(a_1,a_2,a_3,a_4)
CALL RK(b_1,b_2,b_3,b_4)
CALL RK(c_1,c_2,c_3,c_4)
CALL RK(d_1,d_2,d_3,d_4)
SET DRAW mode explicit
END IF
WAIT DELAY 0 !ノートパソコン等の消費電力を押える。
MOUSE POLL mx,my,mlb,mrb
LOOP UNTIL mrb=1
PICTURE PDL1X2(θ1,θ2)
DRAW circle WITH SCALE(0.3)
DRAW PDLM(L1,r1)
DRAW PDLM(L2,r2) WITH ROTATE(θ2-θ1)*SHIFT(0,-L1)
END PICTURE
PICTURE PDLM(L,r)
PLOT LINES: 0,0;0,-L
DRAW disk WITH SCALE(r)*SHIFT(0,-L)
END PICTURE
リラックス
投稿者：SECOND 投稿日：2008年11月27日(木)19時09分8秒!２重振子Chaos
LET g= 9.8 !m/s^2
LET m1=0.1 !kg
LET m2=0.1 !kg
LET L1= 5 !m
LET L2= 5 !m
!
LET dt=0.05 !sec. 演算ピッチ。
!
LET μ2=m2/(m1+m2)
LET L21=L2/L1
DEF ss1(w2,θ1,θ2)=-g/L1*SIN(θ1) -μ2*L21*w2^2*SIN(θ1-θ2)
DEF ss2(w1,θ1,θ2)=-g/L2*SIN(θ2) +w1^2*SIN(θ1-θ2)/L21
DEF D(θ1,θ2)=1-μ2*COS(θ1-θ2)^2
DEF α1(w1,w2,θ1,θ2)=( ss1(w2,θ1,θ2) -L21*μ2*COS(θ1-θ2)*ss2(w1,θ1,θ2) )/D(θ1,θ2)
DEF α2(w1,w2,θ1,θ2)=(-ss1(w2,θ1,θ2)*COS(θ1-θ2)/L21 +ss2(w1,θ1,θ2) )/D(θ1,θ2)
SUB RK(θ1,θ2,w1,w2)
LET w11=w1
LET w12=w2
LET α11=α1(w1,w2,θ1,θ2)
LET α12=α2(w1,w2,θ1,θ2)
!
LET w21=w1+α11*dt/2
LET w22=w2+α12*dt/2
LET α21=α1(w21,w22,θ1+w11*dt/2,θ2+w12*dt/2)
LET α22=α2(w21,w22,θ1+w11*dt/2,θ2+w12*dt/2)
!
LET w31=w1+α21*dt/2
LET w32=w2+α22*dt/2
LET α31=α1(w31,w32,θ1+w21*dt/2,θ2+w22*dt/2)
LET α32=α2(w31,w32,θ1+w21*dt/2,θ2+w22*dt/2)
!
LET w41=w1+α31*dt
LET w42=w2+α32*dt
LET α41=α1(w41,w42,θ1+w31*dt,θ2+w32*dt)
LET α42=α2(w41,w42,θ1+w31*dt,θ2+w32*dt)
!
LET θ1=θ1+(w11+2*w21+2*w31+w41)*dt/6
LET θ2=θ2+(w12+2*w22+2*w32+w42)*dt/6
LET w1=w1+(α11+2*α21+2*α31+α41)*dt/6
LET w2=w2+(α12+2*α22+2*α32+α42)*dt/6
END SUB
!----init.
LET a_1=PI*0.8 !初期角度１
LET a_2=PI*0.9 ! ～ ２
LET a_3=0 ! 初期角速度１
LET a_4=0 ! ～ ２
!
LET b_1=-a_1+0.001
LET b_2=-a_2
LET b_3=0
LET b_4=0
!
LET c_1=a_1
LET c_2=a_2+0.002
LET c_3=0
LET c_4=0
!
LET d_1=-a_1
LET d_2=-a_2+0.003
LET d_3=0
LET d_4=0
!
!----run
LET w=14
SET WINDOW -w,w,-w,w
SET LINE width 2 !4
SET LINE COLOR 2 !43
LET r1=SQR(m1)
LET r2=SQR(m2)
LET t0=TIME
DO
LET t=TIME
IF dt=< ABS(t-t0) THEN
SET DRAW mode hidden
CLEAR
PLOT TEXT,AT 0.25*w,0.9*w:"マウス 右ボタンで、終了。"
PLOT TEXT,AT -0.98*w,0.93*w,USING"演算ピッチ=#.### 秒":dt
PLOT TEXT,AT -0.98*w,0.87*w,USING"描画ピッチ=#.### 秒":t-t0
LET t0=t
SET AREA COLOR 15
DRAW disk WITH SCALE(3.86,4.67)
SET AREA COLOR 1
DRAW PDL1X2(a_1,a_2) WITH ROTATE(a_1)*SHIFT(-3,3)
DRAW PDL1X2(b_1,b_2) WITH ROTATE(b_1)*SHIFT(3,3)
DRAW PDL1X2(c_1,c_2) WITH ROTATE(c_1)*SHIFT(-3,-3)
DRAW PDL1X2(d_1,d_2) WITH ROTATE(d_1)*SHIFT(3,-3)
CALL RK(a_1,a_2,a_3,a_4)
CALL RK(b_1,b_2,b_3,b_4)
CALL RK(c_1,c_2,c_3,c_4)
CALL RK(d_1,d_2,d_3,d_4)
SET DRAW mode explicit
END IF
WAIT DELAY 0 !ノートパソコン等の消費電力を押える。
MOUSE POLL mx,my,mlb,mrb
LOOP UNTIL mrb=1
PICTURE PDL1X2(θ1,θ2)
DRAW circle WITH SCALE(0.3)
DRAW PDLM(L1,r1)
DRAW PDLM(L2,r2) WITH ROTATE(θ2-θ1)*SHIFT(0,-L1)
END PICTURE
PICTURE PDLM(L,r)
PLOT LINES: 0,0;0,-L
DRAW disk WITH SCALE(r)*SHIFT(0,-L)
END PICTURE
END