- Timestamp:
- 2011-09-15T08:41:58+02:00 (13 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2011/dev_r2802_MERCATOR9_floats/NEMOGCM/NEMO/OPA_SRC/FLO/flodom.F90
r2528 r2839 43 43 !! the longitude (degree) and the depth (m). 44 44 !!---------------------------------------------------------------------- 45 LOGICAL :: llinmesh 46 INTEGER :: ji, jj, jk ! DO loop index on 3 directions 47 INTEGER :: jfl, jfl1 ! number of floats 48 INTEGER :: inum ! logical unit for file read 45 CHARACTER (len=21) :: clname ! floats initialisation filename 46 LOGICAL :: llinmesh 47 INTEGER :: ji, jj, jk ! DO loop index on 3 directions 48 INTEGER :: jfl, jfl1 ! number of floats 49 INTEGER :: inum ! logical unit for file read 50 INTEGER :: jtrash ! trash var for reading 51 INTEGER :: ierr 49 52 INTEGER, DIMENSION(jpnfl) :: iimfl, ijmfl, ikmfl ! index mesh of floats 50 53 INTEGER, DIMENSION(jpnfl) :: idomfl, ivtest, ihtest ! - - … … 66 69 67 70 ! read of the restart file 68 READ(inum )( tpifl (jfl), jfl=1, jpnrstflo), &71 READ(inum,*) ( tpifl (jfl), jfl=1, jpnrstflo), & 69 72 ( tpjfl (jfl), jfl=1, jpnrstflo), & 70 73 ( tpkfl (jfl), jfl=1, jpnrstflo), & … … 208 211 ENDIF 209 212 ELSE 210 IF(lwp) WRITE(numout,*) ' init_float read ' 213 214 IF( ln_ariane )THEN 215 216 IF(lwp) WRITE(numout,*) ' init_float read with ariane convention (mesh indexes)' 217 218 ! First initialisation of floats with ariane convention 219 ! 220 ! The indexes are read directly from file (warning ariane 221 ! convention, are refered to 222 ! U,V,W grids - and not T-) 223 ! The isobar advection is managed with the sign of tpkfl ( >0 -> 3D 224 ! advection, <0 -> 2D) 225 ! Some variables are not read, as - gl : time index; 4th 226 ! column 227 ! - transport : transport ; 5th 228 ! column 229 ! and paste in the jtrash var 230 ! At the end, ones need to replace the indexes on T grid 231 ! RMQ : there is no float groups identification ! 232 233 clname='init_float_ariane' 234 235 nisobfl = 1 ! we assume that by default we want 3D advection 236 237 ! we check that the number of floats in the init_file are consistant 238 ! with the namelist 239 IF( lwp ) THEN 240 jfl1=0 241 OPEN( unit=inum, file=clname,status='old',access='sequential',form='formatted') 242 DO WHILE (ierr .GE. 0) 243 jfl1=jfl1+1 244 READ (inum,*, iostat=ierr) 245 END DO 246 CLOSE(inum) 247 IF( (jfl1-1) .NE. jpnfl )THEN 248 WRITE (numout,*) ' STOP the number of floats in' ,clname,' = ',jfl1 249 WRITE (numout,*) ' is not equal to jfl= ',jpnfl 250 STOP 251 ENDIF 252 ENDIF 253 254 ! we get the init values 255 CALL ctl_opn( inum, clname, 'OLD', 'FORMATTED', 'SEQUENTIAL', & 256 & 1, numout, .TRUE., 1 ) 257 DO jfl = 1, jpnfl 258 READ(inum,*) tpifl(jfl),tpjfl(jfl),tpkfl(jfl),jtrash, jtrash 259 if(lwp)write(numout,*)"read : ",tpifl(jfl),tpjfl(jfl),tpkfl(jfl),jtrash, jtrash ; call flush(numout) 260 261 IF ( tpkfl(jfl) .LT. 0. ) nisobfl(jfl) = 0 !set the 2D advection according to init_float 262 ngrpfl(jfl)=jfl 263 END DO 264 265 ! conversion from ariane index to T grid index 266 tpkfl = abs(tpkfl)-0.5 ! reversed vertical axis 267 tpifl = tpifl+0.5 268 tpjfl = tpjfl+0.5 269 270 ! verif of non land point initialisation : no need if correct init 271 272 ELSE 273 IF(lwp) WRITE(numout,*) ' init_float read ' 211 274 212 ! First initialisation of floats213 ! the initials positions of floats are written in a file214 ! with a variable to know if it is a isobar float a number215 ! to identified who want the trajectories of this float and216 ! an index for the number of the float217 ! open the init file218 CALL ctl_opn( inum, 'init_float', 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp )219 READ(inum) (flxx(jfl) , jfl=1, jpnfl), &220 (flyy(jfl) , jfl=1, jpnfl), &221 (flzz(jfl) , jfl=1, jpnfl), &222 (nisobfl(jfl), jfl=1, jpnfl), &223 (ngrpfl(jfl) , jfl=1, jpnfl)224 CLOSE(inum)225 226 ! Test to find the grid point coordonate with the geographical position227 DO jfl = 1, jpnfl228 ihtest(jfl) = 0229 ivtest(jfl) = 0230 ikmfl(jfl) = 0275 ! First initialisation of floats 276 ! the initials positions of floats are written in a file 277 ! with a variable to know if it is a isobar float a number 278 ! to identified who want the trajectories of this float and 279 ! an index for the number of the float 280 ! open the init file 281 CALL ctl_opn( inum, 'init_float', 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp ) 282 READ(inum,*) (flxx(jfl) , jfl=1, jpnfl), & 283 (flyy(jfl) , jfl=1, jpnfl), & 284 (flzz(jfl) , jfl=1, jpnfl), & 285 (nisobfl(jfl), jfl=1, jpnfl), & 286 (ngrpfl(jfl) , jfl=1, jpnfl) 287 CLOSE(inum) 288 289 ! Test to find the grid point coordonate with the geographical position 290 DO jfl = 1, jpnfl 291 ihtest(jfl) = 0 292 ivtest(jfl) = 0 293 ikmfl(jfl) = 0 231 294 # if defined key_mpp_mpi 232 DO ji = MAX(nldi,2), nlei233 DO jj = MAX(nldj,2), nlej ! NO vector opt.234 # else 235 DO ji = 2, jpi236 DO jj = 2, jpj ! NO vector opt.295 DO ji = MAX(nldi,2), nlei 296 DO jj = MAX(nldj,2), nlej ! NO vector opt. 297 # else 298 DO ji = 2, jpi 299 DO jj = 2, jpj ! NO vector opt. 237 300 # endif 238 ! for each float we find the indexes of the mesh301 ! for each float we find the indexes of the mesh 239 302 240 CALL findmesh(glamf(ji-1,jj-1),gphif(ji-1,jj-1), & 241 glamf(ji-1,jj ),gphif(ji-1,jj ), & 242 glamf(ji ,jj ),gphif(ji ,jj ), & 243 glamf(ji ,jj-1),gphif(ji ,jj-1), & 244 flxx(jfl) ,flyy(jfl) , & 245 glamt(ji ,jj ),gphit(ji ,jj ), llinmesh) 246 IF(llinmesh) THEN 247 iimfl(jfl) = ji 248 ijmfl(jfl) = jj 249 ihtest(jfl) = ihtest(jfl)+1 250 DO jk = 1, jpk-1 251 IF( (fsdepw(ji,jj,jk) <= flzz(jfl)) .AND. (fsdepw(ji,jj,jk+1) > flzz(jfl)) ) THEN 252 ikmfl(jfl) = jk 253 ivtest(jfl) = ivtest(jfl) + 1 254 ENDIF 255 END DO 256 ENDIF 303 CALL findmesh(glamf(ji-1,jj-1),gphif(ji-1,jj-1), & 304 glamf(ji-1,jj ),gphif(ji-1,jj ), & 305 glamf(ji ,jj ),gphif(ji ,jj ), & 306 glamf(ji ,jj-1),gphif(ji ,jj-1), & 307 flxx(jfl) ,flyy(jfl) , & 308 glamt(ji ,jj ),gphit(ji ,jj ), llinmesh) 309 IF(llinmesh) THEN 310 iimfl(jfl) = ji 311 ijmfl(jfl) = jj 312 ihtest(jfl) = ihtest(jfl)+1 313 DO jk = 1, jpk-1 314 IF( (fsdepw(ji,jj,jk) <= flzz(jfl)) .AND. (fsdepw(ji,jj,jk+1) > flzz(jfl)) ) THEN 315 ikmfl(jfl) = jk 316 ivtest(jfl) = ivtest(jfl) + 1 317 ENDIF 318 END DO 319 ENDIF 320 END DO 257 321 END DO 258 END DO 259 260 ! If the float is in a mesh computed by an other processor we put iimfl=ijmfl=-1 261 IF( ihtest(jfl) == 0 ) THEN 262 iimfl(jfl) = -1 263 ijmfl(jfl) = -1 264 ENDIF 265 END DO 322 323 ! If the float is in a mesh computed by an other processor we put iimfl=ijmfl=-1 324 IF( ihtest(jfl) == 0 ) THEN 325 iimfl(jfl) = -1 326 ijmfl(jfl) = -1 327 ENDIF 328 END DO 266 329 267 ! A zero in the sum of the arrays "ihtest" and "ivtest"268 IF( lk_mpp ) CALL mpp_sum(ihtest,jpnfl) ! sums over the global domain269 IF( lk_mpp ) CALL mpp_sum(ivtest,jpnfl)270 271 DO jfl = 1, jpnfl272 IF( (ihtest(jfl) > 1 ) .OR. ( ivtest(jfl) > 1 )) THEN273 IF(lwp) WRITE(numout,*) 'THE FLOAT',jfl,' IS NOT IN ONLY ONE MESH'274 ENDIF275 IF( ihtest(jfl) == 0 ) THEN276 IF(lwp) WRITE(numout,*)'THE FLOAT',jfl,' IS IN NO MESH'277 ENDIF278 END DO330 ! A zero in the sum of the arrays "ihtest" and "ivtest" 331 IF( lk_mpp ) CALL mpp_sum(ihtest,jpnfl) ! sums over the global domain 332 IF( lk_mpp ) CALL mpp_sum(ivtest,jpnfl) 333 334 DO jfl = 1, jpnfl 335 IF( (ihtest(jfl) > 1 ) .OR. ( ivtest(jfl) > 1 )) THEN 336 IF(lwp) WRITE(numout,*) 'THE FLOAT',jfl,' IS NOT IN ONLY ONE MESH' 337 ENDIF 338 IF( ihtest(jfl) == 0 ) THEN 339 IF(lwp) WRITE(numout,*)'THE FLOAT',jfl,' IS IN NO MESH' 340 ENDIF 341 END DO 279 342 280 ! We compute the distance between the float and the face of the mesh281 DO jfl = 1, jpnfl282 ! Made only if the float is in the domain of the processor283 IF( (iimfl(jfl) >= 0 ) .AND. ( ijmfl(jfl) >= 0 ) ) THEN284 285 ! TEST TO KNOW IF THE FLOAT IS NOT INITIALISED IN THE COAST286 287 idomfl(jfl) = 0288 IF( tmask(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)) == 0. ) idomfl(jfl)=1289 290 ! Computation of the distance between the float291 ! and the faces of the mesh292 ! zdxab293 ! .294 ! B----.---------C295 ! | . |296 ! |<------>flo |297 ! | ^ |298 ! | |.....|....zdyad299 ! | | |300 ! A--------|-----D301 302 zdxab = dstnce(flxx(jfl),flyy(jfl),glamf(iimfl(jfl)-1,ijmfl(jfl)-1),flyy(jfl))303 zdyad = dstnce(flxx(jfl),flyy(jfl),flxx(jfl),gphif(iimfl(jfl)-1,ijmfl(jfl)-1))304 305 ! Translation of this distances (in meter) in indexes306 307 tpifl(jfl) = (iimfl(jfl)-0.5)+zdxab/ e1u(iimfl(jfl)-1,ijmfl(jfl))+(mig(1)-jpizoom)308 tpjfl(jfl) = (ijmfl(jfl)-0.5)+zdyad/ e2v(iimfl(jfl),ijmfl(jfl)-1)+(mjg(1)-jpjzoom)309 tpkfl(jfl) = (fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)+1) - flzz(jfl))*(ikmfl(jfl)) &310 / (fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)+1) - fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl))) &311 + (flzz(jfl) - fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)))*(ikmfl(jfl)+1) &312 / (fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)+1) - fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)))313 ELSE314 tpifl (jfl) = 0.e0315 tpjfl (jfl) = 0.e0316 tpkfl (jfl) = 0.e0317 idomfl(jfl) = 0318 ENDIF319 END DO343 ! We compute the distance between the float and the face of the mesh 344 DO jfl = 1, jpnfl 345 ! Made only if the float is in the domain of the processor 346 IF( (iimfl(jfl) >= 0 ) .AND. ( ijmfl(jfl) >= 0 ) ) THEN 347 348 ! TEST TO KNOW IF THE FLOAT IS NOT INITIALISED IN THE COAST 349 350 idomfl(jfl) = 0 351 IF( tmask(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)) == 0. ) idomfl(jfl)=1 352 353 ! Computation of the distance between the float 354 ! and the faces of the mesh 355 ! zdxab 356 ! . 357 ! B----.---------C 358 ! | . | 359 ! |<------>flo | 360 ! | ^ | 361 ! | |.....|....zdyad 362 ! | | | 363 ! A--------|-----D 364 365 zdxab = dstnce(flxx(jfl),flyy(jfl),glamf(iimfl(jfl)-1,ijmfl(jfl)-1),flyy(jfl)) 366 zdyad = dstnce(flxx(jfl),flyy(jfl),flxx(jfl),gphif(iimfl(jfl)-1,ijmfl(jfl)-1)) 367 368 ! Translation of this distances (in meter) in indexes 369 370 tpifl(jfl) = (iimfl(jfl)-0.5)+zdxab/ e1u(iimfl(jfl)-1,ijmfl(jfl))+(mig(1)-jpizoom) 371 tpjfl(jfl) = (ijmfl(jfl)-0.5)+zdyad/ e2v(iimfl(jfl),ijmfl(jfl)-1)+(mjg(1)-jpjzoom) 372 tpkfl(jfl) = (fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)+1) - flzz(jfl))*(ikmfl(jfl)) & 373 / (fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)+1) - fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl))) & 374 + (flzz(jfl) - fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)))*(ikmfl(jfl)+1) & 375 / (fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)+1) - fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl))) 376 ELSE 377 tpifl (jfl) = 0.e0 378 tpjfl (jfl) = 0.e0 379 tpkfl (jfl) = 0.e0 380 idomfl(jfl) = 0 381 ENDIF 382 END DO 320 383 321 ! The sum of all the arrays tpifl, tpjfl, tpkfl give 3 arrays with the positions of all the floats. 322 IF( lk_mpp ) CALL mpp_sum( tpifl , jpnfl ) ! sums over the global domain 323 IF( lk_mpp ) CALL mpp_sum( tpjfl , jpnfl ) 324 IF( lk_mpp ) CALL mpp_sum( tpkfl , jpnfl ) 325 IF( lk_mpp ) CALL mpp_sum( idomfl, jpnfl ) 384 ! The sum of all the arrays tpifl, tpjfl, tpkfl give 3 arrays with the positions of all the floats. 385 IF( lk_mpp ) CALL mpp_sum( tpifl , jpnfl ) ! sums over the global domain 386 IF( lk_mpp ) CALL mpp_sum( tpjfl , jpnfl ) 387 IF( lk_mpp ) CALL mpp_sum( tpkfl , jpnfl ) 388 IF( lk_mpp ) CALL mpp_sum( idomfl, jpnfl ) 389 ENDIF 390 326 391 ENDIF 327 392
Note: See TracChangeset
for help on using the changeset viewer.
