Laplacian operator

Hello every one,

I want to compute the Laplacian operator of a vector u in Vh (Th, P1) finite element space.

I use two different approaches:

  1. macro Laplacian(u) (dxx(u) + dyy(u)) //
    Vh u;
    Vh lap_U = Laplacian(u);

with this I always obtain a vector of zeros.

  1. macro div(u1,u2) (dx(u1) + dy(u2))//
    Vh u;
    Vh u1=dx(u);
    Vh u2=dy(u);
    Vh lap_U = div(u1,u2) ;

wirh this I obtain a vector different of zeros.

Is there any recommandations for the computation of the Laplacian operator.

Best regards.

These two methods are very much not the same.

Since P1 finite elements are linear, the second derivative should be zero as you have computed it the first way. If you want a nonzero answer here, I think you will need to define your function on a higher-order FE space.

The second way you’ve programmed it, you are doing something different. First, you compute the gradient of u. Then you INTERPOLATE that gradient (which is piecewise constant in each component) onto new P1 FE functions. Finally, you compute the divergence of those functions to get your answer. I don’t think this is what you want to do…

1 Like

Thank you very much dear Chris for your response.

If I understood your answer correctly, I have to convert my vector u from P1 to P2 finite element space, then compute the Laplacian operator using the first approach.

Code :

Vh1 (Th,P1);
Vh1 u ; //My vector
Vh2 (Th,P2);
Vh2 w;
w=v; //Conversion of the vector u from P1 to P2
Vh2 lap_U=dxx(w)+dyy(w);

If you actually need the value of the Laplacian, it would be better to solve the whole problem with P2 or higher-order elements and compute the Laplacian directly. Here, your code interpolates a bilinear function onto a biquadratic function, which could lead to spurious wiggles.

Thank you very much Chris for your response.

1 Like