geodynamics / aspect

A parallel, extensible finite element code to simulate convection in both 2D and 3D models.
https://aspect.geodynamics.org/
Other
227 stars 237 forks source link

GMG with DC Picard diverges #4563

Open MFraters opened 2 years ago

MFraters commented 2 years ago

I have run into an issue with the GMG and DC Picard (and by extension the Newton solver). In short, single Advection, iterated Stokes with gmg works, single Advection, iterated defect correction Stokes without gmg works, but single Advection, iterated defect correction Stokes doesn't work (see output below).

I have tried to simplify the problem as much as I could and I have tried to compare the normal stokes gmg solver and assemblers with what happens for the Newton solver, but I couldn't find any obvious reasons why it would go wrong. I also simplified the boundary conditions and the compositional fields. The only way I could get the DC Picard + GMG to converge was if I simplified the temperature structure as well such that there wasn't really a nonlinear problem left.

Here are the aspect .prm file and the world builder .wb file. I haven't been able to reproduce the problem without a world builder temperature field, because with the function parser they either become way to easy or way to hard for the solver. This does require the current main branch of the world builder to be linked to aspect.

@tjhei @zjiaqi2018 any ideas on what could cause this issue?

single Advection, iterated Stokes with gmg:

-----------------------------------------------------------------------------
-- This is ASPECT, the Advanced Solver for Problems in Earth's ConvecTion.
--     . version 2.4.0-pre (main, c0287fce1)
--     . using deal.II 9.3.3
--     .       with 64 bit indices and vectorization level 2 (256 bits)
--     . using Trilinos 12.18.1
--     . using p4est 2.2.0
--     . running in OPTIMIZED mode
--     . running with 128 MPI processes
-----------------------------------------------------------------------------

Vectorization over 4 doubles = 256 bits (AVX), VECTORIZATION_LEVEL=2
-----------------------------------------------------------------------------
-- For information on how to cite ASPECT, see:
--   https://aspect.geodynamics.org/citing.html?ver=2.4.0-pre&GWB=1&mf=1&sha=c0287fce1&src=code
-----------------------------------------------------------------------------
Number of active cells: 8,192 (on 1 levels)
Number of degrees of freedom: 656,226 (215,475+9,801+71,825+71,825+71,825+71,825+71,825+71,825)

*** Timestep 0:  t=0 years, dt=0 years
   Solving temperature system... 0 iterations.
   Skipping NA_crust composition solve because RHS is zero.
   Skipping PC_crust composition solve because RHS is zero.
   Skipping spharz composition solve because RHS is zero.
   Skipping total_strain composition solve because RHS is zero.
   Skipping water composition solve because RHS is zero.
   Solving Stokes system... 1492+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 1: 0.0293025

   Solving Stokes system... 2523+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 2: 0.0127178

   Solving Stokes system... 2584+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 3: 0.000660424

   Solving Stokes system... 2989+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 4: 0.000254496

   Solving Stokes system... 2092+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 5: 8.86275e-05

   Solving Stokes system... 2266+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 6: 2.81538e-05

single Advection, iterated defect correction Stokes, with gmg:

-----------------------------------------------------------------------------
-- This is ASPECT, the Advanced Solver for Problems in Earth's ConvecTion.
--     . version 2.4.0-pre (main, c0287fce1)
--     . using deal.II 9.3.3
--     .       with 64 bit indices and vectorization level 2 (256 bits)
--     . using Trilinos 12.18.1
--     . using p4est 2.2.0
--     . running in OPTIMIZED mode
--     . running with 128 MPI processes
-----------------------------------------------------------------------------

Vectorization over 4 doubles = 256 bits (AVX), VECTORIZATION_LEVEL=2
-----------------------------------------------------------------------------
-- For information on how to cite ASPECT, see:
--   https://aspect.geodynamics.org/citing.html?ver=2.4.0-pre&GWB=1&NewtonSolver=1&mf=1&sha=c0287fce1&src=code
-----------------------------------------------------------------------------
Number of active cells: 8,192 (on 1 levels)
Number of degrees of freedom: 656,226 (215,475+9,801+71,825+71,825+71,825+71,825+71,825+71,825)

*** Timestep 0:  t=0 years, dt=0 years
   Solving temperature system... 0 iterations.
   Skipping NA_crust composition solve because RHS is zero.
   Skipping PC_crust composition solve because RHS is zero.
   Skipping spharz composition solve because RHS is zero.
   Skipping total_strain composition solve because RHS is zero.
   Skipping water composition solve because RHS is zero.
   Initial Newton Stokes residual = 2.12495e+19, v = 2.12495e+19, p = 0

   Solving Stokes system... 1492+0 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 1: 1, norm of the rhs: 2.12495e+19

   Solving Stokes system... 2861+0 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 2: 1.0865, norm of the rhs: 2.30875e+19, newton_derivative_scaling_factor: 0

   The linear solver tolerance is set to 0.001. 
   Solving Stokes system... 708+0 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 3: 3.36531, norm of the rhs: 7.15111e+19, newton_derivative_scaling_factor: 0

   The linear solver tolerance is set to 0.001. 
   Solving Stokes system... 378+0 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 4: 5.78038, norm of the rhs: 1.2283e+20, newton_derivative_scaling_factor: 0

   The linear solver tolerance is set to 0.001. 
   Solving Stokes system... 994+0 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 5: 6.55888, norm of the rhs: 1.39373e+20, newton_derivative_scaling_factor: 0

   The linear solver tolerance is set to 0.001. 
   Solving Stokes system... 564+0 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 6: 8.29118, norm of the rhs: 1.76183e+20, newton_derivative_scaling_factor: 0

   The linear solver tolerance is set to 0.001. 

single Advection, iterated defect correction Stokes, with gmg:

-----------------------------------------------------------------------------
-- This is ASPECT, the Advanced Solver for Problems in Earth's ConvecTion.
--     . version 2.4.0-pre (main, c0287fce1)
--     . using deal.II 9.3.3
--     .       with 64 bit indices and vectorization level 2 (256 bits)
--     . using Trilinos 12.18.1
--     . using p4est 2.2.0
--     . running in OPTIMIZED mode
--     . running with 128 MPI processes
-----------------------------------------------------------------------------

Vectorization over 4 doubles = 256 bits (AVX), VECTORIZATION_LEVEL=2
-----------------------------------------------------------------------------
-- For information on how to cite ASPECT, see:
--   https://aspect.geodynamics.org/citing.html?ver=2.4.0-pre&GWB=1&NewtonSolver=1&mf=1&sha=c0287fce1&src=code
-----------------------------------------------------------------------------
Number of active cells: 8,192 (on 1 levels)
Number of degrees of freedom: 656,226 (215,475+9,801+71,825+71,825+71,825+71,825+71,825+71,825)

*** Timestep 0:  t=0 years, dt=0 years
   Solving temperature system... 0 iterations.
   Skipping NA_crust composition solve because RHS is zero.
   Skipping PC_crust composition solve because RHS is zero.
   Skipping spharz composition solve because RHS is zero.
   Skipping total_strain composition solve because RHS is zero.
   Skipping water composition solve because RHS is zero.
   Initial Newton Stokes residual = 2.12495e+19, v = 2.12495e+19, p = 0

   Solving Stokes system... 1492+0 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 1: 1, norm of the rhs: 2.12495e+19

   Solving Stokes system... 2861+0 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 2: 1.0865, norm of the rhs: 2.30875e+19, newton_derivative_scaling_factor: 0

   The linear solver tolerance is set to 0.001. 
   Solving Stokes system... 708+0 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 3: 3.36531, norm of the rhs: 7.15111e+19, newton_derivative_scaling_factor: 0

   The linear solver tolerance is set to 0.001. 
   Solving Stokes system... 378+0 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 4: 5.78038, norm of the rhs: 1.2283e+20, newton_derivative_scaling_factor: 0

   The linear solver tolerance is set to 0.001. 
   Solving Stokes system... 994+0 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 5: 6.55888, norm of the rhs: 1.39373e+20, newton_derivative_scaling_factor: 0

   The linear solver tolerance is set to 0.001. 
   Solving Stokes system... 564+0 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 6: 8.29118, norm of the rhs: 1.76183e+20, newton_derivative_scaling_factor: 0

   The linear solver tolerance is set to 0.001. 

aspect prm:

set World builder file = world_builder_input.wb
set Output directory = ouput
set Dimension = 3
set CFL number              = 0.5
set Max nonlinear iterations = 75
set End time = 10e7
set Nonlinear solver scheme = single Advection, iterated Stokes 
#set Nonlinear solver scheme = single Advection, iterated defect correction Stokes 
set Timing output frequency = 0
set Max nonlinear iterations in pre-refinement = 0
set Pressure normalization                 = surface
set Nonlinear solver tolerance = 5e-15
set Resume computation = false
set Adiabatic surface temperature = 273.

subsection Solver parameters
  subsection Newton solver parameters
    set Max pre-Newton nonlinear iterations = 100
    set Nonlinear Newton solver switch tolerance = 1e-5
    set Max Newton line search iterations = 0
    set Maximum linear Stokes solver tolerance = 1e-3
    set Use Newton residual scaling method = true
    set Use Newton failsafe = false #true
    #set Stabilization preconditioner = none #SPD
    #set Stabilization velocity block = none #SPD
    set Use Eisenstat Walker method for Picard iterations = true
  end
  subsection Stokes solver parameters
    set Maximum number of expensive Stokes solver steps =  5000    
    set Number of cheap Stokes solver steps = 5000
    set Linear solver tolerance = 5e-7
    set Stokes solver type = block GMG
    set GMRES solver restart length = 100
  end
end

subsection Heating model
  set List of model names = adiabatic heating, shear heating 
end

subsection Boundary composition model
  set List of model names = initial composition
  set Fixed composition boundary indicators = 0,1,2,3,5,6,7,8,9
end

subsection Boundary temperature model
  set Fixed temperature boundary indicators   = 0,1,2,3,4,5,6,7,8,9
  set List of model names = initial temperature 
  subsection Initial temperature
    set Minimal temperature = 273.15
    set Maximal temperature = 4000
  end
end

subsection Discretization
  #set Use locally conservative discretization = true
  #set Use discontinuous composition discretization = true
  #set Use discontinuous temperature discretization = true
  subsection Stabilization parameters
    #set Stabilization method = SUPG
    set Use artificial viscosity smoothing = false 
    #set Use limiter for discontinuous composition solution = true
    #set Use limiter for discontinuous temperature solution = true
    set Global composition maximum = 1.0 #,1.0,1.0,1.0,1.0
    set Global composition minimum = 0.0 #,0.0,0.0,0.0,0.0
    set Global temperature minimum = 273.15
    set alpha = 2
    set beta  = 0.78 #0.117 
   end
end

subsection Boundary velocity model
  set Tangential velocity boundary indicators = #top #, lowernorth, lowereast, lowersouth, lowerwest
  set Zero velocity boundary indicators       = top, bottom, lowernorth, lowereast, lowersouth, lowerwest , uppereast, upperwest, uppersouth, uppernorth
end

subsection Boundary traction model #uppersouth: initial lithostatic pressure, upperwest: initial lithostatic pressure,
  #set Prescribed traction boundary indicators =  lowernorth: initial lithostatic pressure, lowereast: initial lithostatic pressure, lowersouth: initial lithostatic pressure, lowerwest: initial lithostatic pressure #, uppersouth: initial lithostatic pressure, uppernorth: initial lithostatic pressure 
  subsection Initial lithostatic pressure
    #set Representative point = 6371e3,-118.01,39.01 #6371e3,-107.01,33.01
    set Representative point = 6371e3,-106.01,32.01
  end
end

subsection Mesh deformation
  #set Mesh deformation boundary indicators = top: free surface
  subsection Free surface
    set Free surface stabilization theta = 1.0
    set Surface velocity projection = vertical
  end
end

subsection Geometry model
  set Model name = chunk with lithosphere boundary indicators

  subsection Chunk with lithosphere boundary indicators
    set Chunk inner radius = 5571e3
    set Chunk outer radius = 6371e3
    set Chunk middle boundary radius = 6271e3
    set Chunk minimum latitude = 32
    set Chunk maximum latitude = 58
    set Chunk minimum longitude = -138
    set Chunk maximum longitude = -106
    set Longitude repetitions = 32
    set Latitude repetitions = 32
    set Inner chunk radius repetitions = 7
  end
end

subsection Material model
  set Material averaging = project to Q1 only viscosity #harmonic average only viscosity  #harmonic average #project to Q1 only viscosity
  set Model name = visco plastic
  subsection Visco Plastic
    set Reference viscosity            = 1.0e21

    set Grain size = 0.01
    set Viscosity averaging scheme = harmonic 
    set Minimum viscosity         = 5e18
    set Maximum viscosity         = 1e24
    set Reference temperature = 293
    #                                              dry oli,  wet q (s),   wet ano (w), wet ano (w), 1e24
    #                                              mante,   weak crust, Pasific crust,  weak crust, weak lith, lower mantle
    # the density of the oceanic crust is usally taken as 3000, but the layer in this model is way too thick because of the resolution
    set Densities                                =     3300, 3100, 3300, 3300, 1e5, 1e5 #,     3300,      3300,       3000,       3200,     3300
    set Thermal expansivities                    =     2e-5
    #                                               Up M.,    NA crust, PC crust, weak crust, weak Lit, lower_mantle #,  water
    set Prefactors for dislocation creep          = 6.51e-15, 8.57e-28, 8.57e-28,    6.51e-16    , 1e-16  , 1e-16
    set Stress exponents for dislocation creep    =      3.5,      4.0,      4.0,         3.5    , 1  , 1
    set Activation energies for dislocation creep =   530.e3,   223.e3,   223.e3,      530.e3    , 1  , 1
    set Activation volumes for dislocation creep  =    18e-6,    18e-6,    18e-6,       18e-6    , 1e-6 , 1e-6
    set Prefactors for diffusion creep            = 8.88e-15, 8.88e-15, 8.88e-15,    8.88e-15    , 1e-15  , 1e-15
    set Activation energies for diffusion creep   =    335e3,    375e3,     375e3,      355e3    , 1  , 1
    set Activation volumes for diffusion creep    =  6.0e-6,    6.0e-6,    6.0e-6,     6.0e-6    , 1e-6 , 1e-6
    set Constant viscosity prefactors             =      1.0,      1.,        0.001,       100.0 , 1.0  , 1.0

    #                                 Up M., NA C., CARBC, Wk C., Wk L.,  Lower Mantle
    set Angles of internal friction =    15,     15,    15,  0  ,0, 0
    set Cohesions                   = 20.e6,   20e6, 10.e4,  20.e6 ,1, 1

  end
end
subsection Compositional fields
  set Number of fields = 5
  set Names of fields =  NA_crust, PC_crust, spharz , total_strain, water
  set Types of fields = chemical composition, chemical composition, chemical composition, generic, generic
  #set Compositional field methods = particles, particles, particles, field 
  set List of normalized fields = 0,1,2
end

subsection Initial temperature model
    set Model name = world builder
end

subsection Initial composition model
    set Model name = function #world builder
  subsection Function
    set Coordinate system = spherical #depth
    set Variable names      = r,long,lat
    set Function expression = 0;0;0;0;0 #if(6371-r<200e3 && lat > 0.69813170079773,1,0);if(6371-r<300e3 && lat <0.69813170079773,1,0);if(6371-r>200e3 && 6371-r>500e3,1,0);0;0
  end
end

subsection Gravity model
  set Model name = radial constant
end

subsection Mesh refinement
  set Initial global refinement = 0
  set Initial adaptive refinement = 0
  set Time steps between mesh refinement = 5
  set Strategy = isosurfaces, thermal energy density #, minimum refinement function, maximum refinement function
   subsection Isosurfaces
                          #minref maxref  mintemp maxtemp
        set Isosurfaces = \
                          max,    max,   Temperature:   0 | 600; \ 
                          max,    max,   PC_crust: 0.5 | 2; \ 
                          max-1,  max-1, Temperature:   600 | 1300; \ 
                          max-1,  max-1, spharz: 0.5 | 2; \ 
                          max-2,  max-2,  Temperature:1300 | 1651; \
                          min,    min,   Temperature:1651 | 3000; \

   end
end

subsection Postprocess
  #set Run postprocessors on nonlinear iterations = true
  set List of postprocessors = visualization, boundary pressures, pressure statistics, velocity statistics, spherical velocity statistics, temperature statistics, composition statistics, load balance statistics, memory statistics 
  subsection Visualization
    set List of output variables = depth, density, viscosity, strain rate,stress second invariant, #, partition #, density #, thermal conductivity
    set Time between graphical output = 25000 #50000
    set Write higher order output = false #true
    set Interpolate output = false
    set Write in background thread = true
  end
end

world_builder_input.txt

tjhei commented 2 years ago

Can you post the output for "single Advection, iterated defect correction Stokes, without gmg" as well?

The only thing I can think of right now is that it has to do with boundary conditions. Are they somehow contributing to the nonlinearity here? I see that we update the constraints only in StokesMatrixFree::setup_dofs(), which is not called between nonlinear steps, obviously.

MFraters commented 2 years ago

Thanks for the reply. Ah, sorry. I accidentally posted single Advection, iterated defect correction Stokes, with gmg twice, luckely I still have the data from single Advection, iterated defect correction Stokes, without gmg. That case does converge.

single Advection, iterated defect correction Stokes, without gmg:

-----------------------------------------------------------------------------
-- This is ASPECT, the Advanced Solver for Problems in Earth's ConvecTion.
--     . version 2.4.0-pre (main, c0287fce1)
--     . using deal.II 9.3.3
--     .       with 64 bit indices and vectorization level 2 (256 bits)
--     . using Trilinos 12.18.1
--     . using p4est 2.2.0
--     . running in OPTIMIZED mode
--     . running with 128 MPI processes
-----------------------------------------------------------------------------

-----------------------------------------------------------------------------
-- For information on how to cite ASPECT, see:
--   https://aspect.geodynamics.org/citing.html?ver=2.4.0-pre&GWB=1&NewtonSolver=1&sha=c0287fce1&src=code
-----------------------------------------------------------------------------
Number of active cells: 8,192 (on 1 levels)
Number of degrees of freedom: 656,226 (215,475+9,801+71,825+71,825+71,825+71,825+71,825+71,825)

*** Timestep 0:  t=0 years, dt=0 years
   Solving temperature system... 0 iterations.
   Skipping NA_crust composition solve because RHS is zero.
   Skipping PC_crust composition solve because RHS is zero.
   Skipping spharz composition solve because RHS is zero.
   Skipping total_strain composition solve because RHS is zero.
   Skipping water composition solve because RHS is zero.
   Initial Newton Stokes residual = 2.12495e+19, v = 2.12495e+19, p = 0

   Rebuilding Stokes preconditioner...
   Solving Stokes system... 160+0 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 1: 1, norm of the rhs: 2.12495e+19

   Rebuilding Stokes preconditioner...
   Solving Stokes system... 121+0 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 2: 0.0214211, norm of the rhs: 4.55187e+17, newton_derivative_scaling_factor: 0

   The linear solver tolerance is set to 0.001. 
   Rebuilding Stokes preconditioner...
   Solving Stokes system... 57+0 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 3: 0.00619948, norm of the rhs: 1.31736e+17, newton_derivative_scaling_factor: 0

   The linear solver tolerance is set to 0.001. 
   Rebuilding Stokes preconditioner...
   Solving Stokes system... 55+0 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 4: 0.00181908, norm of the rhs: 3.86545e+16, newton_derivative_scaling_factor: 0

   The linear solver tolerance is set to 0.001. 
   Rebuilding Stokes preconditioner...
   Solving Stokes system... 55+0 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 5: 0.000703964, norm of the rhs: 1.49589e+16, newton_derivative_scaling_factor: 0

   The linear solver tolerance is set to 0.001. 
   Rebuilding Stokes preconditioner...
   Solving Stokes system... 54+0 iterations.
      Relative nonlinear residual (total Newton system) after nonlinear iteration 6: 0.000293089, norm of the rhs: 6.22799e+15, newton_derivative_scaling_factor: 0
MFraters commented 2 years ago

Can you post the output for "single Advection, iterated defect correction Stokes, without gmg" as well?

The only thing I can think of right now is that it has to do with boundary conditions. Are they somehow contributing to the nonlinearity here? I see that we update the constraints only in StokesMatrixFree::setup_dofs(), which is not called between nonlinear steps, obviously.

hmm, the boundary conditions in this simplified model are all zero velocity, so I don't see how that would be an issue. In the full model I have velocity and traction boundary condition, but I don't think they change with nonlinear iterations, only with advection (which is done online once before the stokes iterations in these schemes).

MFraters commented 2 years ago

I also looked at the boundaries in paraview and the boundaries do indeed seem to be zero velocity, but the velocities inside are way too large (used the output after 5 nonlinear iterations). That shows I think that it also not just a residual computation problem, but a actual problem with the internal velocities.

Doing a bit more testing on that and just looking at the difference for a full stokes solver with and without gmg gives me significant different velocities, which I find a bit worrying. It is a extremely low resolution problem, but still, the residual in both cases are very small.

single Advection, iterated Stokes without gmg:

   Solving Stokes system... 160+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 1: 1

   Rebuilding Stokes preconditioner...
   Solving Stokes system... 144+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 2: 0.362054

   Rebuilding Stokes preconditioner...
   Solving Stokes system... 111+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 3: 0.0214082

   Rebuilding Stokes preconditioner...
   Solving Stokes system... 97+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 4: 0.00620233

   Rebuilding Stokes preconditioner...
   Solving Stokes system... 90+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 5: 0.00181854

   Rebuilding Stokes preconditioner...
   Solving Stokes system... 81+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 6: 0.000703951

   Rebuilding Stokes preconditioner...
   Solving Stokes system... 73+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 7: 0.000293004

   Rebuilding Stokes preconditioner...
   Solving Stokes system... 65+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 8: 0.000125605

   Rebuilding Stokes preconditioner...
   Solving Stokes system... 57+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 9: 5.53917e-05

   Rebuilding Stokes preconditioner...
   Solving Stokes system... 46+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 10: 2.52867e-05

   Postprocessing:
     Writing graphical output:                       cascadia-v01.17.tiny-v01-snogmg/solution/solution-00000
     Pressure at top/bottom of domain:               -0.0001385 Pa, 2.516e+10 Pa
     Pressure min/avg/max:                           -9.944e+07 Pa, 1.203e+10 Pa, 2.516e+10 Pa
     RMS, max velocity:                              0.00114 m/year, 0.0184 m/year
     Radial RMS, tangential RMS, total RMS velocity: 0.000781 m/yr, 0.000834 m/yr, 0.00114 m/yr
     Temperature min/avg/max:                        274.5 K, 1685 K, 2509 K
     Compositions min/max/mass:                      0/0/0 // 0/0/0 // 0/0/0 // 0/0/0 // 0/0/0
     Cells per process min/max/avg:                  64/64/64
     System matrix memory consumption:               5.85 MB

single Advection, iterated Stokes with gmg:

   Solving Stokes system... 1492+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 1: 0.0293025

   Solving Stokes system... 2523+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 2: 0.0127178

   Solving Stokes system... 2584+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 3: 0.000660424

   Solving Stokes system... 2989+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 4: 0.000254496

   Solving Stokes system... 2092+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 5: 8.86275e-05

   Solving Stokes system... 2266+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 6: 2.81538e-05

   Solving Stokes system... 1596+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 7: 1.10608e-05

   Solving Stokes system... 1534+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 8: 4.84137e-06

   Solving Stokes system... 1120+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 9: 2.31725e-06

   Solving Stokes system... 722+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 10: 1.21692e-06

   Postprocessing:
     Writing graphical output:                       cascadia-v01.17.tiny-v01-sgmg/solution/solution-00000
     Pressure at top/bottom of domain:               -0.0001106 Pa, 2.515e+10 Pa
     Pressure min/avg/max:                           -1.369e+08 Pa, 1.203e+10 Pa, 2.515e+10 Pa
     RMS, max velocity:                              0.00203 m/year, 0.0501 m/year
     Radial RMS, tangential RMS, total RMS velocity: 0.00135 m/yr, 0.00152 m/yr, 0.00203 m/yr
     Temperature min/avg/max:                        274.5 K, 1685 K, 2509 K
     Compositions min/max/mass:                      0/0/0 // 0/0/0 // 0/0/0 // 0/0/0 // 0/0/0
     Cells per process min/max/avg:                  64/64/64
     System matrix memory consumption:               1.34 MB

