convection-box-particles.prm

# A modification of the convection box cookbook using physical units
# and particle output. See the manual for more information.

# At the top, we define the number of space dimensions we would like to
# work in:
set Dimension                              = 2

# There are several global variables that have to do with what
# time system we want to work in and what the end time is. We
# also designate an output directory.
set Use years in output instead of seconds = true
set End time                               = 3e10
set Output directory                       = output-tutorial

# Then there are variables that describe how the pressure should
# be normalized. Here, 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

# Then come a number of sections that deal with the setup
# of the problem to solve. The first one deals with the
# geometry of the domain within which we want to solve.
# The sections that follow all have the same basic setup
# where we select the name of a particular model (here,
# the box geometry) and then, in a further subsection,
# set the parameters that are specific to this particular
# model.
subsection Geometry model
  set Model name = box

  subsection Box
    set X extent = 4.2e6
    set Y extent = 3e6
  end
end

# The next section deals with the initial conditions for the
# temperature. Note that there are no initial conditions for the
# velocity variable since the velocity is assumed to always
# be in a static equilibrium with the temperature field.
# There are a number of models with the 'function' model
# a generic one that allows us to enter the actual initial
# conditions in the form of a formula that can contain
# constants. We choose a linear temperature profile that
# matches the boundary conditions defined below plus
# a small perturbation. The variables in this equation are
# described below, and it is important to note that in many
# cases the values correspond to other model parameters
# defined elsewhere. As such, if these model parameters are
# changed, the values below will also need to be adjusted.
#   L - Model length/width
#   D - Model depth/height
#   T_top - Temperature at the top boundary
#   T_bottom - Temperature at the bottom boundary
#   p, k - values related to the small temperature perturbation

subsection Initial temperature model
  set Model name = function

  subsection Function
    set Variable names      = x,y
    set Function constants  = p=-0.01, L=4.2e6, D=3e6, pi=3.1415926536, k=1, T_top=273, T_bottom=3600
    set Function expression = T_top + (T_bottom-T_top)*(1-(y/D) - p*cos(k*pi*x/L)*sin(pi*y/D))
  end
end

# Then follows a section that describes the boundary conditions
# for the temperature. The model we choose is called 'box' and
# allows to set a constant temperature on each of the four sides
# of the box geometry. In our case, we choose something that is
# heated from below and cooled from above, 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).
subsection Boundary temperature model
  set Fixed temperature boundary indicators = bottom, top
  set List of model names = box

  subsection Box
    set Bottom temperature = 3300
    set Top temperature    = 300
  end
end

# The next parameters then describe on which parts of the
# boundary we prescribe a zero or nonzero velocity and
# on which parts the flow is allowed to be tangential.
# Here, all four sides of the box allow tangential
# unrestricted flow but with a zero normal component:
subsection Boundary velocity model
  set Tangential velocity boundary indicators = left, right, bottom, top
end

# The following two sections describe first the
# direction (vertical) and magnitude of gravity and the
# material model (i.e., density, viscosity, etc).
subsection Gravity model
  set Model name = vertical

  subsection Vertical
    set Magnitude = 10.0
  end
end

subsection Material model
  set Model name = simple

  subsection Simple model
    set Viscosity                     = 3e24
    set Reference density             = 3000
    set Thermal conductivity          = 4.86
    set Thermal expansion coefficient = 2e-5
    set Reference specific heat       = 1000
    set Reference temperature         = 0
  end
end

# We also have to specify that we want to use the Boussinesq
# approximation (assuming the density in the temperature
# equation to be constant, and incompressibility).
subsection Formulation
  set Formulation = Boussinesq approximation
end

# 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 global refinement                = 3
  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 1e7 years.
subsection Postprocess
  set List of postprocessors = velocity statistics, temperature statistics, heat flux statistics, visualization, particles, basic statistics

  subsection Visualization
    set Time between graphical output = 1e7
    set Output format                 = vtu
    set List of output variables      = material properties
  end

  subsection Particles
    set Number of particles = 1000
    set Time between data output = 1e7
    set Data output format = vtu
  end
end

subsection Solver parameters
  set Temperature solver tolerance = 1e-10
end