- Timestamp:
- 2011-11-07T16:03:29+01:00 (12 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2011/dev_MERCATOR_2011_MERGE/NEMOGCM/NEMO/OPA_SRC/FLO/flodom.F90
r2528 r3051 4 4 !! Ocean floats : domain 5 5 !!====================================================================== 6 !! History : OPA ! 1998-07 (Y.Drillet, CLIPPER) Original code 6 !! History : OPA ! 1998-07 (Y.Drillet, CLIPPER) Original code 7 !! NEMO_3.3.1 ! 2011-09 (C.Bricaud,S.Law-Chune Mercator-Ocean): 8 ! add Ariane convention, Comsecitc changes 7 9 !!---------------------------------------------------------------------- 8 10 #if defined key_floats || defined key_esopa … … 10 12 !! 'key_floats' float trajectories 11 13 !!---------------------------------------------------------------------- 12 !! flo_dom : initialization of floats 13 !! findmesh : compute index of position 14 !! dstnce : compute distance between face mesh and floats 14 !! flo_dom : initialization of floats 15 !! add_new_floats : add new floats (long/lat/depth) 16 !! add_new_ariane_floats : add new floats with araine convention (i/j/k) 17 !! findmesh : compute index of position 18 !! dstnce : compute distance between face mesh and floats 15 19 !!---------------------------------------------------------------------- 16 20 USE oce ! ocean dynamics and tracers … … 23 27 PRIVATE 24 28 25 PUBLIC flo_dom ! routine called by floats.F90 29 PUBLIC flo_dom ! routine called by floats.F90 30 PUBLIC flo_dom_alloc ! Routine called in floats.F90 31 32 CHARACTER (len=21) :: clname1 = 'init_float' ! floats initialisation filename 33 CHARACTER (len=21) :: clname2 = 'init_float_ariane' ! ariane floats initialisation filename 34 35 36 INTEGER , ALLOCATABLE, DIMENSION(:) :: iimfl, ijmfl, ikmfl ! index mesh of floats 37 INTEGER , ALLOCATABLE, DIMENSION(:) :: idomfl, ivtest, ihtest ! - 38 REAL(wp), ALLOCATABLE, DIMENSION(:) :: zgifl, zgjfl, zgkfl ! distances in indexes 26 39 27 40 !! * Substitutions … … 43 56 !! the longitude (degree) and the depth (m). 44 57 !!---------------------------------------------------------------------- 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 49 INTEGER, DIMENSION(jpnfl) :: iimfl, ijmfl, ikmfl ! index mesh of floats 50 INTEGER, DIMENSION(jpnfl) :: idomfl, ivtest, ihtest ! - - 51 REAL(wp) :: zdxab, zdyad 52 REAL(wp), DIMENSION(jpnnewflo+1) :: zgifl, zgjfl, zgkfl 58 INTEGER :: jfl ! dummy loop 59 INTEGER :: inum ! logical unit for file read 53 60 !!--------------------------------------------------------------------- 54 61 … … 59 66 IF(lwp) WRITE(numout,*) ' jpnfl = ',jpnfl 60 67 61 IF(ln_rstflo) THEN 68 !-------------------------! 69 ! FLOAT RESTART FILE READ ! 70 !-------------------------! 71 IF( ln_rstflo )THEN 72 62 73 IF(lwp) WRITE(numout,*) ' float restart file read' 63 74 64 75 ! open the restart file 76 !---------------------- 65 77 CALL ctl_opn( inum, 'restart_float', 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp ) 66 78 67 79 ! read of the restart file 68 READ(inum )( tpifl (jfl), jfl=1, jpnrstflo), &80 READ(inum,*) ( tpifl (jfl), jfl=1, jpnrstflo), & 69 81 ( tpjfl (jfl), jfl=1, jpnrstflo), & 70 82 ( tpkfl (jfl), jfl=1, jpnrstflo), & … … 74 86 75 87 ! if we want a surface drift ( like PROVOR floats ) 76 IF( ln_argo ) THEN 77 DO jfl = 1, jpnrstflo 78 nisobfl(jfl) = 0 79 END DO 80 ENDIF 81 82 IF(lwp) WRITE(numout,*)' flo_dom: END of florstlec' 88 IF( ln_argo ) nisobfl(1:jpnrstflo) = 0 83 89 84 90 ! It is possible to add new floats. 85 IF(lwp) WRITE(numout,*)' flo_dom:jpnfl jpnrstflo ',jpnfl,jpnrstflo 86 IF( jpnfl > jpnrstflo ) THEN 87 ! open the init file 88 CALL ctl_opn( inum, 'init_float', 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp ) 89 DO jfl = jpnrstflo+1, jpnfl 90 READ(inum,*) flxx(jfl),flyy(jfl),flzz(jfl), nisobfl(jfl),ngrpfl(jfl),jfl1 91 END DO 92 CLOSE(inum) 93 IF(lwp) WRITE(numout,*)' flodom: END reading init_float file' 91 !--------------------------------- 92 IF( jpnfl > jpnrstflo )THEN 93 94 IF(lwp) WRITE(numout,*) ' add new floats' 95 96 IF( ln_ariane )THEN !Add new floats with ariane convention 97 CALL flo_add_new_ariane_floats(jpnrstflo+1,jpnfl) 98 ELSE !Add new floats with long/lat convention 99 CALL flo_add_new_floats(jpnrstflo+1,jpnfl) 100 ENDIF 101 ENDIF 102 103 !--------------------------------------! 104 ! FLOAT INITILISATION: NO RESTART FILE ! 105 !--------------------------------------! 106 ELSE !ln_rstflo 107 108 IF( ln_ariane )THEN !Add new floats with ariane convention 109 CALL flo_add_new_ariane_floats(1,jpnfl) 110 ELSE !Add new floats with long/lat convention 111 CALL flo_add_new_floats(1,jpnfl) 112 ENDIF 113 114 ENDIF 94 115 95 ! Test to find the grid point coordonate with the geographical position 96 DO jfl = jpnrstflo+1, jpnfl 97 ihtest(jfl) = 0 98 ivtest(jfl) = 0 99 ikmfl(jfl) = 0 116 END SUBROUTINE flo_dom 117 118 SUBROUTINE flo_add_new_floats(kfl_start, kfl_end) 119 !! ------------------------------------------------------------- 120 !! *** SUBROUTINE add_new_arianefloats *** 121 !! 122 !! ** Purpose : 123 !! 124 !! First initialisation of floats 125 !! the initials positions of floats are written in a file 126 !! with a variable to know if it is a isobar float a number 127 !! to identified who want the trajectories of this float and 128 !! an index for the number of the float 129 !! open the init file 130 !! 131 !! ** Method : 132 !!---------------------------------------------------------------------- 133 INTEGER, INTENT(in) :: kfl_start, kfl_end 134 !! 135 INTEGER :: inum ! file unit 136 INTEGER :: jfl,ji, jj, jk ! dummy loop indices 137 INTEGER :: itrash ! trash var for reading 138 INTEGER :: ifl ! number of floats to read 139 REAL(wp) :: zdxab, zdyad 140 LOGICAL :: llinmesh 141 CHARACTER(len=80) :: cltmp 142 !!--------------------------------------------------------------------- 143 ifl = kfl_end-kfl_start+1 144 145 ! we get the init values 146 !----------------------- 147 CALL ctl_opn( inum , clname1, 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp ) 148 DO jfl = kfl_start,kfl_end 149 READ(inum,*) flxx(jfl),flyy(jfl),flzz(jfl), nisobfl(jfl),ngrpfl(jfl),itrash 150 if(lwp)write(numout,*)'read:',jfl,flxx(jfl),flyy(jfl),flzz(jfl), nisobfl(jfl),ngrpfl(jfl),itrash ; call flush(numout) 151 END DO 152 CLOSE(inum) 153 154 ! Test to find the grid point coordonate with the geographical position 155 !---------------------------------------------------------------------- 156 DO jfl = kfl_start,kfl_end 157 ihtest(jfl) = 0 158 ivtest(jfl) = 0 159 ikmfl(jfl) = 0 100 160 # if defined key_mpp_mpi 101 102 103 # else 104 105 161 DO ji = MAX(nldi,2), nlei 162 DO jj = MAX(nldj,2), nlej ! NO vector opt. 163 # else 164 DO ji = 2, jpi 165 DO jj = 2, jpj ! NO vector opt. 106 166 # endif 107 ! For each float we find the indexes of the mesh 108 CALL findmesh(glamf(ji-1,jj-1),gphif(ji-1,jj-1), & 109 glamf(ji-1,jj ),gphif(ji-1,jj ), & 110 glamf(ji ,jj ),gphif(ji ,jj ), & 111 glamf(ji ,jj-1),gphif(ji ,jj-1), & 112 flxx(jfl) ,flyy(jfl) , & 113 glamt(ji ,jj ),gphit(ji ,jj ), llinmesh) 114 IF(llinmesh) THEN 115 iimfl(jfl) = ji 116 ijmfl(jfl) = jj 117 ihtest(jfl) = ihtest(jfl)+1 118 DO jk = 1, jpk-1 119 IF( (fsdepw(ji,jj,jk) <= flzz(jfl)) .AND. (fsdepw(ji,jj,jk+1) > flzz(jfl)) ) THEN 120 ikmfl(jfl) = jk 121 ivtest(jfl) = ivtest(jfl) + 1 122 ENDIF 123 END DO 167 ! For each float we find the indexes of the mesh 168 CALL flo_findmesh(glamf(ji-1,jj-1),gphif(ji-1,jj-1), & 169 glamf(ji-1,jj ),gphif(ji-1,jj ), & 170 glamf(ji ,jj ),gphif(ji ,jj ), & 171 glamf(ji ,jj-1),gphif(ji ,jj-1), & 172 flxx(jfl) ,flyy(jfl) , & 173 glamt(ji ,jj ),gphit(ji ,jj ), llinmesh) 174 IF( llinmesh )THEN 175 iimfl(jfl) = ji 176 ijmfl(jfl) = jj 177 ihtest(jfl) = ihtest(jfl)+1 178 DO jk = 1, jpk-1 179 IF( (fsdepw(ji,jj,jk) <= flzz(jfl)) .AND. (fsdepw(ji,jj,jk+1) > flzz(jfl)) ) THEN 180 ikmfl(jfl) = jk 181 ivtest(jfl) = ivtest(jfl) + 1 124 182 ENDIF 125 183 END DO 126 END DO127 IF(lwp) WRITE(numout,*)' flo_dom: END findmesh'128 129 ! If the float is in a mesh computed by an other processor we put iimfl=ijmfl=-1130 IF( ihtest(jfl) == 0 ) THEN131 iimfl(jfl) = -1132 ijmfl(jfl) = -1133 184 ENDIF 134 185 END DO 186 END DO 187 188 ! If the float is in a mesh computed by an other processor we put iimfl=ijmfl=-1 189 IF( ihtest(jfl) == 0 ) THEN 190 iimfl(jfl) = -1 191 ijmfl(jfl) = -1 192 ENDIF 193 END DO 194 195 !Test if each float is in one and only one proc 196 !---------------------------------------------- 197 IF( lk_mpp ) THEN 198 CALL mpp_sum(ihtest,jpnfl) 199 CALL mpp_sum(ivtest,jpnfl) 200 ENDIF 201 DO jfl = kfl_start,kfl_end 202 203 IF( (ihtest(jfl) > 1 ) .OR. ( ivtest(jfl) > 1) ) THEN 204 WRITE(cltmp,'(A10,i4.4,A20)' )'THE FLOAT',jfl,' IS NOT IN ONLY ONE MESH' 205 CALL ctl_stop('STOP',TRIM(cltmp) ) 206 ENDIF 207 IF( (ihtest(jfl) == 0) ) THEN 208 WRITE(cltmp,'(A10,i4.4,A20)' )'THE FLOAT',jfl,' IS IN NO MESH' 209 CALL ctl_stop('STOP',TRIM(cltmp) ) 210 ENDIF 211 END DO 212 213 ! We compute the distance between the float and the face of the mesh 214 !------------------------------------------------------------------- 215 DO jfl = kfl_start,kfl_end 216 217 ! Made only if the float is in the domain of the processor 218 IF( (iimfl(jfl) >= 0) .AND. (ijmfl(jfl) >= 0) ) THEN 219 220 ! TEST TO KNOW IF THE FLOAT IS NOT INITIALISED IN THE COAST 221 idomfl(jfl) = 0 222 IF( tmask(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)) == 0. ) idomfl(jfl) = 1 223 224 ! Computation of the distance between the float and the faces of the mesh 225 ! zdxab 226 ! . 227 ! B----.---------C 228 ! | . | 229 ! |<------>flo | 230 ! | ^ | 231 ! | |.....|....zdyad 232 ! | | | 233 ! A--------|-----D 234 ! 235 zdxab = flo_dstnce( flxx(jfl), flyy(jfl), glamf(iimfl(jfl)-1,ijmfl(jfl)-1), flyy(jfl) ) 236 zdyad = flo_dstnce( flxx(jfl), flyy(jfl), flxx(jfl), gphif(iimfl(jfl)-1,ijmfl(jfl)-1) ) 237 238 ! Translation of this distances (in meter) in indexes 239 zgifl(jfl)= (iimfl(jfl)-0.5) + zdxab/e1u(iimfl(jfl)-1,ijmfl(jfl)) + (mig(1)-jpizoom) 240 zgjfl(jfl)= (ijmfl(jfl)-0.5) + zdyad/e2v(iimfl(jfl),ijmfl(jfl)-1) + (mjg(1)-jpjzoom) 241 zgkfl(jfl) = (( fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)+1) - flzz(jfl) )* ikmfl(jfl)) & 242 & / ( fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)+1) & 243 & - fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl) ) ) & 244 & + (( flzz(jfl)-fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)) ) *(ikmfl(jfl)+1)) & 245 & / ( fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)+1) & 246 & - fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)) ) 247 ELSE 248 zgifl(jfl) = 0.e0 249 zgjfl(jfl) = 0.e0 250 zgkfl(jfl) = 0.e0 251 ENDIF 252 253 END DO 254 255 ! The sum of all the arrays zgifl, zgjfl, zgkfl give 3 arrays with the positions of all the floats. 256 IF( lk_mpp ) THEN 257 CALL mpp_sum( zgjfl, ifl ) ! sums over the global domain 258 CALL mpp_sum( zgkfl, ifl ) 259 ENDIF 135 260 136 ! A zero in the sum of the arrays "ihtest" and "ivtest" 137 # if defined key_mpp_mpi 138 CALL mpp_sum(ihtest,jpnfl) 139 CALL mpp_sum(ivtest,jpnfl) 140 # endif 141 DO jfl = jpnrstflo+1, jpnfl 142 IF( (ihtest(jfl) > 1 ) .OR. ( ivtest(jfl) > 1) ) THEN 143 IF(lwp) WRITE(numout,*) 'THE FLOAT',jfl,' IS NOT IN ONLY ONE MESH' 144 STOP 145 ENDIF 146 IF( (ihtest(jfl) == 0) ) THEN 147 IF(lwp) WRITE(numout,*)'THE FLOAT',jfl,' IS IN NO MESH' 148 STOP 149 ENDIF 150 END DO 151 152 ! We compute the distance between the float and the face of the mesh 153 DO jfl = jpnrstflo+1, jpnfl 154 ! Made only if the float is in the domain of the processor 155 IF( (iimfl(jfl) >= 0) .AND. (ijmfl(jfl) >= 0) ) THEN 156 157 ! TEST TO KNOW IF THE FLOAT IS NOT INITIALISED IN THE COAST 158 159 idomfl(jfl) = 0 160 IF( tmask(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)) == 0. ) idomfl(jfl) = 1 161 162 ! Computation of the distance between the float and the faces of the mesh 163 ! zdxab 164 ! . 165 ! B----.---------C 166 ! | . | 167 ! |<------>flo | 168 ! | ^ | 169 ! | |.....|....zdyad 170 ! | | | 171 ! A--------|-----D 172 ! 173 174 zdxab = dstnce( flxx(jfl), flyy(jfl), glamf(iimfl(jfl)-1,ijmfl(jfl)-1), flyy(jfl) ) 175 zdyad = dstnce( flxx(jfl), flyy(jfl), flxx(jfl), gphif(iimfl(jfl)-1,ijmfl(jfl)-1) ) 176 177 ! Translation of this distances (in meter) in indexes 178 179 zgifl(jfl-jpnrstflo)= (iimfl(jfl)-0.5) + zdxab/e1u(iimfl(jfl)-1,ijmfl(jfl)) + (mig(1)-jpizoom) 180 zgjfl(jfl-jpnrstflo)= (ijmfl(jfl)-0.5) + zdyad/e2v(iimfl(jfl),ijmfl(jfl)-1) + (mjg(1)-jpjzoom) 181 zgkfl(jfl-jpnrstflo) = (( fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)+1) - flzz(jfl) )* ikmfl(jfl)) & 182 & / ( fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)+1) & 183 & - fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl) ) ) & 184 & + (( flzz(jfl)-fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)) ) *(ikmfl(jfl)+1)) & 185 & / ( fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)+1) & 186 & - fsdepw(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)) ) 187 ELSE 188 zgifl(jfl-jpnrstflo) = 0.e0 189 zgjfl(jfl-jpnrstflo) = 0.e0 190 zgkfl(jfl-jpnrstflo) = 0.e0 191 ENDIF 192 END DO 193 194 ! The sum of all the arrays zgifl, zgjfl, zgkfl give 3 arrays with the positions of all the floats. 195 IF( lk_mpp ) THEN 196 CALL mpp_sum( zgjfl, jpnnewflo ) ! sums over the global domain 197 CALL mpp_sum( zgkfl, jpnnewflo ) 198 IF(lwp) WRITE(numout,*) (zgifl(jfl),jfl=1,jpnnewflo) 199 IF(lwp) WRITE(numout,*) (zgjfl(jfl),jfl=1,jpnnewflo) 200 IF(lwp) WRITE(numout,*) (zgkfl(jfl),jfl=1,jpnnewflo) 201 ENDIF 202 203 DO jfl = jpnrstflo+1, jpnfl 204 tpifl(jfl) = zgifl(jfl-jpnrstflo) 205 tpjfl(jfl) = zgjfl(jfl-jpnrstflo) 206 tpkfl(jfl) = zgkfl(jfl-jpnrstflo) 207 END DO 208 ENDIF 209 ELSE 210 IF(lwp) WRITE(numout,*) ' init_float read ' 211 212 ! First initialisation of floats 213 ! the initials positions of floats are written in a file 214 ! with a variable to know if it is a isobar float a number 215 ! to identified who want the trajectories of this float and 216 ! an index for the number of the float 217 ! open the init file 218 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 position 227 DO jfl = 1, jpnfl 228 ihtest(jfl) = 0 229 ivtest(jfl) = 0 230 ikmfl(jfl) = 0 231 # if defined key_mpp_mpi 232 DO ji = MAX(nldi,2), nlei 233 DO jj = MAX(nldj,2), nlej ! NO vector opt. 234 # else 235 DO ji = 2, jpi 236 DO jj = 2, jpj ! NO vector opt. 237 # endif 238 ! for each float we find the indexes of the mesh 239 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 257 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 266 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 domain 269 IF( lk_mpp ) CALL mpp_sum(ivtest,jpnfl) 270 271 DO jfl = 1, jpnfl 272 IF( (ihtest(jfl) > 1 ) .OR. ( ivtest(jfl) > 1 )) THEN 273 IF(lwp) WRITE(numout,*) 'THE FLOAT',jfl,' IS NOT IN ONLY ONE MESH' 274 ENDIF 275 IF( ihtest(jfl) == 0 ) THEN 276 IF(lwp) WRITE(numout,*)'THE FLOAT',jfl,' IS IN NO MESH' 277 ENDIF 278 END DO 279 280 ! We compute the distance between the float and the face of the mesh 281 DO jfl = 1, jpnfl 282 ! Made only if the float is in the domain of the processor 283 IF( (iimfl(jfl) >= 0 ) .AND. ( ijmfl(jfl) >= 0 ) ) THEN 284 285 ! TEST TO KNOW IF THE FLOAT IS NOT INITIALISED IN THE COAST 286 287 idomfl(jfl) = 0 288 IF( tmask(iimfl(jfl),ijmfl(jfl),ikmfl(jfl)) == 0. ) idomfl(jfl)=1 289 290 ! Computation of the distance between the float 291 ! and the faces of the mesh 292 ! zdxab 293 ! . 294 ! B----.---------C 295 ! | . | 296 ! |<------>flo | 297 ! | ^ | 298 ! | |.....|....zdyad 299 ! | | | 300 ! A--------|-----D 301 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 indexes 306 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 ELSE 314 tpifl (jfl) = 0.e0 315 tpjfl (jfl) = 0.e0 316 tpkfl (jfl) = 0.e0 317 idomfl(jfl) = 0 318 ENDIF 319 END DO 320 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 ) 326 ENDIF 327 328 ! Print the initial positions of the floats 261 DO jfl = kfl_start,kfl_end 262 tpifl(jfl) = zgifl(jfl) 263 tpjfl(jfl) = zgjfl(jfl) 264 tpkfl(jfl) = zgkfl(jfl) 265 END DO 266 267 ! WARNING : initial position not in the sea 329 268 IF( .NOT. ln_rstflo ) THEN 330 ! WARNING : initial position not in the sea 331 DO jfl = 1, jpnfl 269 DO jfl = kfl_start,kfl_end 332 270 IF( idomfl(jfl) == 1 ) THEN 333 271 IF(lwp) WRITE(numout,*)'*****************************' … … 341 279 ENDIF 342 280 343 END SUBROUTINE flo_dom 344 345 346 SUBROUTINE findmesh( pax, pay, pbx, pby, & 347 pcx, pcy, pdx, pdy, & 348 px ,py ,ptx, pty, ldinmesh ) 281 END SUBROUTINE flo_add_new_floats 282 283 SUBROUTINE flo_add_new_ariane_floats(kfl_start, kfl_end) 284 !! ------------------------------------------------------------- 285 !! *** SUBROUTINE add_new_arianefloats *** 286 !! 287 !! ** Purpose : 288 !! First initialisation of floats with ariane convention 289 !! 290 !! The indexes are read directly from file (warning ariane 291 !! convention, are refered to 292 !! U,V,W grids - and not T-) 293 !! The isobar advection is managed with the sign of tpkfl ( >0 -> 3D 294 !! advection, <0 -> 2D) 295 !! Some variables are not read, as - gl : time index; 4th 296 !! column 297 !! - transport : transport ; 5th 298 !! column 299 !! and paste in the jtrash var 300 !! At the end, ones need to replace the indexes on T grid 301 !! RMQ : there is no float groups identification ! 302 !! 303 !! 304 !! ** Method : 305 !!---------------------------------------------------------------------- 306 INTEGER, INTENT(in) :: kfl_start, kfl_end 307 !! 308 INTEGER :: inum ! file unit 309 INTEGER :: ierr, ifl 310 INTEGER :: jfl, jfl1 ! dummy loop indices 311 INTEGER :: itrash ! trash var for reading 312 CHARACTER(len=80) :: cltmp 313 314 !!---------------------------------------------------------------------- 315 nisobfl(kfl_start:kfl_end) = 1 ! we assume that by default we want 3D advection 316 317 ifl = kfl_end - kfl_start + 1 ! number of floats to read 318 319 ! we check that the number of floats in the init_file are consistant with the namelist 320 IF( lwp ) THEN 321 322 jfl1=0 323 ierr=0 324 CALL ctl_opn( inum, clname2, 'OLD', 'FORMATTED', 'SEQUENTIAL', 1, numout, .TRUE., 1 ) 325 DO WHILE (ierr .EQ. 0) 326 jfl1=jfl1+1 327 READ(inum,*, iostat=ierr) 328 END DO 329 CLOSE(inum) 330 IF( (jfl1-1) .NE. ifl )THEN 331 WRITE(cltmp,'(A25,A20,A3,i4.4,A10,i4.4)')"the number of floats in ",TRIM(clname2), & 332 " = ",jfl1," is not equal to jfl= ",ifl 333 CALL ctl_stop('STOP',TRIM(cltmp) ) 334 ENDIF 335 336 ENDIF 337 338 ! we get the init values 339 CALL ctl_opn( inum, clname2, 'OLD', 'FORMATTED', 'SEQUENTIAL', 1, numout, .TRUE., 1 ) 340 DO jfl = kfl_start, kfl_end 341 READ(inum,*) tpifl(jfl),tpjfl(jfl),tpkfl(jfl),itrash, itrash 342 343 IF ( tpkfl(jfl) .LT. 0. ) nisobfl(jfl) = 0 !set the 2D advection according to init_float 344 ngrpfl(jfl)=jfl 345 END DO 346 347 ! conversion from ariane index to T grid index 348 tpkfl(kfl_start:kfl_end) = abs(tpkfl)-0.5 ! reversed vertical axis 349 tpifl(kfl_start:kfl_end) = tpifl+0.5 350 tpjfl(kfl_start:kfl_end) = tpjfl+0.5 351 352 353 END SUBROUTINE flo_add_new_ariane_floats 354 355 356 SUBROUTINE flo_findmesh( pax, pay, pbx, pby, & 357 pcx, pcy, pdx, pdy, & 358 px ,py ,ptx, pty, ldinmesh ) 349 359 !! ------------------------------------------------------------- 350 360 !! *** ROUTINE findmesh *** … … 402 412 ENDIF 403 413 ! 404 END SUBROUTINE f indmesh405 406 407 FUNCTION dstnce( pla1, phi1, pla2, phi2 )414 END SUBROUTINE flo_findmesh 415 416 417 FUNCTION flo_dstnce( pla1, phi1, pla2, phi2 ) 408 418 !! ------------------------------------------------------------- 409 419 !! *** Function dstnce *** … … 415 425 REAL(wp), INTENT(in) :: pla1, phi1, pla2, phi2 ! ??? 416 426 !! 417 REAL(wp) :: 418 REAL(wp) :: 427 REAL(wp) :: dly1, dly2, dlx1, dlx2, dlx, dls, dld, dpi 428 REAL(wp) :: flo_dstnce 419 429 !!--------------------------------------------------------------------- 420 430 ! 421 dpi = 2. * ASIN(1.)422 dls = dpi / 180. 431 dpi = 2._wp * ASIN(1._wp) 432 dls = dpi / 180._wp 423 433 dly1 = phi1 * dls 424 434 dly2 = phi2 * dls … … 428 438 dlx = SIN(dly1) * SIN(dly2) + COS(dly1) * COS(dly2) * COS(dlx2-dlx1) 429 439 ! 430 IF( ABS(dlx) > 1.0 ) dlx = 1.0 431 ! 432 dld = ATAN(DSQRT( ( 1-dlx )/( 1+dlx ) )) * 222.24 / dls 433 dstnce = dld * 1000. 434 ! 435 END FUNCTION dstnce 436 437 # else 440 IF( ABS(dlx) > 1.0_wp ) dlx = 1.0_wp 441 ! 442 dld = ATAN(DSQRT( 1._wp * ( 1._wp-dlx )/( 1._wp+dlx ) )) * 222.24_wp / dls 443 flo_dstnce = dld * 1000._wp 444 ! 445 END FUNCTION flo_dstnce 446 447 INTEGER FUNCTION flo_dom_alloc() 448 !!---------------------------------------------------------------------- 449 !! *** FUNCTION flo_dom_alloc *** 450 !!---------------------------------------------------------------------- 451 452 ALLOCATE( iimfl(jpnfl) , ijmfl(jpnfl) , ikmfl(jpnfl) , & 453 idomfl(jpnfl), ivtest(jpnfl), ihtest(jpnfl), & 454 zgifl(jpnfl) , zgjfl(jpnfl) , zgkfl(jpnfl) , STAT=flo_dom_alloc ) 455 ! 456 IF( lk_mpp ) CALL mpp_sum ( flo_dom_alloc ) 457 IF( flo_dom_alloc /= 0 ) CALL ctl_warn('flo_dom_alloc: failed to allocate arrays') 458 END FUNCTION flo_dom_alloc 459 460 461 #else 438 462 !!---------------------------------------------------------------------- 439 463 !! Default option Empty module … … 441 465 CONTAINS 442 466 SUBROUTINE flo_dom ! Empty routine 467 WRITE(*,*) 'flo_dom: : You should not have seen this print! error?' 443 468 END SUBROUTINE flo_dom 444 469 #endif
Note: See TracChangeset
for help on using the changeset viewer.