How to create a bottom velocity bc with no horizontal slip and only vertical flow

ASPECT
1

Hi all,
Hope you are doing well. I wonder if it is feasible to implement a bottom velocity boundary condition with no horizontal slip and only vertical flow in a 2D box model. The top boundary has a prescribed velocity boundary condition. The left and right boundaries have open free-slip boundary conditions.
I also tried setting Zero velocity boundary indicators = bottom x and Tangential velocity boundary indicators = bottom y. However, this doesn’t work.
I would greatly appreciate it if someone could help me set a bottom velocity bc with no horizontal slip and only vertical flow. Thank you very much!

Cheers
ziqi

2

Hi Ziqi,

I wonder if it is feasible to implement a bottom velocity boundary condition with no horizontal slip and only vertical flow in a 2D box model.

Yes, you can do this by prescribing a function for the velocity boundary conditions. Below is an example of how this would work, with free slip on the top, left, and right walls and a function for the bottom boundary with vertical inflow and no horizontal flow.

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

  subsection Function
    set Variable names      = x,y
    set Function constants  = vx=0, vy=1
    set Function expression = vx ; vy
  end
end

There are also lots of examples in the ASPECT manual and tests folder. Please have a careful look at all of these.

Cheers,
John

3

Hi John,
Thanks for your timely and helpful response.
The thing is, we want to implement a prescribed velocity boundary at the top to add a plate velocity of 5cm/yr and your method works 100%.
subsection Boundary velocity model
set Prescribed velocity boundary indicators = top:function
subsection Function
set Variable names = x,y
set Function constants = u0=0.05
set Function expression = if(y=0,u0,0);0
end
end
now we also want to add a bottom boundary condition with only vertical inflow and no horizontal slip (vx=0, vy=1) at 1320 km depth, which is a different function from the top velocity boundary condition. I tried the following, but the later function just covers the previous function.subsection Boundary velocity model
set Prescribed velocity boundary indicators = top:function
subsection Function
set Variable names = x,y
set Function constants = u0=0.05,depth=1320e3,width=2640e3
set Function expression = if(y=0,u0,0);0
end

set Prescribed velocity boundary indicators = bottom:function
subsection Function
set Variable names = x,y
set Function constants = u1=0.01, width=2640e3,depth=1320e3
set Function expression = 0;if(y=depth,u1,0)
end
end

Is there any way we can set a prescribed velocity at the top (vx=0.05, vy=0) and a different prescribed velocity at the bottom (vx=0,vy=1)?

Thank you very much for your time!

Cheers
Ziqi

4

Hi Ziqi,

Is there any way we can set a prescribed velocity at the top (vx=0.05, vy=0) and a different prescribed velocity at the bottom (vx=0,vy=1)?

Yes, you just need to include all the relevant conditional statements in your equation and set the prescribed boundary indicators in a single list.

I think the following should work?

set Prescribed velocity boundary indicators = top: function, bottom: function
  subsection Function
    set Variable names = x,y
    set Function constants = vx_top=0.05, vx_bottom=0.0, vy_top=0.0, vy_bottom=1.0
    set Function expression = if (y=0, vx_bottom, vx_top); if (y=0, vy_bottom, vy_top)
  end
end

Note that you can include multiple nested if statements in your equations if you need a more complicated setup.

Cheers,
John

5

Dear John,
Thank you for your response! I thought that your solution makes a lot of sense and should work. However, when I tried your method as follows:

subsection Boundary velocity model
set Prescribed velocity boundary indicators = top: function, bottom: function
  subsection Function
    set Variable names = x,y
    set Function constants = vx_top=0.05, vx_bottom=0.0, vy_top=0.0, vy_bottom=1.0
    set Function expression = if (y=0, vx_bottom, vx_top); if (y=0, vy_bottom, vy_top)
  end
end

It yielded the following results.

velocity_bc_test
velocity_bc_test1188×615 274 KB

However, there’s still no vertical velocity at the bottom boundary. And the code didn’t show any mistakes during the run. I don’t know where went wrong.

Ziqi

6

Hi Ziqi,

I’m also not sure what the issue is. The snippet I posted looks correct at second glance.

Questions:

  • What are the boundary conditions on the left and right walls?
    Are the velocity or traction BCs?
  • If they are traction BCs, what pressure constraints are you applying and do you remove the nullspace?
  • Are the velocities in the plot above consistent with what you would expect? It does not look you have a smooth gradient in the x-velocity from top to bottom, as one would expect?
  • John
7

Hi John,
Thank you for your reply, and sorry for my late response. I didn’t add any velocity or traction boundary conditions to the left or right boundaries. I think it means that the left and right sides are left open, am I right?
You are right. The velocity gradient is not consistent with what I would expect. It’s not a smooth gradient from top to bottom. And I don’t get a vertical velocity at the bottom at all, which is strange.
I attached my prm file, and I would greatly appreciate it if you could help me have a look when you are available. Thank you very much!
test5.prm (4.5 KB)

8

Hi Ziqi,

I won’t have time to look in detail at the PRM, but below are some suggestions for how to proceed:

  1. Modify the function to be a function of the x and y position, such that x and y velocities also vary as a function of depth (linearly?) from top to bottom. Once that is done, also apply the function to the left and right walls. The output velocity should be identical to what is defined by your function.

  2. See this test for an example of how to correctly apply lithostatic boundary conditions and remove the null space.

Cheers,
John

9

Hi John,
Thanks for your response. We look into the model result and found no matter how we change the vx and vy at the top and bottom, the bottom always has the same velocity boundary conditions as the top. So I want to ask how ASPECT define the axis grid of the box. Is y zero at the bottom, and y=‘Y extent’ at the top, right?
Thank you for your time and help!

Cheers
Ziqi

10

Hi John,
No worries. We have figured out why the top and bottom boundaries have the same velocities. It’s because I made a rookie mistake. After I changed if sentence to if (y==0), it worked. Thank you for your help!
Another question: is it possible that we only set the horizontal velocity at the bottom to zero and left the vertical velocity flow freely without prescribing any vertical velocity?
I tried setting up the bottom x velocity boundary indicator to zero but didn’t work.

subsection Boundary velocity model
  set Prescribed velocity boundary indicators   = top:function,left:function
  set Zero velocity boundary indicators         = bottom x
  subsection Function
    set Variable names          = x,y
    set Function constants      = u0=0.101, u1=0.001, y1=660e3, depth=1320e3
    set Function expression     = if(x<0,0,if(y<0,0,if(y<y1,y/y1*u1,u1+(y-y1)/(depth-y1)*(u0-u1))));0
  end
end

I would greatly appreciate your help!

Cheers
Ziqi