Quasi-Neutral Initialization

Section author: Axel Huebl

In order to initialize the electro-magnetic fields self-consistently, one needs to fulfill Gauss’s law \(\vec \nabla \cdot \vec E = \frac{\rho}{\epsilon_0}\) (and \(\vec \nabla \cdot \vec B = 0\)). The trivial solution to this equation is to start field neutral by microscopically placing a charge-compensating amount of free electrons on the same position as according ions.

Fully Ionized Ions

For fully ionized ions, just use ManipulateDerive in speciesInitialization.param and derive macro-electrons \(1:1\) from macro-ions but increase their weighting by \(1:Z\) of the ion.

using InitPipeline = pmacc::mp_list<
    /* density profile from density.param and
     *     start position from particle.param */
    CreateDensity<
        densityProfiles::YourSelectedProfile,
        startPosition::YourStartPosition,
        Carbon
    >,
    /* create a macro electron for each macro carbon but increase its
     *     weighting by the ion's proton number so it represents all its
     *     electrons after an instantanous ionization */
    ManipulateDerive<
        manipulators::ProtonTimesWeighting,
        Carbon,
        Electrons
    >
>;

If the Carbon species in this example has an attribute boundElectrons (optional, see speciesAttributes.param and speciesDefinition.param) and its value is not manipulated the default value is used (zero bound electrons, fully ionized). If the attribute boundElectrons is not added to the Carbon species the charge state is considered constant and taken from the chargeRatio< ... > particle flag.

Partly Ionized Ions

For partial pre-ionization, the FoilLCT example shows a detailed setup. First, define a functor that manipulates the number of bound electrons in particle.param, e.g. to twice pre-ionized.

#include "picongpu/particles/traits/GetAtomicNumbers.hpp"
// ...

namespace manipulators
{
    //! ionize ions twice
    struct TwiceIonizedImpl
    {
        template< typename T_Particle >
        DINLINE void operator()(
            T_Particle& particle
        )
        {
            constexpr float_X protonNumber =
                GetAtomicNumbers< T_Particle >::type::numberOfProtons;
            particle[ boundElectrons_ ] = protonNumber - float_X( 2. );
        }
    };

    //! definition of TwiceIonizedImpl manipulator
    using TwiceIonized = generic::Free< TwiceIonizedImpl >;

} // namespace manipulators

Then again in speciesInitialization.param set your initialization routines to:

using InitPipeline = pmacc::mp_list<
    /* density profile from density.param and
     *     start position from particle.param */
    CreateDensity<
        densityProfiles::YourSelectedProfile,
        startPosition::YourStartPosition,
        Carbon
    >,
    /* partially pre-ionize the carbons by manipulating the carbon's
     *     `boundElectrons` attribute,
     *     functor defined in particle.param: set to C2+ */
    Manipulate<
        manipulators::TwiceIonized,
        Carbon
    >,
    /* does also manipulate the weighting x2 while deriving the electrons
     *     ("twice pre-ionized") since we set carbon as C2+ */
    ManipulateDerive<
        manipulators::binary::UnboundElectronsTimesWeighting,
        Carbon,
        Electrons
    >
>;