0
$\begingroup$

Continuing this topic, once the Dini-Neumann problem is solved:

$$ \begin{cases} \Delta f = 0 & \text{on} \; \Omega \\ \frac{\text{d}}{\text{d} n} f = y\,n_x - x\,n_y & \text{on} \; \partial\Omega \end{cases} $$

with $n = (n_x,\,n_y)$ exterior normal vector to $\partial\Omega$, the goal is the calculation of:

$$ \color{blue}{J \equiv \iint_{\Omega} \left(f_x\,y - f_y\,x\right)\text{d}x\,\text{d}y}\,. $$

On the other hand, we also have:

$$ J = \iint_{\Omega} \left[\left(f\,y\right)_x - \left(f\,x\right)_y\right]\text{d}x\,\text{d}y $$

which by the divergence theorem is equivalent to:

$$ J = \int_{\partial\Omega} f\left(y\,n_x - x\,n_y\right)\text{d}s $$

and for the boundary condition:

$$ J = \int_{\partial\Omega} f\left(\frac{\text{d}}{\text{d}n}f\right)\text{d}s = \int_{\partial\Omega} \left(f\,f_x\,n_x + f\,f_y\,n_y\right)\text{d}s\,. $$

This done, again by the divergence theorem:

$$ J = \iint_{\Omega} \left[(f\,f_x)_x + (f\,f_y)_y\right]\text{d}x\,\text{d}y = \iint_{\Omega} \left(f_x^2 + f_y^2 + f\,\Delta f\right)\text{d}x\,\text{d}y $$

that is, since $f$ is an harmonic function on $\Omega$:

$$ \color{\red}{J = \iint_{\Omega} \left(f_x^2 + f_y^2\right)\text{d}x\,\text{d}y}\,, $$

known in the literature as the Dirichlet integral.


Therefore, after determining the harmonic function:

Needs["NDSolve`FEM`"];

Ω = ToElementMesh[ImplicitRegion[x^4 + y^4 <= 1, {x, y}], "IncludePoints" -> {{0, 0}}];

f = NDSolveValue[{Laplacian[g[x, y], {x, y}] == NeumannValue[y (4 x^3) - x (4 y^3), True], 
                  DirichletCondition[g[x, y] == 0, x == 0 && y == 0]}, g, {x, y} ∈ Ω];

and verified that it's accurate enough:

NIntegrate[f[x, y], {x, y} ∈ Ω]
NIntegrate[Laplacian[f[x, y], {x, y}], {x, y} ∈ Ω]

-0.00000293157

0.00292529

I tried to calculate the Dirichlet integral in its two equivalent versions:

NIntegrate[D[f[x, y], x] y - D[f[x, y], y] x, {x, y} ∈ Ω]
NIntegrate[D[f[x, y], x]^2 + D[f[x, y], y]^2, {x, y} ∈ Ω]

-0.48313

1.84564

but it's evident that something went wrong, but I don't know what! Ideas? Thank you!


EDIT: Thanks to the precious comment of xzczd I understood that:

  • $(n_x,n_y)$ must necessarily have unitary length;

  • before the Laplace operator it must absolutely be a minus sign.

In fact, writing:

c[x_, y_] = x^4 + y^4;
{nx, ny} = Normalize[Grad[c[x, y], {x, y}]];

Needs["NDSolve`FEM`"];
Ω = ToElementMesh[ImplicitRegion[c[x, y] <= 1, {x, y}],
                  "IncludePoints" -> {{0, 0}}];

f = NDSolveValue[{-Laplacian[g[x, y], {x, y}] == NeumannValue[y nx - x ny,
                  True], DirichletCondition[g[x, y] == 0, x == 0 && y == 0]},
                  g, {x, y} ∈ Ω];

NIntegrate[f[x, y], {x, y} ∈ Ω]
NIntegrate[Laplacian[f[x, y], {x, y}], {x, y} ∈ Ω]

NIntegrate[D[f[x, y], x] y - D[f[x, y], y] x, {x, y} ∈ Ω]
NIntegrate[D[f[x, y], x]^2 + D[f[x, y], y]^2, {x, y} ∈ Ω]

we get:

0.000000504282

-0.000799624

0.126589

0.126595

which is approximately what I want. The only thing that would remain to understand is whether it were possible to obtain a better approximation of those integrals. Thanks again!

$\endgroup$
2
  • 2
    $\begingroup$ Can you add a reference for Dirichlet integral? It doesn't seem to belong to those mentioned in wiki page: en.wikipedia.org/wiki/Dirichlet_integral $\endgroup$
    – xzczd
    Dec 5, 2022 at 2:07
  • 3
    $\begingroup$ For divergence theorem, $(n_x,n_y)$ should be unit vector, isn't it? $\endgroup$
    – xzczd
    Dec 5, 2022 at 8:34

1 Answer 1

2
$\begingroup$

Let's make the mesh denser by setting MaxCellMeasure option:

Ω = 
  ToElementMesh[ImplicitRegion[c[x, y] <= 1, {x, y}], "IncludePoints" -> {{0, 0}}, 
   MaxCellMeasure -> 0.001];

Then the output of

NIntegrate[D[f[x, y], x] y - D[f[x, y], y] x, {x, y} ∈ Ω]
NIntegrate[D[f[x, y], x]^2 + D[f[x, y], y]^2, {x, y} ∈ Ω]

is

0.126603
0.126603
$\endgroup$
0

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.