pynx.scattering
: scattering calculations from vectors of xyz atomic positions and hkl values¶
Description¶
This module aims to help computing scattering (X-ray or neutrons) for atomic structures, especially if they are distorted or disordered.
The library uses GPU computing (although parallel CPU computing is also available as a fall-back), with the following platforms:
nVidia’s CUDA toolkit and the pyCUDA library
OpenCL language, along with pyOpenCL library.
Using GPU computing, PyNX provides fast parallel computation of scattering from large assemblies of atoms (>>1000 atoms) and 1D, 2D or 3D coordinates (>>1000) in reciprocal space.
Typical computing speeds on GPUs more than 10^11 reflections.atoms/s on nVidia cards (3.5x10^11 on 2xTitan, 2x10^11 on a GTX 690, 5x10^10 on a GTX295, 2.5x10^10 on a GTX465), more than 2 orders of magnitude faster than on a CPU.
Note that the main interest of pynx.scattering is the ability to compute scattering from any assembly of atoms (not regularly-spaced) to any set of points in reciprocal space. While a FFT will always be faster, it is much more restrictive since the FFT imposes a strict relation between the sampling in real (atomic positions) and reciprocal (hkl coordinates) space.
API Reference¶
This is the reference for scattering calculations. Note that this part of the code is older than the cdi, ptycho and wavefront modules, so that the API is much more basic.
Structure factor calculations¶
- pynx.scattering.fhkl.Fhkl_thread(h, k, l, x, y, z, occ=None, verbose=False, gpu_name='GTX 295', nbCPUthread=None, sz_imag=None, language='OpenCL', cl_platform='')¶
Compute F(hkl)=SUM_i exp(-2j*pi*(h*x_i + k*y_i + l*z_i)) or equivalently: F(k) =SUM_i exp(-2j*pi*(sx*x_i + sy*y_i + sz*z_i))
nbCPUthread can be used only for CPU computing, when no GPU is available. nbCPUthread can be set to the number of cores available. Using None (the default) makes the program recognize the number of available processors/cores
If sz_imag is not equal to None, then it means that the scattering vector s has an imaginary component (absorption), e.g. because of grazing incidence condition.
- Parameters:
h – numpy vector of h coordinates (reciprocal lattice units)
k – numpy vector of k coordinates (reciprocal lattice units)
l – numpy vector of l coordinates (reciprocal lattice units)
x – numpy vector of x coordinates (fractional coordinates in crystal space)
y – numpy vector of x coordinates (fractional coordinates in crystal space)
z – numpy vector of x coordinates (fractional coordinates in crystal space)
occ – numpy vector of occupancy. Can be None
verbose – if True, will show some computing/threading messages
gpu_name – name of the GPu to be used
nbCPUthread – if using a CPU, number of threads to use
sz_imag – see grazing incidence example
language – either ‘cuda’, ‘opencl’ or ‘cpu’
cl_platform – the opencl platform to use
- Returns:
a tuple with (calculated structure factor as a numpy complex64 array, elapsed time in seconds)
- pynx.scattering.fhkl.Fhkl_gold(h, k, l, x, y, z, occ=None, verbose=False, dtype=<class 'numpy.float64'>)¶
- Compute reference value (using CPU) of:
F(hkl)=SUM_i exp(2j*pi*(h*x_i + k*y_i + l*z_i))
Grazing incidence diffraction¶
This module requires the cctbx library.
Thomson scattering factor¶
- pynx.scattering.fthomson.f_thomson(s, element)¶
Calculate the Thomson scattering factor for one atom or ion, for a given value of sin(theta)/lambda.
- Parameters:
s – \(sin(\vartheta)/\lambda\)
element – string of the element name (‘H’, ‘Fe’, ‘Cu1+’,…)
- Returns: