OK, thank you! I am getting much closer. Elmer does run with that change, but the output does not have the glacier geometry, only the bedrock. I can get results with the correct glacier geometry using the following sif:

- Code: Select all
`Header`

Mesh DB "." "testglacier2"

Include Path ""

Results Directory ""

End

Simulation

Max Output Level = 4

Coordinate System = "Cartesian 2D"

Coordinate Mapping(3) = 1 2 3

Simulation Type = "Transient"

Steady State Max Iterations = 1

Timestepping Method = "BDF"

BDF Order = 1

Timestep Sizes = 1

Timestep Intervals = 3

Output Intervals = 1

Output File = "scarp_thirdtest.result"

Post File = "scarp_thirdtest.vtu" ! use .ep suffix for leagcy format

Initialize Dirichlet Conditions = Logical False

End

Constants

Gravity(4) = 0 -1 0 3.711

Stefan Boltzmann = 5.67e-08

End

Body 1

Name = "Glacier"

Body Force = 1

Equation = 1

Material = 1

Initial Condition = 1

End

Body 2

Name = "free surface"

Equation = 2

Material = 2

Body Force = 2

Initial Condition = 2

End

Equation 1

Name = "Equation1"

Convection = "computed"

Flow Solution Name = String "Flow Solution"

Active Solvers(5) = 1 2 3 4 6

End

Equation 2

Active Solvers(1) = 5

Flow Solution Name = String "Flow Solution"

Convection = Computed

End

Initial Condition 1

Velocity 1 = 0.0

Velocity 2 = 0.0

Pressure = 0.0

Depth = Real 0.0

End

Initial Condition 2

Zs = Variable Coordinate 2

Real cubic

include "icescarp70degground.dat"

End

ReferenceZs = Equals Coordinate 2

Mesh Update 2 = Real 0.0

Mesh Update 1 = Real 0.0

End

Solver 1

Exec Solver = "before simulation"

Equation = "MapCoordinate"

Procedure = "StructuredMeshMapper" "StructuredMeshMapper"

Active Coordinate = Integer 2 ! the mesh-update is y-direction

Mesh Velocity Variable = String "Mesh Velocity 2"

Mesh Velocity First Zero = Logical True

Dot Product Tolerance = Real 0.01

End

! FlowDepth on generally unstructured mesh

! (could be replaced by strucuterd version)

Solver 2

Equation = "Flowdepth"

Exec Solver = "Before Simulation"

Procedure = File "ElmerIceSolvers" "FlowDepthSolver"

Variable = String "Depth"

Variable DOFs = 1

Linear System Solver = "Direct"

Linear System Direct Method = "UMFPACK"

Linear System Max Iterations = 200

Linear System Convergence Tolerance = 1.0E-09

Linear System Abort Not Converged = False

Linear System Preconditioning = "ILU0"

Linear System Residual Output = 1

Gradient = Real -1.0E00

Exported Variable 1 = -dofs 3 "Mesh Velocity"

Calc Free Surface = Logical True

! the name for the exported (if not existing) added variable

! the gradients will be stored in variables with the base

! name given and "Grad1" and (in 3 dimensions) "Grad2" added,

! so in our case "FreeSurfGrad1" and "FreeSurfGrad2"

! again, if those variables did not exist, they will be

! automatically created

!-----------------------------------------------------------

Freesurf Name = String "FreeSurf"

End

Solver 3

Equation = "Navier-Stokes"

Optimize Bandwidth = Logical True

Linear System Solver = Direct

Linear System Direct Method = "UMFPACK"

Linear System Max Iterations = 5000

Linear System Convergence Tolerance = 1.0E-06

Linear System Abort Not Converged = False

Linear System Preconditioning = "ILU1"

Linear System Residual Output = 1

Steady State Convergence Tolerance = 1.0E-05

Stabilization Method = Stabilized

Nonlinear System Convergence Tolerance = 1.0E-04

Nonlinear System Convergence Measure = Solution

Nonlinear System Max Iterations = 50

Nonlinear System Newton After Iterations = 3

Nonlinear System Newton After Tolerance = 1.0E-01

End

Solver 4

Equation = Sij

Procedure = "ElmerIceSolvers" "ComputeDevStress"

Variable = -nooutput "Sij"

Variable DOFs = 1

Exported Variable 1 = -dofs 4 Stress

Stress Variable Name = String "Stress"

Flow Solver Name = String "Flow Solution"

