simulations.xrays

Contents

simulations.xrays#

Classes:

Xray(S, force_recalc, **kwargs)

Base class for all X-ray scattering simulations.

XrayKin(S, force_recalc, **kwargs)

Kinetic X-ray scattering simulations.

XrayDyn(S, force_recalc, **kwargs)

Dynamical X-ray scattering simulations.

XrayDynMag(S, force_recalc, **kwargs)

Dynamical magnetic X-ray scattering simulations.

class udkm1Dsim.simulations.xrays.Xray(S, force_recalc, **kwargs)[source]#

Bases: Simulation

Base class for all X-ray scattering simulations.

Parameters:
  • S (Structure) – sample to do simulations with.

  • force_recalc (boolean) – force recalculation of results.

Keyword Arguments:
  • save_data (boolean) – true to save simulation results.

  • cache_dir (str) – path to cached data.

  • disp_messages (boolean) – true to display messages from within the simulations.

  • progress_bar (boolean) – enable tqdm progress bar.

Attributes:
  • S (Structure) – sample structure to calculate simulations on.

  • force_recalc (boolean) – force recalculation of results.

  • save_data (boolean) – true to save simulation results.

  • cache_dir (str) – path to cached data.

  • disp_messages (boolean) – true to display messages from within the simulations.

  • progress_bar (boolean) – enable tqdm progress bar.

  • energy (ndarray[float]) – photon energies \(E\) of scattering light

  • wl (ndarray[float]) – wavelengths \(\lambda\) of scattering light

  • k (ndarray[float]) – wavenumber \(k\) of scattering light

  • theta (ndarray[float]) – incidence angles \(\theta\) of scattering light

  • qz (ndarray[float]) – scattering vector \(q_z\) of scattering light

  • polarizations (dict) – polarization states and according names.

  • pol_in_state (int) – incoming polarization state as defined in polarizations dict.

  • pol_out_state (int) – outgoing polarization state as defined in polarizations dict.

  • pol_in (float) – incoming polarization factor (can be a complex ndarray).

  • pol_out (float) – outgoing polarization factor (can be a complex ndarray).

Methods:

set_incoming_polarization(pol_in_state)

Must be overwritten by child classes.

set_outgoing_polarization(pol_out_state)

Must be overwritten by child classes.

set_polarization(pol_in_state, pol_out_state)

Sets the incoming and analyzer (outgoing) polarization.

get_hash(strain_vectors, **kwargs)

Calculates an unique hash given by the energy \(E\), \(q_z\) range, polarization states and the strain_vectors as well as the sample structure hash for relevant x-ray parameters.

get_polarization_factor(theta)

Calculates the polarization factor \(P(\vartheta)\) for a given incident angle \(\vartheta\) for the case of s-polarization (pol = 0), or p-polarization (pol = 1), or unpolarized X-rays (pol = 0.5):

update_experiment(caller)

Recalculate energy, wavelength, and wavevector as well as theta and the scattering vector in case any of these has changed.

conv_with_function(y, x, handle)

Convolutes the array \(y(x)\) with a function given by the handle on the argument array \(x\).

disp_message(message)

Wrapper to display messages for that class.

save(full_filename, data, *args)

Save data to file.

set_incoming_polarization(pol_in_state)[source]#

Must be overwritten by child classes.

Parameters:

pol_in_state (int) – incoming polarization state id.

set_outgoing_polarization(pol_out_state)[source]#

Must be overwritten by child classes.

Parameters:

pol_out_state (int) – outgoing polarization state id.

set_polarization(pol_in_state, pol_out_state)[source]#

Sets the incoming and analyzer (outgoing) polarization.

Parameters:
  • pol_in_state (int) – incoming polarization state id.

  • pol_out_state (int) – outgoing polarization state id.

get_hash(strain_vectors, **kwargs)[source]#

Calculates an unique hash given by the energy \(E\), \(q_z\) range, polarization states and the strain_vectors as well as the sample structure hash for relevant x-ray parameters. Optionally, part of the strain_map is used.

Parameters:
  • strain_vectors (dict{ndarray[float]}) – reduced strains per unique layer.

  • **kwargs (ndarray[float]) – spatio-temporal strain profile.

Returns:

hash (str) – unique hash.

get_polarization_factor(theta)[source]#

Calculates the polarization factor \(P(\vartheta)\) for a given incident angle \(\vartheta\) for the case of s-polarization (pol = 0), or p-polarization (pol = 1), or unpolarized X-rays (pol = 0.5):

\[P(\vartheta) = \sqrt{(1-\mbox{pol}) + \mbox{pol} \cdot \cos(2\vartheta)}\]
Parameters:

theta (ndarray[float]) – incidence angle.

Returns:

P (ndarray[float]) – polarization factor.

update_experiment(caller)[source]#

Recalculate energy, wavelength, and wavevector as well as theta and the scattering vector in case any of these has changed.

\[\begin{split}\lambda & = \frac{hc}{E} \\ E & = \frac{hc}{\lambda} \\ k & = \frac{2\pi}{\lambda} \\ \vartheta & = \arcsin{\frac{\lambda q_z}{4\pi}} \\ q_z & = 2k \sin{\vartheta}\end{split}\]
Parameters:

caller (str) – name of calling method.

static conv_with_function(y, x, handle)#

Convolutes the array \(y(x)\) with a function given by the handle on the argument array \(x\).

Parameters:
  • y (ndarray[float]) – y data.

  • x (ndarray[float]) – x data.

  • handle (@lamdba) – convolution function.

Returns:

y_conv (ndarray[float]) – convoluted data.

disp_message(message)#

Wrapper to display messages for that class.

Parameters:

message (str) – message to display.

save(full_filename, data, *args)#

Save data to file. The variable name can be handed as variable argument.

Parameters:
  • full_filename (str) – full file name to data file.

  • data (ndarray) – actual data to save.

  • *args (str, optional) – variable name within the data file.

class udkm1Dsim.simulations.xrays.XrayKin(S, force_recalc, **kwargs)[source]#

Bases: Xray

Kinetic X-ray scattering simulations.

Parameters:
  • S (Structure) – sample to do simulations with.

  • force_recalc (boolean) – force recalculation of results.

Keyword Arguments:
  • save_data (boolean) – true to save simulation results.

  • cache_dir (str) – path to cached data.

  • disp_messages (boolean) – true to display messages from within the simulations.

  • progress_bar (boolean) – enable tqdm progress bar.

Attributes:
  • S (Structure) – sample structure to calculate simulations on.

  • force_recalc (boolean) – force recalculation of results.

  • save_data (boolean) – true to save simulation results.

  • cache_dir (str) – path to cached data.

  • disp_messages (boolean) – true to display messages from within the simulations.

  • progress_bar (boolean) – enable tqdm progress bar.

  • energy (ndarray[float]) – photon energies \(E\) of scattering light

  • wl (ndarray[float]) – wavelengths \(\lambda\) of scattering light

  • k (ndarray[float]) – wavenumber \(k\) of scattering light

  • theta (ndarray[float]) – incidence angles \(\theta\) of scattering light

  • qz (ndarray[float]) – scattering vector \(q_z\) of scattering light

  • polarizations (dict) – polarization states and according names.

  • pol_in_state (int) – incoming polarization state as defined in polarizations dict.

  • pol_out_state (int) – outgoing polarization state as defined in polarizations dict.

  • pol_in (float) – incoming polarization factor (can be a complex ndarray).

  • pol_out (float) – outgoing polarization factor (can be a complex ndarray).

