Implement `icestream` in AD

[AdrienWehrle]

Now that iceslab AD converges, let's apply this setup to icestream!

[AdrienWehrle]

This is now proposed in https://github.com/AdrienWehrle/diuca/pull/20.

This new setup converges too, velocities in the ice stream are of the right order. The only issue is at the glacier front when I apply the hydrostatic pressure through ocean_pressure (which is needed for my final setup)

The system is pretty much blowing up with velocities > 0.1 m.s-1 which is not possible. Independently of the magnitude, one wouldn't expect a max velocity close to the base where the hydrostatic pressure is the highest and in the vicinity of the no slip boundary condition. Expected results would rather be max velocity at the surface at the front of the glacier, and much lower at the base.

Trying to understand what's happening there, I set the same boundary condition on the upstream and downstream sides (u=0.1mh-1) and get results of the right magnitude (but this is not the setup I need).

It hence feels like it's pointing at the pressure BC that I'm not applying properly. Any idea @GiudGiud ? Thanks a lot in advance for your help!

Leaving that downstream side free of any boundary condition results in a much lower velocity (~2x slower) which doesn't make sense since the ocean pressure should stabilize the whole system a lot.

So I thought I was applying it in the wrong direction, but reversing it still results in high velocities but this time one element row back from the front, and some negative velocities at the very front...

I'm a bit lost I must say, it really looks like I'm doing something wrong with this pressure BC!

The input file I'm working on is here.

[GiudGiud]

Is there another way of setting up the simulation? Imposing the water pressure in the ice is definitely pulling on the ice. Do the results get better with mesh refinement?

[AdrienWehrle]

Is there another way of setting up the simulation?

In which way? You mean different BCs?

Imposing the water pressure in the ice is definitely pulling on the ice.

Yes while it should actually push on it right...

Do the results get better with mesh refinement?

Mesh refinement only allows to show how unstable the system is at the glacier front, velocity keeps changing quite a lot back and forth as the mesh is refined/coarsened but it doesn't look better, just changing a lot through those very high unrealistic values...

[GiudGiud]

Results should converge with mesh refinement. This is a telltale sign that something is not right with the setup.

If you impose the ice pressure instead of the water pressure, do you get smoothly converging results

[AdrienWehrle]

do you get smoothly converging results

If I apply ice weight just like

[ice_weight]
    type = ParsedFunction
    expression = 'if(z > 0, 917 * 9.81 * (100-z), -917 * 9.81 * (z-100))'
  []

then no I don't. Which maybe makes sense since ice weight (917 9.81 ice_column) is a pretty similar pressure to the hydrostatic one (1028 9.81 water_column) since there's only 100m above water, and ~1500m below.

[AdrienWehrle]

Hence, I get very similar results than with ocean pressure

[AdrienWehrle]

Do you think I should go into trying to set BCs at a different time in the transient @GiudGiud? Like adding the ocean pressure at the second timestep and not from the beginning? That potentially sounds like a bad idea but at the same time I don't really see other solutions for now..

[GiudGiud]

No I dont think this will change much.

[ice_weight]
    type = ParsedFunction
    expression = 'if(z > 0, 917 * 9.81 * (100-z), -917 * 9.81 * (z-100))'
  []

this does not seem right. Gravity does not change sign around z = 100. You should have an expression that is valid for the entire z domain, and that has higher pressures at the bottom. Use a FunctionAux to check that the imposed pressure function is behaving as expected

[AdrienWehrle]

Yes sorry it's actually much simpler like this (100 is the surface height in meters)

[ice_weight]
    type = ParsedFunction
    expression = '917 * 9.81 * (100 - z)'
  []

(but the above was also right)

[AdrienWehrle]

What I realize is weird is that the pressure at the surface and front of the glacier (so at the corner between surface and front) is p=0, while there's no reason for it to be.

As if the ocean_pressure or ice_weight I'm applying was not added to the background pressure but directly applied.

Although it should be the pressure resulting from the driving stress (created by the surface slope) minus the ocean pressure. But not the ocean pressure alone.

I think that is why the system is struggling so much at the front.

[AdrienWehrle]

Here is the illustration of what I mean above (p=0 at the front surface):

[AdrienWehrle]

I suppose that means what I need is not a DirichletBC since this imposes the variable value. I'd need a BC that adds up to the solution to create the right stress balance on that boundary.

[AdrienWehrle]

More research now shows I might have to use an ADPressure BC directly. But I remember we discarded it earlier in the development process, yet I can't remember the reason.

[AdrienWehrle]

Yet I can't find any example for ADPressure used with INS, it's used most on the time for mechanics... Which makes me think I should find another way but keep that logic of a Pressure BC instead of a Dirichlet one. Do you agree with this @GiudGiud ?

[GiudGiud]

PressureBC is for solid mechanics. It s the equivalent of a neumann BC in those equations. But for the momentum equation, it is not. It would not work here. Above the surface, did you set a free flow BC?

