Now let us assume that the airfoil is hot and that air is there to cool it. Much like in the previous section the heat equation for the temperature \(v\) is $$ \partial_t v -\nabla\cdot(\kappa\nabla v) + u\cdot\nabla v =0, v(t=0)=v_0, \frac{\partial v}{\partial n}|_C=0 $$ But now the domain is outside AND inside \(S\) and \(\kappa\) takes a different value in air and in steel. Furthermore there is convection of heat by the flow, hence the term \(u\cdot\nabla v\) above. Consider the following, to be plugged at the end of the previous program:

... border D(t=0,2){x=1+t;y=0;} // added to have a fine mesh at trail mesh Sh = buildmesh(C(25)+Splus(-90)+Sminus(-90)+D(200)); fespace Wh(Sh,P1); Wh v,vv; int steel=Sh(0.5,0).region, air=Sh(-1,0).region; fespace W0(Sh,P0); W0 k=0.01*(region==air)+0.1*(region==steel); W0 u1=dy(psi)*(region==air), u2=-dx(psi)*(region==air); Wh vold = 120*(region==steel); real dt=0.05, nbT=50; int i; problem thermic(v,vv,init=i,solver=LU)= int2d(Sh)(v*vv/dt + k*(dx(v) * dx(vv) + dy(v) * dy(vv)) + 10*(u1*dx(v)+u2*dy(v))*vv)- int2d(Sh)(vold*vv/dt); for(i=0;i<nbT;i++){ v=vold; thermic; plot(v); }

Notice here

- how steel and air are identified by the mesh parameter region which is defined when buildmesh is

called and takes an integer value corresponding to each connected component of \(\Omega\);

- how the convection terms are added without upwinding. Upwinding is necessary when the

Pecley number \(|u|L/\kappa\) is large (here is a typical length scale), The factor 10 in front of the convection terms is a quick way of multiplying the velocity by 10 (else it is too slow to see something).

- The solver is Gauss' LU factorization and when \tt init\(\neq 0\) the LU decomposition is reused so it

is much faster after the first iteration.