References

Methods:

set_incoming_polarization(pol_in_state)

Sets the incoming polarization factor for sigma, pi, and unpolarized polarization.

set_outgoing_polarization(pol_out_state)

For kinematical X-ray simulation only "no analyzer polarization" is allowed.

get_uc_atomic_form_factors(energy, qz, uc)

Returns the energy- and angle-dependent atomic form factors :math: f(q_z, E) of all atoms in the unit cell as a vector.

get_uc_structure_factor(energy, qz, uc[, strain])

Calculates the energy-, angle-, and strain-dependent structure factor .

homogeneous_reflectivity([strains])

Calculates the reflectivity \(R = E_p^t\,(E_p^t)^*\) of a homogeneous sample structure as well as the reflected field \(E_p^N\) of all substructures.

homogeneous_reflected_field(S, energy, qz, theta)

Calculates the reflected field \(E_p^t\) of the whole sample structure as well as for each sub-structure (\(E_p^N\)).

get_interference_function(qz, z, N)

Calculates the interference function for \(N\) repetitions of the structure with the length \(z\):

get_Ep(energy, qz, theta, uc, strain)

Calculates the reflected field \(E_p\) for one unit cell with a given strain \(\epsilon\):

conv_with_function(y, x, handle)

Convolutes the array \(y(x)\) with a function given by the handle on the argument array \(x\).

disp_message(message)

Wrapper to display messages for that class.

get_hash(strain_vectors, **kwargs)

Calculates an unique hash given by the energy \(E\), \(q_z\) range, polarization states and the strain_vectors as well as the sample structure hash for relevant x-ray parameters.

get_polarization_factor(theta)

Calculates the polarization factor \(P(\vartheta)\) for a given incident angle \(\vartheta\) for the case of s-polarization (pol = 0), or p-polarization (pol = 1), or unpolarized X-rays (pol = 0.5):

save(full_filename, data, *args)

Save data to file.

set_polarization(pol_in_state, pol_out_state)

Sets the incoming and analyzer (outgoing) polarization.

update_experiment(caller)

Recalculate energy, wavelength, and wavevector as well as theta and the scattering vector in case any of these has changed.

set_incoming_polarization(pol_in_state)[source]#

Sets the incoming polarization factor for sigma, pi, and unpolarized polarization.

Parameters:

pol_in_state (int) – incoming polarization state id.

set_outgoing_polarization(pol_out_state)[source]#

For kinematical X-ray simulation only “no analyzer polarization” is allowed.

Parameters:

pol_out_state (int) – outgoing polarization state id.

get_uc_atomic_form_factors(energy, qz, uc)[source]#

Returns the energy- and angle-dependent atomic form factors :math: f(q_z, E) of all atoms in the unit cell as a vector.

Parameters:
  • energy (float, Quantity) – photon energy.

  • qz (ndarray[float, Quantity]) – scattering vectors.

  • uc (UnitCell) – unit cell object.

Returns:

f (ndarray[complex]) – unit cell atomic form factors.

get_uc_structure_factor(energy, qz, uc, strain=0)[source]#

Calculates the energy-, angle-, and strain-dependent structure factor .. math: S(E,q_z,epsilon) of the unit cell:

\[S(E,q_z,\epsilon) = \sum_i^N f_i \, \exp(-i q_z z_i(\epsilon))\]
Parameters:
  • energy (float, Quantity) – photon energy.

  • qz (ndarray[float, Quantity]) – scattering vectors.

  • uc (UnitCell) – unit cell object.

  • strain (float, optional) – strain of the unit cell 0 .. 1. Defaults to 0.

Returns:

S (ndarray[complex]) – unit cell structure factor.

homogeneous_reflectivity(strains=0)[source]#

Calculates the reflectivity \(R = E_p^t\,(E_p^t)^*\) of a homogeneous sample structure as well as the reflected field \(E_p^N\) of all substructures.

Parameters:

strains (ndarray[float], optional) – strains of each sub-structure 0 .. 1. Defaults to 0.

Returns:

(tuple)

  • R (ndarray[complex]) - homogeneous reflectivity.

  • A (ndarray[complex]) - reflected fields of sub-structures.

homogeneous_reflected_field(S, energy, qz, theta, strains=0)[source]#

Calculates the reflected field \(E_p^t\) of the whole sample structure as well as for each sub-structure (\(E_p^N\)). The reflected wave field \(E_p\) from a single layer of unit cells at the detector is calculated according to Ref. [9]:

\[E_p = \frac{i}{\varepsilon_0}\frac{e^2}{m_e c_0^2} \frac{P(\vartheta) S(E,q_z,\epsilon)}{A q_z}\]

For the case of \(N\) similar planes of unit cells one can write:

\[E_p^N = \sum_{n=0}^{N-1} E_p \exp(i q_z z n )\]

where \(z\) is the distance between the planes (c-axis). The above equation can be simplified to:

\[E_p^N = E_p \psi(q_z,z,N)\]

introducing the interference function

\[\begin{split}\psi(q_z,z,N) & = \sum_{n=0}^{N-1} \exp(i q_z z n) \\ & = \frac{1- \exp(i q_z z N)}{1- \exp(i q_z z)}\end{split}\]

The total reflected wave field of all \(i = 1\ldots M\) homogeneous layers (\(E_p^t\)) is the phase-correct summation of all individual \(E_p^{N,i}\):

\[E_p^t = \sum_{i=1}^M E_p^{N,i} \exp(i q_z Z_i)\]

where \(Z_i = \sum_{j=1}^{i-1} N_j z_j\) is the distance of the \(i\)-th layer from the surface.

Parameters:
  • S (Structure, UnitCell) – structure or sub-structure to calculate on.

  • energy (float, Quantity) – photon energy.

  • qz (ndarray[float, Quantity]) – scattering vectors.

  • theta (ndarray[float, Quantity]) – scattering incidence angle.

  • strains (ndarray[float], optional) – strains of each sub-structure 0 .. 1. Defaults to 0.

Returns:

(tuple)

  • Ept (ndarray[complex]) - reflected field.

  • A (ndarray[complex]) - reflected fields of substructures.

get_interference_function(qz, z, N)[source]#

Calculates the interference function for \(N\) repetitions of the structure with the length \(z\):

\[\begin{split}\psi(q_z,z,N) & = \sum_{n=0}^{N-1} \exp(i q_z z n) \\ & = \frac{1- \exp(i q_z z N)}{1- \exp(i q_z z)}\end{split}\]
Parameters:
  • qz (ndarray[float, Quantity]) – scattering vectors.

  • z (float) – thickness/length of the structure.

  • N (int) – repetitions of the structure.

Returns:

psi (ndarray[complex]) – interference function.

get_Ep(energy, qz, theta, uc, strain)[source]#

Calculates the reflected field \(E_p\) for one unit cell with a given strain \(\epsilon\):

\[E_p = \frac{i}{\varepsilon_0} \frac{e^2}{m_e c_0^2} \frac{P S(E,q_z,\epsilon)}{A q_z}\]

