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