Also what are the FV results here?

[AdrienWehrle]

With the free flow BC it doesn't converge anymore

  mu
are added for automatic output by MaterialOutputAction.
Finished Setting Up                                                                      [  0.08 s] [  239 MB]
Framework Information:
MOOSE Version:           git commit 778cc52c3e on 2024-03-02
LibMesh Version:         
PETSc Version:           3.20.3
SLEPc Version:           3.20.1
Current Time:            Thu Mar  7 20:50:45 2024
Executable Timestamp:    Sun Mar  3 15:16:50 2024

Parallelism:
  Num Processors:          1
  Num Threads:             15

Mesh: 
  Parallel Type:           replicated
  Mesh Dimension:          3
  Spatial Dimension:       3
  Nodes:                   2244
  Elems:                   1680
  Num Subdomains:          2

Nonlinear System:
  Num DOFs:                8976
  Num Local DOFs:          8976
  Variables:               "velocity" "p" 
  Finite Element Types:    "LAGRANGE_VEC" "LAGRANGE" 
  Approximation Orders:    "FIRST" "FIRST" 

Auxiliary System:
  Num DOFs:                8412
  Num Local DOFs:          8412
  Variables:               { "vel_x" "vel_y" "vel_z" } "mu" 
  Finite Element Types:    "LAGRANGE" "MONOMIAL" 
  Approximation Orders:    "FIRST" "CONSTANT" 

Execution Information:
  Executioner:             Transient
  TimeStepper:             ConstantDT
  TimeIntegrator:          ImplicitEuler
  Solver Mode:             NEWTON
  PETSc Preconditioner:    lu 
  MOOSE Preconditioner:    SMP

Time Step 0, time = 0

Time Step 1, time = 3.1536e+06, dt = 3.1536e+06
 0 Nonlinear |R| = 3.759854e+05
      0 Linear |R| = 3.759854e+05
      1 Linear |R| = 2.130447e-06
 1 Nonlinear |R| = 3.326504e+05
      0 Linear |R| = 3.326504e+05
      1 Linear |R| = 6.828680e-06
 2 Nonlinear |R| = 7.154014e+04
      0 Linear |R| = 7.154014e+04
      1 Linear |R| = 6.178348e-08
 3 Nonlinear |R| = 7.044624e+04
      0 Linear |R| = 7.044624e+04
      1 Linear |R| = 3.729972e-08
 4 Nonlinear |R| = 1.096641e+05
      0 Linear |R| = 1.096641e+05
      1 Linear |R| = 5.609479e-09
 5 Nonlinear |R| = 3.102064e+04
      0 Linear |R| = 3.102064e+04
      1 Linear |R| = 1.941509e-07
 6 Nonlinear |R| = 1.182982e+06
      0 Linear |R| = 1.182982e+06
      1 Linear |R| = 4.285694e-08
 7 Nonlinear |R| = 5.920402e+04
      0 Linear |R| = 5.920402e+04
      1 Linear |R| = 4.929877e-06
 8 Nonlinear |R| = 3.162043e+07
      0 Linear |R| = 3.162043e+07
      1 Linear |R| = 6.450174e-05
 9 Nonlinear |R| = 9.407620e+04
      0 Linear |R| = 9.407620e+04
      1 Linear |R| = 1.254753e-07
10 Nonlinear |R| = 3.850653e+04
      0 Linear |R| = 3.850653e+04
      1 Linear |R| = 4.257988e-09
11 Nonlinear |R| = 2.359778e+04
      0 Linear |R| = 2.359778e+04
      1 Linear |R| = 1.383340e-08
12 Nonlinear |R| = 4.046347e+04
      0 Linear |R| = 4.046347e+04
      1 Linear |R| = 6.988253e-09
13 Nonlinear |R| = 6.686665e+04
      0 Linear |R| = 6.686665e+04
      1 Linear |R| = 5.053045e-09
14 Nonlinear |R| = 4.639239e+04
      0 Linear |R| = 4.639239e+04
      1 Linear |R| = 4.313262e-09
15 Nonlinear |R| = 1.849833e+04
      0 Linear |R| = 1.849833e+04
      1 Linear |R| = 5.011132e-08
16 Nonlinear |R| = 3.509970e+05
      0 Linear |R| = 3.509970e+05
      1 Linear |R| = 2.451396e-08
17 Nonlinear |R| = 6.757336e+04
      0 Linear |R| = 6.757336e+04
      1 Linear |R| = 9.952722e-09
18 Nonlinear |R| = 1.923171e+04
      0 Linear |R| = 1.923171e+04
      1 Linear |R| = 3.444923e-08
19 Nonlinear |R| = 5.629754e+05
      0 Linear |R| = 5.629754e+05
      1 Linear |R| = 2.101116e-08
20 Nonlinear |R| = 6.738814e+04
      0 Linear |R| = 6.738814e+04
      1 Linear |R| = 3.834287e-08
