Skip to content
Snippets Groups Projects
Commit f1ccca4b authored by wenqing's avatar wenqing
Browse files

[Coupling] Removed the coupling calculation from the classes of single process

parent fe86513b
No related branches found
No related tags found
No related merge requests found
Showing
with 26 additions and 563 deletions
Property of solid phase of liquid flow process.
Biot's constant for the liquid thermal expansion computation.
/**
* \copyright
* Copyright (c) 2012-2017, OpenGeoSys Community (http://www.opengeosys.org)
* Distributed under a Modified BSD License.
* See accompanying file LICENSE.txt or
* http://www.opengeosys.org/project/license
*
* \file HeatConductionFEM-impl.h
* Created on January 17, 2017, 3:41 PM
*/
#pragma once
#include "HeatConductionFEM.h"
#include "NumLib/Function/Interpolation.h"
namespace ProcessLib
{
namespace HeatConduction
{
template <typename ShapeFunction, typename IntegrationMethod,
unsigned GlobalDim>
template <typename LiquidFlowVelocityCalculator>
void LocalAssemblerData<ShapeFunction, IntegrationMethod, GlobalDim>::
assembleHeatTransportLiquidFlow(
double const t, int const material_id, SpatialPosition& pos,
int const gravitational_axis_id,
double const gravitational_acceleration,
Eigen::MatrixXd const& permeability,
ProcessLib::LiquidFlow::LiquidFlowMaterialProperties const&
liquid_flow_properties,
std::vector<double> const& local_x, std::vector<double> const& local_p,
std::vector<double>& local_M_data, std::vector<double>& local_K_data)
{
auto const local_matrix_size = local_x.size();
// This assertion is valid only if all nodal d.o.f. use the same shape
// matrices.
assert(local_matrix_size == ShapeFunction::NPOINTS * NUM_NODAL_DOF);
auto local_M = MathLib::createZeroedMatrix<NodalMatrixType>(
local_M_data, local_matrix_size, local_matrix_size);
auto local_K = MathLib::createZeroedMatrix<NodalMatrixType>(
local_K_data, local_matrix_size, local_matrix_size);
const auto local_p_vec =
MathLib::toVector<NodalVectorType>(local_p, local_matrix_size);
unsigned const n_integration_points =
_integration_method.getNumberOfPoints();
const double porosity_variable = 0.;
for (unsigned ip = 0; ip < n_integration_points; ip++)
{
pos.setIntegrationPoint(ip);
auto const& sm = _shape_matrices[ip];
auto const& wp = _integration_method.getWeightedPoint(ip);
double p = 0.;
NumLib::shapeFunctionInterpolate(local_p, sm.N, p);
double T = 0.;
NumLib::shapeFunctionInterpolate(local_x, sm.N, T);
// Thermal conductivity of solid phase.
auto const k_s = _process_data.thermal_conductivity(t, pos)[0];
// Specific heat capacity of solid phase.
auto const cp_s = _process_data.heat_capacity(t, pos)[0];
// Density of solid phase.
auto const rho_s = _process_data.density(t, pos)[0];
// Thermal conductivity of liquid.
double const k_f = liquid_flow_properties.getThermalConductivity(p, T);
// Specific heat capacity of liquid.
double const cp_f = liquid_flow_properties.getHeatCapacity(p, T);
// Density of liquid.
double const rho_f = liquid_flow_properties.getLiquidDensity(p, T);
// Porosity of porous media.
double const n = liquid_flow_properties.getPorosity(
material_id, t, pos, porosity_variable, T);
// Effective specific heat capacity.
double const effective_cp = (1.0 - n) * cp_s * rho_s + n * cp_f * rho_f;
// Effective thermal conductivity.
double const effective_K = (1.0 - n) * k_s + n * k_f;
double const integration_factor =
sm.detJ * wp.getWeight() * sm.integralMeasure;
local_M.noalias() +=
effective_cp * sm.N.transpose() * sm.N * integration_factor;
local_K.noalias() +=
sm.dNdx.transpose() * effective_K * sm.dNdx * integration_factor;
// Compute fluid flow velocity:
double const mu =
liquid_flow_properties.getViscosity(p, T); // Viscosity
GlobalDimVectorType const& darcy_velocity =
LiquidFlowVelocityCalculator::calculateVelocity(
local_p_vec, sm, permeability, mu,
rho_f * gravitational_acceleration, gravitational_axis_id);
// Add the advection term
local_K.noalias() += cp_f * rho_f * sm.N.transpose() *
(darcy_velocity.transpose() * sm.dNdx) *
integration_factor;
}
}
template <typename ShapeFunction, typename IntegrationMethod,
unsigned GlobalDim>
void LocalAssemblerData<ShapeFunction, IntegrationMethod, GlobalDim>::
assembleWithCoupledTerm(double const t, std::vector<double> const& local_x,
std::vector<double>& local_M_data,
std::vector<double>& local_K_data,
std::vector<double>& /*local_b_data*/,
LocalCoupledSolutions const& coupled_term)
{
/*
for (auto const& coupled_process_pair : coupled_term.coupled_processes)
{
if (coupled_process_pair.first ==
std::type_index(typeid(ProcessLib::LiquidFlow::LiquidFlowProcess)))
{
assert(
dynamic_cast<const ProcessLib::LiquidFlow::LiquidFlowProcess*>(
&(coupled_process_pair.second)) != nullptr);
auto const& pcs =
static_cast<ProcessLib::LiquidFlow::LiquidFlowProcess const&>(
coupled_process_pair.second);
const auto liquid_flow_prop = pcs.getLiquidFlowMaterialProperties();
const auto local_p =
coupled_term.local_coupled_xs.at(coupled_process_pair.first);
SpatialPosition pos;
pos.setElementID(_element.getID());
const int material_id = liquid_flow_prop->getMaterialID(pos);
const Eigen::MatrixXd& K = liquid_flow_prop->getPermeability(
material_id, t, pos, _element.getDimension());
if (K.size() == 1)
{
assembleHeatTransportLiquidFlow<
IsotropicLiquidFlowVelocityCalculator>(
t, material_id, pos, pcs.getGravitationalAxisID(),
pcs.getGravitationalAcceleration(), K, *liquid_flow_prop,
local_x, local_p, local_M_data, local_K_data);
}
else
{
assembleHeatTransportLiquidFlow<
AnisotropicLiquidFlowVelocityCalculator>(
t, material_id, pos, pcs.getGravitationalAxisID(),
pcs.getGravitationalAcceleration(), K, *liquid_flow_prop,
local_x, local_p, local_M_data, local_K_data);
}
}
else
{
OGS_FATAL(
"This coupled process is not presented for "
"HeatConduction process");
}
}*/
}
} // namespace HeatConduction
} // namespace ProcessLib
...@@ -18,6 +18,7 @@ ...@@ -18,6 +18,7 @@
#include "MeshLib/PropertyVector.h" #include "MeshLib/PropertyVector.h"
#include "MaterialLib/Fluid/FluidProperty.h" #include "MaterialLib/Fluid/FluidProperty.h"
#include "MaterialLib/PorousMedium/Permeability/Permeability.h"
#include "MaterialLib/PorousMedium/Porosity/Porosity.h" #include "MaterialLib/PorousMedium/Porosity/Porosity.h"
#include "MaterialLib/PorousMedium/Storage/Storage.h" #include "MaterialLib/PorousMedium/Storage/Storage.h"
#include "MaterialLib/Fluid/FluidProperties/CreateFluidProperties.h" #include "MaterialLib/Fluid/FluidProperties/CreateFluidProperties.h"
...@@ -26,7 +27,6 @@ ...@@ -26,7 +27,6 @@
#include "MaterialLib/PorousMedium/PorousPropertyHeaders.h" #include "MaterialLib/PorousMedium/PorousPropertyHeaders.h"
#include "ProcessLib/Utils/ProcessUtils.h" #include "ProcessLib/Utils/ProcessUtils.h"
#include "ProcessLib/Parameter/ConstantParameter.h"
#include "LiquidFlowMaterialProperties.h" #include "LiquidFlowMaterialProperties.h"
...@@ -97,32 +97,10 @@ createLiquidFlowMaterialProperties( ...@@ -97,32 +97,10 @@ createLiquidFlowMaterialProperties(
BaseLib::reorderVector(porosity_models, mat_ids); BaseLib::reorderVector(porosity_models, mat_ids);
BaseLib::reorderVector(storage_models, mat_ids); BaseLib::reorderVector(storage_models, mat_ids);
//! \ogs_file_param{prj__processes__process__LIQUID_FLOW__material_property__solid}
auto const solid_config = config.getConfigSubtreeOptional("solid");
if (solid_config)
{
auto& solid_thermal_expansion = findParameter<double>(
//! \ogs_file_param_special{prj__processes__process__LIQUID_FLOW__material_property__solid__thermal_expansion}
*solid_config, "thermal_expansion", parameters, 1);
DBUG("Use \'%s\' as solid thermal expansion.",
solid_thermal_expansion.name.c_str());
auto& biot_constant = findParameter<double>(
//! \ogs_file_param_special{prj__processes__process__LIQUID_FLOW__material_property__solid__biot_constant}
*solid_config, "biot_constant", parameters, 1);
return std::make_unique<LiquidFlowMaterialProperties>(
std::move(fluid_properties),
std::move(intrinsic_permeability_models),
std::move(porosity_models), std::move(storage_models),
has_material_ids, material_ids, solid_thermal_expansion,
biot_constant);
}
ConstantParameter<double> void_parameter("void_solid_thermal_expansion",
0.);
return std::make_unique<LiquidFlowMaterialProperties>( return std::make_unique<LiquidFlowMaterialProperties>(
std::move(fluid_properties), std::move(intrinsic_permeability_models), std::move(fluid_properties), std::move(intrinsic_permeability_models),
std::move(porosity_models), std::move(storage_models), has_material_ids, std::move(porosity_models), std::move(storage_models), has_material_ids,
material_ids, void_parameter, void_parameter); material_ids);
} }
} // end of namespace } // end of namespace
......
...@@ -16,10 +16,6 @@ ...@@ -16,10 +16,6 @@
#include "NumLib/Function/Interpolation.h" #include "NumLib/Function/Interpolation.h"
#include "ProcessLib/CoupledSolutionsForStaggeredScheme.h"
#include "ProcessLib/HeatConduction/HeatConductionProcess.h"
namespace ProcessLib namespace ProcessLib
{ {
namespace LiquidFlow namespace LiquidFlow
...@@ -53,60 +49,6 @@ void LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>:: ...@@ -53,60 +49,6 @@ void LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>::
pos, permeability); pos, permeability);
} }
/*
template <typename ShapeFunction, typename IntegrationMethod,
unsigned GlobalDim>
void LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>::
assembleWithCoupledTerm(double const t, std::vector<double> const& local_x,
std::vector<double>& local_M_data,
std::vector<double>& local_K_data,
std::vector<double>& local_b_data,
LocalCouplingTerm const& coupled_term)
{
SpatialPosition pos;
pos.setElementID(_element.getID());
const int material_id = _material_properties.getMaterialID(pos);
const Eigen::MatrixXd& permeability = _material_properties.getPermeability(
material_id, t, pos, _element.getDimension());
// Note: For Inclined 1D in 2D/3D or 2D element in 3D, the first item in
// the assert must be changed to permeability.rows() ==
// _element->getDimension()
assert(permeability.rows() == GlobalDim || permeability.rows() == 1);
const double dt = coupled_term.dt;
for (auto const& coupled_process_pair : coupled_term.coupled_processes)
{
if (coupled_process_pair.first ==
std::type_index(
typeid(ProcessLib::HeatConduction::HeatConductionProcess)))
{
const auto local_T0 =
coupled_term.local_coupled_xs0.at(coupled_process_pair.first);
const auto local_T1 =
coupled_term.local_coupled_xs.at(coupled_process_pair.first);
if (permeability.size() == 1) // isotropic or 1D problem.
assembleWithCoupledWithHeatTransport<IsotropicCalculator>(
material_id, t, dt, local_x, local_T0, local_T1,
local_M_data, local_K_data, local_b_data, pos,
permeability);
else
assembleWithCoupledWithHeatTransport<AnisotropicCalculator>(
material_id, t, dt, local_x, local_T0, local_T1,
local_M_data, local_K_data, local_b_data, pos,
permeability);
}
else
{
OGS_FATAL(
"This coupled process is not presented for "
"LiquidFlow process");
}
}
}
*/
template <typename ShapeFunction, typename IntegrationMethod, template <typename ShapeFunction, typename IntegrationMethod,
unsigned GlobalDim> unsigned GlobalDim>
template <typename LaplacianGravityVelocityCalculator> template <typename LaplacianGravityVelocityCalculator>
...@@ -168,83 +110,6 @@ void LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>:: ...@@ -168,83 +110,6 @@ void LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>::
} }
} }
/*
template <typename ShapeFunction, typename IntegrationMethod,
unsigned GlobalDim>
template <typename LaplacianGravityVelocityCalculator>
void LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>::
assembleWithCoupledWithHeatTransport(
const int material_id, double const t, double const dt,
std::vector<double> const& local_x, std::vector<double> const& local_T0,
std::vector<double> const& local_T1, std::vector<double>& local_M_data,
std::vector<double>& local_K_data, std::vector<double>& local_b_data,
SpatialPosition const& pos, Eigen::MatrixXd const& permeability)
{
auto const local_matrix_size = local_x.size();
assert(local_matrix_size == ShapeFunction::NPOINTS);
auto local_M = MathLib::createZeroedMatrix<NodalMatrixType>(
local_M_data, local_matrix_size, local_matrix_size);
auto local_K = MathLib::createZeroedMatrix<NodalMatrixType>(
local_K_data, local_matrix_size, local_matrix_size);
auto local_b = MathLib::createZeroedVector<NodalVectorType>(
local_b_data, local_matrix_size);
unsigned const n_integration_points =
_integration_method.getNumberOfPoints();
// TODO: The following two variables should be calculated inside the
// the integration loop for non-constant porosity and storage models.
double porosity_variable = 0.;
double storage_variable = 0.;
for (unsigned ip = 0; ip < n_integration_points; ip++)
{
auto const& sm = _shape_matrices[ip];
auto const& wp = _integration_method.getWeightedPoint(ip);
double p = 0.;
NumLib::shapeFunctionInterpolate(local_x, sm.N, p);
double T0 = 0.;
NumLib::shapeFunctionInterpolate(local_T0, sm.N, T0);
double T = 0.;
NumLib::shapeFunctionInterpolate(local_T1, sm.N, T);
const double integration_factor =
sm.integralMeasure * sm.detJ * wp.getWeight();
// Assemble mass matrix, M
local_M.noalias() += _material_properties.getMassCoefficient(
material_id, t, pos, porosity_variable,
storage_variable, p, T) *
sm.N.transpose() * sm.N * integration_factor;
// Compute density:
const double rho = _material_properties.getLiquidDensity(p, T);
const double rho_g = rho * _gravitational_acceleration;
// Compute viscosity:
const double mu = _material_properties.getViscosity(p, T);
// Assemble Laplacian, K, and RHS by the gravitational term
LaplacianGravityVelocityCalculator::calculateLaplacianAndGravityTerm(
local_K, local_b, sm, permeability, integration_factor, mu, rho_g,
_gravitational_axis_id);
// Add the thermal expansion term
auto const solid_thermal_expansion =
_material_properties.getSolidThermalExpansion(t, pos);
auto const biot_constant = _material_properties.getBiotConstant(t, pos);
auto const porosity = _material_properties.getPorosity(
material_id, t, pos, porosity_variable, T);
const double eff_thermal_expansion =
3.0 * (biot_constant - porosity) * solid_thermal_expansion -
porosity * _material_properties.getdLiquidDensity_dT(p, T) / rho;
local_b.noalias() +=
eff_thermal_expansion * (T - T0) * integration_factor * sm.N / dt;
}
}
*/
template <typename ShapeFunction, typename IntegrationMethod, template <typename ShapeFunction, typename IntegrationMethod,
unsigned GlobalDim> unsigned GlobalDim>
std::vector<double> const& std::vector<double> const&
...@@ -272,41 +137,13 @@ LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>:: ...@@ -272,41 +137,13 @@ LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>::
// the assert must be changed to perm.rows() == _element->getDimension() // the assert must be changed to perm.rows() == _element->getDimension()
assert(permeability.rows() == GlobalDim || permeability.rows() == 1); assert(permeability.rows() == GlobalDim || permeability.rows() == 1);
//if (!_coupling_term)
{
computeDarcyVelocity(permeability, local_x, veloctiy_cache_vectors);
}
/*
else
{
auto const local_coupled_xs =
getCurrentLocalSolutionsOfCoupledProcesses(
_coupling_term->coupled_xs, indices);
if (local_coupled_xs.empty())
computeDarcyVelocity(permeability, local_x, veloctiy_cache_vectors);
else
computeDarcyVelocityWithCoupling(permeability, local_x,
local_coupled_xs,
veloctiy_cache_vectors);
}*/
return veloctiy_cache;
}
template <typename ShapeFunction, typename IntegrationMethod,
unsigned GlobalDim>
void LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>::
computeDarcyVelocity(
Eigen::MatrixXd const& permeability,
std::vector<double> const& local_x,
MatrixOfVelocityAtIntegrationPoints& darcy_velocity_at_ips) const
{
if (permeability.size() == 1) // isotropic or 1D problem. if (permeability.size() == 1) // isotropic or 1D problem.
computeDarcyVelocityLocal<IsotropicCalculator>(local_x, permeability, computeDarcyVelocityLocal<IsotropicCalculator>(local_x, permeability,
darcy_velocity_at_ips); veloctiy_cache_vectors);
else else
computeDarcyVelocityLocal<AnisotropicCalculator>(local_x, permeability, computeDarcyVelocityLocal<AnisotropicCalculator>(
darcy_velocity_at_ips); local_x, permeability, veloctiy_cache_vectors);
return veloctiy_cache;
} }
template <typename ShapeFunction, typename IntegrationMethod, template <typename ShapeFunction, typename IntegrationMethod,
...@@ -346,67 +183,6 @@ void LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>:: ...@@ -346,67 +183,6 @@ void LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>::
} }
} }
/*
template <typename ShapeFunction, typename IntegrationMethod,
unsigned GlobalDim>
void LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>::
computeDarcyVelocityWithCoupling(
Eigen::MatrixXd const& permeability, std::vector<double> const& local_x,
std::unordered_map<std::type_index, const std::vector<double>> const&
coupled_local_solutions,
MatrixOfVelocityAtIntegrationPoints& darcy_velocity_at_ips) const
{
const auto local_T = coupled_local_solutions.at(std::type_index(
typeid(ProcessLib::HeatConduction::HeatConductionProcess)));
if (permeability.size() == 1) // isotropic or 1D problem.
computeDarcyVelocityCoupledWithHeatTransportLocal<IsotropicCalculator>(
local_x, local_T, permeability, darcy_velocity_at_ips);
else
computeDarcyVelocityCoupledWithHeatTransportLocal<
AnisotropicCalculator>(local_x, local_T, permeability,
darcy_velocity_at_ips);
}
template <typename ShapeFunction, typename IntegrationMethod,
unsigned GlobalDim>
template <typename LaplacianGravityVelocityCalculator>
void LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>::
computeDarcyVelocityCoupledWithHeatTransportLocal(
std::vector<double> const& local_x,
std::vector<double> const& local_T,
Eigen::MatrixXd const& permeability,
MatrixOfVelocityAtIntegrationPoints& darcy_velocity_at_ips) const
{
auto const local_matrix_size = local_x.size();
assert(local_matrix_size == ShapeFunction::NPOINTS);
const auto local_p_vec =
MathLib::toVector<NodalVectorType>(local_x, local_matrix_size);
unsigned const n_integration_points =
_integration_method.getNumberOfPoints();
for (unsigned ip = 0; ip < n_integration_points; ip++)
{
auto const& sm = _shape_matrices[ip];
double p = 0.;
NumLib::shapeFunctionInterpolate(local_x, sm.N, p);
double T = 0.;
NumLib::shapeFunctionInterpolate(local_T, sm.N, T);
const double rho_g = _material_properties.getLiquidDensity(p, T) *
_gravitational_acceleration;
// Compute viscosity:
const double mu = _material_properties.getViscosity(p, T);
LaplacianGravityVelocityCalculator::calculateVelocity(
ip, local_p_vec, sm, permeability, mu, rho_g,
_gravitational_axis_id, darcy_velocity_at_ips);
}
}
*/
template <typename ShapeFunction, typename IntegrationMethod, template <typename ShapeFunction, typename IntegrationMethod,
unsigned GlobalDim> unsigned GlobalDim>
void LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>:: void LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>::
......
...@@ -31,7 +31,6 @@ ...@@ -31,7 +31,6 @@
namespace ProcessLib namespace ProcessLib
{ {
struct CoupledSolutionsForStaggeredScheme;
namespace LiquidFlow namespace LiquidFlow
{ {
...@@ -42,24 +41,11 @@ class LiquidFlowLocalAssemblerInterface ...@@ -42,24 +41,11 @@ class LiquidFlowLocalAssemblerInterface
public NumLib::ExtrapolatableElement public NumLib::ExtrapolatableElement
{ {
public: public:
LiquidFlowLocalAssemblerInterface() = default;
void setCoupledSolutionsForStaggeredScheme(
std::size_t const /*mesh_item_id*/,
CoupledSolutionsForStaggeredScheme* const coupled_solutions)
{
_coupled_solutions = coupled_solutions;
}
virtual std::vector<double> const& getIntPtDarcyVelocity( virtual std::vector<double> const& getIntPtDarcyVelocity(
const double /*t*/, const double /*t*/,
GlobalVector const& /*current_solution*/, GlobalVector const& /*current_solution*/,
NumLib::LocalToGlobalIndexMap const& /*dof_table*/, NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
std::vector<double>& /*cache*/) const = 0; std::vector<double>& /*cache*/) const = 0;
protected:
/// Pointer that is set from a Process class.
CoupledSolutionsForStaggeredScheme* _coupled_solutions;
}; };
template <typename ShapeFunction, typename IntegrationMethod, template <typename ShapeFunction, typename IntegrationMethod,
...@@ -89,8 +75,7 @@ public: ...@@ -89,8 +75,7 @@ public:
double const gravitational_acceleration, double const gravitational_acceleration,
double const reference_temperature, double const reference_temperature,
LiquidFlowMaterialProperties const& material_propertries) LiquidFlowMaterialProperties const& material_propertries)
: LiquidFlowLocalAssemblerInterface(), : _element(element),
_element(element),
_integration_method(integration_order), _integration_method(integration_order),
_shape_matrices(initShapeMatrices<ShapeFunction, ShapeMatricesType, _shape_matrices(initShapeMatrices<ShapeFunction, ShapeMatricesType,
IntegrationMethod, GlobalDim>( IntegrationMethod, GlobalDim>(
...@@ -107,13 +92,6 @@ public: ...@@ -107,13 +92,6 @@ public:
std::vector<double>& local_K_data, std::vector<double>& local_K_data,
std::vector<double>& local_b_data) override; std::vector<double>& local_b_data) override;
/*
void assembleWithCoupledTerm(
double const t, std::vector<double> const& local_x,
std::vector<double>& local_M_data, std::vector<double>& local_K_data,
std::vector<double>& local_b_data,
LocalCouplingTerm const& coupled_term) override;*/
Eigen::Map<const Eigen::RowVectorXd> getShapeMatrix( Eigen::Map<const Eigen::RowVectorXd> getShapeMatrix(
const unsigned integration_point) const override const unsigned integration_point) const override
{ {
...@@ -187,41 +165,12 @@ private: ...@@ -187,41 +165,12 @@ private:
SpatialPosition const& pos, SpatialPosition const& pos,
Eigen::MatrixXd const& permeability); Eigen::MatrixXd const& permeability);
/*
template <typename LaplacianGravityVelocityCalculator>
void assembleWithCoupledWithHeatTransport(
const int material_id, double const t, double const dt,
std::vector<double> const& local_x, std::vector<double> const& local_T0,
std::vector<double> const& local_T1, std::vector<double>& local_M_data,
std::vector<double>& local_K_data, std::vector<double>& local_b_data,
SpatialPosition const& pos, Eigen::MatrixXd const& permeability);*/
void computeDarcyVelocity(
Eigen::MatrixXd const& permeability,
std::vector<double> const& local_x,
MatrixOfVelocityAtIntegrationPoints& darcy_velocity_at_ips) const;
/*
void computeDarcyVelocityWithCoupling(
Eigen::MatrixXd const& permeability, std::vector<double> const& local_x,
std::unordered_map<std::type_index, const std::vector<double>> const&
coupled_local_solutions,
MatrixOfVelocityAtIntegrationPoints& darcy_velocity_at_ips) const;*/
template <typename LaplacianGravityVelocityCalculator> template <typename LaplacianGravityVelocityCalculator>
void computeDarcyVelocityLocal( void computeDarcyVelocityLocal(
std::vector<double> const& local_x, std::vector<double> const& local_x,
Eigen::MatrixXd const& permeability, Eigen::MatrixXd const& permeability,
MatrixOfVelocityAtIntegrationPoints& darcy_velocity_at_ips) const; MatrixOfVelocityAtIntegrationPoints& darcy_velocity_at_ips) const;
/*
template <typename LaplacianGravityVelocityCalculator>
void computeDarcyVelocityCoupledWithHeatTransportLocal(
std::vector<double> const& local_x,
std::vector<double> const& local_T,
Eigen::MatrixXd const& permeability,
MatrixOfVelocityAtIntegrationPoints& darcy_velocity_at_ips) const;*/
const int _gravitational_axis_id; const int _gravitational_axis_id;
const double _gravitational_acceleration; const double _gravitational_acceleration;
const double _reference_temperature; const double _reference_temperature;
......
...@@ -18,10 +18,9 @@ ...@@ -18,10 +18,9 @@
#include "MeshLib/PropertyVector.h" #include "MeshLib/PropertyVector.h"
#include "ProcessLib/Parameter/Parameter.h"
#include "ProcessLib/Parameter/SpatialPosition.h" #include "ProcessLib/Parameter/SpatialPosition.h"
#include "MaterialLib/Fluid/FluidProperty.h" #include "MaterialLib/PorousMedium/Permeability/Permeability.h"
#include "MaterialLib/PorousMedium/Porosity/Porosity.h" #include "MaterialLib/PorousMedium/Porosity/Porosity.h"
#include "MaterialLib/PorousMedium/Storage/Storage.h" #include "MaterialLib/PorousMedium/Storage/Storage.h"
...@@ -54,17 +53,6 @@ double LiquidFlowMaterialProperties::getLiquidDensity(const double p, ...@@ -54,17 +53,6 @@ double LiquidFlowMaterialProperties::getLiquidDensity(const double p,
MaterialLib::Fluid::FluidPropertyType::Density, vars); MaterialLib::Fluid::FluidPropertyType::Density, vars);
} }
double LiquidFlowMaterialProperties::getdLiquidDensity_dT(const double p,
const double T) const
{
ArrayType vars;
vars[static_cast<int>(MaterialLib::Fluid::PropertyVariableType::T)] = T;
vars[static_cast<int>(MaterialLib::Fluid::PropertyVariableType::p)] = p;
return _fluid_properties->getdValue(
MaterialLib::Fluid::FluidPropertyType::Density, vars,
MaterialLib::Fluid::PropertyVariableType::T);
}
double LiquidFlowMaterialProperties::getViscosity(const double p, double LiquidFlowMaterialProperties::getViscosity(const double p,
const double T) const const double T) const
{ {
...@@ -75,24 +63,14 @@ double LiquidFlowMaterialProperties::getViscosity(const double p, ...@@ -75,24 +63,14 @@ double LiquidFlowMaterialProperties::getViscosity(const double p,
MaterialLib::Fluid::FluidPropertyType::Viscosity, vars); MaterialLib::Fluid::FluidPropertyType::Viscosity, vars);
} }
double LiquidFlowMaterialProperties::getHeatCapacity(const double p, double LiquidFlowMaterialProperties::getPorosity(const int material_id,
const double T) const const double t,
{ const SpatialPosition& pos,
ArrayType vars; const double porosity_variable,
vars[static_cast<int>(MaterialLib::Fluid::PropertyVariableType::T)] = T; const double T) const
vars[static_cast<int>(MaterialLib::Fluid::PropertyVariableType::p)] = p;
return _fluid_properties->getValue(
MaterialLib::Fluid::FluidPropertyType::HeatCapacity, vars);
}
double LiquidFlowMaterialProperties::getThermalConductivity(
const double p, const double T) const
{ {
ArrayType vars; return _porosity_models[material_id]->getValue(t, pos, porosity_variable,
vars[static_cast<int>(MaterialLib::Fluid::PropertyVariableType::T)] = T; T);
vars[static_cast<int>(MaterialLib::Fluid::PropertyVariableType::p)] = p;
return _fluid_properties->getValue(
MaterialLib::Fluid::FluidPropertyType::ThermalConductivity, vars);
} }
double LiquidFlowMaterialProperties::getMassCoefficient( double LiquidFlowMaterialProperties::getMassCoefficient(
...@@ -125,17 +103,5 @@ Eigen::MatrixXd const& LiquidFlowMaterialProperties::getPermeability( ...@@ -125,17 +103,5 @@ Eigen::MatrixXd const& LiquidFlowMaterialProperties::getPermeability(
0.0); 0.0);
} }
double LiquidFlowMaterialProperties::getSolidThermalExpansion(
const double t, const SpatialPosition& pos) const
{
return _solid_thermal_expansion(t, pos)[0];
}
double LiquidFlowMaterialProperties::getBiotConstant(
const double t, const SpatialPosition& pos) const
{
return _biot_constant(t, pos)[0];
}
} // end of namespace } // end of namespace
} // end of namespace } // end of namespace
...@@ -19,16 +19,19 @@ ...@@ -19,16 +19,19 @@
#include "MaterialLib/Fluid/FluidProperty.h" #include "MaterialLib/Fluid/FluidProperty.h"
#include "MaterialLib/Fluid/FluidProperties/FluidProperties.h" #include "MaterialLib/Fluid/FluidProperties/FluidProperties.h"
#include "MaterialLib/PorousMedium/Porosity/Porosity.h"
#include "MaterialLib/PorousMedium/Permeability/Permeability.h"
#include "MaterialLib/PorousMedium/Storage/Storage.h"
namespace MaterialLib namespace MaterialLib
{ {
namespace Fluid namespace Fluid
{ {
class FluidProperties; class FluidProperties;
} }
namespace PorousMedium
{
class Permeability;
class Porosity;
class Storage;
}
} }
namespace BaseLib namespace BaseLib
...@@ -44,9 +47,6 @@ class PropertyVector; ...@@ -44,9 +47,6 @@ class PropertyVector;
namespace ProcessLib namespace ProcessLib
{ {
template <typename T>
struct Parameter;
class SpatialPosition; class SpatialPosition;
namespace LiquidFlow namespace LiquidFlow
...@@ -69,18 +69,14 @@ public: ...@@ -69,18 +69,14 @@ public:
std::vector<std::unique_ptr<MaterialLib::PorousMedium::Storage>>&& std::vector<std::unique_ptr<MaterialLib::PorousMedium::Storage>>&&
storage_models, storage_models,
bool const has_material_ids, bool const has_material_ids,
MeshLib::PropertyVector<int> const& material_ids, MeshLib::PropertyVector<int> const& material_ids)
Parameter<double> const& solid_thermal_expansion,
Parameter<double> const& biot_constant)
: _has_material_ids(has_material_ids), : _has_material_ids(has_material_ids),
_material_ids(material_ids), _material_ids(material_ids),
_fluid_properties(std::move(fluid_properties)), _fluid_properties(std::move(fluid_properties)),
_intrinsic_permeability_models( _intrinsic_permeability_models(
std::move(intrinsic_permeability_models)), std::move(intrinsic_permeability_models)),
_porosity_models(std::move(porosity_models)), _porosity_models(std::move(porosity_models)),
_storage_models(std::move(storage_models)), _storage_models(std::move(storage_models))
_solid_thermal_expansion(solid_thermal_expansion),
_biot_constant(biot_constant)
{ {
} }
...@@ -115,26 +111,11 @@ public: ...@@ -115,26 +111,11 @@ public:
double getLiquidDensity(const double p, const double T) const; double getLiquidDensity(const double p, const double T) const;
double getdLiquidDensity_dT(const double p, const double T) const;
double getViscosity(const double p, const double T) const; double getViscosity(const double p, const double T) const;
double getHeatCapacity(const double p, const double T) const;
double getThermalConductivity(const double p, const double T) const;
double getPorosity(const int material_id, const double t, double getPorosity(const int material_id, const double t,
const SpatialPosition& pos, const SpatialPosition& pos,
const double porosity_variable, const double T) const const double porosity_variable, const double T) const;
{
return _porosity_models[material_id]->getValue(t, pos,
porosity_variable, T);
}
double getSolidThermalExpansion(const double t,
const SpatialPosition& pos) const;
double getBiotConstant(const double t, const SpatialPosition& pos) const;
private: private:
/// A flag to indicate whether the reference member, _material_ids, /// A flag to indicate whether the reference member, _material_ids,
...@@ -155,8 +136,6 @@ private: ...@@ -155,8 +136,6 @@ private:
const std::vector<std::unique_ptr<MaterialLib::PorousMedium::Storage>> const std::vector<std::unique_ptr<MaterialLib::PorousMedium::Storage>>
_storage_models; _storage_models;
Parameter<double> const& _solid_thermal_expansion;
Parameter<double> const& _biot_constant;
// Note: For the statistical data of porous media, they could be read from // Note: For the statistical data of porous media, they could be read from
// vtu files directly. This can be done by using property vectors directly. // vtu files directly. This can be done by using property vectors directly.
// Such property vectors will be added here if they are needed. // Such property vectors will be added here if they are needed.
......
...@@ -20,6 +20,7 @@ ...@@ -20,6 +20,7 @@
#include "LiquidFlowLocalAssembler.h" #include "LiquidFlowLocalAssembler.h"
#include "LiquidFlowMaterialProperties.h" #include "LiquidFlowMaterialProperties.h"
#include "MaterialLib/PorousMedium/Permeability/Permeability.h"
#include "MaterialLib/PorousMedium/Porosity/Porosity.h" #include "MaterialLib/PorousMedium/Porosity/Porosity.h"
#include "MaterialLib/PorousMedium/Storage/Storage.h" #include "MaterialLib/PorousMedium/Storage/Storage.h"
...@@ -112,14 +113,5 @@ void LiquidFlowProcess::computeSecondaryVariableConcrete(const double t, ...@@ -112,14 +113,5 @@ void LiquidFlowProcess::computeSecondaryVariableConcrete(const double t,
_coupled_solutions); _coupled_solutions);
} }
void LiquidFlowProcess::setStaggeredCouplingTermToLocalAssemblers()
{
DBUG("Compute the velocity for LiquidFlowProcess.");
/*
GlobalExecutor::executeMemberOnDereferenced(
&LiquidFlowLocalAssemblerInterface::setStaggeredCouplingTerm,
_local_assemblers, _coupled_solutions);*/
}
} // end of namespace } // end of namespace
} // end of namespace } // end of namespace
...@@ -83,13 +83,6 @@ public: ...@@ -83,13 +83,6 @@ public:
return _gravitational_acceleration; return _gravitational_acceleration;
} }
LiquidFlowMaterialProperties* getLiquidFlowMaterialProperties() const
{
return _material_properties.get();
}
void setStaggeredCouplingTermToLocalAssemblers() override;
private: private:
void initializeConcreteProcess( void initializeConcreteProcess(
NumLib::LocalToGlobalIndexMap const& dof_table, NumLib::LocalToGlobalIndexMap const& dof_table,
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment