1!------------------------------------------------------------------------------- 2 3! This file is part of Code_Saturne, a general-purpose CFD tool. 4! 5! Copyright (C) 1998-2021 EDF S.A. 6! 7! This program is free software; you can redistribute it and/or modify it under 8! the terms of the GNU General Public License as published by the Free Software 9! Foundation; either version 2 of the License, or (at your option) any later 10! version. 11! 12! This program is distributed in the hope that it will be useful, but WITHOUT 13! ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS 14! FOR A PARTICULAR PURPOSE. See the GNU General Public License for more 15! details. 16! 17! You should have received a copy of the GNU General Public License along with 18! this program; if not, write to the Free Software Foundation, Inc., 51 Franklin 19! Street, Fifth Floor, Boston, MA 02110-1301, USA. 20 21!------------------------------------------------------------------------------- 22 23!=============================================================================== 24! Function : 25! -------- 26 27!> \file vericl.f90 28!> 29!> \brief Check boundary condition code. 30!> 31!------------------------------------------------------------------------------- 32 33!------------------------------------------------------------------------------- 34! Arguments 35!______________________________________________________________________________. 36! mode name role ! 37!______________________________________________________________________________! 38!> \param[in] nvar total number of variables 39!> \param[in] nscal total number of scalars 40!> \param[in,out] itypfb face boundary condition type 41!> \param[in,out] icodcl face 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!> - 13 Dirichlet for the advection operator and 54!> Neumann for the diffusion operator 55!_______________________________________________________________________________ 56 57subroutine vericl & 58 ( nvar , nscal , & 59 itypfb , icodcl ) 60 61!=============================================================================== 62 63!=============================================================================== 64! Module files 65!=============================================================================== 66 67use paramx 68use numvar 69use optcal 70use cstnum 71use cstphy 72use entsor 73use albase 74use ppppar 75use ppthch 76use ppincl 77use parall 78use mesh 79use field 80use cs_c_bindings 81 82!=============================================================================== 83 84implicit none 85 86! Arguments 87 88integer nvar, nscal 89 90integer itypfb(nfabor) 91integer icodcl(nfabor,nvar) 92 93! Local variables 94 95character(len=80) :: chaine 96 97double precision grav2 98integer ifac, ivar, icode, f_dim, f_type 99integer nstoni, nstvit, nstopp 100integer nstoke, nstosc, nstovf 101integer nstuvw, nstoup, nstuke 102integer nstrij, nsurij, nstov2 103integer nstuv2, nstokw, nstukw 104integer nstunu, nstonu 105integer nstusc 106integer f_id_rough, f_id_t_rough 107integer iis, icodcu, icodcv, icodcw, icodck, icodce 108integer icodcn 109integer icodcp, icodcf, icodca, icodom 110integer icor11, icor22, icor33, icor12, icor13, icor23 111integer ipp, iokcod, iok 112integer icodni(2), icodvi(3), icodpp(2), icodtb(8), icodsc(2) 113integer icodvf(2), icoduv(3), icodct(11), icodus(4) 114 115double precision, dimension(:), pointer :: bpro_roughness 116double precision, dimension(:), pointer :: bpro_roughness_t 117 118!=============================================================================== 119 120!=============================================================================== 121! 1. Initialisations 122!=============================================================================== 123 124! Initialize variables to avoid compiler warnings 125 126do ipp = 1, 2 127 icodni(ipp) = -1 128 icodpp(ipp) = -1 129 icodsc(ipp) = -1 130 icodvf(ipp) = -1 131enddo 132do ipp = 1, 3 133 icodvi(ipp) = -1 134 icoduv(ipp) = -1 135enddo 136do ipp = 1, 8 137 icodtb(ipp) = -1 138enddo 139do ipp = 1, 11 140 icodct(ipp) = -1 141enddo 142do ipp = 1, 4 143 icodus(ipp) = -1 144enddo 145 146call field_get_id_try("boundary_roughness", f_id_rough) 147if (f_id_rough.ge.0) call field_get_val_s(f_id_rough, bpro_roughness) 148call field_get_id_try("boundary_thermal_roughness", f_id_t_rough) 149if (f_id_t_rough.ge.0) call field_get_val_s(f_id_t_rough, bpro_roughness_t) 150 151!=============================================================================== 152! 2. VERIFICATION DE LA CONSISTANCE DES CL 153!=============================================================================== 154 155! 2.1 INITIALISATION 156! ==================== 157 158! Dans cs_user_boundary_conditions, on se donne une grande liberte 159! pour la specif. des c.l sur les variables. 160! Neanmoins pour limiter la plage des tests, on se 161! donne, pour l'instant, les contraintes suivantes : 162 163! - meme type de c.l pour les 3 composantes de vitesse 164! - pas de conditions de frottement sur la pression 165! - coherence entre les c.l vitesses et pression 166! - coherence entre les c.l vitesses et turbulence 167 168nstoni = 0 169nstosc = 0 170nstovf = 0 171nstusc = 0 172nstvit = 0 173nstopp = 0 174nstoke = 0 175nstrij = 0 176nstov2 = 0 177nstokw = 0 178nstuvw = 0 179nstoup = 0 180nstuke = 0 181nsurij = 0 182nstuv2 = 0 183nstukw = 0 184nstunu = 0 185nstonu = 0 186 187 188! 2.2 VERIFICATIONS QUE TOUTES LES CL SONT INITIALISEES 189! ====================================================== 190 191! --- Premiere boucle rapide 192iokcod = 0 193do ivar = 1, nvar 194 call field_get_type(ivarfl(ivar), f_type) 195 ! Is this field not managed by CDO? 196 if (iand(f_type, FIELD_CDO)/=FIELD_CDO) then 197 198 do ifac = 1, nfabor 199 icode = icodcl(ifac,ivar) 200 if(icode.eq. 0) then 201 iokcod = 1 202 endif 203 enddo 204 205 endif ! CDO ? 206enddo 207 208! --- Seconde boucle lente si pb plus haut 209if (iokcod.ne.0) then 210 do ivar = 1, nvar 211 call field_get_type(ivarfl(ivar), f_type) 212 ! Is this field not managed by CDO? 213 if (iand(f_type, FIELD_CDO)/=FIELD_CDO) then 214 215 do ifac = 1, nfabor 216 icode = icodcl(ifac,ivar) 217 if(icode.eq. 0) then 218 if (itypfb(ifac).gt.0) then 219 itypfb(ifac) = -itypfb(ifac) 220 endif 221 icodni(1) = ivar 222 icodni(2) = ifac 223 nstoni = nstoni + 1 224 endif 225 enddo 226 227 endif ! CDO ? 228 enddo 229endif 230 231 232! 2.3 VERIFICATIONS DE L'ADMISSIBILITE DES CONDITIONS 233! ==================================================== 234 235! --- Conditions admissibles pour les composantes de vitesse 236do ifac = 1, nfabor 237 238 icodcu = icodcl(ifac,iu) 239 icodcv = icodcl(ifac,iv) 240 icodcw = icodcl(ifac,iw) 241 242 if(icodcu.ne. 1.and.icodcu.ne. 2.and.icodcu.ne. 3.and. & 243 icodcu.ne. 4.and.icodcu.ne. 5.and.icodcu.ne. 6.and. & 244 icodcu.ne.11.and. & 245 icodcu.ne. 9.and.icodcu.ne.13.and.icodcu.ne.14) then 246 if (itypfb(ifac).gt.0) then 247 itypfb(ifac) = -itypfb(ifac) 248 endif 249 icodvi(1) = iu 250 icodvi(2) = icodcl(ifac,iu) 251 icodvi(3) = -1 252 nstvit = nstvit + 1 253 endif 254 if(icodcv.ne. 1.and.icodcv.ne. 2.and.icodcv.ne. 3.and. & 255 icodcv.ne. 4.and.icodcv.ne. 5.and.icodcv.ne. 6.and. & 256 icodcv.ne.11.and. & 257 icodcv.ne. 9.and.icodcv.ne.13.and.icodcv.ne.14) then 258 if (itypfb(ifac).gt.0) then 259 itypfb(ifac) = -itypfb(ifac) 260 endif 261 icodvi(1) = iv 262 icodvi(2) = icodcl(ifac,iv) 263 icodvi(3) = -1 264 nstvit = nstvit + 1 265 endif 266 if(icodcw.ne. 1.and.icodcw.ne. 2.and.icodcw.ne. 3.and. & 267 icodcw.ne. 4.and.icodcw.ne. 5.and.icodcw.ne. 6.and. & 268 icodcw.ne.11.and. & 269 icodcw.ne. 9.and.icodcw.ne.13.and.icodcw.ne.14) then 270 if (itypfb(ifac).gt.0) then 271 itypfb(ifac) = -itypfb(ifac) 272 endif 273 icodvi(1) = iw 274 icodvi(2) = icodcl(ifac,iw) 275 icodvi(3) = -1 276 nstvit = nstvit + 1 277 endif 278 279 ! --- Check roughness if rough wall function 280 if (icodcu.eq.6) then 281 282 iok = 0 283 if ((f_id_rough.eq.-1)) then 284 iok = 1 285 else if (bpro_roughness(ifac).le.0.d0) then 286 iok = 1 287 endif 288 if (iok .ge. 1) then 289 if (itypfb(ifac).gt.0) then 290 itypfb(ifac) = -itypfb(ifac) 291 endif 292 icodvi(1) = -6 293 icodvi(2) = icodcl(ifac,iw) 294 icodvi(3) = -1 295 nstvit = nstvit + 1 296 endif 297 endif 298 299enddo 300 301! --- Conditions admissibles pour la pression 302do ifac = 1, nfabor 303 304 if (icodcl(ifac,ipr).ne.1 .and. icodcl(ifac,ipr).ne.2 .and. & 305 icodcl(ifac,ipr).ne.3 .and. icodcl(ifac,ipr).ne.11.and. & 306 icodcl(ifac,ipr).ne.12.and. icodcl(ifac,ipr).ne.13.and. & 307 icodcl(ifac,ipr).ne.15) then 308 if (itypfb(ifac).gt.0) then 309 itypfb(ifac) = -itypfb(ifac) 310 endif 311 icodpp(1) = ipr 312 icodpp(2) = icodcl(ifac,ipr) 313 nstopp = nstopp + 1 314 endif 315 316enddo 317 318! --- Conditions admissibles pour k et epsilon 319if (itytur.eq.2) then 320 321 do ifac = 1, nfabor 322 323 if((icodcl(ifac,ik ).ne. 1.and. & 324 icodcl(ifac,ik ).ne. 2.and. & 325 icodcl(ifac,ik ).ne. 3.and. & 326 icodcl(ifac,ik ).ne. 5.and. & 327 icodcl(ifac,ik ).ne.13.and. & 328 icodcl(ifac,ik ).ne. 6 ).or. & 329 (icodcl(ifac,iep).ne. 1.and. & 330 icodcl(ifac,iep).ne. 2.and. & 331 icodcl(ifac,iep).ne. 3.and. & 332 icodcl(ifac,iep).ne. 5.and. & 333 icodcl(ifac,iep).ne.13.and. & 334 icodcl(ifac,iep).ne. 6 ) )then 335 if (itypfb(ifac).gt.0) then 336 itypfb(ifac) = -itypfb(ifac) 337 endif 338 icodtb(1) = ik 339 icodtb(2) = icodcl(ifac,ik) 340 icodtb(3) = iep 341 icodtb(4) = icodcl(ifac,iep) 342 nstoke = nstoke + 1 343 endif 344 345 enddo 346 347 ! --- Conditions admissibles pour Rij et epsilon 348elseif(itytur.eq.3) then 349 350 ivar = ir11 351 do ifac = 1, nfabor 352 icode = icodcl(ifac,ivar) 353 if (icode.ne. 1.and.icode.ne. 2.and. icode.ne. 3.and. & 354 icode.ne. 4.and.icode.ne. 5.and.icode.ne. 6 .and.icode.ne.13) then 355 if (itypfb(ifac).gt.0) then 356 itypfb(ifac) = -itypfb(ifac) 357 endif 358 icodtb(1) = ivar 359 icodtb(2) = icode 360 nstrij = nstrij + 1 361 endif 362 enddo 363 364 ivar = ir22 365 do ifac = 1, nfabor 366 icode = icodcl(ifac,ivar) 367 if (icode.ne. 1.and.icode.ne. 2.and. icode.ne. 3.and. & 368 icode.ne. 4.and.icode.ne. 5.and.icode.ne. 6 .and.icode.ne.13) then 369 if (itypfb(ifac).gt.0) then 370 itypfb(ifac) = -itypfb(ifac) 371 endif 372 icodtb(1) = ivar 373 icodtb(2) = icode 374 nstrij = nstrij + 1 375 endif 376 enddo 377 378 ivar = ir33 379 do ifac = 1, nfabor 380 icode = icodcl(ifac,ivar) 381 if (icode.ne. 1.and.icode.ne. 2.and. icode.ne. 3.and. & 382 icode.ne. 4.and.icode.ne. 5.and.icode.ne. 6 .and.icode.ne.13) then 383 if (itypfb(ifac).gt.0) then 384 itypfb(ifac) = -itypfb(ifac) 385 endif 386 icodtb(1) = ivar 387 icodtb(2) = icode 388 nstrij = nstrij + 1 389 endif 390 enddo 391 392 ivar = ir12 393 do ifac = 1, nfabor 394 icode = icodcl(ifac,ivar) 395 if (icode.ne. 1.and.icode.ne. 2.and. icode.ne. 3.and. & 396 icode.ne. 4.and.icode.ne. 5.and.icode.ne. 6 .and.icode.ne.13) then 397 if (itypfb(ifac).gt.0) then 398 itypfb(ifac) = -itypfb(ifac) 399 endif 400 icodtb(1) = ivar 401 icodtb(2) = icode 402 nstrij = nstrij + 1 403 endif 404 enddo 405 406 ivar = ir23 407 do ifac = 1, nfabor 408 icode = icodcl(ifac,ivar) 409 if (icode.ne. 1.and.icode.ne. 2.and. icode.ne. 3.and. & 410 icode.ne. 4.and.icode.ne. 5.and.icode.ne. 6 .and.icode.ne.13) then 411 if (itypfb(ifac).gt.0) then 412 itypfb(ifac) = -itypfb(ifac) 413 endif 414 icodtb(1) = ivar 415 icodtb(2) = icode 416 nstrij = nstrij + 1 417 endif 418 enddo 419 420 ivar = ir13 421 do ifac = 1, nfabor 422 icode = icodcl(ifac,ivar) 423 if (icode.ne. 1.and.icode.ne. 2.and. icode.ne. 3.and. & 424 icode.ne. 4.and.icode.ne. 5.and.icode.ne. 6 .and.icode.ne.13) then 425 if (itypfb(ifac).gt.0) then 426 itypfb(ifac) = -itypfb(ifac) 427 endif 428 icodtb(1) = ivar 429 icodtb(2) = icode 430 nstrij = nstrij + 1 431 endif 432 enddo 433 434 do ifac = 1, nfabor 435 icode = icodcl(ifac,iep) 436 if (icode.ne. 1.and.icode.ne. 2.and. icode.ne. 3.and. & 437 icode.ne. 5.and.icode.ne. 6 .and.icode.ne.13) then 438 if (itypfb(ifac).gt.0) then 439 itypfb(ifac) = -itypfb(ifac) 440 endif 441 icodtb(1) = iep 442 icodtb(2) = icode 443 icodtb(3) = -1 444 nstrij = nstrij + 1 445 endif 446 enddo 447 448 if (iturb.eq.32) then 449 do ifac = 1, nfabor 450 icode = icodcl(ifac,ial) 451 if (icode.ne. 1.and.icode.ne. 2.and. icode.ne. 3.and. & 452 icode.ne. 5.and.icode.ne. 6 .and.icode.ne.13) then 453 if (itypfb(ifac).gt.0) then 454 itypfb(ifac) = -itypfb(ifac) 455 endif 456 icodtb(1) = ial 457 icodtb(2) = icode 458 icodtb(3) = -1 459 nstrij = nstrij + 1 460 endif 461 ! No rough wall with EBRSM 462 do ivar = 1, nvar 463 if (icodcl(ifac,ivar).eq.6) then 464 if (itypfb(ifac).gt.0) then 465 itypfb(ifac) = -itypfb(ifac) 466 endif 467 icodtb(1) = ial 468 icodtb(2) = icode 469 icodtb(3) = ivar 470 nstrij = nstrij + 1 471 endif 472 enddo 473 enddo 474 endif 475 476 ! --- Conditions admissibles pour k, epsilon, phi et f_barre 477elseif (iturb.eq.50) then 478 479 do ifac = 1, nfabor 480 481 if((icodcl(ifac,ik ).ne. 1.and. & 482 icodcl(ifac,ik ).ne. 2.and. & 483 icodcl(ifac,ik ).ne. 3.and. & 484 icodcl(ifac,ik ).ne. 5.and. & 485 icodcl(ifac,ik ).ne.13.and. & 486 icodcl(ifac,ik ).ne. 6 ).or. & 487 (icodcl(ifac,iep).ne. 1.and. & 488 icodcl(ifac,iep).ne. 2.and. & 489 icodcl(ifac,iep).ne. 3.and. & 490 icodcl(ifac,iep).ne. 5.and. & 491 icodcl(ifac,iep).ne. 6 ).or. & 492 (icodcl(ifac,iphi).ne. 1.and. & 493 icodcl(ifac,iphi).ne. 2.and. & 494 icodcl(ifac,iphi).ne. 3.and. & 495 icodcl(ifac,iphi).ne. 5.and. & 496 icodcl(ifac,iphi).ne.13.and. & 497 icodcl(ifac,iphi).ne. 6 ).or. & 498 (icodcl(ifac,ifb).ne. 1.and. & 499 icodcl(ifac,ifb).ne. 2.and. & 500 icodcl(ifac,ifb).ne. 3.and. & 501 icodcl(ifac,ifb).ne. 5.and. & 502 icodcl(ifac,ifb).ne.13.and. & 503 icodcl(ifac,ifb).ne. 6 ) )then 504 505 if (itypfb(ifac).gt.0) then 506 itypfb(ifac) = -itypfb(ifac) 507 endif 508 icodtb(1) = ik 509 icodtb(2) = icodcl(ifac,ik) 510 icodtb(3) = iep 511 icodtb(4) = icodcl(ifac,iep) 512 icodtb(5) = iphi 513 icodtb(6) = icodcl(ifac,iphi) 514 icodtb(7) = ifb 515 icodtb(8) = icodcl(ifac,ifb) 516 nstov2 = nstov2 + 1 517 518 endif 519 520 enddo 521 522! --- Conditions admissibles pour k, epsilon, phi et alpha 523elseif (iturb.eq.51) then 524 525 do ifac = 1, nfabor 526 527 if((icodcl(ifac,ik ).ne. 1.and. & 528 icodcl(ifac,ik ).ne. 2.and. & 529 icodcl(ifac,ik ).ne. 3.and. & 530 icodcl(ifac,ik ).ne. 5.and. & 531 icodcl(ifac,ik ).ne.13.and. & 532 icodcl(ifac,ik ).ne. 6 ).or. & 533 (icodcl(ifac,iep).ne. 1.and. & 534 icodcl(ifac,iep).ne. 2.and. & 535 icodcl(ifac,iep).ne. 3.and. & 536 icodcl(ifac,iep).ne. 5.and. & 537 icodcl(ifac,iep).ne.13.and. & 538 icodcl(ifac,iep).ne. 6 ).or. & 539 (icodcl(ifac,iphi).ne. 1.and. & 540 icodcl(ifac,iphi).ne. 2.and. & 541 icodcl(ifac,iphi).ne. 3.and. & 542 icodcl(ifac,iphi).ne. 5.and. & 543 icodcl(ifac,iphi).ne.13.and. & 544 icodcl(ifac,iphi).ne. 6 ).or. & 545 (icodcl(ifac,ial).ne. 1.and. & 546 icodcl(ifac,ial).ne. 2.and. & 547 icodcl(ifac,ial).ne. 3.and. & 548 icodcl(ifac,ial).ne. 5.and. & 549 icodcl(ifac,ial).ne.13.and. & 550 icodcl(ifac,ial).ne. 6 ) )then 551 552 if (itypfb(ifac).gt.0) then 553 itypfb(ifac) = -itypfb(ifac) 554 endif 555 icodtb(1) = ik 556 icodtb(2) = icodcl(ifac,ik) 557 icodtb(3) = iep 558 icodtb(4) = icodcl(ifac,iep) 559 icodtb(5) = iphi 560 icodtb(6) = icodcl(ifac,iphi) 561 icodtb(7) = ial 562 icodtb(8) = icodcl(ifac,ial) 563 nstov2 = nstov2 + 1 564 565 endif 566 567 enddo 568 569! --- Conditions admissibles pour k et omega 570elseif (iturb.eq.60) then 571 572 do ifac = 1, nfabor 573 574 if((icodcl(ifac,ik ).ne. 1.and. & 575 icodcl(ifac,ik ).ne. 2.and. & 576 icodcl(ifac,ik ).ne. 3.and. & 577 icodcl(ifac,ik ).ne. 5.and. & 578 icodcl(ifac,ik ).ne.13.and. & 579 icodcl(ifac,ik ).ne. 6 ).or. & 580 (icodcl(ifac,iomg).ne. 1.and. & 581 icodcl(ifac,iomg).ne. 2.and. & 582 icodcl(ifac,iomg).ne. 3.and. & 583 icodcl(ifac,iomg).ne. 5.and. & 584 icodcl(ifac,iomg).ne.13.and. & 585 icodcl(ifac,iomg).ne. 6 ) )then 586 587 if (itypfb(ifac).gt.0) then 588 itypfb(ifac) = -itypfb(ifac) 589 endif 590 icodtb(1) = ik 591 icodtb(2) = icodcl(ifac,ik) 592 icodtb(3) = iomg 593 icodtb(4) = icodcl(ifac,iomg) 594 nstokw = nstokw + 1 595 596 endif 597 598 enddo 599 600 ! --- Conditions admissibles pour Spalart-Allmaras 601elseif (iturb.eq.70) then 602 603 do ifac = 1, nfabor 604 605 if(icodcl(ifac,inusa ).ne. 1.and. & 606 icodcl(ifac,inusa ).ne. 2.and. & 607 icodcl(ifac,inusa ).ne. 3.and. & 608 icodcl(ifac,inusa ).ne. 5.and. & 609 icodcl(ifac,inusa ).ne.13.and. & 610 icodcl(ifac,inusa ).ne. 6 )then 611 612 if (itypfb(ifac).gt.0) then 613 itypfb(ifac) = -itypfb(ifac) 614 endif 615 icodtb(1) = inusa 616 icodtb(2) = icodcl(ifac,inusa) 617 nstonu = nstonu + 1 618 619 endif 620 621 enddo 622 623endif 624 625! Admissible conditions for additional transporter variables 626! (scalars or vectors) 627if (nscal.ge.1) then 628 do iis = 1,nscal 629 ivar = isca(iis) 630 call field_get_dim(ivarfl(ivar), f_dim) 631 do ifac = 1, nfabor 632 if ((icodcl(ifac,ivar).ne. 1.and. & 633 icodcl(ifac,ivar).ne. 2.and. & 634 icodcl(ifac,ivar).ne. 3.and. & 635 icodcl(ifac,ivar).ne. 4.and. & 636 icodcl(ifac,ivar).ne. 5.and. & 637 icodcl(ifac,ivar).ne. 6.and. & 638 icodcl(ifac,ivar).ne.11.and. & 639 icodcl(ifac,ivar).ne.12.and. & 640 icodcl(ifac,ivar).ne.13.and. & 641 icodcl(ifac,ivar).ne.14.and. & 642 icodcl(ifac,ivar).ne.15) & 643 ! Only for scalars 644 .or. ( icodcl(ifac,ivar).eq.12.and. f_dim.gt.1) & 645 ! Only for vectors 646 .or. ( icodcl(ifac,ivar).eq. 4.and. f_dim.eq.1) & 647 .or. ( icodcl(ifac,ivar).eq.11.and. f_dim.eq.1) & 648 .or. ( icodcl(ifac,ivar).eq.14.and. f_dim.eq.1)) then 649 if (itypfb(ifac).gt.0) then 650 itypfb(ifac) = -itypfb(ifac) 651 endif 652 icodsc(1) = ivar 653 icodsc(2) = icodcl(ifac,ivar) 654 nstosc = nstosc + 1 655 endif 656 if((icodcl(ifac,ivar).eq. 5.or.icodcl(ifac,ivar).eq.15).and. & 657 iscavr(iis).gt.0 ) then 658 if (itypfb(ifac).gt.0) then 659 itypfb(ifac) = -itypfb(ifac) 660 endif 661 icodvf(1) = ivar 662 icodvf(2) = icodcl(ifac,ivar) 663 nstovf = nstovf + 1 664 endif 665 if(icodcl(ifac,ivar).eq. 6.and. & 666 iscavr(iis).gt.0 ) then 667 if (itypfb(ifac).gt.0) then 668 itypfb(ifac) = -itypfb(ifac) 669 endif 670 icodvf(1) = ivar 671 icodvf(2) = icodcl(ifac,ivar) 672 nstovf = nstovf + 1 673 endif 674 675 ! --- Check roughness if rough wall function 676 if(icodcl(ifac,ivar).eq.6)then 677 678 iok = 0 679 if ((f_id_t_rough.eq.-1)) then 680 iok = 1 681 else if (bpro_roughness_t(ifac).le.0.d0) then 682 iok = 1 683 endif 684 if (iok .ge. 1) then 685 686 if (itypfb(ifac).gt.0) then 687 itypfb(ifac) = -itypfb(ifac) 688 endif 689 icodsc(1) = -6 690 icodsc(2) = icodcl(ifac,ivar) 691 nstosc = nstosc + 1 692 endif 693 endif 694 enddo 695 enddo 696endif 697 698! --- Check if the gravity is non zero in case of free-surface 699iok = 0 700if (iale.ge.1) then 701 grav2 = gx**2 + gy**2 + gz**2 702 do ifac = 1, nfabor 703 if (ialtyb(ifac).eq.ifresf.and.grav2.le.epzero**2) then 704 iok = 1 705 endif 706 enddo 707endif 708 709if (irangp.ge.0) call parcmx(iok) 710 711if (iok.ne.0) then 712 write(nfecra,2001) 713 call csexit (1) 714 !========== 715endif 716 717! 2.4 VERIFICATIONS DES COHERENCES INTER VARIABLES 718! ================================================ 719 720! --- Coherence pour les composantes de vitesse 721do ifac = 1, nfabor 722 723 icodcu = icodcl(ifac,iu) 724 icodcv = icodcl(ifac,iv) 725 icodcw = icodcl(ifac,iw) 726 727 if(icodcu.eq.4.or.icodcu.eq.5.or.icodcu.eq.6.or. & 728 icodcu.eq.9.or. & 729 icodcv.eq.4.or.icodcv.eq.5.or.icodcv.eq.6.or. & 730 icodcv.eq.9.or. & 731 icodcw.eq.4.or.icodcw.eq.5.or.icodcw.eq.6.or. & 732 icodcw.eq.9 )then 733 734 if (icodcu.ne.icodcv .or. icodcu.ne.icodcw .or. & 735 icodcv.ne.icodcw ) then 736 if (itypfb(ifac).gt.0) then 737 itypfb(ifac) = -itypfb(ifac) 738 endif 739 icoduv(1) = icodcu 740 icoduv(2) = icodcv 741 icoduv(3) = icodcw 742 nstuvw = nstuvw + 1 743 endif 744 endif 745 746 ! --- Coherence vitesse pression 747 748 ! Remarques : 749 ! Pas de regle stricte de coherence vitesse/pression. 750 ! Avant on imposait un Dirichlet sur la pression pour en 751 ! entree/sortie, mais cela ne semble pas imperatif. Le test 752 ! est laisse en commentaire pour etre recupere si necessaire. 753 754 ! if (icodcu.eq.9 .or. icodcv.eq.9 .or. icodcw.eq.9) THEN 755 ! if (icodcl(ifac,ipriph).ne.1) then 756 ! nstoup = nstoup + 1 757 ! endif 758 ! endif 759 760enddo 761 762! --- Coherence vitesse turbulence 763 764if(itytur.eq.2) then 765 766 do ifac = 1, nfabor 767 768 icodcu = icodcl(ifac,iu) 769 icodcv = icodcl(ifac,iv) 770 icodcw = icodcl(ifac,iw) 771 icodck = icodcl(ifac,ik) 772 icodce = icodcl(ifac,iep) 773 774 if( (icodcu.eq.5 .or. icodcv.eq.5 .or. icodcw.eq.5 .or. & 775 icodck.eq.5 .or. icodce.eq.5) .and. & 776 (icodcu.ne.5 .or. icodcv.ne.5 .or. icodcw.ne.5 .or. & 777 icodck.ne.5 .or. icodce.ne.5) ) then 778 if (itypfb(ifac).gt.0) then 779 itypfb(ifac) = -itypfb(ifac) 780 endif 781 icodct(1) = ik 782 icodct(2) = icodcl(ifac,ik) 783 icodct(3) = iep 784 icodct(4) = icodcl(ifac,iep) 785 icodct(5) = icodcu 786 icodct(6) = icodcv 787 icodct(7) = icodcw 788 nstuke = nstuke + 1 789 endif 790 791 if( (icodcu.eq.6 .or. icodcv.eq.6 .or. icodcw.eq.6 .or. & 792 icodck.eq.6 .or. icodce.eq.6) .and. & 793 (icodcu.ne.6 .or. icodcv.ne.6 .or. icodcw.ne.6 .or. & 794 icodck.ne.6 .or. icodce.ne.6) ) then 795 if (itypfb(ifac).gt.0) then 796 itypfb(ifac) = -itypfb(ifac) 797 endif 798 icodct(1) = ik 799 icodct(2) = icodcl(ifac,ik) 800 icodct(3) = iep 801 icodct(4) = icodcl(ifac,iep) 802 icodct(5) = icodcu 803 icodct(6) = icodcv 804 icodct(7) = icodcw 805 nstuke = nstuke + 1 806 endif 807 808 enddo 809 810elseif(iturb.eq.30.or.iturb.eq.31) then 811 812 do ifac = 1, nfabor 813 814 icodcu = icodcl(ifac,iu) 815 icodcv = icodcl(ifac,iv) 816 icodcw = icodcl(ifac,iw) 817 icor11 = icodcl(ifac,ir11) 818 icor22 = icodcl(ifac,ir22) 819 icor33 = icodcl(ifac,ir33) 820 icor12 = icodcl(ifac,ir12) 821 icor13 = icodcl(ifac,ir13) 822 icor23 = icodcl(ifac,ir23) 823 icodce = icodcl(ifac,iep) 824 825 if( (icodcu.eq.5 .or. icodcv.eq.5 .or. icodcw.eq.5 .or. & 826 icor11.eq.5 .or. icor22.eq.5 .or. & 827 icor33.eq.5 .or. icor12.eq.5 .or. & 828 icor13.eq.5 .or. icor23.eq.5 .or. & 829 icodce.eq.5 ) .and. & 830 (icodcu.ne.5 .or. icodcv.ne.5 .or. icodcw.ne.5 .or. & 831 icor11.ne.5 .or. icor22.ne.5 .or. & 832 icor33.ne.5 .or. icor12.ne.5 .or. & 833 icor13.ne.5 .or. icor23.ne.5 .or. & 834 icodce.ne.5 ) ) then 835 if (itypfb(ifac).gt.0) then 836 itypfb(ifac) = -itypfb(ifac) 837 endif 838 icodct(1) = icor11 839 icodct(2) = icor22 840 icodct(3) = icor33 841 icodct(4) = icor12 842 icodct(5) = icor13 843 icodct(6) = icor23 844 icodct(7) = icodce 845 icodct(8) = icodcu 846 icodct(9) = icodcv 847 icodct(10) = icodcw 848 nsurij = nsurij + 1 849 endif 850 851 if( (icodcu.eq.6 .or. icodcv.eq.6 .or. icodcw.eq.6 .or. & 852 icor11.eq.6 .or. icor22.eq.6 .or. & 853 icor33.eq.6 .or. icor12.eq.6 .or. & 854 icor13.eq.6 .or. icor23.eq.6 .or. & 855 icodce.eq.6 ) .and. & 856 (icodcu.ne.6 .or. icodcv.ne.6 .or. icodcw.ne.6 .or. & 857 icor11.ne.6 .or. icor22.ne.6 .or. & 858 icor33.ne.6 .or. icor12.ne.6 .or. & 859 icor13.ne.6 .or. icor23.ne.6 .or. & 860 icodce.ne.6 ) ) then 861 if (itypfb(ifac).gt.0) then 862 itypfb(ifac) = -itypfb(ifac) 863 endif 864 icodct(1) = icor11 865 icodct(2) = icor22 866 icodct(3) = icor33 867 icodct(4) = icor12 868 icodct(5) = icor13 869 icodct(6) = icor23 870 icodct(7) = icodce 871 icodct(8) = icodcu 872 icodct(9) = icodcv 873 icodct(10) = icodcw 874 nsurij = nsurij + 1 875 endif 876 877 if( (icodcu.eq.4 .or. icodcv.eq.4 .or. icodcw.eq.4 .or. & 878 icor11.eq.4 .or. icor22.eq.4 .or. & 879 icor33.eq.4 .or. icor12.eq.4 .or. & 880 icor13.eq.4 .or. icor23.eq.4 & 881 ) .and. & 882 (icodcu.ne.4 .or. icodcv.ne.4 .or. icodcw.ne.4 .or. & 883 icor11.ne.4 .or. icor22.ne.4 .or. & 884 icor33.ne.4 .or. icor12.ne.4 .or. & 885 icor13.ne.4 .or. icor23.ne.4 .or. & 886 icodce.ne.3) ) then 887 if (itypfb(ifac).gt.0) then 888 itypfb(ifac) = -itypfb(ifac) 889 endif 890 icodct(1) = icor11 891 icodct(2) = icor22 892 icodct(3) = icor33 893 icodct(4) = icor12 894 icodct(5) = icor13 895 icodct(6) = icor23 896 icodct(7) = icodce 897 icodct(8) = icodcu 898 icodct(9) = icodcv 899 icodct(10) = icodcw 900 nsurij = nsurij + 1 901 endif 902 903 enddo 904 905elseif (iturb.eq.32) then 906 do ifac = 1, nfabor 907 908 icodcu = icodcl(ifac,iu) 909 icodcv = icodcl(ifac,iv) 910 icodcw = icodcl(ifac,iw) 911 icor11 = icodcl(ifac,ir11) 912 icor22 = icodcl(ifac,ir22) 913 icor33 = icodcl(ifac,ir33) 914 icor12 = icodcl(ifac,ir12) 915 icor13 = icodcl(ifac,ir13) 916 icor23 = icodcl(ifac,ir23) 917 icodce = icodcl(ifac,iep) 918 icodca = icodcl(ifac,ial) 919 920 if ( (icodcu.eq.5 .or. icodcv.eq.5 .or. icodcw.eq.5 .or. & 921 icor11.eq.5 .or. icor22.eq.5 .or. & 922 icor33.eq.5 .or. icor12.eq.5 .or. & 923 icor13.eq.5 .or. icor23.eq.5 .or. & 924 icodce.eq.5 .or. icodca.eq.5 ) .and. & 925 (icodcu.ne.5 .or. icodcv.ne.5 .or. icodcw.ne.5 .or. & 926 icor11.ne.5 .or. icor22.ne.5 .or. & 927 icor33.ne.5 .or. icor12.ne.5 .or. & 928 icor13.ne.5 .or. icor23.ne.5 .or. & 929 icodce.ne.5 .or. icodca.ne.5 ) ) then 930 if (itypfb(ifac).gt.0) then 931 itypfb(ifac) = -itypfb(ifac) 932 endif 933 icodct(1) = icor11 934 icodct(2) = icor22 935 icodct(3) = icor33 936 icodct(4) = icor12 937 icodct(5) = icor13 938 icodct(6) = icor23 939 icodct(7) = icodce 940 icodct(8) = icodcu 941 icodct(9) = icodca 942 icodct(10) = icodcv 943 icodct(11) = icodcw 944 nsurij = nsurij + 1 945 endif 946 947 if ( (icodcu.eq.6 .or. icodcv.eq.6 .or. icodcw.eq.6 .or. & 948 icor11.eq.6 .or. icor22.eq.6 .or. & 949 icor33.eq.6 .or. icor12.eq.6 .or. & 950 icor13.eq.6 .or. icor23.eq.6 .or. & 951 icodce.eq.6 .or. icodca.eq.6 ) .and. & 952 (icodcu.ne.6 .or. icodcv.ne.6 .or. icodcw.ne.6 .or. & 953 icor11.ne.6 .or. icor22.ne.6 .or. & 954 icor33.ne.6 .or. icor12.ne.6 .or. & 955 icor13.ne.6 .or. icor23.ne.6 .or. & 956 icodce.ne.6 .or. icodca.ne.6 ) ) then 957 if (itypfb(ifac).gt.0) then 958 itypfb(ifac) = -itypfb(ifac) 959 endif 960 icodct(1) = icor11 961 icodct(2) = icor22 962 icodct(3) = icor33 963 icodct(4) = icor12 964 icodct(5) = icor13 965 icodct(6) = icor23 966 icodct(7) = icodce 967 icodct(8) = icodcu 968 icodct(9) = icodca 969 icodct(10) = icodcv 970 icodct(11) = icodcw 971 nsurij = nsurij + 1 972 endif 973 974 if ( (icodcu.eq.4 .or. icodcv.eq.4 .or. icodcw.eq.4 .or. & 975 icor11.eq.4 .or. icor22.eq.4 .or. & 976 icor33.eq.4 .or. icor12.eq.4 .or. & 977 icor13.eq.4 .or. icor23.eq.4 & 978 ) .and. & 979 (icodcu.ne.4 .or. icodcv.ne.4 .or. icodcw.ne.4 .or. & 980 icor11.ne.4 .or. icor22.ne.4 .or. & 981 icor33.ne.4 .or. icor12.ne.4 .or. & 982 icor13.ne.4 .or. icor23.ne.4 .or. & 983 icodce.ne.3 & 984 .or.icodca.ne.3 ) ) then 985 if (itypfb(ifac).gt.0) then 986 itypfb(ifac) = -itypfb(ifac) 987 endif 988 icodct(1) = icor11 989 icodct(2) = icor22 990 icodct(3) = icor33 991 icodct(4) = icor12 992 icodct(5) = icor13 993 icodct(6) = icor23 994 icodct(7) = icodce 995 icodct(8) = icodcu 996 icodct(9) = icodca 997 icodct(10) = icodcv 998 icodct(11) = icodcw 999 nsurij = nsurij + 1 1000 endif 1001 1002 enddo 1003 1004elseif(iturb.eq.50 ) then 1005 1006 do ifac = 1, nfabor 1007 1008 icodcu = icodcl(ifac,iu) 1009 icodcv = icodcl(ifac,iv) 1010 icodcw = icodcl(ifac,iw) 1011 icodck = icodcl(ifac,ik) 1012 icodce = icodcl(ifac,iep) 1013 icodcp = icodcl(ifac,iphi) 1014 icodcf = icodcl(ifac,ifb) 1015 1016 if( (icodcu.eq.5 .or. icodcv.eq.5 .or. icodcw.eq.5 .or. & 1017 icodck.eq.5 .or. icodce.eq.5 .or. icodcp.eq.5 .or. & 1018 icodcf.eq.5 ) .and. & 1019 (icodcu.ne.5 .or. icodcv.ne.5 .or. icodcw.ne.5 .or. & 1020 icodck.ne.5 .or. icodce.ne.5 .or. icodcp.ne.5 .or. & 1021 icodcf.ne.5 ) ) then 1022 1023 if (itypfb(ifac).gt.0) then 1024 itypfb(ifac) = -itypfb(ifac) 1025 endif 1026 icodct(1) = ik 1027 icodct(2) = icodcl(ifac,ik) 1028 icodct(3) = iep 1029 icodct(4) = icodcl(ifac,iep) 1030 icodct(5) = iphi 1031 icodct(6) = icodcl(ifac,iphi) 1032 icodct(7) = ifb 1033 icodct(8) = icodcl(ifac,ifb) 1034 icodct(9) = icodcu 1035 icodct(10) = icodcv 1036 icodct(11) = icodcw 1037 nstuv2 = nstuv2 + 1 1038 1039 if( (icodcu.eq.6 .or. icodcv.eq.6 .or. icodcw.eq.6 .or. & 1040 icodck.eq.6 .or. icodce.eq.6 .or. icodcp.eq.6 .or. & 1041 icodcf.eq.6 ) .and. & 1042 (icodcu.ne.6 .or. icodcv.ne.6 .or. icodcw.ne.6 .or. & 1043 icodck.ne.6 .or. icodce.ne.6 .or. icodcp.ne.6 .or. & 1044 icodcf.ne.6 ) ) then 1045 1046 if (itypfb(ifac).gt.0) then 1047 itypfb(ifac) = -itypfb(ifac) 1048 endif 1049 icodct(1) = ik 1050 icodct(2) = icodcl(ifac,ik) 1051 icodct(3) = iep 1052 icodct(4) = icodcl(ifac,iep) 1053 icodct(5) = iphi 1054 icodct(6) = icodcl(ifac,iphi) 1055 icodct(7) = ifb 1056 icodct(8) = icodcl(ifac,ifb) 1057 icodct(9) = icodcu 1058 icodct(10) = icodcv 1059 icodct(11) = icodcw 1060 nstuv2 = nstuv2 + 1 1061 1062 endif 1063 1064 endif 1065 1066 enddo 1067 1068elseif(iturb.eq.51 ) then 1069 1070 do ifac = 1, nfabor 1071 1072 icodcu = icodcl(ifac,iu) 1073 icodcv = icodcl(ifac,iv) 1074 icodcw = icodcl(ifac,iw) 1075 icodck = icodcl(ifac,ik) 1076 icodce = icodcl(ifac,iep) 1077 icodcp = icodcl(ifac,iphi) 1078 icodca = icodcl(ifac,ial) 1079 1080 if( (icodcu.eq.5 .or. icodcv.eq.5 .or. icodcw.eq.5 .or. & 1081 icodck.eq.5 .or. icodce.eq.5 .or. icodcp.eq.5 .or. & 1082 icodca.eq.5 ) .and. & 1083 (icodcu.ne.5 .or. icodcv.ne.5 .or. icodcw.ne.5 .or. & 1084 icodck.ne.5 .or. icodce.ne.5 .or. icodcp.ne.5 .or. & 1085 icodca.ne.5 ) ) then 1086 1087 if (itypfb(ifac).gt.0) then 1088 itypfb(ifac) = -itypfb(ifac) 1089 endif 1090 icodct(1) = ik 1091 icodct(2) = icodcl(ifac,ik) 1092 icodct(3) = iep 1093 icodct(4) = icodcl(ifac,iep) 1094 icodct(5) = iphi 1095 icodct(6) = icodcl(ifac,iphi) 1096 icodct(7) = ial 1097 icodct(8) = icodcl(ifac,ial) 1098 icodct(9) = icodcu 1099 icodct(10) = icodcv 1100 icodct(11) = icodcw 1101 nstuv2 = nstuv2 + 1 1102 1103 if( (icodcu.eq.6 .or. icodcv.eq.6 .or. icodcw.eq.6 .or. & 1104 icodck.eq.6 .or. icodce.eq.6 .or. icodcp.eq.6 .or. & 1105 icodca.eq.6 ) .and. & 1106 (icodcu.ne.6 .or. icodcv.ne.6 .or. icodcw.ne.6 .or. & 1107 icodck.ne.6 .or. icodce.ne.6 .or. icodcp.ne.6 .or. & 1108 icodca.ne.6 ) ) then 1109 1110 if (itypfb(ifac).gt.0) then 1111 itypfb(ifac) = -itypfb(ifac) 1112 endif 1113 icodct(1) = ik 1114 icodct(2) = icodcl(ifac,ik) 1115 icodct(3) = iep 1116 icodct(4) = icodcl(ifac,iep) 1117 icodct(5) = iphi 1118 icodct(6) = icodcl(ifac,iphi) 1119 icodct(7) = ial 1120 icodct(8) = icodcl(ifac,ial) 1121 icodct(9) = icodcu 1122 icodct(10) = icodcv 1123 icodct(11) = icodcw 1124 nstuv2 = nstuv2 + 1 1125 1126 endif 1127 1128 endif 1129 1130 enddo 1131 1132elseif(iturb.eq.60 ) then 1133 1134 do ifac = 1, nfabor 1135 1136 icodcu = icodcl(ifac,iu) 1137 icodcv = icodcl(ifac,iv) 1138 icodcw = icodcl(ifac,iw) 1139 icodck = icodcl(ifac,ik) 1140 icodom = icodcl(ifac,iomg) 1141 1142 if( (icodcu.eq.5 .or. icodcv.eq.5 .or. icodcw.eq.5 .or. & 1143 icodck.eq.5 .or. icodom.eq.5 ) .and. & 1144 (icodcu.ne.5 .or. icodcv.ne.5 .or. icodcw.ne.5 .or. & 1145 icodck.ne.5 .or. icodom.ne.5 ) ) then 1146 if (itypfb(ifac).gt.0) then 1147 itypfb(ifac) = -itypfb(ifac) 1148 endif 1149 icodct(1) = ik 1150 icodct(2) = icodcl(ifac,ik) 1151 icodct(3) = iomg 1152 icodct(4) = icodcl(ifac,iomg) 1153 icodct(5) = icodcu 1154 icodct(6) = icodcv 1155 icodct(7) = icodcw 1156 nstukw = nstukw + 1 1157 endif 1158 1159 if( (icodcu.eq.6 .or. icodcv.eq.6 .or. icodcw.eq.6 .or. & 1160 icodck.eq.6 .or. icodom.eq.6 ) .and. & 1161 (icodcu.ne.6 .or. icodcv.ne.6 .or. icodcw.ne.6 .or. & 1162 icodck.ne.6 .or. icodom.ne.6 ) ) then 1163 if (itypfb(ifac).gt.0) then 1164 itypfb(ifac) = -itypfb(ifac) 1165 endif 1166 icodct(1) = ik 1167 icodct(2) = icodcl(ifac,ik) 1168 icodct(3) = iomg 1169 icodct(4) = icodcl(ifac,iomg) 1170 icodct(5) = icodcu 1171 icodct(6) = icodcv 1172 icodct(7) = icodcw 1173 nstukw = nstukw + 1 1174 endif 1175 1176 enddo 1177 1178elseif(iturb.eq.70 ) then 1179 1180 do ifac = 1, nfabor 1181 1182 icodcu = icodcl(ifac,iu) 1183 icodcv = icodcl(ifac,iv) 1184 icodcw = icodcl(ifac,iw) 1185 icodcn = icodcl(ifac,inusa) 1186 1187 if( (icodcu.eq.5 .or. icodcv.eq.5 .or. icodcw.eq.5 .or. & 1188 icodcn.eq.5 ) .and. & 1189 (icodcu.ne.5 .or. icodcv.ne.5 .or. icodcw.ne.5 .or. & 1190 icodcn.ne.5 ) ) then 1191 if (itypfb(ifac).gt.0) then 1192 itypfb(ifac) = -itypfb(ifac) 1193 endif 1194 icodct(1) = inusa 1195 icodct(2) = icodcl(ifac,inusa) 1196 icodct(3) = icodcu 1197 icodct(4) = icodcv 1198 icodct(5) = icodcw 1199 nstunu = nstunu + 1 1200 endif 1201 1202 if( (icodcu.eq.6 .or. icodcv.eq.6 .or. icodcw.eq.6 .or. & 1203 icodcn.eq.6 ) .and. & 1204 (icodcu.ne.6 .or. icodcv.ne.6 .or. icodcw.ne.6 .or. & 1205 icodcn.ne.6 ) ) then 1206 if (itypfb(ifac).gt.0) then 1207 itypfb(ifac) = -itypfb(ifac) 1208 endif 1209 icodct(1) = inusa 1210 icodct(2) = icodcl(ifac,inusa) 1211 icodct(3) = icodcu 1212 icodct(4) = icodcv 1213 icodct(5) = icodcw 1214 nstunu = nstunu + 1 1215 endif 1216 1217 enddo 1218endif 1219 1220! 2.5 VERIFICATIONS DES COHERENCES INTER VARIABLES 1221! ================================================ 1222 1223! --- Coherence vitesse scalaires 1224 1225if( nscal.ge.1 ) then 1226 do iis = 1, nscal 1227 if(itytur.eq.2.or.itytur.eq.3) then 1228 ivar = isca(iis) 1229 do ifac = 1, nfabor 1230 icodcu = icodcl(ifac,iu) 1231 if(icodcl(ifac,ivar).eq.5.and.icodcu.ne.5) then 1232 if (itypfb(ifac).gt.0) then 1233 itypfb(ifac) = -itypfb(ifac) 1234 endif 1235 icodus(1) = ivar 1236 icodus(2) = iis 1237 icodus(3) = icodcl(ifac,ivar) 1238 icodus(4) = icodcu 1239 nstusc = nstusc + 1 1240 endif 1241 enddo 1242 endif 1243 enddo 1244endif 1245 1246!=============================================================================== 1247! 3. IMPRESSIONS RECAPITULATIVES 1248!=============================================================================== 1249 1250iok = 0 1251 1252if (nstoni.gt.0 .or. nstosc.gt.0 .or. nstovf.gt.0 .or. nstusc.gt.0 ) then 1253 iok = 1 1254endif 1255 1256if (nstvit.gt.0 .or. nstopp.gt.0 .or. nstoke.gt.0 .or. nstrij.gt.0 .or. & 1257 nstov2.gt.0 .or. nstonu.gt.0 .or. nstuvw.gt.0 .or. nstoup.gt.0 .or. & 1258 nstuke.gt.0 .or. nsurij.gt.0 .or. nstuv2.gt.0 .or. nstunu.gt.0 ) then 1259 iok = 1 1260endif 1261 1262if (irangp.ge.0) call parcmx(iok) 1263 1264if(iok.ne.0) then 1265 1266 call sync_bc_err(nstoni, 2, icodni) 1267 if (nstoni.ne.0) then 1268 call field_get_label(ivarfl(icodni(1)), chaine) 1269 write(nfecra,1000) nstoni, chaine, icodni(2) 1270 endif 1271 1272 call sync_bc_err(nstvit, 3, icodvi) 1273 if (nstvit.ne.0) then 1274 if (icodvi(1) .eq. -6) then 1275 chaine = 'roughness' 1276 else 1277 call field_get_label(ivarfl(icodvi(1)), chaine) 1278 endif 1279 if (icodvi(3).eq.-1) then 1280 write(nfecra,1010) nstvit, chaine, icodvi(2) 1281 endif 1282 endif 1283 1284 if (itytur.eq.2) then 1285 1286 call sync_bc_err(nstoke, 4, icodtb) 1287 if (nstoke.ne.0) then 1288 call field_get_label(ivarfl (icodtb(1)), chaine) 1289 write(nfecra,1010) nstoke, chaine, icodtb(2) 1290 call field_get_label(ivarfl (icodtb(3)), chaine) 1291 write(nfecra,1010) nstoke, chaine, icodtb(4) 1292 endif 1293 1294 elseif (itytur.eq.3) then 1295 1296 call sync_bc_err(nstrij, 3, icodtb) 1297 if (nstrij.ne.0) then 1298 call field_get_label(ivarfl (icodtb(1)), chaine) 1299 write(nfecra,1010) nstrij, chaine, icodtb(2) 1300 if (iturb.eq.32 .and. icodtb(3) .ge. 0) then 1301 write(nfecra,2000) 1302 endif 1303 endif 1304 1305 elseif (itytur.eq.5) then 1306 1307 call sync_bc_err(nstov2, 8, icodtb) 1308 if (nstov2.ne.0) then 1309 call field_get_label(ivarfl (icodtb(1)), chaine) 1310 write(nfecra,1010) nstov2, chaine, icodtb(2) 1311 call field_get_label(ivarfl (icodtb(3)), chaine) 1312 write(nfecra,1010) nstov2, chaine, icodtb(4) 1313 call field_get_label(ivarfl (icodtb(5)), chaine) 1314 write(nfecra,1010) nstov2, chaine, icodtb(6) 1315 call field_get_label(ivarfl (icodtb(7)), chaine) 1316 write(nfecra,1010) nstov2, chaine, icodtb(8) 1317 endif 1318 1319 elseif (itytur.eq.6) then 1320 1321 call sync_bc_err(nstokw, 4, icodtb) 1322 if (nstokw.ne.0) then 1323 call field_get_label(ivarfl (icodtb(1)), chaine) 1324 write(nfecra,1010) nstokw, chaine, icodtb(2) 1325 call field_get_label(ivarfl (icodtb(3)), chaine) 1326 write(nfecra,1010) nstokw, chaine, icodtb(4) 1327 endif 1328 1329 elseif (itytur.eq.7) then 1330 1331 call sync_bc_err(nstonu, 2, icodtb) 1332 if (nstonu.ne.0) then 1333 call field_get_label(ivarfl (icodtb(1)), chaine) 1334 write(nfecra,1000) nstonu, chaine, icodtb(2) 1335 endif 1336 1337 endif 1338 1339 call sync_bc_err(nstosc, 2, icodsc) 1340 if (nstosc.ne.0) then 1341 if (icodsc(1) .eq. -6) then 1342 chaine = 'roughness' 1343 else 1344 call field_get_label(ivarfl (icodsc(1)), chaine) 1345 endif 1346 write(nfecra,1010) nstosc, chaine, icodsc(2) 1347 endif 1348 1349 call sync_bc_err(nstovf, 2, icodvf) 1350 if (nstovf.ne.0) then 1351 call field_get_label(ivarfl (icodvf(1)), chaine) 1352 write(nfecra,1010) nstovf, chaine, icodvf(2) 1353 endif 1354 1355 call sync_bc_err(nstuvw, 3, icoduv) 1356 if (nstuvw.ne.0) then 1357 write(nfecra,1020) nstuvw, icoduv(1), icoduv(2), icoduv(3) 1358 endif 1359 1360 if (itytur.eq.2) then 1361 1362 call sync_bc_err(nstuke, 7, icodct) 1363 if (nstuke.ne.0) then 1364 call field_get_label(ivarfl (icodct(1)), chaine) 1365 write(nfecra,1030) nstuke, chaine, icodct(2), & 1366 icodct(5), icodct(6), icodct(7) 1367 call field_get_label(ivarfl (icodct(3)), chaine) 1368 write(nfecra,1030) nstuke, chaine, icodct(4), & 1369 icodct(5), icodct(6), icodct(7) 1370 endif 1371 1372 elseif (iturb.eq.30 .or. iturb.eq.31) then 1373 1374 call sync_bc_err(nsurij, 10, icodct) 1375 if (nsurij.ne.0) then 1376 write(nfecra,1040) nsurij, icodct(1), icodct(2), & 1377 icodct(3), icodct(4), icodct(5), & 1378 icodct(6), icodct(7), icodct(8), & 1379 icodct(9), icodct(10) 1380 endif 1381 1382 elseif (iturb.eq.32) then 1383 1384 call sync_bc_err(nsurij, 11, icodct) 1385 if (nsurij.ne.0) then 1386 write(nfecra,1041) nsurij, icodct(1), icodct(2), & 1387 icodct(3), icodct(4), icodct(5), & 1388 icodct(6), icodct(7), icodct(8), & 1389 icodct(9), icodct(10), icodct(11) 1390 endif 1391 1392 elseif (itytur.eq.5) then 1393 1394 call sync_bc_err(nstuv2, 11, icodct) 1395 if (nstuv2.ne.0) then 1396 call field_get_label(ivarfl (icodct(1)), chaine) 1397 write(nfecra,1030) nstuv2, chaine, icodct(2), & 1398 icodct(9), icodct(10), icodct(11) 1399 call field_get_label(ivarfl (icodct(3)), chaine) 1400 write(nfecra,1030) nstuv2, chaine, icodct(4), & 1401 icodct(9), icodct(10), icodct(11) 1402 call field_get_label(ivarfl (icodct(5)), chaine) 1403 write(nfecra,1030) nstuv2, chaine, icodct(6), & 1404 icodct(9), icodct(10), icodct(11) 1405 call field_get_label(ivarfl (icodct(7)), chaine) 1406 write(nfecra,1030) nstuv2, chaine, icodct(8), & 1407 icodct(9), icodct(10), icodct(11) 1408 endif 1409 1410 elseif (itytur.eq.6) then 1411 1412 call sync_bc_err(nstukw, 7, icodct) 1413 if (nstukw.ne.0) then 1414 call field_get_label(ivarfl (icodct(1)), chaine) 1415 write(nfecra,1030) nstukw, chaine, icodct(2), & 1416 icodct(5), icodct(6), icodct(7) 1417 call field_get_label(ivarfl (icodct(3)), chaine) 1418 write(nfecra,1030) nstukw, chaine, icodct(4), & 1419 icodct(5), icodct(6), icodct(7) 1420 endif 1421 1422 elseif (itytur.eq.7) then 1423 1424 call sync_bc_err(nstunu, 5, icodct) 1425 if (nstunu.ne.0) then 1426 call field_get_label(ivarfl (icodct(1)), chaine) 1427 write(nfecra,1030) nstunu, chaine, icodct(2), & 1428 icodct(3), icodct(4), icodct(5) 1429 endif 1430 endif 1431 1432 call sync_bc_err(nstusc, 4, icodus) 1433 if (nstusc.ne.0) then 1434 call field_get_label(ivarfl (icodus(1)), chaine) 1435 write(nfecra,1050) nstusc, chaine, & 1436 icodus(2), icodus(3), icodus(4) 1437 endif 1438 1439 if (nstoni.gt.0 .or. nstosc.gt.0 .or. nstovf.gt.0 .or. nstusc.gt.0 ) then 1440 write (nfecra,1901) nstoni, nstosc, nstovf, nstusc 1441 endif 1442 1443 if (nstvit.gt.0 .or. nstopp.gt.0 .or. nstoke.gt.0 .or. nstrij.gt.0 .or. & 1444 nstov2.gt.0 .or. nstonu.gt.0 .or. nstuvw.gt.0 .or. nstoup.gt.0 .or. & 1445 nstuke.gt.0 .or. nsurij.gt.0 .or. nstuv2.gt.0 .or. nstunu.gt.0 ) then 1446 write (nfecra,1902) nstvit,nstopp, nstoke,nstrij, & 1447 nstov2,nstonu, nstuvw,nstoup, & 1448 nstuke,nsurij, nstuv2,nstunu 1449 endif 1450 1451 call boundary_conditions_error(itypfb) 1452 1453endif 1454 1455!=============================================================================== 1456! 3. FORMATS 1457!=============================================================================== 1458 1459 1000 format( & 1460'@ ',/,& 1461'@ UNINITIALIZED BOUNDARY CONDITIONS ',/,& 1462'@ Number of boundary faces:',i10 ,'; variable:',a16 ,/,& 1463'@ icodcl variable last face ', i10 ,/,& 1464'@ ' ) 1465 1010 format( & 1466'@ ',/,& 1467'@ UNEXPECTED BOUNDARY CONDITIONS ',/,& 1468'@ Number of boundary faces:',i10 ,'; variable:',a16 ,/,& 1469'@ icodcl variable last face:', i10 ,/,& 1470'@ ' ) 1471 1020 format( & 1472'@ ',/,& 1473'@ INCOHERENT BOUNDARY CONDITIONS VELOCITY COMPONENT ',/,& 1474'@ Number of boundary faces:',i10 ,/,& 1475'@ icodcl last face:', 3i10 ,/,& 1476'@ ' ) 1477 1030 format( & 1478'@ ',/,& 1479'@ INCOHERENCY BOUNDARY CONDITIONS VELOCITY-VARIABLE ',/,& 1480'@ Number of boundary faces:',i10 ,'; variable:',a16 ,/,& 1481'@ icodcl variable last face:', i10 ,/,& 1482'@ icodcl velocity:',3i10 ,/,& 1483'@ ' ) 1484 1040 format( & 1485'@ ',/,& 1486'@ INCOHERENT BOUNDARY CONDITIONS VELOCITY-RIJ-EPSILON ',/,& 1487'@ Number of boundary faces:',i10 ,/,& 1488'@ icodcl last face Rij-eps: ', 7i5 ,/,& 1489'@ icodcl velocity:', 3i5 ,/,& 1490'@ ' ) 1491 1041 format( & 1492'@ ',/,& 1493'@ INCOHERENT BOUNDARY CONDITIONS VELOCITY-RIJ-EPSILON EBRSM ',/,& 1494'@ Number of boundary faces:',i10 ,/,& 1495'@ icodcl last face Rij-eps: ', 8i5 ,/,& 1496'@ icodcl velocity:', 3i5 ,/,& 1497'@ ' ) 1498 1050 format( & 1499'@ ',/,& 1500'@ INCOHERENT BOUNDARY CONDITIONS VELOCITY-SCALAR ',/,& 1501'@ Number of boundary faces:',i10 ,'; variable:',a16 ,/,& 1502'@ last face:' ,/,& 1503'@ scalar number:',i10 ,/,& 1504'@ icodcl scalar:',i10 ,'; icodcl velocity:',i10 ,/,& 1505'@ ' ) 1506 1901 format( & 1507'@ ',/,& 1508'@ ',/,& 1509'@ ',/,& 1510'@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@',/,& 1511'@ ',/,& 1512'@ @@ WARNING: ABORT DURING THE BOUNDARY CONDITIONS CHECK ',/,& 1513'@ ======== ',/,& 1514'@ ',/,& 1515'@ Uninitialized boundary conditions : ',I10 ,/,& 1516'@ Unexpected boundary conditions: ',/,& 1517'@ on the scalars : ',I10 ,/,& 1518'@ on the scalars representing ',/,& 1519'@ a variance : ',I10 ,/,& 1520'@ Incoherencies: ',/,& 1521'@ between velocity and scalars : ',I10 ,/,& 1522'@ ',/,& 1523'@ The calculation will not be run. ',/,& 1524'@ ',/,& 1525'@ Verify the parameters given via the interface or ',/,& 1526'@ cs_user_boundary_conditions. ',/,& 1527'@ ',/) 1528 1902 format( & 1529'@ ',/,& 1530'@ ',/,& 1531'@ ',/,& 1532'@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@',/,& 1533'@ ',/,& 1534'@ @@ WARNING: ABORT DURING THE BOUNDARY CONDITIONS CHECK ',/,& 1535'@ ======== ',/,& 1536'@ ',/,& 1537'@ Unexpected boundary conditions: ',/,& 1538'@ on the velocity : ',I10 ,/,& 1539'@ on the pressure : ',I10 ,/,& 1540'@ on k and epsilon : ',I10 ,/,& 1541'@ on Rij and epsilon : ',I10 ,/,& 1542'@ on k, epsilon, phi and f_barre : ',I10 ,/,& 1543'@ on nu of Spalart Allmaras model : ',I10 ,/,& 1544'@ Incoherencies: ',/,& 1545'@ between the velocity components : ',I10 ,/,& 1546'@ between velocity and pressure : ',I10 ,/,& 1547'@ between velocity and k-epsilon : ',I10 ,/,& 1548'@ between velocity and Rij-epsilon : ',I10 ,/,& 1549'@ between velocity and v2f : ',I10 ,/,& 1550'@ between velocity and nu : ',I10 ,/,& 1551'@ ',/,& 1552'@ The calculation will not be run. ',/,& 1553'@ ',/,& 1554'@ Verify the parameters given via the interface or ',/,& 1555'@ cs_user_boundary_conditions. ',/,& 1556'@ ',/,& 1557'@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@',/,& 1558'@ ',/) 1559 2000 format( & 1560'@ ',/,& 1561'@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@',/,& 1562'@ ',/,& 1563'@ @@ WARNING: ABORT DURING THE BOUNDARY CONDITIONS CHECK ',/,& 1564'@ ========= ',/,& 1565'@ ROUGH WALL BOUNDARY CONDITIONS INCOMPATIBLE ',/,& 1566'@ WITH EBRSM MODEL ',/,& 1567'@ ',/,& 1568'@ The calculation will not be run. ',/,& 1569'@ ',/,& 1570'@ Verify the parameters given via the interface or ',/,& 1571'@ cs_user_boundary_conditions. ',/,& 1572'@ ',/,& 1573'@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@',/,& 1574'@ ',/) 15752001 format( & 1576'@ ',/,& 1577'@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@',/,& 1578'@ ',/,& 1579'@ @@ WARNING: ABORT DURING THE BOUNDARY CONDITIONS CHECK ',/,& 1580'@ ======= ',/,& 1581'@ FREE-SURFACE CONDITION IN ALE ',/,& 1582'@ MUST BE COMBINED WITH A NON ZERO GRAVITY TERM! ',/,& 1583'@ ',/,& 1584'@ The calculation will not be run. ',/,& 1585'@ ',/,& 1586'@ Verify the parameters or modify the gravity. ',/,& 1587'@ ',/,& 1588'@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@',/,& 1589'@ ',/) 1590 1591return 1592end subroutine vericl 1593 1594!=============================================================================== 1595! Local functions 1596!=============================================================================== 1597 1598!> \brief synchronize boundary condition error logging across MPI ranks. 1599 1600!------------------------------------------------------------------------------- 1601! Arguments 1602!______________________________________________________________________________. 1603! mode name role ! 1604!______________________________________________________________________________! 1605!> \param[in, out] nerloc number of errors (local rank in, global out) 1606!> \param[in] nerrcd number of codes saved at error faces 1607!> \param[in, out] errcod codes saved at one error face (local in, 1608! broadcast out) 1609!_______________________________________________________________________________ 1610 1611subroutine sync_bc_err & 1612 ( nerloc , nerrcd , errcod ) 1613 1614!=============================================================================== 1615! Module files 1616!=============================================================================== 1617 1618use parall 1619 1620!=============================================================================== 1621 1622implicit none 1623 1624! Arguments 1625 1626integer nerloc, nerrcd 1627integer errcod(nerrcd) 1628 1629! Local variables 1630 1631integer irkerr 1632 1633!=============================================================================== 1634 1635if (irangp.ge.0) then 1636 irkerr = -1 1637 if (nerloc.gt.0) irkerr = irangp 1638 call parcpt(nerloc) 1639 if (nerloc .ne. 0) then 1640 call parcmx(irkerr) 1641 call parbci(irkerr, nerrcd, errcod) 1642 endif 1643endif 1644 1645return 1646end subroutine sync_bc_err 1647