I was performing a transient simulation for 1 year on an Italian glacier using Elmer. I attacched my .sif file below. Basically my simulation performes the following tasks:

- It considers a timestep of about 1 month (30/365.25) and it divides the year into 2 periods: the winter with a zero velocity (BasalWinterVelocity) on the bedrock and the summer with a fixed velocity of 50 mm/d (BasalSummerVelocity)
- It applies during the winter an accumulation forcing and during the summer an ablation flux. Both values are retrieved from DAT files (balance_winter.dat and balance_summer.dat) and interpolated over the domain using the getAccumulation function.
- The mesh is initially extruded and deformed vertically with the StructuredMeshMapper solver. The FreeSurfaceSolver computes the elevation of the free surface.

I was using the Arrhenius factor as tuning parameter in order to match the measured surface velocity field of the glacier. So I performed 2 transient simulations using the same .sif file (reported below) and changing the Arrhenius factor only.

I computed the initial volume of the mesh using ParaView (Data Analysis -> Integrated Variables).

The strange fact I have noticed is that the initial mesh volume of the two simulations is different of a magnitude of 5 even if I haven't change anything except the Arrhenius factor. How is that possible if the assigned values of Zs for each node are the same in both the simulations? Analysing the ElmerPost file it seems that both meshes have the same elevation of the bedrock but not for the surface. The differences between Zs of the simulations are reported below for each nodes:

I would expect that the initial volume of the meshes is the same. Someone could give me some hints?

If you need any further files, just ask me. I appreciate your helps!

- Code: Select all
`!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!`

!! !!

!! simC !!

!! !!

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

! Simulation title

$title = "simC1_a"

! Conversion rate

$yearInSec = 365.25*86400

Header

Mesh DB "." "mesh"

End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

Constants

Water Density = Real $ 1000.0/(yearInSec^2)

! Sliding on the bedrock - 50 mm/d - 18.26 m/y

BasalSummerVelocity = Real $ 50*1.0e-3*365.25

! No sliding

BasalWinterVelocity = Real 0.0

! winter accumulation file

WinterAccumulationFile = String balance_winter.dat

! summer accumulation file

SummerAccumulationFile = String balance_summer.dat

End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

Simulation

Coordinate System = "Cartesian 3D"

Simulation Type = Transient

Timestepping Method = "bdf"

BDF Order = 1

TimeStep Intervals = 12

Timestep Sizes = Real $30/365.25

Output Intervals = 1

Steady State Min Iterations = 1

Steady State Max Iterations = 100

Output File = "$title".result"

Post File = "$title".ep"

Max Output Level = 3

Extruded Mesh Levels = Integer 10

End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

Body 1

Equation = 1

Material = 1

Body Force = 1

Initial condition = 1

End

Body 2

Equation = 2

Body Force = 2

Material = 1

Initial Condition = 2

End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

Body Force 1

Flow BodyForce 1 = Real 0.0

Flow BodyForce 2 = Real 0.0

Flow Bodyforce 3 = Real $ -9.81*yearInSec^2

End

Body Force 2

Zs Accumulation = Variable Coordinate 1

Real Procedure "./utility" "getAccumulation"

End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

Initial Condition 1

Pressure = Real 0.00

Velocity 1 = Real 0.00

Velocity 2 = Real 0.00

Velocity 3 = Real 0.00

End

Initial Condition 2

Zs = Variable Coordinate 1

Real Procedure "./utility" "getTopSurface"

End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

Boundary Condition 1 ! lateral side

Target Boundaries = 1

Save Scalars = Logical True

End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

Boundary Condition 2 ! bedrock

Normal-Tangential Velocity = Logical True

Velocity 1 = 0.0

Velocity 2 = Real Procedure "./utility" "getBottomVelocity"

Velocity 3 = Real Procedure "./utility" "getBottomVelocity"

Bottom Surface = Variable Coordinate 1

Real Procedure "./utility" "getBedrock"

Save Scalars = Logical True

End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

Boundary Condition 3 ! upper surface

Body Id = 2

Save Scalars = Logical True

End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

Material 1

Density = Real $ 917.0/(yearInSec^2)

Viscosity Model = String "glen"

Viscosity = 1.0 ! Avoid warning output. This value will not be used.

Limit Temperature = Real -10.0

Rate Factor 1 = Real $ 3.9855e-13*yearInSec

Rate Factor 2 = Real $ 1.916e3*yearInSec

Activation Energy 1 = Real 60e3

Activation Energy 2 = Real 139e3

Glen Enhancement Factor = Real 1.0

Critical Shear Rate = Real $ 0.5e-6*yearInSec

Constant Temperature = Real 0.0

! Bed condition

Min Zs = Variable Coordinate 1

Real Procedure "./utility" "getMinTopSurface"

Max Zs = Real +1.0e10

End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

Solver 1

Equation = "MapCoordinate"

Procedure = "StructuredMeshMapper" "StructuredMeshMapper"

Active Coordinate = Integer 3

Mesh Velocity Variable = String "dSdt"

Mesh Update Variable = String "dS"

Mesh Velocity First Zero = Logical True

Dot Product Tolerance = Real 1.0e-6

Top Surface Variable Name = String "Zs"

End

Solver 2

Equation = Navier-Stokes

Stabilization Method = String Stabilized

Flow Model = Stokes

Exported Variable 1 = -dofs 1 "dSdt"

Exported Variable 2 = -dofs 1 "dS"

Linear System Solver = Direct

Linear System Direct Method = umfpack

Nonlinear System Max Iterations = 50

Nonlinear System Convergence Tolerance = 1.0e-5

! first 5 iterations or convergence < 1.e-2 before moving from Picard to Newton

Nonlinear System Newton After Iterations = 5

Nonlinear System Newton After Tolerance = 1.0e-02

! reset newton to false each new time step

Nonlinear System Reset Newton = Logical True

Steady State Convergence Tolerance = Real 1.0e-3

End

Solver 3

Equation = "Free Surface"

Variable = String Zs

Variable DOFs = 1

Exported Variable 1 = String "Zs Residual"

Exported Variable 1 DOFs = 1

Procedure = "FreeSurfaceSolver" "FreeSurfaceSolver"

Before Linsolve = "EliminateDirichlet" "EliminateDirichlet"

Linear System Solver = Iterative

Linear System Max Iterations = 1500

Linear System Iterative Method = BiCGStab

Linear System Preconditioning = ILU1

Linear System Convergence Tolerance = Real 1.0e-5

Linear System Abort Not Converged = True

Linear System Residual Output = 3

Nonlinear System Max Iterations = 20

Nonlinear System Convergence Tolerance = 1.0e-6

Nonlinear System Relaxation Factor = 1.0

Steady State Convergence Tolerance = 1.0e-03

Stabilization Method = Bubbles

Apply Dirichlet = Logical True

! How much the free surface is relaxed

Relaxation Factor = Real 1.0

! Is there a maximum step-size for the displacement use/or not accumulation

Use Accumulation = Logical True

! take accumulation to be given normal to surface/as vector

Normal Flux = Logical False

End

Solver 4

Exec Solver = After TimeStep

Exec Interval = 1

Equation = "result output"

Procedure = "ResultOutputSolve" "ResultOutputSolver"

Output File Name = File "$title".vtu"

Output Format = String vtu

End

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

Equation 1

Name = "Navier-Stokes"

Active Solvers(3) = 1 2 4

End

Equation 2

Active Solvers(1) = 3

Flow Solution Name = String "Flow Solution"

Convection = String Computed

End