OGS
TH2MProcess.cpp
Go to the documentation of this file.
1
11#include "TH2MProcess.h"
12
13#include <cassert>
14
17#include "ProcessLib/Process.h"
19#include "TH2MFEM.h"
20#include "TH2MProcessData.h"
21
22namespace ProcessLib
23{
24namespace TH2M
25{
26template <int DisplacementDim>
28 std::string name,
29 MeshLib::Mesh& mesh,
30 std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&& jacobian_assembler,
31 std::vector<std::unique_ptr<ParameterLib::ParameterBase>> const& parameters,
32 unsigned const integration_order,
33 std::vector<std::vector<std::reference_wrapper<ProcessVariable>>>&&
34 process_variables,
36 SecondaryVariableCollection&& secondary_variables,
37 bool const use_monolithic_scheme)
38 : Process(std::move(name), mesh, std::move(jacobian_assembler), parameters,
39 integration_order, std::move(process_variables),
40 std::move(secondary_variables), use_monolithic_scheme),
41 _process_data(std::move(process_data))
42{
43 _nodal_forces = MeshLib::getOrCreateMeshProperty<double>(
44 mesh, "NodalForces", MeshLib::MeshItemType::Node, DisplacementDim);
45
46 _hydraulic_flow = MeshLib::getOrCreateMeshProperty<double>(
47 mesh, "MassFlowRate", MeshLib::MeshItemType::Node, 1);
48
49 // TODO (naumov) remove ip suffix. Probably needs modification of the mesh
50 // properties, s.t. there is no "overlapping" with cell/point data.
51 // See getOrCreateMeshProperty.
52 _integration_point_writer.emplace_back(
53 std::make_unique<IntegrationPointWriter>(
54 "sigma_ip",
55 static_cast<int>(mesh.getDimension() == 2 ? 4 : 6) /*n components*/,
56 integration_order, _local_assemblers,
58
59 _integration_point_writer.emplace_back(
60 std::make_unique<IntegrationPointWriter>(
61 "swelling_stress_ip",
62 static_cast<int>(mesh.getDimension() == 2 ? 4 : 6) /*n components*/,
63 integration_order, _local_assemblers,
65
66 _integration_point_writer.emplace_back(
67 std::make_unique<IntegrationPointWriter>(
68 "saturation_ip", 1 /*n components*/, integration_order,
70
71 _integration_point_writer.emplace_back(
72 std::make_unique<IntegrationPointWriter>(
73 "epsilon_ip",
74 static_cast<int>(mesh.getDimension() == 2 ? 4 : 6) /*n components*/,
75 integration_order, _local_assemblers,
77}
78
79template <int DisplacementDim>
81{
82 return false;
83}
84
85template <int DisplacementDim>
88 const int process_id) const
89{
90 // For the monolithic scheme or the M process (deformation) in the staggered
91 // scheme.
92 if (_use_monolithic_scheme || process_id == deformation_process_id)
93 {
94 auto const& l = *_local_to_global_index_map;
95 return {l.dofSizeWithoutGhosts(), l.dofSizeWithoutGhosts(),
96 &l.getGhostIndices(), &this->_sparsity_pattern};
97 }
98
99 // For staggered scheme and T or H process (pressure).
100 auto const& l = *_local_to_global_index_map_with_base_nodes;
101 return {l.dofSizeWithoutGhosts(), l.dofSizeWithoutGhosts(),
102 &l.getGhostIndices(), &_sparsity_pattern_with_linear_element};
103}
104
105template <int DisplacementDim>
107{
108 // Create single component dof in every of the mesh's nodes.
109 _mesh_subset_all_nodes = std::make_unique<MeshLib::MeshSubset>(
110 _mesh, _mesh.getNodes(), _process_data.use_TaylorHood_elements);
111 // Create single component dof in the mesh's base nodes.
112 _base_nodes = MeshLib::getBaseNodes(_mesh.getElements());
113 _mesh_subset_base_nodes = std::make_unique<MeshLib::MeshSubset>(
114 _mesh, _base_nodes, _process_data.use_TaylorHood_elements);
115
116 // TODO move the two data members somewhere else.
117 // for extrapolation of secondary variables of stress or strain
118 std::vector<MeshLib::MeshSubset> all_mesh_subsets_single_component{
119 *_mesh_subset_all_nodes};
120 _local_to_global_index_map_single_component =
121 std::make_unique<NumLib::LocalToGlobalIndexMap>(
122 std::move(all_mesh_subsets_single_component),
123 // by location order is needed for output
125
126 if (_use_monolithic_scheme)
127 {
128 // For gas pressure, which is the first
129 std::vector<MeshLib::MeshSubset> all_mesh_subsets{
130 *_mesh_subset_base_nodes};
131
132 // For capillary pressure, which is the second
133 all_mesh_subsets.push_back(*_mesh_subset_base_nodes);
134
135 // For temperature, which is the third
136 all_mesh_subsets.push_back(*_mesh_subset_base_nodes);
137
138 // For displacement.
139 std::generate_n(
140 std::back_inserter(all_mesh_subsets),
141 getProcessVariables(monolithic_process_id)[deformation_process_id]
142 .get()
143 .getNumberOfGlobalComponents(),
144 [&]() { return *_mesh_subset_all_nodes; });
145
146 std::vector<int> const vec_n_components{
147 n_gas_pressure_components, n_capillary_pressure_components,
148 n_temperature_components, n_displacement_components};
149
150 _local_to_global_index_map =
151 std::make_unique<NumLib::LocalToGlobalIndexMap>(
152 std::move(all_mesh_subsets), vec_n_components,
154 assert(_local_to_global_index_map);
155 }
156 else
157 {
158 OGS_FATAL("A Staggered version of TH2M is not implemented.");
159 }
160}
161
162template <int DisplacementDim>
164 NumLib::LocalToGlobalIndexMap const& dof_table,
165 MeshLib::Mesh const& mesh,
166 unsigned const integration_order)
167{
168 ProcessLib::createLocalAssemblersHM<DisplacementDim, TH2MLocalAssembler>(
169 mesh.getElements(), dof_table, _local_assemblers,
170 mesh.isAxiallySymmetric(), integration_order, _process_data);
171
172 auto add_secondary_variable = [&](std::string const& name,
173 int const num_components,
174 auto get_ip_values_function)
175 {
176 _secondary_variables.addSecondaryVariable(
177 name,
178 makeExtrapolator(num_components, getExtrapolator(),
179 _local_assemblers,
180 std::move(get_ip_values_function)));
181 };
182
183 add_secondary_variable("sigma",
185 DisplacementDim>::RowsAtCompileTime,
187 add_secondary_variable("swelling_stress",
189 DisplacementDim>::RowsAtCompileTime,
191 add_secondary_variable("epsilon",
193 DisplacementDim>::RowsAtCompileTime,
195 add_secondary_variable("velocity_gas", mesh.getDimension(),
197 add_secondary_variable(
198 "velocity_liquid", mesh.getDimension(),
200 add_secondary_variable("saturation", 1,
202 add_secondary_variable("vapour_pressure", 1,
204 add_secondary_variable("porosity", 1,
206 add_secondary_variable("gas_density", 1,
208 add_secondary_variable("solid_density", 1,
210 add_secondary_variable("liquid_density", 1,
212 add_secondary_variable("mole_fraction_gas", 1,
214 add_secondary_variable("mass_fraction_gas", 1,
216 add_secondary_variable(
217 "mass_fraction_liquid", 1,
219
220 add_secondary_variable(
221 "relative_permeability_gas", 1,
223 add_secondary_variable(
224 "relative_permeability_liquid", 1,
226
227 add_secondary_variable("enthalpy_gas", 1,
229
230 add_secondary_variable("enthalpy_liquid", 1,
232
233 add_secondary_variable("enthalpy_solid", 1,
235
236 _process_data.element_saturation = MeshLib::getOrCreateMeshProperty<double>(
237 const_cast<MeshLib::Mesh&>(mesh), "saturation_avg",
239
240 _process_data.gas_pressure_interpolated =
241 MeshLib::getOrCreateMeshProperty<double>(
242 const_cast<MeshLib::Mesh&>(mesh), "gas_pressure_interpolated",
244
245 _process_data.capillary_pressure_interpolated =
246 MeshLib::getOrCreateMeshProperty<double>(
247 const_cast<MeshLib::Mesh&>(mesh), "capillary_pressure_interpolated",
249
250 _process_data.liquid_pressure_interpolated =
251 MeshLib::getOrCreateMeshProperty<double>(
252 const_cast<MeshLib::Mesh&>(mesh), "liquid_pressure_interpolated",
254
255 _process_data.temperature_interpolated =
256 MeshLib::getOrCreateMeshProperty<double>(
257 const_cast<MeshLib::Mesh&>(mesh), "temperature_interpolated",
259
260 // Set initial conditions for integration point data.
261 for (auto const& ip_writer : _integration_point_writer)
262 {
263 // Find the mesh property with integration point writer's name.
264 auto const& name = ip_writer->name();
265 if (!mesh.getProperties().existsPropertyVector<double>(name))
266 {
267 continue;
268 }
269 auto const& _meshproperty =
270 *mesh.getProperties().template getPropertyVector<double>(name);
271
272 // The mesh property must be defined on integration points.
273 if (_meshproperty.getMeshItemType() !=
275 {
276 continue;
277 }
278
279 auto const ip_meta_data =
281
282 // Check the number of components.
283 if (ip_meta_data.n_components !=
284 _meshproperty.getNumberOfGlobalComponents())
285 {
286 OGS_FATAL(
287 "Different number of components in meta data ({:d}) than in "
288 "the integration point field data for '{:s}': {:d}.",
289 ip_meta_data.n_components, name,
290 _meshproperty.getNumberOfGlobalComponents());
291 }
292
293 // Now we have a properly named vtk's field data array and the
294 // corresponding meta data.
295 std::size_t position = 0;
296 for (auto& local_asm : _local_assemblers)
297 {
298 std::size_t const integration_points_read =
299 local_asm->setIPDataInitialConditions(
300 name, &_meshproperty[position],
301 ip_meta_data.integration_order);
302 if (integration_points_read == 0)
303 {
304 OGS_FATAL(
305 "No integration points read in the integration point "
306 "initial conditions set function.");
307 }
308 position += integration_points_read * ip_meta_data.n_components;
309 }
310 }
311
312 // Initialize local assemblers after all variables have been set.
314 &LocalAssemblerInterface::initialize, _local_assemblers,
315 *_local_to_global_index_map);
316}
317
318template <int DisplacementDim>
320{
321 if (_use_monolithic_scheme)
322 {
323 initializeProcessBoundaryConditionsAndSourceTerms(
324 *_local_to_global_index_map, monolithic_process_id);
325 return;
326 }
327
328 // Staggered scheme:
329 OGS_FATAL("A Staggered version of TH2M is not implemented.");
330}
331
332template <int DisplacementDim>
334 std::vector<GlobalVector*>& x, double const t, int const process_id)
335{
336 if (process_id != 0)
337 {
338 return;
339 }
340
341 DBUG("Set initial conditions of TH2MProcess.");
342
345 getDOFTable(process_id), *x[process_id], t, _use_monolithic_scheme,
346 process_id);
347}
348
349template <int DisplacementDim>
351 const double t, double const dt, std::vector<GlobalVector*> const& x,
352 std::vector<GlobalVector*> const& xdot, int const process_id,
354{
355 DBUG("Assemble the equations for TH2M");
356
357 std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
358 dof_table = {std::ref(*_local_to_global_index_map)};
359
360 ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
361
363 _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers,
364 pv.getActiveElementIDs(), dof_table, t, dt, x, xdot, process_id, M, K,
365 b);
366}
367
368template <int DisplacementDim>
370 const double t, double const dt, std::vector<GlobalVector*> const& x,
371 std::vector<GlobalVector*> const& xdot, int const process_id,
373{
374 std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
375 dof_tables;
376 // For the monolithic scheme
377 if (_use_monolithic_scheme)
378 {
379 DBUG("Assemble the Jacobian of TH2M for the monolithic scheme.");
380 dof_tables.emplace_back(*_local_to_global_index_map);
381 }
382 else
383 {
384 // For the staggered scheme
385 OGS_FATAL("A Staggered version of TH2M is not implemented.");
386 }
387
388 ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
389
392 _local_assemblers, pv.getActiveElementIDs(), dof_tables, t, dt, x, xdot,
393 process_id, M, K, b, Jac);
394
395 auto copyRhs = [&](int const variable_id, auto& output_vector)
396 {
397 if (_use_monolithic_scheme)
398 {
399 transformVariableFromGlobalVector(b, variable_id, dof_tables[0],
400 output_vector,
401 std::negate<double>());
402 }
403 else
404 {
405 transformVariableFromGlobalVector(b, 0, dof_tables[process_id],
406 output_vector,
407 std::negate<double>());
408 }
409 };
410 if (_use_monolithic_scheme || process_id == 1)
411 {
412 copyRhs(0, *_hydraulic_flow);
413 }
414 if (_use_monolithic_scheme || process_id == 2)
415 {
416 copyRhs(1, *_nodal_forces);
417 }
418}
419
420template <int DisplacementDim>
422 std::vector<GlobalVector*> const& x, double const t, double const dt,
423 const int process_id)
424{
425 DBUG("PreTimestep TH2MProcess.");
426
427 if (hasMechanicalProcess(process_id))
428 {
430 getProcessVariables(process_id)[0];
431
433 &LocalAssemblerInterface::preTimestep, _local_assemblers,
434 pv.getActiveElementIDs(), *_local_to_global_index_map,
435 *x[process_id], t, dt);
436 }
437}
438
439template <int DisplacementDim>
441 std::vector<GlobalVector*> const& x, double const t, double const dt,
442 const int process_id)
443{
444 DBUG("PostTimestep TH2MProcess.");
445 auto const dof_tables = getDOFTables(x.size());
446 ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
448 &LocalAssemblerInterface::postTimestep, _local_assemblers,
449 pv.getActiveElementIDs(), dof_tables, x, t, dt);
450}
451
452template <int DisplacementDim>
454 double const t, double const dt, std::vector<GlobalVector*> const& x,
455 GlobalVector const& x_dot, const int process_id)
456{
457 if (process_id != 0)
458 {
459 return;
460 }
461
462 DBUG("Compute the secondary variables for TH2MProcess.");
463 auto const dof_tables = getDOFTables(x.size());
464
465 ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
468 pv.getActiveElementIDs(), dof_tables, t, dt, x, x_dot, process_id);
469}
470
471template <int DisplacementDim>
472std::tuple<NumLib::LocalToGlobalIndexMap*, bool>
474{
475 const bool manage_storage = false;
476 return std::make_tuple(_local_to_global_index_map_single_component.get(),
477 manage_storage);
478}
479
480template <int DisplacementDim>
482 const int process_id) const
483{
484 if (hasMechanicalProcess(process_id))
485 {
486 return *_local_to_global_index_map;
487 }
488
489 // For the equation of pressure
490 return *_local_to_global_index_map_with_base_nodes;
491}
492
493template <int DisplacementDim>
494std::vector<NumLib::LocalToGlobalIndexMap const*>
495TH2MProcess<DisplacementDim>::getDOFTables(int const number_of_processes) const
496{
497 std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
498 dof_tables.reserve(number_of_processes);
499 std::generate_n(std::back_inserter(dof_tables), number_of_processes,
500 [&]() { return &getDOFTable(dof_tables.size()); });
501 return dof_tables;
502}
503template class TH2MProcess<2>;
504template class TH2MProcess<3>;
505
506} // namespace TH2M
507} // namespace ProcessLib
#define OGS_FATAL(...)
Definition: Error.h:26
void DBUG(char const *fmt, Args const &... args)
Definition: Logging.h:29
Global vector based on Eigen vector.
Definition: EigenVector.h:28
bool isAxiallySymmetric() const
Definition: Mesh.h:131
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition: Mesh.h:103
unsigned getDimension() const
Returns the dimension of the mesh (determined by the maximum dimension over all elements).
Definition: Mesh.h:76
Properties & getProperties()
Definition: Mesh.h:128
bool existsPropertyVector(std::string const &name) const
virtual void postTimestep(std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, double const t, double const dt)
virtual void preTimestep(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, double const t, double const delta_t)
virtual void initialize(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table)
void setInitialConditions(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, double const t, bool const use_monolithic_scheme, int const process_id)
virtual void computeSecondaryVariable(std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_dot, int const process_id)
std::vector< std::size_t > const & getActiveElementIDs() const
std::vector< std::unique_ptr< IntegrationPointWriter > > _integration_point_writer
Definition: Process.h:352
Handles configuration of several secondary variables from the project file.
void assembleConcreteProcess(const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &xdot, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b) override
bool isLinear() const override
Definition: TH2MProcess.cpp:80
std::vector< std::unique_ptr< LocalAssemblerInterface > > _local_assemblers
Definition: TH2MProcess.h:109
std::tuple< NumLib::LocalToGlobalIndexMap *, bool > getDOFTableForExtrapolatorData() const override
void setInitialConditionsConcreteProcess(std::vector< GlobalVector * > &x, double const t, int const process_id) override
MeshLib::PropertyVector< double > * _nodal_forces
Definition: TH2MProcess.h:141
TH2MProcess(std::string name, MeshLib::Mesh &mesh, std::unique_ptr< ProcessLib::AbstractJacobianAssembler > &&jacobian_assembler, std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const &parameters, unsigned const integration_order, std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > &&process_variables, TH2MProcessData< DisplacementDim > &&process_data, SecondaryVariableCollection &&secondary_variables, bool const use_monolithic_scheme)
Definition: TH2MProcess.cpp:27
void assembleWithJacobianConcreteProcess(const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &xdot, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac) override
void preTimestepConcreteProcess(std::vector< GlobalVector * > const &x, double const t, double const dt, const int process_id) override
void postTimestepConcreteProcess(std::vector< GlobalVector * > const &x, const double t, const double dt, int const) override
MathLib::MatrixSpecifications getMatrixSpecifications(const int process_id) const override
Definition: TH2MProcess.cpp:87
MeshLib::PropertyVector< double > * _hydraulic_flow
Definition: TH2MProcess.h:142
void constructDofTable() override
void computeSecondaryVariableConcrete(double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_dot, const int process_id) override
NumLib::LocalToGlobalIndexMap const & getDOFTable(const int process_id) const override
void initializeConcreteProcess(NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order) override
Process specific initialization called by initialize().
std::vector< NumLib::LocalToGlobalIndexMap const * > getDOFTables(int const number_of_processes) const
void initializeBoundaryConditions() override
void assemble(std::size_t const mesh_item_id, LocalAssemblerInterface &local_assembler, std::vector< std::reference_wrapper< NumLib::LocalToGlobalIndexMap > > const &dof_tables, double const t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &xdot, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b)
void assembleWithJacobian(std::size_t const mesh_item_id, LocalAssemblerInterface &local_assembler, std::vector< std::reference_wrapper< NumLib::LocalToGlobalIndexMap > > const &dof_tables, const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &xdot, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac)
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
Definition: KelvinVector.h:48
static const double t
std::vector< Node * > getBaseNodes(std::vector< Element * > const &elements)
Definition: Utils.h:26
@ BY_LOCATION
Ordering data by spatial location.
void transformVariableFromGlobalVector(GlobalVector const &input_vector, int const variable_id, NumLib::LocalToGlobalIndexMap const &local_to_global_index_map, MeshLib::PropertyVector< double > &output_vector, Functor mapFunction)
Definition: DOFTableUtil.h:59
IntegrationPointMetaData getIntegrationPointMetaData(MeshLib::Properties const &properties, std::string const &name)
SecondaryVariableFunctions makeExtrapolator(const unsigned num_components, NumLib::Extrapolator &extrapolator, LocalAssemblerCollection const &local_assemblers, typename NumLib::ExtrapolatableLocalAssemblerCollection< LocalAssemblerCollection >::IntegrationPointValuesMethod integration_point_values_method)
static void executeSelectedMemberOnDereferenced(Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)
static void executeSelectedMemberDereferenced(Object &object, Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)
static void executeMemberOnDereferenced(Method method, Container const &container, Args &&... args)
virtual std::vector< double > const & getIntPtLiquidDensity(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > getEpsilon() const =0
virtual std::vector< double > const & getIntPtRelativePermeabilityLiquid(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > const & getIntPtEpsilon(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > const & getIntPtEnthalpyGas(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > const & getIntPtVapourPressure(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > const & getIntPtRelativePermeabilityGas(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > const & getIntPtSolidDensity(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > const & getIntPtPorosity(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > const & getIntPtMoleFractionGas(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > const & getIntPtMassFractionGas(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > const & getIntPtSwellingStress(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > const & getIntPtEnthalpyLiquid(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > getSaturation() const =0
virtual std::vector< double > const & getIntPtSaturation(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > const & getIntPtEnthalpySolid(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > const & getIntPtDarcyVelocityGas(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > const & getIntPtSigma(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > const & getIntPtDarcyVelocityLiquid(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > getSwellingStress() const =0
virtual std::vector< double > const & getIntPtGasDensity(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > const & getIntPtMassFractionLiquid(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > getSigma() const =0