Solution of Poisson equation with two regions

Question

I'm trying to figure out how is possible to solve a Poisson equation

$\nabla\cdot[d(x,y)\nabla u]+1=0$

where $d(x,y)$ equals 1 in one region and 2 in another one. Let say I have homogeneous Dirichlet boundary conditions.

The two regions should be physically distinct, I mean do not use just

d[x_,y_]:=If[x<0.5,1,2]

if x=0.5 is the edge between the two regions.

Thanks for the suggestion(s) F

I've tried the "If" solution, and also merging two regions. In the last case there is no internal boundary, further I do not know how to address $d(x,y)$ for the two regions. — Fabio, Nov 13 '15 at 16:56
I am not sure I understand, could you not solve the equation once per region with the appropriate coefficient - maybe you could clarify what you mean by "distinct" - a code example would be best. — user21, Nov 13 '15 at 17:47
Can this be a case for "matched asymptotics" (perturbation methods)?.. or am I over-complicating the problem? — dearN, Nov 13 '15 at 17:51
No, the problem I've illustrated implies that across the interface of the two domains both the function $u$ and the normal flux $d\nabla u\cdot n$ must be continuous. Setting a system of PDE by requiring these constraints would be ok, too. But how to do that? If I had the code I've solved the problem. — Fabio, Nov 13 '15 at 17:53
I've found something in here https://reference.wolfram.com/language/FEMDocumentation/tutorial/SolvingPDEwithFEM.html#509267359 — Fabio, Nov 13 '15 at 18:15
@user21 I think this is related to your answer here. In Fabio's case the automatic feature detection seems to fail and no internal boundary is created by NDSolve. His solution is to use NDSolve`FEM directly to gain better control over the mesh, which I guess is what you'd recommend? — Simon Woods, Nov 14 '15 at 13:55
@SimonWoods, well observed. If one changes the region to {x, y} \[Element] ImplicitRegion[0 <= x <= 1 && 0 <= y <= 1, {x, y}] then the automatic feature detection works. With time more and more regions and triggers should work. — user21, Nov 15 '15 at 07:13

score 4 · Answer 1 · answered Nov 13 '15 at 18:46

I've found

<< NDSolve`FEM`
bndmesh = 
  ToBoundaryMesh[
   "Coordinates" -> {{0, 0}, {0.5, 0}, {0.5, 1}, {0, 1}, {1, 0}, {1, 
      1}},
   "BoundaryElements" -> {LineElement[{{1, 2}, {2, 3}, {3, 4}, {4, 
        1}, {2, 5}, {5, 6}, {6, 3}, {3, 2}}]}];
mesh = ToElementMesh[bndmesh];
d = If[x <= 0.5, 1, 10];
usol = NDSolve[{Div[d Grad[u[x, y], {x, y}], {x, y}] + 1 == 0, 
     DirichletCondition[u[x, y] == 0, 
      x == 0 || x == 1 || y == 0 || y == 1]}, 
    u, {x, y} \[Element] mesh][[1]];
{bndmesh["Wireframe"], mesh["Wireframe"], 
 Plot3D[u[x, y] /. usol, {x, y} \[Element] mesh]}

as a possible solution.

Any enhancement/advice to that?

Solution of Poisson equation with two regions

1 Answers1