Skip to content

Commit

Permalink
Merge pull request #2288 from tjhei/release-2.0.1-picks1
Browse files Browse the repository at this point in the history
[2.0.1] pick #2258 #2253
  • Loading branch information
gassmoeller committed Jun 14, 2018
2 parents 3373d2f + c652fed commit fc33a30
Show file tree
Hide file tree
Showing 11 changed files with 813 additions and 27 deletions.
13 changes: 13 additions & 0 deletions include/aspect/material_model/depth_dependent.h
Original file line number Diff line number Diff line change
Expand Up @@ -45,6 +45,19 @@ namespace aspect
class DepthDependent : public MaterialModel::Interface<dim>, public ::aspect::SimulatorAccess<dim>
{
public:
/**
* Initialize the base model at the beginning of the run.
*/
virtual
void initialize();

/**
* Update the base model and viscosity function at the beginning of
* each timestep.
*/
virtual
void update();

/**
* Function to compute the material properties in @p out given the
* inputs in @p in.
Expand Down
6 changes: 3 additions & 3 deletions include/aspect/simulator/assemblers/interface.h
Original file line number Diff line number Diff line change
Expand Up @@ -216,9 +216,9 @@ namespace aspect

FEValues<dim> finite_element_values;

std_cxx11::shared_ptr<FEFaceValues<dim> > face_finite_element_values;
std_cxx11::shared_ptr<FEFaceValues<dim> > neighbor_face_finite_element_values;
std_cxx11::shared_ptr<FESubfaceValues<dim> > subface_finite_element_values;
std_cxx11::unique_ptr<FEFaceValues<dim> > face_finite_element_values;
std_cxx11::unique_ptr<FEFaceValues<dim> > neighbor_face_finite_element_values;
std_cxx11::unique_ptr<FESubfaceValues<dim> > subface_finite_element_values;

std::vector<types::global_dof_index> local_dof_indices;

Expand Down
25 changes: 25 additions & 0 deletions source/material_model/depth_dependent.cc
Original file line number Diff line number Diff line change
Expand Up @@ -31,6 +31,31 @@ namespace aspect
{
namespace MaterialModel
{
template <int dim>
void
DepthDependent<dim>::initialize()
{
base_model->initialize();
}



template <int dim>
void
DepthDependent<dim>::update()
{
base_model->update();

// we get time passed as seconds (always) but may want
// to reinterpret it in years
if (this->convert_output_to_years())
viscosity_function.set_time (this->get_time() / year_in_seconds);
else
viscosity_function.set_time (this->get_time());
}



template <int dim>
void
DepthDependent<dim>::read_viscosity_file(const std::string &filename,
Expand Down
69 changes: 45 additions & 24 deletions source/simulator/assemblers/interface.cc
Original file line number Diff line number Diff line change
Expand Up @@ -176,27 +176,27 @@ namespace aspect
finite_element_values (mapping,
finite_element, quadrature,
update_flags),
face_finite_element_values ((face_quadrature.size() > 0
?
new FEFaceValues<dim> (mapping,
finite_element, face_quadrature,
face_update_flags)
:
NULL)),
neighbor_face_finite_element_values ((face_quadrature.size() > 0
?
new FEFaceValues<dim> (mapping,
finite_element, face_quadrature,
face_update_flags)
:
NULL)),
subface_finite_element_values ((face_quadrature.size() > 0
?
new FESubfaceValues<dim> (mapping,
finite_element, face_quadrature,
face_update_flags)
:
NULL)),
face_finite_element_values (face_quadrature.size() > 0
?
new FEFaceValues<dim> (mapping,
finite_element, face_quadrature,
face_update_flags)
:
NULL),
neighbor_face_finite_element_values (face_quadrature.size() > 0
?
new FEFaceValues<dim> (mapping,
finite_element, face_quadrature,
face_update_flags)
:
NULL),
subface_finite_element_values (face_quadrature.size() > 0
?
new FESubfaceValues<dim> (mapping,
finite_element, face_quadrature,
face_update_flags)
:
NULL),
local_dof_indices (finite_element.dofs_per_cell),

phi_field (advection_element.dofs_per_cell, numbers::signaling_nan<double>()),
Expand Down Expand Up @@ -263,9 +263,30 @@ namespace aspect
scratch.finite_element_values.get_fe(),
scratch.finite_element_values.get_quadrature(),
scratch.finite_element_values.get_update_flags()),
face_finite_element_values (scratch.face_finite_element_values),
neighbor_face_finite_element_values (scratch.neighbor_face_finite_element_values),
subface_finite_element_values (scratch.subface_finite_element_values),
face_finite_element_values (scratch.face_finite_element_values.get()
?
new FEFaceValues<dim> (scratch.face_finite_element_values->get_mapping(),
scratch.face_finite_element_values->get_fe(),
scratch.face_finite_element_values->get_quadrature(),
scratch.face_finite_element_values->get_update_flags())
:
NULL),
neighbor_face_finite_element_values (scratch.neighbor_face_finite_element_values.get()
?
new FEFaceValues<dim> (scratch.neighbor_face_finite_element_values->get_mapping(),
scratch.neighbor_face_finite_element_values->get_fe(),
scratch.neighbor_face_finite_element_values->get_quadrature(),
scratch.neighbor_face_finite_element_values->get_update_flags())
:
NULL),
subface_finite_element_values (scratch.subface_finite_element_values.get()
?
new FESubfaceValues<dim> (scratch.subface_finite_element_values->get_mapping(),
scratch.subface_finite_element_values->get_fe(),
scratch.subface_finite_element_values->get_quadrature(),
scratch.subface_finite_element_values->get_update_flags())
:
NULL),
local_dof_indices (scratch.finite_element_values.get_fe().dofs_per_cell),

phi_field (scratch.phi_field),
Expand Down
174 changes: 174 additions & 0 deletions tests/depth_dependent_box_function_time_dependent.prm
Original file line number Diff line number Diff line change
@@ -0,0 +1,174 @@
# Test depth-dependent material model using simple material model as
# base model with depth-dependence specified by user-provided function.
# This is an extended copy from depth_dependent_box_function that
# uses a time dependent function to prescribe depth dependent viscosity.

# 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 = 1.0e8
set Output directory = output-depth-dependent-box-function-time-dependent

# Then there are variables that describe the tolerance of
# the linear solver as well as 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 = 3.0e6
set Y extent = 3.0e6
end
end


# The next section deals with the initial conditions for the
# temperature (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:
subsection Initial temperature model
set Model name = function

subsection Function
set Variable names = x,z
set Function constants = p=10.0, L=3.0e6, pi=3.1415926536, k=1
set Function expression = 2773.0 - z/L*(2500.0) + p*cos(k*pi*x/L)*sin(k*pi*z/L/2.0)
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. (As will be seen
# in the next section, the actual temperature prescribed here
# at the left and right does not matter.)
subsection Boundary temperature model
set List of model names = box

subsection Box
set Bottom temperature = 2773
set Left temperature = 0
set Right temperature = 0
set Top temperature = 273.0
end
end


# We then also have to prescribe several other parts of the model such as
# which boundaries actually 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 parameters below this comment were created by the update script
# as replacement for the old 'Model settings' subsection. They can be
# safely merged with any existing subsections with the same name.

subsection Boundary temperature model
set Fixed temperature boundary indicators = bottom, top
end

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). We have
# discussed the settings used here in the introduction to
# this cookbook in the manual already.
subsection Gravity model
set Model name = vertical

subsection Vertical
set Magnitude = 2e-5 # = Ra / Thermal expansion coefficient
end
end


subsection Material model
set Model name = depth dependent
subsection Depth dependent model
set Base model = simple
set Depth dependence method = Function
subsection Viscosity depth function
set Variable names = d,t
set Function expression = if(d>6.70e5,3.0e22*(1.0+t/1e8),1.0e21*(1.0+t/1e8));
end
end
subsection Simple model
set Reference density = 3300.0
set Reference specific heat = 1250.0
set Reference temperature = 0
set Thermal conductivity = 4.0
set Thermal expansion coefficient = 2e-5
set Viscosity = 1e25
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 global refinement = 4
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.01. Given
# our start time (zero) and final time (0.5) this means that
# we will obtain 50 output files.
subsection Postprocess
set List of postprocessors = velocity statistics, temperature statistics, heat flux statistics, depth average

subsection Depth average
set Time between graphical output = 0
end
end

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

subsection Stokes solver parameters
set Linear solver tolerance = 1e-7
end
end

0 comments on commit fc33a30

Please sign in to comment.