Changeset 13204 for utils/tools/DOMAINcfg/src/agrif_user.F90
- Timestamp:
- 2020-07-02T10:38:35+02:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
utils/tools/DOMAINcfg/src/agrif_user.F90
r12414 r13204 1 1 #if defined key_agrif 2 subroutine agrif_initworkspace() 2 SUBROUTINE agrif_user() 3 !!---------------------------------------------------------------------- 4 !! *** ROUTINE agrif_user *** 5 !!---------------------------------------------------------------------- 6 END SUBROUTINE agrif_user 7 8 SUBROUTINE agrif_initworkspace() 3 9 !!---------------------------------------------------------------------- 4 10 !! *** ROUTINE Agrif_InitWorkspace *** 5 11 !!---------------------------------------------------------------------- 6 end subroutine agrif_initworkspace 7 SUBROUTINE Agrif_InitValues 12 END SUBROUTINE agrif_initworkspace 13 14 SUBROUTINE Agrif_InitValues 8 15 !!---------------------------------------------------------------------- 9 16 !! *** ROUTINE Agrif_InitValues *** … … 11 18 !! ** Purpose :: Declaration of variables to be interpolated 12 19 !!---------------------------------------------------------------------- 13 USE Agrif_Util 14 USE dom_oce 15 USE nemogcm 16 USE domain 17 !! 18 IMPLICIT NONE 19 20 21 CALL nemo_init !* Initializations of each fine grid 22 23 CALL dom_nam 24 CALL cfg_write ! create the configuration file 25 26 END SUBROUTINE Agrif_InitValues 27 28 SUBROUTINE Agrif_InitValues_cont 29 30 use dom_oce 31 integer :: irafx, irafy 32 logical :: ln_perio 33 integer nx,ny 34 35 irafx = agrif_irhox() 36 irafy = agrif_irhoy() 37 38 nx=nlci ; ny=nlcj 20 USE Agrif_Util 21 USE dom_oce 22 USE nemogcm 23 USE domain 24 !! 25 IMPLICIT NONE 26 27 ! No temporal refinement 28 CALL Agrif_Set_coeffreft(1) 29 30 CALL nemo_init !* Initializations of each fine grid 31 32 CALL dom_nam 33 34 END SUBROUTINE Agrif_InitValues 35 36 SUBROUTINE Agrif_InitValues_cont 37 !!---------------------------------------------------------------------- 38 !! *** ROUTINE Agrif_InitValues_cont *** 39 !! 40 !! ** Purpose :: Initialisation of variables to be interpolated 41 !!---------------------------------------------------------------------- 42 USE dom_oce 43 USE lbclnk 44 !! 45 IMPLICIT NONE 46 ! 47 INTEGER :: nx, ny 48 INTEGER :: irafx, irafy 49 LOGICAL :: ln_perio 50 ! 51 irafx = agrif_irhox() 52 irafy = agrif_irhoy() 53 54 nx = nlci ; ny = nlcj 39 55 40 56 ! IF(jperio /=1 .AND. jperio/=4 .AND. jperio/=6 ) THEN … … 44 60 ! nbghostcellsfine_tot_y= nbghostcellsfine_y +1 45 61 ! ELSE 46 ! nx = (nbcellsx)+2*nbghostcellsfine_x 62 ! nx = (nbcellsx)+2*nbghostcellsfine_x 47 63 ! ny = (nbcellsy)+2*nbghostcellsfine+2 48 64 ! nbghostcellsfine_tot_x= 1 … … 53 69 ! nx = nbcellsx+irafx 54 70 ! ny = nbcellsy+irafy 55 56 WRITE(*,*) ' ' 57 WRITE(*,*)'Size of the High resolution grid: ',nx,' x ',ny 58 WRITE(*,*) ' ' 59 60 call agrif_init_lonlat() 61 ln_perio=.FALSE. 62 if( jperio ==1 .OR. jperio==2 .OR. jperio==4) ln_perio=.TRUE. 63 64 where (glamt < -180) glamt = glamt +360. 65 if (ln_perio) then 66 glamt(1,:)=glamt(nx-1,:) 67 glamt(nx,:)=glamt(2,:) 68 endif 69 70 where (glamu < -180) glamu = glamu +360. 71 if (ln_perio) then 72 glamu(1,:)=glamu(nx-1,:) 73 glamu(nx,:)=glamu(2,:) 74 endif 75 76 where (glamv < -180) glamv = glamv +360. 77 if (ln_perio) then 78 glamv(1,:)=glamv(nx-1,:) 79 glamv(nx,:)=glamv(2,:) 80 endif 81 82 where (glamf < -180) glamf = glamf +360. 83 if (ln_perio) then 84 glamf(1,:)=glamf(nx-1,:) 85 glamf(nx,:)=glamf(2,:) 86 endif 87 88 call agrif_init_scales() 89 90 91 END SUBROUTINE Agrif_InitValues_cont 92 93 94 subroutine agrif_declare_var() 95 use par_oce 96 use dom_oce 97 use agrif_profiles 98 use agrif_parameters 99 100 IMPLICIT NONE 101 102 integer :: ind1, ind2, ind3 103 integer nx,ny 104 integer nbghostcellsfine_tot_x, nbghostcellsfine_tot_y 105 INTEGER :: irafx 106 !!---------------------------------------------------------------------- 107 108 ! 1. Declaration of the type of variable which have to be interpolated 109 !--------------------------------------------------------------------- 110 nx=nlci ; ny=nlcj 111 112 !ind2 = nbghostcellsfine_tot_x + 1 113 !ind3 = nbghostcellsfine_tot_y + 1 114 ind2 = 2 + nbghostcells 115 ind3 = ind2 116 nbghostcellsfine_tot_x=nbghostcells+1 117 nbghostcellsfine_tot_y=nbghostcells+1 118 119 irafx = Agrif_irhox() 120 121 CALL agrif_nemo_init ! specific namelist part if needed 122 123 CALL agrif_declare_variable((/2,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),glamt_id) 124 CALL agrif_declare_variable((/1,2/),(/ind2-1,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),glamu_id) 125 CALL agrif_declare_variable((/2,1/),(/ind2,ind3-1/),(/'x','y'/),(/1,1/),(/nx,ny/),glamv_id) 126 CALL agrif_declare_variable((/1,1/),(/ind2-1,ind3-1/),(/'x','y'/),(/1,1/),(/nx,ny/),glamf_id) 127 128 CALL agrif_declare_variable((/2,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),gphit_id) 129 CALL agrif_declare_variable((/1,2/),(/ind2-1,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),gphiu_id) 130 CALL agrif_declare_variable((/2,1/),(/ind2,ind3-1/),(/'x','y'/),(/1,1/),(/nx,ny/),gphiv_id) 131 CALL agrif_declare_variable((/1,1/),(/ind2-1,ind3-1/),(/'x','y'/),(/1,1/),(/nx,ny/),gphif_id) 132 133 ! Horizontal scale factors 134 135 CALL agrif_declare_variable((/2,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),e1t_id) 136 CALL agrif_declare_variable((/1,2/),(/ind2-1,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),e1u_id) 137 CALL agrif_declare_variable((/2,1/),(/ind2,ind3-1/),(/'x','y'/),(/1,1/),(/nx,ny/),e1v_id) 138 CALL agrif_declare_variable((/1,1/),(/ind2-1,ind3-1/),(/'x','y'/),(/1,1/),(/nx,ny/),e1f_id) 139 140 CALL agrif_declare_variable((/2,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),e2t_id) 141 CALL agrif_declare_variable((/1,2/),(/ind2-1,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),e2u_id) 142 CALL agrif_declare_variable((/2,1/),(/ind2,ind3-1/),(/'x','y'/),(/1,1/),(/nx,ny/),e2v_id) 143 CALL agrif_declare_variable((/1,1/),(/ind2-1,ind3-1/),(/'x','y'/),(/1,1/),(/nx,ny/),e2f_id) 144 145 ! Bathymetry 146 147 CALL agrif_declare_variable((/2,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),bathy_id) 148 149 ! Vertical scale factors 150 CALL agrif_declare_variable((/2,2,0/),(/ind2,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nx,ny,jpk/),e3t_id) 151 CALL agrif_declare_variable((/2,2,0/),(/ind2,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nx,ny,jpk/),e3t_copy_id) 152 CALL agrif_declare_variable((/2,2,0/),(/ind2,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nx,ny,jpk+1/),e3t_connect_id) 153 154 CALL agrif_declare_variable((/1,2,0/),(/ind2-1,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nx,ny,jpk/),e3u_id) 155 CALL agrif_declare_variable((/2,1,0/),(/ind2,ind3-1,0/),(/'x','y','N'/),(/1,1,1/),(/nx,ny,jpk/),e3v_id) 156 157 ! Bottom level 158 159 CALL agrif_declare_variable((/2,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),bottom_level_id) 160 161 CALL Agrif_Set_bcinterp(glamt_id,interp=AGRIF_linear) 162 CALL Agrif_Set_interp(glamt_id,interp=AGRIF_linear) 163 CALL Agrif_Set_bc( glamt_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 164 165 CALL Agrif_Set_bcinterp(glamu_id,interp=AGRIF_linear) 166 CALL Agrif_Set_interp(glamu_id,interp=AGRIF_linear) 167 CALL Agrif_Set_bc( glamu_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 168 169 CALL Agrif_Set_bcinterp(glamv_id,interp=AGRIF_linear) 170 CALL Agrif_Set_interp(glamv_id,interp=AGRIF_linear) 171 CALL Agrif_Set_bc( glamv_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 172 173 CALL Agrif_Set_bcinterp(glamf_id,interp=AGRIF_linear) 174 CALL Agrif_Set_interp(glamf_id,interp=AGRIF_linear) 175 CALL Agrif_Set_bc( glamf_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 176 177 CALL Agrif_Set_bcinterp(gphit_id,interp=AGRIF_linear) 178 CALL Agrif_Set_interp(gphit_id,interp=AGRIF_linear) 179 CALL Agrif_Set_bc( gphit_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 180 181 CALL Agrif_Set_bcinterp(gphiu_id,interp=AGRIF_linear) 182 CALL Agrif_Set_interp(gphiu_id,interp=AGRIF_linear) 183 CALL Agrif_Set_bc( gphiu_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 184 185 CALL Agrif_Set_bcinterp(gphiv_id,interp=AGRIF_linear) 186 CALL Agrif_Set_interp(gphiv_id,interp=AGRIF_linear) 187 CALL Agrif_Set_bc( gphiv_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 188 189 CALL Agrif_Set_bcinterp(gphif_id,interp=AGRIF_linear) 190 CALL Agrif_Set_interp(gphif_id,interp=AGRIF_linear) 191 CALL Agrif_Set_bc( gphif_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 192 193 ! 194 195 CALL Agrif_Set_bcinterp(e1t_id,interp=AGRIF_ppm) 196 CALL Agrif_Set_interp(e1t_id,interp=AGRIF_ppm) 197 CALL Agrif_Set_bc( e1t_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 198 199 CALL Agrif_Set_bcinterp(e1u_id, interp1=Agrif_linear, interp2=AGRIF_ppm) 200 CALL Agrif_Set_interp(e1u_id, interp1=Agrif_linear, interp2=AGRIF_ppm) 201 CALL Agrif_Set_bc( e1u_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 202 203 CALL Agrif_Set_bcinterp(e1v_id,interp1=AGRIF_ppm, interp2=Agrif_linear) 204 CALL Agrif_Set_interp(e1v_id, interp1=AGRIF_ppm, interp2=Agrif_linear) 205 CALL Agrif_Set_bc( e1v_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 206 207 CALL Agrif_Set_bcinterp(e1f_id,interp=AGRIF_linear) 208 CALL Agrif_Set_interp(e1f_id,interp=AGRIF_linear) 209 CALL Agrif_Set_bc( e1f_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 210 211 CALL Agrif_Set_bcinterp(e2t_id,interp=AGRIF_ppm) 212 CALL Agrif_Set_interp(e2t_id,interp=AGRIF_ppm) 213 CALL Agrif_Set_bc( e2t_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 214 215 CALL Agrif_Set_bcinterp(e2u_id,interp1=Agrif_linear, interp2=AGRIF_ppm) 216 CALL Agrif_Set_interp(e2u_id,interp1=Agrif_linear, interp2=AGRIF_ppm) 217 CALL Agrif_Set_bc( e2u_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 218 219 CALL Agrif_Set_bcinterp(e2v_id,interp1=AGRIF_ppm, interp2=Agrif_linear) 220 CALL Agrif_Set_interp(e2v_id,interp1=AGRIF_ppm, interp2=Agrif_linear) 221 CALL Agrif_Set_bc( e2v_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 222 223 CALL Agrif_Set_bcinterp(e2f_id,interp=AGRIF_linear) 224 CALL Agrif_Set_interp(e2f_id,interp=AGRIF_linear) 225 CALL Agrif_Set_bc( e2f_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 226 227 CALL Agrif_Set_bcinterp(bathy_id,interp=AGRIF_linear) 228 CALL Agrif_Set_interp(bathy_id,interp=AGRIF_linear) 229 CALL Agrif_Set_bc( bathy_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 230 231 ! Vertical scale factors 232 CALL Agrif_Set_bcinterp(e3t_id,interp=AGRIF_ppm) 233 CALL Agrif_Set_interp(e3t_id,interp=AGRIF_ppm) 234 CALL Agrif_Set_bc( e3t_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 235 CALL Agrif_Set_Updatetype( e3t_id, update = AGRIF_Update_Average) 236 237 CALL Agrif_Set_bcinterp(e3t_copy_id,interp=AGRIF_constant) 238 CALL Agrif_Set_interp(e3t_copy_id,interp=AGRIF_constant) 239 CALL Agrif_Set_bc( e3t_copy_id, (/-npt_copy*irafx-1,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/)) 240 241 CALL Agrif_Set_bcinterp(e3t_connect_id,interp=AGRIF_ppm) 242 CALL Agrif_Set_interp(e3t_connect_id,interp=AGRIF_ppm) 243 CALL Agrif_Set_bc( e3t_connect_id, (/-(npt_copy+npt_connect)*irafx-1,-npt_copy*irafx-2/)) 244 245 CALL Agrif_Set_bcinterp(e3u_id, interp1=Agrif_linear, interp2=AGRIF_ppm) 246 CALL Agrif_Set_interp(e3u_id, interp1=Agrif_linear, interp2=AGRIF_ppm) 247 CALL Agrif_Set_bc( e3u_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 248 CALL Agrif_Set_Updatetype(e3u_id,update1 = Agrif_Update_Copy, update2 = Agrif_Update_Average) 249 250 CALL Agrif_Set_bcinterp(e3v_id,interp1=AGRIF_ppm, interp2=Agrif_linear) 251 CALL Agrif_Set_interp(e3v_id, interp1=AGRIF_ppm, interp2=Agrif_linear) 252 CALL Agrif_Set_bc( e3v_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 253 CALL Agrif_Set_Updatetype(e3v_id,update1 = Agrif_Update_Average, update2 = Agrif_Update_Copy) 254 255 ! Bottom level 256 CALL Agrif_Set_bcinterp(bottom_level_id,interp=AGRIF_constant) 257 CALL Agrif_Set_interp(bottom_level_id,interp=AGRIF_constant) 258 CALL Agrif_Set_bc( bottom_level_id, (/-npt_copy*irafx-1,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/)) 259 CALL Agrif_Set_Updatetype( bottom_level_id, update = AGRIF_Update_Max) 260 261 end subroutine agrif_declare_var 262 263 264 subroutine agrif_init_lonlat() 265 use agrif_profiles 266 use agrif_util 267 external :: init_glamt, init_glamu, init_glamv, init_glamf 268 external :: init_gphit, init_gphiu, init_gphiv, init_gphif 269 270 call Agrif_Init_variable(glamt_id, procname = init_glamt) 271 call Agrif_Init_variable(glamu_id, procname = init_glamu) 272 call Agrif_Init_variable(glamv_id, procname = init_glamv) 273 call Agrif_Init_variable(glamf_id, procname = init_glamf) 274 275 call Agrif_Init_variable(gphit_id, procname = init_gphit) 276 call Agrif_Init_variable(gphiu_id, procname = init_gphiu) 277 call Agrif_Init_variable(gphiv_id, procname = init_gphiv) 278 call Agrif_Init_variable(gphif_id, procname = init_gphif) 279 280 end subroutine agrif_init_lonlat 281 282 subroutine agrif_init_scales() 283 use agrif_profiles 284 use agrif_util 285 external :: init_e1t, init_e1u, init_e1v, init_e1f 286 external :: init_e2t, init_e2u, init_e2v, init_e2f 287 288 call Agrif_Init_variable(e1t_id, procname = init_e1t) 289 call Agrif_Init_variable(e1u_id, procname = init_e1u) 290 call Agrif_Init_variable(e1v_id, procname = init_e1v) 291 call Agrif_Init_variable(e1f_id, procname = init_e1f) 292 293 call Agrif_Init_variable(e2t_id, procname = init_e2t) 294 call Agrif_Init_variable(e2u_id, procname = init_e2u) 295 call Agrif_Init_variable(e2v_id, procname = init_e2v) 296 call Agrif_Init_variable(e2f_id, procname = init_e2f) 297 298 end subroutine agrif_init_scales 299 300 301 302 SUBROUTINE init_glamt( ptab, i1, i2, j1, j2, before, nb,ndir) 303 use dom_oce 71 72 WRITE(*,*) ' ' 73 WRITE(*,*)'Size of the High resolution grid: ',nx,' x ',ny 74 WRITE(*,*) ' ' 75 76 CALL agrif_init_lonlat() 77 ln_perio = .FALSE. 78 IF( jperio == 1 .OR. jperio == 2 .OR. jperio == 4 ) ln_perio=.TRUE. 79 80 WHERE (glamt < -180) glamt = glamt +360. 81 WHERE (glamt > +180) glamt = glamt -360. 82 83 CALL lbc_lnk( 'glamt', glamt, 'T', 1._wp) 84 CALL lbc_lnk( 'gphit', gphit, 'T', 1._wp) 85 86 WHERE (glamu < -180) glamu = glamu +360. 87 WHERE (glamu > +180) glamu = glamu -360. 88 CALL lbc_lnk( 'glamu', glamu, 'U', 1._wp) 89 CALL lbc_lnk( 'gphiu', gphiu, 'U', 1._wp) 90 91 WHERE (glamv < -180) glamv = glamv +360. 92 WHERE (glamv > +180) glamv = glamv -360. 93 CALL lbc_lnk( 'glamv', glamv, 'V', 1._wp) 94 CALL lbc_lnk( 'gphiv', gphiv, 'V', 1._wp) 95 96 WHERE (glamf < -180) glamf = glamf +360. 97 WHERE (glamf > +180) glamf = glamf -360. 98 CALL lbc_lnk( 'glamf', glamf, 'F', 1._wp) 99 CALL lbc_lnk( 'gphif', gphif, 'F', 1._wp) 100 101 ! Correct South and North 102 IF ((nbondj == -1).OR.(nbondj == 2)) THEN 103 glamt(:,1) = glamt(:,2) 104 gphit(:,1) = gphit(:,2) 105 glamu(:,1) = glamu(:,2) 106 gphiu(:,1) = gphiu(:,2) 107 glamv(:,1) = glamv(:,2) 108 gphiv(:,1) = gphiv(:,2) 109 ENDIF 110 IF ((nbondj == 1).OR.(nbondj == 2)) THEN 111 glamt(:,jpj) = glamt(:,jpj-1) 112 gphit(:,jpj) = gphit(:,jpj-1) 113 glamu(:,jpj) = glamu(:,jpj-1) 114 gphiu(:,jpj) = gphiu(:,jpj-1) 115 glamv(:,jpj) = glamv(:,jpj-1) 116 gphiv(:,jpj) = gphiv(:,jpj-1) 117 glamf(:,jpj) = glamf(:,jpj-1) 118 gphif(:,jpj) = gphif(:,jpj-1) 119 ENDIF 120 121 ! Correct West and East 122 IF( jperio /= 1 ) THEN 123 IF((nbondi == -1) .OR. (nbondi == 2) ) THEN 124 glamt(1,:) = glamt(2,:) 125 gphit(1,:) = gphit(2,:) 126 glamu(1,:) = glamu(2,:) 127 gphiu(1,:) = gphiu(2,:) 128 glamv(1,:) = glamv(2,:) 129 gphiv(1,:) = gphiv(2,:) 130 ENDIF 131 IF( (nbondi == 1) .OR. (nbondi == 2) ) THEN 132 glamt(jpi,:) = glamt(jpi-1,:) 133 gphit(jpi,:) = gphit(jpi-1,:) 134 glamu(jpi,:) = glamu(jpi-1,:) 135 gphiu(jpi,:) = gphiu(jpi-1,:) 136 glamv(jpi,:) = glamv(jpi-1,:) 137 gphiv(jpi,:) = gphiv(jpi-1,:) 138 glamf(jpi,:) = glamf(jpi-1,:) 139 gphif(jpi,:) = gphif(jpi-1,:) 140 ENDIF 141 ENDIF 142 143 CALL agrif_init_scales() 144 145 ! Correct South and North 146 IF( (nbondj == -1) .OR. (nbondj == 2) ) THEN 147 e1t(:,1) = e1t(:,2) 148 e2t(:,1) = e2t(:,2) 149 e1u(:,1) = e1u(:,2) 150 e2u(:,1) = e2u(:,2) 151 e1v(:,1) = e1v(:,2) 152 e2v(:,1) = e2v(:,2) 153 ENDIF 154 IF( (nbondj == 1) .OR. (nbondj == 2) ) THEN 155 e1t(:,jpj) = e1t(:,jpj-1) 156 e2t(:,jpj) = e2t(:,jpj-1) 157 e1u(:,jpj) = e1u(:,jpj-1) 158 e2u(:,jpj) = e2u(:,jpj-1) 159 e1v(:,jpj) = e1v(:,jpj-1) 160 e2v(:,jpj) = e2v(:,jpj-1) 161 e1f(:,jpj) = e1f(:,jpj-1) 162 e2f(:,jpj) = e2f(:,jpj-1) 163 ENDIF 164 165 ! Correct West and East 166 IF( jperio /= 1 ) THEN 167 IF( (nbondj == -1) .OR. (nbondj == 2) ) THEN 168 e1t(1,:) = e1t(2,:) 169 e2t(1,:) = e2t(2,:) 170 e1u(1,:) = e1u(2,:) 171 e2u(1,:) = e2u(2,:) 172 e1v(1,:) = e1v(2,:) 173 e2v(1,:) = e2v(2,:) 174 ENDIF 175 IF( (nbondj == 1) .OR. (nbondj == 2) ) THEN 176 e1t(jpi,:) = e1t(jpi-1,:) 177 e2t(jpi,:) = e2t(jpi-1,:) 178 e1u(jpi,:) = e1u(jpi-1,:) 179 e2u(jpi,:) = e2u(jpi-1,:) 180 e1v(jpi,:) = e1v(jpi-1,:) 181 e2v(jpi,:) = e2v(jpi-1,:) 182 e1f(jpi,:) = e1f(jpi-1,:) 183 e2f(jpi,:) = e2f(jpi-1,:) 184 ENDIF 185 ENDIF 186 187 END SUBROUTINE Agrif_InitValues_cont 188 189 190 SUBROUTINE agrif_declare_var() 191 !!---------------------------------------------------------------------- 192 !! *** ROUTINE Agrif_InitValues_cont *** 193 !! 194 !! ** Purpose :: Declaration of variables to be interpolated 195 !!---------------------------------------------------------------------- 196 USE par_oce 197 USE dom_oce 198 USE agrif_profiles 199 USE agrif_parameters 200 201 IMPLICIT NONE 202 203 INTEGER :: ind1, ind2, ind3 204 INTEGER :: nx, ny 205 INTEGER ::nbghostcellsfine_tot_x, nbghostcellsfine_tot_y 206 INTEGER :: irafx 207 208 EXTERNAL :: nemo_mapping 209 210 ! 1. Declaration of the type of variable which have to be interpolated 211 !--------------------------------------------------------------------- 212 213 nx=nlci ; ny=nlcj 214 215 ind2 = 2 + nbghostcells_x 216 ind3 = 2 + nbghostcells_y_s 217 nbghostcellsfine_tot_x=nbghostcells_x+1 218 nbghostcellsfine_tot_y=max(nbghostcells_y_s,nbghostcells_y_n)+1 219 220 irafx = Agrif_irhox() 221 222 ! In case of East-West periodicity, prevent AGRIF interpolation at east and west boundaries 223 ! The procnames will not be CALLed at these boundaries 224 if (jperio == 1) THEN 225 CALL Agrif_Set_NearCommonBorderX(.TRUE.) 226 CALL Agrif_Set_DistantCommonBorderX(.TRUE.) 227 endif 228 if (.not.lk_south) THEN 229 CALL Agrif_Set_NearCommonBorderY(.TRUE.) 230 endif 231 232 CALL agrif_declare_variable((/2,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),glamt_id) 233 CALL agrif_declare_variable((/1,2/),(/ind2-1,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),glamu_id) 234 CALL agrif_declare_variable((/2,1/),(/ind2,ind3-1/),(/'x','y'/),(/1,1/),(/nx,ny/),glamv_id) 235 CALL agrif_declare_variable((/1,1/),(/ind2-1,ind3-1/),(/'x','y'/),(/1,1/),(/nx,ny/),glamf_id) 236 237 CALL agrif_declare_variable((/2,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),gphit_id) 238 CALL agrif_declare_variable((/1,2/),(/ind2-1,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),gphiu_id) 239 CALL agrif_declare_variable((/2,1/),(/ind2,ind3-1/),(/'x','y'/),(/1,1/),(/nx,ny/),gphiv_id) 240 CALL agrif_declare_variable((/1,1/),(/ind2-1,ind3-1/),(/'x','y'/),(/1,1/),(/nx,ny/),gphif_id) 241 242 ! Horizontal scale factors 243 244 CALL agrif_declare_variable((/2,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),e1t_id) 245 CALL agrif_declare_variable((/1,2/),(/ind2-1,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),e1u_id) 246 CALL agrif_declare_variable((/2,1/),(/ind2,ind3-1/),(/'x','y'/),(/1,1/),(/nx,ny/),e1v_id) 247 CALL agrif_declare_variable((/1,1/),(/ind2-1,ind3-1/),(/'x','y'/),(/1,1/),(/nx,ny/),e1f_id) 248 249 CALL agrif_declare_variable((/2,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),e2t_id) 250 CALL agrif_declare_variable((/1,2/),(/ind2-1,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),e2u_id) 251 CALL agrif_declare_variable((/2,1/),(/ind2,ind3-1/),(/'x','y'/),(/1,1/),(/nx,ny/),e2v_id) 252 CALL agrif_declare_variable((/1,1/),(/ind2-1,ind3-1/),(/'x','y'/),(/1,1/),(/nx,ny/),e2f_id) 253 254 ! Bathymetry 255 256 CALL agrif_declare_variable((/2,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),bathy_id) 257 258 ! Vertical scale factors 259 CALL agrif_declare_variable((/2,2,0/),(/ind2,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nx,ny,jpk/),e3t_id) 260 CALL agrif_declare_variable((/2,2,0/),(/ind2,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nx,ny,jpk/),e3t_copy_id) 261 CALL agrif_declare_variable((/2,2,0/),(/ind2,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nx,ny,jpk+1/),e3t_connect_id) 262 263 CALL agrif_declare_variable((/1,2,0/),(/ind2-1,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nx,ny,jpk/),e3u_id) 264 CALL agrif_declare_variable((/2,1,0/),(/ind2,ind3-1,0/),(/'x','y','N'/),(/1,1,1/),(/nx,ny,jpk/),e3v_id) 265 266 ! Bottom level 267 268 CALL agrif_declare_variable((/2,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nx,ny/),bottom_level_id) 269 270 CALL Agrif_Set_bcinterp(glamt_id,interp=AGRIF_linear) 271 CALL Agrif_Set_interp(glamt_id,interp=AGRIF_linear) 272 CALL Agrif_Set_bc( glamt_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 273 274 CALL Agrif_Set_bcinterp(glamu_id,interp=AGRIF_linear) 275 CALL Agrif_Set_interp(glamu_id,interp=AGRIF_linear) 276 CALL Agrif_Set_bc( glamu_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)/) ) 277 278 CALL Agrif_Set_bcinterp(glamv_id,interp=AGRIF_linear) 279 CALL Agrif_Set_interp(glamv_id,interp=AGRIF_linear) 280 CALL Agrif_Set_bc( glamv_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)/) ) 281 282 CALL Agrif_Set_bcinterp(glamf_id,interp=AGRIF_linear) 283 CALL Agrif_Set_interp(glamf_id,interp=AGRIF_linear) 284 CALL Agrif_Set_bc( glamf_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)/) ) 285 286 CALL Agrif_Set_bcinterp(gphit_id,interp=AGRIF_linear) 287 CALL Agrif_Set_interp(gphit_id,interp=AGRIF_linear) 288 CALL Agrif_Set_bc( gphit_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 289 290 CALL Agrif_Set_bcinterp(gphiu_id,interp=AGRIF_linear) 291 CALL Agrif_Set_interp(gphiu_id,interp=AGRIF_linear) 292 CALL Agrif_Set_bc( gphiu_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)/) ) 293 294 CALL Agrif_Set_bcinterp(gphiv_id,interp=AGRIF_linear) 295 CALL Agrif_Set_interp(gphiv_id,interp=AGRIF_linear) 296 CALL Agrif_Set_bc( gphiv_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)/) ) 297 298 CALL Agrif_Set_bcinterp(gphif_id,interp=AGRIF_linear) 299 CALL Agrif_Set_interp(gphif_id,interp=AGRIF_linear) 300 CALL Agrif_Set_bc( gphif_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)/) ) 301 302 ! 303 304 CALL Agrif_Set_bcinterp(e1t_id,interp=AGRIF_ppm) 305 CALL Agrif_Set_interp(e1t_id,interp=AGRIF_ppm) 306 CALL Agrif_Set_bc( e1t_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 307 308 CALL Agrif_Set_bcinterp(e1u_id, interp1=Agrif_linear, interp2=AGRIF_ppm) 309 CALL Agrif_Set_interp(e1u_id, interp1=Agrif_linear, interp2=AGRIF_ppm) 310 CALL Agrif_Set_bc( e1u_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)/) ) 311 312 CALL Agrif_Set_bcinterp(e1v_id,interp1=AGRIF_ppm, interp2=Agrif_linear) 313 CALL Agrif_Set_interp(e1v_id, interp1=AGRIF_ppm, interp2=Agrif_linear) 314 CALL Agrif_Set_bc( e1v_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)/) ) 315 316 CALL Agrif_Set_bcinterp(e1f_id,interp=AGRIF_linear) 317 CALL Agrif_Set_interp(e1f_id,interp=AGRIF_linear) 318 CALL Agrif_Set_bc( e1f_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)/) ) 319 320 CALL Agrif_Set_bcinterp(e2t_id,interp=AGRIF_ppm) 321 CALL Agrif_Set_interp(e2t_id,interp=AGRIF_ppm) 322 CALL Agrif_Set_bc( e2t_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 323 324 CALL Agrif_Set_bcinterp(e2u_id,interp1=Agrif_linear, interp2=AGRIF_ppm) 325 CALL Agrif_Set_interp(e2u_id,interp1=Agrif_linear, interp2=AGRIF_ppm) 326 CALL Agrif_Set_bc( e2u_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)/) ) 327 328 CALL Agrif_Set_bcinterp(e2v_id,interp1=AGRIF_ppm, interp2=Agrif_linear) 329 CALL Agrif_Set_interp(e2v_id,interp1=AGRIF_ppm, interp2=Agrif_linear) 330 CALL Agrif_Set_bc( e2v_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)/) ) 331 332 CALL Agrif_Set_bcinterp(e2f_id,interp=AGRIF_linear) 333 CALL Agrif_Set_interp(e2f_id,interp=AGRIF_linear) 334 CALL Agrif_Set_bc( e2f_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)/) ) 335 336 CALL Agrif_Set_bcinterp(bathy_id,interp=AGRIF_linear) 337 CALL Agrif_Set_interp(bathy_id,interp=AGRIF_linear) 338 CALL Agrif_Set_bc( bathy_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 339 340 ! Vertical scale factors 341 CALL Agrif_Set_bcinterp(e3t_id,interp=AGRIF_ppm) 342 CALL Agrif_Set_interp(e3t_id,interp=AGRIF_ppm) 343 CALL Agrif_Set_bc( e3t_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 344 CALL Agrif_Set_Updatetype( e3t_id, update = AGRIF_Update_Average) 345 346 CALL Agrif_Set_bcinterp(e3t_copy_id,interp=AGRIF_constant) 347 CALL Agrif_Set_interp(e3t_copy_id,interp=AGRIF_constant) 348 CALL Agrif_Set_bc( e3t_copy_id, (/-npt_copy*irafx,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/)) 349 350 CALL Agrif_Set_bcinterp(e3t_connect_id,interp=AGRIF_ppm) 351 CALL Agrif_Set_interp(e3t_connect_id,interp=AGRIF_ppm) 352 CALL Agrif_Set_bc( e3t_connect_id, (/-(npt_copy+npt_connect)*irafx-1,-npt_copy*irafx/)) 353 354 CALL Agrif_Set_bcinterp(e3u_id, interp1=Agrif_linear, interp2=AGRIF_ppm) 355 CALL Agrif_Set_interp(e3u_id, interp1=Agrif_linear, interp2=AGRIF_ppm) 356 CALL Agrif_Set_bc( e3u_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 357 CALL Agrif_Set_Updatetype(e3u_id,update1 = Agrif_Update_Copy, update2 = Agrif_Update_Average) 358 359 CALL Agrif_Set_bcinterp(e3v_id,interp1=AGRIF_ppm, interp2=Agrif_linear) 360 CALL Agrif_Set_interp(e3v_id, interp1=AGRIF_ppm, interp2=Agrif_linear) 361 CALL Agrif_Set_bc( e3v_id, (/0,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/) ) 362 CALL Agrif_Set_Updatetype(e3v_id,update1 = Agrif_Update_Average, update2 = Agrif_Update_Copy) 363 364 ! Bottom level 365 CALL Agrif_Set_bcinterp(bottom_level_id,interp=AGRIF_constant) 366 CALL Agrif_Set_interp(bottom_level_id,interp=AGRIF_constant) 367 CALL Agrif_Set_bc( bottom_level_id, (/-npt_copy*irafx,max(nbghostcellsfine_tot_x,nbghostcellsfine_tot_y)-1/)) 368 CALL Agrif_Set_Updatetype( bottom_level_id, update = AGRIF_Update_Max) 369 370 CALL Agrif_Set_ExternalMapping(nemo_mapping) 371 372 END SUBROUTINE agrif_declare_var 373 374 SUBROUTINE nemo_mapping(ndim,ptx,pty,bounds,bounds_chunks,correction_required,nb_chunks) 375 USE dom_oce 376 INTEGER :: ndim 377 INTEGER :: ptx, pty 378 INTEGER,DIMENSION(ndim,2,2) :: bounds 379 INTEGER,DIMENSION(:,:,:,:),allocatable :: bounds_chunks 380 LOGICAL,DIMENSION(:),allocatable :: correction_required 381 INTEGER :: nb_chunks 382 INTEGER :: i 383 384 IF (agrif_debug_interp) THEN 385 DO i = 1, ndim 386 print *,'direction = ',i,bounds(i,1,2),bounds(i,2,2) 387 END DO 388 ENDIF 389 390 IF( bounds(2,2,2) > jpjglo ) THEN 391 IF( bounds(2,1,2) <= jpjglo ) THEN 392 nb_chunks = 2 393 ALLOCATE(bounds_chunks(nb_chunks,ndim,2,2)) 394 ALLOCATE(correction_required(nb_chunks)) 395 DO i = 1, nb_chunks 396 bounds_chunks(i,:,:,:) = bounds 397 END DO 398 399 ! FIRST CHUNCK (for j<=jpjglo) 400 401 ! Original indices 402 bounds_chunks(1,1,1,1) = bounds(1,1,2) 403 bounds_chunks(1,1,2,1) = bounds(1,2,2) 404 bounds_chunks(1,2,1,1) = bounds(2,1,2) 405 bounds_chunks(1,2,2,1) = jpjglo 406 407 bounds_chunks(1,1,1,2) = bounds(1,1,2) 408 bounds_chunks(1,1,2,2) = bounds(1,2,2) 409 bounds_chunks(1,2,1,2) = bounds(2,1,2) 410 bounds_chunks(1,2,2,2) = jpjglo 411 412 ! Correction required or not 413 correction_required(1)=.FALSE. 414 415 ! SECOND CHUNCK (for j>jpjglo) 416 417 !Original indices 418 bounds_chunks(2,1,1,1) = bounds(1,1,2) 419 bounds_chunks(2,1,2,1) = bounds(1,2,2) 420 bounds_chunks(2,2,1,1) = jpjglo-2 421 bounds_chunks(2,2,2,1) = bounds(2,2,2) 422 423 ! Where to find them 424 ! We use the relation TAB(ji,jj)=TAB(jpiglo-ji+2,jpjglo-2-(jj-jpjglo)) 425 426 IF (ptx == 2) THEN ! T, V points 427 bounds_chunks(2,1,1,2) = jpiglo-bounds(1,2,2)+2 428 bounds_chunks(2,1,2,2) = jpiglo-bounds(1,1,2)+2 429 ELSE ! U, F points 430 bounds_chunks(2,1,1,2) = jpiglo-bounds(1,2,2)+1 431 bounds_chunks(2,1,2,2) = jpiglo-bounds(1,1,2)+1 432 ENDIF 433 434 IF (pty == 2) THEN ! T, U points 435 bounds_chunks(2,2,1,2) = jpjglo-2-(bounds(2,2,2) -jpjglo) 436 bounds_chunks(2,2,2,2) = jpjglo-2-(jpjglo-2 -jpjglo) 437 ELSE ! V, F points 438 bounds_chunks(2,2,1,2) = jpjglo-3-(bounds(2,2,2) -jpjglo) 439 bounds_chunks(2,2,2,2) = jpjglo-3-(jpjglo-2 -jpjglo) 440 ENDIF 441 442 ! Correction required or not 443 correction_required(2)=.TRUE. 444 445 ELSE 446 447 nb_chunks = 1 448 ALLOCATE(bounds_chunks(nb_chunks,ndim,2,2)) 449 ALLOCATE(correction_required(nb_chunks)) 450 DO i=1,nb_chunks 451 bounds_chunks(i,:,:,:) = bounds 452 END DO 453 454 bounds_chunks(1,1,1,1) = bounds(1,1,2) 455 bounds_chunks(1,1,2,1) = bounds(1,2,2) 456 bounds_chunks(1,2,1,1) = bounds(2,1,2) 457 bounds_chunks(1,2,2,1) = bounds(2,2,2) 458 459 bounds_chunks(1,1,1,2) = jpiglo-bounds(1,2,2)+2 460 bounds_chunks(1,1,2,2) = jpiglo-bounds(1,1,2)+2 461 462 bounds_chunks(1,2,1,2) = jpjglo-2-(bounds(2,2,2)-jpjglo) 463 bounds_chunks(1,2,2,2) = jpjglo-2-(bounds(2,1,2)-jpjglo) 464 465 IF (ptx == 2) THEN ! T, V points 466 bounds_chunks(1,1,1,2) = jpiglo-bounds(1,2,2)+2 467 bounds_chunks(1,1,2,2) = jpiglo-bounds(1,1,2)+2 468 ELSE ! U, F points 469 bounds_chunks(1,1,1,2) = jpiglo-bounds(1,2,2)+1 470 bounds_chunks(1,1,2,2) = jpiglo-bounds(1,1,2)+1 471 ENDIF 472 473 IF (pty == 2) THEN ! T, U points 474 bounds_chunks(1,2,1,2) = jpjglo-2-(bounds(2,2,2) -jpjglo) 475 bounds_chunks(1,2,2,2) = jpjglo-2-(bounds(2,1,2) -jpjglo) 476 ELSE ! V, F points 477 bounds_chunks(1,2,1,2) = jpjglo-3-(bounds(2,2,2) -jpjglo) 478 bounds_chunks(1,2,2,2) = jpjglo-3-(bounds(2,1,2) -jpjglo) 479 ENDIF 480 481 correction_required(1)=.TRUE. 482 483 ENDIF ! bounds(2,1,2) <= jpjglo 484 485 ELSE IF (bounds(1,1,2) < 1) THEN 486 487 IF (bounds(1,2,2) > 0) THEN 488 nb_chunks = 2 489 ALLOCATE(correction_required(nb_chunks)) 490 correction_required=.FALSE. 491 ALLOCATE(bounds_chunks(nb_chunks,ndim,2,2)) 492 DO i=1,nb_chunks 493 bounds_chunks(i,:,:,:) = bounds 494 END DO 495 496 bounds_chunks(1,1,1,2) = bounds(1,1,2)+jpiglo-2 497 bounds_chunks(1,1,2,2) = 1+jpiglo-2 498 499 bounds_chunks(1,1,1,1) = bounds(1,1,2) 500 bounds_chunks(1,1,2,1) = 1 501 502 bounds_chunks(2,1,1,2) = 2 503 bounds_chunks(2,1,2,2) = bounds(1,2,2) 504 505 bounds_chunks(2,1,1,1) = 2 506 bounds_chunks(2,1,2,1) = bounds(1,2,2) 507 ELSE 508 nb_chunks = 1 509 ALLOCATE(correction_required(nb_chunks)) 510 correction_required=.FALSE. 511 ALLOCATE(bounds_chunks(nb_chunks,ndim,2,2)) 512 DO i=1,nb_chunks 513 bounds_chunks(i,:,:,:) = bounds 514 END DO 515 bounds_chunks(1,1,1,2) = bounds(1,1,2)+jpiglo-2 516 bounds_chunks(1,1,2,2) = bounds(1,2,2)+jpiglo-2 517 518 bounds_chunks(1,1,1,1) = bounds(1,1,2) 519 bounds_chunks(1,1,2,1) = bounds(1,2,2) 520 ENDIF 521 522 ELSE 523 524 nb_chunks=1 525 ALLOCATE(correction_required(nb_chunks)) 526 correction_required=.FALSE. 527 ALLOCATE(bounds_chunks(nb_chunks,ndim,2,2)) 528 DO i=1,nb_chunks 529 bounds_chunks(i,:,:,:) = bounds 530 END DO 531 bounds_chunks(1,1,1,2) = bounds(1,1,2) 532 bounds_chunks(1,1,2,2) = bounds(1,2,2) 533 bounds_chunks(1,2,1,2) = bounds(2,1,2) 534 bounds_chunks(1,2,2,2) = bounds(2,2,2) 535 536 bounds_chunks(1,1,1,1) = bounds(1,1,2) 537 bounds_chunks(1,1,2,1) = bounds(1,2,2) 538 bounds_chunks(1,2,1,1) = bounds(2,1,2) 539 bounds_chunks(1,2,2,1) = bounds(2,2,2) 540 541 ENDIF 542 543 END SUBROUTINE nemo_mapping 544 545 FUNCTION agrif_external_switch_index(ptx,pty,i1,isens) 546 USE dom_oce 547 INTEGER :: ptx, pty, i1, isens 548 INTEGER :: agrif_external_switch_index 549 550 IF( isens == 1 ) THEN 551 IF( ptx == 2 ) THEN ! T, V points 552 agrif_external_switch_index = jpiglo-i1+2 553 ELSE ! U, F points 554 agrif_external_switch_index = jpiglo-i1+1 555 ENDIF 556 ELSE IF (isens ==2) THEN 557 IF (pty == 2) THEN ! T, U points 558 agrif_external_switch_index = jpjglo-2-(i1 -jpjglo) 559 ELSE ! V, F points 560 agrif_external_switch_index = jpjglo-3-(i1 -jpjglo) 561 ENDIF 562 ENDIF 563 564 END FUNCTION agrif_external_switch_index 565 566 SUBROUTINE correct_field(tab2d,i1,i2,j1,j2) 567 USE dom_oce 568 INTEGER :: i1,i2,j1,j2 569 REAL,DIMENSION(i1:i2,j1:j2) :: tab2d 570 571 INTEGER :: i,j 572 REAL,DIMENSION(i1:i2,j1:j2) :: tab2dtemp 573 574 tab2dtemp = tab2d 575 576 DO j=j1,j2 577 DO i=i1,i2 578 tab2d(i,j)=tab2dtemp(i2-(i-i1),j2-(j-j1)) 579 END DO 580 ENDDO 581 582 END SUBROUTINE correct_field 583 584 SUBROUTINE agrif_init_lonlat() 585 USE agrif_profiles 586 USE agrif_util 587 USE dom_oce 588 589 EXTERNAL :: init_glamt, init_glamu, init_glamv, init_glamf 590 EXTERNAL :: init_gphit, init_gphiu, init_gphiv, init_gphif 591 EXTERNAL :: longitude_linear_interp 592 593 INTEGER :: ji,jj,i1,i2,j1,j2 594 REAL, DIMENSION(jpi,jpj) :: tab2dtemp 595 INTEGER :: ind2,ind3 596 INTEGER :: irhox, irhoy 597 598 irhox = agrif_irhox() 599 irhoy = agrif_irhoy() 600 CALL Agrif_Set_external_linear_interp(longitude_linear_interp) 601 602 CALL Agrif_Init_variable(glamt_id, procname = init_glamt) 603 CALL Agrif_Init_variable(glamu_id, procname = init_glamu) 604 CALL Agrif_Init_variable(glamv_id, procname = init_glamv) 605 CALL Agrif_Init_variable(glamf_id, procname = init_glamf) 606 CALL Agrif_UnSet_external_linear_interp() 607 608 CALL Agrif_Init_variable(gphit_id, procname = init_gphit) 609 CALL Agrif_Init_variable(gphiu_id, procname = init_gphiu) 610 CALL Agrif_Init_variable(gphiv_id, procname = init_gphiv) 611 CALL Agrif_Init_variable(gphif_id, procname = init_gphif) 612 613 END SUBROUTINE agrif_init_lonlat 614 615 REAL FUNCTION longitude_linear_interp(x1,x2,coeff) 616 REAL :: x1, x2, coeff 617 REAL :: val_interp 618 619 IF( (x1*x2 <= -50*50) ) THEN 620 IF( x1 < 0 ) THEN 621 val_interp = coeff *(x1+360.) + (1.-coeff) *x2 622 ELSE 623 val_interp = coeff *x1 + (1.-coeff) *(x2+360.) 624 ENDIF 625 IF ((val_interp) >=180.) val_interp = val_interp - 360. 626 ELSE 627 val_interp = coeff * x1 + (1.-coeff) * x2 628 ENDIF 629 longitude_linear_interp = val_interp 630 631 END FUNCTION longitude_linear_interp 632 633 SUBROUTINE agrif_init_scales() 634 USE agrif_profiles 635 USE agrif_util 636 USE dom_oce 637 USE lbclnk 638 LOGICAL :: ln_perio 639 INTEGER nx,ny 640 641 EXTERNAL :: init_e1t, init_e1u, init_e1v, init_e1f 642 EXTERNAL :: init_e2t, init_e2u, init_e2v, init_e2f 643 644 nx = nlci ; ny = nlcj 645 646 ln_perio=.FALSE. 647 if( jperio ==1 .OR. jperio==2 .OR. jperio==4) ln_perio=.TRUE. 648 649 CALL Agrif_Init_variable(e1t_id, procname = init_e1t) 650 CALL Agrif_Init_variable(e1u_id, procname = init_e1u) 651 CALL Agrif_Init_variable(e1v_id, procname = init_e1v) 652 CALL Agrif_Init_variable(e1f_id, procname = init_e1f) 653 654 CALL Agrif_Init_variable(e2t_id, procname = init_e2t) 655 CALL Agrif_Init_variable(e2u_id, procname = init_e2u) 656 CALL Agrif_Init_variable(e2v_id, procname = init_e2v) 657 CALL Agrif_Init_variable(e2f_id, procname = init_e2f) 658 659 CALL lbc_lnk( 'e1t', e1t, 'T', 1._wp) 660 CALL lbc_lnk( 'e2t', e2t, 'T', 1._wp) 661 CALL lbc_lnk( 'e1u', e1u, 'U', 1._wp) 662 CALL lbc_lnk( 'e2u', e2u, 'U', 1._wp) 663 CALL lbc_lnk( 'e1v', e1v, 'V', 1._wp) 664 CALL lbc_lnk( 'e2v', e2v, 'V', 1._wp) 665 CALL lbc_lnk( 'e1f', e1f, 'F', 1._wp) 666 CALL lbc_lnk( 'e2f', e2f, 'F', 1._wp) 667 668 END SUBROUTINE agrif_init_scales 669 670 SUBROUTINE init_glamt( ptab, i1, i2, j1, j2, before, nb,ndir) 671 USE dom_oce 304 672 !!---------------------------------------------------------------------- 305 673 !! *** ROUTINE interpsshn *** 306 !!---------------------------------------------------------------------- 674 !!---------------------------------------------------------------------- 675 INTEGER , INTENT(in ) :: i1, i2, j1, j2 676 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 677 LOGICAL , INTENT(in ) :: before 678 INTEGER , INTENT(in ) :: nb , ndir 679 680 ! 681 !!---------------------------------------------------------------------- 682 ! 683 IF( before) THEN 684 ptab(i1:i2,j1:j2) = glamt(i1:i2,j1:j2) 685 ELSE 686 glamt(i1:i2,j1:j2) = ptab(i1:i2,j1:j2) 687 ENDIF 688 ! 689 END SUBROUTINE init_glamt 690 691 SUBROUTINE init_glamu( ptab, i1, i2, j1, j2, before, nb,ndir) 692 USE dom_oce 693 !!---------------------------------------------------------------------- 694 !! *** ROUTINE interpsshn *** 695 !!---------------------------------------------------------------------- 307 696 INTEGER , INTENT(in ) :: i1, i2, j1, j2 308 697 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab … … 311 700 LOGICAL :: western_side, eastern_side,northern_side,southern_side 312 701 ! 313 !!---------------------------------------------------------------------- 314 ! 315 western_side = (nb == 1).AND.(ndir == 1) 316 eastern_side = (nb == 1).AND.(ndir == 2) 317 southern_side = (nb == 2).AND.(ndir == 1) 318 northern_side = (nb == 2).AND.(ndir == 2) 319 IF( before) THEN 320 ptab(i1:i2,j1:j2) = glamt(i1:i2,j1:j2) 321 ELSE 322 glamt(i1:i2,j1:j2)=ptab 323 ENDIF 324 ! 325 END SUBROUTINE init_glamt 326 327 SUBROUTINE init_glamu( ptab, i1, i2, j1, j2, before, nb,ndir) 328 use dom_oce 702 !!---------------------------------------------------------------------- 703 ! 704 IF( before) THEN 705 ptab(i1:i2,j1:j2) = glamu(i1:i2,j1:j2) 706 ELSE 707 glamu(i1:i2,j1:j2) = ptab(i1:i2,j1:j2) 708 ENDIF 709 ! 710 END SUBROUTINE init_glamu 711 712 SUBROUTINE init_glamv( ptab, i1, i2, j1, j2, before, nb,ndir) 713 USE dom_oce 329 714 !!---------------------------------------------------------------------- 330 715 !! *** ROUTINE interpsshn *** 331 !!---------------------------------------------------------------------- 332 INTEGER , INTENT(in ) :: i1, i2, j1, j2 333 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 334 LOGICAL , INTENT(in ) :: before 335 INTEGER , INTENT(in ) :: nb , ndir 336 LOGICAL :: western_side, eastern_side,northern_side,southern_side 337 ! 338 !!---------------------------------------------------------------------- 339 ! 340 western_side = (nb == 1).AND.(ndir == 1) 341 eastern_side = (nb == 1).AND.(ndir == 2) 342 southern_side = (nb == 2).AND.(ndir == 1) 343 northern_side = (nb == 2).AND.(ndir == 2) 344 IF( before) THEN 345 ptab(i1:i2,j1:j2) = glamu(i1:i2,j1:j2) 346 ELSE 347 glamu(i1:i2,j1:j2)=ptab 348 ENDIF 349 ! 350 END SUBROUTINE init_glamu 351 352 SUBROUTINE init_glamv( ptab, i1, i2, j1, j2, before, nb,ndir) 353 use dom_oce 716 !!---------------------------------------------------------------------- 717 INTEGER , INTENT(in ) :: i1, i2, j1, j2 718 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 719 LOGICAL , INTENT(in ) :: before 720 INTEGER , INTENT(in ) :: nb , ndir 721 ! 722 !!---------------------------------------------------------------------- 723 ! 724 IF( before) THEN 725 ptab(i1:i2,j1:j2) = glamv(i1:i2,j1:j2) 726 ELSE 727 glamv(i1:i2,j1:j2) = ptab(i1:i2,j1:j2) 728 ENDIF 729 ! 730 END SUBROUTINE init_glamv 731 732 SUBROUTINE init_glamf( ptab, i1, i2, j1, j2, before, nb,ndir) 733 USE dom_oce 734 !!---------------------------------------------------------------------- 735 !! *** ROUTINE init_glamf *** 736 !!---------------------------------------------------------------------- 737 INTEGER , INTENT(in ) :: i1, i2, j1, j2 738 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 739 LOGICAL , INTENT(in ) :: before 740 INTEGER , INTENT(in ) :: nb , ndir 741 ! 742 !!---------------------------------------------------------------------- 743 ! 744 IF( before) THEN 745 ptab(i1:i2,j1:j2) = glamf(i1:i2,j1:j2) 746 ELSE 747 glamf(i1:i2,j1:j2) = ptab(i1:i2,j1:j2) 748 ENDIF 749 ! 750 END SUBROUTINE init_glamf 751 752 SUBROUTINE init_gphit( ptab, i1, i2, j1, j2, before, nb,ndir) 753 USE dom_oce 754 !!---------------------------------------------------------------------- 755 !! *** ROUTINE init_gphit *** 756 !!---------------------------------------------------------------------- 757 INTEGER , INTENT(in ) :: i1, i2, j1, j2 758 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 759 LOGICAL , INTENT(in ) :: before 760 INTEGER , INTENT(in ) :: nb , ndir 761 ! 762 !!---------------------------------------------------------------------- 763 ! 764 IF( before) THEN 765 ptab(i1:i2,j1:j2) = gphit(i1:i2,j1:j2) 766 ELSE 767 gphit(i1:i2,j1:j2)=ptab 768 ENDIF 769 ! 770 END SUBROUTINE init_gphit 771 772 SUBROUTINE init_gphiu( ptab, i1, i2, j1, j2, before, nb,ndir) 773 USE dom_oce 774 !!---------------------------------------------------------------------- 775 !! *** ROUTINE init_gphiu *** 776 !!---------------------------------------------------------------------- 777 INTEGER , INTENT(in ) :: i1, i2, j1, j2 778 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 779 LOGICAL , INTENT(in ) :: before 780 INTEGER , INTENT(in ) :: nb , ndir 781 ! 782 !!---------------------------------------------------------------------- 783 ! 784 IF( before) THEN 785 ptab(i1:i2,j1:j2) = gphiu(i1:i2,j1:j2) 786 ELSE 787 gphiu(i1:i2,j1:j2)=ptab 788 ENDIF 789 ! 790 END SUBROUTINE init_gphiu 791 792 SUBROUTINE init_gphiv( ptab, i1, i2, j1, j2, before, nb,ndir) 793 USE dom_oce 794 !!---------------------------------------------------------------------- 795 !! *** ROUTINE init_gphiv *** 796 !!---------------------------------------------------------------------- 797 INTEGER , INTENT(in ) :: i1, i2, j1, j2 798 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 799 LOGICAL , INTENT(in ) :: before 800 INTEGER , INTENT(in ) :: nb , ndir 801 ! 802 !!---------------------------------------------------------------------- 803 804 IF( before) THEN 805 ptab(i1:i2,j1:j2) = gphiv(i1:i2,j1:j2) 806 ELSE 807 gphiv(i1:i2,j1:j2)=ptab 808 ENDIF 809 ! 810 END SUBROUTINE init_gphiv 811 812 813 SUBROUTINE init_gphif( ptab, i1, i2, j1, j2, before, nb,ndir) 814 USE dom_oce 815 !!---------------------------------------------------------------------- 816 !! *** ROUTINE init_gphif *** 817 !!---------------------------------------------------------------------- 818 INTEGER , INTENT(in ) :: i1, i2, j1, j2 819 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 820 LOGICAL , INTENT(in ) :: before 821 INTEGER , INTENT(in ) :: nb , ndir 822 ! 823 !!---------------------------------------------------------------------- 824 ! 825 IF( before) THEN 826 ptab(i1:i2,j1:j2) = gphif(i1:i2,j1:j2) 827 ELSE 828 gphif(i1:i2,j1:j2)=ptab 829 ENDIF 830 ! 831 END SUBROUTINE init_gphif 832 833 834 SUBROUTINE init_e1t( ptab, i1, i2, j1, j2, before, nb,ndir) 835 USE dom_oce 836 USE agrif_parameters 837 !!---------------------------------------------------------------------- 838 !! *** ROUTINE init_e1t *** 839 !!---------------------------------------------------------------------- 840 INTEGER , INTENT(in ) :: i1, i2, j1, j2 841 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 842 LOGICAL , INTENT(in ) :: before 843 INTEGER , INTENT(in ) :: nb , ndir 844 ! 845 !!---------------------------------------------------------------------- 846 ! 847 INTEGER :: jj 848 849 IF( before) THEN 850 ! May need to extend at south boundary 851 IF (j1<1) THEN 852 IF (.NOT.agrif_child(lk_south)) THEN 853 IF ((nbondj == -1).OR.(nbondj == 2)) THEN 854 DO jj=1,j2 855 ptab(i1:i2,jj)=e1t(i1:i2,jj) 856 ENDDO 857 DO jj=j1,0 858 ptab(i1:i2,jj)=e1t(i1:i2,1) 859 ENDDO 860 ENDIF 861 ELSE 862 stop "OUT OF BOUNDS" 863 ENDIF 864 ELSE 865 ptab(i1:i2,j1:j2) = e1t(i1:i2,j1:j2) 866 ENDIF 867 ELSE 868 e1t(i1:i2,j1:j2)=ptab/Agrif_rhoy() 869 ENDIF 870 ! 871 END SUBROUTINE init_e1t 872 873 SUBROUTINE init_e1u( ptab, i1, i2, j1, j2, before, nb,ndir) 874 USE dom_oce 875 USE agrif_parameters 876 !!---------------------------------------------------------------------- 877 !! *** ROUTINE init_e1u *** 878 !!---------------------------------------------------------------------- 879 INTEGER , INTENT(in ) :: i1, i2, j1, j2 880 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 881 LOGICAL , INTENT(in ) :: before 882 INTEGER , INTENT(in ) :: nb , ndir 883 ! 884 !!---------------------------------------------------------------------- 885 ! 886 INTEGER :: jj 887 888 IF( before) THEN 889 IF (j1<1) THEN 890 IF (.NOT.agrif_child(lk_south)) THEN 891 IF ((nbondj == -1).OR.(nbondj == 2)) THEN 892 DO jj=1,j2 893 ptab(i1:i2,jj)=e1u(i1:i2,jj) 894 ENDDO 895 DO jj=j1,0 896 ptab(i1:i2,jj)=e1u(i1:i2,1) 897 ENDDO 898 ENDIF 899 ELSE 900 stop "OUT OF BOUNDS" 901 ENDIF 902 ELSE 903 ptab(i1:i2,j1:j2) = e1u(i1:i2,j1:j2) 904 ENDIF 905 ELSE 906 e1u(i1:i2,j1:j2)=ptab/Agrif_rhoy() 907 ENDIF 908 ! 909 END SUBROUTINE init_e1u 910 911 SUBROUTINE init_e1v( ptab, i1, i2, j1, j2, before, nb,ndir) 912 USE dom_oce 913 !!---------------------------------------------------------------------- 914 !! *** ROUTINE init_e1v *** 915 !!---------------------------------------------------------------------- 916 INTEGER , INTENT(in ) :: i1, i2, j1, j2 917 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 918 LOGICAL , INTENT(in ) :: before 919 INTEGER , INTENT(in ) :: nb , ndir 920 ! 921 !!---------------------------------------------------------------------- 922 ! 923 IF( before) THEN 924 ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) 925 ELSE 926 e1v(i1:i2,j1:j2)=ptab/Agrif_rhoy() 927 ENDIF 928 ! 929 END SUBROUTINE init_e1v 930 931 SUBROUTINE init_e1f( ptab, i1, i2, j1, j2, before, nb,ndir) 932 USE dom_oce 933 !!---------------------------------------------------------------------- 934 !! *** ROUTINE init_e1f *** 935 !!---------------------------------------------------------------------- 936 INTEGER , INTENT(in ) :: i1, i2, j1, j2 937 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 938 LOGICAL , INTENT(in ) :: before 939 INTEGER , INTENT(in ) :: nb , ndir 940 ! 941 !!---------------------------------------------------------------------- 942 ! 943 IF( before) THEN 944 ptab(i1:i2,j1:j2) = e1f(i1:i2,j1:j2) 945 ELSE 946 e1f(i1:i2,j1:j2)=ptab/Agrif_rhoy() 947 ENDIF 948 ! 949 END SUBROUTINE init_e1f 950 951 SUBROUTINE init_e2t( ptab, i1, i2, j1, j2, before, nb,ndir) 952 USE dom_oce 953 USE agrif_parameters 954 !!---------------------------------------------------------------------- 955 !! *** ROUTINE init_e2t *** 956 !!---------------------------------------------------------------------- 957 INTEGER , INTENT(in ) :: i1, i2, j1, j2 958 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 959 LOGICAL , INTENT(in ) :: before 960 INTEGER , INTENT(in ) :: nb , ndir 961 ! 962 !!---------------------------------------------------------------------- 963 ! 964 INTEGER :: jj 965 966 IF( before) THEN 967 IF (j1<1) THEN 968 IF (.NOT.agrif_child(lk_south)) THEN 969 IF ((nbondj == -1).OR.(nbondj == 2)) THEN 970 DO jj=1,j2 971 ptab(i1:i2,jj)=e2t(i1:i2,jj) 972 ENDDO 973 DO jj=j1,0 974 ptab(i1:i2,jj)=e2t(i1:i2,1) 975 ENDDO 976 ENDIF 977 ELSE 978 stop "OUT OF BOUNDS" 979 ENDIF 980 ELSE 981 ptab(i1:i2,j1:j2) = e2t(i1:i2,j1:j2) 982 ENDIF 983 ELSE 984 e2t(i1:i2,j1:j2)=ptab/Agrif_rhoy() 985 ENDIF 986 ! 987 END SUBROUTINE init_e2t 988 989 SUBROUTINE init_e2u( ptab, i1, i2, j1, j2, before, nb,ndir) 990 USE dom_oce 991 USE agrif_parameters 354 992 !!---------------------------------------------------------------------- 355 993 !! *** ROUTINE interpsshn *** 356 !!---------------------------------------------------------------------- 357 INTEGER , INTENT(in ) :: i1, i2, j1, j2 358 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 359 LOGICAL , INTENT(in ) :: before 360 INTEGER , INTENT(in ) :: nb , ndir 361 LOGICAL :: western_side, eastern_side,northern_side,southern_side 362 ! 363 !!---------------------------------------------------------------------- 364 ! 365 western_side = (nb == 1).AND.(ndir == 1) 366 eastern_side = (nb == 1).AND.(ndir == 2) 367 southern_side = (nb == 2).AND.(ndir == 1) 368 northern_side = (nb == 2).AND.(ndir == 2) 369 IF( before) THEN 370 ptab(i1:i2,j1:j2) = glamv(i1:i2,j1:j2) 371 ELSE 372 glamv(i1:i2,j1:j2)=ptab 373 ENDIF 374 ! 375 END SUBROUTINE init_glamv 376 377 SUBROUTINE init_glamf( ptab, i1, i2, j1, j2, before, nb,ndir) 378 use dom_oce 994 !!---------------------------------------------------------------------- 995 INTEGER , INTENT(in ) :: i1, i2, j1, j2 996 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 997 LOGICAL , INTENT(in ) :: before 998 INTEGER , INTENT(in ) :: nb , ndir 999 ! 1000 !!---------------------------------------------------------------------- 1001 ! 1002 INTEGER :: jj 1003 1004 IF( before) THEN 1005 IF (j1<1) THEN 1006 IF (.NOT.agrif_child(lk_south)) THEN 1007 IF ((nbondj == -1).OR.(nbondj == 2)) THEN 1008 DO jj=1,j2 1009 ptab(i1:i2,jj)=e2u(i1:i2,jj) 1010 ENDDO 1011 DO jj=j1,0 1012 ptab(i1:i2,jj)=e2u(i1:i2,1) 1013 ENDDO 1014 ENDIF 1015 ELSE 1016 stop "OUT OF BOUNDS" 1017 ENDIF 1018 ELSE 1019 ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) 1020 ENDIF 1021 ELSE 1022 e2u(i1:i2,j1:j2)=ptab/Agrif_rhoy() 1023 ENDIF 1024 ! 1025 END SUBROUTINE init_e2u 1026 1027 SUBROUTINE init_e2v( ptab, i1, i2, j1, j2, before, nb,ndir) 1028 USE dom_oce 379 1029 !!---------------------------------------------------------------------- 380 1030 !! *** ROUTINE interpsshn *** 381 !!---------------------------------------------------------------------- 382 INTEGER , INTENT(in ) :: i1, i2, j1, j2 383 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 384 LOGICAL , INTENT(in ) :: before 385 INTEGER , INTENT(in ) :: nb , ndir 386 LOGICAL :: western_side, eastern_side,northern_side,southern_side 387 ! 388 !!---------------------------------------------------------------------- 389 ! 390 western_side = (nb == 1).AND.(ndir == 1) 391 eastern_side = (nb == 1).AND.(ndir == 2) 392 southern_side = (nb == 2).AND.(ndir == 1) 393 northern_side = (nb == 2).AND.(ndir == 2) 394 IF( before) THEN 395 ptab(i1:i2,j1:j2) = glamf(i1:i2,j1:j2) 396 ELSE 397 glamf(i1:i2,j1:j2)=ptab 398 ENDIF 399 ! 400 END SUBROUTINE init_glamf 401 402 SUBROUTINE init_gphit( ptab, i1, i2, j1, j2, before, nb,ndir) 403 use dom_oce 1031 !!---------------------------------------------------------------------- 1032 INTEGER , INTENT(in ) :: i1, i2, j1, j2 1033 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 1034 LOGICAL , INTENT(in ) :: before 1035 INTEGER , INTENT(in ) :: nb , ndir 1036 ! 1037 !!---------------------------------------------------------------------- 1038 IF( before) THEN 1039 ptab(i1:i2,j1:j2) = e2v(i1:i2,j1:j2) 1040 ELSE 1041 e2v(i1:i2,j1:j2)=ptab/Agrif_rhoy() 1042 ENDIF 1043 ! 1044 END SUBROUTINE init_e2v 1045 1046 SUBROUTINE init_e2f( ptab, i1, i2, j1, j2, before, nb,ndir) 1047 USE dom_oce 404 1048 !!---------------------------------------------------------------------- 405 1049 !! *** ROUTINE interpsshn *** 406 !!---------------------------------------------------------------------- 407 INTEGER , INTENT(in ) :: i1, i2, j1, j2 408 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 409 LOGICAL , INTENT(in ) :: before 410 INTEGER , INTENT(in ) :: nb , ndir 411 LOGICAL :: western_side, eastern_side,northern_side,southern_side 412 ! 413 !!---------------------------------------------------------------------- 414 ! 415 western_side = (nb == 1).AND.(ndir == 1) 416 eastern_side = (nb == 1).AND.(ndir == 2) 417 southern_side = (nb == 2).AND.(ndir == 1) 418 northern_side = (nb == 2).AND.(ndir == 2) 419 IF( before) THEN 420 ptab(i1:i2,j1:j2) = gphit(i1:i2,j1:j2) 421 ELSE 422 gphit(i1:i2,j1:j2)=ptab 423 ENDIF 424 ! 425 END SUBROUTINE init_gphit 426 427 SUBROUTINE init_gphiu( ptab, i1, i2, j1, j2, before, nb,ndir) 428 use dom_oce 429 !!---------------------------------------------------------------------- 430 !! *** ROUTINE interpsshn *** 431 !!---------------------------------------------------------------------- 432 INTEGER , INTENT(in ) :: i1, i2, j1, j2 433 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 434 LOGICAL , INTENT(in ) :: before 435 INTEGER , INTENT(in ) :: nb , ndir 436 LOGICAL :: western_side, eastern_side,northern_side,southern_side 437 ! 438 !!---------------------------------------------------------------------- 439 ! 440 western_side = (nb == 1).AND.(ndir == 1) 441 eastern_side = (nb == 1).AND.(ndir == 2) 442 southern_side = (nb == 2).AND.(ndir == 1) 443 northern_side = (nb == 2).AND.(ndir == 2) 444 IF( before) THEN 445 ptab(i1:i2,j1:j2) = gphiu(i1:i2,j1:j2) 446 ELSE 447 gphiu(i1:i2,j1:j2)=ptab 448 ENDIF 449 ! 450 END SUBROUTINE init_gphiu 451 452 SUBROUTINE init_gphiv( ptab, i1, i2, j1, j2, before, nb,ndir) 453 use dom_oce 454 !!---------------------------------------------------------------------- 455 !! *** ROUTINE interpsshn *** 456 !!---------------------------------------------------------------------- 457 INTEGER , INTENT(in ) :: i1, i2, j1, j2 458 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 459 LOGICAL , INTENT(in ) :: before 460 INTEGER , INTENT(in ) :: nb , ndir 461 LOGICAL :: western_side, eastern_side,northern_side,southern_side 462 ! 463 !!---------------------------------------------------------------------- 464 ! 465 western_side = (nb == 1).AND.(ndir == 1) 466 eastern_side = (nb == 1).AND.(ndir == 2) 467 southern_side = (nb == 2).AND.(ndir == 1) 468 northern_side = (nb == 2).AND.(ndir == 2) 469 IF( before) THEN 470 ptab(i1:i2,j1:j2) = gphiv(i1:i2,j1:j2) 471 ELSE 472 gphiv(i1:i2,j1:j2)=ptab 473 ENDIF 474 ! 475 END SUBROUTINE init_gphiv 476 477 478 SUBROUTINE init_gphif( ptab, i1, i2, j1, j2, before, nb,ndir) 479 use dom_oce 480 !!---------------------------------------------------------------------- 481 !! *** ROUTINE interpsshn *** 482 !!---------------------------------------------------------------------- 483 INTEGER , INTENT(in ) :: i1, i2, j1, j2 484 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 485 LOGICAL , INTENT(in ) :: before 486 INTEGER , INTENT(in ) :: nb , ndir 487 LOGICAL :: western_side, eastern_side,northern_side,southern_side 488 ! 489 !!---------------------------------------------------------------------- 490 ! 491 western_side = (nb == 1).AND.(ndir == 1) 492 eastern_side = (nb == 1).AND.(ndir == 2) 493 southern_side = (nb == 2).AND.(ndir == 1) 494 northern_side = (nb == 2).AND.(ndir == 2) 495 IF( before) THEN 496 ptab(i1:i2,j1:j2) = gphif(i1:i2,j1:j2) 497 ELSE 498 gphif(i1:i2,j1:j2)=ptab 499 ENDIF 500 ! 501 END SUBROUTINE init_gphif 502 503 504 SUBROUTINE init_e1t( ptab, i1, i2, j1, j2, before, nb,ndir) 505 use dom_oce 506 !!---------------------------------------------------------------------- 507 !! *** ROUTINE interpsshn *** 508 !!---------------------------------------------------------------------- 509 INTEGER , INTENT(in ) :: i1, i2, j1, j2 510 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 511 LOGICAL , INTENT(in ) :: before 512 INTEGER , INTENT(in ) :: nb , ndir 513 LOGICAL :: western_side, eastern_side,northern_side,southern_side 514 ! 515 !!---------------------------------------------------------------------- 516 ! 517 western_side = (nb == 1).AND.(ndir == 1) 518 eastern_side = (nb == 1).AND.(ndir == 2) 519 southern_side = (nb == 2).AND.(ndir == 1) 520 northern_side = (nb == 2).AND.(ndir == 2) 521 IF( before) THEN 522 ptab(i1:i2,j1:j2) = e1t(i1:i2,j1:j2) 523 ELSE 524 e1t(i1:i2,j1:j2)=ptab/Agrif_rhoy() 525 ENDIF 526 ! 527 END SUBROUTINE init_e1t 528 529 SUBROUTINE init_e1u( ptab, i1, i2, j1, j2, before, nb,ndir) 530 use dom_oce 531 !!---------------------------------------------------------------------- 532 !! *** ROUTINE interpsshn *** 533 !!---------------------------------------------------------------------- 534 INTEGER , INTENT(in ) :: i1, i2, j1, j2 535 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 536 LOGICAL , INTENT(in ) :: before 537 INTEGER , INTENT(in ) :: nb , ndir 538 LOGICAL :: western_side, eastern_side,northern_side,southern_side 539 ! 540 !!---------------------------------------------------------------------- 541 ! 542 western_side = (nb == 1).AND.(ndir == 1) 543 eastern_side = (nb == 1).AND.(ndir == 2) 544 southern_side = (nb == 2).AND.(ndir == 1) 545 northern_side = (nb == 2).AND.(ndir == 2) 546 IF( before) THEN 547 ptab(i1:i2,j1:j2) = e1u(i1:i2,j1:j2) 548 ELSE 549 e1u(i1:i2,j1:j2)=ptab/Agrif_rhoy() 550 ENDIF 551 ! 552 END SUBROUTINE init_e1u 553 554 SUBROUTINE init_e1v( ptab, i1, i2, j1, j2, before, nb,ndir) 555 use dom_oce 556 !!---------------------------------------------------------------------- 557 !! *** ROUTINE interpsshn *** 558 !!---------------------------------------------------------------------- 559 INTEGER , INTENT(in ) :: i1, i2, j1, j2 560 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 561 LOGICAL , INTENT(in ) :: before 562 INTEGER , INTENT(in ) :: nb , ndir 563 LOGICAL :: western_side, eastern_side,northern_side,southern_side 564 ! 565 !!---------------------------------------------------------------------- 566 ! 567 western_side = (nb == 1).AND.(ndir == 1) 568 eastern_side = (nb == 1).AND.(ndir == 2) 569 southern_side = (nb == 2).AND.(ndir == 1) 570 northern_side = (nb == 2).AND.(ndir == 2) 571 IF( before) THEN 572 ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) 573 ELSE 574 e1v(i1:i2,j1:j2)=ptab/Agrif_rhoy() 575 ENDIF 576 ! 577 END SUBROUTINE init_e1v 578 579 SUBROUTINE init_e1f( ptab, i1, i2, j1, j2, before, nb,ndir) 580 use dom_oce 581 !!---------------------------------------------------------------------- 582 !! *** ROUTINE interpsshn *** 583 !!---------------------------------------------------------------------- 584 INTEGER , INTENT(in ) :: i1, i2, j1, j2 585 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 586 LOGICAL , INTENT(in ) :: before 587 INTEGER , INTENT(in ) :: nb , ndir 588 LOGICAL :: western_side, eastern_side,northern_side,southern_side 589 ! 590 !!---------------------------------------------------------------------- 591 ! 592 western_side = (nb == 1).AND.(ndir == 1) 593 eastern_side = (nb == 1).AND.(ndir == 2) 594 southern_side = (nb == 2).AND.(ndir == 1) 595 northern_side = (nb == 2).AND.(ndir == 2) 596 IF( before) THEN 597 ptab(i1:i2,j1:j2) = e1f(i1:i2,j1:j2) 598 ELSE 599 e1f(i1:i2,j1:j2)=ptab/Agrif_rhoy() 600 ENDIF 601 ! 602 END SUBROUTINE init_e1f 603 604 SUBROUTINE init_e2t( ptab, i1, i2, j1, j2, before, nb,ndir) 605 use dom_oce 606 !!---------------------------------------------------------------------- 607 !! *** ROUTINE interpsshn *** 608 !!---------------------------------------------------------------------- 609 INTEGER , INTENT(in ) :: i1, i2, j1, j2 610 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 611 LOGICAL , INTENT(in ) :: before 612 INTEGER , INTENT(in ) :: nb , ndir 613 LOGICAL :: western_side, eastern_side,northern_side,southern_side 614 ! 615 !!---------------------------------------------------------------------- 616 ! 617 western_side = (nb == 1).AND.(ndir == 1) 618 eastern_side = (nb == 1).AND.(ndir == 2) 619 southern_side = (nb == 2).AND.(ndir == 1) 620 northern_side = (nb == 2).AND.(ndir == 2) 621 IF( before) THEN 622 ptab(i1:i2,j1:j2) = e2t(i1:i2,j1:j2) 623 ELSE 624 e2t(i1:i2,j1:j2)=ptab/Agrif_rhoy() 625 ENDIF 626 ! 627 END SUBROUTINE init_e2t 628 629 SUBROUTINE init_e2u( ptab, i1, i2, j1, j2, before, nb,ndir) 630 use dom_oce 631 !!---------------------------------------------------------------------- 632 !! *** ROUTINE interpsshn *** 633 !!---------------------------------------------------------------------- 634 INTEGER , INTENT(in ) :: i1, i2, j1, j2 635 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 636 LOGICAL , INTENT(in ) :: before 637 INTEGER , INTENT(in ) :: nb , ndir 638 LOGICAL :: western_side, eastern_side,northern_side,southern_side 639 ! 640 !!---------------------------------------------------------------------- 641 ! 642 western_side = (nb == 1).AND.(ndir == 1) 643 eastern_side = (nb == 1).AND.(ndir == 2) 644 southern_side = (nb == 2).AND.(ndir == 1) 645 northern_side = (nb == 2).AND.(ndir == 2) 646 IF( before) THEN 647 ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) 648 ELSE 649 e2u(i1:i2,j1:j2)=ptab/Agrif_rhoy() 650 ENDIF 651 ! 652 END SUBROUTINE init_e2u 653 654 SUBROUTINE init_e2v( ptab, i1, i2, j1, j2, before, nb,ndir) 655 use dom_oce 656 !!---------------------------------------------------------------------- 657 !! *** ROUTINE interpsshn *** 658 !!---------------------------------------------------------------------- 659 INTEGER , INTENT(in ) :: i1, i2, j1, j2 660 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 661 LOGICAL , INTENT(in ) :: before 662 INTEGER , INTENT(in ) :: nb , ndir 663 LOGICAL :: western_side, eastern_side,northern_side,southern_side 664 ! 665 !!---------------------------------------------------------------------- 666 ! 667 western_side = (nb == 1).AND.(ndir == 1) 668 eastern_side = (nb == 1).AND.(ndir == 2) 669 southern_side = (nb == 2).AND.(ndir == 1) 670 northern_side = (nb == 2).AND.(ndir == 2) 671 IF( before) THEN 672 ptab(i1:i2,j1:j2) = e2v(i1:i2,j1:j2) 673 ELSE 674 e2v(i1:i2,j1:j2)=ptab/Agrif_rhoy() 675 ENDIF 676 ! 677 END SUBROUTINE init_e2v 678 679 SUBROUTINE init_e2f( ptab, i1, i2, j1, j2, before, nb,ndir) 680 use dom_oce 681 !!---------------------------------------------------------------------- 682 !! *** ROUTINE interpsshn *** 683 !!---------------------------------------------------------------------- 684 INTEGER , INTENT(in ) :: i1, i2, j1, j2 685 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 686 LOGICAL , INTENT(in ) :: before 687 INTEGER , INTENT(in ) :: nb , ndir 688 LOGICAL :: western_side, eastern_side,northern_side,southern_side 689 ! 690 !!---------------------------------------------------------------------- 691 ! 692 western_side = (nb == 1).AND.(ndir == 1) 693 eastern_side = (nb == 1).AND.(ndir == 2) 694 southern_side = (nb == 2).AND.(ndir == 1) 695 northern_side = (nb == 2).AND.(ndir == 2) 1050 !!---------------------------------------------------------------------- 1051 INTEGER , INTENT(in ) :: i1, i2, j1, j2 1052 REAL, DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 1053 LOGICAL , INTENT(in ) :: before 1054 INTEGER , INTENT(in ) :: nb , ndir 1055 ! 1056 !!---------------------------------------------------------------------- 1057 ! 696 1058 IF( before) THEN 697 1059 ptab(i1:i2,j1:j2) = e2f(i1:i2,j1:j2) … … 703 1065 704 1066 705 SUBROUTINE agrif_nemo_init706 USE agrif_parameters707 USE in_out_manager 708 USE lib_mpp 709 710 711 !!712 IMPLICIT NONE 713 714 INTEGER :: ios 715 716 NAMELIST/namagrif/ nn_cln_update,ln_spc_dyn,rn_sponge_tra,rn_sponge_dyn,ln_chk_bathy,npt_connect,npt_copy1067 SUBROUTINE agrif_nemo_init 1068 USE agrif_parameters 1069 USE dom_oce 1070 USE in_out_manager 1071 USE lib_mpp 1072 !! 1073 IMPLICIT NONE 1074 1075 INTEGER :: ios 1076 1077 NAMELIST/namagrif/ nn_cln_update,ln_spc_dyn,rn_sponge_tra,rn_sponge_dyn,ln_chk_bathy,npt_connect, & 1078 & npt_copy 717 1079 718 1080 REWIND( numnam_ref ) ! Namelist namagrif in reference namelist : nesting parameters … … 730 1092 WRITE(numout,*) '~~~~~~~' 731 1093 WRITE(numout,*) ' Namelist namagrif : set nesting parameters' 732 WRITE(numout,*) ' npt_copy = ', npt_copy 733 WRITE(numout,*) ' npt_connect = ', npt_connect 734 ENDIF 735 736 END SUBROUTINE agrif_nemo_init 737 738 739 SUBROUTINE Agrif_detect( kg, ksizex ) 1094 WRITE(numout,*) ' npt_copy = ', npt_copy 1095 WRITE(numout,*) ' npt_connect = ', npt_connect 1096 ENDIF 1097 1098 ! Set the number of ghost cells according to periodicity 1099 1100 nbghostcells_x = nbghostcells 1101 nbghostcells_y_s = nbghostcells 1102 nbghostcells_y_n = nbghostcells 1103 1104 lk_west = .NOT. ( Agrif_Ix() == 1 ) 1105 lk_east = .NOT. ( Agrif_Ix() + nbcellsx/AGRIF_Irhox() == Agrif_Parent(jpiglo) -1 ) 1106 lk_south = .NOT. ( Agrif_Iy() == 1 ) 1107 lk_north = .NOT. ( Agrif_Iy() + nbcellsy/AGRIF_Irhoy() == Agrif_Parent(jpjglo) -1 ) 1108 1109 IF (.not.agrif_root()) THEN 1110 IF (jperio == 1) THEN 1111 nbghostcells_x = 0 1112 ENDIF 1113 IF (.NOT.lk_south) THEN 1114 nbghostcells_y_s = 0 1115 ENDIF 1116 ENDIF 1117 1118 END SUBROUTINE agrif_nemo_init 1119 1120 SUBROUTINE Agrif_detect( kg, ksizex ) 740 1121 !!---------------------------------------------------------------------- 741 1122 !! *** ROUTINE Agrif_detect *** 742 1123 !!---------------------------------------------------------------------- 743 INTEGER, DIMENSION(2) :: ksizex 744 INTEGER, DIMENSION(ksizex(1),ksizex(2)) :: kg 745 !!---------------------------------------------------------------------- 746 ! 747 RETURN 748 ! 749 END SUBROUTINE Agrif_detect 750 SUBROUTINE agrif_before_regridding 751 END SUBROUTINE agrif_before_regridding 752 753 SUBROUTINE Agrif_InvLoc( indloc, nprocloc, i, indglob ) 754 !!---------------------------------------------------------------------- 755 !! *** ROUTINE Agrif_InvLoc *** 756 !!---------------------------------------------------------------------- 757 USE dom_oce 758 !! 759 IMPLICIT NONE 760 ! 761 INTEGER :: indglob, indloc, nprocloc, i 762 !!---------------------------------------------------------------------- 763 ! 764 SELECT CASE( i ) 765 CASE(1) ; indglob = indloc + nimppt(nprocloc+1) - 1 766 CASE(2) ; indglob = indloc + njmppt(nprocloc+1) - 1 767 CASE DEFAULT 768 indglob = indloc 769 END SELECT 770 ! 771 END SUBROUTINE Agrif_InvLoc 772 773 SUBROUTINE Agrif_get_proc_info( imin, imax, jmin, jmax ) 1124 INTEGER, DIMENSION(2) :: ksizex 1125 INTEGER, DIMENSION(ksizex(1),ksizex(2)) :: kg 1126 !!---------------------------------------------------------------------- 1127 ! 1128 RETURN 1129 ! 1130 END SUBROUTINE Agrif_detect 1131 1132 SUBROUTINE agrif_before_regridding 1133 END SUBROUTINE agrif_before_regridding 1134 1135 # if defined key_mpp_mpi 1136 SUBROUTINE Agrif_InvLoc( indloc, nprocloc, i, indglob ) 1137 !!---------------------------------------------------------------------- 1138 !! *** ROUTINE Agrif_InvLoc *** 1139 !!---------------------------------------------------------------------- 1140 USE dom_oce 1141 !! 1142 IMPLICIT NONE 1143 ! 1144 INTEGER :: indglob, indloc, nprocloc, i 1145 !!---------------------------------------------------------------------- 1146 ! 1147 SELECT CASE( i ) 1148 CASE(1) ; indglob = indloc + nimppt(nprocloc+1) - 1 1149 CASE(2) ; indglob = indloc + njmppt(nprocloc+1) - 1 1150 CASE DEFAULT 1151 indglob = indloc 1152 END SELECT 1153 ! 1154 END SUBROUTINE Agrif_InvLoc 1155 1156 SUBROUTINE Agrif_get_proc_info( imin, imax, jmin, jmax ) 774 1157 !!---------------------------------------------------------------------- 775 1158 !! *** ROUTINE Agrif_get_proc_info *** 776 1159 !!---------------------------------------------------------------------- 777 USE par_oce778 USE dom_oce779 !!780 IMPLICIT NONE781 !782 INTEGER, INTENT(out) :: imin, imax783 INTEGER, INTENT(out) :: jmin, jmax784 !!----------------------------------------------------------------------785 !786 imin = nimppt(Agrif_Procrank+1) ! ?????787 jmin = njmppt(Agrif_Procrank+1) ! ?????788 imax = imin + jpi - 1789 jmax = jmin + jpj - 1790 !791 END SUBROUTINE Agrif_get_proc_info792 793 SUBROUTINE Agrif_estimate_parallel_cost(imin, imax,jmin, jmax, nbprocs, grid_cost)1160 USE par_oce 1161 USE dom_oce 1162 !! 1163 IMPLICIT NONE 1164 ! 1165 INTEGER, INTENT(out) :: imin, imax 1166 INTEGER, INTENT(out) :: jmin, jmax 1167 !!---------------------------------------------------------------------- 1168 ! 1169 imin = nimppt(Agrif_Procrank+1) ! ????? 1170 jmin = njmppt(Agrif_Procrank+1) ! ????? 1171 imax = imin + jpi - 1 1172 jmax = jmin + jpj - 1 1173 ! 1174 END SUBROUTINE Agrif_get_proc_info 1175 1176 SUBROUTINE Agrif_estimate_parallel_cost(imin, imax,jmin, jmax, nbprocs, grid_cost) 794 1177 !!---------------------------------------------------------------------- 795 1178 !! *** ROUTINE Agrif_estimate_parallel_cost *** 796 1179 !!---------------------------------------------------------------------- 797 USE par_oce798 !!799 IMPLICIT NONE800 !801 INTEGER, INTENT(in) :: imin, imax802 INTEGER, INTENT(in) :: jmin, jmax803 INTEGER, INTENT(in) :: nbprocs804 REAL(wp), INTENT(out) :: grid_cost805 !!----------------------------------------------------------------------806 !807 grid_cost = REAL(imax-imin+1,wp)*REAL(jmax-jmin+1,wp) / REAL(nbprocs,wp)808 !809 END SUBROUTINE Agrif_estimate_parallel_cost810 1180 USE par_oce 1181 !! 1182 IMPLICIT NONE 1183 ! 1184 INTEGER, INTENT(in) :: imin, imax 1185 INTEGER, INTENT(in) :: jmin, jmax 1186 INTEGER, INTENT(in) :: nbprocs 1187 REAL(wp), INTENT(out) :: grid_cost 1188 !!---------------------------------------------------------------------- 1189 ! 1190 grid_cost = REAL(imax-imin+1,wp)*REAL(jmax-jmin+1,wp) / REAL(nbprocs,wp) 1191 ! 1192 END SUBROUTINE Agrif_estimate_parallel_cost 1193 # endif 811 1194 #endif
Note: See TracChangeset
for help on using the changeset viewer.