"platforms/cuda/vscode:/vscode.git/clone" did not exist on "9cd18aeb5140adc10f6309d938b11f61078a13f5"
HelloArgon.cpp 3.39 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
// -----------------------------------------------------------------------------
//           OpenMM(tm) HelloArgon example in C++ (June 2009)
// -----------------------------------------------------------------------------
// This program demonstrates a simple molecular simulation using the OpenMM
// API for GPU-accelerated molecular dynamics simulation.  The system
// represents a small collection of argon atoms in a vaccuum.
// 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.
// -----------------------------------------------------------------------------

11
12
13
#include "OpenMM.h"
#include <cstdio>

Christopher Bruns's avatar
Christopher Bruns committed
14
15
// Forward declaration of writePdb() subroutine for printing atomic 
// coordinates, defined later in this source file.
16
void writePdb(const OpenMM::OpenMMContext& context);
17

18
void simulateArgon()
19
20
{
    // Load any shared libraries containing GPU implementations
21
22
    OpenMM::Platform::loadPluginsFromDirectory(
        OpenMM::Platform::getDefaultPluginsDirectory());
23

24
    // Create a system with nonbonded forces
25
26
    OpenMM::System system;
    OpenMM::NonbondedForce* nonbond = new OpenMM::NonbondedForce(); 
27
28
    system.addForce(nonbond);

29
    // Create three atoms
30
    std::vector<OpenMM::Vec3> initialPositions(3);
31
32
    for (int a = 0; a < 3; ++a) 
    {
Christopher Bruns's avatar
Christopher Bruns committed
33
        system.addParticle(39.95); // mass, grams per mole
34
35

        // charge, sigma, well depth
36
        nonbond->addParticle(0.0, 0.3350, 0.996); 
37

Christopher Bruns's avatar
Christopher Bruns committed
38
        initialPositions[a] = OpenMM::Vec3(0.5*a,0,0); // location, nanometers
39
40
    }

41
    OpenMM::VerletIntegrator integrator(0.004); // step size in picoseconds
42

43
    // Let OpenMM Context choose best platform.
44
    OpenMM::OpenMMContext context(system, integrator);
45
46
    printf( "REMARK  Using OpenMM platform %s\n", 
        context.getPlatform().getName().c_str() );
47

48
49
    // Set the starting positions of the atoms. Velocities will be zero.
    context.setPositions(initialPositions);
50

51
    // Simulate
52
    while(context.getTime() < 10.0) { // picoseconds
53
54
55
56
57
58
59
60
61
62
63
64
        writePdb(context); // output coordinates
        // Run 100 steps at a time, for efficient use of OpenMM
        integrator.step(100);
    }
    writePdb(context); // output final coordinates
}

int main() 
{
    try {
        simulateArgon();
        return 0; // success!
65
    }
66
67
68
69
70
71
72
73
    // Catch and report usage and runtime errors detected by OpenMM and fail.
    catch(const std::exception& e) {
        printf("EXCEPTION: %s\n", e.what());
        return 1; // failure!
    }
}

// writePdb() subroutine for quick-and-dirty trajectory output.
74
void writePdb(const OpenMM::OpenMMContext& context) 
75
76
{
    // Request atomic positions from OpenMM
77
78
    const OpenMM::State state = context.getState(OpenMM::State::Positions);
    const std::vector<OpenMM::Vec3>& pos = state.getPositions();
79

80
81
82
83
84
85
86
87
88
    // write out in PDB format

    // Use PDB MODEL cards to number trajectory frames
    static int modelFrameNumber = 0; 
    modelFrameNumber++;
    printf("MODEL     %d\n", modelFrameNumber); // start of frame
    for (int a = 0; a < context.getSystem().getNumParticles(); ++a)
    {
        printf("ATOM  %5d AR    AR     1    ", a+1); // atom number
89
        printf("%8.3f%8.3f%8.3f  1.00  0.00          AR\n", // coordinates
Christopher Bruns's avatar
Christopher Bruns committed
90
            // "*10" converts nanometers to Angstroms
91
92
93
            pos[a][0]*10, pos[a][1]*10, pos[a][2]*10);
    }
    printf("ENDMDL\n"); // end of frame
94
}