with \(e\) as electron charge, \(m_e\) as electron mass, \(c_0\) as vacuum light velocity, \(\varepsilon_0\) as vacuum permittivity, \(P\) as polarization factor and \(S(E,q_z,\sigma)\) as energy-, angle-, and strain-dependent unit cell structure factor.

Parameters:
  • energy (float, Quantity) – photon energy.

  • qz (ndarray[float, Quantity]) – scattering vectors.

  • theta (ndarray[float, Quantity]) – scattering incidence angle.

  • uc (UnitCell) – unit cell object.

  • strain (float, optional) – strain of the unit cell 0 .. 1. Defaults to 0.

Returns:

Ep (ndarray[complex]) – reflected field.

static conv_with_function(y, x, handle)#

Convolutes the array \(y(x)\) with a function given by the handle on the argument array \(x\).

Parameters:
  • y (ndarray[float]) – y data.

  • x (ndarray[float]) – x data.

  • handle (@lamdba) – convolution function.

Returns:

y_conv (ndarray[float]) – convoluted data.

disp_message(message)#

Wrapper to display messages for that class.

Parameters:

message (str) – message to display.

get_hash(strain_vectors, **kwargs)#

Calculates an unique hash given by the energy \(E\), \(q_z\) range, polarization states and the strain_vectors as well as the sample structure hash for relevant x-ray parameters. Optionally, part of the strain_map is used.

Parameters:
  • strain_vectors (dict{ndarray[float]}) – reduced strains per unique layer.

  • **kwargs (ndarray[float]) – spatio-temporal strain profile.

Returns:

hash (str) – unique hash.

get_polarization_factor(theta)#

Calculates the polarization factor \(P(\vartheta)\) for a given incident angle \(\vartheta\) for the case of s-polarization (pol = 0), or p-polarization (pol = 1), or unpolarized X-rays (pol = 0.5):

\[P(\vartheta) = \sqrt{(1-\mbox{pol}) + \mbox{pol} \cdot \cos(2\vartheta)}\]
Parameters:

theta (ndarray[float]) – incidence angle.

Returns:

P (ndarray[float]) – polarization factor.

save(full_filename, data, *args)#

Save data to file. The variable name can be handed as variable argument.

Parameters:
  • full_filename (str) – full file name to data file.

  • data (ndarray) – actual data to save.

  • *args (str, optional) – variable name within the data file.

set_polarization(pol_in_state, pol_out_state)#

Sets the incoming and analyzer (outgoing) polarization.

Parameters:
  • pol_in_state (int) – incoming polarization state id.

  • pol_out_state (int) – outgoing polarization state id.

update_experiment(caller)#

Recalculate energy, wavelength, and wavevector as well as theta and the scattering vector in case any of these has changed.

\[\begin{split}\lambda & = \frac{hc}{E} \\ E & = \frac{hc}{\lambda} \\ k & = \frac{2\pi}{\lambda} \\ \vartheta & = \arcsin{\frac{\lambda q_z}{4\pi}} \\ q_z & = 2k \sin{\vartheta}\end{split}\]
Parameters:

caller (str) – name of calling method.

class udkm1Dsim.simulations.xrays.XrayDyn(S, force_recalc, **kwargs)[source]#

Bases: Xray

Dynamical X-ray scattering simulations.

Parameters:
  • S (Structure) – sample to do simulations with.

  • force_recalc (boolean) – force recalculation of results.

Keyword Arguments:
  • save_data (boolean) – true to save simulation results.

  • cache_dir (str) – path to cached data.

  • disp_messages (boolean) – true to display messages from within the simulations.

  • progress_bar (boolean) – enable tqdm progress bar.

Attributes:
  • S (Structure) – sample structure to calculate simulations on.

  • force_recalc (boolean) – force recalculation of results.

  • save_data (boolean) – true to save simulation results.

  • cache_dir (str) – path to cached data.

  • disp_messages (boolean) – true to display messages from within the simulations.

  • progress_bar (boolean) – enable tqdm progress bar.

  • energy (ndarray[float]) – photon energies \(E\) of scattering light

  • wl (ndarray[float]) – wavelengths \(\lambda\) of scattering light

  • k (ndarray[float]) – wavenumber \(k\) of scattering light

  • theta (ndarray[float]) – incidence angles \(\theta\) of scattering light

  • qz (ndarray[float]) – scattering vector \(q_z\) of scattering light

  • polarizations (dict) – polarization states and according names.

  • pol_in_state (int) – incoming polarization state as defined in polarizations dict.

  • pol_out_state (int) – outgoing polarization state as defined in polarizations dict.

  • pol_in (float) – incoming polarization factor (can be a complex ndarray).

  • pol_out (float) – outgoing polarization factor (can be a complex ndarray).

  • last_atom_ref_trans_matrices (list) – remember last result of atom ref_trans_matrices to speed up calculation.

Methods:

set_incoming_polarization(pol_in_state)

Sets the incoming polarization factor for sigma, pi, and unpolarized polarization.

set_outgoing_polarization(pol_out_state)

For dynamical X-ray simulation only "no analyzer polarization" is allowed.

homogeneous_reflectivity(*args)

Calculates the reflectivity \(R\) of the whole sample structure and the reflectivity-transmission matrices \(M_{RT}\) for each substructure.

homogeneous_ref_trans_matrix(S, *args)

Calculates the reflectivity-transmission matrices \(M_{RT}\) of the whole sample structure as well as for each sub-structure.

inhomogeneous_reflectivity(strain_map, ...)

Returns the reflectivity of an inhomogeneously strained sample structure for a given strain_map in position and time, as well as for a given set of possible strains for each unit cell in the sample structure (strain_vectors).

sequential_inhomogeneous_reflectivity(...)

Returns the reflectivity of an inhomogeneously strained sample structure for a given strain_map in position and time, as well as for a given set of possible strains for each unit cell in the sample structure (strain_vectors).

parallel_inhomogeneous_reflectivity(...)

Returns the reflectivity of an inhomogeneously strained sample structure for a given strain_map in position and time, as well as for a given set of possible strains for each unit cell in the sample structure (strain_vectors).

distributed_inhomogeneous_reflectivity(...)

This is a stub.

calc_inhomogeneous_reflectivity(strains, ...)

Calculates the reflectivity of a inhomogeneous sample structure for given strain_vectors for a single time step.

calc_inhomogeneous_ref_trans_matrix(...)

Sub-function of calc_inhomogeneous_reflectivity() and for parallel computing (needs to be static) only for calculating the total reflection-transmission matrix \(M_{RT}^t\):

get_all_ref_trans_matrices(*args)

Returns a list of all reflection-transmission matrices for each unique unit cell in the sample structure for a given set of applied strains for each unique unit cell given by the strain_vectors input.

calc_all_ref_trans_matrices(*args)

Calculates a list of all reflection-transmission matrices for each unique unit cell in the sample structure for a given set of applied strains to each unique unit cell given by the strain_vectors input.

get_uc_ref_trans_matrix(uc, *args)

Returns the reflection-transmission matrix of a unit cell:

get_atom_ref_trans_matrix(atom, area, ...)

