# Heaviside function integration

Dear all,

I try to integrate (H o u(x))*v where H is Heaviside function in a 2D mesh by using:

Uh(Th, P1);
varf a(uu,vv) = int2d(Th)((uu == 0. ? 0.5 : (uu>0. ? 1. : 0.))*vv);
matrix Ah = a(Uh,Uh);

However, Freefem raise the following error. I also get the same error with a macro : Heaviside(u) ((u == 0.) ? 0.5 : ((u>0.) ? 1. : 0.)).

varf tesths(uu,vv) = int2d(Th)(((uu == 0.)
error operator == <10LinearCombI7MGauche4C_F0E>, List of choices ( : , ) ( : <10GlgElementIN5Fem2D5MeshLEE>, <10GlgElementIN5Fem2D5MeshLEE> ) ( : <10GlgElementIN5Fem2D5MeshSEE>, <10GlgElementIN5Fem2D5MeshSEE> ) ( : <10GlgElementIN5Fem2D5Mesh3EE>, <10GlgElementIN5Fem2D5Mesh3EE> ) ( : <N12_GLOBAL__N_19lgElementE>, <N12_GLOBAL__N_19lgElementE> ) ( : <NSt7__cxx1112basic_stringIcSt11char_traitsIcESaIcEEE>, <NSt7__cxx1112basic_stringIcSt11char_traitsIcESaIcEEE> ) ( : , ) ( : , )

Does anyone know what is wrong here? Thank you.

So directly what is wrong is that the == operator is not defined for the combination of types you are asking it to take as input. What exactly are you trying to do? There may be an easier way to go about it.

I am trying to evaluate

where H is the function of Heaviside. It is motivated by a contact problem using Nitsche’s method.

More precisely, I have

I understand that Freefem is trying to evaluate u == 0 and it is a nonsense but u(x,y) == 0 makes sense.

So with int2d you are integrating over \Omega not \partial \Omega .
My understanding is that integral will always be 0 as v(x,y) (The test function) evaluated at the boundary is always 0.

Remark, your problem is cleary nolinear and freefem++ is not magic,
You must linearize you problem to solve it !!! with a a fixe point method for exemple
for exemple you can do something like:

set
macro H(u) real (u>0)//
Vh Hu;

Hu = H(u) ;

// you can use Hu

compute Hu and update u form Hu …

Ok I see the problem, thank you very much for you answer.