21 Nonlinear |R| = 1.572175e+04
      0 Linear |R| = 1.572175e+04
      1 Linear |R| = 3.032296e-08
22 Nonlinear |R| = 7.524961e+05
      0 Linear |R| = 7.524961e+05
      1 Linear |R| = 4.202473e-08
23 Nonlinear |R| = 6.856843e+04
      0 Linear |R| = 6.856843e+04
      1 Linear |R| = 1.242913e-08
24 Nonlinear |R| = 1.780686e+04
      0 Linear |R| = 1.780686e+04
      1 Linear |R| = 5.244571e-08
25 Nonlinear |R| = 8.106327e+05

but it does if I apply the pressure BC in the wrong direction (pulling) and shows a much better pattern, see screenshot:

But it still has the problem that I shouldn't use a Dirichlet BC, since the pressure is still 0 at the edge

[AdrienWehrle]

Also what are the FV results here?

I will share those in a minute!

[AdrienWehrle]

Also what are the FV results here?

icestream FV doesn't converge either for now. The input file i used is this one.

Framework Information:
MOOSE Version:           git commit 778cc52c3e on 2024-03-02
LibMesh Version:         
PETSc Version:           3.20.3
SLEPc Version:           3.20.1
Current Time:            Thu Mar  7 21:17:56 2024
Executable Timestamp:    Thu Mar  7 21:07:16 2024

Parallelism:
  Num Processors:          1
  Num Threads:             15

Mesh: 
  Parallel Type:           replicated
  Mesh Dimension:          3
  Spatial Dimension:       3
  Nodes:                   2244
  Elems:                   1680
  Num Subdomains:          2

Nonlinear System:
  Num DOFs:                6720
  Num Local DOFs:          6720
  Variables:               { "vel_x" "vel_y" "vel_z" "pressure" } 
  Finite Element Types:    "MONOMIAL" 
  Approximation Orders:    "CONSTANT" 

Execution Information:
  Executioner:             Transient
  TimeStepper:             ConstantDT
  TimeIntegrator:          ImplicitEuler
  Solver Mode:             Preconditioned JFNK
  PETSc Preconditioner:    lu 

Time Step 0, time = 0

