16#ifndef dealii_integrators_laplace_h
17#define dealii_integrators_laplace_h
54 const double factor = 1.)
56 const unsigned int n_dofs = fe.dofs_per_cell;
57 const unsigned int n_components = fe.get_fe().n_components();
59 for (
unsigned int k = 0;
k < fe.n_quadrature_points; ++
k)
61 const double dx = fe.JxW(
k) * factor;
62 for (
unsigned int i = 0; i < n_dofs; ++i)
65 for (
unsigned int d = 0; d < n_components; ++d)
66 Mii += dx * (fe.shape_grad_component(i,
k, d) *
67 fe.shape_grad_component(i,
k, d));
71 for (
unsigned int j = i + 1;
j < n_dofs; ++
j)
74 for (
unsigned int d = 0; d < n_components; ++d)
75 Mij += dx * (fe.shape_grad_component(
j,
k, d) *
76 fe.shape_grad_component(i,
k, d));
97 const unsigned int nq = fe.n_quadrature_points;
98 const unsigned int n_dofs = fe.dofs_per_cell;
103 for (
unsigned int k = 0;
k <
nq; ++
k)
105 const double dx = factor * fe.JxW(
k);
106 for (
unsigned int i = 0; i < n_dofs; ++i)
107 result(i) += dx * (input[
k] * fe.shape_grad(i,
k));
124 const unsigned int nq = fe.n_quadrature_points;
125 const unsigned int n_dofs = fe.dofs_per_cell;
126 const unsigned int n_comp = fe.get_fe().n_components();
132 for (
unsigned int k = 0;
k <
nq; ++
k)
134 const double dx = factor * fe.JxW(
k);
135 for (
unsigned int i = 0; i < n_dofs; ++i)
136 for (
unsigned int d = 0; d <
n_comp; ++d)
139 dx * (input[d][
k] * fe.shape_grad_component(i,
k, d));
162 const unsigned int n_dofs = fe.dofs_per_cell;
163 const unsigned int n_comp = fe.get_fe().n_components();
168 for (
unsigned int k = 0;
k < fe.n_quadrature_points; ++
k)
170 const double dx = fe.JxW(
k) * factor;
172 for (
unsigned int i = 0; i < n_dofs; ++i)
173 for (
unsigned int j = 0;
j < n_dofs; ++
j)
174 for (
unsigned int d = 0; d <
n_comp; ++d)
175 M(i,
j) += dx * (2. * fe.shape_value_component(i,
k, d) *
176 penalty * fe.shape_value_component(
j,
k, d) -
177 (n * fe.shape_grad_component(i,
k, d)) *
178 fe.shape_value_component(
j,
k, d) -
179 (n * fe.shape_grad_component(
j,
k, d)) *
180 fe.shape_value_component(i,
k, d));
203 const unsigned int n_dofs = fe.dofs_per_cell;
208 for (
unsigned int k = 0;
k < fe.n_quadrature_points; ++
k)
210 const double dx = fe.JxW(
k) * factor;
212 for (
unsigned int i = 0; i < n_dofs; ++i)
213 for (
unsigned int j = 0;
j < n_dofs; ++
j)
220 for (
unsigned int d = 0; d < dim; ++d)
222 udotn += n[d] * fe.shape_value_component(
j,
k, d);
223 vdotn += n[d] * fe.shape_value_component(i,
k, d);
224 ngradun += n * fe.shape_grad_component(
j,
k, d) * n[d];
225 ngradvn += n * fe.shape_grad_component(i,
k, d) * n[d];
228 for (
unsigned int d = 0; d < dim; ++d)
231 fe.shape_value_component(i,
k, d) -
vdotn * n[d];
233 n * fe.shape_grad_component(i,
k, d) -
ngradvn * n[d];
235 fe.shape_value_component(
j,
k, d) -
udotn * n[d];
237 n * fe.shape_grad_component(
j,
k, d) -
ngradun * n[d];
263 const std::vector<double> & input,
265 const std::vector<double> & data,
269 const unsigned int n_dofs = fe.dofs_per_cell;
274 for (
unsigned int k = 0;
k < fe.n_quadrature_points; ++
k)
276 const double dx = factor * fe.JxW(
k);
278 for (
unsigned int i = 0; i < n_dofs; ++i)
280 const double dnv = fe.shape_grad(i,
k) * n;
282 const double v = fe.shape_value(i,
k);
283 const double u = input[
k];
284 const double g = data[
k];
310 const ArrayView<
const std::vector<double>> & input,
312 const ArrayView<
const std::vector<double>> & data,
316 const unsigned int n_dofs = fe.dofs_per_cell;
317 const unsigned int n_comp = fe.get_fe().n_components();
322 for (
unsigned int k = 0;
k < fe.n_quadrature_points; ++
k)
324 const double dx = factor * fe.JxW(
k);
326 for (
unsigned int i = 0; i < n_dofs; ++i)
327 for (
unsigned int d = 0; d <
n_comp; ++d)
329 const double dnv = fe.shape_grad_component(i,
k, d) * n;
331 const double v = fe.shape_value_component(i,
k, d);
332 const double u = input[d][
k];
333 const double g = data[d][
k];
369 const unsigned int n_dofs =
fe1.dofs_per_cell;
383 for (
unsigned int k = 0;
k <
fe1.n_quadrature_points; ++
k)
385 const double dx =
fe1.JxW(
k);
387 for (
unsigned int d = 0; d <
fe1.get_fe().n_components(); ++d)
389 for (
unsigned int i = 0; i < n_dofs; ++i)
391 for (
unsigned int j = 0;
j < n_dofs; ++
j)
393 const double vi =
fe1.shape_value_component(i,
k, d);
394 const double dnvi = n *
fe1.shape_grad_component(i,
k, d);
395 const double ve =
fe2.shape_value_component(i,
k, d);
396 const double dnve = n *
fe2.shape_grad_component(i,
k, d);
397 const double ui =
fe1.shape_value_component(
j,
k, d);
398 const double dnui = n *
fe1.shape_grad_component(
j,
k, d);
399 const double ue =
fe2.shape_value_component(
j,
k, d);
400 const double dnue = n *
fe2.shape_grad_component(
j,
k, d);
442 const unsigned int n_dofs =
fe1.dofs_per_cell;
458 for (
unsigned int k = 0;
k <
fe1.n_quadrature_points; ++
k)
460 const double dx =
fe1.JxW(
k);
462 for (
unsigned int i = 0; i < n_dofs; ++i)
467 for (
unsigned int j = 0;
j < n_dofs; ++
j)
479 for (
unsigned int d = 0; d < dim; ++d)
481 u1dotn += n[d] *
fe1.shape_value_component(
j,
k, d);
482 v1dotn += n[d] *
fe1.shape_value_component(i,
k, d);
483 u2dotn += n[d] *
fe2.shape_value_component(
j,
k, d);
484 v2dotn += n[d] *
fe2.shape_value_component(i,
k, d);
487 ngradv1n += n *
fe1.shape_grad_component(i,
k, d) * n[d];
489 ngradv2n += n *
fe2.shape_grad_component(i,
k, d) * n[d];
495 for (
unsigned int d = 0; d < dim; ++d)
498 fe1.shape_value_component(i,
k, d) -
v1dotn * n[d];
500 n *
fe1.shape_grad_component(i,
k, d) -
ngradv1n * n[d];
503 fe2.shape_value_component(i,
k, d) -
v2dotn * n[d];
505 n *
fe2.shape_grad_component(i,
k, d) -
ngradv2n * n[d];
508 fe1.shape_value_component(
j,
k, d) -
u1dotn * n[d];
513 fe2.shape_value_component(
j,
k, d) -
u2dotn * n[d];
548 const std::vector<double> &
input1,
550 const std::vector<double> &
input2,
565 const unsigned int n_dofs =
fe1.dofs_per_cell;
567 for (
unsigned int k = 0;
k <
fe1.n_quadrature_points; ++
k)
569 const double dx =
fe1.JxW(
k);
572 for (
unsigned int i = 0; i < n_dofs; ++i)
574 const double vi =
fe1.shape_value(i,
k);
577 const double ve =
fe2.shape_value(i,
k);
623 const unsigned int n_comp =
fe1.get_fe().n_components();
624 const unsigned int n1 =
fe1.dofs_per_cell;
636 for (
unsigned int k = 0;
k <
fe1.n_quadrature_points; ++
k)
638 const double dx =
fe1.JxW(
k);
641 for (
unsigned int i = 0; i <
n1; ++i)
642 for (
unsigned int d = 0; d <
n_comp; ++d)
644 const double vi =
fe1.shape_value_component(i,
k, d);
647 const double ve =
fe2.shape_value_component(i,
k, d);
683 template <
int dim,
int spacedim,
typename number>
699 if (
dinfo1.cell->has_children() ^
dinfo2.cell->has_children())
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define Assert(cond, exc)
#define AssertVectorVectorDimension(VEC, DIM1, DIM2)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
void ip_tangential_matrix(FullMatrix< double > &M11, FullMatrix< double > &M12, FullMatrix< double > &M21, FullMatrix< double > &M22, const FEValuesBase< dim > &fe1, const FEValuesBase< dim > &fe2, double penalty, double factor1=1., double factor2=-1.)
void cell_residual(Vector< double > &result, const FEValuesBase< dim > &fe, const std::vector< Tensor< 1, dim > > &input, double factor=1.)
void cell_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const double factor=1.)
double compute_penalty(const MeshWorker::DoFInfo< dim, spacedim, number > &dinfo1, const MeshWorker::DoFInfo< dim, spacedim, number > &dinfo2, unsigned int deg1, unsigned int deg2)
void nitsche_residual(Vector< double > &result, const FEValuesBase< dim > &fe, const std::vector< double > &input, const std::vector< Tensor< 1, dim > > &Dinput, const std::vector< double > &data, double penalty, double factor=1.)
void ip_matrix(FullMatrix< double > &M11, FullMatrix< double > &M12, FullMatrix< double > &M21, FullMatrix< double > &M22, const FEValuesBase< dim > &fe1, const FEValuesBase< dim > &fe2, double penalty, double factor1=1., double factor2=-1.)
void nitsche_tangential_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, double penalty, double factor=1.)
void ip_residual(Vector< double > &result1, Vector< double > &result2, const FEValuesBase< dim > &fe1, const FEValuesBase< dim > &fe2, const std::vector< double > &input1, const std::vector< Tensor< 1, dim > > &Dinput1, const std::vector< double > &input2, const std::vector< Tensor< 1, dim > > &Dinput2, double pen, double int_factor=1., double ext_factor=-1.)
void nitsche_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, double penalty, double factor=1.)
Library of integrals over cells and faces.