]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add FETools::compute_projection_from_quadrature_points.
authorAndrew McBride <mcbride.andrew@gmail.com>
Thu, 29 Jul 2010 07:38:28 +0000 (07:38 +0000)
committerAndrew McBride <mcbride.andrew@gmail.com>
Thu, 29 Jul 2010 07:38:28 +0000 (07:38 +0000)
git-svn-id: https://svn.dealii.org/trunk@21585 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_tools.h
deal.II/deal.II/source/fe/fe_tools.cc

index cb20e06da56b1e7e78e25c65d402bc21a8a52c53..0e2f5c866a98e50855d056e003eb69410b9161ff 100644 (file)
@@ -25,6 +25,8 @@ DEAL_II_NAMESPACE_OPEN
 
 template <typename number> class FullMatrix;
 template <typename number> class Vector;
+template <int rank, int dim> class SymmetricTensor;
+template <int rank, int dim> class Tensor;
 template <int dim> class Quadrature;
 template <int dim, int spacedim> class FiniteElement;
 template <int dim, int spacedim> class DoFHandler;
@@ -666,6 +668,54 @@ class FETools
                                                        const Quadrature<dim>    &quadrature,
                                                        FullMatrix<double>       &I_q);
 
+                                        /**
+                                         * Computes the projection of tensorial
+                                         * (first-order tensor)
+                                         * data stored at the quadrature points
+                                         * @p vector_of_tensors_at_qp
+                                         * to data @p vector_of_tensors_at_nodes
+                                         * at the support points of the cell.
+                                         * The data in
+                                         * @p vector_of_tensors_at_qp
+                                         * is ordered sequentially following the
+                                         * quadrature point numbering.
+                                         * The size of
+                                         * @p vector_of_tensors_at_qp
+                                         * must correspond to the number of columns
+                                         * of @p projection_matrix.
+                                         * The size of @p vector_of_tensors_at_nodes
+                                         * must correspond to the number of rows of
+                                         * @p vector_of_tensors_at_nodes .
+                                         * The projection matrix
+                                         * @p projection_matrix desribes the
+                                         * projection of scalar data from the
+                                         * quadrature points and can be obtained
+                                         * from the
+                                         * FETools::compute_projection_from_quadrature_points_matrix
+                                         * function.
+                                         */
+        template <int dim>
+        static
+        void
+        compute_projection_from_quadrature_points(
+                        const FullMatrix<double>    &projection_matrix,
+                        const std::vector< Tensor<1, dim > >    &vector_of_tensors_at_qp,
+                        std::vector< Tensor<1, dim > >          &vector_of_tensors_at_nodes);
+
+
+
+                                        /**
+                                         * same as last function but for a
+                                         * @p SymmetricTensor .
+                                         */
+        template <int dim>
+        static
+        void
+        compute_projection_from_quadrature_points(
+                        const FullMatrix<double>    &projection_matrix,
+                        const std::vector< SymmetricTensor<2, dim > >   &vector_of_tensors_at_qp,
+                        std::vector< SymmetricTensor<2, dim > >         &vector_of_tensors_at_nodes);
+
 
 
 
