Changeset 636 for trunk/NEMO/NST_SRC/agrif_opa_interp.F90
- Timestamp:
- 2007-03-07T14:28:16+01:00 (17 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/NST_SRC/agrif_opa_interp.F90
r469 r636 1 ! 2 Module agrif_opa_interp 1 MODULE agrif_opa_interp 3 2 #if defined key_agrif 4 USE par_oce 5 USE oce 6 USE dom_oce 7 USE sol_oce 8 9 CONTAINS 10 SUBROUTINE Agrif_tra( kt ) 11 12 Implicit none 13 14 !! * Substitutions 3 USE par_oce 4 USE oce 5 USE dom_oce 6 USE sol_oce 7 8 IMPLICIT NONE 9 PRIVATE 10 11 PUBLIC Agrif_tra, Agrif_dyn, interpu, interpv 12 13 CONTAINS 14 15 SUBROUTINE Agrif_tra( kt ) 16 !!--------------------------------------------- 17 !! *** ROUTINE Agrif_Tra *** 18 !!--------------------------------------------- 15 19 # include "domzgr_substitute.h90" 16 20 # include "vectopt_loop_substitute.h90" 17 ! 18 INTEGER :: kt19 REAL(wp) tatemp(jpi,jpj,jpk) , satemp(jpi,jpj,jpk) 21 22 INTEGER, INTENT(in) :: kt 23 20 24 INTEGER :: ji,jj,jk 21 REAL(wp) :: rhox25 REAL(wp) :: zrhox 22 26 REAL(wp) :: alpha1, alpha2, alpha3, alpha4 23 27 REAL(wp) :: alpha5, alpha6, alpha7 24 ! 25 IF (Agrif_Root()) RETURN 26 27 Agrif_SpecialValue=0. 28 Agrif_UseSpecialValue = .TRUE. 29 tatemp = 0. 30 satemp = 0. 31 32 Call Agrif_Bc_variable(tatemp,tn) 33 Call Agrif_Bc_variable(satemp,sn) 34 Agrif_UseSpecialValue = .FALSE. 35 36 rhox = Agrif_Rhox() 37 38 alpha1 = (rhox-1.)/2. 39 alpha2 = 1.-alpha1 40 41 alpha3 = (rhox-1)/(rhox+1) 42 alpha4 = 1.-alpha3 43 44 alpha6 = 2.*(rhox-1.)/(rhox+1.) 45 alpha7 = -(rhox-1)/(rhox+3) 46 alpha5 = 1. - alpha6 - alpha7 47 48 ! 49 If ((nbondi == 1).OR.(nbondi == 2)) THEN 28 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zta, zsa 29 ! 30 IF(Agrif_Root()) RETURN 31 32 Agrif_SpecialValue=0. 33 Agrif_UseSpecialValue = .TRUE. 34 zta = 0.e0 35 zsa = 0.e0 36 37 CALL Agrif_Bc_variable(zta,tn) 38 CALL Agrif_Bc_variable(zsa,sn) 39 Agrif_UseSpecialValue = .FALSE. 40 41 zrhox = Agrif_Rhox() 42 43 alpha1 = (zrhox-1.)/2. 44 alpha2 = 1.-alpha1 45 46 alpha3 = (zrhox-1)/(zrhox+1) 47 alpha4 = 1.-alpha3 48 49 alpha6 = 2.*(zrhox-1.)/(zrhox+1.) 50 alpha7 = -(zrhox-1)/(zrhox+3) 51 alpha5 = 1. - alpha6 - alpha7 52 53 IF((nbondi == 1).OR.(nbondi == 2)) THEN 54 55 ta(nlci,:,:) = alpha1 * zta(nlci,:,:) + alpha2 * zta(nlci-1,:,:) 56 sa(nlci,:,:) = alpha1 * zsa(nlci,:,:) + alpha2 * zsa(nlci-1,:,:) 57 58 DO jk=1,jpk 59 DO jj=1,jpj 60 IF (umask(nlci-2,jj,jk).EQ.0.) THEN 61 ta(nlci-1,jj,jk) = ta(nlci,jj,jk) * tmask(nlci-1,jj,jk) 62 sa(nlci-1,jj,jk) = sa(nlci,jj,jk) * tmask(nlci-1,jj,jk) 63 ELSE 64 ta(nlci-1,jj,jk)=(alpha4*ta(nlci,jj,jk)+alpha3*ta(nlci-2,jj,jk))*tmask(nlci-1,jj,jk) 65 sa(nlci-1,jj,jk)=(alpha4*sa(nlci,jj,jk)+alpha3*sa(nlci-2,jj,jk))*tmask(nlci-1,jj,jk) 66 IF (un(nlci-2,jj,jk).GT.0.) THEN 67 ta(nlci-1,jj,jk)=( alpha6*ta(nlci-2,jj,jk)+alpha5*ta(nlci,jj,jk) & 68 + alpha7*ta(nlci-3,jj,jk) ) * tmask(nlci-1,jj,jk) 69 sa(nlci-1,jj,jk)=( alpha6*sa(nlci-2,jj,jk)+alpha5*sa(nlci,jj,jk) & 70 + alpha7*sa(nlci-3,jj,jk) ) * tmask(nlci-1,jj,jk) 71 ENDIF 72 ENDIF 73 END DO 74 END DO 75 ENDIF 76 77 IF((nbondj == 1).OR.(nbondj == 2)) THEN 78 79 ta(:,nlcj,:) = alpha1 * zta(:,nlcj,:) + alpha2 * zta(:,nlcj-1,:) 80 sa(:,nlcj,:) = alpha1 * zsa(:,nlcj,:) + alpha2 * zsa(:,nlcj-1,:) 81 82 DO jk=1,jpk 83 DO ji=1,jpi 84 IF (vmask(ji,nlcj-2,jk).EQ.0.) THEN 85 ta(ji,nlcj-1,jk) = ta(ji,nlcj,jk) * tmask(ji,nlcj-1,jk) 86 sa(ji,nlcj-1,jk) = sa(ji,nlcj,jk) * tmask(ji,nlcj-1,jk) 87 ELSE 88 ta(ji,nlcj-1,jk)=(alpha4*ta(ji,nlcj,jk)+alpha3*ta(ji,nlcj-2,jk))*tmask(ji,nlcj-1,jk) 89 sa(ji,nlcj-1,jk)=(alpha4*sa(ji,nlcj,jk)+alpha3*sa(ji,nlcj-2,jk))*tmask(ji,nlcj-1,jk) 90 IF (vn(ji,nlcj-2,jk) .GT. 0.) THEN 91 ta(ji,nlcj-1,jk)=( alpha6*ta(ji,nlcj-2,jk)+alpha5*ta(ji,nlcj,jk) & 92 + alpha7*ta(ji,nlcj-3,jk) ) * tmask(ji,nlcj-1,jk) 93 sa(ji,nlcj-1,jk)=( alpha6*sa(ji,nlcj-2,jk)+alpha5*sa(ji,nlcj,jk) & 94 + alpha7*sa(ji,nlcj-3,jk))*tmask(ji,nlcj-1,jk) 95 ENDIF 96 ENDIF 97 END DO 98 END DO 99 ENDIF 100 101 IF((nbondi == -1).OR.(nbondi == 2)) THEN 102 ta(1,:,:) = alpha1 * zta(1,:,:) + alpha2 * zta(2,:,:) 103 sa(1,:,:) = alpha1 * zsa(1,:,:) + alpha2 * zsa(2,:,:) 104 DO jk=1,jpk 105 DO jj=1,jpj 106 IF (umask(2,jj,jk).EQ.0.) THEN 107 ta(2,jj,jk) = ta(1,jj,jk) * tmask(2,jj,jk) 108 sa(2,jj,jk) = sa(1,jj,jk) * tmask(2,jj,jk) 109 ELSE 110 ta(2,jj,jk)=(alpha4*ta(1,jj,jk)+alpha3*ta(3,jj,jk))*tmask(2,jj,jk) 111 sa(2,jj,jk)=(alpha4*sa(1,jj,jk)+alpha3*sa(3,jj,jk))*tmask(2,jj,jk) 112 IF (un(2,jj,jk).LT.0.) THEN 113 ta(2,jj,jk)=(alpha6*ta(3,jj,jk)+alpha5*ta(1,jj,jk)+alpha7*ta(4,jj,jk))*tmask(2,jj,jk) 114 sa(2,jj,jk)=(alpha6*sa(3,jj,jk)+alpha5*sa(1,jj,jk)+alpha7*sa(4,jj,jk))*tmask(2,jj,jk) 115 ENDIF 116 ENDIF 117 END DO 118 END DO 119 ENDIF 120 121 IF((nbondj == -1).OR.(nbondj == 2)) THEN 122 ta(:,1,:) = alpha1 * zta(:,1,:) + alpha2 * zta(:,2,:) 123 sa(:,1,:) = alpha1 * zsa(:,1,:) + alpha2 * zsa(:,2,:) 124 DO jk=1,jpk 125 DO ji=1,jpi 126 IF (vmask(ji,2,jk).EQ.0.) THEN 127 ta(ji,2,jk)=ta(ji,1,jk) * tmask(ji,2,jk) 128 sa(ji,2,jk)=sa(ji,1,jk) * tmask(ji,2,jk) 129 ELSE 130 ta(ji,2,jk)=(alpha4*ta(ji,1,jk)+alpha3*ta(ji,3,jk))*tmask(ji,2,jk) 131 sa(ji,2,jk)=(alpha4*sa(ji,1,jk)+alpha3*sa(ji,3,jk))*tmask(ji,2,jk) 132 IF (vn(ji,2,jk) .LT. 0.) THEN 133 ta(ji,2,jk)=(alpha6*ta(ji,3,jk)+alpha5*ta(ji,1,jk)+alpha7*ta(ji,4,jk))*tmask(ji,2,jk) 134 sa(ji,2,jk)=(alpha6*sa(ji,3,jk)+alpha5*sa(ji,1,jk)+alpha7*sa(ji,4,jk))*tmask(ji,2,jk) 135 ENDIF 136 ENDIF 137 END DO 138 END DO 139 ENDIF 140 141 END SUBROUTINE Agrif_tra 142 143 SUBROUTINE Agrif_dyn( kt ) 144 !!--------------------------------------------- 145 !! *** ROUTINE Agrif_DYN *** 146 !!--------------------------------------------- 147 USE phycst 148 USE in_out_manager 149 150 # include "domzgr_substitute.h90" 50 151 51 ta(nlci,:,:) = alpha1 * tatemp(nlci,:,:) + alpha2 * tatemp(nlci-1,:,:) 52 sa(nlci,:,:) = alpha1 * satemp(nlci,:,:) + alpha2 * satemp(nlci-1,:,:) 53 54 Do jk=1,jpk 55 Do jj=1,jpj 56 IF (umask(nlci-2,jj,jk).EQ.0.) THEN 57 ta(nlci-1,jj,jk) = ta(nlci,jj,jk) * tmask(nlci-1,jj,jk) 58 sa(nlci-1,jj,jk) = sa(nlci,jj,jk) * tmask(nlci-1,jj,jk) 59 ELSE 60 ta(nlci-1,jj,jk)=(alpha4*ta(nlci,jj,jk)+alpha3*ta(nlci-2,jj,jk))*tmask(nlci-1,jj,jk) 61 sa(nlci-1,jj,jk)=(alpha4*sa(nlci,jj,jk)+alpha3*sa(nlci-2,jj,jk))*tmask(nlci-1,jj,jk) 62 IF (un(nlci-2,jj,jk).GT.0.) THEN 63 ta(nlci-1,jj,jk)=(alpha6*ta(nlci-2,jj,jk)+alpha5*ta(nlci,jj,jk)+alpha7*ta(nlci-3,jj,jk))*tmask(nlci-1,jj,jk) 64 sa(nlci-1,jj,jk)=(alpha6*sa(nlci-2,jj,jk)+alpha5*sa(nlci,jj,jk)+alpha7*sa(nlci-3,jj,jk))*tmask(nlci-1,jj,jk) 65 ENDIF 66 ENDIF 67 End Do 68 enddo 69 ENDIF 70 71 If ((nbondj == 1).OR.(nbondj == 2)) THEN 72 73 ta(:,nlcj,:) = alpha1 * tatemp(:,nlcj,:) + alpha2 * tatemp(:,nlcj-1,:) 74 sa(:,nlcj,:) = alpha1 * satemp(:,nlcj,:) + alpha2 * satemp(:,nlcj-1,:) 75 76 Do jk=1,jpk 77 Do ji=1,jpi 78 IF (vmask(ji,nlcj-2,jk).EQ.0.) THEN 79 ta(ji,nlcj-1,jk) = ta(ji,nlcj,jk) * tmask(ji,nlcj-1,jk) 80 sa(ji,nlcj-1,jk) = sa(ji,nlcj,jk) * tmask(ji,nlcj-1,jk) 81 ELSE 82 ta(ji,nlcj-1,jk)=(alpha4*ta(ji,nlcj,jk)+alpha3*ta(ji,nlcj-2,jk))*tmask(ji,nlcj-1,jk) 83 sa(ji,nlcj-1,jk)=(alpha4*sa(ji,nlcj,jk)+alpha3*sa(ji,nlcj-2,jk))*tmask(ji,nlcj-1,jk) 84 IF (vn(ji,nlcj-2,jk) .GT. 0.) THEN 85 ta(ji,nlcj-1,jk)=(alpha6*ta(ji,nlcj-2,jk)+alpha5*ta(ji,nlcj,jk)+alpha7*ta(ji,nlcj-3,jk))*tmask(ji,nlcj-1,jk) 86 sa(ji,nlcj-1,jk)=(alpha6*sa(ji,nlcj-2,jk)+alpha5*sa(ji,nlcj,jk)+alpha7*sa(ji,nlcj-3,jk))*tmask(ji,nlcj-1,jk) 87 ENDIF 88 ENDIF 89 End Do 90 enddo 91 ENDIF 92 93 IF ((nbondi == -1).OR.(nbondi == 2)) THEN 94 95 ta(1,:,:) = alpha1 * tatemp(1,:,:) + alpha2 * tatemp(2,:,:) 96 sa(1,:,:) = alpha1 * satemp(1,:,:) + alpha2 * satemp(2,:,:) 97 98 Do jk=1,jpk 99 Do jj=1,jpj 100 IF (umask(2,jj,jk).EQ.0.) THEN 101 ta(2,jj,jk) = ta(1,jj,jk) * tmask(2,jj,jk) 102 sa(2,jj,jk) = sa(1,jj,jk) * tmask(2,jj,jk) 103 ELSE 104 ta(2,jj,jk)=(alpha4*ta(1,jj,jk)+alpha3*ta(3,jj,jk))*tmask(2,jj,jk) 105 sa(2,jj,jk)=(alpha4*sa(1,jj,jk)+alpha3*sa(3,jj,jk))*tmask(2,jj,jk) 106 IF (un(2,jj,jk).LT.0.) THEN 107 ta(2,jj,jk)=(alpha6*ta(3,jj,jk)+alpha5*ta(1,jj,jk)+alpha7*ta(4,jj,jk))*tmask(2,jj,jk) 108 sa(2,jj,jk)=(alpha6*sa(3,jj,jk)+alpha5*sa(1,jj,jk)+alpha7*sa(4,jj,jk))*tmask(2,jj,jk) 109 ENDIF 110 ENDIF 111 End Do 112 enddo 113 ENDIF 114 115 IF ((nbondj == -1).OR.(nbondj == 2)) THEN 116 117 ta(:,1,:) = alpha1 * tatemp(:,1,:) + alpha2 * tatemp(:,2,:) 118 sa(:,1,:) = alpha1 * satemp(:,1,:) + alpha2 * satemp(:,2,:) 119 120 Do jk=1,jpk 121 Do ji=1,jpi 122 IF (vmask(ji,2,jk).EQ.0.) THEN 123 ta(ji,2,jk)=ta(ji,1,jk) * tmask(ji,2,jk) 124 sa(ji,2,jk)=sa(ji,1,jk) * tmask(ji,2,jk) 125 ELSE 126 ta(ji,2,jk)=(alpha4*ta(ji,1,jk)+alpha3*ta(ji,3,jk))*tmask(ji,2,jk) 127 sa(ji,2,jk)=(alpha4*sa(ji,1,jk)+alpha3*sa(ji,3,jk))*tmask(ji,2,jk) 128 IF (vn(ji,2,jk) .LT. 0.) THEN 129 ta(ji,2,jk)=(alpha6*ta(ji,3,jk)+alpha5*ta(ji,1,jk)+alpha7*ta(ji,4,jk))*tmask(ji,2,jk) 130 sa(ji,2,jk)=(alpha6*sa(ji,3,jk)+alpha5*sa(ji,1,jk)+alpha7*sa(ji,4,jk))*tmask(ji,2,jk) 131 ENDIF 132 ENDIF 133 End Do 134 enddo 135 ENDIF 136 137 Return 138 End Subroutine Agrif_tra 139 ! 140 ! 141 SUBROUTINE Agrif_dyn(kt) 142 ! 143 USE phycst 144 USE sol_oce 145 USE in_out_manager 146 147 implicit none 148 # include "domzgr_substitute.h90" 149 ! 150 REAL(wp) uatemp(jpi,jpj,jpk) , vatemp(jpi,jpj,jpk) 152 INTEGER, INTENT(in) :: kt 153 154 REAL(wp) :: timeref 155 REAL(wp) :: z2dt, znugdt 156 REAL(wp) :: zrhox, rhoy 157 REAL(wp), DIMENSION(jpi,jpj) :: zua2d, zva2d 158 REAL(wp), DIMENSION(jpi,jpj) :: spgu1,spgv1 159 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zua, zva 151 160 INTEGER :: ji,jj,jk 152 INTEGER kt153 REAL(wp) :: z2dt, znugdt154 REAL(wp), DIMENSION(jpi,jpj) :: uatemp2D, vatemp2D155 REAL(wp) :: timeref156 REAL(wp), DIMENSION(jpi,jpj) :: spgu1,spgv1157 REAL(wp) :: rhox, rhoy158 161 159 162 IF (Agrif_Root()) RETURN 160 163 161 rhox = Agrif_Rhox()164 zrhox = Agrif_Rhox() 162 165 rhoy = Agrif_Rhoy() 163 166 … … 171 174 znugdt = rnu * grav * z2dt 172 175 173 174 175 uatemp= 0.176 vatemp= 0.177 Call Agrif_Bc_variable(uatemp,un,procname=interpu)178 Call Agrif_Bc_variable(vatemp,vn,procname=interpv)179 uatemp2d = 0.180 vatemp2d = 0.181 182 183 184 Call Agrif_Bc_variable(uatemp2d,e1u,calledweight=1.,procname=interpu2d)185 Call Agrif_Bc_variable(vatemp2d,e2v,calledweight=1.,procname=interpv2d)186 187 188 189 If((nbondi == -1).OR.(nbondi == 2)) THEN190 191 DO jj=1,jpj192 laplacu(2,jj) = timeref * (uatemp2d(2,jj)/(rhoy*e2u(2,jj)))*umask(2,jj,1)193 ENDDO194 195 Dojk=1,jpkm1196 DO jj=1,jpj197 ua(1:2,jj,jk) = (uatemp(1:2,jj,jk)/(rhoy*e2u(1:2,jj)))198 #if ! defined key_zco 199 ua(1:2,jj,jk) = ua(1:2,jj,jk) / fse3u(1:2,jj,jk)200 #endif 201 ENDDO202 ENDDO203 204 Dojk=1,jpkm1205 DO jj=1,jpj206 ua(2,jj,jk) = (ua(2,jj,jk) - z2dt * znugdt * laplacu(2,jj))*umask(2,jj,jk)207 ENDDO208 ENDDO209 210 spgu(2,:)=0.211 212 dojk=1,jpkm1213 dojj=1,jpj214 spgu(2,jj)=spgu(2,jj)+fse3u(2,jj,jk)*ua(2,jj,jk)215 enddo216 enddo217 218 DO jj=1,jpj219 IF (umask(2,jj,1).NE.0.) THEN220 spgu(2,jj)=spgu(2,jj)/hu(2,jj)221 ENDIF222 enddo223 224 Dojk=1,jpkm1225 DO jj=1,jpj226 ua(2,jj,jk) = 0.25*(ua(1,jj,jk)+2.*ua(2,jj,jk)+ua(3,jj,jk))227 ua(2,jj,jk) = ua(2,jj,jk) * umask(2,jj,jk)228 ENDDO229 ENDDO230 231 spgu1(2,:)=0.232 233 dojk=1,jpkm1234 dojj=1,jpj235 spgu1(2,jj)=spgu1(2,jj)+fse3u(2,jj,jk)*ua(2,jj,jk)236 enddo237 enddo238 239 DO jj=1,jpj240 IF (umask(2,jj,1).NE.0.) THEN241 spgu1(2,jj)=spgu1(2,jj)/hu(2,jj)242 ENDIF243 enddo244 245 DO jk=1,jpkm1246 DO jj=1,jpj247 ua(2,jj,jk) = (ua(2,jj,jk)+spgu(2,jj)-spgu1(2,jj))*umask(2,jj,jk)248 ENDDO249 ENDDO250 251 Dojk=1,jpkm1252 Dojj=1,jpj253 va(2,jj,jk) = (vatemp(2,jj,jk)/(rhox*e1v(2,jj)))*vmask(2,jj,jk)254 #if ! defined key_zco 255 va(2,jj,jk) = va(2,jj,jk) / fse3v(2,jj,jk)176 Agrif_SpecialValue=0. 177 Agrif_UseSpecialValue = .TRUE. 178 zua = 0. 179 zva = 0. 180 CALL Agrif_Bc_variable(zua,un,procname=interpu) 181 CALL Agrif_Bc_variable(zva,vn,procname=interpv) 182 zua2d = 0. 183 zva2d = 0. 184 185 Agrif_SpecialValue=0. 186 Agrif_UseSpecialValue = .TRUE. 187 CALL Agrif_Bc_variable(zua2d,e1u,calledweight=1.,procname=interpu2d) 188 CALL Agrif_Bc_variable(zva2d,e2v,calledweight=1.,procname=interpv2d) 189 Agrif_UseSpecialValue = .FALSE. 190 191 192 IF((nbondi == -1).OR.(nbondi == 2)) THEN 193 194 DO jj=1,jpj 195 laplacu(2,jj) = timeref * (zua2d(2,jj)/(rhoy*e2u(2,jj)))*umask(2,jj,1) 196 END DO 197 198 DO jk=1,jpkm1 199 DO jj=1,jpj 200 ua(1:2,jj,jk) = (zua(1:2,jj,jk)/(rhoy*e2u(1:2,jj))) 201 #if ! defined key_zco 202 ua(1:2,jj,jk) = ua(1:2,jj,jk) / fse3u(1:2,jj,jk) 203 #endif 204 END DO 205 END DO 206 207 DO jk=1,jpkm1 208 DO jj=1,jpj 209 ua(2,jj,jk) = (ua(2,jj,jk) - z2dt * znugdt * laplacu(2,jj))*umask(2,jj,jk) 210 END DO 211 END DO 212 213 spgu(2,:)=0. 214 215 DO jk=1,jpkm1 216 DO jj=1,jpj 217 spgu(2,jj)=spgu(2,jj)+fse3u(2,jj,jk)*ua(2,jj,jk) 218 END DO 219 END DO 220 221 DO jj=1,jpj 222 IF (umask(2,jj,1).NE.0.) THEN 223 spgu(2,jj)=spgu(2,jj)/hu(2,jj) 224 ENDIF 225 END DO 226 227 DO jk=1,jpkm1 228 DO jj=1,jpj 229 ua(2,jj,jk) = 0.25*(ua(1,jj,jk)+2.*ua(2,jj,jk)+ua(3,jj,jk)) 230 ua(2,jj,jk) = ua(2,jj,jk) * umask(2,jj,jk) 231 END DO 232 END DO 233 234 spgu1(2,:)=0. 235 236 DO jk=1,jpkm1 237 DO jj=1,jpj 238 spgu1(2,jj)=spgu1(2,jj)+fse3u(2,jj,jk)*ua(2,jj,jk) 239 END DO 240 END DO 241 242 DO jj=1,jpj 243 IF (umask(2,jj,1).NE.0.) THEN 244 spgu1(2,jj)=spgu1(2,jj)/hu(2,jj) 245 ENDIF 246 END DO 247 248 DO jk=1,jpkm1 249 DO jj=1,jpj 250 ua(2,jj,jk) = (ua(2,jj,jk)+spgu(2,jj)-spgu1(2,jj))*umask(2,jj,jk) 251 END DO 252 END DO 253 254 DO jk=1,jpkm1 255 DO jj=1,jpj 256 va(2,jj,jk) = (zva(2,jj,jk)/(zrhox*e1v(2,jj)))*vmask(2,jj,jk) 257 #if ! defined key_zco 258 va(2,jj,jk) = va(2,jj,jk) / fse3v(2,jj,jk) 256 259 #endif 257 End Do258 End Do259 260 sshn(2,:)=sshn(3,:)261 sshb(2,:)=sshb(3,:)262 263 264 265 If((nbondi == 1).OR.(nbondi == 2)) THEN266 267 DO jj=1,jpj268 laplacu(nlci-2,jj) = timeref * (uatemp2d(nlci-2,jj)/(rhoy*e2u(nlci-2,jj)))269 ENDDO270 271 Dojk=1,jpkm1272 DO jj=1,jpj273 ua(nlci-2:nlci-1,jj,jk) = (uatemp(nlci-2:nlci-1,jj,jk)/(rhoy*e2u(nlci-2:nlci-1,jj)))274 275 #if ! defined key_zco 276 ua(nlci-2:nlci-1,jj,jk) = ua(nlci-2:nlci-1,jj,jk) / fse3u(nlci-2:nlci-1,jj,jk)277 #endif 278 279 ENDDO280 ENDDO281 282 Dojk=1,jpkm1283 DO jj=1,jpj284 ua(nlci-2,jj,jk) = (ua(nlci-2,jj,jk)- z2dt * znugdt * laplacu(nlci-2,jj))*umask(nlci-2,jj,jk)285 ENDDO286 ENDDO287 288 289 spgu(nlci-2,:)=0.290 291 do jk=1,jpkm1292 do jj=1,jpj293 spgu(nlci-2,jj)=spgu(nlci-2,jj)+fse3u(nlci-2,jj,jk)*ua(nlci-2,jj,jk)294 enddo295 enddo296 297 DO jj=1,jpj298 IF (umask(nlci-2,jj,1).NE.0.) THEN299 spgu(nlci-2,jj)=spgu(nlci-2,jj)/hu(nlci-2,jj)300 ENDIF301 enddo302 303 Dojk=1,jpkm1304 DO jj=1,jpj305 ua(nlci-2,jj,jk) = 0.25*(ua(nlci-3,jj,jk)+2.*ua(nlci-2,jj,jk)+ua(nlci-1,jj,jk))306 307 ua(nlci-2,jj,jk) = ua(nlci-2,jj,jk) * umask(nlci-2,jj,jk)308 309 ENDDO310 ENDDO311 312 spgu1(nlci-2,:)=0.313 314 dojk=1,jpkm1315 dojj=1,jpj316 spgu1(nlci-2,jj)=spgu1(nlci-2,jj)+fse3u(nlci-2,jj,jk)*ua(nlci-2,jj,jk)*umask(nlci-2,jj,jk)317 enddo318 enddo319 320 DO jj=1,jpj321 IF (umask(nlci-2,jj,1).NE.0.) THEN322 spgu1(nlci-2,jj)=spgu1(nlci-2,jj)/hu(nlci-2,jj)323 ENDIF324 enddo325 326 DO jk=1,jpkm1327 DO jj=1,jpj328 ua(nlci-2,jj,jk) = (ua(nlci-2,jj,jk)+spgu(nlci-2,jj)-spgu1(nlci-2,jj))*umask(nlci-2,jj,jk)329 ENDDO330 ENDDO331 332 Dojk=1,jpkm1333 Dojj=1,jpj-1334 va(nlci-1,jj,jk) = (vatemp(nlci-1,jj,jk)/(rhox*e1v(nlci-1,jj)))*vmask(nlci-1,jj,jk)335 #if ! defined key_zco 336 va(nlci-1,jj,jk) = va(nlci-1,jj,jk) / fse3v(nlci-1,jj,jk)337 #endif 338 End Do339 End Do340 341 sshn(nlci-1,:)=sshn(nlci-2,:)342 sshb(nlci-1,:)=sshb(nlci-2,:)343 344 345 If((nbondj == -1).OR.(nbondj == 2)) THEN346 347 DO ji=1,jpi348 laplacv(ji,2) = timeref * (vatemp2d(ji,2)/(rhox*e1v(ji,2)))349 ENDDO350 351 DO jk=1,jpkm1352 DO ji=1,jpi353 va(ji,1:2,jk) = (vatemp(ji,1:2,jk)/(rhox*e1v(ji,1:2)))354 #if ! defined key_zco 355 va(ji,1:2,jk) = va(ji,1:2,jk) / fse3v(ji,1:2,jk)356 #endif 357 ENDDO358 ENDDO359 360 DO jk=1,jpkm1361 DO ji=1,jpi362 va(ji,2,jk) = (va(ji,2,jk) - z2dt * znugdt * laplacv(ji,2))*vmask(ji,2,jk)363 ENDDO364 ENDDO365 366 spgv(:,2)=0.367 368 dojk=1,jpkm1369 doji=1,jpi370 spgv(ji,2)=spgv(ji,2)+fse3v(ji,2,jk)*va(ji,2,jk)371 enddo372 enddo373 374 DO ji=1,jpi375 IF (vmask(ji,2,1).NE.0.) THEN376 spgv(ji,2)=spgv(ji,2)/hv(ji,2)377 ENDIF378 enddo379 380 DO jk=1,jpkm1381 DO ji=1,jpi382 va(ji,2,jk)=0.25*(va(ji,1,jk)+2.*va(ji,2,jk)+va(ji,3,jk))383 va(ji,2,jk)=va(ji,2,jk)*vmask(ji,2,jk)384 ENDDO385 ENDDO386 387 spgv1(:,2)=0.388 389 dojk=1,jpkm1390 doji=1,jpi391 spgv1(ji,2)=spgv1(ji,2)+fse3v(ji,2,jk)*va(ji,2,jk)*vmask(ji,2,jk)392 enddo393 enddo394 395 DO ji=1,jpi396 IF (vmask(ji,2,1).NE.0.) THEN397 spgv1(ji,2)=spgv1(ji,2)/hv(ji,2)398 ENDIF399 enddo400 401 DO jk=1,jpkm1402 DO ji=1,jpi403 va(ji,2,jk) = (va(ji,2,jk)+spgv(ji,2)-spgv1(ji,2))*vmask(ji,2,jk)404 ENDDO405 ENDDO406 407 DO jk=1,jpkm1408 DO ji=1,jpi409 ua(ji,2,jk) = (uatemp(ji,2,jk)/(rhoy*e2u(ji,2)))*umask(ji,2,jk)410 #if ! defined key_zco 411 ua(ji,2,jk) = ua(ji,2,jk) / fse3u(ji,2,jk)260 END DO 261 END DO 262 263 sshn(2,:)=sshn(3,:) 264 sshb(2,:)=sshb(3,:) 265 266 ENDIF 267 268 IF((nbondi == 1).OR.(nbondi == 2)) THEN 269 270 DO jj=1,jpj 271 laplacu(nlci-2,jj) = timeref * (zua2d(nlci-2,jj)/(rhoy*e2u(nlci-2,jj))) 272 END DO 273 274 DO jk=1,jpkm1 275 DO jj=1,jpj 276 ua(nlci-2:nlci-1,jj,jk) = (zua(nlci-2:nlci-1,jj,jk)/(rhoy*e2u(nlci-2:nlci-1,jj))) 277 278 #if ! defined key_zco 279 ua(nlci-2:nlci-1,jj,jk) = ua(nlci-2:nlci-1,jj,jk) / fse3u(nlci-2:nlci-1,jj,jk) 280 #endif 281 282 END DO 283 END DO 284 285 DO jk=1,jpkm1 286 DO jj=1,jpj 287 ua(nlci-2,jj,jk) = (ua(nlci-2,jj,jk)- z2dt * znugdt * laplacu(nlci-2,jj))*umask(nlci-2,jj,jk) 288 END DO 289 END DO 290 291 292 spgu(nlci-2,:)=0. 293 294 do jk=1,jpkm1 295 do jj=1,jpj 296 spgu(nlci-2,jj)=spgu(nlci-2,jj)+fse3u(nlci-2,jj,jk)*ua(nlci-2,jj,jk) 297 enddo 298 enddo 299 300 DO jj=1,jpj 301 IF (umask(nlci-2,jj,1).NE.0.) THEN 302 spgu(nlci-2,jj)=spgu(nlci-2,jj)/hu(nlci-2,jj) 303 ENDIF 304 END DO 305 306 DO jk=1,jpkm1 307 DO jj=1,jpj 308 ua(nlci-2,jj,jk) = 0.25*(ua(nlci-3,jj,jk)+2.*ua(nlci-2,jj,jk)+ua(nlci-1,jj,jk)) 309 310 ua(nlci-2,jj,jk) = ua(nlci-2,jj,jk) * umask(nlci-2,jj,jk) 311 312 END DO 313 END DO 314 315 spgu1(nlci-2,:)=0. 316 317 DO jk=1,jpkm1 318 DO jj=1,jpj 319 spgu1(nlci-2,jj)=spgu1(nlci-2,jj)+fse3u(nlci-2,jj,jk)*ua(nlci-2,jj,jk)*umask(nlci-2,jj,jk) 320 END DO 321 END DO 322 323 DO jj=1,jpj 324 IF (umask(nlci-2,jj,1).NE.0.) THEN 325 spgu1(nlci-2,jj)=spgu1(nlci-2,jj)/hu(nlci-2,jj) 326 ENDIF 327 END DO 328 329 DO jk=1,jpkm1 330 DO jj=1,jpj 331 ua(nlci-2,jj,jk) = (ua(nlci-2,jj,jk)+spgu(nlci-2,jj)-spgu1(nlci-2,jj))*umask(nlci-2,jj,jk) 332 END DO 333 END DO 334 335 DO jk=1,jpkm1 336 DO jj=1,jpj-1 337 va(nlci-1,jj,jk) = (zva(nlci-1,jj,jk)/(zrhox*e1v(nlci-1,jj)))*vmask(nlci-1,jj,jk) 338 #if ! defined key_zco 339 va(nlci-1,jj,jk) = va(nlci-1,jj,jk) / fse3v(nlci-1,jj,jk) 340 #endif 341 END DO 342 END DO 343 344 sshn(nlci-1,:)=sshn(nlci-2,:) 345 sshb(nlci-1,:)=sshb(nlci-2,:) 346 ENDIF 347 348 IF((nbondj == -1).OR.(nbondj == 2)) THEN 349 350 DO ji=1,jpi 351 laplacv(ji,2) = timeref * (zva2d(ji,2)/(zrhox*e1v(ji,2))) 352 END DO 353 354 DO jk=1,jpkm1 355 DO ji=1,jpi 356 va(ji,1:2,jk) = (zva(ji,1:2,jk)/(zrhox*e1v(ji,1:2))) 357 #if ! defined key_zco 358 va(ji,1:2,jk) = va(ji,1:2,jk) / fse3v(ji,1:2,jk) 359 #endif 360 END DO 361 END DO 362 363 DO jk=1,jpkm1 364 DO ji=1,jpi 365 va(ji,2,jk) = (va(ji,2,jk) - z2dt * znugdt * laplacv(ji,2))*vmask(ji,2,jk) 366 END DO 367 END DO 368 369 spgv(:,2)=0. 370 371 DO jk=1,jpkm1 372 DO ji=1,jpi 373 spgv(ji,2)=spgv(ji,2)+fse3v(ji,2,jk)*va(ji,2,jk) 374 END DO 375 END DO 376 377 DO ji=1,jpi 378 IF (vmask(ji,2,1).NE.0.) THEN 379 spgv(ji,2)=spgv(ji,2)/hv(ji,2) 380 ENDIF 381 END DO 382 383 DO jk=1,jpkm1 384 DO ji=1,jpi 385 va(ji,2,jk)=0.25*(va(ji,1,jk)+2.*va(ji,2,jk)+va(ji,3,jk)) 386 va(ji,2,jk)=va(ji,2,jk)*vmask(ji,2,jk) 387 END DO 388 END DO 389 390 spgv1(:,2)=0. 391 392 DO jk=1,jpkm1 393 DO ji=1,jpi 394 spgv1(ji,2)=spgv1(ji,2)+fse3v(ji,2,jk)*va(ji,2,jk)*vmask(ji,2,jk) 395 END DO 396 END DO 397 398 DO ji=1,jpi 399 IF (vmask(ji,2,1).NE.0.) THEN 400 spgv1(ji,2)=spgv1(ji,2)/hv(ji,2) 401 ENDIF 402 END DO 403 404 DO jk=1,jpkm1 405 DO ji=1,jpi 406 va(ji,2,jk) = (va(ji,2,jk)+spgv(ji,2)-spgv1(ji,2))*vmask(ji,2,jk) 407 END DO 408 END DO 409 410 DO jk=1,jpkm1 411 DO ji=1,jpi 412 ua(ji,2,jk) = (zua(ji,2,jk)/(rhoy*e2u(ji,2)))*umask(ji,2,jk) 413 #if ! defined key_zco 414 ua(ji,2,jk) = ua(ji,2,jk) / fse3u(ji,2,jk) 412 415 #endif 413 ENDDO414 ENDDO415 416 sshn(:,2)=sshn(:,3)417 sshb(:,2)=sshb(:,3)418 419 420 If((nbondj == 1).OR.(nbondj == 2)) THEN421 422 DO ji=1,jpi423 laplacv(ji,nlcj-2) = timeref * (vatemp2d(ji,nlcj-2)/(rhox*e1v(ji,nlcj-2)))424 ENDDO425 426 DO jk=1,jpkm1427 DO ji=1,jpi428 va(ji,nlcj-2:nlcj-1,jk) = (vatemp(ji,nlcj-2:nlcj-1,jk)/(rhox*e1v(ji,nlcj-2:nlcj-1)))429 #if ! defined key_zco 430 va(ji,nlcj-2:nlcj-1,jk) = va(ji,nlcj-2:nlcj-1,jk) / fse3v(ji,nlcj-2:nlcj-1,jk)431 #endif 432 ENDDO433 ENDDO434 435 DO jk=1,jpkm1436 DO ji=1,jpi437 va(ji,nlcj-2,jk) = (va(ji,nlcj-2,jk)-z2dt * znugdt * laplacv(ji,nlcj-2))*vmask(ji,nlcj-2,jk)438 ENDDO439 ENDDO440 441 442 spgv(:,nlcj-2)=0.443 444 dojk=1,jpkm1445 doji=1,jpi446 spgv(ji,nlcj-2)=spgv(ji,nlcj-2)+fse3v(ji,nlcj-2,jk)*va(ji,nlcj-2,jk)447 enddo448 enddo449 450 DO ji=1,jpi451 IF (vmask(ji,nlcj-2,1).NE.0.) THEN452 spgv(ji,nlcj-2)=spgv(ji,nlcj-2)/hv(ji,nlcj-2)453 ENDIF454 enddo455 456 DO jk=1,jpkm1457 DO ji=1,jpi458 va(ji,nlcj-2,jk)=0.25*(va(ji,nlcj-3,jk)+2.*va(ji,nlcj-2,jk)+va(ji,nlcj-1,jk))459 va(ji,nlcj-2,jk) = va(ji,nlcj-2,jk) * vmask(ji,nlcj-2,jk)460 ENDDO461 ENDDO462 463 spgv1(:,nlcj-2)=0.464 465 dojk=1,jpkm1466 doji=1,jpi467 spgv1(ji,nlcj-2)=spgv1(ji,nlcj-2)+fse3v(ji,nlcj-2,jk)*va(ji,nlcj-2,jk)468 enddo469 enddo470 471 DO ji=1,jpi472 IF (vmask(ji,nlcj-2,1).NE.0.) THEN473 spgv1(ji,nlcj-2)=spgv1(ji,nlcj-2)/hv(ji,nlcj-2)474 ENDIF475 enddo476 477 DO jk=1,jpkm1478 DO ji=1,jpi479 va(ji,nlcj-2,jk) = (va(ji,nlcj-2,jk)+spgv(ji,nlcj-2)-spgv1(ji,nlcj-2))*vmask(ji,nlcj-2,jk)480 ENDDO481 ENDDO482 483 DO jk=1,jpkm1484 DO ji=1,jpi485 ua(ji,nlcj-1,jk) = (uatemp(ji,nlcj-1,jk)/(rhoy*e2u(ji,nlcj-1)))*umask(ji,nlcj-1,jk)486 #if ! defined key_zco 487 ua(ji,nlcj-1,jk) = ua(ji,nlcj-1,jk) / fse3u(ji,nlcj-1,jk)416 END DO 417 END DO 418 419 sshn(:,2)=sshn(:,3) 420 sshb(:,2)=sshb(:,3) 421 ENDIF 422 423 IF((nbondj == 1).OR.(nbondj == 2)) THEN 424 425 DO ji=1,jpi 426 laplacv(ji,nlcj-2) = timeref * (zva2d(ji,nlcj-2)/(zrhox*e1v(ji,nlcj-2))) 427 END DO 428 429 DO jk=1,jpkm1 430 DO ji=1,jpi 431 va(ji,nlcj-2:nlcj-1,jk) = (zva(ji,nlcj-2:nlcj-1,jk)/(zrhox*e1v(ji,nlcj-2:nlcj-1))) 432 #if ! defined key_zco 433 va(ji,nlcj-2:nlcj-1,jk) = va(ji,nlcj-2:nlcj-1,jk) / fse3v(ji,nlcj-2:nlcj-1,jk) 434 #endif 435 END DO 436 END DO 437 438 DO jk=1,jpkm1 439 DO ji=1,jpi 440 va(ji,nlcj-2,jk) = (va(ji,nlcj-2,jk)-z2dt * znugdt * laplacv(ji,nlcj-2))*vmask(ji,nlcj-2,jk) 441 END DO 442 END DO 443 444 445 spgv(:,nlcj-2)=0. 446 447 DO jk=1,jpkm1 448 DO ji=1,jpi 449 spgv(ji,nlcj-2)=spgv(ji,nlcj-2)+fse3v(ji,nlcj-2,jk)*va(ji,nlcj-2,jk) 450 END DO 451 END DO 452 453 DO ji=1,jpi 454 IF (vmask(ji,nlcj-2,1).NE.0.) THEN 455 spgv(ji,nlcj-2)=spgv(ji,nlcj-2)/hv(ji,nlcj-2) 456 ENDIF 457 END DO 458 459 DO jk=1,jpkm1 460 DO ji=1,jpi 461 va(ji,nlcj-2,jk)=0.25*(va(ji,nlcj-3,jk)+2.*va(ji,nlcj-2,jk)+va(ji,nlcj-1,jk)) 462 va(ji,nlcj-2,jk) = va(ji,nlcj-2,jk) * vmask(ji,nlcj-2,jk) 463 END DO 464 END DO 465 466 spgv1(:,nlcj-2)=0. 467 468 DO jk=1,jpkm1 469 DO ji=1,jpi 470 spgv1(ji,nlcj-2)=spgv1(ji,nlcj-2)+fse3v(ji,nlcj-2,jk)*va(ji,nlcj-2,jk) 471 END DO 472 END DO 473 474 DO ji=1,jpi 475 IF (vmask(ji,nlcj-2,1).NE.0.) THEN 476 spgv1(ji,nlcj-2)=spgv1(ji,nlcj-2)/hv(ji,nlcj-2) 477 ENDIF 478 END DO 479 480 DO jk=1,jpkm1 481 DO ji=1,jpi 482 va(ji,nlcj-2,jk) = (va(ji,nlcj-2,jk)+spgv(ji,nlcj-2)-spgv1(ji,nlcj-2))*vmask(ji,nlcj-2,jk) 483 END DO 484 END DO 485 486 DO jk=1,jpkm1 487 DO ji=1,jpi 488 ua(ji,nlcj-1,jk) = (zua(ji,nlcj-1,jk)/(rhoy*e2u(ji,nlcj-1)))*umask(ji,nlcj-1,jk) 489 #if ! defined key_zco 490 ua(ji,nlcj-1,jk) = ua(ji,nlcj-1,jk) / fse3u(ji,nlcj-1,jk) 488 491 #endif 489 ENDDO 490 ENDDO 491 492 sshn(:,nlcj-1)=sshn(:,nlcj-2) 493 sshb(:,nlcj-1)=sshb(:,nlcj-2) 494 ENDIF 495 496 ! 497 Return 498 End Subroutine Agrif_dyn 499 500 501 subroutine interpu(tabres,i1,i2,j1,j2,k1,k2) 502 Implicit none 503 # include "domzgr_substitute.h90" 504 integer i1,i2,j1,j2,k1,k2 505 integer ji,jj,jk 506 real,dimension(i1:i2,j1:j2,k1:k2) :: tabres 507 508 do jk=k1,k2 509 DO jj=j1,j2 510 DO ji=i1,i2 511 tabres(ji,jj,jk) = e2u(ji,jj) * un(ji,jj,jk) 512 #if ! defined key_zco 513 tabres(ji,jj,jk) = tabres(ji,jj,jk) * fse3u(ji,jj,jk) 514 #endif 515 ENDDO 516 ENDDO 517 ENDDO 518 end subroutine interpu 519 520 subroutine interpu2d(tabres,i1,i2,j1,j2) 521 Implicit none 522 integer i1,i2,j1,j2 523 integer ji,jj 524 real,dimension(i1:i2,j1:j2) :: tabres 525 526 DO jj=j1,j2 527 DO ji=i1,i2 528 tabres(ji,jj) = e2u(ji,jj) * ((gcx(ji+1,jj) - gcx(ji,jj))/e1u(ji,jj)) & 529 *umask(ji,jj,1) 530 ENDDO 531 ENDDO 532 end subroutine interpu2d 533 534 subroutine interpv(tabres,i1,i2,j1,j2,k1,k2) 535 Implicit none 536 # include "domzgr_substitute.h90" 537 integer i1,i2,j1,j2,k1,k2 538 integer ji,jj,jk 539 real,dimension(i1:i2,j1:j2,k1:k2) :: tabres 540 541 do jk=k1,k2 542 DO jj=j1,j2 543 DO ji=i1,i2 544 tabres(ji,jj,jk) = e1v(ji,jj) * vn(ji,jj,jk) 545 #if ! defined key_zco 546 tabres(ji,jj,jk) = tabres(ji,jj,jk) * fse3v(ji,jj,jk) 492 END DO 493 END DO 494 495 sshn(:,nlcj-1)=sshn(:,nlcj-2) 496 sshb(:,nlcj-1)=sshb(:,nlcj-2) 497 ENDIF 498 499 END SUBROUTINE Agrif_dyn 500 501 SUBROUTINE interpu(tabres,i1,i2,j1,j2,k1,k2) 502 !!--------------------------------------------- 503 !! *** ROUTINE interpu *** 504 !!--------------------------------------------- 505 # include "domzgr_substitute.h90" 506 507 INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2 508 REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres 509 510 INTEGER :: ji,jj,jk 511 512 DO jk=k1,k2 513 DO jj=j1,j2 514 DO ji=i1,i2 515 tabres(ji,jj,jk) = e2u(ji,jj) * un(ji,jj,jk) 516 #if ! defined key_zco 517 tabres(ji,jj,jk) = tabres(ji,jj,jk) * fse3u(ji,jj,jk) 518 #endif 519 END DO 520 END DO 521 END DO 522 END SUBROUTINE interpu 523 524 SUBROUTINE interpu2d(tabres,i1,i2,j1,j2) 525 !!--------------------------------------------- 526 !! *** ROUTINE interpu2d *** 527 !!--------------------------------------------- 528 529 INTEGER, INTENT(in) :: i1,i2,j1,j2 530 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres 531 532 INTEGER :: ji,jj 533 534 DO jj=j1,j2 535 DO ji=i1,i2 536 tabres(ji,jj) = e2u(ji,jj) * ((gcx(ji+1,jj) - gcx(ji,jj))/e1u(ji,jj)) & 537 * umask(ji,jj,1) 538 END DO 539 END DO 540 541 END SUBROUTINE interpu2d 542 543 SUBROUTINE interpv(tabres,i1,i2,j1,j2,k1,k2) 544 !!--------------------------------------------- 545 !! *** ROUTINE interpv *** 546 !!--------------------------------------------- 547 # include "domzgr_substitute.h90" 548 549 INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2 550 REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres 551 552 INTEGER :: ji, jj, jk 553 554 DO jk=k1,k2 555 DO jj=j1,j2 556 DO ji=i1,i2 557 tabres(ji,jj,jk) = e1v(ji,jj) * vn(ji,jj,jk) 558 #if ! defined key_zco 559 tabres(ji,jj,jk) = tabres(ji,jj,jk) * fse3v(ji,jj,jk) 547 560 #endif 548 ENDDO 549 ENDDO 550 ENDDO 551 end subroutine interpv 552 553 subroutine interpv2d(tabres,i1,i2,j1,j2) 554 Implicit none 555 integer i1,i2,j1,j2 556 integer ji,jj 557 real,dimension(i1:i2,j1:j2) :: tabres 558 559 DO jj=j1,j2 560 DO ji=i1,i2 561 tabres(ji,jj) = e1v(ji,jj) * ((gcx(ji,jj+1) - gcx(ji,jj))/e2v(ji,jj)) & 562 * vmask(ji,jj,1) 563 ENDDO 564 ENDDO 565 end subroutine interpv2d 561 END DO 562 END DO 563 END DO 564 565 END SUBROUTINE interpv 566 567 SUBROUTINE interpv2d(tabres,i1,i2,j1,j2) 568 !!--------------------------------------------- 569 !! *** ROUTINE interpv2d *** 570 !!--------------------------------------------- 571 572 INTEGER, INTENT(in) :: i1,i2,j1,j2 573 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres 574 575 INTEGER :: ji,jj 576 577 DO jj=j1,j2 578 DO ji=i1,i2 579 tabres(ji,jj) = e1v(ji,jj) * ((gcx(ji,jj+1) - gcx(ji,jj))/e2v(ji,jj)) & 580 * vmask(ji,jj,1) 581 END DO 582 END DO 583 584 END SUBROUTINE interpv2d 566 585 567 586 #else 568 CONTAINS 569 subroutine Agrif_OPA_Interp_empty 570 571 end subroutine Agrif_OPA_Interp_empty 572 #endif 573 End Module agrif_opa_interp 574 587 CONTAINS 588 589 SUBROUTINE Agrif_OPA_Interp_empty 590 !!--------------------------------------------- 591 !! *** ROUTINE agrif_OPA_Interp_empty *** 592 !!--------------------------------------------- 593 WRITE(*,*) 'agrif_opa_interp : You should not have seen this print! error?' 594 END SUBROUTINE Agrif_OPA_Interp_empty 595 #endif 596 END MODULE agrif_opa_interp 597
Note: See TracChangeset
for help on using the changeset viewer.