Linear System Solver = Direct

Linear System Direct Method = umfpack

End

Solver 5

Exec Solver = "After TimeStep"

Equation = String "Free Surface Evolution"

Variable = "Zs"

Variable DOFs = 1

Procedure = "FreeSurfaceSolver" "FreeSurfaceSolver"

Apply Dirichlet = Logical True

Linear System Solver = Iterative

Linear System Iterative Method = BiCGStab

Linear System Max Iterations = 100

Linear System Preconditioning = ILU1

Linear System Convergence Tolerance = 1.0e-8

Nonlinear System Max Iterations = 100

Nonlinear System Min Iterations = 2

Nonlinear System Convergence Tolerance = 1.0e-6

Steady State Convergence Tolerance = 1.0e-4

Stabilization Method = Bubbles

Flow Solution Name = String "Flow Solution"

Exported Variable 1 = Zs Residual

Exported Variable 1 DOFS = 1

Exported Variable 2 = ReferenceZs

Exported Variable 2 DOFS = 1

End

Solver 6

Equation = "Mesh Update"

Exec Solver = "After TimeStep"

Linear System Solver = Iterative

Linear System Iterative Method = BiCGStab

Linear System Max Iterations = 500

Linear System Preconditioning = ILU1

Linear System Convergence Tolerance = 1.0e-6

Nonlinear System Max Iterations = 1

Nonlinear System Convergence Tolerance = 1.0e-06

End

! we use m-yr-MPa system 1 yr = 31556926.0 sec

Material 1

Name = "ice-ice-baby"

Density = Real $910.0*1.0E-06*(31556926.0)^(-2.0)

Viscosity Model = String "Glen"

Viscosity = Real 1.0

Glen Exponent = Real 3.0

Critical Shear Rate = Real 1.0e-10

Rate Factor 1 = Real 1.258e13

Rate Factor 2 = Real 6.046e28

Activation Energy 1 = Real 60e3

Activation Energy 2 = Real 139e3

Glen Enhancement Factor = Real 1.0

Mesh Elastic Modulus = 1.0

Mesh Poisson Ratio = 0.3

! the variable taken to evaluate the Arrhenius law

! in general this should be the temperature relative

! to pressure melting point. The suggestion below plugs

! in the correct value obtained with TemperateIceSolver

! Temperature Field Variable = String "Temp Homologous"

! the temperature to switch between the

! two regimes in the flow law

Limit Temperature = Real 0.0

Constant Temperature = Real MATC "220-273.15"

Cauchy = Logical False

End

Material 2

Density = Real $910.0*1.0E-06*(31556926.0)^(-2.0)

! Min Zs = Real 0.0

! Max Zs = Real 900.0

End

Body Force 1

Name = "BodyForce1"

Heat Source = 1

Flow BodyForce 1 = Real 0.0

Flow BodyForce 2 = Real -3.711e15 !MPa - a - m

End

Body Force 2

Zs Accumulation Flux 1 = Real 0.0

Zs Accumulation Flux 2 = Real 0.0

End

Boundary Condition 1

Name = "bedrock"

Target Boundaries = 1

Compute Normals = Logical True

Bottom Surface = Variable Coordinate 1

Real cubic

include "groundwide.dat"

End

Velocity 1 = Real 0.0e0

Velocity 2 = Real 0.0e0

Mesh Update 1 = Real 0.0

Mesh Update 2 = Real 0.0

End

Boundary Condition 2

Name = "sides"

Target Boundaries(2) = 3 4 ! combine left and right boundary

Velocity 1 = 0.0

Mesh Update 1 = Real 0.0

Mesh Update 2 = Real 0.0

End

Boundary Condition 3

Name = "surface"

Target Boundaries = 2

Body ID = 2

Top Surface = Variable Coordinate 1

Real cubic

include "icescarp70degground.dat"

End

! Top Surface = Equals Zs

Depth = Real 0.0

Mesh Update 1 = Real 0.0

Mesh Update 2 = Variable Zs, ReferenceZs

Real MATC "tx(0) - tx(1)"

End

However, I don't think this is quite correct, because as you pointed out before the surface Boundary Condition is pointing toward what I want to be the initial state. I have played around a lot with when the StructuredMeshMapper is executed like you suggested, and with a Mesh Update solver, but think I am probably using them incorrectly - is it necessary to have Mesh Update terms in the BC, and would it make sense to have the top surface term inside one of the solvers instead of the BC?