From: Andrew McBride Date: Thu, 29 Jul 2010 07:38:28 +0000 (+0000) Subject: Add FETools::compute_projection_from_quadrature_points. X-Git-Tag: v8.0.0~5759 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=fa18a6ecdbd3a4b506a2505f52e5dbbf3994feb6;p=dealii.git Add FETools::compute_projection_from_quadrature_points. git-svn-id: https://svn.dealii.org/trunk@21585 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe_tools.h b/deal.II/deal.II/include/fe/fe_tools.h index cb20e06da5..0e2f5c866a 100644 --- a/deal.II/deal.II/include/fe/fe_tools.h +++ b/deal.II/deal.II/include/fe/fe_tools.h @@ -25,6 +25,8 @@ DEAL_II_NAMESPACE_OPEN template class FullMatrix; template class Vector; +template class SymmetricTensor; +template class Tensor; template class Quadrature; template class FiniteElement; template class DoFHandler; @@ -666,6 +668,54 @@ class FETools const Quadrature &quadrature, FullMatrix &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 + static + void + compute_projection_from_quadrature_points( + const FullMatrix &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 + static + void + compute_projection_from_quadrature_points( + const FullMatrix &projection_matrix, + const std::vector< SymmetricTensor<2, dim > > &vector_of_tensors_at_qp, + std::vector< SymmetricTensor<2, dim > > &vector_of_tensors_at_nodes); + diff --git a/deal.II/deal.II/source/fe/fe_tools.cc b/deal.II/deal.II/source/fe/fe_tools.cc index 316962bedd..f6567234fe 100644 --- a/deal.II/deal.II/source/fe/fe_tools.cc +++ b/deal.II/deal.II/source/fe/fe_tools.cc @@ -1973,6 +1973,130 @@ compute_interpolation_to_quadrature_points_matrix (const FiniteElement +void +FETools::compute_projection_from_quadrature_points( + const FullMatrix &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 component_at_node(n_support_points); + // component at the quadrature point + Vector 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 +void +FETools::compute_projection_from_quadrature_points( + const FullMatrix &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 component_at_node(n_support_points); + // component at the quadrature point + Vector 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 void @@ -2487,6 +2611,22 @@ compute_projection_from_quadrature_points_matrix (const FiniteElement &rhs_quadrature, FullMatrix &X); +template +void +FETools:: +compute_projection_from_quadrature_points( + const FullMatrix &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 &projection_matrix, + const std::vector > &vector_of_tensors_at_qp, + std::vector > &vector_of_tensors_at_nodes); + + template void FETools::