diff --git a/Docs/sphinx/Utility.rst b/Docs/sphinx/Utility.rst index 703609a03..cd4ca4dc6 100644 --- a/Docs/sphinx/Utility.rst +++ b/Docs/sphinx/Utility.rst @@ -99,7 +99,66 @@ This code contains data structures used to handle data read from plt files that Diagnostics =========== -Placeholder. Once the porting of diagnostics from PeleLMeX to PelePhysics/PeleC is complete, documentation can be added here. +Analysing the data a-posteriori can become extremely cumbersome when dealing with extreme datasets. +PelePhysics supports shared capability for diagnostics available at runtime during PeleC and PeleLMeX simulation +and more are under development. +Currently, the list of diagnostic contains: + +* ``DiagFramePlane`` : extract a plane aligned in the 'x','y' or 'z' direction across the AMR hierarchy, writing + a 2D plotfile compatible with Amrvis, Paraview or yt. Only available for 3D simulations. +* ``DiagPDF`` : extract the PDF of a given variable and write it to an ASCII file. +* ``DiagConditional`` : extract statistics (average and standard deviation, integral or sum) of a + set of variables conditioned on the value of given variable and write it to an ASCII file. + +When using `DiagPDF` or `DiagConditional`, it is possible to narrow down the diagnostic to a region of interest +by specifying a set of filters, defining a range of interest for a variable. Note also that for these two diagnostics, +fine-covered regions are masked. An arbitrary number of diagnostics can be named in a list with the ``diagnostics`` +keyword, then each diagnostic name in the list becomes a keyword allowing the type and specific options for +that diagnostic to be specified. +The following provide examples for each diagnostic in PeleLMeX (in PeleC, all diagnostics would be prefixed with +`pelec` instead of `peleLM`: + +:: + + #--------------------------DIAGNOSTICS------------------------ + + peleLM.diagnostics = xnormP condT pdfTest + + peleLM.xnormP.type = DiagFramePlane # Diagnostic type + peleLM.xnormP.file = xNorm5mm # Output file prefix + peleLM.xnormP.normal = 0 # Plane normal (0, 1 or 2 for x, y or z) + peleLM.xnormP.center = 0.005 # Coordinate in the normal direction + peleLM.xnormP.int = 5 # Frequency (as step #) for performing the diagnostic + peleLM.xnormP.interpolation = Linear # [OPT, DEF=Linear] Interpolation type : Linear or Quadratic + peleLM.xnormP.field_names = x_velocity mag_vort density # List of variables outputted to the 2D pltfile + peleLM.xnormP.n_files = 2 # [OPT, DEF="min(256,NProcs)"] Number of files to write per level + peleLM.xnormP.dump_ghost_if_OOB = 1 # [OPT, DEF=false] if the specified coordinate is out-of-bounds, a plane of ghost cells in that direction will be dumped (for debugging purposes). If false, an error is raised if the requested plane is OOB. + + peleLM.condT.type = DiagConditional # Diagnostic type + peleLM.condT.file = condTest # Output file prefix + peleLM.condT.int = 5 # Frequency (as step #) for performing the diagnostic + peleLM.condT.filters = xHigh stoich # [OPT, DEF=None] List of filters + peleLM.condT.xHigh.field_name = x # Filter field + peleLM.condT.xHigh.value_greater = 0.006 # Filter definition : value_greater, value_less, value_inrange + peleLM.condT.stoich.field_name = mixture_fraction # Filter field + peleLM.condT.stoich.value_inrange = 0.053 0.055 # Filter definition : value_greater, value_less, value_inrange + peleLM.condT.conditional_type = Average # Conditional type : Average, Integral or Sum + peleLM.condT.nBins = 50 # Number of bins for the conditioning variable + peleLM.condT.condition_field_name = temp # Conditioning variable name + peleLM.condT.field_names = HeatRelease I_R(CH4) I_R(H2) # List of variables to be treated + + peleLM.pdfTest.type = DiagPDF # Diagnostic type + peleLM.pdfTest.file = PDFTest # Output file prefix + peleLM.pdfTest.int = 5 # Frequency (as step #) for performing the diagnostic + peleLM.pdfTest.filters = innerFlame # [OPT, DEF=None] List of filters + peleLM.pdfTest.innerFlame.field_name = temp # Filter field + peleLM.pdfTest.innerFlame.value_inrange = 450.0 1500.0 # Filter definition : value_greater, value_less, value_inrange + peleLM.pdfTest.nBins = 50 # Number of bins for the PDF + peleLM.pdfTest.normalized = 1 # [OPT, DEF=1] PDF is normalized (i.e. integral is unity) ? + peleLM.pdfTest.volume_weighted = 1 # [OPT, DEF=1] Computation of the PDF is volume weighted ? + peleLM.pdfTest.range = 0.0 2.0 # [OPT, DEF=data min/max] Specify the range of the PDF + peleLM.pdfTest.field_name = x_velocity # Variable of interest + Filter ====== @@ -169,4 +228,4 @@ For example, users can call ``eos.PYT2R`` as follows: // Calculate density in CGS units, then convert to MKS eos.PYT2R(m2c::P(P_mean), massfrac, Temp, rho_cgs); - amrex::Real rho = c2m::Rho(rho_cgs); // Convert eos density to MKS \ No newline at end of file + amrex::Real rho = c2m::Rho(rho_cgs); // Convert eos density to MKS diff --git a/Source/Utility/Diagnostics/DiagFramePlane.H b/Source/Utility/Diagnostics/DiagFramePlane.H index 2b6c5cf81..b79ca612d 100644 --- a/Source/Utility/Diagnostics/DiagFramePlane.H +++ b/Source/Utility/Diagnostics/DiagFramePlane.H @@ -84,6 +84,7 @@ private: // Plane definition int m_normal; amrex::GpuArray m_center; + bool m_dump_ghost_if_OOB; // Interpolation data InterpType m_interpType; diff --git a/Source/Utility/Diagnostics/DiagFramePlane.cpp b/Source/Utility/Diagnostics/DiagFramePlane.cpp index 03dc9142c..79839111d 100644 --- a/Source/Utility/Diagnostics/DiagFramePlane.cpp +++ b/Source/Utility/Diagnostics/DiagFramePlane.cpp @@ -76,6 +76,8 @@ DiagFramePlane::init(const std::string& a_prefix, std::string_view a_diagName) } else if (center.size() == 1) { m_center[m_normal] = center[0]; } + m_dump_ghost_if_OOB = 0; + pp.query("dump_ghost_if_OOB", m_dump_ghost_if_OOB); // Interpolation std::string intType = "Quadratic"; @@ -155,13 +157,11 @@ DiagFramePlane::prepare( dx[m_normal]; int k0 = static_cast(std::round(dist)); dist -= static_cast(k0); - m_k0[lev] = k0; if (m_interpType == Quadratic) { // Quadratic interp. weights on k0-1, k0, k0+1 m_intwgt[lev][0] = 0.5 * (dist - 1.0) * (dist - 2.0); m_intwgt[lev][1] = dist * (2.0 - dist); m_intwgt[lev][2] = 0.5 * dist * (dist - 1.0); - ; } else if (m_interpType == Linear) { // linear interp. weights on k0-1, k0, k0+1 if (dist > 0.0) { @@ -178,6 +178,41 @@ DiagFramePlane::prepare( m_intwgt[lev][2] = 0.0; } } + + if (k0 < a_geoms[lev].Domain().smallEnd(m_normal)) { + if (!m_dump_ghost_if_OOB) { + amrex::Error( + "DiagFramePlane: Requested location is out of bounds (" + m_diagfile + + ")"); + } else if (lev == 0) { + amrex::Print() << "DiagFramePlane: Dumping ghost cells because " + "requested location is " + "of bounds (" + + m_diagfile + ")" + << std::endl; + } + k0 = a_geoms[lev].Domain().smallEnd(m_normal); + m_intwgt[lev][0] = 1.0; + m_intwgt[lev][1] = 0.0; + m_intwgt[lev][2] = 0.0; + } else if (k0 > a_geoms[lev].Domain().bigEnd(m_normal)) { + if (!m_dump_ghost_if_OOB) { + amrex::Error( + "DiagFramePlane: Requested location is out of bounds (" + m_diagfile + + ")"); + } else if (lev == 0) { + amrex::Print() << "DiagFramePlane: Dumping ghost cells because " + "requested location is " + "of bounds (" + + m_diagfile + ")" + << std::endl; + } + k0 = a_geoms[lev].Domain().bigEnd(m_normal); + m_intwgt[lev][0] = 0.0; + m_intwgt[lev][1] = 0.0; + m_intwgt[lev][2] = 1.0; + } + m_k0[lev] = k0; } // Assemble the 2D slice boxArray