Tutorial

Installation and Interface

Installation is described in the README of the base directory of FEASST. Or visit the website at https://pages.nist.gov/feasst

Python or C++?

Although this tutorial will focus upon the python interface, FEASST may also be utilized as a C++ library. For a C++ example, see the file tutorial/tutorial.cpp which is virtually identical to the Python tutorial. Thus Python tutorials are sufficient for learning the C++ library. While the majority of FEASST users prefer the Python interface, FEASST is written almost entirely in C++ for speed. Thus, both interfaces will be supported in the long term.

Some HPC clusters may not have the required Python libraries. If that is the case, do not hesitate to give the C++ interface a try, even if you have never written C++ before. For example, the only minor differences between tutorial/tutorial.cpp and tutorial/tutorial.py are the argument parsing syntax, semi-colons at the end of every line, and compiling any binary with int main().

Canonical ensemble Lennard-Jones Monte Carlo

The following simulation based on tutorial/tutorial.py demonstrates the basics of FEASST.

To begin, FEASST is imported, and it is recommended to log the exact version used for the simulation.

[1]:
import feasst
print("version:", feasst.version())
version: v0.7.0-90-g6cb364c4ae hwh/branch

Then a MonteCarlo object is created and the random number generator is initialized using the C++ implementation of the Mersenne Twister seeded by the time and date.

[2]:
monte_carlo = feasst.MonteCarlo()
monte_carlo.set(feasst.MakeRandomMT19937(feasst.args({"seed" : "time"})))

Care must be taken not to run two identical simulations at the same second on an HPC node using the time and date, or they will have the same seed and may be equivalent. Instead, consider using a thread safe random number generator to seed the simulations.

FEASST standard output and error goes to the Jupyter notebook terminal. Thus, you should see output like “# Info 0 [plugin/math/src/random.cpp:30] time(seed): 1572362164” but with a different seed. This seed is provided in case you wanted to reproduce your simulation exactly, which we will do now. Note that your simulation may still be different than presented here, possibly because of different compiler implementations.

[3]:
monte_carlo.set(feasst.MakeRandomMT19937(feasst.args({"seed" : "1572362164"})))

The second step is to add a Configuration. In this example, a simple cubic periodic box of length 8 is defined with a single type of particle as described in forcefield/data.lj. See Forcefield for more information about the format of the data file, which is a LAMMPS-inspired file with some major differences.

[4]:
monte_carlo.add(feasst.Configuration(feasst.MakeDomain(feasst.args({"cubic_box_length": "8"})),
     feasst.args({"particle_type": feasst.install_dir() + "/forcefield/data.lj"})))

FEASST arguments are input as a dictionary of strings with limited type checking. Thus, care must be taken to input strings which follow the documentation.

Next, initializing the Potential proceeds as follows:

[5]:
monte_carlo.add(feasst.Potential(feasst.MakeLennardJones()))
monte_carlo.add(feasst.Potential(feasst.MakeLongRangeCorrections()))

In this example, we introduce both the pair-wise Lennard-Jones (LJ) model, and also long-range corrections, which approximately account for the cut off of the LJ potential by assuming a pair-wise radial distance distribution function of unity. A FEASST convention is to use a helper function which appends the word Make onto the class name when creating pointers to FEASST derived class objects. This serves two purposes involving C++11 smart pointers and brace enclosed initializer lists.

Initialize ThermoParams, such as temperature, and the acceptance Criteria.

[6]:
monte_carlo.set(feasst.MakeThermoParams(feasst.args({"beta": "1.5"})))
monte_carlo.set(feasst.MakeMetropolis())

A TrialTranslate is then introduced which attempts to translate a random particle by a random distance which is bound in each dimension by a tunable_param. This parameter may be adjusted to obtain a desired acceptance ratio, tunable_target_acceptance, with the help of Tuner.

[7]:
monte_carlo.add(feasst.MakeTrialTranslate(feasst.args(
    {"tunable_param": "2.", "tunable_target_acceptance": "0.2"})))
steps_per = int(1e3)
monte_carlo.add(feasst.MakeTuner(feasst.args({"steps_per" : str(steps_per)})))

With the help of TrialTranslate, we can now initialize the number of particles.

[8]:
feasst.SeekNumParticles(50)\
    .with_thermo_params(feasst.args({"beta": "0.1", "chemical_potential": "10"}))\
    .with_metropolis()\
    .with_trial_add().run(monte_carlo)

A grand canonical simulation is performed here by utilizing a temporary TrialAdd, which is why chemical potential was input to Criteria.

Additional Analyze or Modify may be added at any time to perform some task contingent upon the number of attempted trials.

[9]:
monte_carlo.add(feasst.MakeLog(feasst.args({"steps_per" : str(steps_per),
                                            "file_name" : "log.txt",
                                            "clear_file" : "true"})))
monte_carlo.add(feasst.MakeMovie(feasst.args(
    {"steps_per" : str(steps_per), "file_name" : "movie.xyz"})))
monte_carlo.add(feasst.MakeCheckEnergy(feasst.args(
    {"steps_per" : str(steps_per), "tolerance" : "1e-8"})))

In this example, Log outputs the current status of the trials, Movie outputs the configuration, and CheckEnergy asserts that the optimized energy calculations match the unoptimized ones.

The simulation is finally run for a number of trial attempts.

[10]:
%%time
monte_carlo.attempt(int(1e5))
CPU times: user 797 ms, sys: 3.79 ms, total: 801 ms
Wall time: 798 ms

Now we can analyze the simulation by, for example, plotting the instantaneous energy as a function of the number of attempts.

[11]:
%matplotlib inline
%config InlineBackend.figure_format = 'svg'
import pandas as pd
pd.read_csv("log.txt").plot('attempt', 'energy')
[11]:
<matplotlib.axes._subplots.AxesSubplot at 0x7fd1f117e7d0>
../_images/tutorial_tutorial_21_1.svg

You should also find the movie.xyz trajectory file with an automatically-generated movie.xyz.vmd file for use with VMD (e.g., vmd -e movie.xyz.vmd).

Did this tutorial work as expected? Did you find any inconsistencies or have any comments? Please contact us. Any feedback is appreciated!

Additional tutorials

After completing this basic tutorial, check out the tutorials specific to each Plugin. For example, see the tutorials of the System, MonteCarlo and flat histogram plugins.