Skip to content

Conversation

@n01r
Copy link
Member

@n01r n01r commented May 3, 2017

Note: This also affects 0.2.4, we might only fix it in 0.3.0

Fixes the calculation of the effective principal quantum number nEff for
the implemented ADK ionization model. The quantity nEff arises from scaling
an ionic potential of an arbitrary configuration and charge state Z to the
hydrogenic Bohr model. In the Bohr model:

n = Z / sqrt(2 * E_n)

where E_n denotes the energy of the n-th level/shell.
In order to extrapolate that to arbitrary ions one has to do the following
replacements:

Z   (nuclear charge)    --> Z_i (ionic residual charge)
E_n (n-th shell energy) --> E_i (ionization energy of the state)

Here atomic residual charge means charge AFTER ionization since this is the
charge that the electron-to-be-ionized 'feels'.

Previously it was wrongly implemented since the bare nuclear charge was used to calculate the effective principal quantum number

In the following comparison one can see the ionization rates for the transition from C II to C III (meaning 1+ to 2+).
comparison_carbon_ii
While previously (orange curve ADK C H-like) with the old nEff values this carbon state was ionized earlier than hydrogen that is not the case anymore (green curve ADK C (fixed)).

state ionization energy (eV)
H 13.6
C I 11.26
C II 24.36

While in the hydrogen case the effective charge is just 1 for C II it should be either 2 (proton number minus number of bound electrons in ionized state) or 3.136 when accounting for shielding of the orbitals via SCF methods
(see: Clementi, E.; Raimondi, D. L. (1963). "Atomic Screening Constants from SCF Functions". J. Chem. Phys. 38 (11): 2686–2689.)

This affects all previous simulations done with the ADK ionization model.
Especially low charge states are affected, as can be seen in the example of values of nEff for carbon ions.

state C I C II C III C IV C V C VI
nEff (fixed) 1.099 1.494 1.599 1.837 0.931 1.000
nEff (old) 6.594 4.483 3.197 2.755 1.117 1.000

That previously led to an overestimation of ionization for low charge states.

@n01r n01r added this to the 0.3.0: Release Candidate milestone May 3, 2017
@n01r n01r requested a review from ax3l May 3, 2017 12:51
@n01r n01r added affects latest release a bug that affects the latest stable release bug a bug in the project's code labels May 3, 2017
@n01r n01r force-pushed the fix-ionChargeInADK branch from ba7a30f to a59a4fd Compare May 3, 2017 12:59
@ax3l ax3l self-assigned this May 3, 2017
@ax3l
Copy link
Member

ax3l commented May 3, 2017

offline discussion: the H-like models will be removed as they only made sense for the last bound e- and are error prune (usage) plus lead to the same result as "fixed" for the last e-.

The Zeff model in BSI (and an added follow-up for ADK) are interesting for users to tune their rates but are not necessarily better, as we can e.g. see for the BSI/ADK thresholds where a 24eV C2+ e- is estimated with the same threshold as a 13.6 eV Hydrogen e-, which is obviously weird and caused by the common approximation to increase Zeff only per closed shell. Luckily, we allow the user to change that at will, but we need to document that Zeff does not necessarily mean its "better" as in the example here.

for this PR, we will add a description for the sphinx docs with the models and the image above (without the two orange lines) and a follow-up PR will include add ADK with Zeff, although it might be as (consistently) bad as BSI with Zeff, but at least it's configurable.

@n01r n01r force-pushed the fix-ionChargeInADK branch from a59a4fd to 195184c Compare May 3, 2017 14:37
@ax3l ax3l mentioned this pull request May 3, 2017
3 tasks
@n01r n01r force-pushed the fix-ionChargeInADK branch from 555b542 to accd5a8 Compare May 8, 2017 11:54
Copy link
Member Author

@n01r n01r left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Still some kernel optimizations to do.

plt.draw()
plt.show()

plt.savefig("field_ionization_effective_potentials.svg") No newline at end of file
Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

missing EOF newline

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

fixed

E. Clementi and D. Raimondi.
*Atomic Screening Constant from SCF Functions. II. Atoms with 37 to 86 Electrons*,
The Journal of Chemical Physics 47, 1300-1307 (1967)
https://dx.doi.org/10.1063/1.1712084 No newline at end of file
Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

missing EOF newline

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

fixed

/* the charge that attracts the electron that is to be ionized:
* equals `protonNumber - #allInnerElectrons`
*/
float_X effectiveCharge = chargeState + float_X(1.0);
Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

float_X const

as well as above where:

const float_X protonNumber = GetAtomicNumbers::type::numberOfProtons;
float_X chargeState = attribute::getChargeState(parentIon);

