From: Peter Munch Date: Fri, 22 Oct 2021 15:46:40 +0000 (+0200) Subject: MatrixCreator: make SparseMatrixType a template parameter X-Git-Tag: v9.4.0-rc1~880^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=0c38c71f43bb413c869831f7419f99b73a9fe888;p=dealii.git MatrixCreator: make SparseMatrixType a template parameter --- diff --git a/include/deal.II/numerics/matrix_creator.templates.h b/include/deal.II/numerics/matrix_creator.templates.h index ff9f2fac5c..b1a9a8f54f 100644 --- a/include/deal.II/numerics/matrix_creator.templates.h +++ b/include/deal.II/numerics/matrix_creator.templates.h @@ -150,6 +150,7 @@ namespace MatrixCreator FullMatrix cell_matrix; dealii::Vector cell_rhs; const AffineConstraints * constraints; + bool is_locally_owned; }; } // namespace AssemblerData @@ -162,6 +163,10 @@ namespace MatrixCreator & data, MatrixCreator::internal::AssemblerData::CopyData ©_data) { + copy_data.is_locally_owned = cell->is_locally_owned(); + if (copy_data.is_locally_owned == false) + return; + data.x_fe_values.reinit(cell); const FEValues &fe_values = data.x_fe_values.get_present_fe_values(); @@ -361,14 +366,18 @@ namespace MatrixCreator - template + template void laplace_assembler( const CellIterator &cell, - MatrixCreator::internal::AssemblerData::Scratch + MatrixCreator::internal::AssemblerData::Scratch & data, - MatrixCreator::internal::AssemblerData::CopyData ©_data) + MatrixCreator::internal::AssemblerData::CopyData ©_data) { + copy_data.is_locally_owned = cell->is_locally_owned(); + if (copy_data.is_locally_owned == false) + return; + data.x_fe_values.reinit(cell); const FEValues &fe_values = data.x_fe_values.get_present_fe_values(); @@ -405,7 +414,7 @@ namespace MatrixCreator else { data.rhs_vector_values.resize( - n_q_points, dealii::Vector(n_components)); + n_q_points, dealii::Vector(n_components)); data.rhs_function->vector_value_list( fe_values.get_quadrature_points(), data.rhs_vector_values); } @@ -423,7 +432,7 @@ namespace MatrixCreator else { data.coefficient_vector_values.resize( - n_q_points, dealii::Vector(n_components)); + n_q_points, dealii::Vector(n_components)); data.coefficient->vector_value_list( fe_values.get_quadrature_points(), data.coefficient_vector_values); @@ -432,7 +441,7 @@ namespace MatrixCreator const std::vector &JxW = fe_values.get_JxW_values(); - double add_data; + Number add_data; for (unsigned int i = 0; i < dofs_per_cell; ++i) if (fe.is_primitive()) { @@ -569,6 +578,9 @@ namespace MatrixCreator MatrixType * matrix, VectorType * right_hand_side) { + if (data.is_locally_owned == false) + return; + const unsigned int dofs_per_cell = data.dof_indices.size(); (void)dofs_per_cell; @@ -626,14 +638,16 @@ namespace MatrixCreator namespace MatrixCreator { - template + template void - create_mass_matrix(const Mapping & mapping, - const DoFHandler & dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function *const coefficient, - const AffineConstraints & constraints) + create_mass_matrix( + const Mapping & mapping, + const DoFHandler &dof, + const Quadrature & q, + SparseMatrixType & matrix, + const Function + *const coefficient, + const AffineConstraints &constraints) { Assert(matrix.m() == dof.n_dofs(), ExcDimensionMismatch(matrix.m(), dof.n_dofs())); @@ -643,17 +657,20 @@ namespace MatrixCreator const auto & fe_collection = dof.get_fe_collection(); hp::QCollection q_collection(q); hp::MappingCollection mapping_collection(mapping); - MatrixCreator::internal::AssemblerData::Scratch - assembler_data(fe_collection, - update_values | update_JxW_values | - (coefficient != nullptr ? update_quadrature_points : - UpdateFlags(0)), - coefficient, - /*rhs_function=*/nullptr, - q_collection, - mapping_collection); - - MatrixCreator::internal::AssemblerData::CopyData copy_data; + MatrixCreator::internal::AssemblerData:: + Scratch + assembler_data(fe_collection, + update_values | update_JxW_values | + (coefficient != nullptr ? update_quadrature_points : + UpdateFlags(0)), + coefficient, + /*rhs_function=*/nullptr, + q_collection, + mapping_collection); + + MatrixCreator::internal::AssemblerData::CopyData< + typename SparseMatrixType::value_type> + copy_data; copy_data.cell_matrix.reinit( assembler_data.fe_collection.max_dofs_per_cell(), assembler_data.fe_collection.max_dofs_per_cell()); @@ -670,24 +687,32 @@ namespace MatrixCreator dim, spacedim, typename DoFHandler::active_cell_iterator, - number>, - [&matrix](const internal::AssemblerData::CopyData &data) { + typename SparseMatrixType::value_type>, + [&matrix](const internal::AssemblerData::CopyData< + typename SparseMatrixType::value_type> &data) { MatrixCreator::internal::copy_local_to_global( - data, &matrix, static_cast *>(nullptr)); + data, + &matrix, + static_cast *>( + nullptr)); }, assembler_data, copy_data); + + matrix.compress(VectorOperation::add); } - template + template void - create_mass_matrix(const DoFHandler & dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function *const coefficient, - const AffineConstraints & constraints) + create_mass_matrix( + const DoFHandler &dof, + const Quadrature & q, + SparseMatrixType & matrix, + const Function + *const coefficient, + const AffineConstraints &constraints) { create_mass_matrix(get_default_linear_mapping(dof.get_triangulation()), dof, @@ -699,16 +724,18 @@ namespace MatrixCreator - template + template void - create_mass_matrix(const Mapping & mapping, - const DoFHandler & dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const coefficient, - const AffineConstraints & constraints) + create_mass_matrix( + const Mapping & mapping, + const DoFHandler & dof, + const Quadrature & q, + SparseMatrixType & matrix, + const Function &rhs, + Vector & rhs_vector, + const Function + *const coefficient, + const AffineConstraints &constraints) { Assert(matrix.m() == dof.n_dofs(), ExcDimensionMismatch(matrix.m(), dof.n_dofs())); @@ -718,15 +745,18 @@ namespace MatrixCreator const auto & fe_collection = dof.get_fe_collection(); hp::QCollection q_collection(q); hp::MappingCollection mapping_collection(mapping); - MatrixCreator::internal::AssemblerData::Scratch - assembler_data(fe_collection, - update_values | update_JxW_values | - update_quadrature_points, - coefficient, - &rhs, - q_collection, - mapping_collection); - MatrixCreator::internal::AssemblerData::CopyData copy_data; + MatrixCreator::internal::AssemblerData:: + Scratch + assembler_data(fe_collection, + update_values | update_JxW_values | + update_quadrature_points, + coefficient, + &rhs, + q_collection, + mapping_collection); + MatrixCreator::internal::AssemblerData::CopyData< + typename SparseMatrixType::value_type> + copy_data; copy_data.cell_matrix.reinit( assembler_data.fe_collection.max_dofs_per_cell(), assembler_data.fe_collection.max_dofs_per_cell()); @@ -743,28 +773,32 @@ namespace MatrixCreator dim, spacedim, typename DoFHandler::active_cell_iterator, - number>, - [&matrix, - &rhs_vector](const internal::AssemblerData::CopyData &data) { + typename SparseMatrixType::value_type>, + [&matrix, &rhs_vector](const internal::AssemblerData::CopyData< + typename SparseMatrixType::value_type> &data) { MatrixCreator::internal::copy_local_to_global(data, &matrix, &rhs_vector); }, assembler_data, copy_data); + + matrix.compress(VectorOperation::add); } - template + template void - create_mass_matrix(const DoFHandler & dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const coefficient, - const AffineConstraints & constraints) + create_mass_matrix( + const DoFHandler & dof, + const Quadrature & q, + SparseMatrixType & matrix, + const Function &rhs, + Vector & rhs_vector, + const Function + *const coefficient, + const AffineConstraints &constraints) { create_mass_matrix(get_default_linear_mapping(dof.get_triangulation()), dof, @@ -778,30 +812,35 @@ namespace MatrixCreator - template + template void - create_mass_matrix(const hp::MappingCollection &mapping, - const DoFHandler & dof, - const hp::QCollection & q, - SparseMatrix & matrix, - const Function *const coefficient, - const AffineConstraints & constraints) + create_mass_matrix( + const hp::MappingCollection &mapping, + const DoFHandler & dof, + const hp::QCollection & q, + SparseMatrixType & matrix, + const Function + *const coefficient, + const AffineConstraints &constraints) { Assert(matrix.m() == dof.n_dofs(), ExcDimensionMismatch(matrix.m(), dof.n_dofs())); Assert(matrix.n() == dof.n_dofs(), ExcDimensionMismatch(matrix.n(), dof.n_dofs())); - MatrixCreator::internal::AssemblerData::Scratch - assembler_data(dof.get_fe_collection(), - update_values | update_JxW_values | - (coefficient != nullptr ? update_quadrature_points : - UpdateFlags(0)), - coefficient, - /*rhs_function=*/nullptr, - q, - mapping); - MatrixCreator::internal::AssemblerData::CopyData copy_data; + MatrixCreator::internal::AssemblerData:: + Scratch + assembler_data(dof.get_fe_collection(), + update_values | update_JxW_values | + (coefficient != nullptr ? update_quadrature_points : + UpdateFlags(0)), + coefficient, + /*rhs_function=*/nullptr, + q, + mapping); + MatrixCreator::internal::AssemblerData::CopyData< + typename SparseMatrixType::value_type> + copy_data; copy_data.cell_matrix.reinit( assembler_data.fe_collection.max_dofs_per_cell(), assembler_data.fe_collection.max_dofs_per_cell()); @@ -818,24 +857,32 @@ namespace MatrixCreator dim, spacedim, typename DoFHandler::active_cell_iterator, - number>, - [&matrix](const internal::AssemblerData::CopyData &data) { + typename SparseMatrixType::value_type>, + [&matrix](const internal::AssemblerData::CopyData< + typename SparseMatrixType::value_type> &data) { MatrixCreator::internal::copy_local_to_global( - data, &matrix, static_cast *>(nullptr)); + data, + &matrix, + static_cast *>( + nullptr)); }, assembler_data, copy_data); + + matrix.compress(VectorOperation::add); } - template + template void - create_mass_matrix(const DoFHandler & dof, - const hp::QCollection & q, - SparseMatrix & matrix, - const Function *const coefficient, - const AffineConstraints & constraints) + create_mass_matrix( + const DoFHandler &dof, + const hp::QCollection & q, + SparseMatrixType & matrix, + const Function + *const coefficient, + const AffineConstraints &constraints) { create_mass_matrix(hp::StaticMappingQ1::mapping_collection, dof, @@ -847,31 +894,36 @@ namespace MatrixCreator - template + template void - create_mass_matrix(const hp::MappingCollection &mapping, - const DoFHandler & dof, - const hp::QCollection & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const coefficient, - const AffineConstraints & constraints) + create_mass_matrix( + const hp::MappingCollection & mapping, + const DoFHandler & dof, + const hp::QCollection & q, + SparseMatrixType & matrix, + const Function &rhs, + Vector & rhs_vector, + const Function + *const coefficient, + const AffineConstraints &constraints) { Assert(matrix.m() == dof.n_dofs(), ExcDimensionMismatch(matrix.m(), dof.n_dofs())); Assert(matrix.n() == dof.n_dofs(), ExcDimensionMismatch(matrix.n(), dof.n_dofs())); - MatrixCreator::internal::AssemblerData::Scratch - assembler_data(dof.get_fe_collection(), - update_values | update_JxW_values | - update_quadrature_points, - coefficient, - &rhs, - q, - mapping); - MatrixCreator::internal::AssemblerData::CopyData copy_data; + MatrixCreator::internal::AssemblerData:: + Scratch + assembler_data(dof.get_fe_collection(), + update_values | update_JxW_values | + update_quadrature_points, + coefficient, + &rhs, + q, + mapping); + MatrixCreator::internal::AssemblerData::CopyData< + typename SparseMatrixType::value_type> + copy_data; copy_data.cell_matrix.reinit( assembler_data.fe_collection.max_dofs_per_cell(), assembler_data.fe_collection.max_dofs_per_cell()); @@ -888,28 +940,32 @@ namespace MatrixCreator dim, spacedim, typename DoFHandler::active_cell_iterator, - number>, - [&matrix, - &rhs_vector](const internal::AssemblerData::CopyData &data) { + typename SparseMatrixType::value_type>, + [&matrix, &rhs_vector](const internal::AssemblerData::CopyData< + typename SparseMatrixType::value_type> &data) { MatrixCreator::internal::copy_local_to_global(data, &matrix, &rhs_vector); }, assembler_data, copy_data); + + matrix.compress(VectorOperation::add); } - template + template void - create_mass_matrix(const DoFHandler & dof, - const hp::QCollection & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const coefficient, - const AffineConstraints & constraints) + create_mass_matrix( + const DoFHandler & dof, + const hp::QCollection & q, + SparseMatrixType & matrix, + const Function &rhs, + Vector & rhs_vector, + const Function + *const coefficient, + const AffineConstraints &constraints) { create_mass_matrix(hp::StaticMappingQ1::mapping_collection, dof, @@ -1857,14 +1913,16 @@ namespace MatrixCreator - template + template void - create_laplace_matrix(const Mapping & mapping, - const DoFHandler &dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function *const coefficient, - const AffineConstraints &constraints) + create_laplace_matrix( + const Mapping & mapping, + const DoFHandler &dof, + const Quadrature & q, + SparseMatrixType & matrix, + const Function + *const coefficient, + const AffineConstraints &constraints) { Assert(matrix.m() == dof.n_dofs(), ExcDimensionMismatch(matrix.m(), dof.n_dofs())); @@ -1874,16 +1932,19 @@ namespace MatrixCreator const auto & fe_collection = dof.get_fe_collection(); hp::QCollection q_collection(q); hp::MappingCollection mapping_collection(mapping); - MatrixCreator::internal::AssemblerData::Scratch - assembler_data(fe_collection, - update_gradients | update_JxW_values | - (coefficient != nullptr ? update_quadrature_points : - UpdateFlags(0)), - coefficient, - /*rhs_function=*/nullptr, - q_collection, - mapping_collection); - MatrixCreator::internal::AssemblerData::CopyData copy_data; + MatrixCreator::internal::AssemblerData:: + Scratch + assembler_data(fe_collection, + update_gradients | update_JxW_values | + (coefficient != nullptr ? update_quadrature_points : + UpdateFlags(0)), + coefficient, + /*rhs_function=*/nullptr, + q_collection, + mapping_collection); + MatrixCreator::internal::AssemblerData::CopyData< + typename SparseMatrixType::value_type> + copy_data; copy_data.cell_matrix.reinit( assembler_data.fe_collection.max_dofs_per_cell(), assembler_data.fe_collection.max_dofs_per_cell()); @@ -1899,24 +1960,33 @@ namespace MatrixCreator &MatrixCreator::internal::laplace_assembler< dim, spacedim, - typename DoFHandler::active_cell_iterator>, - [&matrix](const internal::AssemblerData::CopyData &data) { + typename DoFHandler::active_cell_iterator, + typename SparseMatrixType::value_type>, + [&matrix](const internal::AssemblerData::CopyData< + typename SparseMatrixType::value_type> &data) { MatrixCreator::internal::copy_local_to_global( - data, &matrix, static_cast *>(nullptr)); + data, + &matrix, + static_cast *>( + nullptr)); }, assembler_data, copy_data); + + matrix.compress(VectorOperation::add); } - template + template void - create_laplace_matrix(const DoFHandler &dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function *const coefficient, - const AffineConstraints &constraints) + create_laplace_matrix( + const DoFHandler &dof, + const Quadrature & q, + SparseMatrixType & matrix, + const Function + *const coefficient, + const AffineConstraints &constraints) { create_laplace_matrix(get_default_linear_mapping(dof.get_triangulation()), dof, @@ -1928,16 +1998,18 @@ namespace MatrixCreator - template + template void - create_laplace_matrix(const Mapping & mapping, - const DoFHandler &dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const coefficient, - const AffineConstraints &constraints) + create_laplace_matrix( + const Mapping & mapping, + const DoFHandler & dof, + const Quadrature & q, + SparseMatrixType & matrix, + const Function &rhs, + Vector & rhs_vector, + const Function + *const coefficient, + const AffineConstraints &constraints) { Assert(matrix.m() == dof.n_dofs(), ExcDimensionMismatch(matrix.m(), dof.n_dofs())); @@ -1947,15 +2019,18 @@ namespace MatrixCreator const auto & fe_collection = dof.get_fe_collection(); hp::QCollection q_collection(q); hp::MappingCollection mapping_collection(mapping); - MatrixCreator::internal::AssemblerData::Scratch - assembler_data(fe_collection, - update_gradients | update_values | update_JxW_values | - update_quadrature_points, - coefficient, - &rhs, - q_collection, - mapping_collection); - MatrixCreator::internal::AssemblerData::CopyData copy_data; + MatrixCreator::internal::AssemblerData:: + Scratch + assembler_data(fe_collection, + update_gradients | update_values | update_JxW_values | + update_quadrature_points, + coefficient, + &rhs, + q_collection, + mapping_collection); + MatrixCreator::internal::AssemblerData::CopyData< + typename SparseMatrixType::value_type> + copy_data; copy_data.cell_matrix.reinit( assembler_data.fe_collection.max_dofs_per_cell(), assembler_data.fe_collection.max_dofs_per_cell()); @@ -1971,28 +2046,33 @@ namespace MatrixCreator &MatrixCreator::internal::laplace_assembler< dim, spacedim, - typename DoFHandler::active_cell_iterator>, - [&matrix, - &rhs_vector](const internal::AssemblerData::CopyData &data) { + typename DoFHandler::active_cell_iterator, + typename SparseMatrixType::value_type>, + [&matrix, &rhs_vector](const internal::AssemblerData::CopyData< + typename SparseMatrixType::value_type> &data) { MatrixCreator::internal::copy_local_to_global(data, &matrix, &rhs_vector); }, assembler_data, copy_data); + + matrix.compress(VectorOperation::add); } - template + template void - create_laplace_matrix(const DoFHandler &dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const coefficient, - const AffineConstraints &constraints) + create_laplace_matrix( + const DoFHandler & dof, + const Quadrature & q, + SparseMatrixType & matrix, + const Function &rhs, + Vector & rhs_vector, + const Function + *const coefficient, + const AffineConstraints &constraints) { create_laplace_matrix(get_default_linear_mapping(dof.get_triangulation()), dof, @@ -2006,30 +2086,35 @@ namespace MatrixCreator - template + template void - create_laplace_matrix(const hp::MappingCollection &mapping, - const DoFHandler & dof, - const hp::QCollection & q, - SparseMatrix & matrix, - const Function *const coefficient, - const AffineConstraints & constraints) + create_laplace_matrix( + const hp::MappingCollection &mapping, + const DoFHandler & dof, + const hp::QCollection & q, + SparseMatrixType & matrix, + const Function + *const coefficient, + const AffineConstraints &constraints) { Assert(matrix.m() == dof.n_dofs(), ExcDimensionMismatch(matrix.m(), dof.n_dofs())); Assert(matrix.n() == dof.n_dofs(), ExcDimensionMismatch(matrix.n(), dof.n_dofs())); - MatrixCreator::internal::AssemblerData::Scratch - assembler_data(dof.get_fe_collection(), - update_gradients | update_JxW_values | - (coefficient != nullptr ? update_quadrature_points : - UpdateFlags(0)), - coefficient, - /*rhs_function=*/nullptr, - q, - mapping); - MatrixCreator::internal::AssemblerData::CopyData copy_data; + MatrixCreator::internal::AssemblerData:: + Scratch + assembler_data(dof.get_fe_collection(), + update_gradients | update_JxW_values | + (coefficient != nullptr ? update_quadrature_points : + UpdateFlags(0)), + coefficient, + /*rhs_function=*/nullptr, + q, + mapping); + MatrixCreator::internal::AssemblerData::CopyData< + typename SparseMatrixType::value_type> + copy_data; copy_data.cell_matrix.reinit( assembler_data.fe_collection.max_dofs_per_cell(), assembler_data.fe_collection.max_dofs_per_cell()); @@ -2045,24 +2130,33 @@ namespace MatrixCreator &MatrixCreator::internal::laplace_assembler< dim, spacedim, - typename DoFHandler::active_cell_iterator>, - [&matrix](const internal::AssemblerData::CopyData &data) { + typename DoFHandler::active_cell_iterator, + typename SparseMatrixType::value_type>, + [&matrix](const internal::AssemblerData::CopyData< + typename SparseMatrixType::value_type> &data) { MatrixCreator::internal::copy_local_to_global( - data, &matrix, static_cast *>(nullptr)); + data, + &matrix, + static_cast *>( + nullptr)); }, assembler_data, copy_data); + + matrix.compress(VectorOperation::add); } - template + template void - create_laplace_matrix(const DoFHandler &dof, - const hp::QCollection & q, - SparseMatrix & matrix, - const Function *const coefficient, - const AffineConstraints &constraints) + create_laplace_matrix( + const DoFHandler &dof, + const hp::QCollection & q, + SparseMatrixType & matrix, + const Function + *const coefficient, + const AffineConstraints &constraints) { create_laplace_matrix( hp::StaticMappingQ1::mapping_collection, @@ -2075,31 +2169,36 @@ namespace MatrixCreator - template + template void - create_laplace_matrix(const hp::MappingCollection &mapping, - const DoFHandler & dof, - const hp::QCollection & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const coefficient, - const AffineConstraints & constraints) + create_laplace_matrix( + const hp::MappingCollection & mapping, + const DoFHandler & dof, + const hp::QCollection & q, + SparseMatrixType & matrix, + const Function &rhs, + Vector & rhs_vector, + const Function + *const coefficient, + const AffineConstraints &constraints) { Assert(matrix.m() == dof.n_dofs(), ExcDimensionMismatch(matrix.m(), dof.n_dofs())); Assert(matrix.n() == dof.n_dofs(), ExcDimensionMismatch(matrix.n(), dof.n_dofs())); - MatrixCreator::internal::AssemblerData::Scratch - assembler_data(dof.get_fe_collection(), - update_gradients | update_values | update_JxW_values | - update_quadrature_points, - coefficient, - &rhs, - q, - mapping); - MatrixCreator::internal::AssemblerData::CopyData copy_data; + MatrixCreator::internal::AssemblerData:: + Scratch + assembler_data(dof.get_fe_collection(), + update_gradients | update_values | update_JxW_values | + update_quadrature_points, + coefficient, + &rhs, + q, + mapping); + MatrixCreator::internal::AssemblerData::CopyData< + typename SparseMatrixType::value_type> + copy_data; copy_data.cell_matrix.reinit( assembler_data.fe_collection.max_dofs_per_cell(), assembler_data.fe_collection.max_dofs_per_cell()); @@ -2115,28 +2214,33 @@ namespace MatrixCreator &MatrixCreator::internal::laplace_assembler< dim, spacedim, - typename DoFHandler::active_cell_iterator>, - [&matrix, - &rhs_vector](const internal::AssemblerData::CopyData &data) { + typename DoFHandler::active_cell_iterator, + typename SparseMatrixType::value_type>, + [&matrix, &rhs_vector](const internal::AssemblerData::CopyData< + typename SparseMatrixType::value_type> &data) { MatrixCreator::internal::copy_local_to_global(data, &matrix, &rhs_vector); }, assembler_data, copy_data); + + matrix.compress(VectorOperation::add); } - template + template void - create_laplace_matrix(const DoFHandler &dof, - const hp::QCollection & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const coefficient, - const AffineConstraints &constraints) + create_laplace_matrix( + const DoFHandler & dof, + const hp::QCollection & q, + SparseMatrixType & matrix, + const Function &rhs, + Vector & rhs_vector, + const Function + *const coefficient, + const AffineConstraints &constraints) { create_laplace_matrix( hp::StaticMappingQ1::mapping_collection, diff --git a/include/deal.II/numerics/matrix_tools.h b/include/deal.II/numerics/matrix_tools.h index 21a0c2d25c..745c5cc4ac 100644 --- a/include/deal.II/numerics/matrix_tools.h +++ b/include/deal.II/numerics/matrix_tools.h @@ -246,28 +246,32 @@ namespace MatrixCreator * * See the general documentation of this namespace for more information. */ - template + template void create_mass_matrix( - const Mapping & mapping, - const DoFHandler & dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function *const a = nullptr, - const AffineConstraints &constraints = AffineConstraints()); + const Mapping & mapping, + const DoFHandler &dof, + const Quadrature & q, + SparseMatrixType & matrix, + const Function *const a = + nullptr, + const AffineConstraints & + constraints = AffineConstraints()); /** * Call the create_mass_matrix() function, see above, with * mapping=MappingQ@(1). */ - template + template void create_mass_matrix( - const DoFHandler & dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function *const a = nullptr, - const AffineConstraints &constraints = AffineConstraints()); + const DoFHandler &dof, + const Quadrature & q, + SparseMatrixType & matrix, + const Function *const a = + nullptr, + const AffineConstraints & + constraints = AffineConstraints()); /** * Assemble the mass matrix and a right hand side vector. If no coefficient @@ -288,86 +292,98 @@ namespace MatrixCreator * * See the general documentation of this namespace for more information. */ - template + template void create_mass_matrix( - const Mapping & mapping, - const DoFHandler & dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const a = nullptr, - const AffineConstraints &constraints = AffineConstraints()); + const Mapping & mapping, + const DoFHandler & dof, + const Quadrature & q, + SparseMatrixType & matrix, + const Function &rhs, + Vector & rhs_vector, + const Function *const a = + nullptr, + const AffineConstraints & + constraints = AffineConstraints()); /** * Call the create_mass_matrix() function, see above, with * mapping=MappingQ@(1). */ - template + template void create_mass_matrix( - const DoFHandler & dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const a = nullptr, - const AffineConstraints &constraints = AffineConstraints()); + const DoFHandler & dof, + const Quadrature & q, + SparseMatrixType & matrix, + const Function &rhs, + Vector & rhs_vector, + const Function *const a = + nullptr, + const AffineConstraints & + constraints = AffineConstraints()); /** * Same function as above, but for hp-objects. */ - template + template void create_mass_matrix( const hp::MappingCollection &mapping, const DoFHandler & dof, const hp::QCollection & q, - SparseMatrix & matrix, - const Function *const a = nullptr, - const AffineConstraints &constraints = AffineConstraints()); + SparseMatrixType & matrix, + const Function *const a = + nullptr, + const AffineConstraints & + constraints = AffineConstraints()); /** * Same function as above, but for hp-objects. */ - template + template void create_mass_matrix( - const DoFHandler & dof, - const hp::QCollection & q, - SparseMatrix & matrix, - const Function *const a = nullptr, - const AffineConstraints &constraints = AffineConstraints()); + const DoFHandler &dof, + const hp::QCollection & q, + SparseMatrixType & matrix, + const Function *const a = + nullptr, + const AffineConstraints & + constraints = AffineConstraints()); /** * Same function as above, but for hp-objects. */ - template + template void create_mass_matrix( - const hp::MappingCollection &mapping, - const DoFHandler & dof, - const hp::QCollection & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const a = nullptr, - const AffineConstraints &constraints = AffineConstraints()); + const hp::MappingCollection & mapping, + const DoFHandler & dof, + const hp::QCollection & q, + SparseMatrixType & matrix, + const Function &rhs, + Vector & rhs_vector, + const Function *const a = + nullptr, + const AffineConstraints & + constraints = AffineConstraints()); /** * Same function as above, but for hp-objects. */ - template + template void create_mass_matrix( - const DoFHandler & dof, - const hp::QCollection & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const a = nullptr, - const AffineConstraints &constraints = AffineConstraints()); + const DoFHandler & dof, + const hp::QCollection & q, + SparseMatrixType & matrix, + const Function &rhs, + Vector & rhs_vector, + const Function *const a = + nullptr, + const AffineConstraints & + constraints = AffineConstraints()); /** @@ -479,28 +495,32 @@ namespace MatrixCreator * * See the general documentation of this namespace for more information. */ - template + template void create_laplace_matrix( const Mapping & mapping, const DoFHandler &dof, const Quadrature & q, - SparseMatrix & matrix, - const Function *const a = nullptr, - const AffineConstraints &constraints = AffineConstraints()); + SparseMatrixType & matrix, + const Function *const a = + nullptr, + const AffineConstraints & + constraints = AffineConstraints()); /** * Call the create_laplace_matrix() function, see above, with * mapping=MappingQ@(1). */ - template + template void create_laplace_matrix( const DoFHandler &dof, const Quadrature & q, - SparseMatrix & matrix, - const Function *const a = nullptr, - const AffineConstraints &constraints = AffineConstraints()); + SparseMatrixType & matrix, + const Function *const a = + nullptr, + const AffineConstraints & + constraints = AffineConstraints()); /** * Assemble the Laplace matrix and a right hand side vector. If no @@ -520,86 +540,98 @@ namespace MatrixCreator * * See the general documentation of this namespace for more information. */ - template + template void create_laplace_matrix( - const Mapping & mapping, - const DoFHandler &dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const a = nullptr, - const AffineConstraints &constraints = AffineConstraints()); + const Mapping & mapping, + const DoFHandler & dof, + const Quadrature & q, + SparseMatrixType & matrix, + const Function &rhs, + Vector & rhs_vector, + const Function *const a = + nullptr, + const AffineConstraints & + constraints = AffineConstraints()); /** * Call the create_laplace_matrix() function, see above, with * mapping=MappingQ@(1). */ - template + template void create_laplace_matrix( - const DoFHandler &dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const a = nullptr, - const AffineConstraints &constraints = AffineConstraints()); + const DoFHandler & dof, + const Quadrature & q, + SparseMatrixType & matrix, + const Function &rhs, + Vector & rhs_vector, + const Function *const a = + nullptr, + const AffineConstraints & + constraints = AffineConstraints()); /** * Like the functions above, but for hp-objects. */ - template + template void create_laplace_matrix( const hp::MappingCollection &mapping, const DoFHandler & dof, const hp::QCollection & q, - SparseMatrix & matrix, - const Function *const a = nullptr, - const AffineConstraints &constraints = AffineConstraints()); + SparseMatrixType & matrix, + const Function *const a = + nullptr, + const AffineConstraints & + constraints = AffineConstraints()); /** * Like the functions above, but for hp-objects. */ - template + template void create_laplace_matrix( const DoFHandler &dof, const hp::QCollection & q, - SparseMatrix & matrix, - const Function *const a = nullptr, - const AffineConstraints &constraints = AffineConstraints()); + SparseMatrixType & matrix, + const Function *const a = + nullptr, + const AffineConstraints & + constraints = AffineConstraints()); /** * Like the functions above, but for hp-objects. */ - template + template void create_laplace_matrix( - const hp::MappingCollection &mapping, - const DoFHandler & dof, - const hp::QCollection & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const a = nullptr, - const AffineConstraints &constraints = AffineConstraints()); + const hp::MappingCollection & mapping, + const DoFHandler & dof, + const hp::QCollection & q, + SparseMatrixType & matrix, + const Function &rhs, + Vector & rhs_vector, + const Function *const a = + nullptr, + const AffineConstraints & + constraints = AffineConstraints()); /** * Like the functions above, but for hp-objects. */ - template + template void create_laplace_matrix( - const DoFHandler &dof, - const hp::QCollection & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const a = nullptr, - const AffineConstraints &constraints = AffineConstraints()); + const DoFHandler & dof, + const hp::QCollection & q, + SparseMatrixType & matrix, + const Function &rhs, + Vector & rhs_vector, + const Function *const a = + nullptr, + const AffineConstraints & + constraints = AffineConstraints()); /** * Exception diff --git a/source/numerics/matrix_creator.inst.in b/source/numerics/matrix_creator.inst.in index 0eb54b5389..a13918d1a4 100644 --- a/source/numerics/matrix_creator.inst.in +++ b/source/numerics/matrix_creator.inst.in @@ -17,45 +17,62 @@ // MatrixCreator::create_mass_matrix for the matrix only --> only real-valued // scalars -for (scalar : REAL_SCALARS; deal_II_dimension : DIMENSIONS; +for (deal_II_sparse_matrix_type : SPARSE_MATRICES; + deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) { #if deal_II_dimension <= deal_II_space_dimension // non-hp-version of create_mass_matrix - template void MatrixCreator:: - create_mass_matrix( - const Mapping & mapping, - const DoFHandler &dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function *const coefficient, - const AffineConstraints & constraints); - template void MatrixCreator:: - create_mass_matrix( - const DoFHandler &dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function *const coefficient, - const AffineConstraints & constraints); + template void MatrixCreator::create_mass_matrix( + const Mapping & mapping, + const DoFHandler &dof, + const Quadrature & q, + deal_II_sparse_matrix_type & matrix, + const Function + *const coefficient, + const AffineConstraints + &constraints); + template void MatrixCreator::create_mass_matrix( + const DoFHandler &dof, + const Quadrature & q, + deal_II_sparse_matrix_type & matrix, + const Function + *const coefficient, + const AffineConstraints + &constraints); // hp-versions of functions - template void MatrixCreator:: - create_mass_matrix( - const hp::MappingCollection - & mapping, - const DoFHandler &dof, - const hp::QCollection & q, - SparseMatrix & matrix, - const Function *const coefficient, - const AffineConstraints & constraints); + template void MatrixCreator::create_mass_matrix( + const hp::MappingCollection + & mapping, + const DoFHandler &dof, + const hp::QCollection & q, + deal_II_sparse_matrix_type & matrix, + const Function + *const coefficient, + const AffineConstraints + &constraints); - template void MatrixCreator:: - create_mass_matrix( - const DoFHandler &dof, - const hp::QCollection & q, - SparseMatrix & matrix, - const Function *const coefficient, - const AffineConstraints & constraints); + template void MatrixCreator::create_mass_matrix( + const DoFHandler &dof, + const hp::QCollection & q, + deal_II_sparse_matrix_type & matrix, + const Function + *const coefficient, + const AffineConstraints + &constraints); #endif } @@ -67,48 +84,52 @@ for (scalar : REAL_AND_COMPLEX_SCALARS; deal_II_dimension : DIMENSIONS; { #if deal_II_dimension <= deal_II_space_dimension // non-hp-version of create_mass_matrix - template void MatrixCreator:: - create_mass_matrix( - const Mapping & mapping, - const DoFHandler &dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const coefficient, - const AffineConstraints & constraints); - template void MatrixCreator:: - create_mass_matrix( - const DoFHandler &dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const coefficient, - const AffineConstraints & constraints); + template void MatrixCreator::create_mass_matrix>( + const Mapping & mapping, + const DoFHandler &dof, + const Quadrature & q, + SparseMatrix & matrix, + const Function & rhs, + Vector & rhs_vector, + const Function *const coefficient, + const AffineConstraints &constraints); + template void MatrixCreator::create_mass_matrix>( + const DoFHandler &dof, + const Quadrature & q, + SparseMatrix & matrix, + const Function & rhs, + Vector & rhs_vector, + const Function *const coefficient, + const AffineConstraints &constraints); // hp-version - template void MatrixCreator:: - create_mass_matrix( - const hp::MappingCollection - & mapping, - const DoFHandler &dof, - const hp::QCollection & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const coefficient, - const AffineConstraints & constraints); + template void MatrixCreator::create_mass_matrix>( + const hp::MappingCollection + & mapping, + const DoFHandler &dof, + const hp::QCollection & q, + SparseMatrix & matrix, + const Function & rhs, + Vector & rhs_vector, + const Function *const coefficient, + const AffineConstraints &constraints); - template void MatrixCreator:: - create_mass_matrix( - const DoFHandler &dof, - const hp::QCollection & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const coefficient, - const AffineConstraints & constraints); + template void MatrixCreator::create_mass_matrix>( + const DoFHandler &dof, + const hp::QCollection & q, + SparseMatrix & matrix, + const Function & rhs, + Vector & rhs_vector, + const Function *const coefficient, + const AffineConstraints &constraints); #endif } @@ -179,114 +200,192 @@ for (scalar : REAL_AND_COMPLEX_SCALARS; deal_II_dimension : DIMENSIONS; -for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) +for (deal_II_sparse_matrix_type : SPARSE_MATRICES; + deal_II_dimension : DIMENSIONS; + deal_II_space_dimension : SPACE_DIMENSIONS) { #if deal_II_dimension == deal_II_space_dimension // non-hp-versions of create_laplace_matrix - template void MatrixCreator::create_laplace_matrix( - const DoFHandler & dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function *const coefficient, - const AffineConstraints & constraints); - template void MatrixCreator::create_laplace_matrix( - const Mapping & mapping, - const DoFHandler & dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function *const coefficient, - const AffineConstraints & constraints); - template void MatrixCreator::create_laplace_matrix( - const Mapping & mapping, - const DoFHandler & dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const coefficient, - const AffineConstraints & constraints); - template void MatrixCreator::create_laplace_matrix( - const DoFHandler & dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const coefficient, - const AffineConstraints & constraints); + template void + MatrixCreator::create_laplace_matrix( + const DoFHandler &dof, + const Quadrature &q, + deal_II_sparse_matrix_type & matrix, + const Function + *const coefficient, + const AffineConstraints + &constraints); + template void + MatrixCreator::create_laplace_matrix( + const Mapping & mapping, + const DoFHandler &dof, + const Quadrature &q, + deal_II_sparse_matrix_type & matrix, + const Function + *const coefficient, + const AffineConstraints + &constraints); + template void + MatrixCreator::create_laplace_matrix( + const Mapping & mapping, + const DoFHandler & dof, + const Quadrature & q, + deal_II_sparse_matrix_type & matrix, + const Function &rhs, + Vector &rhs_vector, + const Function + *const coefficient, + const AffineConstraints + &constraints); + template void + MatrixCreator::create_laplace_matrix( + const DoFHandler & dof, + const Quadrature & q, + deal_II_sparse_matrix_type & matrix, + const Function &rhs, + Vector &rhs_vector, + const Function + *const coefficient, + const AffineConstraints + &constraints); // hp-versions of create_laplace_matrix - template void MatrixCreator::create_laplace_matrix( + template void + MatrixCreator::create_laplace_matrix( const DoFHandler & dof, const hp::QCollection &q, - SparseMatrix & matrix, - const Function *const coefficient, - const AffineConstraints & constraints); - template void MatrixCreator::create_laplace_matrix( - const hp::MappingCollection &mapping, - const DoFHandler & dof, - const hp::QCollection & q, - SparseMatrix & matrix, - const Function *const coefficient, - const AffineConstraints & constraints); - template void MatrixCreator::create_laplace_matrix( + deal_II_sparse_matrix_type & matrix, + const Function + *const coefficient, + const AffineConstraints + &constraints); + template void + MatrixCreator::create_laplace_matrix( const hp::MappingCollection &mapping, const DoFHandler & dof, const hp::QCollection & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const coefficient, - const AffineConstraints & constraints); - template void MatrixCreator::create_laplace_matrix( - const DoFHandler & dof, - const hp::QCollection &q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const coefficient, - const AffineConstraints & constraints); + deal_II_sparse_matrix_type & matrix, + const Function + *const coefficient, + const AffineConstraints + &constraints); + template void + MatrixCreator::create_laplace_matrix( + const hp::MappingCollection & mapping, + const DoFHandler & dof, + const hp::QCollection & q, + deal_II_sparse_matrix_type & matrix, + const Function &rhs, + Vector &rhs_vector, + const Function + *const coefficient, + const AffineConstraints + &constraints); + template void + MatrixCreator::create_laplace_matrix( + const DoFHandler & dof, + const hp::QCollection & q, + deal_II_sparse_matrix_type & matrix, + const Function &rhs, + Vector &rhs_vector, + const Function + *const coefficient, + const AffineConstraints + &constraints); #endif } -for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) +for (deal_II_sparse_matrix_type : SPARSE_MATRICES; + deal_II_dimension : DIMENSIONS; + deal_II_space_dimension : SPACE_DIMENSIONS) { #if deal_II_dimension < deal_II_space_dimension // non-hp-version of create_laplace_matrix - template void MatrixCreator::create_laplace_matrix( - const Mapping & mapping, - const DoFHandler &dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function *const coefficient, - const AffineConstraints &constraints); - template void MatrixCreator::create_laplace_matrix( - const DoFHandler &dof, - const Quadrature & q, - SparseMatrix & matrix, - const Function *const coefficient, - const AffineConstraints &constraints); - template void MatrixCreator::create_laplace_matrix( + template void + MatrixCreator::create_laplace_matrix( const Mapping & mapping, const DoFHandler &dof, const Quadrature & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const coefficient, - const AffineConstraints &constraints); - template void MatrixCreator::create_laplace_matrix( + deal_II_sparse_matrix_type & matrix, + const Function + *const coefficient, + const AffineConstraints + &constraints); + template void + MatrixCreator::create_laplace_matrix( const DoFHandler &dof, const Quadrature & q, - SparseMatrix & matrix, - const Function & rhs, - Vector & rhs_vector, - const Function *const coefficient, - const AffineConstraints &constraints); + deal_II_sparse_matrix_type & matrix, + const Function + *const coefficient, + const AffineConstraints + &constraints); + template void + MatrixCreator::create_laplace_matrix( + const Mapping & mapping, + const DoFHandler & dof, + const Quadrature & q, + deal_II_sparse_matrix_type & matrix, + const Function &rhs, + Vector &rhs_vector, + const Function + *const coefficient, + const AffineConstraints + &constraints); + template void + MatrixCreator::create_laplace_matrix( + const DoFHandler & dof, + const Quadrature & q, + deal_II_sparse_matrix_type & matrix, + const Function &rhs, + Vector &rhs_vector, + const Function + *const coefficient, + const AffineConstraints + &constraints); #endif } diff --git a/tests/numerics/matrix_creator_01.cc b/tests/numerics/matrix_creator_01.cc new file mode 100644 index 0000000000..07f893cd41 --- /dev/null +++ b/tests/numerics/matrix_creator_01.cc @@ -0,0 +1,97 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2021 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + +// Test MatrixCreator::create_mass_matrix() and +// MatrixCreator::create_laplace_matrix() for Trilinos and PETSc matrices in +// parallel. + +#include +#include + +#include + +#include + +#include + +#include + +#include +#include + +#include + +#include "../tests.h" + +using namespace dealii; + +template +void +test() +{ + parallel::distributed::Triangulation tria(MPI_COMM_WORLD); + GridGenerator::hyper_cube(tria); + tria.refine_global(2); + + FE_Q fe(2); + QGauss quad(3); + + DoFHandler dof_handler(tria); + dof_handler.distribute_dofs(fe); + + const auto mpi_comm = dof_handler.get_communicator(); + const IndexSet &owned_dofs = dof_handler.locally_owned_dofs(); + + IndexSet relevant_dofs; + DoFTools::extract_locally_relevant_dofs(dof_handler, relevant_dofs); + + AffineConstraints constraints; + + DynamicSparsityPattern sparsity_pattern(relevant_dofs); + DoFTools::make_sparsity_pattern(dof_handler, sparsity_pattern, constraints); + SparsityTools::distribute_sparsity_pattern(sparsity_pattern, + owned_dofs, + mpi_comm, + relevant_dofs); + + MatrixType matrix; + matrix.reinit(owned_dofs, owned_dofs, sparsity_pattern, mpi_comm); + + MatrixCreator::create_laplace_matrix(dof_handler, quad, matrix); + deallog << "Laplace matrix: " << matrix.frobenius_norm() << std::endl; + + matrix = 0.0; + MatrixCreator::create_mass_matrix(dof_handler, quad, matrix); + deallog << "Mass matrix: " << matrix.frobenius_norm() << std::endl; +} + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + mpi_initlog(); + + { + deallog.push("Trilinos"); + test<2, TrilinosWrappers::SparseMatrix>(); + deallog.pop(); + } + { + deallog.push("PETSc"); + test<2, PETScWrappers::MPI::SparseMatrix>(); + deallog.pop(); + } +} diff --git a/tests/numerics/matrix_creator_01.with_p4est=true.mpirun=1.with_petsc=on.with_trilinos=on.output b/tests/numerics/matrix_creator_01.with_p4est=true.mpirun=1.with_petsc=on.with_trilinos=on.output new file mode 100644 index 0000000000..b78c103684 --- /dev/null +++ b/tests/numerics/matrix_creator_01.with_p4est=true.mpirun=1.with_petsc=on.with_trilinos=on.output @@ -0,0 +1,5 @@ + +DEAL:Trilinos::Laplace matrix: 34.9481 +DEAL:Trilinos::Mass matrix: 0.0916667 +DEAL:PETSc::Laplace matrix: 34.9481 +DEAL:PETSc::Mass matrix: 0.0916667 diff --git a/tests/numerics/matrix_creator_01.with_p4est=true.mpirun=3.with_petsc=on.with_trilinos=on.output b/tests/numerics/matrix_creator_01.with_p4est=true.mpirun=3.with_petsc=on.with_trilinos=on.output new file mode 100644 index 0000000000..b78c103684 --- /dev/null +++ b/tests/numerics/matrix_creator_01.with_p4est=true.mpirun=3.with_petsc=on.with_trilinos=on.output @@ -0,0 +1,5 @@ + +DEAL:Trilinos::Laplace matrix: 34.9481 +DEAL:Trilinos::Mass matrix: 0.0916667 +DEAL:PETSc::Laplace matrix: 34.9481 +DEAL:PETSc::Mass matrix: 0.0916667