7 #ifndef __readiff_matrix_h 8 #define __readiff_matrix_h 10 #include <amandus/integrator.h> 11 #include <deal.II/integrators/divergence.h> 12 #include <deal.II/integrators/l2.h> 13 #include <deal.II/integrators/laplace.h> 14 #include <deal.II/meshworker/integration_info.h> 15 #include <readiff/parameters.h> 42 Matrix(
const Parameters& par);
44 virtual void cell(MeshWorker::DoFInfo<dim>& dinfo, MeshWorker::IntegrationInfo<dim>& info)
const;
45 virtual void boundary(MeshWorker::DoFInfo<dim>& dinfo,
46 MeshWorker::IntegrationInfo<dim>& info)
const;
47 virtual void face(MeshWorker::DoFInfo<dim>& dinfo1, MeshWorker::DoFInfo<dim>& dinfo2,
48 MeshWorker::IntegrationInfo<dim>& info1,
49 MeshWorker::IntegrationInfo<dim>& info2)
const;
52 SmartPointer<const Parameters, class Matrix<dim>>
parameters;
59 this->use_boundary =
false;
64 Matrix<dim>::cell(MeshWorker::DoFInfo<dim>& dinfo, MeshWorker::IntegrationInfo<dim>& info)
const 66 AssertDimension(dinfo.n_matrices(), 4);
72 if (info.values.size() > 0)
74 AssertDimension(info.values[0].size(), 2);
75 AssertDimension(info.values[0][0].size(), info.fe_values(0).n_quadrature_points);
76 AssertDimension(info.values[0][1].size(), info.fe_values(0).n_quadrature_points);
77 std::vector<double> Du_ru(info.fe_values(0).n_quadrature_points);
78 std::vector<double> Dv_ru(info.fe_values(0).n_quadrature_points);
79 std::vector<double> Dv_rv(info.fe_values(0).n_quadrature_points);
80 std::vector<double> Du_rv(info.fe_values(0).n_quadrature_points);
81 for (
unsigned int k = 0; k < Du_ru.size(); ++k)
83 const double u = info.values[0][0][k];
84 const double v = info.values[0][1][k];
85 Du_ru[k] = p.B1 + p.D1 * v + 2. * p.E1 * u + 2. * p.G1 * u * v + p.H1 * v * v;
86 Dv_ru[k] = p.C1 + p.D1 * u + 2. * p.F1 * v + 2. * p.H1 * u * v + p.G1 * u * u;
87 Dv_rv[k] = p.C2 + p.D2 * u + 2. * p.F2 * v + 2. * p.H2 * u * v + p.G2 * u * u;
88 Du_rv[k] = p.B2 + p.D2 * v + 2. * p.E2 * u + 2. * p.G2 * u * v + p.H2 * v * v;
106 MeshWorker::IntegrationInfo<dim>& info1,
107 MeshWorker::IntegrationInfo<dim>& info2)
const 109 const unsigned int deg = info1.fe_values(0).get_fe().tensor_degree();
111 Laplace::ip_matrix(dinfo1.matrix(0,
false).matrix,
112 dinfo1.matrix(0,
true).matrix,
113 dinfo2.matrix(0,
true).matrix,
114 dinfo2.matrix(0,
false).matrix,
117 Laplace::compute_penalty(dinfo1, dinfo2, deg, deg),
119 Laplace::ip_matrix(dinfo1.matrix(3,
false).matrix,
120 dinfo1.matrix(3,
true).matrix,
121 dinfo2.matrix(3,
true).matrix,
122 dinfo2.matrix(3,
false).matrix,
125 Laplace::compute_penalty(dinfo1, dinfo2, deg, deg),
void weighted_mass_matrix(dealii::FullMatrix< double > &M, const dealii::FEValuesBase< dim > &fe, const dealii::TensorFunction< 2, dim > &weight)
Definition: darcy/integrators.h:41
Definition: readiff/matrix.h:39
void cell_matrix(dealii::FullMatrix< double > &M, const dealii::FEValuesBase< dim > &fe, const dealii::VectorSlice< const std::vector< std::vector< dealii::Tensor< 1, dim >>>> &input, double lambda=0., double mu=1.)
Definition: matrix_integrators.h:23
virtual void cell(MeshWorker::DoFInfo< dim > &dinfo, MeshWorker::IntegrationInfo< dim > &info) const
Definition: readiff/matrix.h:64
virtual void face(MeshWorker::DoFInfo< dim > &dinfo1, MeshWorker::DoFInfo< dim > &dinfo2, MeshWorker::IntegrationInfo< dim > &info1, MeshWorker::IntegrationInfo< dim > &info2) const
Definition: readiff/matrix.h:105
Definition: readiff/matrix.h:31
virtual void boundary(MeshWorker::DoFInfo< dim > &dinfo, MeshWorker::IntegrationInfo< dim > &info) const
Definition: readiff/matrix.h:99
SmartPointer< const Parameters, class Matrix< dim > > parameters
Definition: readiff/matrix.h:52
Definition: integrator.h:29