Pulsar remove managed memory allocation

Examples/Physics_applications/pulsar/inputs.corotating.3d.PM

@@ -13,16 +13,18 @@
 #################################
 ####### GENERAL PARAMETERS ######
 #################################
-max_step = 50000
+max_step = 1
 amr.n_cell = 128 128 128
 amr.max_grid_size = 128
 amr.blocking_factor = 128
 amr.max_level = 0
 geometry.coord_sys   = 0                  # 0: Cartesian
-geometry.is_periodic = 0 0 0     # Is periodic?
-geometry.prob_lo     = 0.0    0.0    0.0
+geometry.prob_lo     = 0 0 0
 geometry.prob_hi     = 180000 180000 180000
 
+boundary.field_lo = pml pml pml
+boundary.field_hi = pml pml pml
+
 #################################
 ############ NUMERICS ###########
 #################################
@@ -31,28 +33,25 @@ warpx.verbose = 1
 warpx.do_dive_cleaning = 0
 warpx.use_filter = 1
 warpx.cfl = .95
-warpx.do_pml = 1
 my_constants.pi    = 3.141592653589793
 my_constants.dens  = 5.544e6
 my_constants.xc  = 90000
 my_constants.yc  = 90000
 my_constants.zc  = 90000
-my_constants.r_star  = 12000
+my_constants.r_star  = 12032
 my_constants.omega = 6245.676
 my_constants.c = 299792458.
 my_constants.B_star = 8.0323e-6                  # magnetic field of NS (T)
-my_constants.dR = 1400
+my_constants.dR = 1406.25
 my_constants.to = 2.e-4
-interpolation.nox = 3
-interpolation.noy = 3
-interpolation.noz = 3
+algo.particle_shape = 3
 
 
 #################################
 ############ PLASMA #############
 #################################
-particles.nspecies = 2
-particles.species_names = plasma_e plasma_p
+#particles.nspecies = 2
+#particles.species_names = plasma_e plasma_p
 
 plasma_e.charge = -q_e
 plasma_e.mass = m_e
@@ -61,7 +60,6 @@ plasma_e.profile = parse_density_function
 plasma_e.density_function(x,y,z) = "( ((( (z-zc)*(z-zc) + (y-yc)*(y-yc) + (x-xc)*(x-xc) )^(0.5))<=(r_star)) * ((( (z-zc)*(z-zc) + (y-yc)*(y-yc) + (x-xc)*(x-xc) )^(0.5))>=(r_star-dR)) )*dens"
 plasma_e.num_particles_per_cell_each_dim = 4 4 4
 plasma_e.momentum_distribution_type = parse_momentum_function
-
 ## Inject stationary pairs
 plasma_e.momentum_function_ux(x,y,z) = "0.0"
 plasma_e.momentum_function_uy(x,y,z) = "0.0"
@@ -77,7 +75,6 @@ plasma_p.profile = parse_density_function
 plasma_p.density_function(x,y,z) = "( ((( (z-zc)*(z-zc) + (y-yc)*(y-yc) + (x-xc)*(x-xc) )^(0.5))<=(r_star)) * ((( (z-zc)*(z-zc) + (y-yc)*(y-yc) + (x-xc)*(x-xc) )^(0.5))>=(r_star-dR)) )*dens"
 plasma_p.num_particles_per_cell_each_dim = 4 4 4
 plasma_p.momentum_distribution_type = parse_momentum_function
-
 ## Inject stationary pairs
 plasma_p.momentum_function_ux(x,y,z) = "0.0"
 plasma_p.momentum_function_uy(x,y,z) = "0.0"
@@ -88,8 +85,8 @@ plasma_p.density_min = 1E0
 
 diagnostics.diags_names = plt chk
 plt.diag_type = Full
-plt.intervals = 20
-plt.fields_to_plot = Ex Ey Ez Bx By Bz jx jy jz part_per_cell rho divE
+plt.intervals = 1
+plt.fields_to_plot = Ex Ey Ez Bx By Bz jx jy jz part_per_cell rho divE Er Etheta Ephi Br Btheta Bphi EdotB Sx Sy Sz
 
 chk.format = checkpoint
 chk.diag_type = Full
@@ -106,11 +103,11 @@ pulsar.ramp_omega_time = 0.0             # time over which to ramp up NS angular
 pulsar.center_star = 90000 90000 90000
 pulsar.R_star = 12032                  # radius of NS (m)
 pulsar.B_star = 8.0323e-6                  # magnetic field of NS (T)
-pulsar.dR = 1400                       # thickness on the surface of the pulsar
+pulsar.dR = 1406.25                       # thickness on the surface of the pulsar
                                        # consistency requires dR = my_constants.dR
 pulsar.verbose = 0                     # [0/1]: turn on verbosity for debugging print statements
-pulsar.EB_external = 1                 # [0/1]: to apply external E and B
-pulsar.E_external_monopole = 0                  # [0/1]
+pulsar.EB_external = 0                 # [0/1]: to apply external E and B
+pulsar.E_external_monopole = 1                  # [0/1]
 pulsar.EB_corotating_maxradius = 12032 # The radius where E-field changes from corotating
                                        # inside the star to quadrapole outside.
                                        # Default is R_star. (r<=cor_radius) i.e. the includes
@@ -126,7 +123,7 @@ pulsar.max_particle_absorption_radius = 12032   # Maximum radius within which pa
                                                 # deleted every timestep.
                                                 # Default is Rstar
 pulsar.rhoGJ_scale = 1e0              # scaling down of rho_GJ
-pulsar.damp_EB_internal = 1             # Damp internal electric field
+pulsar.damp_EB_internal = 0             # Damp internal electric field
 pulsar.damp_EB_radius = 12032         # The radius within which EB should be damped.
                                       # default is r_star (damping will include this radius)
 pulsar.damping_scale = 1000.0            # Damping scale factor for internal electric field
@@ -137,4 +134,20 @@ pulsar.max_nodepos_radius = 0.        # radius within which particle deposition
 pulsar.turnoff_plasmaEB_gather = 0    # (0/1) 0 is default. always gather
                                       # 1 for no gather of self-consistent fields
 pulsar.max_nogather_radius = 0.       # radius within which self-consistent field gather
-                                      # is turned off
+pulsar.init_dipoleBfield = 1          # default is 1
+pulsar.init_corotatingEfield = 1      # default
+pulsar.corotatingE_maxradius = 16250.75  # default is Rstar
+pulsar.enforceDipoleB_maxradius = 12032 # default is Rstar
+pulsar.enforceCorotatingE = 1        # default 1
+pulsar.enforceDipoleB = 1        # default 1
+pulsar.singleParticleTest = 0    # default 0, if 1, then a particle pair will be introduced.
+                                 # Single particle position/momentum will need to be specified
+pulsar.DampBDipoleInRing = 0
+pulsar.Bdamping_scale = 10000
+pulsar.LimitDipoleBInit = 0   # if 1 then Bdipole is iniialization only in a smaller region
+pulsar.DipoleB_init_maxradius = 12032 # read if LimitBdipole is 1
+pulsar.use_theoreticalEB = 0
+pulsar.theory_max_rstar = 100000
+pulsar.EnforceTheoreticalEBInGrid = 1
+pulsar.AddExternalMonopoleOnly = 1
+

Examples/Physics_applications/pulsar/inputs.corotating.singleParticle.3d.PM

