epsproc.MFBLM module

ePSproc MFBLM functions

23/09/19 Testing function caching…

19/09/19 Working & verified basic version (slow loop).

05/09/19 v1 Initial python version.
Based on original Matlab code ePS_MFBLM.m

Formalism

Should be something like this, with possible substitutions or phase swaps.

\[\begin{split}\begin{eqnarray} \beta_{L,-M}^{\mu_{i},\mu_{f}} & = & \sum_{l,m,\mu}\sum_{l',m',\mu'}(-1)^{M}(-1)^{m}(-1)^{(\mu'-\mu_{0})}\left(\frac{(2l+1)(2l'+1)(2L+1)}{4\pi}\right)^{1/2}\left(\begin{array}{ccc} l & l' & L\\ 0 & 0 & 0 \end{array}\right)\left(\begin{array}{ccc} l & l' & L\\ -m & m' & -M \end{array}\right)\nonumber \\ & \times & \sum_{P,R',R}(2P+1)(-1)^{(R'-R)}\left(\begin{array}{ccc} 1 & 1 & P\\ \mu & -\mu' & R' \end{array}\right)\left(\begin{array}{ccc} 1 & 1 & P\\ \mu_{0} & -\mu_{0} & R \end{array}\right)D_{-R',-R}^{P}(R_{\hat{n}})I_{l,m,\mu}^{p_{i}\mu_{i},p_{f}\mu_{f}}(E)I_{l',m',\mu'}^{p_{i}\mu_{i},p_{f}\mu_{f}*}(E) \end{eqnarray}\end{split}\]

Exact numerics may vary.

epsproc.MFBLM.MFBLMCalcLoop(matE, eAngs=[0, 0, 0], thres=1e-06, p=0, R=0, verbose=1)[source]

Calculate inner loop for MFBLMs, based on passed set of matrix elements (Xarray).

Loop code based on Matlab code ePSproc_MFBLM.m

This works, but not very clean or efficient - should be possible to parallelize & vectorise… but this is OK for testing, benchmarking and verification purposes.

Parameters:
  • matE (Xarray) – Contains one set of matrix elements to use for calculation. Currently assumes these are a 1D list, with associated (l,m,mu) parameters, as set by epsproc.MFBLM.mfblm().
  • eAngs ([phi,theta,chi], optional, default = [0,0,0]) – Single set of Euler angles defining polarization geometry.
  • thres (float, optional, default = 1e-4) – Threshold value for testing significance of terms. Terms < thres will be dropped.
  • p (int, optional, default p = 0) – LF polarization term. Currently only valid for p = 0
  • R (int, optional, default R = 0) – LF polarization term (from tensor contraction). Currently only valid for p = 0
  • verbose (int, optional, default 1) –

    Verbosity level:

    • 0: Silent run.
    • 1: Print basic info.
    • 2: Print intermediate C parameter array to terminal when running.
