/
HydroMechanicsLocalAssemblerFracture.h
125 lines (104 loc) · 4.29 KB
/
HydroMechanicsLocalAssemblerFracture.h
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
/**
* \file
* \copyright
* Copyright (c) 2012-2019, OpenGeoSys Community (http://www.opengeosys.org)
* Distributed under a Modified BSD License.
* See accompanying file LICENSE.txt or
* http://www.opengeosys.org/project/license
*
*/
#pragma once
#include <vector>
#include "NumLib/Fem/ShapeMatrixPolicy.h"
#include "ProcessLib/LIE/Common/HMatrixUtils.h"
#include "ProcessLib/LIE/HydroMechanics/HydroMechanicsProcessData.h"
#include "HydroMechanicsLocalAssemblerInterface.h"
#include "IntegrationPointDataFracture.h"
namespace ProcessLib
{
namespace LIE
{
namespace HydroMechanics
{
template <typename ShapeFunctionDisplacement,
typename ShapeFunctionPressure,
typename IntegrationMethod,
int GlobalDim>
class HydroMechanicsLocalAssemblerFracture
: public HydroMechanicsLocalAssemblerInterface
{
public:
HydroMechanicsLocalAssemblerFracture(
HydroMechanicsLocalAssemblerFracture const&) = delete;
HydroMechanicsLocalAssemblerFracture(
HydroMechanicsLocalAssemblerFracture&&) = delete;
HydroMechanicsLocalAssemblerFracture(
MeshLib::Element const& e,
std::size_t const local_matrix_size,
std::vector<unsigned> const& dofIndex_to_localIndex,
bool const is_axially_symmetric,
unsigned const integration_order,
HydroMechanicsProcessData<GlobalDim>& process_data);
void preTimestepConcrete(std::vector<double> const& /*local_x*/,
double const /*t*/,
double const /*delta_t*/) override
{
for (auto& data : _ip_data)
{
data.pushBackState();
}
}
void postTimestepConcreteWithVector(
const double t, double const dt,
Eigen::VectorXd const& local_x) override;
Eigen::Map<const Eigen::RowVectorXd> getShapeMatrix(
const unsigned integration_point) const override
{
auto const& N = _ip_data[integration_point].N_p;
// assumes N is stored contiguously in memory
return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
}
private:
void assembleWithJacobianConcrete(double const t, double const dt,
Eigen::VectorXd const& local_x,
Eigen::VectorXd const& local_xdot,
Eigen::VectorXd& local_b,
Eigen::MatrixXd& local_J) override;
void assembleBlockMatricesWithJacobian(
double const t, double const dt,
Eigen::Ref<const Eigen::VectorXd> const& p,
Eigen::Ref<const Eigen::VectorXd> const& p_dot,
Eigen::Ref<const Eigen::VectorXd> const& g,
Eigen::Ref<const Eigen::VectorXd> const& g_dot,
Eigen::Ref<Eigen::VectorXd> rhs_p, Eigen::Ref<Eigen::VectorXd> rhs_g,
Eigen::Ref<Eigen::MatrixXd> J_pp, Eigen::Ref<Eigen::MatrixXd> J_pg,
Eigen::Ref<Eigen::MatrixXd> J_gg, Eigen::Ref<Eigen::MatrixXd> J_gp);
// Types for displacement.
using ShapeMatricesTypeDisplacement =
ShapeMatrixPolicyType<ShapeFunctionDisplacement, GlobalDim>;
using HMatricesType =
HMatrixPolicyType<ShapeFunctionDisplacement, GlobalDim>;
using HMatrixType = typename HMatricesType::HMatrixType;
// Types for pressure.
using ShapeMatricesTypePressure =
ShapeMatrixPolicyType<ShapeFunctionPressure, GlobalDim>;
// Types for the integration point data
using IntegrationPointDataType =
IntegrationPointDataFracture<HMatricesType,
ShapeMatricesTypeDisplacement,
ShapeMatricesTypePressure,
GlobalDim>;
HydroMechanicsProcessData<GlobalDim>& _process_data;
std::vector<IntegrationPointDataType,
Eigen::aligned_allocator<IntegrationPointDataType>>
_ip_data;
static const int pressure_index = 0;
static const int pressure_size = ShapeFunctionPressure::NPOINTS;
static const int displacement_index = ShapeFunctionPressure::NPOINTS;
static const int displacement_size =
ShapeFunctionDisplacement::NPOINTS * GlobalDim;
};
} // namespace HydroMechanics
} // namespace LIE
} // namespace ProcessLib
#include "HydroMechanicsLocalAssemblerFracture-impl.h"