Atomic and Molecular Modeling Language (textM)

The Atomic and Molecular Modeling Language (textM) is designed as a modular extension of the Scientific Computing Language (textS).

The language provides builtin types, expressions and functions specific for atomistic modeling and simulation.

Structure

An atomic structure can be defined using the Structure literal. It has the following syntax

Structure [optional name] <Table>

The <Table> parameter is a Table literal with specific constraints (such as specific column names, types etc.) listed in the following table:

Name	Description	Type	Data type	Shape	Dimensionality
atoms	individual atoms structures	Series	Table	columns: symbols, x, y, z, px, py, pz, tags, masses	None
cell	unit cell	Series	FloatArray	(3, 3)	[length]
pbc	periodic boundary conditions	Series	BoolArray	(3)	None

Every tuple (row, record) in the top-level Table corresponds to an individual atomic, or molecular, structure. The atoms column is a special case: every element of the atoms column, is also a Table whose own columns hold the properties of chemical species such as: atomic symbols, positions, masses, momenta, etc. So the Table parameter of the Structure literal contains the atoms Table nested within. The columns of the atoms sub-table have the following constraints:

Name	Description	Type	Data type	Shape	Dimensionality
symbols	chemical symbol	Series	String	scalar	None
x	Cartesian x position	Series	Float	scalar	[length]
y	Cartesian y position	Series	Float	scalar	[length]
z	Cartesian z position	Series	Float	scalar	[length]
px	Cartesian momentum along x	Series	Float	scalar	[length mass time-1]
py	Cartesian momentum along y	Series	Float	scalar	[length mass time-1]
pz	Cartesian momentum along z	Series	Float	scalar	[length mass time-1]
masses	atomic mass	Series	Float	scalar	[mass]
tags	tag	Series	Integer	scalar	[dimensionless]

Example:

water = Structure H2O (
         (atoms: ((symbols: 'O', 'H', 'H'),
                  (x: 0., 0., 0.) [nm],
                  (y: 0., 0.763239, -0.763239) [angstrom],
                  (z: 0.119262, -0.477047, -0.477047) [angstrom],
                  (tags: 1, 0, 0),
                  (masses: 16., 1., 1.) [amu]
                 )
         ),
         (cell: [[5., 0., 0.], [0., 5., 0.], [0., 0., 5.]] [bohr]),
         (pbc: [false, false, false])
       )

Loading a structure from a file

An atomic structure can be loaded from a file, e.g.

water = Structure from file 'h2o.xyz'

Note that, when loading from file, there would be no water.name parameter.

Creating a structure by using the builders library

An atomic structure can be created by using builder functions from a specialized built-in library. A builder is a function taking a Table as argument and returning a Structure. For example, to build and view a periodic fcc-Cu system one can use the statements:

use bulk from virtmat.builders
fcc_cu = bulk(((name: 'Cu'), (crystalstructure: 'fcc')))
view structure (fcc_cu,)

Builders with no mandatory parameters must be called without arguments:

use mx2 from virtmat.builders
m = mx2()

By default, the structure is named after the builder function. Optionally, a custom name can be given by passing a second argument to the builder:

use graphene from virtmat.builders
g = graphene(((formula: 'BN')), 'BN')
print(g.name)  # 'BN'

A list of all implemented builder functions can be obtained with

use BUILDERS from virtmat.builders
print(BUILDERS)

Full descriptions of the builders with their parameters can be found in the ASE documentation.

Accessing the structure attributes

The parameters contained in a Structure parameter can be accessed using the common subscripting or slice syntax:

print(water.name) # name of structure
print(water.atoms) # the series with atoms tables
print(water.atoms[0]) # atoms table of the first structure
print(water.atoms[0:1]) # series with atoms inluding only the first structure
print(water.pbc) # series of boolean arrays
print(water.cell) # series of float arrays
print(water.cell[0]) # the cell of the first structure as float array
print(water[0]) # a tuple (atoms, cell, pbc) of the first structure
print(water[0:1]) # a table containing the first structure

Computed parameters available in Structure

All parameters of the structure do not depend on external parameters but some of them depend on other Structure parameters. These dependent parameters can be computed without further input. For example, the kinetic energy depends on the momenta and the masses. These parameters can be accessed in the same way as the parameters in the Structure literal, for example the kinetic energy of the nuclei:

print(water.kinetic_energy)

The computed parameters are summarized in the following table.

Name	Description	Type	Data type	Shape	Dimensionality
kinetic_energy	nuclear kinetic energy	Series	float	scalar	[length2 mass time-2]
temperature	temperature	Series	float	scalar	[temperature]
distance_matrix	matrix of inter-atomic distances	Series	FloatArray	(N, N)	[length]
chemical_formula	brute formula / composition	Series	String	scalar	None
number_of_atoms	number of atoms	Series	Integer	scalar	None
cell_volume	volume of the cell	Series	Float	scalar	[length]3
center_of_mass	center of mass	Series	FloatArray	(3, )	[length]
radius_of_gyration	radius of gyration	Series	Float	scalar	[length]
moments_of_inertia	moments of inertia along the principal axes	Series	FloatArray	(3, )	[mass length2]
angular_momentum	the total angular momentum with respect to the center of mass	Series	FloatArray	(3, )	[mass length2 time-1]

Using build tools to process a structure

The build tools take a structure as input and return a new structure as a result of their operation. The build tools include object-mutating methods of the ase.Atoms class. These are the functions repeat, center, translate, rotate, euler_rotate, rattle and wrap. Furthermore, the build tools include functions from the ase.build module that process one or more atoms objects as inputs. These include, e.g., the functions cut, stack, sort, minimize_tilt, niggli_reduce, minimize_rotation_and_translation.

The full list of build tools provided in textM can be seen by printing the virtmat.builders.TOOLS variable:

use TOOLS from virtmat.builders
print(TOOLS)

The use of the tools is similar to the builder functions. The main difference is that the build tools require a Structure parameter and a Table parameter as first and second argument, respectively.

Calculator

A Calculator parameter describes operations to be performed on a Structure parameter typically employing computational tools external to textM, e.g. molecular dynamics or quantum chemistry packages. It is commonly used within the Property parameter. The Calculator has the following syntax:

Calculator <name> [ >|<|==|>=|<= <version>] <Table>[, task: <task>]
Calculator <name> [ >|<|==|>=|<= <version>] (default)[, task: <task>]
Calculator <name> [ >|<|==|>=|<= <version>] ()[, task: <task>]

The specification of name is mandatory. Currently name must be one of vasp, turbomole, lj, emt, free_electrons and elk.

The version specification is optional. The calculator version must be semantic version and is interpreted for codes are not part of ASE, such as VASP, Turbomole and Elk.

The actual calculator parameters can be specified as series in the top-level table. The parameter names, type and units depend on the calculator. The table below lists the parameters of the currently supported calculators. If only default parameters have to be used then no Table parameter is needed and the syntax is (default) or simply ().

Example for a VASP calculator with one parameter set:

calc = Calculator vasp == 5.4.4 (
        (algo: 'Fast'),
        (ediff: 1e-06) [eV],
        (ediffg: -0.01) [eV/angstrom],
        (encut: 400.0) [eV],
        (ibrion: 2),
        (icharg: 2),
        (isif: 2),
        (ismear: 0),
        (ispin: 2),
        (istart: 0),
        (kpts: [5, 5, 1]),
        (lcharg: false),
        (lreal: 'Auto'),
        (lwave: false),
        (nelm: 250),
        (nsw: 1500),
        (potim: 0.1),
        (prec: 'Normal'),
        (sigma: 0.1) [eV],
        (xc: 'PBE')
), task: local minimum

NOTE: The textM language can be extended to support further calculators. This is fairly straightforward because the current implementation of textM interpreter is based on the Atomic Simulation Environment (ASE). Particularly, the Calculator parameter can reuse any existing ASE Calculator class. The process of integration of further calculators into textM is basically a registration of the name in the grammar and of function to check calculator parameter types and units. See this section for more details.

The parameters stored in a calculator can be accessed using subscripting and slice syntax:

print(calc.parameters) # top-level table of parameters
print(new.name) # calculator's name
print(calc.pinning, calc.version) # if not set `null` will be printed
print(calc[0]) # a tuple of the first parameter set
print(new_calc[0:1]) # table with the first parameter set

NOTE: The calculator parameter names and types are identical with those in the Atomic Simulation Environment (ASE). For this information we refer to the ASE documentation of the calculators. The physical units (or their dimension) of the calculator parameters are not always fully and systematically documented. While units / dimensionality checks are already performed independently by the textM interpreter (with error messages suggesting relevant corrections), this information is not yet provided here. For now, the textM user should consult the manuals of the original codes.

The special parameter task can be currently one of the following strings: single point, local minimum, global minimum, transition state, normal modes, micro-canonical,canonical, isothermal-isobaric, grand-canonical. Currently, the task parameter is only used to perform consistency checks.

NOTE: The Elk calculator is currently in alpha phase: 1) There is an open issue with the ASE Elk calculator that prevents using the calculator in some use cases; Until the issue is fixed, as a work-around, ASE must be installed from this branch; 2) the units definitions in the interpreter spec are not complete and not available for some input blocks.

Constraint

The Constraint parameter describes constraints to be applied to atomic positions, momenta and forces. A set of constraints is pertaining to (but not an integral part of) a specific Structure parameter.

c1 = FixedAtoms where <fixed>
c2 = FixedLine collinear to <direction> where <fixed>
c3 = FixedPlane normal to <direction> where <fixed>
c4 = FixedCOM  # fix the center of mass of the system

The <fixed> parameter is a Series of Boolean data type, whose length should be equal to the number of atoms in the structure to which the constraint will be applied. It indicates which atoms are included in the constraint (true elements) and which are not (false elements).

The <direction> parameter is an Array of length 3 and Integer data type, and represents a vector in Cartesian space [x, y, z].

Example:

h2o_plane = FixedPlane normal to [0, 0, 1] where (fix: true, true, true) # allow changes in xy plane
h2_line = FixedLine collinear to [0, 0, 1] where (fix: true, true) # allow changes along z axis
h2_fix0 = FixedAtoms where (fix: true, false) # fix first atom in H2

NOTE: Extensions of textM with further types of constraints is possible. See this section for more details.

Algorithm

The Algorithm parameter describes operations that are either available as algorithms in the Atomic Simulation Environment (ASE) or are user-defined routines. For example, the BFGS structure optimization is an algorithm. Further available algorithms are summarized in the table below.

The general syntax of Algorithm is:

Algorithm <name>[, many_to_one] <Table>
Algorithm <name>[, many_to_one] (default)
Algorithm <name>[, many_to_one] ()

The <name> is a placeholder for the algorithm name which is one of the listed above. If the optional flag many_to_one is not specified then the algorithm is applied per geometry in a structure, similar to how map is applied to series. If many_to_one is specified then the algorithm is applied to all geometries in a structure, similar to how reduce is applied to series.

The <Table> is a Table that contains input parameters for the algorithm, specified in Series that are the columns of the table. As for the calculators, the numerical parameters of algorithms must have explicit units. The parameter names, data types, default values, dimensionality and usual units are summarized for every algorithm in the table below. If (default) or () is specified then only default input parameters are used.

Algorithm	Parameter name	Parameter data type	Default value	Dimensionality	Sample unit
BFGS	steps	Integer	100000000	[dimensionless]	—
	maxstep	Float	0.2 [Å]	[length]	[Å]
	fmax	Float	0.05 [eV Å-1]	[length mass time2]	[eV Å-1]
	logfile	Boolean	null	—	—
	trajectory	Boolean	null	—	—
	interval	Integer	1	[dimensionless]	—
	alpha	Float	70.0 [eV Å-2]	[mass time-2]	[eV Å-2]
BFGSLineSearch	steps	Integer	100000000	[dimensionless]	—
	maxstep	Float	0.2 [Å]	[length]	[Å]
	fmax	Float	0.05 [eV Å-1]	[length mass time-2]	[eV Å-1]
	logfile	Boolean	null	—	—
	trajectory	Boolean	null	—	—
	interval	Integer	1	[dimensionless]	—
LBFGS	steps	Integer	100000000	[dimensionless}	—
	maxstep	Float	0.2 [Å]	[length]	[Å]
	fmax	Float	0.05 [eV Å-1]	[length mass time-2]	[eV Å-1]
	logfile	Boolean	null	—	—
	trajectory	Boolean	null	—	—
	interval	Integer	1	[dimensionless]	—
	memory	Integer	100	[dimensionless]	—
	damping	Float	1.0	[dimensionless]	—
	alpha	Float	70.0 [eV Å-2]	[mass time-2]	[eV Å-2]
LBFGSLineSearch	steps	Integer	100000000	[dimensionless]	—
	maxstep	Float	0.2 [Å]	[length]	[Å]
	fmax	Float	0.05 [eV Å-1]	[length mass time-2]	[eV Å-1]
	logfile	Boolean	null	—	—
	trajectory	Boolean	null	—	—
	interval	Integer	1	[dimensionless]	—
	memory	Integer	100	[dimensionless]	—
	damping	Float	1.0	[dimensionless]	—
	alpha	Float	70.0 [eV Å-2]	[mass time-2]	[eV Å-2]
QuasiNewton	steps	Integer	100000000	[dimensionless]	—
	maxstep	Float	0.2 [Å]	[length]	[Å]
	fmax	Float	0.05 [eV Å-1]	[length mass time-2]	[eV Å-1]
	logfile	Boolean	null	—	—
	trajectory	Boolean	null	—	—
	interval	Integer	1	[dimensionless]	—
GPMin	steps	Integer	100000000	[dimensionless]	—
	fmax	Float	0.05 [eV Å-1]	[length mass time-2]	[eV Å-1]
	logfile	Boolean	null	—	—
	trajectory	Boolean	null	—	—
	interval	Integer	1	[dimensionless]	—
	prior	Object	null	—	—
	kernel	Object	null	—	—
	update_prior_strategy	String	‘maximum’	—	—
	update_hyperparams	Boolean	false	—	—
	noise	Float	0.005	[dimensionless]	—
	weight	Float	1.0	[dimensionless]	—
	scale	Float	0.4	[dimensionless]	—
FIRE	steps	Integer	100000000	[dimensionless]	—
	maxstep	Float	0.2 [Å]	[length]	[Å]
	fmax	Float	0.05 [eV Å-1]	[length mass time-2]	[eV Å-1]
	logfile	Boolean	null	—	—
	trajectory	Boolean	null	—	—
	interval	Integer	1	[dimensionless]	—
	dt	Float	0.1 [Å (amu / eV)1/2]	[time]	[ps]
	dtmax	Float	1.0 [Å (amu / eV)1/2]	[time]	[ps]
	downhill_check	Boolean	false	—	—
MDMin	steps	Integer	100000000	[dimensionless]	—
	maxstep	Float	0.2 [Å]	[length]	[Å]
	fmax	Float	0.05 [eV Å]-1	[length mass time-2]	[eV Å-1]
	logfile	Boolean	null	—	—
	trajectory	Boolean	null	—	—
	interval	Integer	1	[dimensionless]	—
	dt	Float	0.2 [Å (amu / eV)1/2]	[time]	[ps]
RMSD	reference	Structure	—	—	—
	adjust	Boolean	true	—	—
RDF	rmax	Float	null [Å]	[length]	[Å]
	nbins	Integer	40	[dimensionless]	—
	neighborlist	Table	null	—	—
	elements	String	null	—	—
VelocityDistribution	distribution	String	‘maxwell-boltzmann’	—	—
	temperature_K	Float	—	[temperature]	[K]
	force_temp	Boolean	false	—	—
	stationary	Boolean	true	—	—
	zero_rotation	Boolean	true	—	—
	preserve_temperature	Boolean	true	—	—
	nbins	Integer	40	[dimensionless]	—
EquationOfState	energies	FloatArray	—	[length2 mass time-2]	[eV]
	eos	String	‘sjeos’	—	—
	filename	String	null	—	—
DensityOfState	npts	Integer	401	[dimensionless]	—
	width	Float	0.1 [eV]	[length2 mass time-2]	[eV]
	window	FloatArray	null	[length2 mass time-2]	[eV]
	spin	Integer	null	[dimensionless]	—
BandStructure	filename	String	null	—	—
	emin	Float	null	[length2 mass time-2]	[eV]
	emax	Float	null	[length2 mass time-2]	[eV]
VelocityVerlet	steps	Integer	50	[dimensionless]	—
	logfile	Boolean	null	—	—
	trajectory	Boolean	null	—	—
	interval	Integer	1	[dimensionless]	—
	timestep	Float	—	[time]	[ps]
Langevin	steps	Integer	50	[dimensionless]	—
	logfile	Boolean	null	—	—
	trajectory	Boolean	null	—	—
	interval	Integer	1	[dimensionless]	—
	timestep	Float	—	[time]	[ps]
	temperature_K	Float	—	[temperature]	[K]
	friction	Float	—	[time-1]	[ps-1]
	fixcm	Boolean	false	—	—
MelchionnaNPT	steps	Integer	50	[dimensionless]	—
	logfile	Boolean	null	—	—
	trajectory	Boolean	null	—	—
	interval	Integer	1	[dimensionless]	—
	timestep	Float	—	[time]	[ps]
	temperature_K	Float	—	[temperature]	[K]
	externalstress	Float	—	[mass length-1 time-2]	[N m-2], [Pa]
	ttime	Float	—	[time]	[ps]
	pfactor	Float	—	[mass length-1]	[amu Å-1]
Andersen	steps	Integer	50	[dimensionless]	—
	logfile	Boolean	null	—	—
	trajectory	Boolean	null	—	—
	interval	Integer	1	[dimensionless]	—
	timestep	Float	—	[time]	[ps]
	temperature_K	Float	—	[temperature]	[K]
	andersen_prob	Float	—	[dimensionless]	—
	fixcm	Boolean	true	—	—
NVTBerendsen	steps	Integer	50	[dimensionless]	—
	logfile	Boolean	null	—	—
	trajectory	Boolean	null	—	—
	interval	Integer	1	[dimensionless]	—
	timestep	Float	—	[time]	[ps]
	temperature_K	Float	—	[temperature]	[K]
	taut	Float	—	[time]	[ps]
	fixcm	Boolean	true	—	—
