Source documentation

muesr – Main package

muesr.core – Core functionalities

muesr.core.sample – The Sample definition

class muesr.core.sample.Sample[source]

Bases: object

This object contains all the information about the sample under investigation. It includes a definition of the lattice parameters and the atomic positions, the symmetry, the positions of the muons and the magnetic structure.

In each sample there can only be a single lattice structure. Therefore, also a single symmetry is defined.

Multiple magnetic structures can be defined instead. Every time a magnetic structure is defined, it automatically becomes the current magnetic structure. To select a different magnetic structure the user can use the property current_mm_idx.

>>> yoursample = Sample()  #initialize the sample
>>> yoursample.name = "Test"
add_muon(position, cartesian=False)[source]

Adds a muon position.

Parameters:
  • position (list) – list of three float or numpy array of shape (3,) describing the position of the muon.
  • cartesian (bool) – if True, the position os assumed to be in Cartesian coordinates. If False (default) the position is assumed to be in fractional coordinates.
Returns:

None

Return type:

None

Raises:

MuonError

cell

Returns the atomic structure definition, i.e. an Atoms object.

Getter:returns a Atoms Object.
Setter:Sets the atomic structure definition from a Atoms object.
Type:int
Raises:TypeError, CellError
check_status(cell=False, magdefs=False, muon=False, sym=False)[source]
current_mm_idx

Index of the currently selected magnetic model (starting from 0)

Getter:returns the index of the currently selected Magnetic Model
Setter:Sets the current Magnetic Model from the index.
Type:int
Raises:MagDefError
mm

Current Magnetic Model.

Getter:Returns the current Magnetic Model (a MM object)
Setter:Adds and select a new Magnetic Model (MM object)
Type:MM object
Raises:TypeError, MagDefError
mm_count

Count the loaded Magnetic Model(s).

Getter:Returns the current Magnetic Model (a MM object)
Setter:Read-only
Type:MM object
Raises:None
muons

Muon positions in fractional coordinates.

Getter:Returns a list of numpy array of shape (3,).
name

The sample name.

Getter:Returns the sample name
Setter:Sets the sample name
Type:str
new_mm()[source]

Creates a new empty magnetic model (everything is set to zero)

Returns:None
Return type:None
new_smm(symbolic_parameters)[source]

Creates a new empty Symbolic Magnetic Model (everything is set to zero)

Parameters:symbolic_parameters (str) – String containing a list of symbolic parameters. Example “x,y,z”
Returns:None
Return type:None
Raises:ImportError if sympy is not installed. CellError if lattice is not defined.
sym

Symmetry of the system, i.e. a Spacegroup Object.

Getter:returns a Spacegroup Object.
Setter:Sets the symmetry of the system from a Spacegroup Object.
Type:int
Raises:TypeError, SymError

muesr.core.sampleErrors – Sample related exceptions

exception muesr.core.sampleErrors.CellError[source]

Bases: muesr.core.sampleErrors.SampleException

Lattice structure needed but not defined.

exception muesr.core.sampleErrors.MagDefError[source]

Bases: muesr.core.sampleErrors.SampleException

Magnetic structure not defined or incompatible.

exception muesr.core.sampleErrors.MuonError[source]

Bases: muesr.core.sampleErrors.SampleException

Muon position needed but not defined.

exception muesr.core.sampleErrors.SampleException[source]

Bases: Exception

Exceptions connected with problems with the Sample definition.

exception muesr.core.sampleErrors.SymmetryError[source]

Bases: muesr.core.sampleErrors.SampleException

Symmetry needed but not defined.

muesr.core.atoms – Defined the unit cell

class muesr.core.atoms.Atoms(symbols=None, positions=None, numbers=None, masses=None, magmoms=None, scaled_positions=None, cell=None, pbc=None)[source]

Bases: object

Atoms class compatible with the ASE Atoms class Only stuff needed by muesr is implemented.

>>> a = 4.05  # Gold lattice constant
>>> b = a / 2
>>> fcc = Atoms('Au',
...             cell=[(0, b, b), (b, 0, b), (b, b, 0)],
...             pbc=True)
del_atom(i)[source]

Removes one atom

