LifeV
|
OneDFSIData - Class which read and holds all the data for the One Dimensional Model Solver. More...
#include <OneDFSIData.hpp>
Type definitions | |
enum | OneD_distributionLaw { uniform, linear, pointwise } |
typedef TimeData | time_Type |
typedef std::shared_ptr< time_Type > | timePtr_Type |
typedef RegionMesh< LinearLine > | mesh_Type |
typedef std::shared_ptr< mesh_Type > | meshPtr_Type |
typedef ublas::vector< Real > | scalarVector_Type |
typedef std::array< Real, 2 > | container2D_Type |
Constructors & Destructor | |
OneDFSIData () | |
Empty constructor. More... | |
virtual | ~OneDFSIData () |
Destructor. More... | |
Methods | |
void | setup (const GetPot &dataFile, const std::string §ion="1D_Model") |
Setup method. More... | |
void | LIFEV_DEPRECATED (oldStyleSetup(const GetPot &dataFile, const std::string §ion="1dnetwork")) |
Deprecated setup method. (DO NOT USE THIS) More... | |
void | updateCoefficients () |
Update all the physical coefficients. More... | |
template<typename VectorType > | |
Real | computeSpatialDerivativeAtNode (const VectorType &vector, const UInt &iNode, const UInt &bcFiniteDifferenceOrder=2) |
Compute the spatial derivative of a quantity at a node. More... | |
void | showMe (std::ostream &output=std::cout) const |
Initialize linear parameters (NOT WORKING) More... | |
Set Methods | |
void | setPostprocessingDirectory (const std::string &directory) |
Set the post-processing directory. More... | |
void | setPostprocessingFile (const std::string &file) |
Set the post-processing file name. More... | |
void | setJacobianPerturbationArea (const Real &jacobianPerturbationArea) |
Set the area perturbation parameter to compute the Jacobian matrix (in the Multiscale framework) More... | |
void | setJacobianPerturbationFlowRate (const Real &jacobianPerturbationFlowRate) |
Set the flow rate perturbation parameter to compute the Jacobian matrix (in the Multiscale framework) More... | |
void | setJacobianPerturbationStress (const Real &jacobianPerturbationStress) |
Set the stress perturbation parameter to compute the Jacobian matrix (in the Multiscale framework) More... | |
void | setTimeData (const timePtr_Type timeDataPtr) |
Set data time container. More... | |
void | setDensity (const Real &density) |
Set the fluid density. More... | |
void | setViscosity (const Real &viscosity) |
Set the fluid viscosity. More... | |
void | setDensityWall (const Real &densityWall) |
Set the wall density. More... | |
void | setThickness (const Real &thickness, const UInt &i) |
Set the wall thickness. More... | |
void | setYoung (const Real &young) |
Set the wall Young modulus. More... | |
void | setPoisson (const Real &poisson) |
Set the wall Poisson number. More... | |
void | setExternalPressure (const Real &externalPressure) |
Set the wall external pressure. More... | |
void | setVenousPressure (const Real &venousPressure) |
Set the venous pressure at the terminals. More... | |
void | setArea0 (const Real &area0, const UInt &i) |
Set the reference area. More... | |
void | setBeta0 (const Real &beta0, const UInt &i) |
Set the wall . More... | |
void | setdBeta0dz (const Real &dBeta0dz, const UInt &i) |
Set the wall . More... | |
Get Methods | |
const OneDFSI::physicsType_Type & | physicsType () const |
Get the physics type. More... | |
const OneDFSI::fluxTerm_Type & | fluxType () const |
Get the flux type. More... | |
const OneDFSI::sourceTerm_Type & | sourceType () const |
Get the source type. More... | |
timePtr_Type | dataTime () const |
Get data time container. More... | |
meshPtr_Type | mesh () const |
Get the mesh container. More... | |
Real | length () const |
Get the length of the 1D segment. More... | |
UInt | numberOfElements () const |
Get the number of elements in the 1D segment. More... | |
UInt | numberOfNodes () const |
Get the number of nodes in the 1D segment. More... | |
const bool & | viscoelasticWall () const |
Get the flag identifying if the wall is viscoelastic. More... | |
const Real & | viscoelasticCoefficient (const UInt &i) const |
Get the viscoelastic coefficient . More... | |
const bool & | inertialWall () const |
Get the flag identifying if the wall has inertia. More... | |
const Real & | densityWall () const |
Get the density of the wall. More... | |
const Real & | inertialModulus () const |
Get the inertial coefficient (to be defined) More... | |
const bool & | longitudinalWall () const |
Get the flag identifying if the wall has a longitudinal pre-stress. More... | |
const std::string & | postprocessingDirectory () const |
Get the post-processing directory. More... | |
const std::string & | postprocessingFile () const |
Get the post-processing file. More... | |
const Real & | CFLmax () const |
Get the imposed CFL condition. More... | |
const Real & | jacobianPerturbationArea () const |
Get the area perturbation parameter to compute the Jacobian matrix (in the Multiscale framework) More... | |
const Real & | jacobianPerturbationFlowRate () const |
Get the flow rate perturbation parameter to compute the Jacobian matrix (in the Multiscale framework) More... | |
const Real & | jacobianPerturbationStress () const |
Get the stress perturbation parameter to compute the Jacobian matrix (in the Multiscale framework) More... | |
const Real & | densityRho () const |
Get the fluid density. More... | |
const Real & | viscosity () const |
Get the fluid viscosity. More... | |
const Real & | young () const |
Get the wall Young modulus. More... | |
const Real & | poisson () const |
Get the wall Poisson number. More... | |
const Real & | externalPressure () const |
Get the wall external pressure. More... | |
const Real & | venousPressure () const |
Get the venous pressure. More... | |
const Real & | robertsonCorrection () const |
Get the Robertson correction coefficient (Not tested: maybe wrong in the code) More... | |
const Real & | thickness (const UInt &i) const |
Get the wall thickness. More... | |
const Real & | friction () const |
Get the wall friction coefficient . More... | |
const Real & | area0 (const UInt &i) const |
Get the reference area . More... | |
const Real & | alpha (const UInt &i) const |
Get the Coriolis coefficient . More... | |
const Real & | beta0 (const UInt &i) const |
Get the coefficient. More... | |
const Real & | beta1 (const UInt &i) const |
Get the coefficient. More... | |
const Real & | dArea0dz (const UInt &i) const |
Get . More... | |
const Real & | dAlphadz (const UInt &i) const |
Get . More... | |
const Real & | dBeta0dz (const UInt &i) const |
Get . More... | |
const Real & | dBeta1dz (const UInt &i) const |
Get . More... | |
const Real & | flux11 (const UInt &i) const |
Get the flux coefficient (used only in the linear problem) More... | |
const Real & | flux12 (const UInt &i) const |
Get the flux coefficient (used only in the linear problem) More... | |
const Real & | flux21 (const UInt &i) const |
Get the flux coefficient (used only in the linear problem) More... | |
const Real & | flux22 (const UInt &i) const |
Get the flux coefficient (used only in the linear problem) More... | |
const Real & | celerity1 (const UInt &i) const |
Get the first eigenvector (used only in the linear problem) More... | |
const Real & | celerity2 (const UInt &i) const |
Get the second eigenvector (used only in the linear problem) More... | |
const Real & | leftEigenVector11 (const UInt &i) const |
Get the left eigenvector coefficient (used only in the linear problem) More... | |
const Real & | leftEigenVector12 (const UInt &i) const |
Get the left eigenvector coefficient (used only in the linear problem) More... | |
const Real & | leftEigenVector21 (const UInt &i) const |
Get the left eigenvector coefficient (used only in the linear problem) More... | |
const Real & | leftEigenVector22 (const UInt &i) const |
Get the left eigenvector coefficient (used only in the linear problem) More... | |
const Real & | source10 (const UInt &i) const |
Get the source coefficient (used only in the linear problem) More... | |
const Real & | source20 (const UInt &i) const |
Get the source coefficient (used only in the linear problem) More... | |
const Real & | source11 (const UInt &i) const |
Get the source coefficient (used only in the linear problem) More... | |
const Real & | source12 (const UInt &i) const |
Get the source coefficient (used only in the linear problem) More... | |
const Real & | source21 (const UInt &i) const |
Get the source coefficient (used only in the linear problem) More... | |
const Real & | source22 (const UInt &i) const |
Get the source coefficient (used only in the linear problem) More... | |
Private Methods | |
void | linearInterpolation (scalarVector_Type &vector, const GetPot &dataFile, const std::string &quantity, const Real &defaultValue, const bool &isArea=false) |
Compute the linear interpolation of a general quantity. More... | |
void | computeSpatialDerivatives () |
Compute the derivatives of alpha, area0, beta0, and beta1 using centered differences. More... | |
void | resetContainers () |
Reset all the containers. More... | |
OneDFSIData - Class which read and holds all the data for the One Dimensional Model Solver.
Physical Parameters
Main parameters: .
Euler equations:
with
Linear Parameters
Parameters:
Equations:
The flux matrix has the eigenvalues .
Definition at line 106 of file OneDFSIData.hpp.
Definition at line 113 of file OneDFSIData.hpp.
typedef std::shared_ptr< time_Type > timePtr_Type |
Definition at line 114 of file OneDFSIData.hpp.
typedef RegionMesh< LinearLine > mesh_Type |
Definition at line 116 of file OneDFSIData.hpp.
typedef std::shared_ptr< mesh_Type > meshPtr_Type |
Definition at line 117 of file OneDFSIData.hpp.
typedef ublas::vector< Real > scalarVector_Type |
Definition at line 120 of file OneDFSIData.hpp.
typedef std::array< Real, 2 > container2D_Type |
Definition at line 121 of file OneDFSIData.hpp.
enum OneD_distributionLaw |
Enumerator | |
---|---|
uniform | |
linear | |
pointwise |
Definition at line 123 of file OneDFSIData.hpp.
|
explicit |
Empty constructor.
Definition at line 51 of file OneDFSIData.cpp.
|
inlinevirtual |
Destructor.
Definition at line 140 of file OneDFSIData.hpp.
void setup | ( | const GetPot & | dataFile, |
const std::string & | section = "1D_Model" |
||
) |
Setup method.
dataFile | GetPot dataFile |
section | section in the dataFile |
Definition at line 114 of file OneDFSIData.cpp.
void LIFEV_DEPRECATED | ( | oldStyleSetup(const GetPot &dataFile, const std::string §ion="1dnetwork") | ) |
Deprecated setup method. (DO NOT USE THIS)
This method has been implemented for compatibility reason with some old applications. Please, for any new application use the setup method in place of this one. Future compatibility for this method is not guaranteed.
dataFile | GetPot dataFile |
void updateCoefficients | ( | ) |
Update all the physical coefficients.
This method can be called after any update in the physical coefficient. It recomputes the main coefficients
Definition at line 416 of file OneDFSIData.cpp.
|
inline |
Compute the spatial derivative of a quantity at a node.
Note: works only for homogeneous discretizations.
vector | the quantity vector |
iNode | node |
Definition at line 1022 of file OneDFSIData.hpp.
void showMe | ( | std::ostream & | output = std::cout | ) | const |
Initialize linear parameters (NOT WORKING)
The linearization of the Euler model yields
the vessel stiffer on the left side of interval [xl, xr] These routines change the elastic modulus along the vessel
When x < alpha - delta/2, the Young modulus is E * factor When x > alpha + delta/2, the Young modulus is E When alpha - delta/2 < x < alpha + delta/2, the Young modulus changes smoothly from the larger to the smaller value, according to a polynomial law of order n.
The grid size can be adapted (yesadaptive=1) in the nieghborhood of alpha, where the spatial derivative of the parameter will be maximum. However, the grid size is not allowed to be smaller than min_deltaxMake the vessel stiffer on the right side of interval [xl, xr]
output | Stream where the informations must be printed |
Definition at line 691 of file OneDFSIData.cpp.
|
inline |
Set the post-processing directory.
directory | post-processing directory |
Definition at line 245 of file OneDFSIData.hpp.
|
inline |
Set the post-processing file name.
file | post-processing file name |
Definition at line 254 of file OneDFSIData.hpp.
|
inline |
Set the area perturbation parameter to compute the Jacobian matrix (in the Multiscale framework)
jacobianPerturbationArea | area perturbation parameter |
Definition at line 263 of file OneDFSIData.hpp.
|
inline |
Set the flow rate perturbation parameter to compute the Jacobian matrix (in the Multiscale framework)
jacobianPerturbationFlowRate | flow rate perturbation parameter |
Definition at line 272 of file OneDFSIData.hpp.
|
inline |
Set the stress perturbation parameter to compute the Jacobian matrix (in the Multiscale framework)
jacobianPerturbationStress | stress perturbation parameter |
Definition at line 281 of file OneDFSIData.hpp.
|
inline |
Set data time container.
timeDataPtr | shared_ptr to TimeData container |
Definition at line 290 of file OneDFSIData.hpp.
|
inline |
Set the fluid density.
density | fluid density |
Definition at line 299 of file OneDFSIData.hpp.
|
inline |
Set the fluid viscosity.
viscosity | fluid viscosity |
Definition at line 308 of file OneDFSIData.hpp.
|
inline |
Set the wall density.
densityWall | wall density |
Definition at line 317 of file OneDFSIData.hpp.
Set the wall thickness.
i | node id |
thickness | wall thickness |
Definition at line 327 of file OneDFSIData.hpp.
|
inline |
Set the wall Young modulus.
young | wall Young modulus |
Definition at line 336 of file OneDFSIData.hpp.
|
inline |
Set the wall Poisson number.
poisson | wall Poisson number |
Definition at line 345 of file OneDFSIData.hpp.
|
inline |
Set the wall external pressure.
externalPressure | wall external pressure |
Definition at line 354 of file OneDFSIData.hpp.
|
inline |
Set the venous pressure at the terminals.
venousPressure | venous pressure at the terminals |
Definition at line 363 of file OneDFSIData.hpp.
Set the reference area.
i | node id |
area0 | reference area |
Definition at line 373 of file OneDFSIData.hpp.
|
inline |
|
inline |
|
inline |
|
inline |
Get data time container.
Definition at line 435 of file OneDFSIData.hpp.
|
inline |
Get the mesh container.
Definition at line 444 of file OneDFSIData.hpp.
|
inline |
Get the length of the 1D segment.
Definition at line 453 of file OneDFSIData.hpp.
|
inline |
Get the number of elements in the 1D segment.
Definition at line 462 of file OneDFSIData.hpp.
|
inline |
Get the number of nodes in the 1D segment.
Definition at line 471 of file OneDFSIData.hpp.
|
inline |
Get the flag identifying if the wall is viscoelastic.
Definition at line 480 of file OneDFSIData.hpp.
Get the viscoelastic coefficient .
Definition at line 489 of file OneDFSIData.hpp.
|
inline |
Get the flag identifying if the wall has inertia.
Definition at line 498 of file OneDFSIData.hpp.
|
inline |
Get the density of the wall.
Definition at line 507 of file OneDFSIData.hpp.
|
inline |
Get the inertial coefficient (to be defined)
Definition at line 516 of file OneDFSIData.hpp.
|
inline |
Get the flag identifying if the wall has a longitudinal pre-stress.
Definition at line 525 of file OneDFSIData.hpp.
|
inline |
Get the post-processing directory.
Definition at line 534 of file OneDFSIData.hpp.
|
inline |
Get the post-processing file.
Definition at line 543 of file OneDFSIData.hpp.
|
inline |
Get the imposed CFL condition.
Definition at line 552 of file OneDFSIData.hpp.
|
inline |
Get the area perturbation parameter to compute the Jacobian matrix (in the Multiscale framework)
Definition at line 561 of file OneDFSIData.hpp.
|
inline |
Get the flow rate perturbation parameter to compute the Jacobian matrix (in the Multiscale framework)
Definition at line 570 of file OneDFSIData.hpp.
|
inline |
Get the stress perturbation parameter to compute the Jacobian matrix (in the Multiscale framework)
Definition at line 579 of file OneDFSIData.hpp.
|
inline |
|
inline |
|
inline |
Get the wall Young modulus.
Definition at line 606 of file OneDFSIData.hpp.
|
inline |
Get the wall Poisson number.
Definition at line 615 of file OneDFSIData.hpp.
|
inline |
Get the wall external pressure.
Definition at line 624 of file OneDFSIData.hpp.
|
inline |
const Real & robertsonCorrection | ( | ) | const |
Get the Robertson correction coefficient (Not tested: maybe wrong in the code)
Definition at line 782 of file OneDFSIData.cpp.
|
inline |
Get the wall friction coefficient .
Definition at line 657 of file OneDFSIData.hpp.
Get the reference area .
i | node id |
Definition at line 667 of file OneDFSIData.hpp.
Get the Coriolis coefficient .
i | node id |
Definition at line 677 of file OneDFSIData.hpp.
Get the coefficient.
i | node id |
Definition at line 687 of file OneDFSIData.hpp.
Get the coefficient.
i | node id |
Definition at line 697 of file OneDFSIData.hpp.
Get the flux coefficient (used only in the linear problem)
i | node id |
Definition at line 747 of file OneDFSIData.hpp.
Get the flux coefficient (used only in the linear problem)
i | node id |
Definition at line 757 of file OneDFSIData.hpp.
Get the flux coefficient (used only in the linear problem)
i | node id |
Definition at line 767 of file OneDFSIData.hpp.
Get the flux coefficient (used only in the linear problem)
i | node id |
Definition at line 777 of file OneDFSIData.hpp.
Get the first eigenvector (used only in the linear problem)
i | node id |
Definition at line 787 of file OneDFSIData.hpp.
Get the second eigenvector (used only in the linear problem)
i | node id |
Definition at line 797 of file OneDFSIData.hpp.
Get the left eigenvector coefficient (used only in the linear problem)
i | node id |
Definition at line 807 of file OneDFSIData.hpp.
Get the left eigenvector coefficient (used only in the linear problem)
i | node id |
Definition at line 817 of file OneDFSIData.hpp.
Get the left eigenvector coefficient (used only in the linear problem)
i | node id |
Definition at line 827 of file OneDFSIData.hpp.
Get the left eigenvector coefficient (used only in the linear problem)
i | node id |
Definition at line 837 of file OneDFSIData.hpp.
Get the source coefficient (used only in the linear problem)
i | node id |
Definition at line 847 of file OneDFSIData.hpp.
Get the source coefficient (used only in the linear problem)
i | node id |
Definition at line 857 of file OneDFSIData.hpp.
Get the source coefficient (used only in the linear problem)
i | node id |
Definition at line 867 of file OneDFSIData.hpp.
Get the source coefficient (used only in the linear problem)
i | node id |
Definition at line 877 of file OneDFSIData.hpp.
Get the source coefficient (used only in the linear problem)
i | node id |
Definition at line 887 of file OneDFSIData.hpp.
Get the source coefficient (used only in the linear problem)
i | node id |
Definition at line 897 of file OneDFSIData.hpp.
|
private |
Compute the linear interpolation of a general quantity.
Very useful for tapering.
vector | interpolated vector |
dataFile | data file |
quantity | quantity |
defaultValue | default value |
isArea | flag identifying if the vector is an area (the linear interpolation is done on the radius) |
Definition at line 795 of file OneDFSIData.cpp.
|
private |
Compute the derivatives of alpha, area0, beta0, and beta1 using centered differences.
Definition at line 828 of file OneDFSIData.cpp.
|
private |
Reset all the containers.
Definition at line 840 of file OneDFSIData.cpp.
|
private |
Model.
Definition at line 933 of file OneDFSIData.hpp.
|
private |
Definition at line 934 of file OneDFSIData.hpp.
|
private |
Definition at line 935 of file OneDFSIData.hpp.
|
private |
Data containers for time and mesh.
Definition at line 938 of file OneDFSIData.hpp.
|
private |
Definition at line 939 of file OneDFSIData.hpp.
|
private |
Physical Wall Model.
Definition at line 942 of file OneDFSIData.hpp.
|
private |
Definition at line 943 of file OneDFSIData.hpp.
|
private |
Definition at line 944 of file OneDFSIData.hpp.
|
private |
Definition at line 945 of file OneDFSIData.hpp.
|
private |
Definition at line 947 of file OneDFSIData.hpp.
|
private |
Definition at line 948 of file OneDFSIData.hpp.
|
private |
Definition at line 949 of file OneDFSIData.hpp.
|
private |
Definition at line 950 of file OneDFSIData.hpp.
|
private |
Miscellaneous.
Definition at line 953 of file OneDFSIData.hpp.
|
private |
full directory name (including path)
Definition at line 954 of file OneDFSIData.hpp.
|
private |
output file name
Definition at line 955 of file OneDFSIData.hpp.
|
private |
Jacobian perturbation.
Definition at line 958 of file OneDFSIData.hpp.
|
private |
Definition at line 959 of file OneDFSIData.hpp.
|
private |
Definition at line 960 of file OneDFSIData.hpp.
|
private |
Physical Parameters.
Definition at line 963 of file OneDFSIData.hpp.
|
private |
Definition at line 964 of file OneDFSIData.hpp.
|
private |
Definition at line 966 of file OneDFSIData.hpp.
|
private |
Definition at line 967 of file OneDFSIData.hpp.
|
private |
Definition at line 969 of file OneDFSIData.hpp.
|
private |
Definition at line 970 of file OneDFSIData.hpp.
|
private |
Definition at line 971 of file OneDFSIData.hpp.
|
private |
Definition at line 973 of file OneDFSIData.hpp.
|
private |
Definition at line 974 of file OneDFSIData.hpp.
|
private |
Definition at line 975 of file OneDFSIData.hpp.
|
private |
Definition at line 977 of file OneDFSIData.hpp.
|
private |
Definition at line 978 of file OneDFSIData.hpp.
|
private |
Definition at line 980 of file OneDFSIData.hpp.
|
private |
Definition at line 981 of file OneDFSIData.hpp.
|
private |
Definition at line 982 of file OneDFSIData.hpp.
|
private |
Definition at line 983 of file OneDFSIData.hpp.
|
private |
Derivatives of main coefficients.
Definition at line 986 of file OneDFSIData.hpp.
|
private |
Definition at line 987 of file OneDFSIData.hpp.
|
private |
Definition at line 988 of file OneDFSIData.hpp.
|
private |
Definition at line 989 of file OneDFSIData.hpp.
|
private |
Flux matrix.
Definition at line 992 of file OneDFSIData.hpp.
|
private |
Definition at line 993 of file OneDFSIData.hpp.
|
private |
Definition at line 994 of file OneDFSIData.hpp.
|
private |
Definition at line 995 of file OneDFSIData.hpp.
|
private |
Celerities of the linear problem (eigenvalues of the flux matrix)
Definition at line 998 of file OneDFSIData.hpp.
|
private |
Definition at line 999 of file OneDFSIData.hpp.
|
private |
Eigenvector for first and second eigenvalue.
Definition at line 1002 of file OneDFSIData.hpp.
|
private |
Definition at line 1003 of file OneDFSIData.hpp.
|
private |
Definition at line 1004 of file OneDFSIData.hpp.
|
private |
Definition at line 1005 of file OneDFSIData.hpp.
|
private |
Source matrix.
Definition at line 1008 of file OneDFSIData.hpp.
|
private |
Definition at line 1009 of file OneDFSIData.hpp.
|
private |
Definition at line 1010 of file OneDFSIData.hpp.
|
private |
Definition at line 1011 of file OneDFSIData.hpp.
|
private |
Definition at line 1012 of file OneDFSIData.hpp.
|
private |
Definition at line 1013 of file OneDFSIData.hpp.