NPTBerendsen	steps	Integer	50	[dimensionless]	—
	logfile	Boolean	null	—	—
	trajectory	Boolean	null	—	—
	interval	Integer	1	[dimensionless]	—
	timestep	Float	—	[time]	[ps]
	temperature_K	Float	—	[temperature]	[K]
	pressure_au	Float	—	[mass length-1 time-2]	[N m-2], [Pa]
	taut	Float	0.5 [ps]	[time]	[ps]
	taup	Float	1.0 [ps]	[time]	[ps]
	compressibility	Float	—	[mass length-1 time]-2]	[N m-2], [Pa]
	fixcm	Boolean	true	—	—
IsotropicMTKNPT	steps	Integer	50	[dimensionless]	—
	logfile	Boolean	null	—	—
	trajectory	Boolean	null	—	—
	interval	Integer	1	[dimensionless]	—
	timestep	Float	—	[time]	[fs]
	temperature_K	Float	—	[temperature]	[K]
	pressure_au	Float	—	[mass length-1 time-2]	[N m-2], [Pa]
	tchain	Integer	1	[dimensionless]	—
	pchain	Integer	1	[dimensionless]	—
	tloop	Integer	—	[dimensionless]	—
	ploop	Integer	—	[dimensionless]	—
MTKNPT	steps	Integer	50	[dimensionless]	—
	logfile	Boolean	null	—	—
	trajectory	Boolean	null	—	—
	interval	Integer	1	[dimensionless]	—
	timestep	Float	—	[time]	[fs]
	temperature_K	Float	—	[temperature]	[K]
	pressure_au	Float	—	[mass length-1 time-2]	[N m-2], [Pa]
	tchain	Integer	3	[dimensionless]	—
	pchain	Integer	3	[dimensionless]	—
	tloop	Integer	1	[dimensionless]	—
	ploop	Integer	1	[dimensionless]	—
LangevinBAOAB	steps	Integer	50	[dimensionless]	—
	logfile	Boolean	null	—	—
	trajectory	Boolean	null	—	—
	interval	Integer	1	[dimensionless]	—
	timestep	Float	—	[time]	[fs]
	temperature_K	Float	—	[temperature]	[K]
	externalstress	Float	—	[stress]	[eV Å-3]
	hydrostatic	Boolean	False	—	—
	T_tau	Float	—	[time]	[fs]
	P_tau	Float	—	[time]	[fs]
	P_mass	Float	—	[mass]	[amu]
	P_mass_factor	Float	1.0	[dimensionless]	—
	disable_cell_langevin	Boolean	False	—	—
	rng		—	—	—
NEB	number_of_images	Integer	3	[dimensionless]	—
	trajectory	Boolean	true	—	—
	logfile	Boolean	null	—	—
	interpolate_method	String	‘linear’	—	—
	interpolate_mic	Boolean	false	—	—
	dynamic_relaxation	Boolean	false	—	—
	fmax	Float	0.05 [eV Å-1]	[length mass time-2]	[eV Å-1]
	scale_fmax	Float	0.0 [Å-1]	[length-1]	[Å-1]
	k	Float	0.1 [eV Å-1]	[length mass time-2]	[eV Å-1]
	climb	Boolean	false	—	—
	remove_rotation_and_translation	Boolean	false	—	—
	method	String	‘aseneb’	—	—
	optimizer	Table	—	—	—
	fit	Boolean	false	—	—
	filename	String	null	—	—
Dimer	target	Structure	null	—	—
	fmax	Float	0.05 [eV Å-1]	[length mass time-2]	[eV Å-1]
	trajectory	Boolean	null	—	—
	logfile	Boolean	null	—	—
	eigenmode_logfile	String	null	—	—
	eigenmode_method	String	‘dimer’	—	—
	f_rot_min	Float	0.1 [eV Å-1]	[length mass time-2]	[eV Å-1]
	f_rot_max	Float	1.0 [eV Å-1]	[length mass time-2]	[eV Å-1]
	max_num_rot	Integer	1	[dimensionless]	—
	trial_angle	Float	0.25 $\pi$ [radians]	[angle]	[radians]
	trial_trans_step	Float	0.001 [Å]	[length]	[Å]
	maximum_translation	Float	0.1 [Å]	[length]	[Å]
	cg_translation	Boolean	true	—	—
	use_central_forces	Boolean	true	—	—
	dimer_separation	Float	0.0001 [Å]	[length]	[Å]
	initial_eigenmode_method	String	‘gauss’	—	—
	extrapolate_forces	Boolean	False	—	—
	displacement_method	String	‘gauss’	—	—
	gauss_std	String	0.1 [Å]	[length]	[Å]
	order	Integer	1	[dimensionless]	—
	mask	BoolArray	null	—	—
	displacement_center	FloatArray	null	[length]	[Å]
	displacement_radius	Float	null	[length]	[Å]
	number_of_displacement_atoms	Integer	null	[dimensionless]	—
