Boundary traction model - initial lithostatic pressure implementation

Hi,
We are trying to implement open side boundary conditions in my 2D subduction model, but I am currently testing implementing them in a very simple 2D sinker model with this setup:

subsection Boundary velocity model
  set Tangential velocity boundary indicators = bottom, top
end

subsection Boundary traction model
  set Prescribed traction boundary indicators = left:initial lithostatic pressure, right:initial lithostatic pressure
  subsection Initial lithostatic pressure
    set Representative point= 0,250e3
  end
end

I ran a comparison test of a 2D box that has free-slip conditions on all four boundaries (sinker_closed.prm) and one that has free-slip top and bottom boundaries, but open left and right boundaries (sinker_open.prm). The outputs of running sinker_open.prm shows that the sinker drops down to the left, rather than sinking straight down as it does in the closed model. I am wondering if I am not understanding the open boundaries correctly or the implementation is missing something/has something wrong causing the sinker to fall to the side rather than straight down. I have included the two very simple prm files below. The images show three timestep comparisons of the closed model (left), and the open side boundaries model (right). The black vectors are the unscaled velocity field. There is a horizontal flow across the box from right to left apparent at Time 0 in the open model.



sinker_open.prm (2.6 KB)
sinker_closed.prm (2.4 KB)

Are you removing the net x translation mode? This would be done under the ‘Nullspace removal’ subsection in the .prm file. With free-slip on the top and bottom and stress-free side boundary conditions, I think that the solution is only obtained up to an arbitrary x-translation.

I agree with Max, removing the net x translation will likely fix the issue (assuming this is not already in the .prm file).

From past experience and the results of other groups, even when removing the null space anomalies of sorts can appear when one of the velocity directions is not constrained on any boundary. Such issues were most common when using a free surface (and mesh deformation).

For your subduction models, it may be worth trying to only have open boundaries on one side or setting a very deep no-slip boundary on the base.

Ah thank you Max and John! I was not removing the net x translation. After adding

subsection Nullspace removal
 set Remove nullspace = net x translation
end

The sinker now drops down vertically and the horizontal flow across the box is gone. I will add this to my subduction model trying to have the side open boundaries next. Thank you!

When implementing this in my subduction model which is a 2D spherical shell using the Chunk geometry, does ASPECT convert “net x translation” to spherical coordinates? or is “net rotation” the analogous implementation for models using spherical coordinates?

Hi Becky,

Yes, the removing the net rotation is the equivalent for the spherical models. There are similar options, which involve removing net translation or momentum.

The parameter documentation has the following notes about these options:

For most problems there should not be a significant difference between the momentum and rotation/translation versions of nullspace removal, although the momentum versions are more physically motivated. They are equivalent for constant density simulations, and approximately equivalent when the density variations are small. Note that while more than one operation can be selected it only makes sense to pick one rotational and one translational operation.

So, all should have similar effects but perhaps worth testing a few of them. I have no experience running models with these features on, but perhaps others with experience will chime in.

Let us know how it goes!

FYI, the ASPECT website has a searchable list of all the parameters + documentation, which I find really convenient for quickly looking up the documentation if I don’t have a parameters.prm file handy! - https://aspect.geodynamics.org/doc/parameter_view/parameters.xml

Cheers,
John