Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
l2.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2010 - 2023 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_l2_h
17#define dealii_integrators_l2_h
18
19
20#include <deal.II/base/config.h>
21
24
26#include <deal.II/fe/mapping.h>
27
29
31
33
34namespace LocalIntegrators
35{
41 namespace L2
42 {
56 template <int dim>
57 void
59 const FEValuesBase<dim> &fe,
60 const double factor = 1.)
61 {
62 const unsigned int n_dofs = fe.dofs_per_cell;
63 const unsigned int n_components = fe.get_fe().n_components();
64
65 for (unsigned int k = 0; k < fe.n_quadrature_points; ++k)
66 {
67 const double dx = fe.JxW(k) * factor;
68 for (unsigned int i = 0; i < n_dofs; ++i)
69 {
70 double Mii = 0.0;
71 for (unsigned int d = 0; d < n_components; ++d)
72 Mii += dx * fe.shape_value_component(i, k, d) *
73 fe.shape_value_component(i, k, d);
74
75 M(i, i) += Mii;
76
77 for (unsigned int j = i + 1; j < n_dofs; ++j)
78 {
79 double Mij = 0.0;
80 for (unsigned int d = 0; d < n_components; ++d)
81 Mij += dx * fe.shape_value_component(j, k, d) *
82 fe.shape_value_component(i, k, d);
83
84 M(i, j) += Mij;
85 M(j, i) += Mij;
86 }
87 }
88 }
89 }
90
107 template <int dim>
108 void
110 const FEValuesBase<dim> & fe,
111 const std::vector<double> &weights)
112 {
113 const unsigned int n_dofs = fe.dofs_per_cell;
114 const unsigned int n_components = fe.get_fe().n_components();
115 AssertDimension(M.m(), n_dofs);
116 AssertDimension(M.n(), n_dofs);
117 AssertDimension(weights.size(), fe.n_quadrature_points);
118
119 for (unsigned int k = 0; k < fe.n_quadrature_points; ++k)
120 {
121 const double dx = fe.JxW(k) * weights[k];
122 for (unsigned int i = 0; i < n_dofs; ++i)
123 {
124 double Mii = 0.0;
125 for (unsigned int d = 0; d < n_components; ++d)
126 Mii += dx * fe.shape_value_component(i, k, d) *
127 fe.shape_value_component(i, k, d);
128
129 M(i, i) += Mii;
130
131 for (unsigned int j = i + 1; j < n_dofs; ++j)
132 {
133 double Mij = 0.0;
134 for (unsigned int d = 0; d < n_components; ++d)
135 Mij += dx * fe.shape_value_component(j, k, d) *
136 fe.shape_value_component(i, k, d);
137
138 M(i, j) += Mij;
139 M(j, i) += Mij;
140 }
141 }
142 }
143 }
144
158 template <int dim, typename number>
159 void
161 const FEValuesBase<dim> & fe,
162 const std::vector<double> &input,
163 const double factor = 1.)
164 {
165 const unsigned int n_dofs = fe.dofs_per_cell;
166 AssertDimension(result.size(), n_dofs);
168 AssertDimension(input.size(), fe.n_quadrature_points);
169
170 for (unsigned int k = 0; k < fe.n_quadrature_points; ++k)
171 for (unsigned int i = 0; i < n_dofs; ++i)
172 result(i) += fe.JxW(k) * factor * input[k] * fe.shape_value(i, k);
173 }
174
188 template <int dim, typename number>
189 void
191 const FEValuesBase<dim> & fe,
192 const ArrayView<const std::vector<double>> &input,
193 const double factor = 1.)
194 {
195 const unsigned int n_dofs = fe.dofs_per_cell;
196 const unsigned int n_components = input.size();
197
198 AssertDimension(result.size(), n_dofs);
199 AssertDimension(input.size(), fe.get_fe().n_components());
200
201 for (unsigned int k = 0; k < fe.n_quadrature_points; ++k)
202 for (unsigned int i = 0; i < n_dofs; ++i)
203 for (unsigned int d = 0; d < n_components; ++d)
204 result(i) += fe.JxW(k) * factor *
205 fe.shape_value_component(i, k, d) * input[d][k];
206 }
207
236 template <int dim>
237 void
239 FullMatrix<double> & M12,
240 FullMatrix<double> & M21,
241 FullMatrix<double> & M22,
242 const FEValuesBase<dim> &fe1,
243 const FEValuesBase<dim> &fe2,
244 const double factor1 = 1.,
245 const double factor2 = 1.)
246 {
247 const unsigned int n1_dofs = fe1.n_dofs_per_cell();
248 const unsigned int n2_dofs = fe2.n_dofs_per_cell();
249 const unsigned int n_components = fe1.get_fe().n_components();
250
251 Assert(n1_dofs == n2_dofs, ExcNotImplemented());
252 (void)n2_dofs;
253 AssertDimension(n_components, fe2.get_fe().n_components());
254 AssertDimension(M11.m(), n1_dofs);
255 AssertDimension(M12.m(), n1_dofs);
256 AssertDimension(M21.m(), n2_dofs);
257 AssertDimension(M22.m(), n2_dofs);
258 AssertDimension(M11.n(), n1_dofs);
259 AssertDimension(M12.n(), n2_dofs);
260 AssertDimension(M21.n(), n1_dofs);
261 AssertDimension(M22.n(), n2_dofs);
262
263 for (unsigned int k = 0; k < fe1.n_quadrature_points; ++k)
264 {
265 const double dx = fe1.JxW(k);
266
267 for (unsigned int i = 0; i < n1_dofs; ++i)
268 for (unsigned int j = 0; j < n1_dofs; ++j)
269 for (unsigned int d = 0; d < n_components; ++d)
270 {
271 const double u1 =
272 factor1 * fe1.shape_value_component(j, k, d);
273 const double u2 =
274 -factor2 * fe2.shape_value_component(j, k, d);
275 const double v1 =
276 factor1 * fe1.shape_value_component(i, k, d);
277 const double v2 =
278 -factor2 * fe2.shape_value_component(i, k, d);
279
280 M11(i, j) += dx * u1 * v1;
281 M12(i, j) += dx * u2 * v1;
282 M21(i, j) += dx * u1 * v2;
283 M22(i, j) += dx * u2 * v2;
284 }
285 }
286 }
287 } // namespace L2
288} // namespace LocalIntegrators
289
291
292#endif
const unsigned int dofs_per_cell
Definition fe_values.h:2451
double shape_value_component(const unsigned int i, const unsigned int q_point, const unsigned int component) const
const unsigned int n_quadrature_points
Definition fe_values.h:2433
const FiniteElement< dim, spacedim > & get_fe() const
double JxW(const unsigned int q_point) const
const double & shape_value(const unsigned int i, const unsigned int q_point) const
unsigned int n_components() const
size_type n() const
size_type m() const
size_type size() const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
const unsigned int v1
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
void L2(Vector< number > &result, const FEValuesBase< dim > &fe, const std::vector< double > &input, const double factor=1.)
Definition l2.h:160
void mass_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const double factor=1.)
Definition l2.h:58
void jump_matrix(FullMatrix< double > &M11, FullMatrix< double > &M12, FullMatrix< double > &M21, FullMatrix< double > &M22, const FEValuesBase< dim > &fe1, const FEValuesBase< dim > &fe2, const double factor1=1., const double factor2=1.)
Definition l2.h:238
void weighted_mass_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const std::vector< double > &weights)
Definition l2.h:109
Library of integrals over cells and faces.
Definition advection.h:35