Problems defining a non-linear form

General Discussion
1

Hello everyone,

I have a running code that solves two linear equations successively: first I solve for a polarity vector and then I solve for the velocity vector in a force balance equation with some forces depending on the (previously solved) polarity. The second equation has a bi-linear term representing the friction with the form

\int_\Omega \xi \,\vec{v}\cdot\vec{u},

With \xi being the friction coeficcient, \vec{v}, the velocity vector and \vec{u}, the vectorial test function. I want to modify this model so that the fricition depends on the orientation between the polarity and velocity vectors. With that the friction form becomes

\int_\Omega \xi \,\left(1-\frac{\vec{p}\cdot\vec{v}}{||\vec{p}|| ||\vec{v}||}\right)\,\vec{v}\cdot\vec{u},

which adds a non-linear term to the first bi-linear form. At the time that I am solving for the velocity, the polarity \vec{p} is known, so this is just a quadratic form in \vec{v}. I am aware that I should implement a Newton-Raphson to solve this non-linear form, which is not a problem, but I am unable to properly define it in FreeFem.

My code goes as

macro p [px,py]  //
macro v [vx,vy]  //
macro u [ux,uy]  //

macro normInner(p,v) ((p'*v)/(sqrt(p'*p)*sqrt(v'*v))) // normalized inner product

Then I define the FEM spaces, solve for p and try to define the friction form as

// p-v alignment-dependent friction
varf aLinFric(v,u) = int2d(Th)(xi*(v'*u));
varf aNLinFric(v,u) = int2d(Th)(xi*normInner(p,v)*(v'*u));

resulting in the error

>> error operator *  <10LinearCombI7MGauche4C_F0E>, <10LinearCombI7MGauche4C_F0E>

at the line where I define macro normInner(p,v). I’ve also tried with defining the macro more explicitly as

macro Inner(px,py,vx,vy) (px*vx + vy*py) //
macro normInner(px,py,vx,vy) (Inner(px,py,vx,vy)/(sqrt(Inner(px,py,px,py))*sqrt(Inner(vx,vy,vx,vy)))) // normalized inner product

but I end up obtaining the same error. I understand that I am not properly combining the FEM solutions, but I do not manage to understand how should I define a non-liner term depending on vectors like that.

Thank you so much,

Joan Térmens

2

This is non linear probem now a FreeFem solve directly only linear problem.

macro V(p)  [p#x,p#y] // 
macro Norm(p)  dist(#x,p#y) //def: || V(p)|| 

macro pdot(u,v) (V(u)'*V(v))// u.v 
macro nolinear(p,vp) (1.- pdot(p,vp)/(Norm(p)*Norm(vp))// vp previous v 

mesh Th= .....;
fespace Vh(Th,[P1,P2]); 

Vh V(p),V(u),V(vp);

non linear loop
{
   vpx[]=vx[]; // copie in previous vector
  compute p here ...

 compute  the new v depend ending of vp and p. here ..

vpx[]-=vx[];  the difference
 if( vpx[].linfty < epsilon )break ; // converge 
 

}
3

This way it works,

Thank you.

Joan Térmens