Hi all,
I have a 3D model Cartesian box for which I want to prescribe traction boundary conditions of the kind of pressure on two of the side walls (left and right) and also at the bottom.
I tried using ‘initial lithostatic pressure’, but it didn’t work, it complained about my initial topography (which is a loaded data). We could discuss that eventually.
After getting tired of not understanding, and also thinking that in my model I prefer to know exactly what pressure profile is being forced on the walls, plus, my model is laterally heterogeneous, thus getting a good pressure profile (rho.g.depth) is not so direct; ..after that, I decided to set up for traction boundaries as functions which I can control. So, I am setting a unique expression that specifies traction vector functions over three different walls.
I did:
subsection Boundary traction model
set Prescribed traction boundary indicators = bottom:function , left:function , right:function
subsection Function
set Coordinate system = cartesian
set Variable names = x,y,z
set Function constants = Pbott=8.2e9, dens=3400, g=10, H=2.52e5
set Function expression = if( x<=0, dens*g*(H-z), if(x>=4e5,-dens*g*(H-z),0) ) ; 0 ; if(z<=0,Pbott,0)
set Function expression = if( x<=0, densg(H-z), if(x<4e5,0,-densg(H-z)) ) ; 0 ; if(z<=0,Pbott,0) ###### I also tried this one.
end
end
#########################
I also tried with that second expression (# set Function expression), which is slightly different, perhaps more correct…
For both of those trial functions:
That vector function with if statements will prescribe signed Px(z) on the lateral walls and Pbott at the bottom.
(I am planning to set more complex expressions for Px(z), regarding the different densities of layers, so this was just a starting point).
For both of them the compiler parses ! Aspects runs !
But, … The iterative Stokes solver did not converge 
Among other messages it reads: AztecOO::Iterate error code -2: numerical breakdown
The density I used is the maximum of the system, and H is slightly (2 km) larger than the model box height (cause after learning that trying the original value didn’t work, I decided to add the topography height to it); …So there will be an over pressure on the side walls, a pressure profile that is slightly larger than the theoretical lithostatic pressure…
I know that imbalanced forces on the walls may lead to …global acceleration, and a likely breakdown in the Solver. Here I have exactly opposite function forces on the lateral walls, but they are not exactly equal to the lithostatic pressure function, so I would expect some 3D global deformation adjustment in early times…I also tried increasing the bottom pressure…to no avail.
What do you think I should do?
You can find my prm file, the initial topography file, and the logo.txt log file with the run error message, here:
https://github.com/UfrancoU/longtrain
cheers,
Felipe