search

pohlan / SheetModel.jl

2 stars 1 forks source link

Negative hydraulic potential #3

Closed pohlan closed 3 years ago

pohlan commented 3 years ago

https://github.com/pohlan/SheetModel.jl/blob/fa03e42fbe84dfcedd56ad35c8ef01b5121b083b/examples/run_SHMIP.jl#L1-L4

If it is run like this, the hydraulic potential (phi) is negative in a few lines next to the x=0 boundary, leading to a glacier-upwards flux (in positive x direction) since the boundary condition at x=0 is pw = 0 (hence phi=0 for zb=0).

pohlan commented 3 years ago

It is not the case for other initial conditions in ϕ (here it was ϕ0 = 0 everywhere) but something is still wrong since the final ϕ seems to depend a lot on the initial condition.

luraess commented 3 years ago

https://github.com/pohlan/SheetModel.jl/blob/fa03e42fbe84dfcedd56ad35c8ef01b5121b083b/examples/run_SHMIP.jl#L1-L4

If it is run like this, the hydraulic potential (phi) is negative in a few lines next to the x=0 boundary, leading to a glacier-upwards flux (in positive x direction) since the boundary condition at x=0 is pw = 0 (hence phi=0 for zb=0).

I just ran the "A1" case and it does not produce negative ϕ values. I get ϕ[1,:]=0 but as I understand your comment it is expected.

pohlan commented 3 years ago

I just ran the "A1" case and it does not produce negative ϕ values.

That was because I changed the initial conditions in between, see my previous comment.

pohlan commented 3 years ago

In the end the problem was that the damping parameters and scaling factors above (which are quite extreme), in combination with the definition of phi over phi_previous = phi_current, allowed for super quick convergence (~ 50 iterations), but somehow the solution stayed very close to the initial condition and did not produce the correct result. I'm not sure why though. Now, with more reasonable parameters, the solution is as expected, but the number of iterations is in the order of ~1'000.

https://github.com/pohlan/SheetModel.jl/blob/334a61f6411da19d47c226252cf4dd4ae7a3e976/examples/run_SHMIP.jl#L1-L6

mauro3 commented 3 years ago

Could it be that the updates were so small that the solver thought it converged? The 1000 iterations are to reach steady state?

@luraess: 1000 iterations, is that reasonable?

luraess commented 3 years ago

1000 iterations, is that reasonable?

Depends. For steady state may be fine. Could check if one can make it smaller maybe - order 500 or so may be ok. maybe doing the physical time implicit could help, as discussed.