@@ -0,0 +1,151 @@
+# Description:
+#
+# This inputs file sets up a NS with these properties:
+# - E = -(omega*R)[cross]B inside the NS
+# - external E and B outside the star
+#
+# See the "Pulsar Setup" section at the end for the options
+#
+# This initializes the electrons and positrons with a corotating momentum function.
+# Based on the pulsar_type = "active" or "dead", particles are injected continuously or
+#                             until rho_GJ is reached
+
+#################################
+####### GENERAL PARAMETERS ######
+#################################
+max_step = 20000
+amr.n_cell = 256 256 256
+amr.max_grid_size = 128
+amr.blocking_factor = 128
+amr.max_level = 0
+geometry.coord_sys   = 0                  # 0: Cartesian
+geometry.is_periodic = 0 0 0     # Is periodic?
+geometry.prob_lo     = 0.0    0.0    0.0
+geometry.prob_hi     = 180000 180000 180000
+
+#################################
+############ NUMERICS ###########
+#################################
+#algo.maxwell_fdtd_solver = yee
+warpx.verbose = 1
+warpx.do_dive_cleaning = 0
+warpx.use_filter = 1
+warpx.cfl = .95
+warpx.do_pml = 1
+my_constants.pi    = 3.141592653589793
+my_constants.dens  = 5.544e6
+my_constants.xc  = 90000
+my_constants.yc  = 90000
+my_constants.zc  = 90000
+my_constants.r_star  = 12000
+my_constants.omega = 6245.676
+my_constants.c = 299792458.
+my_constants.B_star = 8.0323e-6                  # magnetic field of NS (T)
+my_constants.dR = 1400
+my_constants.to = 2.e-4
+interpolation.nox = 3
+interpolation.noy = 3
+interpolation.noz = 3
+
+
+#################################
+############ PLASMA #############
+#################################
+particles.nspecies = 2
+particles.species_names = plasma_e plasma_p
+
+plasma_e.charge = -q_e
+plasma_e.mass = m_e
+plasma_e.injection_style = "singleparticle"
+plasma_e.single_particle_pos = 90000. 103000. 90000.
+plasma_e.single_particle_vel = 0. 0. 0.
+plasma_e.single_particle_weight = 1.
+#plasma_e.profile = parse_density_function
+#plasma_e.density_function(x,y,z) = "( ((( (z-zc)*(z-zc) + (y-yc)*(y-yc) + (x-xc)*(x-xc) )^(0.5))<=(r_star)) * ((( (z-zc)*(z-zc) + (y-yc)*(y-yc) + (x-xc)*(x-xc) )^(0.5))>=(r_star-dR)) )*dens"
+#plasma_e.num_particles_per_cell_each_dim = 4 4 4
+#plasma_e.momentum_distribution_type = parse_momentum_function
+#
+### Inject stationary pairs
+#plasma_e.momentum_function_ux(x,y,z) = "0.0"
+#plasma_e.momentum_function_uy(x,y,z) = "0.0"
+### uz is always 0 for the injection methods above
+#plasma_e.momentum_function_uz(x,y,z) = "0.0"
+#plasma_e.do_continuous_injection = 0
+#plasma_e.density_min = 1E0
+
+plasma_p.charge = q_e
+plasma_p.mass = m_e
+plasma_p.injection_style = "singleparticle"
+plasma_p.single_particle_pos = 90000. 103000. 90000.
+plasma_p.single_particle_vel = 0. 0. 0.
+plasma_p.single_particle_weight = 1.
+#plasma_p.profile = parse_density_function
+#plasma_p.density_function(x,y,z) = "( ((( (z-zc)*(z-zc) + (y-yc)*(y-yc) + (x-xc)*(x-xc) )^(0.5))<=(r_star)) * ((( (z-zc)*(z-zc) + (y-yc)*(y-yc) + (x-xc)*(x-xc) )^(0.5))>=(r_star-dR)) )*dens"
+#plasma_p.num_particles_per_cell_each_dim = 4 4 4
+#plasma_p.momentum_distribution_type = parse_momentum_function
+#
+### Inject stationary pairs
+#plasma_p.momentum_function_ux(x,y,z) = "0.0"
+#plasma_p.momentum_function_uy(x,y,z) = "0.0"
+### uz is always 0 for the injection methods above
+#plasma_p.momentum_function_uz(x,y,z) = "0.0"
+#plasma_p.do_continuous_injection = 0
+#plasma_p.density_min = 1E0
+
+diagnostics.diags_names = plt chk
+plt.diag_type = Full
+plt.intervals = 50
+plt.fields_to_plot = Ex Ey Ez Bx By Bz jx jy jz part_per_cell rho divE rho_plasma_e rho_plasma_p
+
+chk.format = checkpoint
+chk.diag_type = Full
+chk.intervals = 1000
+
+#################################
+######### PULSAR SETUP ##########
+#################################
+pulsar.pulsarType = "dead"             # [dead/active]: sets particle injection type
+pulsar.omega_star = 6245.676               # angular velocity of NS (rad/s)
+pulsar.ramp_omega_time = 0.0             # time over which to ramp up NS angular velocity (s)
+                                           # if ramp_omega_time < 0, omega = omega_star for any t
+                                           # consistency requires ramp_omega_time = my_constants.to
+pulsar.center_star = 90000 90000 90000
+pulsar.R_star = 12032                  # radius of NS (m)
+pulsar.B_star = 8.0323e-6                  # magnetic field of NS (T)
+pulsar.dR = 1400                       # thickness on the surface of the pulsar
+                                       # consistency requires dR = my_constants.dR
+pulsar.verbose = 0                     # [0/1]: turn on verbosity for debugging print statements
+pulsar.EB_external = 0                 # [0/1]: to apply external E and B
+pulsar.E_external_monopole = 0                  # [0/1]
+pulsar.EB_corotating_maxradius = 12032 # The radius where E-field changes from corotating
+                                       # inside the star to quadrapole outside.
+                                       # Default is R_star. (r<=cor_radius) i.e. the includes
+                                       # the r specified in the input
+pulsar.max_ndens = 5.54e6              # max ndens == ndens used in initializing density
+pulsar.Ninj_fraction = 0.2            # fraction of sigma injected
+pulsar.ModifyParticleWeight = 1       # (0/1) the fractional injection is achieved
+                                      # by modifying particle weight if 1
+                                      # by modifying num_ppc if 0
+pulsar.particle_inj_rmin = 10632      # default is Rstar-dR (consistent with density function above)
+pulsar.particle_inj_rmax = 12032      # default is Rstar (consistent with density function)
+pulsar.max_particle_absorption_radius = 12032   # Maximum radius within which particles are
+                                                # deleted every timestep.
+                                                # Default is Rstar
+pulsar.rhoGJ_scale = 1e0              # scaling down of rho_GJ
+pulsar.damp_EB_internal = 0             # Damp internal electric field
+pulsar.damp_EB_radius = 12032         # The radius within which EB should be damped.
+                                      # default is r_star (damping will include this radius)
+pulsar.damping_scale = 1000.0            # Damping scale factor for internal electric field
+pulsar.turnoffdeposition = 0         # (0/1) 0 for depositing (default)
+                                      #       1 for no deposition
+pulsar.max_nodepos_radius = 0.        # radius within which particle deposition for j/rho
+                                      # is turned off. (r<=max_radius)
+pulsar.turnoff_plasmaEB_gather = 0    # (0/1) 0 is default. always gather
+                                      # 1 for no gather of self-consistent fields
+pulsar.max_nogather_radius = 0.       # radius within which self-consistent field gather
+pulsar.init_dipoleBfield = 1          # default is 1
+pulsar.init_corotatingEfield = 1      # default
+pulsar.corotatingE_maxradius = 12032  # default is Rstar
+pulsar.enforceDipoleB_maxradius = 12032 # default is Rstar
+pulsar.enforceCorotatingE = 1        # default 1
+pulsar.enforceDipoleB = 1        # default 1