Time Step 1, time = 3.1536e+06, dt = 3.1536e+06
Currently Executing
    Computing Initial Residual
      Computing Residual
        Finished Computing Rhie-Chow coefficients                                        [  0.04 s] [  388 MB]
      Finished Computing Residual                                                        [  0.13 s] [  442 MB]
    Finished Computing Initial Residual                                                  [  0.13 s] [  442 MB]
    |residual|_2 of individual variables:
                   vel_x:    65797.3
                   vel_y:    6.44697e-10
                   vel_z:    421052
                   pressure: 49077
 0 Nonlinear |R| = 4.289785e+05
      0 Linear |R| = 4.289785e+05
      1 Linear |R| = 4.289001e+05
      2 Linear |R| = 4.287488e+05
      3 Linear |R| = 4.280157e+05
      4 Linear |R| = 4.278638e+05
      5 Linear |R| = 4.277476e+05
      6 Linear |R| = 4.274709e+05
      7 Linear |R| = 4.268316e+05
      8 Linear |R| = 4.258378e+05
      9 Linear |R| = 4.257111e+05
     10 Linear |R| = 4.256695e+05
     11 Linear |R| = 4.254529e+05
     12 Linear |R| = 4.248987e+05
     13 Linear |R| = 4.238929e+05
     14 Linear |R| = 4.237199e+05
     15 Linear |R| = 4.236885e+05
     16 Linear |R| = 4.236827e+05
     17 Linear |R| = 4.236827e+05
     18 Linear |R| = 4.236768e+05
     19 Linear |R| = 4.236610e+05
     20 Linear |R| = 4.236331e+05
     21 Linear |R| = 4.221271e+05
     22 Linear |R| = 4.188363e+05
     23 Linear |R| = 4.158506e+05
     24 Linear |R| = 4.125485e+05
     25 Linear |R| = 4.093579e+05
     26 Linear |R| = 4.062735e+05
     27 Linear |R| = 4.032371e+05
     28 Linear |R| = 4.002689e+05
     29 Linear |R| = 3.973503e+05
  Linear solve did not converge due to DIVERGED_BREAKDOWN iterations 30
    |residual|_2 of individual variables:
                   vel_x:    4.56126e+24
                   vel_y:    1.12701e+24
                   vel_z:    2.81028e+24
                   pressure: 7.53948e+20
 1 Nonlinear |R| = 5.474754e+24
      0 Linear |R| = 5.474754e+24
      1 Linear |R| = 5.344840e+24
      2 Linear |R| = 5.342711e+24
      3 Linear |R| = 5.341671e+24
      4 Linear |R| = 5.048252e+24
      5 Linear |R| = 5.044845e+24
      6 Linear |R| = 5.040589e+24
      7 Linear |R| = 5.040498e+24
      8 Linear |R| = 5.040356e+24
      9 Linear |R| = 5.039270e+24
     10 Linear |R| = 5.038068e+24
     11 Linear |R| = 5.038010e+24
     12 Linear |R| = 5.037988e+24
     13 Linear |R| = 5.037597e+24
     14 Linear |R| = 5.037595e+24
     15 Linear |R| = 5.037576e+24
     16 Linear |R| = 5.037396e+24
     17 Linear |R| = 5.036935e+24
     18 Linear |R| = 5.036525e+24
     19 Linear |R| = 5.035927e+24
     20 Linear |R| = 5.035927e+24
     21 Linear |R| = 5.035885e+24
     22 Linear |R| = 5.035776e+24
     23 Linear |R| = 5.035588e+24
     24 Linear |R| = 5.035585e+24
     25 Linear |R| = 5.035573e+24
     26 Linear |R| = 5.035525e+24
     27 Linear |R| = 5.034845e+24
     28 Linear |R| = 5.034815e+24
     29 Linear |R| = 5.034708e+24
     30 Linear |R| = 5.508211e+24
     31 Linear |R| = 5.506706e+24
     32 Linear |R| = 5.483555e+24
     33 Linear |R| = 5.426422e+24
     34 Linear |R| = 5.426125e+24
     35 Linear |R| = 5.424704e+24
     36 Linear |R| = 5.412411e+24
     37 Linear |R| = 5.412299e+24
     38 Linear |R| = 5.405521e+24
     39 Linear |R| = 5.401950e+24
     40 Linear |R| = 5.400867e+24
     41 Linear |R| = 5.339503e+24
     42 Linear |R| = 5.339116e+24
     43 Linear |R| = 5.338785e+24
     44 Linear |R| = 5.338544e+24
     45 Linear |R| = 5.338498e+24
     46 Linear |R| = 5.338440e+24
     47 Linear |R| = 5.335158e+24
     48 Linear |R| = 5.334539e+24
     49 Linear |R| = 5.334335e+24
     50 Linear |R| = 5.332619e+24
     51 Linear |R| = 5.332594e+24
     52 Linear |R| = 5.332581e+24
     53 Linear |R| = 5.332560e+24
     54 Linear |R| = 5.332488e+24
     55 Linear |R| = 5.332448e+24
     56 Linear |R| = 5.331607e+24
     57 Linear |R| = 5.331607e+24
     58 Linear |R| = 5.331604e+24
     59 Linear |R| = 5.331184e+24
     60 Linear |R| = 5.608439e+24
     61 Linear |R| = 5.608033e+24
     62 Linear |R| = 5.580976e+24
     63 Linear |R| = 5.549287e+24
     64 Linear |R| = 5.529411e+24
     65 Linear |R| = 5.505868e+24
     66 Linear |R| = 5.494408e+24
     67 Linear |R| = 5.441813e+24
     68 Linear |R| = 5.441808e+24
^C

[AdrienWehrle]

But I think the Dirichlet Pressure BC explains the big difference between FE and FV on the iceslab problem, where the pressure at the front is applied properly in FV but not in FE.

PressureBC is for solid mechanics. It s the equivalent of a neumann BC in those equations. But for the momentum equation, it is not. It would not work here.

Do you think there's an object I could write to make it work? It seems like that's the main remaining issue we have to get icestream AD to converge.

[GiudGiud]

At the outlet, pressure should be fine. We should first focus on making that work

But I think the Dirichlet Pressure BC explains the big difference between FE and FV on the iceslab problem, where the pressure at the front is applied properly in FV but not in FE.

the pressure should be applied property in both FV and FE. Inlet velocity, bottom no-slip, top free-slip, outlet pressure

[GiudGiud]

did you use the MeshDiagnosticsGenerator on this new mesh?

[AdrienWehrle]

At the outlet, pressure should be fine.

But a Dirichlet BC on Pressure will constrain the pressure on the side to the prescribed function and not add it to the solution, no?

the pressure should be applied property in both FV and FE. Inlet velocity, bottom no-slip, top free-slip, outlet pressure

To me that's the case at the moment in both FV and FE.

[AdrienWehrle]

did you use the MeshDiagnosticsGenerator on this new mesh?

No, but I just did now:

Sideset left (2) is consistently oriented with regards to the blocks it neighbors

*** Info ***
Sideset left (2) does not appear to have side-to-neighbor-side orientation flips. All neighbor sides normal differ by less than pi/2

*** Info ***
Sideset right (3) is consistently oriented with regards to the blocks it neighbors

*** Info ***
Sideset right (3) does not appear to have side-to-neighbor-side orientation flips. All neighbor sides normal differ by less than pi/2

*** Info ***
Sideset upstream (4) is consistently oriented with regards to the blocks it neighbors

*** Info ***
Sideset upstream (4) does not appear to have side-to-neighbor-side orientation flips. All neighbor sides normal differ by less than pi/2

*** Info ***
Sideset downstream (5) is consistently oriented with regards to the blocks it neighbors