/* nameless variable for convenience dFromADK*/
float_X dBase = float_X(4.0) * util::cube(protonNumber) / (eInAU * util::quad(nEff)) ;
float_X dBase = float_X(4.0) * util::cube(effectiveCharge) / (eInAU * util::quad(nEff)) ;
float_X dFromADK = math::pow(dBase,nEff);
Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

float_X const

float_X nEff = effectiveCharge / math::sqrt(float_X(2.0) * iEnergy );
/* nameless variable for convenience dFromADK*/
float_X dBase = float_X(4.0) * util::cube(protonNumber) / (eInAU * util::quad(nEff)) ;
float_X dBase = float_X(4.0) * util::cube(effectiveCharge) / (eInAU * util::quad(nEff)) ;
Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

float_X const

float_X effectiveCharge = chargeState + float_X(1.0);
/* effective principal quantum number (unitless) */
float_X nEff = protonNumber / math::sqrt(float_X(2.0) * iEnergy );
float_X nEff = effectiveCharge / math::sqrt(float_X(2.0) * iEnergy );
Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

float_X const


const float_X pi = precisionCast<float_X>(M_PI);
/* electric field in atomic units - only absolute value */
float_X eInAU = math::abs(eField) / ATOMIC_UNIT_EFIELD;
Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

float_X const

@ax3l ax3l force-pushed the fix-ionChargeInADK branch from 79824fc to e3126e5 Compare May 8, 2017 13:51
Fixes the calculation of the effective principal quantum number `nEff` for
the implemented ADK ionization model. The quantity `nEff` arises from scaling
an ionic potential of an arbitrary configuration and charge state Z to the
hydrogenic Bohr model. In the Bohr model:

    n = Z / sqrt(2 * E_n)

where E_n denotes the energy of the n-th level/shell.
In order to extrapolate that to arbitrary ions one has to do the following
replacements:

    Z   (nuclear charge)    --> Z_i (ionic residual charge)
    E_n (n-th shell energy) --> E_i (ionization energy of the state)

Here atomic residual charge means charge AFTER ionization since this is the
charge that the electron to be ionized 'feels'.
@ax3l ax3l force-pushed the fix-ionChargeInADK branch 2 times, most recently from 55b2306 to 6f75a2c Compare May 8, 2017 14:07
Copy link
Member

@ax3l ax3l left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

pushed my changed directly to the branch

@ax3l ax3l requested a review from PrometheusPi May 8, 2017 14:08
@n01r
Copy link
Member Author

n01r commented May 8, 2017

Thanks! The script refactoring is very helpful!

Copy link
Member

@PrometheusPi PrometheusPi left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

good job - some typos and missing python comments should be fixed


.. moduleauthor:: Marco Garten

Implemented LTE Model: Thomas-Fermi Ionization according to [More1985]_
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

OPTIONAL:
Please remove capitalization model, ionization


Get started here https://github.com/ComputationalRadiationPhysics/picongpu/wiki/Ionization-in-PIConGPU

NLTE Models
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

OPTIONAL:
either write model or capitalize line 6

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

this header is fine, line 6 needs capitalization


.. note::

Most of the calculations and formulae in this section of the docs are done in the system of **Atomic Units (AU)**.
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

OPTIONAL:
**Atomiv Units (AU)** system?

================ =====
length :math:`5.292 \cdot 10^{-11}\,\mathrm{m}`
time :math:`2.419 \cdot 10^{-17}\,\mathrm{s}`
energy :math:`4.360 \cdot 10^{-19}\,\mathrm{J}\quad` (= 27.21 eV = 1 Rydberg)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

REQUIRED:
The Hartree energy is 4.3 10^-18 J

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

but time is correct

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

thanks for catching that typo!

\Gamma_\mathrm{K} = \frac{\left(6 \pi\right)^{1/2}}{2^{5/4}} E_\mathrm{ip} \left( \frac{F}{(2 E_\mathrm{ip})^{3/2}} \right)^{1/2} \exp\left(-\frac{2 \left(2 E_\mathrm{ip}\right)^{3/2}}{3 F}\right)
According to the equation (9) in [BauerMulser1999]_ the Keldysh ionization rate has been implemented. See also [Keldysh]_ for the original work.
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

OPTIONAL:
Better turn parts around:
The Keldysh ionization rate has been implemented according to the equation (9) in [BauerMulser1999]_.

# Atomic units
E_AU = 27.2 # eV
F_AU = 5.1422e11 # V/m
I_AU = 3.5095e16 # W/cm^2
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

QUESTION:
If this is in A.U. why do you give a SI unit?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

that is the unit ^^

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Okay - I was confused by the relation but know it makes sense

plt.ylabel("ionization rate $\Gamma$ [Hz]")
plt.xlabel("field strength $F$ [AU = 5.1422$\cdot 10^{11}$ V/m]")
plt.legend(loc="best")

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

OPTIONAL:
you could add plt.tight_layout()

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

offline: @ax3l says no

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Okay

def V_eff(x, Z_eff, F):
"""
Effective radially symmetric nuclear potential under the influence of a
homogenous external electric field.
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

REQUIRED:
typo homogeneous

Effective radially symmetric nuclear potential under the influence of a
homogenous external electric field.
Unit: 1 AU = 27.2 eV
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

REQUIRED:
Please add brief comment on parameters

plt.ylim([-3., 1.])
plt.ylabel("$V_\mathrm{eff}$ [AU = Rydberg]")
plt.xlabel("$x$ [AU = Bohr radii]")

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

OPTIONAL:
add plt.tight_layout()

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

as discussed offline: will not implement as it should rather be set manually

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Okay

Collisional Ionization
======================

LTE models
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Model needs capitalization

@ax3l
Copy link
Member

ax3l commented May 8, 2017

@n01r pls ammend the changes (are all docs-only) to your last commit when finished (careful: check the branch out fresh first since I force-pushed)

@n01r
Copy link
Member Author

n01r commented May 8, 2017

check the branch out fresh first

Yep, did that immediately.

Thanks for the feedback! I hope it can be merged now 👍

@n01r n01r force-pushed the fix-ionChargeInADK branch from 6f75a2c to af6368d Compare May 8, 2017 16:48
- Split up `ionization.rst` into `field_ionization.rst` and
  `collisionalIonization.rst`
- Added equations for the implemented tunneling models
- Added images showing a comparison of ionization rates as well as effective
  potentials
- Added Python scripts that produce these images
- Updated the references
@n01r n01r force-pushed the fix-ionChargeInADK branch from af6368d to baca9d8 Compare May 8, 2017 22:12
@ax3l
Copy link
Member

ax3l commented May 9, 2017

@PrometheusPi pls update your review :)

Copy link
Member

@PrometheusPi PrometheusPi left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

looks good - thank you for the fixes

@PrometheusPi PrometheusPi merged commit 0d43681 into ComputationalRadiationPhysics:dev May 9, 2017
ax3l added a commit to ax3l/picongpu that referenced this pull request May 11, 2017
Adds the following changes to the `release-0.3.0` branch:

- Multiple ionizers per species:
  `ionizer< ... >` became a sequence of `ionizers< ... >` ComputationalRadiationPhysics#1999
- Thomas-Fermi:
  - fixes in domain decomposition ComputationalRadiationPhysics#2007
  - fixes compile issue (missing include) ComputationalRadiationPhysics#2003
- BSI models restructured ComputationalRadiationPhysics#2013
- ADK: fix effective principal quantum number `nEff` ComputationalRadiationPhysics#2011
- multiple ionization algorithms can be applied per species,
  e.g. cut-off barrier suppression ionization (BSI),
  probabilistic field ionization (ADK) and collisional ionization ComputationalRadiationPhysics#1999
- FieldTmp: Gather support to fill GUARD ComputationalRadiationPhysics#2009
- Fix undefined device variable in tilted laser pulse ComputationalRadiationPhysics#2017
- CompileTime Accessor: Type (::type) ComputationalRadiationPhysics#1998
- TBG Macros: Outdated Comment ComputationalRadiationPhysics#2004

Increases the state of the `0.3.0` release candidate to `rc3`.
ax3l added a commit to ax3l/picongpu that referenced this pull request May 11, 2017
Adds the following changes to the `release-0.3.0` branch:

- Multiple ionizers per species:
  `ionizer< ... >` became a sequence of `ionizers< ... >` ComputationalRadiationPhysics#1999
- Thomas-Fermi:
  - fixes in domain decomposition ComputationalRadiationPhysics#2007
  - fixes compile issue (missing include) ComputationalRadiationPhysics#2003
- BSI models restructured ComputationalRadiationPhysics#2013
- ADK: fix effective principal quantum number `nEff` ComputationalRadiationPhysics#2011
- multiple ionization algorithms can be applied per species,
  e.g. cut-off barrier suppression ionization (BSI),
  probabilistic field ionization (ADK) and collisional ionization ComputationalRadiationPhysics#1999
- FieldTmp: Gather support to fill GUARD ComputationalRadiationPhysics#2009
- Fix undefined device variable in tilted laser pulse ComputationalRadiationPhysics#2017
- CompileTime Accessor: Type (::type) ComputationalRadiationPhysics#1998
- TBG Macros: Outdated Comment ComputationalRadiationPhysics#2004

Increases the state of the `0.3.0` release candidate to `rc3`.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

affects latest release a bug that affects the latest stable release bug a bug in the project's code

Projects

None yet

Development

Successfully merging this pull request may close these issues.

3 participants