Convergence error during simulation of plate thickening process

ASPECT
1

Hello everyone,

I’m currently struggling to learn this software and have encountered several issues. I’m aiming to create a model where the surface experiences compressive velocities, ultimately leading to the formation of a mountain belt. Unfortunately, my model isn’t running as expected.

In previous attempts, I adjusted parameters such as convergence velocity, adiabatic surface temperature, and refinement function, finding that the model could run at lower velocities or higher grid heights. However, this increased costs and produced unsatisfactory results. After balancing and simplifying, I arrived at the following .prm file, but the model stops with an error at step 52. I can’t pinpoint the exact cause. My suspected reasons are: 1. Unreasonable temperature settings, 2. Unrealistic rheological properties for each layer of material, 3. Improper velocity function settings, 4. Too low grid resolution leading to iteration errors, 5. Plugin errors. I think the problem lies more in the initial setting of the component field after seeing the result of step 51. I would greatly appreciate any help in resolving this problem that has been bothering me for a long time. My files are attached here.

base.prm (31.7 KB)

*** Timestep 51:  t=231452 years, dt=920.169 years
   Calling FastScape... 5 timesteps of 184.034 years.
      FastScape runtime... 3.708s
   Solving mesh velocity system... 7 iterations.
   Solving temperature system... 16 iterations.
   Solving sediment_age system ... 21 iterations.
   Solving noninitial_plastic_strain system ... 20 iterations.
   Solving plastic_strain system ... 19 iterations.
   Solving viscous_strain system ... 20 iterations.
   Solving sediment_1 system ... 21 iterations.
   Skipping sediment_2 composition solve because RHS is zero.
   Solving upper system ... 20 iterations.
   Solving lower system ... 21 iterations.
   Solving mantle_upper system ... 21 iterations.
   Skipping mantle_middle composition solve because RHS is zero.
   Skipping mantle_lower composition solve because RHS is zero.
   Solving water_content system ... 19 iterations.
   Skipping max_melt composition solve because RHS is zero.
   Skipping melt_record composition solve because RHS is zero.
   Skipping melt_fraction composition solve because RHS is zero.
   Solving rift_grain_size system ... 19 iterations.
   Skipping new_crust composition solve because RHS is zero.
   Initial Newton Stokes residual = 1.13293e+16, v = 1.13292e+16, p = 5.68073e+13

   Rebuilding Stokes preconditioner...
   Solving Stokes system... 0+11 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 1: 0.00336951, norm of the rhs: 3.81742e+13

   Rebuilding Stokes preconditioner...
   Solving Stokes system... 0+19 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 2: 0.00129852, norm of the rhs: 1.47114e+13, newton_derivative_scaling_factor: 0

   The linear solver tolerance is set to 0.1. 
   Rebuilding Stokes preconditioner...
   Solving Stokes system... 0+1 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 3: 0.000444685, norm of the rhs: 5.03797e+12, newton_derivative_scaling_factor: 0


   Postprocessing:
     RMS, max velocity:                  0.0657 m/year, 0.144 m/year
     Temperature min/avg/max:            273 K, 1608 K, 1852 K
     Heat fluxes through boundary parts: -1.887e+06 W, -1.9e+06 W, 4.365e+06 W, 1.119e+06 W
     Mass fluxes through boundary parts: -3.099e+07 kg/yr, -3.092e+07 kg/yr, 6.199e+07 kg/yr, 5.201e+06 kg/yr
     Compositions min/max/mass:          -0.06828/0.2315/5.752e+08 // -0.001584/0.03256/2.928e+08 // -0.0005882/0.4315/6.315e+09 // -0.007389/0.07323/2.432e+10 // -0.2959/1/2.756e+09 // 0/0/0 // -0.1392/1.296/2.509e+10 // -0.09803/1.107/1.811e+10 // -0.1072/1.086/9.838e+10 // 0/0/0 // 0/0/0 // 0/62.6/3.955e+13 // 0/0/0 // 0/0/0 // 0/0/0 // 0/0.001252/7.909e+08 // 0/0/0
     Topography min/max:                 -709 m, 1187 m
     Writing stokes residuals            20250211-wy-cut-ai/stokes_residuals.txt

*** Timestep 52:  t=236452 years, dt=5000 years
   Calling FastScape... 5 timesteps of 1000 years.
      FastScape runtime... 3.703s
   Solving mesh velocity system... 7 iterations.
   Solving temperature system... retrying linear solve with different preconditioner...
---------------------------------------------------------
TimerOutput objects finalize timed values printed to the
screen by communicating over MPI in their destructors.
Since an exception is currently uncaught, this
synchronization (and subsequent output) will be skipped
to avoid a possible deadlock.
---------------------------------------------------------


----------------------------------------------------
Exception 'ExcMessage ("The " + solver_name + " did not converge. It reported the following error:\n\n" + exc.what() + "\n The required residual for convergence is: " + std::to_string(solver_controls.front().tolerance()) + ".\n See " + output_filename + " for convergence history.")' on rank 0 on processing: 

--------------------------------------------------------
An error occurred in line <405> of file </fs2/home/liuzhonglan/wy/aspect-model-wy/source/simulator/solver.cc> in function
    void aspect::{anonymous}::linear_solver_failed(const string&, const string&, const std::vector<dealii::SolverControl>&, const std::exception&)