*** Info ***
Sideset downstream (5) does not appear to have side-to-neighbor-side orientation flips. All neighbor sides normal differ by less than pi/2

*** Info ***
Sideset bottom (6) is consistently oriented with regards to the blocks it neighbors

*** Info ***
Sideset bottom (6) does not appear to have side-to-neighbor-side orientation flips. All neighbor sides normal differ by less than pi/2

*** Info ***
Sideset surface (7) is consistently oriented with regards to the blocks it neighbors

*** Info ***
Sideset surface (7) does not appear to have side-to-neighbor-side orientation flips. All neighbor sides normal differ by less than pi/2

*** Info ***
Sideset sediment (11) is consistently oriented with regards to the blocks it neighbors

*** Info ***
Sideset sediment (11) does not appear to have side-to-neighbor-side orientation flips. All neighbor sides normal differ by less than pi/2

*** Info ***
Sideset left_sediment (12) is consistently oriented with regards to the blocks it neighbors

*** Info ***
Sideset left_sediment (12) does not appear to have side-to-neighbor-side orientation flips. All neighbor sides normal differ by less than pi/2

*** Info ***
Sideset right_sediment (13) is consistently oriented with regards to the blocks it neighbors

*** Info ***
Sideset right_sediment (13) does not appear to have side-to-neighbor-side orientation flips. All neighbor sides normal differ by less than pi/2

*** Info ***
Sideset upstream_sediment (14) is consistently oriented with regards to the blocks it neighbors

*** Info ***
Sideset upstream_sediment (14) does not appear to have side-to-neighbor-side orientation flips. All neighbor sides normal differ by less than pi/2

*** Info ***
Sideset downstream_sediment (15) is consistently oriented with regards to the blocks it neighbors

*** Info ***
Sideset downstream_sediment (15) does not appear to have side-to-neighbor-side orientation flips. All neighbor sides normal differ by less than pi/2

*** Info ***
Number of elements below prescribed minimum volume : 0

*** Info ***
Number of elements with negative volume : 0

*** Info ***
Number of elements above prescribed maximum volume : 0
Element type in subdomain eleblock1 (1) : HEX8
Element type in subdomain eleblock2 (2) : HEX8
Element type in subdomain eleblock3 (3) : HEX8

*** Info ***
Number of elements overlapping (node-based heuristics): 0

*** Info ***
Number of elements overlapping (centroid-based heuristics): 0
Nonplanar side detected for side 5 of element :  Elem Information
   id()=0, unique_id()=2618, subdomain_id()=1, processor_id()=0
   type()=HEX8
   dim()=3
   n_nodes()=8
   mapping=LAGRANGE_MAP
    0  Node id()=0, processor_id()=0, Point=(x,y,z)=(       0,        0, -168.385)
    DoFs=
    1  Node id()=1, processor_id()=0, Point=(x,y,z)=( 933.333,        0, -168.385)
    DoFs=
    2  Node id()=2, processor_id()=0, Point=(x,y,z)=( 933.333,      625,  -151.75)
    DoFs=
    3  Node id()=3, processor_id()=0, Point=(x,y,z)=(       0,      625,  -151.75)
    DoFs=
    4  Node id()=4, processor_id()=0, Point=(x,y,z)=(       0,        0,  -48.068)
    DoFs=
    5  Node id()=5, processor_id()=0, Point=(x,y,z)=( 933.333,        0, -51.2413)
    DoFs=
    6  Node id()=6, processor_id()=0, Point=(x,y,z)=( 933.333,      625, -37.9333)
    DoFs=
    7  Node id()=7, processor_id()=0, Point=(x,y,z)=(       0,      625, -34.7599)
    DoFs=
   n_sides()=6
    neighbor(0)=1680
    neighbor(1)=nullptr
    neighbor(2)=1
    neighbor(3)=21
    neighbor(4)=nullptr
    neighbor(5)=336
   hmin()=113.817, hmax()=1131.19
   volume()=6.8289e+07
   active()=1, ancestor()=0, subactive()=0, has_children()=0
   parent()=nullptr
   level()=0, p_level()=0
   refinement_flag()=DO_NOTHING
   p_refinement_flag()=DO_NOTHING
   DoFs=
