Page 1 of 1

Glacier volume

Posted: 17 May 2014, 11:14
by sstefano
Hello

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:
  1. 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)
  2. 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.
  3. 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:

Image

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

Re: Glacier volume

Posted: 17 May 2014, 11:39
by raback
Hi Stefano

Cannot really comment on the physics but for computing the volumes there is also the "volume" operator in SaveScalars. Then you know that you get a consistant integral volume.

-Peter