Reference documentation for deal.II version 8.1.0
l2.h
1 // ---------------------------------------------------------------------
2 // @f$Id: l2.h 30036 2013-07-18 16:55:32Z maier @f$
3 //
4 // Copyright (C) 2010 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__integrators_l2_h
18 #define __deal2__integrators_l2_h
19 
20 
21 #include <deal.II/base/config.h>
23 #include <deal.II/base/quadrature.h>
24 #include <deal.II/lac/full_matrix.h>
25 #include <deal.II/fe/mapping.h>
26 #include <deal.II/fe/fe_values.h>
27 #include <deal.II/meshworker/dof_info.h>
28 
30 
31 namespace LocalIntegrators
32 {
40  namespace L2
41  {
56  template <int dim>
57  void mass_matrix (
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 
69  for (unsigned int i=0; i<n_dofs; ++i)
70  for (unsigned int j=0; j<n_dofs; ++j)
71  for (unsigned int d=0; d<n_components; ++d)
72  M(i,j) += dx
73  * fe.shape_value_component(j,k,d)
74  * fe.shape_value_component(i,k,d);
75  }
76  }
77 
88  template <int dim, typename number>
89  void L2 (
90  Vector<number> &result,
91  const FEValuesBase<dim> &fe,
92  const std::vector<double> &input,
93  const double factor = 1.)
94  {
95  const unsigned int n_dofs = fe.dofs_per_cell;
96  AssertDimension(result.size(), n_dofs);
97  AssertDimension(fe.get_fe().n_components(), 1);
98  AssertDimension(input.size(), fe.n_quadrature_points);
99 
100  for (unsigned int k=0; k<fe.n_quadrature_points; ++k)
101  for (unsigned int i=0; i<n_dofs; ++i)
102  result(i) += fe.JxW(k) * factor * input[k] * fe.shape_value(i,k);
103  }
104 
117  template <int dim, typename number>
118  void L2 (
119  Vector<number> &result,
120  const FEValuesBase<dim> &fe,
121  const VectorSlice<const std::vector<std::vector<double> > > &input,
122  const double factor = 1.)
123  {
124  const unsigned int n_dofs = fe.dofs_per_cell;
125  const unsigned int fe_components = fe.get_fe().n_components();
126  const unsigned int n_components = input.size();
127 
128  AssertDimension(result.size(), n_dofs);
129  AssertDimension(input.size(), fe_components);
130 
131  for (unsigned int k=0; k<fe.n_quadrature_points; ++k)
132  for (unsigned int i=0; i<n_dofs; ++i)
133  for (unsigned int d=0; d<n_components; ++d)
134  result(i) += fe.JxW(k) * factor * fe.shape_value_component(i,k,d) * input[d][k];
135  }
136 
153  template <int dim>
154  void jump_matrix (
155  FullMatrix<double> &M11,
156  FullMatrix<double> &M12,
157  FullMatrix<double> &M21,
158  FullMatrix<double> &M22,
159  const FEValuesBase<dim> &fe1,
160  const FEValuesBase<dim> &fe2,
161  const double factor1 = 1.,
162  const double factor2 = 1.)
163  {
164  const unsigned int n1_dofs = fe1.dofs_per_cell;
165  const unsigned int n2_dofs = fe2.dofs_per_cell;
166  const unsigned int n_components = fe1.get_fe().n_components();
167 
168  Assert(n1_dofs == n2_dofs, ExcNotImplemented());
169  AssertDimension(n_components, fe2.get_fe().n_components());
170  AssertDimension(M11.m(), n1_dofs);
171  AssertDimension(M12.m(), n1_dofs);
172  AssertDimension(M21.m(), n2_dofs);
173  AssertDimension(M22.m(), n2_dofs);
174  AssertDimension(M11.n(), n1_dofs);
175  AssertDimension(M12.n(), n2_dofs);
176  AssertDimension(M21.n(), n1_dofs);
177  AssertDimension(M22.n(), n2_dofs);
178 
179  for (unsigned int k=0; k<fe1.n_quadrature_points; ++k)
180  {
181  const double dx = fe1.JxW(k);
182 
183  for (unsigned int i=0; i<n1_dofs; ++i)
184  for (unsigned int j=0; j<n1_dofs; ++j)
185  for (unsigned int d=0; d<n_components; ++d)
186  {
187  const double u1 = factor1*fe1.shape_value_component(j,k,d);
188  const double u2 =-factor2*fe2.shape_value_component(j,k,d);
189  const double v1 = factor1*fe1.shape_value_component(i,k,d);
190  const double v2 =-factor2*fe2.shape_value_component(i,k,d);
191 
192  M11(i,j) += dx * u1*v1;
193  M12(i,j) += dx * u2*v1;
194  M21(i,j) += dx * u1*v2;
195  M22(i,j) += dx * u2*v2;
196  }
197  }
198  }
199  }
200 }
201 
202 DEAL_II_NAMESPACE_CLOSE
203 
204 #endif
double shape_value_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:858
const unsigned int dofs_per_cell
Definition: fe_values.h:1418
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:154
Library of integrals over cells and faces.
Definition: maxwell.h:31
size_type n() const
void mass_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const double factor=1.)
Definition: l2.h:57
#define Assert(cond, exc)
Definition: exceptions.h:299
std::size_t size() const
size_type m() const
const unsigned int n_quadrature_points
Definition: fe_values.h:1411
const double & shape_value(const unsigned int function_no, const unsigned int point_no) const
::ExceptionBase & ExcNotImplemented()
const FiniteElement< dim, spacedim > & get_fe() const
void L2(Vector< number > &result, const FEValuesBase< dim > &fe, const std::vector< double > &input, const double factor=1.)
Definition: l2.h:89
double JxW(const unsigned int quadrature_point) const