Source/Diagnostics/ComputeDiagFunctors/CMakeLists.txt

@@ -7,4 +7,7 @@ target_sources(WarpX
     PartPerCellFunctor.cpp
     PartPerGridFunctor.cpp
     BackTransformFunctor.cpp
+    SphericalComponentFunctor.cpp
+    EdotBFunctor.cpp
+    PoyntingVectorFunctor.cpp
 )

Source/Diagnostics/ComputeDiagFunctors/EdotBFunctor.H

@@ -0,0 +1,31 @@
+#ifndef WARPX_EDOTBFUNCTOR_H_
+#define WARPX_EDOTBFUNCTOR_H_
+
+#include "ComputeDiagFunctor.H"
+
+class
+EdotBFunctor final : public ComputeDiagFunctor
+{
+public:
+    EdotBFunctor ( const amrex::MultiFab * const Ex_src,
+                   const amrex::MultiFab * const Ey_src,
+                   const amrex::MultiFab * const EZ_src,
+                   const amrex::MultiFab * const Bx_src,
+                   const amrex::MultiFab * const By_src,
+                   const amrex::MultiFab * const Bz_src,
+                   const int lev,
+                   const amrex::IntVect crse_ratio,
+                   int ncomp = 1);
+    virtual void operator() (amrex::MultiFab& mf_dst, int dcomp, const int /*i_buffer=0*/) const override;
+    void ComputeEdotB(amrex::MultiFab& mf_dst, int dcomp) const;
+private:
+    amrex::MultiFab const * const m_Ex_src = nullptr;
+    amrex::MultiFab const * const m_Ey_src = nullptr;
+    amrex::MultiFab const * const m_Ez_src = nullptr;
+    amrex::MultiFab const * const m_Bx_src = nullptr;
+    amrex::MultiFab const * const m_By_src = nullptr;
+    amrex::MultiFab const * const m_Bz_src = nullptr;
+    int m_lev;
+};
+
+#endif