Derive on a specific border for a particular calculus

General Discussion
1

Dear users,

In order to obtain a partial reflective B.C. condition for a Helmholtz equation, I’d like to calculate the angle of incidence of a waterwave on a border.
I call u the solution of the equation, and the researched angle is $\gamma$.

Here is the formula (Steward - Panchang, 1999) :

$f(\gamma) = tan(\gamma) - \frac{\frac{1}{k} \frac{\partial arg(u)}{\partial x_s}}{ \frac{1}{k} \frac{\partial arg(u)}{\partial x_n} + \frac{2K_r(cos(k\beta)+K_r)}{1+2K_rcos(k\beta)+K_r^2}cos(\gamma)} $

$u$ is obtained by solving the equation a first time, with $\gamma= 0$.
I use a P1 triangulation.
Here I’d like to determine $\gamma$, knowing k, K_r, $\beta$, and $u$, buy solving $f(\gamma) = 0$ (with bracketing and disection).
x_n is the normal, and x_s is the tangent to the border, and arg(u) is the argument of u.

My problem is to calculate the derivatives of u.

I saw that I could obtain the tangent to the border by rotating the normal by pi/2.

Then on another topic, I read that I could collect the values of u on the border : getting data of specific points

So I tried to collect those values of u on the border, put them in a FE variable, take the argument of it, and apply the “dx” and “dy” operators to it. But I had an error.

Then my question is :
Do you have any idea on how to obtain $darg(u)/dxs$ or $darg(u)/dxn$ on a specific border ?

Here is a simple code I made to test the method :

//param
real k = 0.002; // wavenumber

real Kr = 0.4; //reflection coeff

real theta = 11pi/8 ;
real k1 = kcos(theta);
real k2 = k*sin(theta);

real internRadius = 2000;
real externRadius = 5000;

// interior disc
border Cint(t=0,2pi){x= internRadiuscos(t) ; y = internRadius*sin(t) ;}

// exte disc

border CF(t=5/4.pi, 3/2.pi){x = externRadiuscos(t) ; y = externRadiussin(t);}
border Cext(t=-pi/2, 5/4.pi){x = externRadiuscos(t) ; y = externRadius*sin(t);}

mesh Th = buildmesh(Cint(-90) + CF(45) + Cext(90));

//Finite element space

fespace Vh(Th, P1);
Vh uh, vh, uinc = 0.5exp(-1i(k1x+k2y)); // uinc = incoming wave

Vh u,ur, uim; (for the real/imag solution)

Vh argu; // dNargu; // argu = argu of uh, dNargu = normal derivative ?

// Helmoltz
solve Helmoltz(uh,vh) = -int2d(Th)(dx(uh)dx(vh) + dy(uh)dy(vh)
-k^2uhvh)
-int1d(Th,Cint)(1ik((1-Kr)/(1+Kr))uhvh)
+on(CF, uh = -uinc);

//solution

// Obtain uh on the interior circle

// boundary dof
int ndof;
{ int1d(Th,Cint,qfe=qf1pE)((ndof++)*1.); }

// indexes on the boundary
int[int] i1(ndof);
{
int[int] a(2ndof);
Vh index=0;
for(int i=0;i<Vh.ndof;i++) index[][i]=i;
int n=0;
int1d(Th,Cint,qfe=qf1pElump)((a(n++)=floor(index+0.5))1.);
Vh next=0;
for(int i=0;i<ndof;i++) next[][a(2i+1)]=a(2i);
i1(0)=a(1);
for(int i=1;i<ndof;i++) i1(i)=next[][i1(i-1)];
}

//display the value of uh on the boundary Cint
//for (int i=0;i<ndof;i++){
//cout<<i1(i)<<" uh("<<Th(i1(i)).x<<","<<Th(i1(i)).y<<")="<<uh[][i1(i)]<<endl;
//;
//}

//How to get darg(uh)/dn on Cint??

// plot the solution
u = real(uh);

plot(u);

Thank you all for reading,

Best

2

I have a problem the function argument is only defined with complex number ur+ ui i
and it is atan2(ui, ur) = imag( log( (ur + uii)/sqrt(ur^2+ui^2)) )
then you can derive with maple in function of ui and ur , so not difficulty.

otherwise what is function argument.

3

Hello,

Thank you for your answer.

Indeed, the solution u will be a complex number so I will be able to take its argument easily.
But my problem is that I can’t find how to get all the values of grad(argu(u))*n, where n is the normal to the border on only a single border of my domain.

I’d like to get those values because I’ll have to do calculations with them.

I’m not used to maple, i’ll check this right away.

Thanks again,

Regards

4

you can compute formerly

func dxargu = dx(argut(u))
func dyargu =dy(argu(u)),

so you can do [dxargu, dyargu]’ *[N.x,N.y]

5

I will try using this, thanks a lot !