Points on the interface

Question

We consider the problem

$\left\{\begin{matrix} k(x)\Delta u(x)=f(x) & \text{ in } \Omega\\ u=0 & \text{ in } \Gamma \end{matrix}\right.$

where $\Omega \subset \mathbb{R}^2$ open and bounded with smooth boundary $\Gamma=\partial{\Omega}$ and $k(x) \geq a, \forall x \in \Omega$.

We consider that $\Omega= \Omega_1 \cup \Omega_2$ and that $k$ is always a constant $k(x)=\left\{\begin{matrix} k_1 & \text{ in } \Omega_1 \\ k_2 & \text{ in } \Omega_2 \end{matrix}\right.$

$\Sigma$ is the boundary between $\Omega_1$ and $\Omega_2$, $\Sigma= \overline{\Omega_1} \cap \overline{\Omega_2}, \Gamma_j=\Gamma \cap \partial{\Omega_j}, j=1,2$.

We suppose that $\Omega=[-1,1]^2$ and $\Sigma=[-1,1] \times \{ y=0 \}, \Omega_1=\Omega \cap \{ y>0 \}$ and $\Omega_2=\Omega \cap \{ y<0 \}$.

We construct a grid $\Omega$ so that the interface $\Sigma$ coincides with the sides of the squares.

We write the problem as $AU=f$ in order to find the approximation of the solution.

For the approximations that correspond to the points that are over $\Sigma$ we pick $k1$ and for the approximations that correspond to the points that are under $\Sigma$ we pick $k2$.

What do we do for the approimations that correspond to the points that are on $\Sigma$ ?

EDIT:Finally, we don't want the approximation of the solution to be continuous. We may take $k=\frac{k_1+k_2}{2}$ for the points on the interface. So do we have to take cases for $y$ in order to give a value to $k$?

Answer 1

Rewrite your problem as $$\begin{cases} \Delta u(x)=\frac{f(x)}{k_1}& \text{ in } \Omega_1\\ \Delta u(x)=\frac{f(x)}{k_2}& \text{ in } \Omega_2\\ u=0 & \text{ in } \Gamma \end{cases}$$ you have to add a condition on an interface, for example a jump condition, e.g. $$u^1(x)-u^2(x)=\alpha(x)\quad\text{and}\quad u^1_n(x)-u^2_n(y)=\beta(x) ,\quad x\in\Sigma$$ or smoothness (e.g. $\alpha=0=\beta$) of the solution and its normal derivative across the interface. In the above $u^i$ is the solution in the domain $\Omega_i$, $i\in\{1,2\}$.

One of the best ways to solve such problem would be Difference Potentials Method (based on good theory, need time to learn, worth the time you spend on it), another is Immersed Interface method(easy to understand, implementation is involved and changes with the problem, less accurate). You find a little more under the link provided.

Regularization approach

If the low order of accuracy is enough for you, you may want to try something very simple. You take a continuous approximation of step function and solve a little bit different, but continuous problem, e.g.: $$\begin{cases} \Delta u(x)=\frac{f(x)}{g_p(x)}& \text{ in } \Omega_1\cup \Omega_2\\ u=0 & \text{ in } \Gamma \end{cases}$$ where (you can choose another approximation of step function) $$g_p=\frac{k_2-k_1}{1+e^{-M(x-p)}} + k_1. $$ $M$ supposed to be a big number, the larger the $M$ the sharper the "jump", i.e. the closest to what you need. You can play with different values to see the effect. Note, when $x<p$ the denominator is a big number and so the fraction become close to zero. when $x>p$ the exponential in the denominator is negligible and therefore the fraction is almost the nominator.

Thus, this way you solve a problem without singularities or jumps. However, you solve a little bit different problem, which is asymptotically (as $M\to\infty$) tends to the original one.