Calculates the reflection-transmission matrix of an atom from dynamical x-ray theory:

get_atom_reflection_factor(atom, area, ...)

Calculates the reflection factor from dynamical x-ray theory:

get_atom_transmission_factor(atom, area, ...)

Calculates the transmission factor from dynamical x-ray theory:

get_atom_phase_matrix(distance)

Calculates the phase matrix from dynamical x-ray theory:

get_atom_phase_factor(distance)

Calculates the phase factor \(\phi\) for a distance \(d\) from dynamical x-ray theory:

calc_reflectivity_from_matrix(M)

Calculates the reflectivity from an \(2\times2\) matrix of transmission and reflectivity factors:

conv_with_function(y, x, handle)

Convolutes the array \(y(x)\) with a function given by the handle on the argument array \(x\).

disp_message(message)

Wrapper to display messages for that class.

get_hash(strain_vectors, **kwargs)

Calculates an unique hash given by the energy \(E\), \(q_z\) range, polarization states and the strain_vectors as well as the sample structure hash for relevant x-ray parameters.

get_polarization_factor(theta)

Calculates the polarization factor \(P(\vartheta)\) for a given incident angle \(\vartheta\) for the case of s-polarization (pol = 0), or p-polarization (pol = 1), or unpolarized X-rays (pol = 0.5):

save(full_filename, data, *args)

Save data to file.

set_polarization(pol_in_state, pol_out_state)

Sets the incoming and analyzer (outgoing) polarization.

update_experiment(caller)

Recalculate energy, wavelength, and wavevector as well as theta and the scattering vector in case any of these has changed.

set_incoming_polarization(pol_in_state)[source]#

Sets the incoming polarization factor for sigma, pi, and unpolarized polarization.

Parameters:

pol_in_state (int) – incoming polarization state id.

set_outgoing_polarization(pol_out_state)[source]#

For dynamical X-ray simulation only “no analyzer polarization” is allowed.

Parameters:

pol_out_state (int) – outgoing polarization state id.

homogeneous_reflectivity(*args)[source]#

Calculates the reflectivity \(R\) of the whole sample structure and the reflectivity-transmission matrices \(M_{RT}\) for each substructure. The reflectivity of the \(2\times 2\) matrices for each \(q_z\) is calculates as follow:

\[R = \left|M_{RT}^t(0,1)/M_{RT}^t(1,1)\right|^2\]
Parameters:

*args (ndarray[float], optional) – strains for each substructure.

Returns:

(tuple)

  • R (ndarray[float]) - homogeneous reflectivity.

  • A (ndarray[complex]) - reflectivity-transmission matrices of sub-structures.

homogeneous_ref_trans_matrix(S, *args)[source]#

Calculates the reflectivity-transmission matrices \(M_{RT}\) of the whole sample structure as well as for each sub-structure. The reflectivity-transmission matrix of a single unit cell is calculated from the reflection-transmission matrices \(H_i\) of each atom and the phase matrices between the atoms \(L_i\):

\[M_{RT} = \prod_i H_i \ L_i\]

For \(N\) similar layers of unit cells one can calculate the \(N\)-th power of the unit cell \(\left(M_{RT}\right)^N\). The reflection-transmission matrix for the whole sample \(M_{RT}^t\) consisting of \(j = 1\ldots M\) sub-structures is then again:

\[M_{RT}^t = \prod_{j=1}^M \left(M_{RT^,j}\right)^{N_j}\]
Parameters:
  • S (Structure, UnitCell) – structure or sub-structure to calculate on.

  • *args (ndarray[float], optional) – strains for each substructure.

Returns:

(tuple)

  • RT (ndarray[complex]) - reflectivity-transmission matrix.

  • A (ndarray[complex]) - reflectivity-transmission matrices of sub-structures.

inhomogeneous_reflectivity(strain_map, strain_vectors, **kwargs)[source]#

Returns the reflectivity of an inhomogeneously strained sample structure for a given strain_map in position and time, as well as for a given set of possible strains for each unit cell in the sample structure (strain_vectors). If no reflectivity is saved in the cache it is caluclated. Providing the calc_type for the calculation the corresponding sub-routines for the reflectivity computation are called:

  • parallel parallelization over the time steps utilizing Dask

  • distributed not implemented in Python, but should be possible with Dask as well

  • sequential no parallelization at all

Parameters:
  • strain_map (ndarray[float]) – spatio-temporal strain profile.

  • strain_vectors (list[ndarray[float]]) – reduced strains per unique layer.

  • **kwargs

    • calc_type (str) - type of calculation.

    • dask_client (Dask.Client) - Dask client.

    • job (Dask.job) - Dask job.

    • num_workers (int) - Dask number of workers.

Returns:

R (ndarray[float]) – inhomogeneous reflectivity.

sequential_inhomogeneous_reflectivity(strain_map, strain_vectors, RTM)[source]#

Returns the reflectivity of an inhomogeneously strained sample structure for a given strain_map in position and time, as well as for a given set of possible strains for each unit cell in the sample structure (strain_vectors). The function calculates the results sequentially without parallelization.

Parameters:
  • strain_map (ndarray[float]) – spatio-temporal strain profile.

  • strain_vectors (list[ndarray[float]]) – reduced strains per unique layer.

  • RTM (list[ndarray[complex]]) – reflection-transmission matrices for all given strains per unique layer.

Returns:

R (ndarray[float]) – inhomogeneous reflectivity.

parallel_inhomogeneous_reflectivity(strain_map, strain_vectors, RTM, dask_client)[source]#

Returns the reflectivity of an inhomogeneously strained sample structure for a given strain_map in position and time, as well as for a given set of possible strains for each unit cell in the sample structure (strain_vectors). The function parallelizes the calculation over the time steps, since the results do not depend on each other.

Parameters:
  • strain_map (ndarray[float]) – spatio-temporal strain profile.

  • strain_vectors (list[ndarray[float]]) – reduced strains per unique layer.

  • RTM (list[ndarray[complex]]) – reflection-transmission matrices for all given strains per unique layer.

  • dask_client (Dask.Client) – Dask client.

Returns:

R (ndarray[float]) – inhomogeneous reflectivity.

distributed_inhomogeneous_reflectivity(strain_map, strain_vectors, RTM, job, num_worker)[source]#

This is a stub. Not yet implemented in python.

Parameters:
  • strain_map (ndarray[float]) – spatio-temporal strain profile.

  • strain_vectors (list[ndarray[float]]) – reduced strains per unique layer.

  • RTM (list[ndarray[complex]]) – reflection-transmission matrices for all given strains per unique layer.

  • job (Dask.job) – Dask job.

  • num_workers (int) – Dask number of workers.

Returns:

R (ndarray[float]) – inhomogeneous reflectivity.

calc_inhomogeneous_reflectivity(strains, strain_vectors, RTM)[source]#

Calculates the reflectivity of a inhomogeneous sample structure for given strain_vectors for a single time step. Similar to the homogeneous sample structure, the reflectivity of an unit cell is calculated from the reflection-transmission matrices \(H_i\) of each atom and the phase matrices between the atoms \(L_i\) in the unit cell:

\[M_{RT} = \prod_i H_i \ L_i\]

