======== Tutorial ======== This tutorial is also available as `jupyter `_ and `org-mode `_ notebook. Quick start ----------- `Berni `_ provides a database of **interaction models** and **trajectory samples** to quickly start molecular dynamics or Monte Carlo simulations. The focus is on simple models to study liquids and glasses. The package strives to follow FAIR-data practises. The package provides an API to browse and retrieve models and samples: 1) ``models``: database of all available models 2) ``models.get("lennard_jones")``: get a specific model by its qualified name 3) ``samples``: database all available samples 4) ``samples.get("lennard_jones-13ce47602b259f7802e89e23ffd57f19.xyz")``: get a copy of a sample by its qualified name Both ``models`` and ``samples`` are slightly customized ``TinyDB`` databases, which you can query using the ``TinyDB`` methods. The interaction models can be exported in formats suitable for simulation with a range of backends, such as ``atooms``, ``RUMD``, and ``LAMMPS``. Interaction models ------------------ Let's see how to browse the models database. Available models ~~~~~~~~~~~~~~~~ The metadata of the available models can browsed by iterating over ``berni.models`` .. code:: python import berni for model in berni.models: print(model['name']) :: bernu_hiwatari_hansen_pastore coslovich_pastore dellavalle_gazzillo_frattini_pastore gaussian_core grigera_cavagna_giardina_parisi harmonic_spheres kob_andersen kob_andersen_2 lennard_jones roux_barrat_hansen wahnstrom These are the fields (or columns) of the database .. code:: python berni.models.fields() :: ['_model', 'description', 'doi', 'name', 'notes', 'reference'] You can pretty print the database like this .. code:: python berni.models.pprint(['name', 'description']) :: name description -------------------------------------------------------------------------------------------------------- bernu_hiwatari_hansen_pastore Binary soft-sphere mixture with size ratio of 1.4 coslovich_pastore Short-ranged pairwise-additive model for silica dellavalle_gazzillo_frattini_pastore Binary Lennard-Jones mixture model for NiY alloys gaussian_core One-component Gaussian core model with long-range cutoff grigera_cavagna_giardina_parisi Binary soft-sphere mixture with size ratio of 1.2 and smooth cutoff harmonic_spheres Binary mixture of harmonic spheres with size ratio of 1.4 kob_andersen Binary Kob-Andersen Lennard-Jones mixture kob_andersen_2 Ternary Kob-Andersen Lennard-Jones mixture lennard_jones One-component Lennard-Jones model roux_barrat_hansen Binary soft-sphere mixture with size ratio of 1.2 wahnstrom Binary Lennard-Jones mixture with size ratio of 1.2 Since ``berni.models`` inherits from ``TinyDB``, you can use all the ``TinyDB`` methods to browse and search the database, see the `TinyDB docs `_. The ``name`` field is a readable but "fully qualified" name that you can use to retrieve the specifications of the interaction model Get the one-component Lennard-Jones model and inspect the parameters .. code:: python model = berni.models.get('lennard_jones') pprint(berni.models.get("lennard_jones")) :: {'cutoff': [{'parameters': {'rcut': [[2.5]]}, 'type': 'cut_shift'}], 'potential': [{'parameters': {'epsilon': [[1.0]], 'sigma': [[1.0]]}, 'type': 'lennard_jones'}]} The model specifications are in a dict / json format, with two possibile schemas (1 and 2). The default schema requires the ``potential`` and ``cutoff`` fields. Additional fields are stored in the ``metadata`` field, but are not returned by the ``get`` method. Here is how the default schema looks like .. code:: python pprint(berni.schemas[1]) :: {'$schema': 'https://json-schema.org/draft/2020-12/schema', 'properties': {'cutoff': {'items': {'properties': {'parameters': {'type': 'object'}, 'type': {'type': 'string'}}, 'required': ['type', 'parameters'], 'type': 'object'}, 'type': 'array'}, 'metadata': {'properties': {'doi': {'type': 'string'}, 'name': {'type': 'string'}, 'notes': {'type': 'string'}, 'reference': {'type': 'string'}}, 'required': ['name'], 'type': 'object'}, 'potential': {'items': {'properties': {'parameters': {'type': 'object'}, 'type': {'type': 'string'}}, 'required': ['type', 'parameters'], 'type': 'object'}, 'type': 'array'}}, 'required': ['potential', 'cutoff'], 'type': 'object'} Schema version 2 stores the interaction parameters by pairs of species, instead of a single array. .. code:: python pprint(berni.models.get("kob_andersen", schema_version=2)) :: {'potential': [{'cutoff': {'parameters': {'1-1': {'rcut': 2.5}, '1-2': {'rcut': 2.0}, '2-2': {'rcut': 2.2}}, 'type': 'cut_shift'}, 'parameters': {'1-1': {'epsilon': 1.0, 'sigma': 1.0}, '1-2': {'epsilon': 1.5, 'sigma': 0.8}, '2-2': {'epsilon': 0.5, 'sigma': 0.88}}, 'type': 'lennard_jones'}]} This is useful to generate interaction potentials for some third-party backends. To choose or build your model, we need some more details about the schema layout, so keep reading. Potentials and cutoffs ~~~~~~~~~~~~~~~~~~~~~~ Each model is associated to a list of interaction potentials along with their parameters. If we limit ourselves to two-body potentials, we can write .. math:: u_{\alpha\beta}(r) = u_{\alpha\beta}^{(1)}(r) + u_{\alpha\beta}^{(2)}(r) + \dots where :math:`\alpha` and :math:`\beta` are chemical species in a mixture. Turning back to the LJ model .. code:: python model['potential'] :: [{'type': 'lennard_jones', 'parameters': {'epsilon': [[1.0]], 'sigma': [[1.0]]}}] Notice how the ``epsilon`` and ``sigma`` parameters are both 2d arrays, as they depend on the chemical species of the interacting particles. There is a cutoff field as well, since the potentials are typically cut off at some distance. .. code:: python model['cutoff'] :: [{'type': 'cut_shift', 'parameters': {'rcut': [[2.5]]}}] The potentials and cutoffs have standardized names .. code:: python for potential in sorted(berni.potentials): print(potential) :: fene gaussian harmonic inverse_power lennard_jones sum_inverse_power yukawa The cutoffs have standardized names and some aliases too .. code:: python for cutoff in sorted(berni.cutoffs): print(cutoff) :: cubic_spline cut cut_shift linear linear_cut_shift quadratic quadratic_cut_shift shift Trajectory samples ------------------ ``berni`` provides a database of trajectory samples to make it easy to - start a new simulation - check the implementation of interactions in third-party code - provide a little repository for further analysis Each sample is associated to one the models defined above. Of course, you are free to use these samples for other interaction models. Browsing the database ~~~~~~~~~~~~~~~~~~~~~ Browse the database of trajectory samples. .. code:: python berni.samples.pprint(['name', 'model'], sort_by='model') :: name model ---------------------------------------------------------------------------------------------------- bernu_hiwatari_hansen_pastore-f61d7e58b9656cf9640f6e5754441930.xyz bernu_hiwatari_hansen_pastore coslovich_pastore-488db481cdac35e599922a26129c3e35.xyz coslovich_pastore grigera_cavagna_giardina_parisi-0ac97fa8c69c320e48bd1fca80855e8a.xyz grigera_cavagna_giardina_parisi kob_andersen-8f4a9fe755e5c1966c10b50c9a53e6bf.xyz kob_andersen lennard_jones-5cc3b80bc415fa5c262e83410ca65779.xyz lennard_jones Get a local copy of a Lennard-Jones fluid sample. .. code:: python local_file = berni.samples.get("lennard_jones-13ce47602b259f7802e89e23ffd57f19.xyz") By default, ``local_file`` will be stored in a temporary directory. You can use it to start a simulation, benchmarking your code or further analysis using other packages. We rely on ``TinyDB`` to perform queries on the metadata of the available samples. For instance, here we look for samples of the Lennard-Jones models with a unit density .. code:: python query = berni.Query() samples = berni.samples.search((query.model == 'lennard_jones') & (query.density == 1.0)) pprint(samples) :: [{'density': 1.0, 'format': 'xyz', 'md5_hash': '13ce47602b259f7802e89e23ffd57f19', 'model': 'lennard_jones', 'name': 'lennard_jones-13ce47602b259f7802e89e23ffd57f19.xyz', 'number_of_particles': 256, 'potential_energy_per_particle': -3.8079776291909284, 'url': 'https://framagit.org/coslo/berni/-/raw/master/berni/samples/lennard_jones-13ce47602b259f7802e89e23ffd57f19.xyz'}, {'density': 1.0, 'format': 'xyz', 'md5_hash': '5cc3b80bc415fa5c262e83410ca65779', 'model': 'lennard_jones', 'name': 'lennard_jones-5cc3b80bc415fa5c262e83410ca65779.xyz', 'number_of_particles': 108, 'potential_energy_per_particle': 0.0, 'temperature': 0.0}] We can tidy up the output with ``pprint()`` .. code:: python berni.samples.pprint(['name', 'density', 'number_of_particles'], cond=(query.model == 'lennard_jones') & (query.density == 1.0)) :: name density number_of_particles ------------------------------------------------------------------------------ lennard_jones-13ce47602b259f7802e89e23ffd57f19.xyz 1.0 256 lennard_jones-5cc3b80bc415fa5c262e83410ca65779.xyz 1.0 108 Here is another way to inspect the full metadata of a given sample, identified this time by the file ``name`` .. code:: python query = berni.Query() pprint(berni.samples.search(query.name == "lennard_jones-5cc3b80bc415fa5c262e83410ca65779.xyz")) :: [{'density': 1.0, 'format': 'xyz', 'md5_hash': '5cc3b80bc415fa5c262e83410ca65779', 'model': 'lennard_jones', 'name': 'lennard_jones-5cc3b80bc415fa5c262e83410ca65779.xyz', 'number_of_particles': 108, 'potential_energy_per_particle': 0.0, 'temperature': 0.0}] For more advanced queries, check out what ``TinyDB`` `can do! `_ The database may provide some provenance information via the ``notes`` optional field. We do not aim at "full reproducibility" here. Exporting models for backends ----------------------------- In addition to specifying the interaction models, ``berni`` can generate potentials objects for different interaction "backends". In particular, there are implementations for the native ``f90`` atooms backend and for ``RUMD`` - and more to come! Atooms ~~~~~~ The model specifications of ``berni`` are compatible with the native ``f90`` atooms backend. Here we compute the potential energy of the configuration stored in ``local_file`` using ``lennard_jones`` model. .. code:: python from atooms.backends.f90.trajectory import Trajectory from atooms.backends.f90 import Interaction with Trajectory(local_file) as th: system = th[0] model = berni.get("lennard_jones") system.interaction = Interaction(model) print(system.potential_energy(per_particle=True)) :: -3.807977629191845 We can compare it to the reference value stored in the samples database and check they agree. This is useful if you want to test a new implementation of a backend. .. code:: python query = berni.Query() sample = berni.samples.search(query.name == "lennard_jones-13ce47602b259f7802e89e23ffd57f19.xyz")[0] sample["potential_energy_per_particle"] :: -3.8079776291909284 RUMD ~~~~ In addition to the native ``f90`` backend, there is currently support for ``RUMD``. It is possible to generate potential objects for RUMD like this .. code:: python potentials = berni.rumd.export("kob_andersen") The ``potentials`` list can then be passed to an ``Simulation`` instance defined in RUMD to carry out a simulation. LAMMPS ~~~~~~ .. code:: conf atom_style atomic pair_style lj/cut 2.5 pair_coeff 1 1 1.0 1.0 pair_coeff 2 2 0.5 3.0 pair_coeff 1 2 0.707 1.732