index 316962beddccd88f436cb837ae7494374760cca9..f6567234fe0a3314b01600b0d63c31712cdd3606 100644 (file)
@@ -1973,6 +1973,130 @@ compute_interpolation_to_quadrature_points_matrix (const FiniteElement<dim,space
 
 
 
+template <int dim>
+void
+FETools::compute_projection_from_quadrature_points(
+                const FullMatrix<double>                &projection_matrix,
+                const std::vector< Tensor<1, dim > >    &vector_of_tensors_at_qp,
+                std::vector< Tensor<1, dim > >          &vector_of_tensors_at_nodes)
+{
+
+                // check that the number columns of the projection_matrix
+                // matches the size of the vector_of_tensors_at_qp
+    Assert(projection_matrix.n_cols() == vector_of_tensors_at_qp.size(),
+    ExcDimensionMismatch(projection_matrix.n_cols(),
+                    vector_of_tensors_at_qp.size()));
+
+                // check that the number rows of the projection_matrix
+                // matches the size of the vector_of_tensors_at_nodes
+    Assert(projection_matrix.n_rows() == vector_of_tensors_at_nodes.size(),
+    ExcDimensionMismatch(projection_matrix.n_rows(),
+                    vector_of_tensors_at_nodes.size()));
+
+                // number of support points (nodes) to project to
+    const unsigned int n_support_points = projection_matrix.n_rows();
+                // number of quadrature points to project from
+    const unsigned int n_quad_points = projection_matrix.n_cols();
+
+                 // component projected to the nodes
+    Vector<double> component_at_node(n_support_points);
+                 // component at the quadrature point
+    Vector<double> component_at_qp(n_quad_points);
+
+    for (unsigned int ii = 0; ii < dim; ++ii) {
+
+        component_at_qp = 0;
+
+                // populate the vector of components at the qps
+                // from vector_of_tensors_at_qp
+                // vector_of_tensors_at_qp data is in form:
+                //      columns:        0, 1, ...,  dim
+                //      rows:           0,1,....,  n_quad_points
+                // so extract the ii'th column of vector_of_tensors_at_qp
+        for (unsigned int q = 0; q < n_quad_points; ++q) {
+            component_at_qp(q) = vector_of_tensors_at_qp[q][ii];
+        }
+
+                // project from the qps -> nodes
+                // component_at_node = projection_matrix_u * component_at_qp
+        projection_matrix.vmult(component_at_node, component_at_qp);
+
+                // rewrite the projection of the components
+                // back into the vector of tensors
+        for (unsigned int nn =0; nn <n_support_points; ++nn) {
+            vector_of_tensors_at_nodes[nn][ii] = component_at_node(nn);
+        }
+    }
+}
+
+
+
+template <int dim>
+void
+FETools::compute_projection_from_quadrature_points(
+                const FullMatrix<double>                        &projection_matrix,
+                const std::vector< SymmetricTensor<2, dim > >   &vector_of_tensors_at_qp,
+                std::vector< SymmetricTensor<2, dim > >         &vector_of_tensors_at_nodes)
+{
+
+                // check that the number columns of the projection_matrix
+                // matches the size of the vector_of_tensors_at_qp
+    Assert(projection_matrix.n_cols() == vector_of_tensors_at_qp.size(),
+    ExcDimensionMismatch(projection_matrix.n_cols(),
+                    vector_of_tensors_at_qp.size()));
+
+                // check that the number rows of the projection_matrix
+                // matches the size of the vector_of_tensors_at_nodes
+    Assert(projection_matrix.n_rows() == vector_of_tensors_at_nodes.size(),
+    ExcDimensionMismatch(projection_matrix.n_rows(),
+                    vector_of_tensors_at_nodes.size()));
+
+                // number of support points (nodes)
+    const unsigned int n_support_points = projection_matrix.n_rows();
+                // number of quadrature points to project from
+    const unsigned int n_quad_points = projection_matrix.n_cols();
+
+                // number of unique entries in a symmetric second-order tensor
+    const unsigned int n_independent_components =
+            SymmetricTensor<2, dim >::n_independent_components;
+
+                // component projected to the nodes
+    Vector<double> component_at_node(n_support_points);
+                // component at the quadrature point
+    Vector<double> component_at_qp(n_quad_points);
+
+                // loop over the number of unique dimensions of the tensor
+    for (unsigned int ii = 0; ii < n_independent_components; ++ii) {
+
+        component_at_qp = 0;
+
+                // row-column entry of tensor corresponding the unrolled index
+        TableIndices<2>  row_column_index = SymmetricTensor< 2, dim >::unrolled_to_component_indices(ii);
+        const unsigned int row = row_column_index[0];
+        const unsigned int column = row_column_index[1];
+
+        //  populate the vector of components at the qps
+        //  from vector_of_tensors_at_qp
+        //  vector_of_tensors_at_qp is in form:
+        //      columns:       0, 1, ..., n_independent_components
+        //      rows:           0,1,....,  n_quad_points
+        //  so extract the ii'th column of vector_of_tensors_at_qp
+        for (unsigned int q = 0; q < n_quad_points; ++q) {
+            component_at_qp(q) = (vector_of_tensors_at_qp[q])[row][column];
+        }
+
+            // project from the qps -> nodes
+            // component_at_node = projection_matrix_u * component_at_qp
+        projection_matrix.vmult(component_at_node, component_at_qp);
+
+            // rewrite the projection of the components back into the vector of tensors
+        for (unsigned int nn =0; nn <n_support_points; ++nn) {
+            (vector_of_tensors_at_nodes[nn])[row][column] = component_at_node(nn);
+        }
+    }
+}
+
+
 
 template <int dim, int spacedim>
 void
@@ -2487,6 +2611,22 @@ compute_projection_from_quadrature_points_matrix (const FiniteElement<deal_II_di
                                                   const Quadrature<deal_II_dimension>    &rhs_quadrature,
                                                   FullMatrix<double>       &X);
 
+template
+void
+FETools::
+compute_projection_from_quadrature_points(
+                const FullMatrix<double>                &projection_matrix,
+                const std::vector< Tensor<1, deal_II_dimension > >    &vector_of_tensors_at_qp,
+                std::vector< Tensor<1, deal_II_dimension > >          &vector_of_tensors_at_nodes);
+
+template
+void
+FETools::compute_projection_from_quadrature_points(
+               const FullMatrix<double>                      &projection_matrix,
+               const std::vector<SymmetricTensor<2, deal_II_dimension> > &vector_of_tensors_at_qp,
+               std::vector<SymmetricTensor<2, deal_II_dimension> >       &vector_of_tensors_at_nodes);
+
+
 template
 void
 FETools::

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.