Hello, I am currently trying to do the linear stability around a sphere flow at low Reynolds number (following the work of P. Meliga, J.-M.Chomaz and D.Sippa ( Unsteadiness in the wake of disks and spheres: Instability, receptivity and control using direct and adjoint global stability analyses). I have implemented the incompressible linearized Navier-Stokes equations in cylindrical coordinates and provided the base flow based on a simulations from Freefem++. My main question concerns the boundary conditions to be applied. These are the following one based on the domain notations in attached: [domain] [boundary_cond] The way to impose Dirichlet and Neumann boundary conditions seems relatively clear for me but I have more difficulites to understand how to impose the no-stress boundary condition on the outlet boundary. Would it be possible to head me to some existing scripts or page that could help to impose such kind of boundary condition. I Thank you very much, Maxime Fiore

Hello Maxime, Usually the variational formulation for the no-stress condition is obtained by taking any test function without any boundary condition nor boundary integral. For the free slip condition one has to impose \vec u \cdot n = 0 and take any test function satisfying the same normality condition (still without any boundary integral). I think that the formal writing of the free slip condition should include \partial_r v - v/r instead of \partial_r v alone, because of the free rotating solution v=r that should be allowed. François.

You can see some example implementations here. If you integrate the weak pressure gradient and viscous terms by parts, you simply neglect writing the boundary integral to enforce the free outflow condition. Note that this BC is NOT a no-stress condition, even though it is often incorrectly called that in the literature.

Hello François and Chris, I thank you very much for your answers and the explanations, it seems clear to me now. I will have use the example proposed to implement it. Thank you once again, Maxime

Long Time ago I have write slide on this. Subject BC-STOKES.pdf 159.69 KB The example are in the url directory

Hello, I thank you for the document.

You may look at this post for the implementation of the slip condition on a curved boundary. [image] Boundary condition General Discussion If you want to avoid using a Lagrange multiplier (because it enlarges the matrix to invert), some other ways to penalize the normal boundary conditions are possible. In 3D you can use for example the boundary vertices. You have to define first a continuous approximate normal on the boundary (label 1 is for boundary) fespace Vh1(Th,[P2,P2,P2]); varf bdryn([w1,w2,w3],[v1,v2,v3])=int2d(Th,1)(v1*N.x+v2*N.y+v3*N.z); Vh1 [normalappx,normalappy,normalappz]; normalappx[]=bdryn(0,Vh1); Then you can ad…

Hello François, I thank you for your answer, so if I understand well, the difficulty is not really on the outlet boundary condition but rather on inlet and sphere boundaries. Would it be possible to share with you my code because for inlet and sphere boudaries, I tought impose zero velosity on these conditions was enough but it doesn’t seem to be the case if I understand well: include "getARGV.idp" // read the parameters real m = 1; // helical mode m=1 // data from a laminar simulation using Freefem++ string meshname = getARGV("-mesh", "BASE/ffem-mesh.msh"); string urname = getARGV("-ur", "BASE/ffem-ur.dat"); //ur string utname = getARGV("-ut", "BASE/ffem-ut.dat"); // utheta string uxname = getARGV("-ux", "BASE/ffem-ux.dat"); // ux string nuname = getARGV("-nu", "BASE/ffem-nu.dat"); //viscosity // import mesh mesh Th = readmesh(meshname); // prepare the vector spaces fespace Xh(Th, P2); fespace Mh(Th, P1); // import meanfields on P1 space Mh ur, ut, ux, nu; { ifstream fid(urname); fid >> ur[]; } { ifstream fid(utname); fid >> ut[]; } { ifstream fid(uxname); fid >> ux[]; } { ifstream fid(nuname); fid >> nu[]; } // interpolate implicitely on P2 space Xh U=ur, V=ut, W=ux, Visco=nu; Mh PrMean; Xh dUdr=dr(U), dVdr=dr(V), dWdr=dr(W); Xh dUdx=dx(U), dVdx=dx(V), dWdx=dx(W); fespace Vh(Th, [P2, P2, P2, P1]); Vh [F, G, H, Pr], [X1, X2, X3, X4]; // [ur_fluct,utheta_fluct,ux_fluct,p_fluct] Xh Re = 1.0/Visco; // LHS matrix function (sorted by independent to m, linear to m and quadratic to m), y represents r, momentum equations are multiplied by i * y^2 to prevent issues on the axis and have directly the problem LHS * X= omega * RHS * X varf opLHS([F, G, H, Pr], [X1, X2, X3, X4]) = int2d(Th)(1i*( - y^2*X1*(W*dx(F) + H*dUdx + U*dy(F) + F*dUdr) + y*X1*2*G*V - y^2*X2*(W*dx(G) + H*dVdx + U*dy(G) + F*dVdr) - y*X2*U*G - y*X2*V*F - y^2*X3*(W*dx(H) + H*dWdx + U*dy(H) + F*dWdr) - 1/Re*( y^2*dx(H)*dx(X3) + y^2*dy(H)*dy(X3) + y^2*dx(G)*dx(X2) + y^2*dy(G)*dy(X2) + G*X2 + y^2*dx(F)*dx(X1) + y^2*dy(F)*dy(X1) + F*X1 ) + Pr*(y^2*dx(X3) + y*X1 + y^2*dy(X1)) ) -X4*(y*dx(H) + F + y*dy(F)) +1i*(-y*X1*1i*m*V*F -y*X2*1i*m*V*G -y*X3*1i*m*V*H -1/Re*(+2i*m*X1*G - 2i*m*X2*F) +Pr*y*(-1i*m*X2)) -X4*(+1i*m*G) +1i*(-1/Re*(m^2*H*X3 + m^2*G*X2 + m^2*F*X1)) ) +on(1, F=0.0, G=0.0, H=0.0) +on(2, H=0.0, Pr=0.0) +on(3, F=0.0) +on(4, F=0.0, G=0.0, H=0.0) ; // RHS matrix function varf opRHS([F, G, H, Pr], [X1, X2, X3, X4]) = int2d(Th)(y^2*(F*X1 + G*X2 + H*X3)) +on(1, F=0.0, G=0.0, H=0.0) +on(2, H=0.0, Pr=0.0) +on(3, F=0.0) +on(4, F=0.0, G=0.0, H=0.0) ; // The boundaries are defined this way: // o—————————3————————— o // | | // 4 5 // | ╭—1—╮ | // o-----2-----o 0 o————2 ————— o // 1: sphere, 2: axis, 3: top condition, 4: inlet, 5: outlet I thank you for your kind help.

Did you use tgv=-1 or other for the matrix LHS and RHS ? Maybe you should not put boundary condition in the RHS matrix if you use penalty to impose the boundary condition (ex. tgv =1e30) . When you use EPSsolve to solve the eigenvalue problem, the shift-invert methods might be done by the EPSsolve and RHS^{-1} or something like that will be involved. You could refer to the “note” in the document: https://doc.freefem.org/models/eigen-value-problems.html

Besides, the outlet boundary condition might be called " free outflow" more properly. It is not the strict “no stress”. But for the incompressible NS equations, the results should be almost the same if the viscous terms (with “grad of velocity” or “rate of deformation tensor”) in your variational form is written correctly and computational domain is long enough.

Maxime, My message was rather about the external boundary where you impose a slip condition. If I understand well this boundary is curved (in the sense that the normal is not constant on this boundary). On the boundaries where you impose homogeneous Dirichlet I see no particular problem. However you have to be sure (from you modeling situation) of the boundary conditions that you want.

Boundary condition in eigenvalue problem

General Discussion

fb77 (François Bouchut) May 3, 2024, 9:44am 7

You may look at this post for the implementation of the slip condition on a curved boundary.

Topic		Replies	Views
Boundary conditions in eigenvalue problem General Discussion	15	991	February 10, 2024
Boundary condition General Discussion	19	1212	April 25, 2024
On the boundary condition resembling Robin boundary condition General Discussion	13	881	December 14, 2022
Normal boundary condition General Discussion	11	3741	March 20, 2024
How to implement Stress-Free BC General Discussion	1	586	July 23, 2020

Boundary condition in eigenvalue problem

Related topics