Since all layers are generally inhomogeneously strained we have to traverse all individual unit cells (\(j = 1\ldots M\)) in the sample to calculate the total reflection-transmission matrix \(M_{RT}^t\):

\[M_{RT}^t = \prod_{j=1}^M M_{RT,j}\]

The reflectivity of the \(2\times 2\) matrices for each \(q_z\) is calculates as follow:

\[R = \left|M_{RT}^t(1,2)/M_{RT}^t(2,2)\right|^2\]
Parameters:
  • strain_map (ndarray[float]) – spatio-temporal strain profile.

  • strain_vectors (list[ndarray[float]]) – reduced strains per unique layer.

  • RTM (list[ndarray[complex]]) – reflection-transmission matrices for all given strains per unique layer.

Returns:

R (ndarray[float]) – inhomogeneous reflectivity.

static calc_inhomogeneous_ref_trans_matrix(uc_indices, RT, strains, strain_vectors, RTM)[source]#

Sub-function of calc_inhomogeneous_reflectivity() and for parallel computing (needs to be static) only for calculating the total reflection-transmission matrix \(M_{RT}^t\):

\[M_{RT}^t = \prod_{j=1}^M M_{RT,j}\]
Parameters:
  • uc_indices (ndarray[float]) – unit cell indices.

  • RT (ndarray[complex]) – reflection-transmission matrix.

  • strains (ndarray[float]) – spatial strain profile for single time step.

  • strain_vectors (list[ndarray[float]]) – reduced strains per unique layer.

  • RTM (list[ndarray[complex]]) – reflection-transmission matrices for all given strains per unique layer.

Returns:

RT (ndarray[complex]) – reflection-transmission matrix.

get_all_ref_trans_matrices(*args)[source]#

Returns a list of all reflection-transmission matrices for each unique unit cell in the sample structure for a given set of applied strains for each unique unit cell given by the strain_vectors input. If this data was saved on disk before, it is loaded, otherwise it is calculated.

Parameters:

args (list[ndarray[float]], optional) – reduced strains per unique layer.

Returns:

RTM (list[ndarray[complex]]) – reflection-transmission matrices for all given strains per unique layer.

calc_all_ref_trans_matrices(*args)[source]#

Calculates a list of all reflection-transmission matrices for each unique unit cell in the sample structure for a given set of applied strains to each unique unit cell given by the strain_vectors input.

Args::
args (list[ndarray[float]], optional): reduced strains per unique

layer.

Returns:

RTM (list[ndarray[complex]]) – reflection-transmission matrices for all given strains per unique layer.

get_uc_ref_trans_matrix(uc, *args)[source]#

Returns the reflection-transmission matrix of a unit cell:

\[M_{RT} = \prod_i H_i \ L_i\]

where \(H_i\) and \(L_i\) are the atomic reflection- transmission matrix and the phase matrix for the atomic distances, respectively.

Parameters:
  • uc (UnitCell) – unit cell object.

  • args (float, optional) – strain of unit cell.

Returns:

RTM (list[ndarray[complex]])

reflection-transmission matrices for

all given strains per unique layer.

get_atom_ref_trans_matrix(atom, area, deb_wal_fac)[source]#

Calculates the reflection-transmission matrix of an atom from dynamical x-ray theory:

\[\begin{split}H = \frac{1}{\tau} \begin{bmatrix} \left(\tau^2 - \rho^2\right) & \rho \\ -\rho & 1 \end{bmatrix}\end{split}\]
Parameters:
  • atom (Atom, AtomMixed) – atom or mixed atom

  • area (float) – area of the unit cell [m²]

  • deb_wal_fac (float) – Debye-Waller factor for unit cell

Returns:

H (ndarray[complex]) – reflection-transmission matrix

get_atom_reflection_factor(atom, area, deb_wal_fac)[source]#

Calculates the reflection factor from dynamical x-ray theory:

\[\rho = \frac{-i 4 \pi \ r_e \ f(E,q_z) \ P(\theta) \exp(-M)}{q_z \ A}\]
  • \(r_e\) is the electron radius

  • \(f(E,q_z)\) is the energy and angle dispersive atomic form factor

  • \(P(q_z)\) is the polarization factor

  • \(A\) is the area in \(x-y\) plane on which the atom is placed

  • \(M = 0.5(\mbox{dbf} \ q_z)^2)\) where \(\mbox{dbf}^2 = \langle u^2\rangle\) is the average thermal vibration of the atoms - Debye-Waller factor

Parameters:
  • atom (Atom, AtomMixed) – atom or mixed atom

  • area (float) – area of the unit cell [m²]

  • deb_wal_fac (float) – Debye-Waller factor for unit cell

Returns:

rho (complex) – reflection factor

get_atom_transmission_factor(atom, area, deb_wal_fac)[source]#

Calculates the transmission factor from dynamical x-ray theory:

\[\tau = 1 - \frac{i 4 \pi r_e f(E,0) \exp(-M)}{q_z A}\]
  • \(r_e\) is the electron radius

  • \(f(E,0)\) is the energy dispersive atomic form factor (no angle correction)

  • \(A\) is the area in \(x-y\) plane on which the atom is placed

  • \(M = 0.5(\mbox{dbf} \ q_z)^2\) where \(\mbox{dbf}^2 = \langle u^2\rangle\) is the average thermal vibration of the atoms - Debye-Waller factor

Parameters:
  • atom (Atom, AtomMixed) – atom or mixed atom

  • area (float) – area of the unit cell [m²]

  • deb_wal_fac (float) – Debye-Waller factor for unit cell

Returns:

tau (complex) – transmission factor

get_atom_phase_matrix(distance)[source]#

Calculates the phase matrix from dynamical x-ray theory:

\[\begin{split}L = \begin{bmatrix} \exp(i \phi) & 0 \\ 0 & \exp(-i \phi) \end{bmatrix}\end{split}\]
Parameters:

distance (float) – distance between atomic planes

Returns:

L (ndarray[complex]) – phase matrix

get_atom_phase_factor(distance)[source]#

Calculates the phase factor \(\phi\) for a distance \(d\) from dynamical x-ray theory:

\[\phi = \frac{d \ q_z}{2}\]
Parameters:

distance (float) – distance between atomic planes

Returns:

phi (float) – phase factor

static calc_reflectivity_from_matrix(M)[source]#

Calculates the reflectivity from an \(2\times2\) matrix of transmission and reflectivity factors:

\[R = \left|M(0,1)/M(1,1)\right|^2\]
Parameters:

M (ndarray[complex]) – reflection-transmission matrix

Returns:

R (ndarray[float]) – reflectivity

static conv_with_function(y, x, handle)#

Convolutes the array \(y(x)\) with a function given by the handle on the argument array \(x\).

Parameters:
  • y (ndarray[float]) – y data.

  • x (ndarray[float]) – x data.

  • handle (@lamdba) – convolution function.

Returns:

y_conv (ndarray[float]) – convoluted data.

disp_message(message)#

Wrapper to display messages for that class.

Parameters:

message (str) – message to display.

get_hash(strain_vectors, **kwargs)#

Calculates an unique hash given by the energy \(E\), \(q_z\) range, polarization states and the strain_vectors as well as the sample structure hash for relevant x-ray parameters. Optionally, part of the strain_map is used.

