Davies et al benchmark test

[heronphi]

The Davies et al benchmark test looks to be a good test to see whether the solvers are working correctly. There is a table in the manual (Table 7) with values that can be compared to your model run, to make sure the executable is behaving correctly (and indeed aspect have been compiled properly).

However, the case-.prm examples in the benchmark folder are not set up to reproduce Table 7 (they require global refinement 7 instead of 5 as given). When I run a global refinement level higher than 5 with these caseprm files I get an error:

An error occurred in line <516> of file </ddn/data/klmz32/deal.ii-candi/tmp/unpack/deal.II-v8.5.0/source/lac/trilinos_solver.cc> in function void dealii::TrilinosWrappers::SolverBase::do_solve(const Preconditioner&) [with Preconditioner = dealii::TrilinosWrappers::PreconditionBase] The violated condition was: false Additional information: AztecOO::Iterate error code -3: loss of precision

Is this due to the compilation setup of aspect, or that the prm files cannot be used by simply increasing the global refinement level as specified in the manual?

Does there need to be more information on editing case-*.prm in the benchmark davies_et_al folder to obtain Table 7?

[gassmoeller]

Hi Phil, I never ran those benchmarks (maybe someone else can comment on what to change in the input files), but I can not reproduce your error message, using different setups (recent master of aspect, deal.II developer version and deal.II 8.5, case-1-1.prm, global refinement 7, 16 or 40 processes, debug or release mode). At least the model runs through the first 50 timesteps or so. Could you specify the setup you used, and maybe attach a parameter file and log file? Does the error occur in the first timestep?

[heronphi]

Hi Rene,

Did you try case-2-1.prm? case-1-1.prm actually ran (with global refinement 7 on 32 proc), but did not obtain the same Nu values as in Table 7. Could you (or anyone) run case-2-1.prm with global refinement 7 to see if it fails on the first tilmestep? I suspect my aspect build is not correct:

module load module-git gcc/4.9.1 cmake/3.6.2 blas/gcc/64/1 lapack/gcc/3.1.1 zlib/gcc/1.2.7 sge/current openmpi/gcc/2.1.1

-- This is ASPECT, the Advanced Solver for Problems in Earth's ConvecTion. -- . version 2.0.0-pre -- . running in DEBUG mode -- . running with 32 MPI processes -- . using Trilinos

Number of active cells: 196,608 (on 8 levels) Number of degrees of freedom: 2,566,656 (1,579,008+198,144+789,504)

*** Timestep 0: t=0 seconds Solving temperature system... 0 iterations. Rebuilding Stokes preconditioner... Solving Stokes system... 46+0 iterations.

Postprocessing: Writing graphical output: output/solution/solution-00000 RMS, max velocity: 10.1 m/s, 17.2 m/s Temperature min/avg/max: 0 K, 0.4516 K, 1 K Heat fluxes through boundary parts: -7.665 W, 13.95 W Writing depth average: output/depth_average.gnuplot Max and min velocity along boundary parts: 14.66 m/s, 0.02345 m/s, 17.24 m/s, 0.02523 m/s Radial RMS, tangential RMS, total RMS velocity: 5.97 m/s, 8.09 m/s, 10.1 m/s

*** Timestep 1: t=0.00169151 seconds Solving temperature system... 206 iterations. Solving Stokes system...

[gassmoeller]

Hi Phil, I was able to reproduce the issue in case-2.1.prm. The benchmark still uses our old method of scaling gravity vs thermal expansivity and it seems its scaling is too extreme (very large static pressures compared to dynamic pressures). I am not sure why it worked some time in the past, but reducing the gravity by 3 orders of magnitude and increasing the thermal expansivity by the same amount works for me. It might give slightly less accurate results (the density differences in the model are larger, and therefore less boussinesq-like), but it allows you to experiment a bit. It should be possible to change the model setup to using the actual boussinesq approximation and the nondimensional material model (implemented after this benchmark was added), if you want to give that a try that would be helpful.