Nonplanar side detected for side 5 of element :  Elem Information
   id()=1, unique_id()=2619, subdomain_id()=1, processor_id()=0
   type()=HEX8
   dim()=3
   n_nodes()=8
   mapping=LAGRANGE_MAP
    0  Node id()=1, processor_id()=0, Point=(x,y,z)=( 933.333,        0, -168.385)
    DoFs=
    1  Node id()=8, processor_id()=0, Point=(x,y,z)=( 1866.67,        0, -168.385)
    DoFs=
    2  Node id()=9, processor_id()=0, Point=(x,y,z)=( 1866.67,      625,  -151.75)
    DoFs=
    3  Node id()=2, processor_id()=0, Point=(x,y,z)=( 933.333,      625,  -151.75)
    DoFs=
    4  Node id()=5, processor_id()=0, Point=(x,y,z)=( 933.333,        0, -51.2413)
    DoFs=
    5  Node id()=10, processor_id()=0, Point=(x,y,z)=( 1866.67,        0, -54.4146)
    DoFs=
    6  Node id()=11, processor_id()=0, Point=(x,y,z)=( 1866.67,      625, -41.1066)
    DoFs=
    7  Node id()=6, processor_id()=0, Point=(x,y,z)=( 933.333,      625, -37.9333)
    DoFs=
   n_sides()=6
    neighbor(0)=1681
    neighbor(1)=nullptr
    neighbor(2)=2
    neighbor(3)=22
    neighbor(4)=0
    neighbor(5)=337
   hmin()=110.643, hmax()=1130.82
   volume()=6.64379e+07
   active()=1, ancestor()=0, subactive()=0, has_children()=0
   parent()=nullptr
   level()=0, p_level()=0
   refinement_flag()=DO_NOTHING
   p_refinement_flag()=DO_NOTHING
   DoFs=
Nonplanar side detected for side 5 of element :  Elem Information
   id()=2, unique_id()=2620, subdomain_id()=1, processor_id()=0
   type()=HEX8
   dim()=3
   n_nodes()=8
   mapping=LAGRANGE_MAP
    0  Node id()=8, processor_id()=0, Point=(x,y,z)=( 1866.67,        0, -168.385)
    DoFs=
    1  Node id()=12, processor_id()=0, Point=(x,y,z)=(    2800,        0, -168.385)
    DoFs=
    2  Node id()=13, processor_id()=0, Point=(x,y,z)=(    2800,      625,  -151.75)
    DoFs=
    3  Node id()=9, processor_id()=0, Point=(x,y,z)=( 1866.67,      625,  -151.75)
    DoFs=
    4  Node id()=10, processor_id()=0, Point=(x,y,z)=( 1866.67,        0, -54.4146)
    DoFs=
    5  Node id()=14, processor_id()=0, Point=(x,y,z)=(    2800,        0,  -57.588)
    DoFs=
    6  Node id()=15, processor_id()=0, Point=(x,y,z)=(    2800,      625, -44.2799)
    DoFs=
    7  Node id()=11, processor_id()=0, Point=(x,y,z)=( 1866.67,      625, -41.1066)
    DoFs=
   n_sides()=6
    neighbor(0)=1682
    neighbor(1)=nullptr
    neighbor(2)=3
    neighbor(3)=23
    neighbor(4)=1
    neighbor(5)=338
   hmin()=107.47, hmax()=1130.46
   volume()=6.45868e+07
   active()=1, ancestor()=0, subactive()=0, has_children()=0
   parent()=nullptr
   level()=0, p_level()=0
   refinement_flag()=DO_NOTHING
   p_refinement_flag()=DO_NOTHING
   DoFs=
Nonplanar side detected for side 5 of element :  Elem Information
   id()=3, unique_id()=2621, subdomain_id()=1, processor_id()=0
   type()=HEX8
   dim()=3
   n_nodes()=8
   mapping=LAGRANGE_MAP
    0  Node id()=12, processor_id()=0, Point=(x,y,z)=(    2800,        0, -168.385)
    DoFs=
    1  Node id()=16, processor_id()=0, Point=(x,y,z)=( 3733.33,        0, -168.385)
    DoFs=
    2  Node id()=17, processor_id()=0, Point=(x,y,z)=( 3733.33,      625,  -151.75)
    DoFs=
    3  Node id()=13, processor_id()=0, Point=(x,y,z)=(    2800,      625,  -151.75)
    DoFs=
    4  Node id()=14, processor_id()=0, Point=(x,y,z)=(    2800,        0,  -57.588)
    DoFs=
    5  Node id()=18, processor_id()=0, Point=(x,y,z)=( 3733.33,        0, -60.7613)
    DoFs=
    6  Node id()=19, processor_id()=0, Point=(x,y,z)=( 3733.33,      625, -47.4533)
    DoFs=
    7  Node id()=15, processor_id()=0, Point=(x,y,z)=(    2800,      625, -44.2799)
    DoFs=
   n_sides()=6
    neighbor(0)=1683
    neighbor(1)=nullptr
    neighbor(2)=4
    neighbor(3)=24
    neighbor(4)=2
    neighbor(5)=339
   hmin()=104.297, hmax()=1130.11
   volume()=6.27356e+07
   active()=1, ancestor()=0, subactive()=0, has_children()=0
   parent()=nullptr
   level()=0, p_level()=0
   refinement_flag()=DO_NOTHING
   p_refinement_flag()=DO_NOTHING
   DoFs=
