Changeset 5705 for branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP
- Timestamp:
- 2015-08-25T14:17:45+02:00 (9 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/fldread.F90
r5700 r5705 106 106 INTEGER :: nxt_wgt = 1 ! point to next available space in ref_wgts array 107 107 REAL(wp), PARAMETER :: undeff_lsm = -999.00_wp 108 # include "domzgr_substitute.h90" 108 109 109 110 !$AGRIF_END_DO_NOT_TREAT … … 839 840 USE bdy_oce, ONLY: idx_bdy ! indexing for map <-> ij transformation 840 841 841 REAL(wp), POINTER, DIMENSION(:,:,:), INTENT(in ) :: dta_read ! work space for global data842 REAL(wp), POINTER, DIMENSION(:,:,:), INTENT(in ) :: dta_read_z ! work space for global data843 REAL(wp), POINTER, DIMENSION(:,:,:), INTENT(in ) :: dta_read_dz ! work space for global data844 REAL(wp) , DIMENSION(:,:,:), INTENT(out) :: dta ! output field on model grid (2 dimensional)845 TYPE(MAP_POINTER) , INTENT(in ) :: map ! global-to-local mapping indices846 LOGICAL , INTENT(in), OPTIONAL :: fvl ! grid type, set number and number of vertical levels in the bdy data847 INTEGER , INTENT(in) :: igrd, ibdy, jpk_bdy ! number of levels inbdy data848 INTEGER :: jpkm1_bdy ! number of levels in bdy data minus 1849 REAL(wp) , INTENT(in) :: fv ! fillvalue and alternative -ABS(fv)850 !! 851 INTEGER :: ipi ! length of boundary data on local process852 INTEGER :: ipj ! length of dummy dimension ( = 1 )853 INTEGER :: ipk ! number of vertical levels of dta ( 2D: ipk=1 ; 3D: ipk=jpk )854 INTEGER , INTENT(in) :: ilendta ! length of data in file855 INTEGER :: ib, ik, ikk ! loop counters856 INTEGER :: ji, jj ! loop counters857 REAL(wp) :: zl, zi, zh , zz, zdz! tmp variable for current depth and interpolation factor842 REAL(wp), POINTER, DIMENSION(:,:,:), INTENT(in ) :: dta_read ! work space for global data 843 REAL(wp), POINTER, DIMENSION(:,:,:), INTENT(in ) :: dta_read_z ! work space for global data 844 REAL(wp), POINTER, DIMENSION(:,:,:), INTENT(in ) :: dta_read_dz ! work space for global data 845 REAL(wp) , INTENT(in) :: fv ! fillvalue and alternative -ABS(fv) 846 REAL(wp), DIMENSION(:,:,:), INTENT(out) :: dta ! output field on model grid (2 dimensional) 847 TYPE(MAP_POINTER) , INTENT(in ) :: map ! global-to-local mapping indices 848 LOGICAL , INTENT(in), OPTIONAL :: fvl ! grid type, set number and number of vertical levels in the bdy data 849 INTEGER , INTENT(in) :: igrd, ibdy, jpk_bdy ! number of levels in bdy data 850 INTEGER , INTENT(in) :: ilendta ! length of data in file 851 !! 852 INTEGER :: ipi ! length of boundary data on local process 853 INTEGER :: ipj ! length of dummy dimension ( = 1 ) 854 INTEGER :: ipk ! number of vertical levels of dta ( 2D: ipk=1 ; 3D: ipk=jpk ) 855 INTEGER :: jpkm1_bdy ! number of levels in bdy data minus 1 856 INTEGER :: ib, ik, ikk ! loop counters 857 INTEGER :: ji, jj, zij, zjj ! temporary indices 858 REAL(wp) :: zl, zi, zh ! tmp variable for current depth and interpolation factor 858 859 REAL(wp) :: fv_alt, ztrans, ztrans_new ! fillvalue and alternative -ABS(fv) 860 CHARACTER (LEN=10) :: ibstr 859 861 !!--------------------------------------------------------------------- 860 862 … … 868 870 869 871 ! 870 IF ( map%ll_unstruc ) THEN ! unstructured open boundary data file872 IF ( map%ll_unstruc ) THEN ! unstructured open boundary data file 871 873 872 874 DO ib = 1, ipi 873 875 DO ik = 1, jpk_bdy 874 876 IF( ( dta_read(map%ptr(ib),1,ik) == fv ) ) THEN 875 dta_read_z(map%ptr(ib),1,ik) = fv_alt ! safety: put fillvalue into external depth field so consistent with data876 dta_read_dz(map%ptr(ib),1,ik) = 0._wp ! safety: put 0._wp into external thickness factors to ensure transport is correct877 dta_read_z(map%ptr(ib),1,ik) = fv_alt ! safety: put fillvalue into external depth field so consistent with data 878 dta_read_dz(map%ptr(ib),1,ik) = 0._wp ! safety: put 0._wp into external thickness factors to ensure transport is correct 877 879 ENDIF 878 880 ENDDO … … 880 882 881 883 DO ib = 1, ipi 884 zij = idx_bdy(ibdy)%nbi(ib,igrd) 885 zjj = idx_bdy(ibdy)%nbj(ib,igrd) 886 zh = SUM(dta_read_dz(map%ptr(ib),1,:) ) 882 887 DO ik = 1, ipk 883 zl = gdept_0(idx_bdy(ibdy)%nbi(ib,igrd),idx_bdy(ibdy)%nbj(ib,igrd),ik) ! if using in step could use fsdept instead of gdept_0? 888 SELECT CASE( igrd ) 889 CASE(1) 890 zl = gdept_0(zij,zjj,ik) ! if using in step could use fsdept instead of gdept_0? 891 CASE(2) 892 IF(ln_sco) THEN 893 zl = ( gdept_0(zij,zjj,ik) + gdept_0(zij+1,zjj,ik) ) * 0.5_wp ! if using in step could use fsdept instead of gdept_0? 894 ELSE 895 zl = MIN( gdept_0(zij,zjj,ik), gdept_0(zij+1,zjj,ik) ) 896 ENDIF 897 CASE(3) 898 IF(ln_sco) THEN 899 zl = ( gdept_0(zij,zjj,ik) + gdept_0(zij,zjj+1,ik) ) * 0.5_wp ! if using in step could use fsdept instead of gdept_0? 900 ELSE 901 zl = MIN( gdept_0(zij,zjj,ik), gdept_0(zij,zjj+1,ik) ) 902 ENDIF 903 END SELECT 884 904 IF( zl < dta_read_z(map%ptr(ib),1,1) ) THEN ! above the first level of external data 885 905 dta(ib,1,ik) = dta_read(map%ptr(ib),1,1) 886 906 ELSEIF( zl > MAXVAL(dta_read_z(map%ptr(ib),1,:),1) ) THEN ! below the last level of external data 887 907 dta(ib,1,ik) = dta_read(map%ptr(ib),1,MAXLOC(dta_read_z(map%ptr(ib),1,:),1)) 888 ELSE ! inbetween : vertical interpolation between ikk & ikk+1889 DO ikk = 1, jpkm1_bdy ! when gdept_0(ikk) < zl < gdept_0(ikk+1)890 IF( ( (zl-dta_read_z(map%ptr(ib),1,ikk)) * (zl-dta_read_z(map%ptr(ib),1,ikk+1)) <= 0._wp) 908 ELSE ! inbetween : vertical interpolation between ikk & ikk+1 909 DO ikk = 1, jpkm1_bdy ! when gdept_0(ikk) < zl < gdept_0(ikk+1) 910 IF( ( (zl-dta_read_z(map%ptr(ib),1,ikk)) * (zl-dta_read_z(map%ptr(ib),1,ikk+1)) <= 0._wp) & 891 911 & .AND. (dta_read_z(map%ptr(ib),1,ikk+1) /= fv_alt)) THEN 892 zi = ( zl - dta_read_z(map%ptr(ib),1,ikk) ) / (dta_read_z(map%ptr(ib),1,ikk+1)-dta_read_z(map%ptr(ib),1,ikk)) 912 zi = ( zl - dta_read_z(map%ptr(ib),1,ikk) ) / & 913 & ( dta_read_z(map%ptr(ib),1,ikk+1) - dta_read_z(map%ptr(ib),1,ikk) ) 893 914 dta(ib,1,ik) = dta_read(map%ptr(ib),1,ikk) + & 894 & ( dta_read(map%ptr(ib),1,ikk+1) -dta_read(map%ptr(ib),1,ikk) ) * zi915 & ( dta_read(map%ptr(ib),1,ikk+1) - dta_read(map%ptr(ib),1,ikk) ) * zi 895 916 ENDIF 896 917 END DO … … 899 920 END DO 900 921 901 IF(igrd == 2) THEN ! do we need to adjust the transport term?922 IF(igrd == 2) THEN ! do we need to adjust the transport term? 902 923 DO ib = 1, ipi 903 zh = SUM(dta_read_dz(map%ptr(ib),1,:) ) 924 zij = idx_bdy(ibdy)%nbi(ib,igrd) 925 zjj = idx_bdy(ibdy)%nbj(ib,igrd) 926 zh = SUM(dta_read_dz(map%ptr(ib),1,:) ) 904 927 ztrans = 0._wp 905 928 ztrans_new = 0._wp 906 DO ik = 1, jpk_bdy 907 ztrans = ztrans+ dta_read(map%ptr(ib),1,ik) * dta_read_dz(map%ptr(ib),1,ik)929 DO ik = 1, jpk_bdy ! calculate transport on input grid 930 ztrans = ztrans + dta_read(map%ptr(ib),1,ik) * dta_read_dz(map%ptr(ib),1,ik) 908 931 ENDDO 909 DO ik = 1, ipk 910 zdz = e3u_0(idx_bdy(ibdy)%nbi(ib,igrd),idx_bdy(ibdy)%nbj(ib,igrd),ik) 911 ztrans_new = ztrans_new + dta(ib,1,ik) * zdz 932 DO ik = 1, ipk ! calculate transport on model grid 933 ztrans_new = ztrans_new + dta(ib,1,ik) * e3u_0(zij,zjj,ik) * umask(zij,zjj,ik) 912 934 ENDDO 913 DO ik = 1, ipk 914 zdz = e3u_0(idx_bdy(ibdy)%nbi(ib,igrd),idx_bdy(ibdy)%nbj(ib,igrd),ik) 915 zz = hur(idx_bdy(ibdy)%nbi(ib,igrd),idx_bdy(ibdy)%nbj(ib,igrd)) ! Is array HUR the correct one to use???? 935 DO ik = 1, ipk ! make transport correction 916 936 IF(fvl) THEN ! bdy data are total velocity so adjust bt transport term to match input data 917 dta(ib,1,ik) = dta(ib,1,ik) + ( ztrans - ztrans_new ) * ( zdz * zz)937 dta(ib,1,ik) = ( dta(ib,1,ik) + ( ztrans - ztrans_new ) / hu_0(zij,zjj) ) * umask(zij,zjj,ik) 918 938 ELSE ! we're just dealing with bc velocity so bt transport term should sum to zero 919 dta(ib,1,ik) = dta(ib,1,ik) + ( 0._wp - ztrans_new ) * ( zdz * zz ) 939 IF( ABS(ztrans / hu_0(zij,zjj)) > 0.01_wp ) & 940 & CALL ctl_warn('fld_bdy_interp: barotropic component of > 0.01 ms-1 found in baroclinic velocities at') 941 dta(ib,1,ik) = dta(ib,1,ik) + ( 0._wp - ztrans_new ) / hu_0(zij,zjj) * umask(zij,zjj,ik) 920 942 ENDIF 921 943 ENDDO … … 923 945 ENDIF 924 946 925 IF(igrd == 3) THEN ! do we need to adjust the transport term?947 IF(igrd == 3) THEN ! do we need to adjust the transport term? 926 948 DO ib = 1, ipi 927 zh = SUM(dta_read_dz(map%ptr(ib),1,:) ) 949 zij = idx_bdy(ibdy)%nbi(ib,igrd) 950 zjj = idx_bdy(ibdy)%nbj(ib,igrd) 951 zh = SUM(dta_read_dz(map%ptr(ib),1,:) ) 928 952 ztrans = 0._wp 929 953 ztrans_new = 0._wp 930 DO ik = 1, jpk_bdy 931 ztrans = ztrans+ dta_read(map%ptr(ib),1,ik) * dta_read_dz(map%ptr(ib),1,ik)954 DO ik = 1, jpk_bdy ! calculate transport on input grid 955 ztrans = ztrans + dta_read(map%ptr(ib),1,ik) * dta_read_dz(map%ptr(ib),1,ik) 932 956 ENDDO 933 DO ik = 1, ipk 934 zdz = e3v_0(idx_bdy(ibdy)%nbi(ib,igrd),idx_bdy(ibdy)%nbj(ib,igrd),ik) 935 ztrans_new = ztrans_new + dta(ib,1,ik) * zdz 957 DO ik = 1, ipk ! calculate transport on model grid 958 ztrans_new = ztrans_new + dta(ib,1,ik) * e3v_0(zij,zjj,ik) * vmask(zij,zjj,ik) 936 959 ENDDO 937 DO ik = 1, ipk 938 zdz = e3v_0(idx_bdy(ibdy)%nbi(ib,igrd),idx_bdy(ibdy)%nbj(ib,igrd),ik) 939 zz = hvr(idx_bdy(ibdy)%nbi(ib,igrd),idx_bdy(ibdy)%nbj(ib,igrd)) 960 DO ik = 1, ipk ! make transport correction 940 961 IF(fvl) THEN ! bdy data are total velocity so adjust bt transport term to match input data 941 dta(ib,1,ik) = dta(ib,1,ik) + ( ztrans - ztrans_new ) * ( zdz * zz)962 dta(ib,1,ik) = ( dta(ib,1,ik) + ( ztrans - ztrans_new ) / hv_0(zij,zjj) ) * vmask(zij,zjj,ik) 942 963 ELSE ! we're just dealing with bc velocity so bt transport term should sum to zero 943 dta(ib,1,ik) = dta(ib,1,ik) + ( 0._wp - ztrans_new ) * ( zdz * zz)964 dta(ib,1,ik) = dta(ib,1,ik) + ( 0._wp - ztrans_new ) / hv_0(zij,zjj) * vmask(zij,zjj,ik) 944 965 ENDIF 945 966 ENDDO … … 947 968 ENDIF 948 969 949 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!950 ! At this point write out a single velocity profile/dz/H for model and input data !951 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!952 953 970 ELSE ! structured open boundary file 954 971 … … 958 975 DO ik = 1, jpk_bdy 959 976 IF( ( dta_read(ji,jj,ik) == fv ) ) THEN 960 dta_read_z(ji,jj,ik) = fv_alt ! safety: put fillvalue into external depth field so consistent with data961 dta_read_dz(ji,jj,ik) = 0._wp ! safety: put 0._wp into external thickness factors to ensure transport is correct977 dta_read_z(ji,jj,ik) = fv_alt ! safety: put fillvalue into external depth field so consistent with data 978 dta_read_dz(ji,jj,ik) = 0._wp ! safety: put 0._wp into external thickness factors to ensure transport is correct 962 979 ENDIF 963 980 ENDDO … … 968 985 jj=1+floor(REAL(map%ptr(ib)-1)/REAL(ilendta)) 969 986 ji=map%ptr(ib)-(jj-1)*ilendta 987 zij = idx_bdy(ibdy)%nbi(ib,igrd) 988 zjj = idx_bdy(ibdy)%nbj(ib,igrd) 989 zh = SUM(dta_read_dz(ji,jj,:) ) 970 990 DO ik = 1, ipk 971 zl = gdept_0(idx_bdy(ibdy)%nbi(ib,igrd),idx_bdy(ibdy)%nbj(ib,igrd),ik) ! if using in step could use fsdept instead of gdept_0? 972 IF( zl < dta_read_z(ji,jj,1) ) THEN ! above the first level of external data 991 SELECT CASE( igrd ) ! coded for sco - need zco and zps option using min 992 CASE(1) 993 zl = gdept_0(zij,zjj,ik) ! if using in step could use fsdept instead of gdept_0? 994 CASE(2) 995 IF(ln_sco) THEN 996 zl = ( gdept_0(zij,zjj,ik) + gdept_0(zij+1,zjj,ik) ) * 0.5_wp ! if using in step could use fsdept instead of gdept_0? 997 ELSE 998 zl = MIN( gdept_0(zij,zjj,ik), gdept_0(zij+1,zjj,ik) ) 999 ENDIF 1000 CASE(3) 1001 IF(ln_sco) THEN 1002 zl = ( gdept_0(zij,zjj,ik) + gdept_0(zij,zjj+1,ik) ) * 0.5_wp ! if using in step could use fsdept instead of gdept_0? 1003 ELSE 1004 zl = MIN( gdept_0(zij,zjj,ik), gdept_0(zij,zjj+1,ik) ) 1005 ENDIF 1006 END SELECT 1007 IF( zl < dta_read_z(ji,jj,1) ) THEN ! above the first level of external data 973 1008 dta(ib,1,ik) = dta_read(ji,jj,1) 974 ELSEIF( zl > MAXVAL(dta_read_z(ji,jj,:),1) ) THEN 1009 ELSEIF( zl > MAXVAL(dta_read_z(ji,jj,:),1) ) THEN ! below the last level of external data 975 1010 dta(ib,1,ik) = dta_read(ji,jj,MAXLOC(dta_read_z(ji,jj,:),1)) 976 ELSE 977 DO ikk = 1, jpkm1_bdy 978 IF( ( (zl-dta_read_z(ji,jj,ikk)) * (zl-dta_read_z(ji,jj,ikk+1)) <= 0._wp) 1011 ELSE ! inbetween : vertical interpolation between ikk & ikk+1 1012 DO ikk = 1, jpkm1_bdy ! when gdept_0(ikk) < zl < gdept_0(ikk+1) 1013 IF( ( (zl-dta_read_z(ji,jj,ikk)) * (zl-dta_read_z(ji,jj,ikk+1)) <= 0._wp) & 979 1014 & .AND. (dta_read_z(ji,jj,ikk+1) /= fv_alt)) THEN 980 zi = ( zl - dta_read_z(ji,jj,ikk) ) / (dta_read_z(ji,jj,ikk+1)-dta_read_z(ji,jj,ikk)) 1015 zi = ( zl - dta_read_z(ji,jj,ikk) ) / & 1016 & ( dta_read_z(ji,jj,ikk+1) - dta_read_z(ji,jj,ikk) ) 981 1017 dta(ib,1,ik) = dta_read(ji,jj,ikk) + & 982 & ( dta_read(ji,jj,ikk+1) -dta_read(ji,jj,ikk) ) * zi1018 & ( dta_read(ji,jj,ikk+1) - dta_read(ji,jj,ikk) ) * zi 983 1019 ENDIF 984 1020 END DO … … 987 1023 END DO 988 1024 989 IF(igrd == 2) THEN ! do we need to adjust the transport term?1025 IF(igrd == 2) THEN ! do we need to adjust the transport term? 990 1026 jj=1+floor(REAL(map%ptr(ib)-1)/REAL(ilendta)) 991 1027 ji=map%ptr(ib)-(jj-1)*ilendta 1028 zij = idx_bdy(ibdy)%nbi(ib,igrd) 1029 zjj = idx_bdy(ibdy)%nbj(ib,igrd) 992 1030 DO ib = 1, ipi 993 1031 zh = SUM(dta_read_dz(ji,jj,:) ) 994 1032 ztrans = 0._wp 995 1033 ztrans_new = 0._wp 996 DO ik = 1, jpk_bdy 1034 DO ik = 1, jpk_bdy ! calculate transport on input grid 997 1035 ztrans = ztrans + dta_read(ji,jj,ik) * dta_read_dz(ji,jj,ik) 998 1036 ENDDO 999 DO ik = 1, ipk 1000 zdz = e3u_0(idx_bdy(ibdy)%nbi(ib,igrd),idx_bdy(ibdy)%nbj(ib,igrd),ik) 1001 ztrans_new = ztrans_new + dta(ib,1,ik) * zdz 1037 DO ik = 1, ipk ! calculate transport on model grid 1038 ztrans_new = ztrans_new + dta(ib,1,ik) * e3u_0(zij,zjj,ik) * umask(zij,zjj,ik) 1002 1039 ENDDO 1003 DO ik = 1, ipk 1004 zdz = e3u_0(idx_bdy(ibdy)%nbi(ib,igrd),idx_bdy(ibdy)%nbj(ib,igrd),ik) 1005 zz = hur(idx_bdy(ibdy)%nbi(ib,igrd),idx_bdy(ibdy)%nbj(ib,igrd)) 1040 DO ik = 1, ipk ! make transport correction 1006 1041 IF(fvl) THEN ! bdy data are total velocity so adjust bt transport term to match input data 1007 dta(ib,1,ik) = dta(ib,1,ik) + ( ztrans - ztrans_new ) * ( zdz * zz)1042 dta(ib,1,ik) = ( dta(ib,1,ik) + ( ztrans - ztrans_new ) * hur(zij,zjj) ) * umask(zij,zjj,ik) 1008 1043 ELSE ! we're just dealing with bc velocity so bt transport term should sum to zero 1009 dta(ib,1,ik) = dta(ib,1,ik) + ( 0._wp - ztrans_new ) * ( zdz * zz)1044 dta(ib,1,ik) = ( dta(ib,1,ik) + ( 0._wp - ztrans_new ) * hur(zij,zjj) ) * umask(zij,zjj,ik) 1010 1045 ENDIF 1011 1046 ENDDO … … 1013 1048 ENDIF 1014 1049 1015 IF(igrd == 3) THEN ! do we need to adjust the transport term? 1016 jj=1+floor(REAL(map%ptr(ib)-1)/REAL(ilendta)) 1017 ji=map%ptr(ib)-(jj-1)*ilendta 1050 IF(igrd == 3) THEN ! do we need to adjust the transport term? 1018 1051 DO ib = 1, ipi 1019 zh = SUM(dta_read_dz(ji,jj,:) ) 1052 jj = 1+floor(REAL(map%ptr(ib)-1)/REAL(ilendta)) 1053 ji = map%ptr(ib)-(jj-1)*ilendta 1054 zij = idx_bdy(ibdy)%nbi(ib,igrd) 1055 zjj = idx_bdy(ibdy)%nbj(ib,igrd) 1056 zh = SUM(dta_read_dz(ji,jj,:) ) 1020 1057 ztrans = 0._wp 1021 1058 ztrans_new = 0._wp 1022 DO ik = 1, jpk_bdy 1023 ztrans = ztrans+ dta_read(ji,jj,ik) * dta_read_dz(ji,jj,ik)1059 DO ik = 1, jpk_bdy ! calculate transport on input grid 1060 ztrans = ztrans + dta_read(ji,jj,ik) * dta_read_dz(ji,jj,ik) 1024 1061 ENDDO 1025 DO ik = 1, ipk 1026 zdz = e3v_0(idx_bdy(ibdy)%nbi(ib,igrd),idx_bdy(ibdy)%nbj(ib,igrd),ik) 1027 ztrans_new = ztrans_new + dta(ib,1,ik) * zdz 1062 DO ik = 1, ipk ! calculate transport on model grid 1063 ztrans_new = ztrans_new + dta(ib,1,ik) * e3v_0(zij,zjj,ik) * vmask(zij,zjj,ik) 1028 1064 ENDDO 1029 DO ik = 1, ipk 1030 zdz = e3v_0(idx_bdy(ibdy)%nbi(ib,igrd),idx_bdy(ibdy)%nbj(ib,igrd),ik) 1031 zz = hvr(idx_bdy(ibdy)%nbi(ib,igrd),idx_bdy(ibdy)%nbj(ib,igrd)) 1065 DO ik = 1, ipk ! make transport correction 1032 1066 IF(fvl) THEN ! bdy data are total velocity so adjust bt transport term to match input data 1033 dta(ib,1,ik) = dta(ib,1,ik) + ( ztrans - ztrans_new ) * ( zdz * zz)1067 dta(ib,1,ik) = ( dta(ib,1,ik) + ( ztrans - ztrans_new ) * hvr(zij,zjj) ) * vmask(zij,zjj,ik) 1034 1068 ELSE ! we're just dealing with bc velocity so bt transport term should sum to zero 1035 dta(ib,1,ik) = dta(ib,1,ik) + ( 0._wp - ztrans_new ) * ( zdz * zz)1069 dta(ib,1,ik) = ( dta(ib,1,ik) + ( 0._wp - ztrans_new ) * hvr(zij,zjj) ) * vmask(zij,zjj,ik) 1036 1070 ENDIF 1037 1071 ENDDO
Note: See TracChangeset
for help on using the changeset viewer.