Changeset 13409
- Timestamp:
- 2020-08-17T15:28:54+02:00 (3 years ago)
- Location:
- NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE
- Files:
-
- 36 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/ASM/asminc.F90
r12958 r13409 518 518 INTEGER :: it 519 519 REAL(wp) :: zincwgt ! IAU weight for current time step 520 REAL (wp), DIMENSION(jpi,jpj,jpk) :: fzptnz ! 3d freezing point values520 REAL(wp), DIMENSION(A2D,jpk) :: fzptnz ! 3d freezing point values 521 521 !!---------------------------------------------------------------------- 522 522 ! 523 523 ! freezing point calculation taken from oc_fz_pt (but calculated for all depths) 524 524 ! used to prevent the applied increments taking the temperature below the local freezing point 525 DO jk = 1, jpkm1 526 CALL eos_fzp( pts(:,:,jk,jp_sal,Kmm), fzptnz(:,:,jk), gdept(:,:,jk,Kmm) ) 527 END DO 525 ! TODO: NOT TESTED- logical is forced to False 526 IF( ln_temnofreeze ) THEN 527 DO jk = 1, jpkm1 528 CALL eos_fzp( pts(:,:,jk,jp_sal,Kmm), fzptnz(:,:,jk), gdept(:,:,jk,Kmm) ) 529 END DO 530 ENDIF 528 531 ! 529 532 ! !-------------------------------------- … … 536 539 zincwgt = wgtiau(it) / rn_Dt ! IAU weight for the current time step 537 540 ! 538 IF(lwp) THEN 539 WRITE(numout,*) 540 WRITE(numout,*) 'tra_asm_inc : Tracer IAU at time step = ', kt,' with IAU weight = ', wgtiau(it) 541 WRITE(numout,*) '~~~~~~~~~~~~' 541 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 542 IF(lwp) THEN 543 WRITE(numout,*) 544 WRITE(numout,*) 'tra_asm_inc : Tracer IAU at time step = ', kt,' with IAU weight = ', wgtiau(it) 545 WRITE(numout,*) '~~~~~~~~~~~~' 546 ENDIF 542 547 ENDIF 543 548 ! 544 549 ! Update the tracer tendencies 550 ! TODO: NOT TESTED- logical is forced to False 545 551 DO jk = 1, jpkm1 546 552 IF (ln_temnofreeze) THEN 547 553 ! Do not apply negative increments if the temperature will fall below freezing 548 WHERE(t_bkginc( :,:,jk) > 0.0_wp .OR. &549 & pts( :,:,jk,jp_tem,Kmm) + pts(:,:,jk,jp_tem,Krhs) + t_bkginc(:,:,jk) * wgtiau(it) > fzptnz(:,:,jk) )550 pts( :,:,jk,jp_tem,Krhs) = pts(:,:,jk,jp_tem,Krhs) + t_bkginc(:,:,jk) * zincwgt554 WHERE(t_bkginc(A2D0,jk) > 0.0_wp .OR. & 555 & pts(A2D0,jk,jp_tem,Kmm) + pts(A2D0,jk,jp_tem,Krhs) + t_bkginc(A2D0,jk) * wgtiau(it) > fzptnz(:,:,jk) ) 556 pts(A2D0,jk,jp_tem,Krhs) = pts(A2D0,jk,jp_tem,Krhs) + t_bkginc(A2D0,jk) * zincwgt 551 557 END WHERE 552 558 ELSE 553 pts(:,:,jk,jp_tem,Krhs) = pts(:,:,jk,jp_tem,Krhs) + t_bkginc(:,:,jk) * zincwgt 559 DO_2D_00_00 560 pts(ji,jj,jk,jp_tem,Krhs) = pts(ji,jj,jk,jp_tem,Krhs) + t_bkginc(ji,jj,jk) * zincwgt 561 END_2D 554 562 ENDIF 555 563 IF (ln_salfix) THEN 556 564 ! Do not apply negative increments if the salinity will fall below a specified 557 565 ! minimum value salfixmin 558 WHERE(s_bkginc( :,:,jk) > 0.0_wp .OR. &559 & pts( :,:,jk,jp_sal,Kmm) + pts(:,:,jk,jp_sal,Krhs) + s_bkginc(:,:,jk) * wgtiau(it) > salfixmin )560 pts( :,:,jk,jp_sal,Krhs) = pts(:,:,jk,jp_sal,Krhs) + s_bkginc(:,:,jk) * zincwgt566 WHERE(s_bkginc(A2D0,jk) > 0.0_wp .OR. & 567 & pts(A2D0,jk,jp_sal,Kmm) + pts(A2D0,jk,jp_sal,Krhs) + s_bkginc(A2D0,jk) * wgtiau(it) > salfixmin ) 568 pts(A2D0,jk,jp_sal,Krhs) = pts(A2D0,jk,jp_sal,Krhs) + s_bkginc(A2D0,jk) * zincwgt 561 569 END WHERE 562 570 ELSE 563 pts(:,:,jk,jp_sal,Krhs) = pts(:,:,jk,jp_sal,Krhs) + s_bkginc(:,:,jk) * zincwgt 571 DO_2D_00_00 572 pts(ji,jj,jk,jp_sal,Krhs) = pts(ji,jj,jk,jp_sal,Krhs) + s_bkginc(ji,jj,jk) * zincwgt 573 END_2D 564 574 ENDIF 565 575 END DO … … 567 577 ENDIF 568 578 ! 569 IF ( kt == nitiaufin_r + 1 ) THEN ! For bias crcn to work 570 DEALLOCATE( t_bkginc ) 571 DEALLOCATE( s_bkginc ) 579 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only on the last tile 580 IF ( kt == nitiaufin_r + 1 ) THEN ! For bias crcn to work 581 DEALLOCATE( t_bkginc ) 582 DEALLOCATE( s_bkginc ) 583 ENDIF 572 584 ENDIF 573 585 ! !-------------------------------------- … … 580 592 ! 581 593 ! Initialize the now fields with the background + increment 594 ! TODO: NOT TESTED- logical is forced to False 582 595 IF (ln_temnofreeze) THEN 583 596 ! Do not apply negative increments if the temperature will fall below freezing 584 WHERE( t_bkginc( :,:,:) > 0.0_wp .OR. pts(:,:,:,jp_tem,Kmm) + t_bkginc(:,:,:) > fzptnz(:,:,:) )585 pts( :,:,:,jp_tem,Kmm) = t_bkg(:,:,:) + t_bkginc(:,:,:)597 WHERE( t_bkginc(A2D0,:) > 0.0_wp .OR. pts(A2D0,:,jp_tem,Kmm) + t_bkginc(A2D0,:) > fzptnz(:,:,:) ) 598 pts(A2D0,:,jp_tem,Kmm) = t_bkg(A2D0,:) + t_bkginc(A2D0,:) 586 599 END WHERE 587 600 ELSE 588 pts(:,:,:,jp_tem,Kmm) = t_bkg(:,:,:) + t_bkginc(:,:,:) 601 DO_3D_00_00( 1, jpk ) 602 pts(ji,jj,jk,jp_tem,Kmm) = t_bkg(ji,jj,jk) + t_bkginc(ji,jj,jk) 603 END_3D 589 604 ENDIF 590 605 IF (ln_salfix) THEN 591 606 ! Do not apply negative increments if the salinity will fall below a specified 592 607 ! minimum value salfixmin 593 WHERE( s_bkginc( :,:,:) > 0.0_wp .OR. pts(:,:,:,jp_sal,Kmm) + s_bkginc(:,:,:) > salfixmin )594 pts( :,:,:,jp_sal,Kmm) = s_bkg(:,:,:) + s_bkginc(:,:,:)608 WHERE( s_bkginc(A2D0,:) > 0.0_wp .OR. pts(A2D0,:,jp_sal,Kmm) + s_bkginc(A2D0,:) > salfixmin ) 609 pts(A2D0,:,jp_sal,Kmm) = s_bkg(A2D0,:) + s_bkginc(A2D0,:) 595 610 END WHERE 596 611 ELSE 597 pts(:,:,:,jp_sal,Kmm) = s_bkg(:,:,:) + s_bkginc(:,:,:) 598 ENDIF 599 600 pts(:,:,:,:,Kbb) = pts(:,:,:,:,Kmm) ! Update before fields 612 DO_3D_00_00( 1, jpk ) 613 pts(ji,jj,jk,jp_sal,Kmm) = s_bkg(ji,jj,jk) + s_bkginc(ji,jj,jk) 614 END_3D 615 ENDIF 616 617 DO_3D_00_00( 1, jpk ) 618 pts(ji,jj,jk,:,Kbb) = pts(ji,jj,jk,:,Kmm) ! Update before fields 619 END_3D 601 620 602 621 CALL eos( pts(:,:,:,:,Kbb), rhd, rhop, gdept_0(:,:,:) ) ! Before potential and in situ densities … … 605 624 !!gm 606 625 607 IF( ln_zps .AND. .NOT. lk_c1d .AND. .NOT. ln_isfcav) & 608 & CALL zps_hde ( kt, Kmm, jpts, pts(:,:,:,:,Kbb), gtsu, gtsv, & ! Partial steps: before horizontal gradient 609 & rhd, gru , grv ) ! of t, s, rd at the last ocean level 610 IF( ln_zps .AND. .NOT. lk_c1d .AND. ln_isfcav) & 611 & CALL zps_hde_isf( nit000, Kmm, jpts, pts(:,:,:,:,Kbb), gtsu, gtsv, gtui, gtvi, & ! Partial steps for top cell (ISF) 612 & rhd, gru , grv , grui, grvi ) ! of t, s, rd at the last ocean level 613 614 DEALLOCATE( t_bkginc ) 615 DEALLOCATE( s_bkginc ) 616 DEALLOCATE( t_bkg ) 617 DEALLOCATE( s_bkg ) 626 ! TEMP: This change not necessary after extra haloes development (lbc_lnk removed from zps_hde*) 627 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 628 IF( ln_zps .AND. .NOT. lk_c1d .AND. .NOT. ln_isfcav) & 629 & CALL zps_hde ( kt, Kmm, jpts, pts(:,:,:,:,Kbb), gtsu, gtsv, & ! Partial steps: before horizontal gradient 630 & rhd, gru , grv ) ! of t, s, rd at the last ocean level 631 IF( ln_zps .AND. .NOT. lk_c1d .AND. ln_isfcav) & 632 & CALL zps_hde_isf( nit000, Kmm, jpts, pts(:,:,:,:,Kbb), gtsu, gtsv, gtui, gtvi, & ! Partial steps for top cell (ISF) 633 & rhd, gru , grv , grui, grvi ) ! of t, s, rd at the last ocean level 634 ENDIF 635 636 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only on the last tile 637 DEALLOCATE( t_bkginc ) 638 DEALLOCATE( s_bkginc ) 639 DEALLOCATE( t_bkg ) 640 DEALLOCATE( s_bkg ) 641 ENDIF 642 ! 618 643 ENDIF 619 644 ! 620 645 ENDIF 646 ! TODO: NOT TESTED- logical is forced to False 621 647 ! Perhaps the following call should be in step 622 648 IF ( ln_seaiceinc ) CALL seaice_asm_inc ( kt ) ! apply sea ice concentration increment … … 825 851 INTEGER, INTENT(in), OPTIONAL :: kindic ! flag for disabling the deallocation 826 852 ! 853 INTEGER :: ji, jj 827 854 INTEGER :: it 828 855 REAL(wp) :: zincwgt ! IAU weight for current time step 829 856 #if defined key_si3 830 REAL(wp), DIMENSION( jpi,jpj) :: zofrld, zohicif, zseaicendg, zhicifinc857 REAL(wp), DIMENSION(A2D) :: zofrld, zohicif, zseaicendg, zhicifinc 831 858 REAL(wp) :: zhicifmin = 0.5_wp ! ice minimum depth in metres 832 859 #endif … … 843 870 ! note this is not a tendency so should not be divided by rn_Dt (as with the tracer and other increments) 844 871 ! 845 IF(lwp) THEN 846 WRITE(numout,*) 847 WRITE(numout,*) 'seaice_asm_inc : sea ice conc IAU at time step = ', kt,' with IAU weight = ', wgtiau(it) 848 WRITE(numout,*) '~~~~~~~~~~~~' 872 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 873 IF(lwp) THEN 874 WRITE(numout,*) 875 WRITE(numout,*) 'seaice_asm_inc : sea ice conc IAU at time step = ', kt,' with IAU weight = ', wgtiau(it) 876 WRITE(numout,*) '~~~~~~~~~~~~' 877 ENDIF 849 878 ENDIF 850 879 ! … … 852 881 ! 853 882 #if defined key_si3 854 zofrld (:,:) = 1._wp - at_i(:,:) 855 zohicif(:,:) = hm_i(:,:) 856 ! 857 at_i (:,:) = 1. - MIN( MAX( 1.-at_i (:,:) - seaice_bkginc(:,:) * zincwgt, 0.0_wp), 1.0_wp) 858 at_i_b(:,:) = 1. - MIN( MAX( 1.-at_i_b(:,:) - seaice_bkginc(:,:) * zincwgt, 0.0_wp), 1.0_wp) 859 fr_i(:,:) = at_i(:,:) ! adjust ice fraction 860 ! 861 zseaicendg(:,:) = zofrld(:,:) - (1. - at_i(:,:)) ! find out actual sea ice nudge applied 883 DO_2D_00_00 884 zofrld (ji,jj) = 1._wp - at_i(ji,jj) 885 zohicif(ji,jj) = hm_i(ji,jj) 886 ! 887 at_i (ji,jj) = 1. - MIN( MAX( 1.-at_i (ji,jj) - seaice_bkginc(ji,jj) * zincwgt, 0.0_wp), 1.0_wp) 888 at_i_b(ji,jj) = 1. - MIN( MAX( 1.-at_i_b(ji,jj) - seaice_bkginc(ji,jj) * zincwgt, 0.0_wp), 1.0_wp) 889 fr_i(ji,jj) = at_i(ji,jj) ! adjust ice fraction 890 ! 891 zseaicendg(ji,jj) = zofrld(ji,jj) - (1. - at_i(ji,jj)) ! find out actual sea ice nudge applied 892 END_2D 862 893 ! 863 894 ! Nudge sea ice depth to bring it up to a required minimum depth 864 WHERE( zseaicendg(:,:) > 0.0_wp .AND. hm_i( :,:) < zhicifmin )865 zhicifinc(:,:) = (zhicifmin - hm_i( :,:)) * zincwgt895 WHERE( zseaicendg(:,:) > 0.0_wp .AND. hm_i(A2D0) < zhicifmin ) 896 zhicifinc(:,:) = (zhicifmin - hm_i(A2D0)) * zincwgt 866 897 ELSEWHERE 867 898 zhicifinc(:,:) = 0.0_wp … … 869 900 ! 870 901 ! nudge ice depth 871 hm_i (:,:) = hm_i (:,:) + zhicifinc(:,:) 902 DO_2D_00_00 903 hm_i (ji,jj) = hm_i (ji,jj) + zhicifinc(ji,jj) 904 END_2D 872 905 ! 873 906 ! seaice salinity balancing (to add) … … 876 909 #if defined key_cice && defined key_asminc 877 910 ! Sea-ice : CICE case. Pass ice increment tendency into CICE 878 ndaice_da(:,:) = seaice_bkginc(:,:) * zincwgt / rn_Dt 879 #endif 880 ! 881 IF ( kt == nitiaufin_r ) THEN 882 DEALLOCATE( seaice_bkginc ) 911 DO_2D_00_00 912 ndaice_da(ji,jj) = seaice_bkginc(ji,jj) * zincwgt / rn_Dt 913 END_2D 914 #endif 915 ! 916 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only on the last tile 917 IF ( kt == nitiaufin_r ) THEN 918 DEALLOCATE( seaice_bkginc ) 919 ENDIF 883 920 ENDIF 884 921 ! … … 886 923 ! 887 924 #if defined key_cice && defined key_asminc 888 ndaice_da(:,:) = 0._wp ! Sea-ice : CICE case. Zero ice increment tendency into CICE 925 DO_2D_00_00 926 ndaice_da(ji,jj) = 0._wp ! Sea-ice : CICE case. Zero ice increment tendency into CICE 927 END_2D 889 928 #endif 890 929 ! … … 901 940 ! 902 941 #if defined key_si3 903 zofrld (:,:) = 1._wp - at_i(:,:) 904 zohicif(:,:) = hm_i(:,:) 905 ! 906 ! Initialize the now fields the background + increment 907 at_i(:,:) = 1. - MIN( MAX( 1.-at_i(:,:) - seaice_bkginc(:,:), 0.0_wp), 1.0_wp) 908 at_i_b(:,:) = at_i(:,:) 909 fr_i(:,:) = at_i(:,:) ! adjust ice fraction 910 ! 911 zseaicendg(:,:) = zofrld(:,:) - (1. - at_i(:,:)) ! find out actual sea ice nudge applied 942 DO_2D_00_00 943 zofrld (ji,jj) = 1._wp - at_i(ji,jj) 944 zohicif(ji,jj) = hm_i(ji,jj) 945 ! 946 ! Initialize the now fields the background + increment 947 at_i(ji,jj) = 1. - MIN( MAX( 1.-at_i(ji,jj) - seaice_bkginc(ji,jj), 0.0_wp), 1.0_wp) 948 at_i_b(ji,jj) = at_i(ji,jj) 949 fr_i(ji,jj) = at_i(ji,jj) ! adjust ice fraction 950 ! 951 zseaicendg(ji,jj) = zofrld(ji,jj) - (1. - at_i(ji,jj)) ! find out actual sea ice nudge applied 952 END_2D 912 953 ! 913 954 ! Nudge sea ice depth to bring it up to a required minimum depth 914 WHERE( zseaicendg(:,:) > 0.0_wp .AND. hm_i( :,:) < zhicifmin )915 zhicifinc(:,:) = zhicifmin - hm_i( :,:)955 WHERE( zseaicendg(:,:) > 0.0_wp .AND. hm_i(A2D0) < zhicifmin ) 956 zhicifinc(:,:) = zhicifmin - hm_i(A2D0) 916 957 ELSEWHERE 917 958 zhicifinc(:,:) = 0.0_wp … … 919 960 ! 920 961 ! nudge ice depth 921 hm_i (:,:) = hm_i (:,:) + zhicifinc(:,:) 962 DO_2D_00_00 963 hm_i(ji,jj) = hm_i (ji,jj) + zhicifinc(ji,jj) 964 END_2D 922 965 ! 923 966 ! seaice salinity balancing (to add) … … 926 969 #if defined key_cice && defined key_asminc 927 970 ! Sea-ice : CICE case. Pass ice increment tendency into CICE 928 ndaice_da(:,:) = seaice_bkginc(:,:) / rn_Dt 929 #endif 930 IF ( .NOT. PRESENT(kindic) ) THEN 931 DEALLOCATE( seaice_bkginc ) 932 END IF 971 DO_2D_00_00 972 ndaice_da(ji,jj) = seaice_bkginc(ji,jj) / rn_Dt 973 END_2D 974 #endif 975 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only on the last tile 976 IF ( .NOT. PRESENT(kindic) ) THEN 977 DEALLOCATE( seaice_bkginc ) 978 END IF 979 ENDIF 933 980 ! 934 981 ELSE 935 982 ! 936 983 #if defined key_cice && defined key_asminc 937 ndaice_da(:,:) = 0._wp ! Sea-ice : CICE case. Zero ice increment tendency into CICE 984 DO_2D_00_00 985 ndaice_da(ji,jj) = 0._wp ! Sea-ice : CICE case. Zero ice increment tendency into CICE 986 END_2D 938 987 #endif 939 988 ! -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/BDY/bdytra.F90
r12377 r13409 13 13 !!---------------------------------------------------------------------- 14 14 USE oce ! ocean dynamics and tracers variables 15 USE dom_oce ! ocean space and time domain variables 15 USE dom_oce ! ocean space and time domain variables 16 16 USE bdy_oce ! ocean open boundary conditions 17 17 USE bdylib ! for orlanski library routines … … 157 157 INTEGER :: ib_bdy ! Loop index 158 158 !!---------------------------------------------------------------------- 159 ! TODO: TO BE TILED 160 ! TODO: NOT TESTED- requires bdy 161 ! NOTE: Tiling these BDY loops is nontrivial; IF statements to check whether a point is in the current tile won't work (will be for every ib, every tile). The idx_bdy structure might require modifying to include a %nblen and list of ib indices for the current tile. 162 IF( ntile /= 0 .AND. ntile /= 1 ) RETURN ! Do only for the full domain 159 163 ! 160 164 IF( ln_timing ) CALL timing_start('bdy_tra_dmp') -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/DIA/diaar5.F90
r13054 r13409 303 303 END SUBROUTINE dia_ar5 304 304 305 ! T ODO: These changes and lbc_lnk not necessary if using XIOS (subdomain support, will not output haloes)305 ! TEMP: These changes and lbc_lnk not necessary if using XIOS (subdomain support, will not output haloes) 306 306 SUBROUTINE dia_ar5_hst( ktra, cptr, puflx, pvflx ) 307 307 !!---------------------------------------------------------------------- -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/DIA/diaptr.F90
r13054 r13409 22 22 USE oce ! ocean dynamics and active tracers 23 23 USE dom_oce ! ocean space and time domain 24 ! TEMP: Possibly not necessary if using XIOS (if cumulative axis operations are possible) 24 25 USE domain, ONLY : dom_tile 25 26 USE phycst ! physical constants … … 74 75 CONTAINS 75 76 76 ! T ODO: Most changes and some code in this module not necessary if using XIOS (subdomain support, axis operations)77 ! TEMP: Most changes and some code in this module not necessary if using XIOS (subdomain support, axis operations) 77 78 SUBROUTINE dia_ptr( kt, Kmm, pvtr ) 78 79 !!---------------------------------------------------------------------- … … 81 82 INTEGER , INTENT(in) :: kt ! ocean time-step index 82 83 INTEGER , INTENT(in) :: Kmm ! time level index 83 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL :: pvtr ! j-effective transport 84 REAL(wp), DIMENSION(A2D,jpk) , INTENT(in), OPTIONAL :: pvtr ! j-effective transport 85 !!---------------------------------------------------------------------- 86 ! 87 IF( ln_timing ) CALL timing_start('dia_ptr') 88 89 IF( kt == nit000 .AND. ll_init ) CALL dia_ptr_init 90 ! 91 IF( l_diaptr ) THEN 92 ! Calculate zonal integrals 93 IF( PRESENT( pvtr ) ) THEN 94 CALL dia_ptr_zint( Kmm, pvtr ) 95 ELSE 96 CALL dia_ptr_zint( Kmm ) 97 ENDIF 98 99 ! Calculate diagnostics only when zonal integrals have finished 100 IF( ntile == 0 .OR. ntile == nijtile ) CALL dia_ptr_iom(kt, Kmm, pvtr) 101 ENDIF 102 103 IF( ln_timing ) CALL timing_stop('dia_ptr') 104 ! 105 END SUBROUTINE dia_ptr 106 107 108 SUBROUTINE dia_ptr_iom( kt, Kmm, pvtr ) 109 !!---------------------------------------------------------------------- 110 !! *** ROUTINE dia_ptr_iom *** 111 !!---------------------------------------------------------------------- 112 !! ** Purpose : Calculate diagnostics and send to XIOS 113 !!---------------------------------------------------------------------- 114 INTEGER , INTENT(in) :: kt ! ocean time-step index 115 INTEGER , INTENT(in) :: Kmm ! time level index 116 REAL(wp), DIMENSION(A2D,jpk) , INTENT(in), OPTIONAL :: pvtr ! j-effective transport 84 117 ! 85 118 INTEGER :: ji, jj, jk, jn ! dummy loop indices 86 INTEGER :: itile87 119 REAL(wp), DIMENSION(jpi,jpj) :: z2d ! 2D workspace 88 120 REAL(wp), DIMENSION(jpj) :: zvsum, ztsum, zssum ! 1D workspace … … 95 127 REAL(wp), DIMENSION(jpi,jpj,nptr) :: z3dtr ! i-mean T and S, j-Stream-Function 96 128 !!---------------------------------------------------------------------- 97 !98 IF( ln_timing ) CALL timing_start('dia_ptr')99 100 IF( kt == nit000 .AND. ll_init ) CALL dia_ptr_init101 !102 IF( .NOT. l_diaptr ) RETURN103 104 ! Calculate zonal integrals105 IF( PRESENT( pvtr ) ) THEN106 CALL dia_ptr_zint( Kmm, pvtr )107 ELSE108 CALL dia_ptr_zint( Kmm )109 ENDIF110 111 ! Calculate diagnostics only when zonal integrals have finished112 IF( ntile /= 0 .AND. ntile /= nijtile ) RETURN113 129 114 130 IF( PRESENT( pvtr ) ) THEN … … 155 171 156 172 IF( iom_use( 'sopstbtr' ) .OR. iom_use( 'sophtbtr' ) ) THEN 157 ! Calculate barotropic heat and salt transport here 173 ! Calculate barotropic heat and salt transport here 158 174 DO jn = 1, nptr 159 175 sjk(:,1,jn) = SUM( pvtr_int(:,:,jp_msk,jn), 2 ) … … 182 198 ENDDO 183 199 CALL iom_put( 'sopstbtr', z3dtr ) 184 ENDIF 200 ENDIF 185 201 ! 186 202 hstr_ove(:,:,:) = 0._wp ! Zero before next timestep … … 217 233 ! 218 234 ! ! Advective and diffusive heat and salt transport 219 IF( iom_use( 'sophtadv' ) .OR. iom_use( 'sopstadv' ) ) THEN 220 ! 235 IF( iom_use( 'sophtadv' ) .OR. iom_use( 'sopstadv' ) ) THEN 236 ! 221 237 DO jn = 1, nptr 222 238 z3dtr(1,:,jn) = hstr_adv(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) … … 235 251 ENDIF 236 252 ! 237 IF( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) THEN 238 ! 253 IF( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) THEN 254 ! 239 255 DO jn = 1, nptr 240 256 z3dtr(1,:,jn) = hstr_ldf(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) … … 253 269 ENDIF 254 270 ! 255 IF( iom_use( 'sophteiv' ) .OR. iom_use( 'sopsteiv' ) ) THEN 256 ! 271 IF( iom_use( 'sophteiv' ) .OR. iom_use( 'sopsteiv' ) ) THEN 272 ! 257 273 DO jn = 1, nptr 258 274 z3dtr(1,:,jn) = hstr_eiv(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) … … 288 304 ENDIF 289 305 ! 290 ! TODO: Possibly not necessary if using XIOS (if cumulative axis operations are possible) 306 ! TEMP: Possibly not necessary if using XIOS (if cumulative axis operations are possible) 307 ! TODO: NOT TESTED- hangs on iom_get_var 291 308 IF( iom_use( 'uocetr_vsum_cumul' ) ) THEN 292 itile = ntile 293 CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 309 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 294 310 CALL iom_get_var( 'uocetr_vsum_op', z2d ) ! get uocetr_vsum_op from xml 295 z2d(:,:) = ptr_ci_2d( z2d(:,:) ) 311 z2d(:,:) = ptr_ci_2d( z2d(:,:) ) 296 312 CALL iom_put( 'uocetr_vsum_cumul', z2d ) 297 IF( ntile /= itile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = itile ) ! Revert to tile domain313 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = nijtile ) ! Revert to tile domain 298 314 ENDIF 299 315 ! … … 304 320 pzon_int(:,:,:,:) = 0._wp 305 321 ENDIF 306 ! 307 IF( ln_timing ) CALL timing_stop('dia_ptr') 308 ! 309 END SUBROUTINE dia_ptr 322 END SUBROUTINE dia_ptr_iom 310 323 311 324 … … 330 343 REAL(wp) :: zsfc, zvfc ! i-k surface area 331 344 INTEGER :: ji, jj, jk, jn ! dummy loop indices 345 !!---------------------------------------------------------------------- 332 346 333 347 IF( PRESENT( pvtr ) ) THEN … … 442 456 IF( lk_mpp ) CALL mpp_ini_znl( numout ) ! Define MPI communicator for zonal sum 443 457 444 btmsk(:,:,1) = tmask_i(:,:) 458 btmsk(:,:,:) = 0._wp 459 btmsk(:,:,1) = tmask_i(:,:) 445 460 CALL iom_open( 'subbasins', inum, ldstop = .FALSE. ) 446 461 CALL iom_get( inum, jpdom_global, 'atlmsk', btmsk(:,:,2) ) ! Atlantic basin -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/DOM/domain.F90
r12979 r13409 295 295 !! - nijtile : total number of tiles 296 296 !!---------------------------------------------------------------------- 297 INTEGER, INTENT(out) :: ktsi, ktsj, ktei, ktej ! Tile domain indices297 TYPE(tiledom), INTENT(out) :: ktsi, ktsj, ktei, ktej ! Tile domain indices 298 298 INTEGER, INTENT(in), OPTIONAL :: ktile ! Tile number 299 299 INTEGER :: jt ! dummy loop argument … … 302 302 IF( PRESENT(ktile) .AND. ln_tile ) THEN 303 303 ntile = ktile ! Set domain indices for tile 304 ktsi = ntsi_a(ktile)305 ktsj = ntsj_a(ktile)306 ktei = ntei_a(ktile)307 ktej = ntej_a(ktile)304 ktsi%h0 = ntsi_a(ktile) 305 ktsj%h0 = ntsj_a(ktile) 306 ktei%h0 = ntei_a(ktile) 307 ktej%h0 = ntej_a(ktile) 308 308 ELSE 309 309 ntile = 0 ! Initialise to full domain 310 310 nijtile = 1 311 ktsi = Nis0312 ktsj = Njs0313 ktei = Nie0314 ktej = Nje0311 ktsi%h0 = Nis0 312 ktsj%h0 = Njs0 313 ktei%h0 = Nie0 314 ktej%h0 = Nje0 315 315 316 316 IF( ln_tile ) THEN ! Calculate tile domain indices … … 323 323 ALLOCATE( ntsi_a(0:nijtile), ntsj_a(0:nijtile), ntei_a(0:nijtile), ntej_a(0:nijtile) ) 324 324 325 ntsi_a(0) = ktsi ! Full domain326 ntsj_a(0) = ktsj 327 ntei_a(0) = ktei 328 ntej_a(0) = ktej 325 ntsi_a(0) = ktsi%h0 ! Full domain 326 ntsj_a(0) = ktsj%h0 327 ntei_a(0) = ktei%h0 328 ntej_a(0) = ktej%h0 329 329 330 330 DO jt = 1, nijtile ! Tile domains … … 353 353 ELSE 354 354 WRITE(numout,*) 'No domain tiling' 355 WRITE(numout,*) ' i indices =', ktsi , ':', ktei356 WRITE(numout,*) ' j indices =', ktsj , ':', ktej355 WRITE(numout,*) ' i indices =', ktsi%h0, ':', ktei%h0 356 WRITE(numout,*) ' j indices =', ktsj%h0, ':', ktej%h0 357 357 ENDIF 358 358 ENDIF 359 359 ENDIF 360 361 ktsi%h1 = ktsi%h0 - 1 362 ktsj%h1 = ktsj%h0 - 1 363 ktei%h1 = ktei%h0 + 1 364 ktej%h1 = ktej%h0 + 1 365 ktsi%h2 = ktsi%h0 - nn_hls 366 ktsj%h2 = ktsj%h0 - nn_hls 367 ktei%h2 = ktei%h0 + nn_hls 368 ktej%h2 = ktej%h0 + nn_hls 360 369 END SUBROUTINE dom_tile 361 370 -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/DOM/domutl.F90
r12807 r13409 21 21 PRIVATE 22 22 23 INTERFACE is_tile 24 MODULE PROCEDURE is_tile_2d, is_tile_3d, is_tile_4d 25 END INTERFACE is_tile 26 23 27 PUBLIC dom_ngb ! routine called in iom.F90 module 24 28 PUBLIC dom_uniq ! Called by dommsk and domwri 29 PUBLIC is_tile 25 30 26 31 !!---------------------------------------------------------------------- … … 115 120 ! 116 121 END SUBROUTINE dom_uniq 117 122 123 124 PURE FUNCTION is_tile_2d( pt ) 125 !! 126 REAL(wp), DIMENSION(:,:), INTENT(in) :: pt 127 INTEGER :: is_tile_2d 128 !! 129 IF( ln_tile .AND. SIZE(pt, 1) < jpi ) THEN 130 is_tile_2d = 1 131 ELSE 132 is_tile_2d = 0 133 ENDIF 134 END FUNCTION is_tile_2d 135 136 137 PURE FUNCTION is_tile_3d( pt ) 138 !! 139 REAL(wp), DIMENSION(:,:,:), INTENT(in) :: pt 140 INTEGER :: is_tile_3d 141 !! 142 IF( ln_tile .AND. SIZE(pt, 1) < jpi ) THEN 143 is_tile_3d = 1 144 ELSE 145 is_tile_3d = 0 146 ENDIF 147 END FUNCTION is_tile_3d 148 149 150 PURE FUNCTION is_tile_4d( pt ) 151 !! 152 REAL(wp), DIMENSION(:,:,:,:), INTENT(in) :: pt 153 INTEGER :: is_tile_4d 154 !! 155 IF( ln_tile .AND. SIZE(pt, 1) < jpi ) THEN 156 is_tile_4d = 1 157 ELSE 158 is_tile_4d = 0 159 ENDIF 160 END FUNCTION is_tile_4d 161 118 162 !!====================================================================== 119 163 END MODULE domutl -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/DOM/dtatsd.F90
r12962 r13409 18 18 USE phycst ! physical constants 19 19 USE dom_oce ! ocean space and time domain 20 USE domain, ONLY : dom_tile 20 21 USE fldread ! read input fields 21 22 ! … … 135 136 !! ** Action : ptsd T-S data on medl mesh and interpolated at time-step kt 136 137 !!---------------------------------------------------------------------- 137 INTEGER 138 REAL(wp), DIMENSION( jpi,jpj,jpk,jpts), INTENT( out) :: ptsd ! T & S data138 INTEGER , INTENT(in ) :: kt ! ocean time-step 139 REAL(wp), DIMENSION(A2D,jpk,jpts), INTENT( out) :: ptsd ! T & S data 139 140 ! 140 141 INTEGER :: ji, jj, jk, jl, jkk ! dummy loop indicies … … 144 145 !!---------------------------------------------------------------------- 145 146 ! 146 CALL fld_read( kt, 1, sf_tsd ) !== read T & S data at kt time step ==! 147 ! 148 ! 147 ! NOTE: We can read only part of the global field (kstart, kcount in iom_get & some other mods for tiling) but the lbc_lnk in iom_get can't be overcome. Opening & reading from the file multiple times is also inefficient. Therefore, we call fld_read once for the full domain (since sf_tsd is never deallocated) and do the remaining steps for the tile domain. 148 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only for the full domain 149 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 150 CALL fld_read( kt, 1, sf_tsd ) !== read T & S data at kt time step ==! 151 ! 152 ! 153 ! NOTE: Do this for the full domain too, since it is modifying sf_tsd in-place and only for ORCA2 149 154 !!gm This should be removed from the code ===>>>> T & S files has to be changed 150 ! 151 ! !== ORCA_R2 configuration and T & S damping ==! 152 IF( cn_cfg == "orca" .OR. cn_cfg == "ORCA" ) THEN 153 IF( nn_cfg == 2 .AND. ln_tsd_dmp ) THEN ! some hand made alterations 154 ! 155 ij0 = 101 + nn_hls ; ij1 = 109 + nn_hls ! Reduced T & S in the Alboran Sea 156 ii0 = 141 + nn_hls - 1 ; ii1 = 155 + nn_hls - 1 157 DO jj = mj0(ij0), mj1(ij1) 158 DO ji = mi0(ii0), mi1(ii1) 159 sf_tsd(jp_tem)%fnow(ji,jj,13:13) = sf_tsd(jp_tem)%fnow(ji,jj,13:13) - 0.20_wp 160 sf_tsd(jp_tem)%fnow(ji,jj,14:15) = sf_tsd(jp_tem)%fnow(ji,jj,14:15) - 0.35_wp 161 sf_tsd(jp_tem)%fnow(ji,jj,16:25) = sf_tsd(jp_tem)%fnow(ji,jj,16:25) - 0.40_wp 162 ! 163 sf_tsd(jp_sal)%fnow(ji,jj,13:13) = sf_tsd(jp_sal)%fnow(ji,jj,13:13) - 0.15_wp 164 sf_tsd(jp_sal)%fnow(ji,jj,14:15) = sf_tsd(jp_sal)%fnow(ji,jj,14:15) - 0.25_wp 165 sf_tsd(jp_sal)%fnow(ji,jj,16:17) = sf_tsd(jp_sal)%fnow(ji,jj,16:17) - 0.30_wp 166 sf_tsd(jp_sal)%fnow(ji,jj,18:25) = sf_tsd(jp_sal)%fnow(ji,jj,18:25) - 0.35_wp 155 ! 156 ! !== ORCA_R2 configuration and T & S damping ==! 157 ! TODO: NOT TESTED- requires orca2 158 IF( cn_cfg == "orca" .OR. cn_cfg == "ORCA" ) THEN 159 IF( nn_cfg == 2 .AND. ln_tsd_dmp ) THEN ! some hand made alterations 160 ! 161 ij0 = 101 + nn_hls ; ij1 = 109 + nn_hls ! Reduced T & S in the Alboran Sea 162 ii0 = 141 + nn_hls - 1 ; ii1 = 155 + nn_hls - 1 163 DO jj = mj0(ij0), mj1(ij1) 164 DO ji = mi0(ii0), mi1(ii1) 165 sf_tsd(jp_tem)%fnow(ji,jj,13:13) = sf_tsd(jp_tem)%fnow(ji,jj,13:13) - 0.20_wp 166 sf_tsd(jp_tem)%fnow(ji,jj,14:15) = sf_tsd(jp_tem)%fnow(ji,jj,14:15) - 0.35_wp 167 sf_tsd(jp_tem)%fnow(ji,jj,16:25) = sf_tsd(jp_tem)%fnow(ji,jj,16:25) - 0.40_wp 168 ! 169 sf_tsd(jp_sal)%fnow(ji,jj,13:13) = sf_tsd(jp_sal)%fnow(ji,jj,13:13) - 0.15_wp 170 sf_tsd(jp_sal)%fnow(ji,jj,14:15) = sf_tsd(jp_sal)%fnow(ji,jj,14:15) - 0.25_wp 171 sf_tsd(jp_sal)%fnow(ji,jj,16:17) = sf_tsd(jp_sal)%fnow(ji,jj,16:17) - 0.30_wp 172 sf_tsd(jp_sal)%fnow(ji,jj,18:25) = sf_tsd(jp_sal)%fnow(ji,jj,18:25) - 0.35_wp 173 END DO 167 174 END DO 168 END DO 169 ij0 = 87 + nn_hls ; ij1 = 96 + nn_hls ! Reduced temperature in Red Sea 170 ii0 = 148 + nn_hls - 1 ; ii1 = 160 + nn_hls - 1 171 sf_tsd(jp_tem)%fnow( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 4:10 ) = 7.0_wp 172 sf_tsd(jp_tem)%fnow( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 11:13 ) = 6.5_wp 173 sf_tsd(jp_tem)%fnow( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 14:20 ) = 6.0_wp 174 ENDIF 175 ENDIF 175 ij0 = 87 + nn_hls ; ij1 = 96 + nn_hls ! Reduced temperature in Red Sea 176 ii0 = 148 + nn_hls - 1 ; ii1 = 160 + nn_hls - 1 177 sf_tsd(jp_tem)%fnow( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 4:10 ) = 7.0_wp 178 sf_tsd(jp_tem)%fnow( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 11:13 ) = 6.5_wp 179 sf_tsd(jp_tem)%fnow( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 14:20 ) = 6.0_wp 180 ENDIF 181 ENDIF 176 182 !!gm end 177 ! 178 ptsd(:,:,:,jp_tem) = sf_tsd(jp_tem)%fnow(:,:,:) ! NO mask 179 ptsd(:,:,:,jp_sal) = sf_tsd(jp_sal)%fnow(:,:,:) 183 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 1 ) ! Revert to tile domain 184 ENDIF 185 ! 186 DO_3D_00_00( 1, jpk ) 187 ptsd(ji,jj,jk,jp_tem) = sf_tsd(jp_tem)%fnow(ji,jj,jk) ! NO mask 188 ptsd(ji,jj,jk,jp_sal) = sf_tsd(jp_sal)%fnow(ji,jj,jk) 189 END_3D 180 190 ! 181 191 IF( ln_sco ) THEN !== s- or mixed s-zps-coordinate ==! 182 192 ! 183 IF( kt == nit000 .AND. lwp )THEN 184 WRITE(numout,*) 185 WRITE(numout,*) 'dta_tsd: interpolates T & S data onto the s- or mixed s-z-coordinate mesh' 193 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 194 IF( kt == nit000 .AND. lwp )THEN 195 WRITE(numout,*) 196 WRITE(numout,*) 'dta_tsd: interpolates T & S data onto the s- or mixed s-z-coordinate mesh' 197 ENDIF 186 198 ENDIF 187 199 ! … … 215 227 ELSE !== z- or zps- coordinate ==! 216 228 ! 217 ptsd(:,:,:,jp_tem) = ptsd(:,:,:,jp_tem) * tmask(:,:,:) ! Mask 218 ptsd(:,:,:,jp_sal) = ptsd(:,:,:,jp_sal) * tmask(:,:,:) 219 ! 229 DO_3D_00_00( 1, jpk ) 230 ptsd(ji,jj,jk,jp_tem) = ptsd(ji,jj,jk,jp_tem) * tmask(ji,jj,jk) ! Mask 231 ptsd(ji,jj,jk,jp_sal) = ptsd(ji,jj,jk,jp_sal) * tmask(ji,jj,jk) 232 END_3D 233 ! 234 ! TODO: NOT TESTED- requires zps 220 235 IF( ln_zps ) THEN ! zps-coordinate (partial steps) interpolation at the last ocean level 221 236 DO_2D_11_11 -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/DYN/dynhpg.F90
r12377 r13409 300 300 INTEGER :: iku, ikv ! temporary integers 301 301 REAL(wp) :: zcoef0, zcoef1, zcoef2, zcoef3 ! temporary scalars 302 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zhpi, zhpj 303 REAL(wp), DIMENSION(jpi,jpj) :: zgtsu, zgtsv, zgru, zgrv 302 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zhpi, zhpj 303 REAL(wp), DIMENSION(jpi,jpj,1) :: zgtsu, zgtsv 304 REAL(wp), DIMENSION(jpi,jpj) :: zgru, zgrv 304 305 !!---------------------------------------------------------------------- 305 306 ! -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/IOM/prtctl.F90
r12942 r13409 42 42 43 43 SUBROUTINE prt_ctl (tab2d_1, tab3d_1, mask1, clinfo1, tab2d_2, tab3d_2, & 44 & mask2, clinfo2, kdim, clinfo3 )44 & mask2, clinfo2, kdim, clinfo3, psum1, psum2 ) 45 45 !!---------------------------------------------------------------------- 46 46 !! *** ROUTINE prt_ctl *** … … 75 75 !! clinfo2 : information about the tab[23]d_2 array 76 76 !! kdim : k- direction for 3D arrays 77 !! clinfo3 : additional information 77 !! clinfo3 : additional information 78 !! psum1 : sum of tab[23]d_1 over tiles (required if ln_tile=T) 79 !! psum2 : sum of tab[23]d_2 over tiles (required if ln_tile=T) 78 80 !!---------------------------------------------------------------------- 79 81 REAL(wp), DIMENSION(:,:) , INTENT(in), OPTIONAL :: tab2d_1 … … 87 89 INTEGER , INTENT(in), OPTIONAL :: kdim 88 90 CHARACTER (len=*) , INTENT(in), OPTIONAL :: clinfo3 91 REAL(wp) , INTENT(inout), OPTIONAL :: psum1, psum2 89 92 ! 90 93 CHARACTER (len=15) :: cl2 91 94 INTEGER :: jn, sind, eind, kdir,j_id 95 INTEGER :: iictls, iictle, ijctls, ijctle 92 96 REAL(wp) :: zsum1, zsum2, zvctl1, zvctl2 93 97 REAL(wp), DIMENSION(jpi,jpj) :: ztab2d_1, ztab2d_2 94 98 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zmask1, zmask2, ztab3d_1, ztab3d_2 95 99 !!---------------------------------------------------------------------- 96 ! TODO: TO BE TILED97 IF( ntile /= 0 .AND. ntile /= nijtile ) RETURN98 100 99 101 ! Arrays, scalars initialization … … 120 122 IF( PRESENT(mask1) ) zmask1 (:,:,:) = mask1 (:,:,:) 121 123 IF( PRESENT(mask2) ) zmask2 (:,:,:) = mask2 (:,:,:) 124 125 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 126 IF( PRESENT(psum1) ) psum1 = 0.e0 127 IF( PRESENT(psum2) ) psum2 = 0.e0 128 ENDIF 129 130 ! psum1 and psum2 must be used to sum over tiles 131 IF( ntile /= 0 ) THEN 132 IF( (PRESENT(tab2d_1) .OR. PRESENT(tab3d_1)) .AND. .NOT. PRESENT(psum1) ) & 133 & CALL ctl_stop('prt_ctl must be called with psum1 when ln_tile = T (clinfo1 = "'//TRIM(clinfo1)//'")') 134 IF( (PRESENT(tab2d_2) .OR. PRESENT(tab3d_2)) .AND. .NOT. PRESENT(psum2) ) & 135 & CALL ctl_stop('prt_ctl must be called with psum2 when ln_tile = T (clinfo2 = "'//TRIM(cl2)//'")') 136 ENDIF 122 137 123 138 IF( lk_mpp .AND. jpnij > 1 ) THEN ! processor number … … 158 173 ENDIF 159 174 160 IF( PRESENT(clinfo3)) THEN 161 IF ( clinfo3 == 'tra' ) THEN 162 zvctl1 = t_ctll(jn) 163 zvctl2 = s_ctll(jn) 164 ELSEIF ( clinfo3 == 'dyn' ) THEN 165 zvctl1 = u_ctll(jn) 166 zvctl2 = v_ctll(jn) 175 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only on the last tile 176 IF( PRESENT(clinfo3)) THEN 177 IF ( clinfo3 == 'tra' ) THEN 178 zvctl1 = t_ctll(jn) 179 zvctl2 = s_ctll(jn) 180 ELSEIF ( clinfo3 == 'dyn' ) THEN 181 zvctl1 = u_ctll(jn) 182 zvctl2 = v_ctll(jn) 183 ENDIF 167 184 ENDIF 168 185 ENDIF 169 186 170 ! Compute the sum control 171 ! 2D arrays 172 IF( PRESENT(tab2d_1) ) THEN 173 zsum1 = SUM( ztab2d_1(nictls:nictle,njctls:njctle)*zmask1(nictls:nictle,njctls:njctle,1) ) 174 zsum2 = SUM( ztab2d_2(nictls:nictle,njctls:njctle)*zmask2(nictls:nictle,njctls:njctle,1) ) 187 iictls = MAX( nictls, ntsi%h0 ) 188 iictle = MIN( nictle, ntei%h0 ) 189 ijctls = MAX( njctls, ntsj%h0 ) 190 ijctle = MIN( njctle, ntej%h0 ) 191 192 ! Compute the sum control only where the tile domain and control print area overlap 193 IF( iictle >= iictls .AND. ijctle >= ijctls ) THEN 194 ! 2D arrays 195 IF( PRESENT(tab2d_1) ) THEN 196 zsum1 = SUM( ztab2d_1(iictls:iictle,ijctls:ijctle)*zmask1(iictls:iictle,ijctls:ijctle,1) ) 197 zsum2 = SUM( ztab2d_2(iictls:iictle,ijctls:ijctle)*zmask2(iictls:iictle,ijctls:ijctle,1) ) 198 ENDIF 199 200 ! 3D arrays 201 IF( PRESENT(tab3d_1) ) THEN 202 zsum1 = SUM( ztab3d_1(iictls:iictle,ijctls:ijctle,1:kdir)*zmask1(iictls:iictle,ijctls:ijctle,1:kdir) ) 203 zsum2 = SUM( ztab3d_2(iictls:iictle,ijctls:ijctle,1:kdir)*zmask2(iictls:iictle,ijctls:ijctle,1:kdir) ) 204 ENDIF 205 206 ! Sum over tiles 207 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only on the last tile 208 IF( PRESENT(psum1) ) zsum1 = psum1 + zsum1 209 IF( PRESENT(psum2) ) zsum2 = psum2 + zsum2 210 ELSE 211 IF( PRESENT(psum1) ) psum1 = psum1 + zsum1 212 IF( PRESENT(psum2) ) psum2 = psum2 + zsum2 213 ENDIF 175 214 ENDIF 176 215 177 ! 3D arrays178 IF( PRESENT(tab3d_1) ) THEN179 zsum1 = SUM( ztab3d_1(nictls:nictle,njctls:njctle,1:kdir)*zmask1(nictls:nictle,njctls:njctle,1:kdir) )180 zsum2 = SUM( ztab3d_2(nictls:nictle,njctls:njctle,1:kdir)*zmask2(nictls:nictle,njctls:njctle,1:kdir) )181 ENDIF182 183 216 ! Print the result 184 IF( PRESENT(clinfo3) ) THEN 185 WRITE(j_id,FMT='(a,D23.16,3x,a,D23.16)')clinfo1, zsum1-zvctl1, cl2, zsum2-zvctl2 186 SELECT CASE( clinfo3 ) 187 CASE ( 'tra-ta' ) 188 t_ctll(jn) = zsum1 189 CASE ( 'tra' ) 190 t_ctll(jn) = zsum1 191 s_ctll(jn) = zsum2 192 CASE ( 'dyn' ) 193 u_ctll(jn) = zsum1 194 v_ctll(jn) = zsum2 195 END SELECT 196 ELSEIF ( PRESENT(clinfo2) .OR. PRESENT(tab2d_2) .OR. PRESENT(tab3d_2) ) THEN 197 WRITE(j_id,FMT='(a,D23.16,3x,a,D23.16)')clinfo1, zsum1, cl2, zsum2 198 ELSE 199 WRITE(j_id,FMT='(a,D23.16)')clinfo1, zsum1 217 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only on the last tile 218 IF( PRESENT(clinfo3) ) THEN 219 WRITE(j_id,FMT='(a,D23.16,3x,a,D23.16)')clinfo1, zsum1-zvctl1, cl2, zsum2-zvctl2 220 SELECT CASE( clinfo3 ) 221 CASE ( 'tra-ta' ) 222 t_ctll(jn) = zsum1 223 CASE ( 'tra' ) 224 t_ctll(jn) = zsum1 225 s_ctll(jn) = zsum2 226 CASE ( 'dyn' ) 227 u_ctll(jn) = zsum1 228 v_ctll(jn) = zsum2 229 END SELECT 230 ELSEIF ( PRESENT(clinfo2) .OR. PRESENT(tab2d_2) .OR. PRESENT(tab3d_2) ) THEN 231 WRITE(j_id,FMT='(a,D23.16,3x,a,D23.16)')clinfo1, zsum1, cl2, zsum2 232 ELSE 233 WRITE(j_id,FMT='(a,D23.16)')clinfo1, zsum1 234 ENDIF 200 235 ENDIF 201 236 -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/LDF/ldftra.F90
r12738 r13409 725 725 !! ** Action : pu, pv increased by the eiv transport 726 726 !!---------------------------------------------------------------------- 727 INTEGER 728 INTEGER 729 INTEGER 730 CHARACTER(len=3) 731 REAL(wp), DIMENSION( jpi,jpj,jpk), INTENT(inout) :: pu! in : 3 ocean transport components [m3/s]732 REAL(wp), DIMENSION( jpi,jpj,jpk), INTENT(inout) :: pv! out: 3 ocean transport components [m3/s]733 REAL(wp), DIMENSION( jpi,jpj,jpk), INTENT(inout) :: pw! increased by the eiv [m3/s]727 INTEGER , INTENT(in ) :: kt ! ocean time-step index 728 INTEGER , INTENT(in ) :: kit000 ! first time step index 729 INTEGER , INTENT(in ) :: Kmm, Krhs ! ocean time level indices 730 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 731 REAL(wp), DIMENSION(A2D,jpk), INTENT(inout) :: pu ! in : 3 ocean transport components [m3/s] 732 REAL(wp), DIMENSION(A2D,jpk), INTENT(inout) :: pv ! out: 3 ocean transport components [m3/s] 733 REAL(wp), DIMENSION(A2D,jpk), INTENT(inout) :: pw ! increased by the eiv [m3/s] 734 734 !! 735 735 INTEGER :: ji, jj, jk ! dummy loop indices 736 736 REAL(wp) :: zuwk, zuwk1, zuwi, zuwi1 ! local scalars 737 737 REAL(wp) :: zvwk, zvwk1, zvwj, zvwj1 ! - - 738 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpsi_uw, zpsi_vw 739 !!---------------------------------------------------------------------- 740 ! 741 IF( kt == kit000 ) THEN 742 IF(lwp) WRITE(numout,*) 743 IF(lwp) WRITE(numout,*) 'ldf_eiv_trp : eddy induced advection on ', cdtype,' :' 744 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ add to velocity fields the eiv component' 738 REAL(wp), DIMENSION(A2D,jpk) :: zpsi_uw, zpsi_vw 739 !!---------------------------------------------------------------------- 740 ! 741 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 742 IF( kt == kit000 ) THEN 743 IF(lwp) WRITE(numout,*) 744 IF(lwp) WRITE(numout,*) 'ldf_eiv_trp : eddy induced advection on ', cdtype,' :' 745 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ add to velocity fields the eiv component' 746 ENDIF 745 747 ENDIF 746 748 … … 781 783 !! 782 784 !!---------------------------------------------------------------------- 783 REAL(wp), DIMENSION( jpi,jpj,jpk), INTENT(inout) :: psi_uw, psi_vw ! streamfunction [m3/s]784 INTEGER 785 REAL(wp), DIMENSION(A2D,jpk), INTENT(inout) :: psi_uw, psi_vw ! streamfunction [m3/s] 786 INTEGER , INTENT(in ) :: Kmm ! ocean time level indices 785 787 ! 786 788 INTEGER :: ji, jj, jk ! dummy loop indices 787 789 REAL(wp) :: zztmp ! local scalar 788 REAL(wp), DIMENSION( jpi,jpj) :: zw2d ! 2D workspace789 REAL(wp), DIMENSION( jpi,jpj,jpk) :: zw3d ! 3D workspace790 REAL(wp), DIMENSION(A2D) :: zw2d ! 2D workspace 791 REAL(wp), DIMENSION(A2D,jpk) :: zw3d ! 3D workspace 790 792 !!---------------------------------------------------------------------- 791 793 ! … … 802 804 zw3d(:,:,jpk) = 0._wp ! bottom value always 0 803 805 ! 804 DO jk = 1, jpkm1! e2u e3u u_eiv = -dk[psi_uw]805 zw3d( :,:,jk) = ( psi_uw(:,:,jk+1) - psi_uw(:,:,jk) ) / ( e2u(:,:) * e3u(:,:,jk,Kmm) )806 END DO806 DO_3D_11_11( 1, jpkm1 ) ! e2u e3u u_eiv = -dk[psi_uw] 807 zw3d(ji,jj,jk) = ( psi_uw(ji,jj,jk+1) - psi_uw(ji,jj,jk) ) / ( e2u(ji,jj) * e3u(ji,jj,jk,Kmm) ) 808 END_3D 807 809 CALL iom_put( "uoce_eiv", zw3d ) 808 810 ! 809 DO jk = 1, jpkm1! e1v e3v v_eiv = -dk[psi_vw]810 zw3d( :,:,jk) = ( psi_vw(:,:,jk+1) - psi_vw(:,:,jk) ) / ( e1v(:,:) * e3v(:,:,jk,Kmm) )811 END DO811 DO_3D_11_11( 1, jpkm1 ) ! e1v e3v v_eiv = -dk[psi_vw] 812 zw3d(ji,jj,jk) = ( psi_vw(ji,jj,jk+1) - psi_vw(ji,jj,jk) ) / ( e1v(ji,jj) * e3v(ji,jj,jk,Kmm) ) 813 END_3D 812 814 CALL iom_put( "voce_eiv", zw3d ) 813 815 ! … … 820 822 ! 821 823 IF( iom_use('weiv_masstr') ) THEN ! vertical mass transport & its square value 822 zw2d(:,:) = rho0 * e1e2t(:,:) 824 DO_2D_11_11 825 zw2d(ji,jj) = rho0 * e1e2t(ji,jj) 826 END_2D 823 827 DO jk = 1, jpk 824 828 zw3d(:,:,jk) = zw3d(:,:,jk) * zw2d(:,:) … … 866 870 END_3D 867 871 CALL lbc_lnk( 'ldftra', zw2d, 'V', -1. ) 868 CALL iom_put( "veiv_heattr", zztmp * zw2d ) ! heat transport in j-direction 869 CALL iom_put( "veiv_heattr", zztmp * zw3d ) ! heat transport in j-direction 872 CALL lbc_lnk( 'ldftra', zw3d, 'V', -1. ) 873 CALL iom_put( "veiv_heattr" , zztmp * zw2d ) ! heat transport in j-direction 874 CALL iom_put( "veiv_heattr3d", zztmp * zw3d ) ! heat transport in j-direction 870 875 ! 871 876 IF( iom_use( 'sophteiv' ) ) CALL dia_ptr_hst( jp_tem, 'eiv', 0.5 * zw3d ) … … 893 898 END_3D 894 899 CALL lbc_lnk( 'ldftra', zw2d, 'V', -1. ) 895 CALL iom_put( "veiv_salttr", zztmp * zw2d ) ! salt transport in j-direction 896 CALL iom_put( "veiv_salttr", zztmp * zw3d ) ! salt transport in j-direction 900 CALL lbc_lnk( 'ldftra', zw3d, 'V', -1. ) 901 CALL iom_put( "veiv_salttr" , zztmp * zw2d ) ! salt transport in j-direction 902 CALL iom_put( "veiv_salttr3d", zztmp * zw3d ) ! salt transport in j-direction 897 903 ! 898 904 IF( iom_use( 'sopsteiv' ) ) CALL dia_ptr_hst( jp_sal, 'eiv', 0.5 * zw3d ) -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/eosbn2.F90
r12489 r13409 39 39 !!---------------------------------------------------------------------- 40 40 USE dom_oce ! ocean space and time domain 41 USE domutl, ONLY : is_tile 41 42 USE phycst ! physical constants 42 43 USE stopar ! Stochastic T/S fluctuations … … 186 187 !!---------------------------------------------------------------------- 187 188 CONTAINS 188 189 189 SUBROUTINE eos_insitu( pts, prd, pdep ) 190 !! 191 REAL(wp), DIMENSION(:,:,:,:), INTENT(in ) :: pts ! 1 : potential temperature [Celsius] 192 ! ! 2 : salinity [psu] 193 REAL(wp), DIMENSION(:,:,:) , INTENT( out) :: prd ! in situ density [-] 194 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pdep ! depth [m] 195 !! 196 CALL eos_insitu_t( pts, is_tile(pts), prd, is_tile(prd), pdep, is_tile(pdep) ) 197 END SUBROUTINE eos_insitu 198 199 SUBROUTINE eos_insitu_t( pts, ktts, prd, ktrd, pdep, ktdep ) 190 200 !!---------------------------------------------------------------------- 191 201 !! *** ROUTINE eos_insitu *** … … 221 231 !! TEOS-10 Manual, 2010 222 232 !!---------------------------------------------------------------------- 223 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT(in ) :: pts ! 1 : potential temperature [Celsius] 224 ! ! 2 : salinity [psu] 225 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT( out) :: prd ! in situ density [-] 226 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pdep ! depth [m] 233 INTEGER , INTENT(in ) :: ktts, ktrd, ktdep 234 REAL(wp), DIMENSION(T2D(ktts) ,jpk,jpts), INTENT(in ) :: pts ! 1 : potential temperature [Celsius] 235 ! ! 2 : salinity [psu] 236 REAL(wp), DIMENSION(T2D(ktrd) ,jpk ), INTENT( out) :: prd ! in situ density [-] 237 REAL(wp), DIMENSION(T2D(ktdep),jpk ), INTENT(in ) :: pdep ! depth [m] 227 238 ! 228 239 INTEGER :: ji, jj, jk ! dummy loop indices 229 240 REAL(wp) :: zt , zh , zs , ztm ! local scalars 230 241 REAL(wp) :: zn , zn0, zn1, zn2, zn3 ! - - 242 REAL(wp), SAVE :: zsum1 231 243 !!---------------------------------------------------------------------- 232 244 ! … … 288 300 END SELECT 289 301 ! 290 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=prd, clinfo1=' eos-insitu : ', kdim=jpk )302 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=prd, clinfo1=' eos-insitu : ', kdim=jpk, psum1=zsum1 ) 291 303 ! 292 304 IF( ln_timing ) CALL timing_stop('eos-insitu') 293 305 ! 294 END SUBROUTINE eos_insitu 306 END SUBROUTINE eos_insitu_t 295 307 296 308 297 309 SUBROUTINE eos_insitu_pot( pts, prd, prhop, pdep ) 310 !! 311 REAL(wp), DIMENSION(:,:,:,:), INTENT(in ) :: pts ! 1 : potential temperature [Celsius] 312 ! ! 2 : salinity [psu] 313 REAL(wp), DIMENSION(:,:,:) , INTENT( out) :: prd ! in situ density [-] 314 REAL(wp), DIMENSION(:,:,:) , INTENT( out) :: prhop ! potential density (surface referenced) 315 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pdep ! depth [m] 316 !! 317 CALL eos_insitu_pot_t( pts, is_tile(pts), prd, is_tile(prd), prhop, is_tile(prhop), pdep, is_tile(pdep) ) 318 END SUBROUTINE eos_insitu_pot 319 320 321 SUBROUTINE eos_insitu_pot_t( pts, ktts, prd, ktrd, prhop, ktrhop, pdep, ktdep ) 298 322 !!---------------------------------------------------------------------- 299 323 !! *** ROUTINE eos_insitu_pot *** … … 308 332 !! 309 333 !!---------------------------------------------------------------------- 310 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT(in ) :: pts ! 1 : potential temperature [Celsius] 311 ! ! 2 : salinity [psu] 312 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT( out) :: prd ! in situ density [-] 313 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT( out) :: prhop ! potential density (surface referenced) 314 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pdep ! depth [m] 334 INTEGER , INTENT(in ) :: ktts, ktrd, ktrhop, ktdep 335 REAL(wp), DIMENSION(T2D(ktts) ,jpk,jpts), INTENT(in ) :: pts ! 1 : potential temperature [Celsius] 336 ! ! 2 : salinity [psu] 337 REAL(wp), DIMENSION(T2D(ktrd) ,jpk ), INTENT( out) :: prd ! in situ density [-] 338 REAL(wp), DIMENSION(T2D(ktrhop),jpk ), INTENT( out) :: prhop ! potential density (surface referenced) 339 REAL(wp), DIMENSION(T2D(ktdep) ,jpk ), INTENT(in ) :: pdep ! depth [m] 315 340 ! 316 341 INTEGER :: ji, jj, jk, jsmp ! dummy loop indices … … 318 343 REAL(wp) :: zt , zh , zstemp, zs , ztm ! local scalars 319 344 REAL(wp) :: zn , zn0, zn1, zn2, zn3 ! - - 345 REAL(wp), SAVE :: zsum1, zsum2 320 346 REAL(wp), DIMENSION(:), ALLOCATABLE :: zn0_sto, zn_sto, zsign ! local vectors 321 347 !!---------------------------------------------------------------------- … … 443 469 END SELECT 444 470 ! 445 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=prd, clinfo1=' eos-pot: ', tab3d_2=prhop, clinfo2=' pot : ', kdim=jpk ) 471 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=prd, clinfo1=' eos-pot: ', psum1=zsum1, & 472 & tab3d_2=prhop, clinfo2=' pot : ', psum2=zsum2, kdim=jpk ) 446 473 ! 447 474 IF( ln_timing ) CALL timing_stop('eos-pot') 448 475 ! 449 END SUBROUTINE eos_insitu_pot 476 END SUBROUTINE eos_insitu_pot_t 450 477 451 478 452 479 SUBROUTINE eos_insitu_2d( pts, pdep, prd ) 480 !! 481 REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: pts ! 1 : potential temperature [Celsius] 482 ! ! 2 : salinity [psu] 483 REAL(wp), DIMENSION(:,:) , INTENT(in ) :: pdep ! depth [m] 484 REAL(wp), DIMENSION(:,:) , INTENT( out) :: prd ! in situ density 485 !! 486 CALL eos_insitu_2d_t( pts, is_tile(pts), pdep, is_tile(pdep), prd, is_tile(prd) ) 487 END SUBROUTINE eos_insitu_2d 488 489 490 SUBROUTINE eos_insitu_2d_t( pts, ktts, pdep, ktdep, prd, ktrd ) 453 491 !!---------------------------------------------------------------------- 454 492 !! *** ROUTINE eos_insitu_2d *** … … 461 499 !! 462 500 !!---------------------------------------------------------------------- 463 REAL(wp), DIMENSION(jpi,jpj,jpts), INTENT(in ) :: pts ! 1 : potential temperature [Celsius] 464 ! ! 2 : salinity [psu] 465 REAL(wp), DIMENSION(jpi,jpj) , INTENT(in ) :: pdep ! depth [m] 466 REAL(wp), DIMENSION(jpi,jpj) , INTENT( out) :: prd ! in situ density 501 INTEGER , INTENT(in ) :: ktts, ktdep, ktrd 502 REAL(wp), DIMENSION(T2D(ktts),jpts), INTENT(in ) :: pts ! 1 : potential temperature [Celsius] 503 ! ! 2 : salinity [psu] 504 REAL(wp), DIMENSION(T2D(ktdep) ), INTENT(in ) :: pdep ! depth [m] 505 REAL(wp), DIMENSION(T2D(ktrd) ), INTENT( out) :: prd ! in situ density 467 506 ! 468 507 INTEGER :: ji, jj, jk ! dummy loop indices 469 508 REAL(wp) :: zt , zh , zs ! local scalars 470 509 REAL(wp) :: zn , zn0, zn1, zn2, zn3 ! - - 510 REAL(wp), SAVE :: zsum1 471 511 !!---------------------------------------------------------------------- 472 512 ! … … 530 570 END SELECT 531 571 ! 532 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab2d_1=prd, clinfo1=' eos2d: ' )572 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab2d_1=prd, clinfo1=' eos2d: ', psum1=zsum1 ) 533 573 ! 534 574 IF( ln_timing ) CALL timing_stop('eos2d') 535 575 ! 536 END SUBROUTINE eos_insitu_2d 576 END SUBROUTINE eos_insitu_2d_t 537 577 538 578 539 579 SUBROUTINE rab_3d( pts, pab, Kmm ) 580 !! 581 INTEGER , INTENT(in ) :: Kmm ! time level index 582 REAL(wp), DIMENSION(:,:,:,:), INTENT(in ) :: pts ! pot. temperature & salinity 583 REAL(wp), DIMENSION(:,:,:,:), INTENT( out) :: pab ! thermal/haline expansion ratio 584 !! 585 CALL rab_3d_t( pts, is_tile(pts), pab, is_tile(pab), Kmm ) 586 END SUBROUTINE rab_3d 587 588 589 SUBROUTINE rab_3d_t( pts, ktts, pab, ktab, Kmm ) 540 590 !!---------------------------------------------------------------------- 541 591 !! *** ROUTINE rab_3d *** … … 547 597 !! ** Action : - pab : thermal/haline expansion ratio at T-points 548 598 !!---------------------------------------------------------------------- 549 INTEGER , INTENT(in ) :: Kmm ! time level index 550 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT(in ) :: pts ! pot. temperature & salinity 551 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT( out) :: pab ! thermal/haline expansion ratio 599 INTEGER , INTENT(in ) :: Kmm ! time level index 600 INTEGER , INTENT(in ) :: ktts, ktab 601 REAL(wp), DIMENSION(T2D(ktts),jpk,jpts), INTENT(in ) :: pts ! pot. temperature & salinity 602 REAL(wp), DIMENSION(T2D(ktab),jpk,jpts), INTENT( out) :: pab ! thermal/haline expansion ratio 552 603 ! 553 604 INTEGER :: ji, jj, jk ! dummy loop indices 554 605 REAL(wp) :: zt , zh , zs , ztm ! local scalars 555 606 REAL(wp) :: zn , zn0, zn1, zn2, zn3 ! - - 607 REAL(wp), SAVE :: zsum1, zsum2 556 608 !!---------------------------------------------------------------------- 557 609 ! … … 635 687 END SELECT 636 688 ! 637 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=pab(:,:,:,jp_tem), clinfo1=' rab_3d_t: ', &638 & tab3d_2=pab(:,:,:,jp_sal), clinfo2=' rab_3d_s : ', kdim=jpk )689 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=pab(:,:,:,jp_tem), clinfo1=' rab_3d_t: ', psum1=zsum1, & 690 & tab3d_2=pab(:,:,:,jp_sal), clinfo2=' rab_3d_s : ', psum2=zsum2, kdim=jpk ) 639 691 ! 640 692 IF( ln_timing ) CALL timing_stop('rab_3d') 641 693 ! 642 END SUBROUTINE rab_3d 694 END SUBROUTINE rab_3d_t 643 695 644 696 645 697 SUBROUTINE rab_2d( pts, pdep, pab, Kmm ) 698 !! 699 INTEGER , INTENT(in ) :: Kmm ! time level index 700 REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: pts ! pot. temperature & salinity 701 REAL(wp), DIMENSION(:,:) , INTENT(in ) :: pdep ! depth [m] 702 REAL(wp), DIMENSION(:,:,:), INTENT( out) :: pab ! thermal/haline expansion ratio 703 !! 704 CALL rab_2d_t(pts, is_tile(pts), pdep, is_tile(pdep), pab, is_tile(pab), Kmm) 705 END SUBROUTINE rab_2d 706 707 708 SUBROUTINE rab_2d_t( pts, ktts, pdep, ktdep, pab, ktab, Kmm ) 646 709 !!---------------------------------------------------------------------- 647 710 !! *** ROUTINE rab_2d *** … … 651 714 !! ** Action : - pab : thermal/haline expansion ratio at T-points 652 715 !!---------------------------------------------------------------------- 653 INTEGER , INTENT(in ) :: Kmm ! time level index 654 REAL(wp), DIMENSION(jpi,jpj,jpts) , INTENT(in ) :: pts ! pot. temperature & salinity 655 REAL(wp), DIMENSION(jpi,jpj) , INTENT(in ) :: pdep ! depth [m] 656 REAL(wp), DIMENSION(jpi,jpj,jpts) , INTENT( out) :: pab ! thermal/haline expansion ratio 716 INTEGER , INTENT(in ) :: Kmm ! time level index 717 INTEGER , INTENT(in ) :: ktts, ktdep, ktab 718 REAL(wp), DIMENSION(T2D(ktts),jpts), INTENT(in ) :: pts ! pot. temperature & salinity 719 REAL(wp), DIMENSION(T2D(ktdep) ), INTENT(in ) :: pdep ! depth [m] 720 REAL(wp), DIMENSION(T2D(ktab),jpts), INTENT( out) :: pab ! thermal/haline expansion ratio 657 721 ! 658 722 INTEGER :: ji, jj, jk ! dummy loop indices 659 723 REAL(wp) :: zt , zh , zs ! local scalars 660 724 REAL(wp) :: zn , zn0, zn1, zn2, zn3 ! - - 725 REAL(wp), SAVE :: zsum1, zsum2 661 726 !!---------------------------------------------------------------------- 662 727 ! … … 742 807 END SELECT 743 808 ! 744 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab2d_1=pab(:,:,jp_tem), clinfo1=' rab_2d_t: ', &745 & tab2d_2=pab(:,:,jp_sal), clinfo2=' rab_2d_s : ' )809 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab2d_1=pab(:,:,jp_tem), clinfo1=' rab_2d_t: ', psum1=zsum1, & 810 & tab2d_2=pab(:,:,jp_sal), clinfo2=' rab_2d_s : ', psum2=zsum2 ) 746 811 ! 747 812 IF( ln_timing ) CALL timing_stop('rab_2d') 748 813 ! 749 END SUBROUTINE rab_2d 814 END SUBROUTINE rab_2d_t 750 815 751 816 … … 848 913 849 914 SUBROUTINE bn2( pts, pab, pn2, Kmm ) 915 !! 916 INTEGER , INTENT(in ) :: Kmm ! time level index 917 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT(in ) :: pts ! pot. temperature and salinity [Celsius,psu] 918 REAL(wp), DIMENSION(:,:,:,:) , INTENT(in ) :: pab ! thermal/haline expansion coef. [Celsius-1,psu-1] 919 REAL(wp), DIMENSION(:,:,:) , INTENT( out) :: pn2 ! Brunt-Vaisala frequency squared [1/s^2] 920 !! 921 CALL bn2_t( pts, pab, is_tile(pab), pn2, is_tile(pn2), Kmm ) 922 END SUBROUTINE bn2 923 924 925 SUBROUTINE bn2_t( pts, pab, ktab, pn2, ktn2, Kmm ) 850 926 !!---------------------------------------------------------------------- 851 927 !! *** ROUTINE bn2 *** … … 861 937 !! 862 938 !!---------------------------------------------------------------------- 863 INTEGER , INTENT(in ) :: Kmm ! time level index 864 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT(in ) :: pts ! pot. temperature and salinity [Celsius,psu] 865 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT(in ) :: pab ! thermal/haline expansion coef. [Celsius-1,psu-1] 866 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT( out) :: pn2 ! Brunt-Vaisala frequency squared [1/s^2] 939 INTEGER , INTENT(in ) :: Kmm ! time level index 940 INTEGER , INTENT(in ) :: ktab, ktn2 941 REAL(wp), DIMENSION(jpi,jpj, jpk,jpts), INTENT(in ) :: pts ! pot. temperature and salinity [Celsius,psu] 942 REAL(wp), DIMENSION(T2D(ktab),jpk,jpts), INTENT(in ) :: pab ! thermal/haline expansion coef. [Celsius-1,psu-1] 943 REAL(wp), DIMENSION(T2D(ktn2),jpk ), INTENT( out) :: pn2 ! Brunt-Vaisala frequency squared [1/s^2] 867 944 ! 868 945 INTEGER :: ji, jj, jk ! dummy loop indices 869 946 REAL(wp) :: zaw, zbw, zrw ! local scalars 947 REAL(wp), SAVE :: zsum1 870 948 !!---------------------------------------------------------------------- 871 949 ! … … 884 962 END_3D 885 963 ! 886 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=pn2, clinfo1=' bn2 : ', kdim=jpk )964 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=pn2, clinfo1=' bn2 : ', kdim=jpk, psum1=zsum1 ) 887 965 ! 888 966 IF( ln_timing ) CALL timing_stop('bn2') 889 967 ! 890 END SUBROUTINE bn2 968 END SUBROUTINE bn2_t 891 969 892 970 … … 948 1026 949 1027 950 SUBROUTINE eos_fzp_2d( psal, ptf, pdep ) 1028 SUBROUTINE eos_fzp_2d( psal, ptf, pdep ) 1029 !! 1030 REAL(wp), DIMENSION(jpi,jpj), INTENT(in ) :: psal ! salinity [psu] 1031 REAL(wp), DIMENSION(jpi,jpj), INTENT(in ), OPTIONAL :: pdep ! depth [m] 1032 REAL(wp), DIMENSION(:,:) , INTENT(out ) :: ptf ! freezing temperature [Celsius] 1033 !! 1034 CALL eos_fzp_2d_t( psal, ptf, is_tile(ptf), pdep ) 1035 END SUBROUTINE eos_fzp_2d 1036 1037 1038 SUBROUTINE eos_fzp_2d_t( psal, ptf, kttf, pdep ) 951 1039 !!---------------------------------------------------------------------- 952 1040 !! *** ROUTINE eos_fzp *** … … 960 1048 !! Reference : UNESCO tech. papers in the marine science no. 28. 1978 961 1049 !!---------------------------------------------------------------------- 962 REAL(wp), DIMENSION(jpi,jpj), INTENT(in ) :: psal ! salinity [psu] 963 REAL(wp), DIMENSION(jpi,jpj), INTENT(in ), OPTIONAL :: pdep ! depth [m] 964 REAL(wp), DIMENSION(jpi,jpj), INTENT(out ) :: ptf ! freezing temperature [Celsius] 1050 INTEGER , INTENT(in ) :: kttf 1051 REAL(wp), DIMENSION(jpi,jpj) , INTENT(in ) :: psal ! salinity [psu] 1052 REAL(wp), DIMENSION(jpi,jpj) , INTENT(in ), OPTIONAL :: pdep ! depth [m] 1053 REAL(wp), DIMENSION(T2D(kttf)), INTENT(out ) :: ptf ! freezing temperature [Celsius] 965 1054 ! 966 1055 INTEGER :: ji, jj ! dummy loop indices … … 995 1084 END SELECT 996 1085 ! 997 END SUBROUTINE eos_fzp_2d 1086 END SUBROUTINE eos_fzp_2d_t 998 1087 999 1088 -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/traadv.F90
r12810 r13409 18 18 USE oce ! ocean dynamics and active tracers 19 19 USE dom_oce ! ocean space and time domain 20 ! TEMP: This change not necessary after trd_tra is tiled and extended haloes development 21 USE domain, ONLY : dom_tile 20 22 USE domvvl ! variable vertical scale factors 21 23 USE sbcwave ! wave module … … 65 67 INTEGER, PARAMETER :: np_UBS = 4 ! 3rd order Upstream Biased Scheme 66 68 INTEGER, PARAMETER :: np_QCK = 5 ! QUICK scheme 67 69 70 !! * Substitutions 71 # include "do_loop_substitute.h90" 68 72 !!---------------------------------------------------------------------- 69 73 !! NEMO/OCE 4.0 , NEMO Consortium (2018) … … 85 89 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts,jpt), INTENT(inout) :: pts ! active tracers and RHS of tracer equation 86 90 ! 87 INTEGER :: jk ! dummy loop index 88 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zuu, zvv, zww ! 3D workspace 89 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ztrdt, ztrds 91 ! TEMP: This change not necessary after trd_tra is tiled 92 INTEGER :: itile 93 INTEGER :: ji, jj, jk ! dummy loop index 94 REAL(wp), SAVE :: zsum1, zsum2 95 ! TEMP: This change not necessary and can be A2D if using XIOS (subdomain support) 96 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, SAVE :: zuu, zvv, zww ! 3D workspace 97 ! TEMP: This change not necessary after trd_tra is tiled 98 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, SAVE :: ztrdt, ztrds 99 ! TEMP: This change not necessary after extra haloes development 100 LOGICAL :: lskip 90 101 !!---------------------------------------------------------------------- 91 102 ! 92 103 IF( ln_timing ) CALL timing_start('tra_adv') 93 104 ! 94 ! !== effective transport ==! 95 zuu(:,:,jpk) = 0._wp 96 zvv(:,:,jpk) = 0._wp 97 zww(:,:,jpk) = 0._wp 98 IF( ln_wave .AND. ln_sdw ) THEN 99 DO jk = 1, jpkm1 ! eulerian transport + Stokes Drift 100 zuu(:,:,jk) = e2u (:,:) * e3u(:,:,jk,Kmm) * ( uu(:,:,jk,Kmm) + usd(:,:,jk) ) 101 zvv(:,:,jk) = e1v (:,:) * e3v(:,:,jk,Kmm) * ( vv(:,:,jk,Kmm) + vsd(:,:,jk) ) 102 zww(:,:,jk) = e1e2t(:,:) * ( ww(:,:,jk) + wsd(:,:,jk) ) 103 END DO 104 ELSE 105 DO jk = 1, jpkm1 106 zuu(:,:,jk) = e2u (:,:) * e3u(:,:,jk,Kmm) * uu(:,:,jk,Kmm) ! eulerian transport only 107 zvv(:,:,jk) = e1v (:,:) * e3v(:,:,jk,Kmm) * vv(:,:,jk,Kmm) 108 zww(:,:,jk) = e1e2t(:,:) * ww(:,:,jk) 109 END DO 110 ENDIF 111 ! 112 IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN ! add z-tilde and/or vvl corrections 113 zuu(:,:,:) = zuu(:,:,:) + un_td(:,:,:) 114 zvv(:,:,:) = zvv(:,:,:) + vn_td(:,:,:) 115 ENDIF 116 ! 117 zuu(:,:,jpk) = 0._wp ! no transport trough the bottom 118 zvv(:,:,jpk) = 0._wp 119 zww(:,:,jpk) = 0._wp 120 ! 121 IF( ln_ldfeiv .AND. .NOT. ln_traldf_triad ) & 122 & CALL ldf_eiv_trp( kt, nit000, zuu, zvv, zww, 'TRA', Kmm, Krhs ) ! add the eiv transport (if necessary) 123 ! 124 IF( ln_mle ) CALL tra_mle_trp( kt, nit000, zuu, zvv, zww, 'TRA', Kmm ) ! add the mle transport (if necessary) 125 ! 126 CALL iom_put( "uocetr_eff", zuu ) ! output effective transport 127 CALL iom_put( "vocetr_eff", zvv ) 128 CALL iom_put( "wocetr_eff", zww ) 129 ! 130 !!gm ??? 131 CALL dia_ptr( kt, Kmm, zvv ) ! diagnose the effective MSF 132 !!gm ??? 133 ! 134 135 IF( l_trdtra ) THEN !* Save ta and sa trends 136 ALLOCATE( ztrdt(jpi,jpj,jpk), ztrds(jpi,jpj,jpk) ) 137 ztrdt(:,:,:) = pts(:,:,:,jp_tem,Krhs) 138 ztrds(:,:,:) = pts(:,:,:,jp_sal,Krhs) 139 ENDIF 140 ! 141 SELECT CASE ( nadv ) !== compute advection trend and add it to general trend ==! 142 ! 143 CASE ( np_CEN ) ! Centered scheme : 2nd / 4th order 144 CALL tra_adv_cen ( kt, nit000, 'TRA', zuu, zvv, zww, Kmm, pts, jpts, Krhs, nn_cen_h, nn_cen_v ) 145 CASE ( np_FCT ) ! FCT scheme : 2nd / 4th order 146 CALL tra_adv_fct ( kt, nit000, 'TRA', rDt, zuu, zvv, zww, Kbb, Kmm, pts, jpts, Krhs, nn_fct_h, nn_fct_v ) 147 CASE ( np_MUS ) ! MUSCL 148 CALL tra_adv_mus ( kt, nit000, 'TRA', rDt, zuu, zvv, zww, Kbb, Kmm, pts, jpts, Krhs, ln_mus_ups ) 149 CASE ( np_UBS ) ! UBS 150 CALL tra_adv_ubs ( kt, nit000, 'TRA', rDt, zuu, zvv, zww, Kbb, Kmm, pts, jpts, Krhs, nn_ubs_v ) 151 CASE ( np_QCK ) ! QUICKEST 152 CALL tra_adv_qck ( kt, nit000, 'TRA', rDt, zuu, zvv, zww, Kbb, Kmm, pts, jpts, Krhs ) 153 ! 154 END SELECT 155 ! 156 IF( l_trdtra ) THEN ! save the advective trends for further diagnostics 157 DO jk = 1, jpkm1 158 ztrdt(:,:,jk) = pts(:,:,jk,jp_tem,Krhs) - ztrdt(:,:,jk) 159 ztrds(:,:,jk) = pts(:,:,jk,jp_sal,Krhs) - ztrds(:,:,jk) 160 END DO 161 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_totad, ztrdt ) 162 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_totad, ztrds ) 163 DEALLOCATE( ztrdt, ztrds ) 105 lskip = .FALSE. 106 107 ! TEMP: These changes not necessary if using XIOS (subdomain support) 108 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 109 ALLOCATE( zuu(jpi,jpj,jpk), zvv(jpi,jpj,jpk), zww(jpi,jpj,jpk) ) 110 ! TEMP: This can be A2D after trd_tra is tiled 111 IF( l_trdtra ) ALLOCATE( ztrdt(jpi,jpj,jpk), ztrds(jpi,jpj,jpk) ) 112 ENDIF 113 114 ! TEMP: These changes not necessary after extra haloes development (lbc_lnk removed from tra_adv_*, ldf_eiv_dia) 115 IF( nadv /= np_CEN .OR. (nadv == np_CEN .AND. nn_cen_h == 4) .OR. ln_ldfeiv_dia ) THEN 116 IF( ln_tile ) THEN 117 IF( ntile == 1 ) THEN 118 CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) 119 ELSE 120 lskip = .TRUE. 121 ENDIF 122 ENDIF 123 ENDIF 124 IF( .NOT. lskip ) THEN 125 126 ! TEMP: This change not necessary after trd_tra is tiled 127 itile = ntile 128 ! !== effective transport ==! 129 ! TODO: NOT TESTED- requires waves 130 IF( ln_wave .AND. ln_sdw ) THEN 131 DO_3D_11_11( 1, jpkm1 ) 132 zuu(ji,jj,jk) = e2u (ji,jj) * e3u(ji,jj,jk,Kmm) * ( uu(ji,jj,jk,Kmm) + usd(ji,jj,jk) ) 133 zvv(ji,jj,jk) = e1v (ji,jj) * e3v(ji,jj,jk,Kmm) * ( vv(ji,jj,jk,Kmm) + vsd(ji,jj,jk) ) 134 zww(ji,jj,jk) = e1e2t(ji,jj) * ( ww(ji,jj,jk) + wsd(ji,jj,jk) ) 135 END_3D 136 ELSE 137 DO_3D_11_11( 1, jpkm1 ) 138 zuu(ji,jj,jk) = e2u (ji,jj) * e3u(ji,jj,jk,Kmm) * uu(ji,jj,jk,Kmm) ! eulerian transport only 139 zvv(ji,jj,jk) = e1v (ji,jj) * e3v(ji,jj,jk,Kmm) * vv(ji,jj,jk,Kmm) 140 zww(ji,jj,jk) = e1e2t(ji,jj) * ww(ji,jj,jk) 141 END_3D 142 ENDIF 143 ! 144 ! TODO: NOT TESTED- requires ztilde 145 IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN ! add z-tilde and/or vvl corrections 146 DO_3D_11_11( 1, jpkm1 ) 147 zuu(ji,jj,jk) = zuu(ji,jj,jk) + un_td(ji,jj,jk) 148 zvv(ji,jj,jk) = zvv(ji,jj,jk) + vn_td(ji,jj,jk) 149 END_3D 150 ENDIF 151 ! 152 DO_2D_11_11 153 zuu(ji,jj,jpk) = 0._wp ! no transport trough the bottom 154 zvv(ji,jj,jpk) = 0._wp 155 zww(ji,jj,jpk) = 0._wp 156 END_2D 157 ! 158 ! TEMP: These changes not necessary if using XIOS (subdomain support) 159 IF( ln_ldfeiv .AND. .NOT. ln_traldf_triad ) & 160 & CALL ldf_eiv_trp( kt, nit000, zuu(A2D,:), zvv(A2D,:), zww(A2D,:), & 161 & 'TRA', Kmm, Krhs ) ! add the eiv transport (if necessary) 162 ! 163 IF( ln_mle ) CALL tra_mle_trp( kt, nit000, zuu(A2D,:), zvv(A2D,:), zww(A2D,:), & 164 & 'TRA', Kmm ) ! add the mle transport (if necessary) 165 ! 166 ! TEMP: This change not necessary if using XIOS (subdomain support) 167 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only on the last tile 168 CALL iom_put( "uocetr_eff", zuu ) ! output effective transport 169 CALL iom_put( "vocetr_eff", zvv ) 170 CALL iom_put( "wocetr_eff", zww ) 171 ENDIF 172 ! 173 !!gm ??? 174 ! TEMP: This change not necessary if using XIOS (subdomain support) 175 CALL dia_ptr( kt, Kmm, zvv(A2D,:) ) ! diagnose the effective MSF 176 !!gm ??? 177 ! 178 179 IF( l_trdtra ) THEN !* Save ta and sa trends 180 DO_3D_00_00( 1, jpkm1 ) 181 ztrdt(ji,jj,jk) = pts(ji,jj,jk,jp_tem,Krhs) 182 ztrds(ji,jj,jk) = pts(ji,jj,jk,jp_sal,Krhs) 183 END_3D 184 ENDIF 185 ! 186 SELECT CASE ( nadv ) !== compute advection trend and add it to general trend ==! 187 ! 188 CASE ( np_CEN ) ! Centered scheme : 2nd / 4th order 189 CALL tra_adv_cen ( kt, nit000, 'TRA', zuu, zvv, zww, Kmm, pts, jpts, Krhs, nn_cen_h, nn_cen_v ) 190 CASE ( np_FCT ) ! FCT scheme : 2nd / 4th order 191 CALL tra_adv_fct ( kt, nit000, 'TRA', rDt, zuu, zvv, zww, Kbb, Kmm, pts, jpts, Krhs, nn_fct_h, nn_fct_v ) 192 CASE ( np_MUS ) ! MUSCL 193 CALL tra_adv_mus ( kt, nit000, 'TRA', rDt, zuu, zvv, zww, Kbb, Kmm, pts, jpts, Krhs, ln_mus_ups ) 194 CASE ( np_UBS ) ! UBS 195 CALL tra_adv_ubs ( kt, nit000, 'TRA', rDt, zuu, zvv, zww, Kbb, Kmm, pts, jpts, Krhs, nn_ubs_v ) 196 CASE ( np_QCK ) ! QUICKEST 197 CALL tra_adv_qck ( kt, nit000, 'TRA', rDt, zuu, zvv, zww, Kbb, Kmm, pts, jpts, Krhs ) 198 ! 199 END SELECT 200 ! 201 ! TEMP: These changes not necessary after trd_tra is tiled 202 IF( l_trdtra ) THEN ! save the advective trends for further diagnostics 203 DO_3D_00_00( 1, jpkm1 ) 204 ztrdt(ji,jj,jk) = pts(ji,jj,jk,jp_tem,Krhs) - ztrdt(ji,jj,jk) 205 ztrds(ji,jj,jk) = pts(ji,jj,jk,jp_sal,Krhs) - ztrds(ji,jj,jk) 206 END_3D 207 208 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 209 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 210 211 ! TODO: TO BE TILED- trd_tra 212 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_totad, ztrdt ) 213 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_totad, ztrds ) 214 DEALLOCATE( ztrdt, ztrds ) 215 216 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = itile ) ! Revert to tile domain 217 ENDIF 218 ENDIF 219 220 ! TEMP: This change not necessary after extra haloes development (lbc_lnk removed from tra_adv_*, ldf_eiv_dia) 221 IF( ln_tile .AND. ntile == 0 ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 1 ) 222 164 223 ENDIF 165 224 ! ! print mean trends (used for debugging) 166 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Krhs), clinfo1=' adv - Ta: ', mask1=tmask, & 167 & tab3d_2=pts(:,:,:,jp_sal,Krhs), clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra' ) 225 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Krhs), clinfo1=' adv - Ta: ', mask1=tmask, psum1=zsum1, & 226 & tab3d_2=pts(:,:,:,jp_sal,Krhs), clinfo2= ' Sa: ', mask2=tmask, psum2=zsum2, & 227 & clinfo3='tra' ) 228 229 ! TEMP: This change not necessary if using XIOS (subdomain support) 230 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 231 DEALLOCATE( zuu, zvv, zww ) 232 ENDIF 168 233 ! 169 234 IF( ln_timing ) CALL timing_stop( 'tra_adv' ) -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/traadv_cen.F90
r12377 r13409 12 12 !!---------------------------------------------------------------------- 13 13 USE dom_oce ! ocean space and time domain 14 ! TEMP: This change not necessary after trd_tra is tiled 15 USE domain, ONLY : dom_tile 14 16 USE eosbn2 ! equation of state 15 17 USE traadv_fct ! acces to routine interp_4th_cpt … … 70 72 INTEGER , INTENT(in ) :: kn_cen_h ! =2/4 (2nd or 4th order scheme) 71 73 INTEGER , INTENT(in ) :: kn_cen_v ! =2/4 (2nd or 4th order scheme) 74 ! TEMP: This can be A2D after trd_tra is tiled 72 75 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pU, pV, pW ! 3 ocean volume flux components 73 76 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) :: pt ! tracers and RHS of tracer equation … … 75 78 INTEGER :: ji, jj, jk, jn ! dummy loop indices 76 79 INTEGER :: ierr ! local integer 77 REAL(wp) :: zC2t_u, zC4t_u ! local scalars 78 REAL(wp) :: zC2t_v, zC4t_v ! - - 79 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwx, zwy, zwz, ztu, ztv, ztw 80 ! TEMP: This change not necessary after trd_tra is tiled 81 INTEGER :: itile 82 REAL(wp) :: zC2t_u, zC2t_v ! local scalars 83 REAL(wp), DIMENSION(A2D,jpk) :: zwx, zwy, zwz, ztu, ztv, ztw, zltu, zltv 84 ! TEMP: This change not necessary after trd_tra is tiled 85 REAL(wp), DIMENSION(:,:,:), SAVE, ALLOCATABLE :: ztrdx, ztrdy, ztrdz 80 86 !!---------------------------------------------------------------------- 81 ! 82 IF( kt == kit000 ) THEN 83 IF(lwp) WRITE(numout,*) 84 IF(lwp) WRITE(numout,*) 'tra_adv_cen : centered advection scheme on ', cdtype, ' order h/v =', kn_cen_h,'/', kn_cen_v 85 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~ ' 87 ! TEMP: This change not necessary after trd_tra is tiled 88 itile = ntile 89 ! 90 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 91 IF( kt == kit000 ) THEN 92 IF(lwp) WRITE(numout,*) 93 IF(lwp) WRITE(numout,*) 'tra_adv_cen : centered advection scheme on ', cdtype, ' order h/v =', kn_cen_h,'/', kn_cen_v 94 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~ ' 95 ENDIF 96 ! ! set local switches 97 l_trd = .FALSE. 98 l_hst = .FALSE. 99 l_ptr = .FALSE. 100 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE. 101 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) ) ) l_ptr = .TRUE. 102 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 103 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE. 104 105 ! TEMP: This can be A2D after trd_tra is tiled 106 IF( kt == kit000 .AND. l_trd ) THEN 107 ALLOCATE( ztrdx(jpi,jpj,jpk), ztrdy(jpi,jpj,jpk), ztrdz(jpi,jpj,jpk) ) 108 ENDIF 86 109 ENDIF 87 ! ! set local switches88 l_trd = .FALSE.89 l_hst = .FALSE.90 l_ptr = .FALSE.91 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE.92 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) ) ) l_ptr = .TRUE.93 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. &94 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE.95 110 ! 96 111 ! … … 109 124 ! 110 125 CASE( 4 ) !* 4th order centered 111 ztu(:,:,jpk) = 0._wp ! Bottom value : flux set to zero 112 ztv(:,:,jpk) = 0._wp 113 DO_3D_00_00( 1, jpkm1 ) 114 ztu(ji,jj,jk) = ( pt(ji+1,jj ,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * umask(ji,jj,jk) 115 ztv(ji,jj,jk) = ( pt(ji ,jj+1,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * vmask(ji,jj,jk) 116 END_3D 117 CALL lbc_lnk_multi( 'traadv_cen', ztu, 'U', -1. , ztv, 'V', -1. ) ! Lateral boundary cond. 118 ! 119 DO_3D_00_10( 1, jpkm1 ) 126 zltu(:,:,jpk) = 0._wp ! Bottom value : flux set to zero 127 zltv(:,:,jpk) = 0._wp 128 DO jk = 1, jpkm1 ! Laplacian 129 DO_2D_10_10 130 ztu(ji,jj,jk) = ( pt(ji+1,jj ,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * umask(ji,jj,jk) 131 ztv(ji,jj,jk) = ( pt(ji ,jj+1,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * vmask(ji,jj,jk) 132 END_2D 133 DO_2D_00_00 134 zltu(ji,jj,jk) = ztu(ji,jj,jk) + ztu(ji-1,jj,jk) 135 zltv(ji,jj,jk) = ztv(ji,jj,jk) + ztv(ji,jj-1,jk) 136 END_2D 137 END DO 138 CALL lbc_lnk_multi( 'traadv_cen', zltu, 'T', 1. , zltv, 'T', 1. ) ! Lateral boundary cond. (unchanged sgn) 139 ! 140 DO_3D_10_10( 1, jpkm1 ) 120 141 zC2t_u = pt(ji,jj,jk,jn,Kmm) + pt(ji+1,jj ,jk,jn,Kmm) ! C2 interpolation of T at u- & v-points (x2) 121 142 zC2t_v = pt(ji,jj,jk,jn,Kmm) + pt(ji ,jj+1,jk,jn,Kmm) 122 ! ! C4 interpolation of T at u- & v-points (x2)123 zC4t_u = zC2t_u + r1_6 * ( ztu(ji-1,jj,jk) - ztu(ji+1,jj,jk) )124 zC4t_v = zC2t_v + r1_6 * ( ztv(ji,jj-1,jk) - ztv(ji,jj+1,jk) )125 143 ! ! C4 fluxes 126 zwx(ji,jj,jk) = 0.5_wp * pU(ji,jj,jk) * zC4t_u127 zwy(ji,jj,jk) = 0.5_wp * pV(ji,jj,jk) * zC4t_v144 zwx(ji,jj,jk) = 0.5_wp * pU(ji,jj,jk) * ( zC2t_u + r1_6 * (zltu(ji,jj,jk) - zltu(ji+1,jj,jk)) ) 145 zwy(ji,jj,jk) = 0.5_wp * pV(ji,jj,jk) * ( zC2t_v + r1_6 * (zltv(ji,jj,jk) - zltv(ji,jj+1,jk)) ) 128 146 END_3D 129 147 ! … … 148 166 ! 149 167 IF( ln_linssh ) THEN !* top value (linear free surf. only as zwz is multiplied by wmask) 168 ! TODO: NOT TESTED- requires isf 150 169 IF( ln_isfcav ) THEN ! ice-shelf cavities (top of the ocean) 151 170 DO_2D_11_11 … … 153 172 END_2D 154 173 ELSE ! no ice-shelf cavities (only ocean surface) 155 zwz(:,:,1) = pW(:,:,1) * pt(:,:,1,jn,Kmm) 174 DO_2D_11_11 175 zwz(ji,jj,1) = pW(ji,jj,1) * pt(ji,jj,1,jn,Kmm) 176 END_2D 156 177 ENDIF 157 178 ENDIF … … 164 185 END_3D 165 186 ! ! trend diagnostics 187 ! TEMP: These changes not necessary after trd_tra is tiled 166 188 IF( l_trd ) THEN 167 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, zwx, pU, pt(:,:,:,jn,Kmm) ) 168 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, zwy, pV, pt(:,:,:,jn,Kmm) ) 169 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, zwz, pW, pt(:,:,:,jn,Kmm) ) 170 END IF 171 ! ! "Poleward" heat and salt transports 189 DO_3D_11_11( 1, jpk ) 190 ztrdx(ji,jj,jk) = zwx(ji,jj,jk) 191 ztrdy(ji,jj,jk) = zwy(ji,jj,jk) 192 ztrdz(ji,jj,jk) = zwz(ji,jj,jk) 193 END_3D 194 195 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 196 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 197 198 ! TODO: TO BE TILED- trd_tra 199 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, ztrdx, pU, pt(:,:,:,jn,Kmm) ) 200 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, ztrdy, pV, pt(:,:,:,jn,Kmm) ) 201 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, ztrdz, pW, pt(:,:,:,jn,Kmm) ) 202 203 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = itile ) ! Revert to tile domain 204 ENDIF 205 ENDIF 206 ! ! "Poleward" heat and salt transports 172 207 IF( l_ptr ) CALL dia_ptr_hst( jn, 'adv', zwy(:,:,:) ) 173 208 ! ! heat and salt transport -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/traadv_fct.F90
r12866 r13409 15 15 USE oce ! ocean dynamics and active tracers 16 16 USE dom_oce ! ocean space and time domain 17 ! TEMP: This change not necessary after trd_tra is tiled 18 USE domain, ONLY : dom_tile 17 19 USE trc_oce ! share passive tracers/Ocean variables 18 20 USE trd_oce ! trends: ocean variables … … 78 80 INTEGER , INTENT(in ) :: kn_fct_v ! order of the FCT scheme (=2 or 4) 79 81 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 82 ! TEMP: This can be A2D after trd_tra is tiled 80 83 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pU, pV, pW ! 3 ocean volume flux components 81 84 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) :: pt ! tracers and RHS of tracer equation 82 85 ! 83 86 INTEGER :: ji, jj, jk, jn ! dummy loop indices 87 ! TEMP: This change not necessary after trd_tra is tiled 88 INTEGER :: itile 84 89 REAL(wp) :: ztra ! local scalar 85 90 REAL(wp) :: zfp_ui, zfp_vj, zfp_wk, zC2t_u, zC4t_u ! - - 86 91 REAL(wp) :: zfm_ui, zfm_vj, zfm_wk, zC2t_v, zC4t_v ! - - 87 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwi, zwx, zwy, zwz, ztu, ztv, zltu, zltv, ztw 88 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ztrdx, ztrdy, ztrdz, zptry 89 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zwinf, zwdia, zwsup 92 REAL(wp), DIMENSION(A2D,jpk) :: zwi, zwx, zwy, zwz, ztu, ztv, zltu, zltv, ztw 93 ! TEMP: This change not necessary after trd_tra is tiled 94 REAL(wp), DIMENSION(:,:,:), SAVE, ALLOCATABLE :: ztrdx, ztrdy, ztrdz 95 REAL(wp), DIMENSION(:,:,:) , ALLOCATABLE :: zptry 96 REAL(wp), DIMENSION(:,:,:) , ALLOCATABLE :: zwinf, zwdia, zwsup 90 97 LOGICAL :: ll_zAimp ! flag to apply adaptive implicit vertical advection 91 98 !!---------------------------------------------------------------------- 92 ! 93 IF( kt == kit000 ) THEN 94 IF(lwp) WRITE(numout,*) 95 IF(lwp) WRITE(numout,*) 'tra_adv_fct : FCT advection scheme on ', cdtype 96 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 99 ! TEMP: This change not necessary after trd_tra is tiled 100 itile = ntile 101 ! 102 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 103 IF( kt == kit000 ) THEN 104 IF(lwp) WRITE(numout,*) 105 IF(lwp) WRITE(numout,*) 'tra_adv_fct : FCT advection scheme on ', cdtype 106 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 107 ENDIF 108 ! 109 l_trd = .FALSE. ! set local switches 110 l_hst = .FALSE. 111 l_ptr = .FALSE. 112 ll_zAimp = .FALSE. 113 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype =='TRC' .AND. l_trdtrc ) ) l_trd = .TRUE. 114 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) ) ) l_ptr = .TRUE. 115 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 116 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE. 117 ! 118 ! TEMP: This can be A2D after trd_tra is tiled 119 IF( kt == kit000 .AND. (l_trd .OR. l_hst) ) THEN 120 ALLOCATE( ztrdx(jpi,jpj,jpk), ztrdy(jpi,jpj,jpk), ztrdz(jpi,jpj,jpk) ) 121 ENDIF 97 122 ENDIF 98 123 ! 99 l_trd = .FALSE. ! set local switches100 l_hst = .FALSE.101 l_ptr = .FALSE.102 ll_zAimp = .FALSE.103 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype =='TRC' .AND. l_trdtrc ) ) l_trd = .TRUE.104 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) ) ) l_ptr = .TRUE.105 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. &106 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE.107 !108 IF( l_trd .OR. l_hst ) THEN109 ALLOCATE( ztrdx(jpi,jpj,jpk), ztrdy(jpi,jpj,jpk), ztrdz(jpi,jpj,jpk) )110 ztrdx(:,:,:) = 0._wp ; ztrdy(:,:,:) = 0._wp ; ztrdz(:,:,:) = 0._wp111 ENDIF112 !113 124 IF( l_ptr ) THEN 114 ALLOCATE( zptry( jpi,jpj,jpk) )125 ALLOCATE( zptry(A2D,jpk) ) 115 126 zptry(:,:,:) = 0._wp 116 127 ENDIF … … 123 134 ! If adaptive vertical advection, check if it is needed on this PE at this time 124 135 IF( ln_zad_Aimp ) THEN 125 IF( MAXVAL( ABS( wi( :,:,:) ) ) > 0._wp ) ll_zAimp = .TRUE.136 IF( MAXVAL( ABS( wi(A2D,:) ) ) > 0._wp ) ll_zAimp = .TRUE. 126 137 END IF 127 138 ! If active adaptive vertical advection, build tridiagonal matrix 128 139 IF( ll_zAimp ) THEN 129 ALLOCATE(zwdia( jpi,jpj,jpk), zwinf(jpi,jpj,jpk),zwsup(jpi,jpj,jpk))140 ALLOCATE(zwdia(A2D,jpk), zwinf(A2D,jpk), zwsup(A2D,jpk)) 130 141 DO_3D_00_00( 1, jpkm1 ) 131 142 zwdia(ji,jj,jk) = 1._wp + p2dt * ( MAX( wi(ji,jj,jk ) , 0._wp ) - MIN( wi(ji,jj,jk+1) , 0._wp ) ) / e3t(ji,jj,jk,Krhs) … … 155 166 END_3D 156 167 IF( ln_linssh ) THEN ! top ocean value (only in linear free surface as zwz has been w-masked) 168 ! TODO: NOT TESTED- requires isf 157 169 IF( ln_isfcav ) THEN ! top of the ice-shelf cavities and at the ocean surface 158 170 DO_2D_11_11 … … 193 205 END IF 194 206 ! 207 ! TEMP: This change not necessary after trd_tra is tiled 195 208 IF( l_trd .OR. l_hst ) THEN ! trend diagnostics (contribution of upstream fluxes) 196 ztrdx(:,:,:) = zwx(:,:,:) ; ztrdy(:,:,:) = zwy(:,:,:) ; ztrdz(:,:,:) = zwz(:,:,:) 209 DO_3D_11_11( 1, jpk ) 210 ztrdx(ji,jj,jk) = zwx(ji,jj,jk) ; ztrdy(ji,jj,jk) = zwy(ji,jj,jk) ; ztrdz(ji,jj,jk) = zwz(ji,jj,jk) 211 END_3D 197 212 END IF 198 213 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) … … 322 337 END IF 323 338 ! 339 ! TEMP: These changes not necessary after trd_tra is tiled 324 340 IF( l_trd .OR. l_hst ) THEN ! trend diagnostics // heat/salt transport 325 ztrdx(:,:,:) = ztrdx(:,:,:) + zwx(:,:,:) ! <<< add anti-diffusive fluxes 326 ztrdy(:,:,:) = ztrdy(:,:,:) + zwy(:,:,:) ! to upstream fluxes 327 ztrdz(:,:,:) = ztrdz(:,:,:) + zwz(:,:,:) ! 328 ! 329 IF( l_trd ) THEN ! trend diagnostics 330 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, ztrdx, pU, pt(:,:,:,jn,Kmm) ) 331 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, ztrdy, pV, pt(:,:,:,jn,Kmm) ) 332 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, ztrdz, pW, pt(:,:,:,jn,Kmm) ) 341 DO_3D_11_11( 1, jpk ) 342 ztrdx(ji,jj,jk) = ztrdx(ji,jj,jk) + zwx(ji,jj,jk) ! <<< add anti-diffusive fluxes 343 ztrdy(ji,jj,jk) = ztrdy(ji,jj,jk) + zwy(ji,jj,jk) ! to upstream fluxes 344 ztrdz(ji,jj,jk) = ztrdz(ji,jj,jk) + zwz(ji,jj,jk) ! 345 END_3D 346 ! 347 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 348 IF( l_trd ) THEN ! trend diagnostics 349 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 350 351 ! TODO: TO BE TILED- trd_tra 352 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, ztrdx, pU, pt(:,:,:,jn,Kmm) ) 353 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, ztrdy, pV, pt(:,:,:,jn,Kmm) ) 354 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, ztrdz, pW, pt(:,:,:,jn,Kmm) ) 355 356 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = itile ) ! Revert to tile domain 357 ENDIF 333 358 ENDIF 334 359 ! ! heat/salt transport 335 IF( l_hst ) CALL dia_ar5_hst( jn, 'adv', ztrdx( :,:,:), ztrdy(:,:,:) )360 IF( l_hst ) CALL dia_ar5_hst( jn, 'adv', ztrdx(A2D,:), ztrdy(A2D,:) ) 336 361 ! 337 362 ENDIF … … 346 371 DEALLOCATE( zwdia, zwinf, zwsup ) 347 372 ENDIF 348 IF( l_trd .OR. l_hst ) THEN 349 DEALLOCATE( ztrdx, ztrdy, ztrdz ) 350 ENDIF 373 ! TEMP: These changes not necessary after trd_tra is tiled 374 ! IF( l_trd .OR. l_hst ) THEN 375 ! DEALLOCATE( ztrdx, ztrdy, ztrdz ) 376 ! ENDIF 351 377 IF( l_ptr ) THEN 352 378 DEALLOCATE( zptry ) … … 369 395 !! in-space based differencing for fluid 370 396 !!---------------------------------------------------------------------- 371 INTEGER , INTENT(in ) :: Kmm ! time level index 372 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 373 REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(in ) :: pbef, paft ! before & after field 374 REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(inout) :: paa, pbb, pcc ! monotonic fluxes in the 3 directions 397 INTEGER , INTENT(in ) :: Kmm ! time level index 398 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 399 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in ) :: pbef ! before field 400 REAL(wp), DIMENSION(A2D ,jpk), INTENT(in ) :: paft ! after field 401 REAL(wp), DIMENSION(A2D ,jpk), INTENT(inout) :: paa, pbb, pcc ! monotonic fluxes in the 3 directions 375 402 ! 376 403 INTEGER :: ji, jj, jk ! dummy loop indices … … 378 405 REAL(wp) :: zpos, zneg, zbt, za, zb, zc, zbig, zrtrn ! local scalars 379 406 REAL(wp) :: zau, zbu, zcu, zav, zbv, zcv, zup, zdo ! - - 380 REAL(wp), DIMENSION( jpi,jpj,jpk) :: zbetup, zbetdo, zbup, zbdo407 REAL(wp), DIMENSION(A2D,jpk) :: zbetup, zbetdo, zbup, zbdo 381 408 !!---------------------------------------------------------------------- 382 409 ! … … 388 415 ! -------------------- 389 416 ! max/min of pbef & paft with large negative/positive value (-/+zbig) inside land 390 zbup = MAX( pbef * tmask - zbig * ( 1._wp - tmask ), & 391 & paft * tmask - zbig * ( 1._wp - tmask ) ) 392 zbdo = MIN( pbef * tmask + zbig * ( 1._wp - tmask ), & 393 & paft * tmask + zbig * ( 1._wp - tmask ) ) 417 DO_3D_11_11( 1, jpk ) 418 zbup(ji,jj,jk) = MAX( pbef(ji,jj,jk) * tmask(ji,jj,jk) - zbig * ( 1._wp - tmask(ji,jj,jk) ), & 419 & paft(ji,jj,jk) * tmask(ji,jj,jk) - zbig * ( 1._wp - tmask(ji,jj,jk) ) ) 420 zbdo(ji,jj,jk) = MIN( pbef(ji,jj,jk) * tmask(ji,jj,jk) + zbig * ( 1._wp - tmask(ji,jj,jk) ), & 421 & paft(ji,jj,jk) * tmask(ji,jj,jk) + zbig * ( 1._wp - tmask(ji,jj,jk) ) ) 422 END_3D 394 423 395 424 DO jk = 1, jpkm1 … … 523 552 !!---------------------------------------------------------------------- 524 553 REAL(wp),DIMENSION(jpi,jpj,jpk), INTENT(in ) :: pt_in ! field at t-point 525 REAL(wp),DIMENSION( jpi,jpj,jpk), INTENT( out) :: pt_out ! field interpolated at w-point554 REAL(wp),DIMENSION(A2D ,jpk), INTENT( out) :: pt_out ! field interpolated at w-point 526 555 ! 527 556 INTEGER :: ji, jj, jk ! dummy loop integers 528 557 INTEGER :: ikt, ikb ! local integers 529 REAL(wp),DIMENSION( jpi,jpj,jpk) :: zwd, zwi, zws, zwrm, zwt558 REAL(wp),DIMENSION(A2D,jpk) :: zwd, zwi, zws, zwrm, zwt 530 559 !!---------------------------------------------------------------------- 531 560 ! … … 547 576 !!gm 548 577 ! 578 ! TODO: NOT TESTED- requires isf 549 579 IF ( ln_isfcav ) THEN ! set level two values which may not be set in ISF case 550 580 zwd(:,:,2) = 1._wp ; zwi(:,:,2) = 0._wp ; zws(:,:,2) = 0._wp ; zwrm(:,:,2) = 0._wp … … 612 642 !! The 3d array zwt is used as a work space array. 613 643 !!---------------------------------------------------------------------- 614 REAL(wp),DIMENSION( :,:,:), INTENT(in ) :: pD, pU, PL ! 3-diagonal matrix615 REAL(wp),DIMENSION( :,:,:), INTENT(in ) :: pRHS ! Right-Hand-Side616 REAL(wp),DIMENSION( :,:,:), INTENT( out) :: pt_out !!gm field at level=F(klev)617 INTEGER , INTENT(in ) :: klev ! =1 pt_out at w-level618 ! ! =0 pt at t-level644 REAL(wp),DIMENSION(A2D,jpk), INTENT(in ) :: pD, pU, PL ! 3-diagonal matrix 645 REAL(wp),DIMENSION(A2D,jpk), INTENT(in ) :: pRHS ! Right-Hand-Side 646 REAL(wp),DIMENSION(A2D,jpk), INTENT( out) :: pt_out !!gm field at level=F(klev) 647 INTEGER , INTENT(in ) :: klev ! =1 pt_out at w-level 648 ! ! =0 pt at t-level 619 649 INTEGER :: ji, jj, jk ! dummy loop integers 620 650 INTEGER :: kstart ! local indices 621 REAL(wp),DIMENSION( jpi,jpj,jpk) :: zwt ! 3D work array651 REAL(wp),DIMENSION(A2D,jpk) :: zwt ! 3D work array 622 652 !!---------------------------------------------------------------------- 623 653 ! -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/traadv_mus.F90
r12810 r13409 19 19 USE trc_oce ! share passive tracers/Ocean variables 20 20 USE dom_oce ! ocean space and time domain 21 ! TEMP: This change not necessary after trd_tra is tiled 22 USE domain, ONLY : dom_tile 21 23 USE trd_oce ! trends: ocean variables 22 24 USE trdtra ! tracers trends manager … … 80 82 LOGICAL , INTENT(in ) :: ld_msc_ups ! use upstream scheme within muscl 81 83 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 84 ! TEMP: This can be A2D after trd_tra is tiled 82 85 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pU, pV, pW ! 3 ocean volume flux components 83 86 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) :: pt ! tracers and RHS of tracer equation 84 87 ! 88 ! TEMP: This change not necessary after trd_tra is tiled 89 INTEGER :: itile 85 90 INTEGER :: ji, jj, jk, jn ! dummy loop indices 86 91 INTEGER :: ierr ! local integer 87 92 REAL(wp) :: zu, z0u, zzwx, zw , zalpha ! local scalars 88 93 REAL(wp) :: zv, z0v, zzwy, z0w ! - - 89 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwx, zslpx ! 3D workspace 90 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwy, zslpy ! - - 94 REAL(wp), DIMENSION(A2D,jpk) :: zwx, zslpx ! 3D workspace 95 REAL(wp), DIMENSION(A2D,jpk) :: zwy, zslpy ! - - 96 ! TEMP: This change not necessary after trd_tra is tiled 97 REAL(wp), DIMENSION(:,:,:), SAVE, ALLOCATABLE :: ztrdx, ztrdy, ztrdz 91 98 !!---------------------------------------------------------------------- 92 ! 93 IF( kt == kit000 ) THEN 94 IF(lwp) WRITE(numout,*) 95 IF(lwp) WRITE(numout,*) 'tra_adv : MUSCL advection scheme on ', cdtype 96 IF(lwp) WRITE(numout,*) ' : mixed up-stream ', ld_msc_ups 97 IF(lwp) WRITE(numout,*) '~~~~~~~' 98 IF(lwp) WRITE(numout,*) 99 ! 100 ! Upstream / MUSCL scheme indicator 101 ! 102 ALLOCATE( xind(jpi,jpj,jpk), STAT=ierr ) 103 xind(:,:,:) = 1._wp ! set equal to 1 where up-stream is not needed 104 ! 105 IF( ld_msc_ups ) THEN ! define the upstream indicator (if asked) 106 ALLOCATE( upsmsk(jpi,jpj), STAT=ierr ) 107 upsmsk(:,:) = 0._wp ! not upstream by default 108 ! 109 DO jk = 1, jpkm1 110 xind(:,:,jk) = 1._wp & ! =>1 where up-stream is not needed 111 & - MAX ( rnfmsk(:,:) * rnfmsk_z(jk), & ! =>0 near runoff mouths (& closed sea outflows) 112 & upsmsk(:,:) ) * tmask(:,:,jk) ! =>0 in some user defined area 113 END DO 114 ENDIF 115 ! 116 ENDIF 117 ! 118 l_trd = .FALSE. 119 l_hst = .FALSE. 120 l_ptr = .FALSE. 121 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE. 122 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) ) ) l_ptr = .TRUE. 123 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 124 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE. 99 ! TEMP: This change not necessary after trd_tra is tiled 100 itile = ntile 101 ! 102 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 103 IF( kt == kit000 ) THEN 104 IF(lwp) WRITE(numout,*) 105 IF(lwp) WRITE(numout,*) 'tra_adv : MUSCL advection scheme on ', cdtype 106 IF(lwp) WRITE(numout,*) ' : mixed up-stream ', ld_msc_ups 107 IF(lwp) WRITE(numout,*) '~~~~~~~' 108 IF(lwp) WRITE(numout,*) 109 ! 110 ! Upstream / MUSCL scheme indicator 111 ! 112 ALLOCATE( xind(jpi,jpj,jpk), STAT=ierr ) 113 xind(:,:,:) = 1._wp ! set equal to 1 where up-stream is not needed 114 ! 115 IF( ld_msc_ups ) THEN ! define the upstream indicator (if asked) 116 ALLOCATE( upsmsk(jpi,jpj), STAT=ierr ) 117 upsmsk(:,:) = 0._wp ! not upstream by default 118 ! 119 DO jk = 1, jpkm1 120 xind(:,:,jk) = 1._wp & ! =>1 where up-stream is not needed 121 & - MAX ( rnfmsk(:,:) * rnfmsk_z(jk), & ! =>0 near runoff mouths (& closed sea outflows) 122 & upsmsk(:,:) ) * tmask(:,:,jk) ! =>0 in some user defined area 123 END DO 124 ENDIF 125 ! 126 ENDIF 127 ! 128 l_trd = .FALSE. 129 l_hst = .FALSE. 130 l_ptr = .FALSE. 131 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE. 132 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) ) ) l_ptr = .TRUE. 133 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 134 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE. 135 136 ! TEMP: This can be A2D after trd_tra is tiled 137 IF( kt == kit000 .AND. l_trd ) THEN 138 ALLOCATE( ztrdx(jpi,jpj,jpk), ztrdy(jpi,jpj,jpk), ztrdz(jpi,jpj,jpk) ) 139 ENDIF 140 ENDIF 125 141 ! 126 142 DO jn = 1, kjpt !== loop over the tracers ==! … … 180 196 END_3D 181 197 ! ! trend diagnostics 198 ! TEMP: These changes not necessary after trd_tra is tiled 182 199 IF( l_trd ) THEN 183 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, zwx, pU, pt(:,:,:,jn,Kbb) ) 184 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, zwy, pV, pt(:,:,:,jn,Kbb) ) 200 DO_3D_11_11( 1, jpk ) 201 ztrdx(ji,jj,jk) = zwx(ji,jj,jk) 202 ztrdy(ji,jj,jk) = zwy(ji,jj,jk) 203 END_3D 204 205 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 206 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 207 208 ! TODO: TO BE TILED- trd_tra 209 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, ztrdx, pU, pt(:,:,:,jn,Kbb) ) 210 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, ztrdy, pV, pt(:,:,:,jn,Kbb) ) 211 212 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = nijtile ) ! Revert to tile domain 213 ENDIF 185 214 END IF 186 215 ! ! "Poleward" heat and salt transports … … 194 223 zwx(:,:, 1 ) = 0._wp ! surface & bottom boundary conditions 195 224 zwx(:,:,jpk) = 0._wp 196 DO jk = 2, jpkm1! interior values197 zwx( :,:,jk) = tmask(:,:,jk) * ( pt(:,:,jk-1,jn,Kbb) - pt(:,:,jk,jn,Kbb) )198 END DO225 DO_3D_11_11( 2, jpkm1 ) ! interior values 226 zwx(ji,jj,jk) = tmask(ji,jj,jk) * ( pt(ji,jj,jk-1,jn,Kbb) - pt(ji,jj,jk,jn,Kbb) ) 227 END_3D 199 228 ! !-- Slopes of tracer 200 229 zslpx(:,:,1) = 0._wp ! surface values … … 217 246 END_3D 218 247 IF( ln_linssh ) THEN ! top values, linear free surface only 248 ! TODO: NOT TESTED- requires isf 219 249 IF( ln_isfcav ) THEN ! ice-shelf cavities (top of the ocean) 220 250 DO_2D_11_11 … … 222 252 END_2D 223 253 ELSE ! no cavities: only at the ocean surface 224 zwx(:,:,1) = pW(:,:,1) * pt(:,:,1,jn,Kbb) 254 DO_2D_11_11 255 zwx(ji,jj,1) = pW(ji,jj,1) * pt(ji,jj,1,jn,Kbb) 256 END_2D 225 257 ENDIF 226 258 ENDIF … … 230 262 END_3D 231 263 ! ! send trends for diagnostic 232 IF( l_trd ) CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, zwx, pW, pt(:,:,:,jn,Kbb) ) 264 ! TEMP: These changes not necessary after trd_tra is tiled 265 IF( l_trd ) THEN 266 DO_3D_11_11( 1, jpk ) 267 ztrdz(ji,jj,jk) = zwx(ji,jj,jk) 268 END_3D 269 270 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 271 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 272 273 ! TODO: TO BE TILED- trd_tra 274 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, ztrdz, pW, pt(:,:,:,jn,Kbb) ) 275 276 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = itile ) ! Revert to tile domain 277 ENDIF 278 ENDIF 233 279 ! 234 280 END DO ! end of tracer loop -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/traadv_qck.F90
r12377 r13409 17 17 USE oce ! ocean dynamics and active tracers 18 18 USE dom_oce ! ocean space and time domain 19 ! TEMP: This change not necessary after trd_tra is tiled 20 USE domain, ONLY : dom_tile 19 21 USE trc_oce ! share passive tracers/Ocean variables 20 22 USE trd_oce ! trends: ocean variables … … 90 92 INTEGER , INTENT(in ) :: kjpt ! number of tracers 91 93 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 94 ! TEMP: This can be A2D after trd_tra is tiled 92 95 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pU, pV, pW ! 3 ocean volume transport components 93 96 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) :: pt ! tracers and RHS of tracer equation 94 97 !!---------------------------------------------------------------------- 95 98 ! 96 IF( kt == kit000 ) THEN 97 IF(lwp) WRITE(numout,*) 98 IF(lwp) WRITE(numout,*) 'tra_adv_qck : 3rd order quickest advection scheme on ', cdtype 99 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~' 100 IF(lwp) WRITE(numout,*) 99 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 100 IF( kt == kit000 ) THEN 101 IF(lwp) WRITE(numout,*) 102 IF(lwp) WRITE(numout,*) 'tra_adv_qck : 3rd order quickest advection scheme on ', cdtype 103 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~' 104 IF(lwp) WRITE(numout,*) 105 ENDIF 106 ! 107 l_trd = .FALSE. 108 l_ptr = .FALSE. 109 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE. 110 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) ) ) l_ptr = .TRUE. 101 111 ENDIF 102 !103 l_trd = .FALSE.104 l_ptr = .FALSE.105 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE.106 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) ) ) l_ptr = .TRUE.107 112 ! 108 113 ! … … 126 131 INTEGER , INTENT(in ) :: kjpt ! number of tracers 127 132 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 133 ! TEMP: This can be A2D after trd_tra is tiled 128 134 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pU ! i-velocity components 129 135 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) :: pt ! active tracers and RHS of tracer equation 130 136 !! 137 ! TEMP: This change not necessary after trd_tra is tiled 138 INTEGER :: itile 131 139 INTEGER :: ji, jj, jk, jn ! dummy loop indices 132 140 REAL(wp) :: ztra, zbtr, zdir, zdx, zmsk ! local scalars 133 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwx, zfu, zfc, zfd 141 REAL(wp), DIMENSION(A2D,jpk) :: zwx, zfu, zfc, zfd 142 ! TEMP: This change not necessary after trd_tra is tiled 143 REAL(wp), DIMENSION(:,:,:), SAVE, ALLOCATABLE :: ztrdx 134 144 !---------------------------------------------------------------------- 135 ! 145 ! TEMP: This change not necessary after trd_tra is tiled 146 itile = ntile 147 ! 148 ! TEMP: This can be A2D after trd_tra is tiled 149 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 150 IF( kt == nit000 .AND. l_trd ) THEN 151 ALLOCATE( ztrdx(jpi,jpj,jpk) ) 152 ENDIF 153 ENDIF 136 154 ! ! =========== 137 155 DO jn = 1, kjpt ! tracer loop … … 199 217 END_3D 200 218 ! ! trend diagnostics 201 IF( l_trd ) CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, zwx, pU, pt(:,:,:,jn,Kmm) ) 219 ! TEMP: These changes not necessary after trd_tra is tiled 220 IF( l_trd ) THEN 221 DO_3D_11_11( 1, jpk ) 222 ztrdx(ji,jj,jk) = zwx(ji,jj,jk) 223 END_3D 224 225 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 226 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 227 228 ! TODO: TO BE TILED- trd_tra 229 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, ztrdx, pU, pt(:,:,:,jn,Kmm) ) 230 231 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = itile ) ! Revert to tile domain 232 ENDIF 233 ENDIF 202 234 ! 203 235 END DO … … 215 247 INTEGER , INTENT(in ) :: kjpt ! number of tracers 216 248 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 249 ! TEMP: This can be A2D after trd_tra is tiled 217 250 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pV ! j-velocity components 218 251 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) :: pt ! active tracers and RHS of tracer equation 219 252 !! 253 ! TEMP: This change not necessary after trd_tra is tiled 254 INTEGER :: itile 220 255 INTEGER :: ji, jj, jk, jn ! dummy loop indices 221 256 REAL(wp) :: ztra, zbtr, zdir, zdx, zmsk ! local scalars 222 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwy, zfu, zfc, zfd ! 3D workspace 257 REAL(wp), DIMENSION(A2D,jpk) :: zwy, zfu, zfc, zfd ! 3D workspace 258 ! TEMP: This change not necessary after trd_tra is tiled 259 REAL(wp), DIMENSION(:,:,:), SAVE, ALLOCATABLE :: ztrdy 223 260 !---------------------------------------------------------------------- 224 ! 261 ! TEMP: This change not necessary after trd_tra is tiled 262 itile = ntile 263 ! 264 ! TEMP: This can be A2D after trd_tra is tiled 265 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 266 IF( kt == nit000 .AND. l_trd ) THEN 267 ALLOCATE( ztrdy(jpi,jpj,jpk) ) 268 ENDIF 269 ENDIF 225 270 ! ! =========== 226 271 DO jn = 1, kjpt ! tracer loop … … 295 340 END_3D 296 341 ! ! trend diagnostics 297 IF( l_trd ) CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, zwy, pV, pt(:,:,:,jn,Kmm) ) 342 ! TEMP: These changes not necessary after trd_tra is tiled 343 IF( l_trd ) THEN 344 DO_3D_11_11( 1, jpk ) 345 ztrdy(ji,jj,jk) = zwy(ji,jj,jk) 346 END_3D 347 348 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 349 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 350 351 ! TODO: TO BE TILED- trd_tra 352 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, ztrdy, pV, pt(:,:,:,jn,Kmm) ) 353 354 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = itile ) ! Revert to tile domain 355 ENDIF 356 ENDIF 298 357 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) 299 358 IF( l_ptr ) CALL dia_ptr_hst( jn, 'adv', zwy(:,:,:) ) … … 312 371 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 313 372 INTEGER , INTENT(in ) :: kjpt ! number of tracers 314 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pW ! vertical velocity 373 ! TEMP: This can be A2D after trd_tra is tiled 374 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pW ! vertical velocity 315 375 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) :: pt ! active tracers and RHS of tracer equation 316 376 ! 377 ! TEMP: This change not necessary after trd_tra is tiled 378 INTEGER :: itile 317 379 INTEGER :: ji, jj, jk, jn ! dummy loop indices 318 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwz ! 3D workspace 319 !!---------------------------------------------------------------------- 320 ! 380 REAL(wp), DIMENSION(A2D,jpk) :: zwz ! 3D workspace 381 ! TEMP: This change not necessary after trd_tra is tiled 382 REAL(wp), DIMENSION(:,:,:), SAVE, ALLOCATABLE :: ztrdz 383 !!---------------------------------------------------------------------- 384 ! TEMP: This change not necessary after trd_tra is tiled 385 itile = ntile 386 ! 387 ! TEMP: This can be A2D after trd_tra is tiled 388 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 389 IF( kt == nit000 .AND. l_trd ) THEN 390 ALLOCATE( ztrdz(jpi,jpj,jpk) ) 391 ENDIF 392 ENDIF 393 321 394 zwz(:,:, 1 ) = 0._wp ! surface & bottom values set to zero for all tracers 322 395 zwz(:,:,jpk) = 0._wp … … 330 403 END_3D 331 404 IF( ln_linssh ) THEN !* top value (only in linear free surf. as zwz is multiplied by wmask) 405 ! TODO: NOT TESTED- requires isf 332 406 IF( ln_isfcav ) THEN ! ice-shelf cavities (top of the ocean) 333 407 DO_2D_11_11 … … 335 409 END_2D 336 410 ELSE ! no ocean cavities (only ocean surface) 337 zwz(:,:,1) = pW(:,:,1) * pt(:,:,1,jn,Kmm) 411 DO_2D_11_11 412 zwz(ji,jj,1) = pW(ji,jj,1) * pt(ji,jj,1,jn,Kmm) 413 END_2D 338 414 ENDIF 339 415 ENDIF … … 344 420 END_3D 345 421 ! ! Send trends for diagnostic 346 IF( l_trd ) CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, zwz, pW, pt(:,:,:,jn,Kmm) ) 422 ! TEMP: These changes not necessary after trd_tra is tiled 423 IF( l_trd ) THEN 424 DO_3D_11_11( 1, jpk ) 425 ztrdz(ji,jj,jk) = zwz(ji,jj,jk) 426 END_3D 427 428 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 429 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 430 431 ! TODO: TO BE TILED- trd_tra 432 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, ztrdz, pW, pt(:,:,:,jn,Kmm) ) 433 434 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = itile ) ! Revert to tile domain 435 ENDIF 436 ENDIF 347 437 ! 348 438 END DO … … 358 448 !! ** Method : 359 449 !!---------------------------------------------------------------------- 360 REAL(wp), DIMENSION( jpi,jpj,jpk), INTENT(in ) :: pfu ! second upwind point361 REAL(wp), DIMENSION( jpi,jpj,jpk), INTENT(in ) :: pfd ! first douwning point362 REAL(wp), DIMENSION( jpi,jpj,jpk), INTENT(in ) :: pfc ! the central point (or the first upwind point)363 REAL(wp), DIMENSION( jpi,jpj,jpk), INTENT(inout) :: puc ! input as Courant number ; output as flux450 REAL(wp), DIMENSION(A2D,jpk), INTENT(in ) :: pfu ! second upwind point 451 REAL(wp), DIMENSION(A2D,jpk), INTENT(in ) :: pfd ! first douwning point 452 REAL(wp), DIMENSION(A2D,jpk), INTENT(in ) :: pfc ! the central point (or the first upwind point) 453 REAL(wp), DIMENSION(A2D,jpk), INTENT(inout) :: puc ! input as Courant number ; output as flux 364 454 !! 365 455 INTEGER :: ji, jj, jk ! dummy loop indices -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/traadv_ubs.F90
r12377 r13409 14 14 USE oce ! ocean dynamics and active tracers 15 15 USE dom_oce ! ocean space and time domain 16 ! TEMP: This change not necessary after trd_tra is tiled 17 USE domain, ONLY : dom_tile 16 18 USE trc_oce ! share passive tracers/Ocean variables 17 19 USE trd_oce ! trends: ocean variables … … 91 93 INTEGER , INTENT(in ) :: kn_ubs_v ! number of tracers 92 94 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 95 ! TEMP: This can be A2D after trd_tra is tiled 93 96 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pU, pV, pW ! 3 ocean volume transport components 94 97 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) :: pt ! tracers and RHS of tracer equation 95 98 ! 99 ! TEMP: This change not necessary after trd_tra is tiled 100 INTEGER :: itile 96 101 INTEGER :: ji, jj, jk, jn ! dummy loop indices 97 102 REAL(wp) :: ztra, zbtr, zcoef ! local scalars 98 103 REAL(wp) :: zfp_ui, zfm_ui, zcenut, ztak, zfp_wk, zfm_wk ! - - 99 104 REAL(wp) :: zfp_vj, zfm_vj, zcenvt, zeeu, zeev, z_hdivn ! - - 100 REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztu, ztv, zltu, zltv, zti, ztw ! 3D workspace 101 !!---------------------------------------------------------------------- 102 ! 103 IF( kt == kit000 ) THEN 104 IF(lwp) WRITE(numout,*) 105 IF(lwp) WRITE(numout,*) 'tra_adv_ubs : horizontal UBS advection scheme on ', cdtype 106 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~' 105 REAL(wp), DIMENSION(A2D,jpk) :: ztu, ztv, zltu, zltv, zti, ztw ! 3D workspace 106 ! TEMP: This change not necessary after trd_tra is tiled 107 REAL(wp), DIMENSION(:,:,:), SAVE, ALLOCATABLE :: ztrdx, ztrdy, ztrdz 108 !!---------------------------------------------------------------------- 109 ! TEMP: This change not necessary after trd_tra is tiled 110 itile = ntile 111 ! 112 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 113 IF( kt == kit000 ) THEN 114 IF(lwp) WRITE(numout,*) 115 IF(lwp) WRITE(numout,*) 'tra_adv_ubs : horizontal UBS advection scheme on ', cdtype 116 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~' 117 ENDIF 118 ! 119 l_trd = .FALSE. 120 l_hst = .FALSE. 121 l_ptr = .FALSE. 122 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE. 123 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) ) ) l_ptr = .TRUE. 124 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 125 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE. 126 127 ! TEMP: This can be A2D after trd_tra is tiled 128 IF( kt == kit000 .AND. l_trd ) THEN 129 ALLOCATE( ztrdx(jpi,jpj,jpk), ztrdy(jpi,jpj,jpk), ztrdz(jpi,jpj,jpk) ) 130 ENDIF 107 131 ENDIF 108 !109 l_trd = .FALSE.110 l_hst = .FALSE.111 l_ptr = .FALSE.112 IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE.113 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) ) ) l_ptr = .TRUE.114 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. &115 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE.116 132 ! 117 133 ztw (:,:, 1 ) = 0._wp ! surface & bottom value : set to zero for all tracers … … 152 168 END_3D 153 169 ! 154 zltu(:,:,:) = pt(:,:,:,jn,Krhs) ! store the initial trends before its update 170 DO_3D_11_11( 1, jpk ) 171 zltu(ji,jj,jk) = pt(ji,jj,jk,jn,Krhs) ! store the initial trends before its update 172 END_3D 155 173 ! 156 174 DO jk = 1, jpkm1 !== add the horizontal advective trend ==! … … 163 181 END DO 164 182 ! 165 zltu(:,:,:) = pt(:,:,:,jn,Krhs) - zltu(:,:,:) ! Horizontal advective trend used in vertical 2nd order FCT case 166 ! ! and/or in trend diagnostic (l_trd=T) 167 ! 183 DO_3D_11_11( 1, jpk ) 184 zltu(ji,jj,jk) = pt(ji,jj,jk,jn,Krhs) - zltu(ji,jj,jk) ! Horizontal advective trend used in vertical 2nd order FCT case 185 END_3D ! and/or in trend diagnostic (l_trd=T) 186 ! 187 ! TEMP: These changes not necessary after trd_tra is tiled 168 188 IF( l_trd ) THEN ! trend diagnostics 169 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, ztu, pU, pt(:,:,:,jn,Kmm) ) 170 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, ztv, pV, pt(:,:,:,jn,Kmm) ) 189 DO_3D_11_11( 1, jpk ) 190 ztrdx(ji,jj,jk) = ztu(ji,jj,jk) 191 ztrdy(ji,jj,jk) = ztv(ji,jj,jk) 192 END_3D 193 194 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 195 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 196 197 ! TODO: TO BE TILED- trd_tra 198 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, ztrdx, pU, pt(:,:,:,jn,Kmm) ) 199 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, ztrdy, pV, pt(:,:,:,jn,Kmm) ) 200 201 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = itile ) ! Revert to tile domain 202 ENDIF 171 203 END IF 172 204 ! … … 183 215 CASE( 2 ) ! 2nd order FCT 184 216 ! 185 IF( l_trd ) zltv(:,:,:) = pt(:,:,:,jn,Krhs) ! store pt(:,:,:,:,Krhs) if trend diag. 217 IF( l_trd ) THEN 218 DO_3D_11_11( 1, jpk ) 219 zltv(ji,jj,jk) = pt(ji,jj,jk,jn,Krhs) ! store pt(:,:,:,:,Krhs) if trend diag. 220 END_3D 221 ENDIF 186 222 ! 187 223 ! !* upstream advection with initial mass fluxes & intermediate update ==! … … 192 228 END_3D 193 229 IF( ln_linssh ) THEN ! top ocean value (only in linear free surface as ztw has been w-masked) 230 ! TODO: NOT TESTED- requires isf 194 231 IF( ln_isfcav ) THEN ! top of the ice-shelf cavities and at the ocean surface 195 232 DO_2D_11_11 … … 197 234 END_2D 198 235 ELSE ! no cavities: only at the ocean surface 199 ztw(:,:,1) = pW(:,:,1) * pt(:,:,1,jn,Kbb) 236 DO_2D_11_11 237 ztw(ji,jj,1) = pW(ji,jj,1) * pt(ji,jj,1,jn,Kbb) 238 END_2D 200 239 ENDIF 201 240 ENDIF … … 223 262 ztw(ji,jj,jk) = pW(ji,jj,jk) * ztw(ji,jj,jk) * wmask(ji,jj,jk) 224 263 END_3D 225 IF( ln_linssh ) ztw(:,:, 1 ) = pW(:,:,1) * pt(:,:,1,jn,Kmm) !!gm ISF & 4th COMPACT doesn't work 264 IF( ln_linssh ) THEN 265 DO_2D_11_11 266 ztw(ji,jj,1) = pW(ji,jj,1) * pt(ji,jj,1,jn,Kmm) !!gm ISF & 4th COMPACT doesn't work 267 END_2D 268 ENDIF 226 269 ! 227 270 END SELECT … … 231 274 END_3D 232 275 ! 276 ! TEMP: These changes not necessary after trd_tra is tiled 233 277 IF( l_trd ) THEN ! vertical advective trend diagnostics 234 278 DO_3D_00_00( 1, jpkm1 ) 235 z ltv(ji,jj,jk) = pt(ji,jj,jk,jn,Krhs) - zltv(ji,jj,jk) &279 ztrdz(ji,jj,jk) = pt(ji,jj,jk,jn,Krhs) - zltv(ji,jj,jk) & 236 280 & + pt(ji,jj,jk,jn,Kmm) * ( pW(ji,jj,jk) - pW(ji,jj,jk+1) ) & 237 281 & * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 238 282 END_3D 239 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, zltv ) 283 284 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 285 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 286 287 ! TODO: TO BE TILED- trd_tra 288 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, ztrdz ) 289 290 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = nijtile ) ! Revert to tile domain 291 ENDIF 240 292 ENDIF 241 293 ! … … 258 310 !! in-space based differencing for fluid 259 311 !!---------------------------------------------------------------------- 260 INTEGER , INTENT(in ) 261 REAL(wp), INTENT(in ) 262 REAL(wp), DIMENSION 263 REAL(wp), INTENT(inout), DIMENSION (jpi,jpj,jpk) :: paft ! after field264 REAL(wp), INTENT(inout), DIMENSION (jpi,jpj,jpk) :: pcc ! monotonic flux in the k direction312 INTEGER , INTENT(in ) :: Kmm ! time level index 313 REAL(wp), INTENT(in ) :: p2dt ! tracer time-step 314 REAL(wp), DIMENSION(jpi,jpj,jpk) :: pbef ! before field 315 REAL(wp), INTENT(inout), DIMENSION(A2D ,jpk) :: paft ! after field 316 REAL(wp), INTENT(inout), DIMENSION(A2D ,jpk) :: pcc ! monotonic flux in the k direction 265 317 ! 266 318 INTEGER :: ji, jj, jk ! dummy loop indices 267 319 INTEGER :: ikm1 ! local integer 268 320 REAL(wp) :: zpos, zneg, zbt, za, zb, zc, zbig, zrtrn ! local scalars 269 REAL(wp), DIMENSION( jpi,jpj,jpk) :: zbetup, zbetdo! 3D workspace321 REAL(wp), DIMENSION(A2D,jpk) :: zbetup, zbetdo ! 3D workspace 270 322 !!---------------------------------------------------------------------- 271 323 ! … … 277 329 ! -------------------- 278 330 ! ! large negative value (-zbig) inside land 279 pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) - zbig * ( 1.e0 - tmask(:,:,:) ) 280 paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) - zbig * ( 1.e0 - tmask(:,:,:) ) 331 DO_3D_00_00( 1, jpk ) 332 pbef(ji,jj,jk) = pbef(ji,jj,jk) * tmask(ji,jj,jk) - zbig * ( 1.e0 - tmask(ji,jj,jk) ) 333 paft(ji,jj,jk) = paft(ji,jj,jk) * tmask(ji,jj,jk) - zbig * ( 1.e0 - tmask(ji,jj,jk) ) 334 END_3D 281 335 ! 282 336 DO jk = 1, jpkm1 ! search maximum in neighbourhood … … 289 343 END DO 290 344 ! ! large positive value (+zbig) inside land 291 pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) + zbig * ( 1.e0 - tmask(:,:,:) ) 292 paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) + zbig * ( 1.e0 - tmask(:,:,:) ) 345 DO_3D_00_00( 1, jpk ) 346 pbef(ji,jj,jk) = pbef(ji,jj,jk) * tmask(ji,jj,jk) + zbig * ( 1.e0 - tmask(ji,jj,jk) ) 347 paft(ji,jj,jk) = paft(ji,jj,jk) * tmask(ji,jj,jk) + zbig * ( 1.e0 - tmask(ji,jj,jk) ) 348 END_3D 293 349 ! 294 350 DO jk = 1, jpkm1 ! search minimum in neighbourhood … … 301 357 END DO 302 358 ! ! restore masked values to zero 303 pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) 304 paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) 359 DO_3D_00_00( 1, jpk ) 360 pbef(ji,jj,jk) = pbef(ji,jj,jk) * tmask(ji,jj,jk) 361 paft(ji,jj,jk) = paft(ji,jj,jk) * tmask(ji,jj,jk) 362 END_3D 305 363 ! 306 364 ! Positive and negative part of fluxes and beta terms -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/trabbc.F90
r12489 r13409 17 17 USE oce ! ocean variables 18 18 USE dom_oce ! domain: ocean 19 ! TEMP: This change not necessary after trd_tra is tiled 20 USE domain, ONLY : dom_tile 19 21 USE phycst ! physical constants 20 22 USE trd_oce ! trends: ocean variables … … 79 81 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts,jpt), INTENT(inout) :: pts ! active tracers and RHS of tracer equation 80 82 ! 81 INTEGER :: ji, jj ! dummy loop indices 82 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdt ! 3D workspace 83 INTEGER :: ji, jj, jk ! dummy loop indices 84 REAL(wp), SAVE :: zsum1 85 ! TEMP: This change not necessary after trd_tra is tiled 86 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ztrdt ! 3D workspace 83 87 !!---------------------------------------------------------------------- 84 88 ! 85 89 IF( ln_timing ) CALL timing_start('tra_bbc') 86 90 ! 87 IF( l_trdtra ) THEN ! Save the input temperature trend 88 ALLOCATE( ztrdt(jpi,jpj,jpk) ) 89 ztrdt(:,:,:) = pts(:,:,:,jp_tem,Krhs) 91 IF( l_trdtra ) THEN ! Save the input temperature trend 92 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 93 ! TEMP: This can be A2D after trd_tra is tiled 94 ALLOCATE( ztrdt(jpi,jpj,jpk) ) 95 ENDIF 96 97 DO_3D_00_00( 1, jpkm1 ) 98 ztrdt(ji,jj,jk) = pts(ji,jj,jk,jp_tem,Krhs) 99 END_3D 90 100 ENDIF 91 101 ! ! Add the geothermal trend on temperature … … 94 104 END_2D 95 105 ! 96 CALL lbc_lnk( 'trabbc', pts(:,:,:,jp_tem,Krhs) , 'T', 1. ) 97 ! 98 IF( l_trdtra ) THEN ! Send the trend for diagnostics 99 ztrdt(:,:,:) = pts(:,:,:,jp_tem,Krhs) - ztrdt(:,:,:) 100 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_bbc, ztrdt ) 101 DEALLOCATE( ztrdt ) 102 ENDIF 103 ! 104 CALL iom_put ( "hfgeou" , rho0_rcp * qgh_trd0(:,:) ) 105 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Krhs), clinfo1=' bbc - Ta: ', mask1=tmask, clinfo3='tra-ta' ) 106 ! TEMP: These changes not necessary after trd_tra is tiled, lbc_lnk not necessary if using XIOS (subdomain support, will not output haloes) 107 IF( l_trdtra ) THEN 108 DO_3D_00_00( 1, jpkm1 ) 109 ztrdt(ji,jj,jk) = pts(ji,jj,jk,jp_tem,Krhs) - ztrdt(ji,jj,jk) 110 END_3D 111 ENDIF 112 113 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 114 ! NOTE: I don't think pts needs to be the input here? 115 CALL lbc_lnk( 'trabbc', ztrdt(:,:,:) , 'T', 1. ) 116 ! 117 IF( l_trdtra ) THEN ! Send the trend for diagnostics 118 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 119 120 ! TODO: TO BE TILED- trd_tra 121 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_bbc, ztrdt ) 122 DEALLOCATE( ztrdt ) 123 124 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = nijtile ) ! Revert to tile domain 125 ENDIF 126 ! 127 CALL iom_put ( "hfgeou" , rho0_rcp * qgh_trd0(:,:) ) 128 ENDIF 129 130 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Krhs), clinfo1=' bbc - Ta: ', mask1=tmask, psum1=zsum1, & 131 & clinfo3='tra-ta' ) 106 132 ! 107 133 IF( ln_timing ) CALL timing_stop('tra_bbc') -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/trabbl.F90
r12377 r13409 26 26 USE oce ! ocean dynamics and active tracers 27 27 USE dom_oce ! ocean space and time domain 28 ! TEMP: This change not necessary after trd_tra is tiled 29 USE domain, ONLY : dom_tile 28 30 USE phycst ! physical constant 29 31 USE eosbn2 ! equation of state … … 105 107 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts,jpt), INTENT(inout) :: pts ! active tracers and RHS of tracer equation 106 108 ! 107 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdt, ztrds 109 INTEGER :: ji, jj, jk ! Dummy loop indices 110 REAL(wp), SAVE :: zsum1, zsum2, zsum3, zsum4 111 ! TEMP: This change not necessary after trd_tra is tiled 112 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ztrdt, ztrds 108 113 !!---------------------------------------------------------------------- 109 114 ! … … 111 116 ! 112 117 IF( l_trdtra ) THEN !* Save the T-S input trends 113 ALLOCATE( ztrdt(jpi,jpj,jpk) , ztrds(jpi,jpj,jpk) ) 114 ztrdt(:,:,:) = pts(:,:,:,jp_tem,Krhs) 115 ztrds(:,:,:) = pts(:,:,:,jp_sal,Krhs) 118 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 119 ! TEMP: This can be A2D after trd_tra is tiled 120 ALLOCATE( ztrdt(jpi,jpj,jpk), ztrds(jpi,jpj,jpk) ) 121 ENDIF 122 123 DO_3D_00_00( 1, jpkm1 ) 124 ztrdt(ji,jj,jk) = pts(ji,jj,jk,jp_tem,Krhs) 125 ztrds(ji,jj,jk) = pts(ji,jj,jk,jp_sal,Krhs) 126 END_3D 116 127 ENDIF 117 128 … … 122 133 CALL tra_bbl_dif( pts(:,:,:,:,Kbb), pts(:,:,:,:,Krhs), jpts, Kmm ) 123 134 IF( sn_cfctl%l_prtctl ) & 124 CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Krhs), clinfo1=' bbl_ldf - Ta: ', mask1=tmask, & 125 & tab3d_2=pts(:,:,:,jp_sal,Krhs), clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra' ) 126 ! lateral boundary conditions ; just need for outputs 127 CALL lbc_lnk_multi( 'trabbl', ahu_bbl, 'U', 1. , ahv_bbl, 'V', 1. ) 128 CALL iom_put( "ahu_bbl", ahu_bbl ) ! bbl diffusive flux i-coef 129 CALL iom_put( "ahv_bbl", ahv_bbl ) ! bbl diffusive flux j-coef 135 CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Krhs), clinfo1=' bbl_ldf - Ta: ', mask1=tmask, psum1=zsum1, & 136 & tab3d_2=pts(:,:,:,jp_sal,Krhs), clinfo2= ' Sa: ', mask2=tmask, psum2=zsum2, & 137 & clinfo3='tra' ) 138 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only on the last tile 139 ! lateral boundary conditions ; just need for outputs 140 CALL lbc_lnk_multi( 'trabbl', ahu_bbl, 'U', 1. , ahv_bbl, 'V', 1. ) 141 CALL iom_put( "ahu_bbl", ahu_bbl ) ! bbl diffusive flux i-coef 142 CALL iom_put( "ahv_bbl", ahv_bbl ) ! bbl diffusive flux j-coef 143 ENDIF 130 144 ! 131 145 ENDIF … … 135 149 CALL tra_bbl_adv( pts(:,:,:,:,Kbb), pts(:,:,:,:,Krhs), jpts, Kmm ) 136 150 IF(sn_cfctl%l_prtctl) & 137 CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Krhs), clinfo1=' bbl_adv - Ta: ', mask1=tmask, & 138 & tab3d_2=pts(:,:,:,jp_sal,Krhs), clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra' ) 139 ! lateral boundary conditions ; just need for outputs 140 CALL lbc_lnk_multi( 'trabbl', utr_bbl, 'U', 1. , vtr_bbl, 'V', 1. ) 141 CALL iom_put( "uoce_bbl", utr_bbl ) ! bbl i-transport 142 CALL iom_put( "voce_bbl", vtr_bbl ) ! bbl j-transport 143 ! 144 ENDIF 145 151 CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Krhs), clinfo1=' bbl_adv - Ta: ', mask1=tmask, psum1=zsum3, & 152 & tab3d_2=pts(:,:,:,jp_sal,Krhs), clinfo2= ' Sa: ', mask2=tmask, psum2=zsum4, & 153 & clinfo3='tra' ) 154 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only on the last tile 155 ! lateral boundary conditions ; just need for outputs 156 CALL lbc_lnk_multi( 'trabbl', utr_bbl, 'U', 1. , vtr_bbl, 'V', 1. ) 157 CALL iom_put( "uoce_bbl", utr_bbl ) ! bbl i-transport 158 CALL iom_put( "voce_bbl", vtr_bbl ) ! bbl j-transport 159 ENDIF 160 ! 161 ENDIF 162 163 ! TEMP: These changes not necessary after trd_tra is tiled 146 164 IF( l_trdtra ) THEN ! send the trends for further diagnostics 147 ztrdt(:,:,:) = pts(:,:,:,jp_tem,Krhs) - ztrdt(:,:,:) 148 ztrds(:,:,:) = pts(:,:,:,jp_sal,Krhs) - ztrds(:,:,:) 149 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_bbl, ztrdt ) 150 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_bbl, ztrds ) 151 DEALLOCATE( ztrdt, ztrds ) 165 DO_3D_00_00( 1, jpkm1 ) 166 ztrdt(ji,jj,jk) = pts(ji,jj,jk,jp_tem,Krhs) - ztrdt(ji,jj,jk) 167 ztrds(ji,jj,jk) = pts(ji,jj,jk,jp_sal,Krhs) - ztrds(ji,jj,jk) 168 END_3D 169 170 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 171 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 172 173 ! TODO: TO BE TILED- trd_tra 174 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_bbl, ztrdt ) 175 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_bbl, ztrds ) 176 DEALLOCATE( ztrdt, ztrds ) 177 178 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = nijtile ) ! Revert to tile domain 179 ENDIF 152 180 ENDIF 153 181 ! … … 186 214 INTEGER :: ik ! local integers 187 215 REAL(wp) :: zbtr ! local scalars 188 REAL(wp), DIMENSION( jpi,jpj) :: zptb ! workspace216 REAL(wp), DIMENSION(A2D) :: zptb ! workspace 189 217 !!---------------------------------------------------------------------- 190 218 ! … … 241 269 DO jn = 1, kjpt ! tracer loop 242 270 ! ! =========== 243 DO jj = 1, jpjm1 244 DO ji = 1, jpim1 ! CAUTION start from i=1 to update i=2 when cyclic east-west 245 IF( utr_bbl(ji,jj) /= 0.e0 ) THEN ! non-zero i-direction bbl advection 246 ! down-slope i/k-indices (deep) & up-slope i/k indices (shelf) 247 iid = ji + MAX( 0, mgrhu(ji,jj) ) ; iis = ji + 1 - MAX( 0, mgrhu(ji,jj) ) 248 ikud = mbku_d(ji,jj) ; ikus = mbku(ji,jj) 249 zu_bbl = ABS( utr_bbl(ji,jj) ) 250 ! 251 ! ! up -slope T-point (shelf bottom point) 252 zbtr = r1_e1e2t(iis,jj) / e3t(iis,jj,ikus,Kmm) 253 ztra = zu_bbl * ( pt(iid,jj,ikus,jn) - pt(iis,jj,ikus,jn) ) * zbtr 254 pt_rhs(iis,jj,ikus,jn) = pt_rhs(iis,jj,ikus,jn) + ztra 255 ! 256 DO jk = ikus, ikud-1 ! down-slope upper to down T-point (deep column) 257 zbtr = r1_e1e2t(iid,jj) / e3t(iid,jj,jk,Kmm) 258 ztra = zu_bbl * ( pt(iid,jj,jk+1,jn) - pt(iid,jj,jk,jn) ) * zbtr 259 pt_rhs(iid,jj,jk,jn) = pt_rhs(iid,jj,jk,jn) + ztra 260 END DO 261 ! 262 zbtr = r1_e1e2t(iid,jj) / e3t(iid,jj,ikud,Kmm) 263 ztra = zu_bbl * ( pt(iis,jj,ikus,jn) - pt(iid,jj,ikud,jn) ) * zbtr 264 pt_rhs(iid,jj,ikud,jn) = pt_rhs(iid,jj,ikud,jn) + ztra 265 ENDIF 266 ! 267 IF( vtr_bbl(ji,jj) /= 0.e0 ) THEN ! non-zero j-direction bbl advection 268 ! down-slope j/k-indices (deep) & up-slope j/k indices (shelf) 269 ijd = jj + MAX( 0, mgrhv(ji,jj) ) ; ijs = jj + 1 - MAX( 0, mgrhv(ji,jj) ) 270 ikvd = mbkv_d(ji,jj) ; ikvs = mbkv(ji,jj) 271 zv_bbl = ABS( vtr_bbl(ji,jj) ) 272 ! 273 ! up -slope T-point (shelf bottom point) 274 zbtr = r1_e1e2t(ji,ijs) / e3t(ji,ijs,ikvs,Kmm) 275 ztra = zv_bbl * ( pt(ji,ijd,ikvs,jn) - pt(ji,ijs,ikvs,jn) ) * zbtr 276 pt_rhs(ji,ijs,ikvs,jn) = pt_rhs(ji,ijs,ikvs,jn) + ztra 277 ! 278 DO jk = ikvs, ikvd-1 ! down-slope upper to down T-point (deep column) 279 zbtr = r1_e1e2t(ji,ijd) / e3t(ji,ijd,jk,Kmm) 280 ztra = zv_bbl * ( pt(ji,ijd,jk+1,jn) - pt(ji,ijd,jk,jn) ) * zbtr 281 pt_rhs(ji,ijd,jk,jn) = pt_rhs(ji,ijd,jk,jn) + ztra 282 END DO 283 ! ! down-slope T-point (deep bottom point) 284 zbtr = r1_e1e2t(ji,ijd) / e3t(ji,ijd,ikvd,Kmm) 285 ztra = zv_bbl * ( pt(ji,ijs,ikvs,jn) - pt(ji,ijd,ikvd,jn) ) * zbtr 286 pt_rhs(ji,ijd,ikvd,jn) = pt_rhs(ji,ijd,ikvd,jn) + ztra 287 ENDIF 288 END DO 271 DO_2D_10_10 ! CAUTION start from i=1 to update i=2 when cyclic east-west 272 IF( utr_bbl(ji,jj) /= 0.e0 ) THEN ! non-zero i-direction bbl advection 273 ! down-slope i/k-indices (deep) & up-slope i/k indices (shelf) 274 iid = ji + MAX( 0, mgrhu(ji,jj) ) ; iis = ji + 1 - MAX( 0, mgrhu(ji,jj) ) 275 ikud = mbku_d(ji,jj) ; ikus = mbku(ji,jj) 276 zu_bbl = ABS( utr_bbl(ji,jj) ) 277 ! 278 ! ! up -slope T-point (shelf bottom point) 279 zbtr = r1_e1e2t(iis,jj) / e3t(iis,jj,ikus,Kmm) 280 ztra = zu_bbl * ( pt(iid,jj,ikus,jn) - pt(iis,jj,ikus,jn) ) * zbtr 281 pt_rhs(iis,jj,ikus,jn) = pt_rhs(iis,jj,ikus,jn) + ztra 282 ! 283 DO jk = ikus, ikud-1 ! down-slope upper to down T-point (deep column) 284 zbtr = r1_e1e2t(iid,jj) / e3t(iid,jj,jk,Kmm) 285 ztra = zu_bbl * ( pt(iid,jj,jk+1,jn) - pt(iid,jj,jk,jn) ) * zbtr 286 pt_rhs(iid,jj,jk,jn) = pt_rhs(iid,jj,jk,jn) + ztra 287 END DO 288 ! 289 zbtr = r1_e1e2t(iid,jj) / e3t(iid,jj,ikud,Kmm) 290 ztra = zu_bbl * ( pt(iis,jj,ikus,jn) - pt(iid,jj,ikud,jn) ) * zbtr 291 pt_rhs(iid,jj,ikud,jn) = pt_rhs(iid,jj,ikud,jn) + ztra 292 ENDIF 289 293 ! 290 END DO 291 ! ! =========== 292 END DO ! end tracer 293 ! ! =========== 294 IF( vtr_bbl(ji,jj) /= 0.e0 ) THEN ! non-zero j-direction bbl advection 295 ! down-slope j/k-indices (deep) & up-slope j/k indices (shelf) 296 ijd = jj + MAX( 0, mgrhv(ji,jj) ) ; ijs = jj + 1 - MAX( 0, mgrhv(ji,jj) ) 297 ikvd = mbkv_d(ji,jj) ; ikvs = mbkv(ji,jj) 298 zv_bbl = ABS( vtr_bbl(ji,jj) ) 299 ! 300 ! up -slope T-point (shelf bottom point) 301 zbtr = r1_e1e2t(ji,ijs) / e3t(ji,ijs,ikvs,Kmm) 302 ztra = zv_bbl * ( pt(ji,ijd,ikvs,jn) - pt(ji,ijs,ikvs,jn) ) * zbtr 303 pt_rhs(ji,ijs,ikvs,jn) = pt_rhs(ji,ijs,ikvs,jn) + ztra 304 ! 305 DO jk = ikvs, ikvd-1 ! down-slope upper to down T-point (deep column) 306 zbtr = r1_e1e2t(ji,ijd) / e3t(ji,ijd,jk,Kmm) 307 ztra = zv_bbl * ( pt(ji,ijd,jk+1,jn) - pt(ji,ijd,jk,jn) ) * zbtr 308 pt_rhs(ji,ijd,jk,jn) = pt_rhs(ji,ijd,jk,jn) + ztra 309 END DO 310 ! ! down-slope T-point (deep bottom point) 311 zbtr = r1_e1e2t(ji,ijd) / e3t(ji,ijd,ikvd,Kmm) 312 ztra = zv_bbl * ( pt(ji,ijs,ikvs,jn) - pt(ji,ijd,ikvd,jn) ) * zbtr 313 pt_rhs(ji,ijd,ikvd,jn) = pt_rhs(ji,ijd,ikvd,jn) + ztra 314 ENDIF 315 END_2D 316 ! ! =========== 317 END DO ! end tracer 318 ! ! =========== 294 319 END SUBROUTINE tra_bbl_adv 295 320 … … 332 357 REAL(wp) :: za, zb, zgdrho ! local scalars 333 358 REAL(wp) :: zsign, zsigna, zgbbl ! - - 334 REAL(wp), DIMENSION(jpi,jpj,jpts) :: zts, zab ! 3D workspace 335 REAL(wp), DIMENSION(jpi,jpj) :: zub, zvb, zdep ! 2D workspace 336 !!---------------------------------------------------------------------- 337 ! 338 IF( kt == kit000 ) THEN 339 IF(lwp) WRITE(numout,*) 340 IF(lwp) WRITE(numout,*) 'trabbl:bbl : Compute bbl velocities and diffusive coefficients in ', cdtype 341 IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 359 REAL(wp), DIMENSION(A2D,jpts) :: zts, zab ! 3D workspace 360 REAL(wp), DIMENSION(A2D) :: zub, zvb, zdep ! 2D workspace 361 !!---------------------------------------------------------------------- 362 ! 363 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 364 IF( kt == kit000 ) THEN 365 IF(lwp) WRITE(numout,*) 366 IF(lwp) WRITE(numout,*) 'trabbl:bbl : Compute bbl velocities and diffusive coefficients in ', cdtype 367 IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 368 ENDIF 342 369 ENDIF 343 370 ! !* bottom variables (T, S, alpha, beta, depth, velocity) -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/tradmp.F90
r12738 r13409 24 24 USE oce ! ocean: variables 25 25 USE dom_oce ! ocean: domain variables 26 ! TEMP: This change not necessary after trd_tra is tiled 27 USE domain, ONLY : dom_tile 26 28 USE c1d ! 1D vertical configuration 27 29 USE trd_oce ! trends: ocean variables … … 95 97 ! 96 98 INTEGER :: ji, jj, jk, jn ! dummy loop indices 97 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts) :: zts_dta 98 REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE :: ztrdts 99 REAL(wp), DIMENSION(A2D,jpk,jpts) :: zts_dta 100 REAL(wp), SAVE :: zsum1, zsum2 101 ! TEMP: This change not necessary after trd_tra is tiled 102 REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE, SAVE :: ztrdts 99 103 !!---------------------------------------------------------------------- 100 104 ! … … 102 106 ! 103 107 IF( l_trdtra ) THEN !* Save ta and sa trends 104 ALLOCATE( ztrdts(jpi,jpj,jpk,jpts) ) 105 ztrdts(:,:,:,:) = pts(:,:,:,:,Krhs) 108 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 109 ! TEMP: This can be A2D after trd_tra is tiled 110 ALLOCATE( ztrdts(jpi,jpj,jpk,jpts) ) 111 ENDIF 112 113 DO_3D_00_00( 1, jpkm1 ) 114 ztrdts(ji,jj,jk,:) = pts(ji,jj,jk,:,Krhs) 115 END_3D 106 116 ENDIF 107 117 ! !== input T-S data at kt ==! … … 140 150 END SELECT 141 151 ! 152 ! TEMP: These changes not necessary after trd_tra is tiled 142 153 IF( l_trdtra ) THEN ! trend diagnostic 143 ztrdts(:,:,:,:) = pts(:,:,:,:,Krhs) - ztrdts(:,:,:,:) 144 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_dmp, ztrdts(:,:,:,jp_tem) ) 145 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_dmp, ztrdts(:,:,:,jp_sal) ) 146 DEALLOCATE( ztrdts ) 154 DO_3D_00_00( 1, jpkm1 ) 155 ztrdts(ji,jj,jk,:) = pts(ji,jj,jk,:,Krhs) - ztrdts(ji,jj,jk,:) 156 END_3D 157 158 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 159 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 160 161 ! TODO: TO BE TILED- trd_tra 162 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_dmp, ztrdts(:,:,:,jp_tem) ) 163 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_dmp, ztrdts(:,:,:,jp_sal) ) 164 DEALLOCATE( ztrdts ) 165 166 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = nijtile ) ! Revert to tile domain 167 ENDIF 147 168 ENDIF 148 169 ! ! Control print 149 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Krhs), clinfo1=' dmp - Ta: ', mask1=tmask, & 150 & tab3d_2=pts(:,:,:,jp_sal,Krhs), clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra' ) 170 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Krhs), clinfo1=' dmp - Ta: ', mask1=tmask, psum1=zsum1, & 171 & tab3d_2=pts(:,:,:,jp_sal,Krhs), clinfo2= ' Sa: ', mask2=tmask, psum2=zsum2, & 172 & clinfo3='tra' ) 151 173 ! 152 174 IF( ln_timing ) CALL timing_stop('tra_dmp') -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/traisf.F90
r12377 r13409 11 11 !!---------------------------------------------------------------------- 12 12 USE isf_oce ! Ice shelf variables 13 USE par_oce , ONLY : nijtile, ntile, ntsi, ntei, ntsj, ntej 13 14 USE dom_oce , ONLY : e3t, r1_e1e2t ! ocean space domain variables 14 15 USE isfutils, ONLY : debug ! debug option … … 30 31 CONTAINS 31 32 33 ! TODO: NOT TESTED- requires isf 32 34 SUBROUTINE tra_isf ( kt, Kmm, pts, Krhs ) 33 35 !!---------------------------------------------------------------------- … … 45 47 IF( ln_timing ) CALL timing_start('tra_isf') 46 48 ! 47 IF( kt == nit000 ) THEN 48 IF(lwp) WRITE(numout,*) 49 IF(lwp) WRITE(numout,*) 'tra_isf : Ice shelf heat fluxes' 50 IF(lwp) WRITE(numout,*) '~~~~~~~ ' 49 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 50 IF( kt == nit000 ) THEN 51 IF(lwp) WRITE(numout,*) 52 IF(lwp) WRITE(numout,*) 'tra_isf : Ice shelf heat fluxes' 53 IF(lwp) WRITE(numout,*) '~~~~~~~ ' 54 ENDIF 51 55 ENDIF 52 56 ! … … 75 79 ! 76 80 IF ( ln_isfdebug ) THEN 77 CALL debug('tra_isf: pts(:,:,:,:,Krhs) T', pts(:,:,:,1,Krhs)) 78 CALL debug('tra_isf: pts(:,:,:,:,Krhs) S', pts(:,:,:,2,Krhs)) 81 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 82 CALL debug('tra_isf: pts(:,:,:,:,Krhs) T', pts(:,:,:,1,Krhs)) 83 CALL debug('tra_isf: pts(:,:,:,:,Krhs) S', pts(:,:,:,2,Krhs)) 84 ENDIF 79 85 END IF 80 86 ! … … 83 89 END SUBROUTINE tra_isf 84 90 ! 91 ! TODO: NOT TESTED- requires isf 85 92 SUBROUTINE tra_isf_mlt(ktop, kbot, phtbl, pfrac, ptsc, ptsc_b, pts) 86 93 !!---------------------------------------------------------------------- … … 100 107 INTEGER :: ji,jj,jk ! loop index 101 108 INTEGER :: ikt, ikb ! top and bottom level of the tbl 102 REAL(wp), DIMENSION( jpi,jpj):: ztc ! total ice shelf tracer trend109 REAL(wp), DIMENSION(A2D) :: ztc ! total ice shelf tracer trend 103 110 !!---------------------------------------------------------------------- 104 111 ! 105 112 ! compute 2d total trend due to isf 106 ztc(:,:) = 0.5_wp * ( ptsc(:,:,jp_tem) + ptsc_b(:,:,jp_tem) ) / phtbl(:,:) 113 DO_2D_11_11 114 ztc(ji,jj) = 0.5_wp * ( ptsc(ji,jj,jp_tem) + ptsc_b(ji,jj,jp_tem) ) / phtbl(ji,jj) 115 END_2D 107 116 ! 108 117 ! update pts(:,:,:,:,Krhs) … … 124 133 END SUBROUTINE tra_isf_mlt 125 134 ! 135 ! TODO: NOT TESTED- requires isf 126 136 SUBROUTINE tra_isf_cpl( Kmm, ptsc, ptsa ) 127 137 !!---------------------------------------------------------------------- … … 136 146 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT(in ) :: ptsc 137 147 !!---------------------------------------------------------------------- 138 INTEGER :: j k148 INTEGER :: ji, jj, jk 139 149 !!---------------------------------------------------------------------- 140 150 ! 141 DO jk = 1,jpk142 ptsa( :,:,jk,jp_tem) = ptsa(:,:,jk,jp_tem) + ptsc(:,:,jk,jp_tem) * r1_e1e2t(:,:) / e3t(:,:,jk,Kmm)143 ptsa( :,:,jk,jp_sal) = ptsa(:,:,jk,jp_sal) + ptsc(:,:,jk,jp_sal) * r1_e1e2t(:,:) / e3t(:,:,jk,Kmm)144 END DO151 DO_3D_00_00( 1, jpk ) 152 ptsa(ji,jj,jk,jp_tem) = ptsa(ji,jj,jk,jp_tem) + ptsc(ji,jj,jk,jp_tem) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 153 ptsa(ji,jj,jk,jp_sal) = ptsa(ji,jj,jk,jp_sal) + ptsc(ji,jj,jk,jp_sal) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 154 END_3D 145 155 ! 146 156 END SUBROUTINE tra_isf_cpl -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/traldf.F90
r12942 r13409 17 17 USE oce ! ocean dynamics and tracers 18 18 USE dom_oce ! ocean space and time domain 19 ! TEMP: This change not necessary after trd_tra is tiled 20 USE domain, ONLY : dom_tile 19 21 USE phycst ! physical constants 20 22 USE ldftra ! lateral diffusion: eddy diffusivity & EIV coeff. … … 37 39 PUBLIC tra_ldf ! called by step.F90 38 40 PUBLIC tra_ldf_init ! called by nemogcm.F90 39 41 42 !! * Substitutions 43 # include "do_loop_substitute.h90" 40 44 !!---------------------------------------------------------------------- 41 45 !! NEMO/OCE 4.0 , NEMO Consortium (2018) … … 55 59 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts,jpt), INTENT(inout) :: pts ! active tracers and RHS of tracer equation 56 60 !! 57 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdt, ztrds 61 ! TEMP: This change not necessary after trd_tra is tiled 62 INTEGER :: itile 63 INTEGER :: ji, jj, jk ! dummy loop indices 64 REAL(wp), SAVE :: zsum1, zsum2 65 ! TEMP: This change not necessary after trd_tra is tiled 66 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ztrdt, ztrds 67 ! TEMP: This change not necessary after extra haloes development 68 LOGICAL :: lskip 58 69 !!---------------------------------------------------------------------- 59 70 ! 60 71 IF( ln_timing ) CALL timing_start('tra_ldf') 61 72 ! 62 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only after all tiles finish63 IF( l_trdtra ) THEN !* Save ta and sa trends 64 ! TODO: TO BE TILED65 ALLOCATE( ztrdt(jpi,jpj,jpk) , ztrds(jpi,jpj,jpk) )66 ztrdt(:,:,:) = pts(:,:,:,jp_tem,Krhs)67 ztrds(:,:,:) = pts(:,:,:,jp_sal,Krhs)73 lskip = .FALSE. 74 75 IF( l_trdtra ) THEN 76 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 77 ! TEMP: This can be A2D after trd_tra is tiled 78 ALLOCATE( ztrdt(jpi,jpj,jpk), ztrds(jpi,jpj,jpk) ) 68 79 ENDIF 69 80 ENDIF 70 ! 71 SELECT CASE ( nldf_tra ) !* compute lateral mixing trend and add it to the general trend 72 CASE ( np_lap ) ! laplacian: iso-level operator 73 CALL tra_ldf_lap ( kt, Kmm, nit000,'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, pts(:,:,:,:,Kbb), pts(:,:,:,:,Krhs), jpts, 1 ) 74 CASE ( np_lap_i ) ! laplacian: standard iso-neutral operator (Madec) 75 CALL tra_ldf_iso ( kt, Kmm, nit000,'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, pts(:,:,:,:,Kbb), pts(:,:,:,:,Kbb), pts(:,:,:,:,Krhs), jpts, 1 ) 76 CASE ( np_lap_it ) ! laplacian: triad iso-neutral operator (griffies) 77 CALL tra_ldf_triad( kt, Kmm, nit000,'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, pts(:,:,:,:,Kbb), pts(:,:,:,:,Kbb), pts(:,:,:,:,Krhs), jpts, 1 ) 78 CASE ( np_blp , np_blp_i , np_blp_it ) ! bilaplacian: iso-level & iso-neutral operators 79 CALL tra_ldf_blp ( kt, Kmm, nit000,'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, pts(:,:,:,:,Kbb), pts(:,:,:,:,Krhs), jpts, nldf_tra ) 80 END SELECT 81 ! 82 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only after all tiles finish 83 IF( l_trdtra ) THEN !* save the horizontal diffusive trends for further diagnostics 84 ! TODO: TO BE TILED 85 ztrdt(:,:,:) = pts(:,:,:,jp_tem,Krhs) - ztrdt(:,:,:) 86 ztrds(:,:,:) = pts(:,:,:,jp_sal,Krhs) - ztrds(:,:,:) 87 ! TODO: TO BE TILED 88 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_ldf, ztrdt ) 89 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_ldf, ztrds ) 90 DEALLOCATE( ztrdt, ztrds ) 81 82 ! TEMP: These changes not necessary after extra haloes development (lbc_lnk removed from tra_ldf_blp, zps_hde*) 83 IF( nldf_tra == np_blp .OR. nldf_tra == np_blp_i .OR. nldf_tra == np_blp_it ) THEN 84 IF( ln_tile ) THEN 85 IF( ntile == 1 ) THEN 86 CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) 87 ELSE 88 lskip = .TRUE. 89 ENDIF 91 90 ENDIF 92 91 ENDIF 92 IF( .NOT. lskip ) THEN 93 94 ! TEMP: This change not necessary after trd_tra is tiled 95 itile = ntile 96 97 IF( l_trdtra ) THEN !* Save ta and sa trends 98 DO_3D_00_00( 1, jpkm1 ) 99 ztrdt(ji,jj,jk) = pts(ji,jj,jk,jp_tem,Krhs) 100 ztrds(ji,jj,jk) = pts(ji,jj,jk,jp_sal,Krhs) 101 END_3D 102 ENDIF 103 ! 104 SELECT CASE ( nldf_tra ) !* compute lateral mixing trend and add it to the general trend 105 CASE ( np_lap ) ! laplacian: iso-level operator 106 CALL tra_ldf_lap ( kt, Kmm, nit000,'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, pts(:,:,:,:,Kbb), pts(:,:,:,:,Krhs), jpts, 1 ) 107 CASE ( np_lap_i ) ! laplacian: standard iso-neutral operator (Madec) 108 CALL tra_ldf_iso ( kt, Kmm, nit000,'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, pts(:,:,:,:,Kbb), pts(:,:,:,:,Kbb), pts(:,:,:,:,Krhs), jpts, 1 ) 109 CASE ( np_lap_it ) ! laplacian: triad iso-neutral operator (griffies) 110 CALL tra_ldf_triad( kt, Kmm, nit000,'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, pts(:,:,:,:,Kbb), pts(:,:,:,:,Kbb), pts(:,:,:,:,Krhs), jpts, 1 ) 111 CASE ( np_blp , np_blp_i , np_blp_it ) ! bilaplacian: iso-level & iso-neutral operators 112 CALL tra_ldf_blp ( kt, Kmm, nit000,'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, pts(:,:,:,:,Kbb), pts(:,:,:,:,Krhs), jpts, nldf_tra ) 113 END SELECT 114 ! 115 ! TEMP: These changes not necessary after trd_tra is tiled 116 IF( l_trdtra ) THEN !* save the horizontal diffusive trends for further diagnostics 117 DO_3D_00_00( 1, jpkm1 ) 118 ztrdt(ji,jj,jk) = pts(ji,jj,jk,jp_tem,Krhs) - ztrdt(ji,jj,jk) 119 ztrds(ji,jj,jk) = pts(ji,jj,jk,jp_sal,Krhs) - ztrds(ji,jj,jk) 120 END_3D 121 122 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 123 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 124 125 ! TODO: TO BE TILED- trd_tra 126 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_ldf, ztrdt ) 127 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_ldf, ztrds ) 128 DEALLOCATE( ztrdt, ztrds ) 129 130 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = itile ) ! Revert to tile domain 131 ENDIF 132 ENDIF 133 134 ! TEMP: This change not necessary after extra haloes development (lbc_lnk removed from tra_ldf_blp, zps_hde*) 135 IF( ln_tile .AND. ntile == 0 ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 1 ) 136 137 ENDIF 93 138 ! !* print mean trends (used for debugging) 94 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Krhs), clinfo1=' ldf - Ta: ', mask1=tmask, & 95 & tab3d_2=pts(:,:,:,jp_sal,Krhs), clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra' ) 139 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Krhs), clinfo1=' ldf - Ta: ', mask1=tmask, psum1=zsum1, & 140 & tab3d_2=pts(:,:,:,jp_sal,Krhs), clinfo2= ' Sa: ', mask2=tmask, psum2=zsum2, & 141 & clinfo3='tra' ) 96 142 ! 97 143 IF( ln_timing ) CALL timing_stop('tra_ldf') -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/traldf_iso.F90
r12942 r13409 19 19 USE oce ! ocean dynamics and active tracers 20 20 USE dom_oce ! ocean space and time domain 21 USE domutl, ONLY : is_tile 21 22 USE trc_oce ! share passive tracers/Ocean variables 22 23 USE zdf_oce ! ocean vertical physics … … 36 37 PUBLIC tra_ldf_iso ! routine called by step.F90 37 38 39 LOGICAL :: l_ptr ! flag to compute poleward transport 40 LOGICAL :: l_hst ! flag to compute heat transport 41 38 42 !! * Substitutions 39 43 # include "do_loop_substitute.h90" … … 45 49 CONTAINS 46 50 47 SUBROUTINE tra_ldf_iso( kt, Kmm, kit000, cdtype, pahu, pahv, & 48 & pgu , pgv , pgui, pgvi, & 49 & pt , pt2 , pt_rhs , kjpt , kpass ) 51 SUBROUTINE tra_ldf_iso( kt, Kmm, kit000, cdtype, pahu, pahv, & 52 & pgu , pgv , pgui, pgvi, & 53 & pt, pt2, pt_rhs, kjpt, kpass ) 54 !! 55 INTEGER , INTENT(in ) :: kt ! ocean time-step index 56 INTEGER , INTENT(in ) :: kit000 ! first time step index 57 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 58 INTEGER , INTENT(in ) :: kjpt ! number of tracers 59 INTEGER , INTENT(in ) :: kpass ! =1/2 first or second passage 60 INTEGER , INTENT(in ) :: Kmm ! ocean time level index 61 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pahu, pahv ! eddy diffusivity at u- and v-points [m2/s] 62 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pgu, pgv ! tracer gradient at pstep levels 63 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pgui, pgvi ! tracer gradient at top levels 64 REAL(wp), DIMENSION(:,:,:,:), INTENT(in ) :: pt ! tracer (kpass=1) or laplacian of tracer (kpass=2) 65 REAL(wp), DIMENSION(:,:,:,:), INTENT(in ) :: pt2 ! tracer (only used in kpass=2) 66 REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) :: pt_rhs ! tracer trend 67 !! 68 CALL tra_ldf_iso_t( kt, Kmm, kit000, cdtype, pahu, pahv, is_tile(pahu), & 69 & pgu , pgv , is_tile(pgu) , pgui, pgvi, is_tile(pgui), & 70 & pt, is_tile(pt), pt2, is_tile(pt2), pt_rhs, is_tile(pt_rhs), kjpt, kpass ) 71 END SUBROUTINE tra_ldf_iso 72 73 74 SUBROUTINE tra_ldf_iso_t( kt, Kmm, kit000, cdtype, pahu, pahv, ktah, & 75 & pgu , pgv , ktg , pgui, pgvi, ktgi, & 76 & pt, ktt, pt2, ktt2, pt_rhs, ktt_rhs, kjpt, kpass ) 50 77 !!---------------------------------------------------------------------- 51 78 !! *** ROUTINE tra_ldf_iso *** … … 88 115 !! ** Action : Update pt_rhs arrays with the before rotated diffusion 89 116 !!---------------------------------------------------------------------- 90 INTEGER , INTENT(in ) :: kt ! ocean time-step index 91 INTEGER , INTENT(in ) :: kit000 ! first time step index 92 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 93 INTEGER , INTENT(in ) :: kjpt ! number of tracers 94 INTEGER , INTENT(in ) :: kpass ! =1/2 first or second passage 95 INTEGER , INTENT(in ) :: Kmm ! ocean time level index 96 REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT(in ) :: pahu, pahv ! eddy diffusivity at u- and v-points [m2/s] 97 REAL(wp), DIMENSION(jpi,jpj ,kjpt), INTENT(in ) :: pgu, pgv ! tracer gradient at pstep levels 98 REAL(wp), DIMENSION(jpi,jpj, kjpt), INTENT(in ) :: pgui, pgvi ! tracer gradient at top levels 99 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt ! tracer (kpass=1) or laplacian of tracer (kpass=2) 100 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt2 ! tracer (only used in kpass=2) 101 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 102 ! 103 LOGICAL :: l_ptr ! flag to compute poleward transport 104 LOGICAL :: l_hst ! flag to compute heat transport 117 INTEGER , INTENT(in ) :: kt ! ocean time-step index 118 INTEGER , INTENT(in ) :: kit000 ! first time step index 119 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 120 INTEGER , INTENT(in ) :: kjpt ! number of tracers 121 INTEGER , INTENT(in ) :: kpass ! =1/2 first or second passage 122 INTEGER , INTENT(in ) :: Kmm ! ocean time level index 123 INTEGER , INTENT(in ) :: ktah, ktg, ktgi, ktt, ktt2, ktt_rhs 124 REAL(wp), DIMENSION(T2D(ktah) ,jpk) , INTENT(in ) :: pahu, pahv ! eddy diffusivity at u- and v-points [m2/s] 125 REAL(wp), DIMENSION(T2D(ktg) ,kjpt), INTENT(in ) :: pgu, pgv ! tracer gradient at pstep levels 126 REAL(wp), DIMENSION(T2D(ktgi) ,kjpt), INTENT(in ) :: pgui, pgvi ! tracer gradient at top levels 127 REAL(wp), DIMENSION(T2D(ktt) ,jpk,kjpt), INTENT(in ) :: pt ! tracer (kpass=1) or laplacian of tracer (kpass=2) 128 REAL(wp), DIMENSION(T2D(ktt2) ,jpk,kjpt), INTENT(in ) :: pt2 ! tracer (only used in kpass=2) 129 REAL(wp), DIMENSION(T2D(ktt_rhs),jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 130 ! 105 131 INTEGER :: ji, jj, jk, jn ! dummy loop indices 106 132 INTEGER :: ikt … … 115 141 IF( kpass == 1 .AND. kt == kit000 ) THEN 116 142 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 117 ! TODO: TO BE TILED118 143 IF(lwp) WRITE(numout,*) 119 144 IF(lwp) WRITE(numout,*) 'tra_ldf_iso : rotated laplacian diffusion operator on ', cdtype … … 121 146 ENDIF 122 147 ! 123 DO_3D_ 11_11( 1, jpk )148 DO_3D_00_00( 1, jpk ) 124 149 akz (ji,jj,jk) = 0._wp 125 150 ah_wslp2(ji,jj,jk) = 0._wp 126 151 END_3D 127 152 ENDIF 128 ! 129 l_hst = .FALSE. 130 l_ptr = .FALSE. 131 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) ) l_ptr = .TRUE. 132 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 133 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE. 153 ! 154 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 155 l_hst = .FALSE. 156 l_ptr = .FALSE. 157 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) ) l_ptr = .TRUE. 158 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 159 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE. 160 ENDIF 134 161 ! 135 162 ! … … 170 197 ! 171 198 IF( ln_traldf_blp ) THEN ! bilaplacian operator 172 DO_3D_ 10_10( 2, jpkm1 )199 DO_3D_00_00( 2, jpkm1 ) 173 200 akz(ji,jj,jk) = 16._wp * ah_wslp2(ji,jj,jk) & 174 201 & * ( akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ( e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) ) ) 175 202 END_3D 176 203 ELSEIF( ln_traldf_lap ) THEN ! laplacian operator 177 DO_3D_ 10_10( 2, jpkm1 )204 DO_3D_00_00( 2, jpkm1 ) 178 205 ze3w_2 = e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) 179 206 zcoef0 = rDt * ( akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ze3w_2 ) … … 183 210 ! 184 211 ELSE ! 33 flux set to zero with akz=ah_wslp2 ==>> computed in full implicit 185 DO_3D_ 11_11( 1, jpk )212 DO_3D_00_00( 1, jpk ) 186 213 akz(ji,jj,jk) = ah_wslp2(ji,jj,jk) 187 214 END_3D … … 206 233 zdjt(ji,jj,jk) = ( pt(ji ,jj+1,jk,jn) - pt(ji,jj,jk,jn) ) * vmask(ji,jj,jk) 207 234 END_3D 235 ! TODO: NOT TESTED- requires zps 208 236 IF( ln_zps ) THEN ! botton and surface ocean correction of the horizontal gradient 209 237 DO_2D_10_10 … … 211 239 zdjt(ji,jj,mbkv(ji,jj)) = pgv(ji,jj,jn) 212 240 END_2D 241 ! TODO: NOT TESTED- requires isf 213 242 IF( ln_isfcav ) THEN ! first wet level beneath a cavity 214 243 DO_2D_10_10 … … 326 355 ! ! "Poleward" diffusive heat or salt transports (T-S case only) 327 356 ! note sign is reversed to give down-gradient diffusive transports ) 328 ! TODO: TO BE TILED329 357 IF( l_ptr ) CALL dia_ptr_hst( jn, 'ldf', -zftv(:,:,:) ) 330 358 ! ! Diffusive heat transports 331 ! TODO: TO BE TILED332 359 IF( l_hst ) CALL dia_ar5_hst( jn, 'ldf', -zftu(:,:,:), -zftv(:,:,:) ) 333 360 ! … … 337 364 END DO ! end tracer loop 338 365 ! 339 END SUBROUTINE tra_ldf_iso 366 END SUBROUTINE tra_ldf_iso_t 340 367 341 368 !!============================================================================== -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/traldf_lap_blp.F90
r12377 r13409 13 13 USE oce ! ocean dynamics and active tracers 14 14 USE dom_oce ! ocean space and time domain 15 USE domutl, ONLY : is_tile 15 16 USE ldftra ! lateral physics: eddy diffusivity 16 17 USE traldf_iso ! iso-neutral lateral diffusion (standard operator) (tra_ldf_iso routine) … … 45 46 CONTAINS 46 47 47 SUBROUTINE tra_ldf_lap( kt, Kmm, kit000, cdtype, pahu, pahv , & 48 & pgu , pgv , pgui, pgvi, & 49 & pt , pt_rhs, kjpt, kpass ) 48 SUBROUTINE tra_ldf_lap( kt, Kmm, kit000, cdtype, pahu, pahv, & 49 & pgu , pgv , pgui, pgvi, & 50 & pt, pt_rhs, kjpt, kpass ) 51 !! 52 INTEGER , INTENT(in ) :: kt ! ocean time-step index 53 INTEGER , INTENT(in ) :: kit000 ! first time step index 54 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 55 INTEGER , INTENT(in ) :: kjpt ! number of tracers 56 INTEGER , INTENT(in ) :: kpass ! =1/2 first or second passage 57 INTEGER , INTENT(in ) :: Kmm ! ocean time level index 58 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pahu, pahv ! eddy diffusivity at u- and v-points [m2/s] 59 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pgu, pgv ! tracer gradient at pstep levels 60 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pgui, pgvi ! tracer gradient at top levels 61 REAL(wp), DIMENSION(:,:,:,:), INTENT(in ) :: pt ! before tracer fields 62 REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) :: pt_rhs ! tracer trend 63 !! 64 CALL tra_ldf_lap_t( kt, Kmm, kit000, cdtype, pahu, pahv, is_tile(pahu), & 65 & pgu , pgv , is_tile(pgu) , pgui, pgvi, is_tile(pgui), & 66 & pt, is_tile(pt), pt_rhs, is_tile(pt_rhs), kjpt, kpass ) 67 END SUBROUTINE tra_ldf_lap 68 69 70 SUBROUTINE tra_ldf_lap_t( kt, Kmm, kit000, cdtype, pahu, pahv, ktah, & 71 & pgu , pgv , ktg , pgui, pgvi, ktgi, & 72 & pt, ktt, pt_rhs, ktt_rhs, kjpt, kpass ) 50 73 !!---------------------------------------------------------------------- 51 74 !! *** ROUTINE tra_ldf_lap *** … … 71 94 INTEGER , INTENT(in ) :: kpass ! =1/2 first or second passage 72 95 INTEGER , INTENT(in ) :: Kmm ! ocean time level index 73 REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT(in ) :: pahu, pahv ! eddy diffusivity at u- and v-points [m2/s] 74 REAL(wp), DIMENSION(jpi,jpj ,kjpt), INTENT(in ) :: pgu, pgv ! tracer gradient at pstep levels 75 REAL(wp), DIMENSION(jpi,jpj, kjpt), INTENT(in ) :: pgui, pgvi ! tracer gradient at top levels 76 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt ! before tracer fields 77 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 96 INTEGER , INTENT(in ) :: ktah, ktg, ktgi, ktt, ktt_rhs 97 REAL(wp), DIMENSION(T2D(ktah), jpk) , INTENT(in ) :: pahu, pahv ! eddy diffusivity at u- and v-points [m2/s] 98 REAL(wp), DIMENSION(T2D(ktg), kjpt), INTENT(in ) :: pgu, pgv ! tracer gradient at pstep levels 99 REAL(wp), DIMENSION(T2D(ktgi), kjpt), INTENT(in ) :: pgui, pgvi ! tracer gradient at top levels 100 REAL(wp), DIMENSION(T2D(ktt), jpk,kjpt), INTENT(in ) :: pt ! before tracer fields 101 REAL(wp), DIMENSION(T2D(ktt_rhs),jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 78 102 ! 79 103 INTEGER :: ji, jj, jk, jn ! dummy loop indices 80 104 REAL(wp) :: zsign ! local scalars 81 REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztu, ztv, zaheeu, zaheev 82 !!---------------------------------------------------------------------- 83 ! 84 IF( kt == nit000 .AND. lwp ) THEN 85 WRITE(numout,*) 86 WRITE(numout,*) 'tra_ldf_lap : iso-level laplacian diffusion on ', cdtype, ', pass=', kpass 87 WRITE(numout,*) '~~~~~~~~~~~ ' 88 ENDIF 89 ! 90 l_hst = .FALSE. 91 l_ptr = .FALSE. 92 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) ) l_ptr = .TRUE. 93 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 94 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE. 105 REAL(wp), DIMENSION(A2D,jpk) :: ztu, ztv, zaheeu, zaheev 106 !!---------------------------------------------------------------------- 107 ! 108 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 109 IF( kt == nit000 .AND. lwp ) THEN 110 WRITE(numout,*) 111 WRITE(numout,*) 'tra_ldf_lap : iso-level laplacian diffusion on ', cdtype, ', pass=', kpass 112 WRITE(numout,*) '~~~~~~~~~~~ ' 113 ENDIF 114 ! 115 l_hst = .FALSE. 116 l_ptr = .FALSE. 117 IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) ) l_ptr = .TRUE. 118 IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 119 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) ) l_hst = .TRUE. 120 ENDIF 95 121 ! 96 122 ! !== Initialization of metric arrays used for all tracers ==! … … 111 137 ztv(ji,jj,jk) = zaheev(ji,jj,jk) * ( pt(ji ,jj+1,jk,jn) - pt(ji,jj,jk,jn) ) 112 138 END_3D 139 ! TODO: NOT TESTED- requires zps 113 140 IF( ln_zps ) THEN ! set gradient at bottom/top ocean level 114 141 DO_2D_10_10 … … 116 143 ztv(ji,jj,mbkv(ji,jj)) = zaheev(ji,jj,mbkv(ji,jj)) * pgv(ji,jj,jn) 117 144 END_2D 145 ! TODO: NOT TESTED- requires isf 118 146 IF( ln_isfcav ) THEN ! top in ocean cavities only 119 147 DO_2D_10_10 … … 141 169 ! ! ================== 142 170 ! 143 END SUBROUTINE tra_ldf_lap 171 END SUBROUTINE tra_ldf_lap_t 144 172 145 173 … … 172 200 ! 173 201 INTEGER :: ji, jj, jk, jn ! dummy loop indices 174 REAL(wp), DIMENSION( jpi,jpj,jpk,kjpt) :: zlap ! laplacian at t-point175 REAL(wp), DIMENSION( jpi,jpj, kjpt) :: zglu, zglv ! bottom GRADh of the laplacian (u- and v-points)176 REAL(wp), DIMENSION( jpi,jpj, kjpt) :: zgui, zgvi ! top GRADh of the laplacian (u- and v-points)202 REAL(wp), DIMENSION(A2D,jpk,kjpt) :: zlap ! laplacian at t-point 203 REAL(wp), DIMENSION(A2D, kjpt) :: zglu, zglv ! bottom GRADh of the laplacian (u- and v-points) 204 REAL(wp), DIMENSION(A2D, kjpt) :: zgui, zgvi ! top GRADh of the laplacian (u- and v-points) 177 205 !!--------------------------------------------------------------------- 178 206 ! 179 IF( kt == kit000 .AND. lwp ) THEN 180 WRITE(numout,*) 181 SELECT CASE ( kldf ) 182 CASE ( np_blp ) ; WRITE(numout,*) 'tra_ldf_blp : iso-level bilaplacian operator on ', cdtype 183 CASE ( np_blp_i ) ; WRITE(numout,*) 'tra_ldf_blp : iso-neutral bilaplacian operator on ', cdtype, ' (Standard)' 184 CASE ( np_blp_it ) ; WRITE(numout,*) 'tra_ldf_blp : iso-neutral bilaplacian operator on ', cdtype, ' (triad)' 185 END SELECT 186 WRITE(numout,*) '~~~~~~~~~~~' 207 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 208 IF( kt == kit000 .AND. lwp ) THEN 209 WRITE(numout,*) 210 SELECT CASE ( kldf ) 211 CASE ( np_blp ) ; WRITE(numout,*) 'tra_ldf_blp : iso-level bilaplacian operator on ', cdtype 212 CASE ( np_blp_i ) ; WRITE(numout,*) 'tra_ldf_blp : iso-neutral bilaplacian operator on ', cdtype, ' (Standard)' 213 CASE ( np_blp_it ) ; WRITE(numout,*) 'tra_ldf_blp : iso-neutral bilaplacian operator on ', cdtype, ' (triad)' 214 END SELECT 215 WRITE(numout,*) '~~~~~~~~~~~' 216 ENDIF 187 217 ENDIF 188 218 … … 201 231 CALL lbc_lnk( 'traldf_lap_blp', zlap(:,:,:,:) , 'T', 1. ) ! Lateral boundary conditions (unchanged sign) 202 232 ! ! Partial top/bottom cell: GRADh( zlap ) 233 ! TODO: NOT TESTED- requires zps and isf 203 234 IF( ln_isfcav .AND. ln_zps ) THEN ; CALL zps_hde_isf( kt, Kmm, kjpt, zlap, zglu, zglv, zgui, zgvi ) ! both top & bottom 204 235 ELSEIF( ln_zps ) THEN ; CALL zps_hde ( kt, Kmm, kjpt, zlap, zglu, zglv ) ! only bottom -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/traldf_triad.F90
r12489 r13409 13 13 USE oce ! ocean dynamics and active tracers 14 14 USE dom_oce ! ocean space and time domain 15 ! TEMP: This change not necessary if lbc_lnk is removed from ldf_eiv_dia and XIOS has subdomain support 16 USE domain, ONLY : dom_tile 17 USE domutl, ONLY : is_tile 15 18 USE phycst ! physical constants 16 19 USE trc_oce ! share passive tracers/Ocean variables … … 33 36 PUBLIC tra_ldf_triad ! routine called by traldf.F90 34 37 35 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, SAVE :: zdkt3d !: vertical tracer gradient at 2 levels36 37 38 LOGICAL :: l_ptr ! flag to compute poleward transport 38 39 LOGICAL :: l_hst ! flag to compute heat transport … … 48 49 CONTAINS 49 50 50 SUBROUTINE tra_ldf_triad( kt, Kmm, kit000, cdtype, pahu, pahv, & 51 & pgu , pgv , pgui, pgvi , & 52 & pt , pt2, pt_rhs, kjpt, kpass ) 51 SUBROUTINE tra_ldf_triad( kt, Kmm, kit000, cdtype, pahu, pahv, & 52 & pgu , pgv , pgui, pgvi, & 53 & pt, pt2, pt_rhs, kjpt, kpass ) 54 !! 55 INTEGER , INTENT(in ) :: kt ! ocean time-step index 56 INTEGER , INTENT(in ) :: kit000 ! first time step index 57 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 58 INTEGER , INTENT(in ) :: kjpt ! number of tracers 59 INTEGER , INTENT(in ) :: kpass ! =1/2 first or second passage 60 INTEGER , INTENT(in ) :: Kmm ! ocean time level indices 61 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pahu, pahv ! eddy diffusivity at u- and v-points [m2/s] 62 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pgu , pgv ! tracer gradient at pstep levels 63 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pgui, pgvi ! tracer gradient at top levels 64 REAL(wp), DIMENSION(:,:,:,:), INTENT(in ) :: pt ! tracer (kpass=1) or laplacian of tracer (kpass=2) 65 REAL(wp), DIMENSION(:,:,:,:), INTENT(in ) :: pt2 ! tracer (only used in kpass=2) 66 REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) :: pt_rhs ! tracer trend 67 !! 68 CALL tra_ldf_triad_t( kt, Kmm, kit000, cdtype, pahu, pahv, is_tile(pahu), & 69 & pgu , pgv , is_tile(pgu) , pgui, pgvi, is_tile(pgui), & 70 & pt, is_tile(pt), pt2, is_tile(pt2), pt_rhs, is_tile(pt_rhs), kjpt, kpass ) 71 END SUBROUTINE tra_ldf_triad 72 73 74 SUBROUTINE tra_ldf_triad_t( kt, Kmm, kit000, cdtype, pahu, pahv, ktah, & 75 & pgu , pgv , ktg , pgui, pgvi, ktgi, & 76 & pt, ktt, pt2, ktt2, pt_rhs, ktt_rhs, kjpt, kpass ) 53 77 !!---------------------------------------------------------------------- 54 78 !! *** ROUTINE tra_ldf_triad *** … … 76 100 INTEGER , INTENT(in ) :: kpass ! =1/2 first or second passage 77 101 INTEGER , INTENT(in) :: Kmm ! ocean time level indices 78 REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT(in ) :: pahu, pahv ! eddy diffusivity at u- and v-points [m2/s] 79 REAL(wp), DIMENSION(jpi,jpj ,kjpt), INTENT(in ) :: pgu , pgv ! tracer gradient at pstep levels 80 REAL(wp), DIMENSION(jpi,jpj, kjpt), INTENT(in ) :: pgui, pgvi ! tracer gradient at top levels 81 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt ! tracer (kpass=1) or laplacian of tracer (kpass=2) 82 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt2 ! tracer (only used in kpass=2) 83 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 102 INTEGER , INTENT(in ) :: ktah, ktg, ktgi, ktt, ktt2, ktt_rhs 103 REAL(wp), DIMENSION(T2D(ktah), jpk) , INTENT(in ) :: pahu, pahv ! eddy diffusivity at u- and v-points [m2/s] 104 REAL(wp), DIMENSION(T2D(ktg), kjpt), INTENT(in ) :: pgu , pgv ! tracer gradient at pstep levels 105 REAL(wp), DIMENSION(T2D(ktgi), kjpt), INTENT(in ) :: pgui, pgvi ! tracer gradient at top levels 106 REAL(wp), DIMENSION(T2D(ktt), jpk,kjpt), INTENT(in ) :: pt ! tracer (kpass=1) or laplacian of tracer (kpass=2) 107 REAL(wp), DIMENSION(T2D(ktt2), jpk,kjpt), INTENT(in ) :: pt2 ! tracer (only used in kpass=2) 108 REAL(wp), DIMENSION(T2D(ktt_rhs),jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 84 109 ! 85 110 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 93 118 REAL(wp) :: ze1ur, ze2vr, ze3wr, zdxt, zdyt, zdzt 94 119 REAL(wp) :: zah, zah_slp, zaei_slp 95 REAL(wp), DIMENSION(jpi,jpj ) :: z2d ! 2D workspace 96 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zdit, zdjt, zftu, zftv, ztfw, zpsi_uw, zpsi_vw ! 3D - 120 REAL(wp), DIMENSION(A2D,0:1) :: zdkt3d ! vertical tracer gradient at 2 levels 121 REAL(wp), DIMENSION(A2D ) :: z2d ! 2D workspace 122 REAL(wp), DIMENSION(A2D ,jpk) :: zdit, zdjt, zftu, zftv, ztfw ! 3D - 123 ! TEMP: This can be A2D if lbc_lnk is removed from ldf_eiv_dia and XIOS has subdomain support 124 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpsi_uw, zpsi_vw 97 125 !!---------------------------------------------------------------------- 98 126 ! 99 IF( .NOT.ALLOCATED(zdkt3d) ) THEN 100 ALLOCATE( zdkt3d(jpi,jpj,0:1) , STAT=ierr ) 101 CALL mpp_sum ( 'traldf_triad', ierr ) 102 IF( ierr > 0 ) CALL ctl_stop('STOP', 'tra_ldf_triad: unable to allocate arrays') 103 ENDIF 104 ! 105 IF( kpass == 1 .AND. kt == kit000 ) THEN 106 IF(lwp) WRITE(numout,*) 107 IF(lwp) WRITE(numout,*) 'tra_ldf_triad : rotated laplacian diffusion operator on ', cdtype 108 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~' 109 ENDIF 110 ! 111 l_hst = .FALSE. 112 l_ptr = .FALSE. 113 IF( cdtype == 'TRA' ) THEN 114 IF( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf') ) l_ptr = .TRUE. 115 IF( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 116 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) l_hst = .TRUE. 127 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 128 IF( kpass == 1 .AND. kt == kit000 ) THEN 129 IF(lwp) WRITE(numout,*) 130 IF(lwp) WRITE(numout,*) 'tra_ldf_triad : rotated laplacian diffusion operator on ', cdtype 131 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~' 132 ENDIF 133 ! 134 l_hst = .FALSE. 135 l_ptr = .FALSE. 136 IF( cdtype == 'TRA' ) THEN 137 IF( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf') ) l_ptr = .TRUE. 138 IF( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 139 & iom_use("uadv_salttr") .OR. iom_use("vadv_salttr") ) l_hst = .TRUE. 140 ENDIF 117 141 ENDIF 118 142 ! … … 127 151 IF( kpass == 1 ) THEN !== first pass only and whatever the tracer is ==! 128 152 ! 129 akz (:,:,:) = 0._wp 130 ah_wslp2(:,:,:) = 0._wp 131 IF( ln_ldfeiv_dia ) THEN 132 zpsi_uw(:,:,:) = 0._wp 133 zpsi_vw(:,:,:) = 0._wp 134 ENDIF 153 DO_3D_00_00( 1, jpk ) 154 akz (ji,jj,jk) = 0._wp 155 ah_wslp2(ji,jj,jk) = 0._wp 156 END_3D 135 157 ! 136 158 DO ip = 0, 1 ! i-k triads 137 159 DO kp = 0, 1 138 DO_3D_ 10_10( 1, jpkm1 )139 ze3wr = 1._wp / e3w(ji +ip,jj,jk+kp,Kmm)140 zbu = e1e2u(ji ,jj) * e3u(ji,jj,jk,Kmm)141 zah = 0.25_wp * pahu(ji ,jj,jk)142 zslope_skew = triadi_g(ji +ip,jj,jk,1-ip,kp)160 DO_3D_00_00( 1, jpkm1 ) 161 ze3wr = 1._wp / e3w(ji,jj,jk+kp,Kmm) 162 zbu = e1e2u(ji-ip,jj) * e3u(ji-ip,jj,jk,Kmm) 163 zah = 0.25_wp * pahu(ji-ip,jj,jk) 164 zslope_skew = triadi_g(ji,jj,jk,1-ip,kp) 143 165 ! Subtract s-coordinate slope at t-points to give slope rel to s-surfaces (do this by *adding* gradient of depth) 144 zslope2 = zslope_skew + ( gdept(ji +1,jj,jk,Kmm) - gdept(ji,jj,jk,Kmm) ) * r1_e1u(ji,jj) * umask(ji,jj,jk+kp)166 zslope2 = zslope_skew + ( gdept(ji-ip+1,jj,jk,Kmm) - gdept(ji-ip,jj,jk,Kmm) ) * r1_e1u(ji-ip,jj) * umask(ji-ip,jj,jk+kp) 145 167 zslope2 = zslope2 *zslope2 146 ah_wslp2(ji +ip,jj,jk+kp) = ah_wslp2(ji+ip,jj,jk+kp) + zah * zbu * ze3wr * r1_e1e2t(ji+ip,jj) * zslope2147 akz (ji +ip,jj,jk+kp) = akz (ji+ip,jj,jk+kp) + zah * r1_e1u(ji,jj) &148 & * r1_e1u(ji ,jj) * umask(ji,jj,jk+kp)168 ah_wslp2(ji,jj,jk+kp) = ah_wslp2(ji,jj,jk+kp) + zah * zbu * ze3wr * r1_e1e2t(ji,jj) * zslope2 169 akz (ji,jj,jk+kp) = akz (ji,jj,jk+kp) + zah * r1_e1u(ji-ip,jj) & 170 & * r1_e1u(ji-ip,jj) * umask(ji-ip,jj,jk+kp) 149 171 ! 150 IF( ln_ldfeiv_dia ) zpsi_uw(ji,jj,jk+kp) = zpsi_uw(ji,jj,jk+kp) &151 & + 0.25_wp * aeiu(ji,jj,jk) * e2u(ji,jj) * zslope_skew152 172 END_3D 153 173 END DO … … 156 176 DO jp = 0, 1 ! j-k triads 157 177 DO kp = 0, 1 158 DO_3D_ 10_10( 1, jpkm1 )159 ze3wr = 1.0_wp / e3w(ji,jj +jp,jk+kp,Kmm)160 zbv = e1e2v(ji,jj ) * e3v(ji,jj,jk,Kmm)161 zah = 0.25_wp * pahv(ji,jj ,jk)162 zslope_skew = triadj_g(ji,jj +jp,jk,1-jp,kp)178 DO_3D_00_00( 1, jpkm1 ) 179 ze3wr = 1.0_wp / e3w(ji,jj,jk+kp,Kmm) 180 zbv = e1e2v(ji,jj-jp) * e3v(ji,jj-jp,jk,Kmm) 181 zah = 0.25_wp * pahv(ji,jj-jp,jk) 182 zslope_skew = triadj_g(ji,jj,jk,1-jp,kp) 163 183 ! Subtract s-coordinate slope at t-points to give slope rel to s surfaces 164 184 ! (do this by *adding* gradient of depth) 165 zslope2 = zslope_skew + ( gdept(ji,jj +1,jk,Kmm) - gdept(ji,jj,jk,Kmm) ) * r1_e2v(ji,jj) * vmask(ji,jj,jk+kp)185 zslope2 = zslope_skew + ( gdept(ji,jj-jp+1,jk,Kmm) - gdept(ji,jj-jp,jk,Kmm) ) * r1_e2v(ji,jj-jp) * vmask(ji,jj-jp,jk+kp) 166 186 zslope2 = zslope2 * zslope2 167 ah_wslp2(ji,jj +jp,jk+kp) = ah_wslp2(ji,jj+jp,jk+kp) + zah * zbv * ze3wr * r1_e1e2t(ji,jj+jp) * zslope2168 akz (ji,jj +jp,jk+kp) = akz (ji,jj+jp,jk+kp) + zah * r1_e2v(ji,jj) &169 & * r1_e2v(ji,jj ) * vmask(ji,jj,jk+kp)187 ah_wslp2(ji,jj,jk+kp) = ah_wslp2(ji,jj,jk+kp) + zah * zbv * ze3wr * r1_e1e2t(ji,jj) * zslope2 188 akz (ji,jj,jk+kp) = akz (ji,jj,jk+kp) + zah * r1_e2v(ji,jj-jp) & 189 & * r1_e2v(ji,jj-jp) * vmask(ji,jj-jp,jk+kp) 170 190 ! 171 IF( ln_ldfeiv_dia ) zpsi_vw(ji,jj,jk+kp) = zpsi_vw(ji,jj,jk+kp) &172 & + 0.25 * aeiv(ji,jj,jk) * e1v(ji,jj) * zslope_skew173 191 END_3D 174 192 END DO … … 178 196 ! 179 197 IF( ln_traldf_blp ) THEN ! bilaplacian operator 180 DO_3D_ 10_10( 2, jpkm1 )198 DO_3D_00_00( 2, jpkm1 ) 181 199 akz(ji,jj,jk) = 16._wp * ah_wslp2(ji,jj,jk) & 182 200 & * ( akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ( e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) ) ) 183 201 END_3D 184 202 ELSEIF( ln_traldf_lap ) THEN ! laplacian operator 185 DO_3D_ 10_10( 2, jpkm1 )203 DO_3D_00_00( 2, jpkm1 ) 186 204 ze3w_2 = e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) 187 205 zcoef0 = rDt * ( akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ze3w_2 ) … … 191 209 ! 192 210 ELSE ! 33 flux set to zero with akz=ah_wslp2 ==>> computed in full implicit 193 akz(:,:,:) = ah_wslp2(:,:,:) 194 ENDIF 195 ! 196 IF( ln_ldfeiv_dia .AND. cdtype == 'TRA' ) CALL ldf_eiv_dia( zpsi_uw, zpsi_vw, Kmm ) 211 DO_3D_00_00( 1, jpk ) 212 akz(ji,jj,jk) = ah_wslp2(ji,jj,jk) 213 END_3D 214 ENDIF 215 ! 216 ! TEMP: These changes not necessary if lbc_lnk is removed from ldf_eiv_dia and XIOS has subdomain support 217 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 218 IF( ln_ldfeiv_dia .AND. cdtype == 'TRA' ) THEN 219 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) 220 221 zpsi_uw(:,:,:) = 0._wp 222 zpsi_vw(:,:,:) = 0._wp 223 224 DO jp = 0, 1 225 DO kp = 0, 1 226 DO_3D_10_10( 1, jpkm1 ) 227 zpsi_uw(ji,jj,jk+kp) = zpsi_uw(ji,jj,jk+kp) & 228 & + 0.25_wp * aeiu(ji,jj,jk) * e2u(ji,jj) * triadi_g(ji+jp,jj,jk,1-jp,kp) 229 zpsi_vw(ji,jj,jk+kp) = zpsi_vw(ji,jj,jk+kp) & 230 & + 0.25_wp * aeiv(ji,jj,jk) * e1v(ji,jj) * triadj_g(ji,jj+jp,jk,1-jp,kp) 231 END_3D 232 END DO 233 END DO 234 CALL ldf_eiv_dia( zpsi_uw, zpsi_vw, Kmm ) 235 236 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = nijtile ) 237 ENDIF 238 ENDIF 197 239 ! 198 240 ENDIF !== end 1st pass only ==! … … 211 253 zdjt(ji,jj,jk) = ( pt(ji ,jj+1,jk,jn) - pt(ji,jj,jk,jn) ) * vmask(ji,jj,jk) 212 254 END_3D 255 ! TODO: NOT TESTED- requires zps 213 256 IF( ln_zps .AND. l_grad_zps ) THEN ! partial steps: correction at top/bottom ocean level 214 257 DO_2D_10_10 … … 216 259 zdjt(ji,jj,mbkv(ji,jj)) = pgv(ji,jj,jn) 217 260 END_2D 261 ! TODO: NOT TESTED- requires isf 218 262 IF( ln_isfcav ) THEN ! top level (ocean cavities only) 219 263 DO_2D_10_10 … … 230 274 DO jk = 1, jpkm1 231 275 ! !== Vertical tracer gradient at level jk and jk+1 232 zdkt3d(:,:,1) = ( pt(:,:,jk,jn) - pt(:,:,jk+1,jn) ) * tmask(:,:,jk+1) 276 DO_2D_11_11 277 zdkt3d(ji,jj,1) = ( pt(ji,jj,jk,jn) - pt(ji,jj,jk+1,jn) ) * tmask(ji,jj,jk+1) 278 END_2D 233 279 ! 234 280 ! ! surface boundary condition: zdkt3d(jk=0)=zdkt3d(jk=1) 235 281 IF( jk == 1 ) THEN ; zdkt3d(:,:,0) = zdkt3d(:,:,1) 236 ELSE ; zdkt3d(:,:,0) = ( pt(:,:,jk-1,jn) - pt(:,:,jk,jn) ) * tmask(:,:,jk) 282 ELSE 283 DO_2D_11_11 284 zdkt3d(ji,jj,0) = ( pt(ji,jj,jk-1,jn) - pt(ji,jj,jk,jn) ) * tmask(ji,jj,jk) 285 END_2D 237 286 ENDIF 238 287 ! … … 374 423 END DO ! end tracer loop 375 424 ! ! =============== 376 END SUBROUTINE tra_ldf_triad 425 END SUBROUTINE tra_ldf_triad_t 377 426 378 427 !!============================================================================== -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/tramle.F90
r12489 r13409 78 78 !! Fox-Kemper and Ferrari, JPO, 38, 1166-1179, 2008 79 79 !!---------------------------------------------------------------------- 80 INTEGER 81 INTEGER 82 INTEGER 83 CHARACTER(len=3) 84 REAL(wp), DIMENSION( jpi,jpj,jpk), INTENT(inout) :: pu ! in : 3 ocean transport components85 REAL(wp), DIMENSION( jpi,jpj,jpk), INTENT(inout) :: pv ! out: same 3 transport components86 REAL(wp), DIMENSION( jpi,jpj,jpk), INTENT(inout) :: pw ! increased by the MLE induced transport80 INTEGER , INTENT(in ) :: kt ! ocean time-step index 81 INTEGER , INTENT(in ) :: kit000 ! first time step index 82 INTEGER , INTENT(in ) :: Kmm ! ocean time level index 83 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 84 REAL(wp), DIMENSION(A2D,jpk), INTENT(inout) :: pu ! in : 3 ocean transport components 85 REAL(wp), DIMENSION(A2D,jpk), INTENT(inout) :: pv ! out: same 3 transport components 86 REAL(wp), DIMENSION(A2D,jpk), INTENT(inout) :: pw ! increased by the MLE induced transport 87 87 ! 88 88 INTEGER :: ji, jj, jk ! dummy loop indices … … 90 90 REAL(wp) :: zcuw, zmuw, zc ! local scalar 91 91 REAL(wp) :: zcvw, zmvw ! - - 92 INTEGER , DIMENSION(jpi,jpj) :: inml_mle 93 REAL(wp), DIMENSION(jpi,jpj) :: zpsim_u, zpsim_v, zmld, zbm, zhu, zhv, zn2, zLf_NH, zLf_MH 94 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpsi_uw, zpsi_vw 92 INTEGER , DIMENSION(A2D) :: inml_mle 93 REAL(wp), DIMENSION(A2D) :: zpsim_u, zpsim_v, zmld, zbm, zhu, zhv, zn2, zLf_MH 94 REAL(wp), DIMENSION(A2D,jpk) :: zpsi_uw, zpsi_vw 95 ! TEMP: These changes not necessary if using XIOS (subdomain support) 96 REAL(wp), DIMENSION(:,:), ALLOCATABLE, SAVE :: zLf_NH 97 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, SAVE :: zpsiu_mle, zpsiv_mle 95 98 !!---------------------------------------------------------------------- 96 99 ! 97 100 ! !== MLD used for MLE ==! 98 101 ! ! compute from the 10m density to deal with the diurnal cycle 99 inml_mle(:,:) = mbkt(:,:) + 1 ! init. to number of ocean w-level (T-level + 1) 102 DO_2D_11_11 103 inml_mle(ji,jj) = mbkt(ji,jj) + 1 ! init. to number of ocean w-level (T-level + 1) 104 END_2D 100 105 IF ( nla10 > 0 ) THEN ! avoid case where first level is thicker than 10m 101 106 DO_3DS_11_11( jpkm1, nlb10, -1 ) … … 134 139 END SELECT 135 140 ! ! convert density into buoyancy 136 zbm(:,:) = + grav * zbm(:,:) / MAX( e3t(:,:,1,Kmm), zmld(:,:) ) 141 DO_2D_11_11 142 zbm(ji,jj) = + grav * zbm(ji,jj) / MAX( e3t(ji,jj,1,Kmm), zmld(ji,jj) ) 143 END_2D 137 144 ! 138 145 ! … … 205 212 END DO 206 213 214 ! TEMP: These changes not necessary if using XIOS (subdomain support) 207 215 IF( cdtype == 'TRA') THEN !== outputs ==! 208 ! 209 zLf_NH(:,:) = SQRT( rb_c * zmld(:,:) ) * r1_ft(:,:) ! Lf = N H / f 210 CALL iom_put( "Lf_NHpf" , zLf_NH ) ! Lf = N H / f 216 IF( kt == nit000 .AND. (ntile == 0 .OR. ntile == 1) ) THEN ! Do only on the first tile and timestep 217 ALLOCATE( zLf_NH(jpi,jpj), zpsiu_mle(jpi,jpj,jpk), zpsiv_mle(jpi,jpj,jpk) ) 218 ENDIF 219 ! 220 DO_2D_11_11 221 zLf_NH(ji,jj) = SQRT( rb_c * zmld(ji,jj) ) * r1_ft(ji,jj) ! Lf = N H / f 222 END_2D 211 223 ! 212 224 ! divide by cross distance to give streamfunction with dimensions m^2/s 213 DO jk = 1, ikmax+1 214 zpsi_uw(:,:,jk) = zpsi_uw(:,:,jk) * r1_e2u(:,:) 215 zpsi_vw(:,:,jk) = zpsi_vw(:,:,jk) * r1_e1v(:,:) 216 END DO 217 CALL iom_put( "psiu_mle", zpsi_uw ) ! i-mle streamfunction 218 CALL iom_put( "psiv_mle", zpsi_vw ) ! j-mle streamfunction 225 DO_3D_11_11( 1, ikmax+1 ) 226 zpsiu_mle(ji,jj,jk) = zpsi_uw(ji,jj,jk) * r1_e2u(ji,jj) 227 zpsiv_mle(ji,jj,jk) = zpsi_vw(ji,jj,jk) * r1_e1v(ji,jj) 228 END_3D 229 230 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only on the last tile 231 CALL iom_put( "Lf_NHpf" , zLf_NH ) ! Lf = N H / f 232 CALL iom_put( "psiu_mle", zpsiu_mle ) ! i-mle streamfunction 233 CALL iom_put( "psiv_mle", zpsiv_mle ) ! j-mle streamfunction 234 ENDIF 219 235 ENDIF 220 236 ! -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/tranpc.F90
r12489 r13409 17 17 USE oce ! ocean dynamics and active tracers 18 18 USE dom_oce ! ocean space and time domain 19 ! TEMP: This change not necessary after trd_tra is tiled and extra haloes development (lbc_lnk removed) 20 USE domain, ONLY : dom_tile 19 21 USE phycst ! physical constants 20 22 USE zdf_oce ! ocean vertical physics … … 33 35 PUBLIC tra_npc ! routine called by step.F90 34 36 37 INTEGER :: nnpcc ! number of statically instable water column 38 35 39 !! * Substitutions 36 40 # include "do_loop_substitute.h90" … … 63 67 ! 64 68 INTEGER :: ji, jj, jk ! dummy loop indices 65 INTEGER :: inpcc ! number of statically instable water column66 69 INTEGER :: jiter, ikbot, ikp, ikup, ikdown, ilayer, ik_low ! local integers 67 70 LOGICAL :: l_bottom_reached, l_column_treated … … 69 72 REAL(wp) :: zsa, zbeta, zsum_sali, zsum_beta, zbw, zrw, z1_rDt 70 73 REAL(wp), PARAMETER :: zn2_zero = 1.e-14_wp ! acceptance criteria for neutrality (N2==0) 71 REAL(wp), DIMENSION( jpk ) :: zvn2 ! vertical profile of N2 at 1 given point... 72 REAL(wp), DIMENSION( jpk,jpts) :: zvts, zvab ! vertical profile of T & S , and alpha & betaat 1 given point 73 REAL(wp), DIMENSION(jpi,jpj,jpk ) :: zn2 ! N^2 74 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts) :: zab ! alpha and beta 75 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdt, ztrds ! 3D workspace 74 REAL(wp), DIMENSION( jpk ) :: zvn2 ! vertical profile of N2 at 1 given point... 75 REAL(wp), DIMENSION( jpk,jpts) :: zvts, zvab ! vertical profile of T & S , and alpha & betaat 1 given point 76 REAL(wp), DIMENSION(A2D,jpk ) :: zn2 ! N^2 77 REAL(wp), DIMENSION(A2D,jpk,jpts) :: zab ! alpha and beta 78 ! TEMP: This change not necessary after trd_tra is tiled 79 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ztrdt, ztrds ! 3D workspace 76 80 ! 77 81 LOGICAL, PARAMETER :: l_LB_debug = .FALSE. ! set to true if you want to follow what is … … 81 85 ! 82 86 IF( ln_timing ) CALL timing_start('tra_npc') 87 88 IF( l_trdtra ) THEN 89 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 90 ! TEMP: This can be A2D after trd_tra is tiled 91 ALLOCATE( ztrdt(jpi,jpj,jpk), ztrds(jpi,jpj,jpk) ) 92 ENDIF 93 ENDIF 83 94 ! 84 95 IF( MOD( kt, nn_npc ) == 0 ) THEN 85 96 ! 86 97 IF( l_trdtra ) THEN !* Save initial after fields 87 ALLOCATE( ztrdt(jpi,jpj,jpk) , ztrds(jpi,jpj,jpk) ) 88 ztrdt(:,:,:) = pts(:,:,:,jp_tem,Kaa) 89 ztrds(:,:,:) = pts(:,:,:,jp_sal,Kaa) 90 ENDIF 91 ! 98 DO_3D_00_00( 1, jpkm1 ) 99 ztrdt(ji,jj,jk) = pts(ji,jj,jk,jp_tem,Kaa) 100 ztrds(ji,jj,jk) = pts(ji,jj,jk,jp_sal,Kaa) 101 END_3D 102 ENDIF 103 ! 104 ! TODO: NOT TESTED- requires ORCA2 92 105 IF( l_LB_debug ) THEN 93 106 ! Location of 1 known convection site to follow what's happening in the water column … … 100 113 CALL bn2 ( pts(:,:,:,:,Kaa), zab, zn2, Kmm ) ! after Brunt-Vaisala (given on W-points) 101 114 ! 102 inpcc = 0115 IF( ntile == 0 .OR. ntile == 1 ) nnpcc = 0 ! Do only on the first tile 103 116 ! 104 117 DO_2D_00_00 … … 159 172 ENDIF 160 173 ! 161 IF( jiter == 1 ) inpcc = inpcc + 1174 IF( jiter == 1 ) nnpcc = nnpcc + 1 162 175 ! 163 176 IF( lp_monitor_point ) WRITE(numout, *) 'Negative N2 at ikp =',ikp,' for layer #', ilayer … … 300 313 END_2D 301 314 ! 302 IF( l_trdtra ) THEN ! send the Non penetrative mixing trends for diagnostic 315 ! TEMP: These changes not necessary after trd_tra is tiled and extra haloes development (lbc_lnk removed) 316 IF( l_trdtra ) THEN 303 317 z1_rDt = 1._wp / (2._wp * rn_Dt) 304 ztrdt(:,:,:) = ( pts(:,:,:,jp_tem,Kaa) - ztrdt(:,:,:) ) * z1_rDt 305 ztrds(:,:,:) = ( pts(:,:,:,jp_sal,Kaa) - ztrds(:,:,:) ) * z1_rDt 306 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_npc, ztrdt ) 307 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_npc, ztrds ) 308 DEALLOCATE( ztrdt, ztrds ) 309 ENDIF 310 ! 311 CALL lbc_lnk_multi( 'tranpc', pts(:,:,:,jp_tem,Kaa), 'T', 1., pts(:,:,:,jp_sal,Kaa), 'T', 1. ) 312 ! 313 IF( lwp .AND. l_LB_debug ) THEN 314 WRITE(numout,*) 'Exiting tra_npc , kt = ',kt,', => numb. of statically instable water-columns: ', inpcc 315 WRITE(numout,*) 318 319 DO_3D_00_00( 1, jpkm1 ) 320 ztrdt(ji,jj,jk) = ( pts(ji,jj,jk,jp_tem,Kaa) - ztrdt(ji,jj,jk) ) * z1_rDt 321 ztrds(ji,jj,jk) = ( pts(ji,jj,jk,jp_sal,Kaa) - ztrds(ji,jj,jk) ) * z1_rDt 322 END_3D 323 ENDIF 324 325 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 326 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 327 328 IF( l_trdtra ) THEN ! send the Non penetrative mixing trends for diagnostic 329 ! TODO: TO BE TILED- trd_tra 330 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_npc, ztrdt ) 331 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_npc, ztrds ) 332 DEALLOCATE( ztrdt, ztrds ) 333 ENDIF 334 335 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = nijtile ) ! Revert to tile domain 336 ! 337 CALL lbc_lnk_multi( 'tranpc', pts(:,:,:,jp_tem,Kaa), 'T', 1., pts(:,:,:,jp_sal,Kaa), 'T', 1. ) 338 ! 339 IF( lwp .AND. l_LB_debug ) THEN 340 WRITE(numout,*) 'Exiting tra_npc , kt = ',kt,', => numb. of statically instable water-columns: ', nnpcc 341 WRITE(numout,*) 342 ENDIF 316 343 ENDIF 317 344 ! -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/traqsr.F90
r12738 r13409 22 22 USE phycst ! physical constants 23 23 USE dom_oce ! ocean space and time domain 24 USE domain, ONLY : dom_tile 24 25 USE sbc_oce ! surface boundary condition: ocean 25 26 USE trc_oce ! share SMS/Ocean variables … … 112 113 REAL(wp) :: zz0 , zz1 ! - - 113 114 REAL(wp) :: zCb, zCmax, zze, zpsi, zpsimax, zdelpsi, zCtot, zCze 114 REAL(wp) :: zlogc, zlogc2, zlogc3 115 REAL(wp) :: zlogc, zlogc2, zlogc3 116 REAL(wp), SAVE :: zsum1 115 117 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zekb, zekg, zekr 116 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ze0, ze1, ze2, ze3, zea, ztrdt 118 ! TEMP: These changes not necessary after trd_tra is tiled 119 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ztrdt 120 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ze0, ze1, ze2, ze3, zea 117 121 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zetot, zchl3d 118 122 !!---------------------------------------------------------------------- … … 120 124 IF( ln_timing ) CALL timing_start('tra_qsr') 121 125 ! 122 IF( kt == nit000 ) THEN 123 IF(lwp) WRITE(numout,*) 124 IF(lwp) WRITE(numout,*) 'tra_qsr : penetration of the surface solar radiation' 125 IF(lwp) WRITE(numout,*) '~~~~~~~' 126 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 127 IF( kt == nit000 ) THEN 128 IF(lwp) WRITE(numout,*) 129 IF(lwp) WRITE(numout,*) 'tra_qsr : penetration of the surface solar radiation' 130 IF(lwp) WRITE(numout,*) '~~~~~~~' 131 ENDIF 126 132 ENDIF 127 133 ! 128 134 IF( l_trdtra ) THEN ! trends diagnostic: save the input temperature trend 129 ALLOCATE( ztrdt(jpi,jpj,jpk) ) 130 ztrdt(:,:,:) = pts(:,:,:,jp_tem,Krhs) 135 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 136 ! TEMP: This can be A2D after trd_tra is tiled 137 ALLOCATE( ztrdt(jpi,jpj,jpk) ) 138 ENDIF 139 140 DO_3D_00_00( 1, jpkm1 ) 141 ztrdt(ji,jj,jk) = pts(ji,jj,jk,jp_tem,Krhs) 142 END_3D 131 143 ENDIF 132 144 ! … … 136 148 IF( kt == nit000 ) THEN !== 1st time step ==! 137 149 IF( ln_rstart .AND. iom_varid( numror, 'qsr_hc_b', ldstop = .FALSE. ) > 0 .AND. .NOT.l_1st_euler ) THEN ! read in restart 138 IF(lwp) WRITE(numout,*) ' nit000-1 qsr tracer content forcing field read in the restart file'139 150 z1_2 = 0.5_wp 140 CALL iom_get( numror, jpdom_auto, 'qsr_hc_b', qsr_hc_b, ldxios = lrxios ) ! before heat content trend due to Qsr flux 151 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 152 IF(lwp) WRITE(numout,*) ' nit000-1 qsr tracer content forcing field read in the restart file' 153 CALL iom_get( numror, jpdom_auto, 'qsr_hc_b', qsr_hc_b, ldxios = lrxios ) ! before heat content trend due to Qsr flux 154 ENDIF 141 155 ELSE ! No restart or restart not found: Euler forward time stepping 142 156 z1_2 = 1._wp 143 qsr_hc_b(:,:,:) = 0._wp 157 DO_3D_00_00( 1, jpk ) 158 qsr_hc_b(ji,jj,jk) = 0._wp 159 END_3D 144 160 ENDIF 145 161 ELSE !== Swap of qsr heat content ==! 146 162 z1_2 = 0.5_wp 147 qsr_hc_b(:,:,:) = qsr_hc(:,:,:) 163 DO_3D_00_00( 1, jpk ) 164 qsr_hc_b(ji,jj,jk) = qsr_hc(ji,jj,jk) 165 END_3D 148 166 ENDIF 149 167 ! … … 154 172 CASE( np_BIO ) !== bio-model fluxes ==! 155 173 ! 156 DO jk = 1, nksr157 qsr_hc( :,:,jk) = r1_rho0_rcp * ( etot3(:,:,jk) - etot3(:,:,jk+1) )158 END DO174 DO_3D_00_00( 1, nksr ) 175 qsr_hc(ji,jj,jk) = r1_rho0_rcp * ( etot3(ji,jj,jk) - etot3(ji,jj,jk+1) ) 176 END_3D 159 177 ! 160 178 CASE( np_RGB , np_RGBc ) !== R-G-B fluxes ==! 161 179 ! 162 ALLOCATE( zekb( jpi,jpj) , zekg(jpi,jpj) , zekr (jpi,jpj) , &163 & ze0 ( jpi,jpj,jpk) , ze1 (jpi,jpj,jpk) , ze2 (jpi,jpj,jpk) , &164 & ze3 ( jpi,jpj,jpk) , zea (jpi,jpj,jpk) , zchl3d(jpi,jpj,jpk) )180 ALLOCATE( zekb(A2D) , zekg(A2D) , zekr (A2D) , & 181 & ze0 (A2D,jpk) , ze1 (A2D,jpk) , ze2 (A2D,jpk) , & 182 & ze3 (A2D,jpk) , zea (A2D,jpk) , zchl3d(A2D,jpk) ) 165 183 ! 166 184 IF( nqsr == np_RGBc ) THEN !* Variable Chlorophyll 167 CALL fld_read( kt, 1, sf_chl ) ! Read Chl data and provides it at the current time step 168 DO jk = 1, nksr + 1 169 DO jj = 2, jpjm1 ! Separation in R-G-B depending of the surface Chl 170 DO ji = 2, jpim1 171 zchl = MIN( 10. , MAX( 0.03, sf_chl(1)%fnow(ji,jj,1) ) ) 172 zCtot = 40.6 * zchl**0.459 173 zze = 568.2 * zCtot**(-0.746) 174 IF( zze > 102. ) zze = 200.0 * zCtot**(-0.293) 175 zpsi = gdepw(ji,jj,jk,Kmm) / zze 176 ! 177 zlogc = LOG( zchl ) 178 zlogc2 = zlogc * zlogc 179 zlogc3 = zlogc * zlogc * zlogc 180 zCb = 0.768 + 0.087 * zlogc - 0.179 * zlogc2 - 0.025 * zlogc3 181 zCmax = 0.299 - 0.289 * zlogc + 0.579 * zlogc2 182 zpsimax = 0.6 - 0.640 * zlogc + 0.021 * zlogc2 + 0.115 * zlogc3 183 zdelpsi = 0.710 + 0.159 * zlogc + 0.021 * zlogc2 184 zCze = 1.12 * (zchl)**0.803 185 ! 186 zchl3d(ji,jj,jk) = zCze * ( zCb + zCmax * EXP( -( (zpsi - zpsimax) / zdelpsi )**2 ) ) 187 END DO 188 ! 189 END DO 190 END DO 185 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only for the full domain 186 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 187 CALL fld_read( kt, 1, sf_chl ) ! Read Chl data and provides it at the current time step 188 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 1 ) ! Revert to tile domain 189 ENDIF 190 DO_3D_00_00( 1, nksr + 1 ) ! Separation in R-G-B depending of the surface Chl 191 zchl = MIN( 10. , MAX( 0.03, sf_chl(1)%fnow(ji,jj,1) ) ) 192 zCtot = 40.6 * zchl**0.459 193 zze = 568.2 * zCtot**(-0.746) 194 IF( zze > 102. ) zze = 200.0 * zCtot**(-0.293) 195 zpsi = gdepw(ji,jj,jk,Kmm) / zze 196 ! 197 zlogc = LOG( zchl ) 198 zlogc2 = zlogc * zlogc 199 zlogc3 = zlogc * zlogc * zlogc 200 zCb = 0.768 + 0.087 * zlogc - 0.179 * zlogc2 - 0.025 * zlogc3 201 zCmax = 0.299 - 0.289 * zlogc + 0.579 * zlogc2 202 zpsimax = 0.6 - 0.640 * zlogc + 0.021 * zlogc2 + 0.115 * zlogc3 203 zdelpsi = 0.710 + 0.159 * zlogc + 0.021 * zlogc2 204 zCze = 1.12 * (zchl)**0.803 205 ! 206 zchl3d(ji,jj,jk) = zCze * ( zCb + zCmax * EXP( -( (zpsi - zpsimax) / zdelpsi )**2 ) ) 207 END_3D 191 208 ELSE !* constant chrlorophyll 192 DO jk = 1, nksr + 1193 zchl3d(:,:,jk) = 0.05194 END DO209 DO_3D_11_11( 1, nksr + 1 ) 210 zchl3d(ji,jj,jk) = 0.05 211 END_3D 195 212 ENDIF 196 213 ! … … 256 273 ENDIF 257 274 END_2D 258 CALL lbc_lnk( 'traqsr', fraqsr_1lev(:,:), 'T', 1._wp ) 259 ! 260 IF( iom_use('qsr3d') ) THEN ! output the shortwave Radiation distribution 261 ALLOCATE( zetot(jpi,jpj,jpk) ) 262 zetot(:,:,nksr+1:jpk) = 0._wp ! below ~400m set to zero 263 DO jk = nksr, 1, -1 264 zetot(:,:,jk) = zetot(:,:,jk+1) + qsr_hc(:,:,jk) * rho0_rcp 265 END DO 266 CALL iom_put( 'qsr3d', zetot ) ! 3D distribution of shortwave Radiation 267 DEALLOCATE( zetot ) 268 ENDIF 269 ! 270 IF( lrst_oce ) THEN ! write in the ocean restart file 271 IF( lwxios ) CALL iom_swap( cwxios_context ) 272 CALL iom_rstput( kt, nitrst, numrow, 'qsr_hc_b' , qsr_hc , ldxios = lwxios ) 273 CALL iom_rstput( kt, nitrst, numrow, 'fraqsr_1lev', fraqsr_1lev, ldxios = lwxios ) 274 IF( lwxios ) CALL iom_swap( cxios_context ) 275 ENDIF 276 ! 275 ! TEMP: This change not necessary after extra haloes development (lbc_lnk removed) 276 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only on the last tile 277 CALL lbc_lnk( 'traqsr', fraqsr_1lev(:,:), 'T', 1._wp ) 278 ENDIF 279 ! 280 ! TEMP: This change not necessary and working array can use A2D if using XIOS (subdomain support) 281 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 282 IF( iom_use('qsr3d') ) THEN ! output the shortwave Radiation distribution 283 ALLOCATE( zetot(jpi,jpj,jpk) ) 284 zetot(:,:,nksr+1:jpk) = 0._wp ! below ~400m set to zero 285 DO jk = nksr, 1, -1 286 zetot(:,:,jk) = zetot(:,:,jk+1) + qsr_hc(:,:,jk) * rho0_rcp 287 END DO 288 CALL iom_put( 'qsr3d', zetot ) ! 3D distribution of shortwave Radiation 289 DEALLOCATE( zetot ) 290 ENDIF 291 ENDIF 292 ! 293 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only on the last tile 294 IF( lrst_oce ) THEN ! write in the ocean restart file 295 IF( lwxios ) CALL iom_swap( cwxios_context ) 296 CALL iom_rstput( kt, nitrst, numrow, 'qsr_hc_b' , qsr_hc , ldxios = lwxios ) 297 CALL iom_rstput( kt, nitrst, numrow, 'fraqsr_1lev', fraqsr_1lev, ldxios = lwxios ) 298 IF( lwxios ) CALL iom_swap( cxios_context ) 299 ENDIF 300 ENDIF 301 ! 302 ! TEMP: These changes not necessary after trd_tra is tiled 277 303 IF( l_trdtra ) THEN ! qsr tracers trends saved for diagnostics 278 ztrdt(:,:,:) = pts(:,:,:,jp_tem,Krhs) - ztrdt(:,:,:) 279 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_qsr, ztrdt ) 280 DEALLOCATE( ztrdt ) 304 DO_3D_00_00( 1, jpkm1 ) 305 ztrdt(ji,jj,jk) = pts(ji,jj,jk,jp_tem,Krhs) - ztrdt(ji,jj,jk) 306 END_3D 307 308 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only for the full domain 309 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) ! Use full domain 310 311 ! TODO: TO BE TILED- trd_tra 312 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_qsr, ztrdt ) 313 DEALLOCATE( ztrdt ) 314 315 IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = nijtile ) ! Revert to tile domain 316 ENDIF 281 317 ENDIF 282 318 ! ! print mean trends (used for debugging) 283 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Krhs), clinfo1=' qsr - Ta: ', mask1=tmask, clinfo3='tra-ta' ) 319 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Krhs), clinfo1=' qsr - Ta: ', mask1=tmask, psum1=zsum1, & 320 & clinfo3='tra-ta' ) 284 321 ! 285 322 IF( ln_timing ) CALL timing_stop('tra_qsr') -
NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/TRA/trasbc.F90
r12738 r13409 19 19 USE sbc_oce ! surface boundary condition: ocean 20 20 USE dom_oce ! ocean space domain variables 21 ! TEMP: This change not necessary after trd_tra is tiled 22 USE domain, ONLY : dom_tile 21 23 USE phycst ! physical constant 22 24 USE eosbn2 ! Equation Of State … … 78 80 INTEGER :: ikt, ikb ! local integers 79 81 REAL(wp) :: zfact, z1_e3t, zdep, ztim ! local scalar 80 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdt, ztrds 82 REAL(wp), SAVE :: zsum1, zsum2 83 ! TEMP: This change not necessary after trd_tra is tiled 84 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ztrdt, ztrds 81 85 !!---------------------------------------------------------------------- 82 86 ! 83 87 IF( ln_timing ) CALL timing_start('tra_sbc') 84 88 ! 85 IF( kt == nit000 ) THEN 86 IF(lwp) WRITE(numout,*) 87 IF(lwp) WRITE(numout,*) 'tra_sbc : TRAcer Surface Boundary Condition' 88 IF(lwp) WRITE(numout,*) '~~~~~~~ ' 89 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 90 IF( kt == nit000 ) THEN 91 IF(lwp) WRITE(numout,*) 92 IF(lwp) WRITE(numout,*) 'tra_sbc : TRAcer Surface Boundary Condition' 93 IF(lwp) WRITE(numout,*) '~~~~~~~ ' 94 ENDIF 89 95 ENDIF 90 96 ! 91 97 IF( l_trdtra ) THEN !* Save ta and sa trends 92 ALLOCATE( ztrdt(jpi,jpj,jpk) , ztrds(jpi,jpj,jpk) ) 93 ztrdt(:,:,:) = pts(:,:,:,jp_tem,Krhs) 94 ztrds(:,:,:) = pts(:,:,:,jp_sal,Krhs) 98 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 99 ! TEMP: This can be A2D after trd_tra is tiled 100 ALLOCATE( ztrdt(jpi,jpj,jpk), ztrds(jpi,jpj,jpk) ) 101 ENDIF 102 103 DO_3D_00_00( 1, jpkm1 ) 104 ztrdt(ji,jj,jk) = pts(ji,jj,jk,jp_tem,Krhs) 105 ztrds(ji,jj,jk) = pts(ji,jj,jk,jp_sal,Krhs) 106 END_3D 95 107 ENDIF 96 108 ! 97 109 !!gm This should be moved into sbcmod.F90 module ? (especially now that ln_traqsr is read in namsbc namelist) 98 110 IF( .NOT.ln_traqsr ) THEN ! no solar radiation penetration 99 qns(:,:) = qns(:,:) + qsr(:,:) ! total heat flux in qns 100 qsr(:,:) = 0._wp ! qsr set to zero 111 DO_2D_00_00 112 qns(ji,jj) = qns(ji,jj) + qsr(ji,jj) ! total heat flux in qns 113 qsr(ji,jj) = 0._wp ! qsr set to zero 114 END_2D 101 115 ENDIF 102 116 … … 108 122 IF( ln_rstart .AND. & ! Restart: read in restart file 109 123 & iom_varid( numror, 'sbc_hc_b', ldstop = .FALSE. ) > 0 ) THEN 110 IF(lwp) WRITE(numout,*) ' nit000-1 sbc tracer content field read in the restart file'111 124 zfact = 0.5_wp 112 sbc_tsc(:,:,:) = 0._wp 113 CALL iom_get( numror, jpdom_auto, 'sbc_hc_b', sbc_tsc_b(:,:,jp_tem), ldxios = lrxios ) ! before heat content sbc trend 114 CALL iom_get( numror, jpdom_auto, 'sbc_sc_b', sbc_tsc_b(:,:,jp_sal), ldxios = lrxios ) ! before salt content sbc trend 125 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 126 IF(lwp) WRITE(numout,*) ' nit000-1 sbc tracer content field read in the restart file' 127 sbc_tsc(:,:,:) = 0._wp 128 CALL iom_get( numror, jpdom_auto, 'sbc_hc_b', sbc_tsc_b(:,:,jp_tem), ldxios = lrxios ) ! before heat content sbc trend 129 CALL iom_get( numror, jpdom_auto, 'sbc_sc_b', sbc_tsc_b(:,:,jp_sal), ldxios = lrxios ) ! before salt content sbc trend 130 ENDIF 115 131 ELSE ! No restart or restart not found: Euler forward time stepping 116 132 zfact = 1._wp 117 sbc_tsc(:,:,:) = 0._wp 118 sbc_tsc_b(:,:,:) = 0._wp 133 DO_2D_00_00 134 sbc_tsc(ji,jj,:) = 0._wp 135 sbc_tsc_b(ji,jj,:) = 0._wp 136 END_2D 119 137 ENDIF 120 138 ELSE !* other time-steps: swap of forcing fields 121 139 zfact = 0.5_wp 122 sbc_tsc_b(:,:,:) = sbc_tsc(:,:,:) 140 &