Nonplanar side detected for side 5 of element :  Elem Information
   id()=4, unique_id()=2622, subdomain_id()=1, processor_id()=0
   type()=HEX8
   dim()=3
   n_nodes()=8
   mapping=LAGRANGE_MAP
    0  Node id()=16, processor_id()=0, Point=(x,y,z)=( 3733.33,        0, -168.385)
    DoFs=
    1  Node id()=20, processor_id()=0, Point=(x,y,z)=( 4666.67,        0, -168.385)
    DoFs=
    2  Node id()=21, processor_id()=0, Point=(x,y,z)=( 4666.67,      625,  -151.75)
    DoFs=
    3  Node id()=17, processor_id()=0, Point=(x,y,z)=( 3733.33,      625,  -151.75)
    DoFs=
    4  Node id()=18, processor_id()=0, Point=(x,y,z)=( 3733.33,        0, -60.7613)
    DoFs=
    5  Node id()=22, processor_id()=0, Point=(x,y,z)=( 4666.67,        0, -63.9346)
    DoFs=
    6  Node id()=23, processor_id()=0, Point=(x,y,z)=( 4666.67,      625, -50.6266)
    DoFs=
    7  Node id()=19, processor_id()=0, Point=(x,y,z)=( 3733.33,      625, -47.4533)
    DoFs=
   n_sides()=6
    neighbor(0)=1684
    neighbor(1)=nullptr
    neighbor(2)=5
    neighbor(3)=25
    neighbor(4)=3
    neighbor(5)=340
   hmin()=101.123, hmax()=1129.76
   volume()=6.08845e+07
   active()=1, ancestor()=0, subactive()=0, has_children()=0
   parent()=nullptr
   level()=0, p_level()=0
   refinement_flag()=DO_NOTHING
   p_refinement_flag()=DO_NOTHING
   DoFs=
Nonplanar side detected for side 5 of element :  Elem Information
   id()=8, unique_id()=2626, subdomain_id()=1, processor_id()=0
   type()=HEX8
   dim()=3
   n_nodes()=8
   mapping=LAGRANGE_MAP
    0  Node id()=32, processor_id()=0, Point=(x,y,z)=( 7466.67,        0, -168.385)
    DoFs=
    1  Node id()=36, processor_id()=0, Point=(x,y,z)=(    8400,        0, -168.385)
    DoFs=
    2  Node id()=37, processor_id()=0, Point=(x,y,z)=(    8400,      625,  -151.75)
    DoFs=
    3  Node id()=33, processor_id()=0, Point=(x,y,z)=( 7466.67,      625,  -151.75)
    DoFs=
    4  Node id()=34, processor_id()=0, Point=(x,y,z)=( 7466.67,        0, -73.4546)
    DoFs=
    5  Node id()=38, processor_id()=0, Point=(x,y,z)=(    8400,        0,  -76.628)
    DoFs=
    6  Node id()=39, processor_id()=0, Point=(x,y,z)=(    8400,      625, -63.3199)
    DoFs=
    7  Node id()=35, processor_id()=0, Point=(x,y,z)=( 7466.67,      625, -60.1466)
    DoFs=
   n_sides()=6
    neighbor(0)=1688
    neighbor(1)=nullptr
    neighbor(2)=9
    neighbor(3)=29
    neighbor(4)=7
    neighbor(5)=344
   hmin()=88.43, hmax()=1128.47
   volume()=5.34801e+07
   active()=1, ancestor()=0, subactive()=0, has_children()=0
   parent()=nullptr
   level()=0, p_level()=0
   refinement_flag()=DO_NOTHING
   p_refinement_flag()=DO_NOTHING
   DoFs=
Nonplanar side detected for side 5 of element :  Elem Information
   id()=10, unique_id()=2628, subdomain_id()=1, processor_id()=0
   type()=HEX8
   dim()=3
   n_nodes()=8
   mapping=LAGRANGE_MAP
    0  Node id()=40, processor_id()=0, Point=(x,y,z)=( 9333.33,        0, -168.385)
    DoFs=
    1  Node id()=44, processor_id()=0, Point=(x,y,z)=( 10266.7,        0, -168.385)
    DoFs=
    2  Node id()=45, processor_id()=0, Point=(x,y,z)=( 10266.7,      625,  -151.75)
    DoFs=
    3  Node id()=41, processor_id()=0, Point=(x,y,z)=( 9333.33,      625,  -151.75)
    DoFs=
    4  Node id()=42, processor_id()=0, Point=(x,y,z)=( 9333.33,        0, -79.8013)
    DoFs=
    5  Node id()=46, processor_id()=0, Point=(x,y,z)=( 10266.7,        0, -82.9746)
    DoFs=
    6  Node id()=47, processor_id()=0, Point=(x,y,z)=( 10266.7,      625, -69.6666)
    DoFs=
    7  Node id()=43, processor_id()=0, Point=(x,y,z)=( 9333.33,      625, -66.4933)
    DoFs=
   n_sides()=6
    neighbor(0)=1690
    neighbor(1)=nullptr
    neighbor(2)=11
    neighbor(3)=31
    neighbor(4)=9
    neighbor(5)=346
   hmin()=82.0833, hmax()=1127.88
   volume()=4.97779e+07
   active()=1, ancestor()=0, subactive()=0, has_children()=0
   parent()=nullptr
   level()=0, p_level()=0
   refinement_flag()=DO_NOTHING
   p_refinement_flag()=DO_NOTHING
   DoFs=
