OGS
FractureProperty.h
Go to the documentation of this file.
1
11#pragma once
12
13#include <Eigen/Core>
14#include <memory>
15
16#include "BranchProperty.h"
17#include "JunctionProperty.h"
21
22namespace ParameterLib
23{
24template <typename T>
25struct Parameter;
26}
27
28namespace ProcessLib
29{
30namespace LIE
31{
33{
34 int fracture_id = 0;
35 int mat_id = 0;
36 Eigen::Vector3d point_on_fracture;
37 Eigen::Vector3d normal_vector;
39 Eigen::MatrixXd R;
42 std::vector<BranchProperty> branches_master;
43 std::vector<BranchProperty> branches_slave;
44
45 FractureProperty(int const fracture_id_, int const material_id,
46 ParameterLib::Parameter<double> const& initial_aperture)
47 : fracture_id(fracture_id_),
48 mat_id(material_id),
49 aperture0(initial_aperture)
50 {
51 }
52
53 virtual ~FractureProperty() = default;
54};
55
57{
58 FracturePropertyHM(int const fracture_id_, int const material_id,
59 ParameterLib::Parameter<double> const& initial_aperture,
60 ParameterLib::Parameter<double> const& specific_storage_,
61 ParameterLib::Parameter<double> const& biot_coefficient_)
62 : FractureProperty(fracture_id_, material_id, initial_aperture),
63 specific_storage(specific_storage_),
64 biot_coefficient(biot_coefficient_)
65 {
66 }
69
70 std::unique_ptr<MaterialLib::Fracture::Permeability::Permeability>
72};
73
76inline void setFractureProperty(int const dim, MeshLib::Element const& e,
77 FractureProperty& frac_prop)
78{
79 auto& n = frac_prop.normal_vector;
80 // 1st node is used but using other node is also possible, because
81 // a fracture is not curving
82 for (int j = 0; j < 3; j++)
83 {
84 frac_prop.point_on_fracture[j] = getCenterOfGravity(e).data()[j];
85 }
86
87 const MeshLib::ElementCoordinatesMappingLocal ele_local_coord(e, dim);
88
89 // Global to local rotation matrix:
90 Eigen::MatrixXd const global2local_rotation =
91 ele_local_coord.getRotationMatrixToGlobal().transpose();
92 n = global2local_rotation.row(dim - 1);
93
94 frac_prop.R = global2local_rotation.topLeftCorner(dim, dim);
95
96 DBUG("Normal vector of the fracture element {:d}: [{:g}, {:g}, {:g}]",
97 e.getID(), n[0], n[1], n[2]);
98}
99
101 FractureProperty const& master_frac,
102 FractureProperty const& slave_frac)
103{
104 BranchProperty branch{branchNode, master_frac.fracture_id,
105 slave_frac.fracture_id};
106
107 // set a normal vector from the master to the slave fracture
108 Eigen::Vector3d branch_vector =
109 slave_frac.point_on_fracture - branch.coords;
110 double sign = (branch_vector.dot(master_frac.normal_vector) < 0) ? -1 : 1;
111 branch.normal_vector_branch = sign * master_frac.normal_vector;
112 return branch;
113}
114} // namespace LIE
115} // namespace ProcessLib
Definition of the Element class.
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:30
const double * data() const
Definition Point3d.h:60
const RotationMatrix & getRotationMatrixToGlobal() const
return a rotation matrix converting to global coordinates
std::size_t getID() const
Returns the ID of the element.
Definition Element.h:89
void setFractureProperty(int const dim, MeshLib::Element const &e, FractureProperty &frac_prop)
BranchProperty createBranchProperty(MeshLib::Node const &branchNode, FractureProperty const &master_frac, FractureProperty const &slave_frac)
ParameterLib::Parameter< double > const & biot_coefficient
std::unique_ptr< MaterialLib::Fracture::Permeability::Permeability > permeability_model
FracturePropertyHM(int const fracture_id_, int const material_id, ParameterLib::Parameter< double > const &initial_aperture, ParameterLib::Parameter< double > const &specific_storage_, ParameterLib::Parameter< double > const &biot_coefficient_)
ParameterLib::Parameter< double > const & specific_storage
virtual ~FractureProperty()=default
ParameterLib::Parameter< double > const & aperture0
Initial aperture.
std::vector< BranchProperty > branches_slave
FractureProperty(int const fracture_id_, int const material_id, ParameterLib::Parameter< double > const &initial_aperture)
Eigen::MatrixXd R
Rotation matrix from global to local coordinates.
std::vector< BranchProperty > branches_master