The only thing I do here is turning the line with block GMG on and off. The RMS provided by the GMG solver is almost twice the value and the max velocity is more than twice while the residuals are very low. The values in both cases are almost the same for 5 nonlinear iterations, so they are not seem to be converging to each other.

tjhei commented 2 years ago

Something seems broken here. I also don't like that GMG requires 1000+ iterations. That can't be right.

tjhei commented 2 years ago

Can you share the prm for the last test that you ran?

anne-glerum commented 2 years ago

I'm not sure it's the same problem, but I just ran the gmg_mesh_def_ghost_nodes.prm test for #4248 and this also shows significant differences in RMS and max vel and iteration counts with and without block GMG. (It uses the normal Picard solver.)

MFraters commented 2 years ago

Can you share the prm for the last test that you ran?

Like I said above I use the current world builder main branch and the .wb file is above as well. The only change I make comment and uncomment the line with set Stokes solver type = block GMG

set World builder file = iofiles/cascadia-v01/cascadia-v01.17.wb
set Output directory = iofiles/cascadia-v01/cascadia-v01.17.tiny-v01-sgmg
set Dimension = 3
set CFL number              = 0.5
set Max nonlinear iterations = 10
set End time = 0
set Nonlinear solver scheme = single Advection, iterated Stokes 
#set Nonlinear solver scheme = single Advection, iterated defect correction Stokes 
set Timing output frequency = 0
set Max nonlinear iterations in pre-refinement = 0
set Pressure normalization                 = surface
set Nonlinear solver tolerance = 5e-15
set Resume computation = false
set Adiabatic surface temperature = 273.

subsection Solver parameters
  subsection Newton solver parameters
    set Max pre-Newton nonlinear iterations = 100
    set Nonlinear Newton solver switch tolerance = 1e-5
    set Max Newton line search iterations = 0
    set Maximum linear Stokes solver tolerance = 1e-4
    set Use Newton residual scaling method = true
    set Use Newton failsafe = false #true
    #set Stabilization preconditioner = none #SPD
    #set Stabilization velocity block = none #SPD
    set Use Eisenstat Walker method for Picard iterations = true
  end
  subsection Stokes solver parameters
    set Maximum number of expensive Stokes solver steps =  5000    
    set Number of cheap Stokes solver steps = 5000
    set Linear solver tolerance = 5e-7
    #set Stokes solver type = block GMG
    set GMRES solver restart length = 100
  end
end

subsection Heating model
  set List of model names = adiabatic heating, shear heating 
end

subsection Boundary composition model
  set List of model names = initial composition
  set Fixed composition boundary indicators = 0,1,2,3,5,6,7,8,9
end

subsection Boundary temperature model
  set Fixed temperature boundary indicators   = 0,1,2,3,4,5,6,7,8,9
  set List of model names = initial temperature 
  subsection Initial temperature
    set Minimal temperature = 273.15
    set Maximal temperature = 4000
  end
end

subsection Discretization
  #set Use locally conservative discretization = true
  #set Use discontinuous composition discretization = true
  #set Use discontinuous temperature discretization = true
  subsection Stabilization parameters
    #set Stabilization method = SUPG
    set Use artificial viscosity smoothing = false 
    #set Use limiter for discontinuous composition solution = true
    #set Use limiter for discontinuous temperature solution = true
    set Global composition maximum = 1.0 #,1.0,1.0,1.0,1.0
    set Global composition minimum = 0.0 #,0.0,0.0,0.0,0.0
    set Global temperature minimum = 273.15
    set alpha = 2
    set beta  = 0.78 #0.117 
   end
end

subsection Boundary velocity model
  set Tangential velocity boundary indicators = #top #, lowernorth, lowereast, lowersouth, lowerwest
  set Zero velocity boundary indicators       = top, bottom, lowernorth, lowereast, lowersouth, lowerwest , uppereast, upperwest, uppersouth, uppernorth
end

subsection Boundary traction model #uppersouth: initial lithostatic pressure, upperwest: initial lithostatic pressure,
  #set Prescribed traction boundary indicators =  lowernorth: initial lithostatic pressure, lowereast: initial lithostatic pressure, lowersouth: initial lithostatic pressure, lowerwest: initial lithostatic pressure #, uppersouth: initial lithostatic pressure, uppernorth: initial lithostatic pressure 
  subsection Initial lithostatic pressure
    #set Representative point = 6371e3,-118.01,39.01 #6371e3,-107.01,33.01
    set Representative point = 6371e3,-106.01,32.01
  end
end

subsection Mesh deformation
  #set Mesh deformation boundary indicators = top: free surface
  subsection Free surface
    set Free surface stabilization theta = 1.0
    set Surface velocity projection = vertical
  end
end

subsection Geometry model
  set Model name = chunk with lithosphere boundary indicators

  subsection Chunk with lithosphere boundary indicators
    set Chunk inner radius = 5571e3
    set Chunk outer radius = 6371e3
    set Chunk middle boundary radius = 6271e3
    set Chunk minimum latitude = 32
    set Chunk maximum latitude = 58
    set Chunk minimum longitude = -138
    set Chunk maximum longitude = -106
    set Longitude repetitions = 32
    set Latitude repetitions = 32
    set Inner chunk radius repetitions = 7
  end
end

subsection Material model
  set Material averaging = project to Q1 only viscosity #harmonic average only viscosity  #harmonic average #project to Q1 only viscosity
  set Model name = visco plastic
  subsection Visco Plastic
    set Reference viscosity            = 1.0e21

    set Grain size = 0.01
    set Viscosity averaging scheme = harmonic 
    set Minimum viscosity         = 5e18
    set Maximum viscosity         = 1e24
    set Reference temperature = 293
    #                                              dry oli,  wet q (s),   wet ano (w), wet ano (w), 1e24
    #                                              mante,   weak crust, Pasific crust,  weak crust, weak lith, lower mantle
    # the density of the oceanic crust is usally taken as 3000, but the layer in this model is way too thick because of the resolution
    set Densities                                =     3300, 3100, 3300, 3300, 1e5, 1e5 #,     3300,      3300,       3000,       3200,     3300
    set Thermal expansivities                    =     2e-5
    #                                               Up M.,    NA crust, PC crust, weak crust, weak Lit, lower_mantle #,  water
    set Prefactors for dislocation creep          = 6.51e-15, 8.57e-28, 8.57e-28,    6.51e-16    , 1e-16  , 1e-16
    set Stress exponents for dislocation creep    =      3.5,      4.0,      4.0,         3.5    , 1  , 1
    set Activation energies for dislocation creep =   530.e3,   223.e3,   223.e3,      530.e3    , 1  , 1
    set Activation volumes for dislocation creep  =    18e-6,    18e-6,    18e-6,       18e-6    , 1e-6 , 1e-6
    set Prefactors for diffusion creep            = 8.88e-15, 8.88e-15, 8.88e-15,    8.88e-15    , 1e-15  , 1e-15
    set Activation energies for diffusion creep   =    335e3,    375e3,     375e3,      355e3    , 1  , 1
    set Activation volumes for diffusion creep    =  6.0e-6,    6.0e-6,    6.0e-6,     6.0e-6    , 1e-6 , 1e-6
    set Constant viscosity prefactors             =      1.0,      1.,        0.001,       100.0 , 1.0  , 1.0

    #                                 Up M., NA C., CARBC, Wk C., Wk L.,  Lower Mantle
    set Angles of internal friction =    15,     15,    15,  0  ,0, 0
    set Cohesions                   = 20.e6,   20e6, 10.e4,  20.e6 ,1, 1

  end
end
subsection Compositional fields
  set Number of fields = 5
  set Names of fields =  NA_crust, PC_crust, spharz , total_strain, water
  set Types of fields = chemical composition, chemical composition, chemical composition, generic, generic
  #set Compositional field methods = particles, particles, particles, field 
  set List of normalized fields = 0,1,2
end

subsection Initial temperature model
    set Model name = world builder
end

subsection Initial composition model
    set Model name = function #world builder
  subsection Function
    set Coordinate system = spherical #depth
    set Variable names      = r,long,lat
    set Function expression = 0;0;0;0;0 #if(6371-r<200e3 && lat > 0.69813170079773,1,0);if(6371-r<300e3 && lat <0.69813170079773,1,0);if(6371-r>200e3 && 6371-r>500e3,1,0);0;0
  end
end

subsection Gravity model
  set Model name = radial constant
end

subsection Mesh refinement
  set Initial global refinement = 0
  set Initial adaptive refinement = 0
  set Time steps between mesh refinement = 5
  set Strategy = isosurfaces, thermal energy density #, minimum refinement function, maximum refinement function
   subsection Isosurfaces
                          #minref maxref  mintemp maxtemp
        set Isosurfaces = \
                          max,    max,   Temperature:   0 | 600; \ 
                          max,    max,   PC_crust: 0.5 | 2; \ 
                          max-1,  max-1, Temperature:   600 | 1300; \ 
                          max-1,  max-1, spharz: 0.5 | 2; \ 
                          max-2,  max-2,  Temperature:1300 | 1651; \
                          min,    min,   Temperature:1651 | 3000; \

   end
end

subsection Postprocess
  #set Run postprocessors on nonlinear iterations = true
  set List of postprocessors = visualization, boundary pressures, pressure statistics, velocity statistics, spherical velocity statistics, temperature statistics, composition statistics, load balance statistics, memory statistics 
  subsection Visualization
    set List of output variables = depth, density, viscosity, strain rate,stress second invariant, #, partition #, density #, thermal conductivity
    set Time between graphical output = 25000 #50000
    set Write higher order output = false #true
    set Interpolate output = false
    set Write in background thread = true
  end
end
MFraters commented 2 years ago

I'm not sure it's the same problem, but I just ran the gmg_mesh_def_ghost_nodes.prm test for #4248 and this also shows significant differences in RMS and max vel and iteration counts with and without block GMG. (It uses the normal Picard solver.)

Yes, that looks like the same, although the differences are not as large, but still look significant. What was your residual? In summary there seem to be two issues:

  1. Full Picard and dc Picard without gmg give different velocities as compared to full Picard with gmg.
  2. Dc Picard with gmg it doesn't converge for this model (some other models do converge, so it is not always an issue)
anne-glerum commented 2 years ago

Yes, that looks like the same, although the differences are not as large, but still look significant. What was your residual?

The test just does a single Stokes solve, the residual is not reported.

tjhei commented 2 years ago

@zjiaqi2018 Can you take a look at this issue as well? Let us know if you have any suggestions on what to try or what could be the problem here.

zjiaqi2018 commented 2 years ago

I cannot run the prm file:

----------------------------------------------------
Exception on rank 0 on processing: 
AssertThrow `false` failed in ../contrib/world_builder/source/parameters.cc at line 168: Invalid schema: #/properties/coordinate%20system
Invalid keyword: oneOfInvalid schema: #/coordinate%20system
Error document: 
#/coordinate%20system{
    "oneOf": {
        "errors": [
            {
                "enum": {
                    "instanceRef": "#/coordinate%20system/model",
                    "schemaRef": "#/properties/coordinate%20system/oneOf/0/properties/model"
                }
            },
            {
                "enum": {
                    "instanceRef": "#/coordinate%20system/depth%20method",
                    "schemaRef": "#/properties/coordinate%20system/oneOf/1/properties/depth%20method"
                }
            }
        ],
        "instanceRef": "#/coordinate%20system",
        "schemaRef": "#/properties/coordinate%20system"
    }
}

Aborting!
tjhei commented 2 years ago

Menno has been using a newer world builder version. @MFraters can you check if this is still broken with the current master?

elodie-kendall commented 2 years ago

Any news here? Can I help at all? Thanks!

MFraters commented 2 years ago

@tjhei Sorry for the late reply. Due to traveling I didn't get around to it. I check it and the results are still the same.

tjhei commented 2 years ago

@elodie-kendall What would help is having a set of simple(!) examples that currently work or don't work. With that, it would be much easier to find the issue.

tjhei commented 2 years ago

In particular: Is a free surface required to see the problem? Or is it related to boundary conditions?

tjhei commented 2 years ago

Hi Menno, we ran a couple of defect correction with GMG tests in the last few days and all of them seem to work okay. I might have fixed recently an issue with incorrect material properties being requested in the GMG solver.

Can you try the problem you have again? Otherwise, any suggestions on how to narrow in on the issue?

MFraters commented 1 year ago

I had tested it and it still didn't work. I have now taken the time to make a simpler example which runs on the current main branch. It still requires the world builder and a bit of resolution, but I took out everything non-standard and got the problem to reproduce on a lower resolution.

With the line set Stokes solver type = block GMG commented out it works and the results look fine. With the line enabled the gmg solver fails (see the last block).

Let me know if you have any ideas on what it could be or what further to test.

prm file:

set World builder file = iofiles/gmg/gmg-test-v4.wb
set Output directory = iofiles/gmg/gmg-test-gmg-v4
set Dimension = 3
set CFL number              = 0.5
set Max nonlinear iterations = 2 
set End time = 10e7
set Nonlinear solver scheme = single Advection, iterated Stokes 
set Timing output frequency = 0
set Max nonlinear iterations in pre-refinement = 0
set Pressure normalization                 = surface
set Nonlinear solver tolerance =1e-4
set Adiabatic surface temperature = 273.

subsection Solver parameters
  subsection Stokes solver parameters
    set Maximum number of expensive Stokes solver steps =  5000    
    set Number of cheap Stokes solver steps = 1000 #80000
    #set Linear solver tolerance = 1e-2
    set Stokes solver type = block GMG
    #set GMRES solver restart length = 100
  end
end

subsection Heating model
  set List of model names = adiabatic heating, shear heating 
end

subsection Boundary composition model
  set List of model names = initial composition
  set Fixed composition boundary indicators = 0,1,2,3,5,6,7,8,9
end

subsection Boundary temperature model
  set Fixed temperature boundary indicators   = 0,1,2,3,4,5,6,7,8,9
  set List of model names = initial temperature 
  subsection Initial temperature
    set Minimal temperature = 273.15
    set Maximal temperature = 4000
  end
end

subsection Discretization
  subsection Stabilization parameters
    set Use artificial viscosity smoothing = false 
    set Global composition maximum = 1.0 
    set Global composition minimum = 0.0 
    set Global temperature minimum = 273.15
    set alpha = 2
    set beta  = 0.78 #0.117 
   end
end

subsection Boundary velocity model
  set Prescribed velocity boundary indicators = uppereast: function, upperwest: function 
  set Tangential velocity boundary indicators = top
  set Zero velocity boundary indicators       = bottom
  subsection Function
    set Coordinate system = spherical
    set Use spherical unit vectors = true
    set Variable names = r,phi,theta,t
    set Function expression = 0; -0.01;0  
  end
end

subsection Boundary traction model 
  set Prescribed traction boundary indicators =  lowernorth: initial lithostatic pressure, lowereast: initial lithostatic pressure, lowersouth: initial lithostatic pressure, lowerwest: initial lithostatic pressure, uppersouth: initial lithostatic pressure, uppernorth: initial lithostatic pressure 
  subsection Initial lithostatic pressure
    set Representative point = 6371e3,-106.01,32.01
  end
end

subsection Geometry model
  set Model name = chunk with lithosphere boundary indicators

  subsection Chunk with lithosphere boundary indicators
    set Chunk inner radius = 5571e3
    set Chunk outer radius = 6371e3
    set Chunk middle boundary radius = 6271e3
    set Chunk minimum latitude = 32
    set Chunk maximum latitude = 58
    set Chunk minimum longitude = -138
    set Chunk maximum longitude = -106
    set Longitude repetitions = 32
    set Latitude repetitions = 32
    set Inner chunk radius repetitions = 7
  end
end

subsection Material model
  set Material averaging = harmonic average only viscosity  #project to Q1 only viscosity #harmonic average only viscosity  #harmonic average #project to Q1 only viscosity
  set Model name = visco plastic
  subsection Visco Plastic
        set Reference temperature = 273
        set Minimum viscosity = 1e19
        set Maximum viscosity = 1e24
        set Phase transition depths = background:80e3|410e3|520e3|560e3|670e3|670e3|670e3|670e3, PC_crust:80e3|665e3|720e3, spharz: 410e3|520e3|560e3|670e3|670e3|670e3|670e3
        set Phase transition widths = background:5e3|5e3|5e3|5e3|10e3|5e3|5e3|5e3, PC_crust:5e3|5e3|5e3, spharz: 5e3|5e3|5e3|5e3|5e3|5e3|5e3
        set Phase transition temperatures = background:1662.0|1662.0|1662.0|1662.0|1662.0|1662.0|1662.0|1662.0, PC_crust:1173.0|1662.0|1662.0, spharz: 1662.0|1662.0|1662.0|1662.0|1662.0|1662.0|1662.0
        set Phase transition Clapeyron slopes = background:0.0|4e6|4.1e6|4e6|-2e6|4e6|-3.1e6|1.3e6, PC_crust:0.0|4e6|1.3e6, spharz: 4e6|4.1e6|4e6|-2e6|4e6|-3.1e6|1.3e6
        set Thermal diffusivities = 1.0e-6
        set Heat capacities = 1250.0
        set Densities = background: 3300.0|3300.0|3394.4|3442.1|3453.2|3617.6|3691.5|3774.7|3929.1, NA_crust:3000.0, PC_crust:3000.0|3540.0|3613.0|3871.7, spharz: 3235.0|3372.3|3441.7|3441.7|3680.8|3717.8|3759.4|3836.6 #,  total_strain:1e5, water:1e5
        set Thermal expansivities = 3.1e-5
        set Viscosity averaging scheme = harmonic
        set Viscous flow law = composite
        set Yield mechanism = drucker
        set Grain size = 1.0000e-02
        set Prefactors for diffusion creep =            background:8.5700e-28|7.1768e-15|7.1768e-15|7.1768e-15|7.1768e-15|2.4977e-19|2.4977e-19|2.4977e-19|2.4977e-19, NA_crust:8.5700e-28, PC_crust:8.5700e-28|7.1768e-15|2.4977e-19|2.4977e-19, spharz:7.1768e-15|7.1768e-15|7.1768e-15|7.1768e-15|2.4977e-19|2.4977e-19|2.4977e-19|2.4977e-19 #, total_strain:1e-15, water: 1e-15
        set Grain size exponents for diffusion creep =  background:3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00, NA_crust:3.0000e+00, PC_crust:3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00, spharz:3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00 #, total_strain:1, water:1
        set Activation energies for diffusion creep =   background:3.7500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05, NA_crust:3.7500e+05, PC_crust:3.7500e+05|3.2500e+05|3.2500e+05|3.2500e+05, spharz:3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05 #, total_strain:1, water:1
        set Activation volumes for diffusion creep =    background:2.3000e-05|2.3000e-05|2.3000e-05|2.3000e-05|2.3000e-05|3.0000e-06|3.0000e-06|3.0000e-06|3.0000e-06, NA_crust:2.3000e-05, PC_crust:2.3000e-05|2.3000e-05|3.0000e-06|3.0000e-06, spharz:2.3000e-05|2.3000e-05|2.3000e-05|2.3000e-05|3.0000e-06|3.0000e-06|3.0000e-06|3.0000e-06 #, total_strain:1, water:1
        set Prefactors for dislocation creep =          background:8.5700e-28|4.4668e-16|4.4668e-16|4.4668e-16|4.4668e-16|5.0000e-32|5.0000e-32|5.0000e-32|5.0000e-32, NA_crust:8.5700e-28, PC_crust:8.5700e-28|4.4668e-16|5.0000e-32|5.0000e-32, spharz:4.4668e-16|4.4668e-16|4.4668e-16|4.4668e-16|5.0000e-32|5.0000e-32|5.0000e-32|5.0000e-32 #, total_strain:1e-15, water:1e-15
        set Stress exponents for dislocation creep =    background:4.0000e+00|3.5000e+00|3.5000e+00|3.5000e+00|3.5000e+00|1.0000e+00|1.0000e+00|1.0000e+00|1.0000e+00, NA_crust:3.5000e+00, PC_crust:4.0000e+00|3.5000e+00|1.0000e+00|1.0000e+00, spharz:3.5000e+00|3.5000e+00|3.5000e+00|3.5000e+00|1.0000e+00|1.0000e+00|1.0000e+00|1.0000e+00 #, total_strain:1, water:1
        set Activation energies for dislocation creep = background:2.2300e+05|4.8000e+05|4.8000e+05|4.8000e+05|4.8000e+05|0.0000e+00|0.0000e+00|0.0000e+00|0.0000e+00, NA_crust:2.2300e+05, PC_crust:2.2300e+05|4.8000e+05|0.0000e+00|0.0000e+00, spharz:4.8000e+05|4.8000e+05|4.8000e+05|4.8000e+05|0.0000e+00|0.0000e+00|0.0000e+00|0.0000e+00 #, total_strain:1, water:1
        set Activation volumes for dislocation creep =  background:2.4000e-05|2.4000e-05|2.4000e-05|2.4000e-05|2.4000e-05|0.0000e+00|0.0000e+00|0.0000e+00|0.0000e+00, NA_crust:2.4000e-05, PC_crust:2.4000e-05|2.4000e-05|0.0000e+00|0.0000e+00, spharz:2.4000e-05|2.4000e-05|2.4000e-05|2.4000e-05|0.0000e+00|0.0000e+00|0.0000e+00|0.0000e+00 #, total_strain:1, water:1

    set Angles of internal friction = background: 45|45|0|0|0|0|0|0|0, NA_crust:45,  PC_crust:45|45|0|0, spharz:45|45|0|0|0|0|0|0 #, total_strain:0, water:0
    set Cohesions                   = background: 20e8|20.e8|20e8|20e8|20e8|20e8|20e8|20e8|20e8, NA_crust:20e6, PC_crust:10.e4|20e6|20e8|20e8, spharz:20e6|20e8|20e8|20e8|20e8|20e8|20e8|20e8 #, total_strain:1, water:1

    set Prefactor strain weakening factors = 0.01
    set Friction strain weakening factors = 0.1
    set Cohesion strain weakening factors = 0.1
    set Start plasticity strain weakening intervals = 0.1
    set Start prefactor strain weakening intervals = 0.1
    set Strain healing temperature dependent recovery rate =  1e-20
    set Strain healing mechanism                     = temperature dependent
  end
end

subsection Compositional fields
  set Number of fields = 3
  set Names of fields =  NA_crust, PC_crust, spharz #, total_strain, water
  set Types of fields = chemical composition, chemical composition, chemical composition #, generic, generic
  set List of normalized fields = 0,1,2
end

subsection Initial temperature model
    set Model name = world builder
end

subsection Initial composition model
    set Model name = world builder
end

subsection Gravity model
  set Model name = radial constant
end

subsection Mesh refinement
  set Initial global refinement = 2 #3
  set Initial adaptive refinement = 2
  set Time steps between mesh refinement = 5
  set Strategy = isosurfaces, thermal energy density 
   subsection Isosurfaces
                          #minref maxref  mintemp maxtemp
        set Isosurfaces = \
                          max,    max,   Temperature:   0 | 350; \ 
                          max,    max,   Temperature:   0 | 1300, PC_crust: 0.5 | 2; \ 
                          max-1,  max-1, Temperature:   350 | 400; \ 
                          max-2,  max-2, Temperature:   400 | 1300; \ 
                          max-1,    max,   PC_crust: 0.5 | 2; \ 
                          max-1,  max-1, spharz: 0.5 | 2; \ 
                          max-2,  max-2,  Temperature:1300 | 1660; \
                          min,    min,   Temperature:1660 | 3000; \

   end
end

subsection Postprocess
  #set Run postprocessors on nonlinear iterations = true
  set List of postprocessors = visualization, boundary pressures, pressure statistics, velocity statistics, spherical velocity statistics, temperature statistics, composition statistics, load balance statistics, memory statistics 
  subsection Visualization
    set List of output variables = depth, density, viscosity, strain rate,stress second invariant
    set Time between graphical output = 12500 
    set Write higher order output = false 
    set Interpolate output = false
    set Write in background thread = true
  end
end

world builder file:

{
  // changed oceanic plate depth from 95e3 to 100e3 and changed oceanic plate velocity from 0.02 to 0.02
  "version":"0.4",
  "coordinate system":{"model":"spherical", "depth method":"begin segment"},
  //"coordinate system":{"model":"spherical", "depth method":"begin at end segment"},
  //"cross section":[[-128,50.5],[-113,46.5]], //[[-131,47],[-117,43]],
  "cross section":[[-131,47],[-117,43]],
  //"maximum distance between coordinates":0.01,
  "interpolation":"continuous monotone spline",
  "features":
  [
    {"model":"mantle layer", "name":"mantle", "min depth":100e3, "max depth":1200e3, "coordinates":[[-150,20],[-100,20],[-100,80],[-150,80]],
      "composition models":[{"model":"uniform", "min depth":100e3, "max depth":1200e3, "compositions":[4], "fractions":[400]}]
    },
    {"model":"continental plate", "name":"boundary", "min depth":-1e2, "max depth":120e3, "coordinates":[[-150,20],[-100,20],[-100,80],[-150,80]],
     "temperature models":[{"model":"linear", "min depth":-1e2, "max depth":100e3, "top temperature":273.15}]//,
      //"composition models":[{"model":"uniform", "min depth":-1e2, "max depth":15e3, "compositions":[0]}]
    },
    {"model":"continental plate", "name":"inner part", "min depth":-1e2, "max depth":120e3, "coordinates":
    [[-135,53], [-135,36],[-111,36],[-111,53]],
     "temperature models":[{"model":"adiabatic", "min depth":-1e2, "max depth":100e3}]
    },
    {"model":"oceanic plate", "name":"Pasific", "min depth":-1e2,"max depth":100e3, "coordinates":
     [[-119.53587316701834,35.026738764076697],[-119.57383732651698,35.056444555160567],[-119.60939570151095,35.084327995405658],[-119.63417894839137,35.106540426058871],
     [-119.65858565209328,35.125880316256769],[-119.71763954289975,35.171045968622082],[-119.74821603067983,35.197882826954526],[-119.76978766719469,35.218363290625462],
     [-119.79409829064866,35.241824187169698],[-119.82807176026733,35.272843231920774],[-119.86055950731412,35.304971252158566],[-119.90632396186282,35.351792300798877],
     [-119.94692681693385,35.39038972770112],[-119.97804209353507,35.42203524916863],[-119.99889777420293,35.446300194659898],[-120.03159767449893,35.483016977920613],
     [-120.06943818848163,35.521453648690397],[-120.1121292412505,35.563296615814409],[-120.14862523818482,35.59870591463266],[-120.18506317056102,35.633728585922995],
     [-120.20593586300009,35.653689383733365],[-120.22783449400163,35.674173828857363],[-120.24859639654591,35.695779235772477],[-120.28108807142428,35.726943855131026],
     [-120.30506054229693,35.754745248739425],[-120.32579703346602,35.788598829563568],[-120.35500781938225,35.825358236558714],[-120.38094720203441,35.849557704201459],
     [-120.42164544113507,35.884236198652957],[-120.44946117595816,35.910487135624123],[-120.47477761013653,35.933559258904154],[-120.51162105275768,35.964166488899593],
     [-120.54669696564423,35.99339198982949],[-120.58260325960816,36.022430244199256],[-120.61725375997969,36.05478276694754],[-120.65117422382883,36.085785007481036],
     [-120.67630445821555,36.10849694324088],[-120.75261850319879,36.176845516870685],[-120.78691544416824,36.207853475097693],[-120.82941491103179,36.246844661554348],
     [-120.86409649358581,36.278691639242709],[-120.89398206909289,36.309107397971445],[-120.92541506386243,36.33833146896319],[-120.95735784000937,36.369754311823044],
     [-121.01576174330518,36.431157243479788],[-121.07207309013523,36.485893948235002],[-121.10791832828488,36.519877354310154],[-121.14071876667788,36.554153106004833],
     [-121.16371695120017,36.575907615284507],[-121.18670053002239,36.59376134384371],[-121.23916541841322,36.6393384352711],[-121.26368275517387,36.66066564742988],
     [-121.2926215875213,36.683474164488587],[-121.33607120903258,36.715009966943455],[-121.36875476511875,36.739378826299514],[-121.40886245141422,36.764290491092368],
     [-121.45673589161618,36.793529056408488],[-121.49115617325845,36.815708963134853],[-121.52610334369285,36.840114879061275],[-121.55508229475919,36.86255855283656],
     [-121.58944549593645,36.884433460794298],[-121.61260729844474,36.901346802870705],[-121.64048543975076,36.919630706516728],[-121.67950895795178,36.946478562229061],
     [-121.74247097000728,37.001658514256235],[-121.76465348371033,37.020752188503081],[-121.80503194745722,37.049836399722381],[-121.84121748147328,37.072019692455228],
     [-121.88399134347702,37.094449022037679],[-121.92357308919088,37.118855092475712],[-121.95555511187678,37.141996698957144],[-121.98598162856348,37.167992306611893],
     [-122.03105205620977,37.204234039470748],[-122.06084978266296,37.230940161988258],[-122.09556558307531,37.26296665270047],[-122.13081161527418,37.290741949604353],
     [-122.18537362833922,37.335571261058476],[-122.21486053342659,37.365971508253097],[-122.24885806811528,37.404279337111973],[-122.27151125838043,37.428654657266122],
     [-122.2964121118514,37.453237354638418],[-122.33354466861408,37.49091562348093],[-122.38065642027948,37.545925782427048],[-122.45281742259908,37.631268840719542],
     [-122.49843602016591,37.677498177477901],[-122.54652300999368,37.74016580820421],[-122.61834876039984,37.832996147222048],[-122.68624318948224,37.933754160924821],
     [-122.77673386980899,38.047758175794456],[-122.84707112471938,38.138129618891696],[-123.04243147171962,38.354804001848436],[-123.15680798818414,38.465500437349476],
     [-123.23553320492562,38.535877182151239],[-123.35086732408223,38.651910036719414],[-123.45763301882351,38.757692204931345],[-123.60327705800381,38.902664743325033],
     [-123.68955032886964,38.990249403417863],[-123.79457029839853,39.10506637042954],[-123.86454292131367,39.19856257350574],[-123.89959026792752,39.315996211438403],
     [-123.90342677323696,39.400114622094463],[-123.90735483128094,39.496257560632387],[-123.91125237174668,39.577212134567844],[-123.91901693509999,39.66712659439122],
     [-123.94234114273848,39.80768398526817],[-123.96566535037681,39.932953534579156],[-123.99292861035991,39.992593264718494],[-124.19510143902914,40.120597209247762],
     [-124.32344561619641,40.206810567148295],[-124.46425306411533,40.286066178214298],[-124.53798284749558,40.334662538185285],[-124.56169595004297,40.365903678824168],
     [-124.58085657573872,40.400869998897178],[-124.6037466475421,40.461963489125139],[-124.62492811920265,40.535454327553111],[-124.63991615536975,40.609661038047022],
     [-124.66248748380093,40.712361881530057],[-124.68502829465405,40.831881254680354],[-124.70008601946768,40.96270930117629],[-124.72265734789886,41.110147714779771],
     [-124.74519815875209,41.206337164386809],[-124.76025588356561,41.302455930510632],[-124.7828272119969,41.415326020018654],[-124.81291214404587,41.539342548312504],
     [-124.83542243732086,41.674259504659858],[-124.85802428333028,41.797814459606172],[-124.88807869780112,41.921006597274413],[-124.90313642261469,42.027126701748216],
     [-124.91065002623236,42.144380122681355],[-124.9181331122719,42.227948569109287],[-124.93319083708553,42.305843035079135],[-124.94070444070331,42.389213085149095],
     [-124.94824856189911,42.478078166155967],[-124.95055654943178,42.511816891401907],[-124.95543755909887,42.541800212543365],[-124.95576216551683,42.58881926305952],
     [-124.97081989033035,42.693952048585857],[-124.97830297636995,42.798873979195889],[-124.98631479194438,42.876810111869816],[-124.98581657998773,42.903615546829997],
     [-124.99336070118352,43.013764220027099],[-125.00090482237937,43.123570293871126],[-125.00090482237937,43.222335971073676],[-125.00090482237937,43.320869270583898],
     [-125.00090482237937,43.424739199265161],[-125.00090482237937,43.550164741977767],[-125.00090482237937,43.680834262654059],[-125.00090482237937,43.800350067213913],
     [-125.00090482237937,43.946773924627792],[-125.00090482237937,44.044076346615896],[-125.00090482237937,44.173658102565923],[-125.00090482237937,44.356727126113867],
     [-125.00841842599704,44.491089257144154],[-125.01593202961476,44.598242117204791],[-125.01593202961476,44.721246723067111],[-125.01593202961476,44.812068668851964],
     [-125.01593202961476,44.913302675262571],[-125.02341511565436,45.003729477628951],[-125.03095923685021,45.08880602490143],[-125.03847284046788,45.184277606617457],
     [-125.05353056528139,45.306017103338434],[-125.0610136513211,45.411743579181575],[-125.06858829009502,45.496088910970172],[-125.07607137613462,45.596168871404302],
     [-125.09112910094814,45.690786885662931],[-125.09112910094814,45.764275948766283],[-125.09864270456592,45.863903263486748],[-125.11370042937943,45.978971691538902],
     [-125.15129896504618,46.093798696997055],[-125.18886698313474,46.198052810571312],[-125.24152324361501,46.327918979406036],[-125.29417950409515,46.457661619232738],
     [-125.34680524699735,46.592122971339904],[-125.3994615074775,46.721215771239486],[-125.44457364676197,46.844836463286015],[-125.49716887208598,46.978390295793474],
     [-125.57239646099754,47.106472109759466],[-125.58313273818118,47.125105963257056],[-125.64756301475273,47.229191965514246],[-125.70773287885066,47.336314462151108],
     [-125.76793326052677,47.432993137616904],[-125.81304539981119,47.529552543485295],[-125.89572555718422,47.641223234149948],[-125.94083769646852,47.717159604540939],
     [-126.06120794224267,47.924184795803001],[-126.1487682387085,48.03953114645401],[-126.21908517094903,48.120254144163653],[-126.26419731023338,48.175467946903154],
     [-126.35448262395846,48.265665026546912],[-126.41459145290025,48.335776961601425],[-126.48230543819403,48.410625607137149],[-126.63269958086079,48.575129134706515],
     [-126.76806651629209,48.719130573816187],[-126.88834520933187,48.853004335396349],[-127.02374266234136,49.001120051799319],[-127.11396694091007,49.08494958196394],
     [-127.21921842671441,49.173547160076851],[-127.28690189443,49.237403726137131],[-127.35461587972384,49.296325726520649],[-127.41306401972383,49.349829478057643],
     [-127.45238427948851,49.379540060022237],[-127.53733753510443,49.444140695318822],[-127.62531923300844,49.52137975333693],[-127.75314204724395,49.623730752804363],
     //[-127.88853950025339,49.740528222935779],[-127.98624686486187,49.827896997616847],[-128.10658659305773,49.924880982887657],
     [-128.18181418196917,49.982933141910962],
     [-128.30963699620474,50.040858753763359],
     //[-128.46003113887144,50.127742619548656],[-128.61796940273416,50.204882944220287],[-128.84300994873041,50.314043045044059],
     //[-128.98594665527338,50.400043487548885],[-129.08899688720703,50.461957931518668],[-129.31197357177723,50.600048065185717],[-129.38500213623041,50.644979476928768],
     //[-129.58605489946569,50.770968205269241],[-129.69413545125235,50.852225456444955],
     [-130.04454776264879,51.124572457229874],[-130.13697774629304,51.156850961964722],
     [-130.21995544433594,51.200006484985465],[-130.33099365234375,51.260034561157283],[-130.47798538208002,51.401062011718921],[-130.73699951171869,51.513031005859375],
     [-130.86196899414062,51.599987030029354],[-131.04196929931635,51.72603797912609],[-131.16999053955067,51.799993515014762],[-131.35295104980463,51.906002044677905],
     [-131.46995544433582,51.999969482421875],[-131.66001892089844,52.153038024902457],[-131.72998809814442,52.200067520141772],[-131.97396850585937,52.365087509155387],
     [-132.00598907470697,52.40007400512701],[-132.19000244140625,52.600048065185661],[-132.22995758056635,52.644052505493278],[-132.44000244140614,52.800054550170898],
     [-132.65498352050781,52.95908355712902],[-132.71800994873036,53.000030517578182],[-132.92400360107411,53.135065078735522],[-132.97396850585932,53.20003700256359],
     [-133.08198547363276,53.341958999633903],[-133.17196655273438,53.400012969970817],[-133.40899658203119,53.552030563354606],[-133.44698333740229,53.600017547607479],
     [-133.5709609985351,53.75904655456543],[-133.60298156738276,53.799993515014762],[-133.73397064208973,53.966032028198242],[-133.77001190185547,53.999969482421989],
     [-133.95497894287104,54.173048019409237],[-133.99897766113276,54.20009803771984],[-134.23699951171869,54.346004486084041],[-134.29598236083984,54.400074005126953],
     [-134.46300506591791,54.553020477295092],[-134.49595642089832,54.600078582763786],[-134.63549804687494,54.797028541565055],[-134.81088461534318,55.017699237567456],
     [-134.82935989875682,55.043150180419673], 
     [-136,55.043150180419673], [-136,35.026738764076697]],
          "temperature models":[{"model":"plate model", "min depth":-1e2, "max depth":100e3, "spreading velocity":0.04, "ridge coordinates":
     [
     [-130.04499816894531,51.127998352050781],[-130.2239990234375,51],[-130.50300598144531,50.799999237060547],[-130.65899658203125,50.687999725341797],
     [-129.71800231933594,50.187999725341797],[-129.91000366210937,50],[-130.11500549316406,49.799999237060547],[-130.32000732421875,49.599998474121094],[-130.52499389648437,49.400001525878906],[-130.53399658203125,49.390998840332031],
     [-128.79600524902344,48.655998229980469],[-128.822998046875,48.599998474121094],[-128.91799926757812,48.400001525878906],[-129.01300048828125,48.200000762939453],[-129.10800170898437,48],[-129.17999267578125,47.847999572753906],
     [-128.69400024414062,47.7239990234375],[-128.75799560546875,47.599998474121094],[-128.86199951171875,47.400001525878906],[-128.96600341796875,47.200000762939453],[-129.07000732421875,47],[-129.17300415039062,46.799999237060547],[-129.27699279785156,46.599998474121094],[-129.38099670410156,46.400001525878906],[-129.48500061035156,46.200000762939453],[-129.58900451660156,46],[-129.69200134277344,45.799999237060547],[-129.79600524902344,45.599998474121094],[-129.89999389648437,45.400001525878906],[-130.00399780273437,45.200000762939453],[-130.10699462890625,45],[-130.21099853515625,44.799999237060547],[-130.31500244140625,44.599998474121094],
     [-128.4949951171875,43.944000244140625],[-128.61700439453125,43.799999237060547],[-128.66099548339844,43.748001098632812],
     [-126.67299652099609,43.035999298095689],[-126.68399810791016,43.011001586914055],[-126.68900299072266,42.999999999999986],[-126.69499969482422,42.985000610351555],[-126.70600128173828,42.959999084472642],[-126.71700286865234,42.935001373291016],[-126.72899627685547,42.909000396728509],[-126.73999786376953,42.883998870849609],[-126.75099945068359,42.859001159667976],[-126.76200103759766,42.833000183105469],[-126.77300262451172,42.80799865722657],[-126.77700042724609,42.79999923706054],[-126.78399658203125,42.783000946044929],[-126.79499816894531,42.756999969482422],[-126.80599975585937,42.731998443603509],[-126.81700134277344,42.707000732421868],
     [-126.8280029296875,42.680999755859382],[-126.83999633789063,42.655998229980469],[-126.85099792480469,42.631000518798828],[-126.86199951171875,42.604999542236335],[-126.86399841308594,42.599998474121094],[-126.87300109863281,42.58000183105468],[-126.88400268554687,42.555000305175781],[-126.89499664306641,42.528999328613288],[-126.90599822998047,42.504001617431626],[-126.91699981689453,42.47900009155272],[-126.92800140380859,42.452999114990213],[-126.93900299072266,42.428001403808594],[-126.95099639892578,42.40299987792968],[-126.95200347900391,42.400001525878892],[-126.96199798583984,42.377998352050781],[-126.97299957275391,42.352001190185533],
     [-126.98400115966797,42.326999664306641],[-126.99500274658203,42.301998138427734],[-127.00599670410156,42.276000976562493],[-127.01699829101562,42.250999450683601],[-127.02799987792969,42.226001739501953],[-127.03900146484375,42.200000762939453],[-127.05000305175781,42.174999237060547],[-127.06199645996094,42.150001525878906],[-127.072998046875,42.124000549316406],[-127.08399963378906,42.098999023437486],[-127.09500122070312,42.074001312255852],[-127.10600280761719,42.048000335693359],[-127.11699676513672,42.022998809814446],[-127.12699890136719,42],[-127.12799835205078,41.998001098632805],[-127.13899993896484,41.972000122070313],
     [-127.15000152587891,41.946998596191406],[-127.16100311279297,41.922000885009766],[-127.17299652099609,41.895999908447266],[-127.18399810791016,41.870998382568366],[-127.19499969482422,41.846000671386719],[-127.20600128173828,41.819999694824219],[-127.21499633789062,41.799999237060547],[-127.28167774175972,41.63504645876381],[-127.43599080535165,41.634943347434728],[-127.43800354003906,41.610000610351555],[-127.43900299072266,41.59999847412108],[-127.44000244140625,41.583999633789048],[-127.44100189208984,41.558998107910156],[-127.44300079345703,41.533000946044908],[-127.44499969482422,41.507999420166009],[-127.44699859619141,41.483001708984375],
     [-127.447998046875,41.457000732421861],[-127.44999694824219,41.431999206542969],[-127.45200347900391,41.405998229980462],[-127.45200347900391,41.400001525878899],[-127.45400238037109,41.381000518798828],[-127.45500183105469,41.355998992919922],[-127.45700073242187,41.330001831054688],[-127.45899963378906,41.305000305175788],[-127.46099853515625,41.278999328613281],[-127.46199798583984,41.254001617431641],[-127.46399688720703,41.228000640869141],[-127.46600341796875,41.202999114990227],[-127.46600341796875,41.20000076293946],[-127.46800231933594,41.178001403808601],[-127.46900177001953,41.152000427246094],[-127.47100067138672,41.126998901367188],
     [-127.47299957275391,41.101001739501967],[-127.47499847412109,41.076000213623047],[-127.47599792480469,41.050998687744141],[-127.47799682617187,41.025001525878899],[-127.48000335693359,41],[-127.48000335693359,41.00999832153321],[-127.48200225830078,40.9739990234375],[-127.48300170898438,40.949001312255852],[-127.48500061035156,40.923999786376946],[-127.48699951171875,40.897998809814453],[-127.48899841308594,40.873001098632798],[-127.48999786376953,40.847000122070313],[-127.49199676513672,40.821998596191399],[-127.49400329589844,40.79999923706054],[-127.49400329589844,40.797000885009744],[-127.49600219726562,40.770999908447266],[-127.49700164794922,40.745998382568359],
     [-127.49900054931641,40.720001220703125],[-127.50099945068359,40.694999694824219],[-127.50299835205078,40.668998718261705],[-127.50399780273437,40.644001007080071],[-127.50599670410156,40.618999481201165],[-127.50700378417969,40.599998474121094],[-127.50800323486328,40.592998504638665],[-127.51000213623047,40.568000793457024],[-127.51100158691406,40.541999816894524],[-127.51300048828125,40.516998291015625],[-127.51499938964844,40.492000579833984],[-127.51699829101562,40.465999603271477],[-127.51799774169922,40.441001892089844],[-127.51999664306641,40.415000915527344],[-127.52100372314453,40.400001525878913],[-127.52200317382813,40.389999389648438]
     ]}],
     "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,0.5], "min depth":-1e2, "max depth":15e3},
                           {"model":"uniform", "compositions":[2], "min depth":15e3, "max depth":30e3}]
    },

    {"model":"oceanic plate", "name":"Pasific southest (0)", "min depth":-1e2,"max depth":100e3, "coordinates":
     [[-124.58085657573872,40.400869998897178],[-124.56169595004297,40.365903678824168],[-124.53798284749558,40.334662538185285],[-124.46425306411533,40.286066178214298],
     [-124.32344561619641,40.206810567148295],[-124.19510143902914,40.120597209247762],[-123.99292861035991,39.992593264718494],[-123.96566535037681,39.932953534579156],
     [-123.94234114273848,39.80768398526817],[-123.91901693509999,39.66712659439122],[-123.91125237174668,39.577212134567844],[-123.90735483128094,39.496257560632387],
     [-123.90342677323696,39.400114622094463],[-123.89959026792752,39.315996211438403],[-123.86454292131367,39.19856257350574],[-123.79457029839853,39.10506637042954],
     [-123.68955032886964,38.990249403417863],[-123.60327705800381,38.902664743325033],[-123.45763301882351,38.757692204931345],[-123.35086732408223,38.651910036719414],
     [-123.23553320492562,38.535877182151239],[-123.15680798818414,38.465500437349476],[-123.04243147171962,38.354804001848436],[-122.84707112471938,38.138129618891696],
     [-122.77673386980899,38.047758175794456],[-122.68624318948224,37.933754160924821],[-122.61834876039984,37.832996147222048],[-122.54652300999368,37.74016580820421],
     [-122.49843602016591,37.677498177477901],[-122.45281742259908,37.631268840719542],[-122.38065642027948,37.545925782427048],[-122.33354466861408,37.49091562348093],
     [-122.2964121118514,37.453237354638418],[-122.27151125838043,37.428654657266122],[-122.24885806811528,37.404279337111973],[-122.21486053342659,37.365971508253097],
     [-122.18537362833922,37.335571261058476],[-122.13081161527418,37.290741949604353],[-122.09556558307531,37.26296665270047],[-122.06084978266296,37.230940161988258],
     [-122.03105205620977,37.204234039470748],[-121.98598162856348,37.167992306611893],[-121.95555511187678,37.141996698957144],[-121.92357308919088,37.118855092475712],
     [-121.88399134347702,37.094449022037679],[-121.84121748147328,37.072019692455228],[-121.80503194745722,37.049836399722381],[-121.76465348371033,37.020752188503081],
     [-121.74247097000728,37.001658514256235],[-121.67950895795178,36.946478562229061],[-121.64048543975076,36.919630706516728],[-121.61260729844474,36.901346802870705],
     [-121.58944549593645,36.884433460794298],[-121.55508229475919,36.86255855283656],[-121.52610334369285,36.840114879061275],[-121.49115617325845,36.815708963134853],
     [-121.45673589161618,36.793529056408488],[-121.40886245141422,36.764290491092368],[-121.36875476511875,36.739378826299514],[-121.33607120903258,36.715009966943455],
     [-121.2926215875213,36.683474164488587],[-121.26368275517387,36.66066564742988],[-121.23916541841322,36.6393384352711],[-121.18670053002239,36.59376134384371],
     [-121.16371695120017,36.575907615284507],[-121.14071876667788,36.554153106004833],[-121.10791832828488,36.519877354310154],[-121.07207309013523,36.485893948235002],
     [-121.01576174330518,36.431157243479788],[-120.95735784000937,36.369754311823044],[-120.92541506386243,36.33833146896319],[-120.89398206909289,36.309107397971445],
     [-120.86409649358581,36.278691639242709],[-120.82941491103179,36.246844661554348],[-120.78691544416824,36.207853475097693],[-120.75261850319879,36.176845516870685],
     [-120.67630445821555,36.10849694324088],[-120.65117422382883,36.085785007481036],[-120.61725375997969,36.05478276694754],[-120.58260325960816,36.022430244199256],
     [-120.54669696564423,35.99339198982949],[-120.51162105275768,35.964166488899593],[-120.47477761013653,35.933559258904154],[-120.44946117595816,35.910487135624123],
     [-120.42164544113507,35.884236198652957],[-120.38094720203441,35.849557704201459],[-120.35500781938225,35.825358236558714],[-120.32579703346602,35.788598829563568],
     [-120.30506054229693,35.754745248739425],[-120.28108807142428,35.726943855131026],[-120.24859639654591,35.695779235772477],[-120.22783449400163,35.674173828857363],
     [-120.20593586300009,35.653689383733365],[-120.18506317056102,35.633728585922995],[-120.14862523818482,35.59870591463266],[-120.1121292412505,35.563296615814409],
     [-120.06943818848163,35.521453648690397],[-120.03159767449893,35.483016977920613],[-119.99889777420293,35.446300194659898],[-119.97804209353507,35.42203524916863],
     [-119.94692681693385,35.39038972770112],[-119.90632396186282,35.351792300798877],[-119.86055950731412,35.304971252158566],[-119.82807176026733,35.272843231920774],
     [-119.79409829064866,35.241824187169698],[-119.76978766719469,35.218363290625462],[-119.74821603067983,35.197882826954526],[-119.71763954289975,35.171045968622082],
     [-119.65858565209328,35.125880316256769],[-119.63417894839137,35.106540426058871],[-119.60939570151095,35.084327995405658],[-119.57383732651698,35.056444555160567],
     [-119.53587316701834,35.026738764076697],[-136,35.05755615234375],[-136,40.400869998897178]],
     "temperature models":[{"model":"linear", "min depth":-1e2, "max depth":100e3, "bottom temperature":-1}],
     "composition models":[{"model":"uniform", "compositions":[1], "min depth":-1e2, "max depth":15e3},
                           {"model":"uniform", "compositions":[2], "min depth":15e3, "max depth":30e3}]
    },

    {"model":"continental plate", "name":"NA", "min depth":-1e2, "coordinates":
     [[-119.53587316701834,35.026738764076697],[-119.57383732651698,35.056444555160567],[-119.60939570151095,35.084327995405658],[-119.63417894839137,35.106540426058871],
     [-119.65858565209328,35.125880316256769],[-119.71763954289975,35.171045968622082],[-119.74821603067983,35.197882826954526],[-119.76978766719469,35.218363290625462],
     [-119.79409829064866,35.241824187169698],[-119.82807176026733,35.272843231920774],[-119.86055950731412,35.304971252158566],[-119.90632396186282,35.351792300798877],
     [-119.94692681693385,35.39038972770112],[-119.97804209353507,35.42203524916863],[-119.99889777420293,35.446300194659898],[-120.03159767449893,35.483016977920613],
     [-120.06943818848163,35.521453648690397],[-120.1121292412505,35.563296615814409],[-120.14862523818482,35.59870591463266],[-120.18506317056102,35.633728585922995],
     [-120.20593586300009,35.653689383733365],[-120.22783449400163,35.674173828857363],[-120.24859639654591,35.695779235772477],[-120.28108807142428,35.726943855131026],
     [-120.30506054229693,35.754745248739425],[-120.32579703346602,35.788598829563568],[-120.35500781938225,35.825358236558714],[-120.38094720203441,35.849557704201459],
     [-120.42164544113507,35.884236198652957],[-120.44946117595816,35.910487135624123],[-120.47477761013653,35.933559258904154],[-120.51162105275768,35.964166488899593],
     [-120.54669696564423,35.99339198982949],[-120.58260325960816,36.022430244199256],[-120.61725375997969,36.05478276694754],[-120.65117422382883,36.085785007481036],
     [-120.67630445821555,36.10849694324088],[-120.75261850319879,36.176845516870685],[-120.78691544416824,36.207853475097693],[-120.82941491103179,36.246844661554348],
     [-120.86409649358581,36.278691639242709],[-120.89398206909289,36.309107397971445],[-120.92541506386243,36.33833146896319],[-120.95735784000937,36.369754311823044],
     [-121.01576174330518,36.431157243479788],[-121.07207309013523,36.485893948235002],[-121.10791832828488,36.519877354310154],[-121.14071876667788,36.554153106004833],
     [-121.16371695120017,36.575907615284507],[-121.18670053002239,36.59376134384371],[-121.23916541841322,36.6393384352711],[-121.26368275517387,36.66066564742988],
     [-121.2926215875213,36.683474164488587],[-121.33607120903258,36.715009966943455],[-121.36875476511875,36.739378826299514],[-121.40886245141422,36.764290491092368],
     [-121.45673589161618,36.793529056408488],[-121.49115617325845,36.815708963134853],[-121.52610334369285,36.840114879061275],[-121.55508229475919,36.86255855283656],
     [-121.58944549593645,36.884433460794298],[-121.61260729844474,36.901346802870705],[-121.64048543975076,36.919630706516728],[-121.67950895795178,36.946478562229061],
     [-121.74247097000728,37.001658514256235],[-121.76465348371033,37.020752188503081],[-121.80503194745722,37.049836399722381],[-121.84121748147328,37.072019692455228],
     [-121.88399134347702,37.094449022037679],[-121.92357308919088,37.118855092475712],[-121.95555511187678,37.141996698957144],[-121.98598162856348,37.167992306611893],
     [-122.03105205620977,37.204234039470748],[-122.06084978266296,37.230940161988258],[-122.09556558307531,37.26296665270047],[-122.13081161527418,37.290741949604353],
     [-122.18537362833922,37.335571261058476],[-122.21486053342659,37.365971508253097],[-122.24885806811528,37.404279337111973],[-122.27151125838043,37.428654657266122],
     [-122.2964121118514,37.453237354638418],[-122.33354466861408,37.49091562348093],[-122.38065642027948,37.545925782427048],[-122.45281742259908,37.631268840719542],
     [-122.49843602016591,37.677498177477901],[-122.54652300999368,37.74016580820421],[-122.61834876039984,37.832996147222048],[-122.68624318948224,37.933754160924821],
     [-122.77673386980899,38.047758175794456],[-122.84707112471938,38.138129618891696],[-123.04243147171962,38.354804001848436],[-123.15680798818414,38.465500437349476],
     [-123.23553320492562,38.535877182151239],[-123.35086732408223,38.651910036719414],[-123.45763301882351,38.757692204931345],[-123.60327705800381,38.902664743325033],
     [-123.68955032886964,38.990249403417863],[-123.79457029839853,39.10506637042954],[-123.86454292131367,39.19856257350574],[-123.89959026792752,39.315996211438403],
     [-123.90342677323696,39.400114622094463],[-123.90735483128094,39.496257560632387],[-123.91125237174668,39.577212134567844],[-123.91901693509999,39.66712659439122],
     [-123.94234114273848,39.80768398526817],[-123.96566535037681,39.932953534579156],[-123.99292861035991,39.992593264718494],[-124.19510143902914,40.120597209247762],
     [-124.32344561619641,40.206810567148295],[-124.46425306411533,40.286066178214298],[-124.53798284749558,40.334662538185285],[-124.56169595004297,40.365903678824168],
     [-124.58085657573872,40.400869998897178],[-124.6037466475421,40.461963489125139],[-124.62492811920265,40.535454327553111],[-124.63991615536975,40.609661038047022],
     [-124.66248748380093,40.712361881530057],[-124.68502829465405,40.831881254680354],[-124.70008601946768,40.96270930117629],[-124.72265734789886,41.110147714779771],
     [-124.74519815875209,41.206337164386809],[-124.76025588356561,41.302455930510632],[-124.7828272119969,41.415326020018654],[-124.81291214404587,41.539342548312504],
     [-124.83542243732086,41.674259504659858],[-124.85802428333028,41.797814459606172],[-124.88807869780112,41.921006597274413],[-124.90313642261469,42.027126701748216],
     [-124.91065002623236,42.144380122681355],[-124.9181331122719,42.227948569109287],[-124.93319083708553,42.305843035079135],[-124.94070444070331,42.389213085149095],
     [-124.94824856189911,42.478078166155967],[-124.95055654943178,42.511816891401907],[-124.95543755909887,42.541800212543365],[-124.95576216551683,42.58881926305952],
     [-124.97081989033035,42.693952048585857],[-124.97830297636995,42.798873979195889],[-124.98631479194438,42.876810111869816],[-124.98581657998773,42.903615546829997],
     [-124.99336070118352,43.013764220027099],[-125.00090482237937,43.123570293871126],[-125.00090482237937,43.222335971073676],[-125.00090482237937,43.320869270583898],
     [-125.00090482237937,43.424739199265161],[-125.00090482237937,43.550164741977767],[-125.00090482237937,43.680834262654059],[-125.00090482237937,43.800350067213913],
     [-125.00090482237937,43.946773924627792],[-125.00090482237937,44.044076346615896],[-125.00090482237937,44.173658102565923],[-125.00090482237937,44.356727126113867],
     [-125.00841842599704,44.491089257144154],[-125.01593202961476,44.598242117204791],[-125.01593202961476,44.721246723067111],[-125.01593202961476,44.812068668851964],
     [-125.01593202961476,44.913302675262571],[-125.02341511565436,45.003729477628951],[-125.03095923685021,45.08880602490143],[-125.03847284046788,45.184277606617457],
     [-125.05353056528139,45.306017103338434],[-125.0610136513211,45.411743579181575],[-125.06858829009502,45.496088910970172],[-125.07607137613462,45.596168871404302],
     [-125.09112910094814,45.690786885662931],[-125.09112910094814,45.764275948766283],[-125.09864270456592,45.863903263486748],[-125.11370042937943,45.978971691538902],
     [-125.15129896504618,46.093798696997055],[-125.18886698313474,46.198052810571312],[-125.24152324361501,46.327918979406036],[-125.29417950409515,46.457661619232738],
     [-125.34680524699735,46.592122971339904],[-125.3994615074775,46.721215771239486],[-125.44457364676197,46.844836463286015],[-125.49716887208598,46.978390295793474],
     [-125.57239646099754,47.106472109759466],[-125.58313273818118,47.125105963257056],[-125.64756301475273,47.229191965514246],[-125.70773287885066,47.336314462151108],
     [-125.76793326052677,47.432993137616904],[-125.81304539981119,47.529552543485295],[-125.89572555718422,47.641223234149948],[-125.94083769646852,47.717159604540939],
     [-126.06120794224267,47.924184795803001],[-126.1487682387085,48.03953114645401],[-126.21908517094903,48.120254144163653],[-126.26419731023338,48.175467946903154],
     [-126.35448262395846,48.265665026546912],[-126.41459145290025,48.335776961601425],[-126.48230543819403,48.410625607137149],[-126.63269958086079,48.575129134706515],
     [-126.76806651629209,48.719130573816187],[-126.88834520933187,48.853004335396349],[-127.02374266234136,49.001120051799319],[-127.11396694091007,49.08494958196394],
     [-127.21921842671441,49.173547160076851],[-127.28690189443,49.237403726137131],[-127.35461587972384,49.296325726520649],[-127.41306401972383,49.349829478057643],
     [-127.45238427948851,49.379540060022237],[-127.53733753510443,49.444140695318822],[-127.62531923300844,49.52137975333693],[-127.75314204724395,49.623730752804363],
     //[-127.88853950025339,49.740528222935779],[-127.98624686486187,49.827896997616847],[-128.10658659305773,49.924880982887657],
     [-128.18181418196917,49.982933141910962],
     [-128.30963699620474,50.040858753763359],
     //[-128.46003113887144,50.127742619548656],[-128.61796940273416,50.204882944220287],[-128.84300994873041,50.314043045044059],
     //[-128.98594665527338,50.400043487548885],[-129.08899688720703,50.461957931518668],[-129.31197357177723,50.600048065185717],[-129.38500213623041,50.644979476928768],
     //[-129.58605489946569,50.770968205269241],[-129.69413545125235,50.852225456444955],
     [-130.04454776264879,51.124572457229874],[-130.13697774629304,51.156850961964722],
     [-130.21995544433594,51.200006484985465],[-130.33099365234375,51.260034561157283],[-130.47798538208002,51.401062011718921],[-130.73699951171869,51.513031005859375],
     [-130.86196899414062,51.599987030029354],[-131.04196929931635,51.72603797912609],[-131.16999053955067,51.799993515014762],[-131.35295104980463,51.906002044677905],
     [-131.46995544433582,51.999969482421875],[-131.66001892089844,52.153038024902457],[-131.72998809814442,52.200067520141772],[-131.97396850585937,52.365087509155387],
     [-132.00598907470697,52.40007400512701],[-132.19000244140625,52.600048065185661],[-132.22995758056635,52.644052505493278],[-132.44000244140614,52.800054550170898],
     [-132.65498352050781,52.95908355712902],[-132.71800994873036,53.000030517578182],[-132.92400360107411,53.135065078735522],[-132.97396850585932,53.20003700256359],
     [-133.08198547363276,53.341958999633903],[-133.17196655273438,53.400012969970817],[-133.40899658203119,53.552030563354606],[-133.44698333740229,53.600017547607479],
     [-133.5709609985351,53.75904655456543],[-133.60298156738276,53.799993515014762],[-133.73397064208973,53.966032028198242],[-133.77001190185547,53.999969482421989],
     [-133.95497894287104,54.173048019409237],[-133.99897766113276,54.20009803771984],[-134.23699951171869,54.346004486084041],[-134.29598236083984,54.400074005126953],
     [-134.46300506591791,54.553020477295092],[-134.49595642089832,54.600078582763786],[-134.63549804687494,54.797028541565055],[-134.81088461534318,55.017699237567456],
     [-134.82935989875682,55.043150180419673],[-109.98626157568839,54.987757598234737],[-110.03077376778907,35.001783345023512],[-119.53587316701834,35.026738764076697]],
     "max depth":120e3,
     "temperature models":[{"model":"linear", "min depth":-1e2, "max depth":120e3, "bottom temperature":-1}],
     "composition models":[{"model":"uniform","compositions":[0], "min depth":-1e2, "max depth":30e3}//,
     //                           {"model":"uniform","compositions":[7], "min depth":30e3, "max depth":100e3}
    ]},
    {"model":"subducting plate", "name":"Subducting plate", "dip point":[135,45], "coordinates":
     //[[-125.045,40.492],[-124.946,40.719],[-125.306,41.670],[-125.273,42.637],[-125.336,43.130],[-125.363,44.101],[-125.374,44.601],[-125.439,45.094],[-125.550,45.588],[-125.699,46.077],[-125.853,46.565],[-126.017,47.052],[-126.088,47.348],[-126.212,47.529],[-129.995,51.109]],
     [[-125.045,40.492],[-124.946,40.719],[-125.1,41.670],[-125.273,42.637],[-125.336,43.130],[-125.363,44.101],[-125.374,44.601],[-125.439,45.094],[-125.550,45.588],[-125.699,46.077],[-125.853,46.565],[-126.017,47.052],[-126.108,47.348],[-126.212,47.529],[-129.995,51.109]],
     "segments":[
       {"length":46.762e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[0.000,8.244],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=0.000km
       {"length":11.216e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[8.244,8.851],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=5.728km
       {"length":11.246e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[8.851,9.123],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=7.280km
       {"length":11.234e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[9.123,9.219],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=9.052km
       {"length":11.247e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[9.219,9.406],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=10.767km
       {"length":11.256e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[9.406,9.772],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=12.585km
       {"length":11.266e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[9.772,10.362],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=14.474km
       {"length":11.282e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[10.362,11.151],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=16.444km
       {"length":11.315e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[11.151,12.124],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=18.520km
       {"length":11.353e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[12.124,13.246],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=20.788km
       {"length":11.402e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[13.246,14.489],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=23.255km
       {"length":11.465e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[14.489,15.795],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=25.958km
       {"length":11.520e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[15.795,17.008],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=28.933km
       {"length":11.596e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[17.008,17.984],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,1.787], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=32.134km
       {"length":11.629e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[17.984,18.544],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,1.438], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=35.617km
       {"length":11.649e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[18.544,18.683],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,1.077], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=39.227km
       {"length":11.650e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[18.683,18.530],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,0.708], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=42.921km
       {"length":11.634e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[18.530,18.368],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,0.336], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=46.636km
       {"length":11.625e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[18.368,18.614]},  // depth=50.319km
       {"length":11.677e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[18.614,19.741]},  // depth=53.994km
       {"length":11.785e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[19.741,22.164]},  // depth=57.850km
       {"length":12.039e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[22.164,26.115]},  // depth=62.036km
       {"length":12.535e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[26.115,31.500]},  // depth=66.906km
       {"length":13.339e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[31.500,37.724]},  // depth=72.911km
       {"length":14.461e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[37.724,43.851]},  // depth=80.464km
       {"length":15.961e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[43.851,49.173]},  // depth=89.875km
       {"length":17.371e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[49.173,53.194]},  // depth=101.477km
       {"length":18.861e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[53.194,55.894]},  // depth=114.973km
       {"length":19.745e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[55.894,57.319]},  // depth=130.359km
       {"length":20.251e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[57.319,57.642]},  // depth=146.837km
       {"length":20.179e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[57.642,57.046]},  // depth=163.937km
       {"length":19.702e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[57.046,55.798]},  // depth=180.970km
       {"length":19.018e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[55.798,54.310]},  // depth=197.455km
       {"length":18.540e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[54.310,53.068]},  // depth=213.135km
       {"length":17.930e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[53.068,52.553]},  // depth=228.251km
       {"length":18.000e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[52.553,53.002]},  // depth=242.632km
       {"length":18.355e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[53.002,54.341]},  // depth=257.118km
       {"length":19.218e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[54.341,56.310]},  // depth=272.062km
       {"length":20.443e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[56.310,58.597]},  // depth=288.071km
       {"length":21.748e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[58.597,60.907]},  // depth=305.548km
       {"length":23.763e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[60.907,63.085]},  // depth=324.554km
       {"length":25.102e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[63.085,64.903]},  // depth=345.852km
       {"length":26.317e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[64.903,66.236]},
        {"length":0.0, "thickness":[300.0], "top truncation":[-100.0], "angle":[66.236]},
        {"length":0.0, "thickness":[300.0], "top truncation":[-100.0], "angle":[66.236]},
        {"length":0.0, "thickness":[300.0], "top truncation":[-100.0], "angle":[66.236]},
        {"length":0.0, "thickness":[300.0], "top truncation":[-100.0], "angle":[66.236]},
        {"length":0.0, "thickness":[300.0], "top truncation":[-100.0], "angle":[66.236]},
        {"length":0.0, "thickness":[300.0], "top truncation":[-100.0], "angle":[66.236]}
     //]  
     //}
    ],
    "temperature models":[{"model":"plate model","adiabatic heating":true, "density":3300, "plate velocity":0.04, "max distance slab top":50e3//, //,"coupling depth":50e3, "forearc cooling factor":10, "taper distance":50e3, "min distance slab top":-200e3,"max distance slab top":300e3,
    //"temperature models":[{"model":"plate model", "density":3300, "plate velocity":0.04,"max distance slab top":20e3},{"model":"adiabatic", "min distance slab top":20e3}]

     }],
     "composition models":[{"model":"uniform", "compositions":[1,3],  "fractions":[1,0], "max distance slab top":15e3},
                           {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                           {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]
    },
    {"model":"continental plate", "name":"South fault", "min depth":-10e3, "max depth":100e3, "coordinates":
     [[-124.68085657573872,40.400869998897178],[-124.66169595004297,40.365903678824168],[-124.63798284749558,40.334662538185285],[-124.56425306411533,40.286066178214298],
     [-124.42344561619641,40.206810567148295],[-124.29510143902914,40.120597209247762],[-124.09292861035991,39.992593264718494],[-124.06566535037681,39.932953534579156],
     [-124.04234114273848,39.80768398526817],[-124.01901693509999,39.66712659439122],[-124.01125237174668,39.577212134567844],[-124.00735483128094,39.496257560632387],
     [-124.00342677323696,39.400114622094463],[-123.99959026792752,39.315996211438403],[-123.96454292131367,39.19856257350574],[-123.89457029839853,39.10506637042954],
     [-123.78955032886964,38.990249403417863],[-123.70327705800381,38.902664743325033],[-123.55763301882351,38.757692204931345],[-123.45086732408223,38.651910036719414],
     [-123.33553320492562,38.535877182151239],[-123.25680798818414,38.465500437349476],[-123.14243147171962,38.354804001848436],[-122.94707112471938,38.138129618891696],
     [-122.87673386980899,38.047758175794456],[-122.78624318948224,37.933754160924821],[-122.71834876039984,37.832996147222048],[-122.64652300999368,37.74016580820421],
     [-122.59843602016591,37.677498177477901],[-122.55281742259908,37.631268840719542],[-122.48065642027948,37.545925782427048],[-122.43354466861408,37.49091562348093],
     [-122.3964121118514,37.453237354638418],[-122.37151125838043,37.428654657266122],[-122.34885806811528,37.404279337111973],[-122.31486053342659,37.365971508253097],
     [-122.28537362833922,37.335571261058476],[-122.23081161527418,37.290741949604353],[-122.19556558307531,37.26296665270047],[-122.16084978266296,37.230940161988258],
     [-122.13105205620977,37.204234039470748],[-122.08598162856348,37.167992306611893],[-122.05555511187678,37.141996698957144],[-122.02357308919088,37.118855092475712],
     [-121.98399134347702,37.094449022037679],[-121.94121748147328,37.072019692455228],[-121.90503194745722,37.049836399722381],[-121.86465348371033,37.020752188503081],
     [-121.84247097000728,37.001658514256235],[-121.77950895795178,36.946478562229061],[-121.74048543975076,36.919630706516728],[-121.71260729844474,36.901346802870705],
     [-121.68944549593645,36.884433460794298],[-121.65508229475919,36.86255855283656],[-121.62610334369285,36.840114879061275],[-121.59115617325845,36.815708963134853],
     [-121.55673589161618,36.793529056408488],[-121.50886245141422,36.764290491092368],[-121.46875476511875,36.739378826299514],[-121.43607120903258,36.715009966943455],
     [-121.3926215875213,36.683474164488587],[-121.36368275517387,36.66066564742988],[-121.33916541841322,36.6393384352711],[-121.28670053002239,36.59376134384371],
     [-121.26371695120017,36.575907615284507],[-121.24071876667788,36.554153106004833],[-121.20791832828488,36.519877354310154],[-121.17207309013523,36.485893948235002],
     [-121.11576174330518,36.431157243479788],[-121.05735784000937,36.369754311823044],[-121.02541506386243,36.33833146896319],[-120.99398206909289,36.309107397971445],
     [-120.96409649358581,36.278691639242709],[-120.92941491103179,36.246844661554348],[-120.88691544416824,36.207853475097693],[-120.85261850319879,36.176845516870685],
     [-120.77630445821555,36.10849694324088],[-120.75117422382883,36.085785007481036],[-120.71725375997969,36.05478276694754],[-120.68260325960816,36.022430244199256],
     [-120.64669696564423,35.99339198982949],[-120.61162105275768,35.964166488899593],[-120.57477761013653,35.933559258904154],[-120.54946117595816,35.910487135624123],
     [-120.52164544113507,35.884236198652957],[-120.48094720203441,35.849557704201459],[-120.45500781938225,35.825358236558714],[-120.42579703346602,35.788598829563568],
     [-120.40506054229693,35.754745248739425],[-120.38108807142428,35.726943855131026],[-120.34859639654591,35.695779235772477],[-120.32783449400163,35.674173828857363],
     [-120.30593586300009,35.653689383733365],[-120.28506317056102,35.633728585922995],[-120.24862523818482,35.59870591463266],[-120.2121292412505,35.563296615814409],
     [-120.16943818848163,35.521453648690397],[-120.13159767449893,35.483016977920613],[-120.09889777420293,35.446300194659898],[-120.07804209353507,35.42203524916863],
     [-120.04692681693385,35.39038972770112],[-120.00632396186282,35.351792300798877],[-119.96055950731412,35.304971252158566],[-119.92807176026733,35.272843231920774],
     [-119.89409829064866,35.241824187169698],[-119.86978766719469,35.218363290625462],[-119.84821603067983,35.197882826954526],[-119.81763954289975,35.171045968622082],
     [-119.75858565209328,35.125880316256769],[-119.73417894839137,35.106540426058871],[-119.70939570151095,35.084327995405658],[-119.67383732651698,35.056444555160567],
     [-119.63587316701834,35.026738764076697],
     [-124.48085657573872,40.400869998897178],[-124.46169595004297,40.365903678824168],[-124.43798284749558,40.334662538185285],[-124.36425306411533,40.286066178214298],
     [-124.22344561619641,40.206810567148295],[-124.09510143902914,40.120597209247762],[-123.89292861035991,39.992593264718494],[-123.86566535037681,39.932953534579156],
     [-123.84234114273848,39.80768398526817],[-123.81901693509999,39.66712659439122],[-123.81125237174668,39.577212134567844],[-123.80735483128094,39.496257560632387],
     [-123.80342677323696,39.400114622094463],[-123.79959026792752,39.315996211438403],[-123.76454292131367,39.19856257350574],[-123.69457029839853,39.10506637042954],
     [-123.58955032886964,38.990249403417863],[-123.50327705800381,38.902664743325033],[-123.35763301882351,38.757692204931345],[-123.25086732408223,38.651910036719414],
     [-123.13553320492562,38.535877182151239],[-123.05680798818414,38.465500437349476],[-122.94243147171962,38.354804001848436],[-122.74707112471938,38.138129618891696],
     [-122.67673386980899,38.047758175794456],[-122.58624318948224,37.933754160924821],[-122.51834876039984,37.832996147222048],[-122.44652300999368,37.74016580820421],
     [-122.39843602016591,37.677498177477901],[-122.35281742259908,37.631268840719542],[-122.28065642027948,37.545925782427048],[-122.23354466861408,37.49091562348093],
     [-122.1964121118514,37.453237354638418],[-122.17151125838043,37.428654657266122],[-122.14885806811528,37.404279337111973],[-122.11486053342659,37.365971508253097],
     [-122.08537362833922,37.335571261058476],[-122.03081161527418,37.290741949604353],[-121.99556558307531,37.26296665270047],[-121.96084978266296,37.230940161988258],
     [-121.93105205620977,37.204234039470748],[-121.88598162856348,37.167992306611893],[-121.85555511187678,37.141996698957144],[-121.82357308919088,37.118855092475712],
     [-121.78399134347702,37.094449022037679],[-121.74121748147328,37.072019692455228],[-121.70503194745722,37.049836399722381],[-121.66465348371033,37.020752188503081],
     [-121.64247097000728,37.001658514256235],[-121.57950895795178,36.946478562229061],[-121.54048543975076,36.919630706516728],[-121.51260729844474,36.901346802870705],
     [-121.48944549593645,36.884433460794298],[-121.45508229475919,36.86255855283656],[-121.42610334369285,36.840114879061275],[-121.39115617325845,36.815708963134853],
     [-121.35673589161618,36.793529056408488],[-121.30886245141422,36.764290491092368],[-121.26875476511875,36.739378826299514],[-121.23607120903258,36.715009966943455],
     [-121.1926215875213,36.683474164488587],[-121.16368275517387,36.66066564742988],[-121.13916541841322,36.6393384352711],[-121.08670053002239,36.59376134384371],
     [-121.06371695120017,36.575907615284507],[-121.04071876667788,36.554153106004833],[-121.00791832828488,36.519877354310154],[-120.97207309013523,36.485893948235002],
     [-120.91576174330518,36.431157243479788],[-120.85735784000937,36.369754311823044],[-120.82541506386243,36.33833146896319],[-120.79398206909289,36.309107397971445],
     [-120.76409649358581,36.278691639242709],[-120.72941491103179,36.246844661554348],[-120.68691544416824,36.207853475097693],[-120.65261850319879,36.176845516870685],
     [-120.57630445821555,36.10849694324088],[-120.55117422382883,36.085785007481036],[-120.51725375997969,36.05478276694754],[-120.48260325960816,36.022430244199256],
     [-120.44669696564423,35.99339198982949],[-120.41162105275768,35.964166488899593],[-120.37477761013653,35.933559258904154],[-120.34946117595816,35.910487135624123],
     [-120.32164544113507,35.884236198652957],[-120.28094720203441,35.849557704201459],[-120.25500781938225,35.825358236558714],[-120.22579703346602,35.788598829563568],
     [-120.20506054229693,35.754745248739425],[-120.18108807142428,35.726943855131026],[-120.14859639654591,35.695779235772477],[-120.12783449400163,35.674173828857363],
     [-120.10593586300009,35.653689383733365],[-120.08506317056102,35.633728585922995],[-120.04862523818482,35.59870591463266],[-120.0121292412505,35.563296615814409],
     [-119.96943818848163,35.521453648690397],[-119.93159767449893,35.483016977920613],[-119.89889777420293,35.446300194659898],[-119.87804209353507,35.42203524916863],
     [-119.84692681693385,35.39038972770112],[-119.80632396186282,35.351792300798877],[-119.76055950731412,35.304971252158566],[-119.72807176026733,35.272843231920774],
     [-119.69409829064866,35.241824187169698],[-119.66978766719469,35.218363290625462],[-119.64821603067983,35.197882826954526],[-119.61763954289975,35.171045968622082],
     [-119.55858565209328,35.125880316256769],[-119.53417894839137,35.106540426058871],[-119.50939570151095,35.084327995405658],[-119.47383732651698,35.056444555160567],
     [-119.43587316701834,35.026738764076697]
     ],
     "composition models":[{"model":"uniform", "compositions":[3], "fractions":[1]}]
    }
  ]
}

output:

-----------------------------------------------------------------------------
-- This is ASPECT, the Advanced Solver for Problems in Earth's ConvecTion.
--     . version 2.5.0-pre (main, 88736bc93)
--     . using deal.II 9.4.0
--     .       with 64 bit indices and vectorization level 3 (512 bits)
--     . using Trilinos 12.18.1
--     . using p4est 2.3.2
--     . running in OPTIMIZED mode
--     . running with 192 MPI processes
-----------------------------------------------------------------------------

Vectorization over 8 doubles = 512 bits (AVX512), VECTORIZATION_LEVEL=3
-----------------------------------------------------------------------------
-- For information on how to cite ASPECT, see:
--   https://aspect.geodynamics.org/citing.html?ver=2.5.0-pre&GWB=1&mf=1&sha=88736bc93&src=code
-----------------------------------------------------------------------------
Number of active cells: 524,288 (on 3 levels)
Number of degrees of freedom: 30,601,448 (12,879,555+549,153+4,293,185+4,293,185+4,293,185+4,293,185)

*** Timestep 0:  t=0 years, dt=0 years
   Solving temperature system... 0 iterations.
   Solving NA_crust system ... 0 iterations.
   Solving PC_crust system ... 0 iterations.
   Solving spharz system ... 0 iterations.
   Solving Stokes system... 145+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 1: 0.240263

+----------------------------------------------+------------+------------+
| Total wallclock time elapsed since start     |      54.1s |            |
|                                              |            |            |
| Section                          | no. calls |  wall time | % of total |
+----------------------------------+-----------+------------+------------+
| Assemble Stokes system rhs       |         1 |      1.53s |       2.8% |
| Assemble composition system      |         3 |      8.32s |        15% |
| Assemble temperature system      |         1 |      3.18s |       5.9% |
| Build Stokes preconditioner      |         1 |       2.2s |       4.1% |
| Build composition preconditioner |         3 |     0.803s |       1.5% |
| Build temperature preconditioner |         1 |      0.29s |      0.54% |
| Initialization                   |         1 |      1.46s |       2.7% |
| Setup dof systems                |         1 |      1.53s |       2.8% |
| Setup initial conditions         |         1 |      19.2s |        35% |
| Setup matrices                   |         1 |       1.4s |       2.6% |
| Solve Stokes system              |         1 |      9.54s |        18% |
| Solve composition system         |         3 |     0.106s |       0.2% |
| Solve temperature system         |         1 |    0.0506s |         0% |
+----------------------------------+-----------+------------+------------+

-- Total wallclock time elapsed including restarts: 54s
Number of active cells: 342,211 (on 4 levels)
Number of degrees of freedom: 22,353,663 (9,399,546+421,389+3,133,182+3,133,182+3,133,182+3,133,182)

*** Timestep 0:  t=0 years, dt=0 years
   Solving temperature system... 0 iterations.
   Solving NA_crust system ... 0 iterations.
   Solving PC_crust system ... 0 iterations.
   Solving spharz system ... 0 iterations.
   Solving Stokes system... 197+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 1: 0.227578

+----------------------------------------------+------------+------------+
| Total wallclock time elapsed since start     |       112s |            |
|                                              |            |            |
| Section                          | no. calls |  wall time | % of total |
+----------------------------------+-----------+------------+------------+
| Assemble Stokes system rhs       |         2 |      2.56s |       2.3% |
| Assemble composition system      |         6 |      13.8s |        12% |
| Assemble temperature system      |         2 |      5.29s |       4.7% |
| Build Stokes preconditioner      |         2 |      4.54s |       4.1% |
| Build composition preconditioner |         6 |      1.44s |       1.3% |
| Build temperature preconditioner |         2 |     0.453s |       0.4% |
| Initialization                   |         1 |      1.46s |       1.3% |
| Refine mesh structure, part 1    |         1 |      2.69s |       2.4% |
| Refine mesh structure, part 2    |         1 |     0.539s |      0.48% |
| Setup dof systems                |         2 |      2.85s |       2.5% |
| Setup initial conditions         |         2 |      41.2s |        37% |
| Setup matrices                   |         2 |       3.1s |       2.8% |
| Solve Stokes system              |         2 |      21.8s |        19% |
| Solve composition system         |         6 |     0.494s |      0.44% |
| Solve temperature system         |         2 |     0.179s |      0.16% |
+----------------------------------+-----------+------------+------------+

-- Total wallclock time elapsed including restarts: 112s
Number of active cells: 928,594 (on 5 levels)
Number of degrees of freedom: 63,686,356 (26,769,732+1,223,648+8,923,244+8,923,244+8,923,244+8,923,244)

*** Timestep 0:  t=0 years, dt=0 years
   Solving temperature system... 0 iterations.
   Solving NA_crust system ... 0 iterations.
   Solving PC_crust system ... 0 iterations.
   Solving spharz system ... 0 iterations.
   Solving Stokes system... 740+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 1: 0.12594

   Solving Stokes system... 1000+---------------------------------------------------------
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 (exception_message.str())' on rank 0 on processing: 

--------------------------------------------------------
An error occurred in line <2779> of file </work/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/source/utilities.cc> in function
    void aspect::Utilities::linear_solver_failed(const string&, const string&, const std::vector<dealii::SolverControl>&, const std::exception&, const MPI_Comm&, const string&)
The violated condition was: 
    false
Additional information: 
    The iterative Stokes solver in
    StokesMatrixFreeHandlerImplementation::solve did not converge.

    The initial residual was: 4.083311e+18
    The final residual is: 1.443057e+14
    The required residual for convergence is: 4.250178e+13
    See iofiles/gmg/gmg-test-gmg-v3/solver_history.txt for the full
    convergence history.

    The solver reported the following error:

    --------------------------------------------------------
    An error occurred in line <2779> of file
    </work/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/source/utilities.cc>
    in function
    void aspect::Utilities::linear_solver_failed(const string&, const
    string&, const std::vector<dealii::SolverControl>&, const
    std::exception&, const MPI_Comm&, const string&)
    The violated condition was:
    false
    Additional information:
    The iterative (bottom right) solver in
    BlockSchurGMGPreconditioner::vmult did not converge.

    The initial residual was: 1.213044e-01
    The final residual is: nan
    The required residual for convergence is: 1.213044e-07

    The solver reported the following error:

    --------------------------------------------------------
    An error occurred in line <1280> of file
    </work/06923/tg862222/stampede2/software/2022-11-01-dealii-9.4-candi-native-64b/deal.II-v9.4.0/include/deal.II/lac/solver_cg.h>

    in function
    void dealii::SolverCG<VectorType>::solve(const MatrixType&,
    VectorType&, const VectorType&, const PreconditionerType&) [with
    MatrixType = aspect::MatrixFreeStokesOperators::MassMatrixOperator<3,
    1, double>; PreconditionerType = dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double>,
    dealii::MGTransferMatrixFree<3, double> >; VectorType =
    dealii::LinearAlgebra::distributed::Vector<double>]
    The violated condition was:
    solver_state == SolverControl::success
    Additional information:
    Iterative method reported convergence failure in step 1. The residual
    in the last step was nan.

    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.

    Stacktrace:
    -----------
    #0
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    void
    dealii::SolverCG<dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>
    >::solve<aspect::MatrixFreeStokesOperators::MassMatrixOperator<3, 1,
    double>, dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >
    >(aspect::MatrixFreeStokesOperators::MassMatrixOperator<3, 1, double>
    const&, dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>&,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host> const&, dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >
    const&)
    #1
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    ) [0xb52bc9]
    #2
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    void
    dealii::SolverFGMRES<dealii::LinearAlgebra::distributed::BlockVector<double>

    >::solve<aspect::MatrixFreeStokesOperators::StokesOperator<3, 2,
    double>,
    aspect::internal::BlockSchurGMGPreconditioner<aspect::MatrixFreeStokesOperators::StokesOperator<3,

    2, double>, aspect::MatrixFreeStokesOperators::ABlockOperator<3, 2,
    double>, aspect::MatrixFreeStokesOperators::MassMatrixOperator<3, 1,
    double>, dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >,
    dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >
    > >(aspect::MatrixFreeStokesOperators::StokesOperator<3, 2, double>
    const&, dealii::LinearAlgebra::distributed::BlockVector<double>&,
    dealii::LinearAlgebra::distributed::BlockVector<double> const&,
    aspect::internal::BlockSchurGMGPreconditioner<aspect::MatrixFreeStokesOperators::StokesOperator<3,

    2, double>, aspect::MatrixFreeStokesOperators::ABlockOperator<3, 2,
    double>, aspect::MatrixFreeStokesOperators::MassMatrixOperator<3, 1,
    double>, dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >,
    dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >
    > const&)
    #3
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    aspect::StokesMatrixFreeHandlerImplementation<3, 2>::solve()
    #4
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    aspect::Simulator<3>::solve_stokes()
    #5
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    aspect::Simulator<3>::assemble_and_solve_stokes(bool, double*)
    #6
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    aspect::Simulator<3>::solve_single_advection_iterated_stokes()
    #7
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    aspect::Simulator<3>::solve_timestep()
    #8
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    aspect::Simulator<3>::run()
    #9
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    void run_simulator<3>(std::__cxx11::basic_string<char,
    std::char_traits<char>, std::allocator<char> > const&,
    std::__cxx11::basic_string<char, std::char_traits<char>,
    std::allocator<char> > const&, bool, bool, bool)
    #10
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    main
    --------------------------------------------------------

    Stacktrace:
    -----------
    #0
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    aspect::Utilities::linear_solver_failed(std::__cxx11::basic_string<char,
    std::char_traits<char>, std::allocator<char> > const&,
    std::__cxx11::basic_string<char, std::char_traits<char>,
    std::allocator<char> > const&, std::vector<dealii::SolverControl,
    std::allocator<dealii::SolverControl> > const&, std::exception const&,
    int const&, std::__cxx11::basic_string<char, std::char_traits<char>,
    std::allocator<char> > const&)
    #1
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    ) [0xb52d65]
    #2
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    void
    dealii::SolverFGMRES<dealii::LinearAlgebra::distributed::BlockVector<double>
    >::solve<aspect::MatrixFreeStokesOperators::StokesOperator<3, 2,
    double>,
    aspect::internal::BlockSchurGMGPreconditioner<aspect::MatrixFreeStokesOperators::StokesOperator<3,
    2, double>, aspect::MatrixFreeStokesOperators::ABlockOperator<3, 2,
    double>, aspect::MatrixFreeStokesOperators::MassMatrixOperator<3, 1,
    double>, dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >,
    dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >
    > >(aspect::MatrixFreeStokesOperators::StokesOperator<3, 2, double>
    const&, dealii::LinearAlgebra::distributed::BlockVector<double>&,
    dealii::LinearAlgebra::distributed::BlockVector<double> const&,
    aspect::internal::BlockSchurGMGPreconditioner<aspect::MatrixFreeStokesOperators::StokesOperator<3,
    2, double>, aspect::MatrixFreeStokesOperators::ABlockOperator<3, 2,
    double>, aspect::MatrixFreeStokesOperators::MassMatrixOperator<3, 1,
    double>, dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >,
    dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >
    > const&)
    #3
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    aspect::StokesMatrixFreeHandlerImplementation<3, 2>::solve()
    #4
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    aspect::Simulator<3>::solve_stokes()
    #5
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    aspect::Simulator<3>::assemble_and_solve_stokes(bool, double*)
    #6
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    aspect::Simulator<3>::solve_single_advection_iterated_stokes()
    #7
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    aspect::Simulator<3>::solve_timestep()
    #8
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    aspect::Simulator<3>::run()
    #9
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    void run_simulator<3>(std::__cxx11::basic_string<char,
    std::char_traits<char>, std::allocator<char> > const&,
    std::__cxx11::basic_string<char, std::char_traits<char>,
    std::allocator<char> > const&, bool, bool, bool)
    #10
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    main
    --------------------------------------------------------

Stacktrace:
-----------
#0  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect: aspect::Utilities::linear_solver_failed(std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, std::vector<dealii::SolverControl, std::allocator<dealii::SolverControl> > const&, std::exception const&, int const&, std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&)
#1  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect: aspect::StokesMatrixFreeHandlerImplementation<3, 2>::solve()
#2  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect: aspect::Simulator<3>::solve_stokes()
#3  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect: aspect::Simulator<3>::assemble_and_solve_stokes(bool, double*)
#4  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect: aspect::Simulator<3>::solve_single_advection_iterated_stokes()
#5  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect: aspect::Simulator<3>::solve_timestep()
#6  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect: aspect::Simulator<3>::run()
#7  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect: void run_simulator<3>(std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, bool, bool, bool)
#8  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect: main
--------------------------------------------------------

Aborting!
----------------------------------------------------
MFraters commented 1 year ago

Based on @anne-glerum's comment that it might have to do with viscosity contrasts, I explored the issue a bit more. I managed to get the issue to show on a bit under a million DoF's instead of around 63 million. Also, the setup is much simpler (see images below). Again, his is with the current main branch of aspect (using the world builder shipped with aspect) and it works when the line set Stokes solver type = block GMG is commented out.

gmg-stokes_convergence-issue-v5

prm file:

set World builder file = iofiles/gmg/gmg-test-v5.wb
set Output directory = iofiles/gmg/gmg-test-gmg-v5 #-gmg
set Dimension = 3
set CFL number              = 0.5
set Max nonlinear iterations = 4 
set End time = 0
set Nonlinear solver scheme = single Advection, iterated Stokes 
#set Nonlinear solver scheme = single Advection, iterated defect correction Stokes 
set Timing output frequency = 0
set Max nonlinear iterations in pre-refinement = 0
set Pressure normalization                 = surface
set Nonlinear solver tolerance =1e-4
set Resume computation = auto
set Adiabatic surface temperature = 273.

subsection Solver parameters
  subsection Newton solver parameters
    set Max pre-Newton nonlinear iterations = 100
    set Nonlinear Newton solver switch tolerance = 1e-5
    set Max Newton line search iterations = 0
    set Maximum linear Stokes solver tolerance = 1e-2
    set Use Newton residual scaling method = true
    set Use Newton failsafe = false 
    set Use Eisenstat Walker method for Picard iterations = false #true
  end
  subsection Stokes solver parameters
    set Maximum number of expensive Stokes solver steps =  5000    
    #set Number of cheap Stokes solver steps = 600 #80000
    set Number of cheap Stokes solver steps = 1000 #80000
    #set Linear solver tolerance = 1e-2
    set Stokes solver type = block GMG
    #set GMRES solver restart length = 100
  end
end

subsection Heating model
  set List of model names = adiabatic heating, shear heating 
end

subsection Boundary composition model
  set List of model names = initial composition
  set Fixed composition boundary indicators = 0,1,2,3,5,6,7,8,9
end

subsection Boundary temperature model
  set Fixed temperature boundary indicators   = 0,1,2,3,4,5,6,7,8,9
  set List of model names = initial temperature 
  subsection Initial temperature
    set Minimal temperature = 273.15
    set Maximal temperature = 4000
  end
end

subsection Discretization
  subsection Stabilization parameters
    set Use artificial viscosity smoothing = false 
    set Global composition maximum = 1.0 
    set Global composition minimum = 0.0 
    set Global temperature minimum = 273.15
    set alpha = 2
    set beta  = 0.78 #0.117 
   end
end

subsection Boundary velocity model
  set Prescribed velocity boundary indicators = uppereast: function, upperwest: function
  set Tangential velocity boundary indicators = top
  set Zero velocity boundary indicators       = bottom
  subsection Function
    set Coordinate system = spherical
    set Use spherical unit vectors = true
    set Variable names = r,phi,theta,t
    set Function expression = 0; -0.0000001;0 
  end
end

subsection Boundary traction model 
  set Prescribed traction boundary indicators =  lowernorth: initial lithostatic pressure, lowereast: initial lithostatic pressure, lowersouth: initial lithostatic pressure, lowerwest: initial lithostatic pressure, uppersouth: initial lithostatic pressure, uppernorth: initial lithostatic pressure 
  subsection Initial lithostatic pressure
    set Representative point = 6371e3,-114.01,39.01
  end
end

subsection Geometry model
  set Model name = chunk with lithosphere boundary indicators

  subsection Chunk with lithosphere boundary indicators
    set Chunk inner radius = 5571e3
    set Chunk outer radius = 6371e3
    set Chunk middle boundary radius = 6271e3
    set Chunk minimum latitude = 39 
    set Chunk maximum latitude = 46 
    set Chunk minimum longitude = -131 
    set Chunk maximum longitude = -114
    set Longitude repetitions = 15 #20 #32
    set Latitude repetitions = 8 #32
    set Inner chunk radius repetitions = 7
  end
end

subsection Material model
  set Material averaging = harmonic average # only viscosity  #project to Q1 only viscosity #harmonic average only viscosity  #harmonic average #project to Q1 only viscosity
  set Model name = visco plastic
  subsection Visco Plastic
        set Reference temperature = 273
        #set Reference viscosity = 1e20
        set Minimum viscosity = 1e19
        set Maximum viscosity = 1e24
        # the background rheologically becomes oceanic curst at depths lower than 80km, this is to hopefullly fix the boundaries. The densities remain the same for the open boundaries to remain working.
        set Phase transition depths = background:80e3|410e3|520e3|560e3|670e3|670e3|670e3|670e3, PC_crust:80e3|665e3|720e3, spharz: 410e3|520e3|560e3|670e3|670e3|670e3|670e3
        set Phase transition widths = background:5e3|5e3|5e3|5e3|10e3|5e3|5e3|5e3, PC_crust:5e3|5e3|5e3, spharz: 5e3|5e3|5e3|5e3|5e3|5e3|5e3
        set Phase transition temperatures = background:1662.0|1662.0|1662.0|1662.0|1662.0|1662.0|1662.0|1662.0, PC_crust:1173.0|1662.0|1662.0, spharz: 1662.0|1662.0|1662.0|1662.0|1662.0|1662.0|1662.0
        set Phase transition Clapeyron slopes = background:0.0|4e6|4.1e6|4e6|-2e6|4e6|-3.1e6|1.3e6, PC_crust:0.0|4e6|1.3e6, spharz: 4e6|4.1e6|4e6|-2e6|4e6|-3.1e6|1.3e6
        set Thermal diffusivities = 1.0e-6
        set Heat capacities = 1250.0
        set Densities = background: 3300.0|3300.0|3394.4|3442.1|3453.2|3617.6|3691.5|3774.7|3929.1, NA_crust:3000.0, PC_crust:3000.0|3540.0|3613.0|3871.7, spharz: 3235.0|3372.3|3441.7|3441.7|3680.8|3717.8|3759.4|3836.6 #,  total_strain:1e5, water:1e5
        #set Densities = background: 3300.0|3300|3300|3300|3300|3300|3300|3300, NA_crust:3100.0|3100.0|3100.0|3100, PC_crust:3100.0|3100|3100.0|3100, spharz: 3100|3100|3100|3100|3100|3100|3100|3100 , total_strain:1e5, water: 1e5
        set Thermal expansivities = 3.1e-5
        set Viscosity averaging scheme = harmonic
        set Viscous flow law = composite
        set Yield mechanism = drucker
        set Grain size = 1.0000e-02
        set Prefactors for diffusion creep =            background:8.5700e-28|7.1768e-15|7.1768e-15|7.1768e-15|7.1768e-15|2.4977e-19|2.4977e-19|2.4977e-19|2.4977e-19, NA_crust:8.5700e-28, PC_crust:8.5700e-28|7.1768e-15|2.4977e-19|2.4977e-19, spharz:7.1768e-15|7.1768e-15|7.1768e-15|7.1768e-15|2.4977e-19|2.4977e-19|2.4977e-19|2.4977e-19 
        set Grain size exponents for diffusion creep =  background:3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00, NA_crust:3.0000e+00, PC_crust:3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00, spharz:3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00 
        set Activation energies for diffusion creep =   background:3.7500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05, NA_crust:3.7500e+05, PC_crust:3.7500e+05|3.2500e+05|3.2500e+05|3.2500e+05, spharz:3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05 
        set Activation volumes for diffusion creep =    background:2.3000e-05|2.3000e-05|2.3000e-05|2.3000e-05|2.3000e-05|3.0000e-06|3.0000e-06|3.0000e-06|3.0000e-06, NA_crust:2.3000e-05, PC_crust:2.3000e-05|2.3000e-05|3.0000e-06|3.0000e-06, spharz:2.3000e-05|2.3000e-05|2.3000e-05|2.3000e-05|3.0000e-06|3.0000e-06|3.0000e-06|3.0000e-06 
        set Prefactors for dislocation creep =          background:8.5700e-28|4.4668e-16|4.4668e-16|4.4668e-16|4.4668e-16|5.0000e-32|5.0000e-32|5.0000e-32|5.0000e-32, NA_crust:8.5700e-28, PC_crust:8.5700e-28|4.4668e-16|5.0000e-32|5.0000e-32, spharz:4.4668e-16|4.4668e-16|4.4668e-16|4.4668e-16|5.0000e-32|5.0000e-32|5.0000e-32|5.0000e-32 
        set Stress exponents for dislocation creep =    background:4.0000e+00|3.5000e+00|3.5000e+00|3.5000e+00|3.5000e+00|1.0000e+00|1.0000e+00|1.0000e+00|1.0000e+00, NA_crust:3.5000e+00, PC_crust:4.0000e+00|3.5000e+00|1.0000e+00|1.0000e+00, spharz:3.5000e+00|3.5000e+00|3.5000e+00|3.5000e+00|1.0000e+00|1.0000e+00|1.0000e+00|1.0000e+00 
        set Activation energies for dislocation creep = background:2.2300e+05|4.8000e+05|4.8000e+05|4.8000e+05|4.8000e+05|0.0000e+00|0.0000e+00|0.0000e+00|0.0000e+00, NA_crust:2.2300e+05, PC_crust:2.2300e+05|4.8000e+05|0.0000e+00|0.0000e+00, spharz:4.8000e+05|4.8000e+05|4.8000e+05|4.8000e+05|0.0000e+00|0.0000e+00|0.0000e+00|0.0000e+00
        set Activation volumes for dislocation creep =  background:2.4000e-05|2.4000e-05|2.4000e-05|2.4000e-05|2.4000e-05|0.0000e+00|0.0000e+00|0.0000e+00|0.0000e+00, NA_crust:2.4000e-05, PC_crust:2.4000e-05|2.4000e-05|0.0000e+00|0.0000e+00, spharz:2.4000e-05|2.4000e-05|2.4000e-05|2.4000e-05|0.0000e+00|0.0000e+00|0.0000e+00|0.0000e+00 

        #set Manually define phase method crust = background: 0.0, PC_crust:1.3, spharz:0.0  # CDPT
        #set Manually define phase method pyrolite = background:1.0, PC_crust: 0.0, spharz:0.0  # CDPT
        #set Manually define phase method harzburgite = background:0.0, PC_crust: 0.0, spharz:1.0  # CDPT
        #set Constant viscosity prefactors      = background:10.,  NA_crust:1.,  PC_crust:1., spharz:1.,  total_strain:1, water:1 #,       100.0 , 1.0
        #set Constant viscosity prefactors      =2.5, 1., 1., 1., 1, 1 #,       100.0 , 1.0

    set Angles of internal friction = background: 45|45|0|0|0|0|0|0|0, NA_crust:45,  PC_crust:45|45|0|0, spharz:45|45|0|0|0|0|0|0 #, total_strain:0, water:0
    set Cohesions                   = background: 20e8|20.e8|20e8|20e8|20e8|20e8|20e8|20e8|20e8, NA_crust:20e6, PC_crust:10.e4|20e6|20e8|20e8, spharz:20e6|20e8|20e8|20e8|20e8|20e8|20e8|20e8 #, total_strain:1, water:1
  end
end
subsection Compositional fields
  set Number of fields = 3
  set Names of fields =  NA_crust, PC_crust, spharz 
  set Types of fields = chemical composition, chemical composition, chemical composition 
  set List of normalized fields = 0,1,2
end

subsection Initial temperature model
    set Model name = world builder
end

subsection Initial composition model
    set Model name = world builder
end

subsection Gravity model
  set Model name = radial constant
end

subsection Mesh refinement
  set Initial global refinement = 0
  set Initial adaptive refinement = 2
  set Time steps between mesh refinement = 5
  set Strategy = isosurfaces, thermal energy density #, minimum refinement function, maximum refinement function
   subsection Isosurfaces
                          #minref maxref  mintemp maxtemp
        set Isosurfaces = \
                          max,    max,   Temperature:   0 | 350; \ 
                          max,    max,   Temperature:   0 | 1300, PC_crust: 0.5 | 2; \ 
                          max-1,  max-1, Temperature:   350 | 400; \ 
                          max-2,  max-2, Temperature:   400 | 1300; \ 
                          max-1,    max,   PC_crust: 0.5 | 2; \ 
                          max-1,  max-1, spharz: 0.5 | 2; \ 
                          max-3,  max-3,  Temperature:1300 | 1660; \
                          min,    min,   Temperature:1660 | 3000; \

   end
end

subsection Postprocess
  set List of postprocessors = visualization, boundary pressures, pressure statistics, velocity statistics, spherical velocity statistics, temperature statistics, composition statistics, load balance statistics, memory statistics 
  subsection Visualization
    set List of output variables = depth, density, viscosity, strain rate,stress second invariant 
    set Time between graphical output = 12500 
    set Write higher order output = false 
    set Interpolate output = false
    set Write in background thread = true
  end
end

.wb file:

{
  // changed oceanic plate depth from 95e3 to 100e3 and changed oceanic plate velocity from 0.02 to 0.02
  "version":"0.4",
  "coordinate system":{"model":"spherical", "depth method":"begin segment"},
  //"coordinate system":{"model":"spherical", "depth method":"begin at end segment"},
  //"cross section":[[-128,50.5],[-113,46.5]], //[[-131,47],[-117,43]],
  "cross section":[[-131,47],[-117,43]],
  //"maximum distance between coordinates":0.01,
  "interpolation":"continuous monotone spline",
  "features":
  [
    {"model":"mantle layer", "name":"mantle", "min depth":100e3, "max depth":1200e3, "coordinates":[[-150,20],[-100,20],[-100,80],[-150,80]],
      "composition models":[{"model":"uniform", "min depth":100e3, "max depth":1200e3, "compositions":[4], "fractions":[400]}]
    },
    {"model":"continental plate", "name":"boundary", "min depth":-1e2, "max depth":120e3, "coordinates":[[-150,20],[-100,20],[-100,80],[-150,80]],
     "temperature models":[{"model":"linear", "min depth":-1e2, "max depth":100e3, "top temperature":273.15}]//,
      //"composition models":[{"model":"uniform", "min depth":-1e2, "max depth":15e3, "compositions":[0]}]
    },
    {"model":"continental plate", "name":"inner part", "min depth":-1e2, "max depth":120e3, "coordinates":
     [[-135,45], [-135,40],[-115,40],[-115,45]],
     "temperature models":[{"model":"adiabatic", "min depth":-1e2, "max depth":100e3}]
    },
    {"model":"oceanic plate", "name":"Pasific", "min depth":-1e2,"max depth":100e3, "coordinates":
     [[-130,45], [-130,40],[-125,40],[-125,45]],
          "temperature models":[{"model":"plate model", "min depth":-1e2, "max depth":100e3, "spreading velocity":0.01, "ridge coordinates":[[-129,41],[-129,44]]}],
     "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,0.5], "min depth":-1e2, "max depth":15e3},
                           {"model":"uniform", "compositions":[2], "min depth":15e3, "max depth":30e3}]
    },

    {"model":"continental plate", "name":"NA", "min depth":-1e2, "coordinates":
     [[-125,40],[-125,45],[-115,45],[-115,40]],
     "max depth":120e3,
     "temperature models":[{"model":"linear", "min depth":-1e2, "max depth":120e3, "bottom temperature":-1}],
     "composition models":[{"model":"uniform","compositions":[0], "min depth":-1e2, "max depth":30e3}//,
     //                           {"model":"uniform","compositions":[7], "min depth":30e3, "max depth":100e3}
    ]},
    {"model":"subducting plate", "name":"Subducting plate", "dip point":[135,45], "coordinates":
     [[-125,41],[-125,44]],
     "segments":[
       {"length":46.762e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[0.000,8.244],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=0.000km
       {"length":11.216e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[8.244,8.851],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=5.728km
       {"length":11.246e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[8.851,9.123],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=7.280km
       {"length":11.234e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[9.123,9.219],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=9.052km
       {"length":11.247e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[9.219,9.406],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=10.767km
       {"length":11.256e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[9.406,9.772],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=12.585km
       {"length":11.266e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[9.772,10.362],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=14.474km
       {"length":11.282e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[10.362,11.151],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=16.444km
       {"length":11.315e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[11.151,12.124],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=18.520km
       {"length":11.353e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[12.124,13.246],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=20.788km
       {"length":11.402e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[13.246,14.489],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=23.255km
       {"length":11.465e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[14.489,15.795],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=25.958km
       {"length":11.520e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[15.795,17.008],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=28.933km
       {"length":11.596e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[17.008,17.984],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,1.787], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=32.134km
       {"length":11.629e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[17.984,18.544],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,1.438], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=35.617km
       {"length":11.649e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[18.544,18.683],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,1.077], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=39.227km
       {"length":11.650e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[18.683,18.530],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,0.708], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=42.921km
       {"length":11.634e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[18.530,18.368],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,0.336], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=46.636km
       {"length":11.625e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[18.368,18.614]},  // depth=50.319km
       {"length":11.677e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[18.614,19.741]},  // depth=53.994km
       {"length":11.785e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[19.741,22.164]},  // depth=57.850km
       {"length":12.039e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[22.164,26.115]},  // depth=62.036km
       {"length":12.535e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[26.115,31.500]},  // depth=66.906km
       {"length":13.339e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[31.500,37.724]},  // depth=72.911km
       {"length":14.461e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[37.724,43.851]},  // depth=80.464km
       {"length":15.961e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[43.851,49.173]},  // depth=89.875km
       {"length":17.371e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[49.173,53.194]},  // depth=101.477km
       {"length":18.861e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[53.194,55.894]},  // depth=114.973km
       {"length":19.745e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[55.894,57.319]},  // depth=130.359km
       {"length":20.251e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[57.319,57.642]},  // depth=146.837km
       {"length":20.179e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[57.642,57.046]},  // depth=163.937km
       {"length":19.702e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[57.046,55.798]},  // depth=180.970km
       {"length":19.018e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[55.798,54.310]},  // depth=197.455km
       {"length":18.540e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[54.310,53.068]},  // depth=213.135km
       {"length":17.930e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[53.068,52.553]},  // depth=228.251km
       {"length":18.000e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[52.553,53.002]},  // depth=242.632km
       {"length":18.355e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[53.002,54.341]},  // depth=257.118km
       {"length":19.218e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[54.341,56.310]},  // depth=272.062km
       {"length":20.443e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[56.310,58.597]},  // depth=288.071km
       {"length":21.748e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[58.597,60.907]},  // depth=305.548km
       {"length":23.763e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[60.907,63.085]},  // depth=324.554km
       {"length":25.102e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[63.085,64.903]},  // depth=345.852km
       {"length":26.317e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[64.903,66.236]},
        {"length":0.0, "thickness":[300.0], "top truncation":[-100.0], "angle":[66.236]},
        {"length":0.0, "thickness":[300.0], "top truncation":[-100.0], "angle":[66.236]},
        {"length":0.0, "thickness":[300.0], "top truncation":[-100.0], "angle":[66.236]},
        {"length":0.0, "thickness":[300.0], "top truncation":[-100.0], "angle":[66.236]},
        {"length":0.0, "thickness":[300.0], "top truncation":[-100.0], "angle":[66.236]},
        {"length":0.0, "thickness":[300.0], "top truncation":[-100.0], "angle":[66.236]}
    ],
    "temperature models":[{"model":"plate model","adiabatic heating":true, "density":3300, "plate velocity":0.04, "max distance slab top":100e3}],
     "composition models":[{"model":"uniform", "compositions":[1,3],  "fractions":[1,0], "max distance slab top":15e3},
                           {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                           {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]
    }
  ]
}

output:

-----------------------------------------------------------------------------
-- This is ASPECT, the Advanced Solver for Problems in Earth's ConvecTion.
--     . version 2.5.0-pre (main, 88736bc93)
--     . using deal.II 9.4.0
--     .       with 64 bit indices and vectorization level 3 (512 bits)
--     . using Trilinos 12.18.1
--     . using p4est 2.3.2
--     . running in OPTIMIZED mode
--     . running with 192 MPI processes
-----------------------------------------------------------------------------

Vectorization over 8 doubles = 512 bits (AVX512), VECTORIZATION_LEVEL=3
-----------------------------------------------------------------------------
-- For information on how to cite ASPECT, see:
--   https://aspect.geodynamics.org/citing.html?ver=2.5.0-pre&GWB=1&mf=1&sha=88736bc93&src=code
-----------------------------------------------------------------------------
Number of active cells: 960 (on 1 levels)
Number of degrees of freedom: 64,009 (26,877+1,296+8,959+8,959+8,959+8,959)

*** Timestep 0:  t=0 years, dt=0 years
   Solving temperature system... 0 iterations.
   Solving NA_crust system ... 35 iterations.
   Solving PC_crust system ... 0 iterations.
   Solving spharz system ... 0 iterations.
   Solving Stokes system... 94+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 1: 0.0513504

+----------------------------------------------+------------+------------+
| Total wallclock time elapsed since start     |      2.15s |            |
|                                              |            |            |
| Section                          | no. calls |  wall time | % of total |
+----------------------------------+-----------+------------+------------+
| Assemble Stokes system rhs       |         1 |    0.0689s |       3.2% |
| Assemble composition system      |         3 |     0.092s |       4.3% |
| Assemble temperature system      |         1 |    0.0526s |       2.4% |
| Build Stokes preconditioner      |         1 |    0.0105s |      0.49% |
| Build composition preconditioner |         3 |   0.00218s |       0.1% |
| Build temperature preconditioner |         1 |   0.00927s |      0.43% |
| Initialization                   |         1 |     0.584s |        27% |
| Setup dof systems                |         1 |     0.239s |        11% |
| Setup initial conditions         |         1 |     0.257s |        12% |
| Setup matrices                   |         1 |     0.131s |       6.1% |
| Solve Stokes system              |         1 |     0.343s |        16% |
| Solve composition system         |         3 |    0.0545s |       2.5% |
| Solve temperature system         |         1 |    0.0217s |         1% |
+----------------------------------+-----------+------------+------------+

-- Total wallclock time elapsed including restarts: 2s
Number of active cells: 2,815 (on 2 levels)
Number of degrees of freedom: 189,480 (79,608+3,728+26,536+26,536+26,536+26,536)

*** Timestep 0:  t=0 years, dt=0 years
   Solving temperature system... 0 iterations.
   Solving NA_crust system ... 0 iterations.
   Solving PC_crust system ... 0 iterations.
   Solving spharz system ... 0 iterations.
   Solving Stokes system... 201+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 1: 0.051873

+----------------------------------------------+------------+------------+
| Total wallclock time elapsed since start     |      4.67s |            |
|                                              |            |            |
| Section                          | no. calls |  wall time | % of total |
+----------------------------------+-----------+------------+------------+
| Assemble Stokes system rhs       |         2 |     0.106s |       2.3% |
| Assemble composition system      |         6 |     0.226s |       4.8% |
| Assemble temperature system      |         2 |     0.114s |       2.4% |
| Build Stokes preconditioner      |         2 |    0.0373s |       0.8% |
| Build composition preconditioner |         6 |   0.00712s |      0.15% |
| Build temperature preconditioner |         2 |    0.0102s |      0.22% |
| Initialization                   |         1 |     0.584s |        13% |
| Refine mesh structure, part 1    |         1 |    0.0871s |       1.9% |
| Refine mesh structure, part 2    |         1 |    0.0369s |      0.79% |
| Setup dof systems                |         2 |     0.374s |         8% |
| Setup initial conditions         |         2 |     0.656s |        14% |
| Setup matrices                   |         2 |     0.264s |       5.7% |
| Solve Stokes system              |         2 |      1.55s |        33% |
| Solve composition system         |         6 |      0.09s |       1.9% |
| Solve temperature system         |         2 |    0.0334s |      0.71% |
+----------------------------------+-----------+------------+------------+

-- Total wallclock time elapsed including restarts: 5s
Number of active cells: 14,463 (on 3 levels)
Number of degrees of freedom: 931,468 (391,662+17,590+130,554+130,554+130,554+130,554)

*** Timestep 0:  t=0 years, dt=0 years
   Solving temperature system... 0 iterations.
   Solving NA_crust system ... 0 iterations.
   Solving PC_crust system ... 0 iterations.
   Solving spharz system ... 0 iterations.
   Solving Stokes system... 498+0 iterations.
      Relative nonlinear residual (Stokes system) after nonlinear iteration 1: 0.051838

   Solving Stokes system... 1000+---------------------------------------------------------
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 (exception_message.str())' on rank 0 on processing: 

--------------------------------------------------------
An error occurred in line <2779> of file </work/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/source/utilities.cc> in function
    void aspect::Utilities::linear_solver_failed(const string&, const string&, const std::vector<dealii::SolverControl>&, const std::exception&, const MPI_Comm&, const string&)
The violated condition was: 
    false
Additional information: 
    The iterative Stokes solver in
    StokesMatrixFreeHandlerImplementation::solve did not converge.

    The initial residual was: 1.813766e+18
    The final residual is: 6.478978e+12
    The required residual for convergence is: 1.917930e+11
    See iofiles/gmg/gmg-test-gmg-v5/solver_history.txt for the full
    convergence history.

    The solver reported the following error:

    --------------------------------------------------------
    An error occurred in line <2779> of file
    </work/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/source/utilities.cc>
    in function
    void aspect::Utilities::linear_solver_failed(const string&, const
    string&, const std::vector<dealii::SolverControl>&, const
    std::exception&, const MPI_Comm&, const string&)
    The violated condition was:
    false
    Additional information:
    The iterative (bottom right) solver in
    BlockSchurGMGPreconditioner::vmult did not converge.

    The initial residual was: 5.803082e-01
    The final residual is: nan
    The required residual for convergence is: 5.803082e-07

    The solver reported the following error:

    --------------------------------------------------------
    An error occurred in line <1280> of file
    </work/06923/tg862222/stampede2/software/2022-11-01-dealii-9.4-candi-native-64b/deal.II-v9.4.0/include/deal.II/lac/solver_cg.h>

    in function
    void dealii::SolverCG<VectorType>::solve(const MatrixType&,
    VectorType&, const VectorType&, const PreconditionerType&) [with
    MatrixType = aspect::MatrixFreeStokesOperators::MassMatrixOperator<3,
    1, double>; PreconditionerType = dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double>,
    dealii::MGTransferMatrixFree<3, double> >; VectorType =
    dealii::LinearAlgebra::distributed::Vector<double>]
    The violated condition was:
    solver_state == SolverControl::success
    Additional information:
    Iterative method reported convergence failure in step 1. The residual
    in the last step was nan.

    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.

    Stacktrace:
    -----------
    #0
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    void
    dealii::SolverCG<dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>
    >::solve<aspect::MatrixFreeStokesOperators::MassMatrixOperator<3, 1,
    double>, dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >
    >(aspect::MatrixFreeStokesOperators::MassMatrixOperator<3, 1, double>
    const&, dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>&,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host> const&, dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >
    const&)
    #1
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    ) [0xb52bc9]
    #2
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    void
    dealii::SolverFGMRES<dealii::LinearAlgebra::distributed::BlockVector<double>

    >::solve<aspect::MatrixFreeStokesOperators::StokesOperator<3, 2,
    double>,
    aspect::internal::BlockSchurGMGPreconditioner<aspect::MatrixFreeStokesOperators::StokesOperator<3,

    2, double>, aspect::MatrixFreeStokesOperators::ABlockOperator<3, 2,
    double>, aspect::MatrixFreeStokesOperators::MassMatrixOperator<3, 1,
    double>, dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >,
    dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >
    > >(aspect::MatrixFreeStokesOperators::StokesOperator<3, 2, double>
    const&, dealii::LinearAlgebra::distributed::BlockVector<double>&,
    dealii::LinearAlgebra::distributed::BlockVector<double> const&,
    aspect::internal::BlockSchurGMGPreconditioner<aspect::MatrixFreeStokesOperators::StokesOperator<3,

    2, double>, aspect::MatrixFreeStokesOperators::ABlockOperator<3, 2,
    double>, aspect::MatrixFreeStokesOperators::MassMatrixOperator<3, 1,
    double>, dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >,
    dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >
    > const&)
    #3
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    aspect::StokesMatrixFreeHandlerImplementation<3, 2>::solve()
    #4
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    aspect::Simulator<3>::solve_stokes()
    #5
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    aspect::Simulator<3>::assemble_and_solve_stokes(bool, double*)
    #6
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    aspect::Simulator<3>::solve_single_advection_iterated_stokes()
    #7
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    aspect::Simulator<3>::solve_timestep()
    #8
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    aspect::Simulator<3>::run()
    #9
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    void run_simulator<3>(std::__cxx11::basic_string<char,
    std::char_traits<char>, std::allocator<char> > const&,
    std::__cxx11::basic_string<char, std::char_traits<char>,
    std::allocator<char> > const&, bool, bool, bool)
    #10
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:

    main
    --------------------------------------------------------

    Stacktrace:
    -----------
    #0
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    aspect::Utilities::linear_solver_failed(std::__cxx11::basic_string<char,
    std::char_traits<char>, std::allocator<char> > const&,
    std::__cxx11::basic_string<char, std::char_traits<char>,
    std::allocator<char> > const&, std::vector<dealii::SolverControl,
    std::allocator<dealii::SolverControl> > const&, std::exception const&,
    int const&, std::__cxx11::basic_string<char, std::char_traits<char>,
    std::allocator<char> > const&)
    #1
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    ) [0xb52d65]
    #2
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    void
    dealii::SolverFGMRES<dealii::LinearAlgebra::distributed::BlockVector<double>
    >::solve<aspect::MatrixFreeStokesOperators::StokesOperator<3, 2,
    double>,
    aspect::internal::BlockSchurGMGPreconditioner<aspect::MatrixFreeStokesOperators::StokesOperator<3,
    2, double>, aspect::MatrixFreeStokesOperators::ABlockOperator<3, 2,
    double>, aspect::MatrixFreeStokesOperators::MassMatrixOperator<3, 1,
    double>, dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >,
    dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >
    > >(aspect::MatrixFreeStokesOperators::StokesOperator<3, 2, double>
    const&, dealii::LinearAlgebra::distributed::BlockVector<double>&,
    dealii::LinearAlgebra::distributed::BlockVector<double> const&,
    aspect::internal::BlockSchurGMGPreconditioner<aspect::MatrixFreeStokesOperators::StokesOperator<3,
    2, double>, aspect::MatrixFreeStokesOperators::ABlockOperator<3, 2,
    double>, aspect::MatrixFreeStokesOperators::MassMatrixOperator<3, 1,
    double>, dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >,
    dealii::PreconditionMG<3,
    dealii::LinearAlgebra::distributed::Vector<double,
    dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >
    > const&)
    #3
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    aspect::StokesMatrixFreeHandlerImplementation<3, 2>::solve()
    #4
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    aspect::Simulator<3>::solve_stokes()
    #5
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    aspect::Simulator<3>::assemble_and_solve_stokes(bool, double*)
    #6
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    aspect::Simulator<3>::solve_single_advection_iterated_stokes()
    #7
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    aspect::Simulator<3>::solve_timestep()
    #8
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    aspect::Simulator<3>::run()
    #9
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    void run_simulator<3>(std::__cxx11::basic_string<char,
    std::char_traits<char>, std::allocator<char> > const&,
    std::__cxx11::basic_string<char, std::char_traits<char>,
    std::allocator<char> > const&, bool, bool, bool)
    #10
    /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect:
    main
    --------------------------------------------------------

Stacktrace:
-----------
#0  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect: aspect::Utilities::linear_solver_failed(std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, std::vector<dealii::SolverControl, std::allocator<dealii::SolverControl> > const&, std::exception const&, int const&, std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&)
#1  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect: aspect::StokesMatrixFreeHandlerImplementation<3, 2>::solve()
#2  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect: aspect::Simulator<3>::solve_stokes()
#3  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect: aspect::Simulator<3>::assemble_and_solve_stokes(bool, double*)
#4  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect: aspect::Simulator<3>::solve_single_advection_iterated_stokes()
#5  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect: aspect::Simulator<3>::solve_timestep()
#6  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect: aspect::Simulator<3>::run()
#7  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect: void run_simulator<3>(std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, bool, bool, bool)
#8  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-rl-main/aspect: main
--------------------------------------------------------

Aborting!
----------------------------------------------------
tjhei commented 1 year ago

When I run it I get

Assert `add_angle_inner >= 0 && add_angle_inner <= 1` failed in ../contrib/world_builder/source/utilities.cc at line 894: Internal error: The variable add_angle_inner is smaller than zero or larger then one,which causes the std::acos to return nan. If it is only a little bit larger then one, this is probably caused by that begin and end segment are the same and round off error. The value of add_angle_inner = 1

It runs and fails to converge in release mode. I suspect the combination of traction and prescribed flow conditions might be the problem, but I have to play with this some more. Do you think you can change the problem to not use traction boundaries?

MFraters commented 1 year ago

ah, sorry, I forgot to test the simplfied version in debug mode. I only had to remove the 0 length segments and now it works for me in debug mode (see below).

Thanks for looking into this. I will try to make a case without traction boundaries. I might be able to simplify the slab a bit more as well with less segments so that the file is a bit shorter, I will look into that.

{
  // changed oceanic plate depth from 95e3 to 100e3 and changed oceanic plate velocity from 0.02 to 0.02
  "version":"0.4",
  "coordinate system":{"model":"spherical", "depth method":"begin segment"},
  //"coordinate system":{"model":"spherical", "depth method":"begin at end segment"},
  //"cross section":[[-128,50.5],[-113,46.5]], //[[-131,47],[-117,43]],
  "cross section":[[-131,47],[-117,43]],
  //"maximum distance between coordinates":0.01,
  "interpolation":"continuous monotone spline",
  "features":
  [
    {"model":"mantle layer", "name":"mantle", "min depth":100e3, "max depth":1200e3, "coordinates":[[-150,20],[-100,20],[-100,80],[-150,80]],
      "composition models":[{"model":"uniform", "min depth":100e3, "max depth":1200e3, "compositions":[4], "fractions":[400]}]
    },
    {"model":"continental plate", "name":"boundary", "min depth":-1e2, "max depth":120e3, "coordinates":[[-150,20],[-100,20],[-100,80],[-150,80]],
     "temperature models":[{"model":"linear", "min depth":-1e2, "max depth":100e3, "top temperature":273.15}]//,
      //"composition models":[{"model":"uniform", "min depth":-1e2, "max depth":15e3, "compositions":[0]}]
    },
    {"model":"continental plate", "name":"inner part", "min depth":-1e2, "max depth":120e3, "coordinates":
     [[-135,45], [-135,40],[-115,40],[-115,45]],
     "temperature models":[{"model":"adiabatic", "min depth":-1e2, "max depth":100e3}]
    },
    {"model":"oceanic plate", "name":"Pasific", "min depth":-1e2,"max depth":100e3, "coordinates":
     [[-130,45], [-130,40],[-125,40],[-125,45]],
          "temperature models":[{"model":"plate model", "min depth":-1e2, "max depth":100e3, "spreading velocity":0.01, "ridge coordinates":[[-129,41],[-129,44]]}],
     "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,0.5], "min depth":-1e2, "max depth":15e3},
                           {"model":"uniform", "compositions":[2], "min depth":15e3, "max depth":30e3}]
    },

    {"model":"continental plate", "name":"NA", "min depth":-1e2, "coordinates":
     [[-125,40],[-125,45],[-115,45],[-115,40]],
     "max depth":120e3,
     "temperature models":[{"model":"linear", "min depth":-1e2, "max depth":120e3, "bottom temperature":-1}],
     "composition models":[{"model":"uniform","compositions":[0], "min depth":-1e2, "max depth":30e3}//,
     //                           {"model":"uniform","compositions":[7], "min depth":30e3, "max depth":100e3}
    ]},
    {"model":"subducting plate", "name":"Subducting plate", "dip point":[135,45], "coordinates":
     [[-125,41],[-125,44]],
     "segments":[
       {"length":46.762e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[0.000,8.244],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=0.000km
       {"length":11.216e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[8.244,8.851],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=5.728km
       {"length":11.246e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[8.851,9.123],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=7.280km
       {"length":11.234e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[9.123,9.219],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=9.052km
       {"length":11.247e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[9.219,9.406],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=10.767km
       {"length":11.256e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[9.406,9.772],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=12.585km
       {"length":11.266e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[9.772,10.362],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=14.474km
       {"length":11.282e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[10.362,11.151],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=16.444km
       {"length":11.315e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[11.151,12.124],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=18.520km
       {"length":11.353e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[12.124,13.246],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=20.788km
       {"length":11.402e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[13.246,14.489],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=23.255km
       {"length":11.465e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[14.489,15.795],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=25.958km
       {"length":11.520e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[15.795,17.008],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,2.0], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=28.933km
       {"length":11.596e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[17.008,17.984],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,1.787], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=32.134km
       {"length":11.629e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[17.984,18.544],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,1.438], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=35.617km
       {"length":11.649e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[18.544,18.683],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,1.077], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=39.227km
       {"length":11.650e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[18.683,18.530],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,0.708], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=42.921km
       {"length":11.634e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[18.530,18.368],
        "composition models":[{"model":"uniform", "compositions":[1,3], "fractions":[1,0.336], "max distance slab top":15e3},
                              {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                              {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]},  // depth=46.636km
       {"length":11.625e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[18.368,18.614]},  // depth=50.319km
       {"length":11.677e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[18.614,19.741]},  // depth=53.994km
       {"length":11.785e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[19.741,22.164]},  // depth=57.850km
       {"length":12.039e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[22.164,26.115]},  // depth=62.036km
       {"length":12.535e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[26.115,31.500]},  // depth=66.906km
       {"length":13.339e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[31.500,37.724]},  // depth=72.911km
       {"length":14.461e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[37.724,43.851]},  // depth=80.464km
       {"length":15.961e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[43.851,49.173]},  // depth=89.875km
       {"length":17.371e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[49.173,53.194]},  // depth=101.477km
       {"length":18.861e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[53.194,55.894]},  // depth=114.973km
       {"length":19.745e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[55.894,57.319]},  // depth=130.359km
       {"length":20.251e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[57.319,57.642]},  // depth=146.837km
       {"length":20.179e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[57.642,57.046]},  // depth=163.937km
       {"length":19.702e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[57.046,55.798]},  // depth=180.970km
       {"length":19.018e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[55.798,54.310]},  // depth=197.455km
       {"length":18.540e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[54.310,53.068]},  // depth=213.135km
       {"length":17.930e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[53.068,52.553]},  // depth=228.251km
       {"length":18.000e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[52.553,53.002]},  // depth=242.632km
       {"length":18.355e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[53.002,54.341]},  // depth=257.118km
       {"length":19.218e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[54.341,56.310]},  // depth=272.062km
       {"length":20.443e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[56.310,58.597]},  // depth=288.071km
       {"length":21.748e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[58.597,60.907]},  // depth=305.548km
       {"length":23.763e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[60.907,63.085]},  // depth=324.554km
       {"length":25.102e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[63.085,64.903]},  // depth=345.852km
       {"length":26.317e03, "thickness":[300.0e03], "top truncation":[-100.000e03], "angle":[64.903,66.236]}
    ],
    "temperature models":[{"model":"plate model","adiabatic heating":true, "density":3300, "plate velocity":0.04, "max distance slab top":100e3}],
     "composition models":[{"model":"uniform", "compositions":[1,3],  "fractions":[1,0], "max distance slab top":15e3},
                           {"model":"uniform", "compositions":[0,2], "fractions":[0,1], "min distance slab top":15e3, "max distance slab top":30e3},
                           {"model":"uniform", "compositions":[0], "fractions":[0], "min distance slab top":30e3, "max distance slab top":100e3}]
    }
  ]
}
tjhei commented 1 year ago

Don't worry about simplifying the WB file. I am pretty sure there is an issue related to boundary conditions. If the system can not preserve mass due to the boundary conditions, the GMG will fail (while AMG might still work).

MFraters commented 1 year ago

I have simplified the problem a bit more, now it has a chunk instead of a chunk with lithospheric boundary conditions and no traction boundaries anymore. Stokes works on this, but gmg not (in debug mode), with with a different error message (see below)

set World builder file = iofiles/gmg/gmg-test-v6.wb
set Output directory = iofiles/gmg/gmg-test-gmg-v7 #-gmg
set Dimension = 3
set CFL number              = 0.5
set Max nonlinear iterations = 4 
set End time = 0
set Nonlinear solver scheme = single Advection, iterated Stokes 
#set Nonlinear solver scheme = single Advection, iterated defect correction Stokes 
set Timing output frequency = 0
set Max nonlinear iterations in pre-refinement = 0
set Pressure normalization                 = surface
set Nonlinear solver tolerance =1e-4
set Resume computation = auto
set Adiabatic surface temperature = 273.

subsection Solver parameters
  subsection Newton solver parameters
    set Max pre-Newton nonlinear iterations = 100
    set Nonlinear Newton solver switch tolerance = 1e-5
    set Max Newton line search iterations = 0
    set Maximum linear Stokes solver tolerance = 1e-2
    set Use Newton residual scaling method = true
    set Use Newton failsafe = false 
    set Use Eisenstat Walker method for Picard iterations = false #true
  end
  subsection Stokes solver parameters
    set Maximum number of expensive Stokes solver steps =  5000    
    #set Number of cheap Stokes solver steps = 600 #80000
    set Number of cheap Stokes solver steps = 1000 #80000
    #set Linear solver tolerance = 1e-2
    set Stokes solver type = block GMG
    #set GMRES solver restart length = 100
  end
end

subsection Heating model
  set List of model names = adiabatic heating, shear heating 
end

subsection Boundary composition model
  set List of model names = initial composition
  set Fixed composition boundary indicators = 0,1,2,3,5 #,6,7,8,9
end

subsection Boundary temperature model
  set Fixed temperature boundary indicators   = 0,1,2,3,4,5 #,6,7,8,9
  set List of model names = initial temperature 
  subsection Initial temperature
    set Minimal temperature = 273.15
    set Maximal temperature = 4000
  end
end

subsection Discretization
  subsection Stabilization parameters
    set Use artificial viscosity smoothing = false 
    set Global composition maximum = 1.0 
    set Global composition minimum = 0.0 
    set Global temperature minimum = 273.15
    set alpha = 2
    set beta  = 0.78 #0.117 
   end
end

subsection Boundary velocity model
  #set Prescribed velocity boundary indicators = uppereast: function, upperwest: function
  set Tangential velocity boundary indicators = top, north, east, south, west 
  #set Tangential velocity boundary indicators = top, lowernorth, lowereast, lowersouth, lowerwest, uppersouth, uppernorth, uppereast, upperwest
  set Zero velocity boundary indicators       = bottom
  subsection Function
    set Coordinate system = spherical
    set Use spherical unit vectors = true
    set Variable names = r,phi,theta,t
    set Function expression = 0; -0.0000001;0 
  end
end

subsection Boundary traction model 
  #set Prescribed traction boundary indicators =  lowernorth: initial lithostatic pressure, lowereast: initial lithostatic pressure, lowersouth: initial lithostatic pressure, lowerwest: initial lithostatic pressure, uppersouth: initial lithostatic pressure, uppernorth: initial lithostatic pressure 
  subsection Initial lithostatic pressure
    set Representative point = 6371e3,-114.01,39.01
  end
end

subsection Geometry model
  set Model name = chunk #with lithosphere boundary indicators

  subsection Chunk #with lithosphere boundary indicators
    set Chunk inner radius = 5571e3
    set Chunk outer radius = 6371e3
    set Chunk minimum latitude = 39 
    set Chunk maximum latitude = 46 
    set Chunk minimum longitude = -131 
    set Chunk maximum longitude = -114
    set Longitude repetitions = 15 #20 #32
    set Latitude repetitions = 8 #32
    set Radius repetitions = 7
  end
  subsection Chunk with lithosphere boundary indicators
    set Chunk inner radius = 5571e3
    set Chunk outer radius = 6371e3
    set Chunk middle boundary radius = 6271e3
    set Chunk minimum latitude = 39 
    set Chunk maximum latitude = 46 
    set Chunk minimum longitude = -131 
    set Chunk maximum longitude = -114
    set Longitude repetitions = 15 #20 #32
    set Latitude repetitions = 8 #32
    set Inner chunk radius repetitions = 7
  end
end

subsection Material model
  set Material averaging = harmonic average # only viscosity  #project to Q1 only viscosity #harmonic average only viscosity  #harmonic average #project to Q1 only viscosity
  set Model name = visco plastic
  subsection Visco Plastic
        set Reference temperature = 273
        #set Reference viscosity = 1e20
        set Minimum viscosity = 1e19
        set Maximum viscosity = 1e24
        # the background rheologically becomes oceanic curst at depths lower than 80km, this is to hopefullly fix the boundaries. The densities remain the same for the open boundaries to remain working.
        set Phase transition depths = background:80e3|410e3|520e3|560e3|670e3|670e3|670e3|670e3, PC_crust:80e3|665e3|720e3, spharz: 410e3|520e3|560e3|670e3|670e3|670e3|670e3
        set Phase transition widths = background:5e3|5e3|5e3|5e3|10e3|5e3|5e3|5e3, PC_crust:5e3|5e3|5e3, spharz: 5e3|5e3|5e3|5e3|5e3|5e3|5e3
        set Phase transition temperatures = background:1662.0|1662.0|1662.0|1662.0|1662.0|1662.0|1662.0|1662.0, PC_crust:1173.0|1662.0|1662.0, spharz: 1662.0|1662.0|1662.0|1662.0|1662.0|1662.0|1662.0
        set Phase transition Clapeyron slopes = background:0.0|4e6|4.1e6|4e6|-2e6|4e6|-3.1e6|1.3e6, PC_crust:0.0|4e6|1.3e6, spharz: 4e6|4.1e6|4e6|-2e6|4e6|-3.1e6|1.3e6
        set Thermal diffusivities = 1.0e-6
        set Heat capacities = 1250.0
        set Densities = background: 3300.0|3300.0|3394.4|3442.1|3453.2|3617.6|3691.5|3774.7|3929.1, NA_crust:3000.0, PC_crust:3000.0|3540.0|3613.0|3871.7, spharz: 3235.0|3372.3|3441.7|3441.7|3680.8|3717.8|3759.4|3836.6 #,  total_strain:1e5, water:1e5
        #set Densities = background: 3300.0|3300|3300|3300|3300|3300|3300|3300, NA_crust:3100.0|3100.0|3100.0|3100, PC_crust:3100.0|3100|3100.0|3100, spharz: 3100|3100|3100|3100|3100|3100|3100|3100 , total_strain:1e5, water: 1e5
        set Thermal expansivities = 3.1e-5
        set Viscosity averaging scheme = harmonic
        set Viscous flow law = composite
        set Yield mechanism = drucker
        set Grain size = 1.0000e-02
        set Prefactors for diffusion creep =            background:8.5700e-28|7.1768e-15|7.1768e-15|7.1768e-15|7.1768e-15|2.4977e-19|2.4977e-19|2.4977e-19|2.4977e-19, NA_crust:8.5700e-28, PC_crust:8.5700e-28|7.1768e-15|2.4977e-19|2.4977e-19, spharz:7.1768e-15|7.1768e-15|7.1768e-15|7.1768e-15|2.4977e-19|2.4977e-19|2.4977e-19|2.4977e-19 
        set Grain size exponents for diffusion creep =  background:3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00, NA_crust:3.0000e+00, PC_crust:3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00, spharz:3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00 
        set Activation energies for diffusion creep =   background:3.7500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05, NA_crust:3.7500e+05, PC_crust:3.7500e+05|3.2500e+05|3.2500e+05|3.2500e+05, spharz:3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05 
        set Activation volumes for diffusion creep =    background:2.3000e-05|2.3000e-05|2.3000e-05|2.3000e-05|2.3000e-05|3.0000e-06|3.0000e-06|3.0000e-06|3.0000e-06, NA_crust:2.3000e-05, PC_crust:2.3000e-05|2.3000e-05|3.0000e-06|3.0000e-06, spharz:2.3000e-05|2.3000e-05|2.3000e-05|2.3000e-05|3.0000e-06|3.0000e-06|3.0000e-06|3.0000e-06 
        set Prefactors for dislocation creep =          background:8.5700e-28|4.4668e-16|4.4668e-16|4.4668e-16|4.4668e-16|5.0000e-32|5.0000e-32|5.0000e-32|5.0000e-32, NA_crust:8.5700e-28, PC_crust:8.5700e-28|4.4668e-16|5.0000e-32|5.0000e-32, spharz:4.4668e-16|4.4668e-16|4.4668e-16|4.4668e-16|5.0000e-32|5.0000e-32|5.0000e-32|5.0000e-32 
        set Stress exponents for dislocation creep =    background:4.0000e+00|3.5000e+00|3.5000e+00|3.5000e+00|3.5000e+00|1.0000e+00|1.0000e+00|1.0000e+00|1.0000e+00, NA_crust:3.5000e+00, PC_crust:4.0000e+00|3.5000e+00|1.0000e+00|1.0000e+00, spharz:3.5000e+00|3.5000e+00|3.5000e+00|3.5000e+00|1.0000e+00|1.0000e+00|1.0000e+00|1.0000e+00 
        set Activation energies for dislocation creep = background:2.2300e+05|4.8000e+05|4.8000e+05|4.8000e+05|4.8000e+05|0.0000e+00|0.0000e+00|0.0000e+00|0.0000e+00, NA_crust:2.2300e+05, PC_crust:2.2300e+05|4.8000e+05|0.0000e+00|0.0000e+00, spharz:4.8000e+05|4.8000e+05|4.8000e+05|4.8000e+05|0.0000e+00|0.0000e+00|0.0000e+00|0.0000e+00
        set Activation volumes for dislocation creep =  background:2.4000e-05|2.4000e-05|2.4000e-05|2.4000e-05|2.4000e-05|0.0000e+00|0.0000e+00|0.0000e+00|0.0000e+00, NA_crust:2.4000e-05, PC_crust:2.4000e-05|2.4000e-05|0.0000e+00|0.0000e+00, spharz:2.4000e-05|2.4000e-05|2.4000e-05|2.4000e-05|0.0000e+00|0.0000e+00|0.0000e+00|0.0000e+00 

        #set Manually define phase method crust = background: 0.0, PC_crust:1.3, spharz:0.0  # CDPT
        #set Manually define phase method pyrolite = background:1.0, PC_crust: 0.0, spharz:0.0  # CDPT
        #set Manually define phase method harzburgite = background:0.0, PC_crust: 0.0, spharz:1.0  # CDPT
        #set Constant viscosity prefactors      = background:10.,  NA_crust:1.,  PC_crust:1., spharz:1.,  total_strain:1, water:1 #,       100.0 , 1.0
        #set Constant viscosity prefactors      =2.5, 1., 1., 1., 1, 1 #,       100.0 , 1.0

    set Angles of internal friction = background: 45|45|0|0|0|0|0|0|0, NA_crust:45,  PC_crust:45|45|0|0, spharz:45|45|0|0|0|0|0|0 #, total_strain:0, water:0
    set Cohesions                   = background: 20e8|20.e8|20e8|20e8|20e8|20e8|20e8|20e8|20e8, NA_crust:20e6, PC_crust:10.e4|20e6|20e8|20e8, spharz:20e6|20e8|20e8|20e8|20e8|20e8|20e8|20e8 #, total_strain:1, water:1
  end
end
subsection Compositional fields
  set Number of fields = 3
  set Names of fields =  NA_crust, PC_crust, spharz 
  set Types of fields = chemical composition, chemical composition, chemical composition 
  set List of normalized fields = 0,1,2
end

subsection Initial temperature model
    set Model name = world builder
end

subsection Initial composition model
    set Model name = world builder
end

subsection Gravity model
  set Model name = radial constant
end

subsection Mesh refinement
  set Initial global refinement = 0
  set Initial adaptive refinement = 2
  set Time steps between mesh refinement = 5
  set Strategy = isosurfaces, thermal energy density #, minimum refinement function, maximum refinement function
   subsection Isosurfaces
                          #minref maxref  mintemp maxtemp
        set Isosurfaces = \
                          max,    max,   Temperature:   0 | 350; \ 
                          max,    max,   Temperature:   0 | 1300, PC_crust: 0.5 | 2; \ 
                          max-1,  max-1, Temperature:   350 | 400; \ 
                          max-2,  max-2, Temperature:   400 | 1300; \ 
                          max-1,    max,   PC_crust: 0.5 | 2; \ 
                          max-1,  max-1, spharz: 0.5 | 2; \ 
                          max-3,  max-3,  Temperature:1300 | 1660; \
                          min,    min,   Temperature:1660 | 3000; \

   end
end

subsection Postprocess
  set List of postprocessors = visualization, boundary pressures, pressure statistics, velocity statistics, spherical velocity statistics, temperature statistics, composition statistics, load balance statistics, memory statistics 
  subsection Visualization
    set List of output variables = depth, density, viscosity, strain rate,stress second invariant 
    set Time between graphical output = 12500 
    set Write higher order output = false 
    set Interpolate output = false
    set Write in background thread = true
  end
end
--------------------------------------------------------

Calling MPI_Abort now.
To break execution in a GDB session, execute 'break MPI_Abort' before running. You can also put the following into your ~/.gdbinit:
  set breakpoint pending on
  break MPI_Abort
  set breakpoint pending auto
application called MPI_Abort(MPI_COMM_WORLD, 255) - process 39

--------------------------------------------------------
An error occurred in line <1910> of file </work2/06923/tg862222/stampede2/software/2022-11-01-dealii-9.4-candi-native-64b/tmp/unpack/deal.II-v9.4.0/include/deal.II/lac/vector_operations_internal.h> in function
    static Number dealii::internal::VectorOperations::functions<Number, Number2, dealii::MemorySpace::Host>::dot(const std::shared_ptr<dealii::parallel::internal::TBBPartitioner>&, dealii::internal::VectorOperations::size_type, const dealii::MemorySpace::MemorySpaceData<Number2, dealii::MemorySpace::Host>&, dealii::MemorySpace::MemorySpaceData<Number, dealii::MemorySpace::Host>&) [with Number = double; Number2 = double; dealii::internal::VectorOperations::size_type = long unsigned int]
The violated condition was: 
    dealii::numbers::is_finite(sum)
Additional information: 
    In a significant number of places, deal.II checks that some
    intermediate value is a finite number (as opposed to plus or minus
    infinity, or NaN/Not a Number). In the current function, we
    encountered a number that is not finite (its value is (-nan,0) and
    therefore violates the current assertion).

    This may be due to the fact that some operation in this function
    created such a value, or because one of the arguments you passed to
    the function already had this value from some previous operation. In
    the latter case, this function only triggered the error but may not
    actually be responsible for the computation of the number that is not
    finite.

    There are two common cases where this situation happens. First, your
    code (or something in deal.II) divides by zero in a place where this
    should not happen. Or, you are trying to solve a linear system with an
    unsuitable solver (such as an indefinite or non-symmetric linear
    system using a Conjugate Gradient solver); such attempts oftentimes
    yield an operation somewhere that tries to divide by zero or take the
    square root of a negative value.

    In any case, when trying to find the source of the error, recall that
    the location where you are getting this error is simply the first
    place in the program where there is a check that a number (e.g., an
    element of a solution vector) is in fact finite, but that the actual
    error that computed the number may have happened far earlier. To find
    this location, you may want to add checks for finiteness in places of
    your program visited before the place where this error is produced.
    One way to check for finiteness is to use the 'AssertIsFinite' macro.

Stacktrace:
-----------
#0  /work/06923/tg862222/stampede2/software/2022-11-01-dealii-9.4-candi-native-64b/deal.II-v9.4.0/lib/libdeal_II.g.so.9.4.0: dealii::internal::VectorOperations::functions<double, double, dealii::MemorySpace::Host>::dot(std::shared_ptr<dealii::parallel::internal::TBBPartitioner> const&, unsigned long, dealii::MemorySpace::MemorySpaceData<double, dealii::MemorySpace::Host> const&, dealii::MemorySpace::MemorySpaceData<double, dealii::MemorySpace::Host>&)
#1  /work/06923/tg862222/stampede2/software/2022-11-01-dealii-9.4-candi-native-64b/deal.II-v9.4.0/lib/libdeal_II.g.so.9.4.0: double dealii::LinearAlgebra::distributed::Vector<double, dealii::MemorySpace::Host>::inner_product_local<double>(dealii::LinearAlgebra::distributed::Vector<double, dealii::MemorySpace::Host> const&) const
#2  /work/06923/tg862222/stampede2/software/2022-11-01-dealii-9.4-candi-native-64b/deal.II-v9.4.0/lib/libdeal_II.g.so.9.4.0: dealii::LinearAlgebra::distributed::Vector<double, dealii::MemorySpace::Host>::operator*(dealii::LinearAlgebra::VectorSpaceVector<double> const&) const
#3  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-dg-main/aspect: dealii::internal::SolverCG::IterationWorker<dealii::LinearAlgebra::distributed::Vector<double, dealii::MemorySpace::Host>, aspect::MatrixFreeStokesOperators::MassMatrixOperator<3, 1, double>, dealii::PreconditionMG<3, dealii::LinearAlgebra::distributed::Vector<double, dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >, int>::do_iteration(unsigned int)
#4  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-dg-main/aspect: void dealii::SolverCG<dealii::LinearAlgebra::distributed::Vector<double, dealii::MemorySpace::Host> >::solve<aspect::MatrixFreeStokesOperators::MassMatrixOperator<3, 1, double>, dealii::PreconditionMG<3, dealii::LinearAlgebra::distributed::Vector<double, dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> > >(aspect::MatrixFreeStokesOperators::MassMatrixOperator<3, 1, double> const&, dealii::LinearAlgebra::distributed::Vector<double, dealii::MemorySpace::Host>&, dealii::LinearAlgebra::distributed::Vector<double, dealii::MemorySpace::Host> const&, dealii::PreconditionMG<3, dealii::LinearAlgebra::distributed::Vector<double, dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> > const&)
#5  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-dg-main/aspect: aspect::internal::BlockSchurGMGPreconditioner<aspect::MatrixFreeStokesOperators::StokesOperator<3, 2, double>, aspect::MatrixFreeStokesOperators::ABlockOperator<3, 2, double>, aspect::MatrixFreeStokesOperators::MassMatrixOperator<3, 1, double>, dealii::PreconditionMG<3, dealii::LinearAlgebra::distributed::Vector<double, dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >, dealii::PreconditionMG<3, dealii::LinearAlgebra::distributed::Vector<double, dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> > >::vmult(dealii::LinearAlgebra::distributed::BlockVector<double>&, dealii::LinearAlgebra::distributed::BlockVector<double> const&) const
#6  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-dg-main/aspect: void dealii::SolverFGMRES<dealii::LinearAlgebra::distributed::BlockVector<double> >::solve<aspect::MatrixFreeStokesOperators::StokesOperator<3, 2, double>, aspect::internal::BlockSchurGMGPreconditioner<aspect::MatrixFreeStokesOperators::StokesOperator<3, 2, double>, aspect::MatrixFreeStokesOperators::ABlockOperator<3, 2, double>, aspect::MatrixFreeStokesOperators::MassMatrixOperator<3, 1, double>, dealii::PreconditionMG<3, dealii::LinearAlgebra::distributed::Vector<double, dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >, dealii::PreconditionMG<3, dealii::LinearAlgebra::distributed::Vector<double, dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> > > >(aspect::MatrixFreeStokesOperators::StokesOperator<3, 2, double> const&, dealii::LinearAlgebra::distributed::BlockVector<double>&, dealii::LinearAlgebra::distributed::BlockVector<double> const&, aspect::internal::BlockSchurGMGPreconditioner<aspect::MatrixFreeStokesOperators::StokesOperator<3, 2, double>, aspect::MatrixFreeStokesOperators::ABlockOperator<3, 2, double>, aspect::MatrixFreeStokesOperators::MassMatrixOperator<3, 1, double>, dealii::PreconditionMG<3, dealii::LinearAlgebra::distributed::Vector<double, dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> >, dealii::PreconditionMG<3, dealii::LinearAlgebra::distributed::Vector<double, dealii::MemorySpace::Host>, dealii::MGTransferMatrixFree<3, double> > > const&)
#7  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-dg-main/aspect: aspect::StokesMatrixFreeHandlerImplementation<3, 2>::solve()
#8  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-dg-main/aspect: aspect::Simulator<3>::solve_stokes()
#9  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-dg-main/aspect: aspect::Simulator<3>::assemble_and_solve_stokes(bool, double*)
#10  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-dg-main/aspect: aspect::Simulator<3>::solve_single_advection_iterated_stokes()
#11  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-dg-main/aspect: aspect::Simulator<3>::solve_timestep()
#12  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-dg-main/aspect: aspect::Simulator<3>::run()
#13  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-dg-main/aspect: void run_simulator<3>(std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, bool, bool, bool)
#14  /work2/06923/tg862222/stampede2/software/aspect/2023-02-01-aspect-B/aspect/build-dg-main/aspect: main
--------------------------------------------------------
tjhei commented 1 year ago

Thank you. I simplified it further (DoFs, just Dirichlet conditions):

set World builder file = gmg-test-v6.wb
set Output directory = gmg-test-gmg-v7 #-gmg
set Dimension = 3
set CFL number              = 0.5
set Max nonlinear iterations = 3
set End time = 0
set Nonlinear solver scheme = single Advection, iterated Stokes 
#set Nonlinear solver scheme = single Advection, iterated defect correction Stokes 
set Timing output frequency = 0
set Max nonlinear iterations in pre-refinement = 10
set Pressure normalization                 = surface
set Nonlinear solver tolerance =1e-4
set Resume computation = auto
set Adiabatic surface temperature = 273.

subsection Solver parameters
  subsection Newton solver parameters
    set Max pre-Newton nonlinear iterations = 100
    set Nonlinear Newton solver switch tolerance = 1e-5
    set Max Newton line search iterations = 0
    set Maximum linear Stokes solver tolerance = 1e-2
    set Use Newton residual scaling method = true
    set Use Newton failsafe = false 
    set Use Eisenstat Walker method for Picard iterations = false #true
  end
  subsection Stokes solver parameters
    set Maximum number of expensive Stokes solver steps =  5000    
    #set Number of cheap Stokes solver steps = 0 #80000
    set Number of cheap Stokes solver steps = 1000 #80000
    #set Linear solver tolerance = 1e-2
    set Stokes solver type = block GMG
    #set GMRES solver restart length = 100
  end
    subsection Matrix Free
    set Output details         = true
  end

end

subsection Heating model
  set List of model names = adiabatic heating, shear heating 
end

subsection Boundary composition model
  set List of model names = initial composition
  set Fixed composition boundary indicators = 0,1,2,3,5 #,6,7,8,9
end

subsection Boundary temperature model
  set Fixed temperature boundary indicators   = 0,1,2,3,4,5 #,6,7,8,9
  set List of model names = initial temperature 
  subsection Initial temperature
    set Minimal temperature = 273.15
    set Maximal temperature = 4000
  end
end

subsection Discretization
  subsection Stabilization parameters
    set Use artificial viscosity smoothing = false 
    set Global composition maximum = 1.0 
    set Global composition minimum = 0.0 
    set Global temperature minimum = 273.15
    set alpha = 2
    set beta  = 0.78 #0.117 
   end
end

subsection Boundary velocity model
  #set Prescribed velocity boundary indicators = uppereast: function, upperwest: function
  set Tangential velocity boundary indicators = #top, north, east, south, west 
  set Zero velocity boundary indicators       = bottom, top, north, east, south, west 
  subsection Function
    set Coordinate system = spherical
    set Use spherical unit vectors = true
    set Variable names = r,phi,theta,t
    set Function expression = 0; -0.0000001;0 
  end
end

subsection Boundary traction model 
  #set Prescribed traction boundary indicators =  lowernorth: initial lithostatic pressure, lowereast: initial lithostatic pressure, lowersouth: initial lithostatic pressure, lowerwest: initial lithostatic pressure, uppersouth: initial lithostatic pressure, uppernorth: initial lithostatic pressure 
  subsection Initial lithostatic pressure
    set Representative point = 6371e3,-114.01,39.01
  end
end

subsection Geometry model
  set Model name = chunk #with lithosphere boundary indicators

  subsection Chunk #with lithosphere boundary indicators
    set Chunk inner radius = 5571e3
    set Chunk outer radius = 6371e3
    set Chunk minimum latitude = 39 
    set Chunk maximum latitude = 46 
    set Chunk minimum longitude = -131 
    set Chunk maximum longitude = -114
    set Longitude repetitions = 15 #20 #32
    set Latitude repetitions = 8 #32
    set Radius repetitions = 7
  end
  subsection Chunk with lithosphere boundary indicators
    set Chunk inner radius = 5571e3
    set Chunk outer radius = 6371e3
    set Chunk middle boundary radius = 6271e3
    set Chunk minimum latitude = 39 
    set Chunk maximum latitude = 46 
    set Chunk minimum longitude = -131 
    set Chunk maximum longitude = -114
    set Longitude repetitions = 15 #20 #32
    set Latitude repetitions = 8 #32
    set Inner chunk radius repetitions = 7
  end
end

subsection Material model
  set Material averaging = harmonic average # only viscosity  #project to Q1 only viscosity #harmonic average only viscosity  #harmonic average #project to Q1 only viscosity
  set Model name = visco plastic
  subsection Visco Plastic
        set Reference temperature = 273
        #set Reference viscosity = 1e20
        set Minimum viscosity = 1e19
        set Maximum viscosity = 1e24
        # the background rheologically becomes oceanic curst at depths lower than 80km, this is to hopefullly fix the boundaries. The densities remain the same for the open boundaries to remain working.
        set Phase transition depths = background:80e3|410e3|520e3|560e3|670e3|670e3|670e3|670e3, PC_crust:80e3|665e3|720e3, spharz: 410e3|520e3|560e3|670e3|670e3|670e3|670e3
        set Phase transition widths = background:5e3|5e3|5e3|5e3|10e3|5e3|5e3|5e3, PC_crust:5e3|5e3|5e3, spharz: 5e3|5e3|5e3|5e3|5e3|5e3|5e3
        set Phase transition temperatures = background:1662.0|1662.0|1662.0|1662.0|1662.0|1662.0|1662.0|1662.0, PC_crust:1173.0|1662.0|1662.0, spharz: 1662.0|1662.0|1662.0|1662.0|1662.0|1662.0|1662.0
        set Phase transition Clapeyron slopes = background:0.0|4e6|4.1e6|4e6|-2e6|4e6|-3.1e6|1.3e6, PC_crust:0.0|4e6|1.3e6, spharz: 4e6|4.1e6|4e6|-2e6|4e6|-3.1e6|1.3e6
        set Thermal diffusivities = 1.0e-6
        set Heat capacities = 1250.0
        set Densities = background: 3300.0|3300.0|3394.4|3442.1|3453.2|3617.6|3691.5|3774.7|3929.1, NA_crust:3000.0, PC_crust:3000.0|3540.0|3613.0|3871.7, spharz: 3235.0|3372.3|3441.7|3441.7|3680.8|3717.8|3759.4|3836.6 #,  total_strain:1e5, water:1e5
        #set Densities = background: 3300.0|3300|3300|3300|3300|3300|3300|3300, NA_crust:3100.0|3100.0|3100.0|3100, PC_crust:3100.0|3100|3100.0|3100, spharz: 3100|3100|3100|3100|3100|3100|3100|3100 , total_strain:1e5, water: 1e5
        set Thermal expansivities = 3.1e-5
        set Viscosity averaging scheme = harmonic
        set Viscous flow law = composite
        set Yield mechanism = drucker
        set Grain size = 1.0000e-02
        set Prefactors for diffusion creep =            background:8.5700e-28|7.1768e-15|7.1768e-15|7.1768e-15|7.1768e-15|2.4977e-19|2.4977e-19|2.4977e-19|2.4977e-19, NA_crust:8.5700e-28, PC_crust:8.5700e-28|7.1768e-15|2.4977e-19|2.4977e-19, spharz:7.1768e-15|7.1768e-15|7.1768e-15|7.1768e-15|2.4977e-19|2.4977e-19|2.4977e-19|2.4977e-19 
        set Grain size exponents for diffusion creep =  background:3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00, NA_crust:3.0000e+00, PC_crust:3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00, spharz:3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00|3.0000e+00 
        set Activation energies for diffusion creep =   background:3.7500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05, NA_crust:3.7500e+05, PC_crust:3.7500e+05|3.2500e+05|3.2500e+05|3.2500e+05, spharz:3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05|3.2500e+05 
        set Activation volumes for diffusion creep =    background:2.3000e-05|2.3000e-05|2.3000e-05|2.3000e-05|2.3000e-05|3.0000e-06|3.0000e-06|3.0000e-06|3.0000e-06, NA_crust:2.3000e-05, PC_crust:2.3000e-05|2.3000e-05|3.0000e-06|3.0000e-06, spharz:2.3000e-05|2.3000e-05|2.3000e-05|2.3000e-05|3.0000e-06|3.0000e-06|3.0000e-06|3.0000e-06 
        set Prefactors for dislocation creep =          background:8.5700e-28|4.4668e-16|4.4668e-16|4.4668e-16|4.4668e-16|5.0000e-32|5.0000e-32|5.0000e-32|5.0000e-32, NA_crust:8.5700e-28, PC_crust:8.5700e-28|4.4668e-16|5.0000e-32|5.0000e-32, spharz:4.4668e-16|4.4668e-16|4.4668e-16|4.4668e-16|5.0000e-32|5.0000e-32|5.0000e-32|5.0000e-32 
        set Stress exponents for dislocation creep =    background:4.0000e+00|3.5000e+00|3.5000e+00|3.5000e+00|3.5000e+00|1.0000e+00|1.0000e+00|1.0000e+00|1.0000e+00, NA_crust:3.5000e+00, PC_crust:4.0000e+00|3.5000e+00|1.0000e+00|1.0000e+00, spharz:3.5000e+00|3.5000e+00|3.5000e+00|3.5000e+00|1.0000e+00|1.0000e+00|1.0000e+00|1.0000e+00 
        set Activation energies for dislocation creep = background:2.2300e+05|4.8000e+05|4.8000e+05|4.8000e+05|4.8000e+05|0.0000e+00|0.0000e+00|0.0000e+00|0.0000e+00, NA_crust:2.2300e+05, PC_crust:2.2300e+05|4.8000e+05|0.0000e+00|0.0000e+00, spharz:4.8000e+05|4.8000e+05|4.8000e+05|4.8000e+05|0.0000e+00|0.0000e+00|0.0000e+00|0.0000e+00
        set Activation volumes for dislocation creep =  background:2.4000e-05|2.4000e-05|2.4000e-05|2.4000e-05|2.4000e-05|0.0000e+00|0.0000e+00|0.0000e+00|0.0000e+00, NA_crust:2.4000e-05, PC_crust:2.4000e-05|2.4000e-05|0.0000e+00|0.0000e+00, spharz:2.4000e-05|2.4000e-05|2.4000e-05|2.4000e-05|0.0000e+00|0.0000e+00|0.0000e+00|0.0000e+00 

        #set Manually define phase method crust = background: 0.0, PC_crust:1.3, spharz:0.0  # CDPT
        #set Manually define phase method pyrolite = background:1.0, PC_crust: 0.0, spharz:0.0  # CDPT
        #set Manually define phase method harzburgite = background:0.0, PC_crust: 0.0, spharz:1.0  # CDPT
        #set Constant viscosity prefactors      = background:10.,  NA_crust:1.,  PC_crust:1., spharz:1.,  total_strain:1, water:1 #,       100.0 , 1.0
        #set Constant viscosity prefactors      =2.5, 1., 1., 1., 1, 1 #,       100.0 , 1.0

    set Angles of internal friction = background: 45|45|0|0|0|0|0|0|0, NA_crust:45,  PC_crust:45|45|0|0, spharz:45|45|0|0|0|0|0|0 #, total_strain:0, water:0
    set Cohesions                   = background: 20e8|20.e8|20e8|20e8|20e8|20e8|20e8|20e8|20e8, NA_crust:20e6, PC_crust:10.e4|20e6|20e8|20e8, spharz:20e6|20e8|20e8|20e8|20e8|20e8|20e8|20e8 #, total_strain:1, water:1
  end
end
subsection Compositional fields
  set Number of fields = 3
  set Names of fields =  NA_crust, PC_crust, spharz 
  set Types of fields = chemical composition, chemical composition, chemical composition 
  set List of normalized fields = 0,1,2
end

subsection Initial temperature model
    set Model name = world builder
end

subsection Initial composition model
    set Model name = world builder
end

subsection Gravity model
  set Model name = radial constant
end

subsection Mesh refinement
  set Initial global refinement = 1
  set Initial adaptive refinement = 1 
  set Run postprocessors on initial refinement            = true
  set Refinement fraction = 0.1
  set Time steps between mesh refinement = 5
  set Strategy = thermal energy density #isosurfaces, , minimum refinement function, maximum refinement function
   subsection Isosurfaces
                          #minref maxref  mintemp maxtemp
        set Isosurfaces = \
                          max,    max,   Temperature:   0 | 350; \ 
                          max,    max,   Temperature:   0 | 1300, PC_crust: 0.5 | 2; \ 
                          max-1,  max-1, Temperature:   350 | 400; \ 
                          max-2,  max-2, Temperature:   400 | 1300; \ 
                          max-1,    max,   PC_crust: 0.5 | 2; \ 
                          max-1,  max-1, spharz: 0.5 | 2; \ 
                          max-3,  max-3,  Temperature:1300 | 1660; \
                          min,    min,   Temperature:1660 | 3000; \

   end
end

subsection Postprocess
  set List of postprocessors = visualization, boundary pressures, pressure statistics, velocity statistics, spherical velocity statistics, temperature statistics, composition statistics, load balance statistics, memory statistics 
  subsection Visualization
    set List of output variables = depth, density, viscosity, strain rate,stress second invariant 
    set Time between graphical output = 0 
    set Write higher order output = false 
    set Interpolate output = false
  end
end

This is what happens:

MFraters commented 1 year ago

hmm, yes, now realized that I did notice that the crash always happened after all the cheap iterations where done and it hadn't converted yet. I don't know the gmg solver well enough to have ideas what could cause this issue. Is the gmg solver using a different kind of expensive preconditioner than the AMG solver? Let me know if there is anything I can do to help looking further into the problem.

tjhei commented 1 year ago

set Nonlinear solver scheme = single Advection, iterated Stokes

I was reading through the defect correction algorithm but I only now realized that this is just iterated Stokes!

MFraters commented 1 year ago

ah, yes, I thought that would be a good simplification, since that is what most code was designed and tested for in the first place and most people are familiar with. I should have be more clear about that I used that :)

lhy11009 commented 1 year ago

@MFraters, @mibillen I had a temporary fix for this issue (block GMG fails after the cheap stokes in an iterated solver). The fix is to skip the expensive stokes after the cheap stokes. This leads to the solver's failure to reach the linear tolerance in one iteration, but doesn't prevent it from reaching the nonlinear tolerance. The implementation detail in commit 01428009480c6255f45a8c74d71e93e77075bc30

Here is an example's output file: log.txt Timestep 1, failure to converge in iterations 1 and 2, the nonlinear solver still converges to a residual below 1e-3.

not a solution for the issue, research purpose, continue to run with convergence and benefit from the faster speed of the block GMG

lhy11009 commented 1 year ago

I just summarize what we discussed 09272023 (@mibillen )

https://github.com/lhy11009/aspect/tree/skip_expensive_stokes