Nonplanar side detected for side 5 of element :  Elem Information
   id()=11, unique_id()=2629, subdomain_id()=1, processor_id()=0
   type()=HEX8
   dim()=3
   n_nodes()=8
   mapping=LAGRANGE_MAP
    0  Node id()=44, processor_id()=0, Point=(x,y,z)=( 10266.7,        0, -168.385)
    DoFs=
    1  Node id()=48, processor_id()=0, Point=(x,y,z)=(   11200,        0, -168.385)
    DoFs=
    2  Node id()=49, processor_id()=0, Point=(x,y,z)=(   11200,      625,  -151.75)
    DoFs=
    3  Node id()=45, processor_id()=0, Point=(x,y,z)=( 10266.7,      625,  -151.75)
    DoFs=
    4  Node id()=46, processor_id()=0, Point=(x,y,z)=( 10266.7,        0, -82.9746)
    DoFs=
    5  Node id()=50, processor_id()=0, Point=(x,y,z)=(   11200,        0,  -86.148)
    DoFs=
    6  Node id()=51, processor_id()=0, Point=(x,y,z)=(   11200,      625, -72.8399)
    DoFs=
    7  Node id()=47, processor_id()=0, Point=(x,y,z)=( 10266.7,      625, -69.6666)
    DoFs=
   n_sides()=6
    neighbor(0)=1691
    neighbor(1)=nullptr
    neighbor(2)=12
    neighbor(3)=32
    neighbor(4)=10
    neighbor(5)=347
   hmin()=78.91, hmax()=1127.6
   volume()=4.79268e+07
   active()=1, ancestor()=0, subactive()=0, has_children()=0
   parent()=nullptr
   level()=0, p_level()=0
   refinement_flag()=DO_NOTHING
   p_refinement_flag()=DO_NOTHING
   DoFs=
Nonplanar side detected for side 5 of element :  Elem Information
   id()=12, unique_id()=2630, subdomain_id()=1, processor_id()=0
   type()=HEX8
   dim()=3
   n_nodes()=8
   mapping=LAGRANGE_MAP
    0  Node id()=48, processor_id()=0, Point=(x,y,z)=(   11200,        0, -168.385)
    DoFs=
    1  Node id()=52, processor_id()=0, Point=(x,y,z)=( 12133.3,        0, -168.385)
    DoFs=
    2  Node id()=53, processor_id()=0, Point=(x,y,z)=( 12133.3,      625,  -151.75)
    DoFs=
    3  Node id()=49, processor_id()=0, Point=(x,y,z)=(   11200,      625,  -151.75)
    DoFs=
    4  Node id()=50, processor_id()=0, Point=(x,y,z)=(   11200,        0,  -86.148)
    DoFs=
    5  Node id()=54, processor_id()=0, Point=(x,y,z)=( 12133.3,        0, -89.3213)
    DoFs=
    6  Node id()=55, processor_id()=0, Point=(x,y,z)=( 12133.3,      625, -76.0133)
    DoFs=
    7  Node id()=51, processor_id()=0, Point=(x,y,z)=(   11200,      625, -72.8399)
    DoFs=
   n_sides()=6
    neighbor(0)=1692
    neighbor(1)=nullptr
    neighbor(2)=13
    neighbor(3)=33
    neighbor(4)=11
    neighbor(5)=348
   hmin()=75.7366, hmax()=1127.33
   volume()=4.60756e+07
   active()=1, ancestor()=0, subactive()=0, has_children()=0
   parent()=nullptr
   level()=0, p_level()=0
   refinement_flag()=DO_NOTHING
   p_refinement_flag()=DO_NOTHING
   DoFs=
Maximum output reached, log is silenced

*** Info ***
Number of non-planar element sides detected: 2554

*** Info ***
Number of non-conformal nodes: 0

*** Info ***
Number of non-conformal nodes likely due to mesh refinement detected by heuristic: 0

*** Info ***
Number of elements with a negative Jacobian: 0

*** Info ***
Number of element sides with negative Jacobians: 0

[AdrienWehrle]

PressureBC is for solid mechanics. It s the equivalent of a neumann BC in those equations. But for the momentum equation, it is not. It would not work here.

@GiudGiud , does this mean I should try and write a new object based on the PressureBC of solid mechanics but for INSAD? Again, I think what I'm doing here with the ADFunctionDirichletBC on pressure is to prescribe the very pressure value at the front while what I want is to add the ocean pressure to the solution. I don't want p=0 a the surface at the front, but just to not add the hydrostatic pressure there.

[GiudGiud]

No I don’t think so. We should be able to work with a pressure boundary