Dear all,
I would need to write a boundary velocity Function using the Chunk geometry. I used r, x and y as r, phi, theta variable names, and used the colatitude in the coordinates within the function (you can see the function below…). Also changed the component to r, phi theta.
However, the velocity field is not as expected in my domain.
Surely is there something that I am doing wrong. Do you have any clues?
Thanks a lot!
Eleonora
subsection Boundary velocity model
set Tangential velocity boundary indicators = west, east
set Prescribed velocity boundary indicators = top: function, north: zero velocity, south: zero velocity
subsection Function
set Coordinate system = spherical
set Use spherical unit vectors = true
set Variable names = z, x, y
set Function constants = cm=0.01, year=1
set Function expression = 0; if (x<=-16.85 && y>=90.0 | x<=-25.41 && y<=90.0,1.0cm/year,-1.0cm/year); 0
end
end
@nora A good starting point for debugging problems is to be specific about what the problem is, so that one can speculate about solutions. All you say is “Surely is there something that I am doing wrong. Do you have any clues?”, but we really don’t: We don’t even know that what you are doing is wrong because we don’t know what it is you get (and we can’t try ourselves because we don’t know what your complete input file is).
So here are some questions:
What happens when you try with your function?
Does the program provide an error?
If not, does it run to completion?
If it does, how does the solution look like?
In which specific ways does what you see differ from what you expect?
@Nora - As you debug the issue following @bangerth’s suggestions above, it may be helpful to take a look at and try modifying the following test case where a 2D spherical shell has a velocity function applied on the top boundary - boundary_velocity_function_spherical.prm.
Dear Wolfgang and John,
thanks a lot for your replies. You are absolutely right, I was so unspecific about the problem. That’s because I was sure there was something macroscopically wrong in the function I was writing, but still without saying what I was expecting from my function, sorry for that!
The function I wrote should give a transform fault in the XY plane, in the middle of a 3d chunk domain.
The program runs to completion without giving any error, but the velocity field on the top gives a flow in a unique direction, instead of giving a transform fault (see attached figure).
@jbnaliboff thanks for this suggestion, I will try implementing modifications to this .prm
@Nora Separately, looking at your input file, here are a couple of comments:
It is confusing to see the spherical coordinates named x,y,z. (Or, in fact, z,x,y.) Use the correct names that include r,phi,theta in an order you should verify. It is very difficult to debug problems when variable names suggest one thing, but represent another thing.
You use | as the or operator. Should this be ||?
I think that what the expression represents is something where you assume that or binds lower than and (&&). That’s perhaps true, but I’m always confused about this, and so are probably 50% of all programmers. It’s worth using parentheses.
I changed the variable names to make the function clearer and used parentheses as below, also tried with both coordinates (with theta in colatitude) expressed in radians (commented option) and degrees (non commented option), but the result is always the velocity field attached in the previous email, so a unidirectional velocity field. It seems that the function is not using the “else” option.
subsection Function
set Coordinate system = spherical
set Use spherical unit vectors = true
set Variable names = r, phi, theta
set Function constants = cm=0.01, year=1
set Function expression = 0; if (((phi<=-0.44348816 && theta>=1.5708) || (phi<=-0.29408798 && theta<=1.5708)), 1.451cm/year, -1.451cm/year); 0
set Function expression = 0; if ((phi<=-25.41 && theta>=90.0) && (phi<=-16.85 && theta<=90.0), 1.451cm/year, -1.451cm/year); 0
end
phi ranges between 0 and 2*pi (so you should add 2pi to your negative phi values), and both theta and phi are in radians. Sorry for missing this at the beginning, it is also something that I had to look up.
