Merge pull request #1962 from gassmoeller/unify_rebuild_of_stokes_matrix

Unify rebuild of stokes matrix

source/simulator/core.cc

@@ -825,16 +825,6 @@ namespace aspect
         rebuild_stokes_matrix = rebuild_stokes_preconditioner = true;
       }
 
-
-    // TODO: do this in a more efficient way (TH)? we really only need
-    // to make sure that the time dependent velocity boundary conditions
-    // end up in the right hand side in the right way; we currently do
-    // that by re-assembling the entire system
-    if (!boundary_velocity.empty())
-      rebuild_stokes_matrix = rebuild_stokes_preconditioner = assemble_newton_stokes_matrix = true;
-
-
-
     // notify different system components that we started the next time step
     // TODO: implement this for all plugins that might need it at one place.
     // Temperature BC are currently updated in compute_current_constraints

source/simulator/helper_functions.cc

@@ -1236,15 +1236,17 @@ namespace aspect
   bool
   Simulator<dim>::stokes_matrix_depends_on_solution() const
   {
-    // currently, the only coefficient that really appears on the
-    // left hand side of the Stokes equation is the viscosity. note
-    // that our implementation of compressible materials makes sure
-    // that the density does not appear on the LHS.
-    // if melt transport is included in the simulation, we have an
+    // Currently, the only coefficient that really appears on the
+    // left hand side of the Stokes equation is the viscosity and possibly
+    // the density in the case of the implicit reference density profile
+    // approximation.
+    // If melt transport is included in the simulation, we have an
     // additional equation with more coefficients on the left hand
     // side.
 
     return (material_model->get_model_dependence().viscosity != MaterialModel::NonlinearDependence::none)
+           || (parameters.formulation_mass_conservation ==
+               Parameters<dim>::Formulation::MassConservation::implicit_reference_density_profile)
            || parameters.include_melt_transport;
   }
 

source/simulator/solver_schemes.cc

