PetIBM 0.5.4
Toolbox and applications of the immersed-boundary method for distributed-memory architectures
|
Collections of the namespaces containing rarely used functions. More...
Namespaces | |
namespace | petibm::delta |
A namespace of all kinds of discretized delta functions. | |
namespace | petibm::io |
A namespace of I/O-related functions. | |
namespace | petibm::misc |
A namespace holding miscellaneous functions. | |
namespace | petibm::parser |
A collection of YAML node parsers. | |
Classes | |
class | petibm::misc::LinInterpBase |
Abstract Base Class of a linear interpolation object. More... | |
class | petibm::misc::TriLinInterp |
Class to perform a trilinear interpolation at a point in a 3D domain. More... | |
class | petibm::misc::BiLinInterp |
Class to perform a bilinear interpolation at a point in a 3D domain. More... | |
struct | petibm::misc::LoopBound |
A helper struct to make looping function easier. More... | |
class | petibm::misc::ProbeBase |
Abstract Base Class of a probe. More... | |
class | petibm::misc::ProbeVolume |
Probe class to monitor a volume region of the domain. More... | |
class | petibm::misc::ProbePoint |
Probe class to monitor the solution at a single point. More... | |
Typedefs | |
typedef std::shared_ptr< misc::LinInterpBase > | petibm::type::LinInterp |
Type definition of LinInterp. More... | |
typedef std::shared_ptr< misc::ProbeBase > | petibm::type::Probe |
Type definition of Probe. More... | |
Functions | |
PetscReal | petibm::delta::Roma_et_al_1999 (const PetscReal &r, const PetscReal &dr) |
Regularized delta function from Roma et al. (1999). More... | |
PetscReal | petibm::delta::Peskin_2002 (const PetscReal &r, const PetscReal &dr) |
Regularized delta function from Peskin (2002). More... | |
PetscReal | petibm::delta::delta (const std::vector< PetscReal > &source, const std::vector< PetscReal > &target, const std::vector< PetscReal > &widths, const DeltaKernel &kernel) |
Discrete delta function. More... | |
PetscErrorCode | petibm::io::readLagrangianPoints (const std::string &file, PetscInt &nPts, type::RealVec2D &coords) |
Read the number Lagrangian points and their coordinates from a file. More... | |
PetscErrorCode | petibm::io::print (const std::string &info) |
Print information of a parallel object to standard output. More... | |
PetscErrorCode | petibm::io::writeHDF5Vecs (const MPI_Comm comm, const std::string &filePath, const std::string &loc, const std::vector< std::string > &names, const std::vector< Vec > &vecs, const PetscFileMode mode=FILE_MODE_WRITE) |
Write a vector of Vec objects to a HDF5 file. More... | |
PetscErrorCode | petibm::io::writeHDF5Vecs (const MPI_Comm comm, const std::string &filePath, const std::string &loc, const std::vector< std::string > &names, const std::vector< PetscInt > &n, const std::vector< PetscReal * > &vecs, const PetscFileMode mode=FILE_MODE_WRITE) |
Write a vector of raw arrays to a HDF5 file. More... | |
PetscErrorCode | petibm::io::writeHDF5Vecs (const MPI_Comm comm, const std::string &filePath, const std::string &loc, const std::vector< std::string > &names, const type::RealVec2D &vecs, const PetscFileMode mode=FILE_MODE_WRITE) |
Write petibm::type::RealVec2D to a HDF5 file. More... | |
PetscErrorCode | petibm::io::readHDF5Vecs (const MPI_Comm comm, const std::string &filePath, const std::string &loc, const std::vector< std::string > &names, std::vector< Vec > &vecs) |
Read a vector of Vec objects from a HDF5 file. More... | |
PetscErrorCode | petibm::io::writePetscLog (const MPI_Comm comm, const std::string &filePath) |
Write a summary of the PETSc logging into a ASCII file. More... | |
PetscErrorCode | petibm::misc::createLinInterp (const MPI_Comm &comm, const type::RealVec1D &point, const type::Mesh &mesh, const type::Field &field, type::LinInterp &interp) |
Factory function to create a linear interpolation object. More... | |
PetscErrorCode | petibm::misc::checkPeriodicBC (const type::IntVec2D &bcTypes, type::BoolVec2D &periodic) |
Check if there is any periodic boundary condition and if these periodic BCs make sense. More... | |
PetscErrorCode | petibm::misc::checkBoundaryProc (const DM &da, const type::IntVec1D &n, const type::BCLoc &loc, PetscBool &onThisProc) |
Check if a boundary is on this process. More... | |
PetscErrorCode | petibm::misc::getGhostPointList (const type::Mesh &mesh, const type::Field &field, const type::BCLoc &loc, type::GhostPointsList &points) |
Get a list of ghost points on a desired boundary. More... | |
PetscErrorCode | petibm::misc::getPerpendAxes (const PetscInt &self, type::IntVec1D &pAxes) |
An utility to get the perpendicular axes of a desired axis. More... | |
PetscErrorCode | petibm::misc::getGhostTargetStencil (const type::IntVec1D &n, const type::BCLoc &loc, const type::IntVec1D &pIdx, MatStencil &ghost, MatStencil &target) |
Get the stencils of a desired ghost point and its corresponding boundary point. More... | |
PetscErrorCode | petibm::misc::stretchGrid (const PetscReal &bg, const PetscReal &ed, const PetscInt &n, const PetscReal &r, type::RealVec1D &dL) |
Calculate and return cell sizes of stretched grid in one direction. More... | |
PetscErrorCode | petibm::misc::doubleLoops (const LoopBound &bound1, const LoopBound &bound2, const std::function< PetscErrorCode(const PetscInt &, const PetscInt &)> &f) |
A helper function to carry out a double loop on a given function. More... | |
PetscErrorCode | petibm::misc::tripleLoops (const LoopBound &bound1, const LoopBound &bound2, const LoopBound &bound3, const std::function< PetscErrorCode(const PetscInt &, const PetscInt &, const PetscInt &)> &f) |
A helper function to carry out a triple loop on a given function. More... | |
PetscErrorCode | petibm::parser::getSettings (YAML::Node &node) |
Get configuration settings. More... | |
PetscErrorCode | petibm::parser::readYAMLFile (const std::string &filePath, YAML::Node &node) |
Load the content of a YAML file into a YAML node. More... | |
PetscErrorCode | petibm::parser::parseMesh (const YAML::Node &meshNode, PetscInt &dim, type::RealVec1D &bg, type::RealVec1D &ed, type::IntVec1D &nTotal, type::RealVec2D &dL) |
Parse a YAML node of cartesianMesh. More... | |
PetscErrorCode | petibm::parser::parseOneAxis (const YAML::Node &axis, PetscInt &dir, PetscReal &bg, PetscReal &ed, PetscInt &nTotal, type::RealVec1D &dL) |
Parse the info of only one direction from YAML node. More... | |
PetscErrorCode | petibm::parser::parseSubDomains (const YAML::Node &subs, const PetscReal bg, PetscInt &nTotal, PetscReal &ed, type::RealVec1D &dL) |
Parse all sub-domains in a direction. More... | |
PetscErrorCode | petibm::parser::parseOneSubDomain (const YAML::Node &sub, const PetscReal bg, PetscInt &n, PetscReal &ed, type::RealVec1D &dL) |
Parse only one sub-domain. More... | |
PetscErrorCode | petibm::parser::parseBCs (const YAML::Node &node, type::IntVec2D &bcTypes, type::RealVec2D &bcValues) |
Parse boundary conditions from a YAML node. More... | |
PetscErrorCode | petibm::parser::parseICs (const YAML::Node &node, std::vector< SymEngine::LambdaRealDoubleVisitor > &ics) |
Parse initial conditions from a YAML node. More... | |
PetscErrorCode | petibm::misc::createProbe (const MPI_Comm &comm, const YAML::Node &node, const type::Mesh &mesh, type::Probe &probe) |
Factory function to create a probe to monitor the solution. More... | |
Collections of the namespaces containing rarely used functions.
This module contains all other public functions that are rarely used by API users. Four namespaces are categorized into this module: petibm::io, petibm::parser, petibm::delta, and petibm::misc. These functions are sometimes useful at some specific situations.
Also, user-defined YAML converters are categorized into here.
typedef std::shared_ptr<misc::LinInterpBase> petibm::type::LinInterp |
Type definition of LinInterp.
Please use petibm::misc::createLinInterp to create a LinInterp object.
Definition at line 172 of file lininterp.h.
typedef std::shared_ptr<misc::ProbeBase> petibm::type::Probe |
Type definition of Probe.
Please use petibm::misc::createProbe to create a Probe object.
Example usage:
PetscErrorCode petibm::misc::checkBoundaryProc | ( | const DM & | da, |
const type::IntVec1D & | n, | ||
const type::BCLoc & | loc, | ||
PetscBool & | onThisProc | ||
) |
Check if a boundary is on this process.
da | [in] the DMDA object where this boundary belongs to. |
n | [in] a 1D array of the numbers of the cells in x, y and z directions. |
loc | [in] the desired location of this boundary. |
onThisProc | [out] a PetscBool indicating if is on this process. |
PetscErrorCode petibm::misc::checkPeriodicBC | ( | const type::IntVec2D & | bcTypes, |
type::BoolVec2D & | periodic | ||
) |
Check if there is any periodic boundary condition and if these periodic BCs make sense.
bcTypes | [in] a type::IntVec2D obtained from petibm::parser::parseBCs. |
periodic | [out] returned 2D array of PetscBool indicating the status of periodic BCs. |
This function checks if the settings of periodic boundary conditions make sense, if there is any. Also, it returns a 2D array of bools that indicates where a specific direction in a specific velocity field is periodic. For example, assuming p
is the 2D bools returned by this function, then if p[1][1]
is PETSC_TRUE
, this means the y direction of v velocity field is periodic.
PetscErrorCode petibm::misc::createLinInterp | ( | const MPI_Comm & | comm, |
const type::RealVec1D & | point, | ||
const type::Mesh & | mesh, | ||
const type::Field & | field, | ||
type::LinInterp & | interp | ||
) |
Factory function to create a linear interpolation object.
comm | [in] MPI communicator |
point | [in] Coordinates of the point to interpolate |
mesh | [in] Cartesian mesh object |
field | [in] Index of the field to monitor |
interp | [out] Linear interpolation object |
Definition at line 16 of file lininterp.cpp.
PetscErrorCode petibm::misc::createProbe | ( | const MPI_Comm & | comm, |
const YAML::Node & | node, | ||
const type::Mesh & | mesh, | ||
type::Probe & | probe | ||
) |
Factory function to create a probe to monitor the solution.
comm | [in] MPI communicator |
node | [in] YAML configuration node |
mesh | [in] Cartesian mesh object |
probe | [out] Probe |
Definition at line 23 of file probes.cpp.
PetscReal petibm::delta::delta | ( | const std::vector< PetscReal > & | source, |
const std::vector< PetscReal > & | target, | ||
const std::vector< PetscReal > & | widths, | ||
const DeltaKernel & | kernel | ||
) |
PetscErrorCode petibm::misc::getGhostPointList | ( | const type::Mesh & | mesh, |
const type::Field & | field, | ||
const type::BCLoc & | loc, | ||
type::GhostPointsList & | points | ||
) |
Get a list of ghost points on a desired boundary.
mesh | [in] a type::Mesh object. |
field | [in] the desired velocity field. |
loc | [in] the location of the desired boundary. |
points | [out] return type::GhostPointsList. |
PetscErrorCode petibm::misc::getGhostTargetStencil | ( | const type::IntVec1D & | n, |
const type::BCLoc & | loc, | ||
const type::IntVec1D & | pIdx, | ||
MatStencil & | ghost, | ||
MatStencil & | target | ||
) |
Get the stencils of a desired ghost point and its corresponding boundary point.
n | [in] a 1D array of numbers of cells in different directions of desired velocity field. |
loc | [in] the location of the desired boundary. |
pIdx | [in] a 1D array of the indices of the perpendicular axes. |
ghost | [out] the returned stencil of the desired ghost point. |
target | [out] the stencil of the corresponding boundary point. |
For example, in a 5x5x5 3D mesh, if we want to get the stencils of the ghost and corresponding points on the XPLUS
boundary with j
and k
indices equal to 2
and 3
respectively, then the following arguments can be used:
Then, the returned ghost and boundary points will be:
Note, the order of i
, j
, and k
in a MatStencil
is {k, j, i}
.
PetscErrorCode petibm::misc::getPerpendAxes | ( | const PetscInt & | self, |
type::IntVec1D & | pAxes | ||
) |
PetscErrorCode petibm::parser::getSettings | ( | YAML::Node & | node | ) |
Get configuration settings.
The function parses the command-line arguments and reads the various YAML files that contain simulation settings. If a non-empty YAML node is given, sub-nodes might be overwritten.
node | [out] YAML node with simulation settings. |
The function looks for the following optional command-line arguments:
-directory
: working directory (default: current working directory).-output
: output directory (default: folder output
in the working directory).-config
: YAML file path with all settings (default: config.yaml
file in the working directory).-mesh
: YAML file path with mesh settings (default: None, i.e., does nothing).-flow
: YAML file path with flow settings (default: None, i.e. does nothing).-parameters
: YAML file path with simulation parameters (default: None, i.e., does nothing).-bodies
: YAML file path with settings for the immersed bodies (default: None, i.e., does nothing). Definition at line 175 of file parser.cpp.
PetscErrorCode petibm::parser::parseBCs | ( | const YAML::Node & | node, |
type::IntVec2D & | bcTypes, | ||
type::RealVec2D & | bcValues | ||
) |
Parse boundary conditions from a YAML node.
node | [in] YAML configuration node. |
bcTypes | [out] BC types of different fields on different boundaries. |
bcValues | [out] BC values of different fields on different boundaries. |
This function will look into the key boundaryConditions
under the key flow
in the provided YAML node.
Definition at line 359 of file parser.cpp.
PetscErrorCode petibm::parser::parseICs | ( | const YAML::Node & | node, |
std::vector< SymEngine::LambdaRealDoubleVisitor > & | ics | ||
) |
Parse initial conditions from a YAML node.
node | [in] YAML configuration node. |
ics | [out] Lambda functions for calculating IC values of different fields. |
This function will look into the key initialVelocity
under the key flow
in the provided YAML node.
Definition at line 396 of file parser.cpp.
PetscErrorCode petibm::parser::parseMesh | ( | const YAML::Node & | meshNode, |
PetscInt & | dim, | ||
type::RealVec1D & | bg, | ||
type::RealVec1D & | ed, | ||
type::IntVec1D & | nTotal, | ||
type::RealVec2D & | dL | ||
) |
Parse a YAML node of cartesianMesh.
meshNode | [in] YAML configuration node. |
dim | [out] Number of dimensions |
bg | [out] Vector with starting coordinate in each direction. |
ed | [out] Vector with ending coordinate in each direction. |
nTotal | [out] Number of cells in each direction. |
dL | [out] Nested vector with grid cell sizes in each direction. |
This function will retrieve information from and do necessary calculations based on the data provided under the key mesh
in the YAML node.
Definition at line 239 of file parser.cpp.
PetscErrorCode petibm::parser::parseOneAxis | ( | const YAML::Node & | axis, |
PetscInt & | dir, | ||
PetscReal & | bg, | ||
PetscReal & | ed, | ||
PetscInt & | nTotal, | ||
type::RealVec1D & | dL | ||
) |
Parse the info of only one direction from YAML node.
axis | [in] YAML configuration node. |
dir | [out] Label of the direction. |
bg | [out] Starting coordinate. |
ed | [out] Ending coordinate. |
nTotal | [out] Number of cells along the direction. |
dL | [out] 1D vector with cell sizes along the direction. |
This function only parse information and do calculation for only one axis.
Definition at line 275 of file parser.cpp.
PetscErrorCode petibm::parser::parseOneSubDomain | ( | const YAML::Node & | sub, |
const PetscReal | bg, | ||
PetscInt & | n, | ||
PetscReal & | ed, | ||
type::RealVec1D & | dL | ||
) |
Parse only one sub-domain.
sub | [in] YAML configuration node. |
bg | [in] Starting coordinate. |
n | [out] Number of cells along the direction of the sub-domain. |
ed | [out] Ending coordinate. |
dL | [out] Grid cell sizes along the direction of the sub-domain. |
Definition at line 334 of file parser.cpp.
PetscErrorCode petibm::parser::parseSubDomains | ( | const YAML::Node & | subs, |
const PetscReal | bg, | ||
PetscInt & | nTotal, | ||
PetscReal & | ed, | ||
type::RealVec1D & | dL | ||
) |
Parse all sub-domains in a direction.
subs | [in] YAML configuration node. |
bg | [in] Starting coordinate. |
nTotal | [out] Number of cell along the direction. |
ed | [out] Ending coordinate. |
dL | [out] Grid cell sizes along the direction. |
This function parse the subDomains
section of a direction.
Definition at line 298 of file parser.cpp.
PetscReal petibm::delta::Peskin_2002 | ( | const PetscReal & | r, |
const PetscReal & | dr | ||
) |
PetscErrorCode petibm::io::print | ( | const std::string & | info | ) |
PetscErrorCode petibm::io::readHDF5Vecs | ( | const MPI_Comm | comm, |
const std::string & | filePath, | ||
const std::string & | loc, | ||
const std::vector< std::string > & | names, | ||
std::vector< Vec > & | vecs | ||
) |
Read a vector of Vec objects from a HDF5 file.
comm | [in] MPI communicator (should be the same as the one in Vecs). |
filePath | [in] Path of the file to read from. |
loc | [in] Location in the HDF5 file of the data to read. |
names | [in] Vector with the name of each Vec object to read. |
vecs | [out] Vector of Vec objects to fill. |
Note: this function doesn't check the length of the vector of Vec objects and of the vector of names.
PetscErrorCode petibm::io::readLagrangianPoints | ( | const std::string & | file, |
PetscInt & | nPts, | ||
type::RealVec2D & | coords | ||
) |
PetscErrorCode petibm::parser::readYAMLFile | ( | const std::string & | filePath, |
YAML::Node & | node | ||
) |
PetscReal petibm::delta::Roma_et_al_1999 | ( | const PetscReal & | r, |
const PetscReal & | dr | ||
) |
|
inline |
Calculate and return cell sizes of stretched grid in one direction.
bg | [in] start of the stretched region. |
ed | [in] end of the stretched region. |
n | [in] number of total cells in the stretched region. |
r | [in] value of the stretched ratio. |
dL | [out] cell sizes. |
PetscErrorCode petibm::io::writeHDF5Vecs | ( | const MPI_Comm | comm, |
const std::string & | filePath, | ||
const std::string & | loc, | ||
const std::vector< std::string > & | names, | ||
const std::vector< PetscInt > & | n, | ||
const std::vector< PetscReal * > & | vecs, | ||
const PetscFileMode | mode = FILE_MODE_WRITE |
||
) |
Write a vector of raw arrays to a HDF5 file.
comm | [in] MPI communicator. |
filePath | [in] Path of the file to write in. |
loc | [in] Location in the HDF5 file of data to write. |
names | [in] Vector with the name of each raw array to write. |
n | [in] Vector with the length of each raw array to write. |
vecs | [in] Vector of raw arrays to write. |
mode | [in] Either FILE_MODE_WRITE (default) or FILE_MODE_APPEND. |
Note: this function doesn't check the length of the vector of raw arrays and of the vector of names.
PetscErrorCode petibm::io::writeHDF5Vecs | ( | const MPI_Comm | comm, |
const std::string & | filePath, | ||
const std::string & | loc, | ||
const std::vector< std::string > & | names, | ||
const std::vector< Vec > & | vecs, | ||
const PetscFileMode | mode = FILE_MODE_WRITE |
||
) |
Write a vector of Vec objects to a HDF5 file.
comm | [in] MPI communicator (should be the same as the one in Vecs). |
filePath | [in] Path of the file to write in. |
loc | [in] Location in the HDF5 file for data to write. |
names | [in] Vector with the name of each Vec object to write. |
vecs | [in] Vector of Vec objects to write. |
mode | [in] Either FILE_MODE_WRITE (default) or FILE_MODE_APPEND. |
Note: this function doesn't check the length of the vector of Vec objects and of the vector of names.
PetscErrorCode petibm::io::writeHDF5Vecs | ( | const MPI_Comm | comm, |
const std::string & | filePath, | ||
const std::string & | loc, | ||
const std::vector< std::string > & | names, | ||
const type::RealVec2D & | vecs, | ||
const PetscFileMode | mode = FILE_MODE_WRITE |
||
) |
Write petibm::type::RealVec2D to a HDF5 file.
comm | [in] MPI communicator. |
filePath | [in] Path of the file to write in. |
loc | [in] Location in the HDF5 file of the data to write. |
names | [in] Vector with the name of each RealVec1D to write. |
vecs | [in] RealVec2D (i.e., std::vector<type::RealVec1D>) to write. |
mode | [in] Either FILE_MODE_WRITE (default) or FILE_MODE_APPEND. |
Note: this function doesn't check the length of the RealVec2D object and of the vector of names.