tsunami_lab
-
typedef std::size_t tsunami_lab::t_idx
integral type for cell-ids, pointer arithmetic, etc.
-
typedef amrex_real tsunami_lab::t_real
floating point type
-
template<typename T>
T *tsunami_lab::aligned_alloc(T *&rawPtr, size_t size, size_t alignment = alignof(T)) Returns an aligned array zero initalized
- Template Parameters:
T – type of the array
- Parameters:
rawPtr – pointer to the nonaligned array for deletion
size – size of the array to allocated
alignment – the alignment of the array
- Returns:
the aligned array or on fail a nullptr
amr
-
class AMRCoreWavePropagation2d : public amrex::AmrCore
- #include <AMRCoreWavePropagation2d.h>
Two-dimensional adaptiv multi resolution wave propagation
Public Types
Public Functions
-
AMRCoreWavePropagation2d(tsunami_lab::setups::Setup *setup)
Construct Wave Propagation with AMR using the given setup
- Parameters:
setup – the setup to initialize the data with
-
void PrintParameters()
Print some parameters to the console
-
void setTimeStep(amrex::Real timeStep, int level)
Set the time step of the simulation
- Parameters:
timeStep – set the timeStep to step forward
level – the level to set
-
void Evolve()
Advance solution to final time
-
void ErrorEst(int level, amrex::TagBoxArray &tags, amrex::Real time, int ngrow)
Tag cells for refinement. TagBoxArray tags is built on level grids.
- Parameters:
level – the level to estimate and tag
tags – the box that holds the tag values
time – the current time
ngrow – the ghost cells need for estimation calculation
-
void MakeNewLevelFromScratch(int level, amrex::Real time, const amrex::BoxArray &ba, const amrex::DistributionMapping &dm)
Make a new level from scratch using provided BoxArray and DistributionMapping. Only used during initialization.
- Parameters:
level – the level to make from scratch
time – the current time
ba – the provided BoxArray
dm – the provided DistributionMapping
-
void MakeNewLevelFromCoarse(int level, amrex::Real time, const amrex::BoxArray &ba, const amrex::DistributionMapping &dm)
Make a new level using provided BoxArray and DistributionMapping and fill with interpolated coarse level data.
- Parameters:
level – the level to make
time – the current time
ba – the provided BoxArray
dm – the provided DistributionMapping
-
void RemakeLevel(int level, amrex::Real time, const amrex::BoxArray &ba, const amrex::DistributionMapping &dm)
Remake an existing level using provided BoxArray and DistributionMapping and fill with existing fine and coarse data.
- Parameters:
level – the level to Remake
time – the current time
ba – the provided BoxArray
dm – the provided DistributionMapping
-
void ClearLevel(int level)
Delete level data
- Parameters:
level – the level to delete
Private Functions
-
void FillFinePatch(int level, amrex::Real time, amrex::MultiFab &mf)
Fill an entire multifab by interpolating from the coarser level This comes into play when a new level of refinement appears
- Parameters:
level – the level to fill
time – the current time
mf – the multiFab to interpolate to
-
void FixFinePatch(amrex::MultiFab &mf, const amrex::MultiFab &const_mf)
Fix an entire multifab that was interpolating from the coarser level. This comes into play when the fine level was created or updated. This will file the cell near the shore i.e. |bathymetry| < bathymetryMinValue with the values of const_mf instead of using the mf and set any height on the coast to zero.
- Parameters:
mf – the multiFab to interpolate to
const_mf – the multiFab wth piecewise constant bathymetry interpolation
-
void GetData(int level, amrex::Real time, amrex::Vector<amrex::MultiFab*> &data, amrex::Vector<amrex::Real> &datatime)
Utility to copy in data from gridOld and/or gridNew into another multifab
- Parameters:
level – the level to get data from
time – the current time
data – the MutliFab to which the data is written
datatime – the Vector to which the time is written
-
void FillPatch(int level, amrex::Real time, amrex::MultiFab &mf)
Fill a patch with data from the grid
- Parameters:
level – the level to fill
time – the current time
mf – the multifab to fill
-
void timeStepWithSubcycling(int level, amrex::Real time, int iteration)
Advance a level by dt - includes a recursive call for finer levels
- Parameters:
level – the level to sub cycle
time – the current time
iteration – the current iteration step
-
void AdvanceGridAtLevel(int level, amrex::Real time, amrex::Real dtLevel, int iteration, int nCycle)
Advance grid at a single level for a single time step
- Parameters:
level – the level to advance
time – the current time
dtLevel – the time step of this level
iteration – the current iteration
the – current cycle step
-
void WritePlotFile() const
Write plotfile to disk
-
void AverageDownTo(int coarseLevel)
More flexible version of AverageDown() that lets you average down across multiple levels
- Parameters:
coarseLevel – the level to average down to
-
void ReadParameters()
Read in parameters from input file
-
void InitData(int level)
Read the data from the setup into the grid
- Parameters:
level – the level to setup
Private Members
-
const int nComponents = 5
number of components i.e. Height, MomentumX, MomentumY, Bathymetry
-
const int nGhostRow = 1
number of ghost cell around the boundary of the domain
-
amrex::Real bathymetryMinValue = 20
Minimum Height of the bathymetry.
-
amrex::Interpolater *interpolator = &amrex::lincc_interp
interpolator going from coarse to fine
-
amrex::Vector<int> step
which step?
-
amrex::Vector<int> nSubSteps
how many substeps on each level?
-
amrex::Vector<amrex::Real> tNew
keep track of new time step at each level
-
amrex::Vector<amrex::Real> tOld
keep track of old time step at each level
-
amrex::Vector<amrex::Real> dt
keep track of time step at each level
-
amrex::Vector<amrex::MultiFab> gridNew
array of multifabs to store the solution at each level of refinement after advancing a level we use “swap”
-
amrex::Vector<amrex::MultiFab> gridOld
array of multifabs to store the solution at each level of refinement after advancing a level we use “swap”
-
amrex::Vector<amrex::BCRec> physicalBoundary
this is essentially a 2*DIM integer array storing the physical boundary condition types at the lo/hi walls in each direction 4-components: Height, MomentumX, MomentumY, Bathymetry
-
amrex::Real simulationTime = std::numeric_limits<amrex::Real>::max()
time to simulate
-
int regridFrequency = 2
how often each level regrids the higher levels of refinement (after a level advances that many time steps)
-
std::string plotFile = {"plt"}
root name of plot file
-
std::string plotFolder = {"solution"}
folder for plot files
-
int writeFrequency = -1
frequency to write the output
-
amrex::Vector<amrex::Real> gridErr
the error bound for each level
-
tsunami_lab::setups::Setup *setup
the setup to init the data with
-
AMRCoreWavePropagation2d(tsunami_lab::setups::Setup *setup)
io
-
class Csv
- #include <Csv.h>
IO-routines for writing a snapshot as Comma Separated Values (CSV).
Public Static Functions
-
static void write(t_real i_dxy, t_idx i_nx, t_idx i_ny, t_idx i_stride, t_real const *i_h, t_real const *i_hu, t_real const *i_hv, t_real const *i_b, t_real const *i_hTotal, std::ostream &io_stream)
Writes the data as CSV to the given stream.
- Parameters:
i_dxy – cell width in x- and y-direction.
i_nx – number of cells in x-direction.
i_ny – number of cells in y-direction.
i_stride – stride of the data arrays in y-direction (x is assumed to be stride-1).
i_h – water height of the cells; optional: use nullptr if not required.
i_hu – momentum in x-direction of the cells; optional: use nullptr if not required.
i_hv – momentum in y-direction of the cells; optional: use nullptr if not required.
i_b – the bathymetry at x-position; optional: use nullptr if not required.
i_hTotal – the sum of water height and bathymetry; optional: use nullptr if not required.
io_stream – stream to which the CSV-data is written.
-
static bool next_middle_states(std::ifstream &stream, t_real &o_hLeft, t_real &o_hRight, t_real &o_huLeft, t_real &o_huRight, t_real &o_hStar)
gets the next parsed value pair from the middle_state.csv file stream
- Parameters:
stream – file stream of middle_state.csv
o_hLeft – output the height left
o_hRight – output the height right
o_huLeft – output the momentum left
o_huRight – output the momentum right
o_hStar – output the computed height
- Returns:
success
-
static void write(t_real i_dxy, t_idx i_nx, t_idx i_ny, t_idx i_stride, t_real const *i_h, t_real const *i_hu, t_real const *i_hv, t_real const *i_b, t_real const *i_hTotal, std::ostream &io_stream)
-
class NetCdf
- #include <NetCdf.h>
Supports either reading or writing to a netCDF file
Public Types
-
enum VarType
The data type parsed from a netCDF file. The given options contain the same value as the nc_type e.g. VarType::CHAR == NC_CHAR.
Values:
-
enumerator NONE
-
enumerator CHAR
-
enumerator SHORT
-
enumerator INT
-
enumerator FLOAT
-
enumerator DOUBLE
-
enumerator UCHAR
-
enumerator USHORT
-
enumerator UINT
-
enumerator INT64
-
enumerator UINT64
-
enumerator STRING
-
enumerator NONE
Public Functions
-
NetCdf(t_idx writeStep, std::string filePath, t_idx l_nx, t_idx l_ny, t_idx l_k, t_real l_scaleX, t_real l_scaleY, t_idx l_stride, const t_real *bathymetry, bool useSpherical)
Write-Only Constructor of NetCdf for the continuation from a checkpoint.
- Parameters:
writeStep – the writeStep of the checkpoint
filePath – filepath of the netCDF file
l_nx – len in x dimension
l_ny – len in y dimension
l_k – number of cells to average several neighboring cells of the computational grid into one cell
l_scaleX – the scale in x direction in meters
l_scaleY – the scale in y direction in meters
l_stride – the stride of the data-set to write
bathymetry – bathymetry data to write if no bathymetry should be written pass a nullptr
useSpherical – use spherical Longitude & Latitude in degrees for the X and Y Axis instead of meters
-
NetCdf(std::string filePath, t_idx l_nx, t_idx l_ny, t_idx l_k, t_real l_scaleX, t_real l_scaleY, t_idx l_stride, const t_real *bathymetry, bool useSpherical, bool useMomenta)
Write-Only Constructor of NetCdf.
- Parameters:
filePath – filepath of the netCDF file
l_nx – len in x dimension
l_ny – len in y dimension
l_k – number of cells to average several neighboring cells of the computational grid into one cell
l_scaleX – the scale in x direction in meters
l_scaleY – the scale in y direction in meters
l_stride – the stride of the data-set to write
bathymetry – bathymetry data to write if no bathymetry should be written pass a nullptr
useSpherical – use spherical Longitude & Latitude in degrees for the X and Y Axis instead of meters
useMomenta – if true also create variables for momentumX and momentumY and enable writing to these
-
NetCdf(std::string filePath, t_idx l_nx, t_idx l_ny, t_real l_scaleX, t_real l_scaleY, t_idx l_stride, const t_real *bathymetry, const char *commandLine, t_real hMax, const t_real *totalHeight, const t_real *momentumX, const t_real *momentumY, t_real simulationTime, const t_idx writeCount)
Write-Only Constructor of NetCdf for Checkpoints.
- Parameters:
filePath – filepath of the netCDF file
l_nx – len in x dimension
l_ny – len in y dimension
l_scaleX – the scale in x direction in meters
l_scaleY – the scale in y direction in meters
l_stride – the stride of the data-set to write
bathymetry – bathymetry data to write
commandLine – the current input of the commandLine as a string. Is not used if the file is not a checkpoint.
hMax – the current hMax. Is not used if the file is not a checkpoint.
totalHeight – total heights of cells
momentumX – momentum of cells in x direction
momentumY – momentum of cells in y direction
simulationTime – the current simulation time in seconds
writeCount – the current simulation internal write count
-
void averageSeveral(const t_real simulationTime, const t_real *totalHeight, const t_real *momentumX, const t_real *momentumY)
Averages the input arrays so that l_k cells are combined into one cell.
- Parameters:
simulationTime – the current simulation time in seconds
totalHeight – total heights of cells
momentumX – momentum of cells in x direction
momentumY – momentum of cells in y direction
-
void write(const t_real simulationTime, const t_real *totalHeight, const t_real *momentumX, const t_real *momentumY)
Write current time step to a netCDF file.
- Parameters:
simulationTime – the current simulation time in seconds
totalHeight – total heights of cells
momentumX – momentum of cells in x direction
momentumY – momentum of cells in y direction
-
template<size_t N>
inline void read(const char *filepath, const char *(&variableName)[N], VarArray (&outData)[N], long long int timeStep = 0) read the data from the file at the given filepath and parse all given variables. the outData array contains all the needed information (array, type, length, stride) to work with the data. If the timeStep does not matter or should be ignored set it to zero.
IMPORTANT the size of variableName and outData should match N
- Template Parameters:
N – number of variables to parse
- Parameters:
filepath – the path to the netCDF file
variableName – the names of the variables to parse
outData – the parsed data for each variable with the type, length and stride. If the VarArray is destroyed the read data ‘VarArray::array’ is destroyed to.
timeStep – OPTIONAL start of data parsing. The timeStep uses the dimension directly therefore index values are required e.g. 0, 1, 2, … . Can also be negativ to read from the back.
-
void read(const char *filepath, const char *variableName, VarArray &outData, long long int timeStep = 0)
read the data from the file at the given filepath and parse the one given variable. If the timeStep does not matter or should be ignored set it to zero.
- Parameters:
filepath – filepath the path to the netCDF file
variableName – the names of the variable to parse
outData – the parsed data with the type, length and stride. If the VarArray is destroyed the read data ‘VarArray::array’ is destroyed to.
timeStep – OPTIONAL start of data parsing. The timeStep uses the dimension directly therefore index values are required e.g. 0, 1, 2, … . Can also be negativ to read from the back.
Public Static Functions
-
static void checkNcErr(int i_err, std::string text)
Check if command has worked.
- Parameters:
i_err – 0 if everything is ok, else 1
text – text to log
Private Functions
-
void _read(const char *filepath, const char **variableName, VarArray *outData, long long int timeStep, size_t size)
read the data from the file at the given filepath and parse all given variables. the outData array contains all the needed information (array, type, length, stride) to work with the data If the timeStep does not matter or should be ignored set it to zero.
- Parameters:
filepath – the path to the netCDF file
variableName – the names of the variables to parse
outData – the parsed data for each variable with the type, length and stride. If the VarArray is destroyed the read data ‘VarArray::array’ is destroyed to.
timeStep – start of data parsing. The timeStep uses the dimension directly therefore index values are required e.g. 0, 1, 2, … . Can also be negativ to read from the back.
size – the size of outData and variableName
-
void _write(const t_real simulationTime, const t_real *totalHeight, const t_real *momentumX, const t_real *momentumY, const t_idx nx, const t_idx ny, const t_idx stride, const t_idx writeCount)
Write current time step to a netCDF file.
- Parameters:
simulationTime – the current simulation time in seconds
totalHeight – total heights of cells
momentumX – momentum of cells in x direction
momentumY – momentum of cells in y direction
nx – len in x dimension
ny – len in y dimension
stride – the stride of the data-set to write
writeCount – the current simulation internal write count (only used for checkpoint)
-
NetCdf(std::string filePath, t_idx l_nx, t_idx l_ny, t_idx l_k, t_real l_scaleX, t_real l_scaleY, t_idx l_stride, const t_real *bathymetry, bool useSpherical, bool useMomenta, const char *commandLine, t_real hMax)
Write-Only Constructor of NetCdf.
- Parameters:
filePath – filepath of the netCDF file
l_nx – len in x dimension
l_ny – len in y dimension
l_k – number of cells to average several neighboring cells of the computational grid into one cell
l_scaleX – the scale in x direction in meters
l_scaleY – the scale in y direction in meters
l_stride – the stride of the data-set to write
bathymetry – bathymetry data to write if no bathymetry should be written pass a nullptr
useSpherical – use spherical Longitude & Latitude in degrees for the X and Y Axis instead of meters
useMomenta – if true also create variables for momentumX and momentumY and enable writing to these
commandLine – the current input of the commandLine as a string and if empty the writer does not write checkpoint data. Is not used if the file is not a checkpoint.
hMax – the current hMax. Is not used if the file is not a checkpoint.
Private Members
-
bool isReadMode
indicates whether the netCdf is in read or write mode
-
std::string m_filePath = ""
name of the netCDF file
-
t_idx m_k = 1
number of cells to average several neighboring cells of the computational grid into one cell
-
size_t m_time = 0
time of write operation
-
int m_ncId = -1
id of the netCDF file
-
int m_dimTimeId = -1
id of time dimension
-
int m_dimXId = -1
id of x dimension
-
int m_xId = -1
id of longitude
-
int m_dimYId = -1
id of y dimension
-
int m_yId = -1
id of latitude
-
int m_timeId = -1
id of time
-
int m_totalHeightId = -1
id of total height
-
int m_bathymetryId = -1
id of bathymetry
-
int m_momentumXId = -1
id of momentumX
-
int m_momentumYId = -1
id of momentumY
-
int m_writeCountId = -1
id of writeCount
-
int m_hMaxID = -1
id of hMax
-
bool isCheckpoint = false
indicates whether this is a checkpoint
-
const char *commandLine = ""
the current command line only used if this is a checkpoint
-
struct VarArray
- #include <NetCdf.h>
struct that stores the array/data with the corresponding type, length and stride
-
enum VarType
setups
-
class ArtificialTsunami2d : public tsunami_lab::setups::Setup
- #include <ArtificialTsunami2d.h>
artificial tsunami setup. This Setup only provides one hard-coded example.
Public Functions
-
ArtificialTsunami2d()
Default constructor with the default example.
-
virtual t_real getHeight(t_real, t_real) const
Gets the water height at a given point.
- Returns:
height at the given point.
-
virtual t_real getMomentumX(t_real, t_real) const
Gets the momentum in x-direction.
- Returns:
momentum in x-direction.
-
ArtificialTsunami2d()
-
class CircularDamBreak2d : public tsunami_lab::setups::Setup
- #include <CircularDamBreak2d.h>
circular dam break setup.
Public Functions
-
CircularDamBreak2d()
Default constructor with the default example.
-
CircularDamBreak2d(t_real i_heightCenter, t_real i_heightOutside, t_real i_locationCenter[2], t_real i_scaleCenter)
Constructor.
- Parameters:
i_heightCenter – water height inside of the circular dam.
i_heightOutside – water height outside of the circular dam.
i_locationCenter – location (x-coordinate, y-coordinate) of the circular dam.
i_scaleCenter – radius of the circular dam
-
virtual t_real getHeight(t_real i_x, t_real i_y) const
Gets the water height at a given point.
- Parameters:
i_x – x-coordinate of the queried point.
i_y – y-coordinate of the queried point.
- Returns:
height at the given point.
-
virtual t_real getMomentumX(t_real, t_real) const
Gets the momentum in x-direction.
- Returns:
momentum in x-direction.
-
CircularDamBreak2d()
-
class DamBreak1d : public tsunami_lab::setups::Setup
- #include <DamBreak1d.h>
1d dam break setup.
Public Functions
-
DamBreak1d(t_real i_heightLeft, t_real i_heightRight, t_real i_locationDam)
Constructor.
- Parameters:
i_heightLeft – water height on the left side of the dam.
i_heightRight – water height on the right side of the dam.
i_locationDam – location (x-coordinate) of the dam.
-
virtual t_real getHeight(t_real i_x, t_real) const
Gets the water height at a given point.
- Parameters:
i_x – x-coordinate of the queried point.
- Returns:
height at the given point.
-
virtual t_real getMomentumX(t_real, t_real) const
Gets the momentum in x-direction.
- Returns:
momentum in x-direction.
-
DamBreak1d(t_real i_heightLeft, t_real i_heightRight, t_real i_locationDam)
-
class MiddleStates1d : public tsunami_lab::setups::Setup
- #include <MiddleStates1d.h>
1d meddle-states setup.
Public Functions
-
MiddleStates1d(t_real i_heightLeft, t_real i_heightRight, t_real i_momentumLeft, t_real i_momentumRight, t_real i_location)
Construct a new middleStates 1d object.
- Parameters:
i_heightLeft – water height of left side of the rare location.
i_heightRight – momentum of the water of the right side.
i_momentumLeft – momentum of the water of the left side.
i_momentumRight – momentum of the water of the right side.
i_location – location (x-coordinate) of the middle state.
-
inline ~MiddleStates1d()
Destroy the middleStates 1d object.
-
virtual t_real getHeight(t_real i_x, t_real) const
Get the ware height at a given point.
- Parameters:
i_x – x-coordinate of the queried point.
- Returns:
height at the given point.
-
virtual t_real getMomentumX(t_real i_x, t_real) const
Get the momentum in x-direction at a given point.
- Parameters:
i_x – x-coordinate of the queried point.
- Returns:
momentum in x-direction at a given point.
-
MiddleStates1d(t_real i_heightLeft, t_real i_heightRight, t_real i_momentumLeft, t_real i_momentumRight, t_real i_location)
-
class RareRare1d : public tsunami_lab::setups::Setup
- #include <RareRare1d.h>
1d rare-rare setup.
Public Functions
-
RareRare1d(t_real i_heightLeft, t_real i_momentumLeft, t_real i_locationRare)
Construct a new rare-rare 1d object.
- Parameters:
i_heightLeft – water height of left side of the rare location.
i_momentumLeft – momentum of the water of the left side.
i_locationRare – location (x-coordinate) of the rare.
-
inline ~RareRare1d()
Destroy the rare-rare 1d object.
-
virtual t_real getHeight(t_real, t_real) const
Get the ware height at a given point.
- Returns:
height at the given point.
-
virtual t_real getMomentumX(t_real i_x, t_real) const
Get the momentum in x-direction at a given point.
- Parameters:
i_x – x-coordinate of the queried point.
- Returns:
momentum in x-direction at a given point.
-
RareRare1d(t_real i_heightLeft, t_real i_momentumLeft, t_real i_locationRare)
-
class Setup
- #include <Setup.h>
Base setup.
Subclassed by tsunami_lab::setups::ArtificialTsunami2d, tsunami_lab::setups::CircularDamBreak2d, tsunami_lab::setups::DamBreak1d, tsunami_lab::setups::MiddleStates1d, tsunami_lab::setups::RareRare1d, tsunami_lab::setups::ShockShock1d, tsunami_lab::setups::SubcriticalFlow1d, tsunami_lab::setups::SupercriticalFlow1d, tsunami_lab::setups::TsunamiEvent1d, tsunami_lab::setups::TsunamiEvent2d
Public Functions
-
inline virtual ~Setup()
Virtual destructor for base class.
-
virtual t_real getHeight(t_real i_x, t_real i_y) const = 0
Gets the water height at a given point.
- Parameters:
i_x – x-coordinate of the queried point.
i_y – y-coordinate of the queried point.
- Returns:
water height at the given point.
-
virtual t_real getMomentumX(t_real i_x, t_real i_y) const = 0
Gets the momentum in x-direction.
- Parameters:
i_x – x-coordinate of the queried point.
i_y – y-coordinate of the queried point.
- Returns:
momentum in x-direction.
-
inline virtual ~Setup()
-
class ShockShock1d : public tsunami_lab::setups::Setup
- #include <ShockShock1d.h>
1d shock-shock setup.
Public Functions
-
ShockShock1d(t_real i_heightLeft, t_real i_momentumLeft, t_real i_locationShock)
Construct a new shock-shock 1d object.
- Parameters:
i_heightLeft – water height on the left side of the shock location.
i_momentumLeft – momentum of the water of the left side.
i_locationShock – location (x-coordinate) of the shock.
-
inline ~ShockShock1d()
Destroy the shock-shock 1d object.
-
virtual t_real getHeight(t_real, t_real) const
Get the water height at a given point.
- Returns:
height at the given point.
-
virtual t_real getMomentumX(t_real i_x, t_real) const
Get the wate momentum in x-direction at a given point.
- Parameters:
i_x – x-coordinate of the queried point.
- Returns:
momentum in x-direction.
-
ShockShock1d(t_real i_heightLeft, t_real i_momentumLeft, t_real i_locationShock)
-
class SubcriticalFlow1d : public tsunami_lab::setups::Setup
- #include <SubcriticalFlow1d.h>
setup for subcritical flow
Public Functions
-
SubcriticalFlow1d()
Default supercritical flow example
-
SubcriticalFlow1d(t_real momentum, t_real range[2], t_real bathymetryOutRange, t_real (*bathymetryInRange)(t_real))
Customized supercritical flow example
- Parameters:
momentum – the fixed momentum of the water
range – the range were the bathymetry is created by the bathymetryInRange function, range[0] < range[1]
bathymetryOutRange – fixed bathymetry outside of the range
bathymetryInRange – function for creating dynamic bathymetry to first parameter is the x-coordinate of the bathymetry i.e. the values between range[0] and range[1]
-
~SubcriticalFlow1d()
default destructor
-
virtual t_real getHeight(t_real i_x, t_real) const override
Gets the water height at a given point.
- Parameters:
i_x – x-coordinate of the queried point.
- Returns:
water height at the given point.
-
virtual t_real getMomentumX(t_real, t_real) const override
Gets the momentum in x-direction.
- Returns:
momentum in x-direction.
Private Members
-
SubcriticalFlow1d()
-
class SupercriticalFlow1d : public tsunami_lab::setups::Setup
- #include <SupercriticalFlow1d.h>
setup for supercritical flow
Public Functions
-
SupercriticalFlow1d()
Default supercritical flow example
-
SupercriticalFlow1d(t_real momentum, t_real range[2], t_real bathymetryOutRange, t_real (*bathymetryInRange)(t_real))
Customized supercritical flow example
- Parameters:
momentum – the fixed momentum of the water
range – the range were the bathymetry is created by the bathymetryInRange function, range[0] < range[1]
bathymetryOutRange – fixed bathymetry outside of the range
bathymetryInRange – function for creating dynamic bathymetry to first parameter is the x-coordinate of the bathymetry i.e. the values between range[0] and range[1]
-
~SupercriticalFlow1d()
default destructor
-
virtual t_real getHeight(t_real i_x, t_real) const override
Gets the water height at a given point.
- Parameters:
i_x – x-coordinate of the queried point.
- Returns:
water height at the given point.
-
virtual t_real getMomentumX(t_real, t_real) const override
Gets the momentum in x-direction.
- Returns:
momentum in x-direction.
Private Members
-
SupercriticalFlow1d()
-
class TsunamiEvent1d : public tsunami_lab::setups::Setup
- #include <TsunamiEvent1d.h>
TsunamiEvent1d setup.
Public Functions
-
TsunamiEvent1d(std::string i_filePath, t_real i_delta = 20, t_real i_scale = 10)
Constructor.
- Parameters:
i_filePath – string to the bathymetry_profile.csv
i_delta – avoids running into numerical issues due to missing support for wetting and drying in our solver
i_scale – width of grid
-
virtual t_real getHeight(t_real i_x, t_real) const
Gets the water height at a given point.
- Parameters:
i_x – x-coordinate of the queried point.
- Returns:
height at the given point.
-
virtual t_real getMomentumX(t_real, t_real) const
Gets the momentum in x-direction.
- Returns:
momentum in x-direction.
-
virtual t_real getMomentumY(t_real, t_real) const
Gets the momentum in y-direction.
- Returns:
momentum in y-direction.
-
TsunamiEvent1d(std::string i_filePath, t_real i_delta = 20, t_real i_scale = 10)
-
class TsunamiEvent2d : public tsunami_lab::setups::Setup
- #include <TsunamiEvent2d.h>
TsunamiEvent2d setup.
Public Functions
-
TsunamiEvent2d(const char *bathymetryFilePath, const char *(&bathymetryVariable)[3], const char *displacementFilePath, const char *(&displacementVariable)[3], t_real scaleX, t_real scaleY, t_real delta = 20)
Creates the 2d tsunami event that reads two files. The first should be the bathymetry file and the second one should be the displacement of the bathymetry. Variables must be defined for each file, which are then read from the netCDF file in order to read the data successfully.
- Parameters:
bathymetryFilePath – filepath to the bathymetry file
bathymetryVariable – netCDF variables to read from the bathymetry file. The first should be the X dimension, the second the Y dimension and the third the bathymetry data
displacementFilePath – filepath to the displacement file
displacementVariable – netCDF variables to read from the displacement file. The first should be the X dimension, the second the Y dimension and the third the displacement data
scaleX – width of the grid
scaleY – height of the grid
delta – avoids running into numerical issues due to missing support for wetting and drying in our solver
-
inline ~TsunamiEvent2d()
Destroys the 2d tsunami event.
-
virtual t_real getHeight(t_real i_x, t_real i_y) const
Gets the water height at a given point.
- Parameters:
i_x – x-coordinate of the queried point.
i_y – y-coordinate of the queried point.
- Returns:
height at the given point.
-
virtual t_real getMomentumX(t_real, t_real) const
Gets the momentum in x-direction.
- Returns:
momentum in x-direction.
Private Functions
-
t_real getValueAscending(const t_real *const data[3], const t_idx size[3], const t_idx zStride, const t_real x, const t_real y) const
get the closest data point from an ascending sorted array
- Parameters:
data – data arrays where index 0 is x, 1 is y, 2 is z/data
size – the size of the data arrays where index 0 is x, 1 is y, 2 is z/data
zStride – the stride of the third/z array
x – the x-coordinate to obtain the closest value
y – the y-coordinate to obtain the closest value
- Returns:
the closest value found or 0 if no lower not be found or calculated index is out of range
Private Members
-
tsunami_lab::io::NetCdf::VarArray bathymetryData[3]
the bathymetry data from the file where index 0 is x, 1 is y, 2 is z/data
-
t_real *bathymetry[3] = {nullptr, nullptr, nullptr}
the the converted bathymetry data where the array is stored.
-
tsunami_lab::io::NetCdf::VarArray displacementData[3]
the vertical displacement data from the file where index 0 is x, 1 is y, 2 is z/data
-
TsunamiEvent2d(const char *bathymetryFilePath, const char *(&bathymetryVariable)[3], const char *displacementFilePath, const char *(&displacementVariable)[3], t_real scaleX, t_real scaleY, t_real delta = 20)
solvers
-
class FWave
- #include <FWave.h>
The f-wave solver for solving the Initial Value Problem (IVP) for the shallow water equations.
Public Static Functions
-
static void netUpdates(t_real i_hL, t_real i_hR, t_real i_huL, t_real i_huR, t_real o_netUpdate[2][2])
Computes the net-updates without bathymetry.
- Parameters:
i_hL – height of the left side.
i_hR – height of the right side.
i_huL – momentum of the left side.
i_huR – momentum of the right side.
o_netUpdate – will be set to the net-updates for the 0: left 1: right side; 0: height, 1: momentum.
-
static void netUpdates(t_real i_hL, t_real i_hR, t_real i_huL, t_real i_huR, t_real i_bL, t_real i_bR, t_real o_netUpdate[2][2])
Computes the net-updates with bathymetry.
- Parameters:
i_hL – height of the left side.
i_hR – height of the right side.
i_huL – momentum of the left side.
i_huR – momentum of the right side.
i_bL – height of bathymetry of the left side.
i_bR – height of bathymetry of the right side.
o_netUpdate – will be set to the net-updates for the 0: left 1: right side; 0: height, 1: momentum.
Private Static Functions
-
static void computeEigenvalues(t_real i_hL, t_real i_hR, t_real i_uL, t_real i_uR, t_real &o_eigenvalue1, t_real &o_eigenvalue2)
Computes the eigenvalues for the left and right wave.
- Parameters:
i_hL – height of the left side.
i_hR – height of the right side.
i_uL – particle velocity of the left side.
i_uR – particles velocity of the right side.
o_eigenvalue1 – output: Roe eigenvalue one.
o_eigenvalue2 – output: Roe eigenvalue two.
-
static void computeDeltaFlux(t_real i_hL, t_real i_hR, t_real i_uL, t_real i_uR, t_real i_huL, t_real i_huR, t_real o_deltaFlux[2])
Computes difference of the flux function after inserting the left and right quantities.
- Parameters:
i_hL – height of the left side.
i_hR – height of the right side.
i_uL – particle velocity of the left side.
i_uR – particles velocity of the right side.
i_huL – momentum of the left side.
i_huR – momentum of the right side.
o_deltaFlux – output: difference of left and right quantities
-
static void computeEigencoefficients(t_real i_eigenvalue1, t_real i_eigenvalue2, t_real i_deltaFlux[2], t_real &o_eigencoefficient1, t_real &o_eigencoefficient2)
Compute the eigencoefficients for left and right wave.
- Parameters:
i_eigenvalue1 – Roe eigenvalue one.
i_eigenvalue2 – Roe eigenvalue two.
i_deltaFlux – the compute delta Flux with the same inputs used for the eigenvalues
o_eigencoefficient1 – ouput: the eigencoefficient ot eigenvalue one.
o_eigencoefficient2 – ouput: the eigencoefficient ot eigenvalue two.
-
static void computeBathymetryEffects(t_real i_hL, t_real i_hR, t_real i_bL, t_real i_bR, t_real o_bathymetryEffect[2])
- Parameters:
i_hL – height of the left side.
i_hR – height of the right side.
i_bL – height of bathymetry of the left side.
i_bR – height of bathymetry of the right side.
o_bathymetryEffect – output: effect of the bathymetry
-
static void netUpdates(t_real i_hL, t_real i_hR, t_real i_huL, t_real i_huR, t_real o_netUpdate[2][2])