Athena++/Atmosphere
Planetary Atmosphere Simulator
meshblock.cpp
Go to the documentation of this file.
1 //========================================================================================
2 // Athena++ astrophysical MHD code
3 // Copyright(C) 2014 James M. Stone <jmstone@princeton.edu> and other code contributors
4 // Licensed under the 3-clause BSD License, see LICENSE file for details
5 //========================================================================================
7 // \brief implementation of functions in MeshBlock class
8 
9 // C headers
10 
11 // C++ headers
12 #include <algorithm> // sort()
13 #include <cstdlib>
14 #include <cstring> // memcpy()
15 #include <ctime> // clock(), CLOCKS_PER_SEC, clock_t
16 #include <iomanip>
17 #include <iostream>
18 #include <sstream>
19 #include <stdexcept> // runtime_error
20 #include <string> // c_str()
21 
22 // Athena++ headers
23 #include "../athena.hpp"
24 #include "../athena_arrays.hpp"
25 #include "../bvals/bvals.hpp"
26 #include "../coordinates/coordinates.hpp"
27 #include "../eos/eos.hpp"
28 #include "../fft/athena_fft.hpp"
29 #include "../field/field.hpp"
30 #include "../globals.hpp"
31 #include "../gravity/gravity.hpp"
32 #include "../gravity/mg_gravity.hpp"
33 #include "../hydro/hydro.hpp"
34 #include "../parameter_input.hpp"
35 #include "../reconstruct/reconstruction.hpp"
36 #include "../scalars/scalars.hpp"
37 #include "../utils/buffer_utils.hpp"
38 #include "mesh.hpp"
39 #include "mesh_refinement.hpp"
40 #include "meshblock_tree.hpp"
41 #include "../thermodynamics/thermodynamics.hpp"
42 #include "../chemistry/chemistry.hpp"
43 #include "../radiation/radiation.hpp"
44 #include "../physics/physics.hpp"
45 #include "../diagnostics/diagnostics.hpp"
46 #include "../debugger/debugger.hpp"
47 
48 //----------------------------------------------------------------------------------------
49 // MeshBlock constructor: constructs coordinate, boundary condition, hydro, field
50 // and mesh refinement objects.
51 
52 MeshBlock::MeshBlock(int igid, int ilid, LogicalLocation iloc, RegionSize input_block,
53  BoundaryFlag *input_bcs, Mesh *pm, ParameterInput *pin,
54  int igflag, bool ref_flag) :
55  pmy_mesh(pm), loc(iloc), block_size(input_block),
56  gid(igid), lid(ilid), gflag(igflag), nuser_out_var(), prev(nullptr), next(nullptr),
57  new_block_dt_{}, new_block_dt_hyperbolic_{}, new_block_dt_parabolic_{},
58  new_block_dt_user_{},
59  nreal_user_meshblock_data_(), nint_user_meshblock_data_(), cost_(1.0) {
60  // initialize grid indices
61  is = NGHOST;
62  ie = is + block_size.nx1 - 1;
63 
64  ncells1 = block_size.nx1 + 2*NGHOST;
65  ncc1 = block_size.nx1/2 + 2*NGHOST;
66  if (pmy_mesh->f2) {
67  js = NGHOST;
68  je = js + block_size.nx2 - 1;
69  ncells2 = block_size.nx2 + 2*NGHOST;
70  ncc2 = block_size.nx2/2 + 2*NGHOST;
71  } else {
72  js = je = 0;
73  ncells2 = 1;
74  ncc2 = 1;
75  }
76 
77  if (pmy_mesh->f3) {
78  ks = NGHOST;
79  ke = ks + block_size.nx3 - 1;
80  ncells3 = block_size.nx3 + 2*NGHOST;
81  ncc3 = block_size.nx3/2 + 2*NGHOST;
82  } else {
83  ks = ke = 0;
84  ncells3 = 1;
85  ncc3 = 1;
86  }
87 
88  if (pm->multilevel) {
89  cnghost = (NGHOST + 1)/2 + 1;
90  cis = NGHOST; cie = cis + block_size.nx1/2 - 1;
91  cjs = cje = cks = cke = 0;
92  if (pmy_mesh->f2) // 2D or 3D
93  cjs = NGHOST, cje = cjs + block_size.nx2/2 - 1;
94  if (pmy_mesh->f3) // 3D
95  cks = NGHOST, cke = cks + block_size.nx3/2 - 1;
96  }
97 
98  // (probably don't need to preallocate space for references in these vectors)
99  vars_cc_.reserve(3);
100  vars_fc_.reserve(3);
101 
102  // construct objects stored in MeshBlock class. Note in particular that the initial
103  // conditions for the simulation are set in problem generator called from main, not
104  // in the Hydro constructor
105 
106  // mesh-related objects
107  // Boundary
108  pbval = new BoundaryValues(this, input_bcs, pin);
109 
110  // Coordinates
111  if (std::strcmp(COORDINATE_SYSTEM, "cartesian") == 0) {
112  pcoord = new Cartesian(this, pin, false);
113  } else if (std::strcmp(COORDINATE_SYSTEM, "cylindrical") == 0) {
114  pcoord = new Cylindrical(this, pin, false);
115  } else if (std::strcmp(COORDINATE_SYSTEM, "spherical_polar") == 0) {
116  pcoord = new SphericalPolar(this, pin, false);
117  } else if (std::strcmp(COORDINATE_SYSTEM, "minkowski") == 0) {
118  pcoord = new Minkowski(this, pin, false);
119  } else if (std::strcmp(COORDINATE_SYSTEM, "schwarzschild") == 0) {
120  pcoord = new Schwarzschild(this, pin, false);
121  } else if (std::strcmp(COORDINATE_SYSTEM, "kerr-schild") == 0) {
122  pcoord = new KerrSchild(this, pin, false);
123  } else if (std::strcmp(COORDINATE_SYSTEM, "gr_user") == 0) {
124  pcoord = new GRUser(this, pin, false);
125  }
126 
127  // Reconstruction: constructor may implicitly depend on Coordinates, and PPM variable
128  // floors depend on EOS, but EOS isn't needed in Reconstruction constructor-> this is ok
129  precon = new Reconstruction(this, pin);
130 
131  if (pm->multilevel) pmr = new MeshRefinement(this, pin);
132 
133  // physics-related, per-MeshBlock objects: may depend on Coordinates for diffusion
134  // terms, and may enroll quantities in AMR and BoundaryVariable objs. in BoundaryValues
135 
136  // TODO(felker): prepare this section of the MeshBlock ctor to become more complicated
137  // for several extensions:
138  // 1) allow solver to compile without a Hydro class (or with a Hydro class for the
139  // background fluid that is not dynamically evolved)
140  // 2) MPI ranks containing MeshBlocks that solve a subset of the physics, e.g. Gravity
141  // but not Hydro.
142  // 3) MAGNETIC_FIELDS_ENABLED, SELF_GRAVITY_ENABLED, NSCALARS, (future) FLUID_ENABLED,
143  // etc. become runtime switches
144 
145  // if (FLUID_ENABLED) {
146  // if (this->hydro_block)
147  phydro = new Hydro(this, pin);
148  // } else
149  // }
150  // Regardless, advance MeshBlock's local counter (initialized to bvars_next_phys_id=1)
151  // Greedy reservation of phys IDs (only 1 of 2 needed for Hydro if multilevel==false)
152  pbval->AdvanceCounterPhysID(HydroBoundaryVariable::max_phys_id);
153  // }
154  if (MAGNETIC_FIELDS_ENABLED) {
155  // if (this->field_block)
156  pfield = new Field(this, pin);
157  pbval->AdvanceCounterPhysID(FaceCenteredBoundaryVariable::max_phys_id);
158  }
159  if (SELF_GRAVITY_ENABLED) {
160  // if (this->grav_block)
161  pgrav = new Gravity(this, pin);
162  pbval->AdvanceCounterPhysID(CellCenteredBoundaryVariable::max_phys_id);
163  if (SELF_GRAVITY_ENABLED == 2)
164  pmg = new MGGravity(pmy_mesh->pmgrd, this);
165  }
166  if (NSCALARS > 0) {
167  // if (this->scalars_block)
168  pscalars = new PassiveScalars(this, pin);
169  pbval->AdvanceCounterPhysID(CellCenteredBoundaryVariable::max_phys_id);
170  }
171  // KGF: suboptimal solution, since developer must copy/paste BoundaryVariable derived
172  // class type that is used in each PassiveScalars, Gravity, Field, Hydro, ... etc. class
173  // in order to correctly advance the BoundaryValues::bvars_next_phys_id_ local counter.
174 
175  // TODO(felker): check that local counter pbval->bvars_next_phys_id_ agrees with shared
176  // Mesh::next_phys_id_ counter (including non-BoundaryVariable / per-MeshBlock reserved
177  // values). Compare both private member variables via BoundaryValues::CheckCounterPhysID
178 
179  peos = new EquationOfState(this, pin);
180  pthermo = new Thermodynamics(this, pin);
181  pchem = new CHEMISTRY(this, pin);
182  prad = new Radiation(this, pin);
183  pphy = new Physics(this, pin);
184  pdiag = new Diagnostics(this, pin);
185  pdebug = new Debugger(this);
186 
187  // Create user mesh data
189 
190  return;
191 }
192 
193 //----------------------------------------------------------------------------------------
194 // MeshBlock constructor for restarts
195 
196 MeshBlock::MeshBlock(int igid, int ilid, Mesh *pm, ParameterInput *pin,
197  LogicalLocation iloc, RegionSize input_block,
198  BoundaryFlag *input_bcs,
199  double icost, char *mbdata, int igflag) :
200  pmy_mesh(pm), loc(iloc), block_size(input_block),
201  gid(igid), lid(ilid), gflag(igflag), nuser_out_var(), prev(nullptr), next(nullptr),
202  new_block_dt_{}, new_block_dt_hyperbolic_{}, new_block_dt_parabolic_{},
203  new_block_dt_user_{},
204  nreal_user_meshblock_data_(), nint_user_meshblock_data_(), cost_(icost) {
205  // initialize grid indices
206  is = NGHOST;
207  ie = is + block_size.nx1 - 1;
208 
209  ncells1 = block_size.nx1 + 2*NGHOST;
210  ncc1 = block_size.nx1/2 + 2*NGHOST;
211  if (pmy_mesh->f2) {
212  js = NGHOST;
213  je = js + block_size.nx2 - 1;
214  ncells2 = block_size.nx2 + 2*NGHOST;
215  ncc2 = block_size.nx2/2 + 2*NGHOST;
216  } else {
217  js = je = 0;
218  ncells2 = 1;
219  ncc2 = 1;
220  }
221 
222  if (pmy_mesh->f3) {
223  ks = NGHOST;
224  ke = ks + block_size.nx3 - 1;
225  ncells3 = block_size.nx3 + 2*NGHOST;
226  ncc3 = block_size.nx3/2 + 2*NGHOST;
227  } else {
228  ks = ke = 0;
229  ncells3 = 1;
230  ncc3 = 1;
231  }
232 
233  if (pm->multilevel) {
234  cnghost = (NGHOST + 1)/2 + 1;
235  cis = NGHOST; cie = cis + block_size.nx1/2 - 1;
236  cjs = cje = cks = cke = 0;
237  if (pmy_mesh->f2) // 2D or 3D
238  cjs = NGHOST, cje = cjs + block_size.nx2/2 - 1;
239  if (pmy_mesh->f3) // 3D
240  cks = NGHOST, cke = cks + block_size.nx3/2 - 1;
241  }
242 
243  // (re-)create mesh-related objects in MeshBlock
244 
245  // Boundary
246  pbval = new BoundaryValues(this, input_bcs, pin);
247 
248  // Coordinates
249  if (std::strcmp(COORDINATE_SYSTEM, "cartesian") == 0) {
250  pcoord = new Cartesian(this, pin, false);
251  } else if (std::strcmp(COORDINATE_SYSTEM, "cylindrical") == 0) {
252  pcoord = new Cylindrical(this, pin, false);
253  } else if (std::strcmp(COORDINATE_SYSTEM, "spherical_polar") == 0) {
254  pcoord = new SphericalPolar(this, pin, false);
255  } else if (std::strcmp(COORDINATE_SYSTEM, "minkowski") == 0) {
256  pcoord = new Minkowski(this, pin, false);
257  } else if (std::strcmp(COORDINATE_SYSTEM, "schwarzschild") == 0) {
258  pcoord = new Schwarzschild(this, pin, false);
259  } else if (std::strcmp(COORDINATE_SYSTEM, "kerr-schild") == 0) {
260  pcoord = new KerrSchild(this, pin, false);
261  } else if (std::strcmp(COORDINATE_SYSTEM, "gr_user") == 0) {
262  pcoord = new GRUser(this, pin, false);
263  }
264 
265  // Reconstruction (constructor may implicitly depend on Coordinates)
266  precon = new Reconstruction(this, pin);
267 
268  if (pm->multilevel) pmr = new MeshRefinement(this, pin);
269 
270  // (re-)create physics-related objects in MeshBlock
271 
272  // if (FLUID_ENABLED) {
273  // if (this->hydro_block)
274  phydro = new Hydro(this, pin);
275  // } else
276  // }
277  // Regardless, advance MeshBlock's local counter (initialized to bvars_next_phys_id=1)
278  // Greedy reservation of phys IDs (only 1 of 2 needed for Hydro if multilevel==false)
279  pbval->AdvanceCounterPhysID(HydroBoundaryVariable::max_phys_id);
280  // }
281  if (MAGNETIC_FIELDS_ENABLED) {
282  // if (this->field_block)
283  pfield = new Field(this, pin);
284  pbval->AdvanceCounterPhysID(FaceCenteredBoundaryVariable::max_phys_id);
285  }
286  if (SELF_GRAVITY_ENABLED) {
287  // if (this->grav_block)
288  pgrav = new Gravity(this, pin);
289  pbval->AdvanceCounterPhysID(CellCenteredBoundaryVariable::max_phys_id);
290  if (SELF_GRAVITY_ENABLED == 2)
291  pmg = new MGGravity(pmy_mesh->pmgrd, this);
292  }
293 
294  if (NSCALARS > 0) {
295  // if (this->scalars_block)
296  pscalars = new PassiveScalars(this, pin);
297  pbval->AdvanceCounterPhysID(CellCenteredBoundaryVariable::max_phys_id);
298  }
299 
300  peos = new EquationOfState(this, pin);
301  pthermo = new Thermodynamics(this, pin);
302  pchem = new CHEMISTRY(this, pin);
303  prad = new Radiation(this, pin);
304  pphy = new Physics(this, pin);
305  pdiag = new Diagnostics(this, pin);
306  pdebug = new Debugger(this);
307 
309 
310  std::size_t os = 0;
311  // NEW_OUTPUT_TYPES:
312 
313  // load hydro and field data
314  std::memcpy(phydro->u.data(), &(mbdata[os]), phydro->u.GetSizeInBytes());
315  // outflow boundary condition
316  if (pbval->block_bcs[inner_x1] == BoundaryFlag::outflow) {
317  int kl = block_size.nx3 == 1 ? ks : ks-NGHOST;
318  int ku = block_size.nx3 == 1 ? ke : ke+NGHOST;
319  int jl = block_size.nx2 == 1 ? js : js-NGHOST;
320  int ju = block_size.nx2 == 1 ? je : je+NGHOST;
322  phydro->w, pfield->bcc, pcoord,
323  is-NGHOST, is-1, jl, ju, kl, ku);
324  phydro->w1 = phydro->w;
325  }
326  // load it into the other memory register(s) too
327  std::memcpy(phydro->u1.data(), &(mbdata[os]), phydro->u1.GetSizeInBytes());
328  os += phydro->u.GetSizeInBytes();
329  if (GENERAL_RELATIVITY) {
330  std::memcpy(phydro->w.data(), &(mbdata[os]), phydro->w.GetSizeInBytes());
331  os += phydro->w.GetSizeInBytes();
332  std::memcpy(phydro->w1.data(), &(mbdata[os]), phydro->w1.GetSizeInBytes());
333  os += phydro->w1.GetSizeInBytes();
334  }
335  if (MAGNETIC_FIELDS_ENABLED) {
336  std::memcpy(pfield->b.x1f.data(), &(mbdata[os]), pfield->b.x1f.GetSizeInBytes());
337  std::memcpy(pfield->b1.x1f.data(), &(mbdata[os]), pfield->b1.x1f.GetSizeInBytes());
338  os += pfield->b.x1f.GetSizeInBytes();
339  std::memcpy(pfield->b.x2f.data(), &(mbdata[os]), pfield->b.x2f.GetSizeInBytes());
340  std::memcpy(pfield->b1.x2f.data(), &(mbdata[os]), pfield->b1.x2f.GetSizeInBytes());
341  os += pfield->b.x2f.GetSizeInBytes();
342  std::memcpy(pfield->b.x3f.data(), &(mbdata[os]), pfield->b.x3f.GetSizeInBytes());
343  std::memcpy(pfield->b1.x3f.data(), &(mbdata[os]), pfield->b1.x3f.GetSizeInBytes());
344  os += pfield->b.x3f.GetSizeInBytes();
345  }
346 
347  // (conserved variable) Passive scalars:
348  if (NSCALARS > 0) {
349  std::memcpy(pscalars->s.data(), &(mbdata[os]), pscalars->s.GetSizeInBytes());
350  // load it into the other memory register(s) too
351  std::memcpy(pscalars->s1.data(), &(mbdata[os]), pscalars->s1.GetSizeInBytes());
352  os += pscalars->s.GetSizeInBytes();
353  }
354 
355  // load user MeshBlock data
356  for (int n=0; n<nint_user_meshblock_data_; n++) {
357  std::memcpy(iuser_meshblock_data[n].data(), &(mbdata[os]),
358  iuser_meshblock_data[n].GetSizeInBytes());
360  }
361  for (int n=0; n<nreal_user_meshblock_data_; n++) {
362  std::memcpy(ruser_meshblock_data[n].data(), &(mbdata[os]),
363  ruser_meshblock_data[n].GetSizeInBytes());
365  }
366 
367  // load physics data
368  os += pphy->LoadRestartData(&(mbdata[os]));
369 
370  return;
371 }
372 
373 //----------------------------------------------------------------------------------------
374 // MeshBlock destructor
375 
377  if (prev != nullptr) prev->next = next;
378  if (next != nullptr) next->prev = prev;
379 
380  delete pcoord;
381  delete precon;
382  if (pmy_mesh->multilevel) delete pmr;
383 
384  delete phydro;
385  if (MAGNETIC_FIELDS_ENABLED) delete pfield;
386  delete peos;
387  if (SELF_GRAVITY_ENABLED) delete pgrav;
388  if (NSCALARS > 0) delete pscalars;
389 
390  // BoundaryValues should be destructed AFTER all BoundaryVariable objects are destroyed
391  delete pbval;
392  // delete user output variables array
393  if (nuser_out_var > 0) {
394  delete [] user_out_var_names_;
395  delete [] user_out_var_longnames_;
396  delete [] user_out_var_units_;
397  }
398  // delete user MeshBlock data
401 
402  delete pthermo;
403  delete pchem;
404  delete prad;
405  delete pphy;
406 
407  while (pdiag->prev != nullptr)
408  delete pdiag->prev;
409  while (pdiag->next != nullptr)
410  delete pdiag->next;
411  delete pdiag;
412 
413  while (pdebug->prev != nullptr)
414  delete pdebug->prev;
415  while (pdebug->next != nullptr)
416  delete pdebug->next;
417  delete pdebug;
418 }
419 
420 //----------------------------------------------------------------------------------------
422 // \brief Allocate Real AthenaArrays for user-defned data in MeshBlock
423 
425  if (nreal_user_meshblock_data_ != 0) {
426  std::stringstream msg;
427  msg << "### FATAL ERROR in MeshBlock::AllocateRealUserMeshBlockDataField"
428  << std::endl << "User MeshBlock data arrays are already allocated" << std::endl;
429  ATHENA_ERROR(msg);
430  }
433  return;
434 }
435 
436 //----------------------------------------------------------------------------------------
438 // \brief Allocate integer AthenaArrays for user-defned data in MeshBlock
439 
441  if (nint_user_meshblock_data_ != 0) {
442  std::stringstream msg;
443  msg << "### FATAL ERROR in MeshBlock::AllocateIntusermeshblockDataField"
444  << std::endl << "User MeshBlock data arrays are already allocated" << std::endl;
445  ATHENA_ERROR(msg);
446  return;
447  }
450  return;
451 }
452 
453 //----------------------------------------------------------------------------------------
455 // \brief Allocate user-defined output variables
456 
458  if (n <= 0) return;
459  if (nuser_out_var != 0) {
460  std::stringstream msg;
461  msg << "### FATAL ERROR in MeshBlock::AllocateUserOutputVariables"
462  << std::endl << "User output variables are already allocated." << std::endl;
463  ATHENA_ERROR(msg);
464  return;
465  }
466  nuser_out_var = n;
467  user_out_var.NewAthenaArray(nuser_out_var, ncells3, ncells2, ncells1);
468  user_out_var_names_ = new std::string[n];
469  user_out_var_longnames_ = new std::string[n];
470  user_out_var_units_ = new std::string[n];
471  return;
472 }
473 
474 
475 //----------------------------------------------------------------------------------------
477 // \brief set the user-defined output variable name
478 
479 void MeshBlock::SetUserOutputVariableName(int n, const char *name,
480  const char *long_name, const char *units) {
481  if (n >= nuser_out_var) {
482  std::stringstream msg;
483  msg << "### FATAL ERROR in MeshBlock::SetUserOutputVariableName"
484  << std::endl << "User output variable is not allocated." << std::endl;
485  ATHENA_ERROR(msg);
486  return;
487  }
488  user_out_var_names_[n] = name;
489  user_out_var_longnames_[n] = long_name;
490  user_out_var_units_[n] = units;
491  return;
492 }
493 
494 //----------------------------------------------------------------------------------------
496 // \brief Calculate the block data size required for restart.
497 
499  std::size_t size;
500  // NEW_OUTPUT_TYPES:
501  size = phydro->u.GetSizeInBytes();
502  if (GENERAL_RELATIVITY) {
503  size += phydro->w.GetSizeInBytes();
504  size += phydro->w1.GetSizeInBytes();
505  }
506  if (MAGNETIC_FIELDS_ENABLED)
507  size += (pfield->b.x1f.GetSizeInBytes() + pfield->b.x2f.GetSizeInBytes()
508  + pfield->b.x3f.GetSizeInBytes());
509  if (SELF_GRAVITY_ENABLED)
510  size += pgrav->phi.GetSizeInBytes();
511  if (NSCALARS > 0)
512  size += pscalars->s.GetSizeInBytes();
513 
514  // calculate user MeshBlock data size
515  for (int n=0; n<nint_user_meshblock_data_; n++)
516  size += iuser_meshblock_data[n].GetSizeInBytes();
517  for (int n=0; n<nreal_user_meshblock_data_; n++)
518  size += ruser_meshblock_data[n].GetSizeInBytes();
519 
520  // physics data
521  size += pphy->RestartDataSizeInBytes();
522 
523  return size;
524 }
525 
526 //----------------------------------------------------------------------------------------
528 // \brief stop time measurement and accumulate it in the MeshBlock cost
529 
531  if (pmy_mesh->lb_manual_) {
532  cost_ = std::min(cost, TINY_NUMBER);
533  pmy_mesh->lb_flag_ = true;
534  }
535 }
536 
537 //----------------------------------------------------------------------------------------
539 // \brief reset the MeshBlock cost for automatic load balancing
540 
542  if (pmy_mesh->lb_automatic_) cost_ = TINY_NUMBER;
543 }
544 
545 //----------------------------------------------------------------------------------------
547 // \brief start time measurement for automatic load balancing
548 
550  if (pmy_mesh->lb_automatic_) {
551 #ifdef OPENMP_PARALLEL
552  lb_time_ = omp_get_wtime();
553 #else
554  lb_time_ = static_cast<double>(clock());
555 #endif
556  }
557 }
558 
559 //----------------------------------------------------------------------------------------
561 // \brief stop time measurement and accumulate it in the MeshBlock cost
562 
564  if (pmy_mesh->lb_automatic_) {
565 #ifdef OPENMP_PARALLEL
566  lb_time_ = omp_get_wtime() - lb_time_;
567 #else
568  lb_time_ = static_cast<double>(clock()) - lb_time_;
569 #endif
570  cost_ += lb_time_;
571  }
572 }
573 
574 
576  vars_cc_.push_back(pvar_cc);
577  return;
578 }
579 
580 
582  vars_fc_.push_back(pvar_fc);
583  return;
584 }
585 
586 
587 // TODO(felker): consider merging the MeshRefinement::pvars_cc/fc_ into the
588 // MeshBlock::pvars_cc/fc_. Would need to weaken the MeshBlock std::vector to use tuples
589 // of pointers instead of a std::vector of references, so that:
590 // - nullptr can be passed for the second entry if multilevel==false
591 // - we can rebind the pointers to Hydro for GR purposes in bvals_refine.cpp
592 // If GR, etc. in the future requires additional flexiblity from non-refinement load
593 // balancing, we will need to use ptrs instead of references anyways, and add:
594 
595 // void MeshBlock::SetHydroData(HydroBoundaryQuantity hydro_type)
596 // Hydro *ph = pmy_block_->phydro;
597 // // hard-coded assumption that, if multilevel, then Hydro is always present
598 // // and enrolled in mesh refinement in the first pvars_cc_ vector entry
599 // switch (hydro_type) {
600 // case (HydroBoundaryQuantity::cons): {
601 // pvars_cc_.front() = &ph->u;
602 // break;
603 // }
604 // case (HydroBoundaryQuantity::prim): {
605 // pvars_cc_.front() = &ph->w;
606 // break;
607 // }
608 // }
609 // return;
610 // }
std::size_t GetSizeInBytes() const
Diagnostics * next
Definition: diagnostics.hpp:18
Diagnostics * prev
Definition: diagnostics.hpp:18
void ConservedToPrimitive(AthenaArray< Real > &cons, const AthenaArray< Real > &prim_old, const FaceField &b, AthenaArray< Real > &prim, AthenaArray< Real > &bcc, Coordinates *pco, int il, int iu, int jl, int ju, int kl, int ku)
AthenaArray< Real > u1
Definition: hydro.hpp:46
AthenaArray< Real > w
Definition: hydro.hpp:45
AthenaArray< Real > u
Definition: hydro.hpp:45
AthenaArray< Real > w1
Definition: hydro.hpp:46
MeshRefinement * pmr
Definition: mesh.hpp:122
void AllocateUserOutputVariables(int n)
Definition: meshblock.cpp:457
std::string * user_out_var_units_
Definition: mesh.hpp:112
Reconstruction * precon
Definition: mesh.hpp:121
double cost_
Definition: mesh.hpp:187
Physics * pphy
Definition: mesh.hpp:135
void StopTimeMeasurement()
Definition: meshblock.cpp:563
BoundaryValues * pbval
Definition: mesh.hpp:120
Debugger * pdebug
Definition: mesh.hpp:137
int nint_user_meshblock_data_
Definition: mesh.hpp:170
Gravity * pgrav
Definition: mesh.hpp:127
int cie
Definition: mesh.hpp:99
int ncc2
Definition: mesh.hpp:96
friend class BoundaryValues
Definition: mesh.hpp:69
MeshBlock * prev
Definition: mesh.hpp:139
void RegisterMeshBlockData(AthenaArray< Real > &pvar_cc)
Definition: meshblock.cpp:575
friend class Hydro
Definition: mesh.hpp:73
RegionSize block_size
Definition: mesh.hpp:91
Hydro * phydro
Definition: mesh.hpp:125
AthenaArray< Real > * ruser_meshblock_data
Definition: mesh.hpp:115
EquationOfState * peos
Definition: mesh.hpp:130
int cjs
Definition: mesh.hpp:99
Mesh * pmy_mesh
Definition: mesh.hpp:89
MGGravity * pmg
Definition: mesh.hpp:128
int ie
Definition: mesh.hpp:97
int cnghost
Definition: mesh.hpp:99
void AllocateRealUserMeshBlockDataField(int n)
Definition: meshblock.cpp:424
AthenaArray< Real > user_out_var
Definition: mesh.hpp:109
int js
Definition: mesh.hpp:97
std::vector< std::reference_wrapper< AthenaArray< Real > > > vars_cc_
Definition: mesh.hpp:171
int ncells1
Definition: mesh.hpp:94
Diagnostics * pdiag
Definition: mesh.hpp:136
int je
Definition: mesh.hpp:97
void AllocateIntUserMeshBlockDataField(int n)
Definition: meshblock.cpp:440
AthenaArray< int > * iuser_meshblock_data
Definition: mesh.hpp:116
std::string * user_out_var_longnames_
Definition: mesh.hpp:111
int ncc3
Definition: mesh.hpp:96
int ke
Definition: mesh.hpp:97
int cis
Definition: mesh.hpp:99
PassiveScalars * pscalars
Definition: mesh.hpp:129
int ncells3
Definition: mesh.hpp:94
void StartTimeMeasurement()
Definition: meshblock.cpp:549
double lb_time_
Definition: mesh.hpp:187
Radiation * prad
Definition: mesh.hpp:134
int ncc1
Definition: mesh.hpp:96
Chemistry * pchem
Definition: mesh.hpp:133
int is
Definition: mesh.hpp:97
int cks
Definition: mesh.hpp:99
int nuser_out_var
Definition: mesh.hpp:108
void ResetTimeMeasurement()
Definition: meshblock.cpp:541
int ks
Definition: mesh.hpp:97
Coordinates * pcoord
Definition: mesh.hpp:119
int cke
Definition: mesh.hpp:99
int cje
Definition: mesh.hpp:99
int nreal_user_meshblock_data_
Definition: mesh.hpp:170
void InitUserMeshBlockData(ParameterInput *pin)
void SetUserOutputVariableName(int n, const char *name, const char *long_name="", const char *units="")
Definition: meshblock.cpp:479
Thermodynamics * pthermo
Definition: mesh.hpp:132
std::vector< std::reference_wrapper< FaceField > > vars_fc_
Definition: mesh.hpp:172
int ncells2
Definition: mesh.hpp:94
std::size_t GetBlockSizeInBytes()
Definition: meshblock.cpp:498
void SetCostForLoadBalancing(double cost)
Definition: meshblock.cpp:530
MeshBlock(int igid, int ilid, LogicalLocation iloc, RegionSize input_size, BoundaryFlag *input_bcs, Mesh *pm, ParameterInput *pin, int igflag, bool ref_flag=false)
Definition: meshblock.cpp:52
MeshBlock * next
Definition: mesh.hpp:139
std::string * user_out_var_names_
Definition: mesh.hpp:110
Field * pfield
Definition: mesh.hpp:126
Definition: mesh.hpp:197
bool lb_flag_
Definition: mesh.hpp:315
const bool multilevel
Definition: mesh.hpp:240
bool lb_manual_
Definition: mesh.hpp:315
bool lb_automatic_
Definition: mesh.hpp:315
const bool f3
Definition: mesh.hpp:238
MGGravityDriver * pmgrd
Definition: mesh.hpp:257
const bool f2
Definition: mesh.hpp:238
double min(double x1, double x2, double x3)
Definition: core.h:15
int nx3
Definition: athena.hpp:88
int nx2
Definition: athena.hpp:88
int nx1
Definition: athena.hpp:88