Tosi_base.prm

# A description of viscoplastic convection in a 2d box as
# described in Tosi et al. (2015).
# Tosi, N., et al. (2015),
# A community benchmark for viscoplastic thermal
# convection in a 2-D square box,
# Geochem. Geophys. Geosyst., 16, 2175–2196,
# doi:10.1002/2015GC005807.

# We add the necessary plugins through a shared library
set Additional shared libraries            = ./libtosi_benchmark.so

# This is a 2D benchmark.
set Dimension                              = 2

# We run the model from time = 0 to time = 1 and to handle
# the nonlinearity of the viscosity correctly, we use
# the iterated Advection and Stokes nonlinear solver with 200 cheap
# Stokes solves.
set Use years instead of seconds           = false
set End time                               = 1
set Nonlinear solver scheme                = iterated Advection and Stokes
set Max nonlinear iterations               = 30
set Nonlinear solver tolerance             = 1e-5

# We choose a zero average pressure
# at the surface of the domain (for the current geometry, the
# surface is defined as the top boundary).
set Pressure normalization                 = surface
set Surface pressure                       = 0

subsection Solver parameters
  subsection Stokes solver parameters
    set Number of cheap Stokes solver steps    = 200
    set Linear solver tolerance                = 1e-6
  end
end

# Our model domain is the unit square.
subsection Geometry model
  set Model name = box

  subsection Box
    set X extent = 1
    set Y extent = 1
  end
end

# We prescribe a linear temperature profile that
# matches the boundary conditions defined below plus
# a small perturbation with amplitude A.
subsection Initial temperature model
  set Model name = function

  subsection Function
    set Variable names      = x,z
    set Function constants  = A=0.01, pi=3.1415926536
    set Function expression = (1.0-z) + A*cos(pi*x)*sin(pi*z)
  end
end

# The top and bottom boundary carry a prescribed boundary temperature, whereas
# all other parts of the boundary are insulated (i.e., no heat flux through
# these boundaries; this is also often used to specify symmetry boundaries).
# The temperature at the bottom is set to 1; at the top to 0.
subsection Boundary temperature model
  set Fixed temperature boundary indicators   = bottom, top
  set List of model names = box

  subsection Box
    set Bottom temperature = 1
    set Left temperature   = 0
    set Right temperature  = 0
    set Top temperature    = 0
  end
end

subsection Boundary velocity model
  # Here, all four sides of the box allow tangential
  # unrestricted flow but with a zero normal component:
  set Tangential velocity boundary indicators = left, right, bottom, top
end

subsection Formulation
  set Formulation = Boussinesq approximation
end

# Gravity is vertical and of magnitude 1e8.
# The Rayleigh number is defined as
# $Ra = \frac{\rho g \alpha \Delta_T D^3}{\kappa \eta_0}$,
# with $\rho$ the density, g the gravity magnitude, $\alpha$ thermal
# expansivity, $\Delta_T$ the temperature contrast across the mantle and D the domain height.
# $\eta_0$ is the reference viscosity and $\kappa$ the thermal diffusivity.
# With the parameter values in this and the next section, this leads to
# $Ra = \frac{1 * 1e8 * 1e-6 * 1 * 1^3}{\frac{1}{1*1} * 1e-1} = 100$.
# Note that $\eta_0$ is not the same as the Tosi material model "Reference viscosity".
subsection Gravity model
  set Model name = vertical

  subsection Vertical
    set Magnitude = 1e8
  end
end

# Each case has a surface Rayleigh number of 100
# and the parameter controlling the temperature
# dependence of viscosity "Thermal viscosity parameter"
# is always 1e5.
# Between the cases the "Pressure viscosity parameter" ($\Delta_{\eta}_z$),
# the "Yield stress" ($\sigma_y$) and the "Nonlinear viscosity constant" ($\eta^*$)
# are varied, see Eq. (6-8) and Table 1 of the paper.
subsection Material model
  set Model name = TosiMaterial

  subsection Tosi benchmark
    set Reference density             = 1
    set Reference specific heat       = 1
    set Reference temperature         = 0
    set Thermal conductivity          = 1
    set Thermal expansion coefficient = 1e-06
    set Thermal viscosity parameter   = 1e5
    set Initial viscosity        = 1e-1
    set Minimum viscosity        = 1e-6
    set Maximum viscosity        = 10
  end
end

# The settings above all pertain to the description of the
# continuous partial differential equations we want to solve.
# The following section deals with the discretization of
# this problem, namely the kind of mesh we want to compute
# on. We here use a globally refined mesh without
# adaptive mesh refinement.
subsection Mesh refinement
  set Initial adaptive refinement              = 0
  set Time steps between mesh refinement       = 0
end

# The final part is to specify what ASPECT should do with the
# solution once computed at the end of every time step. The
# process of evaluating the solution is called `postprocessing'
# and we choose to compute velocity and temperature statistics,
# statistics about the heat flux through the boundaries of the
# domain, and to generate graphical output files for later
# visualization. These output files are created every time
# a time step crosses time points separated by 0.05. Given
# our start time (0.0) and final time (1.0) this means that
# we will obtain 20 output files.
#
# The statistics file includes the data reported in Tosi et al. (2015):
# - average temperature (temperature statistics)
# - top and bottom Nusselt numbers (heat flux statistics)
# - Vrms over the whole domain (velocity statistics)
# - Vrms along the surface (velocity boundary statistics)
# - Vmax along the surface (velocity boundary statistics)
# - average rate of viscous dissipation (TosiPostprocessor)
# - average rate of work against gravity (TosiPostprocessor)
# - percentage error difference work and dissipation (TosiPostprocessor)
#
# Tosi et al. (2015) also report depth averages; averages for viscosity
# and temperature are produced every 0.05 units of time.
subsection Postprocess
  set List of postprocessors = velocity statistics, temperature statistics, heat flux statistics, \
                               visualization, TosiPostprocessor, depth average, velocity boundary statistics

  subsection Depth average
    set List of output variables = temperature, viscosity
    set Time between graphical output = 0.05
    set Number of zones = 10
  end

  subsection Visualization
    set List of output variables = material properties, strain rate, gravity, heating
    set Time between graphical output = 0.05

    subsection Material properties
      set List of material properties = viscosity, density, specific heat, thermal expansivity
    end
  end
end