Returns:

  • BLMX (Xarray) – Set of B(L,M; eAngs, Eke) terms for supplied matrix elements, in an Xarray. For cases where no values are calculated (below threshold), return an array with B00 = 0 only.
  • Limitations & To Do
  • ——————–
  • * Currently set with (p,R) values passed, but only valid for (0,0) (not full sum over R terms as shown in formalism above.)
  • * Set to accept a single set of matrix elements (single E), assuming looping over E (and other parameters) elsewhere.
  • * Not explicitly parallelized here, should be done by calling function.
  • (Either via Xarray methods, or numba/dask…? http (//xarray.pydata.org/en/stable/computation.html#wrapping-custom-computation))
  • * Coded for ease, not efficiency - there will be lots of repeated ang. mom. calcs. when run over many sets of matrix elements.
  • .. rst-class:: strike
  • * Scale factor currently not propagated.
  • * (SF now propagated via Xarrays and implemented in main calling function )

epsproc.MFBLM.Wigner3jCached[source]

Wrapper for 3j caching with functools.lru_cache

epsproc.MFBLM.Wigner_D_element_Cached[source]

Wrapper for WignerD caching with functools.lru_cache

epsproc.MFBLM.blmXarray(BLM, Eke)[source]

Create Xarray from BLM list, format BLM = [L, M, Beta], at a single energy Eke.

Array sorting only valid for 2D BLM array, for B00=0 case pass BLM = None.

epsproc.MFBLM.mfblm(daIn, selDims={'Type': 'L'}, eAngs=[0, 0, 0], thres=0.0001, sumDims=('l', 'm', 'mu', 'Cont', 'Targ', 'Total', 'it'), SFflag=True, verbose=1)[source]

Calculate MFBLMs for a range of (E, sym) cases. Default is to calculated for all symmetries at each energy.

Parameters:
  • da (Xarray) – Contains matrix elements to use for calculation. Matrix elements will be sorted by energy and BLMs calculated for each set.
  • selDims (dict, optional, default = {'Type':'L'}) – Additional sub-selection to be applied to matrix elements before BLM calculation. Default selects just length gauge results.
  • eAngs ([phi,theta,chi], optional, default = [0,0,0]) – Single set of Euler angles defining polarization geometry.
  • thres (float, optional, default = 1e-4) – Threshold value for testing significance of terms. Terms < thres will be dropped.
  • sumDims (tuple, optional, default = ('l','m','mu','Cont','Targ','Total','it')) – Defines which terms are summed over (coherent) in the MFBLM calculation. (These are used to flatten the Xarray before calculation.) Default includes sum over (l,m), symmetries and degeneracies (but not energies).
  • SFflag (bool, default = True) – Normalise by scale factor to give X-sections (B00) in Mb
  • verbose (int, optional, default 1) –

    Verbosity level:

    • 0: Silent run.
    • 1: Print basic info.
    • 2: Print intermediate C parameter array to terminal when running.
Returns:

  • Xarray – Calculation results BLM, dims (Euler, Eke, l,m). Some global attributes are also appended.
  • Limitations
  • ———–
  • Currently set to loop calcualtions over energy only, and all symmetries.
  • Pass single {‘Cont’ (‘sym’} to calculated for only one symmetry group.)
  • TODO (In future this will be more elegant.)
  • TODO (Setting selDims in output structure needs more thought for netCDF save compatibility.)

epsproc.MFBLM.mfblmEuler(da, selDims={'Type': 'L'}, eAngs=[0, 0, 0], thres=0.0001, sumDims=('l', 'm', 'mu', 'Cont', 'Targ', 'Total', 'it'), SFflag=True, verbose=1)[source]

Wrapper for epsproc.MFBLM.mfblm() for a set of Euler angles. All other parameters are simply passed to mfblm(). Calculate MFBLMs for a range of (E, sym) cases. Default is to calculated for all symmetries at each energy.

Parameters:
  • da (Xarray) – Contains matrix elements to use for calculation. Matrix elements will be sorted by energy and BLMs calculated for each set.
  • selDims (dict, optional, default = {'Type':'L'}) – Additional sub-selection to be applied to matrix elements before BLM calculation. Default selects just length gauge results.
  • eAngs ([phi,theta,chi], optional, default = [0,0,0]) – Set of Euler angles defining polarization geometry. List or np.array, dims(N, 3). Basic Xarray support also in place, as generated by epsproc.sphCalc.setPolGeoms()
  • thres (float, optional, default = 1e-4) – Threshold value for testing significance of terms. Terms < thres will be dropped.
  • sumDims (tuple, optional, default = ('l','m','mu','Cont','Targ','Total','it')) – Defines which terms are summed over (coherent) in the MFBLM calculation. (These are used to flatten the Xarray before calculation.) Default includes sum over (l,m), symmetries and degeneracies (but not energies).
  • SFflag (bool, default = True) – Normalise by scale factor to give X-sections (B00) in Mb
  • verbose (int, optional, default 1) –

    Verbosity level:

    • 0: Silent run.
    • 1: Print basic info.
    • 2: Print intermediate C parameter array to terminal when running.
Returns:

  • Xarray – Calculation results BLM, dims (Euler(P,T,C), Eke, BLM(l,m)) - as per epsproc.util.BLMdimList(). Some global attributes are also appended.
  • Limitations
  • ———–
  • Currently set to loop calcualtions over energy only, and all symmetries.
  • Pass single {‘Cont’ (‘sym’} to calculated for only one symmetry group.)
  • TODO (In future this will be more elegant.)
  • TODO (Setting selDims in output structure needs more thought for netCDF save compatibility.)