two ions and electron species

preview example

Author: psychocoderHPC

Diff for: examples/LaserWakefield/cmakeFlags

@@ -30,16 +30,7 @@
 #   - increase by 1, no gaps
 
 flags[0]="-DCUDA_ARCH=sm_20"
-flags[1]="-DCUDA_ARCH=sm_20 -DPARAM_OVERWRITES:LIST=-DPARAM_FIELDSOLVER=fieldSolverLehe;-DPARAM_PARTICLEPUSHER=Vay"
-flags[2]="-DCUDA_ARCH=sm_20 -DPARAM_OVERWRITES:LIST=-DPARAM_CURRENTSOLVER=Esirkepov;-DPARAM_PARTICLESHAPE=CIC"
-flags[3]="-DCUDA_ARCH=sm_20 -DPARAM_OVERWRITES:LIST=-DPARAM_CURRENTSOLVER=VillaBune;-DPARAM_PARTICLESHAPE=CIC"
-flags[4]="-DCUDA_ARCH=sm_20 -DPARAM_OVERWRITES:LIST=-DPARAM_PRECISION=precision64Bit"
-flags[5]="-DCUDA_ARCH=sm_35 -DPARAM_OVERWRITES:LIST=-DENABLE_CURRENT=0"
-flags[6]="-DCUDA_ARCH=sm_20 -DPARAM_OVERWRITES:LIST=-DPARAM_DIMENSION=DIM2"
-flags[7]="-DCUDA_ARCH=sm_20 -DPARAM_OVERWRITES:LIST=-DPARAM_CURRENTSOLVER=Esirkepov;-DPARAM_PARTICLESHAPE=CIC;-DPARAM_DIMENSION=DIM2"
-flags[8]="-DCUDA_ARCH=sm_20 -DPARAM_OVERWRITES:LIST=-DPARAM_PRECISION=precision64Bit;-DPARAM_DIMENSION=DIM2"
-flags[9]="-DCUDA_ARCH=sm_35 -DPARAM_OVERWRITES:LIST=-DENABLE_CURRENT=0;-DPARAM_DIMENSION=DIM2"
-flags[10]="-DCUDA_ARCH=sm_35 -DPARAM_OVERWRITES:LIST=-DPARAM_IONS=1;-DPARAM_IONIZATION=1"
+
 ################################################################################
 # execution
 

Diff for: examples/LaserWakefield/include/simulation_defines/param/gasConfig.param

