Changeset 15066 for NEMO/branches/2021/dev_r14116_HPC-10_mcastril_Mixed_Precision_implementation/src/OCE/TRA
- Timestamp:
- 2021-07-01T13:17:50+02:00 (3 years ago)
- Location:
- NEMO/branches/2021/dev_r14116_HPC-10_mcastril_Mixed_Precision_implementation/src/OCE/TRA
- Files:
-
- 2 edited
-
eosbn2.F90 (modified) (14 diffs)
-
tranpc.F90 (modified) (1 diff)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2021/dev_r14116_HPC-10_mcastril_Mixed_Precision_implementation/src/OCE/TRA/eosbn2.F90
r14986 r15066 56 56 ! !! * Interface 57 57 INTERFACE eos 58 MODULE PROCEDURE eos_insitu, eos_insitu_pot, eos_insitu_2d, eos_insitu_pot_2d 58 MODULE PROCEDURE eos_insitu_wp, eos_insitu_pot_wp, eos_insitu_2d, eos_insitu_pot_2d 59 MODULE PROCEDURE eos_insitu_mixed, eos_insitu_pot_mixed 59 60 END INTERFACE 60 61 ! 61 62 INTERFACE eos_rab 62 MODULE PROCEDURE rab_3d, rab_2d, rab_0d 63 MODULE PROCEDURE rab_3d_wp, rab_2d, rab_0d 64 MODULE PROCEDURE rab_3d_mixed 63 65 END INTERFACE 64 66 ! … … 190 192 CONTAINS 191 193 192 SUBROUTINE eos_insitu ( pts, prd, pdep )194 SUBROUTINE eos_insitu_wp( pts, prd, pdep ) 193 195 !! 194 196 REAL(wp), DIMENSION(:,:,:,:), INTENT(in ) :: pts ! 1 : potential temperature [Celsius] … … 197 199 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pdep ! depth [m] 198 200 !! 199 CALL eos_insitu_t( pts, is_tile(pts), prd, is_tile(prd), pdep, is_tile(pdep) ) 200 END SUBROUTINE eos_insitu 201 202 SUBROUTINE eos_insitu_t( pts, ktts, prd, ktrd, pdep, ktdep ) 203 !!---------------------------------------------------------------------- 204 !! *** ROUTINE eos_insitu *** 201 CALL eos_insitu_t_wp( pts, is_tile(pts), prd, is_tile(prd), pdep, is_tile(pdep) ) 202 END SUBROUTINE eos_insitu_wp 203 204 SUBROUTINE eos_insitu_mixed( pts, prd, pdep ) 205 !! 206 REAL(dp), DIMENSION(:,:,:,:), INTENT(in ) :: pts ! 1 : potential temperature [Celsius] 207 ! ! 2 : salinity [psu] 208 REAL(sp), DIMENSION(:,:,:) , INTENT( out) :: prd ! in situ density [-] 209 REAL(sp), DIMENSION(:,:,:) , INTENT(in ) :: pdep ! depth [m] 210 !! 211 CALL eos_insitu_t_mixed( pts, is_tile(pts), prd, is_tile(prd), pdep, is_tile(pdep) ) 212 END SUBROUTINE 213 214 SUBROUTINE eos_insitu_t_wp( pts, ktts, prd, ktrd, pdep, ktdep ) 215 !!---------------------------------------------------------------------- 216 !! *** ROUTINE 205 217 !! 206 218 !! ** Purpose : Compute the in situ density (ratio rho/rho0) from … … 306 318 IF( ln_timing ) CALL timing_stop('eos-insitu') 307 319 ! 308 END SUBROUTINE eos_insitu_t 309 310 311 SUBROUTINE eos_insitu_pot( pts, prd, prhop, pdep ) 320 END SUBROUTINE eos_insitu_t_wp 321 SUBROUTINE eos_insitu_t_mixed( pts, ktts, prd, ktrd, pdep, ktdep ) 322 !!---------------------------------------------------------------------- 323 !! *** ROUTINE eos_insitu *** 324 !! 325 !! ** Purpose : Compute the in situ density (ratio rho/rho0) from 326 !! potential temperature and salinity using an equation of state 327 !! selected in the nameos namelist 328 !! 329 !! ** Method : prd(t,s,z) = ( rho(t,s,z) - rho0 ) / rho0 330 !! with prd in situ density anomaly no units 331 !! t TEOS10: CT or EOS80: PT Celsius 332 !! s TEOS10: SA or EOS80: SP TEOS10: g/kg or EOS80: psu 333 !! z depth meters 334 !! rho in situ density kg/m^3 335 !! rho0 reference density kg/m^3 336 !! 337 !! ln_teos10 : polynomial TEOS-10 equation of state is used for rho(t,s,z). 338 !! Check value: rho = 1028.21993233072 kg/m^3 for z=3000 dbar, ct=3 Celsius, sa=35.5 g/kg 339 !! 340 !! ln_eos80 : polynomial EOS-80 equation of state is used for rho(t,s,z). 341 !! Check value: rho = 1028.35011066567 kg/m^3 for z=3000 dbar, pt=3 Celsius, sp=35.5 psu 342 !! 343 !! ln_seos : simplified equation of state 344 !! prd(t,s,z) = ( -a0*(1+lambda/2*(T-T0)+mu*z+nu*(S-S0))*(T-T0) + b0*(S-S0) ) / rho0 345 !! linear case function of T only: rn_alpha<>0, other coefficients = 0 346 !! linear eos function of T and S: rn_alpha and rn_beta<>0, other coefficients=0 347 !! Vallis like equation: use default values of coefficients 348 !! 349 !! ** Action : compute prd , the in situ density (no units) 350 !! 351 !! References : Roquet et al, Ocean Modelling, in preparation (2014) 352 !! Vallis, Atmospheric and Oceanic Fluid Dynamics, 2006 353 !! TEOS-10 Manual, 2010 354 !!---------------------------------------------------------------------- 355 INTEGER , INTENT(in ) :: ktts, ktrd, ktdep 356 REAL(dp), DIMENSION(A2D_T(ktts) ,JPK,JPTS), INTENT(in ) :: pts ! 1 : potential temperature [Celsius] 357 ! ! 2 : salinity [psu] 358 REAL(sp), DIMENSION(A2D_T(ktrd) ,JPK ), INTENT( out) :: prd ! in situ density [-] 359 REAL(sp), DIMENSION(A2D_T(ktdep),JPK ), INTENT(in ) :: pdep ! depth [m] 360 ! 361 INTEGER :: ji, jj, jk ! dummy loop indices 362 REAL(sp) :: zt , zh , zs , ztm ! local scalars 363 REAL(sp) :: zn , zn0, zn1, zn2, zn3 ! - - 364 !!---------------------------------------------------------------------- 365 ! 366 IF( ln_timing ) CALL timing_start('eos-insitu') 367 ! 368 SELECT CASE( neos ) 369 ! 370 CASE( np_teos10, np_eos80 ) !== polynomial TEOS-10 / EOS-80 ==! 371 ! 372 DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 ) 373 ! 374 zh = pdep(ji,jj,jk) * r1_Z0 ! depth 375 zt = pts (ji,jj,jk,jp_tem) * r1_T0 ! temperature 376 zs = SQRT( ABS( pts(ji,jj,jk,jp_sal) + rdeltaS ) * r1_S0 ) ! square root salinity 377 ztm = tmask(ji,jj,jk) ! tmask 378 ! 379 zn3 = EOS013*zt & 380 & + EOS103*zs+EOS003 381 ! 382 zn2 = (EOS022*zt & 383 & + EOS112*zs+EOS012)*zt & 384 & + (EOS202*zs+EOS102)*zs+EOS002 385 ! 386 zn1 = (((EOS041*zt & 387 & + EOS131*zs+EOS031)*zt & 388 & + (EOS221*zs+EOS121)*zs+EOS021)*zt & 389 & + ((EOS311*zs+EOS211)*zs+EOS111)*zs+EOS011)*zt & 390 & + (((EOS401*zs+EOS301)*zs+EOS201)*zs+EOS101)*zs+EOS001 391 ! 392 zn0 = (((((EOS060*zt & 393 & + EOS150*zs+EOS050)*zt & 394 & + (EOS240*zs+EOS140)*zs+EOS040)*zt & 395 & + ((EOS330*zs+EOS230)*zs+EOS130)*zs+EOS030)*zt & 396 & + (((EOS420*zs+EOS320)*zs+EOS220)*zs+EOS120)*zs+EOS020)*zt & 397 & + ((((EOS510*zs+EOS410)*zs+EOS310)*zs+EOS210)*zs+EOS110)*zs+EOS010)*zt & 398 & + (((((EOS600*zs+EOS500)*zs+EOS400)*zs+EOS300)*zs+EOS200)*zs+EOS100)*zs+EOS000 399 ! 400 zn = ( ( zn3 * zh + zn2 ) * zh + zn1 ) * zh + zn0 401 ! 402 prd(ji,jj,jk) = ( zn * r1_rho0 - 1._wp ) * ztm ! density anomaly (masked) 403 ! 404 END_3D 405 ! 406 CASE( np_seos ) !== simplified EOS ==! 407 ! 408 DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 ) 409 zt = pts (ji,jj,jk,jp_tem) - 10._wp 410 zs = pts (ji,jj,jk,jp_sal) - 35._wp 411 zh = pdep (ji,jj,jk) 412 ztm = tmask(ji,jj,jk) 413 ! 414 zn = - rn_a0 * ( 1._wp + 0.5_wp*rn_lambda1*zt + rn_mu1*zh ) * zt & 415 & + rn_b0 * ( 1._wp - 0.5_wp*rn_lambda2*zs - rn_mu2*zh ) * zs & 416 & - rn_nu * zt * zs 417 ! 418 prd(ji,jj,jk) = zn * r1_rho0 * ztm ! density anomaly (masked) 419 END_3D 420 ! 421 END SELECT 422 ! 423 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=REAL(prd, wp), clinfo1=' eos-insitu : ', kdim=jpk ) 424 ! 425 IF( ln_timing ) CALL timing_stop('eos-insitu') 426 ! 427 END SUBROUTINE eos_insitu_t_mixed 428 429 430 SUBROUTINE eos_insitu_pot_wp( pts, prd, prhop, pdep ) 312 431 !! 313 432 REAL(wp), DIMENSION(:,:,:,:), INTENT(in ) :: pts ! 1 : potential temperature [Celsius] … … 317 436 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pdep ! depth [m] 318 437 !! 319 CALL eos_insitu_pot_t( pts, is_tile(pts), prd, is_tile(prd), prhop, is_tile(prhop), pdep, is_tile(pdep) ) 320 END SUBROUTINE eos_insitu_pot 321 322 323 SUBROUTINE eos_insitu_pot_t( pts, ktts, prd, ktrd, prhop, ktrhop, pdep, ktdep ) 438 CALL eos_insitu_pot_t_wp( pts, is_tile(pts), prd, is_tile(prd), prhop, is_tile(prhop), pdep, is_tile(pdep) ) 439 END SUBROUTINE eos_insitu_pot_wp 440 441 SUBROUTINE eos_insitu_pot_mixed( pts, prd, prhop, pdep ) 442 !! 443 REAL(dp), DIMENSION(:,:,:,:), INTENT(in ) :: pts ! 1 : potential temperature [Celsius] 444 ! ! 2 : salinity [psu] 445 REAL(sp), DIMENSION(:,:,:) , INTENT( out) :: prd ! in situ density [-] 446 REAL(dp), DIMENSION(:,:,:) , INTENT( out) :: prhop ! potential density (surface referenced) 447 REAL(sp), DIMENSION(:,:,:) , INTENT(in ) :: pdep ! depth [m] 448 !! 449 CALL eos_insitu_pot_t_mixed( pts, is_tile(pts), prd, is_tile(prd), prhop, is_tile(prhop), pdep, is_tile(pdep) ) 450 END SUBROUTINE eos_insitu_pot_mixed 451 452 453 SUBROUTINE eos_insitu_pot_t_wp( pts, ktts, prd, ktrd, prhop, ktrhop, pdep, ktdep ) 324 454 !!---------------------------------------------------------------------- 325 455 !! *** ROUTINE eos_insitu_pot *** … … 475 605 IF( ln_timing ) CALL timing_stop('eos-pot') 476 606 ! 477 END SUBROUTINE eos_insitu_pot_t 607 END SUBROUTINE eos_insitu_pot_t_wp 608 609 SUBROUTINE eos_insitu_pot_t_mixed( pts, ktts, prd, ktrd, prhop, ktrhop, pdep, ktdep ) 610 !!---------------------------------------------------------------------- 611 !! *** ROUTINE eos_insitu_pot *** 612 !! 613 !! ** Purpose : Compute the in situ density (ratio rho/rho0) and the 614 !! potential volumic mass (Kg/m3) from potential temperature and 615 !! salinity fields using an equation of state selected in the 616 !! namelist. 617 !! 618 !! ** Action : - prd , the in situ density (no units) 619 !! - prhop, the potential volumic mass (Kg/m3) 620 !! 621 !!---------------------------------------------------------------------- 622 INTEGER , INTENT(in ) :: ktts, ktrd, ktrhop, ktdep 623 REAL(dp), DIMENSION(A2D_T(ktts) ,JPK,JPTS), INTENT(in ) :: pts ! 1 : potential temperature [Celsius] 624 ! ! 2 : salinity [psu] 625 REAL(sp), DIMENSION(A2D_T(ktrd) ,JPK ), INTENT( out) :: prd ! in situ density [-] 626 REAL(dp), DIMENSION(A2D_T(ktrhop),JPK ), INTENT( out) :: prhop ! potential density (surface referenced) 627 REAL(sp), DIMENSION(A2D_T(ktdep) ,JPK ), INTENT(in ) :: pdep ! depth [m] 628 ! 629 INTEGER :: ji, jj, jk, jsmp ! dummy loop indices 630 INTEGER :: jdof 631 REAL(wp) :: zt , zh , zstemp, zs , ztm ! local scalars 632 REAL(wp) :: zn , zn0, zn1, zn2, zn3 ! - - 633 REAL(wp), DIMENSION(:), ALLOCATABLE :: zn0_sto, zn_sto, zsign ! local vectors 634 !!---------------------------------------------------------------------- 635 ! 636 IF( ln_timing ) CALL timing_start('eos-pot') 637 ! 638 SELECT CASE ( neos ) 639 ! 640 CASE( np_teos10, np_eos80 ) !== polynomial TEOS-10 / EOS-80 ==! 641 ! 642 ! Stochastic equation of state 643 IF ( ln_sto_eos ) THEN 644 ALLOCATE(zn0_sto(1:2*nn_sto_eos)) 645 ALLOCATE(zn_sto(1:2*nn_sto_eos)) 646 ALLOCATE(zsign(1:2*nn_sto_eos)) 647 DO jsmp = 1, 2*nn_sto_eos, 2 648 zsign(jsmp) = 1._wp 649 zsign(jsmp+1) = -1._wp 650 END DO 651 ! 652 DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 ) 653 ! 654 ! compute density (2*nn_sto_eos) times: 655 ! (1) for t+dt, s+ds (with the random TS fluctutation computed in sto_pts) 656 ! (2) for t-dt, s-ds (with the opposite fluctuation) 657 DO jsmp = 1, nn_sto_eos*2 658 jdof = (jsmp + 1) / 2 659 zh = pdep(ji,jj,jk) * r1_Z0 ! depth 660 zt = (pts (ji,jj,jk,jp_tem) + pts_ran(ji,jj,jk,jp_tem,jdof) * zsign(jsmp)) * r1_T0 ! temperature 661 zstemp = pts (ji,jj,jk,jp_sal) + pts_ran(ji,jj,jk,jp_sal,jdof) * zsign(jsmp) 662 zs = SQRT( ABS( zstemp + rdeltaS ) * r1_S0 ) ! square root salinity 663 ztm = tmask(ji,jj,jk) ! tmask 664 ! 665 zn3 = EOS013*zt & 666 & + EOS103*zs+EOS003 667 ! 668 zn2 = (EOS022*zt & 669 & + EOS112*zs+EOS012)*zt & 670 & + (EOS202*zs+EOS102)*zs+EOS002 671 ! 672 zn1 = (((EOS041*zt & 673 & + EOS131*zs+EOS031)*zt & 674 & + (EOS221*zs+EOS121)*zs+EOS021)*zt & 675 & + ((EOS311*zs+EOS211)*zs+EOS111)*zs+EOS011)*zt & 676 & + (((EOS401*zs+EOS301)*zs+EOS201)*zs+EOS101)*zs+EOS001 677 ! 678 zn0_sto(jsmp) = (((((EOS060*zt & 679 & + EOS150*zs+EOS050)*zt & 680 & + (EOS240*zs+EOS140)*zs+EOS040)*zt & 681 & + ((EOS330*zs+EOS230)*zs+EOS130)*zs+EOS030)*zt & 682 & + (((EOS420*zs+EOS320)*zs+EOS220)*zs+EOS120)*zs+EOS020)*zt & 683 & + ((((EOS510*zs+EOS410)*zs+EOS310)*zs+EOS210)*zs+EOS110)*zs+EOS010)*zt & 684 & + (((((EOS600*zs+EOS500)*zs+EOS400)*zs+EOS300)*zs+EOS200)*zs+EOS100)*zs+EOS000 685 ! 686 zn_sto(jsmp) = ( ( zn3 * zh + zn2 ) * zh + zn1 ) * zh + zn0_sto(jsmp) 687 END DO 688 ! 689 ! compute stochastic density as the mean of the (2*nn_sto_eos) densities 690 prhop(ji,jj,jk) = 0._wp ; prd(ji,jj,jk) = 0._wp 691 DO jsmp = 1, nn_sto_eos*2 692 prhop(ji,jj,jk) = prhop(ji,jj,jk) + zn0_sto(jsmp) ! potential density referenced at the surface 693 ! 694 prd(ji,jj,jk) = prd(ji,jj,jk) + ( zn_sto(jsmp) * r1_rho0 - 1._wp ) ! density anomaly (masked) 695 END DO 696 prhop(ji,jj,jk) = 0.5_wp * prhop(ji,jj,jk) * ztm / nn_sto_eos 697 prd (ji,jj,jk) = 0.5_wp * prd (ji,jj,jk) * ztm / nn_sto_eos 698 END_3D 699 DEALLOCATE(zn0_sto,zn_sto,zsign) 700 ! Non-stochastic equation of state 701 ELSE 702 DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 ) 703 ! 704 zh = pdep(ji,jj,jk) * r1_Z0 ! depth 705 zt = pts (ji,jj,jk,jp_tem) * r1_T0 ! temperature 706 zs = SQRT( ABS( pts(ji,jj,jk,jp_sal) + rdeltaS ) * r1_S0 ) ! square root salinity 707 ztm = tmask(ji,jj,jk) ! tmask 708 ! 709 zn3 = EOS013*zt & 710 & + EOS103*zs+EOS003 711 ! 712 zn2 = (EOS022*zt & 713 & + EOS112*zs+EOS012)*zt & 714 & + (EOS202*zs+EOS102)*zs+EOS002 715 ! 716 zn1 = (((EOS041*zt & 717 & + EOS131*zs+EOS031)*zt & 718 & + (EOS221*zs+EOS121)*zs+EOS021)*zt & 719 & + ((EOS311*zs+EOS211)*zs+EOS111)*zs+EOS011)*zt & 720 & + (((EOS401*zs+EOS301)*zs+EOS201)*zs+EOS101)*zs+EOS001 721 ! 722 zn0 = (((((EOS060*zt & 723 & + EOS150*zs+EOS050)*zt & 724 & + (EOS240*zs+EOS140)*zs+EOS040)*zt & 725 & + ((EOS330*zs+EOS230)*zs+EOS130)*zs+EOS030)*zt & 726 & + (((EOS420*zs+EOS320)*zs+EOS220)*zs+EOS120)*zs+EOS020)*zt & 727 & + ((((EOS510*zs+EOS410)*zs+EOS310)*zs+EOS210)*zs+EOS110)*zs+EOS010)*zt & 728 & + (((((EOS600*zs+EOS500)*zs+EOS400)*zs+EOS300)*zs+EOS200)*zs+EOS100)*zs+EOS000 729 ! 730 zn = ( ( zn3 * zh + zn2 ) * zh + zn1 ) * zh + zn0 731 ! 732 prhop(ji,jj,jk) = zn0 * ztm ! potential density referenced at the surface 733 ! 734 prd(ji,jj,jk) = ( zn * r1_rho0 - 1._wp ) * ztm ! density anomaly (masked) 735 END_3D 736 ENDIF 737 738 CASE( np_seos ) !== simplified EOS ==! 739 ! 740 DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 ) 741 zt = pts (ji,jj,jk,jp_tem) - 10._wp 742 zs = pts (ji,jj,jk,jp_sal) - 35._wp 743 zh = pdep (ji,jj,jk) 744 ztm = tmask(ji,jj,jk) 745 ! ! potential density referenced at the surface 746 zn = - rn_a0 * ( 1._wp + 0.5_wp*rn_lambda1*zt ) * zt & 747 & + rn_b0 * ( 1._wp - 0.5_wp*rn_lambda2*zs ) * zs & 748 & - rn_nu * zt * zs 749 prhop(ji,jj,jk) = ( rho0 + zn ) * ztm 750 ! ! density anomaly (masked) 751 zn = zn - ( rn_a0 * rn_mu1 * zt + rn_b0 * rn_mu2 * zs ) * zh 752 prd(ji,jj,jk) = zn * r1_rho0 * ztm 753 ! 754 END_3D 755 ! 756 END SELECT 757 ! 758 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=REAL(prd, wp), clinfo1=' eos-pot: ', & 759 & tab3d_2=REAL(prhop, wp), clinfo2=' pot : ', kdim=jpk ) 760 ! 761 IF( ln_timing ) CALL timing_stop('eos-pot') 762 ! 763 END SUBROUTINE eos_insitu_pot_t_mixed 478 764 479 765 … … 661 947 662 948 663 SUBROUTINE rab_3d ( pts, pab, Kmm )949 SUBROUTINE rab_3d_wp( pts, pab, Kmm ) 664 950 !! 665 951 INTEGER , INTENT(in ) :: Kmm ! time level index … … 667 953 REAL(wp), DIMENSION(:,:,:,:), INTENT( out) :: pab ! thermal/haline expansion ratio 668 954 !! 669 CALL rab_3d_t( pts, is_tile(pts), pab, is_tile(pab), Kmm ) 670 END SUBROUTINE rab_3d 671 672 673 SUBROUTINE rab_3d_t( pts, ktts, pab, ktab, Kmm ) 955 CALL rab_3d_t_wp( pts, is_tile(pts), pab, is_tile(pab), Kmm ) 956 END SUBROUTINE rab_3d_wp 957 958 SUBROUTINE rab_3d_mixed( pts, pab, Kmm ) 959 !! 960 INTEGER , INTENT(in ) :: Kmm ! time level index 961 REAL(dp), DIMENSION(:,:,:,:), INTENT(in ) :: pts ! pot. temperature & salinity 962 REAL(sp), DIMENSION(:,:,:,:), INTENT( out) :: pab ! thermal/haline expansion ratio 963 !! 964 CALL rab_3d_t_mixed( pts, is_tile(pts), pab, is_tile(pab), Kmm ) 965 END SUBROUTINE rab_3d_mixed 966 967 968 SUBROUTINE rab_3d_t_wp( pts, ktts, pab, ktab, Kmm ) 674 969 !!---------------------------------------------------------------------- 675 970 !! *** ROUTINE rab_3d *** … … 775 1070 IF( ln_timing ) CALL timing_stop('rab_3d') 776 1071 ! 777 END SUBROUTINE rab_3d_t 778 779 780 SUBROUTINE rab_2d( pts, pdep, pab, Kmm ) 781 !! 782 INTEGER , INTENT(in ) :: Kmm ! time level index 783 REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: pts ! pot. temperature & salinity 784 REAL(wp), DIMENSION(:,:) , INTENT(in ) :: pdep ! depth [m] 785 REAL(wp), DIMENSION(:,:,:), INTENT( out) :: pab ! thermal/haline expansion ratio 786 !! 787 CALL rab_2d_t(pts, is_tile(pts), pdep, is_tile(pdep), pab, is_tile(pab), Kmm) 788 END SUBROUTINE rab_2d 789 790 791 SUBROUTINE rab_2d_t( pts, ktts, pdep, ktdep, pab, ktab, Kmm ) 792 !!---------------------------------------------------------------------- 793 !! *** ROUTINE rab_2d *** 794 !! 795 !! ** Purpose : Calculates thermal/haline expansion ratio for a 2d field (unmasked) 1072 END SUBROUTINE rab_3d_t_wp 1073 1074 SUBROUTINE rab_3d_t_mixed( pts, ktts, pab, ktab, Kmm ) 1075 !!---------------------------------------------------------------------- 1076 !! *** ROUTINE rab_3d *** 1077 !! 1078 !! ** Purpose : Calculates thermal/haline expansion ratio at T-points 1079 !! 1080 !! ** Method : calculates alpha / beta at T-points 796 1081 !! 797 1082 !! ** Action : - pab : thermal/haline expansion ratio at T-points 798 1083 !!---------------------------------------------------------------------- 799 INTEGER , INTENT(in ) :: Kmm ! time level index 800 INTEGER , INTENT(in ) :: ktts, ktdep, ktab 801 REAL(wp), DIMENSION(A2D_T(ktts),JPTS), INTENT(in ) :: pts ! pot. temperature & salinity 802 REAL(wp), DIMENSION(A2D_T(ktdep) ), INTENT(in ) :: pdep ! depth [m] 803 REAL(wp), DIMENSION(A2D_T(ktab),JPTS), INTENT( out) :: pab ! thermal/haline expansion ratio 1084 INTEGER , INTENT(in ) :: Kmm ! time level index 1085 INTEGER , INTENT(in ) :: ktts, ktab 1086 REAL(dp), DIMENSION(A2D_T(ktts),JPK,JPTS), INTENT(in ) :: pts ! pot. temperature & salinity 1087 REAL(sp), DIMENSION(A2D_T(ktab),JPK,JPTS), INTENT( out) :: pab ! thermal/haline expansion ratio 804 1088 ! 805 1089 INTEGER :: ji, jj, jk ! dummy loop indices 806 REAL(wp) :: zt , zh , zs 1090 REAL(wp) :: zt , zh , zs , ztm ! local scalars 807 1091 REAL(wp) :: zn , zn0, zn1, zn2, zn3 ! - - 808 1092 !!---------------------------------------------------------------------- 809 1093 ! 810 IF( ln_timing ) CALL timing_start('rab_2d') 811 ! 812 pab(:,:,:) = 0._wp 1094 IF( ln_timing ) CALL timing_start('rab_3d') 813 1095 ! 814 1096 SELECT CASE ( neos ) … … 816 1098 CASE( np_teos10, np_eos80 ) !== polynomial TEOS-10 / EOS-80 ==! 817 1099 ! 818 DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) 819 ! 820 zh = pdep(ji,jj) * r1_Z0 ! depth 821 zt = pts (ji,jj,jp_tem) * r1_T0 ! temperature 822 zs = SQRT( ABS( pts(ji,jj,jp_sal) + rdeltaS ) * r1_S0 ) ! square root salinity 1100 DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 ) 1101 ! 1102 zh = gdept(ji,jj,jk,Kmm) * r1_Z0 ! depth 1103 zt = pts (ji,jj,jk,jp_tem) * r1_T0 ! temperature 1104 zs = SQRT( ABS( pts(ji,jj,jk,jp_sal) + rdeltaS ) * r1_S0 ) ! square root salinity 1105 ztm = tmask(ji,jj,jk) ! tmask 823 1106 ! 824 1107 ! alpha … … 841 1124 zn = ( ( zn3 * zh + zn2 ) * zh + zn1 ) * zh + zn0 842 1125 ! 843 pab(ji,jj,j p_tem) = zn * r1_rho01126 pab(ji,jj,jk,jp_tem) = zn * r1_rho0 * ztm 844 1127 ! 845 1128 ! beta … … 862 1145 zn = ( ( zn3 * zh + zn2 ) * zh + zn1 ) * zh + zn0 863 1146 ! 1147 pab(ji,jj,jk,jp_sal) = zn / zs * r1_rho0 * ztm 1148 ! 1149 END_3D 1150 ! 1151 CASE( np_seos ) !== simplified EOS ==! 1152 ! 1153 DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 ) 1154 zt = pts (ji,jj,jk,jp_tem) - 10._wp ! pot. temperature anomaly (t-T0) 1155 zs = pts (ji,jj,jk,jp_sal) - 35._wp ! abs. salinity anomaly (s-S0) 1156 zh = gdept(ji,jj,jk,Kmm) ! depth in meters at t-point 1157 ztm = tmask(ji,jj,jk) ! land/sea bottom mask = surf. mask 1158 ! 1159 zn = rn_a0 * ( 1._wp + rn_lambda1*zt + rn_mu1*zh ) + rn_nu*zs 1160 pab(ji,jj,jk,jp_tem) = zn * r1_rho0 * ztm ! alpha 1161 ! 1162 zn = rn_b0 * ( 1._wp - rn_lambda2*zs - rn_mu2*zh ) - rn_nu*zt 1163 pab(ji,jj,jk,jp_sal) = zn * r1_rho0 * ztm ! beta 1164 ! 1165 END_3D 1166 ! 1167 CASE DEFAULT 1168 WRITE(ctmp1,*) ' bad flag value for neos = ', neos 1169 CALL ctl_stop( 'rab_3d:', ctmp1 ) 1170 ! 1171 END SELECT 1172 ! 1173 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=REAL(pab(:,:,:,jp_tem), wp), clinfo1=' rab_3d_t: ', & 1174 & tab3d_2=REAL(pab(:,:,:,jp_sal), wp), clinfo2=' rab_3d_s : ', kdim=jpk ) 1175 ! 1176 IF( ln_timing ) CALL timing_stop('rab_3d') 1177 ! 1178 END SUBROUTINE rab_3d_t_mixed 1179 1180 1181 SUBROUTINE rab_2d( pts, pdep, pab, Kmm ) 1182 !! 1183 INTEGER , INTENT(in ) :: Kmm ! time level index 1184 REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: pts ! pot. temperature & salinity 1185 REAL(wp), DIMENSION(:,:) , INTENT(in ) :: pdep ! depth [m] 1186 REAL(wp), DIMENSION(:,:,:), INTENT( out) :: pab ! thermal/haline expansion ratio 1187 !! 1188 CALL rab_2d_t(pts, is_tile(pts), pdep, is_tile(pdep), pab, is_tile(pab), Kmm) 1189 END SUBROUTINE rab_2d 1190 1191 1192 SUBROUTINE rab_2d_t( pts, ktts, pdep, ktdep, pab, ktab, Kmm ) 1193 !!---------------------------------------------------------------------- 1194 !! *** ROUTINE rab_2d *** 1195 !! 1196 !! ** Purpose : Calculates thermal/haline expansion ratio for a 2d field (unmasked) 1197 !! 1198 !! ** Action : - pab : thermal/haline expansion ratio at T-points 1199 !!---------------------------------------------------------------------- 1200 INTEGER , INTENT(in ) :: Kmm ! time level index 1201 INTEGER , INTENT(in ) :: ktts, ktdep, ktab 1202 REAL(wp), DIMENSION(A2D_T(ktts),JPTS), INTENT(in ) :: pts ! pot. temperature & salinity 1203 REAL(wp), DIMENSION(A2D_T(ktdep) ), INTENT(in ) :: pdep ! depth [m] 1204 REAL(wp), DIMENSION(A2D_T(ktab),JPTS), INTENT( out) :: pab ! thermal/haline expansion ratio 1205 ! 1206 INTEGER :: ji, jj, jk ! dummy loop indices 1207 REAL(wp) :: zt , zh , zs ! local scalars 1208 REAL(wp) :: zn , zn0, zn1, zn2, zn3 ! - - 1209 !!---------------------------------------------------------------------- 1210 ! 1211 IF( ln_timing ) CALL timing_start('rab_2d') 1212 ! 1213 pab(:,:,:) = 0._wp 1214 ! 1215 SELECT CASE ( neos ) 1216 ! 1217 CASE( np_teos10, np_eos80 ) !== polynomial TEOS-10 / EOS-80 ==! 1218 ! 1219 DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) 1220 ! 1221 zh = pdep(ji,jj) * r1_Z0 ! depth 1222 zt = pts (ji,jj,jp_tem) * r1_T0 ! temperature 1223 zs = SQRT( ABS( pts(ji,jj,jp_sal) + rdeltaS ) * r1_S0 ) ! square root salinity 1224 ! 1225 ! alpha 1226 zn3 = ALP003 1227 ! 1228 zn2 = ALP012*zt + ALP102*zs+ALP002 1229 ! 1230 zn1 = ((ALP031*zt & 1231 & + ALP121*zs+ALP021)*zt & 1232 & + (ALP211*zs+ALP111)*zs+ALP011)*zt & 1233 & + ((ALP301*zs+ALP201)*zs+ALP101)*zs+ALP001 1234 ! 1235 zn0 = ((((ALP050*zt & 1236 & + ALP140*zs+ALP040)*zt & 1237 & + (ALP230*zs+ALP130)*zs+ALP030)*zt & 1238 & + ((ALP320*zs+ALP220)*zs+ALP120)*zs+ALP020)*zt & 1239 & + (((ALP410*zs+ALP310)*zs+ALP210)*zs+ALP110)*zs+ALP010)*zt & 1240 & + ((((ALP500*zs+ALP400)*zs+ALP300)*zs+ALP200)*zs+ALP100)*zs+ALP000 1241 ! 1242 zn = ( ( zn3 * zh + zn2 ) * zh + zn1 ) * zh + zn0 1243 ! 1244 pab(ji,jj,jp_tem) = zn * r1_rho0 1245 ! 1246 ! beta 1247 zn3 = BET003 1248 ! 1249 zn2 = BET012*zt + BET102*zs+BET002 1250 ! 1251 zn1 = ((BET031*zt & 1252 & + BET121*zs+BET021)*zt & 1253 & + (BET211*zs+BET111)*zs+BET011)*zt & 1254 & + ((BET301*zs+BET201)*zs+BET101)*zs+BET001 1255 ! 1256 zn0 = ((((BET050*zt & 1257 & + BET140*zs+BET040)*zt & 1258 & + (BET230*zs+BET130)*zs+BET030)*zt & 1259 & + ((BET320*zs+BET220)*zs+BET120)*zs+BET020)*zt & 1260 & + (((BET410*zs+BET310)*zs+BET210)*zs+BET110)*zs+BET010)*zt & 1261 & + ((((BET500*zs+BET400)*zs+BET300)*zs+BET200)*zs+BET100)*zs+BET000 1262 ! 1263 zn = ( ( zn3 * zh + zn2 ) * zh + zn1 ) * zh + zn0 1264 ! 864 1265 pab(ji,jj,jp_sal) = zn / zs * r1_rho0 865 1266 ! … … 997 1398 !! 998 1399 INTEGER , INTENT(in ) :: Kmm ! time level index 999 REAL( wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT(in ) :: pts ! pot. temperature and salinity [Celsius,psu]1400 REAL(dp), DIMENSION(jpi,jpj,jpk,jpts), INTENT(in ) :: pts ! pot. temperature and salinity [Celsius,psu] 1000 1401 REAL(wp), DIMENSION(:,:,:,:) , INTENT(in ) :: pab ! thermal/haline expansion coef. [Celsius-1,psu-1] 1001 1402 REAL(wp), DIMENSION(:,:,:) , INTENT( out) :: pn2 ! Brunt-Vaisala frequency squared [1/s^2] … … 1021 1422 INTEGER , INTENT(in ) :: Kmm ! time level index 1022 1423 INTEGER , INTENT(in ) :: ktab, ktn2 1023 REAL( wp), DIMENSION(jpi,jpj, jpk,jpts), INTENT(in ) :: pts ! pot. temperature and salinity [Celsius,psu]1424 REAL(dp), DIMENSION(jpi,jpj, jpk,jpts), INTENT(in ) :: pts ! pot. temperature and salinity [Celsius,psu] 1024 1425 REAL(wp), DIMENSION(A2D_T(ktab),JPK,JPTS), INTENT(in ) :: pab ! thermal/haline expansion coef. [Celsius-1,psu-1] 1025 1426 REAL(wp), DIMENSION(A2D_T(ktn2),JPK ), INTENT( out) :: pn2 ! Brunt-Vaisala frequency squared [1/s^2] -
NEMO/branches/2021/dev_r14116_HPC-10_mcastril_Mixed_Precision_implementation/src/OCE/TRA/tranpc.F90
r14986 r15066 101 101 ! 102 102 CALL eos_rab( CASTWP(pts(:,:,:,:,Kaa)), zab, Kmm ) ! after alpha and beta (given on T-points) 103 CALL bn2 ( CASTWP(pts(:,:,:,:,Kaa)), zab, zn2, Kmm ) ! after Brunt-Vaisala (given on W-points)103 CALL bn2 ( pts(:,:,:,:,Kaa), zab, zn2, Kmm ) ! after Brunt-Vaisala (given on W-points) 104 104 ! 105 105 IF( .NOT. l_istiled .OR. ntile == 1 ) nnpcc = 0 ! Do only on the first tile
Note: See TracChangeset
for help on using the changeset viewer.
