Skip to content

Commit

Permalink
Add benchmark prm
Browse files Browse the repository at this point in the history
  • Loading branch information
@anne-glerum @gassmoeller
anne-glerum authored and gassmoeller committed Nov 11, 2024
1 parent cbcf865 commit 5239c3a
Show file tree
Hide file tree
Showing 8 changed files with 5,750 additions and 0 deletions.
162 changes: 162 additions & 0 deletions benchmarks/diffusion_of_hill/1_sine_zero_flux_ci.prm
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,162 @@
set Additional shared libraries = libanalytical_topography.so
set Dimension = 2
set Use years in output instead of seconds = false
set End time = 0.005
set Maximum time step = 0.0005
set Output directory = output-1_sine_zero_flux_ci_toppresvel_subtractmeshvel/
set Nonlinear solver scheme = single Advection, single Stokes
set Pressure normalization = surface
set Surface pressure = 0

# 1x1 box with an initial hill topography
subsection Geometry model
set Model name = box

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

subsection Initial topography model
set Model name = function

subsection Function
set Function constants = A=0.150, L=1.
set Function expression = \
if(x<0.5,A * sin((x+0.5)*pi),0)
end
end
end

# Temperature effects are ignored
subsection Initial temperature model
set Model name = function

subsection Function
set Function expression = 0
end
end

subsection Boundary temperature model
set Fixed temperature boundary indicators = bottom, top
set List of model names = initial temperature
end

subsection Compositional fields
set Number of fields = 3
set Names of fields = layer_1, layer_2, layer_3
end

subsection Initial composition model
set Model name = function

subsection Function
set Variable names = x,z
set Coordinate system = cartesian #depth
#set Function expression = if(x<.05, 1, 0); if(x>=0.05 && x<0.1, 1, 0); 0
set Function expression = if(z>1.07, 1, 0); 0; 0
end
end

subsection Boundary composition model
set Fixed composition boundary indicators = top
set Allow fixed composition on outflow boundaries = false
set List of model names = function
subsection Function
set Variable names = x,z,t
set Function expression = if(t==0,1,0); 0; if(t==0,0,1)
end
end

# Free slip on all boundaries
subsection Boundary velocity model
set Tangential velocity boundary indicators = left, right, bottom
set Prescribed velocity boundary indicators = top: function
subsection Function
set Variable names = x,z,t
set Function expression = 0;0
end
end

# The mesh will be deformed according to the displacement
# of the surface nodes due to diffusion of the topography.
# The mesh is allowed to move vertical along the left and
# right boundary.
subsection Mesh deformation
set Mesh deformation boundary indicators = top: diffusion
set Additional tangential mesh velocity boundary indicators = left, right

subsection Diffusion
# The diffusivity
set Hillslope transport coefficient = 0.25
end
end

# Vertical gravity
subsection Gravity model
set Model name = vertical

subsection Vertical
set Magnitude = 0
end
end

# One material with unity properties
subsection Material model
set Model name = simple

subsection Simple model
set Reference density = 1
set Reference specific heat = 1
set Reference temperature = 0
set Thermal conductivity = 1
set Thermal expansion coefficient = 1
set Viscosity = 1
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

# We here use a globally refined mesh without
# adaptive mesh refinement.
subsection Mesh refinement
set Initial global refinement = 6
set Initial adaptive refinement = 0
set Time steps between mesh refinement = 0
set Strategy = minimum refinement function
subsection Minimum refinement function
set Variable names = x,z
set Coordinate system = cartesian
set Function expression = if(z>0.87,7,3)
end
end

# We output the computed topography and the analytical topography
# value to file.
subsection Postprocess
set List of postprocessors = velocity statistics, temperature statistics, heat flux statistics, visualization #, analytical topography

subsection Topography
set Output to file = true
set Analytical solution of example = 1
set Diffusivity = 0.25
set Initial sinusoidal topography amplitude = 0.075
end

subsection Visualization
set Time between graphical output = 0.0001
set Output mesh velocity = true
set Interpolate output = false
end
end

subsection Solver parameters
subsection Stokes solver parameters
set Linear solver tolerance = 1e-7
end
end
8 changes: 8 additions & 0 deletions source/simulator/helper_functions.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2306,6 +2306,7 @@ namespace aspect
update_quadrature_points | update_JxW_values);

std::vector<Tensor<1,dim>> face_current_velocity_values (fe_face_values.n_quadrature_points);
std::vector<Tensor<1,dim>> face_current_mesh_velocity_values (fe_face_values.n_quadrature_points);

const auto &tangential_velocity_boundaries =
boundary_velocity_manager.get_tangential_boundary_velocity_indicators();
Expand Down Expand Up @@ -2335,6 +2336,10 @@ namespace aspect
fe_face_values.reinit (cell, face_number);
fe_face_values[introspection.extractors.velocities].get_function_values(current_linearization_point,
face_current_velocity_values);
// get the mesh velocity, as we need to subtract it off of the advection systems
if (parameters.mesh_deformation_enabled)
fe_face_values[introspection.extractors.velocities].get_function_values(mesh_deformation->mesh_velocity,
face_current_mesh_velocity_values);

// ... check if the face is an outflow boundary by integrating the normal velocities
// (flux through the boundary) as: int u*n ds = Sum_q u(x_q)*n(x_q) JxW(x_q)...
Expand All @@ -2349,6 +2354,9 @@ namespace aspect
boundary_velocity = boundary_velocity_manager.boundary_velocity(face->boundary_id(),
fe_face_values.quadrature_point(q));

if (parameters.mesh_deformation_enabled)
boundary_velocity -= face_current_mesh_velocity_values[q];

integrated_flow += (boundary_velocity * fe_face_values.normal_vector(q)) *
fe_face_values.JxW(q);
}
Expand Down
98 changes: 98 additions & 0 deletions tests/hill_diffusion_field_boundary_conditions.prm
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,98 @@
# This setup tests the compositional field boundary conditions
# along a deforming boundary. To determine whether a part of the
# boundary experiences in- or outflow, the mesh velocity is
# subtracted from the material velocity. In this test, the
# prescribed material velocity is zero, and only the erosion
# or deposition (decrease or increase in height) by the diffusion
# applied to the surface topography determines whether boundary
# conditions should be applied (inflow/deposition),
# or not (outflow/erosion).

include $ASPECT_SOURCE_DIR/tests/hill_diffusion.prm

set Dimension = 2
set Nonlinear solver scheme = single Advection, single Stokes

# 1x1 box with an initial hill topography
subsection Geometry model
set Model name = box

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

subsection Initial topography model
set Model name = function

subsection Function
set Function constants = A=0.150, L=1.
set Function expression = \
if(x<0.5,A * sin((x+0.5)*pi),0)
end
end
end

# We here use a globally refined mesh without
# adaptive mesh refinement.
subsection Mesh refinement
set Initial global refinement = 4
end

subsection Compositional fields
set Number of fields = 3
set Names of fields = layer_1, layer_2, layer_3
end

subsection Initial composition model
set Model name = function

subsection Function
set Variable names = x,z
set Coordinate system = cartesian #depth
#set Function expression = if(x<.05, 1, 0); if(x>=0.05 && x<0.1, 1, 0); 0
set Function expression = if(z>1.07, 1, 0); 0; 0
end
end

subsection Boundary composition model
set Fixed composition boundary indicators = top
set Allow fixed composition on outflow boundaries = false
set List of model names = function
subsection Function
set Variable names = x,z,t
set Function expression = if(t==0,1,0); 0; if(t==0,0,1)
end
end

# Free slip on all boundaries except the top boundary.
subsection Boundary velocity model
set Tangential velocity boundary indicators = left, right, bottom
set Prescribed velocity boundary indicators = top: function
subsection Function
set Variable names = x,z,t
set Function expression = 0;0
end
end

# Vertical gravity
subsection Gravity model
set Model name = vertical

subsection Vertical
set Magnitude = 0
end
end

# We output the computed topography
# value to file.
subsection Postprocess
set List of postprocessors = velocity statistics, topography, visualization

subsection Visualization
set Time between graphical output = 0.0001
set Output mesh velocity = true
set Interpolate output = false
set Output format = gnuplot
end
end
63 changes: 63 additions & 0 deletions tests/hill_diffusion_field_boundary_conditions/screen-output
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,63 @@
-----------------------------------------------------------------------------
-----------------------------------------------------------------------------
-----------------------------------------------------------------------------

-----------------------------------------------------------------------------
-----------------------------------------------------------------------------
Number of active cells: 256 (on 5 levels)
Number of degrees of freedom: 6,823 (2,178+289+1,089+1,089+1,089+1,089)

Number of mesh deformation degrees of freedom: 578
Solving mesh displacement system... 0 iterations.
Added initial topography to grid

*** Timestep 0: t=0 seconds, dt=0 seconds
Solving mesh displacement system... 0 iterations.
Skipping temperature solve because RHS is zero.
Solving layer_1 system ... 0 iterations.
Skipping layer_2 composition solve because RHS is zero.
Skipping layer_3 composition solve because RHS is zero.
Solving Stokes system (GMG)... 0+0 iterations.

Postprocessing:
RMS, max velocity: 0 m/s, 0 m/s
Topography min/max: 0 m, 0.15 m
Writing graphical output: output-hill_diffusion_field_boundary_conditions/solution/solution-00000

*** Timestep 1: t=0.0005 seconds, dt=0.0005 seconds
Solving mesh surface diffusion
Solving mesh displacement system... 4 iterations.
Skipping temperature solve because RHS is zero.
Solving layer_1 system ... 11 iterations.
Skipping layer_2 composition solve because RHS is zero.
Solving layer_3 system ... 12 iterations.
Solving Stokes system (GMG)... 0+0 iterations.

Postprocessing:
RMS, max velocity: 0 m/s, 0 m/s
Topography min/max: -0.0002308 m, 0.1498 m
Writing graphical output: output-hill_diffusion_field_boundary_conditions/solution/solution-00001

*** Timestep 2: t=0.001 seconds, dt=0.0005 seconds
Solving mesh surface diffusion
Solving mesh displacement system... 4 iterations.
Skipping temperature solve because RHS is zero.
Solving layer_1 system ... 12 iterations.
Skipping layer_2 composition solve because RHS is zero.
Solving layer_3 system ... 12 iterations.
Solving Stokes system (GMG)... 0+0 iterations.

Postprocessing:
RMS, max velocity: 0 m/s, 0 m/s
Topography min/max: -0.0003541 m, 0.1496 m
Writing graphical output: output-hill_diffusion_field_boundary_conditions/solution/solution-00002

Termination requested by criterion: end time


+----------------------------------------------+------------+------------+
+----------------------------------+-----------+------------+------------+
+----------------------------------+-----------+------------+------------+

-----------------------------------------------------------------------------
-----------------------------------------------------------------------------
Loading

0 comments on commit 5239c3a

Please sign in to comment.