!!--------------------------------------------------------!! ! Flowline template ! !!--------------------------------------------------------!! check keywords warn echo on $namerun = "Experiment2l2d" $timestep = 0.1 ! ! working units are MPa, a, m ! $yearinsec = 365.25*24*60*60 $rhoi = 910.0/(1.0e6*yearinsec^2) $rhow = 1000.0/(1.0e6*yearinsec^2) $A = 3e-25*yearinsec*1.0e18 $Ascaling = 1.0 $Afactor = A*Ascaling $icetemp = -15.0 $critshearrate= 1.0e-10*yearinsec $n = 3.0 $eta = 1.0/(2.0*Afactor)^(1.0/n) $gravity = -9.8*yearinsec^2 $C = 7.624e6/(1.0e6*yearinsec^(1.0/n)) $GLTolerance = 1.0e-4 $H_init = 100.0 ! initial ice thickness !$Afactor = 80.0 $etai = 1.0/(2*Afactor)^(1.0/n) $W= 100.0e3 $Kspring = etai* (n+1)^(1/n) / (rhoi * W^(1+(1/n))) $W2= 300.0e3 $Kspring2= etai* (n+1)^(1/n) / (rhoi * W2^(1+(1/n))) $m = 1.0/3.0 $function bedrock(x) {\ _bedrock = ( 729.0 - 2184.8 * (x/750000.0)^2 + 1031.72 * (x/750000.0)^4 - 151.72 * (x/750000.0)^6 ) \ } ! lower_surf_init - the initial lower surface. Depths are negative downwards from sea level. $function lower_surf_init(x) import rhoi, rhow, H_init { \ _lower_surf_init = 0.0 ;\ bed_depth = ( 729.0 - 2184.8 * (x/750000.0)^2 + 1031.72 * (x/750000.0)^4 - 151.72 * (x/750000.0)^6 ) ;\ flot_depth = -1.0 * H_init * rhoi / rhow;\ if (bed_depth > flot_depth) { \ _lower_surf_init = bed_depth; \ } else { \ _lower_surf_init = flot_depth; \ } \ } ! upper_surf_init - the initial upper surface. Depths are negative downwards from sea level. $function upper_surf_init(x) import rhoi, rhow, H_init { \ bed_depth = ( 729.0 - 2184.8 * (x/750000.0)^2 + 1031.72 * (x/750000.0)^4 - 151.72 * (x/750000.0)^6 ) ;\ flot_depth = -1.0 * H_init * rhoi / rhow ;\ if (bed_depth > flot_depth) {\ _upper_surf_init = bed_depth + H_init;\ } else { \ _upper_surf_init = flot_depth + H_init ;\ }\ } $ function waterpressure(Z) import gravity, rhow {\ waterline = 0.0;\ _waterpressure = 0.0;\ if (Z>waterline) {\ _waterpressure = 0.0;\ }else {\ _waterpressure = 1.0 * rhow * gravity * (waterline - Z);\ }\ } Header Mesh DB "../" "mesh2d" End Constants Water Density = Real $rhow Gas Constant = Real 8.314 !Joule/mol x K End !--------------------------------------------------- !---------------- SIMULATION ----------------------- !--------------------------------------------------- Simulation Coordinate System = Cartesian 2D Simulation Type = transient Timestepping Method = "bdf" BDF Order = 1 Timestep Intervals = 10000 Output Intervals = 20 Timestep Sizes = Real $timestep Initialize Dirichlet Conditions = Logical False Steady State Max Iterations = 1 Steady State Min Iterations = 1 Output File = "$namerun".result" Post File = "$namerun".vtu" Restart File = "Experiment1b.result" Restart Time = Real 19995.0 max output level = 3 End !--------------------------------------------------- !---------------- BODIES --------------------------- !--------------------------------------------------- ! the ice Body 1 Name = "ice" Equation = 1 Body Force = 1 Material = 1 Initial Condition = 1 End ! The upper surface Body 2 Name= "top free surface" Equation = 2 Material = 2 Body Force = 2 Initial Condition = 2 End ! the lower surface Body 3 Name= "free surface sea/ice-shelf" Equation = 3 Material = 3 Body Force = 3 Initial Condition = 3 End !--------------------------------------------------- !---------------- INITIAL CONDITIONS --------------- !--------------------------------------------------- !! for ice Initial Condition 1 Pressure = Real 0.0 Velocity 1 = Real 0.0 Velocity 2 = Real 0.0 Mesh Velocity 1 = Real 0.0 Mesh Velocity 2 = Real 0.0 !Bedrock = Variable Coordinate 1 ! Real MATC "720 - 2184.8*(tx/7.5E5)^2 + 1031.72*(tx/7.5E5)^4 - 151.72*(tx/7.5E5)^6" End !! for top free surface Initial Condition 2 Zs = Variable Coordinate 1 Real MATC "upper_surf_init(tx)" ReferenceZs = Variable Coordinate 1 Real MATC "upper_surf_init(tx)" End !! for free surface sea/ice-shelf Initial Condition 3 Zb = Variable Coordinate 1 Real MATC "lower_surf_init(tx)" ReferenceZb = Variable Coordinate 1 Real MATC "lower_surf_init(tx)" End !--------------------------------------------------- !---------------- BODY FORCES ---------------------- !--------------------------------------------------- Body Force 1 ! Flow BodyForce 1 = Real 0.0 Flow BodyForce 1 = Variable Coordinate 1 Real Procedure "ElmerIceUSF" "LateralFriction_x" ! Flow BodyForce 2 = Real $gravity Flow BodyForce 2 = Variable Coordinate 1 Real Procedure "ElmerIceUSF" "LateralFriction_y" Lateral Friction Gravity 1 = Real 0.0 Lateral Friction Gravity 2 = Real $gravity Lateral Friction Coefficient = variable time Real 0 0 5000 0 10000 $Kspring 15000 $Kspring 20000 $Kspring 20010 $Kspring2 30000 $Kspring2 End Lateral Friction Exponent = Real $1.0/n Flow Solver Name = String Flow Solution End !! accumulation flux in m/year Body Force 2 Zs Accumulation Flux 1 = Real 0.0e0 Zs Accumulation Flux 2 = Real 0.3e0 !m/a End !! no melting/accretion under ice/shelf Body Force 3 Zb Accumulation Flux 1 = Real 0.0e0 Zb Accumulation Flux 2 = Real 0.0e0 End !--------------------------------------------------- !---------------- MATERIALS ------------------------ !--------------------------------------------------- !! ice material properties in MPa - m - a system Material 1 Viscosity Model = String "Glen" Density = Real $rhoi Viscosity = Real 1.0 Glen Exponent = Real 3.0 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 3.0 Limit Temperature = Real $icetemp Constant Temperature = Real $icetemp Critical Shear Rate = Real $critshearrate Sea level = Real 0.0 Cauchy = Logical True End Material 2 Density = Real $rhoi Min Zs = Real 0.0 Max Zs = Variable ReferenceZs Real MATC "tx(0) + 10000.0" End Material 3 Density = Real $rhoi Min Zb = Variable Coordinate 1 Real MATC "bedrock(tx)" Max Zb = Real 1000.0 Min Zs Bottom = Variable Coordinate 1 Real MATC "bedrock(tx)" ! Bedrock variable for grounded solver End !--------------------------------------------------- !---------------- SOLVERS -------------------------- !--------------------------------------------------- Solver 1 !Exec Solver = "never" Exec Solver = Before All Equation = "MapCoordinateIni" Procedure = "StructuredMeshMapper" "StructuredMeshMapper" Active Coordinate = Integer 2 Mesh Velocity Variable = String "dSdt" Mesh Update Variable = String "dS" Mesh Velocity First Zero = Logical True Top Surface Variable Name = String "Zs" Bottom Surface Variable Name = String "Zb" Displacement Mode = Logical False Correct Surface = Logical True Minimum Height = Real 1.0 End Solver 2 !Exec Solver = "never" Exec Solver = Before All Equation = GroundedMaskIni Procedure = "ElmerIceSolvers" "GroundedSolver" Variable = GroundedMask Variable DOFs = 1 Toler = Real 1.0e-3 Exported Variable 1 = -dofs 3 "Mesh Velocity" End Solver 3 !Exec Solver = "never" Equation = "NormalVector" Procedure = "ElmerIceSolvers" "ComputeNormalSolver" Variable = String "Normal Vector" Variable DOFs = 2 ComputeAll = Logical False Optimize Bandwidth = Logical False End Solver 4 !Exec Solver = "never" Equation = Fw Procedure = "ElmerIceSolvers" "GetHydrostaticLoads" Variable = Fw[Fwater:2] Variable DOFs = 2 End Solver 5 !Exec Solver = "never" Equation = "Navier-Stokes" Linear System Solver = Direct Linear System Direct Method = Mumps Nonlinear System Max Iterations = 50 Nonlinear System Convergence Tolerance = 1.0e-5 Nonlinear System Newton After Iterations = 50 Nonlinear System Newton After Tolerance = 1.0e-05 Nonlinear System Relaxation Factor = 1.00 Nonlinear System Reset Newton = Logical True Steady State Convergence Tolerance = Real 1.0e-4 Stabilization Method = String Bubbles Exported Variable 1 = Flow Solution Loads[Stress Vector:2 CEQ Residual:1] Calculate Loads = Logical True Exported Variable 2 = -dofs 1 "dSdt" Exported Variable 3 = -dofs 1 "dS" Exported Variable 4 = -dofs 1 "Bedrock" Flow Model = String "Stokes" End Solver 6 Equation = String "StressSolver" Procedure = File "ElmerIceSolvers" "ComputeDevStress" ! this is just a dummy, hence no output is needed !----------------------------------------------------------------------- Variable = -nooutput "Sij" Variable DOFs = 1 ! the name of the variable containing the flow solution (U,V,W,Pressure) !----------------------------------------------------------------------- Flow Solver Name = String "Flow Solution" ! no default value anymore for "Stress Variable Name" Stress Variable Name = String "Stress" !----------------------------------------------------------------------- Exported Variable 1 = "Stress" ! [Sxx, Syy, Szz, Sxy] in 2D ! [Sxx, Syy, Szz, Sxy, Syz, Szx] in 3D Exported Variable 1 DOFs = 4 ! 4 in 2D, 6 in 3D Linear System Solver = "Iterative" Linear System Iterative Method = "BiCGStab" Linear System Max Iterations = 300 Linear System Convergence Tolerance = 1.0E-09 Linear System Abort Not Converged = True Linear System Preconditioning = "ILU0" Linear System Residual Output = 1 End Solver 7 !Exec Solver = "never" Equation = "HeightDepth" Procedure = "StructuredProjectToPlane" "StructuredProjectToPlane" Active Coordinate = Integer 2 Dot Product Tolerance = Real 1.0e-3 Operator 1 = Depth Operator 2 = Height ! Export Zb on the Upper surface Variable 3 = Zb Operator 3 = Bottom End Solver 8 !Exec Solver = "never" Equation = "Free Surface Top" Procedure = "FreeSurfaceSolver" "FreeSurfaceSolver" Variable = "Zs" Variable DOFs = 1 Exported Variable 1 = "Zs Residual" Exported Variable 1 DOFs = 1 Exported Variable 2 = "ReferenceZs" Exported Variable 2 DOFs = 1 !Before Linsolve = "EliminateDirichlet" "EliminateDirichlet" Linear System Solver = Iterative Linear System Direct Method = Mumps Linear System Max Iterations = 1500 Linear System Iterative Method = BiCGStab Linear System Preconditioning = ILU0 Linear System Convergence Tolerance = Real 1.0e-6 Linear System Abort Not Converged = False Linear System Residual Output = 1 Nonlinear System Max Iterations = 100 Nonlinear System Convergence Tolerance = 1.0e-5 Nonlinear System Relaxation Factor = 1.00 Steady State Convergence Tolerance = 1.0e-03 Stabilization Method = Bubbles Apply Dirichlet = Logical True Relaxation Factor = Real 1.0 End Solver 9 !Exec Solver = "never" Equation = "Free Surface Sea/Shelf" Procedure = "FreeSurfaceSolver" "FreeSurfaceSolver" Variable = "Zb" Variable DOFS = 1 Exported Variable 1 = "Zb Residual" Exported Variable 1 DOFs = 1 Exported Variable 2 = "ReferenceZb" Exported Variable 2 DOFs = 1 !Before Linsolve = "EliminateDirichlet" "EliminateDirichlet" Linear System Solver = Iterative Linear System Max Iterations = 1500 Linear System Direct Method = Mumps Linear System Iterative Method = BiCGStab Linear System Preconditioning = ILU0 Linear System Convergence Tolerance = Real 1.0e-8 Linear System Abort Not Converged = False Linear System Residual Output = 1 Nonlinear System Max Iterations = 100 Nonlinear System Convergence Tolerance = 1.0e-5 Nonlinear System Relaxation Factor = 1.00 Steady State Convergence Tolerance = 1.0e-03 Stabilization Method = Bubbles Apply Dirichlet = Logical True Relaxation Factor = Real 1.0 End Solver 10 Exec Solver = "After timestep" Equation = "MapCoordinate" Procedure = "StructuredMeshMapper" "StructuredMeshMapper" Active Coordinate = Integer 2 Mesh Velocity Variable = String "dSdt" Mesh Update Variable = String "dS" ! Mesh Velocity First Zero = Logical True Top Surface Variable Name = String "Zs" Bottom Surface Variable Name = String "Zb" Displacement Mode = Logical False Correct Surface = Logical True Minimum Height = Real 1.0 End !! Compute the Mask Solver 11 Exec Solver = "After timestep" Equation = GroundedMask Procedure = "ElmerIceSolvers" "GroundedSolver" Variable = GroundedMask Variable DOFs = 1 Toler = Real 1.0e-3 ! Bedrock Variable = String "Bedrock" End Solver 12 Exec Solver = After TimeStep Equation = "Save Scalars" Procedure = File "SaveData" "SaveScalars" Filename = "results.dat" File Append = Logical False Variable 1 = String "Time" Variable 2 = String "flow solution" Operator 2 = String "Volume" Variable 3 = String "Velocity 1" Operator 3 = String "Max" Variable 4 = String "GroundedMask" Operator 4 = String "Sum" Operator 5 = String "cpu time" End Solver 13 !Exec Solver = Never equation = "SubShelfMeltCalculation" procedure = "ElmerIceSolvers" "SubShelfMelt" Variable = meltRate melt function = String bgf_scaled lower surface variable name = String Zb bedrock variable name = String bedrock grounded mask name = string GroundedMask grounding line melt = logical True omega = Real 0.2 far field ocean temperature = Real 1.0 water column scaling = Real 75.0 !ice base scaling depth = Real 100.0 !ice base scaling offset = Real 40.0 iceBaseScaling = logical False End !--------------------------------------------------- !---------------- EQUATIONS ------------------------ !--------------------------------------------------- Equation 1 Active Solvers (7) = 1 3 5 6 7 10 12 End Equation 2 Active Solvers(1) = 8 Flow Solution Name = String "Flow Solution" Convection = String Computed End Equation 3 Active Solvers(5) = 2 4 9 11 13 Flow Solution Name = String "Flow Solution" Convection = String Computed End !--------------------------------------------------- !---------------- BOUNDARY CONDITIONS -------------- !--------------------------------------------------- !! BC Bedrock + Shelf Boundary Condition 1 Name = "bottom" Target Boundaries = 1 Body Id = 3 ComputeNormal = Logical True Height = Real 0.0 Normal-Tangential Velocity = Logical True Flow Force BC = Logical True ComputeNormal Condition = Variable GroundedMask Real MATC "tx + 0.5" ! ! Condition where the bed is stuck ! Zb = Variable Coordinate 1 Real MATC "bedrock(tx)" Zb Condition = Variable GroundedMask Real MATC "tx + 0.5" ! ! Bedrock conditions ! Slip Coefficient 2 = Variable Coordinate 1 Real Procedure "ElmerIceUSF" "SlidCoef_Contact" Sliding Law = String "Coulomb" !! Parameters needed for the Coulomb Friction Law Friction Law Sliding Coefficient = Real 4.1613e5 Friction Law Post-Peak Exponent = Real 1.0 !(q=1) Friction Law Maximum Value = Real 0.1 !(C=1) Friction Law PowerLaw Exponent = Real 3.0 !(m = n = 3 Glen's law) Friction Law Linear Velocity = Real 0.01 !Sliding Law = String "Weertman" !Weertman Friction Coefficient = Real $C !Weertman Exponent = Real $(1.0/n) !Weertman Linear Velocity = Real 1.0 ! Options are 'Last Grounded' (default), 'First Floating' or 'Discontinuous' ! Grounding Line Definition = String "Last Grounded" Grounding Line Definition = String "Discontinuous" ! Grounding Line Definition = String "First Floating" Test Contact Tolerance = real 1.0e-3 Velocity 1 = Real 0.0 Velocity 1 Condition = Variable GroundedMask Real MATC "tx + 0.5" ! ! Shelf conditions ! External Pressure = Variable Coordinate 2 Real Procedure "ElmerIceUSF" "SeaPressure" Slip Coefficient 1 = Variable Coordinate 2 Real Procedure "ElmerIceUSF" "SeaSpring" ComputeNormal Condition = Variable GroundedMask Real MATC "tx + 0.5" Compute Sea Pressure = Logical True Compute Sea Spring = Logical True End !! BC Lateral Ice-Shelf (air or sea contact) Boundary Condition 2 Name = "front" Target Boundaries = 2 Flow Force BC = Logical True External Pressure = Variable Coordinate 2 Real Procedure "ElmerIceUSF" "SeaPressure" Compute Sea Pressure = Logical True ComputeNormal = Logical False End !! BC Free surface Top Boundary Condition 3 Name = "top" Target Boundaries = 3 Body Id = 2 Depth = Real 0.0 End !! Symetry axis Boundary Condition 4 Name = "back" Target Boundaries = 4 Velocity 1 = Real 0.0e0 End RUN