Gradient for surface FEM

General Discussion
1

Hi,

There are examples of surface FEM for a scalar function, such as LapEigenBeltrami.edp, however I haven’t found an example of “expressing the gradient of a vector on a surface mesh ThS”.

For a scalar u, the surface gradient is computed as follows:

meshS ThS=square3(nx,ny,[torex,torey,torez],removeduplicate=true) ;

fespace Vh(ThS,P1);

Vh u.

macro Grad3(u) [dx(u),dy(u),dz(u)]  // EOM

real sigma = -1;
varf aS(u,v) = int2d(Th)(Grad3(u)'Grad3(v)- sigma uv);
varf mS(u,v) = int2d(Th)( uv);

But how to express the gradient of a vector (u1, u2) (intrinsic form: two components in the local coordinates)?

fespace VhS(Th, [P1, P1] );
Vh [u1, u2] .

My week form in the local coordinate evolves the computation of

image

where

image

with g_ij being the metric tensor, g^ij being its inverse, and " ; “ denotes the covariant derivative. Alternatively, the non-symmetric form may also be used if it’s easier:image

How to express these expressions on the surface mesh?

For the following code,

fespace VhS(Th, [P1, P1] );
Vh [u1, u2] .

If we express u1 and u2 using the shape function \phi, the derivative ultimately acts on the shape function \phi (x, y, z). The relevant question is: how does FreeFEM compute dx(\phi)?

In some FEM software, it’s only meaningful to define the derivatives with respect to the local coordinates. Can I say

macro Grad3(u) [dx(u),dy(u),dz(u)]  // EOM

has not been projected to the surface yet and the user has to do the projection?

If my assumption is true, then it seems the weak form in exmaple LapEigenBeltrami.edp has a problem: (Grad3(u)'Grad3(v) includes the normal component, which is wrong?

2

Hi,
For a scalar u defined on a meshS surface, the gradient \nabla u defined as [dx(u),dy(u),dz(u)] is a 3d vector which is tangent to the meshS. In means that interpreted as a element of the tangent space, it is the gradient in the Riemanian sense. In particular, if \bf x is on the surface and \bf y also an the surface but approaching \bf x, one has u(\bf y)\simeq u(\bf x)+\nabla u\cdot (\bf y-\bf x).
This is true indendently of any choice of local coordinates.

In the case of a vector field, I am not sure on how to proceed.

A first approach is to define a 3d vector with the space [P1,P1,P1]. I think that you would like to describe a tangent vector field, thus when defining this 3d vector field you have to built it so that it is tangent to the surface.
Then you can use the operators dx(), dy(), dz() on each of the 3 components of the vector field to build up the gradient. You have to keep in mind that only the 3d vector [dx(),dy(),dz()] is geometrically meaningful (indeed dx() gives the x component of the gradient, dy() the y component, dz() the z component).

A second approach is to consider 2 scalar coordinates by decomposing your vector field in a local basis of the tangent space. But then you have to know what is the gradient of this local basis. It could be computed by the first approach.

3

Thank you a lot @fb77 for your comments!
If the vector is defined by an extrinsic form u = (ux, uy, uz) in the Euclidean space, maybe the projection operator P = I - n x n can be used to build the surface gradient ,

func Pk=[P1,P1,P1];
fespace Vh(ThS, Pk);

Vh [u1, u2, u3],[uh1, uh2, uh3];

macro Nor [N.x, N.y, N.z] //
macro vect(u) [u#1, u#2, u#3]//
macro P(u) (vect(u) - (vect(u)'*Nor)*Nor ) //

Because u has to be projected to the tangent space first, then take the gradient, and then project to the tangent space again, i.e.: something like
P (Grad(Pu))

The macro is a bit messy…. and I haven’t got it right yet.

Moreover, I am looking for an intrinsic form with u = u^i \partial_i being expressed locally. It seems that the integration in the parameter space together with metric tensor is more convenient than the surface FEM.

4

The normal to the surface is [Ns.x,Ns.y,Ns.z].
If your metric is attached to the scalar product inherited from R^3, there is hope that your operator D(u) can be expressed in terms of the freefem gradient applied to the 3 components [u1,u2,u3].
Then if [u1,u2,u3] is your unknown and needs to satisfy u\cdot N_s=0, you could use a Lagrange multiplier formulation.
It leads to look for u1,u2,u3 scalars defined on ThS, and also p scalar defined on ThS such that (for example), denoting v1,v2,v3 and q the test functions
\int_S \nabla[u1,u2,u3]:\nabla[v1,v2,v3]+\int_S p [v1,v2,v3]\cdot N_s+\int_S q[u1,u2,u3]\cdot N_s-\int_S f\cdot [v1,v2,v3]=0

Concerning the expression of your “intrinsic” operators, consider \xi_1,\xi_2 coordinates for your surface and \frac{\partial x}{\partial\xi_i} the associate two vectors (i=1,2), that span the tangent space (x stands for the position in R^3).

As a start consider a scalar P. Then \frac{\partial P}{\partial\xi_i}=\nabla P\cdot\frac{\partial x}{\partial\xi_i} for i=1,2. Denote by G the 3\times 2 matrix defined by G_{ij}=\frac{\partial x_i}{\partial\xi_j}.
Then [\frac{\partial P}{\partial {\xi_1}},\frac{\partial P}{\partial {\xi_2}}]=G^t\nabla P
Next for a tangent vector u=[u1,u2,u3] you can apply the previous formula for each component u1,u2,u3 and you get an expression of the partial derivatives \frac{\partial u_i}{\partial\xi_j} in terms of the (matrix) gradient of u=[u1,u2,u3].
Then you need to express your two local coordinates of u in terms of the vector [u1,u2,u3], hopefully obtained via the scalar product of [u1,u2,u3] with the vectors \frac{\partial x}{\partial\xi_i} for i=1,2, and differentiale the result with respect to \xi_i …