From 237b1d5e70b58053c1852c58723db89ea780ecf3 Mon Sep 17 00:00:00 2001 From: bangerth Date: Wed, 13 Aug 2008 01:04:45 +0000 Subject: [PATCH] No longer build the preconditioner for the pressure mass matrix in the (1,1) block of the system matrix but in the preconditioner_matrix proper. That also means that we can't throw away the preconditioner_matrix at the end of building the preconditioner any more, though. git-svn-id: https://svn.dealii.org/trunk@16529 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/examples/step-31/step-31.cc | 36 +++++++++++++++++------------ 1 file changed, 21 insertions(+), 15 deletions(-) diff --git a/deal.II/examples/step-31/step-31.cc b/deal.II/examples/step-31/step-31.cc index 21b72dbcdf..97deffcca1 100644 --- a/deal.II/examples/step-31/step-31.cc +++ b/deal.II/examples/step-31/step-31.cc @@ -1033,6 +1033,7 @@ void BoussinesqFlowProblem::setup_dofs (const bool setup_matrices) sparsity_pattern.copy_from (csp); system_matrix.reinit (sparsity_pattern); + preconditioner_matrix.reinit (sparsity_pattern); } // As last action in this function, @@ -1088,6 +1089,7 @@ BoussinesqFlowProblem::assemble_preconditioner () FEValues fe_values (fe, quadrature_formula, update_JxW_values | + update_values | update_gradients); const unsigned int dofs_per_cell = fe.dofs_per_cell; @@ -1097,8 +1099,10 @@ BoussinesqFlowProblem::assemble_preconditioner () std::vector local_dof_indices (dofs_per_cell); std::vector > phi_grad_u (dofs_per_cell); + std::vector phi_p (dofs_per_cell); const FEValuesExtractors::Vector velocities (0); + const FEValuesExtractors::Scalar pressure (dim); typename DoFHandler::active_cell_iterator cell = dof_handler.begin_active(), @@ -1111,11 +1115,16 @@ BoussinesqFlowProblem::assemble_preconditioner () for (unsigned int q=0; q::assemble_system () local_matrix(i,j) += (phi_grads_u[i] * phi_grads_u[j] - div_phi_u[i] * phi_p[j] - phi_p[i] * div_phi_u[j] - + phi_p[i] * phi_p[j] + phi_T[i] * phi_T[j]) * fe_values.JxW(q); @@ -1522,7 +1530,6 @@ void BoussinesqFlowProblem::assemble_system () std::cout << " Rebuilding preconditioner..." << std::flush; - preconditioner_matrix.reinit (sparsity_pattern); assemble_preconditioner (); /*A_preconditioner @@ -1537,18 +1544,17 @@ void BoussinesqFlowProblem::assemble_system () preconditioner_matrix.block(0,0), true)); + // TODO: we could throw away the (0,0) + // block here since things have been + // copied over to Trilinos. we need to + // keep the (1,1) block, though + Mp_preconditioner = boost::shared_ptr > (new SparseILU); - Mp_preconditioner->initialize (system_matrix.block(1,1), + Mp_preconditioner->initialize (preconditioner_matrix.block(1,1), SparseILU::AdditionalData()); - // Throw away the preconditioner - // matrix since everything has been - // copied to the Epetra objects - // of the preconditioner. - preconditioner_matrix.clear (); - std::cout << std::endl; rebuild_preconditioner = false; @@ -1932,7 +1938,7 @@ void BoussinesqFlowProblem::solve () // Set up inverse matrix for // pressure mass matrix InverseMatrix,SparseILU > - mp_inverse (system_matrix.block(1,1), *Mp_preconditioner); + mp_inverse (preconditioner_matrix.block(1,1), *Mp_preconditioner); // Set up block Schur preconditioner /*BlockSchurPreconditioner::type, @@ -2138,7 +2144,7 @@ void BoussinesqFlowProblem::run () GridGenerator::hyper_cube (triangulation); - triangulation.refine_global (6); + triangulation.refine_global (4); break; } @@ -2224,7 +2230,7 @@ void BoussinesqFlowProblem::run () if (timestep_number % 10 == 0) refine_mesh (); } - while (time <= 50); + while (time <= 10); } -- 2.39.5