Parameters:
  • strain_vectors (dict{ndarray[float]}) – reduced strains per unique layer.

  • **kwargs (ndarray[float]) – spatio-temporal strain profile.

Returns:

hash (str) – unique hash.

get_polarization_factor(theta)#

Calculates the polarization factor \(P(\vartheta)\) for a given incident angle \(\vartheta\) for the case of s-polarization (pol = 0), or p-polarization (pol = 1), or unpolarized X-rays (pol = 0.5):

\[P(\vartheta) = \sqrt{(1-\mbox{pol}) + \mbox{pol} \cdot \cos(2\vartheta)}\]
Parameters:

theta (ndarray[float]) – incidence angle.

Returns:

P (ndarray[float]) – polarization factor.

save(full_filename, data, *args)#

Save data to file. The variable name can be handed as variable argument.

Parameters:
  • full_filename (str) – full file name to data file.

  • data (ndarray) – actual data to save.

  • *args (str, optional) – variable name within the data file.

set_polarization(pol_in_state, pol_out_state)#

Sets the incoming and analyzer (outgoing) polarization.

Parameters:
  • pol_in_state (int) – incoming polarization state id.

  • pol_out_state (int) – outgoing polarization state id.

update_experiment(caller)#

Recalculate energy, wavelength, and wavevector as well as theta and the scattering vector in case any of these has changed.

\[\begin{split}\lambda & = \frac{hc}{E} \\ E & = \frac{hc}{\lambda} \\ k & = \frac{2\pi}{\lambda} \\ \vartheta & = \arcsin{\frac{\lambda q_z}{4\pi}} \\ q_z & = 2k \sin{\vartheta}\end{split}\]
Parameters:

caller (str) – name of calling method.

class udkm1Dsim.simulations.xrays.XrayDynMag(S, force_recalc, **kwargs)[source]#

Bases: Xray

Dynamical magnetic X-ray scattering simulations.

Adapted from Elzo et.al. [10] and initially realized in Project Dyna.

Original copyright notice:

Copyright Institut Neel, CNRS, Grenoble, France

Project Collaborators:

Questions to:

Parameters:
  • S (Structure) – sample to do simulations with.

  • force_recalc (boolean) – force recalculation of results.

Keyword Arguments:
  • save_data (boolean) – true to save simulation results.

  • cache_dir (str) – path to cached data.

  • disp_messages (boolean) – true to display messages from within the simulations.

  • progress_bar (boolean) – enable tqdm progress bar.

Attributes:
  • S (Structure) – sample structure to calculate simulations on.

  • force_recalc (boolean) – force recalculation of results.

  • save_data (boolean) – true to save simulation results.

  • cache_dir (str) – path to cached data.

  • disp_messages (boolean) – true to display messages from within the simulations.

  • progress_bar (boolean) – enable tqdm progress bar.

  • energy (ndarray[float]) – photon energies \(E\) of scattering light

  • wl (ndarray[float]) – wavelengths \(\lambda\) of scattering light

  • k (ndarray[float]) – wavenumber \(k\) of scattering light

  • theta (ndarray[float]) – incidence angles \(\theta\) of scattering light

  • qz (ndarray[float]) – scattering vector \(q_z\) of scattering light

  • polarizations (dict) – polarization states and according names.

  • pol_in_state (int) – incoming polarization state as defined in polarizations dict.

  • pol_out_state (int) – outgoing polarization state as defined in polarizations dict.

  • pol_in (float) – incoming polarization factor (can be a complex ndarray).

  • pol_out (float) – outgoing polarization factor (can be a complex ndarray).

  • last_atom_ref_trans_matrices (list) – remember last result of atom ref_trans_matrices to speed up calculation.

References

Methods:

get_hash(**kwargs)

Calculates an unique hash given by the energy \(E\), \(q_z\) range, polarization states as well as the sample structure hash for relevant x-ray and magnetic parameters.

set_incoming_polarization(pol_in_state)

Sets the incoming polarization factor for circular +, circular -, sigma, pi, and unpolarized polarization.

set_outgoing_polarization(pol_out_state)

Sets the outgoing polarization factor for circular +, circular -, sigma, pi, and unpolarized polarization.

homogeneous_reflectivity(*args)

Calculates the reflectivity \(R\) of the whole sample structure allowing only for homogeneous strain and magnetization.

calc_homogeneous_matrix(S, last_A, ...)

Calculates the product of all reflection-transmission matrices of the sample structure

inhomogeneous_reflectivity([strain_map, ...])

Returns the reflectivity and transmissivity of an inhomogeneously strained and magnetized sample structure for a given _strain_map_ and _magnetization_map_ in space and time for each unit cell or amorphous layer in the sample structure.

sequential_inhomogeneous_reflectivity(...)

Returns the reflectivity and transmission of an inhomogeneously strained sample structure for a given strain_map and magnetization_map in space and time.

parallel_inhomogeneous_reflectivity(...)

Returns the reflectivity and transmission of an inhomogeneously strained sample structure for a given strain_map and magnetization_map in space and time.

distributed_inhomogeneous_reflectivity(...)

This is a stub.

calc_inhomogeneous_matrix(last_A, ...)

Calculates the product of all reflection-transmission matrices of the sample structure for every atomic layer.

calc_uc_boundary_phase_matrix(uc, last_A, ...)

Calculates the product of all reflection-transmission matrices of a single unit cell for a given strain:

conv_with_function(y, x, handle)

Convolutes the array \(y(x)\) with a function given by the handle on the argument array \(x\).

disp_message(message)

Wrapper to display messages for that class.

get_atom_boundary_phase_matrix(atom, ...[, ...])

Returns the boundary and phase matrices of an atom from Elzo formalism [10].

get_polarization_factor(theta)

Calculates the polarization factor \(P(\vartheta)\) for a given incident angle \(\vartheta\) for the case of s-polarization (pol = 0), or p-polarization (pol = 1), or unpolarized X-rays (pol = 0.5):

save(full_filename, data, *args)

Save data to file.

set_polarization(pol_in_state, pol_out_state)

Sets the incoming and analyzer (outgoing) polarization.

update_experiment(caller)

Recalculate energy, wavelength, and wavevector as well as theta and the scattering vector in case any of these has changed.

calc_atom_boundary_phase_matrix(atom, ...)

Calculates the boundary and phase matrices of an atom from Elzo formalism [10].

calc_reflectivity_transmissivity_from_matrix(RT, ...)

Calculates the actual reflectivity and transmissivity from the reflectivity-transmission matrix for a given incoming and analyzer polarization from Elzo formalism [10].

calc_kerr_effect_from_matrix(RT)

Calculates the Kerr rotation and ellipticity for sigma and pi incident polarization from the reflectivity-transmission matrix independent of the given incoming and analyzer polarization from Elzo formalism [10].

calc_roughness_matrix(roughness, k_z, last_k_z)

Calculates the roughness matrix for an interface with a gaussian roughness for the Elzo formalism [10].

get_hash(**kwargs)[source]#

Calculates an unique hash given by the energy \(E\), \(q_z\) range, polarization states as well as the sample structure hash for relevant x-ray and magnetic parameters. Optionally, part of the strain_map and magnetization_map are used.

Parameters:

**kwargs (ndarray[float]) – spatio-temporal strain and magnetization profile.