@@ -215,6 +215,14 @@ namespace aspect
   double Simulator<dim>::assemble_and_solve_stokes (const bool compute_initial_residual,
                                                     double *initial_residual)
   {
+    // If the Stokes matrix depends on the solution, or the boundary conditions
+    // for the Stokes system have changed rebuild the matrix and preconditioner
+    // before solving.
+    if (stokes_matrix_depends_on_solution()
+        ||
+        (parameters.prescribed_velocity_boundary_indicators.size() > 0))
+      rebuild_stokes_matrix = rebuild_stokes_preconditioner = true;
+
     assemble_stokes_system ();
     build_stokes_preconditioner();
 
@@ -256,14 +264,6 @@ namespace aspect
   {
     assemble_and_solve_temperature();
     assemble_and_solve_composition();
-
-    // the Stokes matrix depends on the viscosity. if the viscosity
-    // depends on other solution variables, then after we need to
-    // update the Stokes matrix in every time step and so need to set
-    // the following flag. if we change the Stokes matrix we also
-    // need to update the Stokes preconditioner.
-    if (stokes_matrix_depends_on_solution() == true)
-      rebuild_stokes_matrix = rebuild_stokes_preconditioner = true;
     assemble_and_solve_stokes();
 
     if (parameters.run_postprocessors_on_nonlinear_iterations)
@@ -288,13 +288,6 @@ namespace aspect
       parameters.max_nonlinear_iterations;
     do
       {
-        if ((stokes_matrix_depends_on_solution() == true)
-            ||
-            (parameters.prescribed_velocity_boundary_indicators.size() > 0))
-          rebuild_stokes_matrix = true;
-        if (stokes_matrix_depends_on_solution() == true)
-          rebuild_stokes_preconditioner = true;
-
         const double relative_stokes_residual =
           assemble_and_solve_stokes(nonlinear_iteration == 0, &initial_stokes_residual);
 
@@ -340,19 +333,6 @@ namespace aspect
         const std::vector<double>  relative_composition_residual =
           assemble_and_solve_composition(nonlinear_iteration == 0, &initial_composition_residual);
 
-        // TODO: This is a leftover from early history and we keep it for now
-        // to keep test results unchanged. Unify the setting of these
-        // variables between all nonlinear solver schemes, by moving them into
-        // assemble_and_solve_stokes().
-        rebuild_stokes_matrix = true;
-        // the Stokes matrix depends on the viscosity. if the viscosity
-        // depends on other solution variables, then after we need to
-        // update the Stokes matrix in every time step and so need to set
-        // the following flag. if we change the Stokes matrix we also
-        // need to update the Stokes preconditioner.
-        if (stokes_matrix_depends_on_solution() == true)
-          rebuild_stokes_matrix = rebuild_stokes_preconditioner = true;
-
         const double relative_stokes_residual =
           assemble_and_solve_stokes(nonlinear_iteration == 0, &initial_stokes_residual);
 
@@ -421,13 +401,6 @@ namespace aspect
 
     do
       {
-        // rebuild the matrix if it actually depends on the solution
-        // of the previous iteration.
-        if ((stokes_matrix_depends_on_solution() == true)
-            ||
-            (parameters.prescribed_velocity_boundary_indicators.size() > 0))
-          rebuild_stokes_matrix = rebuild_stokes_preconditioner = true;
-
         const double relative_stokes_residual =
           assemble_and_solve_stokes(nonlinear_iteration == 0, &initial_stokes_residual);
 

tests/annulus/screen-output

@@ -8,13 +8,14 @@ Number of degrees of freedom: 2,832 (1,728+240+864)
    Rebuilding Stokes preconditioner...
    Solving Stokes system... 0+23 iterations.
       Relative Stokes residual after nonlinear iteration 1: 1
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 0+0 iterations.
-      Relative Stokes residual after nonlinear iteration 2: 2.41392e-13
+      Relative Stokes residual after nonlinear iteration 2: 2.41488e-13
 
    Postprocessing:
      Writing graphical output:               output-annulus/solution/solution-00000
      RMS, max velocity:                      1.08 m/s, 2.15 m/s
-     Pressure at top/bottom of domain:       -4.108e-16 Pa, 1000 Pa
+     Pressure at top/bottom of domain:       -4.401e-16 Pa, 1000 Pa
      Computing dynamic topography            
      Errors u_L1, p_L1, u_L2, p_L2, topo_L2: 3.470600e-02, 5.832529e-01, 1.113925e-02, 3.426625e-01, 5.533544e-02
 

tests/burstedde/screen-output

@@ -8,12 +8,13 @@ Number of degrees of freedom: 3,041 (2,187+125+729)
    Rebuilding Stokes preconditioner...
    Solving Stokes system... 0+XYZ iterations.
       Relative Stokes residual after nonlinear iteration 1: 1
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 0+XYZ iterations.
-      Relative Stokes residual after nonlinear iteration 2: 9.71521e-12
+      Relative Stokes residual after nonlinear iteration 2: 4.81093e-13
 
    Postprocessing:
      RMS, max velocity:             4.3 m/s, 13.2 m/s
-     Errors u_L1, p_L1, u_L2, p_L2: 6.029937e-01, 1.448270e-03, 1.059164e+00, 3.023064e-03
+     Errors u_L1, p_L1, u_L2, p_L2: 6.029938e-01, 1.448270e-03, 1.059164e+00, 3.023064e-03
 
 Termination requested by criterion: end time
 

tests/burstedde_stokes_rhs/screen-output

@@ -8,8 +8,9 @@ Number of degrees of freedom: 3,041 (2,187+125+729)
    Rebuilding Stokes preconditioner...
    Solving Stokes system... 0+XYZ iterations.
       Relative Stokes residual after nonlinear iteration 1: 1
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 0+XYZ iterations.
-      Relative Stokes residual after nonlinear iteration 2: 6.80562e-13
+      Relative Stokes residual after nonlinear iteration 2: 6.95385e-13
 
    Postprocessing:
      RMS, max velocity:             4.3 m/s, 13.2 m/s

tests/composition_reaction_iterated_IMPES/screen-output

@@ -19,6 +19,7 @@ Number of degrees of freedom: 26,439 (8,450+1,089+4,225+4,225+4,225+4,225)
    Solving C_1 system ... 11 iterations.
    Solving C_2 system ... 0 iterations.
    Solving C_3 system ... 11 iterations.
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 40+0 iterations.
       Relative nonlinear residuals: 1.72329e-16, 0, 5.7186e-13, 0, 0.00772709
       Total relative residual after nonlinear iteration 2: 0.00772709
@@ -28,6 +29,7 @@ Number of degrees of freedom: 26,439 (8,450+1,089+4,225+4,225+4,225+4,225)
    Solving C_1 system ... 0 iterations.
    Solving C_2 system ... 0 iterations.
    Solving C_3 system ... 0 iterations.
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 0+0 iterations.
       Relative nonlinear residuals: 1.82797e-16, 0, 5.7186e-13, 0, 8.29714e-08
       Total relative residual after nonlinear iteration 3: 8.29714e-08
@@ -53,6 +55,7 @@ Number of degrees of freedom: 26,439 (8,450+1,089+4,225+4,225+4,225+4,225)
    Solving C_1 system ... 11 iterations.
    Solving C_2 system ... 0 iterations.
    Solving C_3 system ... 11 iterations.
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 40+0 iterations.
       Relative nonlinear residuals: 1.89051e-07, 0.5, 9.41726e-13, 0.5, 0.00783271
       Total relative residual after nonlinear iteration 2: 0.5
@@ -62,6 +65,7 @@ Number of degrees of freedom: 26,439 (8,450+1,089+4,225+4,225+4,225+4,225)
    Solving C_1 system ... 0 iterations.
    Solving C_2 system ... 0 iterations.
    Solving C_3 system ... 0 iterations.
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 0+0 iterations.
       Relative nonlinear residuals: 1.40557e-07, 9.41813e-13, 9.41721e-13, 9.41813e-13, 8.42599e-08
       Total relative residual after nonlinear iteration 3: 1.40557e-07

tests/compressibility_iterated_stokes/screen-output

@@ -8,12 +8,16 @@ Number of degrees of freedom: 13,764 (8,450+1,089+4,225)
    Rebuilding Stokes preconditioner...
    Solving Stokes system... 33+0 iterations.
       Relative Stokes residual after nonlinear iteration 1: 1
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 28+0 iterations.
       Relative Stokes residual after nonlinear iteration 2: 0.0138303
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 19+0 iterations.
       Relative Stokes residual after nonlinear iteration 3: 0.000750098
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 11+0 iterations.
       Relative Stokes residual after nonlinear iteration 4: 2.71234e-05
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 5+0 iterations.
       Relative Stokes residual after nonlinear iteration 5: 6.8004e-07
 

tests/hollow_sphere/screen-output

@@ -8,8 +8,9 @@ Number of degrees of freedom: 28,690 (20,790+970+6,930)
    Rebuilding Stokes preconditioner...
    Solving Stokes system... 0+46 iterations.
       Relative Stokes residual after nonlinear iteration 1: 1
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 0+0 iterations.
-      Relative Stokes residual after nonlinear iteration 2: 4.31947e-13
+      Relative Stokes residual after nonlinear iteration 2: 4.28149e-13
 
    Postprocessing:
      Writing graphical output:              output-hollow_sphere/solution/solution-00000

tests/material_model_dependencies_burstedde/screen-output

@@ -8,8 +8,9 @@ Number of degrees of freedom: 3,041 (2,187+125+729)
    Rebuilding Stokes preconditioner...
    Solving Stokes system... 0+54 iterations.
       Relative Stokes residual after nonlinear iteration 1: 1
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 0+0 iterations.
-      Relative Stokes residual after nonlinear iteration 2: 7.02232e-13
+      Relative Stokes residual after nonlinear iteration 2: 4.81093e-13
 
    Postprocessing:
 viscosity depends on

tests/no_adiabatic_heating/screen-output

@@ -8,12 +8,16 @@ Number of degrees of freedom: 3,556 (2,178+289+1,089)
    Rebuilding Stokes preconditioner...
    Solving Stokes system... 32+0 iterations.
       Relative Stokes residual after nonlinear iteration 1: 1
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 25+0 iterations.
       Relative Stokes residual after nonlinear iteration 2: 0.0245986
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 18+0 iterations.
       Relative Stokes residual after nonlinear iteration 3: 0.00132971
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 12+0 iterations.
       Relative Stokes residual after nonlinear iteration 4: 4.7718e-05
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 5+0 iterations.
       Relative Stokes residual after nonlinear iteration 5: 1.17448e-06
 
@@ -190,6 +194,7 @@ Number of degrees of freedom: 268 (162+25+81)
    Rebuilding Stokes preconditioner...
    Solving Stokes system... 3+0 iterations.
       Relative Stokes residual after nonlinear iteration 1: 2.98496e-05
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 3+0 iterations.
       Relative Stokes residual after nonlinear iteration 2: 1.35729e-07
 

tests/no_adiabatic_heating_02/screen-output

@@ -8,12 +8,16 @@ Number of degrees of freedom: 3,556 (2,178+289+1,089)
    Rebuilding Stokes preconditioner...
    Solving Stokes system... 32+0 iterations.
       Relative Stokes residual after nonlinear iteration 1: 1
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 25+0 iterations.
       Relative Stokes residual after nonlinear iteration 2: 0.0245986
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 18+0 iterations.
       Relative Stokes residual after nonlinear iteration 3: 0.00132971
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 12+0 iterations.
       Relative Stokes residual after nonlinear iteration 4: 4.7718e-05
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 5+0 iterations.
       Relative Stokes residual after nonlinear iteration 5: 1.17448e-06
 
@@ -226,6 +230,7 @@ Number of degrees of freedom: 268 (162+25+81)
    Rebuilding Stokes preconditioner...
    Solving Stokes system... 3+0 iterations.
       Relative Stokes residual after nonlinear iteration 1: 2.98496e-05
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 3+0 iterations.
       Relative Stokes residual after nonlinear iteration 2: 1.35729e-07
 

tests/simple_compressibility_iterated_stokes/screen-output

@@ -8,12 +8,16 @@ Number of degrees of freedom: 13,764 (8,450+1,089+4,225)
    Rebuilding Stokes preconditioner...
    Solving Stokes system... 33+0 iterations.
       Relative Stokes residual after nonlinear iteration 1: 1
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 28+0 iterations.
       Relative Stokes residual after nonlinear iteration 2: 0.0138303
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 19+0 iterations.
       Relative Stokes residual after nonlinear iteration 3: 0.000750098
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 11+0 iterations.
       Relative Stokes residual after nonlinear iteration 4: 2.71234e-05
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 5+0 iterations.
       Relative Stokes residual after nonlinear iteration 5: 6.8004e-07
 

tests/tangurnis_tala_implicit/screen-output

@@ -8,6 +8,7 @@ Number of degrees of freedom: 3,556 (2,178+289+1,089)
    Rebuilding Stokes preconditioner...
    Solving Stokes system... 38+0 iterations.
       Relative Stokes residual after nonlinear iteration 1: 1
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 3+0 iterations.
       Relative Stokes residual after nonlinear iteration 2: 7.84568e-08
 
@@ -17,6 +18,7 @@ Number of degrees of freedom: 3,556 (2,178+289+1,089)
      Heating rate (average/total):  4.234e-05 W/kg, 5.493e-05 W
 
 *** Timestep 1:  t=0.0001 seconds
+   Rebuilding Stokes preconditioner...
    Solving Stokes system... 0+0 iterations.
       Relative Stokes residual after nonlinear iteration 1: 5.0054e-08