Skip to content

Instrumenting CXX Applications

Kevin Huck edited this page Oct 28, 2019 · 6 revisions

This page is an overview of how to use TAU to instrument a C++ application that uses CMake as a configuration/build tool. THe application shown in this example is OpenMM, an open-source molecular dynamics application.

The most important thing to understand about instrumenting C++ applications is that without some selective instrumentation and/or filtering of timers, the application can quickly be overwhelmed with timer calls due to the high number of getter/setter methods that are common in C++ applications. When blindly instrumenting an entire application (especially with compiler-based instrumentation), all of the member functions and template expansions of some standard library objects will also be instrumented, introducing several orders of magnitude of overhead. To prevent that, we use selective instrumentation and the PDT infrastructure for source-to-source instrumentation.

Installing TAU

More detailed instructions for installing TAU are located here.

We will configure PDT with:

wget http://tau.uoregon.edu/pdt.tgz
tar -xzf pdt.tgz
cd pdtoolkit-3.25.1
# Perform in-place installation
./configure -GNU
make -j install

One particular problem with OpenMM is it's use of transient threads. The source code will detect the number of cores on the node, and launch that many pthreads (despite what the documentation says about the OPENMM_CPU_THREADS variable). We will configure TAU with thread recycling (an experimental feature), because for some applciations OpenMM will launch N threads, destroy them, and re-launch them (possibly many times). Currently, TAU has a hard-coded thread count limit of 128-512 (depending on the platform). In this particular example, we will configure TAU with (your path to PDT may vary):

./configure -useropt=-DTAU_RECYCLE_THREADS#-O2#-g -pdt=/home/users/khuck/src/pdtoolkit-3.25.1 -pthread -bfd=download -unwind=download -otf=download

The -useropt=-DTAU_RECYCLE_THREADS#-O2#-g argument will replace the default -O2 -g compiler arguments to build TAU (the # character is used to represent white space, and will be converted by the configuration process).

Building the application

As mentioned earlier, for this example we are building OpenMM. To get the source:

git clone https://github.com/openmm/openmm.git
cd openmm

On this example test system, we are building on an IBM POWER9 node with 40 cores (160 HW threads) and NVIDIA V100 GPUs. These are the modules and environment settings required to build (OpenMM dependencies):

module purge

module load cuda/9.2
module load fftw/3.3.8-openmpi4-gcc6.3
module load openmpi/4.0.1-gcc6.3
module load gcc/6.3
module load cmake/3.15.4

module list

# Add Swig 4 to the path
PATH=/usr/local/packages/swig/4.0.1/bin:$PATH

# Set the TAU settings
export TAU_MAKEFILE=/home/users/khuck/src/tau2/ibm64linux/lib/Makefile.tau-gnu-pthread-pdt
export TAU_OPTIONS="-optQuiet -optRevert -optNoCompInst -optShared -optSelectFile=/home/users/khuck/src/openmm/select.tau"

The meanings of the TAU options are:

  • -optQuiet - parse/instrument/compile silently
  • -optRevert - if the parser fails, try other options, then just compile if all fail
  • -optNoCompInst - don't use compiler-based instrumentation if PDT instrumentation fails
  • -optShared - link with the libTAUsh.so shared library, not libTAU.a
  • -optSelectFile - use this select file to specify which files/functions to instrument

The contents of select.tau are (to instrument only the main function and any functions in the OpenMM namespace):

BEGIN_INCLUDE_LIST
"int main#"
"#OpenMM#"
END_INCLUDE_LIST

These are the settings required to build. However, the CMake configuration for OpenMM uses the compiler to auto-detect hardware and OS configuration settings - including some tests designed to fail. Because we are using the TAU compiler wrappers to instrument the application, we need to tell CMake about the compiler wrapper but disable TAU instrumentation - otherwise the TAU error output will interfere with the auto-detection process. To do that we change the TAU_OPTIONS environment variable for cmake, then change it back:

mkdir build-tau
cd build-tau

# save the TAU options
OLD_TAU_OPTIONS=${TAU_OPTIONS}
TAU_OPTIONS="-optQuiet -optLinkOnly -optRevert -optNoCompInst"

cmake -G "Unix Makefiles" .. \
    -DCMAKE_C_COMPILER=`which tau_cc.sh` \
    -DCMAKE_CXX_COMPILER=`which tau_cxx.sh` \
    -DCMAKE_BUILD_TYPE=RelWithDebInfo

# Restore the TAU options for building
TAU_OPTIONS=${OLD_TAU_OPTIONS}
make -j40

The current (as of October 24, 2019) C++ parser used by PDT cannot support some language features, so several compiler errors will be seen when building, but the -optRevert option will tell tau_cxx.sh to fall back to different settings until either one source parser succeeds or the file is left uninstrumented.

After building, you can validate the build with make test in the build-tau directory (although some tests take a long time to run).

Running an individual test (for example, ./TestReferenceCustomIntegrator in the build directory will write a bunch of profile.* files (one per concurrent thread & process), and they can be summarized on the command line with the pprof TAU program:

$ pprof .
FUNCTION SUMMARY (mean):
---------------------------------------------------------------------------------------
%Time    Exclusive    Inclusive       #Call      #Subrs  Inclusive Name
              msec   total msec                          usec/call 
---------------------------------------------------------------------------------------
 68.9      0.00273           74     1.98758     1.98758      37424 OpenMM::threadBody(void*) 
 68.9      0.00761           74     1.98758     4.19876      37422 void *OpenMM::threadBody(void *) 
 68.9           74           74     3.97516           0      18709 void OpenMM::ThreadPool::syncThreads() 
 54.4           18           58           1     9466.81      58735 .TAU application
  8.8            8            9     499.932     5198.64         19 double OpenMM::NonbondedForceImpl::calcForcesAndEnergy(OpenMM::ContextImpl &, bool, bool, int) 
  1.9            1            2     621.124     2981.39          3 double OpenMM::HarmonicBondForceImpl::calcForcesAndEnergy(OpenMM::ContextImpl &, bool, bool, int) [THROTTLED]
  1.2      0.00284            1    0.136646     4.40373       9748 void OpenMM::ReferenceConstraints::ReferenceConstraints(const OpenMM::System &) 
  0.9        0.425        0.966   0.0124224     1.98758      77779 void OpenMM::ThreadPool::ThreadPool(int) 
  0.5        0.592        0.592     372.671  0.00621118          2 void OpenMM::ReferenceAndersenThermostat::applyThermostat(const std::vector<std::vector<int, std::allocator<int>>, std::allocator<std::vector<int, std::allocator<int>>>> &, std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &, std::vector<double, std::allocator<double>> &, double, double, double) const 
  0.5        0.541        0.541     1.98758           0        272 pthread_create
  0.5       0.0738          0.5    0.310559     5.42236       1611 void OpenMM::MonteCarloBarostatImpl::updateContextState(OpenMM::ContextImpl &, bool &) 
  0.3        0.272         0.35     621.124     347.329          1 double OpenMM::SimTKOpenMMUtilities::getNormallyDistributedRandomNumber() [THROTTLED]
  0.3        0.225        0.286     621.124     273.658          0 double OpenMM::SimTKOpenMMUtilities::getUniformlyDistributedRandomNumber() [THROTTLED]
  0.2        0.187        0.253     621.124     186.335          0 void OpenMM::ReferenceBondForce::calculateForce(int, std::vector<std::vector<int, std::allocator<int>>, std::allocator<std::vector<int, std::allocator<int>>>> &, std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &, std::vector<std::vector<double, std::allocator<double>>, std::allocator<std::vector<double, std::allocator<double>>>> &, std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &, double *, OpenMM::ReferenceBondIxn &) [THROTTLED]
  0.2      0.00131        0.206   0.0124224     1.98758      16583 void OpenMM::ThreadPool::~ThreadPool() 
  0.2        0.205        0.205     1.98758           0        103 pthread_join
  0.2          0.2        0.203     621.124     12.4224          0 void OpenMM::ReferenceHarmonicBondIxn::calculateBondIxn(std::vector<int, std::allocator<int>> &, std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &, std::vector<double, std::allocator<double>> &, std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &, double *, double *) [THROTTLED]
  0.2        0.172        0.172     621.124           0          0 void OpenMM::ReferenceVirtualSites::computePositions(const OpenMM::System &, std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &) [THROTTLED]
  0.1        0.161        0.161     621.124           0          0 void OpenMM::ReferenceForce::getDeltaRPeriodic(const OpenMM::Vec3 &, const OpenMM::Vec3 &, const OpenMM::Vec3 *, double *) [THROTTLED]
  0.1        0.159        0.159     621.124           0          0 void OpenMM::ReferenceForce::ReferenceForce() [THROTTLED]
  0.1        0.158        0.158     621.124           0          0 double OpenMM::Integrator::getStepSize() const [THROTTLED]
  0.1        0.155        0.155     621.124           0          0 void OpenMM::ReferenceHarmonicBondIxn::ReferenceHarmonicBondIxn() [THROTTLED]
  0.1        0.155        0.155   0.0124224           0      12506 void OpenMM::ThreadPool::waitForThreads() 
  0.1        0.151        0.151     621.124           0          0 void OpenMM::ReferenceDynamics::setReferenceConstraintAlgorithm(OpenMM::ReferenceConstraintAlgorithm *) [THROTTLED]
  0.1        0.149        0.149     621.124           0          0 void OpenMM::Integrator::setStepSize(double) [THROTTLED]
  0.1        0.143        0.143     621.124           0          0 void OpenMM::ReferenceBondForce::~ReferenceBondForce() [THROTTLED]
  0.1        0.142        0.142     621.124           0          0 int OpenMM::ReferenceDynamics::incrementTimeStep() [THROTTLED]
  0.1        0.141        0.141     621.124           0          0 void OpenMM::ReferenceBondIxn::~ReferenceBondIxn() [THROTTLED]
  0.1         0.14         0.14     621.124           0          0 void OpenMM::ReferenceForce::getDeltaR(const OpenMM::Vec3 &, const OpenMM::Vec3 &, double *) [THROTTLED]
  0.1         0.14         0.14     621.124           0          0 double OpenMM::Integrator::getConstraintTolerance() const [THROTTLED]
  0.1        0.138        0.138     621.124           0          0 uint32_t OpenMM_SFMT::gen_rand32(OpenMM_SFMT::SFMT &) [THROTTLED]
  0.1        0.138        0.138     621.124           0          0 void OpenMM::ReferenceBondForce::ReferenceBondForce() [THROTTLED]
  0.1        0.138        0.138     621.124           0          0 const std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &OpenMM::State::getPositions() const [THROTTLED]
  0.1        0.137        0.137     621.124           0          0 double OpenMM::ReferenceDynamics::getDeltaT() const [THROTTLED]
  0.1        0.137        0.137     621.124           0          0 void OpenMM::ReferenceForce::~ReferenceForce() [THROTTLED]
  0.1        0.137        0.137     621.124           0          0 void OpenMM::ReferenceHarmonicBondIxn::~ReferenceHarmonicBondIxn() [THROTTLED]
  0.1        0.136        0.136     621.124           0          0 void OpenMM::ReferenceBondIxn::ReferenceBondIxn() [THROTTLED]
  0.1        0.122        0.122     512.255           0          0 void OpenMM::ReferenceVirtualSites::distributeForces(const OpenMM::System &, const std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &, std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &) 
  0.1        0.114        0.114     509.224           0          0 const std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &OpenMM::State::getVelocities() const 
  0.1        0.104        0.104     462.646           0          0 void OpenMM::ReferenceLJCoulomb14::ReferenceLJCoulomb14() 
  0.1        0.103        0.103     462.646           0          0 void OpenMM::ReferenceLJCoulomb14::~ReferenceLJCoulomb14() 
  0.1       0.0422       0.0627     89.6087     89.6087          1 void OpenMM::State::StateBuilder::setPeriodicBoxVectors(const OpenMM::Vec3 &, const OpenMM::Vec3 &, const OpenMM::Vec3 &) 
  0.0       0.0276       0.0424     57.2857     57.2857          1 void OpenMM::State::StateBuilder::setEnergy(double, double) 
  0.0       0.0244       0.0409     50.3168     50.3168          1 void OpenMM::State::StateBuilder::setPositions(const std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &) 
  0.0       0.0265       0.0265     89.6087           0          0 OpenMM::State OpenMM::State::StateBuilder::getState() 
  0.0       0.0208       0.0249     12.4348     12.4224          2 void OpenMM::ReferenceConstraints::applyToVelocities(std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &, std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &, std::vector<double, std::allocator<double>> &, double) 
  0.0       0.0217       0.0217     89.6087           0          0 void OpenMM::State::State(double) 
  0.0        0.016       0.0215     18.7019     12.4224          1 void OpenMM::ReferenceConstraints::apply(std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &, std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &, std::vector<double, std::allocator<double>> &, double) 
  0.0       0.0205       0.0205     89.6087           0          0 void OpenMM::State::setPeriodicBoxVectors(const OpenMM::Vec3 &, const OpenMM::Vec3 &, const OpenMM::Vec3 &) 
  0.0       0.0204       0.0204     93.3727           0          0 void OpenMM::System::getConstraintParameters(int, int &, int &, double &) const 
  0.0       0.0201       0.0201     89.6087           0          0 void OpenMM::State::StateBuilder::StateBuilder(double) 
  0.0       0.0165       0.0165     50.3168           0          0 void OpenMM::State::setPositions(const std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &) 
  0.0       0.0147       0.0147     57.2857           0          0 void OpenMM::State::setEnergy(double, double) 
  0.0       0.0131       0.0131      57.677           0          0 double OpenMM::State::getKineticEnergy() const 
  0.0      0.00654       0.0114     13.0497     13.0497          1 void OpenMM::State::StateBuilder::setVelocities(const std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &) 
  0.0      0.00804      0.00804     31.1242           0          0 OpenMM::ReferenceConstraintAlgorithm *OpenMM::ReferenceDynamics::getReferenceConstraintAlgorithm() const 
  0.0      0.00565      0.00565     25.6398           0          0 double OpenMM::State::getPotentialEnergy() const 
  0.0      0.00546      0.00546     12.4224           0          0 void OpenMM::ReferenceSETTLEAlgorithm::apply(std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &, std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &, std::vector<double, std::allocator<double>> &, double) 
  0.0      0.00486      0.00486     13.0497           0          0 void OpenMM::State::setVelocities(const std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &) 
  0.0      0.00472      0.00476  0.00621118   0.0807453        767 void OpenMM::CustomBondForceImpl::initialize(OpenMM::ContextImpl &) 
  0.0      0.00468      0.00472  0.00621118   0.0807453        760 void OpenMM::CustomAngleForceImpl::initialize(OpenMM::ContextImpl &) 
  0.0      0.00436       0.0044     1.49068    0.136646          3 OpenMM::Kernel &OpenMM::Kernel::operator=(const OpenMM::Kernel &) 
  0.0      0.00405      0.00405     12.4224           0          0 void OpenMM::ReferenceSETTLEAlgorithm::applyToVelocities(std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &, std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &, std::vector<double, std::allocator<double>> &, double) 
  0.0      0.00285      0.00285     12.8137           0          0 double OpenMM::System::getParticleMass(int) const 
  0.0       0.0016       0.0016    0.310559           0          5 void OpenMM::ReferenceMonteCarloBarostat::applyBarostat(std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &, const OpenMM::Vec3 *, double, double, double) 
  0.0      0.00155      0.00155    0.149068           0         10 void OpenMM::SimTKOpenMMUtilities::setRandomNumberSeed(uint32_t) 
  0.0     0.000981      0.00129   0.0559006    0.298137         23 void OpenMM::NonbondedForceImpl::initialize(OpenMM::ContextImpl &) 
  0.0      0.00126      0.00127     2.32298   0.0186335          1 void OpenMM::Kernel::~Kernel() 
  0.0      0.00127      0.00127     6.21118           0          0 double OpenMM::State::getTime() const 
  0.0      0.00111      0.00117   0.0124224    0.236025         94 void OpenMM::CustomExternalForceImpl::initialize(OpenMM::ContextImpl &) 
  0.0      0.00111      0.00115    0.198758    0.198758          6 void OpenMM_SFMT::init_gen_rand(uint32_t, OpenMM_SFMT::SFMT &) 
  0.0     0.000863      0.00111    0.248447    0.745342          4 void OpenMM::SerializationProxy::registerProxy(const std::type_info &, const OpenMM::SerializationProxy *) 
  0.0     0.000155      0.00054  0.00621118  0.00621118         87 Lepton::CustomFunction *createReferenceTabulatedFunction(const OpenMM::TabulatedFunction &) C 
  0.0      0.00046     0.000516   0.0248447    0.180124         21 void OpenMM::HarmonicBondForceImpl::initialize(OpenMM::ContextImpl &) 
  0.0     0.000242     0.000466    0.223602    0.223602          2 Lepton::CustomFunction *OpenMM::ReferenceContinuous1DFunction::clone() const 
  0.0      0.00018     0.000447   0.0124224   0.0124224         36 void OpenMM::ThreadPool::execute(std::function<void (OpenMM::ThreadPool &, int)>) 
  0.0     0.000255     0.000385  0.00621118   0.0124224         62 void OpenMM::ReferenceContinuous1DFunction::ReferenceContinuous1DFunction(const OpenMM::Continuous1DFunction &) 
  0.0     0.000342     0.000342     1.51553           0          0 void OpenMM::Kernel::Kernel() 
  0.0     0.000248     0.000286  0.00621118   0.0372671         46 void OpenMM::MonteCarloBarostatImpl::initialize(OpenMM::ContextImpl &) 
  0.0      0.00028      0.00028    0.807453           0          0 void OpenMM::KernelImpl::KernelImpl(std::string, const OpenMM::Platform &) 
  0.0     0.000255     0.000267    0.136646   0.0372671          2 void OpenMM::ReferenceConstraints::~ReferenceConstraints() 
  0.0     0.000267     0.000267   0.0124224           0         22 void OpenMM::ThreadPool::resumeThreads() 
  0.0     0.000217     0.000236   0.0124224   0.0745342         19 int QUERN_compute_qr(int, int, const int *, const int *, const double *, const int *, int **, int **, double **, int **, int **, double **) C 
  0.0     0.000224      0.00023   0.0186335   0.0372671         12 double OpenMM::CustomAngleForceImpl::calcForcesAndEnergy(OpenMM::ContextImpl &, bool, bool, int) 
  0.0     0.000143     0.000224    0.223602    0.223602          1 void OpenMM::ReferenceContinuous1DFunction::ReferenceContinuous1DFunction(const OpenMM::ReferenceContinuous1DFunction &) 
  0.0     0.000186     0.000211   0.0248447    0.149068          9 double OpenMM::CustomExternalForceImpl::calcForcesAndEnergy(OpenMM::ContextImpl &, bool, bool, int) 
  0.0      0.00018     0.000193   0.0248447   0.0496894          8 OpenMM::ForceImpl *OpenMM::HarmonicBondForce::createImpl() const 
  0.0     0.000174     0.000193   0.0124224   0.0745342         16 std::map<std::string, double, std::less<std::string>, std::allocator<std::pair<const std::string, double>>> OpenMM::CustomAngleForceImpl::getDefaultParameters() 
  0.0     0.000186     0.000186    0.807453           0          0 void OpenMM::Kernel::Kernel(OpenMM::KernelImpl *) 
  0.0     0.000168      0.00018   0.0124224   0.0745342         15 std::map<std::string, double, std::less<std::string>, std::allocator<std::pair<const std::string, double>>> OpenMM::CustomBondForceImpl::getDefaultParameters() 
  0.0     0.000155     0.000161  0.00621118   0.0248447         26 void OpenMM::NonbondedForceImpl::calcEwaldParameters(const OpenMM::System &, const OpenMM::NonbondedForce &, double &, int &, int &, int &) 
  0.0     0.000143     0.000155   0.0372671   0.0372671          4 void OpenMM::State::StateBuilder::setForces(const std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &) 
  0.0      0.00013     0.000137  0.00621118   0.0124224         22 OpenMM::ForceImpl *OpenMM::CustomAngleForce::createImpl() const 
  0.0     0.000124     0.000137  0.00621118   0.0186335         22 OpenMM::ForceImpl *OpenMM::MonteCarloBarostat::createImpl() const 
  0.0     0.000124      0.00013   0.0124224   0.0248447         10 OpenMM::ForceImpl *OpenMM::CustomExternalForce::createImpl() const 
  0.0      0.00013      0.00013   0.0745342           0          2 OpenMM_SFMT::SFMTData *OpenMM_SFMT::createSFMTData() 
  0.0     0.000118      0.00013   0.0248447   0.0248447          5 double OpenMM::ReferenceContinuous1DFunction::evaluate(const double *) const 
  0.0     0.000124      0.00013  0.00621118  0.00621118         21 void OpenMM::SplineFitter::createNaturalSpline(const std::vector<double, std::allocator<double>> &, const std::vector<double, std::allocator<double>> &, std::vector<double, std::allocator<double>> &) 
  0.0     0.000118     0.000124  0.00621118   0.0124224         20 OpenMM::ForceImpl *OpenMM::CustomBondForce::createImpl() const 
  0.0     0.000106     0.000106    0.621118           0          0 double OpenMM::State::getPeriodicBoxVolume() const 
  0.0     9.94E-05     9.94E-05   0.0186335           0          5 double OpenMM::CustomBondForceImpl::calcForcesAndEnergy(OpenMM::ContextImpl &, bool, bool, int) 
  0.0     9.32E-05     9.32E-05    0.248447           0          0 const std::string &OpenMM::SerializationProxy::getTypeName() const 
  0.0     8.07E-05     8.07E-05    0.248447           0          0 std::map<const std::string, const OpenMM::SerializationProxy *, std::less<const std::string>, std::allocator<std::pair<const std::string, const OpenMM::SerializationProxy *>>> &OpenMM::SerializationProxy::getProxiesByName() 
  0.0     8.07E-05     8.07E-05    0.229814           0          0 void OpenMM::Continuous1DFunction::getFunctionParameters(std::vector<double, std::allocator<double>> &, double &, double &) const 
  0.0     8.07E-05     8.07E-05    0.136646           0          1 void OpenMM::ReferenceMonteCarloBarostat::restorePositions(std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &) 
  0.0     7.45E-05     7.45E-05    0.248447           0          0 std::map<const std::string, const OpenMM::SerializationProxy *, std::less<const std::string>, std::allocator<std::pair<const std::string, const OpenMM::SerializationProxy *>>> &OpenMM::SerializationProxy::getProxiesByType() 
  0.0     6.83E-05     6.83E-05   0.0559006           0          1 double OpenMM::NonbondedForceImpl::calcDispersionCorrection(const OpenMM::System &, const OpenMM::NonbondedForce &) 
  0.0     6.21E-05     6.21E-05    0.310559           0          0 int OpenMM::ReferenceContinuous1DFunction::getNumArguments() const 
  0.0     6.21E-05     6.21E-05    0.136646           0          0 void OpenMM::System::~System() 
  0.0     5.59E-05     5.59E-05   0.0993789           0          1 int OpenMM::System::addConstraint(int, int, double) 
  0.0     5.59E-05     5.59E-05   0.0745342           0          1 int QUERN_multiply_with_q_transpose(int, const int *, const int *, const double *, double *) C 
  0.0     4.97E-05     4.97E-05    0.136646           0          0 const OpenMM::Force &OpenMM::System::getForce(int) const 
  0.0     4.97E-05     4.97E-05    0.111801           0          0 std::map<std::string, double, std::less<std::string>, std::allocator<std::pair<const std::string, double>>> OpenMM::NonbondedForceImpl::getDefaultParameters() 
  0.0     4.97E-05     4.97E-05    0.136646           0          0 void OpenMM::System::System() 
  0.0     4.97E-05     4.97E-05     0.15528           0          0 void OpenMM::System::getDefaultPeriodicBoxVectors(OpenMM::Vec3 &, OpenMM::Vec3 &, OpenMM::Vec3 &) const 
  0.0     4.35E-05     4.35E-05   0.0745342           0          1 int QUERN_solve_with_r(int, const int *, const int *, const double *, const double *, double *) C 
  0.0     4.35E-05     4.35E-05    0.149068           0          0 void OpenMM::Integrator::setConstraintTolerance(double) 
  0.0     4.35E-05     4.35E-05    0.136646           0          0 void OpenMM::ReferenceDynamics::ReferenceDynamics(int, double, double) 
  0.0     4.35E-05     4.35E-05    0.136646           0          0 void OpenMM::ReferenceDynamics::~ReferenceDynamics() 
  0.0     3.73E-05     3.73E-05   0.0559006           0          1 int OpenMM::HarmonicBondForce::addBond(int, int, double, double) 
  0.0     3.73E-05     3.73E-05    0.136646           0          0 void OpenMM::Integrator::Integrator() 
  0.0     3.73E-05     3.73E-05    0.248447           0          0 void OpenMM::SerializationProxy::SerializationProxy(const std::string &) 
  0.0     3.73E-05     3.73E-05    0.198758           0          0 void OpenMM_SFMT::period_certification(OpenMM_SFMT::SFMT &) 
  0.0     3.11E-05     3.11E-05    0.149068           0          0 double *OpenMM::ReferenceForce::getVariablePointer(Lepton::CompiledExpression &, const std::string &) 
  0.0     3.11E-05     3.11E-05   0.0745342           0          0 int OpenMM::ThreadPool::getNumThreads() const 
  0.0     3.11E-05     3.11E-05   0.0559006           0          1 std::vector<std::string, std::allocator<std::string>> OpenMM::NonbondedForceImpl::getKernelNames() 
  0.0     3.11E-05     3.11E-05    0.068323           0          0 void OpenMM_SFMT::deleteSFMTData(OpenMM_SFMT::SFMTData *) 
  0.0     2.48E-05     2.48E-05   0.0559006           0          0 const std::string &OpenMM::CustomAngleForce::getGlobalParameterName(int) const 
  0.0     2.48E-05     2.48E-05   0.0559006           0          0 void OpenMM::NonbondedForceImpl::~NonbondedForceImpl() 
  0.0     2.48E-05     2.48E-05    0.149068           0          0 void OpenMM::ReferenceForce::setVariable(double *, double) 
  0.0     6.21E-06     1.86E-05   0.0124224   0.0124224          2 OpenMM::ForceImpl *OpenMM::AndersenThermostat::createImpl() const 
  0.0     1.86E-05     1.86E-05   0.0745342           0          0 bool copy_row(int, const int *, const double *, SparseVector &) 
  0.0     1.86E-05     1.86E-05   0.0372671           0          0 double *OpenMM::SimTKOpenMMUtilities::allocateOneDRealOpenMMArray(int, double *, int, double, const std::string &) 
  0.0     1.86E-05     1.86E-05   0.0124224           0          2 int OpenMM::CustomExternalForce::addParticle(int, const std::vector<double, std::allocator<double>> &) 
  0.0     1.86E-05     1.86E-05   0.0124224           0          2 std::map<std::string, double, std::less<std::string>, std::allocator<std::pair<const std::string, double>>> OpenMM::MonteCarloBarostatImpl::getDefaultParameters() 
  0.0     1.86E-05     1.86E-05   0.0248447           0          1 void QUERN_free_result(int *, int *, double *) C 
  0.0     1.24E-05     1.24E-05   0.0248447           0          0 bool OpenMM::HarmonicBondForce::usesPeriodicBoundaryConditions() const 
  0.0     1.24E-05     1.24E-05   0.0372671           0          0 double OpenMM::CustomBondForce::getGlobalParameterDefaultValue(int) const 
  0.0     1.24E-05     1.24E-05   0.0248447           0          0 double OpenMM::SplineFitter::evaluateSpline(const std::vector<double, std::allocator<double>> &, const std::vector<double, std::allocator<double>> &, const std::vector<double, std::allocator<double>> &, double) 
  0.0     1.24E-05     1.24E-05   0.0186335           0          1 int OpenMM::CustomAngleForce::addGlobalParameter(const std::string &, double) 
  0.0     1.24E-05     1.24E-05  0.00621118           0          2 int OpenMM::CustomBondForce::addBond(int, int, const std::vector<double, std::allocator<double>> &) 
  0.0     1.24E-05     1.24E-05   0.0186335           0          1 int OpenMM::CustomBondForce::addGlobalParameter(const std::string &, double) 
  0.0     1.24E-05     1.24E-05   0.0248447           0          0 std::map<std::string, double, std::less<std::string>, std::allocator<std::pair<const std::string, double>>> OpenMM::CustomExternalForceImpl::getDefaultParameters() 
  0.0     1.24E-05     1.24E-05   0.0248447           0          0 std::vector<std::string, std::allocator<std::string>> OpenMM::HarmonicBondForceImpl::getKernelNames() 
  0.0     1.24E-05     1.24E-05   0.0559006           0          0 void OpenMM::HarmonicBondForce::getBondParameters(int, int &, int &, double &, double &) const 
  0.0     1.24E-05     1.24E-05    0.136646           0          0 void OpenMM::Integrator::~Integrator() 
  0.0     1.24E-05     1.24E-05   0.0372671           0          0 void OpenMM::SimTKOpenMMUtilities::freeOneDRealOpenMMArray(double *, const std::string &) 
  0.0     1.24E-05     1.24E-05   0.0372671           0          0 void OpenMM::State::setForces(const std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &) 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 bool OpenMM::CustomBondForce::usesPeriodicBoundaryConditions() const 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 bool OpenMM::System::usesPeriodicBoundaryConditions() const 
  0.0     6.21E-06     6.21E-06   0.0124224           0          0 const std::string &OpenMM::CustomAngleForce::getEnergyParameterDerivativeName(int) const 
  0.0     6.21E-06     6.21E-06   0.0124224           0          0 const std::string &OpenMM::CustomBondForce::getEnergyParameterDerivativeName(int) const 
  0.0     6.21E-06     6.21E-06   0.0559006           0          0 const std::string &OpenMM::CustomBondForce::getGlobalParameterName(int) const 
  0.0     6.21E-06     6.21E-06   0.0124224           0          0 const std::string &OpenMM::CustomExternalForce::getEnergyFunction() const 
  0.0     6.21E-06     6.21E-06   0.0372671           0          0 double OpenMM::CustomAngleForce::getGlobalParameterDefaultValue(int) const 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 int OpenMM::CustomAngleForce::addAngle(int, int, int, const std::vector<double, std::allocator<double>> &) 
  0.0     6.21E-06     6.21E-06   0.0186335           0          0 int OpenMM::NonbondedForceImpl::findZero(const OpenMM::NonbondedForceImpl::ErrorFunction &, int) 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 std::vector<std::string, std::allocator<std::string>> OpenMM::CustomAngleForceImpl::getKernelNames() 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 std::vector<std::string, std::allocator<std::string>> OpenMM::CustomBondForceImpl::getKernelNames() 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 std::vector<std::string, std::allocator<std::string>> OpenMM::MonteCarloBarostatImpl::getKernelNames() 
  0.0     6.21E-06     6.21E-06   0.0124224           0          0 void OpenMM::AndersenThermostat::AndersenThermostat(double, double) 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 void OpenMM::AndersenThermostatProxy::AndersenThermostatProxy() 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 void OpenMM::Continuous1DFunction::Continuous1DFunction(const std::vector<double, std::allocator<double>> &, double, double) 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 void OpenMM::CustomAngleForce::CustomAngleForce(const std::string &) 
  0.0     6.21E-06     6.21E-06   0.0124224           0          0 void OpenMM::CustomAngleForce::getAngleParameters(int, int &, int &, int &, std::vector<double, std::allocator<double>> &) const 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 void OpenMM::CustomAngleForceImpl::CustomAngleForceImpl(const OpenMM::CustomAngleForce &) 
  0.0     6.21E-06     6.21E-06   0.0124224           0          0 void OpenMM::CustomBondForce::getBondParameters(int, int &, int &, std::vector<double, std::allocator<double>> &) const 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 void OpenMM::CustomBondForceImpl::CustomBondForceImpl(const OpenMM::CustomBondForce &) 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 void OpenMM::CustomBondForceImpl::~CustomBondForceImpl() 
  0.0     6.21E-06     6.21E-06   0.0248447           0          0 void OpenMM::HarmonicBondForce::HarmonicBondForce() 
  0.0     6.21E-06     6.21E-06   0.0248447           0          0 void OpenMM::HarmonicBondForceImpl::HarmonicBondForceImpl(const OpenMM::HarmonicBondForce &) 
  0.0     6.21E-06     6.21E-06   0.0248447           0          0 void OpenMM::HarmonicBondForceImpl::~HarmonicBondForceImpl() 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 void OpenMM::MonteCarloAnisotropicBarostatProxy::MonteCarloAnisotropicBarostatProxy() 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 void OpenMM::MonteCarloBarostatImpl::MonteCarloBarostatImpl(const OpenMM::MonteCarloBarostat &) 
  0.0     6.21E-06     6.21E-06   0.0124224           0          0 void OpenMM::ReferenceAndersenThermostat::~ReferenceAndersenThermostat() 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 void OpenMM::ReferenceMonteCarloBarostat::ReferenceMonteCarloBarostat(int, const std::vector<std::vector<int, std::allocator<int>>, std::allocator<std::vector<int, std::allocator<int>>>> &) 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 void OpenMM::ReferenceMonteCarloBarostat::~ReferenceMonteCarloBarostat() 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 void OpenMM::ReferenceSETTLEAlgorithm::ReferenceSETTLEAlgorithm(const std::vector<int, std::allocator<int>> &, const std::vector<int, std::allocator<int>> &, const std::vector<int, std::allocator<int>> &, const std::vector<double, std::allocator<double>> &, const std::vector<double, std::allocator<double>> &, std::vector<double, std::allocator<double>> &) 
  0.0     6.21E-06     6.21E-06   0.0372671           0          0 void OpenMM::SimTKOpenMMUtilities::crossProductVector3(double *, double *, double *) 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 void OpenMM::SplineFitter::solveTridiagonalMatrix(const std::vector<double, std::allocator<double>> &, const std::vector<double, std::allocator<double>> &, const std::vector<double, std::allocator<double>> &, const std::vector<double, std::allocator<double>> &, std::vector<double, std::allocator<double>> &) 
  0.0     6.21E-06     6.21E-06  0.00621118           0          1 void runPlatformTests() 
  0.0            0            0  0.00621118           0          0 bool OpenMM::CustomAngleForce::usesPeriodicBoundaryConditions() const 
  0.0            0            0  0.00621118           0          0 const std::string &OpenMM::CustomAngleForce::getEnergyFunction() const 
  0.0            0            0  0.00621118           0          0 const std::string &OpenMM::CustomBondForce::getEnergyFunction() const 
  0.0            0            0   0.0186335           0          0 const std::vector<OpenMM::Vec3, std::allocator<OpenMM::Vec3>> &OpenMM::State::getForces() const 
  0.0            0            0   0.0124224           0          0 std::vector<std::string, std::allocator<std::string>> OpenMM::CustomExternalForceImpl::getKernelNames() 
  0.0            0            0  0.00621118           0          0 void OpenMM::CMMotionRemoverProxy::CMMotionRemoverProxy() 
  0.0            0            0   0.0124224           0          0 void OpenMM::CustomAngleForce::addEnergyParameterDerivative(const std::string &) 
  0.0            0            0  0.00621118           0          0 void OpenMM::CustomAngleForceImpl::~CustomAngleForceImpl() 
  0.0            0            0  0.00621118           0          0 void OpenMM::CustomBondForce::CustomBondForce(const std::string &) 
  0.0            0            0   0.0124224           0          0 void OpenMM::CustomBondForce::addEnergyParameterDerivative(const std::string &) 
  0.0            0            0   0.0124224           0          0 void OpenMM::CustomExternalForce::CustomExternalForce(const std::string &) 
  0.0            0            0   0.0248447           0          0 void OpenMM::CustomExternalForce::getParticleParameters(int, int &, std::vector<double, std::allocator<double>> &) const 
  0.0            0            0   0.0124224           0          0 void OpenMM::CustomExternalForceImpl::CustomExternalForceImpl(const OpenMM::CustomExternalForce &) 
  0.0            0            0   0.0124224           0          0 void OpenMM::CustomExternalForceImpl::~CustomExternalForceImpl() 
  0.0            0            0  0.00621118           0          0 void OpenMM::MonteCarloBarostat::MonteCarloBarostat(double, double, int) 
  0.0            0            0  0.00621118           0          0 void OpenMM::MonteCarloBarostatProxy::MonteCarloBarostatProxy() 
  0.0            0            0  0.00621118           0          0 void OpenMM::MonteCarloMembraneBarostatProxy::MonteCarloMembraneBarostatProxy() 
  0.0            0            0   0.0559006           0          0 void OpenMM::NonbondedForceImpl::NonbondedForceImpl(const OpenMM::NonbondedForce &) 
  0.0            0            0   0.0124224           0          0 void OpenMM::ReferenceAndersenThermostat::ReferenceAndersenThermostat() 
  0.0            0            0   0.0248447           0          0 void OpenMM::System::setDefaultPeriodicBoxVectors(const OpenMM::Vec3 &, const OpenMM::Vec3 &, const OpenMM::Vec3 &) 
  0.0            0            0  0.00621118           0          0 void OpenMM::System::setParticleMass(int, double)