function out=TD3_GCH2535() dp=[11:-2:1]'*1e-6; %diamètres des gouttelettes en mètres da=2.0e-3; %diamètre moyen des particules d'alumine en mètres poro=0.395; %porosité du lit d'alumine beta=2; relax=.5; % estimés initiaux A=3.55; B=1.67; C=3.55; K=3; init(2)=A; init(3)=B; init(4)=C; init(1)= K; %fraction de capture par unité d'épaisseur (1/metre) alpha=[1186.46 896.12 653.52 298.41 17.42 1.52 1171.1 972.01 683.4 322.02 20.7 1.67 626.9 445.77 199.86 21.39 6.48 1.99 191.74 77.29 17.18 6.36 5.62 2.22]; % nombre de Reynolds Re=[40.23 40.23 40.23 40.23 40.23 40.23 184.31 184.31 184.31 184.31 184.31 184.31 130.8 130.8 130.8 130.8 130.8 130.8 98.49 98.49 98.49 98.49 98.49 98.49]; % le nombre de Stokes (non effectif!) St=[6.54E-02 4.38E-02 2.65E-02 1.35E-02 4.86E-03 5.40E-04 5.23E-02 3.50E-02 2.12E-02 1.08E-02 3.89E-03 4.32E-04 1.91E-02 1.28E-02 7.75E-03 3.95E-03 1.42E-03 1.58E-04 9.71E-03 6.50E-03 3.93E-03 2.01E-03 7.22E-04 8.03E-05]; % Efficacité de capture d’une particule d’alumine isolée: somme des mécanismes % d'interception (n_int) et de diffusion thermique (n_dif) (fraction sans unité) nautre=[0.041481159 0.034056162 0.026648466 0.019274022 0.011981822 0.005084731 0.041499037 0.034075926 0.026670877 0.019300539 0.012016055 0.005144025 0.041774382 0.034380331 0.02701604 0.019708942 0.012543301 0.006057241 0.042054285 0.034689776 0.027366917 0.020124105 0.013079274 0.006985574]; Ste= St.*(39.18+2.12*(Re).^.5); fun=@(K,A,B,C,Ste,nautre) K*(1-poro)/da*(beta*Ste.^A./(B + Ste.^C)+nautre); icompt=0; critere=1000; while critere> 1e-12; f=fun(K,A,B,C,Ste,nautre) % derivee par rapport à K dfdK=-((nautre + (Ste.^A*beta)./(B + Ste.^C))*(poro - 1))/da; % derivee par rapport à A dfdA=-(K*Ste.^A*beta.*log(Ste)*(poro - 1))./(da*(B + Ste.^C)); % derivee par rapport à B dfdB=(K*Ste.^A*beta*(poro - 1))./(da*(B + Ste.^C).^2); % derivee par rapport à C dfdC=(K*Ste.^A.*Ste.^C*beta.*log(Ste)*(poro - 1))./(da*(B + Ste.^C).^2); matA=[dfdK' dfdA' dfdB' dfdC']; vectb=alpha'-f' + K*dfdK' + A*dfdA' + B*dfdB' + C*dfdC'; condmatA=cond(matA); condvectb=cond(vectb); solav=[K A B C]'; sol=matA\vectb; critere=norm(solav-sol)/norm(solav); K=K+relax*(sol(1)-K); A=A+relax*(sol(2)-A); B=B+relax*(sol(3)-B); C=C+relax*(sol(4)-C); icompt=icompt+1; epsilon(icompt)=critere; end figure semilogy(epsilon,'bx-') ylabel('{\epsilon} (-)','FontSize',16) xlabel('itérations (-)','FontSize',16) s = cellstr({'gd','ko','b^','rs'}); ss = cellstr({'g-','k-','b-','r-'}); figure posi=get(gcf,'Position') cx=posi(1); cy=posi(2); width=posi(3); length=posi(4); erre=0; for j=1:4, x=St(1+(j-1)*6:j*6); y=alpha(1+(j-1)*6:j*6); loglog(x,y,char(s(j)),'MarkerSize',12) hold on dxx=(max(x)-min(x))/1000; xx=min(x):dxx:max(x); nautre_interp=spline(St(1+(j-1)*6:j*6),nautre(1+(j-1)*6:j*6),xx); Re_interp=spline(St(1+(j-1)*6:j*6),Re(1+(j-1)*6:j*6),xx); yy=fun(K,A,B,C,xx.*(39.18+2.12*Re_interp.^.5),nautre_interp); loglog(xx,yy,char(ss(j)),'MarkerSize',12) meanx=mean(x); meany=mean(y); delta_y=abs(max(y)-min(y)); erreur=abs(fun(K,A,B,C,x.*(39.18+2.12*Re(1+(j-1)*6:j*6).^.5),nautre(1+(j-1)*6:j*6))-y); erre=erre+erreur; S(j)=(sum(erreur.^2)/(6-4))^.5; Sb(j)=S(j)/delta_y; Sy(j)=sum((y-meany).^2); r(j)=(1-sum(erreur.^2)/Sy(j))^0.5; end S=S'; Sb=Sb'; r2=(r.^2)'; dat = {colText('Cas 1', 'black'), S(1),Sb(1),r2(1);... colText('Cas 2', 'green'), S(2),Sb(2),r2(2);... colText('Cas 3', 'blue'), S(3),Sb(3),r2(3);... colText('Cas 4', 'red'), S(4),Sb(4),r2(4);}; columnname = {'Conditions', 'S','Ŝ','R²' }; columnformat = {'char', 'numeric', 'numeric','numeric'}; t = uitable('Units','normalized','Position',... [0.15 0.67 0.54 0.225], 'Data', dat,... 'ColumnName', columnname,... 'ColumnFormat', columnformat,... 'RowName',[],'FontName','Arial'); ylabel('{\alpha} (m^{-1})','FontSize',16) xlabel('Stokes (-)','FontSize',16) % Avec la fonction fitnlm de Matlab clear x x(:,1)=Ste'; x(:,2)=nautre'; modelFun = @(b,x)b(1)*(1-poro)/da*(beta*x(:,1).^b(2)./(b(3) + x(:,1).^b(4))+x(:,2)); wnlm = fitnlm(x,alpha',modelFun,init) pvalue=wnlm.Coefficients.pValue(4) b=table2array(wnlm.Coefficients(:,1)) diffsol=max((sol-b)./sol*100) end function outHtml = colText(inText, inColor) % return a HTML string with colored font %http://www.mathworks.com/matlabcentral/answers/92437-how-do-i-change-the-color-of-individual-cells-in-uitable-in-matlab-7-10-r2010a outHtml = ['', ... inText, ... '']; end