From: wolf Date: Tue, 3 Jun 2003 14:50:47 +0000 (+0000) Subject: Add. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=1a62692902388246ae23671599a30062021caa96;p=dealii-svn.git Add. git-svn-id: https://svn.dealii.org/trunk@7729 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/tests/bits/block_sparse_matrix_1.cc b/tests/bits/block_sparse_matrix_1.cc new file mode 100644 index 0000000000..2decc20e67 --- /dev/null +++ b/tests/bits/block_sparse_matrix_1.cc @@ -0,0 +1,124 @@ +//---------------------------- block_sparse_matrix_1.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2002, 2003 by the deal.II authors and Brian Carnes +// +// 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. +// +//---------------------------- block_sparse_matrix_1.cc --------------------------- + +// test by Brian: check some of the scaling operations on matrices + +#include +#include +#include + +#include +#include + +#include +#include + +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include + +#include +#include + + + +int main() +{ + std::ofstream logfile("block_sparse_matrix_1.output"); + deallog.attach(logfile); + deallog.depth_console(0); + + BlockSparseMatrix B; + + Triangulation<2> tria; + GridGenerator::hyper_cube (tria,0,1); + tria.refine_global (1); + + FE_Q_Hierarchical<2> fe (1); + DoFHandler<2> dof_handler (tria); + dof_handler.distribute_dofs (fe); + + BlockSparsityPattern sparsity_pattern; + sparsity_pattern.reinit (2,2); + sparsity_pattern.collect_sizes (); + + sparsity_pattern.block(0,0).reinit (dof_handler.n_dofs(), + dof_handler.n_dofs(), + dof_handler.max_couplings_between_dofs()); + sparsity_pattern.block(0,1).reinit (dof_handler.n_dofs(), + 1, + 1); + sparsity_pattern.block(1,0).reinit (1, + dof_handler.n_dofs(), + dof_handler.n_dofs()); + sparsity_pattern.block(1,1).reinit (1, + 1, + 1); + sparsity_pattern.collect_sizes (); + + DoFTools::make_sparsity_pattern (dof_handler, sparsity_pattern.block(0,0)); + + for (unsigned int j=0; j qr (2); + FEValues<2> fe_values (dof_handler.get_fe(), + qr, UpdateFlags(update_values | + update_gradients | + update_q_points | + update_JxW_values)); + + MatrixTools::create_laplace_matrix (dof_handler, qr, B.block(0,0)); + B.block(1,1).add (0,0,1.); + + B.print_formatted (deallog.get_file_stream (),3,false); + + B.block(0,0) *= 10.; + B.print_formatted (deallog.get_file_stream (),3,false); + + B *= 10.; + B.print_formatted (deallog.get_file_stream (),3,false); + + B /= 10.; + B.print_formatted (deallog.get_file_stream (),3,false); + + B.block(0,0) /= 10.; + B.print_formatted (deallog.get_file_stream (),3,false); + + B.clear(); + + return 0; +}