The violated condition was: 
    false
Additional information: 
    The iterative advection solver did not converge. It reported the following error:


--------------------------------------------------------
An error occurred in line <1092> of file </fs2/home/liuzhonglan/wy/dealii-9.2/include/deal.II/lac/solver_gmres.h> in function
    void dealii::SolverGMRES<VectorType>::solve(const MatrixType&, VectorType&, const VectorType&, const PreconditionerType&) [with MatrixType = dealii::TrilinosWrappers::SparseMatrix; PreconditionerType = dealii::TrilinosWrappers::PreconditionILU; VectorType = dealii::TrilinosWrappers::MPI::Vector]
The violated condition was: 
    iteration_state == SolverControl::success
Additional information: 
Iterative method reported convergence failure in step 1000. The residual in the last step was 1.19688e+18.

This error message can indicate that you have simply not allowed a sufficiently large number of iterations for your iterative solver to converge. This often happens when you increase the size of your problem. In such cases, the last residual will likely still be very small, and you can make the error go away by increasing the allowed number of iterations when setting up the SolverControl object that determines the maximal number of iterations you allow.

The other situation where this error may occur is when your matrix is not invertible (e.g., your matrix has a null-space), or if you try to apply the wrong solver to a matrix (e.g., using CG for a matrix that is not symmetric or not positive definite). In these cases, the residual in the last iteration is likely going to be large.
--------------------------------------------------------

 The required residual for convergence is: 102300091.679163.
 See 20250211-wy-cut-ai/solver_history.txt for convergence history.
--------------------------------------------------------

Aborting!
----------------------------------------------------
yhrun: error: cn5184: task 0: Exited with exit code 1
slurmstepd: error:  mpi/pmix_v3: _errhandler: cn5184 [0]: pmixp_client_v2.c:212: Error handler invoked: status = -25, source = [slurm.pmix.4860408.0:0]
slurmstepd: error: *** STEP 4860408.0 ON cn5184 CANCELLED AT 2025-02-11T17:10:27 ***
yhrun: Job step aborted: Waiting up to 32 seconds for job step to finish.
yhrun: error: cn5188: tasks 224-279: Killed
yhrun: error: cn5186: tasks 112-167: Killed
yhrun: error: cn5187: tasks 168-171,173-223: Killed
yhrun: error: cn5187: task 172: Killed
yhrun: error: cn5185: tasks 56-82,84-111: Killed
yhrun: error: cn5185: task 83: Killed
yhrun: error: cn5184: tasks 1-55: Killed

velocity field
velocity field1705×1146 234 KB

The model result I want to achieve：
image
image875×314 48.2 KB

This is a model that previously produced good results but had incorrect parameters.

Additionally, as a beginner, are there any simpler files I can use for learning purposes? Thank you very much again.

Best,
Noriel Wang

2

Hi Noriel,

Welcome and thank you for posting to the forum!

The Advection solver error you are getting is quite early on in the computation and due to the complexity of the PRM it’s difficult to know exactly what the source of the issue is.

As you noted, it may be easier to start with a simpler model and then build up in complexity.

Here are a few suggestions for models to begin with:

  1. Continental Extension cookbook (note that I am hoping to finish updating this sometime in the next month or so)
  2. A recent paper on rift inversion: Vasey et al. (2024)

This is a model that previously produced good results but had incorrect parameters.

Was that a model run with ASPECT or a different code? If you produced this model with ASPECT, what changes were made that lead to the current issues with the Advection solver?

Cheers,
John

3

Thank you for your response!

I have completed the numerical simulation of the Continental Extension Cookbook and found that using simple models is indeed a very practical approach, especially for beginners. I inherited the ASPECT plugins and code from my predecessors, which has presented some challenges. I greatly appreciate your kind guidance and advice. Apart from the user manual, are there any other resources you would recommend for learning ASPECT, such as video tutorials?

Cheers！
Noriel

4

Hi Noriel,

I inherited the ASPECT plugins and code from my predecessors, which has presented some challenges

I can imagine, I think that is a situation everyone goes through at least once in there careers.

Apart from the user manual, are there any other resources you would recommend for learning ASPECT, such as video tutorials?

There are the ASPECT tectonics modeling tutorial videos on youtube, but I’m not sure how much help that would be if you have already gone through (and feel comfortable with) the relevant manual and cookbooks.

I think the best approach would be to go through relevant recent papers in ASPECT (and other codes) and use those to begin building on your predecessors work. On that note, here is another excellent recent paper using ASPECT to study on orogenic systems and coupled surface processes. that may be of use.

Do you by chance know the version of ASPECT, deal.II, and the related dependencies that your predecessors used to run their models? This information can be found in the log.txt files. If you want to rerun their exact models, you could build those versions of the packages (or use an older ASPECT docker container) and hopefully get the same results.

But again, a more fruitful path may be to using models from the more recent papers and systematically varying them towards a design that is aligned with your project goals.

Cheers,
John

5

Dear John,

Thank you so much for your helpful advice. Your insights on exploring recent papers and checking log files for version information were extremely valuable.

As a student, I greatly appreciate your willingness to share your expertise.

Warm regards, :blush:

Noriel