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!> \file dvvpst.f90 24!> \brief Standard output of variables on post-processing meshes 25!> (called after \ref cs_user_extra_operations). 26!> 27!------------------------------------------------------------------------------ 28 29!------------------------------------------------------------------------------ 30! Arguments 31!------------------------------------------------------------------------------ 32! mode name role 33!------------------------------------------------------------------------------ 34!> \param[in] nummai post-processing mesh number 35!> \param[in] numtyp post-processing type number 36!> - -1: volume 37!> - -2: edge 38!> - default: nummai 39!> \param[in] nvar total number of variables 40!> \param[in] ncelps post-processing mesh cells number 41!> \param[in] nfbrps number of boundary faces 42!> \param[in] lstcel post-processing mesh cell numbers 43!> \param[in] lstfbr post-processing mesh boundary faces numbers 44!> \param[in,out] tracel post processing cell real values 45!> \param[in,out] trafbr post processing boundary faces real values 46!______________________________________________________________________________ 47 48subroutine dvvpst & 49 ( nummai , numtyp , & 50 nvar , & 51 ncelps , nfbrps , & 52 lstcel , lstfbr , & 53 tracel , trafbr ) 54 55!=============================================================================== 56! Module files 57!=============================================================================== 58 59use paramx 60use pointe 61use entsor 62use cstnum 63use cstphy 64use optcal 65use numvar 66use parall 67use period 68use lagran 69use ppppar 70use ppthch 71use ppincl 72use cplsat 73use mesh 74use field 75use field_operator 76use post 77use cs_f_interfaces 78use rotation 79use turbomachinery 80 81!=============================================================================== 82 83implicit none 84 85! Arguments 86 87integer nummai , numtyp 88integer nvar 89integer ncelps , nfbrps 90 91integer lstcel(ncelps), lstfbr(nfbrps) 92 93double precision tracel(ncelps*3) 94double precision trafbr(nfbrps*3) 95 96! Local variables 97 98character(len=80) :: name80 99 100logical ientla, ivarpr 101integer inc , iccocg 102integer ifac , iloc , ivar 103integer ipp , idimt , kk , ll, iel 104integer fldid, fldprv, keycpl, iflcpl 105integer ifcsii, iflpst, itplus, iprev, f_id 106 107double precision rbid(1) 108double precision vr(3) 109double precision visls_0 110 111double precision, allocatable, dimension(:,:) :: grad 112double precision, dimension(:), pointer :: tplusp 113double precision, dimension(:), pointer :: valsp, coefap, coefbp 114double precision, dimension(:,:), pointer :: valvp, cofavp, cofbvp 115double precision, dimension(:,:,:), pointer :: cofbtp 116double precision, dimension(:), pointer :: crom 117double precision, dimension(:,:), pointer :: vel 118double precision, dimension(:), pointer :: cpotr, cpoti, cvisii 119double precision, dimension(:), pointer :: cvar_pr 120 121!=============================================================================== 122 123! Initialize variables to avoid compiler warnings 124 125ipp = 0 126 127!=============================================================================== 128! Fluid domain 129!=============================================================================== 130 131if (numtyp .eq. -1) then 132 133 ! Map field arrays 134 call field_get_val_s(ivarfl(ipr), cvar_pr) 135 call field_get_val_v(ivarfl(iu), vel) 136 137 ! Automatic additional variables 138 ! ------------------------------ 139 140 ! Relative pressure and velocity in case of turbomachinery 141 142 if (iturbo.ne.0) then 143 144 call field_get_val_s(icrom, crom) 145 146 idimt = 1 147 ientla = .true. 148 ivarpr = .false. 149 150 do iloc = 1, ncelps 151 152 iel = lstcel(iloc) 153 if (irotce(iel).gt.0) then 154 call rotation_velocity(irotce(iel), xyzcen(:,iel), vr) 155 else 156 vr(1) = 0 157 vr(2) = 0 158 vr(3) = 0 159 endif 160 161 tracel(iloc) = cvar_pr(iel) & 162 - crom(iel)*0.5d0*(vr(1)**2 + vr(2)**2 + vr(3)**2) 163 164 enddo 165 166 call post_write_var(nummai, 'Rel Pressure', idimt, ientla, ivarpr, & 167 ntcabs, ttcabs, tracel, rbid, rbid) 168 169 idimt = 3 170 ientla = .true. 171 ivarpr = .false. 172 173 do iloc = 1, ncelps 174 175 iel = lstcel(iloc) 176 if (irotce(iel).gt.0) then 177 call rotation_velocity(irotce(iel), xyzcen(:,iel), vr) 178 else 179 vr(1) = 0 180 vr(2) = 0 181 vr(3) = 0 182 endif 183 184 tracel(1 + (iloc-1)*idimt) = vel(1,iel) - vr(1) 185 tracel(2 + (iloc-1)*idimt) = vel(2,iel) - vr(2) 186 tracel(3 + (iloc-1)*idimt) = vel(3,iel) - vr(3) 187 188 enddo 189 190 call post_write_var(nummai, 'Rel Velocity', idimt, ientla, ivarpr, & 191 ntcabs, ttcabs, tracel, rbid, rbid) 192 193 endif 194 195 ! Absolute pressure and velocity in case of relative coordinate system 196 197 if (icorio.eq.1) then 198 199 call field_get_val_s(icrom, crom) 200 201 idimt = 1 202 ientla = .true. 203 ivarpr = .false. 204 205 do iloc = 1, ncelps 206 207 iel = lstcel(iloc) 208 call rotation_velocity(0, xyzcen(:,iel), vr) 209 210 tracel(iloc) = cvar_pr(iel) + & 211 0.5d0*crom(iel)*(vr(1)**2 + vr(2)**2 + vr(3)**2) 212 213 enddo 214 215 call post_write_var(nummai, 'Abs Pressure', idimt, ientla, ivarpr, & 216 ntcabs, ttcabs, tracel, rbid, rbid) 217 218 idimt = 3 219 ientla = .true. 220 ivarpr = .false. 221 222 do iloc = 1, ncelps 223 224 iel = lstcel(iloc) 225 call rotation_velocity(0, xyzcen(:,iel), vr) 226 227 tracel(1 + (iloc-1)*idimt) = vel(1,iel) + vr(1) 228 tracel(2 + (iloc-1)*idimt) = vel(2,iel) + vr(2) 229 tracel(3 + (iloc-1)*idimt) = vel(3,iel) + vr(3) 230 231 enddo 232 233 call post_write_var(nummai, 'Abs Velocity', idimt, ientla, ivarpr, & 234 ntcabs, ttcabs, tracel, rbid, rbid) 235 236 endif 237 238!=============================================================================== 239! Boundary 240!=============================================================================== 241 242else if (numtyp .eq. -2) then 243 244 ! Projection of variables at boundary with no reconstruction 245 ! ---------------------------------------------------------- 246 247 call field_get_key_id('coupled', keycpl) 248 249 fldprv = -1 250 251 do ivar = 1, nvar ! Loop on main cell-based variables 252 253 fldid = ivarfl(ivar) 254 255 if (fldid .eq. fldprv) cycle ! already output for multiple components 256 257 fldprv = fldid 258 259 call field_get_key_int(fldid, keyvis, iflpst) 260 261 if (iand(iflpst, POST_BOUNDARY_NR) .eq. 0) cycle ! nothing to do 262 263 call field_get_dim (fldid, idimt) 264 call field_get_name(fldid, name80(4:80)) 265 name80(1:3) = 'bc_' 266 267 ! Compute non-reconstructed values at boundary faces 268 269 if (idimt.ne.1) then 270 call field_get_key_int(fldid, keycpl, iflcpl) 271 else 272 iflcpl = 0 273 endif 274 275 if (idimt.eq.1) then ! Scalar 276 277 call field_get_val_s(fldid, valsp) 278 call field_get_coefa_s(fldid, coefap) 279 call field_get_coefb_s(fldid, coefbp) 280 281 do iloc = 1, nfbrps 282 283 ifac = lstfbr(iloc) 284 iel = ifabor(ifac) 285 286 trafbr(iloc) = coefap(ifac) + coefbp(ifac)*valsp(iel) 287 288 enddo 289 290 else if (iflcpl.eq.0) then ! Uncoupled vector or tensor 291 292 call field_get_val_v(fldid, valvp) 293 call field_get_coefa_v(fldid, cofavp) 294 call field_get_coefb_uv(fldid, cofbvp) 295 296 do kk = 0, idimt-1 297 298 do iloc = 1, nfbrps 299 300 ifac = lstfbr(iloc) 301 iel = ifabor(ifac) 302 303 trafbr(kk + (iloc-1)*idimt + 1) & 304 = cofavp(kk+1,ifac) & 305 + cofbvp(kk+1,ifac)*valvp(kk+1,iel) 306 307 enddo 308 309 enddo 310 311 else ! Coupled vector or tensor 312 313 call field_get_val_v(fldid, valvp) 314 call field_get_coefa_v(fldid, cofavp) 315 call field_get_coefb_v(fldid, cofbtp) 316 317 do kk = 0, idimt-1 318 319 do iloc = 1, nfbrps 320 321 ifac = lstfbr(iloc) 322 iel = ifabor(ifac) 323 324 trafbr(kk + (iloc-1)*idimt + 1) = cofavp(kk+1,ifac) 325 326 do ll = 1, idimt 327 trafbr(kk + (iloc-1)*idimt + 1) & 328 = trafbr(kk + (iloc-1)*idimt + 1) & 329 + cofbtp(kk+1,ll,ifac)*valvp(ll,iel) 330 enddo 331 332 enddo 333 334 enddo 335 336 endif ! test on field dimension and interleaving 337 338 ientla = .true. ! interleaved result values 339 ivarpr = .false. ! defined on work array 340 341 call post_write_var(nummai, trim(name80), idimt, ientla, ivarpr, & 342 ntcabs, ttcabs, rbid, rbid, trafbr) 343 344 enddo ! End of loop on variables 345 346 ! Handle stresses at boundary 347 ! --------------------------- 348 349 if (iand(ipstdv(ipstfo), 1) .ne. 0) then 350 351 ! Compute variable values on boundary faces 352 353 call post_stress(nfbrps, lstfbr, trafbr) 354 355 idimt = 3 ! variable dimension 356 ientla = .true. ! interleaved values 357 ivarpr = .false. ! defined on work array 358 359 call post_write_var(nummai, 'Stress', idimt, ientla, ivarpr, & 360 ntcabs, ttcabs, rbid, rbid, trafbr) 361 362 endif 363 364 if (iand(ipstdv(ipstfo), 2) .ne. 0) then 365 366 ! Compute variable values on boundary faces 367 368 call post_stress_tangential(nfbrps, lstfbr, trafbr) 369 370 idimt = 3 ! variable dimension 371 ientla = .true. ! interleaved values 372 ivarpr = .false. ! defined on work array 373 374 call post_write_var(nummai, 'Shear Stress', idimt, ientla, ivarpr, & 375 ntcabs, ttcabs, rbid, rbid, trafbr) 376 377 endif 378 379 if (iand(ipstdv(ipstfo), 4) .ne. 0) then 380 381 ! Calcul des valeurs de la variable sur les faces de bord 382 383 call post_stress_normal(nfbrps, lstfbr, trafbr) 384 385 idimt = 1 ! variable dimension 386 ientla = .true. ! interleaved values 387 ivarpr = .false. ! defined on work array 388 389 call post_write_var(nummai, 'Normal Stress', idimt, ientla, ivarpr, & 390 ntcabs, ttcabs, rbid, rbid, trafbr) 391 392 endif 393 394 ! T+ near the boundary 395 ! -------------------- 396 397 if (ipstdv(ipsttp).ne.0) then 398 399 call field_get_id_try('tplus', itplus) 400 401 if (itplus.ge.0) then 402 403 call field_get_val_s(itplus, tplusp) 404 405 idimt = 1 ! variable dimension 406 ientla = .true. ! interleaved values 407 ivarpr = .true. ! defined on parent array 408 409 if (itherm .eq. 1) then 410 name80 = 'Tplus' 411 else if (itherm .eq. 2) then 412 name80 = 'Hplus' 413 else if (itherm .eq. 3) then 414 name80 = 'Eplus' 415 else 416 return 417 endif 418 419 call post_write_var(nummai, name80, idimt, ientla, ivarpr, & 420 ntcabs, ttcabs, rbid, rbid, tplusp) 421 422 endif ! end of test on presence ot T+ 423 424 endif ! end of test on output of y+ 425 426 ! Thermal flux at boundary 427 ! ------------------------ 428 ! If working with enthalpy, compute an enthalpy flux 429 430 if (ipstdv(ipstft).ne.0) then 431 432 if (iscalt.gt.0) then 433 434 call post_boundary_thermal_flux(nfbrps, lstfbr, trafbr) 435 436 idimt = 1 ! variable dimension 437 ientla = .true. ! interleaved values 438 ivarpr = .false. ! defined on work array 439 440 call post_write_var(nummai, 'Input thermal flux', idimt, ientla, ivarpr, & 441 ntcabs, ttcabs, rbid, rbid, trafbr) 442 443 endif 444 445 endif 446 447 ! Nusselt at the boundary 448 ! ----------------------- 449 450 if (ipstdv(ipstnu).ne.0) then 451 452 idimt = 1 ! variable dimension 453 ientla = .true. ! interleaved values 454 ivarpr = .false. ! defined on work array 455 456 ! Compute variable on boundary faces 457 458 call post_boundary_nusselt(nfbrps, lstfbr, trafbr) 459 460 call post_write_var(nummai, 'Dimensionless heat flux', idimt, ientla, ivarpr, & 461 ntcabs, ttcabs, rbid, rbid, trafbr) 462 463 endif ! end of test on output of Nusselt 464 465endif ! end of test on postprocessing mesh number 466 467!=============================================================================== 468! Electric module variables 469!=============================================================================== 470 471if (numtyp.eq.-1) then 472 473 if ( ippmod(ieljou).ge.1 & 474 .or. ippmod(ielarc).ge.1) then 475 476 allocate(grad(3,ncelet)) 477 478 ! For Joule Heating by direct conduction: 479 ! gradient of the imaginary component of the potential 480 481 if (ippmod(ieljou).eq.2 .or. ippmod(ieljou).eq.4) then 482 483 call field_get_id('elec_pot_i', f_id) 484 485 inc = 1 486 iprev = 0 487 iccocg = 1 488 489 call field_gradient_scalar(f_id, iprev, 0, inc, iccocg, grad) 490 491 idimt = 3 492 ientla = .true. 493 ivarpr = .true. 494 495 call post_write_var(nummai, 'Pot_Gradient_Im', idimt, ientla, ivarpr, & 496 ntcabs, ttcabs, grad, rbid, rbid) 497 498 endif 499 500 ! For Joule heating by direct conduction: 501 ! imaginary component of the current density 502 503 if (ippmod(ieljou).eq.2 .or. ippmod(ieljou).eq.4) then 504 505 call field_get_id('elec_pot_i', f_id) 506 507 ! As in elflux 508 509 inc = 1 510 iprev = 0 511 iccocg = 1 512 513 call field_gradient_scalar(f_id, iprev, 0, inc, iccocg, grad) 514 515 call field_get_key_int (f_id, kivisl, ifcsii) 516 if (ifcsii .ge. 0) then 517 call field_get_val_s(ifcsii, cvisii) 518 do iloc = 1, ncelps 519 iel = lstcel(iloc) 520 tracel(1 + (iloc-1)*idimt) = -cvisii(iel)*grad(1,iel) 521 tracel(2 + (iloc-1)*idimt) = -cvisii(iel)*grad(2,iel) 522 tracel(3 + (iloc-1)*idimt) = -cvisii(iel)*grad(3,iel) 523 enddo 524 else 525 call field_get_key_double(f_id, kvisl0, visls_0) 526 do iloc = 1, ncelps 527 iel = lstcel(iloc) 528 tracel(1 + (iloc-1)*idimt) = -visls_0*grad(1,iel) 529 tracel(2 + (iloc-1)*idimt) = -visls_0*grad(2,iel) 530 tracel(3 + (iloc-1)*idimt) = -visls_0*grad(3,iel) 531 enddo 532 endif 533 534 idimt = 3 535 ientla = .true. 536 ivarpr = .false. 537 538 call post_write_var(nummai, 'Current_Im', idimt, ientla, ivarpr, & 539 ntcabs, ttcabs, tracel, rbid, rbid) 540 541 endif 542 543 ! Calculation of Module and Argument of the complex potential if IELJOU = 4 544 545 if (ippmod(ieljou).eq.4) then 546 547 ivar = 0 548 549 call field_get_val_s_by_name('elec_pot_r', cpotr) 550 call field_get_val_s_by_name('elec_pot_i', cpoti) 551 552 do iloc = 1, ncelps 553 iel = lstcel(iloc) 554 tracel(iloc) = sqrt(cpotr(iel)*cpotr(iel) + cpoti(iel)*cpoti(iel)) 555 enddo 556 557 idimt = 1 558 ientla = .true. 559 ivarpr = .false. 560 561 call post_write_var(nummai, 'Pot_Module', idimt, ientla, ivarpr, & 562 ntcabs, ttcabs, tracel, rbid, rbid) 563 564 do iloc = 1, ncelps 565 566 iel = lstcel(iloc) 567 568 if (cpotr(iel) .ne. 0.d0) then 569 if (cpotr(iel) .ge. 0.d0) then 570 tracel(iloc) = atan(cpoti(iel)/cpotr(iel)) 571 else 572 if (cpoti(iel) .gt. 0.d0) then 573 tracel(iloc) = 4.d0*atan(1.d0) & 574 + atan(cpoti(iel) / cpotr(iel)) 575 else 576 tracel(iloc) = -4.d0*atan(1.d0) & 577 + atan(cpoti(iel) / cpotr(iel)) 578 endif 579 endif 580 else 581 tracel(iloc) = 2.d0*atan(1.d0) 582 endif 583 584 if (tracel(iloc) .lt. 0.d0) then 585 tracel(iloc) = tracel(iloc) + 8.d0**atan(1.d0) 586 endif 587 588 enddo 589 590 idimt = 1 591 ientla = .true. 592 ivarpr = .false. 593 594 call post_write_var(nummai, 'Pot_Arg', idimt, ientla, ivarpr, & 595 ntcabs, ttcabs, tracel, rbid, rbid) 596 597 endif 598 599 ! Free memory 600 deallocate(grad) 601 602 endif 603 604endif 605 606!---- 607! End 608!---- 609 610return 611end subroutine 612