Non-Dimensional Equation error

General Discussion
1

Dear @frederichecht, @prj Sir,

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

topic name

relevant problem by prj

relevant problem by frederichecht

I am trying to solve the same problem base on the suggestion and adopted the procedure given by prj Problem link and Research problem by frederichecht.

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);

}