vmd教程 forces-tutorial

1TCLFORCES18

Therestoftheforcecalculationcodeisidenticaltothecodewehadinthelastexample,exceptthatinsteadofcallingaddforcewitheachforcevector,wenowputthevectorintheforceslistvariable:

if{$forcecount==$forcesRecalcFreq}{

print"Recalculatingforcesat$printcount"

loadcoordscoords

setcomsum"000"

settotalmass0

foreachatom$atomsmass$masses{

settmp[vecscale$mass$coords($atom)]

setcomsum[vecadd$comsum$tmp]

settotalmass[expr$totalmass+$mass]

}

setcom[vecscale[expr1.0/$totalmass]$comsum]

print"Centerofmass=$com"

setforces{}

foreachatom$atomsmass$masses{

setlinforce[vecscale$mass$linaccelsetr[vecsub$coords($atom)$com]

setx[lindex$r0]

sety[lindex$r1]

setrho[exprsqrt($x*$x+$y*$y)]

setphi[expratan2($y,$x)+$Mif{$atom==1}{

print"atom$atom:

}phi=$phi"

setangdir"[exprcos($phi)][exprsin($phi)]0.0"

setangmag[expr$angaccel*$rho*$mass]

setangforce[vecscale$angdir]

setforce[vecadd$linforce$angforce]

lappendforces$force

}

7Finally,weprintsomeinformation,callclearconfig,andresetforcecount.

print"Step${printcount}:Recalculated\

[llength$forces]forces"

setforcecount0

clearconfig

}

incrforcecount

incrprintcount

return

}