Lattice="R1x R1y R1z R2x R2y R2z R3x R3y R3z"
where \(R1x R1y R1z\) are the Cartesian x-, y- and z-components of the
first lattice vector (\(\mathbf{a}\)), \(R2x R2y R2z\) those of the second
lattice vector (\(\mathbf{b}\)) and \(R3x R3y R3z\) those of the
third lattice vector (\(\mathbf{c}\)).
The list of properties in the file is described by the \(Properties\)
parameter, which should take the form of a series of colon separated
triplets giving the name, format (\(R\) for real, \(I\) for integer) and
number of columns of each property. For example:
Properties="species:S:1:pos:R:3:vel:R:3:select:I:1"
indicates the first column represents atomic species, the next three
columns represent atomic positions, the next three velcoities, and the
last is an single integer called \(select\). With this property
definition, the line
Si 4.08000000 4.08000000 1.36000000 0.00000000 0.00000000 0.00000000 1
would describe a silicon atom at position (4.08,4.08,1.36) with zero
velocity and the \(select\) property set to 1.
The property names \(pos\), \(Z\), \(mass\), and \(charge\) map to ASE
ase.atoms.Atoms.arrays
entries named
\(positions\), \(numbers\), \(masses\) and \(charges\) respectively.
Additional key-value pairs in the comment line are parsed into the
ase.Atoms.atoms.info
dictionary, with the following conventions
Values can be quoted with \(""\), \(''\), \([]\) or \({}\) (the latter are
included to ease command-line usage as the \({}\) are not treated
specially by the shell)
Quotes within keys or values can be escaped with \(\"\).
Keys with special names \(stress\) or \(virial\) are treated as 3x3 matrices
in Fortran order, as for \(Lattice\) above.
Otherwise, values with multiple elements are treated as 1D arrays, first
assuming integer format and falling back to float if conversion is
unsuccessful.
A missing value defaults to \(True\), e.g. the comment line
\("cutoff=3.4 have_energy"\) leads to
\({'cutoff': 3.4, 'have_energy': True}\) in \(atoms.info\).
Value strings starting with \("_JSON"\) are interpreted as JSON content;
similarly, when writing, anything which does not match the criteria above
is serialised as JSON.
The extended XYZ format is also supported by the
the Ovito visualisation tool
(from v2.4 beta
onwards).
ase.io.extxyz.write_extxyz(fileobj, images, comment='', columns=None, write_info=True, write_results=True, plain=False, vec_cell=False, append=False)
Write output in extended XYZ format
Optionally, specify which columns (arrays) to include in output,
whether to write the contents of the \(atoms.info\) dict to the
XYZ comment line (default is True), the results of any
calculator attached to this Atoms. The \(plain\) argument
can be used to write a simple XYZ file with no additional information.
\(vec_cell\) can be used to write the cell vectors as additional
pseudo-atoms. If \(append\) is set to True, the file is for append (mode \(a\)),
otherwise it is overwritten (mode \(w\)).
See documentation for read_xyz()
for further details of the extended
XYZ file format.
ase.io.gaussian.read_gaussian_in(fd, attach_calculator=False)[source]
Reads a gaussian input file and returns an Atoms object.
Parameters:
fd (file-like) – Contains the contents of a gaussian input file
attach_calculator (bool) – When set to True
, a Gaussian calculator will be
attached to the Atoms object which is returned.
This will mean that additional information is read
from the input file (see below for details).
Returns:
An Atoms object containing the following information that has been
read from the input file: symbols, positions, cell.
Return type:
Atoms
Notes
Gaussian input files can be read where the atoms’ locations are set with
cartesian coordinates or as a z-matrix. Variables may be used in the
z-matrix definition, but currently setting constants for constraining
a geometry optimisation is not supported. Note that the \(alternative\)
z-matrix definition where two bond angles are set instead of a bond angle
and a dihedral angle is not currently supported.
If the parameter attach_calculator
is set to True
, then the Atoms
object is returned with a Gaussian calculator attached.
This Gaussian calculator will contain a parameters dictionary which will
contain the Link 0 commands and the options and keywords set in the route
section of the Gaussian input file, as well as:
The charge and multiplicity
The selected level of output
The method, basis set and (optional) fitting basis set.
Basis file name/definition if set
Nuclear properties
Masses set in the nuclear properties section or in the ReadIsotopes
section (if freq=ReadIso
is set). These will be saved in an array
with keyword isolist
, in the parameters dictionary. For these to
actually be used in calculations and/or written out to input files,
the Atoms object’s masses must be manually updated to these values
(this is not done automatically)
If the Gaussian input file has been generated by ASE’s
write_gaussian_in
method, then the basis set, method and fitting
basis will be saved under the basis
, method
and fitting_basis
keywords in the parameters dictionary. Otherwise they are saved as
keywords in the parameters dict.
Currently this does not support reading of any other sections which would
be found below the molecule specification that have not been mentioned
here (such as the ModRedundant section).
ase.io.gaussian.write_gaussian_in(fd, atoms, properties=['energy'], method=None, basis=None, fitting_basis=None, output_type='P', basisfile=None, basis_set=None, xc=None, charge=None, mult=None, extra=None, ioplist=None, addsec=None, spinlist=None, zefflist=None, qmomlist=None, nmagmlist=None, znuclist=None, radnuclearlist=None, **params)[source]
Generates a Gaussian input file
Parameters:
fd (file-like) – where the Gaussian input file will be written
atoms (Atoms) – Structure to write to the input file
properties (list) – Properties to calculate
method (str) – Level of theory to use, e.g. hf
, ccsd
, mp2
, or b3lyp
.
Overrides xc
(see below).
xc (str) – Level of theory to use. Translates several XC functionals from
their common name (e.g. PBE
) to their internal Gaussian name
(e.g. PBEPBE
).
basis (str) – The basis set to use. If not provided, no basis set will be requested,
which usually results in STO-3G
. Maybe omitted if basisfile is set
(see below).
fitting_basis (str) – The name of the fitting basis set to use.
output_type (str) – Level of output to record in the Gaussian
output file - this may be N
- normal or P
-
additional.
basisfile (str) – The name of the basis file to use. If a value is provided, basis may
be omitted (it will be automatically set to ‘gen’)
basis_set (str) – The basis set definition to use. This is an alternative
to basisfile, and would be the same as the contents
of such a file.
charge (int) – The system charge. If not provided, it will be automatically
determined from the Atoms
object’s initial_charges.
mult (int) – The system multiplicity (spin + 1
). If not provided, it will be
automatically determined from the Atoms
object’s
initial_magnetic_moments
.
extra (str) – Extra lines to be included in the route section verbatim.
It should not be necessary to use this, but it is included for
backwards compatibility.
ioplist (list) – A collection of IOPs definitions to be included in the route line.
addsec (str) – Text to be added after the molecular geometry specification, e.g. for
defining masses with freq=ReadIso
.
spinlist (list) – A list of nuclear spins to be added into the nuclear
propeties section of the molecule specification.
zefflist (list) – A list of effective charges to be added into the nuclear
propeties section of the molecule specification.
qmomlist (list) – A list of nuclear quadropole moments to be added into
the nuclear propeties section of the molecule
specification.
nmagmlist (list) – A list of nuclear magnetic moments to be added into
the nuclear propeties section of the molecule
specification.
znuclist (list) – A list of nuclear charges to be added into the nuclear
propeties section of the molecule specification.
radnuclearlist (list) – A list of nuclear radii to be added into the nuclear
propeties section of the molecule specification.
params (dict) – Contains any extra keywords and values that will be included in either
the link0 section or route section of the gaussian input file.
To be included in the link0 section, the keyword must be one of the
following: mem
, chk
, oldchk
, schk
, rwf
,
oldmatrix
, oldrawmatrix
, int
, d2e
, save
,
nosave
, errorsave
, cpu
, nprocshared
, gpucpu
,
lindaworkers
, usessh
, ssh
, debuglinda
.
Any other keywords will be placed (along with their values) in the
route section.
ase.io.gen.read_gen(fileobj)[source]
Read structure in GEN format (refer to DFTB+ manual).
Multiple snapshot are not allowed.
ase.io.gen.write_gen(fileobj, images: Union[Atoms, Sequence[Atoms]], fractional
: bool = False)[source]
Write structure in GEN format (refer to DFTB+ manual).
Multiple snapshots are not allowed.
ase.io.gpumd.read_gpumd(fd, species=None, isotope_masses=None)[source]
Read Atoms object from a GPUMD structure input file
Parameters:
fd (file | str) – File object or name of file from which to read the Atoms object
species (List[str]) – List with the chemical symbols that correspond to each type, will take
precedence over isotope_masses
isotope_masses (Dict[str, List[float]]) – Dictionary with chemical symbols and lists of the associated atomic
masses, which is used to identify the chemical symbols that correspond
to the types not found in species_types. The default is to find the
closest match ase.data.atomic_masses
.
Returns:
atoms – Atoms object
Return type:
Atoms
Raises:
ValueError – Raised if the list of species is incompatible with the input file
ase.io.gpumd.write_gpumd(fd, atoms, maximum_neighbors=None, cutoff=None, groupings=None, use_triclinic=False, species=None)[source]
Writes atoms into GPUMD input format.
Parameters:
fd (file) – File like object to which the atoms object should be written
atoms (Atoms) – Input structure
maximum_neighbors (int) – Maximum number of neighbors any atom can ever have (not relevant when
using force constant potentials)
cutoff (float) – Initial cutoff distance used for building the neighbor list (not
relevant when using force constant potentials)
groupings (list[list[list[int]]]) – Groups into which the individual atoms should be divided in the form of
a list of list of lists. Specifically, the outer list corresponds to
the grouping methods, of which there can be three at the most, which
contains a list of groups in the form of lists of site indices. The
sum of the lengths of the latter must be the same as the total number
of atoms.
use_triclinic (bool) – Use format for triclinic cells
species (List[str]) – GPUMD uses integers to define atom types. This list allows customized
such definitions (e.g, [‘Pd’, ‘H’] means Pd is type 0 and H type 1).
If None, this list is built by assigning each distinct species to
an integer in the order of appearance in \(atoms\).
Raises:
ValueError – Raised if parameters are incompatible
ase.io.gromacs.read_gromacs(fd)[source]
From:
http://manual.gromacs.org/current/online/gro.html
C format
“%5d%-5s%5s%5d%8.3f%8.3f%8.3f%8.4f%8.4f%8.4f”
python: starting from 0, including first excluding last
0:4 5:10 10:15 15:20 20:28 28:36 36:44 44:52 52:60 60:68
Import gromacs geometry type files (.gro).
Reads atom positions,
velocities(if present) and
simulation cell (if present)
ase.io.gromacs.write_gromacs(fileobj, atoms)[source]
Write gromacs geometry files (.gro).
Writes:
atom positions,
velocities (if present, otherwise 0)
simulation cell (if present)
ase.io.gromos.read_gromos(fileobj)[source]
Read gromos geometry files (.g96).
Reads:
atom positions,
and simulation cell (if present)
tries to set atom types
ase.io.gromos.write_gromos(fileobj, atoms)[source]
Write gromos geometry files (.g96).
Writes:
atom positions,
and simulation cell (if present)
ase.io.lammpsdata.read_lammps_data(fileobj, Z_of_type: dict = None, sort_by_id: bool = True, read_image_flags: bool = True, units: str = 'metal', atom_style: str = None, style: str = None)[source]
Method which reads a LAMMPS data file.
Parameters:
fileobj (file | str) – File from which data should be read.
Z_of_type (dict[int, int], optional) – Mapping from LAMMPS atom types (typically starting from 1) to atomic
numbers. If None, if there is the “Masses” section, atomic numbers are
guessed from the atomic masses. Otherwise, atomic numbers of 1 (H), 2
(He), etc. are assigned to atom types of 1, 2, etc. Default is None.
sort_by_id (bool, optional) – Order the particles according to their id. Might be faster to set it
False. Default is True.
read_image_flags (bool, default True) – If True, the lattice translation vectors derived from image flags are
added to atomic positions.
units (str, optional) – LAMMPS units. Default is
“metal”.
atom_style ({"atomic", "charge", "full"} etc., optional) – LAMMPS atom style.
If None, \(atom_style\) is guessed in the following priority (1) comment
after \(Atoms\) (2) length of fields (valid only \(atomic\) and \(full\)).
Default is None.
ase.io.lammpsdata.write_lammps_data(fd, atoms: Atoms, *
, specorder: list = None, reduce_cell: bool = False, force_skew: bool = False, prismobj: Prism = None, write_image_flags: bool = False, masses: bool = False, velocities: bool = False, units: str = 'metal', atom_style: str = 'atomic')[source]
Write atomic structure data to a LAMMPS data file.
Parameters:
fd (file|str) – File to which the output will be written.
atoms (Atoms) – Atoms to be written.
specorder (list[str], optional) – Chemical symbols in the order of LAMMPS atom types, by default None
force_skew (bool, optional) – Force to write the cell as a
triclinic box,
by default False
reduce_cell (bool, optional) – Whether the cell shape is reduced or not, by default False
prismobj (Prism|None, optional) – Prism, by default None
write_image_flags (bool, default False) – If True, the image flags, i.e., in which images of the periodic
simulation box the atoms are, are written.
masses (bool, optional) – Whether the atomic masses are written or not, by default False
velocities (bool, optional) – Whether the atomic velocities are written or not, by default False
units (str, optional) – LAMMPS units,
by default “metal”
atom_style ({"atomic", "charge", "full"}, optional) – LAMMPS atom style,
by default “atomic”
ase.io.lammpsrun.read_lammps_dump_binary(fileobj, index=-1, colnames=None, intformat='SMALLBIG', **kwargs)[source]
Read binary dump-files (after binary2txt.cpp from lammps/tools)
Parameters:
fileobj – file-stream containing the binary lammps data
index – integer or slice object (default: get the last timestep)
colnames – data is columns and identified by a header
intformat – lammps support different integer size. Parameter set at compile-time and can unfortunately not derived from data file
Returns:
list of Atoms-objects
Return type:
ase.io.lammpsrun.read_lammps_dump_text(fileobj, index=-1, **kwargs)[source]
Process cleartext lammps dumpfiles
Parameters:
fileobj – filestream providing the trajectory data
index – integer or slice object (default: get the last timestep)
Returns:
list of Atoms objects
Return type:
ase.io.magres.write_magres(fd, image)[source]
A writing function for magres files. Two steps: first data are arranged
into structures, then dumped to the actual file
ase.io.mustem.read_mustem(fd)[source]
Import muSTEM input file.
Reads cell, atom positions, etc. from muSTEM xtl file.
The mustem xtl save the root mean square (RMS) displacement, which is
convert to Debye-Waller (in Ų) factor by:
\[B = RMS * 8\pi^2\]
ase.io.mustem.write_mustem(fd,
*args,
**kwargs)[source]
Write muSTEM input file.
Parameters:
atoms: Atoms object
keV: floatEnergy of the electron beam in keV required for the image simulation.
debye_waller_factors: float or dictionary of float with atom type as keyDebye-Waller factor of each atoms. Since the prismatic/computem
software use root means square RMS) displacements, the Debye-Waller
factors (B) needs to be provided in Ų and these values are converted
to RMS displacement by:
\[RMS = \frac{B}{8\pi^2}\]
occupancies: float or dictionary of float with atom type as key (optional)Occupancy of each atoms. Default value is \(1.0\).
comment: str (optional)Comments to be written in the first line of the file. If not
provided, write the total number of atoms and the chemical formula.
fit_cell_to_atoms: bool (optional)If \(True\), fit the cell to the atoms positions. If negative coordinates
are present in the cell, the atoms are translated, so that all
positions are positive. If \(False\) (default), the atoms positions and
the cell are unchanged.
ase.io.nwchem.write_nwchem_in(fd, atoms, properties=None, echo=False, **params)[source]
Writes NWChem input file.
See NWChem
for available params.
Parameters:
fd – file descriptor
atoms – atomic configuration
properties – list of properties to compute; by default only the
calculation of the energy is requested
echo – if True include the \(echo\) keyword at the top of the file,
which causes the content of the input file to be included
in the output file
params – dict of instructions blocks to be included
ase.io.nwchem.read_nwchem_out(fobj, index=-1)[source]
Splits an NWChem output file into chunks corresponding to
individual single point calculations.
ase.io.onetep.read_onetep_in(fd, **kwargs)[source]
Read a single ONETEP input.
This function can be used to visually check ONETEP inputs,
using the ase gui. It can also be used to get access to
the input parameters attached to the ONETEP calculator
returned
The function should work on inputs which contain
‘include_file’ command(s), (possibly recursively
but untested)
The function should work on input which contains
exotic element(s) name(s) if the specie block is
present to map them back to real element(s)
Parameters:
fd (io-object) – File to read.
Returns:
structure – Atoms object with cell and a Onetep calculator
attached which contains the keywords dictionary
Return type:
Atoms
ase.io.onetep.write_onetep_in(fd, atoms, edft=False, xc='PBE', ngwf_count=-1, ngwf_radius=9.0, keywords={}, pseudopotentials={}, pseudo_path='.', pseudo_suffix=None, **kwargs)[source]
Write a single ONETEP input.
This function will be used by ASE to perform
various workflows (Opt, NEB…) or can be used
manually to quickly create ONETEP input file(s).
The function will first write keywords in
alphabetic order in lowercase. Secondly, blocks
will be written in alphabetic order in uppercase.
Two ways to work with the function:
By providing only (simple) keywords present in
the parameters. ngwf_count and ngwf_radius
accept multiple types as described in the Parameters
section.
If the keywords parameters is provided as a dictionary
these keywords will be used to write the input file and
will take priority.
If no pseudopotentials are provided in the parameters and
the function will try to look for suitable pseudopotential
in the pseudo_path.
Parameters:
fd (file) – File to write.
atoms (Atoms) – Atoms including Cell object to write.
edft (Bool) – Activate EDFT.
xc (str) – DFT xc to use e.g (PBE, RPBE, …)
ngwf_count (int|list|dict) –
- Behaviour depends on the type:
int: every species will have this amount
of ngwfs.
list: list of int, will be attributed
alphabetically to species:
dict: keys are species name(s),
value are their number:
ngwf_radius (int|list|dict) –
- Behaviour depends on the type:
float: every species will have this radius.
list: list of float, will be attributed
alphabetically to species:
[10.0, 9.0]
dict: keys are species name(s),
value are their radius:
{‘Na’: 9.0, ‘Cl’: 10.0}
keywords (dict) – Dictionary with ONETEP keywords to write,
keywords with lists as values will be
treated like blocks, with each element
of list being a different line.
pseudopotentials (dict) –
- Behaviour depends on the type:
keys are species name(s) their
value are the pseudopotential file to use:
{‘Na’: ‘Na.usp’, ‘Cl’: ‘Cl.usp’}
pseudo_path (str) – Where to look for pseudopotential, correspond
to the pseudo_path keyword of ONETEP.
pseudo_suffix (str) – Suffix for the pseudopotential filename
to look for, useful if you have multiple sets of
pseudopotentials in pseudo_path.
ase.io.onetep.read_onetep_out(fd, index=-1, improving=False, **kwargs)[source]
Read ONETEP output(s).
This function will be used by ASE when performing
various workflows (Opt, NEB…)
- Parameters:
fd (file) – File to read.
index (slice) – Which atomic configuration to read
improving (Bool) – If the output is a geometry optimisation,
improving = True will keep line search
configuration from BFGS
- Yields:
structure (Atoms|list of Atoms)
ase.io.prismatic.read_prismatic(fd)[source]
Import prismatic and computem xyz input file as an Atoms object.
Reads cell, atom positions, occupancies and Debye Waller factor.
The occupancy values and the Debye Waller factors are obtained using the
\(get_array\) method and the \(occupancies\) and \(debye_waller_factors\) keys,
respectively. The root means square (RMS) values from the
prismatic/computem xyz file are converted to Debye-Waller factors (B) in Ų
\[B = RMS^2 * 8\pi^2\]
ase.io.prismatic.write_prismatic(fd,
*args,
**kwargs)[source]
Write xyz input file for the prismatic and computem software. The cell
needs to be orthorhombric. If the cell contains the \(occupancies\) and
\(debye_waller_factors\) arrays (see the \(set_array\) method to set them),
these array will be written to the file.
If the occupancies is not specified, the default value will be set to 0.
Parameters:
atoms: Atoms object
debye_waller_factors: float or dictionary of float or None (optional).If float, set this value to all
atoms. If dictionary, each atom needs to specified with the symbol as
key and the corresponding Debye-Waller factor as value.
If None, the \(debye_waller_factors\) array of the Atoms object needs to
be set (see the \(set_array\) method to set them), otherwise raise a
ValueError. Since the prismatic/computem software use root means square
(RMS) displacements, the Debye-Waller factors (B) needs to be provided
in Ų and these values are converted to RMS displacement by:
\[RMS = \sqrt{\frac{B}{8\pi^2}}\]
Default is None.
comment: str (optional)Comments to be written in the first line of the file. If not
provided, write the total number of atoms and the chemical formula.
ase.io.proteindatabank.read_proteindatabank(fileobj, index=-1, read_arrays=True)[source]
Read PDB files.
ase.io.proteindatabank.write_proteindatabank(fileobj, images, write_arrays=True)[source]
Write images to PDB-file.
Read data from QBox output file
Inputs:f - str or fileobj, path to file or file object to read from
index - int or slice, which frames to return
Returns:
list of Atoms or atoms, requested frame(s)
ase.io.res.read_res(filename, index=-1)[source]
Read input in SHELX (.res) format
Multiple frames are read if \(filename\) contains a wildcard character,
e.g. \(file_*.res\). \(index\) specifes which frames to retun: default is
last frame only (index=-1).
ase.io.res.write_res(filename, images, write_info=True, write_results=True, significant_figures=6)[source]
Write output in SHELX (.res) format
To write multiple images, include a % format string in filename,
e.g. \(file_%03d.res\).
Optionally include contents of Atoms.info dictionary if \(write_info\)
is True, and/or results from attached calculator if \(write_results\)
is True (only energy results are supported).
ase.io.rmc6f.read_rmc6f(filename, atom_type_map=None)[source]
Parse a RMCProfile rmc6f file into ASE Atoms object
Parameters:
filename (file|str) – A file like object or filename.
atom_type_map (dict{str:str}) –
Map of atom types for conversions. Mainly used if there is
an atom type in the file that is not supported by ASE but
want to map to a supported atom type instead.
Example to map deuterium to hydrogen:
atom_type_map = { ‘D’: ‘H’ }
Returns:
structure – The Atoms object read in from the rmc6f file.
Return type:
Atoms
ase.io.rmc6f.write_rmc6f(filename, atoms, order=None, atom_type_map=None)[source]
Write output in rmc6f format - RMCProfile v6 fractional coordinates
Parameters:
filename (file|str) – A file like object or filename.
atoms (Atoms object) – The Atoms object to be written.
order (list[str]) – If not None, gives a list of atom types for the order
to write out each.
atom_type_map (dict{str:str}) –
Map of atom types for conversions. Mainly used if there is
an atom type in the Atoms object that is a placeholder
for a different atom type. This is used when the atom type
is not supported by ASE but is in RMCProfile.
Example to map hydrogen to deuterium:
atom_type_map = { ‘H’: ‘D’ }
ase.io.turbomole.read_turbomole(fd)[source]
Method to read turbomole coord file
coords in bohr, atom types in lowercase, format:
$coord
x y z atomtype
x y z atomtype f
Above ‘f’ means a fixed atom.
ase.io.v_sim.read_v_sim(fd)[source]
Import V_Sim input file.
Reads cell, atom positions, etc. from v_sim ascii file.
V_sim format is specified here:
https://l_sim.gitlab.io/v_sim/sample.html#sample_ascii
ase.io.v_sim.write_v_sim(fd, atoms)[source]
Write V_Sim input file.
Writes the atom positions and unit cell.
V_sim format is specified here:
https://l_sim.gitlab.io/v_sim/sample.html#sample_ascii
ase.io.vasp.read_vasp(filename='CONTCAR')[source]
Import POSCAR/CONTCAR type file.
Reads unitcell, atom positions and constraints from the POSCAR/CONTCAR
file and tries to read atom types from POSCAR/CONTCAR header, if this fails
the atom types are read from OUTCAR or POTCAR file.
ase.io.vasp.write_vasp(filename, atoms, label=None, direct=False, sort=None, symbol_count=None, long_format=True, vasp5=True, ignore_constraints=False, wrap=False)[source]
Method to write VASP position (POSCAR/CONTCAR) files.
Writes label, scalefactor, unitcell, # of various kinds of atoms,
positions in cartesian or scaled coordinates (Direct), and constraints
to file. Cartesian coordinates is default and default label is the
atomic species, e.g. ‘C N H Cu’.
ase.io.vasp.read_vasp_out(filename='OUTCAR', index=-1)[source]
Import OUTCAR type file.
Reads unitcell, atom positions, energies, and forces from the OUTCAR file
and attempts to read constraints (if any) from CONTCAR/POSCAR, if present.
ase.io.vasp.read_vasp_xdatcar(filename='XDATCAR', index=-1)[source]
Import XDATCAR file.
Parameters:
index (int or slice or str) – Which frame(s) to read. The default is -1 (last frame).
See ase.io.read()
for details.
Notes
Constraints ARE NOT stored in the XDATCAR, and as such, Atoms objects
retrieved from the XDATCAR will not have constraints.
ase.io.vasp.write_vasp_xdatcar(fd, images, label=None)[source]
Write VASP MD trajectory (XDATCAR) file
Only Vasp 5 format is supported (for consistency with read_vasp_xdatcar)
Parameters:
fd (str, fp) – Output file
images (iterable of Atoms) – Atoms images to write. These must have
consistent atom order and lattice vectors - this will not be
checked.
label (str) – Text for first line of file. If empty, default to list of
elements.
ase.io.vasp.read_vasp_xml(filename='vasprun.xml', index=-1)[source]
Parse vasprun.xml file.
Reads unit cell, atom positions, energies, forces, and constraints
from vasprun.xml file
Examples
>>> import ase.io
>>> ase.io.write("out.traj", ase.io.read("vasprun.xml", index=":"))
ase.io.wannier90.read_wout(fileobj: IO[str], include_wannier_function_centers: bool = True) → Atoms[source]
Read atoms and wannier function centers (as symbol X).
ase.io.x3d.write_x3d(fd, atoms, format='X3D', style=None)[source]
Writes to html using X3DOM.
Parameters:
object (filename or output file) –
object –
rendered (atoms - Atoms object to be) –
str (format -) – to be readable by Blender. \(None\) to detect format based on file
extension (‘.html’ -> ‘X3DOM’, ‘.x3d’ -> ‘X3D’)
'X3D' (either 'X3DOM' for web-browser compatibility or) – to be readable by Blender. \(None\) to detect format based on file
extension (‘.html’ -> ‘X3DOM’, ‘.x3d’ -> ‘X3D’)
dict (style -) –
element (css style attributes for the X3D) –
ase.io.xsd.write_xsd(fd, images, connectivity=None)[source]
Takes Atoms object, and write materials studio file
atoms: Atoms object
filename: path of the output file
connectivity: number of atoms by number of atoms matrix for connectivity
between atoms (0 not connected, 1 connected)
note: material studio file cannot use a partial periodic system. If partial
perodic system was inputted, full periodicity was assumed.
ase.io.xtd.read_xtd(filename, index=-1)[source]
Import xtd file (Materials Studio)
Xtd files always come with arc file, and arc file
contains all the relevant information to make atoms
so only Arc file needs to be read
ase.io.xtd.write_xtd(filename, images, connectivity=None, moviespeed=10)[source]
Takes Atoms object, and write materials studio file
atoms: Atoms object
filename: path of the output file
moviespeed: speed of animation. between 0 and 10
note: material studio file cannot use a partial periodic system. If partial
perodic system was inputted, full periodicity was assumed.