# Coefficient that varies on domain

Hello, I’m trying to implement a coefficient value that varies depending on the area of the domain, this is the definition

``````//parameters
func real kappa() {
int ni = 9, nj = 9;
real r1=1./32.*3./ni, k1=1., k3=1.;
real r2=1./16.*3./ni, k2=10000.0;

for(int j=0; j<(nj-1); j++){
real yj=(2.*j+1.)/(2.*nj);
for(int i=0; i<(ni-1); i++){
real xi=(2.*i+1.)/(2.*ni);
if( (x>=(xi-r1)) && (x<=(xi+r1)) && (y>=(yj-r1)) && (y<=(yj+r1)))
return k1;
if ( (x>=(xi-r2)) && (x<=(xi+r2)) && (y>=(yj-r2)) && (y<=(yj+r2)))
return k2;
}
}
return k3;

}
``````

I need this coefficient for the integral for the Darcy Model

``````   varf ahK(etaK, varphi)
``````

but I’m getting the following error

``````varf ahK(etaK, varphi)
289 :     =  int2d(Th[K])(kappa*grad(etaK)   [dx(etaK), dy(etaK)]'* error operator *  <11Polymorphic>, <10LinearCombI7MGauche4C_F0E>
``````

can someone help me?

this is an example of the code I want to implement on FreeFem++

``````local function func_kappa(x,y)
ni=27; nj=27;
-- (not so) aligned to the mesh
r1=1./32.*3./ni ; k1=1.; k3=1.; -- radius and kappa of the inner squared cell
r2=1./16.*3./ni ; k2=10000.0; -- radius and kappa of the outer square cell

for j=0,(nj-1) do
yj=(2.*j+1.)/(2.*nj);
for i=0,(ni-1) do
xi=(2.*i+1.)/(2.*ni);
if ( (x>=(xi-r1)) and (x<=(xi+r1)) and (y>=(yj-r1)) and (y<=(yj+r1)))  then
return k1;
end
if ( (x>=(xi-r2)) and (x<=(xi+r2)) and (y>=(yj-r2)) and (y<=(yj+r2)))  then
return k2;
end
end
end
return k3; -- kappa of the embedding matrix
end

kappa = func_kappa,
``````

The technics is very slow, you function kappa is defined by
with a Gh grid (x_i,y_i) of size ni, nj you want

let us call dist00(x,y) is the minimal distance of (x,y) the the point of the grid in norme infini
now the fonction is just k3+(dist00(x,y) < r1)(k1-k3) + (dist00(x,y) < r2)(k2-k1-k3)

the trick to code dist00(x,y) arrive to morrow!

because

It worked like this:

``````//parameters
func real varkappa(int in, int jn){
int ni = in, nj = jn;
real r1=1./32.*3./ni, k1=1., k3=1.;
real r2=1./16.*3./ni, k2=10000.0;

for(int j=0; j<(nj-1); j++){
real yj=(2.*j+1.)/(2.*nj);
for(int i=0; i<(ni-1); i++){
real xi=(2.*i+1.)/(2.*ni);
if( (x>=(xi-r1)) && (x<=(xi+r1)) && (y>=(yj-r1)) && (y<=(yj+r1))){
return k1;
}
if ( (x>=(xi-r2)) && (x<=(xi+r2)) && (y>=(yj-r2)) && (y<=(yj+r2))){
return k2;
}
}
}
return k3;
}

real kappa = varkappa(9,9);
``````

so my final version (the complexity of function Kappa is O(1) not O(ni*nj)

yyy.edp (1.1 KB)

1 Like