Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
laplace.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2010 - 2019 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_integrators_laplace_h
17 #define dealii_integrators_laplace_h
18 
19 
20 #include <deal.II/base/config.h>
21 
22 #include <deal.II/base/exceptions.h>
23 #include <deal.II/base/quadrature.h>
24 
25 #include <deal.II/fe/fe_values.h>
26 #include <deal.II/fe/mapping.h>
27 
28 #include <deal.II/lac/full_matrix.h>
29 
30 #include <deal.II/meshworker/dof_info.h>
31 
32 DEAL_II_NAMESPACE_OPEN
33 
34 namespace LocalIntegrators
35 {
43  namespace Laplace
44  {
55  template <int dim>
56  void
58  const FEValuesBase<dim> &fe,
59  const double factor = 1.)
60  {
61  const unsigned int n_dofs = fe.dofs_per_cell;
62  const unsigned int n_components = fe.get_fe().n_components();
63 
64  for (unsigned int k = 0; k < fe.n_quadrature_points; ++k)
65  {
66  const double dx = fe.JxW(k) * factor;
67  for (unsigned int i = 0; i < n_dofs; ++i)
68  {
69  double Mii = 0.0;
70  for (unsigned int d = 0; d < n_components; ++d)
71  Mii += dx * (fe.shape_grad_component(i, k, d) *
72  fe.shape_grad_component(i, k, d));
73 
74  M(i, i) += Mii;
75 
76  for (unsigned int j = i + 1; j < n_dofs; ++j)
77  {
78  double Mij = 0.0;
79  for (unsigned int d = 0; d < n_components; ++d)
80  Mij += dx * (fe.shape_grad_component(j, k, d) *
81  fe.shape_grad_component(i, k, d));
82 
83  M(i, j) += Mij;
84  M(j, i) += Mij;
85  }
86  }
87  }
88  }
89 
95  template <int dim>
96  inline void
98  const FEValuesBase<dim> & fe,
99  const std::vector<Tensor<1, dim>> &input,
100  double factor = 1.)
101  {
102  const unsigned int nq = fe.n_quadrature_points;
103  const unsigned int n_dofs = fe.dofs_per_cell;
104  Assert(input.size() == nq, ExcDimensionMismatch(input.size(), nq));
105  Assert(result.size() == n_dofs,
106  ExcDimensionMismatch(result.size(), n_dofs));
107 
108  for (unsigned int k = 0; k < nq; ++k)
109  {
110  const double dx = factor * fe.JxW(k);
111  for (unsigned int i = 0; i < n_dofs; ++i)
112  result(i) += dx * (input[k] * fe.shape_grad(i, k));
113  }
114  }
115 
116 
122  template <int dim>
123  inline void
125  const FEValuesBase<dim> & fe,
126  const ArrayView<const std::vector<Tensor<1, dim>>> &input,
127  double factor = 1.)
128  {
129  const unsigned int nq = fe.n_quadrature_points;
130  const unsigned int n_dofs = fe.dofs_per_cell;
131  const unsigned int n_comp = fe.get_fe().n_components();
132 
134  Assert(result.size() == n_dofs,
135  ExcDimensionMismatch(result.size(), n_dofs));
136 
137  for (unsigned int k = 0; k < nq; ++k)
138  {
139  const double dx = factor * fe.JxW(k);
140  for (unsigned int i = 0; i < n_dofs; ++i)
141  for (unsigned int d = 0; d < n_comp; ++d)
142  {
143  result(i) +=
144  dx * (input[d][k] * fe.shape_grad_component(i, k, d));
145  }
146  }
147  }
148 
149 
163  template <int dim>
164  void
166  const FEValuesBase<dim> &fe,
167  double penalty,
168  double factor = 1.)
169  {
170  const unsigned int n_dofs = fe.dofs_per_cell;
171  const unsigned int n_comp = fe.get_fe().n_components();
172 
173  Assert(M.m() == n_dofs, ExcDimensionMismatch(M.m(), n_dofs));
174  Assert(M.n() == n_dofs, ExcDimensionMismatch(M.n(), n_dofs));
175 
176  for (unsigned int k = 0; k < fe.n_quadrature_points; ++k)
177  {
178  const double dx = fe.JxW(k) * factor;
179  const Tensor<1, dim> n = fe.normal_vector(k);
180  for (unsigned int i = 0; i < n_dofs; ++i)
181  for (unsigned int j = 0; j < n_dofs; ++j)
182  for (unsigned int d = 0; d < n_comp; ++d)
183  M(i, j) += dx * (2. * fe.shape_value_component(i, k, d) *
184  penalty * fe.shape_value_component(j, k, d) -
185  (n * fe.shape_grad_component(i, k, d)) *
186  fe.shape_value_component(j, k, d) -
187  (n * fe.shape_grad_component(j, k, d)) *
188  fe.shape_value_component(i, k, d));
189  }
190  }
191 
207  template <int dim>
208  void
210  const FEValuesBase<dim> &fe,
211  double penalty,
212  double factor = 1.)
213  {
214  const unsigned int n_dofs = fe.dofs_per_cell;
215  AssertDimension(fe.get_fe().n_components(), dim);
216  Assert(M.m() == n_dofs, ExcDimensionMismatch(M.m(), n_dofs));
217  Assert(M.n() == n_dofs, ExcDimensionMismatch(M.n(), n_dofs));
218 
219  for (unsigned int k = 0; k < fe.n_quadrature_points; ++k)
220  {
221  const double dx = fe.JxW(k) * factor;
222  const Tensor<1, dim> n = fe.normal_vector(k);
223  for (unsigned int i = 0; i < n_dofs; ++i)
224  for (unsigned int j = 0; j < n_dofs; ++j)
225  {
226  double udotn = 0.;
227  double vdotn = 0.;
228  double ngradun = 0.;
229  double ngradvn = 0.;
230 
231  for (unsigned int d = 0; d < dim; ++d)
232  {
233  udotn += n[d] * fe.shape_value_component(j, k, d);
234  vdotn += n[d] * fe.shape_value_component(i, k, d);
235  ngradun += n * fe.shape_grad_component(j, k, d) * n[d];
236  ngradvn += n * fe.shape_grad_component(i, k, d) * n[d];
237  }
238 
239  for (unsigned int d = 0; d < dim; ++d)
240  {
241  const double v_t =
242  fe.shape_value_component(i, k, d) - vdotn * n[d];
243  const double dnv_t =
244  n * fe.shape_grad_component(i, k, d) - ngradvn * n[d];
245  const double u_t =
246  fe.shape_value_component(j, k, d) - udotn * n[d];
247  const double dnu_t =
248  n * fe.shape_grad_component(j, k, d) - ngradun * n[d];
249 
250  M(i, j) += dx * (2. * penalty * u_t * v_t - dnu_t * v_t -
251  dnv_t * u_t);
252  }
253  }
254  }
255  }
256 
273  template <int dim>
274  void
276  const FEValuesBase<dim> & fe,
277  const std::vector<double> & input,
278  const std::vector<Tensor<1, dim>> &Dinput,
279  const std::vector<double> & data,
280  double penalty,
281  double factor = 1.)
282  {
283  const unsigned int n_dofs = fe.dofs_per_cell;
284  AssertDimension(input.size(), fe.n_quadrature_points);
285  AssertDimension(Dinput.size(), fe.n_quadrature_points);
286  AssertDimension(data.size(), fe.n_quadrature_points);
287 
288  for (unsigned int k = 0; k < fe.n_quadrature_points; ++k)
289  {
290  const double dx = factor * fe.JxW(k);
291  const Tensor<1, dim> n = fe.normal_vector(k);
292  for (unsigned int i = 0; i < n_dofs; ++i)
293  {
294  const double dnv = fe.shape_grad(i, k) * n;
295  const double dnu = Dinput[k] * n;
296  const double v = fe.shape_value(i, k);
297  const double u = input[k];
298  const double g = data[k];
299 
300  result(i) +=
301  dx * (2. * penalty * (u - g) * v - dnv * (u - g) - dnu * v);
302  }
303  }
304  }
305 
323  template <int dim>
324  void
326  const FEValuesBase<dim> & fe,
327  const ArrayView<const std::vector<double>> & input,
328  const ArrayView<const std::vector<Tensor<1, dim>>> &Dinput,
329  const ArrayView<const std::vector<double>> & data,
330  double penalty,
331  double factor = 1.)
332  {
333  const unsigned int n_dofs = fe.dofs_per_cell;
334  const unsigned int n_comp = fe.get_fe().n_components();
338 
339  for (unsigned int k = 0; k < fe.n_quadrature_points; ++k)
340  {
341  const double dx = factor * fe.JxW(k);
342  const Tensor<1, dim> n = fe.normal_vector(k);
343  for (unsigned int i = 0; i < n_dofs; ++i)
344  for (unsigned int d = 0; d < n_comp; ++d)
345  {
346  const double dnv = fe.shape_grad_component(i, k, d) * n;
347  const double dnu = Dinput[d][k] * n;
348  const double v = fe.shape_value_component(i, k, d);
349  const double u = input[d][k];
350  const double g = data[d][k];
351 
352  result(i) +=
353  dx * (2. * penalty * (u - g) * v - dnv * (u - g) - dnu * v);
354  }
355  }
356  }
357 
377  template <int dim>
378  void
380  FullMatrix<double> & M12,
381  FullMatrix<double> & M21,
382  FullMatrix<double> & M22,
383  const FEValuesBase<dim> &fe1,
384  const FEValuesBase<dim> &fe2,
385  double penalty,
386  double factor1 = 1.,
387  double factor2 = -1.)
388  {
389  const unsigned int n_dofs = fe1.dofs_per_cell;
390  AssertDimension(M11.n(), n_dofs);
391  AssertDimension(M11.m(), n_dofs);
392  AssertDimension(M12.n(), n_dofs);
393  AssertDimension(M12.m(), n_dofs);
394  AssertDimension(M21.n(), n_dofs);
395  AssertDimension(M21.m(), n_dofs);
396  AssertDimension(M22.n(), n_dofs);
397  AssertDimension(M22.m(), n_dofs);
398 
399  const double nui = factor1;
400  const double nue = (factor2 < 0) ? factor1 : factor2;
401  const double nu = .5 * (nui + nue);
402 
403  for (unsigned int k = 0; k < fe1.n_quadrature_points; ++k)
404  {
405  const double dx = fe1.JxW(k);
406  const Tensor<1, dim> n = fe1.normal_vector(k);
407  for (unsigned int d = 0; d < fe1.get_fe().n_components(); ++d)
408  {
409  for (unsigned int i = 0; i < n_dofs; ++i)
410  {
411  for (unsigned int j = 0; j < n_dofs; ++j)
412  {
413  const double vi = fe1.shape_value_component(i, k, d);
414  const double dnvi = n * fe1.shape_grad_component(i, k, d);
415  const double ve = fe2.shape_value_component(i, k, d);
416  const double dnve = n * fe2.shape_grad_component(i, k, d);
417  const double ui = fe1.shape_value_component(j, k, d);
418  const double dnui = n * fe1.shape_grad_component(j, k, d);
419  const double ue = fe2.shape_value_component(j, k, d);
420  const double dnue = n * fe2.shape_grad_component(j, k, d);
421  M11(i, j) +=
422  dx * (-.5 * nui * dnvi * ui - .5 * nui * dnui * vi +
423  nu * penalty * ui * vi);
424  M12(i, j) +=
425  dx * (.5 * nui * dnvi * ue - .5 * nue * dnue * vi -
426  nu * penalty * vi * ue);
427  M21(i, j) +=
428  dx * (-.5 * nue * dnve * ui + .5 * nui * dnui * ve -
429  nu * penalty * ui * ve);
430  M22(i, j) +=
431  dx * (.5 * nue * dnve * ue + .5 * nue * dnue * ve +
432  nu * penalty * ue * ve);
433  }
434  }
435  }
436  }
437  }
438 
453  template <int dim>
454  void
456  FullMatrix<double> & M12,
457  FullMatrix<double> & M21,
458  FullMatrix<double> & M22,
459  const FEValuesBase<dim> &fe1,
460  const FEValuesBase<dim> &fe2,
461  double penalty,
462  double factor1 = 1.,
463  double factor2 = -1.)
464  {
465  const unsigned int n_dofs = fe1.dofs_per_cell;
466  AssertDimension(fe1.get_fe().n_components(), dim);
467  AssertDimension(fe2.get_fe().n_components(), dim);
468  AssertDimension(M11.n(), n_dofs);
469  AssertDimension(M11.m(), n_dofs);
470  AssertDimension(M12.n(), n_dofs);
471  AssertDimension(M12.m(), n_dofs);
472  AssertDimension(M21.n(), n_dofs);
473  AssertDimension(M21.m(), n_dofs);
474  AssertDimension(M22.n(), n_dofs);
475  AssertDimension(M22.m(), n_dofs);
476 
477  const double nui = factor1;
478  const double nue = (factor2 < 0) ? factor1 : factor2;
479  const double nu = .5 * (nui + nue);
480 
481  for (unsigned int k = 0; k < fe1.n_quadrature_points; ++k)
482  {
483  const double dx = fe1.JxW(k);
484  const Tensor<1, dim> n = fe1.normal_vector(k);
485  for (unsigned int i = 0; i < n_dofs; ++i)
486  {
487  // We compute the tangential component by subtracting
488  // the normal component from the total. Thus, compute
489  // the normal component first
490  for (unsigned int j = 0; j < n_dofs; ++j)
491  {
492  double u1dotn = 0.;
493  double v1dotn = 0.;
494  double u2dotn = 0.;
495  double v2dotn = 0.;
496 
497  double ngradu1n = 0.;
498  double ngradv1n = 0.;
499  double ngradu2n = 0.;
500  double ngradv2n = 0.;
501 
502  for (unsigned int d = 0; d < dim; ++d)
503  {
504  u1dotn += n[d] * fe1.shape_value_component(j, k, d);
505  v1dotn += n[d] * fe1.shape_value_component(i, k, d);
506  u2dotn += n[d] * fe2.shape_value_component(j, k, d);
507  v2dotn += n[d] * fe2.shape_value_component(i, k, d);
508 
509  ngradu1n += n * fe1.shape_grad_component(j, k, d) * n[d];
510  ngradv1n += n * fe1.shape_grad_component(i, k, d) * n[d];
511  ngradu2n += n * fe2.shape_grad_component(j, k, d) * n[d];
512  ngradv2n += n * fe2.shape_grad_component(i, k, d) * n[d];
513  }
514 
515  // The following code is equal to ip_matrix() with
516  // the only exception that all variables introduced
517  // below denote tangential components.
518  for (unsigned int d = 0; d < dim; ++d)
519  {
520  const double vi =
521  fe1.shape_value_component(i, k, d) - v1dotn * n[d];
522  const double dnvi =
523  n * fe1.shape_grad_component(i, k, d) - ngradv1n * n[d];
524 
525  const double ve =
526  fe2.shape_value_component(i, k, d) - v2dotn * n[d];
527  const double dnve =
528  n * fe2.shape_grad_component(i, k, d) - ngradv2n * n[d];
529 
530  const double ui =
531  fe1.shape_value_component(j, k, d) - u1dotn * n[d];
532  const double dnui =
533  n * fe1.shape_grad_component(j, k, d) - ngradu1n * n[d];
534 
535  const double ue =
536  fe2.shape_value_component(j, k, d) - u2dotn * n[d];
537  const double dnue =
538  n * fe2.shape_grad_component(j, k, d) - ngradu2n * n[d];
539 
540  M11(i, j) +=
541  dx * (-.5 * nui * dnvi * ui - .5 * nui * dnui * vi +
542  nu * penalty * ui * vi);
543  M12(i, j) +=
544  dx * (.5 * nui * dnvi * ue - .5 * nue * dnue * vi -
545  nu * penalty * vi * ue);
546  M21(i, j) +=
547  dx * (-.5 * nue * dnve * ui + .5 * nui * dnui * ve -
548  nu * penalty * ui * ve);
549  M22(i, j) +=
550  dx * (.5 * nue * dnve * ue + .5 * nue * dnue * ve +
551  nu * penalty * ue * ve);
552  }
553  }
554  }
555  }
556  }
557 
568  template <int dim>
569  void
571  Vector<double> & result2,
572  const FEValuesBase<dim> & fe1,
573  const FEValuesBase<dim> & fe2,
574  const std::vector<double> & input1,
575  const std::vector<Tensor<1, dim>> &Dinput1,
576  const std::vector<double> & input2,
577  const std::vector<Tensor<1, dim>> &Dinput2,
578  double pen,
579  double int_factor = 1.,
580  double ext_factor = -1.)
581  {
582  Assert(fe1.get_fe().n_components() == 1,
583  ExcDimensionMismatch(fe1.get_fe().n_components(), 1));
584  Assert(fe2.get_fe().n_components() == 1,
585  ExcDimensionMismatch(fe2.get_fe().n_components(), 1));
586 
587  const double nui = int_factor;
588  const double nue = (ext_factor < 0) ? int_factor : ext_factor;
589  const double penalty = .5 * pen * (nui + nue);
590 
591  const unsigned int n_dofs = fe1.dofs_per_cell;
592 
593  for (unsigned int k = 0; k < fe1.n_quadrature_points; ++k)
594  {
595  const double dx = fe1.JxW(k);
596  const Tensor<1, dim> n = fe1.normal_vector(k);
597 
598  for (unsigned int i = 0; i < n_dofs; ++i)
599  {
600  const double vi = fe1.shape_value(i, k);
601  const Tensor<1, dim> &Dvi = fe1.shape_grad(i, k);
602  const double dnvi = Dvi * n;
603  const double ve = fe2.shape_value(i, k);
604  const Tensor<1, dim> &Dve = fe2.shape_grad(i, k);
605  const double dnve = Dve * n;
606 
607  const double ui = input1[k];
608  const Tensor<1, dim> &Dui = Dinput1[k];
609  const double dnui = Dui * n;
610  const double ue = input2[k];
611  const Tensor<1, dim> &Due = Dinput2[k];
612  const double dnue = Due * n;
613 
614  result1(i) += dx * (-.5 * nui * dnvi * ui - .5 * nui * dnui * vi +
615  penalty * ui * vi);
616  result1(i) += dx * (.5 * nui * dnvi * ue - .5 * nue * dnue * vi -
617  penalty * vi * ue);
618  result2(i) += dx * (-.5 * nue * dnve * ui + .5 * nui * dnui * ve -
619  penalty * ui * ve);
620  result2(i) += dx * (.5 * nue * dnve * ue + .5 * nue * dnue * ve +
621  penalty * ue * ve);
622  }
623  }
624  }
625 
626 
638  template <int dim>
639  void
641  Vector<double> & result2,
642  const FEValuesBase<dim> & fe1,
643  const FEValuesBase<dim> & fe2,
644  const ArrayView<const std::vector<double>> & input1,
645  const ArrayView<const std::vector<Tensor<1, dim>>> &Dinput1,
646  const ArrayView<const std::vector<double>> & input2,
647  const ArrayView<const std::vector<Tensor<1, dim>>> &Dinput2,
648  double pen,
649  double int_factor = 1.,
650  double ext_factor = -1.)
651  {
652  const unsigned int n_comp = fe1.get_fe().n_components();
653  const unsigned int n1 = fe1.dofs_per_cell;
654 
656  AssertVectorVectorDimension(Dinput1, n_comp, fe1.n_quadrature_points);
658  AssertVectorVectorDimension(Dinput2, n_comp, fe2.n_quadrature_points);
659 
660  const double nui = int_factor;
661  const double nue = (ext_factor < 0) ? int_factor : ext_factor;
662  const double penalty = .5 * pen * (nui + nue);
663 
664 
665  for (unsigned int k = 0; k < fe1.n_quadrature_points; ++k)
666  {
667  const double dx = fe1.JxW(k);
668  const Tensor<1, dim> n = fe1.normal_vector(k);
669 
670  for (unsigned int i = 0; i < n1; ++i)
671  for (unsigned int d = 0; d < n_comp; ++d)
672  {
673  const double vi = fe1.shape_value_component(i, k, d);
674  const Tensor<1, dim> &Dvi = fe1.shape_grad_component(i, k, d);
675  const double dnvi = Dvi * n;
676  const double ve = fe2.shape_value_component(i, k, d);
677  const Tensor<1, dim> &Dve = fe2.shape_grad_component(i, k, d);
678  const double dnve = Dve * n;
679 
680  const double ui = input1[d][k];
681  const Tensor<1, dim> &Dui = Dinput1[d][k];
682  const double dnui = Dui * n;
683  const double ue = input2[d][k];
684  const Tensor<1, dim> &Due = Dinput2[d][k];
685  const double dnue = Due * n;
686 
687  result1(i) += dx * (-.5 * nui * dnvi * ui -
688  .5 * nui * dnui * vi + penalty * ui * vi);
689  result1(i) += dx * (.5 * nui * dnvi * ue -
690  .5 * nue * dnue * vi - penalty * vi * ue);
691  result2(i) += dx * (-.5 * nue * dnve * ui +
692  .5 * nui * dnui * ve - penalty * ui * ve);
693  result2(i) += dx * (.5 * nue * dnve * ue +
694  .5 * nue * dnue * ve + penalty * ue * ve);
695  }
696  }
697  }
698 
699 
700 
715  template <int dim, int spacedim, typename number>
716  double
719  unsigned int deg1,
720  unsigned int deg2)
721  {
722  const unsigned int normal1 =
724  const unsigned int normal2 =
726  const unsigned int deg1sq = (deg1 == 0) ? 1 : deg1 * (deg1 + 1);
727  const unsigned int deg2sq = (deg2 == 0) ? 1 : deg2 * (deg2 + 1);
728 
729  double penalty1 = deg1sq / dinfo1.cell->extent_in_direction(normal1);
730  double penalty2 = deg2sq / dinfo2.cell->extent_in_direction(normal2);
731  if (dinfo1.cell->has_children() ^ dinfo2.cell->has_children())
732  {
733  Assert(dinfo1.face == dinfo2.face, ExcInternalError());
734  Assert(dinfo1.face->has_children(), ExcInternalError());
735  penalty1 *= 2;
736  }
737  const double penalty = 0.5 * (penalty1 + penalty2);
738  return penalty;
739  }
740  } // namespace Laplace
741 } // namespace LocalIntegrators
742 
743 
744 DEAL_II_NAMESPACE_CLOSE
745 
746 #endif
size_type m() const
void cell_residual(Vector< double > &result, const FEValuesBase< dim > &fe, const std::vector< Tensor< 1, dim >> &input, double factor=1.)
Definition: laplace.h:97
Triangulation< dim, spacedim >::cell_iterator cell
The current cell.
Definition: dof_info.h:78
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1567
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.)
Definition: laplace.h:379
void nitsche_tangential_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, double penalty, double factor=1.)
Definition: laplace.h:209
const unsigned int dofs_per_cell
Definition: fe_values.h:2106
void cell_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const double factor=1.)
Definition: laplace.h:57
const FiniteElement< dim, spacedim > & get_fe() const
#define AssertVectorVectorDimension(VEC, DIM1, DIM2)
Definition: exceptions.h:1595
double shape_value_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
Library of integrals over cells and faces.
Definition: advection.h:34
size_type n() const
const double & shape_value(const unsigned int function_no, const unsigned int point_no) const
#define Assert(cond, exc)
Definition: exceptions.h:1407
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
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.)
Definition: laplace.h:570
const unsigned int n_quadrature_points
Definition: fe_values.h:2099
void nitsche_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, double penalty, double factor=1.)
Definition: laplace.h:165
Tensor< 1, spacedim > shape_grad_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
double compute_penalty(const MeshWorker::DoFInfo< dim, spacedim, number > &dinfo1, const MeshWorker::DoFInfo< dim, spacedim, number > &dinfo2, unsigned int deg1, unsigned int deg2)
Definition: laplace.h:717
double JxW(const unsigned int quadrature_point) const
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.)
Definition: laplace.h:275
size_type size() const
unsigned int face_number
Definition: dof_info.h:89
const Tensor< 1, spacedim > & shape_grad(const unsigned int function_no, const unsigned int quadrature_point) const
Triangulation< dim, spacedim >::face_iterator face
The current face.
Definition: dof_info.h:81
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.)
Definition: laplace.h:455
static ::ExceptionBase & ExcInternalError()
const Tensor< 1, spacedim > & normal_vector(const unsigned int i) const