From: Wolfgang Bangerth Date: Mon, 22 Oct 2007 18:59:34 +0000 (+0000) Subject: Only rebuild the matrix when necessary. X-Git-Tag: v8.0.0~9721 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=cf98114f361f2cbcdc6cdb66c790e37216d518bc;p=dealii.git Only rebuild the matrix when necessary. git-svn-id: https://svn.dealii.org/trunk@15364 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/examples/step-22/step-22.cc b/deal.II/examples/step-22/step-22.cc index 2508f19d10..3e9e170600 100644 --- a/deal.II/examples/step-22/step-22.cc +++ b/deal.II/examples/step-22/step-22.cc @@ -92,7 +92,8 @@ class BoussinesqFlowProblem BlockVector system_rhs; boost::shared_ptr A_preconditioner; - + + bool rebuild_matrices; bool rebuild_preconditioner; }; @@ -402,6 +403,7 @@ BoussinesqFlowProblem::BoussinesqFlowProblem (const unsigned int degree) FE_DGQ(degree-1), 1), dof_handler (triangulation), time_step (0), + rebuild_matrices (true), rebuild_preconditioner (true) {} @@ -499,15 +501,24 @@ scalar_product (const Tensor<2,dim> &a, template void BoussinesqFlowProblem::assemble_system () { - system_matrix=0; - system_rhs=0; - + if (rebuild_matrices == true) + { + system_matrix=0; + system_rhs=0; + } + QGauss quadrature_formula(degree+2); QGauss face_quadrature_formula(degree+2); - FEValues fe_values (fe, quadrature_formula, - update_values | update_gradients | - update_quadrature_points | update_JxW_values); + FEValues fe_values (fe, quadrature_formula, + update_values | + update_quadrature_points | + update_JxW_values | + (rebuild_matrices == true + ? + update_gradients + : + UpdateFlags(0))); FEFaceValues fe_face_values (fe, face_quadrature_formula, update_values | update_normal_vectors | update_quadrature_points | update_JxW_values); @@ -551,27 +562,31 @@ void BoussinesqFlowProblem::assemble_system () { const Tensor<1,dim> phi_i_u = extract_u (fe_values, i, q); - const Tensor<2,dim> phi_i_grads_u= extract_grad_s_u (fe_values, i, q); - const double div_phi_i_u = extract_div_u (fe_values, i, q); - const double phi_i_p = extract_p (fe_values, i, q); - const double phi_i_T = extract_T (fe_values, i, q); - const Tensor<1,dim> grad_phi_i_T = extract_grad_T(fe_values, i, q); - - for (unsigned int j=0; j phi_j_grads_u = extract_grad_s_u (fe_values, j, q); - const double div_phi_j_u = extract_div_u (fe_values, j, q); - const double phi_j_p = extract_p (fe_values, j, q); - const double phi_j_T = extract_T (fe_values, j, q); + const Tensor<2,dim> phi_i_grads_u= extract_grad_s_u (fe_values, i, q); + const double div_phi_i_u = extract_div_u (fe_values, i, q); + const double phi_i_p = extract_p (fe_values, i, q); + const double phi_i_T = extract_T (fe_values, i, q); + const Tensor<1,dim> grad_phi_i_T = extract_grad_T(fe_values, i, q); + + for (unsigned int j=0; j phi_j_grads_u = extract_grad_s_u (fe_values, j, q); + const double div_phi_j_u = extract_div_u (fe_values, j, q); + const double phi_j_p = extract_p (fe_values, j, q); + const double phi_j_T = extract_T (fe_values, j, q); - local_matrix(i,j) += (scalar_product(phi_i_grads_u, phi_j_grads_u) - - div_phi_i_u * phi_j_p - - phi_i_p * div_phi_j_u - + phi_i_p * phi_j_p - + phi_i_T * phi_j_T) - * fe_values.JxW(q); + local_matrix(i,j) += (scalar_product(phi_i_grads_u, phi_j_grads_u) + - div_phi_i_u * phi_j_p + - phi_i_p * div_phi_j_u + + phi_i_p * phi_j_p + + phi_i_T * phi_j_T) + * fe_values.JxW(q); + } } - + const Point gravity = fe_values.quadrature_point(q) / fe_values.quadrature_point(q).norm(); @@ -607,36 +622,47 @@ void BoussinesqFlowProblem::assemble_system () } cell->get_dof_indices (local_dof_indices); - for (unsigned int i=0; i component_mask (dim+2, true); - component_mask[dim] = component_mask[dim+1] = false; - std::map boundary_values; - VectorTools::interpolate_boundary_values (dof_handler, - 0, - ZeroFunction(dim+2), - boundary_values, - component_mask); - MatrixTools::apply_boundary_values (boundary_values, - system_matrix, - solution, - system_rhs); - } + if (rebuild_matrices == true) + { + std::vector component_mask (dim+2, true); + component_mask[dim] = component_mask[dim+1] = false; + std::map boundary_values; + VectorTools::interpolate_boundary_values (dof_handler, + 0, + ZeroFunction(dim+2), + boundary_values, + component_mask); + MatrixTools::apply_boundary_values (boundary_values, + system_matrix, + solution, + system_rhs); + } if (rebuild_preconditioner == true) { + Assert (rebuild_matrices == true, + ExcMessage ("There is no point in rebuilding the preconditioner " + "without a rebuilt matrix!")); + std::cout << " Rebuilding preconditioner..." << std::flush; A_preconditioner @@ -647,6 +673,8 @@ void BoussinesqFlowProblem::assemble_system () rebuild_preconditioner = false; } + + rebuild_matrices = false; } @@ -1092,6 +1120,7 @@ BoussinesqFlowProblem::refine_mesh () Vector tmp (dof_handler.n_dofs()); soltrans.interpolate(x_old_solution, tmp); + rebuild_matrices = true; rebuild_preconditioner = true; old_solution = tmp;