]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
transfer more operators
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 19 Oct 2010 20:25:21 +0000 (20:25 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 19 Oct 2010 20:25:21 +0000 (20:25 +0000)
git-svn-id: https://svn.dealii.org/trunk@22389 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/integrators/differential.h [new file with mode: 0644]
deal.II/deal.II/include/integrators/l2.h [new file with mode: 0644]

diff --git a/deal.II/deal.II/include/integrators/differential.h b/deal.II/deal.II/include/integrators/differential.h
new file mode 100644 (file)
index 0000000..5beec75
--- /dev/null
@@ -0,0 +1,192 @@
+//---------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2010 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//---------------------------------------------------------------------------
+#ifndef __deal2__integrators_differential_h
+#define __deal2__integrators_differential_h
+
+
+#include <base/config.h>
+#include <base/exceptions.h>
+#include <base/quadrature.h>
+#include <lac/full_matrix.h>
+#include <fe/mapping.h>
+#include <fe/fe_values.h>
+#include <numerics/mesh_worker_info.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace LocalIntegrators
+{
+/**
+ * @brief Local integrators related to first order differential operators and their traces.
+ *
+ * @author Guido Kanschat
+ * @date 2010
+ */
+  namespace Differential
+  {
+/**
+ * Cell matrix for divergence. The derivative is on the trial function.
+ *
+ * \f[
+ * \int_Z v\nabla \cdot u \,dx
+ * \f]
+ */
+    template <int dim>
+    void divergence_matrix (
+      FullMatrix<double>& M,
+      const FEValuesBase<dim>& fe,
+      const FEValuesBase<dim>& fetest,
+      double factor = 1.)
+    {
+      unsigned int fecomp = fe.get_fe().n_components();
+      const unsigned int n_dofs = fe.dofs_per_cell;
+      const unsigned int t_dofs = fetest.dofs_per_cell;
+      AssertDimension(fecomp, dim);
+      AssertDimension(M.m(), t_dofs);
+      AssertDimension(M.n(), n_dofs);
+       
+      for (unsigned k=0;k<fe.n_quadrature_points;++k)
+       {
+         const double dx = fe.JxW(k) * factor;
+         for (unsigned i=0;i<t_dofs;++i)
+           {
+             const double vv = fetest.shape_value(i,k);
+             for (unsigned int d=0;d<dim;++d)
+               for (unsigned j=0;j<n_dofs;++j)
+                 {
+                   const double du = fe.shape_grad_component(j,k,d)[d];
+                   M(i,j) += dx * du * vv;
+                 }
+           }
+       }
+    }
+/**
+ * Cell matrix for divergence. The derivative is on the test function.
+ *
+ * \f[
+ * \int_Z \nabla u \cdot v\,dx
+ * \f]
+ */
+    template <int dim>
+    void gradient_matrix(
+      FullMatrix<double>& M,
+      const FEValuesBase<dim>& fe,
+      const FEValuesBase<dim>& fetest,
+      double factor = 1.)
+    {
+      unsigned int fecomp = fetest.get_fe().n_components();
+      const unsigned int t_dofs = fetest.dofs_per_cell;
+      const unsigned int n_dofs = fe.dofs_per_cell;
+
+      AssertDimension(fecomp, dim);
+      AssertDimension(fe.get_fe().n_components(), 1);
+      AssertDimension(M.m(), t_dofs);
+      AssertDimension(M.n(), n_dofs);
+       
+      for (unsigned k=0;k<fe.n_quadrature_points;++k)
+       {
+         const double dx = fe.JxW(k) * factor;
+         for (unsigned int d=0;d<dim;++d)
+           for (unsigned i=0;i<t_dofs;++i)
+             {
+               const double vv = fetest.shape_value_component(i,k,d);
+               for (unsigned j=0;j<n_dofs;++j)
+                 {
+                   const Tensor<1,dim>& Du = fe.shape_grad(j,k);
+                   M(i,j) += dx * vv * Du[d];
+                 }
+             }
+       }
+    }
+
+/**
+ * The trace of the divergence
+ * operator, namely the product
+ * of the normal component of the
+ * vector valued trial space and
+ * the test space.
+ * @f[
+ * \int_F (\mathbf u\cdot \mathbf n) v \,ds
+ * @f]
+ */
+    template<int dim>
+    void
+    u_dot_n_matrix (
+      FullMatrix<double>& M,
+      const FEValuesBase<dim>& fe,
+      const FEValuesBase<dim>& fetest,
+      double factor = 1.)
+    {
+      unsigned int fecomp = fe.get_fe().n_components();
+      const unsigned int n_dofs = fe.dofs_per_cell;
+      const unsigned int t_dofs = fetest.dofs_per_cell;
+       
+      AssertDimension(fecomp, dim);
+      AssertDimension(fetest.get_fe().n_components(), 1);
+      AssertDimension(M.m(), t_dofs);
+      AssertDimension(M.n(), n_dofs);
+       
+      for (unsigned k=0;k<fe.n_quadrature_points;++k)
+       {
+         const Tensor<1,dim> ndx = factor * fe.JxW(k) * fe.normal_vector(k);
+         for (unsigned i=0;i<t_dofs;++i)
+           for (unsigned j=0;j<n_dofs;++j)
+             for (unsigned int d=0;d<dim;++d)
+               M(i,j) += ndx[d] * fe.shape_value_component(j,k,d)
+                         * fetest.shape_value(i,k);
+       }
+    }
+
+/**
+ * The trace of the divergence
+ * operator, namely the product
+ * of the normal component of the
+ * vector valued trial space and
+ * the test space.
+ * @f[
+ * \int_F (\mathbf f\cdot \mathbf n) v \,ds
+ * @f]
+ */
+    template<int dim>
+    void
+    u_dot_n_residual (
+      Vector<double>& result,
+      const FEValuesBase<dim>& fe,
+      const FEValuesBase<dim>& fetest,
+      const VectorSlice<const std::vector<std::vector<double> > >& data,
+      double factor = 1.)
+    {
+      unsigned int fecomp = fe.get_fe().n_components();
+      const unsigned int t_dofs = fetest.dofs_per_cell;
+       
+      AssertDimension(fecomp, dim);
+      AssertDimension(fetest.get_fe().n_components(), 1);
+      AssertDimension(result.size(), t_dofs);  
+      AssertVectorVectorDimension (data, dim, fe.n_quadrature_points); 
+       
+      for (unsigned k=0;k<fe.n_quadrature_points;++k)
+       {
+         const Tensor<1,dim> ndx = factor * fe.normal_vector(k) * fe.JxW(k);      
+           
+         for (unsigned i=0;i<t_dofs;++i)
+           for (unsigned int d=0;d<dim;++d)
+             result(i) += ndx[d] * fetest.shape_value(i,k) * data[d][k];
+       }
+    }
+  }
+}
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
diff --git a/deal.II/deal.II/include/integrators/l2.h b/deal.II/deal.II/include/integrators/l2.h
new file mode 100644 (file)
index 0000000..eb5e9f6
--- /dev/null
@@ -0,0 +1,122 @@
+//---------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2010 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//---------------------------------------------------------------------------
+#ifndef __deal2__integrators_l2_h
+#define __deal2__integrators_l2_h
+
+
+#include <base/config.h>
+#include <base/exceptions.h>
+#include <base/quadrature.h>
+#include <lac/full_matrix.h>
+#include <fe/mapping.h>
+#include <fe/fe_values.h>
+#include <numerics/mesh_worker_info.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace LocalIntegrators
+{
+/**
+ * @brief Local integrators related to <i>L<sup>2</sup</i>-inner products.
+ *
+ * @author Guido Kanschat
+ * @date 2010
+ */
+  namespace L2
+  {
+/**
+ * The mass matrix.
+ *
+ * \f[
+ * (a u,v)
+ * \f]
+ */
+    template <int dim>
+    void mass_matrix (
+      FullMatrix<double>& M,
+      const FEValuesBase<dim>& fe,
+      const double factor = 1.)
+    {
+      const unsigned int n_dofs = fe.dofs_per_cell;
+      const unsigned int n_components = fe.get_fe().n_components();
+      
+      for (unsigned k=0;k<fe.n_quadrature_points;++k)
+       {
+         const double dx = fe.JxW(k) * factor;
+         
+         for (unsigned i=0;i<n_dofs;++i)
+           for (unsigned j=0;j<n_dofs;++j)
+             for (unsigned int d=0;d<n_components;++d)
+               M(i,j) += dx
+                         * fe.shape_value_component(j,k,d)
+                         * fe.shape_value_component(i,k,d);
+       }
+    }
+
+/**
+ * <i>L<sup>2</sup</i>-inner product for scalar functions.
+ *
+ * \f[
+ * (f,v)
+ * \f]
+ */
+    template <int dim>
+    void L2 (
+      Vector<double>& result,
+      const FEValuesBase<dim>& fe,
+      const std::vector<double>& input,
+      const double factor = 1.)
+    {
+      const unsigned int n_dofs = fe.dofs_per_cell;
+      AssertDimension(result.size(), n_dofs);
+      AssertDimension(fe.get_fe().n_components(), 1);
+      AssertDimension(input.size(), fe.n_quadrature_points);
+      
+      for (unsigned k=0;k<fe.n_quadrature_points;++k)
+       for (unsigned i=0;i<n_dofs;++i)
+         result(i) += fe.JxW(k) * factor * input[k] * fe.shape_value(i,k);
+    }
+
+/**
+ * <i>L<sup>2</sup</i>-inner product for a slice of a vector valued
+ * right hand side.
+ *
+ * \f[
+ * \int_Z f\cdot v\,dx
+ * \f]
+ */
+    template <int dim>
+    void L2 (
+      Vector<double>& result,
+      const FEValuesBase<dim>& fe,
+      const VectorSlice<const std::vector<std::vector<double> > >& input,
+      const double factor = 1.)
+    {
+      const unsigned int n_dofs = fe.dofs_per_cell;
+      const unsigned int fe_components = fe.get_fe().n_components();
+      const unsigned int n_components = input.size();
+       
+      AssertDimension(result.size(), n_dofs);
+      AssertDimension(input.size(), fe_components);
+       
+      for (unsigned k=0;k<fe.n_quadrature_points;++k)
+       for (unsigned i=0;i<n_dofs;++i)
+         for (unsigned int d=0;d<n_components;++d)
+           result(i) += fe.JxW(k) * factor * fe.shape_value_component(i,k,d) * input[d][k];
+    } 
+  }
+}
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.