1 #ifndef AMREX_PARTICLEINIT_H
2 #define AMREX_PARTICLEINIT_H
3 #include <AMReX_Config.H>
33 template <
typename ParticleType,
int NArrayReal,
int NArrayInt,
34 template<
class>
class Allocator,
class CellAssignor>
36 ParticleContainer_impl<ParticleType, NArrayReal, NArrayInt, Allocator, CellAssignor>
37 ::InitFromAsciiFile (
const std::string& file,
int extradata,
const IntVect* Nrep)
39 BL_PROFILE(
"ParticleContainer<NSR, NSI, NAR, NAI>::InitFromAsciiFile()");
47 int NReaders = MaxReaders();
78 if (Nrep !=
nullptr) {
83 Long how_many_read = 0;
87 if (extradata > NStructReal) { nreals.resize(extradata - NStructReal); }
95 ifs.rdbuf()->pubsetbuf(io_buffer.
dataPtr(), io_buffer.
size());
97 ifs.open(file.c_str(), std::ios::in);
106 ifs >> cnt >> std::ws;
112 if (extradata > NStructReal) {
r.resize(extradata - NStructReal); }
114 const int Chunk = cnt / NReaders;
116 for (
int i = 0; i <
MyProc*Chunk; i++)
124 std::string msg(
"ParticleContainer::InitFromAsciiFile(");
126 msg +=
") failed @ 1";
132 if (
MyProc == (NReaders - 1))
134 MyCnt += cnt % NReaders;
139 for (
int i = 0; i < MyCnt; i++)
145 for (
int n = 0; n < extradata; n++)
153 ifs >>
r[n - NStructReal];
159 std::string msg(
"ParticleContainer::InitFromAsciiFile(");
160 msg += file; msg +=
") failed @ 2";
171 amrex::AllPrint() <<
"BAD PARTICLE ID WOULD BE " << ParticleType::NextID() <<
'\n'
172 <<
"BAD PARTICLE POS "
178 amrex::Abort(
"ParticleContainer::InitFromAsciiFile(): invalid particle");
183 p.id() = ParticleType::NextID();
187 if(nreals.
size() > extradata - NStructReal) {
188 for (
int n = NStructReal; n < extradata; n++)
190 nreals[n-NStructReal].push_back(
r[n-NStructReal]);
197 const Real DomSize[AMREX_SPACEDIM] =
201 int rep[AMREX_SPACEDIM];
203 #if AMREX_SPACEDIM > 2
204 for (rep[2] = 1; rep[2] <= lNrep[2]; rep[2]++)
207 #if AMREX_SPACEDIM > 1
208 for (rep[1] = 1; rep[1] <= lNrep[1]; rep[1]++)
211 for (rep[0] = 1; rep[0] <= lNrep[0]; rep[0]++)
213 if (!(
AMREX_D_TERM( (rep[0] == 1), && (rep[1] == 1), && (rep[2] == 1) ) ) )
216 for (
int d=0; d<AMREX_SPACEDIM; ++d)
218 p_rep.pos(d) =
static_cast<ParticleReal
>(p.pos(d) + Real(rep[d]-1)*DomSize[d]);
222 for (
int n = 0; n < extradata; n++)
226 p_rep.rdata(n) = p.rdata(n);
230 if (!Where(p_rep, pld))
232 PeriodicShift(p_rep);
233 if (!Where(p_rep, pld))
236 amrex::AllPrint() <<
"BAD REPLICATED PARTICLE ID WOULD BE " << ParticleType::NextID() <<
"\n";
238 amrex::Abort(
"ParticleContainer_impl<ParticleType, NArrayReal, NArrayInt>::InitFromAsciiFile(): invalid replicated particle");
242 p_rep.id() = ParticleType::NextID();
247 if (nreals.
size() > extradata - NStructReal) {
248 for (
int n = NStructReal; n < extradata; n++)
250 nreals[n-NStructReal].push_back(
r[n-NStructReal]);
257 #if AMREX_SPACEDIM > 1
260 #if AMREX_SPACEDIM > 2
269 int NRedist_chunk = NReaders / NRedist;
273 for (
int nr = 0; nr < NRedist; nr++)
276 host_particles.reserve(15);
277 host_particles.resize(finestLevel()+1);
280 host_real_attribs.reserve(15);
281 host_real_attribs.resize(finestLevel()+1);
285 << nr*NRedist_chunk <<
" to "
286 << (nr+1)*NRedist_chunk-1 <<
'\n';
289 for (
int which = nr*NRedist_chunk; which < (nr+1)*NRedist_chunk; which++)
293 while (!nparticles.
empty())
300 if (nreals.
size() > extradata - NStructReal && NArrayReal > 0)
302 for (
int n = NStructReal; n < extradata; n++)
304 Real rdata = nreals[n-NStructReal].back();
305 host_real_attribs[pld.
m_lev][std::make_pair(pld.
m_grid, pld.
m_tile)][n-NStructReal].push_back(rdata);
311 if (nreals.
size() > extradata - NStructReal)
313 for (
int n = NStructReal; n < extradata; n++)
315 nreals[n-NStructReal].pop_back();
322 for (
int lev = 0; lev < static_cast<int>(host_particles.
size()); ++lev)
324 for (
auto& kv : host_particles[lev])
327 auto grid = kv.first.first;
328 auto tile = kv.first.second;
329 const auto& src_tile = kv.second;
331 auto& dst_tile = GetParticles(lev)[std::make_pair(grid,tile)];
332 auto old_size = dst_tile.GetArrayOfStructs().
size();
333 auto new_size = old_size + src_tile.size();
334 dst_tile.resize(new_size);
337 dst_tile.GetArrayOfStructs().begin() + old_size);
339 if((host_real_attribs[lev][std::make_pair(grid, tile)]).
size() > (
long unsigned int) NArrayReal) {
340 for (
int i = 0; i < NArrayReal; ++i) {
342 host_real_attribs[lev][std::make_pair(grid,tile)][i].
begin(),
343 host_real_attribs[lev][std::make_pair(grid,tile)][i].
end(),
344 dst_tile.GetStructOfArrays().GetRealData(i).begin() + old_size);
356 if (NRedist*NRedist_chunk < NReaders) {
358 << NRedist*NRedist_chunk <<
" to "
362 for (
int which = NRedist*NRedist_chunk; which < NReaders; which++)
365 host_particles.reserve(15);
366 host_particles.resize(finestLevel()+1);
369 host_real_attribs.reserve(15);
370 host_real_attribs.resize(finestLevel()+1);
374 while (!nparticles.
empty())
380 if((host_real_attribs[pld.
m_lev][std::make_pair(pld.
m_grid, pld.
m_tile)]).size() > (
long unsigned int) (extradata - NStructReal)) {
381 for (
int n = NStructReal; n < extradata; n++)
383 Real rdata = nreals[n-NStructReal].back();
384 host_real_attribs[pld.
m_lev][std::make_pair(pld.
m_grid, pld.
m_tile)][n-NStructReal].push_back(rdata);
390 if (nreals.
size() > extradata - NStructReal) {
391 for (
int n = NStructReal; n < extradata; n++)
393 nreals[n-NStructReal].pop_back();
399 for (
int lev = 0; lev < static_cast<int>(host_particles.
size()); ++lev)
401 for (
auto& kv : host_particles[lev])
403 auto grid = kv.first.first;
404 auto tile = kv.first.second;
405 const auto& src_tile = kv.second;
407 auto& dst_tile = GetParticles(lev)[std::make_pair(grid,tile)];
408 auto old_size = dst_tile.GetArrayOfStructs().
size();
409 auto new_size = old_size + src_tile.size();
410 dst_tile.resize(new_size);
413 dst_tile.GetArrayOfStructs().begin() + old_size);
415 for (
int i = 0; i < NArrayReal; ++i) {
417 host_real_attribs[lev][std::make_pair(grid,tile)][i].
begin(),
418 host_real_attribs[lev][std::make_pair(grid,tile)][i].
end(),
419 dst_tile.GetStructOfArrays().GetRealData(i).begin() + old_size);
433 Long num_particles = how_many;
437 if (
AMREX_D_TERM(lNrep[0] == 1, && lNrep[1] == 1, && lNrep[2] == 1))
439 amrex::Print() <<
"Total number of particles: " << num_particles <<
'\n';
443 Long num_particles_read = how_many_read;
448 <<
AMREX_D_TERM(lNrep[0] <<
" ", << lNrep[1] <<
" ", << lNrep[2]) <<
"\n"
449 <<
"Total number of particles read in : " << num_particles_read <<
'\n'
450 <<
"Total number of particles after replication: " << num_particles <<
'\n';
464 amrex::Print() <<
"InitFromAsciiFile() time: " << runtime <<
'\n';
479 template <
typename ParticleType,
int NArrayReal,
int NArrayInt,
480 template<
class>
class Allocator,
class CellAssignor>
486 BL_PROFILE(
"ParticleContainer<NSR, NSI, NAR, NAI>::InitFromBinaryFile()");
499 const int NReaders = MaxReaders();
505 const Long NPartPerRedist = MaxParticlesPerRead();
509 amrex::Print() <<
"Reading with " << NReaders <<
" readers\n"
510 <<
"Redistributing after every " << NPartPerRedist <<
" particles for each reader\n";
518 tmp_particles.reserve(15);
519 tmp_particles.resize(finestLevel()+1);
537 std::set<int> readers;
548 for (
int i = 0; i <
NProcs; i++) {
575 while (
static_cast<int>(readers.size()) < NReaders);
581 for (
auto it = readers.cbegin(), End = readers.cend();
597 readers.insert(rprocs.begin(), rprocs.end());
604 int RealSizeInFile = 0;
606 if (readers.find(
MyProc) != readers.end())
610 ifs.rdbuf()->pubsetbuf(io_buffer.
dataPtr(), io_buffer.
size());
612 ifs.open(file.c_str(), std::ios::in|std::ios::binary);
618 ifs.read((
char*)&NP,
sizeof(NP));
619 ifs.read((
char*)&DM,
sizeof(DM));
620 ifs.read((
char*)&NX,
sizeof(NX));
625 amrex::Abort(
"ParticleContainer_impl<ParticleType, NArrayReal, NArrayInt>::InitFromBinaryFile(): NP <= 0");
630 if (DM != AMREX_SPACEDIM) {
631 amrex::Abort(
"ParticleContainer_impl<ParticleType, NArrayReal, NArrayInt>::InitFromBinaryFile(): DM != AMREX_SPACEDIM");
636 if (NX < 0 || NX > NStructReal) {
637 amrex::Abort(
"ParticleContainer_impl<ParticleType, NArrayReal, NArrayInt>::InitFromBinaryFile(): NX < 0 || NX > N");
642 if (extradata > NX) {
643 amrex::Abort(
"ParticleContainer_impl<ParticleType, NArrayReal, NArrayInt>::InitFromBinaryFile(): extradata > NX");
650 const std::streamoff CURPOS = ifs.tellg();
658 const std::streamoff ENDPOS = ifs.tellg();
660 RealSizeInFile =
int((ENDPOS - CURPOS) / (NP*(DM+NX)));
662 AMREX_ASSERT(RealSizeInFile ==
sizeof(
float) || RealSizeInFile ==
sizeof(
double));
666 ifs.seekg(CURPOS, std::ios::beg);
671 for ( ;
id < NReaders;
id++) {
672 if (rprocs[
id] ==
MyProc) {
679 const std::streamoff NSKIP =
id * (NP/NReaders) * (DM+NX) * RealSizeInFile;
683 ifs.seekg(NSKIP, std::ios::cur);
688 std::string msg(
"ParticleContainer::InitFromBinaryFile(");
690 msg +=
") failed @ 1";
701 Long MyCnt = NP / NReaders;
703 if (
MyProc == rprocs[0]) {
707 MyCnt += NP % NReaders;
710 Long how_many_redists = NP / (NPartPerRedist*NReaders), how_many_read = 0;
712 if (NP % (NPartPerRedist*NReaders)) { how_many_redists++; }
719 fxtra.resize(extradata);
720 dxtra.resize(extradata);
723 if ((NX-extradata) > 0)
725 fignore.resize(NX-extradata);
726 dignore.resize(NX-extradata);
731 for (
int j = 0; j < how_many_redists; j++)
735 host_particles.reserve(15);
736 host_particles.resize(finestLevel()+1);
738 if (readers.find(
MyProc) != readers.end())
744 const Long NRead =
std::min((MyCnt-how_many_read), NPartPerRedist);
746 for (Long i = 0; i < NRead; i++)
752 if (RealSizeInFile ==
sizeof(
float))
754 float fpos[AMREX_SPACEDIM];
756 ifs.read((
char*)&fpos[0], AMREX_SPACEDIM*
sizeof(
float));
760 p.pos(2) = fpos[2];);
763 else if (RealSizeInFile ==
sizeof(
double))
765 double dpos[AMREX_SPACEDIM];
767 ifs.read((
char*)&dpos[0], AMREX_SPACEDIM*
sizeof(
double));
769 AMREX_D_TERM(p.pos(0) =
static_cast<ParticleReal
>(dpos[0]);,
770 p.pos(1) =
static_cast<ParticleReal
>(dpos[1]);,
771 p.pos(2) =
static_cast<ParticleReal
>(dpos[2]););
779 if (RealSizeInFile ==
sizeof(
float))
781 ifs.read((
char*)fxtra.data(), std::streamsize(extradata*
sizeof(
float)));
783 for (
int ii = 0; ii < extradata; ii++) {
784 p.rdata(ii) =
static_cast<ParticleReal
>(fxtra[ii]);
787 else if (RealSizeInFile ==
sizeof(
double))
789 ifs.read((
char*)dxtra.data(), std::streamsize(extradata*
sizeof(
double)));
791 for (
int ii = 0; ii < extradata; ii++) {
792 p.rdata(ii) =
static_cast<ParticleReal
>(dxtra[ii]);
799 if ((NX-extradata) > 0)
801 if (RealSizeInFile ==
sizeof(
float))
803 ifs.read((
char*)fignore.data(), std::streamsize((NX-extradata)*
sizeof(
float)));
805 else if (RealSizeInFile ==
sizeof(
double))
807 ifs.read((
char*)dignore.data(), std::streamsize((NX-extradata)*
sizeof(
double)));
813 std::string msg(
"ParticleContainer::InitFromBinaryFile(");
815 msg +=
") failed @ 2";
826 amrex::AllPrint() <<
"BAD PARTICLE ID WOULD BE " << ParticleType::NextID() <<
'\n'
827 <<
"BAD PARTICLE POS "
833 amrex::Abort(
"ParticleContainer_impl<ParticleType, NArrayReal, NArrayInt>::InitFromBinaryFile(): invalid particle");
837 p.id() = ParticleType::NextID();
843 how_many_read += NRead;
846 for (
int host_lev = 0; host_lev < static_cast<int>(host_particles.
size()); ++host_lev)
848 for (
auto& kv : host_particles[host_lev]) {
849 auto grid = kv.first.first;
850 auto tile = kv.first.second;
851 const auto& src_tile = kv.second;
853 auto& dst_tile = GetParticles(host_lev)[std::make_pair(grid,tile)];
854 auto old_size = dst_tile.GetArrayOfStructs().
size();
855 auto new_size = old_size + src_tile.size();
856 dst_tile.resize(new_size);
859 dst_tile.GetArrayOfStructs().begin() + old_size);
871 for (
int lev = 0; lev < m_particles.size(); lev++)
873 auto& pmap = m_particles[lev];
874 auto& tmp_pmap = tmp_particles[lev];
876 for (
auto& kv : pmap) {
877 auto& aos = kv.second.GetArrayOfStructs()();
878 auto& tmp_aos = tmp_pmap[kv.first].GetArrayOfStructs()();
880 tmp_aos.insert(tmp_aos.end(), aos.begin(), aos.end());
890 tmp_particles.swap(m_particles);
896 Long num_particles_read = how_many_read;
900 amrex::Print() <<
"\nTotal number of particles: " << num_particles_read <<
'\n';
913 amrex::Print() <<
"InitFromBinaryFile() time: " << runtime <<
'\n';
925 template <
typename ParticleType,
int NArrayReal,
int NArrayInt,
926 template<
class>
class Allocator,
class CellAssignor>
932 BL_PROFILE(
"ParticleContainer<NSR, NSI, NAR, NAI>::InitFromBinaryMetaFile()");
935 std::ifstream ifs(metafile.c_str(), std::ios::in);
941 std::getline(ifs,file);
943 if (!ifs.good()) {
break; }
946 amrex::Print() <<
"InitFromBinaryMetaFile: processing file: " << file <<
'\n';
949 InitFromBinaryFile(file, extradata);
960 amrex::Print() <<
"InitFromBinaryMetaFile() time: " << runtime <<
'\n';
964 template <
typename ParticleType,
int NArrayReal,
int NArrayInt,
965 template<
class>
class Allocator,
class CellAssignor>
974 BL_PROFILE(
"ParticleContainer<NSR, NSI, NAR, NAI>::InitRandom()");
992 if (!containing_bx.
ok()) { containing_bx = geom.
ProbDomain(); }
1001 const Real* xlo = containing_bx.
lo();
1002 const Real* xhi = containing_bx.
hi();
1011 "InitRandom has serialize=true, but this would cause too much "
1012 "particle data to be sent from IOProc. Set serialize=false, "
1013 "or use fewer than " +
1017 amrex::Real(AMREX_SPACEDIM)
1033 for (Long j = 0; j < icount; j++)
1035 for (
int i = 0; i < AMREX_SPACEDIM; i++)
1042 while (
static_cast<ParticleReal
>(
x) <
static_cast<ParticleReal
>(xlo[i]) ||
static_cast<ParticleReal
>(
x) >=
static_cast<ParticleReal
>(xhi[i]));
1044 pos[j*AMREX_SPACEDIM + i] =
static_cast<ParticleReal
>(
x);
1054 host_particles.reserve(15);
1055 host_particles.resize(finestLevel()+1);
1058 host_real_attribs.reserve(15);
1059 host_real_attribs.resize(finestLevel()+1);
1062 host_int_attribs.reserve(15);
1063 host_int_attribs.resize(finestLevel()+1);
1066 host_idcpu.reserve(15);
1067 host_idcpu.resize(finestLevel()+1);
1069 for (Long j = 0; j < icount; j++)
1073 for (
int i = 0; i < AMREX_SPACEDIM; i++) {
1074 ptest.
pos(i) = pos[j*AMREX_SPACEDIM + i];
1077 if (!Where(ptest, pld)) {
1078 amrex::Abort(
"ParticleContainer::InitRandom(): invalid particle");
1084 const int who = ParticleDistributionMap(pld.
m_lev)[pld.
m_grid];
1088 if constexpr(!ParticleType::is_soa_particle)
1091 for (
int i = 0; i < AMREX_SPACEDIM; i++) {
1092 p.pos(i) = pos[j*AMREX_SPACEDIM + i];
1096 p.id() = ParticleType::NextID();
1099 for (
int i = 0; i < NStructInt; i++) {
1103 for (
int i = 0; i < NStructReal; i++) {
1108 host_particles[pld.
m_lev][ind].push_back(p);
1111 for (
int i = 0; i < NArrayReal; i++) {
1112 host_real_attribs[pld.
m_lev][ind][i].push_back(
static_cast<ParticleReal
>(pdata.
real_array_data[i]));
1116 for (
int i = 0; i < NArrayInt; i++) {
1120 for (
int i = 0; i < AMREX_SPACEDIM; i++) {
1121 host_real_attribs[pld.
m_lev][ind][i].push_back(pos[j*AMREX_SPACEDIM+i]);
1124 host_idcpu[pld.
m_lev][ind].push_back(0);
1128 host_particles[pld.
m_lev][ind];
1131 for (
int i = AMREX_SPACEDIM; i < NArrayReal; i++) {
1132 host_real_attribs[pld.
m_lev][ind][i].push_back(
static_cast<ParticleReal
>(pdata.
real_array_data[i]));
1136 for (
int i = 2; i < NArrayInt; i++) {
1143 for (
int host_lev = 0; host_lev < static_cast<int>(host_particles.
size()); ++host_lev)
1145 for (
auto& kv : host_particles[host_lev]) {
1146 auto grid = kv.first.first;
1147 auto tile = kv.first.second;
1148 const auto& src_tile = kv.second;
1150 auto& dst_tile = GetParticles(host_lev)[std::make_pair(grid,tile)];
1151 auto old_size = dst_tile.
size();
1152 auto new_size = old_size;
1153 if constexpr(!ParticleType::is_soa_particle)
1155 new_size += src_tile.size();
1157 new_size += host_real_attribs[host_lev][std::make_pair(grid,tile)][0].
size();
1159 dst_tile.resize(new_size);
1161 if constexpr(!ParticleType::is_soa_particle)
1164 dst_tile.GetArrayOfStructs().begin() + old_size);
1167 host_idcpu[host_lev][std::make_pair(grid,tile)].
begin(),
1168 host_idcpu[host_lev][std::make_pair(grid,tile)].
end(),
1169 dst_tile.GetStructOfArrays().GetIdCPUData().begin() + old_size);
1172 for (
int i = 0; i < NArrayReal; ++i) {
1174 host_real_attribs[host_lev][std::make_pair(grid,tile)][i].
begin(),
1175 host_real_attribs[host_lev][std::make_pair(grid,tile)][i].
end(),
1176 dst_tile.GetStructOfArrays().GetRealData(i).begin() + old_size);
1179 for (
int i = 0; i < NArrayInt; ++i) {
1181 host_int_attribs[host_lev][std::make_pair(grid,tile)][i].
begin(),
1182 host_int_attribs[host_lev][std::make_pair(grid,tile)][i].
end(),
1183 dst_tile.GetStructOfArrays().GetIntData(i).begin() + old_size);
1203 host_particles.reserve(15);
1204 host_particles.resize(finestLevel()+1);
1207 host_real_attribs.reserve(15);
1208 host_real_attribs.resize(finestLevel()+1);
1211 host_int_attribs.reserve(15);
1212 host_int_attribs.resize(finestLevel()+1);
1215 host_idcpu.reserve(15);
1216 host_idcpu.resize(finestLevel()+1);
1218 for (Long icnt = 0; icnt <
M; icnt++) {
1220 for (
int i = 0; i < AMREX_SPACEDIM; i++) {
1225 while (
static_cast<ParticleReal
>(
x) <
static_cast<ParticleReal
>(xlo[i]) ||
static_cast<ParticleReal
>(
x) >=
static_cast<ParticleReal
>(xhi[i]));
1227 ptest.
pos(i) =
static_cast<ParticleReal
>(
x);
1232 ptest.
id() = ParticleType::NextID();
1236 if (!Where(ptest, pld))
1238 amrex::Abort(
"ParticleContainer_impl<ParticleType, NArrayReal, NArrayInt>::InitRandom(): invalid particle");
1243 if constexpr(!ParticleType::is_soa_particle)
1246 for (
int i = 0; i < AMREX_SPACEDIM; i++) {
1247 p.pos(i) = ptest.
pos(i);;
1250 p.id() = ptest.
id();
1251 p.cpu() = ptest.
cpu();
1253 for (
int i = 0; i < NStructReal; i++) {
1257 for (
int i = 0; i < NStructInt; i++) {
1262 host_particles[pld.
m_lev][ind].push_back(p);
1265 for (
int i = 0; i < NArrayReal; i++) {
1266 host_real_attribs[pld.
m_lev][ind][i].push_back(
static_cast<ParticleReal
>(pdata.
real_array_data[i]));
1270 for (
int i = 0; i < NArrayInt; i++) {
1274 for (
int i = 0; i < AMREX_SPACEDIM; i++) {
1275 host_real_attribs[pld.
m_lev][ind][i].push_back(ptest.
pos(i));
1278 host_idcpu[pld.
m_lev][ind].push_back(0);
1282 host_particles[pld.
m_lev][ind];
1285 for (
int i = AMREX_SPACEDIM; i < NArrayReal; i++) {
1286 host_real_attribs[pld.
m_lev][ind][i].push_back(
static_cast<ParticleReal
>(pdata.
real_array_data[i]));
1290 for (
int i = 2; i < NArrayInt; i++) {
1296 for (
int host_lev = 0; host_lev < static_cast<int>(host_particles.
size()); ++host_lev)
1298 for (
auto& kv : host_particles[host_lev]) {
1299 auto grid = kv.first.first;
1300 auto tile = kv.first.second;
1301 const auto& src_tile = kv.second;
1303 auto& dst_tile = GetParticles(host_lev)[std::make_pair(grid,tile)];
1304 auto old_size = dst_tile.
size();
1305 auto new_size = old_size;
1306 if constexpr(!ParticleType::is_soa_particle)
1308 new_size += src_tile.size();
1310 new_size += host_real_attribs[host_lev][std::make_pair(grid,tile)][0].
size();
1312 dst_tile.resize(new_size);
1314 if constexpr(!ParticleType::is_soa_particle)
1317 dst_tile.GetArrayOfStructs().begin() + old_size);
1320 host_idcpu[host_lev][std::make_pair(grid,tile)].
begin(),
1321 host_idcpu[host_lev][std::make_pair(grid,tile)].
end(),
1322 dst_tile.GetStructOfArrays().GetIdCPUData().begin() + old_size);
1325 for (
int i = 0; i < NArrayReal; ++i) {
1327 host_real_attribs[host_lev][std::make_pair(grid,tile)][i].
begin(),
1328 host_real_attribs[host_lev][std::make_pair(grid,tile)][i].
end(),
1329 dst_tile.GetStructOfArrays().GetRealData(i).begin() + old_size);
1332 for (
int i = 0; i < NArrayInt; ++i) {
1334 host_int_attribs[host_lev][std::make_pair(grid,tile)][i].
begin(),
1335 host_int_attribs[host_lev][std::make_pair(grid,tile)][i].
end(),
1336 dst_tile.GetStructOfArrays().GetIntData(i).begin() + old_size);
1352 amrex::Print() <<
"ParticleContainer_impl<ParticleType, NArrayReal, NArrayInt>::InitRandom() time: " << stoptime <<
'\n';
1358 template <
typename ParticleType,
int NArrayReal,
int NArrayInt,
1359 template<
class>
class Allocator,
class CellAssignor>
1366 BL_PROFILE(
"ParticleContainer<NSR, NSI, NAR, NAI>::InitRandomPerBox()");
1383 std::mt19937 mt(iseed);
1384 std::uniform_real_distribution<double> dist(0.0, 1.0);
1386 m_particles.resize(m_gdb->finestLevel()+1);
1388 for (
int lev = 0; lev < m_particles.size(); lev++)
1394 for (
MFIter mfi(*m_dummy_mf[0],
false); mfi.
isValid(); ++mfi)
1396 Box grid = m_gdb->ParticleBoxArray(0)[mfi.index()];
1399 for (Long icnt = 0; icnt < icount_per_box; icnt++) {
1400 for (Long jcnt = 0; jcnt < icount_per_box; jcnt++) {
1401 for (Long kcnt = 0; kcnt < icount_per_box; kcnt++)
1404 p.pos(0) =
static_cast<ParticleReal
>(grid_box.
lo(0) + (dist(mt) +
double(icnt)) /
double(icount_per_box) * grid_box.
length(0));,
1405 p.pos(1) =
static_cast<ParticleReal
>(grid_box.
lo(1) + (dist(mt) +
double(jcnt)) /
double(icount_per_box) * grid_box.
length(1));,
1406 p.pos(2) =
static_cast<ParticleReal
>(grid_box.
lo(2) + (dist(mt) +
double(kcnt)) /
double(icount_per_box) * grid_box.
length(2));
1409 for (
int i = 0; i < AMREX_SPACEDIM; i++) {
1414 for (
int i = 0; i < NStructReal; i++) {
1419 p.id() = ParticleType::NextID();
1422 for (
int i = 0; i < NStructInt; i++) {
1427 if (!Where(p, pld)) {
1428 amrex::Abort(
"ParticleContainer::InitRandomPerBox(): invalid particle");
1434 m_particles[pld.
m_lev][ind].push_back(p);
1437 for (
int i = 0; i < NArrayReal; i++) {
1438 m_particles[pld.
m_lev][ind].push_back_real(i,
static_cast<ParticleReal
>(pdata.
real_array_data[i]));
1442 for (
int i = 0; i < NArrayInt; i++) {
1455 amrex::Print() <<
"ParticleContainer_impl<ParticleType, NArrayReal, NArrayInt>::InitRandomPerBox() time: " << stoptime <<
'\n';
1459 template <
typename ParticleType,
int NArrayReal,
int NArrayInt,
1460 template<
class>
class Allocator,
class CellAssignor>
1467 BL_PROFILE(
"ParticleContainer<NSR, NSI, NAR, NAI>::InitOnePerCell()");
1473 AMREX_ASSERT(x_off >= 0. && y_off >= 0. && z_off >= 0.);
1474 AMREX_ASSERT(x_off <= 1. && y_off <= 1. && z_off <= 1.);
1485 for (
MFIter mfi(*m_dummy_mf[0],
false); mfi.
isValid(); ++mfi) {
1486 Box grid = ParticleBoxArray(0)[mfi.index()];
1487 auto ind = std::make_pair(mfi.index(), mfi.LocalTileIndex());
1493 AMREX_D_TERM(p.pos(0) =
static_cast<ParticleReal
>(grid_box.
lo(0) + (x_off + cell[0]-beg[0])*dx[0]);,
1494 p.pos(1) =
static_cast<ParticleReal
>(grid_box.
lo(1) + (y_off + cell[1]-beg[1])*dx[1]);,
1495 p.pos(2) =
static_cast<ParticleReal
>(grid_box.
lo(2) + (z_off + cell[2]-beg[2])*dx[2]););
1497 for (
int d = 0; d < AMREX_SPACEDIM; ++d) {
1501 for (
int i = 0; i < NStructReal; i++) {
1506 p.id() = ParticleType::NextID();
1509 for (
int i = 0; i < NStructInt; i++) {
1517 for (
int i = 0; i < NArrayReal; i++) {
1522 for (
int i = 0; i < NArrayInt; i++) {
1534 if (m_verbose > 1) {
1539 amrex::Print() <<
"ParticleContainer_impl<ParticleType, NArrayReal, NArrayInt>::InitOnePerCell() time: " << stoptime <<
'\n';
1543 template <
typename ParticleType,
int NArrayReal,
int NArrayInt,
1544 template<
class>
class Allocator,
class CellAssignor>
1549 BL_PROFILE(
"ParticleContainer<NSR, NSI, NAR, NAI>::InitNRandomPerCell()");
1562 for (
int lev = 0; lev < m_particles.size(); lev++) {
1571 for (
MFIter mfi(*m_dummy_mf[0],
false); mfi.
isValid(); ++mfi)
1573 Box grid = ParticleBoxArray(0)[mfi.index()];
1577 host_particles.reserve(15);
1578 host_particles.resize(finestLevel()+1);
1581 host_real_attribs.reserve(15);
1582 host_real_attribs.resize(finestLevel()+1);
1585 host_int_attribs.reserve(15);
1586 host_int_attribs.resize(finestLevel()+1);
1591 for (
int n = 0; n < n_per_cell; n++)
1594 for (
int i = 0; i < AMREX_SPACEDIM; i++) {
1595 constexpr
int max_iter = 10;
1597 while (iter < max_iter) {
1599 p.pos(i) =
static_cast<ParticleReal
>(grid_box.
lo(i) + (
r + Real(cell[i]-beg[i]))*dx[i]);
1600 if (p.pos(i) < grid_box.
hi(i)) {
break; }
1606 for (
int i = 0; i < NStructReal; i++) {
1611 p.id() = ParticleType::NextID();
1614 for (
int i = 0; i < NStructInt; i++) {
1619 if (!Where(p, pld)) {
1620 amrex::Abort(
"ParticleContainer_impl<ParticleType, NArrayReal, NArrayInt>::InitNRandomPerCell(): invalid particle");
1626 host_particles[pld.
m_lev][ind].push_back(p);
1629 for (
int i = 0; i < NArrayReal; i++) {
1634 for (
int i = 0; i < NArrayInt; i++) {
1640 for (
int host_lev = 0; host_lev < static_cast<int>(host_particles.
size()); ++host_lev)
1642 for (
auto& kv : host_particles[host_lev]) {
1643 auto gid = kv.first.first;
1644 auto tid = kv.first.second;
1645 const auto& src_tid = kv.second;
1647 auto& dst_tile = GetParticles(host_lev)[std::make_pair(gid,tid)];
1648 auto old_size = dst_tile.GetArrayOfStructs().
size();
1649 auto new_size = old_size + src_tid.size();
1650 dst_tile.resize(new_size);
1653 dst_tile.GetArrayOfStructs().begin() + old_size);
1655 for (
int i = 0; i < NArrayReal; ++i)
1658 host_real_attribs[host_lev][std::make_pair(gid,tid)][i].
begin(),
1659 host_real_attribs[host_lev][std::make_pair(gid,tid)][i].
end(),
1660 dst_tile.GetStructOfArrays().GetRealData(i).begin() + old_size);
1663 for (
int i = 0; i < NArrayInt; ++i)
1666 host_int_attribs[host_lev][std::make_pair(gid,tid)][i].
begin(),
1667 host_int_attribs[host_lev][std::make_pair(gid,tid)][i].
end(),
1668 dst_tile.GetStructOfArrays().GetIntData(i).begin() + old_size);
1681 amrex::Print() <<
"ParticleContainer_impl<ParticleType, NArrayReal, NArrayInt>::InitNRandomPerCell() time: " << stoptime <<
'\n';
#define BL_PROFILE(a)
Definition: AMReX_BLProfiler.H:551
#define AMREX_ASSERT(EX)
Definition: AMReX_BLassert.H:38
#define AMREX_D_TERM(a, b, c)
Definition: AMReX_SPACE.H:129
#define AMREX_D_DECL(a, b, c)
Definition: AMReX_SPACE.H:104
Print on all processors of the default communicator.
Definition: AMReX_Print.H:117
AMREX_GPU_HOST_DEVICE const IntVectND< dim > & smallEnd() const &noexcept
Get the smallend of the BoxND.
Definition: AMReX_Box.H:105
AMREX_GPU_HOST_DEVICE void next(IntVectND< dim > &) const noexcept
Step through the rectangle. It is a runtime error to give a point not inside rectangle....
Definition: AMReX_Box.H:1059
AMREX_GPU_HOST_DEVICE const IntVectND< dim > & bigEnd() const &noexcept
Get the bigend.
Definition: AMReX_Box.H:116
const Real * CellSize() const noexcept
Returns the cellsize for each coordinate direction.
Definition: AMReX_CoordSys.H:71
Rectangular problem domain geometry.
Definition: AMReX_Geometry.H:73
const Real * ProbLo() const noexcept
Returns the lo end of the problem domain in each dimension.
Definition: AMReX_Geometry.H:178
const RealBox & ProbDomain() const noexcept
Returns the problem domain.
Definition: AMReX_Geometry.H:170
const Real * ProbHi() const noexcept
Returns the hi end of the problem domain in each dimension.
Definition: AMReX_Geometry.H:180
Real ProbLength(int dir) const noexcept
Returns length of problem domain in specified dimension.
Definition: AMReX_Geometry.H:208
Definition: AMReX_MFIter.H:57
bool isValid() const noexcept
Is the iterator valid i.e. is it associated with a FAB?
Definition: AMReX_MFIter.H:141
Definition: AMReX_PODVector.H:246
T & back() noexcept
Definition: AMReX_PODVector.H:589
void pop_back() noexcept
Definition: AMReX_PODVector.H:571
bool empty() const noexcept
Definition: AMReX_PODVector.H:579
void push_back(const T &a_value)
Definition: AMReX_PODVector.H:556
A distributed container for Particles sorted onto the levels, grids, and tiles of a block-structured ...
Definition: AMReX_ParticleContainer.H:145
std::map< std::pair< int, int >, ParticleTileType > ParticleLevel
Definition: AMReX_ParticleContainer.H:186
typename AoS::ParticleVector ParticleVector
Definition: AMReX_ParticleContainer.H:192
T_ParticleType ParticleType
Definition: AMReX_ParticleContainer.H:147
This class provides the user with a few print options.
Definition: AMReX_Print.H:35
A Box with real dimensions. A RealBox is OK iff volume >= 0.
Definition: AMReX_RealBox.H:21
AMREX_GPU_HOST_DEVICE bool ok() const noexcept
Is the RealBox OK; i.e. does it have non-negative volume?
Definition: AMReX_RealBox.H:81
void setLo(const Real *a_lo) noexcept
Sets lo side.
Definition: AMReX_RealBox.H:64
AMREX_GPU_HOST_DEVICE bool contains(const Real *point, Real eps=0.0) const noexcept
Is the specified point contained in the RealBox?
Definition: AMReX_RealBox.H:102
AMREX_GPU_HOST_DEVICE const Real * hi() const &noexcept
Returns hide side.
Definition: AMReX_RealBox.H:51
AMREX_GPU_HOST_DEVICE const Real * lo() const &noexcept
Returns lo side.
Definition: AMReX_RealBox.H:46
void setHi(const Real *a_hi) noexcept
Sets hi side.
Definition: AMReX_RealBox.H:72
AMREX_GPU_HOST_DEVICE Real length(int dir) const noexcept
Returns length in specified direction.
Definition: AMReX_RealBox.H:62
This class is a thin wrapper around std::vector. Unlike vector, Vector::operator[] provides bound che...
Definition: AMReX_Vector.H:27
Long size() const noexcept
Definition: AMReX_Vector.H:50
T * dataPtr() noexcept
get access to the underlying data pointer
Definition: AMReX_Vector.H:46
AMREX_GPU_HOST_DEVICE Long size(T const &b) noexcept
integer version
Definition: AMReX_GpuRange.H:26
void copyAsync(HostToDevice, InIter begin, InIter end, OutIter result) noexcept
A host-to-device copy routine. Note this is just a wrapper around memcpy, so it assumes contiguous st...
Definition: AMReX_GpuContainers.H:233
static constexpr HostToDevice hostToDevice
Definition: AMReX_GpuContainers.H:98
void streamSynchronize() noexcept
Definition: AMReX_GpuDevice.H:237
int NProcs()
Process ID in MPI_COMM_WORLD.
Definition: AMReX_MPMD.cpp:122
int MyProc()
Definition: AMReX_MPMD.cpp:117
void ReduceLongMax(Long &)
Long max reduction.
Definition: AMReX_ParallelDescriptor.cpp:1224
void ReduceLongSum(Long &)
Long sum reduction.
Definition: AMReX_ParallelDescriptor.cpp:1223
void Bcast(void *, int, MPI_Datatype, int, MPI_Comm)
Definition: AMReX_ParallelDescriptor.cpp:1282
int IOProcessorNumber() noexcept
Definition: AMReX_ParallelDescriptor.H:266
bool IOProcessor() noexcept
Is this CPU the I/O Processor? To get the rank number, call IOProcessorNumber()
Definition: AMReX_ParallelDescriptor.H:275
void ReduceRealMax(Vector< std::reference_wrapper< Real > > const &)
Definition: AMReX_ParallelDescriptor.cpp:1215
@ min
Definition: AMReX_ParallelReduce.H:18
@ max
Definition: AMReX_ParallelReduce.H:17
static constexpr int M
Definition: AMReX_OpenBC.H:13
void FileOpenFailed(const std::string &file)
Output a message and abort when couldn't open the file.
Definition: AMReX_Utility.cpp:131
void InitRandom(ULong cpu_seed, int nprocs, ULong gpu_seed)
Set the seed of the random number generator.
Definition: AMReX_Random.cpp:89
void copyParticles(DstTile &dst, const SrcTile &src) noexcept
Copy particles from src to dst. This version copies all the particles, writing them to the beginning ...
Definition: AMReX_ParticleTransformation.H:158
Real Random()
Generate a psuedo-random double from uniform distribution.
Definition: AMReX_Random.cpp:123
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 end(BoxND< dim > const &box) noexcept
Definition: AMReX_Box.H:1890
double second() noexcept
Definition: AMReX_Utility.cpp:922
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void ignore_unused(const Ts &...)
This shuts up the compiler about unused variables.
Definition: AMReX.H:111
void Error(const std::string &msg)
Print out message to cerr and exit via amrex::Abort().
Definition: AMReX.cpp:219
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 begin(BoxND< dim > const &box) noexcept
Definition: AMReX_Box.H:1881
void Abort(const std::string &msg)
Print out message to cerr and exit via abort().
Definition: AMReX.cpp:225
const int[]
Definition: AMReX_BLProfiler.cpp:1664
Definition: AMReX_Particle.H:140
Definition: AMReX_Particle.H:37
A struct used to pass initial data into the various Init methods. This struct is used to pass initial...
Definition: AMReX_ParticleContainer.H:116
std::array< int, NStructInt > int_struct_data
Definition: AMReX_ParticleContainer.H:118
std::array< int, NArrayInt > int_array_data
Definition: AMReX_ParticleContainer.H:120
std::array< double, NArrayReal > real_array_data
Definition: AMReX_ParticleContainer.H:119
std::array< double, NStructReal > real_struct_data
Definition: AMReX_ParticleContainer.H:117
A struct used for storing a particle's position in the AMR hierarchy.
Definition: AMReX_ParticleContainer.H:91
int m_grid
Definition: AMReX_ParticleContainer.H:93
int m_tile
Definition: AMReX_ParticleContainer.H:94
int m_lev
Definition: AMReX_ParticleContainer.H:92
Definition: AMReX_ParticleTile.H:693
void push_back(const ParticleType &p)
Definition: AMReX_ParticleTile.H:914
int numParticles() const
Returns the number of real particles (excluding neighbors)
Definition: AMReX_ParticleTile.H:836
void push_back_int(int comp, int v)
Definition: AMReX_ParticleTile.H:1010
void push_back_real(int comp, ParticleReal v)
Definition: AMReX_ParticleTile.H:958
The struct used to store particles.
Definition: AMReX_Particle.H:295
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ParticleIDWrapper id() &
Definition: AMReX_Particle.H:315
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE RealVect pos() const &
Definition: AMReX_Particle.H:338
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ParticleCPUWrapper cpu() &
Definition: AMReX_Particle.H:312