aUCBLogo Demos and Tests / molecules3dc
be molecules3dc
mol
end
be mol
max_ = 1000
min_ = 1
maxb = 6
dopt = 35
dbind = trunc 2*dopt
radius = dopt/8
ffein = 100
Temperature=300
cfein = 4000000/Temperature
dopf = dopt*ffein*1.07 ;calibration
dopf2 = dopt*ffein*0.8
maxf = ffein*dbind
tE = 1.2
anfE = 0.4
expo = 5
ep = 1
fac = 1
fac2 = fac*0.001
gravV = -0.1*fac
mov = 1*fac
vfac = 50
cmin = 20
deltaTFac = 1.5
sqrDeltaTFac= deltaTFac*deltaTFac
gravity=false
threeDimensional=1
radx=dopt/2
rady=dopt/2
radz=dopt/2
sizeX=2
sizeY=2
c=IntArray max_
r=Array 3
ri=Array 3
v=Array 3
a=Array 3
for [k 1 3]
[ r.k=FloatArray max_
ri.k=IntArray max_
v.k=FloatArray max_
a.k=FloatArray max_
]
banz=IntArray max_
b=Array max_
for [i 1 max_]
[ b.i=IntArray maxb
]
f=(FloatArray maxf+1 0)
f2=(FloatArray maxf+1 0)
sizehx=Int 400/dbind+1
sizex=2*sizehx+1
sizehy=Int 300/dbind+1
sizey=2*sizehy+1
sizehz=Int 300/dbind+1
sizez=2*sizehz+1
sizeh=Array 3
sizeh.1=sizehx
sizeh.2=sizehy
sizeh.3=sizehz
m=(Array sizex -sizehx)
for [mi -sizehx sizehx]
[ mx=(Array sizey -sizehy)
for [mj -sizehy sizehy]
[ my=(Array sizez -sizehz)
for [mk -sizehz sizehz]
[ my.mk=[]
]
mx.mj=my
]
m.mi=mx
]
onePoint=true
tooSlow = false
tooFast = false
topteil = 0
disposalY= 0
clicked=false
col=[]
lineColor=RGB 1 0 1
bindColor=RGB 1 1 1
eye=Array 3
phi=0
theta=30
dtheta=1
center={0 0 0}
upvector={0 1 0}
dphi=1
ddphi=dphi/3
rotatescene_r=800
dr=1.1
setScreenColor 0
refresh
setUpdateGraph false
perspective
setPenSize [0 0]
hideTurtle
PenUp
be init
be initforcetable
; setItems 0 f (rSeqFA 1 0 int maxf/2)^2/2
; setItems int maxf/2 f (rSeqFA 0 1 int maxf/2)^2/100* -1
; stop
for [i 4 maxf]
[ f.i=fac*( -((i/dopf)^(-expo))+(i/dopf)^(-expo-ep))
]
for [i 0 3]
[ f.i=0
]
for [i 4 maxf]
[ f2.i=fac*( -((i/dopf)^(-expo))+(i/dopf)^(-expo-ep))
+fac*0.07*exp -sqr (i/maxf-0.6)*3
]
for [i 0 3]
[ f2.i=0
]
end
be squarehex side x_ y_ angle v_ vangle
angle= angle
vangle=vangle
vx_=v_*Cos vangle
vy_=v_*Sin vangle
vxx= dopt*Cos angle
vxy= dopt*Sin angle
vyx= dopt*Sin angle
vyy=-dopt*Cos angle
kx=side/2+(mod trunc side/2 2)/2
ky=side/2*(Sqrt 3)/2
x_= x_-(vxx*kx+vyx*ky)
y_= y_-(vxy*kx+vyy*ky)
local [i]
i=1
for [yi 1 side]
[ for [xi 1 side]
[ kx=xi+(mod yi 2)/2
ky=yi*(Sqrt 3)/2
r.(1).i= x_+vxx*kx+vyx*ky
r.(2).i= y_+vxy*kx+vyy*ky
r.(3).i= 0
for [k 1 3]
[ ri.k.i=sizeh.k+round r.k.i/dbind
]
v.(1).i=vx_
v.(2).i=vy_
v.(3).i=0
banz.i=0
ifElse i <= max_
[ i=i+1
][ print [Too many parts!]
]
]
]
topteil=i
end
be squaresquare side x_ y_ angle v_ vangle
angle= angle
vangle=vangle
vx_=v_*Cos vangle
vy_=v_*Sin vangle
vxx= dopt*Cos angle
vxy= dopt*Sin angle
vyx= dopt*Sin angle
vyy=-dopt*Cos angle
kx=side/2+(mod trunc side/2 2)/2
ky=side/2*(Sqrt 3)/2
x_= x_-(vxx*kx+vyx*ky)
y_= y_-(vxy*kx+vyy*ky)
local [i]
i=1
for [yi 1 side]
[ for [xi 1 side]
[ kx=xi
ky=yi
r.(1).i= x_+vxx*kx+vyx*ky
r.(2).i= y_+vxy*kx+vyy*ky
r.(3).i= 0
for [k 1 3]
[ ri.k.i=sizeh.k+round r.k.i/dbind
]
v.(1).i=vx_
v.(2).i=vy_
v.(3).i=0
banz.i=0
ifElse i <= max_
[ i=i+1
][ print [Too many parts!]
]
]
]
topteil=i
end
be cubesc side x_ y_ z_ angle v_ vangle
angle= angle
vangle=vangle
vx_=v_*Cos vangle
vy_=v_*Sin vangle
vxx= dopt*Cos angle
vxy= dopt*Sin angle
vyx= dopt*Sin angle
vyy=-dopt*Cos angle
kx=side/2+(mod trunc side/2 2)/2
ky=side/2*(Sqrt 3)/2
x_= x_-(vxx*kx+vyx*ky)
y_= y_-(vxy*kx+vyy*ky)
z_= z_-dopt
local [i]
i=1
for [yi 1 side]
[ for [xi 1 side]
[ for [zi 1 side]
[ kx=xi
ky=yi
kz=zi
r.(1).i= x_+vxx*kx+vyx*ky
r.(2).i= y_+vxy*kx+vyy*ky
r.(3).i= z_+dopt*kz
for [k 1 3]
[ ri.k.i=sizeh.k+round r.k.i/dbind
]
v.(1).i=vx_
v.(2).i=vy_
v.(3).i=0
banz.i=0
ifElse i <= max_
[ i=i+1
][ print [Too many parts!]
]
]
]
]
topteil=i
end
be cubehcp side x_ y_ z_ angle v_ vangle
angle= angle
vangle=vangle
vx_=v_*Cos vangle
vy_=v_*Sin vangle
vxx= dopt*Cos angle
vxy= dopt*Sin angle
vyx= dopt*Sin angle
vyy=-dopt*Cos angle
kx=side/2+(mod trunc side/2 2)/2
ky=side/2*(Sqrt 3)/2
x_= x_-(vxx*kx+vyx*ky)
y_= y_-(vxy*kx+vyy*ky)
z_= z_-dopt
local [i]
i=1
for [yi 1 side]
[ for [xi 1 side]
[ for [zi 1 side]
[ kx=xi+(mod yi 2)/2-(mod zi 2)/2
ky=yi*(Sqrt 3)/2-(mod zi 2)*(Sqrt 3)/4
kz=zi
r.(1).i= x_+vxx*kx+vyx*ky
r.(2).i= y_+vxy*kx+vyy*ky
r.(3).i= z_+dopt*kz
for [k 1 3]
[ ri.k.i=sizeh.k+round r.k.i/dbind
]
v.(1).i=vx_
v.(2).i=vy_
v.(3).i=0
banz.i=0
ifElse i <= max_
[ i=i+1
][ print [Too many parts!]
]
]
]
]
topteil=i
end
ifElse threeDimensional==1
[ ;cubehcp 5 0 0 0 30 0 90
cubesc 5 0 0 0 30 0 90
][ ;squarehex 10 0 0 30 0 90
squaresquare 10 0 0 30 0 90
]
topteil=topteil-1
initforcetable
setXY 0 -270 setH 90
Label [[RETURN]=splines [+]=heat [-]=cool
[G]=gravity [other Key]=cS Mouse: L=pull R=del]
col=loadpalette "TEILE.PAL
repeat count col
[ i=repcount
cc=reRGB col.i
col.i=RGBA cc.1 cc.2 cc.3 0.5
]
end
be movethem
be faster
local [k l]
for [k 1 topteil]
[ for [l 1 3]
[ v.l.k=v.l.k*deltaTFac
]
]
for [k 0 maxf]
[ f.k=f.k*sqrDeltaTFac
]
gravV:= gravV*sqrDeltaTFac
end
be slower
local [k]
for [k 1 topteil]
[ for [l 1 3]
[ v.l.k=v.l.k/deltaTFac
]
]
for [k 0 maxf]
[ f.k=f.k/sqrDeltaTFac
]
gravV:= gravV/sqrDeltaTFac
end
be preparevars
local [i k]
if tooSlow [faster]
if tooFast [slower]
tooSlow=true
tooFast=false
for [i 1 topteil]
[ for [k 1 3]
[ a.k.i=0
]
c.i=0
]
end
be energyloss
local [h k]
h=FloatArray 3
for [k 1 3]
[ h.k=(v.k.i+v.k.j)/2
]
if i >= min_
[ for [k 1 3]
[ v.k.i=h.k
]
]
banz.i=banz.i+1
b.i=fput j b.i
for [k 1 3]
[ v.k.j=h.k
]
banz.j=banz.j+1
b.j=fput i b.j
end
be ionize
if member? j b.i
[ b.i=remove j b.i
banz.i=banz.i-1
b.j=remove i b.j
banz.j=banz.j-1
]
end
be unboundf i j
output not member? j b.i
end
local [i j bi di
d da
fx fy f0 force _c nomml]
unbound=true
da=FloatArray 3
Tag "nomml
preparevars
for [i 1 topteil]
[ for [ix ri.(1).i-1 ri.(1).i+1]
[ for [iy ri.(2).i-1 ri.(2).i+1]
[ for [iz ri.(3).i-1 ri.(3).i+1]
[ l=m.ix.iy.iz
while [not empty? l]
[ j=first l
l=butFirst l
da.1=r.(1).i-r.(1).j
if (abs da.1) > dbind [continueLoop]
da.2=r.(2).i-r.(2).j
if (abs da.2) > dbind [continueLoop]
da.3=r.(3).i-r.(3).j
if (abs da.3) > dbind [continueLoop]
d=Sqrt
(Sqr da.1)+
(Sqr da.2)+
(Sqr da.3)
if d > dbind
[
ionize
]
if d < dbind
[ unbound=unboundf i j
if unbound
and2 (banz.i < maxb)
and2 (banz.j < maxb)
[ ;if not yet bound & free
if (abs d-dopt) < 0.01
[ ;and d around dopt
energyloss ;then "emitt a Photon"
]
]
]
if unbound and2 (d > dopt)
[
continueLoop
]
d=d*ffein
di=Int d
if di >= maxf-1 [di=maxf-1]
f0=f.di
force=f0 ;+(d-Int d)*(f.(di+1)-f0)
for [k 1 3]
[ a.k.i=a.k.i+da.k*force
]
]
]
]
]
]
for [i min_ topteil]
[ c.i=Int (sqrt
(sqr a.(1).i)+
(sqr a.(2).i)+
(sqr a.(3).i))*cfein
if c.i > 250
[
tooFast=true
tooSlow=false
goto "nomml
]
if c.i > cmin
[ tooSlow=false
]
]
for [i 1 topteil]
[ for [k 1 3]
[ v.k.i=v.k.i+a.k.i
]
if gravity
[ v.(2).i=v.(2).i+gravV
]
rxi=ri.(1).i
ryi=ri.(2).i
rzi=ri.(3).i
m.rxi.ryi.rzi=remove i m.rxi.ryi.rzi
for [k 1 3]
[ r.k.i=r.k.i+v.k.i
]
if r.(1).i < -400+radx or2 r.(1).i > 400-radx
[ v.(1).i=-v.(1).i
r.(1).i=r.(1).i+v.(1).i
]
if r.(2).i < -300+rady or2 r.(2).i > 300-rady
[ v.(2).i=-v.(2).i
r.(2).i=r.(2).i+v.(2).i
]
if r.(3).i < -300+radz or2 r.(3).i > 300-radz
[ v.(3).i=-v.(3).i
r.(3).i=r.(3).i+v.(3).i
]
for [k 1 3]
[ ri.k.i=round r.k.i/dopt
]
rxi=ri.(1).i
ryi=ri.(2).i
rzi=ri.(3).i
m.rxi.ryi.rzi=fPut i m.rxi.ryi.rzi
]
end
be drawthem
be draw i
setXYZ
r.(1).i
r.(2).i
r.(3).i
setPC col.(c.i+2)
Sphere radius
setPC lineColor
PenDown
setXYZ
r.(1).i+v.(1).i*vfac
r.(2).i+v.(2).i*vfac
r.(3).i+v.(3).i*vfac
PenUp
end
cs
for [i 1 topteil]
[ draw i
repeat maxb
[ j=b.i.repcount
if j > 0
[ setXYZ
r.(1).i
r.(2).i
r.(3).i
setPC bindColor
PenDown
setXYZ
r.(1).j
r.(2).j
r.(3).j
PenUp
]
]
]
end
be cooling
local [i]
for [i 1 topteil]
[ for [k 1 3]
[ v.k.i=v.k.i/tE
]
]
end
be heating
local [i]
for [i 1 topteil]
[ for [k 1 3]
[ v.k.i=v.k.i*tE
]
]
end
be stopthem
local [i]
for [i 1 topteil]
[ for [k 1 3]
[ v.k.i=0
]
]
end
be findnearest hx hy
local [i j dmin d]
dmin=IntMax
for [i 1 topteil]
[ d=trunc Sqrt (Sqr hx-r.(1).i)+(Sqr hy-r.(2).i)
if d < dmin
[ dmin=d
j=i
]
]
output j
end
be showmark x y
setPC 12
setXYZ x y z
circle dopt/4
setPixelXYZ x y z 0
updateGraph
setPC 0
setXYZ x y z
circle dopt/4
end
be mousepulling
mou=MousePos
if not clicked
[ clicki=findnearest mou.1 mou.2
clicked=true
]
i=clicki
; showmark(ox,oy);
d=((Sqr mou.1-r.(1).i)+Sqr mou.2-r.(2).i)^0.3
for [k 1 3]
[ v.k.i=v.k.i+mov*(mou.k-r.k.i)/d/dbind/te
; r.k.i=r.k.i+mov*(mou.k-r.k.i)/d
]
end
be mousespecials
local [i mx my]
mx=MouseX
my=MouseY
i=findnearest mx my
showmark r.1.i r.2.i r.3.i
for [i 1 3]
[ v.k.i=0
]
r.1.i=radx
r.2.i=rady+dopt*disposalY
r.3.i=0
disposalY= Mod (disposalY+1) 6
while [MouseButtons!=0]
[ dispatchMessages
]
end
be onCharHandler
ch=KeyboardValue
if ch==wxk_escape [OnChar [] running=false]
if ch==wxk_return [onePoint=not onePoint]
if ch==ASCII "- [cooling]
if ch==ASCII "+ [heating]
if ch==ASCII 0 [stopthem]
if (ch==ASCII "g) or2 (ch==ASCII "G) [gravity=not gravity]
if ch==ASCII "\ [clearScreen]
if ch==wxk_right [phi=phi+dphi]
if ch==wxk_left [phi=phi-dphi]
if ch==wxk_up [theta=theta+dtheta]
if ch==wxk_down [theta=theta-dtheta]
if ch==wxk_prior [rotatescene_r=rotatescene_r/dr]
if ch==wxk_next [rotatescene_r=rotatescene_r*dr]
eye.1=rotatescene_r*(cos theta)*sin phi
eye.2=rotatescene_r* sin theta
eye.3=rotatescene_r*(cos theta)*cos phi
setEye eye center upvector
setLightPos {1000 1000 1000}
redraw
updateGraph
end
init
pd
Line [[-400 0][400 0]] 1
setPixelXY rSeqFA -400 400 maxf+1 f*1000 2
setPixelXY rSeqFA -400 400 maxf+1 f2*1000 4
onCharHandler
ignore readChar
OnChar [onCharHandler]
t0=0
running=true
while [running]
[
;ignore[
Molecules3dMoveThem
topteil maxf maxb r v a ri f f2 c banz b m
gravV dbind dopt deltaTFac
ffein cfein cmin gravity
;] movethem
drawthem
setXYZ -350 -250 0
setPC RGB 1 1 1
PenDown
setHeading 90
t=TimeMilli
Label list 1000/(t-t0) "fps
t0=TimeMilli
PenUp
updateGraph
;updateVars
ConsoleSetFocus
dispatchMessages
ifElse MouseButtons==1
[ mousepulling
][ ifElse MouseButtons==2
[ mousespecials
][ clicked=false
]
]
]
end