1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
|
! -----------------------------------------------------------------------------
! OpenMM(tm) HelloArgon example in Fortran 95 (June 2009)
! -----------------------------------------------------------------------------
! This program demonstrates a simple molecular simulation using the OpenMM
! API for GPU-accelerated molecular dynamics simulation. The primary goal is
! to make sure you can compile, link, and run with OpenMM and view the output.
! The example is available in C++, C, and Fortran 95.
!
! The system modeled here is a small number of argon atoms in a vacuum.
! A multi-frame PDB file is written to stdout which can be read by VMD or
! other visualization tool to produce an animation of the resulting trajectory.
! -----------------------------------------------------------------------------
INCLUDE 'OpenMMFortranModule.f90'
PROGRAM HelloArgon
use OpenMM; implicit none
type(OpenMM_System) system
type(OpenMM_VerletIntegrator) verlet
type(OpenMM_Context) context
type(OpenMM_Platform) platform
type(OpenMM_NonbondedForce) nonbond
type(OpenMM_Vec3Array) initPosInNm
type(OpenMM_State) state
type(OpenMM_StringArray) pluginList
real*8 timeInPs
integer*4 a, ix, frameNum
character*10 platformName
character*100 dirName
! Load any shared libraries containing GPU implementations.
call OpenMM_Platform_getDefaultPluginsDirectory(dirName)
call OpenMM_Platform_loadPluginsFromDirectory(dirName, pluginList)
call OpenMM_StringArray_destroy(pluginList)
! Create a system with nonbonded forces. System takes ownership
! of Force; don't destroy it yourself. (We're using transfer here
! to recast the specific NonbondedForce to a general Force.)
call OpenMM_System_create(system)
call OpenMM_NonbondedForce_create(nonbond)
ix = OpenMM_System_addForce(system, transfer(nonbond, OpenMM_Force(0)))
! Create three atoms.
call OpenMM_Vec3Array_create(initPosInNm, 3)
do a=1,3
! Space the atoms out evenly by atom index.
call OpenMM_Vec3Array_set(initPosInNm, a, (/ 0.5d0*(a-1), 0d0, 0d0 /))
ix = OpenMM_System_addParticle(system, 39.95d0) !mass of Ar, grams/mole
! charge, L-J sigma (nm), well depth (kJ) (vdWRad(Ar)=.188 nm)
ix = OpenMM_NonbondedForce_addParticle(nonbond, 0d0, 0.3350d0, 0.996d0)
end do
! Create particular integrator, and recast to generic one.
call OpenMM_VerletIntegrator_create(verlet, 4d-3) !step size in ps
! Let OpenMM Context choose best platform.
call OpenMM_Context_create(context, system, &
transfer(verlet, OpenMM_Integrator(0)))
call OpenMM_Context_getPlatform(context, platform)
call OpenMM_Platform_getName(platform, platformName)
print "('REMARK Using OpenMM platform ', A)", platformName
! Set starting positions of the atoms. Leave time and velocity zero.
call OpenMM_Context_setPositions(context, initPosInNm)
! Simulate.
frameNum = 1
do
! Output current state information.
call OpenMM_Context_getState(context, OpenMM_State_Positions, OpenMM_False, state)
timeInPs = OpenMM_State_getTime(state)
call writePdbFrame(frameNum, state) !output coordinates
call OpenMM_State_destroy(state)
if (timeInPs .ge. 10.) then
exit
end if
! Advance state many steps at a time, for efficient use of OpenMM.
! (use a lot more than 10 normally)
call OpenMM_VerletIntegrator_step(verlet, 10)
frameNum = frameNum + 1
end do
! Free heap space for all the objects created above.
call OpenMM_Vec3Array_destroy(initPosInNm)
call OpenMM_Context_destroy(context)
call OpenMM_VerletIntegrator_destroy(verlet)
call OpenMM_System_destroy(system)
END PROGRAM
! Handy homebrew PDB writer for quick-and-dirty trajectory output.
SUBROUTINE writePDBFrame(frameNum, state)
use OpenMM; implicit none
integer frameNum
type(OpenMM_State) state
type(OpenMM_Vec3Array) allPosInNm
real*8 posInNm(3), posInAng(3)
integer n
! Reference atomic positions in the OpenMM State.
call OpenMM_State_getPositions(state, allPosInNm)
print "('MODEL',5X,I0)", frameNum ! start of frame
do n = 1,OpenMM_Vec3Array_getSize(allPosInNm)
call OpenMM_Vec3Array_get(allPosInNm, n, posInNm)
call OpenMM_Vec3_scale(posInNm, 10d0, posInAng)
print "('ATOM ', I5, ' AR AR 1 ', 3F8.3, ' 1.00 0.00')", &
n, posInAng
end do
print "('ENDMDL')"
END SUBROUTINE
|
