Creating Atoms and Accessing Properties
This tutorial explains how to create atoms and access their relevant properties. This module allows for various isotopes to be specified as well as the natural abundances. To conclude this tutorial there is a short example of creating a water molecule with deuterium instead of hydrogen.
[1]:
# imports used for this tutorial
from MDMC.MD.structures import Atom
from MDMC.trajectory_analysis.observables.fqt import calc_incoherent_scatt_length
from MDMC.MD import *
Now that the Atom class has been imported, we can create an atom. This can be done to create an atom which represents the natural abundance of the element, or a specific isotope.
[2]:
# Creating an oxygen atom using MDMC Atom class
# Standard chemical symbol will give natural abundance
oxygen_na = Atom('O')
# Chemical symbol followed by the mass number will give the specific isotope
oxygen_16 = Atom('O[16]')
# Example for specifying different isotope
oxygen_17 = Atom('O[17]')
Now that the atoms have been created, we can start accessing their associated properties. MDMC uses the Python package periodictable for a database of atomic properties.
[3]:
# Accessing general properties
print("atomic mass: ", oxygen_na.mass)
print("atomic number: ", oxygen_na.element.number)
# Showing the atomic mass for the different oxygen atoms
print("O[16] mass: ", oxygen_16.mass)
print("O[17] mass: ", oxygen_17.mass)
atomic mass: 15.9994 amu
atomic number: 8
O[16] mass: 15.9949146221 amu
O[17] mass: 16.9991315 amu
As well as basic properties like mass and atomic number, we can access scattering-specific properties like scattering lengths and cross-sections.
[4]:
# Accessing scattering specific properties, like scattering lengths and cross sections
print("O[16] coherent cross section: ", oxygen_16.element.neutron.coherent)
print("O[16] incoherent cross section ", oxygen_16.element.neutron.incoherent)
print("O[16] coherent scattering length: ", oxygen_16.element.neutron.b_c)
print("O[16] incoherent scattering length: ", oxygen_16.element.neutron.b_c_i)
print(" ") # for better formatting
# Many incoherent scattering lengths are stored as 0, but with non-zero cross sections.
# The lengths can be approximated in MDMC using 'calc_incoherent_scatt_length'.status
nitrogen_na = Atom('N')
print("Nitrogen incoherent cross section ", nitrogen_na.element.neutron.incoherent)
print("Nitrogen (periodictable) incoherent scattering length: ", nitrogen_na.element.neutron.b_c_i)
print("Nitrogen (MDMC) incoherent scattering length: ", calc_incoherent_scatt_length('N'))
O[16] coherent cross section: 4.232
O[16] incoherent cross section 0.0
O[16] coherent scattering length: 5.805
O[16] incoherent scattering length: None
Nitrogen incoherent cross section 0.5
Nitrogen (periodictable) incoherent scattering length: None
Nitrogen (MDMC) incoherent scattering length: 1.9947114020071635
For further information about key attributes of a periodictable object, we can use help().
[5]:
help(oxygen_na.element.neutron)
Help on Neutron in module periodictable.nsf object:
class Neutron(builtins.object)
| Neutron scattering factors are attached to each element in the periodic
| table for which values are available. If no information is available,
| then the neutron field of the element will be *None*. Even when neutron
| information is available, it may not be complete, so individual fields
| may be *None*.
|
| The following fields are defined:
|
| * b_c (fm)
| Bounds coherent scattering length.
|
| * total (barn)
| Total scattering cross section $\sigma_s$. This does not include the
| absorption cross section. To compute the total collision cross
| section use $\sigma_t = \sigma_s + \sigma_a$
|
| * absorption (barn)
| Absorption cross section $\sigma_a$ at 1.798 |Ang|. Scale to your beam
| by dividing by periodictable.nsf.ABSORPTION_WAVELENGTH and multiplying
| by your wavelength.
|
| * b_c_complex (fm)
| Complex coherent scattering length derived from the tabulated
| values using $b_c - i \sigma_a / (1000 \cdot 2 \lambda)$.
|
| Additional columns not used for calculation include:
|
| * b_c_i (fm)
| Imaginary bound coherent scattering length. This is related to
| absorption cross section by $\sigma_a = 4 \pi \mathrm{Im}(b_c)/k$ where
| $k = 2 \pi/\lambda$ and an additional factor of 1000 for converting
| between |Ang|\ |cdot|\ fm and barns. b_c_i is not available for
| all isotopes for which absorption cross sections have been measured.
|
| * bp, bm (fm)
| Spin-dependent scattering for I+1/2 and I-1/2 (not always available).
| Incoherent scattering arises from the spin-dependent scattering b+
| and b-. The Neutron Data Booklet\ [#Rauch2003]_ gives formulas for
| calculating coherent and incoherent scattering from b+ and b- alone.
|
| * bp_i, bm_i (fm)
| Imaginary portion of bp and bm.
|
| * is_energy_dependent (boolean)
| Do not use this data if scattering is energy dependent.
|
| * coherent (barn)
| Coherent scattering cross section. This is tabulated but not used.
| In theory coherent scattering is related to bound coherent scattering
| by $\sigma_c = 4 \pi |\mathrm{Re}(b_c) + i \mathrm{Im}(b_c)|^2/100$.
| In practice, these values are different, with the following table
| showing the largest relative difference:
|
| ======== ======== ======== ======== ========
| Sc 3% Ti 4% V 34% Mn 1% Cd 2%
| Te 4% Xe 9% Sm 19% Eu 44% Tb 1%
| Ho 11% W 4% Au 7% Hg 2% Ra 3%
| ======== ======== ======== ======== ========
|
| * incoherent (barn)
| Incoherent scattering cross section $\sigma_i$. This is tabulated but
| not used. Instead, the incoherent cross section is computed from the
| total cross section minus the coherent cross section even for single
| atoms so that results from compounds are consistent with results from
| single atoms.
|
| For elements, the scattering cross-sections are based on the natural
| abundance of the individual isotopes. Individual isotopes may have
| the following additional fields
|
| * abundance (%)
| Isotope abundance used to compute the properties of the element in
| natural abundance.
|
| * nuclear_spin (string)
| Spin on the nucleus: '0', '1/2', '3/2', etc.
|
| Each field ``T`` above has a corresponding ``T_units`` attribute with
| the name of the units.
|
| For scattering calculations the scattering length density is the value
| of interest. This is computed from the *number_density* of the individual
| elements, as derived from the element density and atomic mass.
|
| .. Note:: 1 barn = 100 |fm^2|
|
| Methods defined here:
|
| __init__(self)
| Initialize self. See help(type(self)) for accurate signature.
|
| __str__(self)
| Return str(self).
|
| has_sld(self)
| Returns *True* if sld is defined for this element/isotope.
|
| scattering(self, wavelength=1.798)
| Returns neutron scattering information for the element at natural
| abundance and density.
|
| :Parameters:
| *wavelength* \: float(s) | |Ang|
|
| :Returns:
| *sld* \: float(s), float(s), float(s) | |1e-6/Ang^2|
| (*real*, -*imaginary*, *incoherent*) scattering length density
|
| *xs* \: float(s), float(s), float(s) | |1/cm|
| (*coherent*, *absorption*, *incoherent*) cross sections.
| :w
|
| *penetration* \: float(s) | cm
| 1/e penetration length.
|
| Returns (None, None, None) if sld is not known for this element.
|
| See :func:`neutron_scattering` for details.
|
| scattering_by_wavelength(self, wavelength)
| Return scattering length and total cross section for each wavelength.
|
| For rare earth isotopes this returns the energy-dependent
| $\mathrm{Re}(b_c)$ and $\mathrm{Im}(b_c)$ interpolated into the
| scattering length tables. Values are extrapolated with constant
| values at the ends of the table. Total scattering is returned as
| $4\pi/100 |b_c|^2$ with no contribution for bound incoherent
| scattering.
|
| :Parameters:
| *wavelength* \: float(s) | |Ang|
|
| :Returns:
| *b_c* \: complex(s) | fm
|
| *sigma_s* \: float(s) | barn
|
| sld(self, wavelength=1.798)
| Returns scattering length density for the element at natural
| abundance and density.
|
| :Parameters:
| *wavelength* \: float(s) | |Ang|
|
| :Returns:
| *sld* \: float(s), float(s), float(s) | |1e-6/Ang^2|
| (*real*, -*imaginary*, *incoherent*) scattering length density.
|
| Returns (None, None, None) if sld is not known for this element.
|
| See :func:`neutron_scattering` for details.
|
| ----------------------------------------------------------------------
| Data descriptors defined here:
|
| __dict__
| dictionary for instance variables
|
| __weakref__
| list of weak references to the object
|
| ----------------------------------------------------------------------
| Data and other attributes defined here:
|
| absorption = None
|
| absorption_units = 'barn'
|
| abundance = 0.0
|
| abundance_units = '%'
|
| b_c = None
|
| b_c_complex = None
|
| b_c_complex_units = 'fm'
|
| b_c_i = None
|
| b_c_i_units = 'fm'
|
| b_c_units = 'fm'
|
| bm = None
|
| bm_i = None
|
| bm_units = 'fm'
|
| bp = None
|
| bp_i = None
|
| bp_units = 'fm'
|
| coherent = None
|
| coherent_units = 'barn'
|
| incoherent = None
|
| incoherent_units = 'barn'
|
| is_energy_dependent = False
|
| nsf_table = None
|
| total = None
|
| total_units = 'barn'
Now as a final example of this functionality, let’s create a molecule. We will create a heavy water molecule, using deuterium in the place of hydrogen.
[6]:
# Conveniently, deuterium can be specified directly through this artificial chemical symbol 'D'
# It can also be specified using 'H[2]' as expected.
D1 = Atom('D') # or Atom("H[2]")
D2 = Atom('D', position=(0., 1.63298, 0.))
O = Atom('O', position=(0., 0.81649, 0.57736))
D_coulombic = Coulombic(atoms=[D1, D2], cutoff=10.)
O_coulombic = Coulombic(atoms=O, cutoff=10.)
heavy_water_mol = Molecule(position=(0, 0, 0),
velocity=(0, 0, 0),
atoms=[D1, D2, O],
interactions=[Bond((D1, O), (D2, O), constrained=True),
BondAngle(D1, O, D2, constrained=True)],
name='heavy_water')
For a more in depth explanation of molecule and bond creation please see the MDMC tutorial on creating atomic configurations..