Skip to content

Commit

Permalink
add test case
Browse files Browse the repository at this point in the history
  • Loading branch information
@danieldouglas92
danieldouglas92 committed Nov 4, 2024
1 parent 30de055 commit c2f6a2c
Show file tree
Hide file tree
Showing 8 changed files with 3,270 additions and 393 deletions.
390 changes: 0 additions & 390 deletions source/material_model/reaction_model/tian2019_parameterized_reactions.cc
View file

This file was deleted.

7 changes: 4 additions & 3 deletions source/simulator/helper_functions.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2331,7 +2331,7 @@ namespace aspect
MaterialModel::MaterialModelInputs<dim> in(fe_face_values.n_quadrature_points, introspection.n_compositional_fields);
MaterialModel::MaterialModelOutputs<dim> out(fe_face_values.n_quadrature_points, introspection.n_compositional_fields);
MeltHandler<dim>::create_material_model_outputs(out);
MaterialModel::MeltOutputs<dim> *fluid_out = out.template get_additional_output<MaterialModel::MeltOutputs<dim>>();
MaterialModel::MeltOutputs<dim> *fluid_out = nullptr;

const auto &tangential_velocity_boundaries =
boundary_velocity_manager.get_tangential_boundary_velocity_indicators();
Expand Down Expand Up @@ -2366,12 +2366,13 @@ namespace aspect
{
in.reinit(fe_face_values, cell, introspection, solution);
material_model->evaluate(in, out);
fluid_out = out.template get_additional_output<MaterialModel::MeltOutputs<dim>>();
}

if (consider_fluid_velocity)
{
const FEValuesExtractors::Vector ex_u_f = introspection.variable("fluid velocity").extractor_vector();
fe_face_values[ex_u_f].get_function_values(current_linearization_point, face_current_fluid_velocity_values);
fe_face_values[introspection.variable("fluid velocity").extractor_vector()].get_function_values(current_linearization_point,
face_current_fluid_velocity_values);
}

// ... check if the face is an outflow boundary by integrating the normal velocities
Expand Down
188 changes: 188 additions & 0 deletions tests/darcy_outflow.prm
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,188 @@
# This test defines an initial blob of porosity implaced in a solid within a 2D cartesian box.
# Because the porosity blob is less dense, the blob will rise to the surface of the model and
# flow out of the top of the box. This tests the implementation of composition dependent outflow
# boundary conditions for when a composition is advected with the Darcy field advection method.
set Dimension = 2
set Use years in output instead of seconds = true
set End time = 10e3
set Pressure normalization = surface
set Maximum time step = 100
set CFL number = 0.5
set Surface pressure = 0
set Use operator splitting = true
set Output directory = output-darcy
set Nonlinear solver scheme = iterated Advection and Stokes
set Max nonlinear iterations = 3

# Define a function which
subsection Boundary velocity model
set Prescribed velocity boundary indicators = top:function, bottom:function, right:function, left:function

subsection Function
set Variable names = x,y,
set Function expression = 0;-0.1
end
end

subsection Gravity model
set Model name = vertical

subsection Vertical
set Magnitude = 10
end
end

subsection Formulation
set Formulation = Boussinesq approximation
end

subsection Solver parameters
set Temperature solver tolerance = 1e-10
subsection Stokes solver parameters
set Stokes solver type = block GMG
set GMRES solver restart length = 1000
set Number of cheap Stokes solver steps = 20000
set Use full A block as preconditioner = true
set Linear solver tolerance = 1e-7
set Maximum number of expensive Stokes solver steps = 0
end
end

# 10 km x 10 km box
subsection Geometry model
set Model name = box

subsection Box
set X extent = 7.5e3
set Y extent = 7.5e3

set X repetitions = 10
set Y repetitions = 10
end
end

# Uniform temperature of 293 K
subsection Initial temperature model
set Model name = function

subsection Function
set Function expression = 293
end
end

subsection Adiabatic conditions model
set Model name = function

subsection Function
set Variable names = z
set Function expression = 293; 3.3e4*z; 3300
end
end

subsection Boundary temperature model
set Allow fixed temperature on outflow boundaries = false
set List of model names = box
set Fixed temperature boundary indicators = top, bottom, right, left

subsection Box
set Bottom temperature = 293
set Top temperature = 293
set Left temperature = 293
set Right temperature = 293
end
end

subsection Compositional fields
set Number of fields = 3
set Names of fields = porosity, bound_fluid, solid_phase
set Compositional field methods = darcy field, field, field
end

subsection Melt settings
set Include melt transport = false
end

# Initialize a porosity of 1% (0.01) everywhere.
subsection Initial composition model
set Model name = function

subsection Function
set Variable names = x,y,t
set Function constants = pi=3.1415926, x0=3750, y0=3750, c=1000
set Function expression = 0.01 * exp(-((x-x0)*(x-x0)+(y-y0)*(y-y0))/(2*c*c)); 0.0; if(y>=7000, 1, 0)
end
end

subsection Boundary composition model
set Fixed composition boundary indicators = bottom, top
set Allow fixed composition on outflow boundaries = false
set List of model names = initial composition
end

subsection Material model
set Model name = reactive fluid transport
set Material averaging = geometric average only viscosity

subsection Reactive Fluid Transport Model
set Base model = visco plastic
set Reference fluid density = 1000 # density of water
set Shear to bulk viscosity ratio = 1
set Reference fluid viscosity = 10
set Reference permeability = 1e-6
set Exponential fluid weakening factor = 0
set Fluid compressibility = 0
set Fluid reaction time scale for operator splitting = 1e50
set Fluid-solid reaction scheme = zero solubility
end

subsection Visco Plastic
set Viscosity averaging scheme = geometric
set Viscous flow law = diffusion
set Prefactors for diffusion creep = 5e-21
set Stress exponents for diffusion creep = 1.0
set Activation energies for diffusion creep = 0
set Activation volumes for diffusion creep = 0
set Densities = 3000

set Angles of internal friction = 0
set Cohesions = 1e50

set Minimum viscosity = 1e21
set Maximum viscosity = 1e21
set Thermal expansivities = 0
end
end

# Set the global refinement to 0, 1 km x 1 km mesh.
subsection Mesh refinement
set Coarsening fraction = 0.1
set Refinement fraction = 0.9
set Initial adaptive refinement = 2
set Initial global refinement = 0
set Strategy = composition threshold, minimum refinement function
set Time steps between mesh refinement = 1

# Minimum of 4 global refinements
subsection Minimum refinement function
set Function expression = 0
end

# Refine where the porosity is bigger than 1e-6. Other compositions
# are set to 1e50 to ensure that we do not refine based on these
# compositions.
subsection Composition threshold
set Compositional field thresholds = 1e-3, 1e50, 1e50
end
end


# Output the darcy velocity
subsection Postprocess
set List of postprocessors = velocity statistics, visualization

subsection Visualization
set List of output variables = darcy velocity
set Time between graphical output = 250
set Output format = vtu
end
end
Loading

0 comments on commit c2f6a2c

Please sign in to comment.