From 19125bfe0f4dee8e18b6c30db3acac773a02fa85 Mon Sep 17 00:00:00 2001 From: bangerth Date: Tue, 19 Aug 2008 20:18:47 +0000 Subject: [PATCH] Rebuild temperature matrices only when needed. Do not forget to condense the temperature system with hanging nodes. git-svn-id: https://svn.dealii.org/trunk@16597 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/examples/step-31/step-31.cc | 186 ++++++++++++++++++++-------- 1 file changed, 134 insertions(+), 52 deletions(-) diff --git a/deal.II/examples/step-31/step-31.cc b/deal.II/examples/step-31/step-31.cc index b1f21f0e36..37ccb5b971 100644 --- a/deal.II/examples/step-31/step-31.cc +++ b/deal.II/examples/step-31/step-31.cc @@ -323,6 +323,7 @@ class BoussinesqFlowProblem void build_stokes_preconditioner (); void assemble_stokes_system (); void assemble_temperature_system (); + void assemble_temperature_matrix (); double get_maximal_velocity () const; double get_maximal_temperature () const; void solve (); @@ -352,6 +353,8 @@ class BoussinesqFlowProblem ConstraintMatrix temperature_constraints; SparsityPattern temperature_sparsity_pattern; + SparseMatrix temperature_mass_matrix; + SparseMatrix temperature_stiffness_matrix; SparseMatrix temperature_matrix; Vector temperature_solution; @@ -367,8 +370,9 @@ class BoussinesqFlowProblem boost::shared_ptr Amg_preconditioner; boost::shared_ptr > Mp_preconditioner; - bool rebuild_matrices; - bool rebuild_preconditioner; + bool rebuild_stokes_matrix; + bool rebuild_temperature_matrix; + bool rebuild_stokes_preconditioner; }; @@ -818,8 +822,9 @@ BoussinesqFlowProblem::BoussinesqFlowProblem (const unsigned int degree) time_step (0), old_time_step (0), timestep_number (0), - rebuild_matrices (true), - rebuild_preconditioner (true) + rebuild_stokes_matrix (true), + rebuild_temperature_matrix (true), + rebuild_stokes_preconditioner (true) {} @@ -1115,6 +1120,8 @@ void BoussinesqFlowProblem::setup_dofs () } { + temperature_mass_matrix.clear (); + temperature_stiffness_matrix.clear (); temperature_matrix.clear (); CompressedSetSparsityPattern csp (n_T, n_T); @@ -1123,6 +1130,8 @@ void BoussinesqFlowProblem::setup_dofs () temperature_sparsity_pattern.copy_from (csp); temperature_matrix.reinit (temperature_sparsity_pattern); + temperature_mass_matrix.reinit (temperature_sparsity_pattern); + temperature_stiffness_matrix.reinit (temperature_sparsity_pattern); } @@ -1233,6 +1242,9 @@ template void BoussinesqFlowProblem::build_stokes_preconditioner () { + if (rebuild_stokes_preconditioner == false) + return; + std::cout << " Rebuilding Stokes preconditioner..." << std::flush; @@ -1256,13 +1268,6 @@ BoussinesqFlowProblem::build_stokes_preconditioner () // include them using a // shared pointer structure from the // boost library. - // - // When all work is done, we - // change the flags - // rebuild_preconditioner - // and - // rebuild_matrices - // to false. assemble_stokes_preconditioner (); Amg_preconditioner = boost::shared_ptr @@ -1319,7 +1324,7 @@ BoussinesqFlowProblem::build_stokes_preconditioner () std::cout << std::endl; - rebuild_preconditioner = false; + rebuild_stokes_preconditioner = false; } @@ -1397,7 +1402,7 @@ void BoussinesqFlowProblem::assemble_stokes_system () { std::cout << " Assembling..." << std::flush; - if (rebuild_matrices == true) + if (rebuild_stokes_matrix == true) stokes_matrix=0; stokes_rhs=0; @@ -1409,7 +1414,7 @@ void BoussinesqFlowProblem::assemble_stokes_system () update_values | update_quadrature_points | update_JxW_values | - (rebuild_matrices == true + (rebuild_stokes_matrix == true ? update_gradients : @@ -1542,7 +1547,7 @@ void BoussinesqFlowProblem::assemble_stokes_system () for (unsigned int k=0; k::assemble_stokes_system () // define viscosity const double eta = 1; - if (rebuild_matrices) + if (rebuild_stokes_matrix) for (unsigned int i=0; i::assemble_stokes_system () // requested. cell->get_dof_indices (local_dof_indices); - if (rebuild_matrices == true) + if (rebuild_stokes_matrix == true) for (unsigned int i=0; i::assemble_stokes_system () // possibly, to the matrix. stokes_constraints.condense (stokes_rhs); - if (rebuild_matrices == true) + if (rebuild_stokes_matrix == true) { stokes_constraints.condense (stokes_matrix); - + rebuild_stokes_matrix = false; + // std::map boundary_values; // typename DoFHandler::active_cell_iterator @@ -1675,8 +1681,6 @@ void BoussinesqFlowProblem::assemble_stokes_system () // system_rhs); } - rebuild_matrices = false; - std::cout << std::endl; } @@ -1771,12 +1775,113 @@ double compute_viscosity( // and subfaces. The update // flags at face level are the // same as in step-12. +template +void BoussinesqFlowProblem::assemble_temperature_matrix () +{ + if (rebuild_temperature_matrix == false) + return; + + temperature_mass_matrix = 0; + temperature_stiffness_matrix = 0; + + QGauss quadrature_formula(temperature_degree+2); + FEValues temperature_fe_values (temperature_fe, quadrature_formula, + update_values | update_gradients | + update_JxW_values); + + const unsigned int dofs_per_cell = temperature_fe.dofs_per_cell; + const unsigned int n_q_points = quadrature_formula.size(); + + FullMatrix local_mass_matrix (dofs_per_cell, dofs_per_cell); + FullMatrix local_stiffness_matrix (dofs_per_cell, dofs_per_cell); + + std::vector local_dof_indices (dofs_per_cell); + + std::vector gamma_values (n_q_points); + + std::vector phi_T (dofs_per_cell); + std::vector > grad_phi_T (dofs_per_cell); + + // Now, let's start the loop + // over all cells in the + // triangulation. The first + // actions within the loop + // are, 0as usual, the evaluation + // of the FE basis functions + // and the old and present + // solution at the quadrature + // points. + typename DoFHandler::active_cell_iterator + cell = temperature_dof_handler.begin_active(), + endc = temperature_dof_handler.end(); + for (; cell!=endc; ++cell) + { + local_mass_matrix = 0; + local_stiffness_matrix = 0; + + temperature_fe_values.reinit (cell); + + const double kappa = 1e-6; + + for (unsigned int q=0; qget_dof_indices (local_dof_indices); + + for (unsigned int i=0; i void BoussinesqFlowProblem::assemble_temperature_system () { const bool use_bdf2_scheme = (timestep_number != 0); - temperature_matrix = 0; + if (use_bdf2_scheme == true) + { + temperature_matrix.copy_from (temperature_mass_matrix); + temperature_matrix *= (2*time_step + old_time_step) / + (time_step + old_time_step); + temperature_matrix.add (time_step, temperature_stiffness_matrix); + } + else + { + temperature_matrix.copy_from (temperature_mass_matrix); + temperature_matrix.add (time_step, temperature_stiffness_matrix); + } + temperature_rhs = 0; QGauss quadrature_formula(temperature_degree+2); @@ -1848,7 +1953,6 @@ void BoussinesqFlowProblem::assemble_temperature_system () for (; cell!=endc; ++cell, ++stokes_cell) { local_rhs = 0; - local_matrix = 0; temperature_fe_values.reinit (cell); stokes_fe_values.reinit (stokes_cell); @@ -1910,16 +2014,6 @@ void BoussinesqFlowProblem::assemble_temperature_system () if (use_bdf2_scheme == true) { - for (unsigned int i=0; i::assemble_temperature_system () } else { - for (unsigned int i=0; i::assemble_temperature_system () cell->get_dof_indices (local_dof_indices); - for (unsigned int i=0; i::refine_mesh (const unsigned int max_grid_level) temperature_solution = tmp[0]; old_temperature_solution = tmp[1]; - rebuild_matrices = true; - rebuild_preconditioner = true; + rebuild_stokes_matrix = true; + rebuild_temperature_matrix = true; + rebuild_stokes_preconditioner = true; } @@ -2328,9 +2411,8 @@ void BoussinesqFlowProblem::run () << std::endl; assemble_stokes_system (); - - if (rebuild_preconditioner == true) - build_stokes_preconditioner (); + build_stokes_preconditioner (); + assemble_temperature_matrix (); solve (); -- 2.39.5