deal.II version GIT relicensing-1972-g22a7b89abe 2024-10-11 21:20:00+00:00
\(\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// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2010 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_integrators_l2_h
16#define dealii_integrators_l2_h
17
18
19#include <deal.II/base/config.h>
20
23
25#include <deal.II/fe/mapping.h>
26
28
30
32
33namespace LocalIntegrators
34{
40 namespace L2
41 {
55 template <int dim>
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_value_component(i, k, d) *
72 fe.shape_value_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_value_component(j, k, d) *
81 fe.shape_value_component(i, k, d);
82
83 M(i, j) += Mij;
84 M(j, i) += Mij;
85 }
86 }
87 }
88 }
89
106 template <int dim>
109 const FEValuesBase<dim> &fe,
110 const std::vector<double> &weights)
111 {
112 const unsigned int n_dofs = fe.dofs_per_cell;
113 const unsigned int n_components = fe.get_fe().n_components();
114 AssertDimension(M.m(), n_dofs);
115 AssertDimension(M.n(), n_dofs);
116 AssertDimension(weights.size(), fe.n_quadrature_points);
117
118 for (unsigned int k = 0; k < fe.n_quadrature_points; ++k)
119 {
120 const double dx = fe.JxW(k) * weights[k];
121 for (unsigned int i = 0; i < n_dofs; ++i)
122 {
123 double Mii = 0.0;
124 for (unsigned int d = 0; d < n_components; ++d)
125 Mii += dx * fe.shape_value_component(i, k, d) *
126 fe.shape_value_component(i, k, d);
127
128 M(i, i) += Mii;
129
130 for (unsigned int j = i + 1; j < n_dofs; ++j)
131 {
132 double Mij = 0.0;
133 for (unsigned int d = 0; d < n_components; ++d)
134 Mij += dx * fe.shape_value_component(j, k, d) *
135 fe.shape_value_component(i, k, d);
136
137 M(i, j) += Mij;
138 M(j, i) += Mij;
139 }
140 }
141 }
142 }
143
157 template <int dim, typename number>
160 const FEValuesBase<dim> &fe,
161 const std::vector<double> &input,
162 const double factor = 1.)
163 {
164 const unsigned int n_dofs = fe.dofs_per_cell;
165 AssertDimension(result.size(), n_dofs);
167 AssertDimension(input.size(), fe.n_quadrature_points);
168
169 for (unsigned int k = 0; k < fe.n_quadrature_points; ++k)
170 for (unsigned int i = 0; i < n_dofs; ++i)
171 result(i) += fe.JxW(k) * factor * input[k] * fe.shape_value(i, k);
172 }
173
187 template <int dim, typename number>
190 const FEValuesBase<dim> &fe,
191 const ArrayView<const std::vector<double>> &input,
192 const double factor = 1.)
193 {
194 const unsigned int n_dofs = fe.dofs_per_cell;
195 const unsigned int n_components = input.size();
196
197 AssertDimension(result.size(), n_dofs);
198 AssertDimension(input.size(), fe.get_fe().n_components());
199
200 for (unsigned int k = 0; k < fe.n_quadrature_points; ++k)
201 for (unsigned int i = 0; i < n_dofs; ++i)
202 for (unsigned int d = 0; d < n_components; ++d)
203 result(i) += fe.JxW(k) * factor *
204 fe.shape_value_component(i, k, d) * input[d][k];
205 }
206
235 template <int dim>
241 const FEValuesBase<dim> &fe1,
242 const FEValuesBase<dim> &fe2,
243 const double factor1 = 1.,
244 const double factor2 = 1.)
245 {
246 const unsigned int n1_dofs = fe1.n_dofs_per_cell();
247 const unsigned int n2_dofs = fe2.n_dofs_per_cell();
248 const unsigned int n_components = fe1.get_fe().n_components();
249
250 Assert(n1_dofs == n2_dofs, ExcNotImplemented());
251 (void)n2_dofs;
252 AssertDimension(n_components, fe2.get_fe().n_components());
253 AssertDimension(M11.m(), n1_dofs);
254 AssertDimension(M12.m(), n1_dofs);
255 AssertDimension(M21.m(), n2_dofs);
256 AssertDimension(M22.m(), n2_dofs);
257 AssertDimension(M11.n(), n1_dofs);
258 AssertDimension(M12.n(), n2_dofs);
259 AssertDimension(M21.n(), n1_dofs);
260 AssertDimension(M22.n(), n2_dofs);
261
262 for (unsigned int k = 0; k < fe1.n_quadrature_points; ++k)
263 {
264 const double dx = fe1.JxW(k);
265
266 for (unsigned int i = 0; i < n1_dofs; ++i)
267 for (unsigned int j = 0; j < n1_dofs; ++j)
268 for (unsigned int d = 0; d < n_components; ++d)
269 {
270 const double u1 =
271 factor1 * fe1.shape_value_component(j, k, d);
272 const double u2 =
273 -factor2 * fe2.shape_value_component(j, k, d);
274 const double v1 =
275 factor1 * fe1.shape_value_component(i, k, d);
276 const double v2 =
277 -factor2 * fe2.shape_value_component(i, k, d);
278
279 M11(i, j) += dx * u1 * v1;
280 M12(i, j) += dx * u2 * v1;
281 M21(i, j) += dx * u1 * v2;
282 M22(i, j) += dx * u2 * v2;
283 }
284 }
285 }
286 } // namespace L2
287} // namespace LocalIntegrators
288
290
291#endif
const unsigned int dofs_per_cell
double shape_value_component(const unsigned int i, const unsigned int q_point, const unsigned int component) const
const unsigned int n_quadrature_points
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
virtual size_type size() const override
#define DEAL_II_DEPRECATED
Definition config.h:205
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:498
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:499
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:159
void mass_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const double factor=1.)
Definition l2.h:57
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:237
void weighted_mass_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const std::vector< double > &weights)
Definition l2.h:108
Library of integrals over cells and faces.
Definition advection.h:34