edit_atom(i, symbol=None, number=None, magmom=None, mass=None)[source]
extend(symbol=None, position=None, number=None, mass=None, magmom=None, scaled_position=None)[source]

Extends the current atomic structure by one atom.

get_atomic_numbers()[source]
get_cell()[source]
get_chemical_symbols()[source]
get_magnetic_moments()[source]

Get magnetic moments if set. None if nothing set.

get_masses()[source]
get_number_of_atoms()[source]
get_positions()[source]
get_scaled_positions()[source]
get_volume()[source]
numbers_to_symbols()[source]
set_cell(cell)[source]

Set lattice vectors parameters.

Parameters:cell

numpy array of the form

[[v1(1), v1(2), v1(3)],
v2(1), v2(2), v2(3)],
v3(1), v3(2), v3(3) ]]

where v1, v2, v3 are the lattice vectors.

set_chemical_symbols(symbols)[source]
set_magnetic_moments(magmoms=None)[source]
set_masses(masses)[source]
set_positions(cart_positions)[source]
set_scaled_positions(scaled_positions)[source]
symbols_to_masses()[source]
symbols_to_numbers()[source]

muesr.core.cells – Functions for lattice

muesr.core.cells.R2S(x, y, z, H, K, L)[source]

Given reciprocal-space coordinates of a vecotre, calculate its coordinates in the Cartesian space.

muesr.core.cells.S2R(qx, qy, qz, cell)[source]

Given cartesian coordinates of a vector in the S System, calculate its Miller indexes.

muesr.core.cells.get_Delaunay_reduction(lattice, tolerance)[source]
muesr.core.cells.get_angles(lattice)[source]

Get alpha, beta and gamma angles from lattice vectors.

>>> get_angles( np.diag([1,2,3]) )
(90.0, 90.0, 90.0)
muesr.core.cells.get_atom_distance(scaled_pos, center, cell, tolerance=1e-05)[source]

Return the shortest distance to atom from specified position

muesr.core.cells.get_atom_vec(scaled_pos, center, cell, tolerance=1e-05)[source]

Return the shortest distance to atom from specified position

muesr.core.cells.get_cell_matrix(a, b, c, alpha, beta, gamma)[source]
muesr.core.cells.get_cell_parameters(lattice)[source]
muesr.core.cells.get_distance(atoms, a0, a1, tolerance=1e-05)[source]

Return the shortest distance between a pair of atoms in PBC

muesr.core.cells.get_distance_with_center(atoms, center, atom_num, tolerance=1e-05)[source]

Return the shortest distance to atom from specified position

muesr.core.cells.get_reciprocal_lattice(lattice)[source]
muesr.core.cells.get_reduced_bases(cell, tolerance=1e-05)[source]

This is an implementation of Delaunay reduction. Some information is found in International table.

muesr.core.cells.get_shortest_bases_from_extented_bases(extended_bases, tolerance)[source]
muesr.core.cells.get_simple_supercell(sample, multi)[source]

This function creates a simple supercell by expanding the unit cell in a, b and c directions.

muesr.core.cells.gtensor(lattice)[source]

calculates the metric tensor of a lattice

muesr.core.cells.print_cell(cell, mapping=None)[source]

Print lattice structure.

cell : Aroms instance, the atomic structure to be printed mapping : list or None, can be used to show relations between atoms.

muesr.core.cells.reciprocate(x, y, z, cell)[source]

calculate miller indexes of a vector defined by its fractional cell coords

muesr.core.cells.reduce_bases(extended_bases, tolerance)[source]

muesr.core.magmodel – Magnetic Model

class muesr.core.magmodel.MM(cell_size, latt_vects=None)[source]

Bases: object

Magnetic model class.

The magnetic structures are defined in terms of a single propagation vector, complex Fourier components and phases.

Parameters:
  • cell_size (int) – number of Fourier components=number of atoms
  • latt_vects – lattice vectors as numpy 3x3 ndarray. If None, the Fourier componentscan only be specified in cartesian coordinates.
Raises:

TypeError

desc

Description of the magnetic structure

Getter:Returns the description
Setter:Sets the description
Type:str
fc

