Merge pull request ComputationalRadiationPhysics#852 from BeyondEspresso/topic-update-TWTS-polarization-BeyondEspresso

[BG fields] Add in-plane polarization to TWTS laser

Author: ax3l

examples/Bunch/include/simulation_defines/param/fieldBackground.param

@@ -35,20 +35,26 @@ namespace picongpu
     class FieldBackgroundE
     {
     public:
-    
+
         /* Add this additional field for pushing particles */
         static const bool InfluenceParticlePusher = PARAM_INCLUDE_FIELDBACKGROUND;
 
         /* We use this to calculate your SI input back to our unit system */
         PMACC_ALIGN(unitField, const float3_64);
 
         /* TWTS E-fields need to be initialized on host,
-         *  so they can look up global grid dimensions. */
-        PMACC_ALIGN(twtsFieldE, const templates::twts::EField);
-        
+         * so they can look up global grid dimensions.
+         *
+         * Note: No PMACC_ALIGN(...) used, since this *additional* memory alignment would require
+         *       roughly double the number of registers in the corresponding kernel on the device.
+         */
+        const templates::twts::EField twtsFieldE;
+
         /* Constructor is host-only, because of subGrid and halfSimSize initialization */
-        HINLINE FieldBackgroundE( const float3_64 unitField ) : unitField(unitField),
-        
+        HINLINE FieldBackgroundE( const float3_64 unitField ) :
+
+            unitField(unitField),
+
             twtsFieldE(
                 /* focus_y [m], the distance to the laser focus in y-direction */
                 30.0e-6,
@@ -78,10 +84,9 @@ namespace picongpu
         operator()( const DataSpace<simDim>& cellIdx,
                     const uint32_t currentStep ) const
         {
-            /* example 1: TWTS background pulse */
-            /** unit: meter */
+            /* unit: meter */
             const double WAVE_LENGTH_SI = 0.8e-6;
-            
+
             /** UNITCONV */
             /* const double UNITCONV_Intens_to_A0 = SI::ELECTRON_CHARGE_SI
                 * SI::ELECTRON_CHARGE_SI * 2.0 * WAVE_LENGTH_SI * WAVE_LENGTH_SI / (4.0 * PI * PI
@@ -92,24 +97,27 @@ namespace picongpu
                 * SI::ELECTRON_MASS_SI * SI::SPEED_OF_LIGHT_SI
                 * SI::SPEED_OF_LIGHT_SI / SI::ELECTRON_CHARGE_SI;
 
-            /** unit: W / m^2 */
-            /* const double _PEAK_INTENSITY_SI = 3.4e19 * 1.0e4; */
-            /** unit: none */
-            /* const double _A0  = _PEAK_INTENSITY_SI * UNITCONV_Intens_to_A0; */
+            /* unit: W / m^2
+             * const double _PEAK_INTENSITY_SI = 3.4e19 * 1.0e4;
+             * unit: none
+             * const double _A0  = _PEAK_INTENSITY_SI * UNITCONV_Intens_to_A0;
+             */
 
-            /** unit: none */
+            /* unit: none */
             const double _A0  = 1.0;
 
-            /** unit: Volt /meter */
-            /*\todo #738 implement math::vector, native type operations */
+            /* unit: Volt /meter
+             *\todo #738 implement math::vector, native type operations
+             */
             const float3_64 invUnitField = float3_64(1.0 / unitField[0],
                                                      1.0 / unitField[1],
                                                      1.0 / unitField[2] );
-            
+
             /* laser amplitude in picongpu units [ unit: (Volt /meter) / unitField-factor ]
-             * Note: the laser amplitude is included in all field components 
+             * Note: the laser amplitude is included in all field components
              * polarization and other properties are established by the peak amplitude
-             * normalized twtsFieldE(...) */
+             * normalized twtsFieldE(...)
+             */
             const float3_X amplitude = precisionCast<float_X>(
                                 float_64(_A0 * UNITCONV_A0_to_Amplitude_SI) * invUnitField );
 
@@ -125,14 +133,20 @@ namespace picongpu
         static const bool InfluenceParticlePusher = PARAM_INCLUDE_FIELDBACKGROUND;
 
         /* TWTS B-fields need to be initialized on host,
-         * so they can look up global grid dimensions. */
-        PMACC_ALIGN(twtsFieldB, const templates::twts::BField);
-        
+         * so they can look up global grid dimensions.
+         *
+         * Note: No PMACC_ALIGN(...) used, since this *additional* memory alignment would require
+         *       roughly double the number of registers in the corresponding kernel on the device.
+         */
+        const templates::twts::BField twtsFieldB;
+
         /* We use this to calculate your SI input back to our unit system */
         PMACC_ALIGN(unitField, const float3_64);
-        
-        HINLINE FieldBackgroundB( const float3_64 unitField ) : unitField(unitField),
-        
+
+        HINLINE FieldBackgroundB( const float3_64 unitField ) :
+
+            unitField(unitField),
+
             twtsFieldB(
                 /* focus_y [m], the distance to the laser focus in y-direction */
                 30.0e-6,
@@ -162,9 +176,9 @@ namespace picongpu
         operator()( const DataSpace<simDim>& cellIdx,
                     const uint32_t currentStep ) const
         {
-            /** unit: meter */
+            /* unit: meter */
             const double WAVE_LENGTH_SI = 0.8e-6;
-            
+
             /** UNITCONV */
             /* const double UNITCONV_Intens_to_A0 = SI::ELECTRON_CHARGE_SI
                 * SI::ELECTRON_CHARGE_SI * 2.0 * WAVE_LENGTH_SI * WAVE_LENGTH_SI / (4.0 * PI * PI
@@ -175,10 +189,11 @@ namespace picongpu
                 * SI::ELECTRON_MASS_SI * SI::SPEED_OF_LIGHT_SI
                 * SI::SPEED_OF_LIGHT_SI / SI::ELECTRON_CHARGE_SI;
 
-            /** unit: W / m^2 */
-            /* const double _PEAK_INTENSITY_SI = 3.4e19 * 1.0e4; */
-            /** unit: none */
-            /* const double _A0  = _PEAK_INTENSITY_SI * UNITCONV_Intens_to_A0; */
+            /* unit: W / m^2
+             * const double _PEAK_INTENSITY_SI = 3.4e19 * 1.0e4;
+             * unit: none
+             * const double _A0  = _PEAK_INTENSITY_SI * UNITCONV_Intens_to_A0;
+             */
 
             /** unit: none */
             const double _A0  = 1.0;
@@ -187,14 +202,15 @@ namespace picongpu
             const float3_64 invUnitField = float3_64( 1.0 / unitField[0],
                                                       1.0 / unitField[1],
                                                       1.0 / unitField[2] );
-            
+
             /* laser amplitude in picongpu units [ unit: (Volt/meter) / unitField-factor ]
-             * Note: the laser amplitude is included in all field components 
+             * Note: the laser amplitude is included in all field components
              * polarization and other properties are established by the peak amplitude
-             * normalized twtsFieldB(...) */
+             * normalized twtsFieldB(...)
+             */
             const float3_X amplitude = precisionCast<float_X>(
                     float_64(_A0 * UNITCONV_A0_to_Amplitude_SI) * invUnitField );
-            
+
             /* Note: twtsFieldB(...) is normalized, such that peak amplitude equals unity. */
             return amplitude * twtsFieldB( cellIdx, currentStep );
         }
@@ -208,7 +224,7 @@ namespace picongpu
 
         /* We use this to calculate your SI input back to our unit system */
         PMACC_ALIGN(unitField, const float3_64);
-        
+
         HDINLINE FieldBackgroundJ( const float3_64 unitField ) : unitField(unitField)
         {}
 

src/picongpu/include/fields/background/templates/TWTS/BField.hpp

@@ -29,23 +29,39 @@
 
 namespace picongpu
 {
-/** Load pre-defined background field */
+/* Load pre-defined background field */
 namespace templates
 {
-/** Traveling-wave Thomson scattering laser pulse */
+/* Traveling-wave Thomson scattering laser pulse */
 namespace twts
 {
 
 class BField
 {
 public:
     typedef float_X float_T;
-    
+
+    enum PolarizationType
+    {
+        /* The linear polarization of the TWTS laser is defined
+         * relative to the plane of the pulse front tilt (reference plane).
+         *
+         * Polarisation is normal to the reference plane.
+         * Use Ex-fields (and corresponding B-fields) in TWTS laser internal coordinate system.
+         */
+        LINEAR_X = 1u,
+        /* Polarization lies within the reference plane.
+         * Use Ey-fields (and corresponding B-fields) in TWTS laser internal coordinate system.
+         */
+        LINEAR_YZ = 2u,
+    };
+
     /* Center of simulation volume in number of cells */
     PMACC_ALIGN(halfSimSize,DataSpace<simDim>);
     /* y-position of TWTS coordinate origin inside the simulation coordinates [meter]
-       The other origin coordinates (x and z) default to globally centered values
-       with respect to the simulation volume. */
+     * The other origin coordinates (x and z) default to globally centered values
+     * with respect to the simulation volume.
+     */
     const PMACC_ALIGN(focus_y_SI,float_64);
     /* Laser wavelength [meter] */
     const PMACC_ALIGN(wavelength_SI,float_64);
@@ -69,14 +85,17 @@ class BField
     /* TWTS laser time delay */
     PMACC_ALIGN(tdelay,float_64);
     /* Should the TWTS laser time delay be chosen automatically, such that
-    the laser gradually enters the simulation volume? [Default: TRUE] */
+     * the laser gradually enters the simulation volume? [Default: TRUE]
+     */
     const PMACC_ALIGN(auto_tdelay,bool);
-    
+    /* Polarization of TWTS laser */
+    const PMACC_ALIGN(pol,PolarizationType);
+
     /** Magnetic field of the TWTS laser
      *
      * \param focus_y_SI the distance to the laser focus in y-direction [m]
      * \param wavelength_SI central wavelength [m]
-     * \param pulselength_SI sigma of std. gauss for intensity (E^2), 
+     * \param pulselength_SI sigma of std. gauss for intensity (E^2),
      *  pulselength_SI = FWHM_of_Intensity / 2.35482 [seconds (sigma)]
      * \param w_x beam waist: distance from the axis where the pulse electric field
      *  decreases to its 1/e^2-th part at the focus position of the laser [m]
@@ -88,6 +107,8 @@ class BField
      * \param tdelay_user manual time delay if auto_tdelay is false
      * \param auto_tdelay calculate the time delay such that the TWTS pulse is not
      *  inside the simulation volume at simulation start timestep = 0 [default = true]
+     * \param pol determines the TWTS laser polarization, which is either normal or parallel
+     *  to the laser pulse front tilt plane [ default= LINEAR_X , LINEAR_YZ ]
      */
     HINLINE
     BField( const float_64 focus_y_SI,
@@ -98,9 +119,10 @@ class BField
             const float_X phi               = 90.*(PI / 180.),
             const float_X beta_0            = 1.0,
             const float_64 tdelay_user_SI   = 0.0,
-            const bool auto_tdelay          = true );
-    
-    
+            const bool auto_tdelay          = true,
+            const PolarizationType pol      = LINEAR_X );
+
+
     /** Specify your background field B(r,t) here
      *
      * \param cellIdx The total cell id counted from the start at t=0
@@ -109,21 +131,41 @@ class BField
     operator()( const DataSpace<simDim>& cellIdx,
                 const uint32_t currentStep ) const;
 
-    /** Calculate the By(r,t) field here
+    /** Calculate the By(r,t) field, when electric field vector (Ex,0,0)
+     *  is normal to the pulse-front-tilt plane (y,z)
      *
      * \param pos Spatial position of the target field.
      * \param time Absolute time (SI, including all offsets and transformations)
      *  for calculating the field */
     HDINLINE float_T
     calcTWTSBy( const float3_64& pos, const float_64 time ) const;
 
-    /** Calculate the Bz(r,t) field here
+    /** Calculate the Bz(r,t) field, when electric field vector (Ex,0,0)
+     *  is normal to the pulse-front-tilt plane (y,z)
+     *
+     * \param pos Spatial position of the target field.
+     * \param time Absolute time (SI, including all offsets and transformations)
+     *  for calculating the field */
+    HDINLINE float_T
+    calcTWTSBz_Ex( const float3_64& pos, const float_64 time ) const;
+
+    /** Calculate the By(r,t) field, when electric field vector (0,Ey,0)
+     *  lies within the pulse-front-tilt plane (y,z)
+     *
+     * \param pos Spatial position of the target field.
+     * \param time Absolute time (SI, including all offsets and transformations)
+     *  for calculating the field */
+    HDINLINE float_T
+    calcTWTSBx( const float3_64& pos, const float_64 time ) const;
+
+    /** Calculate the Bz(r,t) field here (electric field vector (0,Ey,0)
+     *  lies within the pulse-front-tilt plane (y,z)
      *
      * \param pos Spatial position of the target field.
      * \param time Absolute time (SI, including all offsets and transformations)
      *  for calculating the field */
     HDINLINE float_T
-    calcTWTSBz( const float3_64& pos, const float_64 time ) const;
+    calcTWTSBz_Ey( const float3_64& pos, const float_64 time ) const;
 
     /** Calculate the B-field vector of the TWTS laser in SI units.
      * \tparam T_dim Specializes for the simulation dimension
@@ -134,7 +176,17 @@ class BField
     getTWTSBfield_Normalized(
             const PMacc::math::Vector<floatD_64,detail::numComponents>& eFieldPositions_SI,
             const float_64 time) const;
-    
+
+    /** Calculate the B-field vector of the "in-plane" polarized TWTS laser in SI units.
+     * \tparam T_dim Specializes for the simulation dimension
+     * \param cellIdx The total cell id counted from the start at timestep 0
+     * \return B-field vector of the rotated TWTS field in SI units */
+    template<unsigned T_dim>
+    HDINLINE float3_X
+    getTWTSBfield_Normalized_Ey(
+            const PMacc::math::Vector<floatD_64,detail::numComponents>& eFieldPositions_SI,
+            const float_64 time) const;
+
 };
 
 } /* namespace twts */