* fichier : cl_D_4.dgibi ************************************************************************ * Section : Fluides Permanent ************************************************************************ *C*********************************************************************** *C* * *C* PROJET : Opérateur CLMI * *C* NOM : cl_D_4 .dgibi * *C* DESCRIPTION : Jeu de données pour le calcul de la couche limite * *C* turbulente en utilisant les relations de fermeture de * *C* Cousteix - Cas de la plaque plane avec gradient de * *C* pression * *C* LANGAGE : Gibiane * *C* AUTEUR : Guillaume VENCO - DRN/DMT/SEMT/LTMF * *C* * *C*********************************************************************** *C* * *C* APPELES : Opérateur DETO * *C* * *C*********************************************************************** *C* * *C* * *C*********************************************************************** *C* * *C* VERSION : 02/05/2000 * *C* CREATION : 18/04/2000 * *C* * *C*********************************************************************** *vitesse en puissance m = 0. 'OPTION' dime 2 elem seg2 ; 'OPTION' trace x ; GRAPH=FAUX; *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> *>>>>>>>>>>>>>Maillage et domaine>>>>>>>>>>>>>>>>>>>>>>>>>>> *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> p1 = 0.3 0.; p2 = 0.5 0. ; nb1 = 10; *p1p2 = 'DROIT' p1 p2 'DINI' 0.001 'DFIN' 0.1 ; p1p2 = 'DROIT' p1 p2 10; mt = p1p2 ; * 'TRACER' mt ; *Définition du modèle Mmt = 'CHANGER' mt 'QUAF' ; Mail = 'CHANGER' Mmt 'POI1'; Pt0 = Mail 'ELEM' 1; $mt = 'MODELISER' Mmt 'NAVIER_STOKES' 'LINE' ; *Données Ue0 = 1.; NU = 1.d-6 ; *Valeurs initiales données par la loi de puissance D20=6.26d-4; Q10=5.64d-3; D20=8.67d-4; Q10=7.796d-3; *Pas de temps et nombre de pas DT = 1.d-1; NBIT =1; *Champ de vitesse à l'extérieur *vitesse en puissance K=5; *m=(2.); m=(1.); m1= 'DOMA' $mt 'SOMMET'; X = 'COORDONNEE' 1 m1; X1 = ((X-0.3)**m); UEX = 1.+(K*X1); UE = 'KCHT' $mt scal 'SOMMET' UEX; X2 = (X**(m-1)); DUEX = K*m*X2 ; DUE = 'KCHT' $mt scal 'SOMMET' DUEX; *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> *>>>>>>>>>>>>Création de la table RV et de la table RK>>>>>> *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> *Dans DETO le 1er coef correspond au type de fermeture *le 2ème coef correspond au n°de l'équation à traiter *1=QDM 2=Energie cinétique rv = 'EQEX' 'OMEGA' 1. 'NITER' 1 'ITMA' 0 'FIDT' 1 'OPTI' 'EF' 'IMPL' 'CENTREE' 'BDF2' 'ZONE' $mt 'OPER' 'DFDT' 1. 'D2NM' 'D2N2' DT 'INCO' 'D2' ; * 'OPER' 'DFDT' 1. 'D2NM' DT 'INCO' 'D2' ; rv = 'EQEX' rv 'OPTI' 'EF' 'IMPL' 'SUPG' 'ZONE' $mt 'OPER' 'CLMI' 4 1 UE DUE 'D2NM' 'Q1NM' 'INCO' 'D2' ; rv = 'EQEX' rv 'OPTI' 'EF' 'IMPL' 'CENTREE' 'BDF2' 'ZONE' $mt 'OPER' 'DFDT' 1. 'Q1NM' 'Q1N2' DT 'INCO' 'Q1' ; * 'OPER' 'DFDT' 1. 'D1NM' DT 'INCO' 'D1' ; rv = 'EQEX' rv 'OPTI' 'EF' 'IMPL' 'SUPG' 'ZONE' $mt 'OPER' 'CLMI' 4 3 UE DUE 'Q1NM' 'D2NM' 'INCO' 'Q1' ; rv = 'EQEX' rv 'CLIM' 'D2' 'TIMP' pt0 D20 'CLIM' 'Q1' 'TIMP' pt0 Q10; rv.inco = 'TABLE' inco ; rv.inco.'D2'= 'KCHT' $mt scal sommet D20; rv.inco.'D2NM'= 'KCHT' $mt scal sommet D20; rv.inco.'D2N2'= 'KCHT' $mt scal sommet D20; rv.inco.'Q1'= 'KCHT' $mt scal sommet Q10; rv.inco.'Q1NM'= 'KCHT' $mt scal sommet Q10; rv.inco.'Q1N2'= 'KCHT' $mt scal sommet Q10; ltps = 'PROG' 0. ; lteta = 'PROG' 0. ; TPS = 0. ; nupt = 0 ; mts1 = 'DOMA' $mt 'SOMMET' ; pt1 = 'POIN' mts1 2 ; 'LISTE' pt1 ; *'OPTION' donn 5 ; 'REPETER' BCLT nbit ; nupt=nupt+1; tps=tps+dt ; mess ' Pas de temps DT=' dt ' Temps = ' tps ; ltps=ltps et (prog tps) ; EXEC rv ; teta='EXTRAIRE' (rv.inco.'D2') 'SCAL' pt1 ; mess ' Teta=' teta ; lteta= lteta et ('PROG' teta); YNM = rv.'INCO'.'D2NM'; ZNM = rv.'INCO'.'Q1NM'; rv.'INCO'.'D2N2' = rv.'INCO'.'D2NM'; rv.'INCO'.'D2NM' = rv.'INCO'.'D2'; rv.'INCO'.'Q1N2' = rv.'INCO'.'Q1NM'; rv.'INCO'.'Q1NM' = rv.'INCO'.'Q1'; 'FIN' bclt ; li4 = 'PROG' 'EXPO' 'A' -1. 'B' 0. ltps; evt2 = 'EVOL' 'ROUGE' manu ltps li4 ; evt= evol manu ltps lteta; *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> *>>>>>>>>>>>>>>>>>>>>>>> résultats>>>>>>>>>>>>>>>>>>>>>>>>>> *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> *évolution de la vitesse extérieure et du gradient vit = 'EVOL' 'CHPO' UE mt; gvit = 'EVOL' 'ROUG' 'CHPO' DUE mt; *quantité de mouvement D2S = (rv.inco.'D2'); RD2S = UE*D2S/NU; ev1s = 'EVOL' 'ROUGE' 'CHPO' RD2S mt; *différence entre (d-d1) Q1S = (rv.inco.'Q1'); RQ1S = UE*Q1S/NU; ev2s = 'EVOL' 'ROUGE' 'CHPO' RQ1S mt; *Epaisseur de déplacement D1S = (rv.'INCO'.'D1'); RD1S = UE*D1S/NU; ev3s = 'EVOL' 'ROUGE' 'CHPO' RD1S mt; *évolution du frottement à la paroi CfS = rv.'INCO'.'CF'; ev4s = 'EVOL' 'ROUGE' 'CHPO' CfS mt; *évolution du facteur de forme H final HS = rv.'INCO'.'H'; ev5s = 'EVOL' 'ROUGE' 'CHPO' HS mt; *évolution du facteur H* HHS = Q1S*(D2S**(-1)); ev6s = 'EVOL' 'ROUGE' 'CHPO' HHS mt; *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> *solution théorique cf:SCHLICHTING >>>>>>>>>>>>>>>>>>>>>>>> *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> RX = X*UE*(1./NU); *Quantité de mvt *----Rteta donné par Nikuradze RD2N = (0.01738/1.304)*(RX**(0.861)); ev1n = 'EVOL' 'BLEU' 'CHPO' RD2N mt; *----Rteta par la loi de puissance RD2P = (1./1.e-6)*(0.036*X*(RX**(-1./5)))*UE; ev1p = 'EVOL' 'VERT' 'CHPO' RD2P mt; *différence (d-d1) théorique Q1P=(0.37-(1.2857*0.036))*(X*(RX**(-1./5))); RQ1P=UE*Q1P/NU; ev2p='EVOL' 'VERT' 'CHPO' RQ1P mt; *Epaisseur de déplacement *----Rdelta1 par la Nikuradze RD1N = 1.304*RD2N; ev3n = 'EVOL' 'BLEU' 'CHPO' RD1N mt; *----Rdelta1 par la loi de puissance RD1P = (9./7.)*RD2P; ev3p = 'EVOL' 'VERT' 'CHPO' RD1P mt; *Coefficient de frottement *----Formule par Nikuradze CfN = 0.02296*(RX**(-0.139)); ev4n = 'EVOL' 'BLEU' 'CHPO' CfN mt; *----Formule par loi de puissance CfP = 0.0576*(RX**(-1./5.)); ev4p = 'EVOL' 'VERT' 'CHPO' CfP mt; *facteur de forme *----H par Nikuradze HN = 'KCHT' $mt scal sommet 1.304; ev5n = 'EVOL' 'BLEU' 'CHPO' HN mt; *----H par loi de puissance HP = 'KCHT' $mt scal sommet 1.2857; ev5p = 'EVOL' 'VERT' 'CHPO' HP mt; *coef H* théorique HHP='KCHT' $mt scal sommet 7.89; ev6p= 'EVOL' 'BLEU' 'CHPO' HHP mt; *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> *Courbes>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> *>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> Si GRAPH; 'DESSIN' (evt) 'TITR' 'Convergence de la valeur au dernier noeud' 'TITY' 'Valeurs de Y' 'TITX' 'Temps'; 'DESSIN' (vit 'ET' gvit) 'TITRE' 'Vitesse et gradient de vitesse extérieure'; TAB1 = 'TABLE' ; TAB1.'TITRE' = 'TABLE' ; TAB1.'TITRE' . 1 = MOT ' Simulation'; TAB1.'TITRE' . 2 = MOT ' Nikuradze'; TAB1.'TITRE' . 3 = MOT ' Loi de puissance'; TAB2 = 'TABLE' ; TAB2.'TITRE' = 'TABLE' ; TAB2.'TITRE' . 1 = MOT ' Simulation'; TAB2.'TITRE' . 2 = MOT ' Loi de puissance'; 'DESSIN' (ev1s 'ET' ev1n 'ET' ev1p) 'TITRE' 'Epaisseur de quantité de mouvement' 'TITY' 'Rteta' 'TITX' 'Abscisses' 'GRIL' 'LEGE' TAB1; 'DESSIN' (ev2s 'ET' ev2p) 'TITRE' 'différence entre delta et delta1' 'TITY' 'Rdelta-Rdelta1' 'TITX' 'Abscisses' 'GRIL' 'LEGE' TAB2; 'DESSIN' (ev3s 'ET' ev3n 'ET' ev3p) 'TITRE' 'Epaisseur de déplacement' 'TITY' 'Rdelta1' 'TITX' 'Abscisses' 'GRIL' 'LEGE' TAB1; 'DESSIN' (ev4s 'ET' ev4n 'ET' ev4p) 'TITRE' 'Coefficient de frottement à la paroi' 'TITX' 'Abscisses' 'TITY' 'Cf' 'GRIL' 'LEGE' TAB1; 'DESSIN' (ev5s 'ET' ev5n 'ET' ev5p) 'TITRE' 'Evolution du facteur de forme' 'TITX' 'Abscisses' 'TITY' 'H=delta1/theta' 'LEGE' TAB1; 'DESSIN' (ev6s 'ET' ev6p) 'TITRE' 'Evolution du facteur H*' 'TITX' 'Abscisses' 'TITY' 'H' 'LEGE' TAB1; Finsi ; *calcul de l'erreur au sens de la norme L2 diag = 'DOMA' $mt 'XXDIAGSI'; r = (RD2S-RD2N)**2; r2= (RD2S-RD2P)**2; EL2N = 'SOMT'(diag*r); EL2P = 'SOMT'(diag*r2); 'MESSAGE' 'Erreur sur teta au sens de la norme L2='; 'LISTE' EL2N; 'LISTE' EL2P; r3=(RQ1S-RQ1P)**2; EL3P= 'SOMT'(diag*r3); 'MESSAGE' 'Erreur sur (delta-delta1) au sens de la norme L2='; 'LISTE' EL3P; 'FIN';