1 : * PRODT PROCEDUR JC220346 10/11/18 21:15:01 6795 2 : ************************************************************************ 3 : * NOM : PRODT 4 : * DESCRIPTION : Renvoie une mesure du tenseur gradient des vitesses. 5 : * Cette grandeur intervient notamment dans l'expression 6 : * du terme de production turbulente. 7 : ************************************************************************ 8 : * HISTORIQUE : 14/03/2000 : MAGN : création de la procédure 9 : * HISTORIQUE : 28/10/2010 : JCARDO : généralisation de PRODT 10 : * 1) ajout de 3 autres mesures pour le tenseur gradient, 11 : * en plus du taux de déformation: - taux de rotation 12 : * - tenseur complet 13 : * - expression mixte 14 : * 2) prise en compte du cas où NCO est différent de NDIM 15 : * 3) gestion de composantes de noms quelconques dans UN 16 : * 4) modification de la syntaxe et de la notice 17 : * HISTORIQUE : 18 : * HISTORIQUE : 19 : ************************************************************************ 20 : * Prière de PRENDRE LE TEMPS DE COMPLÉTER LES COMMENTAIRES 21 : * en cas de modification de ce sous-programme afin de faciliter 22 : * la maintenance ! 23 : ************************************************************************ 24 : * Syntaxe : 25 : * _________ 26 : * 27 : * CHPO1 = PRODT $MD (| 'TODEF' |) UN (GB TN) ; 28 : * | 'TOROT' | 29 : * | 'COMPL' | 30 : * | 'MIXTE' | 31 : * 32 : * 33 : ************************************************************************ 34 : 35 : 36 : DEBP PRODT ; 37 : 38 : 39 : * Dimension de l'espace physique 40 : NDIM = VALE 'DIME' ; 41 : 42 : SI (EGA NDIM 1) ; 43 : MESS '******************************************************' ; 44 : MESS ' /!\ ERREUR dans PRODT /!\ ' ; 45 : MESS ' ' ; 46 : MESS ' Cette procédure est inutile pour les problèmes 1D.' ; 47 : MESS '******************************************************' ; 48 : QUIT PRODT ; 49 : FINS ; 50 : 51 : 52 : 53 : * +====================================================================+ 54 : * | | 55 : * | L E C T U R E D E S A R G U M E N T S | 56 : * | | 57 : * +====================================================================+ 58 : 59 : * +========+ 60 : * | MODÈLE | 61 : * +========+ 62 : 63 : ARGU $MD*'MMODEL' ; 64 : 65 : 66 : * +=========+ 67 : * | MOT-CLÉ | 68 : * +=========+ 69 : 70 : ARGU MCLE/'MOT' ; 71 : SI (NON (EXIS MCLE)) ; 72 : MCLE = 'TODEF' ; 73 : FINS ; 74 : 75 : B_TODEF = (EGA MCLE 'TODEF') ; 76 : B_TOROT = (EGA MCLE 'TOROT') ; 77 : B_COMPL = (EGA MCLE 'COMPL') ; 78 : B_MIXTE = (EGA MCLE 'MIXTE') ; 79 : 80 : 81 : * +==================+ 82 : * | CHAMP DE VITESSE | 83 : * +==================+ 84 : 85 : ARGU UN*'CHPOINT' ; 86 : LCO = EXTR UN 'COMP' ; 87 : NCO = DIME LCO ; 88 : 89 : 90 : * +==============================================+ 91 : * | CHAMP DE TEMPÉRATURE ET VECTEUR FLOTTABILITÉ | 92 : * +==============================================+ 93 : 94 : ARGU TN/'CHPOINT' ; 95 : 96 : B_TN = FAUX ; 97 : SI (EXIS TN) ; 98 : 99 : ARGU GB/'POINT' ; 100 : SI (NON (EXIS GB)) ; 101 : MESS '******************************************************' ; 102 : MESS '/!\ ERREUR dans PRODT :' ; 103 : MESS ' Champ de température fourni => on attendait GB' ; 104 : MESS '******************************************************' ; 105 : ERRE 5 ; 106 : QUIT PRODT ; 107 : FINS ; 108 : 109 : B_TN = VRAI ; 110 : FINS ; 111 : 112 : 113 : 114 : 115 : * +====================================================================+ 116 : * | | 117 : * | C A L C U L D E S R É S U L T A T S | 118 : * | | 119 : * +====================================================================+ 120 : 121 : * Extraction des composantes du tenseur gradient 122 : * ============================================== 123 : 124 : UNX = EXCO UN (EXTR LCO 1) ; 125 : GDUX = KOPS 'GRADS' UNX $MD ; 126 : DUDX = EXCO GDUX 'UX' ; 127 : DUDY = EXCO GDUX 'UY' ; 128 : 129 : SI ((NCO > 1) ET (NDIM > 1)) ; 130 : UNY = EXCO UN (EXTR LCO 2) ; 131 : GDUY = KOPS 'GRADS' UNY $MD ; 132 : DVDX = EXCO GDUY 'UX' ; 133 : DVDY = EXCO GDUY 'UY' ; 134 : 135 : SI ((NCO > 2) ET (NDIM EGA 3)) ; 136 : UNZ = EXCO UN (EXTR LCO 3) ; 137 : GDUZ = KOPS 'GRADS' UNZ $MD ; 138 : DWDX = EXCO GDUZ 'UX' ; 139 : DWDY = EXCO GDUZ 'UY' ; 140 : DWDZ = EXCO GDUZ 'UZ' ; 141 : SINON ; 142 : DWDX = 0. ; 143 : DWDY = 0. ; 144 : DWDZ = 0. ; 145 : FINS ; 146 : 147 : SI (NDIM EGA 3) ; 148 : DUDZ = EXCO GDUX 'UZ' ; 149 : DVDZ = EXCO GDUY 'UZ' ; 150 : SINON ; 151 : DUDZ = 0. ; 152 : DVDZ = 0. ; 153 : FINS ; 154 : 155 : SINON ; 156 : SI (NDIM EGA 3) ; 157 : DUDZ = EXCO GDUX 'UZ' ; 158 : SINON ; 159 : DUDZ = 0. ; 160 : FINS ; 161 : DVDX = 0. ; 162 : DVDY = 0. ; 163 : DVDZ = 0. ; 164 : DWDX = 0. ; 165 : DWDY = 0. ; 166 : DWDZ = 0. ; 167 : FINS ; 168 : 169 : 170 : 171 : * Calcul du carré des normes des parties symétrique et anti-symétrique 172 : * ==================================================================== 173 : 174 : A2 = (DUDY*DUDY) + (DUDZ*DUDZ) + (DVDX*DVDX) + (DVDZ*DVDZ) 175 : + (DWDX*DWDX) + (DWDY*DWDY) ; 176 : A3 = (DUDY*DVDX) + (DUDZ*DWDX) + (DVDZ*DWDY) ; 177 : SI (NON B_TOROT) ; 178 : A1 = (DUDX*DUDX) + (DVDY*DVDY) + (DWDZ*DWDZ) ; 179 : FINS ; 180 : 181 : 182 : * Mot-clé 'TODEF': taux de déformation => 2.SijSij (par défaut) 183 : SI (B_TODEF) ; 184 : P = (2. * (A1 + A3)) + A2 ; 185 : FINS ; 186 : 187 : * Mot-clé 'TOROT': taux de rotation => 2.RijRij 188 : SI (B_TOROT) ; 189 : * /!\ Valeur absolue pour éviter les zéros "négatifs" 190 : P = ABS ( A2 + (-2. * A3) ) ; 191 : FINS ; 192 : 193 : * Mot-clé 'COMPL': tenseur gradient complet => UijUij 194 : SI (B_COMPL) ; 195 : P = A1 + A2 ; 196 : FINS ; 197 : 198 : * Mot-clé 'MIXTE': Expression proposée par Daclès-Mariani et al. 199 : SI (B_MIXTE) ; 200 : CPROD = 2. ; 201 : DRR = (ABS (A2 + (-2. * A3))) ** 0.5 ; 202 : DSS = ((2. * (A1 + A3)) + A2) ** 0.5 ; 203 : DMX = KOPS (DSS - DRR) '>|' 0. ; 204 : P = (DRR + (CPROD * DMX)) ** 2 ; 205 : FINS ; 206 : 207 : 208 : 209 : 210 : * Prise en compte de la contribution des forces de flottabilité 211 : * ============================================================= 212 : * (Production turbulente moindre si stratification stable) 213 : 214 : SI (B_TN) ; 215 : * Prandtl turbulent 216 : SGT = 0.7 ; 217 : 218 : * Calcul du produit scalaire gB.Grad(T) 219 : GDT = KOPS 'GRADS' TN $MD ; 220 : DTDX = EXCO GDT 'UX' ; 221 : DTDY = EXCO GDT 'UY' ; 222 : AX = COOR 1 GB ; 223 : AY = COOR 2 GB ; 224 : G = (DTDX*AX) + (DTDY*AY) ; 225 : SI (EGA NDIM 3) ; 226 : DTDZ = EXCO GDT 'UZ' ; 227 : AZ = COOR 3 GB ; 228 : G = G + (DTDZ*AZ) ; 229 : FINS ; 230 : 231 : * G > 0 => dT/dz < 0 => Stratification instable 232 : * G < 0 => dT/dz > 0 => Stratification stable 233 : IG = MASQ G 'INFERIEUR' 0. ; 234 : G = G * IG * (1./SGT) ; 235 : 236 : * Modélisation selon W.Rodi 237 : C3 = 0.8 ; 238 : PG = P + G ; 239 : IPG = MASQ (ABS PG) 'INFERIEUR' 1.E-30 ; 240 : RF = (-1.) * G * (INVE (PG + (IPG*1.E-30))) ; 241 : P = PG * (1. + (C3 * RF)) ; 242 : FINS ; 243 : 244 : 245 : * Empêche une éventuelle division par zéro 246 : P = P + 1.E-30 ; 247 : 248 : RESP P ; 249 : FINP ; 250 : 251 : 252 :
© Cast3M 2003 - All rights reserved.
Disclaimer