# Non-Dimensional Equation error

**URL:** <https://community.freefem.org/t/non-dimensional-equation-error/489>\
**Category:** General Discussion\
**Created:** [July 4, 2020, 5:42pm UTC](https://community.freefem.org/t/non-dimensional-equation-error/489 "2020-07-04T17:42:47Z")\
**Posts on this page:** 1\
**Page:** 1

<div class="post-metadata">

**Author:** ![sumantkr](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/sumantkr/32/2681_2.png) [@sumantkr](https://community.freefem.org/u/sumantkr)\
**Post date:** [July 4, 2020, 5:42pm UTC](https://community.freefem.org/t/non-dimensional-equation-error/489/1 "2020-07-04T17:42:47Z")

</div>

Dear @frederichecht, @prj Sir,

I am trying to formulate the non-linear equation based on your suggestions and link provide before.

[topic name](https://community.freefem.org/t/operation-error-in-problem-definition/466/5)

[relevant problem by prj](https://www.um.es/freefem/ff++/pmwiki.php?n=Main.Newton3D)

[relevant problem by frederichecht](https://www.ljll.math.upmc.fr/hecht/ftp/ff++/2018-XVIII-las-Palmas/edp-XVIII/04Rayleigh-Benard.edp)

I am trying to solve the same problem base on the suggestion and adopted the procedure given by [prj Problem link](https://www.um.es/freefem/ff++/pmwiki.php?n=Main.Newton3D) and [Research problem by frederichecht](https://www.ljll.math.upmc.fr/hecht/ftp/ff++/2018-XVIII-las-Palmas/edp-XVIII/04Rayleigh-Benard.edp).

but, unfortunately, the error at each iteration not decreasing after step passing.

Kindly Suggest me to fix it (Convergence ):

```
// all coefficient used in problem

real Ra = 1e6;
real Pr = 0.71;
real Da = 1e-5;
real pRaPr = Ra*Pr;
real dPrDa = Pr/Da;

real[int] aRa = [1e3];

int nn = 20;
real ccc = 1.5;//2.;
macro fff(x) ((1+tanh(ccc*(x*2-1))/tanh(ccc))/2)//

mesh Th=square(nn,nn,[x,(y*(1-x*0.0))]);
plot(Th,wait=1);

// finite element fespace

fespace Vh(Th,[P2,P2,P1,P2]);

Vh [u1, u2, p, T];
Vh [uw1, uw2, pw, Tw];
Vh [v1, v2, q, TT];

//macro list used in the problem
macro div(u1, u2) ( dx(u1) + dy (u2) ) //
macro grad( u ) [dx(u), dy(u)] //
macro Grad(u1, u2) [grad(u1), grad(u2)] //
macro Ugradv(u1, u2, v) ([u1, u2]'* grad(v)) //
macro UgradV(u1,u2, v1, v2 ) [Ugradv(u1, u2, v1), Ugradv(u1, u2, v2)] //

problem porous( [uw1, uw2, pw, Tw] , [v1, v2, q, TT] ) =

        int2d(Th)( UgradV(uw1, uw2, u1, u2)'*[v1, v2] +
                   UgradV(u1, u2, uw1, uw2)'*[v1, v2]
                   - q*div(uw1, uw2)
                   - pw*div(v1, v2)
                   + Pr* (Grad(uw1, uw2): Grad(v1, v2))
                   + dPrDa * [uw1, uw2]'*[v1, v2]
                   - pRaPr * v2*Tw
                   +Ugradv(uw1, uw2, T)* TT
                   +Ugradv(u1, u2, Tw)* TT
                   + grad(Tw)'*grad(TT)
                   + 1e-8 * pw* q)

        - int2d(Th) ( UgradV(u1,u2,u1,u2)'* [v1, v2]

                    - q*div(u1, u2)
                    - p* div(v1, v2)
                    +Pr*(Grad(u1,u2):Grad(v1,v2))
                    + dPrDa* [u1, u2]'*[v1, v2]
                    -pRaPr*(T*v2)
                    +Ugradv(u1, u2, T)*TT
                    + grad(T)'*grad(TT)
                      + 1e-8*p*q )
                    + on(1, 2,3,4, uw1 =0, uw2 =0, pw =0 )
                    + on(2,3, Tw = 0)
                    + on(1, Tw = 1);

real err;
real tol = 1e-8;
bool Newton = true;
[u1, u2, p, T] =[0, 0, 0, 1];

for(int kkk=0; kkk<aRa.n; ++kkk)
{
for (int step = 0; step <100 ; ++step)
{

  uw1[] = u1[];
  uw2[] = u2[];
  pw[] = p[];
  Tw[] = T[];

// solving problem
  porous;

  u1[] = u1[] + uw1[];
  u2[] = u2[] + uw2[];
  p[] = p[] + pw[];
  T[] = T[] + Tw[];

//plot([u1, u2], p, T, wait =1);

  real Lu1 = u1[].l2, Lu2 = u2[].l2, Lp = p[].l2, LT = T[].l2;

  err = uw1[].l2/Lu1 + uw2[].l2/Lu2 + pw[].l2/Lp + Tw[].l2/LT;

  cout << "err value is"<< err<<endl;

  if (err < tol) break;

  if ( step >3 && err >10 ) break;

}

if(kkk==0) {
   Th=square(nn,nn,[fff(x),fff(y)]);
   // u1,u2,p,T privious mesh
   [uw1,uw2,pw,Tw]=[0,0,0,0];// // up1,up2,pp,Tp
   uw1[] = u1[];
   uw2[] = u2[];
   pw[] = p[];
   Tw[] = T[];
    // new mesh
   [u1,u2,p,T]=[0,0,0,0];// // up1,up2,pp,Tp
   u1[]=uw1[];
   u2[] = uw2[];
   p[] = pw[];
   T[] = Tw[];
 }
plot([u1,u2],T,wait=1,cmm="Ra ="+Ra);

}
```