Vibrations	delta	Float	0.01 [Å]	[lenght]	[Å]
	nfree	Integer	2	[dimensionless]	—
	method	String	‘standard’	—	—
	direction	String	‘central’	—	—
	all_atoms	Boolean	false	—	—
	imag_tol	Float	1e-5 [eV]	[length2 mass time-2]	[eV]
NeighborList	cutoffs	FloatArray	null	[length]	[Å]
	skin	Float	0.3 [Å]	[length]	[Å]
	sorted	Boolean	false	—	—
	self_interaction	Boolean	true	—	—
	bothways	Boolean	false	—	—
	method	String	‘quadratic’	—	—

NOTE: If the default is --- then the parameter must be specified in the input. If the default is null then the parameter is not required in the input and will be set automatically by the algorithm.

NOTE: If trajectory / logfile are set to true a trajectory file / logfile will be created in the datastore location. If logfile is set to false then the step logging will be written on the screen. By default, logfile is set to null, i.e. step logging is turned off completely.

NOTE: For the Dimer algorithm only the textM-specific parameters are listed in the table. All other parameters are identical with the ASE parameters.

NOTE: To fix the center of mass in MD simulations with the Langevin thermostat, do no set (fixcm: true). Instead, use the FixedCOM constraint and add it to Property.

Example: Structure optimization of a structure h2o with calculator calc using the algorithm BFGS.

algo = Algorithm BFGS ((fmax: 1e-4) [hartree/bohr], (steps: 30))
prop = Property energy, forces ((structure: h2o), (calculator: calc), (algorithm: algo))

Example: Compute the radial distribution function for all particles for all geometries in a molecular dynamics trajectory using the algorithm RDF.

prop_samp = ... # sampling with molecular dynamics
struct = prop_samp.trajectory[0].structure
algo = Algorithm RDF, many_to_one ((rmax: 1.7) [angstrom], (nbins: 50))
prop = Property rdf, rdf_distance ((structure: struct), (algorithm: algo))
view lineplot (prop.rdf, prop.rdf_distance [angstrom], (legend: 'RDF'))

Trajectory

Some algorithms and calculators provide trajectories as properties. The trajectory type has these properties:

• description: A table including information about the algorithm that has produced the trajectory. For example, here the description from a Langevin molecular dynamics run:

  ((type: 'molecular-dynamics'),
   (md-type: 'Langevin'),
   (timestep: 0.09822694750253275) [angstrom * unified_atomic_mass_unit ** 0.5 / electron_volt ** 0.5],
   (temperature_K: 300.0) [kelvin],
   (friction: 0.5090252855379708) [electron_volt ** 0.5 / angstrom / unified_atomic_mass_unit ** 0.5],
   (fixcm: true),
   (interval: 1))
  

  Note that the units are the units used by ASE internally.

• structure: A structure containing all configurations of the processed structure.

• properties: A table with the properties stored for each configuration.

• constraints: A tuple of constraints

• filename: The path to the file where the trajectory is stored.

The trajectory can be extracted from the property using the common syntax prop.trajectory which returns a series of trajectory parameters. The individual trajectories can be accessed by subscripting, e.g. the first element can be retrieved with prop.trajectory[0]. The trajectory attributes can be then accessed using the common syntax trajectory[0].description.

Property

Using the Property parameter a set of properties will be computed for a specific Structure with a specific method. The method can be specified either by a Calculator or by an Algorithm, or a combination of Algorithm and Calculator.

Syntax:

prop = Property name1[, name2[...]]
        ((structure: struct),
         (calculator: calc),
         (algorithm: algo),
         (constraints: (c1[, c2[, ...]]))
        )

The names of currently interpreted properties are summarized in the table below. The parameters struct, calc, algo, and c1, … are names of variables of the corresponding types Structure, Calculator, Algorithm and Constraint respectively. The ordering of these parameters is arbitrary. Only the structure parameter is mandatory. The following table provides an overview of currently supported properties with their types and default output units.