Fourier components in Cartesian coordinates. Same as fcCart

fcCart

Get Fourier components in Cartesian coordinates.

Getter:Returns a numpy array of size (size,3) with the fourier components.
Setter:Sets the fourier compinents
Type:numpy ndarray of size (size,3)
fcLattBM

Fourier components in Bohr Magneton units, with x||a, y||b and z||c

Getter:Returns a numpy array of size (size,3) with the fourier components.
Setter:Sets the fourier compinents
Type:numpy ndarray of size (size,3)
fcLattBMA

Get Fourier components in Bohr Magneton/Angstrom units, with x||a, y||b and z||c

Getter:Returns a numpy array of size (size,3) with the fourier components.
Setter:Sets the fourier compinents
Type:numpy ndarray of size (size,3)
fc_get(coord_system=0)[source]

Retrives Fourier components.

Params int coord_system:
 

requested coordinate system.

  • 0 means Cartesian coordinates
  • 1 means Lattice coordinates (values in units of Bohr Magneton/Angstrom)
  • 2 means Bohr magnetons in each lattice direction (values in units Bohr Magneton)

See http://magcryst.org/resources/magnetic-coordinates/

Raises:

ValueError

fc_set(value, coord_system=0)[source]

Sets Fourier components.

Parameters:
  • value – numpy array containing the fourier components for all the atoms.
  • coord_system (int) –

    requested coordinate system.

    • 0 means Cartesian coordinates
    • 1 means Lattice coordinates (values in units of Bohr Magneton/Angstrom)
    • 2 means Bohr magnetons in each lattice direction (values in units Bohr Magneton)

    See http://magcryst.org/resources/magnetic-coordinates/

Raises:

ValueError, TypeError

isSymbolic

True if Fourier components are defined as sympy symbols, False otherwise.

ALWAYS False for MM objects.

k

The propagation vector in reciprocal lattice units.

Getter:Returns a numpy array of shape (3,) with the propagation vector.
Setter:Sets the propagation vector.
Type:numpy ndarray of shape (3,)
lattice_params

Lattice parameters used to convert the Fourier components to the various coordinates systems.

phi

The phase for each fourier component in units of 2 PI.

Getter:Returns a numpy array of size size with the phases.
Setter:Sets the phases. Type can be list of numpy array.
Type:numpy ndarray
size

Number of Fourier components.

muesr.core.spg – Spacegroups (from ASE)

Definition of the Spacegroup class.

This module only depends on NumPy and the space group database.

class muesr.core.spg.Spacegroup(spacegroup, setting=1, datafile=None)[source]

Bases: object

A space group class.

The instances of Spacegroup describes the symmetry operations for the given space group.

Example:

>>> from muesr.core.spg import Spacegroup
>>>
>>> sg = Spacegroup(225)
>>> print('Space group', sg.no, sg.symbol)
Space group 225 F m -3 m
>>> sg.scaled_primitive_cell
array([[ 0. ,  0.5,  0.5],
       [ 0.5,  0. ,  0.5],
       [ 0.5,  0.5,  0. ]])
>>> sites, kinds = sg.equivalent_sites([[0,0,0]])
>>> sites
array([[ 0. ,  0. ,  0. ],
       [ 0. ,  0.5,  0.5],
       [ 0.5,  0. ,  0.5],
       [ 0.5,  0.5,  0. ]])
centrosymmetric

Whether a center of symmetry exists.

equivalent_lattice_points(uvw)[source]

Return all lattice points equivalent to any of the lattice points in uvw with respect to rotations only.

Only equivalent lattice points that conserves the distance to origo are included in the output (making this a kind of real space version of the equivalent_reflections() method).

Example:

>>> from ase.spacegroup import Spacegroup
>>> sg = Spacegroup(225)  # fcc
>>> sg.equivalent_lattice_points([[0, 0, 2]])
array([[ 0,  0, -2],
       [ 0, -2,  0],
       [-2,  0,  0],
       [ 2,  0,  0],
       [ 0,  2,  0],
       [ 0,  0,  2]])
equivalent_reflections(hkl)[source]

Return all equivalent reflections to the list of Miller indices in hkl.

Example:

>>> from muesr.core.spg import Spacegroup
>>> sg = Spacegroup(225)  # fcc
>>> sg.equivalent_reflections([[0, 0, 2]])
array([[ 0,  0, -2],
       [ 0, -2,  0],
       [-2,  0,  0],
       [ 2,  0,  0],
       [ 0,  2,  0],
       [ 0,  0,  2]])
equivalent_sites(scaled_positions, onduplicates='error', symprec=0.001)[source]

Returns the scaled positions and all their equivalent sites.

Parameters:

scaled_positions: list | array
List of non-equivalent sites given in unit cell coordinates.
onduplicates : ‘keep’ | ‘replace’ | ‘warn’ | ‘error’

Action if scaled_positions contain symmetry-equivalent positions:

‘keep’
ignore additional symmetry-equivalent positions
‘replace’
replace
‘warn’
like ‘keep’, but issue an UserWarning
‘error’
raises a SpacegroupValueError
symprec: float
Minimum “distance” betweed two sites in scaled coordinates before they are counted as the same site.

Returns:

sites: array
A NumPy array of equivalent sites.
kinds: list
A list of integer indices specifying which input site is equivalent to the corresponding returned site.

Example:

>>> from muesr.core.spg import Spacegroup
>>> sg = Spacegroup(225)  # fcc
>>> sites, kinds = sg.equivalent_sites([[0, 0, 0], [0.5, 0.0, 0.0]])
>>> sites
array([[ 0. ,  0. ,  0. ],
       [ 0. ,  0.5,  0.5],
       [ 0.5,  0. ,  0.5],
       [ 0.5,  0.5,  0. ],
       [ 0.5,  0. ,  0. ],
       [ 0. ,  0.5,  0. ],
       [ 0. ,  0. ,  0.5],
       [ 0.5,  0.5,  0.5]])
>>> kinds
[0, 0, 0, 0, 1, 1, 1, 1]
get_op()[source]

Returns all symmetry operations (including inversions and subtranslations), but unlike get_symop(), they are returned as two ndarrays.

get_rotations()[source]

Return all rotations, including inversions for centrosymmetric crystals.

get_symop()[source]

Returns all symmetry operations (including inversions and subtranslations) as a sequence of (rotation, translation) tuples.

lattice

Lattice type:

P primitive I body centering, h+k+l=2n F face centering, h,k,l all odd or even A,B,C single face centering, k+l=2n, h+l=2n, h+k=2n R rhombohedral centering, -h+k+l=3n (obverse); h-k+l=3n (reverse)

no

Space group number in International Tables of Crystallography.

nsubtrans

Number of cell-subtranslation vectors.

nsymop

Total number of symmetry operations.

reciprocal_cell

Tree Miller indices that span all kinematically non-forbidden reflections as a matrix with the Miller indices along the rows.

rotations

Symmetry rotation matrices. The invertions are not included for centrosymmetrical crystals.

scaled_primitive_cell

Primitive cell in scaled coordinates as a matrix with the primitive vectors along the rows.

setting

Space group setting. Either one or two.

subtrans

Translations vectors belonging to cell-sub-translations.

symbol

Hermann-Mauguin (or international) symbol for the space group.

symmetry_normalised_reflections(hkl)[source]

Returns an array of same size as hkl, containing the corresponding symmetry-equivalent reflections of lowest indices.

Example:

>>> from muesr.core.spg import Spacegroup
>>> sg = Spacegroup(225)  # fcc
>>> sg.symmetry_normalised_reflections([[2, 0, 0], [0, 2, 0]])
array([[ 0,  0, -2],
       [ 0,  0, -2]])
symmetry_normalised_sites(scaled_positions, map_to_unitcell=True)[source]

Returns an array of same size as scaled_positions, containing the corresponding symmetry-equivalent sites of lowest indices.

If map_to_unitcell is true, the returned positions are all mapped into the unit cell, i.e. lattice translations are included as symmetry operator.

Example:

>>> from muesr.core.spg import Spacegroup
>>> sg = Spacegroup(225)  # fcc
>>> sg.symmetry_normalised_sites([[0.0, 0.5, 0.5], [1.0, 1.0, 0.0]])
array([[ 0.,  0.,  0.],
       [ 0.,  0.,  0.]])
tag_sites(scaled_positions, symprec=0.001)[source]

Returns an integer array of the same length as scaled_positions, tagging all equivalent atoms with the same index.

Example:

>>> from muesr.core.spg import Spacegroup
>>> sg = Spacegroup(225)  # fcc
>>> sg.tag_sites([[0.0, 0.0, 0.0],
...               [0.5, 0.5, 0.0],
...               [1.0, 0.0, 0.0],
...               [0.5, 0.0, 0.0]])
array([0, 0, 0, 1])
todict()[source]
translations

Symmetry translations. The invertions are not included for centrosymmetrical crystals.

unique_reflections(hkl)[source]

Returns a subset hkl containing only the symmetry-unique reflections.

Example:

>>> from muesr.core.spg import Spacegroup
>>> sg = Spacegroup(225)  # fcc
>>> sg.unique_reflections([[ 2,  0,  0],
...                        [ 0, -2,  0],
...                        [ 2,  2,  0],
...                        [ 0, -2, -2]])
array([[2, 0, 0],
       [2, 2, 0]])
unique_sites(scaled_positions, symprec=0.001, output_mask=False, map_to_unitcell=True)[source]

Returns a subset of scaled_positions containing only the symmetry-unique positions. If output_mask is True, a boolean array masking the subset is also returned.

If map_to_unitcell is true, all sites are first mapped into the unit cell making e.g. [0, 0, 0] and [1, 0, 0] equivalent.

Example:

>>> from muesr.core.spg import Spacegroup
>>> sg = Spacegroup(225)  # fcc
>>> sg.unique_sites([[0.0, 0.0, 0.0],
...                  [0.5, 0.5, 0.0],
...                  [1.0, 0.0, 0.0],
...                  [0.5, 0.0, 0.0]])
array([[ 0. ,  0. ,  0. ],
       [ 0.5,  0. ,  0. ]])

muesr.i_o – Input Output functions

muesr.i_o.sampleIO – Input Output functions for the Sample object

muesr.i_o.sampleIO.load_sample(filename='', fileobj=None)[source]

This function load a sample from a file in YAML format.

Parameters:
  • filename (str) – the filename used to store data.
  • fileobj (file) – an optional file object. If specified, this supersede the filename input.
Returns:

a sample object

Return type:

Sample object or None

Raises:

ValueError, FileNotFoundError, IsADirectoryError

muesr.i_o.sampleIO.save_sample(sample, filename='', fileobj=None, overwrite=False)[source]

This function saves the sample provided in YAML format.

Parameters:
  • sample – the sample object
  • filename (str) – the filename used to store data.
  • fileobj (file) – a file object used in place of filename.
  • bool (overwrite) – if selected file should be overwritten.
Returns:

None

Return type:

None

Raises:

TypeError, ValueError, IsADirectoryError

muesr.i_o.exportFPS – Input Output for FullProf Studio

muesr.i_o.exportFPS.export_fpstudio(sample, filename)[source]

Exports magnetic structure in FullProfStudio format.

Parameters:
  • sample – the sample object
  • filename (str) – the filename for the FP Studio file
Returns:

None.

Return type:

None

muesr.i_o.xsf.xsf – Input Output for Xcrysden

muesr.i_o.xsf.xsf.load_xsf(sample, filename)[source]

Loads structural data from Xcrysden Files in sample object.

WARNING: this will reset all current muon positions, magnetic definitions, and the symmetry definition.

Parameters:
  • sample – a sample object.
  • filename (str) – the filename
Returns:

True if succeful, False otherwise

Return type:

bool

Raises:

TypeError

muesr.i_o.xsf.xsf.save_xsf(sample, filename, supercell=[1, 1, 1], addMuon=True)[source]

Export structure to XCrysDen.

Parameters:
  • sample – a sample object.
  • filename (str) – path of the destination file.
  • supercell (list) – a list containig the number of replica along the three lattice parameters
  • addMuon (bool) – if true, adds the muon positions (if any) in the central unit cell
Returns:

True if succeful, False otherwise

Return type:

bool

Raises:

CellError, TypeError

muesr.i_o.cif.cif – Input in CIF and mCIF format

This file was adapted from ASE. Module to read and write atoms in cif file format.

See http://www.iucr.org/resources/cif/spec/version1.1/cifsyntax for a description of the file format. STAR extensions as save frames, global blocks, nested loops and multi-data values are not supported.

muesr.i_o.cif.cif.convert_to_float(frac_str)[source]

This function converts fractions to float, ex. -1/3 = -0.33333…

muesr.i_o.cif.cif.convert_value(value)[source]

Convert CIF value string to corresponding python type.

muesr.i_o.cif.cif.load_cif(sample, filename, reset_muon=True, reset_sym=True)[source]

Loads the structural information from a Cif file.

Parameters:
  • sample (muesr.core.sample.Sample) – the sample object.
  • filename (str) – the cif file path (path + filename).
  • reset_muon (bool) – if true the muon positions is reinitialized.
  • reset_sym (bool) – if true the symmetry is reinitialized. ALWAYS TRUE
Returns:

True if succesfull, false otherwise

Return type:

bool

muesr.i_o.cif.cif.load_mcif(sample, filename, reset_muon=True, reset_sym=True)[source]

Loads both the crystalline structure and the magnetic order from a mcif file. N.B.: This function is EXPERIMENTAL.

Note

Only the first lattice and magnetic structure is loaded. mCif files hosting multiple structures are not supported.

Parameters:
  • sample (muesr.core.sample.Sample) – the sample object.
  • filename (str) – the mcif file path (path + filename).
  • reset_muon (bool) – if true the muon positions is reinitialized.
  • reset_sym (bool) – if true the symmetry is reinitialized.
Returns:

True if succesfull, false otherwise

Return type:

bool

muesr.i_o.cif.cif.parse_block(lines, line)[source]

Parse a CIF data block and return a tuple with the block name and a dict with all tags.

muesr.i_o.cif.cif.parse_cif(fileobj)[source]

Parse a CIF file. Returns a list of blockname and tag pairs. All tag names are converted to lower case.

muesr.i_o.cif.cif.parse_items(lines, line)[source]

Parse a CIF data items and return a dict with all tags.

muesr.i_o.cif.cif.parse_loop(lines)[source]

Parse a CIF loop. Returns a dict with column tag names as keys and a lists of the column content as values.

muesr.i_o.cif.cif.parse_magn_operation_xyz_string(tag)[source]

Parse the symmetry operations of the magnetic part of mcif

muesr.i_o.cif.cif.parse_multiline_string(lines, line)[source]

Parse semicolon-enclosed multiline string and return it.

muesr.i_o.cif.cif.parse_singletag(lines, line)[source]

Parse a CIF tag (entries starting with underscore). Returns a key-value pair.

muesr.i_o.cif.cif.read_cif(fileobj, index, store_tags=False, primitive_cell=False, subtrans_included=True)[source]

Read Atoms object from CIF file. index specifies the data block number or name (if string) to return.

If index is None or a slice object, a list of atoms objects will be returned. In the case of index is None or slice(None), only blocks with valid crystal data will be included.

If store_tags is true, the info attribute of the returned Atoms object will be populated with all tags in the corresponding cif data block.

If primitive_cell is true, the primitive cell will be built instead of the conventional cell.

If subtrans_included is true, sublattice translations are assumed to be included among the symmetry operations listed in the CIF file (seems to be the common behaviour of CIF files). Otherwise the sublattice translations are determined from setting 1 of the extracted space group. A result of setting this flag to true, is that it will not be possible to determine the primitive cell.

muesr.i_o.cif.cif.split_chem_form(comp_name)[source]

Returns e.g. AB2 as [‘A’, ‘1’, ‘B’, ‘2’]

muesr.i_o.cif.cif.tags2atoms(tags, store_tags=False, primitive_cell=False, subtrans_included=True)[source]

Returns an Atoms object from a cif tags dictionary. See read_cif() for a description of the arguments.

muesr.i_o.cif.cif.write_cif(fileobj, images, format='default')[source]

Write images to CIF file.

muesr.engines – Muesr calculators

muesr.engines.clfc – The C powered Local Fields Calculator

muesr.utilities – Various useful functions

muesr.utilities.dft_grid – Generate grid for DFT simulations

muesr.utilities.dft_grid.build_uniform_grid(sample, size, min_distance_from_atoms=1.0)[source]

Generates a grid of symmetry inequivalent interstitial positions with a specified minimum distance from the atoms of the sample. Especially intended for DFT simulations.

Parameters:
  • sample – A sample object.
  • size (int) – The number of steps in the three lattice directions. Only equispaced grids are supported at the moment.
  • min_distance_from_atoms (float) – Minimum distance between a interstitial position and the atoms of the lattice. Units are Angstrom.
Returns:

A list of symmetry inequivalent positions.

Return type:

list

muesr.utilities.symsearch – Identify symmetry (based on spglib)

muesr.utilities.symsearch.symsearch(sample, precision=0.0001)[source]

Identifies symmetry operations of the unit cell using spglib and update the sample definition.

Parameters:
  • sample – A sample object.
  • precision (float) – atoms are assumed equivalent if distances are smaller than precision. In Angstrom.
Returns:

True if successful, False otherwise.

Return type:

bool

muesr.utilities.ms – Functions facilitating magnetic description

muesr.utilities.ms.mago_add(sample, coordinates='b-c', fcs=None, kvalue=None)[source]

Adds a magnetic model (fourier components and K vector). The order is automatically selected if succesfully added.

Parameters:
  • sample – A sample object.
  • coordinates (string) –

    coordinates system and units of the Fourier components. Options are: ‘b-c’, ‘b/a-l’, ‘b-l’.

    b-c : Fourier components in Bohr magnetons and Cartesian coordinates.

    b/a-l: Fourier components in Bohr magnetons/Angstrom and lattice coordinates.

    b-l : Fourier components in Bohr magnetons and in lattice coordinates.

  • fcs (np.complex) – Fourier components in coordinate system (default: Bohr magnetoc/ Cartesian coordinates)
  • kvalue (np.ndarray) – Propagation vector in r.l.u.
Returns:

True if successful, False otherwise.

Return type:

bool

muesr.utilities.ms.mago_set_FC(sample, fcs=None, atoms_types=None, mm=None, inputConvention='b-c')[source]

Defines fourier components for the unit cell.

muesr.utilities.ms.mago_set_k(sample, kvalue=None, mm=None)[source]

Set propagation with respect to the conventional reciprocal cell kavalue : propagation vector (either string or tuple or list), if None, input is prompetd mm : Magnetic Model to be used. If None the current magnetic model is used returns: None

muesr.utilities.muon – Functions facilitating muon position definition

muesr.utilities.muon.muon_find_equiv(sample, eps=0.001)[source]

Given the unit cell symmetry, finds the equivalent muon sites. Magnetic order is NOT considered :params: eps, number of decimal figures for comparison

muesr.utilities.muon.muon_reset(sample)[source]

Removes all previously set muon positions. returns: True

muesr.utilities.muon.muon_set_frac(sample, arg=None)[source]

Set muon position in lattice coordinates.

returns: None

muesr.utilities.visualize – Interfaces to XCrysDen and VESTA

muesr.utilities.visualize.run_vesta(fname, block=True)[source]
muesr.utilities.visualize.run_xcrysden(fname, block=True)[source]
muesr.utilities.visualize.show_structure(sample, supercell=[1, 1, 1], askConfirm=True, block=True, visualizationTool='')[source]

Creates a supercell and shows it in XCrysDen or VESTA. This function only works interactively. WARNING: this can potentially harm your system!

Parameters:
  • sample – a sample object.
  • supercell (list) – a list containig the number of replica along the three lattice parameters
  • askConfirm (bool) – for CLI, ask for confirmation before running XCrysDen
  • block (bool) – the visualization program is forked if block=False, only works on unix
  • visualizationTool (str) – the preferred application for visualization. Valid values are ‘X’, ‘XCrySDen’, ‘V’ or ‘VESTA’. The selection is case insensitive. If missing, any available application will be used.
Returns:

True if succeful, False otherwise

Return type:

bool

Raises:

TypeError