From 821bbb29534081b589c79ca14a240359c85880f6 Mon Sep 17 00:00:00 2001 From: Ralf Hartmann Date: Mon, 27 Aug 2001 13:04:59 +0000 Subject: [PATCH] Extend the create_mass_matrix and create_laplace_matrix functions to use vector valued coefficient functions for system elements. git-svn-id: https://svn.dealii.org/trunk@4919 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/source/numerics/matrices.cc | 233 ++++++++++++++------ 1 file changed, 170 insertions(+), 63 deletions(-) diff --git a/deal.II/deal.II/source/numerics/matrices.cc b/deal.II/deal.II/source/numerics/matrices.cc index f8a70f1f27..163515d7d6 100644 --- a/deal.II/deal.II/source/numerics/matrices.cc +++ b/deal.II/deal.II/source/numerics/matrices.cc @@ -97,7 +97,6 @@ void MatrixCreator::create_mass_matrix (const Mapping &mapping, -// TODO:[RH] extend this function to use vector valued coefficient functions for system elements. template void MatrixCreator::create_mass_matrix_1 (const Mapping &mapping, const DoFHandler &dof, @@ -110,7 +109,7 @@ void MatrixCreator::create_mass_matrix_1 (const Mapping &mapping, UpdateFlags update_flags = UpdateFlags(update_values | update_JxW_values); if (coefficient != 0) update_flags = UpdateFlags (update_flags | update_q_points); - + FEValues fe_values (mapping, dof.get_fe(), q, update_flags); const unsigned int dofs_per_cell = fe_values.dofs_per_cell, @@ -118,8 +117,13 @@ void MatrixCreator::create_mass_matrix_1 (const Mapping &mapping, const FiniteElement &fe = fe_values.get_fe(); const unsigned int n_components = fe.n_components(); + Assert(coefficient->n_components==1 || + coefficient->n_components==n_components, ExcComponentMismatch()); + FullMatrix cell_matrix (dofs_per_cell, dofs_per_cell); std::vector coefficient_values (n_q_points); + std::vector > coefficient_vector_values (n_q_points, + Vector (n_components)); std::vector dof_indices (dofs_per_cell); @@ -136,18 +140,39 @@ void MatrixCreator::create_mass_matrix_1 (const Mapping &mapping, if (coefficient != 0) { - coefficient->value_list (fe_values.get_quadrature_points(), - coefficient_values); - for (unsigned int point=0; pointn_components==1) + { + coefficient->value_list (fe_values.get_quadrature_points(), + coefficient_values); + for (unsigned int point=0; pointvector_value_list (fe_values.get_quadrature_points(), + coefficient_vector_values); + for (unsigned int point=0; point &mapping, -// TODO:[RH] extend this function to use vector valued coefficient functions for system elements. template void MatrixCreator::create_mass_matrix_2 (const Mapping &mapping, @@ -249,10 +273,15 @@ MatrixCreator::create_mass_matrix_2 (const Mapping &mapping, const FiniteElement &fe = fe_values.get_fe(); const unsigned int n_components = fe.n_components(); + Assert(coefficient->n_components==1 || + coefficient->n_components==n_components, ExcComponentMismatch()); + FullMatrix cell_matrix (dofs_per_cell, dofs_per_cell); Vector local_rhs (dofs_per_cell); std::vector rhs_values (fe_values.n_quadrature_points); std::vector coefficient_values (n_q_points); + std::vector > coefficient_vector_values (n_q_points, + Vector (n_components)); std::vector dof_indices (dofs_per_cell); @@ -271,23 +300,47 @@ MatrixCreator::create_mass_matrix_2 (const Mapping &mapping, if (coefficient != 0) { - coefficient->value_list (fe_values.get_quadrature_points(), - coefficient_values); - for (unsigned int point=0; pointn_components==1) + { + coefficient->value_list (fe_values.get_quadrature_points(), + coefficient_values); + for (unsigned int point=0; pointvector_value_list (fe_values.get_quadrature_points(), + coefficient_vector_values); + for (unsigned int point=0; point &mapping, -// TODO:[RH] extend this function to use vector valued coefficient functions for system elements. template void MatrixCreator::create_laplace_matrix_1 (const Mapping &mapping, const DoFHandler &dof, @@ -791,8 +843,13 @@ void MatrixCreator::create_laplace_matrix_1 (const Mapping &mapping, const FiniteElement &fe = fe_values.get_fe(); const unsigned int n_components = fe.n_components(); + Assert(coefficient->n_components==1 || + coefficient->n_components==n_components, ExcComponentMismatch()); + FullMatrix cell_matrix (dofs_per_cell, dofs_per_cell); std::vector coefficient_values (n_q_points); + std::vector > coefficient_vector_values (n_q_points, + Vector (n_components)); std::vector dof_indices (dofs_per_cell); @@ -810,18 +867,40 @@ void MatrixCreator::create_laplace_matrix_1 (const Mapping &mapping, if (coefficient != 0) { - coefficient->value_list (fe_values.get_quadrature_points(), - coefficient_values); - for (unsigned int point=0; pointn_components==1) + { + coefficient->value_list (fe_values.get_quadrature_points(), + coefficient_values); + for (unsigned int point=0; pointvector_value_list (fe_values.get_quadrature_points(), + coefficient_vector_values); + for (unsigned int point=0; point &mapping, -// TODO:[RH] extend this function to use vector valued coefficient functions for system elements. template void MatrixCreator::create_laplace_matrix_2 (const Mapping &mapping, @@ -962,10 +1040,15 @@ MatrixCreator::create_laplace_matrix_2 (const Mapping &mapping, const FiniteElement &fe = fe_values.get_fe(); const unsigned int n_components = fe.n_components(); + Assert(coefficient->n_components==1 || + coefficient->n_components==n_components, ExcComponentMismatch()); + FullMatrix cell_matrix (dofs_per_cell, dofs_per_cell); Vector local_rhs (dofs_per_cell); std::vector rhs_values (fe_values.n_quadrature_points); std::vector coefficient_values (n_q_points); + std::vector > coefficient_vector_values (n_q_points, + Vector (n_components)); std::vector dof_indices (dofs_per_cell); @@ -986,23 +1069,47 @@ MatrixCreator::create_laplace_matrix_2 (const Mapping &mapping, if (coefficient != 0) { - coefficient->value_list (fe_values.get_quadrature_points(), - coefficient_values); - for (unsigned int point=0; pointn_components==1) + { + coefficient->value_list (fe_values.get_quadrature_points(), + coefficient_values); + for (unsigned int point=0; pointvector_value_list (fe_values.get_quadrature_points(), + coefficient_vector_values); + for (unsigned int point=0; point