Dear all,
I met a problem these days. I have to split my domain into 1 and 2, because different equations are solved in Th1 and Th2.( for example, q1=[u,v,p] for the imcompressible NS equations, while q2 = T for the thermal conduction). I have already known how to solve the base-flow. However, I didn’t know how to make the connection condition to a matrix(such as the non-diagonal matrix in the picture), so I could solve the eigenvalue problem of the whole system.
Thanks a lot!
I have thought about that I might be able to introduce a third vector qc (the quantity in the connection boundary C) , so I could solve for the (q1, qc, q2). But I am not sure whether this will work or is there more direct way in the FreeFem++?
This is not quite right certainly not optimal but thought it was an intersting
approach. I guess you could define two meshes but then you need to
fabricate the boundary equation probably with some messy interpolation matrix
although that is easy to generate in FF with “interpolate” command. It is probably
possible to just use one dof and a discontinuous element assigning the boundary
to one or the other variables. However, I’ve never used those and myabe
someone can comment. This solution appears to work although it wastes
a lot of dofs. And note the picture just sums the variables together causing the
boundary “value” to be inflated ( it is not a sum ) but ok for this picture.
I made up BC and relation between vars but this appears to work and i simple
to write. Thoughts? Maybe someone can modify it to use just one dof if you
have the same element types, this example used diffeing ones 
load "medit"
macro buildouter()
border left(t=0, 1){x=-szx/2; y=szy*(t-.5);label=1; }
border top(t=0, 1){x=szx*(t-.5);y=szy/2; label=2; }
border right(t=0, 1){x=szx/2; y=-szy*(t-.5);label=3; }
border bottom(t=0, 1){x=-szx*(t-.5);y=-szy/2; label=4; }
// EOM
real szx=2;
real szy=1;
int nx=20;
int ny=20;
buildouter
border middle(t=0, 1){x=0; y=-szy*(t-.5);label=5; }
mesh Th = buildmesh(left(-ny)+top(-nx)+right(-ny)+bottom(-nx)
+ middle(ny));
//medit("Th",Th);
//plot(Th,wait=1,fill=1);
//medit("Th",Th);
macro Op(xx,vxx) (dx(xx)*dx(vxx)+dy(xx)*dy(vxx)) //
fespace Vh1(Th,P1);
fespace Vh2(Th,P2);
Vh1 a1,va1,b1=0;
Vh2 a2,va2,b2=0;
// https://doc.freefem.org/examples/mesh-generation.html
int rl=Th(-szx/2,0).region;
int rr=Th(szx/2,0).region;;
cout<<" rl="<<rl<<" rr="<<rr<<endl; cout.flush;
problem both([a1,a2],[va1,va2],solver=sparsesolver,tgv=-1)
//problem both([a1,a2],[va1,va2],solver=CG)
= int2d(Th,rl)(Op(a1,va1))+int2d(Th,rl)(b1*va1)
+ int2d(Th,rr)(Op(a2,va2))+int2d(Th,rr)(b2*va2)
// otherwise use CG
+int2d(Th,rr)(a1*va1)+int2d(Th,rl)(a2*va2)
// a1 to a2 boundary
+int1d(Th,5)(a2*va2)-int1d(Th,5)(a1*2.0*va2)
+int1d(Th,5)(a2*va1)-int1d(Th,5)(a1*2.0*va1)
+on(1,a1=2*y)+on(3,a2=-3*y)
;
both;
Vh1 net=a1+a2;
plot(a1,cmm="a1",wait=1,fill=true,value=1);
plot(a2,cmm="a2",wait=1,fill=true,value=1);
plot(net,cmm="net",wait=1,fill=true,value=1);
medit("net",Th,net);
Thanks very much!
I have tried this method in an imcompressible flow to make a test (for example, the wake of the cylinder flow), I found some interesting things~
First, for the steady flow calculation(with Newton iteration), there may exist some numerical oscillation in the ux and p (but I think for a lower Re number, it’s ok, just like the follwing results; however the results may diverge for a higher Re , like Re=47 ).
Sorry, what is the nature of the discontinuity around -10?
If the domains are physically connected as is common in these kinds
of problems the contiguous mesh would seem to make sense.
Boundary conditions are a hassle and if you can define a material property
at each mesh point that may make more sense in all ways. You can eliminate
one of the variables if they are the same “thing” but just in different regions.
This should also work well if you expect the boundary to move.
I haven’t used the eigen stuff enough to comment not even sure if they
are normalized but I take it they look like noise instead of modes.
No idea.
.
This is just an example to test my code… in this example, there is no physical discontinuity.
Actually my target problem is that: in the left part Navier Stokes equations are solved in the low mach number limit and in the right part is the wave equation for the acoustic wave.
Now I have already solved this problem by defining different varfs and the off-diagonal matrixs with interpolating. (I use the overlapping mesh now, and give a source term in the overlapping part to make it change more smoothly. Although I lose some accuracy due to this method, but the results are acceptable…)
Thanks again~
