# Boundary velocity (flow in at bottom) leads to a failure of convergence

**URL:** <https://community.geodynamics.org/t/boundary-velocity-flow-in-at-bottom-leads-to-a-failure-of-convergence/3627>\
**Category:** ASPECT\
**Tags:** convergence\
**Created:** [September 2, 2024, 5:27pm UTC](https://community.geodynamics.org/t/boundary-velocity-flow-in-at-bottom-leads-to-a-failure-of-convergence/3627 "2024-09-02T17:27:31Z")\
**Posts on this page:** 6\
**Page:** 1

<div class="post-metadata">

**Author:** ![optimux](https://avatars.discourse-cdn.com/v4/letter/o/f14d63/32.png) [@optimux](https://community.geodynamics.org/u/optimux)\
**Post date:** [September 2, 2024, 5:27pm UTC](https://community.geodynamics.org/t/boundary-velocity-flow-in-at-bottom-leads-to-a-failure-of-convergence/3627/1 "2024-09-02T17:27:31Z")

</div>

Hi all,

Hope you’re doing well. I’m pretty new to ASPECT tbh and I’m trying to reproduce plume cases in an early paper of Dr. Dannberg (Compressible magma/mantle dynamics: 3-D, adaptive simulations in ASPECT, GJI 2016). The repo could be accessed at [GJI2016 of Dr. Dannberg](https://github.com/jdannberg/melt-transport-data).

The problem is that for subcases plume\_compressible and plume\_hydrous\_compressible, the following boundary condition would lead to a failure of convergence:

```auto
subsection Boundary velocity model
  set Tangential velocity boundary indicators = left, right, top
  set Prescribed velocity boundary indicators = bottom:function
  set Zero velocity boundary indicators = 

  subsection Function
    set Function constants = b=100000, c=20000
    set Variable names = x, y
    set Function expression = 0.0; -0.024995 + 0.1 * exp(-((x-b)*(x-b)+y*y)/(2*c*c))
  end
end

```

And I believe it is the vertical velocity at the bottom (the math expression) that caused the problem. A full prm file is attached  
[plume\_compressible.prm](https://community.geodynamics.org/uploads/short-url/AeVIjJkXMZtJYFTCAOzgrurGs69.prm) (6.6 KB)  
. This prm is modified accordingly since I’m using 2.6. You can see the log  
[log.txt](https://community.geodynamics.org/uploads/short-url/fLbNuOX8Py9DzBGnixg0fBL0yTK.txt) (1.4 KB)  
and solver history  
[solver\_history.txt](https://community.geodynamics.org/uploads/short-url/d28X1BaG51AUZbD9oxLnKG4gf5k.txt) (15.6 KB).

What I’m thinking is that this y component velocity is a source of mass loss or gain of the system, since there is no other boundary condition for mass balance that I can think of for this case. Surprisingly, in Table 3 of Dr. Dannberg’s paper, `u_s` (solid velocity) is omitted, I interpreted it as (0, 0). But in section 4.6, a description of velocity boundary condition is added:  
_the velocity boundary conditions are free slip everywhere except for the bottom boundary layer, where the hydrostatic pressure is applied, but material is allowed to flow in and out._  
So the last piece says y component of us is actually free. So the Table, the boundary condition in the context and boundary condition in prm are not consistent? Which is right?

When I set `u_s` (0, 0), i.e., in the prm `set Function expression = 0.0; 0.0`, it becomes better but again fails during the first step. The log file  
[log.txt](https://community.geodynamics.org/uploads/short-url/8FbMkCU71PlhdHnXqUuhJY4HPRW.txt) (2.3 KB)  
and solver history file  
[solver\_history.txt](https://community.geodynamics.org/uploads/short-url/hr1J1REmDZwrxb9NmzPkUOqQihk.txt) (11.6 KB)  
are attached.

Another interesting thing is that velocity field for step 0 is not as expected:

 ![Screen Shot 2024-09-02 at 12.10.08](https://canada1.discourse-cdn.com/flex036/uploads/geodynamics/original/2X/c/cdfe8574327d3c54f9fa32c30fbec03942af35d7.jpeg)  
`u_s,x`, `u_s,y`, `u_s,z` at bottom boundary are consistent with what is given in the prm. But the velocity over the domain should be (0, 0), isn’t it? Or is this the result of some initial guess/adjustment?

Any thoughts? Really appreciate your help!

Mingming

---

<div class="post-metadata">

**Author:** ![optimux](https://avatars.discourse-cdn.com/v4/letter/o/f14d63/32.png) [@optimux](https://community.geodynamics.org/u/optimux)\
**Post date:** [September 2, 2024, 5:46pm UTC](https://community.geodynamics.org/t/boundary-velocity-flow-in-at-bottom-leads-to-a-failure-of-convergence/3627/2 "2024-09-02T17:46:08Z")

</div>

PS: smaller initial perturbation like `set Amplitude = 50` plus `u_s = (0, 0)` would help.

---

<div class="post-metadata">

**Author:** ![jbnaliboff](https://avatars.discourse-cdn.com/v4/letter/j/58956e/32.png) [@jbnaliboff](https://community.geodynamics.org/u/jbnaliboff)\
**Post date:** [September 5, 2024, 9:16pm UTC](https://community.geodynamics.org/t/boundary-velocity-flow-in-at-bottom-leads-to-a-failure-of-convergence/3627/3 "2024-09-05T21:16:26Z")

</div>

@optimux - Welcome to the community and thank you for posting your question to the forum! Second, apologies it took awhile for a reply.

The issue is that the Boundary velocity model you prescribed in the updated PRM is not actually what was used in the PRM files from [GJI2016 of Dr. Dannberg](https://github.com/jdannberg/melt-transport-data).

In the “plume” series of PRM files, the prescribed velocities constraints are:

```auto
  set Prescribed velocity boundary indicators =

  set Tangential velocity boundary indicators = 0,1,3
  set Zero velocity boundary indicators =

```

There is a prescribed velocity section elsewhere in the PRM file, but the function from that section is not actually applied to the bottom boundary.

Per the code snippet above and what you noted above from the paper, the velocity boundary conditions are free-slip on the left (0), right (1), and top (3) boundaries. A condition is not applied for the bottom (2) boundary, which means the default boundary setting (pressure) is applied.

Can you try re-running the model after adding the following section to produce a lithostatic pressure condition on the bottom boundary?:

```auto
subsection Boundary traction model
  set Prescribed traction boundary indicators = bottom :initial lithostatic pressure

  subsection Initial lithostatic pressure
    set Number of integration points = 1000
    set Representative point = 2000000,0
  end
end

```

Note - the above snippet was taken and modified from [this cookbook](https://github.com/geodynamics/aspect/tree/main/cookbooks/grain_size_ridge).

Cheers,  
John

---

<div class="post-metadata">

**Author:** ![optimux](https://avatars.discourse-cdn.com/v4/letter/o/f14d63/32.png) [@optimux](https://community.geodynamics.org/u/optimux)\
**Post date:** [September 5, 2024, 11:48pm UTC](https://community.geodynamics.org/t/boundary-velocity-flow-in-at-bottom-leads-to-a-failure-of-convergence/3627/4 "2024-09-05T23:48:34Z")

</div>

Hi John,

Thank you so much for the reply, I was so wrong! The syntax is quite misleading.

The original PRM from Juliane2016 has the following:

```auto
subsection Boundary velocity model
  subsection Function
    set Function constants = b=100000, c=20000
    set Variable names = x,y
    set Function expression = 0.0; -0.024995 + 0.1 * exp(-((x-b)*(x-b)+y*y)/(2*c*c))
  end
end

```

and

```auto
subsection Model settings
  ...
  set Prescribed velocity boundary indicators =

  set Tangential velocity boundary indicators = 0,1,3
  set Zero velocity boundary indicators = 
  ...
end

```

So the first question is the long math expression, what is it describing? y component of solid velocity?

The second question is the role that `subsection Boundary velocity model` plays. To which boundary shall this `Boundary velocity model` be applied? It’s like a body without a head (not a good metaphor though). I was thinking this `Boundary velocity model` would be applied to the bottom boundary because the `Prescribed velocity boundary indicators` is left blank/default, i.e., no 2 up there.

The third question is how do you know `Prescribed velocity boundary indicators` is DEFAULT to pressure. Is it in the source code or cookbook? In manual 2.1.0, I don’t see any clues.

The snippet you provided is helpful, it’s working. Thank you !

Cheers,  
Mingming

---

<div class="post-metadata">

**Author:** ![jbnaliboff](https://avatars.discourse-cdn.com/v4/letter/j/58956e/32.png) [@jbnaliboff](https://community.geodynamics.org/u/jbnaliboff)\
**Post date:** [September 9, 2024, 10:38pm UTC](https://community.geodynamics.org/t/boundary-velocity-flow-in-at-bottom-leads-to-a-failure-of-convergence/3627/5 "2024-09-09T22:38:27Z")

</div>

Hi Mingming,

> So the first question is the long math expression, what is it describing? y component of solid velocity?

Yes, it was describing the y-component of velocity. However, in this case that function was not actually used (i.e., it was just extra code from previous testing).

> The second question is the role that `subsection Boundary velocity model` plays. …

I don’t recall the exact syntax for this earlier version of ASPECT, but if one changed the model setting to have `set Prescribed velocity boundary indicators = 2` and indicated that the velocity would be prescribed through a function, then the aforementioned function would have utilized.

Quite some time ago many of the parameters in subsection Model settings were moved to more relevant subsection to make the PRM files easier to ready. For example, all of the parameters related to velocity boundary indicators under `subsection Model settings` are now specified under `subsection Boundary velocity model`.

> The third question is how do you know `Prescribed velocity boundary indicators` is DEFAULT to pressure. Is it in the source code or cookbook? In manual 2.1.0, I don’t see any clues.

That is a great question. I know this anecdotally, but a quick search in the manual did not reveal anything. This could certainly be sorted out by looking at the source code, but it should be in the manual. [I made an issue](https://github.com/geodynamics/aspect/issues/6031) about this on the github page.

> The snippet you provided is helpful, it’s working. Thank you !

Great, glad that it is working!

Cheers,  
John

---

<div class="post-metadata">

**Author:** ![optimux](https://avatars.discourse-cdn.com/v4/letter/o/f14d63/32.png) [@optimux](https://community.geodynamics.org/u/optimux)\
**Post date:** [September 10, 2024, 2:16am UTC](https://community.geodynamics.org/t/boundary-velocity-flow-in-at-bottom-leads-to-a-failure-of-convergence/3627/6 "2024-09-10T02:16:59Z")

</div>

Hi John,

Thank you so much! BTW, in Section 4.6 of Juliane2016 paper, she explained: _the velocity boundary conditions are free slip everywhere except for the bottom boundary layer, where the hydrostatic pressure is applied, but material is allowed to flow in and out._ But unfortunately in her original repo [Juliane2016\_repo](https://github.com/jdannberg/melt-transport-data/blob/master/plume/plume_compressible.prm), there is no subsection of `Initial lithostatic pressure`. May I ask her to update the repo? @jdannberg

Cheers,  
Mingming