Propery name	Data type	Units (of output)	Description
energy_minimum	Boolean	None	true if geometry is minimum, false otherwise; available only with a calculator with task normal modes
transition_state	Boolean	None	true if geometry is transition state, false otherwise; available only with a calculator with task normal modes
energy	Float	[eV]	total energy without nuclear kinetic energy; if the structure has been modified by the algorithm or calculator then this is the energy of the final structure
forces	FloatArray	[eV / Å]	the forces for each atom, i.e. the first derivative of energy with a minus sign; if the structure has been modified by the algorithm or calculator then these are the forces in the final structure
dipole	FloatArray	[e Å]	the total electrostatic dipole moment; if the structure has been modified by the algorithm or calculator then this is the dipole moment of the final structure
magmom	Float	[dimensionless]	Total magnetic moment of the system
magmoms	FloatArray	[dimensionless]	Atomic magnetic moments
vibrational_energies	Series(Float)	[eV]	vibrational energies (frequencies) corresponding to the non-negative eigenvalues of the harmonic normal modes; available only with a calculator with task normal modes
stress	FloatArray	[eV / Å3]	the total stress in Voigt notation
trajectory	Trajectory	None	Trajectory provided by some algorithms and calculators
minimum_energy	Float	[eV]	the energy of the minimum of the fitted equation of state
optimal_volume	Float	[Å3]	the cell volume at minimum energy from the fitted equation of state
bulk_modulus	Float	[eV / Å3]	the bulk modulus, the slope of the fitted equation of state
eos_energy	FloatArray	[eV]	energies fitted to the equation of state
eos_volume	FloatArray	[Å3]	volumes fitted to the equation of state
dos_energy	FloatArray	[eV]	energy array for the density of states
dos	FloatArray	[dimensionless]	density of states
band_structure	Table	None	table holding the band structure
neb_energy	FloatArray	[eV]	energy array of the images from an NEB calculation
activation_energy	Float	[eV]	the activation energy (the energy barrier)
reaction_energy	Float	[eV]	the energy change from the initial to the final state
maximum_force	Float	[eV / Å]	the maximum force
velocity	FloatArray	[(eV / amu)1/2]	velocity array (the bin edges) for the velocity distribution histogram
vdf	FloatArray	[dimensionless]	velocity distribution function (histogram)

The output units can be changed using a modifier.

After initializing a variable with a Property parameter, the individual properties can be accessed using subscripting, i.e. with prop.<property name>.

NOTE: The expression prop.<property name> will always return a Series of the type specified in the table above.

Examples:

props = Property energy, forces ((structure: h2o), (calculator: calc))

print(props) # print the whole property parameter
print(props[0]) # return the first tuple of the property parameter
print(props.structure) # return the initial Structure parameter
print(props.calculator) # return the initial Calculator parameter
print(props.names) # return the property names
print(props.calculator[0:1]) # return a slice of the initial calc
print(props.calculator.name) # return the name of the inital calc
print(props.calculator.parameters) # return a Table of the calc parameters
print(props.calculator.version) # return the calc version
print(props.results) # return a Table with all calculated properties
print(props.output_structure) # return the output Structure shich should contain only one geometry
print(props.output_structure[0:1]) # return only the first geometry (slice) of a Structure parameter containing several geometries
print(props.energy) # return energy as Series for all output structures
print(props.energy[0:1]) # return a slice (Series) with the energy of the first output structure
print(props.energy[0]) # return energy (Float) of the first output structure
print(props.forces) # return forces as series for all output structures
print(props.forces[0]) # return an Array of forces for the first output structure

The Property parameter is evaluated in the following way. The input Structure and Calculator parameters are subscriptable and contain Series of structure geometries and calculator parameters, respectively. The Cartesian product of these two Series is taken to find all combinations and for each combination a calculation is performed. Every calculation is captured as a tuple (a row) in props.results Table. The props.results Table has as columns the names of all requested properties, as well as the corresponding input geometry, calculator parameters and final (output) geometry.

For example, if two geometries of the same molecule have to be processed by the same calculator with three different parameter sets, then the results Table will contain six tuples. The computed constraints and the properties are the same for all structure-calculator-algorithm combinations. The results can be accessed via subscripting the Property parameter, as it is shown in the examples above.

Computational resources for Property

Very often Property parameters require multiple processor cores and multiple computing nodes and long computing times (up to several days) to evaluate. In only a few exceptional cases, e.g. computing small molecules or periodic cells with a few atoms, one can evaluate a property in an interactive statement. To enable batch evaluation one has to add resource requirements at the end of the statement as introduced here.

Storage resources for Property

In many use cases, the evaluation of Property parameters yields large output data. This data does not fit into a single document in the database. Therefore a mechanism is in place that stores data with excessive volumes in a separate datastore. Typical use cases where the mechanism gets activated are storing trajectories from molecular dynamics and Monte Carlo simulations, and charge densities from density functional calculations. See this section to learn more about this mechanism and how to customize it, i.e. to change the threshold volume, type of store etc.

Property with random Algorithm/Calculator

There are plenty of algorithms and calculators whose results are pseudo-random. These results are always deterministic and can be reproduced by setting the exact algorithm for (pseudo) random number generation and a seed, e.g.

