Derivatives w.r.t. electric field (qcip_tools.derivatives_e)

Tools to help the manipulation of derivatives of the energy with respect to a static (F) or dynamic (D) electric field (though the first one is a special case of the second).

\[\newcommand{\sbb}[2]{\beta_{#1}\,\beta_{#2}} \newcommand{\sgg}[2]{\gamma_{#1}\,\gamma_{#2}}\]

Theory

Light can be described as an electromagnetic radiation, with time and space varying components. Generally, only the phenomena associated with its electric components, because dominant, are considered. If the space-dependent term is left constant, the electric part of the field reads

\[\mathbf{f}_i(\omega) = \mathbf{f}^0_i\,(e^{-i\omega t} + e^{i\omega t}),\]

where \(\mathbf{f}^0\) is the amplitude of the the electric field, \(\omega\) is the pulsation (equals to \(2\pi\nu\), with \(\nu\) the frequency in \(s^{-1}\)) and \(t\) the time (in \(s\)). This situation refers to the electric dipole approximation.

A material is an ensemble of particles, negatively (electrons) and positively (nuclei) charged. An oscillating electric field interacts with all particles, though the motion of the electrons is the most significant part (because nuclei have a smaller speed). This leads to an oscillating induced dipole moment, or, at the macroscopic scale, a (change of) polarization. The total dipole of a molecule (sum of the induced and permanent dipole, \(\vec{\mu}_0\)), can be described using a Taylor series. For the \(i\)-components of the dipole moment,

\[\begin{split}\begin{align} \mu_i = \mu^0_i &+ \sum^{xyz}_j \alpha_{ij}(-\omega_{\sigma};\omega_1)\,\mathbf{f}_{j}(\omega_1) \nonumber\\ &+ \frac{1}{2!} \sum^{xyz}_{jk} \beta_{ijk}(-\omega_{\sigma};\omega_1,\omega_2)\,\mathbf{f}_{j}(\omega_1)\,\mathbf{f}_{k}(\omega_2) \nonumber\\ &+ \frac{1}{3!} \sum^{xyz}_{jkl} \gamma_{ijkl}(-\omega_{\sigma};\omega_1,\omega_2,\omega_3)\,\mathbf{f}_{j}(\omega_1)\,\mathbf{f}_{k}(\omega_2)\,\mathbf{f}_{l}(\omega_3)+ \ldots \end{align}\end{split}\]

where \(\alpha\) and \(\beta\) are tensors of rank 2 and 3, called the polarizability and the first hyperpolarizability, respectively (generaly expressed in atomic units, au). \(\omega_1\), \(\omega_2\) and \(\omega_3\) are the pulsations of the fields applied in the \(j\), \(k\) and \(l\) directions, while \(\omega_\sigma = \sum_i \omega_i\). For small electric fields, the polarizability term is dominant, but the hyperpolarizability terms becomes non-negligible when the field gets large.

Sources

  • R.W. Boyd, « Nonlinear Optics (3rd edition) ». Elsevier (2008).

  • F. Castet et al. J. Chem. Phys. 136, 024506 (2012).

  • F. Castet et al. Acc. Chem. Res. 46, 2656 (2013).

API documentation

class qcip_tools.derivatives_e.BaseElectricalDerivativeTensor(tensor=None, input_fields=None, frequency='static')

Base for all derivatives of the energy with respect to the electric field.

rank()

Alias for the order

Return type:

int

response_to_electric_field(field)

Response to an electric field, located at the (0, 0, 0) position and for which the direction is given by field.

\[E_i = \sum_{ij\ldot}^{x,y,z} \chi_{ij\ldots}\,F_i\,F_j\,\ldots\]
Parameters:

field (numpy.array) – The electric field

Return type:

numpy.array

to_name(simplified=False)

Get the name of the tensor from its representation, by using REPRESENTATION.

Return type:

str

to_string(threshold=1e-05, **kwargs)

Rewritten to get a better output with this kind of tensor

class qcip_tools.derivatives_e.ElectricDipole(dipole=None)

Dipole moment, commonly written \(\mu\).

compute_properties()

Function to compute the tensor’s properties

norm()

Norm of the dipole

Return type:

float

to_string(threshold=1e-05, disable_extras=False, **kwargs)

Rewritten to add information

class qcip_tools.derivatives_e.FirstHyperpolarisabilityTensor(tensor=None, input_fields=(1, 1), frequency='static')

First hyperpolarisability tensor, commonly written \(\beta(-\omega_\sigma;\omega_1,\omega_2)\).

beta_hrs()

Hyper-Rayleigh scattering quantity:

\[\beta_{HRS}=\sqrt{\langle\beta^2_{ZZZ}\rangle + \langle\beta^2_{ZXX}\rangle}\]
Return type:

float

beta_kerr(dipole)

Measured quantity in a dc-Kerr \(\beta(-\omega;\omega,0)\) experiment.

\[\beta^K = \frac{3}{2}\,(\beta_{||}-\beta_{\perp})\]
Parameters:

dipole (numpy.ndarray) – dipole moment

Returns:

beta dc-Kerr

Return type:

float

beta_parallel(dipole)

Fetch beta value under parallel polarization of all fields:

\[\begin{split}\begin{align} \beta_{||} &= \frac{1}{5}\,\sum_i \frac{\mu_i}{||\vec{\mu}||}\,\sum_{j} \beta_{ijj} + \beta_{jij} + \beta_{jji} \\ &= \frac{3}{5}\,\sum_i\frac{\mu_i\,\beta_{i}}{||\vec{\mu}||} \end{align}\end{split}\]

where \(\beta_{i}\) is a component of the beta “tensor”.

Parameters:

dipole (numpy.ndarray|qcip_tools.derivatives_e.ElectricDipole) – dipole moment of the molecule

Return type:

float

beta_perpendicular(dipole)

Fetch beta value under perpendicular polarization of optical and static field.

\[\beta_{\perp} = \frac{1}{5}\,\sum_i \frac{\mu_i}{||\vec{\mu}||}\,\sum_{j} 2\beta_{ijj} - 3\beta_{jij} + 2 \beta_{jji}.\]
Parameters:

dipole (numpy.ndarray|qcip_tools.derivatives_e.ElectricDipole) – the dipole moment

Return type:

float

beta_squared_zxx()

Compute \(\langle\beta^2_{ZXX}\rangle\):

\[\begin{split}\begin{align} \langle\beta^2_{ZXX}\rangle &= \frac{1}{105}\, \sum_{ijk}^{xyz} 6\,\beta_{ijk}^2 + 3\,\beta_{ijj} \beta_{ikk} -2\, (\beta_{iij} \beta_{jkk} +\beta_{iij} \beta_{kjk} +\beta_{ijk} \beta_{jik})\\ &\approx \frac{1}{105}\,\sum_{ijk}^{xyz}4\, \beta_{ijk}^2-\beta_{ijj} \beta_{ikk} \end{align}\end{split}\]

(actually computed using the polarization dependant intensity)

Return type:

float

beta_squared_zzz()

Compute \(\langle\beta^2_{ZZZ}\rangle\):

\[\begin{split}\begin{align} \langle\beta^2_{ZZZ}\rangle &= \frac{1}{105}\,\sum_{ijk}^{xyz} 2\, \beta_{ijk}^2 + \beta_{ijj} \beta_{ikk} + 4\, (\beta_{iij} \beta_{jkk} + \beta_{iij} \beta_{kjk} + \beta_{ijk} \beta_{jik})\\ &\approx \frac{1}{35}\,\sum_{ijk}^{xyz} 2\,\beta_{ijk}^2 + 3 \beta_{ijj} \beta_{ikk} \end{align}\end{split}\]

(actually computed using the polarization dependant intensity)

Return type:

float

beta_vector()

return the hyperpolarizability vector

\[\beta_i = \frac{1}{3} \sum_j \beta_{ijj} + \beta_{jij} + \beta_{jji}\]
Return type:

numpy.ndarray

compute_properties(**kwargs)

Function to compute the tensor’s properties

compute_sum(what)

Shorthand to simplify the computation of invariant quantities. what is the list of quantities to compute: beta_tensor.compute_sum([(2/5, 'ijjikk'), (3/5, 'iijkjk']) is equivalent to

\[\sum_{ijk}^{xyz} \frac{2}{5}\,\sbb{ijj}{ijj} + \frac{2}{5}\,\sbb{iij}{kjk}.\]
Parameters:

what (list) – what to compute ? list of tuple of the form (fraction, string of 6 char [either i, j or k])

Return type:

float

depolarization_ratio()

Hyper-Rayleigh depolarization ratio:

\[DR = \frac{\langle\beta^2_{ZZZ}\rangle}{\langle\beta^2_{XZZ}\rangle}\]
dipolar_contribution(old_version=True)
Parameters:

old_version (bool) – Use the previous (with Kleinman’s conditions) version (coherent with betaprint’s g16)

Return type:

float

dipolar_contribution_squared(old_version=True)

Compute the square of the dipolar contribution

Parameters:

old_version (bool) – Use the previous (with Kleinman’s conditions) version

Return type:

float

is_shg()

Check if the tensor correspond to a SHG one

Return type:

bool

nonlinear_anisotropy(old_version=True)

Compute the nonlinear anisotropy:

\[\rho_{3/1}^2 = \frac{|\beta_{J=3}|^2}{|\beta_{J=1}|^2}\]
Parameters:

old_version (bool) – Use the previous (with Kleinman’s conditions) version

Return type:

float

octupolar_contribution(old_version=True)
Parameters:

old_version (bool) – Use the previous (with Kleinman’s conditions) version

Return type:

float

octupolar_contribution_squared(old_version=True)

Compute the square of the octupolar contribution

Parameters:

old_version (bool) – Use the previous (with Kleinman’s conditions) version (coherent with betaprint’s g16)

Return type:

float

polarization_angle_dependant_intensity(angle)

Compute the angle (\(\Psi\)) dependant intensity in the SHS setup.

\[\begin{split}\langle\beta^2\rangle = \frac{1}{105}\,\begin{pmatrix} 4-26\cos^2\Psi+20\cos^4\Psi\\ 4+2\cos^2\Psi-8\cos^4\Psi \\ 1-10\cos^2\Psi+12\cos^4\Psi \\ 2+8\cos^2\Psi-4\cos^4\Psi \\ 4+2\cos^2\Psi-8\cos^4\Psi \end{pmatrix}^T\,\begin{pmatrix} [g\beta^2]_A\\ [g\beta^2]_B\\ [g\beta^2]_C\\ [g\beta^2]_D\\ [g\beta^2]_E\\ \end{pmatrix}.\end{split}\]
Parameters:

angle (float) – angle (in degree)

Return type:

float

spherical_J1_contribution_squared()

Compute the \(|\beta_{J=1}|^2\) contribution

\[\begin{align} |\beta_{J=1}|^2 &= |\beta_{J=1\alpha}|^2 + |\beta_{J=1\beta}|^2 + 2\,|\beta_{J=1\alpha\beta}|^2 &= \frac{1}{5}\,\sum_{ijk}^{xyz} 2\,\sbb{ijj}{ikk} + 3\,\sbb{iij}{kjk} - 2\,\sbb{iij}{jkk} \end{align}\]
Return type:

float

spherical_J1a_contribution_squared()

Compute the \(|\beta_{J=1\alpha}|^2\) contribution

\[|\beta_{J=1\alpha}|^2 = \frac{1}{5}\,\sum_{ijk}^{xyz} 2\,\sbb{ijj}{ikk}\]
Return type:

float

spherical_J1ab_contribution_squared()

Compute the \(|\beta_{J=1\alpha\beta}|^2\) contribution

\[|\beta_{J=1\alpha\beta}|^2 = -\frac{1}{5}\,\sum_{ijk}^{xyz} \sbb{iij}{jkk}\]
Return type:

float

spherical_J1b_contribution_squared()

Compute the \(|\beta_{J=1\beta}|^2\) contribution

\[|\beta_{J=1\beta}|^2 = \frac{1}{5}\,\sum_{ijk}^{xyz} 3\,\sbb{iij}{kjk}\]
Return type:

float

spherical_J2_contribution_squared()

Compute the \(|\beta_{J=2}|^2\) contribution

\[|\beta_{J=2}|^2 = \frac{1}{3}\,\sum_{ijk}^{xyz} 2\beta_{ijk}^2 -2\,\sbb{ijk}{jik}-\sbb{ijj}{ikk} -\,\sbb{iij}{kjk} +2\,\sbb{iij}{jkk}\]
Return type:

float

spherical_J3_contribution_squared()

Compute the \(|\beta_{J=3}|^2\) contribution

\[|\beta_{J=3}|^2 = \frac{1}{15}\,\sum_{ijk}^{xyz} 5\beta_{ijk}^2 + 10\,\sbb{ijk}{jik}-\sbb{ijj}{ikk} -4\,\sbb{iij}{kjk} -4\,\sbb{iij}{jkk}\]
Return type:

float

to_string(threshold=1e-05, disable_extras=False, dipole=None, **kwargs)

Rewritten to add information

exception qcip_tools.derivatives_e.NotSHG

Raised when trying to compute a quantity from SHG experiment on a tensor which is not

exception qcip_tools.derivatives_e.NotTHG

Raised when trying to compute a quantity from THG experiment on a tensor which is not

class qcip_tools.derivatives_e.PolarisabilityTensor(tensor=None, input_fields=(1,), frequency='static')

Polarisability tensor, commonly written \(\alpha(-\omega_\sigma;\omega_1)\).

anisotropic_value()

Ansotropic value:

\[\Delta\alpha= \left[\frac{1}{2}\sum_{ij} 3\,\alpha_{ij}^2 - \alpha_{ii}\,\alpha_{jj}\right]^{1/2}.\]
Return type:

float

compute_properties()

Function to compute the tensor’s properties

isotropic_value()

Isotropic value:

\[\bar{\alpha}=\frac{1}{3}\sum_i \alpha_{ii}.\]
Return type:

float

to_string(threshold=1e-05, disable_extras=False, **kwargs)

Rewritten to add information

qcip_tools.derivatives_e.REPRESENTATIONS = {'alpha(-w;w)': 'dD', 'alpha(0;0)': 'FF', 'beta(-2w;w,w)': 'XDD', 'beta(-w;w,0)': 'dDF', 'beta(0;0,0)': 'FFF', 'beta(0;w,-w)': 'FDd', 'gamma(-2w;w,w,0)': 'XDDF', 'gamma(-3w;w,w,w)': 'XDDD', 'gamma(-w;0,0,w)': 'dFFD', 'gamma(-w;w,0,0)': 'dDFF', 'gamma(-w;w,w,-w)': 'dDDd', 'gamma(0;0,0,0)': 'FFFF', 'mu': 'F'}

Correspondence between a name and a representation

qcip_tools.derivatives_e.SIMPLIFIED_NAMES = {'alpha': 'α(0;0)', 'alpha(-w;w)': 'α(-w;w)', 'beta': 'β(0;0,0)', 'beta(-2w;w,w)': 'β(-2w;w,w)', 'beta(-w;w,0)': 'β(-w;w,0)', 'energy': 'E', 'gamma': 'γ(0;0,0,0)', 'gamma(-2w;w,w,0)': 'γ(-2w;w,w,0)', 'gamma(-3w;w,w,w)': 'γ(-3w;w,w,w)', 'gamma(-w;0,0,w)': 'γ(-w;0,0,w)', 'gamma(-w;w,0,0)': 'γ(-w;0,0,w)', 'gamma(-w;w,w,-w)': 'γ(-w;w,w,-w)', 'mu': 'µ'}

List of all simplified names

class qcip_tools.derivatives_e.SecondHyperpolarizabilityTensor(tensor=None, input_fields=(1, 1, 1), frequency='static')

Second hyperpolarisability tensor, commonly written \(\gamma(-\omega_\sigma;\omega_1,\omega_2,\omega_3)\).

compute_properties(**kwargs)

Function to compoute the tensor’s properties

compute_sum(what)

Shorthand to simplify the computation of invariant quantities. what is the list of quantities to compute: gamma_tensor.compute_sum([(2/5, 'ijklijkl'), (3/5, 'iijjkkll']) is equivalent to

\[\sum_{ijk}^{xyz} \frac{2}{5}\,\sgg{ijkl}{ijkl} + \frac{2}{5}\,\sgg{iijj}{kkll}.\]
Parameters:

what (list) – what to compute ? list of tuple of the form (fraction, string of 8 char [either i, j, k or l])

Return type:

float

depolarization_ratio()

Compute the depolarization ratio:

\[DR = \frac{\langle\gamma^2_{ZZZZ}\rangle}{\langle\gamma^2_{ZXXX}\rangle}\]
Return type:

float

gamma_kerr()

Measured quantity in a dc-Kerr (-w;w,0,0) experiment

\[\gamma^K = \frac{3}{2}\,(\gamma_{||}-\gamma_{\perp})\]
Return type:

float

gamma_parallel()

Fetch parallel isotropic average of gamma

\[\gamma_{||} = \frac{1}{15}\,\sum_{ij} \gamma_{iijj} + \gamma_{ijji} + \gamma_{ijij}\]
Return type:

float

gamma_perpendicular()

Fetch perpendicular isotropic average of gamma

\[\gamma_{\perp} = \frac{1}{15}\,\sum_{ij} 2\gamma_{iijj} - \gamma_{ijij}\]
Return type:

float

gamma_squared_zxxx()

Compute \(\langle\gamma^2_{ZXXX}\rangle\):

\[\begin{split}\langle\gamma_{ZXXX}^2\rangle &= \frac{1}{630}\sum_{ijkl}^{xyz} \left\{\begin{array}{l} 16\,\gamma_{ijkl}^2 +24\,\gamma_{ijjk} \gamma_{ikll} - 12\,\gamma_{iijk} \gamma_{jllk}\\ - 6\,(\gamma_{iijk} \gamma_{ljlk} + \gamma_{ijkl} \gamma_{jikl} ) -3\,(\gamma_{iijj} \gamma_{kllk}+\gamma_{ijjk} \gamma_{kill}) \end{array}\right\}\\ &\approx \frac{1}{630} \sum_{ijkl}^{xyz} 10\,\gamma^2_{ijkl} +3\,\gamma_{iijk}\gamma_{jkll} -3\gamma_{iijj}\gamma_{kkll}\end{split}\]

(actually computed using the polarization dependant intensity)

Return type:

float

gamma_squared_zzzz()

Compute \(\langle\gamma^2_{ZZZZ}\rangle\):

\[\begin{split}\langle\gamma_{ZZZZ}^2\rangle &= \frac{1}{315}\sum_{ijkl}^{xyz} \left\{\begin{array}{l} 2\,\gamma_{ijkl}^2 + 12\,\gamma_{iijk} \gamma_{jllk} + 6\,(\gamma_{iijk} \gamma_{ljlk} + \gamma_{ijkl} \gamma_{jikl} )\\ + 3\,(\gamma_{ijjk} \gamma_{ikll} +\gamma_{iijj} \gamma_{kllk} +\gamma_{ijjk} \gamma_{kill}) \end{array}\right\}\\ &\approx \frac{1}{315}\,\sum_{ijkl}^{xyz} 8\,\gamma^2_{ijkl} + 24\, \gamma_{iijk}\gamma_{jkll}+3\,\gamma_{iijj}\gamma_{kkll}\end{split}\]

(actually computed using the polarization dependant intensity)

Return type:

float

gamma_ths()

Compute \(\gamma_{THS}\), the quantity from a third harmonic scattering experiment.

\[\gamma_{THS} = \sqrt{\langle\gamma^2_{ZZZZ}\rangle + \langle\gamma^2_{ZXXX}\rangle}\]
Return type:

float

hexadecapolar_contribution(old_version=True)

Compute the hexadecapolar contribution

Parameters:

old_version (bool) – Use the previous (with Kleinman’s conditions) version

Return type:

float

hexadecapolar_contribution_squared(old_version=True)

Compute the square of the hexadecapolar (bless you!) contribution

Parameters:

old_version (bool) – Use the previous (with Kleinman’s conditions) version

Return type:

float

is_thg()

Check if a tensor is a THG one

Return type:

bool

isotropic_contribution(old_version=True)

Compute the isotropic contribution

Parameters:

old_version (bool) – Use the previous (with Kleinman’s conditions) version

Return type:

float

isotropic_contribution_squared(old_version=True)

Compute the square of the isotropic contribution

Parameters:

old_version (bool) – Use the previous (with Kleinman’s conditions) version

Return type:

float

polarization_angle_dependant_intensity(angle)

Compute the angle (\(\Psi\)) dependant intensity in the THS setup.

\[\begin{split}\langle\gamma^2\rangle = \frac{1}{630}\, \begin{pmatrix} 6-81\cos^2\Psi+198\cos^4\Psi-126\cos^6\Psi\\ 24-108\cos^2\Psi-72\cos^4\Psi+144\cos^6\Psi\\ 12+54\cos^2\Psi-90\cos^4\Psi+18\cos^6\Psi\\ 6-54\cos^2\Psi+36\cos^4\Psi+36\cos^6\Psi\\ 4+36\cos^2\Psi-12\cos^4\Psi-12\cos^6\Psi \\ 6-81\cos^2\Psi+198\cos^4\Psi-126\cos^6\Psi\\ 12+54\cos^2\Psi-90\cos^4\Psi+18\cos^6\Psi \end{pmatrix}^T\,\begin{pmatrix} [g\gamma^2]_A\\ [g\gamma^2]_B\\ [g\gamma^2]_C\\ [g\gamma^2]_D\\ [g\gamma^2]_E\\ [g\gamma^2]_F\\ [g\gamma^2]_G\\ \end{pmatrix}.\end{split}\]
Parameters:

angle (float) – angle (in degree)

Return type:

float

quadrupolar_contribution(old_version=True)

Compute the quadrupolar contribution

Parameters:

old_version (bool) – Use the previous (with Kleinman’s conditions) version

Return type:

float

quadrupolar_contribution_squared(old_version=True)

Compute the square of the quadrupolar contribution

Parameters:

old_version (bool) – Use the previous (with Kleinman’s conditions) version

Return type:

float

spherical_J0_contribution_squared()

Compute the \(|\gamma_{J=0}|^2\) contribution:

\[|\gamma_{J=0}|^2 = \frac{1}{5} \sum^{xyz}_{ijkl} \sgg{iijj}{kllk}\]
Return type:

float

spherical_J1_contribution_squared()

Compute the \(|\gamma_{J=1}|^2\) contribution:

\[|\gamma_{J=1}|^2 = \frac{1}{10} \sum^{xyz}_{ijkl} 3\,\sgg{ijjk}{ikll} - 3\,\sgg{ijjk}{kill}\]
Return type:

float

spherical_J2_contribution_squared()

Compute the \(|\gamma_{J=2}|^2\) contribution:

\[\begin{split}\begin{align} |\gamma_{J=2}|^2 &= |\gamma_{J=2\alpha}|^2 + |\gamma_{J=2\beta}|^2 + 2\,|\gamma_{J=2\alpha\beta}|^2 \\ &= \frac{1}{14}\,\sum_{ijkl}^{xyz} 10\,\sgg{iijk}{ljlk}+5\,\sgg{ijjk}{kill} +5\,\sgg{ijjk}{ikll}-8\,\sgg{iijk}{jllk}-4\,\sgg{iijj}{kllk} \end{align}\end{split}\]
Return type:

float

spherical_J2a_contribution_squared()

Compute the \(|\gamma_{J=2\alpha}|^2\) contribution:

\[|\gamma_{J=2\alpha}|^2 = \frac{1}{42}\,\sum_{ijkl}^{xyz} 15\,\sgg{ijjk}{kill}+15\,\sgg{ijjk}{ikll}-10\,\sgg{iijj}{kllk}\]
Return type:

float

spherical_J2ab_contribution_squared()

Compute the \(|\gamma_{J=2\alpha\beta}|^2\) contribution:

\[|\gamma_{J=2\alpha\beta}|^2 =\frac{1}{21}\,\sum_{ijkl}^{xyz} 2\,\sgg{iijj}{kllk}-6\,\sgg{iijk}{jllk}\]
Return type:

float

spherical_J2b_contribution_squared()

Compute the \(|\gamma_{J=2\beta}|^2\) contribution:

\[|\gamma_{J=2\beta}|^2 = \frac{1}{21} \sum_{ijkl}^{xyz} 15\,\sgg{iijk}{ljlk}-5\,\sgg{iijj}{kllk}\]
Return type:

float

spherical_J3_contribution_squared()

Compute the \(|\gamma_{J=3}|^2\) contribution:

\[\begin{split}|\gamma_{J=3}|^2 = \frac{1}{20} \sum^{xyz}_{ijkl}\left\{ \begin{array}{l} 15\,\gamma_{ijkl}^2+20\,\sgg{iijk}{jllk}+\sgg{ijjk}{kill}\\ -10\,\sgg{iijk}{ljlk}-15\,\sgg{ijkl}{jikl}-11\sgg{ijjk}{ikll} \end{array} \right\}\end{split}\]
Return type:

float

spherical_J4_contribution_squared()

Compute the \(|\gamma_{J=4}|^2\) contribution:

\[\begin{split}|\gamma_{J=4}|^2 = \frac{1}{140} \sum^{xyz}_{ijkl}\left\{ \begin{array}{l} 35\,\gamma_{ijkl}^2+105\,\sgg{ijkl}{jikl}+12\,\sgg{iijj}{kllk}\\ -60\,\sgg{iijk}{jllk}-30\,\sgg{iijk}{ljlk}-15\sgg{ijjk}{ikll}-15\sgg{ijjk}{kill} \end{array} \right\}\end{split}\]
Return type:

float

to_string(threshold=1e-05, disable_extras=False, dipole=None, **kwargs)

Rewritten to add information

qcip_tools.derivatives_e.convert_energy_to(e, unit)

Convert energy to something else

Parameters:
  • e (float) – the energy (in atomic unit)

  • unit (str) – the unit (nm, ev, cm-1)

qcip_tools.derivatives_e.convert_frequency_from_string(f)

Convert the string value to the value in atomic units.

Parameters:

f (str|float) – the string value

Return type:

float

qcip_tools.derivatives_e.fields_to_str(input_)

Get the string representation of a series of fields and the response, such as “-3w;w,w,w”.

Parameters:

input (tuple) – list of fields

Returns:

representation

Return type:

str