Returns:

hash (str) – unique hash.

set_incoming_polarization(pol_in_state)[source]#

Sets the incoming polarization factor for circular +, circular -, sigma, pi, and unpolarized polarization.

Parameters:

pol_in_state (int) – incoming polarization state id.

set_outgoing_polarization(pol_out_state)[source]#

Sets the outgoing polarization factor for circular +, circular -, sigma, pi, and unpolarized polarization.

Parameters:

pol_out_state (int) – outgoing polarization state id.

homogeneous_reflectivity(*args)[source]#

Calculates the reflectivity \(R\) of the whole sample structure allowing only for homogeneous strain and magnetization.

The reflection-transmission matrices

\[RT = A_f^{-1} \prod_m \left( A_m P_m A_m^{-1} \right) A_0\]

are calculated for every substructure \(m\) before post-processing the incoming and analyzer polarizations and calculating the actual reflectivities as function of energy and \(q_z\).

Parameters:

args (ndarray[float], optional) – strains and magnetization for each sub-structure.

Returns:

(tuple)

  • R (ndarray[float]) - homogeneous reflectivity.

  • R_phi (ndarray[float]) - homogeneous reflectivity for opposite magnetization.

calc_homogeneous_matrix(S, last_A, last_A_phi, last_k_z, *args)[source]#

Calculates the product of all reflection-transmission matrices of the sample structure

\[RT = \prod_m \left(P_m A_m^{-1} A_{m-1} \right)\]

If the sub-structure \(m\) consists of \(N\) unit cells the matrix exponential rule is applied:

\[RT_m = \left( P_{UC} A_{UC}^{-1} A_{UC} \right)^N\]

Roughness is also included by a gaussian width

Parameters:
  • S (Structure, UnitCell, AmorphousLayer) – structure, sub-structure, unit cell or amorphous layer to calculate on.

  • last_A (ndarray[complex]) – last atom boundary matrix.

  • last_A_phi (ndarray[complex]) – last atom boundary matrix for opposite magnetization.

  • last_k_z (ndarray[float]) – last internal wave vector

  • args (ndarray[float], optional) – strains and magnetization for each sub-structure.

Returns:

(tuple)

  • RT (ndarray[complex]) - reflection-transmission matrix.

  • RT_phi (ndarray[complex]) - reflection-transmission matrix for opposite magnetization.

  • A (ndarray[complex]) - atom boundary matrix.

  • A_phi (ndarray[complex]) - atom boundary matrix for opposite magnetization.

  • A_inv (ndarray[complex]) - inverted atom boundary matrix.

  • A_inv_phi (ndarray[complex]) - inverted atom boundary matrix for opposite magnetization.

  • k_z (ndarray[float]) - internal wave vector.

inhomogeneous_reflectivity(strain_map=array([], dtype=float64), magnetization_map=array([], dtype=float64), **kwargs)[source]#

Returns the reflectivity and transmissivity of an inhomogeneously strained and magnetized sample structure for a given _strain_map_ and _magnetization_map_ in space and time for each unit cell or amorphous layer in the sample structure. If no reflectivity is saved in the cache it is caluclated. Providing the calc_type for the calculation the corresponding sub-routines for the reflectivity computation are called:

  • parallel parallelization over the time steps utilizing Dask

  • distributed not implemented in Python, but should be possible with Dask as well

  • sequential no parallelization at all

Parameters:
  • strain_map (ndarray[float], optional) – spatio-temporal strain profile.

  • magnetization_map (ndarray[float], optional) – spatio-temporal magnetization profile.

  • **kwargs

    • calc_type (str) - type of calculation.

    • dask_client (Dask.Client) - Dask client.

    • job (Dask.job) - Dask job.

    • num_workers (int) - Dask number of workers.

Returns:

(tuple)

  • R (ndarray[float]) - inhomogeneous reflectivity.

  • R_phi (ndarray[float]) - inhomogeneous reflectivity for opposite magnetization.

  • T (ndarray[float]) - inhomogeneous transmissivity.

  • T_phi (ndarray[float]) - inhomogeneous transmissivity for opposite magnetization.

sequential_inhomogeneous_reflectivity(strain_map, magnetization_map)[source]#

Returns the reflectivity and transmission of an inhomogeneously strained sample structure for a given strain_map and magnetization_map in space and time. The function calculates the results sequentially for every layer without parallelization.

Parameters:
  • strain_map (ndarray[float]) – spatio-temporal strain profile.

  • magnetization_map (ndarray[float]) – spatio-temporal magnetization profile.

Returns:

(tuple)

  • R (ndarray[float]) - inhomogeneous reflectivity.

  • R_phi (ndarray[float]) - inhomogeneous reflectivity for opposite magnetization.

  • T (ndarray[float]) - inhomogeneous transmission.

  • T_phi (ndarray[float]) - inhomogeneous transmission for opposite magnetization.

parallel_inhomogeneous_reflectivity(strain_map, magnetization_map, dask_client)[source]#

Returns the reflectivity and transmission of an inhomogeneously strained sample structure for a given strain_map and magnetization_map in space and time. The function tries to parallelize the calculation over the time steps, since the results do not depend on each other.

Parameters:
  • strain_map (ndarray[float]) – spatio-temporal strain profile.

  • magnetization_map (ndarray[float]) – spatio-temporal magnetization profile.

  • dask_client (Dask.Client) – Dask client.

Returns:

(tuple)

  • R (ndarray[float]) - inhomogeneous reflectivity.

  • R_phi (ndarray[float]) - inhomogeneous reflectivity for opposite magnetization.

  • T (ndarray[float]) - inhomogeneous transmission.

  • T_phi (ndarray[float]) - inhomogeneous transmission for opposite magnetization.

distributed_inhomogeneous_reflectivity(strain_map, magnetization_map, job, num_worker)[source]#

This is a stub. Not yet implemented in python.

Parameters:
  • strain_map (ndarray[float]) – spatio-temporal strain profile.

  • magnetization_map (ndarray[float]) – spatio-temporal magnetization profile.

  • job (Dask.job) – Dask job.

  • num_workers (int) – Dask number of workers.

Returns:

(tuple)

  • R (ndarray[float]) - inhomogeneous reflectivity.

  • R_phi (ndarray[float]) - inhomogeneous reflectivity for opposite magnetization.

calc_inhomogeneous_matrix(last_A, last_A_phi, last_k_z, strains, magnetizations)[source]#

Calculates the product of all reflection-transmission matrices of the sample structure for every atomic layer.

\[RT = \prod_m \left( P_m A_m^{-1} A_{m-1} \right)\]
Parameters:
  • last_A (ndarray[complex]) – last atom boundary matrix.

  • last_A_phi (ndarray[complex]) – last atom boundary matrix for opposite magnetization.

  • last_k_z (ndarray[float]) – last internal wave vector

  • strains (ndarray[float]) – spatial strain profile for single time step.

  • magnetizations (ndarray[float]) – spatial magnetization profile for single time step.

Returns:

(tuple)

  • RT (ndarray[complex]) - reflection-transmission matrix.

  • RT_phi (ndarray[complex]) - reflection-transmission matrix for opposite magnetization.

  • A (ndarray[complex]) - atom boundary matrix.

  • A_phi (ndarray[complex]) - atom boundary matrix for opposite magnetization.

  • A_inv (ndarray[complex]) - inverted atom boundary matrix.

  • A_inv_phi (ndarray[complex]) - inverted atom boundary matrix for opposite magnetization.

  • k_z (ndarray[float]) - internal wave vector.

calc_uc_boundary_phase_matrix(uc, last_A, last_A_phi, last_k_z, strain, magnetization, force_recalc=False)[source]#

Calculates the product of all reflection-transmission matrices of a single unit cell for a given strain:

\[RT = \prod_m \left( P_m A_m^{-1} A_{m-1}\right)\]

and returns also the last matrices \(A, A^{-1}, k_z\).

Parameters:
  • uc (UnitCell) – unit cell

  • last_A (ndarray[complex]) – last atom boundary matrix.

  • last_A_phi (ndarray[complex]) – last atom boundary matrix for opposite magnetization.

  • last_k_z (ndarray[float]) – last internal wave vector

  • strain (float) – strain of unit cell for a single time step.

  • magnetization (ndarray[float]) – magnetization of unit cell for a single time step.

  • force_recalc (boolean, optional) – force recalculation of boundary phase matrix if True. Defaults to False.

Returns:

(tuple)

  • RT (ndarray[complex]) - reflection-transmission matrix.

  • RT_phi (ndarray[complex]) - reflection-transmission matrix for opposite magnetization.

  • A (ndarray[complex]) - atom boundary matrix.

  • A_phi (ndarray[complex]) - atom boundary matrix for opposite magnetization.

  • A_inv (ndarray[complex]) - inverted atom boundary matrix.

  • A_inv_phi (ndarray[complex]) - inverted atom boundary matrix for opposite magnetization.

  • k_z (ndarray[float]) - internal wave vector.

static conv_with_function(y, x, handle)#

Convolutes the array \(y(x)\) with a function given by the handle on the argument array \(x\).

Parameters:
  • y (ndarray[float]) – y data.

  • x (ndarray[float]) – x data.

  • handle (@lamdba) – convolution function.

Returns:

y_conv (ndarray[float]) – convoluted data.

disp_message(message)#

Wrapper to display messages for that class.

Parameters:

message (str) – message to display.

get_atom_boundary_phase_matrix(atom, density, distance, force_recalc=False, *args)[source]#

Returns the boundary and phase matrices of an atom from Elzo formalism [10]. The results for a given atom, energy, \(q_z\), polarization, and magnetization are stored to RAM to avoid recalculation.

Parameters:
  • atom (Atom, AtomMixed) – atom or mixed atom.

  • density (float) – density around the atom [kg/m³].

  • distance (float) – distance towards the next atomic [m].

  • force_recalc (boolean, optional) – force recalculation of boundary phase matrix if True. Defaults to False.

  • args (ndarray[float]) – magnetization vector.

Returns:

(tuple)

  • A (ndarray[complex]) - atom boundary matrix.

  • A_phi (ndarray[complex]) - atom boundary matrix for opposite magnetization.

  • P (ndarray[complex]) - atom phase matrix.

  • P_phi (ndarray[complex]) - atom phase matrix for opposite magnetization.

  • A_inv (ndarray[complex]) - inverted atom boundary matrix.

  • A_inv_phi (ndarray[complex]) - inverted atom boundary matrix for opposite magnetization.

  • k_z (ndarray[float]) - internal wave vector.

get_polarization_factor(theta)#

Calculates the polarization factor \(P(\vartheta)\) for a given incident angle \(\vartheta\) for the case of s-polarization (pol = 0), or p-polarization (pol = 1), or unpolarized X-rays (pol = 0.5):

\[P(\vartheta) = \sqrt{(1-\mbox{pol}) + \mbox{pol} \cdot \cos(2\vartheta)}\]
Parameters:

theta (ndarray[float]) – incidence angle.

Returns:

P (ndarray[float]) – polarization factor.

save(full_filename, data, *args)#

Save data to file. The variable name can be handed as variable argument.

Parameters:
  • full_filename (str) – full file name to data file.

  • data (ndarray) – actual data to save.

  • *args (str, optional) – variable name within the data file.

set_polarization(pol_in_state, pol_out_state)#

Sets the incoming and analyzer (outgoing) polarization.

Parameters:
  • pol_in_state (int) – incoming polarization state id.

  • pol_out_state (int) – outgoing polarization state id.

update_experiment(caller)#

Recalculate energy, wavelength, and wavevector as well as theta and the scattering vector in case any of these has changed.

\[\begin{split}\lambda & = \frac{hc}{E} \\ E & = \frac{hc}{\lambda} \\ k & = \frac{2\pi}{\lambda} \\ \vartheta & = \arcsin{\frac{\lambda q_z}{4\pi}} \\ q_z & = 2k \sin{\vartheta}\end{split}\]
Parameters:

caller (str) – name of calling method.

calc_atom_boundary_phase_matrix(atom, density, distance, *args)[source]#

Calculates the boundary and phase matrices of an atom from Elzo formalism [10].

Parameters:
  • atom (Atom, AtomMixed) – atom or mixed atom.

  • density (float) – density around the atom [kg/m³].

  • distance (float) – distance towards the next atomic [m].

  • args (ndarray[float]) – magnetization vector.

Returns:

(tuple)

  • A (ndarray[complex]) - atom boundary matrix.

  • A_phi (ndarray[complex]) - atom boundary matrix for opposite magnetization.

  • P (ndarray[complex]) - atom phase matrix.

  • P_phi (ndarray[complex]) - atom phase matrix for opposite magnetization.

  • A_inv (ndarray[complex]) - inverted atom boundary matrix.

  • A_inv_phi (ndarray[complex]) - inverted atom boundary matrix for opposite magnetization.

  • k_z (ndarray[float]) - internal wave vector.

static calc_reflectivity_transmissivity_from_matrix(RT, pol_in, pol_out)[source]#

Calculates the actual reflectivity and transmissivity from the reflectivity-transmission matrix for a given incoming and analyzer polarization from Elzo formalism [10].

Parameters:
  • RT (ndarray[complex]) – reflection-transmission matrix.

  • pol_in (ndarray[complex]) – incoming polarization factor.

  • pol_out (ndarray[complex]) – outgoing polarization factor.

Returns:

(tuple)

  • R (ndarray[float]) - reflectivity.

  • T (ndarray[float]) - transmissivity.

static calc_kerr_effect_from_matrix(RT)[source]#

Calculates the Kerr rotation and ellipticity for sigma and pi incident polarization from the reflectivity-transmission matrix independent of the given incoming and analyzer polarization from Elzo formalism [10].

Parameters:

RT (ndarray[complex]) – reflection-transmission matrix.

Returns:

K (ndarray[float]) – kerr.

static calc_roughness_matrix(roughness, k_z, last_k_z)[source]#

Calculates the roughness matrix for an interface with a gaussian roughness for the Elzo formalism [10].

Parameters:
  • roughness (float) – gaussian roughness of the interface [m].

  • k_z (ndarray[float) – internal wave vector.

  • last_k_z (ndarray[float) – last internal wave vector.

Returns:

W (ndarray[float]) – roughness matrix.