Code indexing in gitaly is broken and leads to code not being visible to the user. We work on the issue with highest priority.

Skip to content
Snippets Groups Projects
Commit 615eca59 authored by adelmann's avatar adelmann :reminder_ribbon:
Browse files

add space charge NOT tested

parent e18474e1
No related branches found
No related tags found
No related merge requests found
...@@ -17,7 +17,7 @@ Value,{OPALVERSION}; ...@@ -17,7 +17,7 @@ Value,{OPALVERSION};
// Global Parameters // Global Parameters
REAL rf_freq = 1.3e9; // RF frequency. (Hz) REAL rf_freq = 1.3e9; // RF frequency. (Hz)
REAL n_particles = 1E6; // Number of particles in simulation. REAL n_particles = 1E5; // Number of particles in simulation.
REAL beam_bunch_charge = 1e-9; // Charge of bunch. (C) REAL beam_bunch_charge = 1e-9; // Charge of bunch. (C)
// ---------------------------------------------------------------------------- // ----------------------------------------------------------------------------
...@@ -32,6 +32,7 @@ REAL P0 = gamma*beta*EMASS; //inital z momentum ...@@ -32,6 +32,7 @@ REAL P0 = gamma*beta*EMASS; //inital z momentum
value, {Edes, P0}; value, {Edes, P0};
D1: DRIFT, L = 1.0, ELEMEDGE = 0.0; D1: DRIFT, L = 1.0, ELEMEDGE = 0.0;
myLine: Line = (D1);
BEAM1: BEAM, PARTICLE = ELECTRON, pc = P0, NPART = n_particles, BEAM1: BEAM, PARTICLE = ELECTRON, pc = P0, NPART = n_particles,
BFREQ = rf_freq, BCURRENT = beam_bunch_charge * rf_freq * 1E6, CHARGE = -1; BFREQ = rf_freq, BCURRENT = beam_bunch_charge * rf_freq * 1E6, CHARGE = -1;
...@@ -41,9 +42,9 @@ FS1: Fieldsolver, NX=32, NY=32, NZ=32, TYPE=FFT, PARFFTX = true, PARFFTY = true, ...@@ -41,9 +42,9 @@ FS1: Fieldsolver, NX=32, NY=32, NZ=32, TYPE=FFT, PARFFTX = true, PARFFTY = true,
BCFFTX=OPEN, BCFFTY=OPEN, BCFFTZ=OPEN, BCFFTX=OPEN, BCFFTY=OPEN, BCFFTZ=OPEN,
BBOXINCR = 1, GREENSF = INTEGRATED; BBOXINCR = 1, GREENSF = INTEGRATED;
Dist1: DISTRIBUTION, TYPE= GAUSS, SIGMAX=1, SIGMAY=1, SIGMAZ=1, SIGMAPX=1, SIGMAPY=1, SIGMAPZ=1; Dist1: DISTRIBUTION, TYPE= GAUSS, SIGMAX=1e-3, SIGMAY=1e-3, SIGMAZ=1e-3, SIGMAPX=1e-6, SIGMAPY=1e-6, SIGMAPZ=1e-6;
myLine: Line = (D1);
TRACK, LINE = myLine, BEAM = BEAM1, MAXSTEPS = 1, DT = {1e-12}, ZSTOP=1.0; TRACK, LINE = myLine, BEAM = BEAM1, MAXSTEPS = 1, DT = {1e-12}, ZSTOP=1.0;
RUN, METHOD = "PARALLEL", BEAM = BEAM1, FIELDSOLVER = FS1, DISTRIBUTION = Dist1; RUN, METHOD = "PARALLEL", BEAM = BEAM1, FIELDSOLVER = FS1, DISTRIBUTION = Dist1;
......
SDDS1
&description
text="Statistics data '../input-files/test-0002.in' 20/06/2024 23:23:03",
contents="stat parameters"
&end
&parameter
name=processors,
type=long,
description="Number of Cores used"
&end
&parameter
name=revision,
type=string,
description="git revision of opal"
&end
&parameter
name=flavor,
type=string,
description="OPAL flavor that wrote file"
&end
&column
name=t,
type=double,
units=ns,
description="1 Time"
&end
&column
name=s,
type=double,
units=m,
description="2 Path length"
&end
&column
name=numParticles,
type=long,
units=1,
description="3 Number of Macro Particles"
&end
&column
name=charge,
type=double,
units=1,
description="4 Bunch Charge"
&end
&column
name=energy,
type=double,
units=MeV,
description="5 Mean Bunch Energy"
&end
&column
name=rms_x,
type=double,
units=m,
description="6 RMS Beamsize in x"
&end
&column
name=rms_y,
type=double,
units=m,
description="7 RMS Beamsize in y"
&end
&column
name=rms_s,
type=double,
units=m,
description="8 RMS Beamsize in s"
&end
&column
name=rms_px,
type=double,
units=1,
description="9 RMS Normalized Momenta in x"
&end
&column
name=rms_py,
type=double,
units=1,
description="10 RMS Normalized Momenta in y"
&end
&column
name=rms_ps,
type=double,
units=1,
description="11 RMS Normalized Momenta in s"
&end
&column
name=emit_x,
type=double,
units=m,
description="12 Normalized Emittance x"
&end
&column
name=emit_y,
type=double,
units=m,
description="13 Normalized Emittance y"
&end
&column
name=emit_s,
type=double,
units=m,
description="14 Normalized Emittance s"
&end
&column
name=mean_x,
type=double,
units=m,
description="15 Mean Beam Position in x"
&end
&column
name=mean_y,
type=double,
units=m,
description="16 Mean Beam Position in y"
&end
&column
name=mean_s,
type=double,
units=m,
description="17 Mean Beam Position in s"
&end
&column
name=ref_x,
type=double,
units=m,
description="18 x coordinate of reference particle in lab cs"
&end
&column
name=ref_y,
type=double,
units=m,
description="19 y coordinate of reference particle in lab cs"
&end
&column
name=ref_z,
type=double,
units=m,
description="20 z coordinate of reference particle in lab cs"
&end
&column
name=ref_px,
type=double,
units=1,
description="21 x momentum of reference particle in lab cs"
&end
&column
name=ref_py,
type=double,
units=1,
description="22 y momentum of reference particle in lab cs"
&end
&column
name=ref_pz,
type=double,
units=1,
description="23 z momentum of reference particle in lab cs"
&end
&column
name=max_x,
type=double,
units=m,
description="24 Max Beamsize in x"
&end
&column
name=max_y,
type=double,
units=m,
description="25 Max Beamsize in y"
&end
&column
name=max_s,
type=double,
units=m,
description="26 Max Beamsize in s"
&end
&column
name=xpx,
type=double,
units=1,
description="27 Correlation xpx"
&end
&column
name=ypy,
type=double,
units=1,
description="28 Correlation ypy"
&end
&column
name=zpz,
type=double,
units=1,
description="29 Correlation zpz"
&end
&column
name=Dx,
type=double,
units=m,
description="30 Dispersion in x"
&end
&column
name=DDx,
type=double,
units=1,
description="31 Derivative of dispersion in x"
&end
&column
name=Dy,
type=double,
units=m,
description="32 Dispersion in y"
&end
&column
name=DDy,
type=double,
units=1,
description="33 Derivative of dispersion in y"
&end
&column
name=Bx_ref,
type=double,
units=T,
description="34 Bx-Field component of ref particle"
&end
&column
name=By_ref,
type=double,
units=T,
description="35 By-Field component of ref particle"
&end
&column
name=Bz_ref,
type=double,
units=T,
description="36 Bz-Field component of ref particle"
&end
&column
name=Ex_ref,
type=double,
units=MV/m,
description="37 Ex-Field component of ref particle"
&end
&column
name=Ey_ref,
type=double,
units=MV/m,
description="38 Ey-Field component of ref particle"
&end
&column
name=Ez_ref,
type=double,
units=MV/m,
description="39 Ez-Field component of ref particle"
&end
&column
name=dE,
type=double,
units=MeV,
description="40 energy spread of the beam"
&end
&column
name=dt,
type=double,
units=ns,
description="41 time step size"
&end
&column
name=partsOutside,
type=double,
units=1,
description="42 outside n*sigma of the beam"
&end
&column
name=DebyeLength,
type=double,
units=m,
description="43 Debye length in the boosted frame"
&end
&column
name=plasmaParameter,
type=double,
units=1,
description="44 Plasma parameter that gives no. of particles in a Debye sphere"
&end
&column
name=temperature,
type=double,
units=K,
description="45 Temperature of the beam"
&end
&column
name=rmsDensity,
type=double,
units=1,
description="46 RMS number density of the beam"
&end
&data
mode=ascii,
no_row_counts=1
&end
1
OPALX 2024.1.00 git rev. #640dbb1f34e2dc23314bad26c94148a0e578fe17
opal-t
1.000000000000000e-03 2.997924188990765e-04 100000 -1.000000000000000e-09 9.999999999999980e+02 9.886771617229491e-04 9.900456978581636e-04 9.861933284523727e-04 9.970719439431197e-07 9.971772595976683e-07 1.003251673042217e-06 9.857814897437524e-10 9.872476063403890e-10 9.893927462964761e-10 1.695475357332832e-20 1.706100538623190e-20 4.746420270707041e-20 1.695475356915250e-20 1.706100538100345e-20 2.997924188990765e-04 -2.727234331922033e-23 -3.414707447737305e-23 1.957950928190183e+03 3.000478311466378e-03 2.999914659247034e-03 2.996972298154323e-03 -1.249744196182289e-03 2.643507480662082e-03 3.857264382145191e-03 -7.000706356679798e-13 4.605817787904672e-15 3.761377802323194e-12 -1.333863982837420e-16 0.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 1.621166691380730e-04 1.000000000000000e-03 0 3.883428884728593e-10 7.942613332657146e-15 1.025293560371121e-09 3.237641375683589e+13
/* /*
Test simple Gaussian distribution in a drift with space charge Test simple Gaussian distribution
*/ */
OPTION, PSDUMPFREQ = 10; // 6d data written every 10th time step (h5). OPTION, PSDUMPFREQ = 10; // 6d data written every 10th time step (h5).
...@@ -9,7 +9,7 @@ OPTION, AUTOPHASE = 4; // Autophase is on, and phase of max energy ...@@ -9,7 +9,7 @@ OPTION, AUTOPHASE = 4; // Autophase is on, and phase of max energy
// gain will be found automatically for cavities. // gain will be found automatically for cavities.
OPTION, VERSION=10900; OPTION, VERSION=10900;
Title, string="Test simple Gaussian distribution "; Title, string="Test simple Gaussian distribution with space charge";
Value,{OPALVERSION}; Value,{OPALVERSION};
...@@ -32,6 +32,7 @@ REAL P0 = gamma*beta*EMASS; //inital z momentum ...@@ -32,6 +32,7 @@ REAL P0 = gamma*beta*EMASS; //inital z momentum
value, {Edes, P0}; value, {Edes, P0};
D1: DRIFT, L = 1.0, ELEMEDGE = 0.0; D1: DRIFT, L = 1.0, ELEMEDGE = 0.0;
myLine: Line = (D1);
BEAM1: BEAM, PARTICLE = ELECTRON, pc = P0, NPART = n_particles, BEAM1: BEAM, PARTICLE = ELECTRON, pc = P0, NPART = n_particles,
BFREQ = rf_freq, BCURRENT = beam_bunch_charge * rf_freq * 1E6, CHARGE = -1; BFREQ = rf_freq, BCURRENT = beam_bunch_charge * rf_freq * 1E6, CHARGE = -1;
...@@ -41,11 +42,11 @@ FS1: Fieldsolver, NX=32, NY=32, NZ=32, TYPE=FFT, PARFFTX = true, PARFFTY = true, ...@@ -41,11 +42,11 @@ FS1: Fieldsolver, NX=32, NY=32, NZ=32, TYPE=FFT, PARFFTX = true, PARFFTY = true,
BCFFTX=OPEN, BCFFTY=OPEN, BCFFTZ=OPEN, BCFFTX=OPEN, BCFFTY=OPEN, BCFFTZ=OPEN,
BBOXINCR = 1, GREENSF = INTEGRATED; BBOXINCR = 1, GREENSF = INTEGRATED;
Dist1: DISTRIBUTION, TYPE= GAUSS, SIGMAX=1, SIGMAY=1, SIGMAZ=1, SIGMAPX=1, SIGMAPY=1, SIGMAPZ=1; Dist1: DISTRIBUTION, TYPE= GAUSS, SIGMAX=1e-3, SIGMAY=1e-3, SIGMAZ=1e-3, SIGMAPX=1e-6, SIGMAPY=1e-6, SIGMAPZ=1e-6;
myLine: Line = (D1);
TRACK, LINE = myLine, BEAM = BEAM1, MAXSTEPS = 10, DT = {1e-12}, ZSTOP=1.0; TRACK, LINE = myLine, BEAM = BEAM1, MAXSTEPS = 2, DT = {1e-12}, ZSTOP=1.0;
RUN, METHOD = "PARALLEL", BEAM = BEAM1, FIELDSOLVER = FS1, DISTRIBUTION = Dist1; RUN, METHOD = "PARALLEL", BEAM = BEAM1, FIELDSOLVER = FS1, DISTRIBUTION = Dist1;
ENDTRACK; ENDTRACK;
QUIT; QUIT;
SDDS1
&description
text="Statistics data '../input-files/test-0003.in' 20/06/2024 23:34:58",
contents="stat parameters"
&end
&parameter
name=processors,
type=long,
description="Number of Cores used"
&end
&parameter
name=revision,
type=string,
description="git revision of opal"
&end
&parameter
name=flavor,
type=string,
description="OPAL flavor that wrote file"
&end
&column
name=t,
type=double,
units=ns,
description="1 Time"
&end
&column
name=s,
type=double,
units=m,
description="2 Path length"
&end
&column
name=numParticles,
type=long,
units=1,
description="3 Number of Macro Particles"
&end
&column
name=charge,
type=double,
units=1,
description="4 Bunch Charge"
&end
&column
name=energy,
type=double,
units=MeV,
description="5 Mean Bunch Energy"
&end
&column
name=rms_x,
type=double,
units=m,
description="6 RMS Beamsize in x"
&end
&column
name=rms_y,
type=double,
units=m,
description="7 RMS Beamsize in y"
&end
&column
name=rms_s,
type=double,
units=m,
description="8 RMS Beamsize in s"
&end
&column
name=rms_px,
type=double,
units=1,
description="9 RMS Normalized Momenta in x"
&end
&column
name=rms_py,
type=double,
units=1,
description="10 RMS Normalized Momenta in y"
&end
&column
name=rms_ps,
type=double,
units=1,
description="11 RMS Normalized Momenta in s"
&end
&column
name=emit_x,
type=double,
units=m,
description="12 Normalized Emittance x"
&end
&column
name=emit_y,
type=double,
units=m,
description="13 Normalized Emittance y"
&end
&column
name=emit_s,
type=double,
units=m,
description="14 Normalized Emittance s"
&end
&column
name=mean_x,
type=double,
units=m,
description="15 Mean Beam Position in x"
&end
&column
name=mean_y,
type=double,
units=m,
description="16 Mean Beam Position in y"
&end
&column
name=mean_s,
type=double,
units=m,
description="17 Mean Beam Position in s"
&end
&column
name=ref_x,
type=double,
units=m,
description="18 x coordinate of reference particle in lab cs"
&end
&column
name=ref_y,
type=double,
units=m,
description="19 y coordinate of reference particle in lab cs"
&end
&column
name=ref_z,
type=double,
units=m,
description="20 z coordinate of reference particle in lab cs"
&end
&column
name=ref_px,
type=double,
units=1,
description="21 x momentum of reference particle in lab cs"
&end
&column
name=ref_py,
type=double,
units=1,
description="22 y momentum of reference particle in lab cs"
&end
&column
name=ref_pz,
type=double,
units=1,
description="23 z momentum of reference particle in lab cs"
&end
&column
name=max_x,
type=double,
units=m,
description="24 Max Beamsize in x"
&end
&column
name=max_y,
type=double,
units=m,
description="25 Max Beamsize in y"
&end
&column
name=max_s,
type=double,
units=m,
description="26 Max Beamsize in s"
&end
&column
name=xpx,
type=double,
units=1,
description="27 Correlation xpx"
&end
&column
name=ypy,
type=double,
units=1,
description="28 Correlation ypy"
&end
&column
name=zpz,
type=double,
units=1,
description="29 Correlation zpz"
&end
&column
name=Dx,
type=double,
units=m,
description="30 Dispersion in x"
&end
&column
name=DDx,
type=double,
units=1,
description="31 Derivative of dispersion in x"
&end
&column
name=Dy,
type=double,
units=m,
description="32 Dispersion in y"
&end
&column
name=DDy,
type=double,
units=1,
description="33 Derivative of dispersion in y"
&end
&column
name=Bx_ref,
type=double,
units=T,
description="34 Bx-Field component of ref particle"
&end
&column
name=By_ref,
type=double,
units=T,
description="35 By-Field component of ref particle"
&end
&column
name=Bz_ref,
type=double,
units=T,
description="36 Bz-Field component of ref particle"
&end
&column
name=Ex_ref,
type=double,
units=MV/m,
description="37 Ex-Field component of ref particle"
&end
&column
name=Ey_ref,
type=double,
units=MV/m,
description="38 Ey-Field component of ref particle"
&end
&column
name=Ez_ref,
type=double,
units=MV/m,
description="39 Ez-Field component of ref particle"
&end
&column
name=dE,
type=double,
units=MeV,
description="40 energy spread of the beam"
&end
&column
name=dt,
type=double,
units=ns,
description="41 time step size"
&end
&column
name=partsOutside,
type=double,
units=1,
description="42 outside n*sigma of the beam"
&end
&column
name=DebyeLength,
type=double,
units=m,
description="43 Debye length in the boosted frame"
&end
&column
name=plasmaParameter,
type=double,
units=1,
description="44 Plasma parameter that gives no. of particles in a Debye sphere"
&end
&column
name=temperature,
type=double,
units=K,
description="45 Temperature of the beam"
&end
&column
name=rmsDensity,
type=double,
units=1,
description="46 RMS number density of the beam"
&end
&data
mode=ascii,
no_row_counts=1
&end
1
OPALX 2024.1.00 git rev. #640dbb1f34e2dc23314bad26c94148a0e578fe17
opal-t
1.000000000000000e-03 2.997924188990765e-04 100000 -1.000000000000000e-09 9.999999999999980e+02 9.886771617229491e-04 9.900456978581636e-04 9.861933284523727e-04 9.970719439431197e-07 9.971772595976683e-07 1.003251673042217e-06 9.857814897437524e-10 9.872476063403890e-10 9.893927462964761e-10 1.695475357332832e-20 1.706100538623190e-20 4.746420270707041e-20 1.695475356915250e-20 1.706100538100345e-20 2.997924188990765e-04 -2.727234331922033e-23 -3.414707447737305e-23 1.957950928190183e+03 3.000478311466378e-03 2.999914659247034e-03 2.996972298154323e-03 -1.249744196182289e-03 2.643507480662082e-03 3.857264382145191e-03 -7.000706356679798e-13 4.605817787904672e-15 3.761377802323194e-12 -1.333863982837420e-16 0.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 1.621166691380730e-04 1.000000000000000e-03 0 3.883428884728593e-10 7.942613332657146e-15 1.025293560371121e-09 3.237641375683589e+13
2.000000000000000e-03 8.993772566972294e-04 100000 -1.000000000000000e-09 9.999999999999980e+02 9.886771617229491e-04 9.900456978581636e-04 9.861933284523727e-04 9.970719439431197e-07 9.971772595976683e-07 1.003251673042217e-06 9.857814897437524e-10 9.872476063403890e-10 9.893927462964761e-10 1.695475357332832e-20 1.706100538623190e-20 4.746420270707041e-20 1.695475356497669e-20 1.706100537577501e-20 5.995848377981529e-04 -2.727234331922033e-23 -3.414707447737305e-23 1.957950928190183e+03 3.000478311466378e-03 2.999914659247034e-03 2.996972298154323e-03 -1.249744196182289e-03 2.643507480662082e-03 3.857264382145191e-03 -7.000706356679798e-13 4.605817787904672e-15 3.761377802323194e-12 -1.333863982837420e-16 0.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 0.000000000000000e+00 1.621166691380730e-04 1.000000000000000e-03 0 3.883407084178780e-10 7.942568744911947e-15 1.025293560371121e-09 3.237677726526380e+13
...@@ -551,6 +551,11 @@ void ParallelTracker::emitParticles(long long step) { ...@@ -551,6 +551,11 @@ void ParallelTracker::emitParticles(long long step) {
} }
void ParallelTracker::computeSpaceChargeFields(unsigned long long step) { void ParallelTracker::computeSpaceChargeFields(unsigned long long step) {
Inform m("INFORM_ALL_NODES");
m << "in ParallelTracker::computeSpaceChargeFields" << endl;
if (!itsBunch_m->hasFieldSolver()) { if (!itsBunch_m->hasFieldSolver()) {
return; return;
} }
......
...@@ -213,54 +213,103 @@ void PartBunch<double,3>::bunchUpdate() { ...@@ -213,54 +213,103 @@ void PartBunch<double,3>::bunchUpdate() {
template <> template <>
void PartBunch<double,3>::computeSelfFields() { void PartBunch<double,3>::computeSelfFields() {
static IpplTimings::TimerRef SolveTimer = IpplTimings::getTimer("SolveTimer"); Inform m("INFORM_ALL_NODES");
IpplTimings::startTimer(SolveTimer); static IpplTimings::TimerRef SolveTimer = IpplTimings::getTimer("SolveTimer");
IpplTimings::startTimer(SolveTimer);
Field_t<3>& rho = this->fcontainer_m->getRho();
ippl::ParticleAttrib<double>& Q = this->pcontainer_m->Q; Field_t<3>& rho = this->fcontainer_m->getRho();
typename Base::particle_position_type& R = this->pcontainer_m->R; ippl::ParticleAttrib<double>& Q = this->pcontainer_m->Q;
typename Base::particle_position_type& R = this->pcontainer_m->R;
rho = 0.0;
Q = Q*this->pcontainer_m->dt; rho = 0.0;
this->qi_m = this->qi_m * getdT(); Q = Q*this->pcontainer_m->dt;
scatter(Q, rho, R); this->qi_m = this->qi_m * getdT();
Q = Q/this->pcontainer_m->dt; scatter(Q, rho, R);
this->qi_m = this->qi_m / getdT(); Q = Q/this->pcontainer_m->dt;
rho = rho/getdT(); this->qi_m = this->qi_m / getdT();
rho = rho/getdT();
double gammaz = this->pcontainer_m->getMeanGammaZ();
double scaleFactor = 1;
// double scaleFactor = Physics::c * getdT(); // some debug output ------------------------------------------------------------
//and get meshspacings in real units [m] //
Vector_t<double, 3> hr_scaled = hr_m * scaleFactor; /*
hr_scaled[2] *= gammaz; Kokkos::parallel_for("print q", ippl::getRangePolicy(Qview),
KOKKOS_LAMBDA(const int i) {
double localDensity2 = 0, Density2=0; if (i<5){
double tmp2 = 1 / hr_scaled[0] * 1 / hr_scaled[1] * 1 / hr_scaled[2]; double myQ = Qview(i);
std::cout << "qi= " << myQ << std::endl;
using index_array_type = typename ippl::RangePolicy<Dim>::index_array_type; }
});
//divide charge by a 'grid-cube' volume to get [C/m^3] */
rho = rho *tmp2;
auto rhoView = rho.getView();
double Npoints = nr_m[0] * nr_m[1] * nr_m[2]; using index_array_type_3d = typename ippl::RangePolicy<Dim>::index_array_type;
const index_array_type_3d& args3d = {8,8,8};
auto rhoView = rho.getView(); auto res3d = ippl::apply(rhoView,args3d);
localDensity2 = 0.; m << "Type of variable: " << typeid(res3d).name() << " rho[8,8,8]= " << res3d << endl;
ippl::parallel_reduce(
"Density stats", ippl::getRangePolicy(rhoView), double qsum = 0.0;
KOKKOS_LAMBDA(const index_array_type& args, double& den) { for (int i=0; i<nr_m[0]; i++) {
double val = ippl::apply(rhoView, args); for (int j=0; j<nr_m[1]; j++) {
den += Kokkos::pow(val, 2); for (int k=0; k<nr_m[2]; k++) {
}, const index_array_type_3d& args3d = {i,j,k};
Kokkos::Sum<double>(localDensity2) ); qsum += ippl::apply(rhoView,args3d);
ippl::Comm->reduce(localDensity2, Density2, 1, std::plus<double>()); }
}
rmsDensity_m = std::sqrt( (1.0 /Npoints) * Density2 / Physics::q_e / Physics::q_e ); }
this->calcDebyeLength(); m << "sum(rho_m)= " << qsum << endl;
IpplTimings::stopTimer(SolveTimer);
auto Qview = this->pcontainer_m->Q.getView();
using index_array_type_1d = typename ippl::RangePolicy<1>::index_array_type;
const index_array_type_1d& args1d = {0};
auto res1d = ippl::apply(Qview,args1d);
m << "Type of variable: " << typeid(res1d).name() << " value " << res1d << endl;
//
// -------------------------------------------------------------------------------
double gammaz = this->pcontainer_m->getMeanGammaZ();
double scaleFactor = 1;
// double scaleFactor = Physics::c * getdT();
//and get meshspacings in real units [m]
Vector_t<double, 3> hr_scaled = hr_m * scaleFactor;
hr_scaled[2] *= gammaz;
double localDensity2 = 0, Density2=0;
double tmp2 = 1 / hr_scaled[0] * 1 / hr_scaled[1] * 1 / hr_scaled[2];
using index_array_type = typename ippl::RangePolicy<Dim>::index_array_type;
//divide charge by a 'grid-cube' volume to get [C/m^3]
rho = rho *tmp2;
double Npoints = nr_m[0] * nr_m[1] * nr_m[2];
// auto rhoView = rho.getView();
localDensity2 = 0.;
ippl::parallel_reduce(
"Density stats", ippl::getRangePolicy(rhoView),
KOKKOS_LAMBDA(const index_array_type& args, double& den) {
double val = ippl::apply(rhoView, args);
den += Kokkos::pow(val, 2);
},
Kokkos::Sum<double>(localDensity2) );
ippl::Comm->reduce(localDensity2, Density2, 1, std::plus<double>());
rmsDensity_m = std::sqrt( (1.0 /Npoints) * Density2 / Physics::q_e / Physics::q_e );
this->calcDebyeLength();
m << "gammaz = " << gammaz << endl;
m << "hr_scaled = " << hr_scaled << endl;
this->fsolver_m->runSolver();
gather(this->pcontainer_m->E, this->fcontainer_m->getE(), this->pcontainer_m->R);
IpplTimings::stopTimer(SolveTimer);
} }
template <> template <>
...@@ -279,10 +328,14 @@ void PartBunch<double,3>::scatterCIC() { ...@@ -279,10 +328,14 @@ void PartBunch<double,3>::scatterCIC() {
scatter(*q, *rho, *R); scatter(*q, *rho, *R);
m << "gammz= " << this->pcontainer_m->getMeanP()[2] << endl;
double relError = std::fabs((Q - (*rho).sum()) / Q); double relError = std::fabs((Q - (*rho).sum()) / Q);
size_type TotalParticles = 0; size_type TotalParticles = 0;
size_type localParticles = this->pcontainer_m->getLocalNum(); size_type localParticles = this->pcontainer_m->getLocalNum();
m << "computeSelfFields sum rho " << (*rho).sum() << endl;
ippl::Comm->reduce(localParticles, TotalParticles, 1, std::plus<size_type>()); ippl::Comm->reduce(localParticles, TotalParticles, 1, std::plus<size_type>());
if (ippl::Comm->rank() == 0) { if (ippl::Comm->rank() == 0) {
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment