Hi, I’m a student and I’m approaching Elmer/Ice, trying to simulate the evolution of an existing glacier over time.
In particular, I use, as initial mesh, the extrusion in 10 layers of the glacier’s footprint, deformed following existing topographic data. As body force on the upper surface I give a summer ablation and a winter accumulation reading from a DAT file, and it works. A 0-velocity condition on the bedrock also works. But…
…I have a problem when I try to give a more realisic sliding condition on the bottom of the glacier. In particular with the help of Elmer courses’s material I tried to give a condition like this
Flow Force BC = Logical True
Normal-Tangential Velocity = Logical True
Velocity 1 = Real 0.0
Slip Coefficient 2 = Real 0.1
Slip Coefficient 3 = Real 0.1
But it just diverges. I just can’t solve this question. Anybody can kindly check my sif file to find anything I missed?
Here it is
Thank you
Bedrock Sliding Boundary Condition
Re: Bedrock Sliding Boundary Condition
Hi,
Have you tried to use only Picard iterations (I see in your Navier-Stokes solver that you switch to Newton after 5 iterations)?
This is the first thing to try when the non linear iterations diverge.
cheers
Fabien
Have you tried to use only Picard iterations (I see in your Navier-Stokes solver that you switch to Newton after 5 iterations)?
This is the first thing to try when the non linear iterations diverge.
cheers
Fabien
Re: Bedrock Sliding Boundary Condition
Thank you for the reply. Following your suggestion I set the Navier-Stokes solver like following in order to avoid newton method
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-02
Nonlinear System Relaxation Factor = 1.00
Nonlinear System Reset Newton = Logical True
But unfortunately this didn't solve the problem. Looking at the terminal output I can see that the calculation diverges AFTER the first Navier-Stokes Solver. I post the terminal output below. It seems that it crashes when using my external function "getBalance". BUT that function works perfectly if I restore the boundary condition on the bedrock to "0-velocity". I attach also the code of my external functions named "utils.f90". I thank you and anybody else that can spend even a couple of minutes with my problem.
Here is the function code And here the output:
ELMER SOLVER (v 7.0) STARTED AT: 2015/08/22 11:53:45
MAIN:
MAIN: =============================================================
MAIN: ElmerSolver finite element software, Welcome!
MAIN: This program is free software licensed under (L)GPL
MAIN: Copyright 1st April 1995 - , CSC - IT Center for Science Ltd.
MAIN: Webpage http://www.csc.fi/elmer, Email elmeradm@csc.fi
MAIN: Library version: 7.0 (Rev: 1b38bc5)
MAIN: HYPRE library linked in.
MAIN: =============================================================
MAIN:
MAIN:
MAIN: -------------------------------------
MAIN: Reading Model: 1997-2007_test_Slip_22-08-2015.sif
Model Input: Unlisted keyword: [water density] in section: [constants]
Model Input: Unlisted keyword: [surfacedem] in section: [constants]
Model Input: Unlisted keyword: [bedrockdem] in section: [constants]
Model Input: Unlisted keyword: [switchperiod] in section: [constants]
Model Input: Unlisted keyword: [winterarrheniusfactor] in section: [constants]
Model Input: Unlisted keyword: [summerarrheniusfactor] in section: [constants]
Model Input: Unlisted keyword: [basalwintervelocity] in section: [constants]
Model Input: Unlisted keyword: [basalsummervelocity] in section: [constants]
Model Input: Unlisted keyword: [winteraccumulationfile] in section: [constants]
Model Input: Unlisted keyword: [summeraccumulationfile] in section: [constants]
Model Input: Unlisted keyword: [minicethickness] in section: [constants]
Model Input: Unlisted keyword: [zs accumulation] in section: [body force 2]
Loading user function library: [./utils_eros_fede]...[getBalance]
Model Input: Unlisted keyword: [zs] in section: [initial condition 2]
Loading user function library: [./utils_eros_fede]...[getSurfaceElevation]
Loading user function library: [./utils_eros_fede]...[getBedrockElevation]
Model Input: Unlisted keyword: [set arrhenius factor] in section: [material 1]
Model Input: Unlisted keyword: [arrhenius factor] in section: [material 1]
Loading user function library: [./utils_eros_fede]...[getArrheniusFactor]
Model Input: Unlisted keyword: [min zs] in section: [material 1]
Loading user function library: [./utils_eros_fede]...[getMinSurfaceElevation]
Model Input: Unlisted keyword: [max zs] in section: [material 1]
Model Input: Unlisted keyword: [mesh velocity variable] in section: [solver 1]
Model Input: Unlisted keyword: [mesh update variable] in section: [solver 1]
Model Input: Unlisted keyword: [mesh velocity first zero] in section: [solver 1]
Model Input: Unlisted keyword: [top surface variable name] in section: [solver 1]
Model Input: Unlisted keyword: [displacement mode] in section: [solver 1]
Model Input: Unlisted keyword: [correct surface] in section: [solver 1]
Model Input: Unlisted keyword: [minimum height] in section: [solver 1]
Model Input: Unlisted keyword: [stress variable name] in section: [solver 4]
Model Input: Unlisted keyword: [flow solver name] in section: [solver 4]
Model Input: Unlisted keyword: [relaxation factor] in section: [solver 6]
MAIN:
MAIN: -------------------------------------
MAIN: Time: 1/24 2592000.0000000000
MAIN: -------------------------------------
MAIN:
SolveEquations: -------------------------------------
SolveEquations: Coupled system iteration: 1
SolveEquations: -------------------------------------
StructuredMeshMapper: >Correct Surface< in case of intersecting upper and lower surface
StructuredMeshMapper: Adjusting upper surface to maintain minimum height to 0.5000E+00
getArrheniusFactor: Set value of Arrhenius Factor to 2.4000000E-24
ComputeChange: NS (ITER=1) (NRM,RELC): ( 0.66460863E-04 2.0000000 ) :: navier-stokes
ComputeChange: NS (ITER=2) (NRM,RELC): ( 0.10736614E-03 0.47064352 ) :: navier-stokes
ComputeChange: NS (ITER=3) (NRM,RELC): ( 0.14503338E-03 0.29847314 ) :: navier-stokes
ComputeChange: NS (ITER=4) (NRM,RELC): ( 0.17635048E-03 0.19488907 ) :: navier-stokes
ComputeChange: NS (ITER=5) (NRM,RELC): ( 0.20055361E-03 0.12843124 ) :: navier-stokes
ComputeChange: NS (ITER=6) (NRM,RELC): ( 0.21831422E-03 0.84802943E-01 ) :: navier-stokes
ComputeChange: NS (ITER=7) (NRM,RELC): ( 0.23089744E-03 0.56023585E-01 ) :: navier-stokes
ComputeChange: NS (ITER=8) (NRM,RELC): ( 0.23960695E-03 0.37021991E-01 ) :: navier-stokes
ComputeChange: NS (ITER=9) (NRM,RELC): ( 0.24554298E-03 0.24470931E-01 ) :: navier-stokes
ComputeChange: NS (ITER=10) (NRM,RELC): ( 0.24954789E-03 0.16178463E-01 ) :: navier-stokes
ComputeChange: NS (ITER=11) (NRM,RELC): ( 0.25223184E-03 0.10697718E-01 ) :: navier-stokes
ComputeChange: NS (ITER=12) (NRM,RELC): ( 0.25402263E-03 0.70747007E-02 ) :: navier-stokes
ComputeChange: NS (ITER=13) (NRM,RELC): ( 0.25726561E-03 0.12685526E-01 ) :: navier-stokes
ComputeChange: NS (ITER=14) (NRM,RELC): ( 0.25753787E-03 0.10577202E-02 ) :: navier-stokes
ComputeChange: NS (ITER=15) (NRM,RELC): ( 0.25755675E-03 0.73293926E-04 ) :: navier-stokes
ComputeChange: NS (ITER=16) (NRM,RELC): ( 0.25755835E-03 0.62131185E-05 ) :: navier-stokes
ComputeChange: SS (ITER=1) (NRM,RELC): ( 0.25755835E-03 2.0000000 ) :: navier-stokes
ComputeChange: NS (ITER=1) (NRM,RELC): ( 394880.79 2.0000000 ) :: sij
ComputeChange: NS (ITER=2) (NRM,RELC): ( 421559.88 0.65354617E-01 ) :: sij
ComputeChange: NS (ITER=3) (NRM,RELC): ( 431568.68 0.23463770E-01 ) :: sij
ComputeChange: NS (ITER=4) (NRM,RELC): ( 309568.48 0.32922437 ) :: sij
ComputeChange: NS (ITER=5) (NRM,RELC): ( 235878.41 0.27020070 ) :: sij
ComputeChange: NS (ITER=6) (NRM,RELC): ( 239879.73 0.16820801E-01 ) :: sij
ComputeChange: SS (ITER=1) (NRM,RELC): ( 239879.73 2.0000000 ) :: sij
ComputeChange: SS (ITER=1) (NRM,RELC): ( 0.0000000 0.0000000 ) :: eigenstresses
getBalance: Set accumulation file to balance_winter.dat
3 0.4927E-01
6 0.4202E+01
9 0.3234E+01
12 0.1999E+01
15 0.2004E+01
18 0.2006E+01
21 0.2006E+01
...
Cutted by me
...
1026 0.1892+139
1029 0.6784+140
1032 0.1688+142
1032 0.1688+142
ERROR:: IterSolve: Failed convergence tolerances.
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-02
Nonlinear System Relaxation Factor = 1.00
Nonlinear System Reset Newton = Logical True
But unfortunately this didn't solve the problem. Looking at the terminal output I can see that the calculation diverges AFTER the first Navier-Stokes Solver. I post the terminal output below. It seems that it crashes when using my external function "getBalance". BUT that function works perfectly if I restore the boundary condition on the bedrock to "0-velocity". I attach also the code of my external functions named "utils.f90". I thank you and anybody else that can spend even a couple of minutes with my problem.
Here is the function code And here the output:
ELMER SOLVER (v 7.0) STARTED AT: 2015/08/22 11:53:45
MAIN:
MAIN: =============================================================
MAIN: ElmerSolver finite element software, Welcome!
MAIN: This program is free software licensed under (L)GPL
MAIN: Copyright 1st April 1995 - , CSC - IT Center for Science Ltd.
MAIN: Webpage http://www.csc.fi/elmer, Email elmeradm@csc.fi
MAIN: Library version: 7.0 (Rev: 1b38bc5)
MAIN: HYPRE library linked in.
MAIN: =============================================================
MAIN:
MAIN:
MAIN: -------------------------------------
MAIN: Reading Model: 1997-2007_test_Slip_22-08-2015.sif
Model Input: Unlisted keyword: [water density] in section: [constants]
Model Input: Unlisted keyword: [surfacedem] in section: [constants]
Model Input: Unlisted keyword: [bedrockdem] in section: [constants]
Model Input: Unlisted keyword: [switchperiod] in section: [constants]
Model Input: Unlisted keyword: [winterarrheniusfactor] in section: [constants]
Model Input: Unlisted keyword: [summerarrheniusfactor] in section: [constants]
Model Input: Unlisted keyword: [basalwintervelocity] in section: [constants]
Model Input: Unlisted keyword: [basalsummervelocity] in section: [constants]
Model Input: Unlisted keyword: [winteraccumulationfile] in section: [constants]
Model Input: Unlisted keyword: [summeraccumulationfile] in section: [constants]
Model Input: Unlisted keyword: [minicethickness] in section: [constants]
Model Input: Unlisted keyword: [zs accumulation] in section: [body force 2]
Loading user function library: [./utils_eros_fede]...[getBalance]
Model Input: Unlisted keyword: [zs] in section: [initial condition 2]
Loading user function library: [./utils_eros_fede]...[getSurfaceElevation]
Loading user function library: [./utils_eros_fede]...[getBedrockElevation]
Model Input: Unlisted keyword: [set arrhenius factor] in section: [material 1]
Model Input: Unlisted keyword: [arrhenius factor] in section: [material 1]
Loading user function library: [./utils_eros_fede]...[getArrheniusFactor]
Model Input: Unlisted keyword: [min zs] in section: [material 1]
Loading user function library: [./utils_eros_fede]...[getMinSurfaceElevation]
Model Input: Unlisted keyword: [max zs] in section: [material 1]
Model Input: Unlisted keyword: [mesh velocity variable] in section: [solver 1]
Model Input: Unlisted keyword: [mesh update variable] in section: [solver 1]
Model Input: Unlisted keyword: [mesh velocity first zero] in section: [solver 1]
Model Input: Unlisted keyword: [top surface variable name] in section: [solver 1]
Model Input: Unlisted keyword: [displacement mode] in section: [solver 1]
Model Input: Unlisted keyword: [correct surface] in section: [solver 1]
Model Input: Unlisted keyword: [minimum height] in section: [solver 1]
Model Input: Unlisted keyword: [stress variable name] in section: [solver 4]
Model Input: Unlisted keyword: [flow solver name] in section: [solver 4]
Model Input: Unlisted keyword: [relaxation factor] in section: [solver 6]
MAIN:
MAIN: -------------------------------------
MAIN: Time: 1/24 2592000.0000000000
MAIN: -------------------------------------
MAIN:
SolveEquations: -------------------------------------
SolveEquations: Coupled system iteration: 1
SolveEquations: -------------------------------------
StructuredMeshMapper: >Correct Surface< in case of intersecting upper and lower surface
StructuredMeshMapper: Adjusting upper surface to maintain minimum height to 0.5000E+00
getArrheniusFactor: Set value of Arrhenius Factor to 2.4000000E-24
ComputeChange: NS (ITER=1) (NRM,RELC): ( 0.66460863E-04 2.0000000 ) :: navier-stokes
ComputeChange: NS (ITER=2) (NRM,RELC): ( 0.10736614E-03 0.47064352 ) :: navier-stokes
ComputeChange: NS (ITER=3) (NRM,RELC): ( 0.14503338E-03 0.29847314 ) :: navier-stokes
ComputeChange: NS (ITER=4) (NRM,RELC): ( 0.17635048E-03 0.19488907 ) :: navier-stokes
ComputeChange: NS (ITER=5) (NRM,RELC): ( 0.20055361E-03 0.12843124 ) :: navier-stokes
ComputeChange: NS (ITER=6) (NRM,RELC): ( 0.21831422E-03 0.84802943E-01 ) :: navier-stokes
ComputeChange: NS (ITER=7) (NRM,RELC): ( 0.23089744E-03 0.56023585E-01 ) :: navier-stokes
ComputeChange: NS (ITER=8) (NRM,RELC): ( 0.23960695E-03 0.37021991E-01 ) :: navier-stokes
ComputeChange: NS (ITER=9) (NRM,RELC): ( 0.24554298E-03 0.24470931E-01 ) :: navier-stokes
ComputeChange: NS (ITER=10) (NRM,RELC): ( 0.24954789E-03 0.16178463E-01 ) :: navier-stokes
ComputeChange: NS (ITER=11) (NRM,RELC): ( 0.25223184E-03 0.10697718E-01 ) :: navier-stokes
ComputeChange: NS (ITER=12) (NRM,RELC): ( 0.25402263E-03 0.70747007E-02 ) :: navier-stokes
ComputeChange: NS (ITER=13) (NRM,RELC): ( 0.25726561E-03 0.12685526E-01 ) :: navier-stokes
ComputeChange: NS (ITER=14) (NRM,RELC): ( 0.25753787E-03 0.10577202E-02 ) :: navier-stokes
ComputeChange: NS (ITER=15) (NRM,RELC): ( 0.25755675E-03 0.73293926E-04 ) :: navier-stokes
ComputeChange: NS (ITER=16) (NRM,RELC): ( 0.25755835E-03 0.62131185E-05 ) :: navier-stokes
ComputeChange: SS (ITER=1) (NRM,RELC): ( 0.25755835E-03 2.0000000 ) :: navier-stokes
ComputeChange: NS (ITER=1) (NRM,RELC): ( 394880.79 2.0000000 ) :: sij
ComputeChange: NS (ITER=2) (NRM,RELC): ( 421559.88 0.65354617E-01 ) :: sij
ComputeChange: NS (ITER=3) (NRM,RELC): ( 431568.68 0.23463770E-01 ) :: sij
ComputeChange: NS (ITER=4) (NRM,RELC): ( 309568.48 0.32922437 ) :: sij
ComputeChange: NS (ITER=5) (NRM,RELC): ( 235878.41 0.27020070 ) :: sij
ComputeChange: NS (ITER=6) (NRM,RELC): ( 239879.73 0.16820801E-01 ) :: sij
ComputeChange: SS (ITER=1) (NRM,RELC): ( 239879.73 2.0000000 ) :: sij
ComputeChange: SS (ITER=1) (NRM,RELC): ( 0.0000000 0.0000000 ) :: eigenstresses
getBalance: Set accumulation file to balance_winter.dat
3 0.4927E-01
6 0.4202E+01
9 0.3234E+01
12 0.1999E+01
15 0.2004E+01
18 0.2006E+01
21 0.2006E+01
...
Cutted by me
...
1026 0.1892+139
1029 0.6784+140
1032 0.1688+142
1032 0.1688+142
ERROR:: IterSolve: Failed convergence tolerances.
Re: Bedrock Sliding Boundary Condition
Hi,
It's the iterative method for the linear system of the free surface (Solver 6) which is not converging.
Increasing the verbosity level (Max Output Level = I in the Simulation section; max I=42) can sometimes help to see where the problem comes from.
Have you checked your velocity field (you can set Exec Solver =never to your Free Surface solver to visualise your results of the NS solver)?
Units for the splip coeff in SI are Pa.s.m-1 (Tau=beta u). 0.1 in SI units is very low; you can strat by putting something like 1MPa.a.m-1 (be carefull my units), this will be very high friction (so you will get a solution close to your 0-velocity); then decrease the value.
In general:
- if iteratives methods for the linear system are not converging you can always switch to a direct solver to get a solution.
- Conditionning of the Free Surface solver is usually good and ietrative methods usually converge. I assume here that you have crazy velocities because of too low friction; and that why it is not converging
- The free surface equatiopn is linear if you are not using min/max (Apply dirichlet = False); so you need only 1 Non linear iteration to get the solution (Nonlinear System Max Iterations = 1, and Nonlinear System Relaxation Factor = 1).
Fabien
It's the iterative method for the linear system of the free surface (Solver 6) which is not converging.
Increasing the verbosity level (Max Output Level = I in the Simulation section; max I=42) can sometimes help to see where the problem comes from.
Have you checked your velocity field (you can set Exec Solver =never to your Free Surface solver to visualise your results of the NS solver)?
Units for the splip coeff in SI are Pa.s.m-1 (Tau=beta u). 0.1 in SI units is very low; you can strat by putting something like 1MPa.a.m-1 (be carefull my units), this will be very high friction (so you will get a solution close to your 0-velocity); then decrease the value.
In general:
- if iteratives methods for the linear system are not converging you can always switch to a direct solver to get a solution.
- Conditionning of the Free Surface solver is usually good and ietrative methods usually converge. I assume here that you have crazy velocities because of too low friction; and that why it is not converging
- The free surface equatiopn is linear if you are not using min/max (Apply dirichlet = False); so you need only 1 Non linear iteration to get the solution (Nonlinear System Max Iterations = 1, and Nonlinear System Relaxation Factor = 1).
Fabien
Re: Bedrock Sliding Boundary Condition
Thank you so much for your answers.
We've solved the problem.
In particular I'm working in SI system, but the Slide Coefficient value of 0.1 was for MPa-m-a system. Calculating the appropriate value of coefficients gets a realistic velocity field and the Free Surface solver doesn't fail anymore.
Thank you again
We've solved the problem.
In particular I'm working in SI system, but the Slide Coefficient value of 0.1 was for MPa-m-a system. Calculating the appropriate value of coefficients gets a realistic velocity field and the Free Surface solver doesn't fail anymore.
Thank you again