@@ -43,39 +43,46 @@ const double GAS_DENSITY_SI = 1.e25;
 //##########################################################################
 namespace gasProfiles
 {
-struct GaussianParameter
-{
-    /** Gas Formular:
-                const float_X exponent = fabs((y - GAS_CENTER) / GAS_SIGMA);
-                const float_X density = __expf(GAS_FACTOR*__powf(exponent, GAS_POWER));
 
-        takes GAS_CENTER_LEFT for y < GAS_CENTER_LEFT,
-              GAS_CENTER_RIGHT for y > GAS_CENTER_RIGHT,
-        and exponent = float_X(0.0) for GAS_CENTER_LEFT < y < GAS_CENTER_RIGHT */
-    static const float_X GAS_FACTOR = -1.0;
-    static const float_X GAS_POWER = 4.0;
+/** Middle of the constant sphere
+ *  unit: meter */
+CONST_VECTOR(float_64, simDim, SphereFlanksParam_center, 8.0e3, 8.0e3, 8.0e3);
 
+struct SphereFlanksParam
+{
     /** height of vacuum area on top border
      *      this vacuum is really important because of the laser initialization,
      *      which is done in the first cell of the simulation
      *  unit: cells */
-    static const uint32_t VACUUM_CELLS_Y = 50;
+    static const uint32_t vacuum_y_cells = 50;
 
     struct SI
     {
-        /** The central position of the gas distribution
+        /** Radius of the constant sphere
          *  unit: meter */
-        static const float_64 GAS_CENTER_LEFT = 8.0e-5;
-        static const float_64 GAS_CENTER_RIGHT = 10.0e-5;
-        /** the distance from GAS_CENTER until the gas density decreases to its 1/e-th part
+        static const float_64 r = 1.0e3;
+
+        /** Inner radius if you want to build a shell/ring
          *  unit: meter */
-        static const float_64 GAS_SIGMA_LEFT = 8.0e-5;
-        static const float_64 GAS_SIGMA_RIGHT = 8.0e-5;
+        static const float_64 ri = 0.0;
+
+        /** Middle of the constant sphere
+         *  unit: meter */
+        const SphereFlanksParam_center_t center;
+
+        /** Parameters for the exponential slope
+         *  For radius > GAS_R_SI:
+         *    let radius' = radius - r
+         *    \rho = exp[ - radius' * exponent ]
+         *  unit: 1/m
+         */
+        static const float_64 exponent = 1.0e-3;
+
     };
 };
 
-/* definition of GasGaussian*/
-typedef GaussianImpl<GaussianParameter> Gaussian;
+/* definition of gas sphere with flanks*/
+typedef SphereFlanksImpl<SphereFlanksParam> SphereFlanks;
 
 }//namespace gasProfiles
 

Diff for: examples/LaserWakefield/include/simulation_defines/param/particleConfig.param

@@ -41,7 +41,24 @@ namespace particles
  *  unit: none */
 const float_X MIN_WEIGHTING = 10.0;
 
-const uint32_t TYPICAL_PARTICLES_PER_CELL = 2;
+const uint32_t TYPICAL_PARTICLES_PER_CELL = 10;
+
+namespace manipulators
+{
+    struct SixTimesWeightingFunctor
+    {
+        template<typename T_Particle>
+        DINLINE void operator()(T_Particle& particle)
+        {
+            particle[weighting_]*=float_X(6.0);
+        }
+    };
+
+    /* definition of SetDrift start*/
+    typedef FreeImpl<SixTimesWeightingFunctor> SixTimesWeighting;
+
+
+} //namespace manipulators
 
 namespace startPosition
 {
@@ -55,22 +72,11 @@ struct RandomParameter
 /* definition of random particle start*/
 typedef RandomImpl<RandomParameter> Random;
 
-struct QuietParameter
-{
-    /** Count of particles per cell per direction at initial state
-     *  unit: none */
-    typedef typename mCT::shrinkTo<mCT::Int<1, TYPICAL_PARTICLES_PER_CELL, 1>, simDim>::type numParticlesPerDimension;
-};
-
-/* definition of random particle start*/
-typedef QuietImpl<QuietParameter> Quiet;
-
-
 } //namespace startPosition
 
 namespace manipulators
 {
-    
+
     struct BoundElectronsWhenNeutral
     {
         /* space for attribute initial parameters */

Diff for: examples/LaserWakefield/include/simulation_defines/param/physicalConstants.param

@@ -0,0 +1,85 @@
+/**
+ * Copyright 2013 Axel Huebl, Heiko Burau, Rene Widera
+ *
+ * This file is part of PIConGPU.
+ *
+ * PIConGPU is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * PIConGPU is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with PIConGPU.
+ * If not, see <http://www.gnu.org/licenses/>.
+ */
+
+
+
+#pragma once
+
+namespace picongpu
+{
+    const float_64 PI = 3.141592653589793238462643383279502884197169399;
+
+    /** Threshold used for calculations that want to separate between
+     *  high-precision formulas for relativistic and non-relativistic
+     *  use-cases, e.g. energy-binning algorithms. */
+    const float_X GAMMA_THRESH = float_X(1.005);
+
+    namespace SI
+    {
+        /** unit: m / s */
+        const float_64 SPEED_OF_LIGHT_SI = 2.99792458e8;
+
+        /** unit: N / A^2 */
+        const float_64 MUE0_SI = PI * 4.e-7;
+        /** unit: C / (V m) */
+        const float_64 EPS0_SI = 1.0 / MUE0_SI / SPEED_OF_LIGHT_SI
+            / SPEED_OF_LIGHT_SI;
+
+        // Electron properties
+        /** unit: kg */
+        const float_64 ELECTRON_MASS_SI = 9.109382e-31;
+        /** unit: C */
+        const float_64 ELECTRON_CHARGE_SI = -1.602176e-19;
+
+        // Ion / Proton / Positron properties
+        /** unit: kg */
+        const float_64 ION_MASS_SI = 1.6726217e-27;
+
+        const float_64 ION_MASS_C_SI = 1.660538921-27;
+    }
+
+    // converts
+    //
+    // UNIT_A to UNIT_B
+    //
+    // CONVENTION: WE DO NOT CONVERT FROM ANY STRANGE UNIT TO UNITLESS UNITS DIRECTLY!
+    //             convert steps: INPUT -> float_64_convert to SI -> float_64_convert to unitless
+    //                                  -> cast to float
+    // WE DO NOT define "UNIT_ENERGY_keV" or something similar! Never!
+    // Stay SI, stay free ;-)
+    //
+    // example:
+    //   // some particle physicist beloved input:
+    //   const float_64 An_Arbitrary_Energy_Input_keV = 30.0; // unit: keV
+    //
+    //   // first convert to SI (because SI stays our standard Unit System!)
+    //   const float_64 An_Arbitrary_Energy_Input_SI = An_Arbitrary_Energy_Input_keV * UNITCONV_keV_to_Joule // unit: Joule
+    //
+    //   // now the "real" convert to our internal unitless system
+    //   const float_X An_Arbitrary_Energy_Input = float_X(An_Arbitrary_Energy_Input_SI / UNIT_ENERGY) // unit: none
+    //
+    // As a convention, we DO NOT use the short track:
+    //   const float_64 An_Arbitrary_Energy_Input_keV = 30.0; // unit: keV
+    //   const float_X An_Arbitrary_Energy_Input = float_X(An_Arbitrary_Energy_Input_SI * UNITCONV_keV_to_Joule / UNIT_ENERGY) // unit: none
+    //
+    const float_64 UNITCONV_keV_to_Joule = 1.60217646e-16;
+    const float_64 UNITCONV_Joule_to_keV = (1.0 / UNITCONV_keV_to_Joule);
+
+}

Diff for: examples/LaserWakefield/include/simulation_defines/param/species.param

@@ -34,17 +34,7 @@
 /*enable (1) or disable (0) electrons*/
 #define ENABLE_ELECTRONS 1
 /*enable (1) or disable (0) ions*/
-#ifndef PARAM_IONS
-#define PARAM_IONS 0
-#endif
-
-#define ENABLE_IONS PARAM_IONS
-
-/*enable (1) or disable (0) ionization*/
-#ifndef PARAM_IONIZATION
-#define PARAM_IONIZATION 0
-#endif
-
+#define ENABLE_IONS 1
 
 namespace picongpu
 {
@@ -76,7 +66,7 @@ typedef FieldToParticleInterpolation<UsedParticleShape, AssignedTrilinearInterpo
  *   without optimization (~4x slower and needs more shared memory)
  */
 #ifndef PARAM_CURRENTSOLVER
-#define PARAM_CURRENTSOLVER Esirkepov
+#define PARAM_CURRENTSOLVER ZigZag
 #endif
 typedef currentSolver::PARAM_CURRENTSOLVER<UsedParticleShape> UsedParticleCurrentSolver;
 

Diff for: examples/LaserWakefield/include/simulation_defines/param/speciesDefinition.param

@@ -75,8 +75,7 @@ typedef
 typename MakeSeq<
 DefaultParticleAttributes,
 AttributMomentum_mt1,
-AttributRadiationFlag,
-boundElectrons
+AttributRadiationFlag
 >::type AttributeSeqIons;
 
 /*########################### end particle attributes ########################*/
@@ -95,13 +94,22 @@ typedef bmpl::vector<
 
 /* define species: electrons */
 typedef Particles<ParticleDescription<
-    bmpl::string<'e'>,
+    bmpl::string<'e','H'>,
     SuperCellSize,
     AttributeSeqElectrons,
     ParticleFlagsElectrons,
     typename MakeSeq<CommunicationId<PAR_ELECTRONS> >::type  >
 > PIC_Electrons;
 
+/* define species: electrons */
+typedef Particles<ParticleDescription<
+    bmpl::string<'e','C'>,
+    SuperCellSize,
+    AttributeSeqElectrons,
+    ParticleFlagsElectrons,
+    typename MakeSeq<CommunicationId<PAR_ELECTRONS+1> >::type  >
+> PIC_ElectronsC;
+
 /*--------------------------- ions -------------------------------------------*/
 
 /*! Specify (chemical) element
@@ -120,41 +128,51 @@ typedef Particles<ParticleDescription<
 struct Hydrogen
 {
     static const float_X numberOfProtons  = 1.0;
-    static const float_X numberOfNeutrons = 0.0;
+    static const float_X numberOfNeutrons = 1.0;
 };
 
-/*! Ionization Model Configuration ----------------------------------------
- *
- * For development purposes: ---------------------------------------------
- * - None : no particle is ionized
- * - BSI : simple barrier suppression ionization
- *
- * Usage: Add a flag to the list of particle flags that has the following structure
- *
- *        ionizer<IonizationModel<Species2BCreated> >
- */
-
-typedef particles::shapes::CIC IonShape;
 typedef bmpl::vector<
     particlePusher<UsedParticlePusher>,
     shape<UsedParticleShape>,
     interpolation<UsedField2Particle>,
     current<UsedParticleCurrentSolver>,
-    #if(PARAM_IONIZATION == 1)
-    ionizer<particles::ionization::BSI<PIC_Electrons> >,
-    #endif
     atomicNumbers<Hydrogen>
 > ParticleFlagsIons;
 
 /* define species: ions */
 typedef Particles<ParticleDescription<
-    bmpl::string<'i'>,
+    bmpl::string<'i','H'>,
     SuperCellSize,
     AttributeSeqIons,
     ParticleFlagsIons,
     typename MakeSeq<CommunicationId<PAR_IONS> >::type >
 > PIC_Ions;
 
+//##### carbon
+
+struct Carbon
+{
+    static const float_X numberOfProtons  = 6.0;
+    static const float_X numberOfNeutrons = 6.0;
+};
+typedef bmpl::vector<
+    particlePusher<UsedParticlePusher>,
+    shape<UsedParticleShape>,
+    interpolation<UsedField2Particle>,
+    current<UsedParticleCurrentSolver>,
+    atomicNumbers<Carbon>
+> ParticleFlagsIonsC;
+
+/* define species: ions */
+typedef Particles<ParticleDescription<
+    bmpl::string<'i','C'>,
+    SuperCellSize,
+    AttributeSeqIons,
+    ParticleFlagsIonsC,
+    typename MakeSeq<CommunicationId<PAR_ELECTRONS+2> >::type >
+> PIC_IonsC;
+
+
 /*########################### end species ####################################*/
 
 
@@ -175,7 +193,9 @@ PIC_Ions
 
 typedef typename MakeSeq<
 Species1,
-Species2
+Species2,
+PIC_IonsC,
+PIC_ElectronsC
 >::type VectorAllSpecies;
 
 

Diff for: examples/LaserWakefield/include/simulation_defines/param/speciesInitialization.param

@@ -32,19 +32,14 @@ namespace picongpu
  * the functors are called in order (from first to last functor)
  */
 typedef mpl::vector<
-#if (PARAM_IONIZATION == 0)
-
-    particles::CreateGas<gasProfiles::Gaussian,particles::startPosition::Random,PIC_Electrons>
-    #if (ENABLE_IONS == 1)
-        ,particles::CloneSpecies<PIC_Electrons,PIC_Ions>
-    #endif
-
-#else
-
-    particles::CreateGas<gasProfiles::Gaussian,particles::startPosition::Random,PIC_Ions>,
-    particles::Manipulate<particles::manipulators::SetBoundElectrons,PIC_Ions>
-
-#endif
+    particles::CreateGas<gasProfiles::SphereFlanks,particles::startPosition::Random,PIC_IonsC>,
+    particles::CloneSpecies<PIC_IonsC,PIC_Ions>,
+    particles::CloneSpecies<PIC_IonsC,PIC_Ions>,
+    particles::CloneSpecies<PIC_IonsC,PIC_Ions>,
+    particles::Manipulate<particles::manipulators::RandomPosImpl<>,PIC_Ions>,
+    particles::CloneSpecies<PIC_Ions,PIC_Electrons>,
+    particles::CloneSpecies<PIC_IonsC,PIC_ElectronsC>,
+    particles::Manipulate<particles::manipulators::SixTimesWeighting,PIC_ElectronsC>
 > InitPipeline;
 
 } //namespace picongpu