v0.14.0
Loading...
Searching...
No Matches
Public Member Functions | Public Attributes | List of all members
UFOperators2D::BlockData Struct Reference

#include <users_modules/softmech/unsaturated_2Dflow/src/Unsaturated2DFlowUDOP.hpp>

Collaboration diagram for UFOperators2D::BlockData:
[legend]

Public Member Functions

 BlockData ()
 
MoFEMErrorCode setOptions ()
 
template<typename T >
T get_waterContent (T head)
 
template<typename T >
T get_conductivity (T head, double head_norm)
 
template<typename T >
T get_relativeConductivity (T head)
 
template<typename T >
T get_Fe (T eff_satu)
 
template<typename T >
T get_eff_satuStar (T eff_satu)
 
template<typename T >
T get_effSaturation (T head)
 
template<typename T >
T get_capacity (T head, double head_norm)
 
double get_diffCapacity (double head)
 
double get_diffConductivity (double head)
 

Public Attributes

int block_id
 
double K_s
 
double ll
 
double theta_r
 
double theta_s
 
double theta_m
 
double nn
 
double alpha
 
double h_s
 
Range block_ents
 

Detailed Description

Definition at line 26 of file Unsaturated2DFlowUDOP.hpp.

Constructor & Destructor Documentation

◆ BlockData()

UFOperators2D::BlockData::BlockData ( )
inline

Member Function Documentation

◆ get_capacity()

template<typename T >
T UFOperators2D::BlockData::get_capacity ( T  head,
double  head_norm 
)
inline

Definition at line 102 of file Unsaturated2DFlowUDOP.hpp.

102 {
103 const T m = 1 - 1.0 / nn;
104 T ret_val;
105 if (head_norm < h_s) {
106 ret_val = alpha * m * nn * (theta_m - theta_r) *
107 pow(-alpha * head, nn - 1) /
108 pow(1 + pow(-alpha * head, nn), m + 1);
109 } else {
110 ret_val = 0;
111 }
112 return ret_val;
113 }
FTensor::Index< 'm', SPACE_DIM > m
const double T

◆ get_conductivity()

template<typename T >
T UFOperators2D::BlockData::get_conductivity ( T  head,
double  head_norm 
)
inline

Definition at line 71 of file Unsaturated2DFlowUDOP.hpp.

71 {
72 T ret_val;
73 if (head_norm < h_s) {
74 ret_val = K_s * get_relativeConductivity(head);
75 } else {
76 ret_val = K_s;
77 }
78 return ret_val;
79 }

◆ get_diffCapacity()

double UFOperators2D::BlockData::get_diffCapacity ( double  head)
inline

Definition at line 115 of file Unsaturated2DFlowUDOP.hpp.

115 {
116 std::complex<double> cx_head = head + eps * 1i;
117 auto capacity = get_capacity(cx_head, head);
118 return capacity.imag() / eps;
119 }
FTensor::Index< 'i', 3 > i
constexpr double eps
T get_capacity(T head, double head_norm)

◆ get_diffConductivity()

double UFOperators2D::BlockData::get_diffConductivity ( double  head)
inline

Definition at line 121 of file Unsaturated2DFlowUDOP.hpp.

121 {
122 std::complex<double> cx_head = head + eps * 1i;
123 auto conductivity = get_conductivity(cx_head, head);
124 return conductivity.imag() / eps;
125 }
T get_conductivity(T head, double head_norm)

◆ get_eff_satuStar()

template<typename T >
T UFOperators2D::BlockData::get_eff_satuStar ( T  eff_satu)
inline

Definition at line 93 of file Unsaturated2DFlowUDOP.hpp.

93 {
94 return (theta_s - theta_r) / (theta_m - theta_r) * eff_satu;
95 }

◆ get_effSaturation()

template<typename T >
T UFOperators2D::BlockData::get_effSaturation ( T  head)
inline

Definition at line 97 of file Unsaturated2DFlowUDOP.hpp.

97 {
98 const T theta = get_waterContent(head);
99 return (theta - theta_r) / (theta_s - theta_r);
100 }

◆ get_Fe()

template<typename T >
T UFOperators2D::BlockData::get_Fe ( T  eff_satu)
inline

