1
$\begingroup$

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

Improve this question
$\endgroup$
9
  • 1
    $\begingroup$ Well, what have you tried so far? $\endgroup$
    – MarcoB
    Commented Nov 13, 2015 at 16:51
  • $\begingroup$ 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. $\endgroup$
    – Fabio
    Commented Nov 13, 2015 at 16:56
  • $\begingroup$ 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. $\endgroup$
    – user21
    Commented Nov 13, 2015 at 17:47
  • $\begingroup$ Can this be a case for "matched asymptotics" (perturbation methods)?.. or am I over-complicating the problem? $\endgroup$
    – dearN
    Commented Nov 13, 2015 at 17:51
  • $\begingroup$ 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. $\endgroup$
    – Fabio
    Commented Nov 13, 2015 at 17:53

1 Answer 1

Reset to default
4
$\begingroup$

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]}

enter image description here

as a possible solution.

Any enhancement/advice to that?

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.