14 #if !defined(KRATOS_TWO_FLUID_NAVIER_STOKES_DATA_H)
15 #define KRATOS_TWO_FLUID_NAVIER_STOKES_DATA_H
32 template<
size_t TDim,
size_t TNumNodes >
132 const Vector& BDFVector = rProcessInfo[BDF_COEFFICIENTS];
147 for (
unsigned int i = 0;
i < TNumNodes;
i++){
177 noalias(this->DN_DXenr) = rDN_DXenr;
185 for (
unsigned int i = 0;
i < TNumNodes;
i++)
218 const double c1 = 2.0*
mu;
219 const double c2 =
mu;
228 if constexpr (TDim == 2)
231 const double volumetric_part = trace/2.0;
238 else if constexpr (TDim == 3)
241 const double volumetric_part = trace/3.0;
259 if constexpr (TDim == 3)
264 this->
StrainRate[3] =
DN(0,0)*
v(0,1) +
DN(0,1)*
v(0,0) +
DN(1,0)*
v(1,1) +
DN(1,1)*
v(1,0) +
DN(2,0)*
v(2,1) +
DN(2,1)*
v(2,0) +
DN(3,0)*
v(3,1) +
DN(3,1)*
v(3,0);
265 this->
StrainRate[4] =
DN(0,1)*
v(0,2) +
DN(0,2)*
v(0,1) +
DN(1,1)*
v(1,2) +
DN(1,2)*
v(1,1) +
DN(2,1)*
v(2,2) +
DN(2,2)*
v(2,1) +
DN(3,1)*
v(3,2) +
DN(3,2)*
v(3,1);
266 this->
StrainRate[5] =
DN(0,0)*
v(0,2) +
DN(0,2)*
v(0,0) +
DN(1,0)*
v(1,2) +
DN(1,2)*
v(1,0) +
DN(2,0)*
v(2,2) +
DN(2,2)*
v(2,0) +
DN(3,0)*
v(3,2) +
DN(3,2)*
v(3,0);
269 else if constexpr (TDim == 2)
273 this->
StrainRate[2] =
DN(0,1)*
v(0,0) +
DN(0,0)*
v(0,1) +
DN(1,1)*
v(1,0) +
DN(1,0)*
v(1,1) +
DN(2,1)*
v(2,0) +
DN(2,0)*
v(2,1);
279 double strain_rate_norm;
281 if constexpr (TDim == 3)
283 strain_rate_norm = std::sqrt(2.*
S[0] *
S[0] + 2.*
S[1] *
S[1] + 2.*
S[2] *
S[2] +
284 S[3] *
S[3] +
S[4] *
S[4] +
S[5] *
S[5]);
287 else if constexpr (TDim == 2)
289 strain_rate_norm = std::sqrt(2.*
S[0] *
S[0] + 2.*
S[1] *
S[1] +
S[2] *
S[2]);
291 return strain_rate_norm;
297 for (
unsigned int i = 0;
i < TNumNodes;
i++)
298 dist += this->
N[
i] * Distance[
i];
302 for (
unsigned int i = 0;
i < TNumNodes;
i++)
317 for (
unsigned int i = 0;
i < TNumNodes;
i++)
318 dist += this->
N[
i] * Distance[
i];
321 double dynamic_viscosity = 0.0;
322 for (
unsigned int i = 0;
i < TNumNodes;
i++)
337 length_scale *= length_scale;
346 for (
size_t i = 0;
i < TNumNodes;
i++) {
347 for (
size_t j = 0;
j < TDim;
j++) {
Base class for all Elements.
Definition: element.h:60
PropertiesType & GetProperties()
Definition: element.h:1024
static double GradientsElementSize(const BoundedMatrix< double, 3, 2 > &rDN_DX)
Element size based on the shape functions gradients. Triangle element version.
Definition: element_size_calculator.cpp:1456
Base class for data containers used within FluidElement and derived types.
Definition: fluid_element_data.h:37
void FillFromProperties(double &rData, const Variable< double > &rVariable, const Properties &rProperties)
Definition: fluid_element_data.cpp:192
void FillFromHistoricalNodalData(NodalScalarData &rData, const Variable< double > &rVariable, const Geometry< Node > &rGeometry)
Definition: fluid_element_data.cpp:65
virtual void Initialize(const Element &rElement, const ProcessInfo &rProcessInfo)
Definition: fluid_element_data.cpp:18
unsigned int IntegrationPointIndex
Definition: fluid_element_data.h:100
ShapeFunctionsType N
Definition: fluid_element_data.h:104
Vector StrainRate
Strain rate (symmetric gradient of velocity) vector in Voigt notation.
Definition: fluid_element_data.h:110
Vector ShearStress
Shear stress vector in Voigt notation.
Definition: fluid_element_data.h:114
virtual void UpdateGeometryValues(unsigned int IntegrationPointIndex, double NewWeight, const MatrixRowType &rN, const ShapeDerivativesType &rDN_DX)
Definition: fluid_element_data.cpp:52
double EffectiveViscosity
Effective viscosity (in dynamic units) produced by the constitutive law.
Definition: fluid_element_data.h:124
ShapeDerivativesType DN_DX
Definition: fluid_element_data.h:106
void FillFromProcessInfo(double &rData, const Variable< double > &rVariable, const ProcessInfo &rProcessInfo)
Definition: fluid_element_data.cpp:153
Matrix C
Constitutive tensor C (expressed as a Matrix).
Definition: fluid_element_data.h:118
static void GetNewtonianConstitutiveMatrix(const double DynamicViscosity, BoundedMatrix< double, VoigtVector2DSize, VoigtVector2DSize > &rConstitutiveMatrix)
Definition: fluid_element_utilities.cpp:58
GeometryType & GetGeometry()
Returns the reference of the geometry.
Definition: geometrical_object.h:158
Geometry base class.
Definition: geometry.h:71
Definition: amatrix_interface.h:41
void resize(std::size_t NewSize1, std::size_t NewSize2, bool preserve=0)
Definition: amatrix_interface.h:224
void clear()
Definition: amatrix_interface.h:284
static double Norm(const Vector &a)
Calculates the norm of vector "a".
Definition: math_utils.h:703
ProcessInfo holds the current value of different solution parameters.
Definition: process_info.h:59
Properties encapsulates data shared by different Elements or Conditions. It can store any type of dat...
Definition: properties.h:69
Definition: two_fluid_navier_stokes_data.h:34
bool IsCut()
Definition: two_fluid_navier_stokes_data.h:197
double DarcyTerm
Definition: two_fluid_navier_stokes_data.h:70
size_t NumNegativeNodes
Definition: two_fluid_navier_stokes_data.h:101
BoundedMatrix< double, TNumNodes, TNumNodes *(TDim+1)> H
Definition: two_fluid_navier_stokes_data.h:80
static int Check(const Element &rElement, const ProcessInfo &rProcessInfo)
Definition: two_fluid_navier_stokes_data.h:181
GeometryType::ShapeFunctionsGradientsType ShapeFunctionsGradientsType
Definition: two_fluid_navier_stokes_data.h:46
typename FluidElementData< TDim, TNumNodes, true >::MatrixRowType MatrixRowType
Definition: two_fluid_navier_stokes_data.h:44
ShapeFunctionsGradientsType DN_DX_pos_side
Definition: two_fluid_navier_stokes_data.h:88
array_1d< double, TNumNodes > rhs_ee
Definition: two_fluid_navier_stokes_data.h:82
double SmagorinskyConstant
Definition: two_fluid_navier_stokes_data.h:67
bool IsAir()
Definition: two_fluid_navier_stokes_data.h:201
double Density
Definition: two_fluid_navier_stokes_data.h:63
unsigned int NumberOfDivisions
Definition: two_fluid_navier_stokes_data.h:102
typename FluidElementData< TDim, TNumNodes, true >::NodalVectorData NodalVectorData
Definition: two_fluid_navier_stokes_data.h:41
Vector w_gauss_neg_side
Definition: two_fluid_navier_stokes_data.h:95
void UpdateGeometryValues(unsigned int IntegrationPointIndex, double NewWeight, const MatrixRowType &rN, const BoundedMatrix< double, TNumNodes, TDim > &rDN_DX, const MatrixRowType &rNenr, const BoundedMatrix< double, TNumNodes, TDim > &rDN_DXenr)
Definition: two_fluid_navier_stokes_data.h:166
ShapeFunctionsGradientsType DN_DX_neg_side
Definition: two_fluid_navier_stokes_data.h:89
void CalculateAirMaterialResponse()
Definition: two_fluid_navier_stokes_data.h:205
void CalculateDensityAtGaussPoint()
Definition: two_fluid_navier_stokes_data.h:294
void ComputeStrain()
Definition: two_fluid_navier_stokes_data.h:252
double ElementSize
Definition: two_fluid_navier_stokes_data.h:84
BoundedMatrix< double, TNumNodes *(TDim+1), TNumNodes > V
Definition: two_fluid_navier_stokes_data.h:79
NodalScalarData Distance
Definition: two_fluid_navier_stokes_data.h:59
Matrix N_pos_side
Definition: two_fluid_navier_stokes_data.h:86
Vector w_gauss_pos_side
Definition: two_fluid_navier_stokes_data.h:94
Matrix N_neg_side
Definition: two_fluid_navier_stokes_data.h:87
double bdf0
Definition: two_fluid_navier_stokes_data.h:72
typename FluidElementData< TDim, TNumNodes, true >::ShapeFunctionsType ShapeFunctionsType
Definition: two_fluid_navier_stokes_data.h:42
NodalScalarData Pressure
Definition: two_fluid_navier_stokes_data.h:58
Geometry< Node > GeometryType
Definition: two_fluid_navier_stokes_data.h:45
BoundedMatrix< double, TNumNodes, TNumNodes > Enr_Neg_Interp
Definition: two_fluid_navier_stokes_data.h:92
NodalScalarData NodalDensity
Definition: two_fluid_navier_stokes_data.h:60
BoundedMatrix< double, TNumNodes *(TDim+1), TNumNodes *(TDim+1)> lhs
Definition: two_fluid_navier_stokes_data.h:77
BoundedMatrix< double, TNumNodes, TNumNodes > Enr_Pos_Interp
Definition: two_fluid_navier_stokes_data.h:91
NodalVectorData MeshVelocity
Definition: two_fluid_navier_stokes_data.h:55
double ComputeStrainNorm()
Definition: two_fluid_navier_stokes_data.h:277
void CalculateEffectiveViscosityAtGaussPoint()
Definition: two_fluid_navier_stokes_data.h:314
double LinearDarcyCoefficient
Definition: two_fluid_navier_stokes_data.h:68
NodalVectorData Velocity_OldStep1
Definition: two_fluid_navier_stokes_data.h:53
void UpdateGeometryValues(unsigned int IntegrationPointIndex, double NewWeight, const MatrixRowType &rN, const BoundedMatrix< double, TNumNodes, TDim > &rDN_DX) override
Definition: two_fluid_navier_stokes_data.h:155
array_1d< double, TNumNodes *(TDim+1)> rhs
Definition: two_fluid_navier_stokes_data.h:78
double DynamicTau
Definition: two_fluid_navier_stokes_data.h:66
BoundedMatrix< double, TNumNodes, TNumNodes > Kee
Definition: two_fluid_navier_stokes_data.h:81
NodalScalarData NodalDynamicViscosity
Definition: two_fluid_navier_stokes_data.h:61
size_t NumPositiveNodes
Definition: two_fluid_navier_stokes_data.h:100
double bdf1
Definition: two_fluid_navier_stokes_data.h:73
NodalVectorData BodyForce
Definition: two_fluid_navier_stokes_data.h:56
NodalVectorData Velocity
Definition: two_fluid_navier_stokes_data.h:52
void Initialize(const Element &rElement, const ProcessInfo &rProcessInfo) override
Definition: two_fluid_navier_stokes_data.h:109
ShapeDerivativesType DN_DXenr
Definition: two_fluid_navier_stokes_data.h:98
NodalVectorData Velocity_OldStep2
Definition: two_fluid_navier_stokes_data.h:54
void ComputeDarcyTerm()
Definition: two_fluid_navier_stokes_data.h:343
typename FluidElementData< TDim, TNumNodes, true >::ShapeDerivativesType ShapeDerivativesType
Definition: two_fluid_navier_stokes_data.h:43
typename FluidElementData< TDim, TNumNodes, true >::NodalScalarData NodalScalarData
Definition: two_fluid_navier_stokes_data.h:40
double NonLinearDarcyCoefficient
Definition: two_fluid_navier_stokes_data.h:69
double VolumeError
Definition: two_fluid_navier_stokes_data.h:71
double DynamicViscosity
Definition: two_fluid_navier_stokes_data.h:64
double bdf2
Definition: two_fluid_navier_stokes_data.h:74
double DeltaTime
Definition: two_fluid_navier_stokes_data.h:65
ShapeFunctionsType Nenr
Definition: two_fluid_navier_stokes_data.h:97
Short class definition.
Definition: array_1d.h:61
#define KRATOS_CHECK_VARIABLE_IN_NODAL_DATA(TheVariable, TheNode)
Definition: checks.h:171
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
KratosZeroVector< double > ZeroVector
Definition: amatrix_interface.h:561
KratosZeroMatrix< double > ZeroMatrix
Definition: amatrix_interface.h:559
T & noalias(T &TheMatrix)
Definition: amatrix_interface.h:484
TABLE_NUMBER_ANGULAR_VELOCITY TABLE_NUMBER_MOMENT I33 BEAM_INERTIA_ROT_UNIT_LENGHT_Y KRATOS_DEFINE_APPLICATION_VARIABLE(DEM_APPLICATION, double, BEAM_INERTIA_ROT_UNIT_LENGHT_Z) typedef std double
Definition: DEM_application_variables.h:182
float dist
Definition: edgebased_PureConvection.py:89
float density
Definition: face_heat.py:56
v
Definition: generate_convection_diffusion_explicit_element.py:114
DN
Definition: generate_convection_diffusion_explicit_element.py:98
stress
Stress vector definition.
Definition: generate_droplet_dynamics.py:82
mu
Definition: generate_frictional_mortar_condition.py:127
int strain_size
Definition: generate_hyper_elastic_simo_taylor_neo_hookean.py:16
strain
HERE WE MUST SUBSTITUTE EXPRESSIONS.
Definition: generate_stokes_twofluid_element.py:104
S
Definition: generate_total_lagrangian_mixed_volumetric_strain_element.py:39
int j
Definition: quadrature.py:648
integer i
Definition: TensorModule.f:17