Definition at line 88 of file Unsaturated2DFlowUDOP.hpp.

88 {
89 const T m = 1 - 1.0 / nn;
90 const T S_eStar = get_eff_satuStar(eff_satu);
91 return pow(1 - pow(S_eStar, 1.0 / m), m);
92 }

◆ get_relativeConductivity()

template<typename T >
T UFOperators2D::BlockData::get_relativeConductivity ( T  head)
inline

Definition at line 81 of file Unsaturated2DFlowUDOP.hpp.

81 {
82 const T S_e = get_effSaturation(head);
83 const T F_e = get_Fe(S_e);
84 const T F_1 = get_Fe(1.0);
85 return pow(S_e, ll) * pow(((1 - F_e) / (1 - F_1)), 2);
86 }

◆ get_waterContent()

template<typename T >
T UFOperators2D::BlockData::get_waterContent ( T  head)
inline

Definition at line 64 of file Unsaturated2DFlowUDOP.hpp.

64 {
65 T ret_val;
66 T m = 1 - 1.0 / nn;
67 return ret_val = theta_r +
68 (theta_m - theta_r) / pow(1 + pow(-alpha * head, nn), m);
69 }

◆ setOptions()

MoFEMErrorCode UFOperators2D::BlockData::setOptions ( )
inline

Definition at line 43 of file Unsaturated2DFlowUDOP.hpp.

43 {
45 ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "", "Material options", "none");
46 CHKERR PetscOptionsScalar("-K_s", "K_s", "", K_s, &K_s, PETSC_NULL);
47 CHKERR PetscOptionsScalar("-saturated_conductivity", "K_s", "", K_s, &K_s,
48 PETSC_NULL);
49 CHKERR PetscOptionsScalar("-l", "l", "", ll, &ll, PETSC_NULL);
50 CHKERR PetscOptionsScalar("-theta_r", "theta_r", "", theta_r, &theta_r,
51 PETSC_NULL);
52 CHKERR PetscOptionsScalar("-theta_s", "theta_s", "", theta_s, &theta_s,
53 PETSC_NULL);
54 CHKERR PetscOptionsScalar("-theta_m", "theta_m", "", theta_m, &theta_m,
55 PETSC_NULL);
56 CHKERR PetscOptionsScalar("-n", "n", "", nn, &nn, PETSC_NULL);
57 CHKERR PetscOptionsScalar("-alpha", "alpha", "", alpha, &alpha, PETSC_NULL);
58 CHKERR PetscOptionsScalar("-h_s", "h_s", "", h_s, &h_s, PETSC_NULL);
59 ierr = PetscOptionsEnd();
62 };
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:483
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76

Member Data Documentation

◆ alpha

double UFOperators2D::BlockData::alpha

Definition at line 34 of file Unsaturated2DFlowUDOP.hpp.

◆ block_ents

Range UFOperators2D::BlockData::block_ents

Definition at line 37 of file Unsaturated2DFlowUDOP.hpp.

◆ block_id

int UFOperators2D::BlockData::block_id

Definition at line 27 of file Unsaturated2DFlowUDOP.hpp.

◆ h_s

double UFOperators2D::BlockData::h_s

Definition at line 35 of file Unsaturated2DFlowUDOP.hpp.

◆ K_s

double UFOperators2D::BlockData::K_s

Definition at line 28 of file Unsaturated2DFlowUDOP.hpp.

◆ ll

double UFOperators2D::BlockData::ll

Definition at line 29 of file Unsaturated2DFlowUDOP.hpp.

◆ nn

double UFOperators2D::BlockData::nn

Definition at line 33 of file Unsaturated2DFlowUDOP.hpp.

◆ theta_m

double UFOperators2D::BlockData::theta_m

Definition at line 32 of file Unsaturated2DFlowUDOP.hpp.

◆ theta_r

double UFOperators2D::BlockData::theta_r

Definition at line 30 of file Unsaturated2DFlowUDOP.hpp.

◆ theta_s

double UFOperators2D::BlockData::theta_s

Definition at line 31 of file Unsaturated2DFlowUDOP.hpp.


The documentation for this struct was generated from the following file: