vmd教程 forces-tutorial(10)

vmd教程 forces-tutorial

1TCLFORCES

setcomsum"000"

settotalmass0

foreachatom$atoms{

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

setcomsum[vecadd$comsum$tmp

settotalmass[expr$totalmass+$masses($atom)]

}

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

print"Centerofmass=$com"

5Thelastpartofthescriptiswheretherealdi erenceslie:

foreachatom$atoms{

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

setx[lindex$r0]

sety[lindex$r1]

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

setphi[expratan2($y,$x)+$MPI/2]

if{$atom==1}{

print"atom$atom:

}phi=$phi"10

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

setangmag[expr$masses($atom)*$angaccelnamd*$rho]

setangforce[vecscale$angmag$angdir]

setforce[vecadd$linforce$angforce]

addforce$atom$force

}

}

The rstsectionabovecalculatesthedistancefromtherotationaxis,rho,andtheanglefromthexaxisatwhichtheforceshouldbeapplied,phi(90 fromtheangleofthepositionvector,hencetheadditionofπ/2radians.)Fig.3showshowthesevariablesarerelated.

Afterprintingtheforceangleassociatedwithoneatomsowecanmonitortherotation,we ndaunitvector(i.e.avectorwithlengthone)pointinginthephidirection.Wethenmultiplythisbythemassoftheatomanditsdistancefromourrotationaxis,givingusthecorrectangularforcefortheangularaccelerationrequestedintheNAMDcon guration le.Finally,weaddthelinearandangularforcestogetthetotalforce.