5 #ifndef _neo_hookean_two_fields_interface_h_
6 #define _neo_hookean_two_fields_interface_h_
9 #include <deal2lkit/parsed_function.h>
26 template <
int dim,
int spacedim>
31 typedef FEValuesCache<dim,spacedim>
Scratch;
33 typedef Assembly::CopyData::piDoMUSSystem<dim,dim>
CopySystem;
46 template<
typename Number>
50 Number &energy)
const;
52 template<
typename Number>
56 Number &energy)
const;
61 const std::vector<shared_ptr<MAT> >,
71 mutable shared_ptr<TrilinosWrappers::PreconditionAMG> Amg_preconditioner;
72 mutable shared_ptr<TrilinosWrappers::PreconditionJacobi> Mp_preconditioner;
73 mutable shared_ptr<TrilinosWrappers::PreconditionJacobi> T_preconditioner;
77 template <
int dim,
int spacedim>
80 "FESystem[FE_Q(2)^d-FE_DGP(1)]",
81 "u,u,u,p",
"1,1; 1,0",
"1,0; 0,1",
87 template <
int dim,
int spacedim>
88 template<
typename Number>
94 Number alpha = this->alpha;
95 this->reinit(alpha, cell, fe_cache);
97 auto &JxW = fe_cache.get_JxW_values();
101 auto &ps = fe_cache.get_values(
"solution",
"p", pressure, alpha);
102 auto &grad_us = fe_cache.get_gradients(
"solution",
"grad_u",displacement, alpha);
104 const unsigned int n_q_points = ps.size();
107 for (
unsigned int q=0; q<n_q_points; ++q)
109 const Number &p = ps[q];
112 energy += (scalar_product(grad_u,grad_u) +
117 template <
int dim,
int spacedim>
118 template<
typename Number>
122 Number &energy)
const
124 Number alpha = this->alpha;
125 this->reinit(alpha, cell, fe_cache);
128 auto &us = fe_cache.get_values(
"solution",
"u", displacement, alpha);
129 auto &us_dot = fe_cache.get_values(
"solution_dot",
"u_dot", displacement, alpha);
130 auto &Fs = fe_cache.get_deformation_gradients(
"solution",
"Fu", displacement, alpha);
133 auto &ps = fe_cache.get_values(
"solution",
"p", pressure, alpha);
135 auto &JxW = fe_cache.get_JxW_values();
137 const unsigned int n_q_points = ps.size();
140 for (
unsigned int q=0; q<n_q_points; ++q)
146 const Number &p = ps[q];
150 Number Ic = trace(C);
151 Number J = determinant(F);
153 Number psi = (
mu/2.)*(Ic-dim) +p*(J-1.);
154 energy += (u*u_dot + psi)*JxW[q];
159 template <
int dim,
int spacedim>
163 this->add_parameter(prm, &E,
"Young's modulus",
"10.0",
Patterns::Double(0.0));
164 this->add_parameter(prm, &nu,
"Poisson's ratio",
"0.3",
Patterns::Double(0.0));
167 template <
int dim,
int spacedim>
171 mu = E/(2.0*(1.+nu));
172 lambda = (E *nu)/((1.+nu)*(1.-2.*nu));
175 template <
int dim,
int spacedim>
180 const std::vector<shared_ptr<MAT> >,
185 std::vector<std::vector<bool> > constant_modes;
200 Mp_preconditioner->initialize (preconditioner_matrix.
block(1,1));
201 Amg_preconditioner->initialize (preconditioner_matrix.
block(0,0),
206 auto A = linear_operator< TrilinosWrappers::MPI::Vector >( matrix.
block(0,0) );
207 auto Bt = linear_operator< TrilinosWrappers::MPI::Vector >( matrix.
block(0,1) );
209 auto B = linear_operator< TrilinosWrappers::MPI::Vector >( matrix.
block(1,0) );
210 auto ZeroP = 0*linear_operator< TrilinosWrappers::MPI::Vector >( matrix.
block(1,1) );
212 auto Mp = linear_operator< TrilinosWrappers::MPI::Vector >( preconditioner_matrix.
block(1,1) );
221 auto P10 = Schur_inv * B * A_inv;
222 auto P11 = -1 * Schur_inv;
225 system_op = block_operator<2, 2, VEC >({{
234 prec_op = block_operator<2, 2, VEC >({{
Definition: conservative.h:20
double aggregation_threshold
std::vector< std::vector< bool > > constant_modes
ActiveSelector::active_cell_iterator active_cell_iterator
TrilinosWrappers::MPI::BlockVector VEC
Definition: neo_hookean_two_fields.h:34
virtual void parse_parameters_call_back()
Definition: neo_hookean_two_fields.h:168
unsigned int smoother_sweeps
NeoHookeanTwoFieldsInterface()
Definition: neo_hookean_two_fields.h:62
virtual void declare_parameters(ParameterHandler &prm)
Definition: neo_hookean_two_fields.h:116
void preconditioner_energy(const typename DoFHandler< dim, spacedim >::active_cell_iterator &, Scratch &, CopyPreconditioner &, Number &energy) const
Definition: neo_hookean_two_fields.h:89
LinearOperator< Range, Domain > null_operator(const LinearOperator< Range, Domain > &op)
void compute_system_operators(const DoFHandler< dim, spacedim > &, const std::vector< shared_ptr< LATrilinos::BlockMatrix > >, LinearOperator< LATrilinos::VectorType > &, LinearOperator< LATrilinos::VectorType > &) const
Definition: neo_hookean_two_fields.h:131
Assembly::CopyData::piDoMUSSystem< dim, dim > CopySystem
Definition: neo_hookean_two_fields.h:33
LinearOperator< typename Solver::vector_type, typename Solver::vector_type > inverse_operator(const LinearOperator< typename Solver::vector_type, typename Solver::vector_type > &op, Solver &solver, const Preconditioner &preconditioner)
FEValuesCache< dim, spacedim > Scratch
Definition: neo_hookean_two_fields.h:31
TrilinosWrappers::BlockSparseMatrix MAT
Definition: neo_hookean_two_fields.h:35
Definition: neo_hookean_two_fields.h:27
BlockType & block(const unsigned int row, const unsigned int column)
Assembly::CopyData::piDoMUSPreconditioner< dim, dim > CopyPreconditioner
Definition: neo_hookean_two_fields.h:32
void system_energy(const typename DoFHandler< dim, spacedim >::active_cell_iterator &, Scratch &, CopySystem &, Number &energy) const
Definition: neo_hookean_two_fields.h:119
const FiniteElement< dim, spacedim > & get_fe() const
bool higher_order_elements