1!------------------------------------------------------------------------------- 2 3!VERS 4 5! This file is part of Code_Saturne, a general-purpose CFD tool. 6! 7! Copyright (C) 1998-2021 EDF S.A. 8! 9! This program is free software; you can redistribute it and/or modify it under 10! the terms of the GNU General Public License as published by the Free Software 11! Foundation; either version 2 of the License, or (at your option) any later 12! version. 13! 14! This program is distributed in the hope that it will be useful, but WITHOUT 15! ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS 16! FOR A PARTICULAR PURPOSE. See the GNU General Public License for more 17! details. 18! 19! You should have received a copy of the GNU General Public License along with 20! this program; if not, write to the Free Software Foundation, Inc., 51 Franklin 21! Street, Fifth Floor, Boston, MA 02110-1301, USA. 22 23!------------------------------------------------------------------------------- 24 25!=============================================================================== 26! Function: 27! --------- 28!> \file cs_user_boundary_conditions-pulverized_coal_lagrangian.f90 29!> \brief Example of cs_f_user_boundary_conditions subroutine.f90 for 30!> lagrangian pulverized coal 31! 32!------------------------------------------------------------------------------- 33 34!------------------------------------------------------------------------------- 35! Arguments 36!______________________________________________________________________________. 37! mode name role ! 38!______________________________________________________________________________! 39!> \param[in] nvar total number of variables 40!> \param[in] nscal total number of scalars 41!> \param[out] icodcl boundary condition code: 42!> - 1 Dirichlet 43!> - 2 Radiative outlet 44!> - 3 Neumann 45!> - 4 sliding and 46!> \f$ \vect{u} \cdot \vect{n} = 0 \f$ 47!> - 5 smooth wall and 48!> \f$ \vect{u} \cdot \vect{n} = 0 \f$ 49!> - 6 rough wall and 50!> \f$ \vect{u} \cdot \vect{n} = 0 \f$ 51!> - 9 free inlet/outlet 52!> (input mass flux blocked to 0) 53!> \param[in] itrifb indirection for boundary faces ordering 54!> \param[in,out] itypfb boundary face types 55!> \param[in,out] izfppp boundary face zone number 56!> \param[in] dt time step (per cell) 57!> \param[in,out] rcodcl boundary condition values: 58!> - rcodcl(1) value of the dirichlet 59!> - rcodcl(2) value of the exterior exchange 60!> coefficient (infinite if no exchange) 61!> - rcodcl(3) value flux density 62!> (negative if gain) in w/m2 63!> -# for the velocity \f$ (\mu+\mu_T) 64!> \gradt \, \vect{u} \cdot \vect{n} \f$ 65!> -# for the pressure \f$ \Delta t 66!> \grad P \cdot \vect{n} \f$ 67!> -# for a scalar \f$ cp \left( K + 68!> \dfrac{K_T}{\sigma_T} \right) 69!> \grad T \cdot \vect{n} \f$ 70!_______________________________________________________________________________ 71 72subroutine cs_f_user_boundary_conditions & 73 ( nvar , nscal , & 74 icodcl , itrifb , itypfb , izfppp , & 75 dt , & 76 rcodcl ) 77 78!=============================================================================== 79 80!=============================================================================== 81! Module files 82!=============================================================================== 83 84use paramx 85use numvar 86use optcal 87use cstphy 88use cstnum 89use entsor 90use parall 91use period 92use ppppar 93use ppthch 94use coincl 95use cpincl 96use ppincl 97use ppcpfu 98use atincl 99use ctincl 100use cs_fuel_incl 101use mesh 102use field 103 104!=============================================================================== 105 106implicit none 107 108! Arguments 109 110integer nvar , nscal 111 112integer icodcl(nfabor,nvar) 113integer itrifb(nfabor), itypfb(nfabor) 114integer izfppp(nfabor) 115 116double precision dt(ncelet) 117double precision rcodcl(nfabor,nvar,3) 118 119! Local variables 120 121!< [loc_var_dec] 122integer ifac, ii 123integer izone 124 125integer ilelt, nlelt 126 127double precision uref2, d2s3 128double precision xkent, xeent 129 130integer, allocatable, dimension(:) :: lstelt 131!< [loc_var_dec] 132 133!=============================================================================== 134 135!=============================================================================== 136! Initialization 137!=============================================================================== 138 139allocate(lstelt(nfabor)) ! temporary array for boundary faces selection 140 141d2s3 = 2.d0/3.d0 142 143!=============================================================================== 144! Assign boundary conditions to boundary faces here 145 146! For each subset: 147! - use selection criteria to filter boundary faces of a given subset 148! - loop on faces from a subset 149! - set the boundary condition for each face 150!=============================================================================== 151 152! ---- Face de type entree correspondant a une entree d'air 153! Par exemple : Air primaire , secondaire ou Air tertiaire 154 155!< [example_1] 156call getfbr('12', nlelt, lstelt) 157!========== 158 159do ilelt = 1, nlelt 160 161 ifac = lstelt(ilelt) 162 163! Type de condition aux limites pour les variables standard 164 itypfb(ifac) = ientre 165 166! Numero de zone (on les numerote de 1 a n) 167 izone = 1 168 169! - Reperage de la zone a laquelle appartient la face 170 izfppp(ifac) = izone 171 172! ------ Pour ces faces d'entree , on est a debit impose 173 174 ientat(izone) = 1 175 iqimp(izone) = 1 176! - Debit en kg/s pour l'air 177 qimpat(izone) = 1.46d-03 178! - Temperature en K pour l'air 179 timpat(izone) = 400.d0 + tkelvi 180 181 182! ------ On impose en couleur 12 une entree a debit impose 183! L'utilisateur donne donc ici uniquement 184! la direction du vecteur vitesse 185 186 rcodcl(ifac,iu,1) = 0.d0 187 rcodcl(ifac,iv,1) = 0.d0 188 rcodcl(ifac,iw,1) = 5.d0 189 190! ------ Traitement de la turbulence 191 192! La turbulence est calculee par defaut si ICALKE different de 0 193! - soit a partir du diametre hydraulique, d'une vitesse 194! de reference adaptes a l'entree courante si ICALKE = 1 195! - soit a partir du diametre hydraulique, d'une vitesse 196! de reference et de l'intensite turvulente 197! adaptes a l'entree courante si ICALKE = 2 198 199! Choix pour le calcul automatique ICALKE = 1 ou 2 200 icalke(izone) = 1 201! Saisie des donnees 202 dh(izone) = 0.1d0 203 xintur(izone) = 0.1d0 204 205 206 207! Exemple de cas ou ICALKE(IZONE) = 0 : DEBUT 208! Eliminer ces lignes pour la clarte si on a fait le choix ICALKE(IZONE) = 1 209 210 if(icalke(izone).eq.0) then 211 212! Calcul de k et epsilon en entree (XKENT et XEENT) a partir 213! l'intensite turbulente et de lois standards en conduite 214! circulaire (leur initialisation est inutile mais plus 215! propre) 216 uref2 = rcodcl(ifac,iu,1)**2 & 217 +rcodcl(ifac,iv,1)**2 & 218 +rcodcl(ifac,iw,1)**2 219 uref2 = max(uref2,1.d-12) 220 xkent = epzero 221 xeent = epzero 222 223! (ITYTUR est un indicateur qui vaut ITURB/10) 224 if (itytur.eq.2) then 225 226 rcodcl(ifac,ik,1) = xkent 227 rcodcl(ifac,iep,1) = xeent 228 229 elseif(itytur.eq.3) then 230 231 rcodcl(ifac,ir11,1) = d2s3*xkent 232 rcodcl(ifac,ir22,1) = d2s3*xkent 233 rcodcl(ifac,ir33,1) = d2s3*xkent 234 rcodcl(ifac,ir12,1) = 0.d0 235 rcodcl(ifac,ir13,1) = 0.d0 236 rcodcl(ifac,ir23,1) = 0.d0 237 rcodcl(ifac,iep,1) = xeent 238 239 elseif (iturb.eq.50) then 240 241 rcodcl(ifac,ik,1) = xkent 242 rcodcl(ifac,iep,1) = xeent 243 rcodcl(ifac,iphi,1) = d2s3 244 rcodcl(ifac,ifb,1) = 0.d0 245 246 elseif (iturb.eq.60) then 247 248 rcodcl(ifac,ik,1) = xkent 249 rcodcl(ifac,iomg,1) = xeent/cmu/xkent 250 251 elseif (iturb.eq.70) then 252 253 rcodcl(ifac,inusa,1) = cmu*xkent**2/xeent 254 255 endif 256 257 endif 258! Exemple de cas ou ICALKE(IZONE) = 0 : FIN 259 260! ------ Traitement des scalaires physiques particulieres 261! Ils sont traites automatiquement 262 263 264! ------ Traitement des scalaires utilisateurs 265 266 if ( (nscal-nscapp).gt.0 ) then 267 do ii = 1, (nscal-nscapp) 268 rcodcl(ifac,isca(ii),1) = 1.d0 269 enddo 270 endif 271 272enddo 273!< [example_1] 274 275! --- On impose en couleur 15 une paroi laterale 276 277!< [example_2] 278call getfbr('15', nlelt, lstelt) 279!========== 280 281do ilelt = 1, nlelt 282 283 ifac = lstelt(ilelt) 284 285! PAROI : DEBIT NUL (FLUX NUL POUR LA PRESSION) 286! FROTTEMENT POUR LES VITESSES (+GRANDEURS TURB) 287! FLUX NUL SUR LES SCALAIRES 288 289! Type de condition aux limites pour les variables standard 290 itypfb(ifac) = iparoi 291 292 293! Numero de zone (on les numerote de 1 a n) 294 izone = 2 295 296! - Reperage de la zone a laquelle appartient la face 297 izfppp(ifac) = izone 298 299enddo 300!< [example_2] 301 302! --- On impose en couleur 19 une sortie 303 304!< [example_3] 305call getfbr('19', nlelt, lstelt) 306!========== 307 308do ilelt = 1, nlelt 309 310 ifac = lstelt(ilelt) 311 312! SORTIE : FLUX NUL VITESSE ET TEMPERATURE, PRESSION IMPOSEE 313! Noter que la pression sera recalee a P0 314! sur la premiere face de sortie libre (ISOLIB) 315 316! Type de condition aux limites pour les variables standard 317 itypfb(ifac) = isolib 318 319! Numero de zone (on les numerote de 1 a n) 320 izone = 3 321 322! - Reperage de la zone a laquelle appartient la face 323 izfppp(ifac) = izone 324 325enddo 326!< [example_3] 327 328 329! --- On impose en couleur 14 et 4 une symetrie 330 331!< [example_4] 332call getfbr('14 or 4', nlelt, lstelt) 333!========== 334 335do ilelt = 1, nlelt 336 337 ifac = lstelt(ilelt) 338 339! SYMETRIES 340 341! Type de condition aux limites pour les variables standard 342 itypfb(ifac) = isymet 343 344! Numero de zone (on les numerote de 1 a n) 345 izone = 4 346 347! - Reperage de la zone a laquelle appartient la face 348 izfppp(ifac) = izone 349 350enddo 351!< [example_4] 352 353!-------- 354! Formats 355!-------- 356 357!---- 358! End 359!---- 360 361deallocate(lstelt) ! temporary array for boundary faces selection 362 363return 364end subroutine cs_f_user_boundary_conditions 365