struct = ...
algo = Algorithm VelocityDistribution ((temperature_K: 298.15) [K])
prop = Property vdf ((algorithm: algo), (structure: struct),
                     (rng: 'MT19937'), (seed: 'aed8c93f1db74abebc32d228a25f7628'))

See the Random parameter in textS for more details about rng and seed.

When rng is not specified then the default bit generator will be selected. After time, another bit generator might be selected by default in numpy (see this link) and if the generator has not been explicitly specified, then the result cannot be reproduced. If seed is not specified, a random seed will be generated and in this case the result of the model cannot be reproduced if evaluated repeatedly.

Whenever reproducibility is a concern, it is recommended to set both rng and seed.

NOTE: In workflow evaluation mode, the values of rng and seed used in the evaluation are stored and can later be found in the value of the Property parameter.

Checkpointing and recovery of Property evaluations

Failed Property evaluation can be partially recovered based on a checkpointing mechanism described in detail here.

Using textM parameters as table-like iterables

As seen above, some textM parameters behave like iterables. All operations defined for the Table type can be applied to the following textM parameters: Structure, Calculator, Algorithm and Property. These operations include, e.g. subscripting, slicing, filter function and expression, map and reduce. These operations base on a table attribute in each of these parameters:

Parameter	Iterable
Structure	Table with atoms, cell, pbc etc.
Algorithm	Table with parameters
Calculator	Table with parameters
Property	Table with results

Examples:

calc = Calculator vasp == 5.4.4 ((lreal: true, false),
                                 (ediff: 1e-06, 1e-05) [eV],
                                 (ediffg: -0.01, -0.05) [eV/angstrom])
calc_0 = calc[0] # return the first parameters tuple
calc_01 = calc[:1] # return a calc with the first row of parameters
calc_lr = calc select ediff, ediffg where column:lreal == false
c_f = filter((x: x.lreal), calc)
c_m = map((x: {lreal: x.lreal, ediff: 2*x.ediff}), calc)
c_r = reduce((x, y: {lreal: (x.lreal and y.lreal)}), calc)

Using the view statement

The view statement can be used to visualize and analyze different textM parameters.

Visualize / analyze atomic structures

Parameters of type Structure with optional geometry constraints can be viewed using this syntax:

view structure (struct, [(constr1, [[constr2], ...])])

Visualize a trajectory

Parameters of type Trajectory can be viewed using the same syntax but the keyword trajectory:

view trajectory (traj)

Visualize a vibrational normal mode

Vibrational normal modes can be visualized by the following statement

view vibration (structure, hessian, ((mode: 4), (constraints: (constr,))))

The hessian is an array returned by the Vibrations algorithm and by some calculators. Only one mode can be viewed at a time. The list of constraints is optional and must be the one that has been used to compute the hessian.

Visualize a nudged elastic band (NEB) simulation

Trajectory that has been produced from an NEB simulation can be evaluated using the view statement:

view neb (traj)

Visualize fits of the equation of state

The equation of state fits can be visualized using this syntax

view eos (<volumes array>, <energies array>, <eos type str>)

For supported values of <eos type str> we refer to the ASE documentation.

Example:

atoms = ((symbols: 'Ag'), (x: 0.) [angstrom], (y: 0.) [angstrom], (z: 0.) [angstrom])
ag = Structure Ag (
                (atoms: atoms, atoms, atoms),
                (pbc: [true, true, true], [true, true, true], [true, true, true]),
                (cell: [[0., 1.9, 1.9], [1.9, 0., 1.9], [1.9, 1.9, 0.]] [angstrom],
                       [[0., 2.0, 2.0], [2.0, 0., 2.0], [2.0, 2.0, 0.]] [angstrom],
                       [[0., 2.1, 2.1], [2.1, 0., 2.1], [2.1, 2.1, 0.]] [angstrom]
                )
)
calc = Calculator emt (), task: single point
prop_en = Property energy ((structure: ag), (calculator: calc))
view eos (ag.cell_volume:array, prop_en.energy:array, 'sj')

Visualize band structure

The band structure can be visualized using this syntax:

view bs (<table of band structure>, <table of plot parameters>)

Example:

cu = Structure ((atoms: ((symbols: 'Cu'), (x: 0.) [angstrom], (y: 0.) [angstrom], (z: 0.) [angstrom])),
                (cell: [[0.0, 1.805, 1.805], [1.805, 0.0, 1.805], [1.805, 1.805, 0.0]] [angstrom]),
                (pbc: [true, true, true]))
calc = Calculator free_electrons ((nvalence: 1), (kpts: ((path: 'GXWLGK'), (npoints: 200))))
algo = Algorithm BandStructure ()
prop = Property band_structure ((structure: cu), (calculator: calc), (algorithm: algo))
view bs (prop.band_structure[0], ((emin: 0) [eV], (emax: 20) [eV]))