|
|
イオン結晶の高速化-分子動力学法プログラム025fasterIonMD2D.basを公開します。（[023]の高速化バージョンです。）
本プログラムは十進BASIC 8.0.0 / macOS 10.7.5 でテストしました。
! ========= molecular dynamics 2D - ion ==========
!
! 025NaClIMD2D.bas
! Copyright(C) 2018 Mitsuru Ikeuchi
! Released under the MIT license ( https://opensource.org/licenses/MIT )
!
! ver 0.0.1 2018.05.11 created
!
OPTION ARITHMETIC NATIVE
DECLARE EXTERNAL SUB imd2d.setInitialCondition, imd2d.moveParticles, imd2d.drawParticles
DECLARE EXTERNAL FUNCTION INKEY$
LET tempMode = 0 !tempMode: 0:adiabatic 1:constant temperature
LET contTemp = 300 !contTemp: controled constant temperature(K)
LET ddTemp = 0 !contTemp = contTemp + ddTemp
LET drawMode = 0 !drawMode: 0:bond 1:circle+bond 2:velocitySpace
LET pauseFlag = 0 !if pauseFlag=0, CALL moveParticles(tempMode,contTemp)
!setInitialCondition(material,boxSizeInNM,contTemp)
CALL setInitialCondition(1,12.0,contTemp) !material = 1:NaCl 2:MgO 3:CaO 4:BaO 5:NaF 6:KF 7:KCl
DO
LET contTemp = contTemp + ddTemp
IF contTemp>3000 THEN LET contTemp = 3000
IF contTemp<10 THEN LET contTemP = 10
IF pauseFlag=0 THEN CALL moveParticles(tempMode,contTemp) ELSE WAIT DELAY 0.05
CALL drawParticles(tempMode,contTemp,drawMode)
LET S$=INKEY$
IF S$="1" OR S$="2" OR S$="3" OR S$="4" OR S$="5" OR S$="6" OR S$="7" THEN
LET material = VAL(S$)
LET tempMode = 0
LET contTemp = 300
CALL setInitialCondition(material,12.0,contTemp)
ELSEIF S$="T" OR S$="t" THEN
LET tempMode = MOD(tempMode+1,2)
IF tempMode=0 THEN LET ddTemp = 0
ELSEIF S$="K" OR S$="k" THEN
LET tempMode = 1
IF ddtemp=0 THEN LET ddTemp = 1 ELSE LET ddTemp = 0
ELSEIF S$="J" OR S$="j" THEN
LET tempMode = 1
IF ddTemp=0 THEN LET ddTemp = -1 ELSE LET ddTemp = 0
ELSEIF S$="D" OR S$="d" THEN
LET drawMode = MOD(drawMode+1,3)
ELSEIF S$=" " THEN
LET pauseFlag = MOD(pauseFlag+1,2)
END IF
LOOP UNTIL S$=CHR$(27)
END
EXTERNAL FUNCTION INKEY$ !--- from decimal BASIC library inkey$.bas
OPTION ARITHMETIC NATIVE
SET ECHO "OFF"
LET S$=""
CHARACTER INPUT NOWAIT: S$
LET INKEY$=S$
END FUNCTION
! ---------- molecular dynamics 2D - ion ----------
!
! method: velocity Verlet ( F=m*d^2r/dt^2 -> r(t+dt)=r(t)+v*dt,v=v+(F/m)*dt )
! (1) vi = vi + (Fi/mi)*(0.5dt)
! (2) ri = ri + vi*dt
! (3) calculation Fi <- {r1,r2,...,rn} Fi=sum(Fij,j=1 to n),Fij=F(ri-rj)
! (4) vi = vi + (Fi/mi)*(0.5dt)
! (6) goto (1)
!
! force: ion f(r) = fc + fr + fa
! fc = eForceConst*zi*zj*(EXP(-r/6.5e-10)/r)*(1.0/r+1.0/6.5e-10) !Debye-screened Coulomb force
! fr = 6.9742e-11*EXP((a-r)/b) !repulsive force
! fa = -6.9742e-21*(c/r^6) !attractive force
!
MODULE imd2d
MODULE OPTION ARITHMETIC NATIVE
PUBLIC SUB setInitialCondition !(molKind,boxSizeInNM,xtalSizeInNM,contTemp)
PUBLIC SUB moveParticles !(tempMode,contTemp)
PUBLIC SUB drawParticles !(drawMode:0:realSpace 1:velocitySpace)
SHARE NUMERIC ionKind1,ionKind2, sysTime, dt, nMolec, xMax, yMax, Nsx,Nsy, rCutoff, hh
SHARE NUMERIC xx(5000),yy(5000) ! (xx(i),yy(i)) : position of i-th particle
SHARE NUMERIC vx(5000),vy(5000) ! (vx(i),vy(i)) : velocity of i-th particle
SHARE NUMERIC ffx(5000),ffy(5000) ! (ffx(i),ffy(i)): total force of i-th particle
SHARE NUMERIC kind(5000),mas(5000) ! kind(i),mas(i) : molec kind, mass of i-th particle
SHARE NUMERIC reg(5000,0 TO 100) ! register near molec reg(i,0):number of near i-th molec
SHARE NUMERIC section(0 TO 101,0 TO 101,0 TO 20) !use pre-registration
SHARE NUMERIC ion_m(0 TO 17),ion_z(0 TO 17) ! ion mass,ion charge
SHARE NUMERIC ion_a(0 TO 17),ion_b(0 TO 17) ! ion force param a,ion force param b
SHARE NUMERIC ion_c(0 TO 17),ion_r(0 TO 17) ! ion force param c, ion radius
SHARE NUMERIC ion_col(0 TO 17) ! ion draw color
SHARE STRING ion_str$(0 TO 17) ! ion string
SHARE NUMERIC forceTable(0 TO 17, 0 TO 17,0 TO 1001) ! force table
LET ionKind1 = 3 ! ion kind1 - 3:Na+
LET ionKind2 = 7 ! ion kind2 - 7:Cl-
LET sysTime = 0.0 ! system time (s) in the module
LET dt = 2.0*1.0e-15 ! time step (s)
LET nMolec = 100 ! number of particles
LET xMax = 6.0E-9 ! x-Box size (m)
LET yMax = 6.0E-9 ! y-Box size (m)
LET Nsx = 100
LET Nsy = 100
LET rCutoff = 1.0e-9 ! force cutoff radius (m)
LET hh = 1e-12 ! force table step (m)
LET ion_z(0)=0.0 ! if ion_z(0)=0.0 then set ion data
! ---------- set initial condition
EXTERNAL SUB setInitialCondition(material,boxSizeInNM,contTemp)
DECLARE EXTERNAL SUB setIonData,ajustVelocity,setForceTable
DECLARE EXTERNAL FUNCTION setNaClTypeBlock
RANDOMIZE
LET sysTime = 0.0
CALL setIonData
CALL setForceTable
LET xMax = boxSizeInNM*1.0e-9
LET yMax = boxSizeInNM*1.0e-9
IF material=1 THEN !NaCl
LET ionKind1 = 3 !Na+
LET ionKind2 = 7 !Cl-
LET lattice = 5.6407e-10
ELSEIF material=2 THEN !MgO
LET ionKind1 = 4 !Mg++
LET ionKind2 = 1 !O--
LET lattice = 4.212e-10*0.94 !0.94: correction factor
ELSEIF material=3 THEN !CaO
LET ionKind1 = 9 !Ca++
LET ionKind2 = 1 !O--
LET lattice = 4.80e-10*0.94 !0.94: correction factor
ELSEIF material=4 THEN !BaO
LET ionKind1 = 12 !Ba++
LET ionKind2 = 1 !O--
LET lattice = 5.536e-10
ELSEIF material=5 THEN !NaF
LET ionKind1 = 3 !Na+
LET ionKind2 = 2 !F-
LET lattice = 4.62e-10
ELSEIF material=6 THEN !KF
LET ionKind1 = 8 !K+
LET ionKind2 = 2 !F-
LET lattice = 5.34e-10
ELSEIF material=7 THEN !KCl
LET ionKind1 = 8 !K+
LET ionKind2 = 7 !Cl-
LET lattice = 6.29e-10
END IF
LET nL = int(0.8*xMax/lattice)
LET s = 0.5*(xMax - lattice*nL)
! setNaClTypeBlock(ii, knd1, knd2, nx, ny, lattice, xPos, yPos)
LET nMolec = setNaClTypeBlock(1, ionKind1, ionKind2, nL, nL, lattice, s, s)
CALL ajustVelocity(contTemp)
! set window
SET WINDOW 0,500,0,500
END SUB
EXTERNAL SUB setIonData
! ion potential data
! 0 mass,1 charge, 2 a , 3 b , 4 c , 5 r-ion, 6 color 7 str$
DATA 10.81, 3.0, 0.720e-10, 0.080e-10, 0.0, 0.23e-10, 13 , "B+++" !0
DATA 16.00, -2.0, 1.626e-10, 0.085e-10, 20.0, 1.40e-10, 2 , "O--" !1
DATA 19.00, -1.0, 1.565e-10, 0.085e-10, 20.0, 1.33e-10, 2 , "F-" !2
DATA 22.99, 1.0, 1.260e-10, 0.080e-10, 20.0, 1.02e-10, 4 , "Na+" !3
DATA 24.31, 2.0, 1.161e-10, 0.080e-10, 10.0, 0.72e-10, 4 , "Mg++" !4
DATA 26.98, 3.0, 1.064e-10, 0.080e-10, 2.0, 0.53e-10, 13 , "Al+++" !5
DATA 28.09, 4.0, 1.012e-10, 0.080e-10, 0.0, 0.40e-10, 7 , "Si++++" !6
DATA 35.45, -1.0, 1.950e-10, 0.090e-10, 30.0, 1.81e-10, 2 , "Cl-" !7
DATA 39.10, 1.0, 1.595e-10, 0.080e-10, 15.0, 1.38e-10, 4 , "K+" !8
DATA 40.08, 2.0, 1.414e-10, 0.080e-10, 10.0, 1.00e-10, 4 , "Ca++" !9
DATA 47.88, 4.0, 1.235e-10, 0.080e-10, 0.0, 0.61e-10, 7 , "Ti++++" !10
DATA 87.62, 2.0, 1.632e-10, 0.080e-10, 15.0, 1.16e-10, 4 , "Sr++" !11
DATA 137.3, 2.0, 1.820e-10, 0.080e-10, 20.0, 1.36e-10, 4 , "Ba++" !12
DATA 4.003, 0.0, 1.200e-10, 0.110e-10, 4.76, 1.28e-10, 3 , "He" !13
DATA 20.18, 0.0, 1.415e-10, 0.112e-10,11.03, 1.37e-10, 3 , "Ne" !14
DATA 39.95, 0.0, 1.878e-10, 0.117e-10,38.53, 1.70e-10, 3 , "Ar" !15
DATA 83.80, 0.0, 2.041e-10, 0.130e-10,55.33, 1.83e-10, 3 , "Kr" !16
DATA 131.3, 0.0, 2.258e-10, 0.145e-10,85.55, 1.99e-10, 3 , "Xe" !17
IF ion_z(0)=0.0 THEN
FOR i=0 TO 17
READ ion_m(i),ion_z(i),ion_a(i),ion_b(i),ion_c(i),ion_r(i),ion_col(i),ion_str$(i)
LET ion_m(i) = ion_m(i)*1.661e-27
LET ion_z(i) = ion_z(i)*1.602e-19
NEXT i
END IF
END SUB
EXTERNAL FUNCTION setNaClTypeBlock(ii, knd1, knd2, nx, ny, lattice, xPos, yPos)
DECLARE EXTERNAL SUB setParticle
LET ipp = ii
LET a = lattice/2.0
FOR i=0 TO 2*nx-1
FOR j=0 TO 2*ny-1
LET x = xPos + a*i
LET y = yPos + a*j
IF MOD(i+j,2)=0 THEN
LET knd = knd1
ELSE
LET knd = knd2
END IF
CALL setParticle(ipp, knd, x, y)
LET ipp = ipp + 1
NEXT j
NEXT i
LET setNaClTypeBlock = ipp - 1
END FUNCTION
EXTERNAL SUB setParticle(i, knd, x, y)
LET xx(i) = x
LET yy(i) = y
LET vx(i) = 200.0*(RND+RND+RND+RND+RND+RND-3)
LET vy(i) = 200.0*(RND+RND+RND+RND+RND+RND-3)
LET ffx(i) = 0.0
LET ffy(i) = 0.0
LET mas(i) = ion_m(knd)
LET kind(i) = knd
END SUB
EXTERNAL SUB setForceTable
DECLARE EXTERNAL FUNCTION cutoff
LET eForceConst = 1.0/(4.0*PI*8.8542e-12) !epsilon0=8.8542e-12
FOR ki=0 TO 17
FOR kj=0 TO 17
LET a = ion_a(ki)+ion_a(kj)
LET b = ion_b(ki)+ion_b(kj)
LET c = ion_c(ki)*ion_c(kj)*1.0e-60
LET zi = ion_z(ki)
LET zj = ion_z(kj)
FOR ir=10 TO 1001
LET r = ir*hh
LET fc = eForceConst*zi*zj*(EXP(-r/6.5e-10)/r)*(1.0/r+1.0/6.5e-10) !Debye-screened Coulomb force
LET fr = 6.9742e-11*EXP((a-r)/b) !repulsive force
LET fa = -6.9742e-21*(c/(r*r*r*r*r*r)) !attractive force
LET forceTable(ki,kj,ir) = cutoff(r)*(fc + fr + fa)
NEXT ir
FOR ir=0 TO 9
LET forceTable(ki,kj,ir) = forceTable(ki,kj,10)
NEXT ir
NEXT kj
NEXT ki
END SUB
EXTERNAL FUNCTION cutoff(r)
IF r>0 AND r<0.8*rCutoff THEN
LET ret = 1
ELSEIF r>=0.8*rCutoff AND r<rCutoff THEN
LET ret = 0.5+0.5*COS(PI*(r-0.8*rCutoff)/(0.2*rCutoff))
else
LET ret = 0
END IF
LET cutoff = ret
END FUNCTION
! ---------- move particles
EXTERNAL SUB moveParticles(tempMode,contTemp) !tempMode 0:adiabatic 1:constant-temp
DECLARE EXTERNAL SUB moveParticlesDT,ajustVelocity,registerNearMolec,registration
IF (tempMode=1) THEN CALL ajustVelocity(contTemp)
CALL registration
!CALL registerNearMolec
FOR i=1 TO 20
CALL moveParticlesDT
NEXT i
END SUB
EXTERNAL SUB moveParticlesDT ! velocity Verlet method
DECLARE EXTERNAL SUB calcForce
FOR i=1 TO nMolec
LET a = 0.5*dt/mas(i)
LET vx(i) = vx(i)+a*ffx(i)
LET vy(i) = vy(i)+a*ffy(i)
LET xx(i) = xx(i)+vx(i)*dt
LET yy(i) = yy(i)+vy(i)*dt
NEXT i
CALL calcForce
FOR i=1 TO nMolec
LET a = 0.5*dt/mas(i)
LET vx(i) = vx(i)+a*ffx(i)
LET vy(i) = vy(i)+a*ffy(i)
NEXT i
LET sysTime=sysTime+dt
END SUB
EXTERNAL SUB calcForce
DECLARE EXTERNAL FUNCTION force,boundaryForce
LET s = 0.5*3.418e-10
FOR i=1 TO nMolec
LET ffx(i) = 0.0
LET ffy(i) = 0.0
NEXT i
FOR i=1 TO nMolec-1
FOR k=1 TO reg(i,0)-1
LET j = reg(i,k)
LET xij = xx(i)-xx(j)
LET yij = yy(i)-yy(j)
LET rij = SQR(xij*xij+yij*yij)
IF rij<rCutoff THEN
LET f = force(rij,kind(i),kind(j))
LET fxij = f*xij/rij
LET fyij = f*yij/rij
LET ffx(i) = ffx(i)+fxij
LET ffy(i) = ffy(i)+fyij
LET ffx(j) = ffx(j)-fxij
LET ffy(j) = ffy(j)-fyij
END IF
NEXT k
NEXT i
FOR i=1 TO nMolec ! boundary force
LET ffx(i) = ffx(i)+boundaryForce(xx(i)+s)+boundaryForce(xx(i)-xMax-s)
LET ffy(i) = ffy(i)+boundaryForce(yy(i)+s)+boundaryForce(yy(i)-yMax-s)
NEXT i
END SUB
EXTERNAL FUNCTION force(r,ki,kj) !force(r) <-- forceTable - linear interporation
LET ir = INT(r/hh)
LET a = r - ir*hh
LET force = ((hh-a)*forceTable(ki,kj,ir) + a*forceTable(ki,kj,ir+1))/hh
END FUNCTION
EXTERNAL FUNCTION boundaryForce(r)
LET r6 = (3.418e-10/r)^6
LET boundaryForce = (24.0*(0.5*1.711e-21)*r6*(2.0*r6-1.0)/r)
END FUNCTION
EXTERNAL SUB registerNearMolec
LET rCut = rCutoff+20*2000*dt
LET rcut2 = rCut*rCut
FOR i=1 TO nMolec-1
LET k = 1
FOR j=i+1 TO nMolec
LET r2 = (xx(i)-xx(j))*(xx(i)-xx(j))+(yy(i)-yy(j))*(yy(i)-yy(j))
IF (r2<rcut2) THEN
LET reg(i,k) = j
LET k = k + 1
END IF
NEXT j
LET reg(i,0) = k
NEXT i
END SUB
EXTERNAL FUNCTION maxNearMolec
LET mx = 0
FOR i=1 TO nMolec-1
IF mx<reg(i,0) THEN LET mx = reg(i,0)
NEXT i
LET maxNearMolec = mx-1
END FUNCTION
EXTERNAL SUB registration !faster registration
DECLARE EXTERNAL SUB preRegistration
CALL preRegistration
LET rreg = rCutoff+20*2000*dt
LET rreg2 = rreg*rreg
FOR ipp=1 TO nMolec-1
LET kp = 1
LET i0 = INT(Nsx*(xx(ipp)-rreg)/xMax)
IF (i0<0) THEN LET i0 = 0
LET i1 = INT(Nsx*(xx(ipp)+rreg)/xMax )
IF (i1>=Nsx) THEN LET i1 = Nsx-1
LET j0 = INT(Nsy*(yy(ipp)-rreg)/yMax )
IF (j0<0) THEN LET j0 = 0
LET j1 = INT(Nsy*(yy(ipp)+rreg)/yMax )
IF (j1>=Nsy) THEN LET j1 = Nsy-1
FOR i=i0 TO i1
FOR j=j0 TO j1
FOR iq=1 TO section(i,j,0)
LET jp = section(i,j,iq)
IF (jp>ipp) THEN
LET r2=(xx(ipp)-xx(jp))*(xx(ipp)-xx(jp))+(yy(ipp)-yy(jp))*(yy(ipp)-yy(jp))
IF (r2<rreg2) THEN
LET reg(ipp,kp) = jp
LET kp = kp + 1
END IF
END IF
NEXT iq
NEXT j
NEXT i
LET reg(ipp,0) = kp
NEXT ipp
END SUB
EXTERNAL SUB preRegistration
FOR i=0 TO Nsx-1
FOR j=0 TO Nsy-1
LET section(i,j,0) = 0
NEXT j
NEXT i
FOR ipp=1 TO nMolec
LET i = INT(Nsx*xx(ipp)/xMax)
IF i>=Nsx THEN LET i = Nsx-1
LET j = INT(Nsy*yy(ipp)/yMax)
IF j>=Nsy THEN LET j = Nsy-1
LET iq = section(i,j,0) + 1
LET section(i,j,0) = iq
LET section(i,j,iq) = ipp
NEXT ipp
END SUB
! ---------- utility
EXTERNAL FUNCTION systemTemprature
LET kB = 1.38e-23 ! Boltzman's constant (J/K)
LET ek= 0.0 !kinetic energy (J)
FOR i=1 TO nMolec
LET ek = ek + 0.5*mas(i)*(vx(i)^2+vy(i)^2)
NEXT i
LET systemTemprature = ek/(nMolec*kB) !2D: E/N=kT, 3D: E/N=(3/2)kT
END FUNCTION
EXTERNAL SUB ajustVelocity(temp)
DECLARE EXTERNAL FUNCTION systemTemprature
LET r = sqr(temp/systemTemprature)
FOR i=1 TO nMolec
LET vx(i) = r*vx(i)
LET vy(i) = r*vy(i)
NEXT i
END SUB
! ---------- draw particles
EXTERNAL SUB drawParticles(tempMode,contTemp,drawMode)
DECLARE EXTERNAL FUNCTION systemTemprature,maxNearMolec
DECLARE EXTERNAL sub realSpace,velocitySpace,plotBond
SET DRAW MODE HIDDEN
CLEAR
IF drawMode=0 OR drawMode=1 THEN !--- 0:circle 1:circle+bond
call plotBond(drawMode)
ELSEIF drawMode=2 THEN
call velocitySpace
END IF
!--- draw caption
SET TEXT HEIGHT 10
SET TEXT COLOR 1 ! black
LET tmp$ = "adiabatic constantTemp "
PLOT TEXT, AT 50, 70 ,USING "time =#####.## (ps) temp =####.## (K)":sysTime*1E12,systemTemprature
PLOT TEXT, AT 50, 55 :ion_str$(ionKind1)&","&ion_str$(ionKind2)
PLOT TEXT, AT 100, 55 ,USING " N =####":nMolec
PLOT TEXT, AT 50, 40 ,USING "tempMode = ## ":tempMode
PLOT TEXT, AT 200, 40 :tmp$(tempMode*12+1:tempMode*12+12)
PLOT TEXT, AT 50, 25 ,USING "controled Temperature =####.# (K)":contTemp
PLOT TEXT, AT 50, 10 :"2-dimensional ion - molecular dynamics"
PLOT TEXT, AT 50,480 :"[esc] exit [space]pause/go [D]change draw mode"
PLOT TEXT, AT 50,460 :"[1]NaCl [2]MgO [3]CaO [4]BaO [5]NaF [6]KF [7]KCl"
PLOT TEXT, AT 50,440 :"[T] toggle 0:adiabatic/1:constant temperature"
PLOT TEXT, AT 50,420 :"temp control -> [J]coolDown/stop [k]heatUp/stop"
SET DRAW MODE EXPLICIT
END SUB
EXTERNAL sub plotBond(drawMode)
LET boxSize = 300
LET mag = boxSize/xMax
LET xp = 100
LET yp = 100
SET LINE COLOR 1 ! black : !--- box
PLOT LINES: xp,yp; boxSize+xp,yp; boxSize+xp,boxSize+yp; xp,boxSize+yp; xp,yp
SET TEXT HEIGHT 6
SET TEXT COLOR 1 ! black
PLOT TEXT, AT xp,boxSize+2+yp ,USING "box size =##.# x ##.# (nm)":xMax*1e9,yMax*1e9
FOR i=1 TO nMolec
SET LINE COLOR ion_col(kind(i))
DRAW circle WITH SCALE(ion_r(kind(i))*mag)*SHIFT(xx(i)*mag+xp,yy(i)*mag+yp)
NEXT i
IF drawMode=1 THEN
FOR i=1 TO nMolec-1
FOR k=1 TO reg(i,0)-1
LET j = reg(i,k)
LET r = SQR((xx(i)-xx(j))*(xx(i)-xx(j))+(yy(i)-yy(j))*(yy(i)-yy(j)))
LET r0 = (ion_r(kind(i))+ion_r(kind(j)))
IF r<1.1*r0 AND kind(i)<>kind(j) THEN
SET LINE COLOR 8 !gray
PLOT LINES: xx(i)*mag+xp,yy(i)*mag+yp;xx(j)*mag+xp,yy(j)*mag+yp
END IF
NEXT k
NEXT i
END IF
END sub
EXTERNAL sub velocitySpace
LET boxSize = 300
LET xp = 100
LET yp = 100
SET LINE COLOR 1 !black : axis
PLOT LINES: xp,boxSize/2+yp; boxSize+xp,boxSize/2+yp !vx-axis
PLOT LINES: boxSize/2+xp,yp; boxSize/2+xp,boxSize+yp !vy-axis
SET TEXT HEIGHT 6
SET TEXT COLOR 1 ! black
PLOT TEXT, AT boxSize+xp,boxSize/2+yp: "vx"
PLOT TEXT, AT boxSize+xp,boxSize/2-12+yp: "2000m/s"
PLOT TEXT, AT boxSize/2-12+xp,boxSize+yp: "vy 2000m/s"
PLOT TEXT, AT boxSize/2-8+xp,boxSize/2-10+yp: "0"
PLOT TEXT, AT xp,boxSize+8+yp: "velocity space (vx,vy)"
LET mag = boxSize/4000
FOR i=1 TO nMolec
IF kind(i)=ionKind1 THEN SET LINE COLOR 4 ELSE SET LINE COLOR 2 !4:red, 2:blue
DRAW circle WITH SCALE(5)*SHIFT(vx(i)*mag+boxSize/2+xp,vy(i)*mag+boxSize/2+yp)
NEXT i
END sub
END MODULE
http://mike1336.web.fc2.com/
|
|