aboutsummaryrefslogtreecommitdiff
path: root/Source
diff options
context:
space:
mode:
Diffstat (limited to 'Source')
-rw-r--r--Source/Evolve/WarpXEvolveEM.cpp12
-rw-r--r--Source/FieldSolver/SpectralSolver/SpectralAlgorithms/PsatdAlgorithm.H20
-rw-r--r--Source/FieldSolver/SpectralSolver/SpectralAlgorithms/PsatdAlgorithm.cpp169
-rw-r--r--Source/FieldSolver/SpectralSolver/SpectralFieldData.H1
-rw-r--r--Source/FieldSolver/SpectralSolver/SpectralFieldData.cpp63
-rw-r--r--Source/FieldSolver/SpectralSolver/SpectralKSpace.cpp6
-rw-r--r--Source/FieldSolver/WarpXFFT.cpp1
-rw-r--r--Source/Initialization/WarpXInitData.cpp13
-rw-r--r--Source/Parallelization/WarpXComm.cpp7
-rw-r--r--Source/Particles/PhysicalParticleContainer.cpp2
-rw-r--r--Source/Particles/WarpXParticleContainer.cpp4
-rw-r--r--Source/WarpX.cpp2
-rw-r--r--Source/main.cpp3
13 files changed, 242 insertions, 61 deletions
diff --git a/Source/Evolve/WarpXEvolveEM.cpp b/Source/Evolve/WarpXEvolveEM.cpp
index e98561be1..02c21529b 100644
--- a/Source/Evolve/WarpXEvolveEM.cpp
+++ b/Source/Evolve/WarpXEvolveEM.cpp
@@ -69,24 +69,30 @@ WarpX::EvolveEM (int numsteps)
// At the beginning, we have B^{n} and E^{n}.
// Particles have p^{n} and x^{n}.
// is_synchronized is true.
+ amrex::Print() << " in evolve before fill boundary \n";
if (is_synchronized) {
FillBoundaryE();
FillBoundaryB();
UpdateAuxilaryData();
// on first step, push p by -0.5*dt
+ amrex::Print() << " in evolve before pushP \n";
for (int lev = 0; lev <= finest_level; ++lev) {
mypc->PushP(lev, -0.5*dt[lev],
*Efield_aux[lev][0],*Efield_aux[lev][1],*Efield_aux[lev][2],
*Bfield_aux[lev][0],*Bfield_aux[lev][1],*Bfield_aux[lev][2]);
}
is_synchronized = false;
+ amrex::Print() << " in evolve after pushP \n";
+
} else {
// Beyond one step, we have E^{n} and B^{n}.
// Particles have p^{n-1/2} and x^{n}.
FillBoundaryE();
FillBoundaryB();
UpdateAuxilaryData();
+
}
+ amrex::Print() << " in evolve after fill boundary \n";
if (do_subcycling == 0 || finest_level == 0) {
OneStep_nosub(cur_time);
@@ -96,6 +102,7 @@ WarpX::EvolveEM (int numsteps)
amrex::Print() << "Error: do_subcycling = " << do_subcycling << std::endl;
amrex::Abort("Unsupported do_subcycling type");
}
+ amrex::Print() << " in evolve after onestep no sub \n";
if (num_mirrors>0){
applyMirrors(cur_time);
@@ -267,7 +274,10 @@ WarpX::OneStep_nosub (Real cur_time)
if (warpx_py_particlescraper) warpx_py_particlescraper();
if (warpx_py_beforedeposition) warpx_py_beforedeposition();
#endif
+ amrex::Print() << " before push \n";
PushParticlesandDepose(cur_time);
+ amrex::Print() << " after push \n";
+
#ifdef WARPX_USE_PY
if (warpx_py_afterdeposition) warpx_py_afterdeposition();
#endif
@@ -280,7 +290,9 @@ WarpX::OneStep_nosub (Real cur_time)
// (And update guard cells immediately afterwards)
#ifdef WARPX_USE_PSATD
PushPSATD(dt[0]);
+ amrex::Print() << " before fill bndry E \n";
FillBoundaryE();
+ amrex::Print() << " before fill bndry B \n";
FillBoundaryB();
#else
EvolveF(0.5*dt[0], DtType::FirstHalf);
diff --git a/Source/FieldSolver/SpectralSolver/SpectralAlgorithms/PsatdAlgorithm.H b/Source/FieldSolver/SpectralSolver/SpectralAlgorithms/PsatdAlgorithm.H
index 0487e5226..34743525e 100644
--- a/Source/FieldSolver/SpectralSolver/SpectralAlgorithms/PsatdAlgorithm.H
+++ b/Source/FieldSolver/SpectralSolver/SpectralAlgorithms/PsatdAlgorithm.H
@@ -1,6 +1,8 @@
#ifndef WARPX_PSATD_ALGORITHM_H_
#define WARPX_PSATD_ALGORITHM_H_
+#include <SpectralKSpace.H>
+#include <SpectralFieldData.H>
#include <SpectralBaseAlgorithm.H>
/* \brief Class that updates the field in spectral space
@@ -8,16 +10,26 @@
*/
class PsatdAlgorithm : public SpectralBaseAlgorithm
{
+ using SpectralCoefficients = amrex::FabArray< amrex::BaseFab <amrex::Real> >;
+
public:
PsatdAlgorithm(const SpectralKSpace& spectral_kspace,
const amrex::DistributionMapping& dm,
const int norder_x, const int norder_y,
- const int norder_z, const bool nodal,
- const amrex::Real dt);
- // Redefine update equation from base class
- virtual void pushSpectralFields(SpectralFieldData& f) const override final;
+ const int norder_z, const bool nodal, const amrex::Real dt);
+ PsatdAlgorithm() = default; // Default constructor
+ PsatdAlgorithm& operator=(PsatdAlgorithm&& algorithm) = default;
+ void pushSpectralFields(SpectralFieldData& f) const;
+ void InitializeCoefficience(const SpectralKSpace& spectral_kspace,
+ const amrex::DistributionMapping& dm,
+ const amrex::Real dt);
private:
+ // Modified finite-order vectors
+ KVectorComponent modified_kx_vec, modified_kz_vec;
+#if (AMREX_SPACEDIM==3)
+ KVectorComponent modified_ky_vec;
+#endif
SpectralCoefficients C_coef, S_ck_coef, X1_coef, X2_coef, X3_coef;
};
diff --git a/Source/FieldSolver/SpectralSolver/SpectralAlgorithms/PsatdAlgorithm.cpp b/Source/FieldSolver/SpectralSolver/SpectralAlgorithms/PsatdAlgorithm.cpp
index 37892d35a..8dd2a830f 100644
--- a/Source/FieldSolver/SpectralSolver/SpectralAlgorithms/PsatdAlgorithm.cpp
+++ b/Source/FieldSolver/SpectralSolver/SpectralAlgorithms/PsatdAlgorithm.cpp
@@ -22,58 +22,60 @@ PsatdAlgorithm::PsatdAlgorithm(const SpectralKSpace& spectral_kspace,
X2_coef = SpectralCoefficients(ba, dm, 1, 0);
X3_coef = SpectralCoefficients(ba, dm, 1, 0);
- // Fill them with the right values:
- // Loop over boxes and allocate the corresponding coefficients
- // for each box owned by the local MPI proc
- for (MFIter mfi(ba, dm); mfi.isValid(); ++mfi){
-
- const Box& bx = ba[mfi];
-
- // Extract pointers for the k vectors
- const Real* modified_kx = modified_kx_vec[mfi].dataPtr();
-#if (AMREX_SPACEDIM==3)
- const Real* modified_ky = modified_ky_vec[mfi].dataPtr();
-#endif
- const Real* modified_kz = modified_kz_vec[mfi].dataPtr();
- // Extract arrays for the coefficients
- Array4<Real> C = C_coef[mfi].array();
- Array4<Real> S_ck = S_ck_coef[mfi].array();
- Array4<Real> X1 = X1_coef[mfi].array();
- Array4<Real> X2 = X2_coef[mfi].array();
- Array4<Real> X3 = X3_coef[mfi].array();
-
- // Loop over indices within one box
- ParallelFor(bx,
- [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept
- {
- // Calculate norm of vector
- const Real k_norm = std::sqrt(
- std::pow(modified_kx[i], 2) +
-#if (AMREX_SPACEDIM==3)
- std::pow(modified_ky[j], 2) +
- std::pow(modified_kz[k], 2));
-#else
- std::pow(modified_kz[j], 2));
-#endif
-
- // Calculate coefficients
- constexpr Real c = PhysConst::c;
- constexpr Real ep0 = PhysConst::ep0;
- if (k_norm != 0){
- C(i,j,k) = std::cos(c*k_norm*dt);
- S_ck(i,j,k) = std::sin(c*k_norm*dt)/(c*k_norm);
- X1(i,j,k) = (1. - C(i,j,k))/(ep0 * c*c * k_norm*k_norm);
- X2(i,j,k) = (1. - S_ck(i,j,k)/dt)/(ep0 * k_norm*k_norm);
- X3(i,j,k) = (C(i,j,k) - S_ck(i,j,k)/dt)/(ep0 * k_norm*k_norm);
- } else { // Handle k_norm = 0, by using the analytical limit
- C(i,j,k) = 1.;
- S_ck(i,j,k) = dt;
- X1(i,j,k) = 0.5 * dt*dt / ep0;
- X2(i,j,k) = c*c * dt*dt / (6.*ep0);
- X3(i,j,k) = - c*c * dt*dt / (3.*ep0);
- }
- });
- }
+ InitializeCoefficience(spectral_kspace, dm, dt);
+// // Fill them with the right values:
+// // Loop over boxes and allocate the corresponding coefficients
+// // for each box owned by the local MPI proc
+// for (MFIter mfi(ba, dm); mfi.isValid(); ++mfi){
+//
+// //const Box& bx = ba[mfi];
+// const Box bx = ba[mfi];
+//
+// // Extract pointers for the k vectors
+// const Real* modified_kx = modified_kx_vec[mfi].dataPtr();
+//#if (AMREX_SPACEDIM==3)
+// const Real* modified_ky = modified_ky_vec[mfi].dataPtr();
+//#endif
+// const Real* modified_kz = modified_kz_vec[mfi].dataPtr();
+// // Extract arrays for the coefficients
+// Array4<Real> C = C_coef[mfi].array();
+// Array4<Real> S_ck = S_ck_coef[mfi].array();
+// Array4<Real> X1 = X1_coef[mfi].array();
+// Array4<Real> X2 = X2_coef[mfi].array();
+// Array4<Real> X3 = X3_coef[mfi].array();
+//
+// // Loop over indices within one box
+// ParallelFor(bx,
+// [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept
+// {
+// // Calculate norm of vector
+// const Real k_norm = std::sqrt(
+// std::pow(modified_kx[i], 2) +
+//#if (AMREX_SPACEDIM==3)
+// std::pow(modified_ky[j], 2) +
+// std::pow(modified_kz[k], 2));
+//#else
+// std::pow(modified_kz[j], 2));
+//#endif
+//
+// // Calculate coefficients
+// constexpr Real c = PhysConst::c;
+// constexpr Real ep0 = PhysConst::ep0;
+// if (k_norm != 0){
+// C(i,j,k) = std::cos(c*k_norm*dt);
+// S_ck(i,j,k) = std::sin(c*k_norm*dt)/(c*k_norm);
+// X1(i,j,k) = (1. - C(i,j,k))/(ep0 * c*c * k_norm*k_norm);
+// X2(i,j,k) = (1. - S_ck(i,j,k)/dt)/(ep0 * k_norm*k_norm);
+// X3(i,j,k) = (C(i,j,k) - S_ck(i,j,k)/dt)/(ep0 * k_norm*k_norm);
+// } else { // Handle k_norm = 0, by using the analytical limit
+// C(i,j,k) = 1.;
+// S_ck(i,j,k) = dt;
+// X1(i,j,k) = 0.5 * dt*dt / ep0;
+// X2(i,j,k) = c*c * dt*dt / (6.*ep0);
+// X3(i,j,k) = - c*c * dt*dt / (3.*ep0);
+// }
+// });
+// }
};
/* Advance the E and B field in spectral space (stored in `f`)
@@ -130,7 +132,7 @@ PsatdAlgorithm::pushSpectralFields(SpectralFieldData& f) const{
#endif
constexpr Real c2 = PhysConst::c*PhysConst::c;
constexpr Real inv_ep0 = 1./PhysConst::ep0;
- constexpr Complex I = Complex{0,1};
+ const Complex I = Complex{0,1};
const Real C = C_arr(i,j,k);
const Real S_ck = S_ck_arr(i,j,k);
const Real X1 = X1_arr(i,j,k);
@@ -160,3 +162,64 @@ PsatdAlgorithm::pushSpectralFields(SpectralFieldData& f) const{
});
}
};
+
+void PsatdAlgorithm::InitializeCoefficience(const SpectralKSpace& spectral_kspace,
+ const amrex::DistributionMapping& dm,
+ const amrex::Real dt)
+{
+ const BoxArray& ba = spectral_kspace.spectralspace_ba;
+ // Fill them with the right values:
+ // Loop over boxes and allocate the corresponding coefficients
+ // for each box owned by the local MPI proc
+ for (MFIter mfi(ba, dm); mfi.isValid(); ++mfi){
+
+ //const Box& bx = ba[mfi];
+ const Box bx = ba[mfi];
+
+ // Extract pointers for the k vectors
+ const Real* modified_kx = modified_kx_vec[mfi].dataPtr();
+#if (AMREX_SPACEDIM==3)
+ const Real* modified_ky = modified_ky_vec[mfi].dataPtr();
+#endif
+ const Real* modified_kz = modified_kz_vec[mfi].dataPtr();
+ // Extract arrays for the coefficients
+ Array4<Real> C = C_coef[mfi].array();
+ Array4<Real> S_ck = S_ck_coef[mfi].array();
+ Array4<Real> X1 = X1_coef[mfi].array();
+ Array4<Real> X2 = X2_coef[mfi].array();
+ Array4<Real> X3 = X3_coef[mfi].array();
+
+ // Loop over indices within one box
+ ParallelFor(bx,
+ [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept
+ {
+ // Calculate norm of vector
+ const Real k_norm = std::sqrt(
+ std::pow(modified_kx[i], 2) +
+#if (AMREX_SPACEDIM==3)
+ std::pow(modified_ky[j], 2) +
+ std::pow(modified_kz[k], 2));
+#else
+ std::pow(modified_kz[j], 2));
+#endif
+
+
+ // Calculate coefficients
+ constexpr Real c = PhysConst::c;
+ constexpr Real ep0 = PhysConst::ep0;
+ if (k_norm != 0){
+ C(i,j,k) = std::cos(c*k_norm*dt);
+ S_ck(i,j,k) = std::sin(c*k_norm*dt)/(c*k_norm);
+ X1(i,j,k) = (1. - C(i,j,k))/(ep0 * c*c * k_norm*k_norm);
+ X2(i,j,k) = (1. - S_ck(i,j,k)/dt)/(ep0 * k_norm*k_norm);
+ X3(i,j,k) = (C(i,j,k) - S_ck(i,j,k)/dt)/(ep0 * k_norm*k_norm);
+ } else { // Handle k_norm = 0, by using the analytical limit
+ C(i,j,k) = 1.;
+ S_ck(i,j,k) = dt;
+ X1(i,j,k) = 0.5 * dt*dt / ep0;
+ X2(i,j,k) = c*c * dt*dt / (6.*ep0);
+ X3(i,j,k) = - c*c * dt*dt / (3.*ep0);
+ }
+ });
+ }
+}
diff --git a/Source/FieldSolver/SpectralSolver/SpectralFieldData.H b/Source/FieldSolver/SpectralSolver/SpectralFieldData.H
index 8e58aa1d8..1d64817ef 100644
--- a/Source/FieldSolver/SpectralSolver/SpectralFieldData.H
+++ b/Source/FieldSolver/SpectralSolver/SpectralFieldData.H
@@ -26,6 +26,7 @@ class SpectralFieldData
// the local MPI rank)
#ifdef AMREX_USE_GPU
// Add cuFFT-specific code
+ using FFTplans = amrex::LayoutData<cufftHandle>;
#else
using FFTplans = amrex::LayoutData<fftw_plan>;
#endif
diff --git a/Source/FieldSolver/SpectralSolver/SpectralFieldData.cpp b/Source/FieldSolver/SpectralSolver/SpectralFieldData.cpp
index 02fa2015f..5998bdd2b 100644
--- a/Source/FieldSolver/SpectralSolver/SpectralFieldData.cpp
+++ b/Source/FieldSolver/SpectralSolver/SpectralFieldData.cpp
@@ -54,6 +54,25 @@ SpectralFieldData::SpectralFieldData( const BoxArray& realspace_ba,
IntVect fft_size = realspace_ba[mfi].length();
#ifdef AMREX_USE_GPU
// Add cuFFT-specific code
+ // Creating 3D plan for real to complex -- double precision
+ cudaDeviceSynchronize();
+ cufftResult result;
+ result = cufftPlan3d( &forward_plan[mfi], fft_size[2],
+ fft_size[1],fft_size[0], CUFFT_D2Z);
+ if ( result != CUFFT_SUCCESS ) {
+ amrex::Print() << " cufftplan3d forward failed! \n";
+ }
+ // Add 2D cuFFT-spacific code for D2Z
+ // Note that D2Z is inherently forward plan
+
+ result = cufftPlan3d( &backward_plan[mfi], fft_size[2],
+ fft_size[1], fft_size[0], CUFFT_Z2D);
+ // Add 2D cuFFT-specific code for Z2D
+ if ( result != CUFFT_SUCCESS ) {
+ amrex::Print() << " cufftplan3d backward failed! \n";
+ }
+ cudaDeviceSynchronize();
+
#else
// Create FFTW plans
forward_plan[mfi] =
@@ -87,6 +106,8 @@ SpectralFieldData::~SpectralFieldData()
for ( MFIter mfi(tmpRealField); mfi.isValid(); ++mfi ){
#ifdef AMREX_USE_GPU
// Add cuFFT-specific code
+ cufftDestroy( forward_plan[mfi] );
+ cufftDestroy( backward_plan[mfi] );
#else
// Destroy FFTW plans
fftw_destroy_plan( forward_plan[mfi] );
@@ -129,14 +150,30 @@ SpectralFieldData::ForwardTransform( const MultiFab& mf,
Array4<Real> tmp_arr = tmpRealField[mfi].array();
ParallelFor( realspace_bx,
[=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept {
- tmp_arr(i,j,k) = mf_arr(i,j,k,i_comp);
+ tmp_arr(i,j,k) = mf_arr(i,j,k,i_comp);
});
+//#ifdef AMREX_USE_GPU
+// cudaDeviceSynchronize();
+//#endif
}
// Perform Fourier transform from `tmpRealField` to `tmpSpectralField`
#ifdef AMREX_USE_GPU
// Add cuFFT-specific code ; make sure that this is done on the same
// GPU stream as the above copy
+ cudaDeviceSynchronize();
+ cufftResult result;
+ cudaStream_t stream = amrex::Cuda::Device::cudaStream();
+ amrex::Print() << " stream is " << stream << "\n";
+ cufftSetStream ( forward_plan[mfi], stream);
+ result = cufftExecD2Z( forward_plan[mfi],
+ tmpRealField[mfi].dataPtr(),
+ reinterpret_cast<cuDoubleComplex*>(
+ tmpSpectralField[mfi].dataPtr()) );
+ if ( result != CUFFT_SUCCESS ) {
+ amrex::Print() << " cufftplan3d execute failed ! \n";
+ }
+ cudaDeviceSynchronize();
#else
fftw_execute( forward_plan[mfi] );
#endif
@@ -227,6 +264,22 @@ SpectralFieldData::BackwardTransform( MultiFab& mf,
#ifdef AMREX_USE_GPU
// Add cuFFT-specific code ; make sure that this is done on the same
// GPU stream as the above copy
+ cudaDeviceSynchronize();
+ cufftResult result;
+ cudaStream_t stream = amrex::Cuda::Device::cudaStream();
+ amrex::Print() << " stream is " << stream << "\n";
+ cufftSetStream ( backward_plan[mfi], stream);
+ result = cufftExecZ2D( backward_plan[mfi],
+ reinterpret_cast<cuDoubleComplex*>(
+ tmpSpectralField[mfi].dataPtr()),
+ tmpRealField[mfi].dataPtr() );
+ if ( result != CUFFT_SUCCESS ) {
+ amrex::Print() << " cufftplan3d execute inverse failed ! \n";
+ }
+ if ( result == CUFFT_SUCCESS ) {
+ amrex::Print() << " created cufft inverse transform\n";
+ }
+ cudaDeviceSynchronize();
#else
fftw_execute( backward_plan[mfi] );
#endif
@@ -245,6 +298,14 @@ SpectralFieldData::BackwardTransform( MultiFab& mf,
// Copy and normalize field
mf_arr(i,j,k,i_comp) = inv_N*tmp_arr(i,j,k);
});
+//#ifdef AMREX_USE_GPU
+// cudaDeviceSynchronize();
+//#endif
+
+#ifdef AMREX_USE_GPU
+ cudaDeviceSynchronize();
+#endif
+ amrex::Print() << " divided by 1/N \n";
}
}
}
diff --git a/Source/FieldSolver/SpectralSolver/SpectralKSpace.cpp b/Source/FieldSolver/SpectralSolver/SpectralKSpace.cpp
index 2fe78cedd..6a88a52a3 100644
--- a/Source/FieldSolver/SpectralSolver/SpectralKSpace.cpp
+++ b/Source/FieldSolver/SpectralSolver/SpectralKSpace.cpp
@@ -142,9 +142,13 @@ SpectralKSpace::getSpectralShiftFactor( const DistributionMapping& dm,
case ShiftType::TransformFromCellCentered: sign = -1.; break;
case ShiftType::TransformToCellCentered: sign = 1.;
}
- constexpr Complex I{0,1};
+ const Complex I{0,1};
for (int i=0; i<k.size(); i++ ){
+#ifdef AMREX_USE_GPU
+ shift[i] = thrust::exp( I*sign*k[i]*0.5*dx[i_dim] );
+#else
shift[i] = std::exp( I*sign*k[i]*0.5*dx[i_dim] );
+#endif
}
}
return shift_factor;
diff --git a/Source/FieldSolver/WarpXFFT.cpp b/Source/FieldSolver/WarpXFFT.cpp
index 1cf5460f2..c4e0461f9 100644
--- a/Source/FieldSolver/WarpXFFT.cpp
+++ b/Source/FieldSolver/WarpXFFT.cpp
@@ -490,4 +490,5 @@ WarpX::PushPSATD (int lev, amrex::Real /* dt */)
{
amrex::Abort("WarpX::PushPSATD: TODO");
}
+ amrex::Print() << " coped data from fft to valid \n";
}
diff --git a/Source/Initialization/WarpXInitData.cpp b/Source/Initialization/WarpXInitData.cpp
index 23637ec97..b9f27f07e 100644
--- a/Source/Initialization/WarpXInitData.cpp
+++ b/Source/Initialization/WarpXInitData.cpp
@@ -318,16 +318,28 @@ WarpX::InitLevelData (int lev, Real time)
void
WarpX::InitLevelDataFFT (int lev, Real time)
{
+
+ amrex::Print() << " print out pointer " << &Efield_fp_fft[lev][0] << "\n";
Efield_fp_fft[lev][0]->setVal(0.0);
+ amrex::Print() << " ex set \n";
Efield_fp_fft[lev][1]->setVal(0.0);
+ amrex::Print() << " ey set \n";
Efield_fp_fft[lev][2]->setVal(0.0);
+ amrex::Print() << " ez set \n";
Bfield_fp_fft[lev][0]->setVal(0.0);
+ amrex::Print() << " bx set \n";
Bfield_fp_fft[lev][1]->setVal(0.0);
+ amrex::Print() << " by set \n";
Bfield_fp_fft[lev][2]->setVal(0.0);
+ amrex::Print() << " bz set \n";
current_fp_fft[lev][0]->setVal(0.0);
+ amrex::Print() << " jx set \n";
current_fp_fft[lev][1]->setVal(0.0);
+ amrex::Print() << " jy set \n";
current_fp_fft[lev][2]->setVal(0.0);
+ amrex::Print() << " jz set \n";
rho_fp_fft[lev]->setVal(0.0);
+ amrex::Print() << " rhofp set \n";
if (lev > 0)
{
@@ -342,6 +354,7 @@ WarpX::InitLevelDataFFT (int lev, Real time)
current_cp_fft[lev][2]->setVal(0.0);
rho_cp_fft[lev]->setVal(0.0);
}
+
}
#endif
diff --git a/Source/Parallelization/WarpXComm.cpp b/Source/Parallelization/WarpXComm.cpp
index 348b31c8b..e9d814a94 100644
--- a/Source/Parallelization/WarpXComm.cpp
+++ b/Source/Parallelization/WarpXComm.cpp
@@ -222,7 +222,9 @@ WarpX::FillBoundaryE ()
{
for (int lev = 0; lev <= finest_level; ++lev)
{
+ amrex::Print() << " lev " << lev << " in fill bndry E \n";
FillBoundaryE(lev);
+ amrex::Print() << " lev " << lev << " after fill bndry E \n";
}
}
@@ -246,7 +248,7 @@ void
WarpX::FillBoundaryE (int lev, PatchType patch_type)
{
if (patch_type == PatchType::fine)
- {
+ {
if (do_pml && pml[lev]->ok())
{
pml[lev]->ExchangeE(patch_type,
@@ -271,9 +273,12 @@ WarpX::FillBoundaryE (int lev, PatchType patch_type)
pml[lev]->FillBoundaryE(patch_type);
}
+ amrex::Print() << " before defining multifab \n";
const auto& cperiod = Geom(lev-1).periodicity();
Vector<MultiFab*> mf{Efield_cp[lev][0].get(),Efield_cp[lev][1].get(),Efield_cp[lev][2].get()};
+ amrex::Print() << " before amrex fill bndry \n";
amrex::FillBoundary(mf, cperiod);
+ amrex::Print() << " after amrex fill bndry \n";
}
}
diff --git a/Source/Particles/PhysicalParticleContainer.cpp b/Source/Particles/PhysicalParticleContainer.cpp
index 2fa39d87d..a1d4e1319 100644
--- a/Source/Particles/PhysicalParticleContainer.cpp
+++ b/Source/Particles/PhysicalParticleContainer.cpp
@@ -1332,7 +1332,9 @@ PhysicalParticleContainer::Evolve (int lev,
pti.GetPosition(m_xp[thread_num], m_yp[thread_num], m_zp[thread_num]);
BL_PROFILE_VAR_STOP(blp_copy);
+ amrex::Print() << " before deposit chage \n";
if (rho) DepositCharge(pti, wp, rho, crho, 0, np_current, np, thread_num, lev);
+ amrex::Print() << " after deposit chage \n";
if (! do_not_push)
{
diff --git a/Source/Particles/WarpXParticleContainer.cpp b/Source/Particles/WarpXParticleContainer.cpp
index 80f3882a0..89e21df1c 100644
--- a/Source/Particles/WarpXParticleContainer.cpp
+++ b/Source/Particles/WarpXParticleContainer.cpp
@@ -608,6 +608,7 @@ WarpXParticleContainer::DepositCharge ( WarpXParIter& pti, RealVector& wp,
const std::array<Real, 3>& xyzmin = xyzmin_tile;
#ifdef AMREX_USE_GPU
+ amrex::Print() << " before icomp data ptr " << icomp << "\n";
data_ptr = (*rhomf)[pti].dataPtr(icomp);
auto rholen = (*rhomf)[pti].length();
#else
@@ -630,6 +631,7 @@ WarpXParticleContainer::DepositCharge ( WarpXParIter& pti, RealVector& wp,
const long nz = rholen[1]-1-2*ngRho;
#endif
BL_PROFILE_VAR_START(blp_pxr_chd);
+ amrex::Print() << " before warpxcharge deposition \n";
warpx_charge_deposition(data_ptr, &np_current,
m_xp[thread_num].dataPtr(),
m_yp[thread_num].dataPtr(),
@@ -665,7 +667,7 @@ WarpXParticleContainer::DepositCharge ( WarpXParIter& pti, RealVector& wp,
const std::array<Real,3>& cxyzmin_tile = WarpX::LowerCorner(ctilebox, lev-1);
#ifdef AMREX_USE_GPU
- data_ptr = (*crhomf)[pti].dataPtr();
+ data_ptr = (*crhomf)[pti].dataPtr(icomp);
auto rholen = (*crhomf)[pti].length();
#else
tile_box = amrex::convert(ctilebox, IntVect::TheUnitVector());
diff --git a/Source/WarpX.cpp b/Source/WarpX.cpp
index a3a24897a..83788d5b9 100644
--- a/Source/WarpX.cpp
+++ b/Source/WarpX.cpp
@@ -556,7 +556,9 @@ WarpX::MakeNewLevelFromScratch (int lev, Real time, const BoxArray& new_grids,
InitLevelData(lev, time);
#ifdef WARPX_USE_PSATD
+ amrex::Print() << " alloc level data fft \n";
AllocLevelDataFFT(lev);
+ amrex::Print() << " init level data fft \n";
InitLevelDataFFT(lev, time);
#endif
}
diff --git a/Source/main.cpp b/Source/main.cpp
index d89e89e47..2e1a5e9cf 100644
--- a/Source/main.cpp
+++ b/Source/main.cpp
@@ -30,10 +30,13 @@ int main(int argc, char* argv[])
const Real strt_total = amrex::second();
{
+ amrex::Print() << " about to construct warpx \n";
WarpX warpx;
+ amrex::Print() << " call warpx init data \n";
warpx.InitData();
+ amrex::Print() << " call warpx evolve \n";
warpx.Evolve();
Real end_total = amrex::second() - strt_total;
generated by cgit v1.2.3 (git 2.39.1) at 2025-08-07 19:43:04 +0000