Changeset 72 for trunk


Ignore:
Timestamp:
02/12/13 23:04:18 (11 years ago)
Author:
smasson
Message:

put bugfixes from NEMO v3.4.1

Location:
trunk/NEMOGCM/NEMO
Files:
24 edited

Legend:

Unmodified
Added
Removed
  • trunk/NEMOGCM/NEMO/OPA_SRC/DIA/diadct.F90

    r7 r72  
    3838  USE dianam          ! build name of file 
    3939  USE lib_mpp         ! distributed memory computing library 
    40 #if defined key_lim2 || defined key_lim3 
    41   USE ice 
     40#if defined key_lim2 
     41  USE ice_2 
     42#endif 
     43#if defined key_lim3 
     44  USE ice_3 
    4245#endif 
    4346  USE domvvl 
     
    362365              WRITE(numout,*)"      List of points in global domain:" 
    363366              DO jpt=1,iptglo 
    364                  WRITE(numout,*)'        # I J ',jpt,coordtemp(jpt) 
     367                 WRITE(numout,*)'        # I J ',jpt,coordtemp(jpt),directemp(jpt) 
    365368              ENDDO                   
    366369           ENDIF 
     
    403406 
    404407              IF(jsec==nn_secdebug .AND. secs(jsec)%nb_point .NE. 0)THEN 
    405               WRITE(narea+200,*)'avant secs(jsec)%nb_point iptloc ',secs(jsec)%nb_point,iptloc 
    406408              DO jpt = 1,iptloc 
    407409                 iiglo = secs(jsec)%listPoint(jpt)%I + jpizoom - 1 + nimpp - 1 
    408410                 ijglo = secs(jsec)%listPoint(jpt)%J + jpjzoom - 1 + njmpp - 1 
    409                  WRITE(narea+200,*)'avant # I J : ',iiglo,ijglo 
    410411              ENDDO 
    411412              ENDIF 
     
    421422           ENDIF 
    422423           IF(jsec==nn_secdebug .AND. secs(jsec)%nb_point .NE. 0)THEN 
    423               WRITE(narea+200,*)'apres secs(jsec)%nb_point iptloc ',secs(jsec)%nb_point,iptloc 
    424424              DO jpt = 1,secs(jsec)%nb_point 
    425425                 iiglo = secs(jsec)%listPoint(jpt)%I + jpizoom - 1 + nimpp - 1 
    426426                 ijglo = secs(jsec)%listPoint(jpt)%J + jpjzoom - 1 + njmpp - 1 
    427                  WRITE(narea+200,*)'apres # I J : ',iiglo,ijglo 
    428427              ENDDO 
    429428           ENDIF 
     
    626625        ELSE                                ; isgnv =  1 
    627626        ENDIF 
    628  
    629         IF( ld_debug )write(numout,*)"isgnu isgnv ",isgnu,isgnv 
     627        IF( sec%slopeSection .GE. 9999. )     isgnv =  1 
     628 
     629        IF( ld_debug )write(numout,*)"sec%slopeSection isgnu isgnv ",sec%slopeSection,isgnu,isgnv 
    630630 
    631631        !--------------------------------------! 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DIA/diahth.F90

    r7 r72  
    313313      ! surface boundary condition 
    314314      IF( lk_vvl ) THEN   ;   zthick(:,:) = 0._wp       ;   htc3(:,:) = 0._wp                                    
    315       ELSE                ;   zthick(:,:) = sshn(:,:)   ;   htc3(:,:) = tsn(:,:,jk,jp_tem) * sshn(:,:) * tmask(:,:,jk)    
     315      ELSE                ;   zthick(:,:) = sshn(:,:)   ;   htc3(:,:) = tsn(:,:,1,jp_tem) * sshn(:,:) * tmask(:,:,1)    
    316316      ENDIF 
    317317      ! integration down to ilevel 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DIA/diaptr.F90

    r46 r72  
    455455      NAMELIST/namptr/ ln_diaptr, ln_diaznl, ln_subbas, ln_ptrcomp, nn_fptr, nn_fwri 
    456456      !!---------------------------------------------------------------------- 
    457       IF( nn_timing == 1 )   CALL timing_start('dia_ptr_init') 
    458457 
    459458      REWIND( numnam )                 ! Read Namelist namptr : poleward transport parameters 
     
    474473       
    475474      IF( ln_diaptr) THEN   
     475      
     476         IF( nn_timing == 1 )   CALL timing_start('dia_ptr_init') 
    476477       
    477478         IF( ln_subbas ) THEN   ;   nptr = 5       ! Global, Atlantic, Pacific, Indian, Indo-Pacific 
     
    528529         nidom_ptr = FLIO_DOM_NONE 
    529530#endif 
     531      IF( nn_timing == 1 )   CALL timing_stop('dia_ptr_init') 
     532      ! 
    530533      ENDIF  
    531       !  
    532       IF( nn_timing == 1 )   CALL timing_stop('dia_ptr_init') 
    533534      !  
    534535   END SUBROUTINE dia_ptr_init 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90

    r55 r72  
    173173         z3d(:,:,jpk) = 0.e0 
    174174         DO jk = 1, jpkm1 
    175             z3d(:,:,jk) = rau0 * un(:,:,jk) * e1u(:,:) * fse3u(:,:,jk) 
     175            z3d(:,:,jk) = rau0 * un(:,:,jk) * e2u(:,:) * fse3u(:,:,jk) 
    176176         END DO 
    177177         CALL iom_put( "u_masstr", z3d )                  ! mass transport in i-direction 
     
    188188         CALL iom_put( "u_heattr", z2d )                  ! heat transport in i-direction 
    189189         DO jk = 1, jpkm1 
    190             z3d(:,:,jk) = rau0 * vn(:,:,jk) * e2v(:,:) * fse3v(:,:,jk) 
     190            z3d(:,:,jk) = rau0 * vn(:,:,jk) * e1v(:,:) * fse3v(:,:,jk) 
    191191         END DO 
    192192         CALL iom_put( "v_masstr", z3d )                  ! mass transport in j-direction 
     
    700700      !!---------------------------------------------------------------------- 
    701701      !  
    702       IF( nn_timing == 1 )   CALL timing_start('dia_wri_state') 
    703  
    704702      ! 0. Initialisation 
    705703      ! ----------------- 
     
    796794      ENDIF 
    797795#endif 
    798         
    799       IF( nn_timing == 1 )   CALL timing_stop('dia_wri_state') 
    800796      !  
    801797 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90

    r7 r72  
    12691269      !                                                     ! ================ ! 
    12701270      ! 
    1271       !                                        ! envelop bathymetry saved in hbatt 
     1271      ! Fill ghost rows with appropriate values to avoid undefined e3 values with some mpp decompositions 
     1272      DO ji = nlci+1, jpi  
     1273         zenv(ji,1:nlcj) = zenv(nlci,1:nlcj) 
     1274      END DO 
     1275      ! 
     1276      DO jj = nlcj+1, jpj 
     1277         zenv(:,jj) = zenv(:,nlcj) 
     1278      END DO 
     1279      ! 
     1280      ! Envelope bathymetry saved in hbatt 
    12721281      hbatt(:,:) = zenv(:,:)  
    12731282      IF( MINVAL( gphit(:,:) ) * MAXVAL( gphit(:,:) ) <= 0._wp ) THEN 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/istate.F90

    r61 r72  
    107107         hdivb(:,:,:) = 0._wp   ;   hdivn(:,:,:) = 0._wp 
    108108         ! 
    109          !                                       ! define e3u_b, e3v_b from e3t_b initialized in domzgr 
    110          CALL dom_vvl_2( nit000, fse3u_b(:,:,:), fse3v_b(:,:,:) ) 
    111          ! 
    112109         IF( cp_cfg == 'eel' ) THEN 
    113110            CALL istate_eel                      ! EEL   configuration : start from pre-defined U,V T-S fields 
     
    134131            ENDDO 
    135132         ENDIF 
     133         !                                       ! define e3u_b, e3v_b from e3t_b initialized in domzgr 
     134         CALL dom_vvl_2( nit000, fse3u_b(:,:,:), fse3v_b(:,:,:) ) 
    136135         !  
    137136      ENDIF 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DYN/dynnept.F90

    r1 r72  
    115115                            ln_neptramp, rn_htrmin, rn_htrmax 
    116116      !!---------------------------------------------------------------------- 
    117       !                                                           ! Dynamically allocate local work arrays 
    118       CALL wrk_alloc( jpi, jpj     , ht, htn, tscale, tsp, hur_n, hvr_n, hu_n, hv_n  )  
    119       CALL wrk_alloc( jpi, jpj, jpk, znmask                                          )  
    120       ! 
    121117      ! Define the (simplified) Neptune parameters 
    122118      ! ========================================== 
    123119 
    124 !!    WRITE(numout,*) ' start dynnept namelist' 
    125 !!    CALL FLUSH(numout) 
    126120      REWIND( numnam )                  ! Read Namelist namdyn_nept:  Simplified Neptune 
    127121      READ  ( numnam, namdyn_nept ) 
    128 !!    WRITE(numout,*) ' dynnept namelist done' 
    129 !!    CALL FLUSH(numout) 
    130122 
    131123      IF(lwp) THEN                      ! Control print 
    132124         WRITE(numout,*) 
    133          WRITE(numout,*) 'dyn_nept_init : Simplified Neptune module enabled' 
     125         WRITE(numout,*) 'dyn_nept_init : Simplified Neptune module' 
    134126         WRITE(numout,*) '~~~~~~~~~~~~~' 
    135127         WRITE(numout,*) ' -->   Reading namelist namdyn_nept parameters:' 
    136128         WRITE(numout,*) '       ln_neptsimp          = ', ln_neptsimp 
    137129         WRITE(numout,*) 
    138          WRITE(numout,*) '       ln_smooth_neptvel    = ', ln_smooth_neptvel 
    139          WRITE(numout,*) '       rn_tslse             = ', rn_tslse 
    140          WRITE(numout,*) '       rn_tslsp             = ', rn_tslsp 
    141          WRITE(numout,*) 
    142          WRITE(numout,*) '       ln_neptramp          = ', ln_neptramp 
    143          WRITE(numout,*) '       rn_htrmin            = ', rn_htrmin 
    144          WRITE(numout,*) '       rn_htrmax            = ', rn_htrmax 
    145          WRITE(numout,*) 
    146          CALL FLUSH(numout) 
    147       ENDIF 
     130         IF( ln_neptsimp ) THEN 
     131            WRITE(numout,*) '       ln_smooth_neptvel    = ', ln_smooth_neptvel 
     132            WRITE(numout,*) '       rn_tslse             = ', rn_tslse 
     133            WRITE(numout,*) '       rn_tslsp             = ', rn_tslsp 
     134            WRITE(numout,*) 
     135            WRITE(numout,*) '       ln_neptramp          = ', ln_neptramp 
     136            WRITE(numout,*) '       rn_htrmin            = ', rn_htrmin 
     137            WRITE(numout,*) '       rn_htrmax            = ', rn_htrmax 
     138            WRITE(numout,*) 
     139         ENDIF 
     140      ENDIF 
     141      ! 
     142      IF( .NOT. ln_neptsimp ) RETURN 
     143      !                                 ! Dynamically allocate local work arrays 
     144      CALL wrk_alloc( jpi, jpj     , ht, htn, tscale, tsp, hur_n, hvr_n, hu_n, hv_n  )  
     145      CALL wrk_alloc( jpi, jpj, jpk, znmask                                          )  
    148146 
    149147      IF( ln_smooth_neptvel ) THEN 
     
    151149      ELSE 
    152150         IF(lwp) WRITE(numout,*) ' -->   neptune velocities will not be smoothed' 
    153       ENDIF 
    154  
    155       IF( ln_neptsimp ) THEN 
    156           IF(lwp) WRITE(numout,*) ' -->   ln_neptsimp enabled, solving for U-UN' 
    157       ELSE 
    158           IF(lwp) WRITE(numout,*) ' -->   ln_neptsimp disabled' 
    159           RETURN 
    160151      ENDIF 
    161152 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DYN/dynvor.F90

    r7 r72  
    595595         IF(lwp) WRITE(numout,*) 'dyn:vor_een : vorticity term: energy and enstrophy conserving scheme' 
    596596         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 
    597          IF( .NOT.lk_vvl ) THEN 
     597#if ! defined key_vvl 
     598         IF( .NOT.ALLOCATED(ze3f) ) THEN 
    598599            ALLOCATE( ze3f(jpi,jpj,jpk) , STAT=ierr ) 
    599600            IF( lk_mpp    )   CALL mpp_sum ( ierr ) 
    600601            IF( ierr /= 0 )   CALL ctl_stop( 'STOP', 'dyn:vor_een : unable to allocate arrays' ) 
    601602         ENDIF 
     603#endif 
    602604      ENDIF 
    603605 
  • trunk/NEMOGCM/NEMO/OPA_SRC/LBC/lib_mpp.F90

    r7 r72  
    28282828 
    28292829   SUBROUTINE mppstop 
    2830       WRITE(*,*) 'mppstop: You should not have seen this print if running in mpp mode! error?...' 
    2831       WRITE(*,*) 'mppstop: ..otherwise this is a stop condition raised by ctl_stop in single processor mode' 
    2832       STOP 
     2830      STOP      ! non MPP case, just stop the run 
    28332831   END SUBROUTINE mppstop 
    28342832 
  • trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90

    r66 r72  
    737737               !                                                       ! (geographical to local grid -> rotate the components) 
    738738               CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), srcv(jpr_otx1)%clgrid, 'en->i', ztx )    
    739                frcv(jpr_otx1)%z3(:,:,1) = ztx(:,:)      ! overwrite 1st component on the 1st grid 
    740739               IF( srcv(jpr_otx2)%laction ) THEN 
    741740                  CALL rot_rep( frcv(jpr_otx2)%z3(:,:,1), frcv(jpr_oty2)%z3(:,:,1), srcv(jpr_otx2)%clgrid, 'en->j', zty )    
     
    743742                  CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), srcv(jpr_otx1)%clgrid, 'en->j', zty )   
    744743               ENDIF 
     744               frcv(jpr_otx1)%z3(:,:,1) = ztx(:,:)      ! overwrite 1st component on the 1st grid 
    745745               frcv(jpr_oty1)%z3(:,:,1) = zty(:,:)      ! overwrite 2nd component on the 2nd grid 
    746746            ENDIF 
     
    965965               !                                                       ! (geographical to local grid -> rotate the components) 
    966966               CALL rot_rep( frcv(jpr_itx1)%z3(:,:,1), frcv(jpr_ity1)%z3(:,:,1), srcv(jpr_itx1)%clgrid, 'en->i', ztx )    
    967                frcv(jpr_itx1)%z3(:,:,1) = ztx(:,:)      ! overwrite 1st component on the 1st grid 
    968967               IF( srcv(jpr_itx2)%laction ) THEN 
    969968                  CALL rot_rep( frcv(jpr_itx2)%z3(:,:,1), frcv(jpr_ity2)%z3(:,:,1), srcv(jpr_itx2)%clgrid, 'en->j', zty )    
     
    971970                  CALL rot_rep( frcv(jpr_itx1)%z3(:,:,1), frcv(jpr_ity1)%z3(:,:,1), srcv(jpr_itx1)%clgrid, 'en->j', zty )   
    972971               ENDIF 
     972               frcv(jpr_itx1)%z3(:,:,1) = ztx(:,:)      ! overwrite 1st component on the 1st grid 
    973973               frcv(jpr_ity1)%z3(:,:,1) = zty(:,:)      ! overwrite 2nd component on the 1st grid 
    974974            ENDIF 
  • trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbcdcy.F90

    r7 r72  
    102102         rcc(:,:) = zconvrad * glamt(:,:) - rpi 
    103103         ! time of midday 
    104          rtmd(:,:) = 0.5 - glamt(:,:) / 360. 
    105          rtmd(:,:) = MOD( (rtmd(:,:) + 1.), 1. ) 
     104         rtmd(:,:) = 0.5_wp - glamt(:,:) / 360._wp 
     105         rtmd(:,:) = MOD( (rtmd(:,:) + 1._wp) , 1._wp) 
    106106      ENDIF 
    107107 
     
    118118         zdsws = REAL(11 + nday_year, wp) 
    119119         ! declination of the earths orbit 
    120          zdecrad = (-23.5 * zconvrad) * COS( zdsws * ztwopi / REAL(nyear_len(1),wp) ) 
     120         zdecrad = (-23.5_wp * zconvrad) * COS( zdsws * ztwopi / REAL(nyear_len(1),wp) ) 
    121121         ! Compute A and B needed to compute the time integral of the diurnal cycle 
    122122         
     
    136136         DO jj = 1, jpj 
    137137            DO ji = 1, jpi 
    138                IF ( ABS(rab(ji,jj)) < 1 ) THEN         ! day duration is less than 24h 
     138               IF ( ABS(rab(ji,jj)) < 1._wp ) THEN         ! day duration is less than 24h 
    139139         ! When is it night? 
    140140                  ztx = zinvtwopi * (ACOS(rab(ji,jj)) - rcc(ji,jj)) 
    141141                  ztest = -rbb(ji,jj) * SIN( rcc(ji,jj) + ztwopi * ztx ) 
    142142         ! is it dawn or dusk? 
    143                   IF ( ztest > 0 ) THEN 
     143                  IF ( ztest > 0._wp ) THEN 
    144144                     rdawn(ji,jj) = ztx 
    145145                     rdusk(ji,jj) = rtmd(ji,jj) + ( rtmd(ji,jj) - rdawn(ji,jj) ) 
     
    149149                  ENDIF 
    150150               ELSE 
    151                   rdawn(ji,jj) = rtmd(ji,jj) + 0.5 
     151                  rdawn(ji,jj) = rtmd(ji,jj) + 0.5_wp 
    152152                  rdusk(ji,jj) = rdawn(ji,jj) 
    153153               ENDIF 
     
    157157         rdusk(:,:) = MOD( (rdusk(:,:) + 1._wp), 1._wp ) 
    158158 
    159          !     2.2 Compute the scalling function: 
    160          !         S* = the inverse of the time integral of the diurnal cycle from dawm to dusk 
     159         !     2.2 Compute the scaling function: 
     160         !         S* = the inverse of the time integral of the diurnal cycle from dawn to dusk 
     161         !         Avoid possible infinite scaling factor, associated with very short daylight 
     162         !         periods, by ignoring periods less than 1/1000th of a day (ticket #1040) 
    161163         DO jj = 1, jpj 
    162164            DO ji = 1, jpi 
    163                IF ( ABS(rab(ji,jj)) < 1 ) THEN         ! day duration is less than 24h 
     165               IF ( ABS(rab(ji,jj)) < 1._wp ) THEN         ! day duration is less than 24h 
     166                  rscal(ji,jj) = 0.0_wp 
    164167                  IF ( rdawn(ji,jj) < rdusk(ji,jj) ) THEN      ! day time in one part 
    165                      rscal(ji,jj) = fintegral(rdawn(ji,jj), rdusk(ji,jj), raa(ji,jj), rbb(ji,jj), rcc(ji,jj))  
    166                      rscal(ji,jj) = 1. / rscal(ji,jj) 
     168                     IF( (rdusk(ji,jj) - rdawn(ji,jj) ) .ge. 0.001_wp ) THEN 
     169                       rscal(ji,jj) = fintegral(rdawn(ji,jj), rdusk(ji,jj), raa(ji,jj), rbb(ji,jj), rcc(ji,jj))  
     170                       rscal(ji,jj) = 1._wp / rscal(ji,jj) 
     171                     ENDIF 
    167172                  ELSE                                         ! day time in two parts 
    168                      rscal(ji,jj) = fintegral(0., rdusk(ji,jj), raa(ji,jj), rbb(ji,jj), rcc(ji,jj))   & 
    169                         &         + fintegral(rdawn(ji,jj), 1., raa(ji,jj), rbb(ji,jj), rcc(ji,jj))  
    170                      rscal(ji,jj) = 1. / rscal(ji,jj) 
     173                     IF( (rdusk(ji,jj) + (1._wp - rdawn(ji,jj)) ) .ge. 0.001_wp ) THEN 
     174                       rscal(ji,jj) = fintegral(0._wp, rdusk(ji,jj), raa(ji,jj), rbb(ji,jj), rcc(ji,jj))   & 
     175                          &         + fintegral(rdawn(ji,jj), 1._wp, raa(ji,jj), rbb(ji,jj), rcc(ji,jj))  
     176                       rscal(ji,jj) = 1. / rscal(ji,jj) 
     177                     ENDIF 
    171178                  ENDIF 
    172179               ELSE 
    173180                  IF ( raa(ji,jj) > rbb(ji,jj) ) THEN         ! 24h day 
    174                      rscal(ji,jj) = fintegral(0., 1., raa(ji,jj), rbb(ji,jj), rcc(ji,jj))  
    175                      rscal(ji,jj) = 1. / rscal(ji,jj) 
     181                     rscal(ji,jj) = fintegral(0._wp, 1._wp, raa(ji,jj), rbb(ji,jj), rcc(ji,jj))  
     182                     rscal(ji,jj) = 1._wp / rscal(ji,jj) 
    176183                  ELSE                                          ! No day 
    177                      rscal(ji,jj) = 0.e0 
     184                     rscal(ji,jj) = 0.0_wp 
    178185                  ENDIF 
    179186               ENDIF 
     
    191198      DO jj = 1, jpj 
    192199         DO ji = 1, jpi 
    193             IF( ABS(rab(ji,jj)) < 1 ) THEN         ! day duration is less than 24h 
     200            IF( ABS(rab(ji,jj)) < 1._wp ) THEN         ! day duration is less than 24h 
    194201               ! 
    195202               IF( rdawn(ji,jj) < rdusk(ji,jj) ) THEN       ! day time in one part 
     
    218225                  ! 
    219226               ELSE                                         ! No day 
    220                   zqsrout(ji,jj) = 0.e0 
     227                  zqsrout(ji,jj) = 0.0_wp 
    221228               ENDIF 
    222229            ENDIF 
  • trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbcrnf.F90

    r43 r72  
    346346            ! 
    347347            nk_rnf(:,:) = 0                               ! set the number of level over which river runoffs are applied 
    348             DO jj = 1, jpj   
    349                DO ji = 1, jpi   
    350                   IF( h_rnf(ji,jj) > 0._wp ) THEN   
    351                      jk = 2   
    352                      DO WHILE ( jk /= mbkt(ji,jj) .AND. fsdept(ji,jj,jk) < h_rnf(ji,jj) ) ;  jk = jk + 1 ;  END DO   
    353                      nk_rnf(ji,jj) = jk   
    354                   ELSEIF( h_rnf(ji,jj) == -1   ) THEN   ;  nk_rnf(ji,jj) = 1   
    355                   ELSEIF( h_rnf(ji,jj) == -999 ) THEN   ;  nk_rnf(ji,jj) = mbkt(ji,jj) 
    356                   ELSEIF( h_rnf(ji,jj) /=  0   ) THEN   
    357                      CALL ctl_stop( 'runoff depth not positive, and not -999 or -1, rnf value in file fort.999'  )   
    358                      WRITE(999,*) 'ji, jj, rnf(ji,jj) :', ji, jj, rnf(ji,jj)   
    359                   ENDIF   
    360                END DO   
    361             END DO   
    362             DO jj = 1, jpj                                ! set the associated depth  
    363                DO ji = 1, jpi  
     348            DO jj = 1, jpj 
     349               DO ji = 1, jpi 
     350                  IF( h_rnf(ji,jj) > 0._wp ) THEN 
     351                     jk = 2 
     352                     DO WHILE ( jk /= mbkt(ji,jj) .AND. fsdept(ji,jj,jk) < h_rnf(ji,jj) ) ;  jk = jk + 1 ;  END DO 
     353                     nk_rnf(ji,jj) = jk 
     354                  ELSEIF( h_rnf(ji,jj) == -1._wp   ) THEN   ;  nk_rnf(ji,jj) = 1 
     355                  ELSEIF( h_rnf(ji,jj) == -999._wp ) THEN   ;  nk_rnf(ji,jj) = mbkt(ji,jj) 
     356                  ELSE 
     357                     CALL ctl_stop( 'runoff depth not positive, and not -999 or -1, rnf value in file fort.999'  ) 
     358                     WRITE(999,*) 'ji, jj, rnf(ji,jj) :', ji, jj, rnf(ji,jj) 
     359                  ENDIF 
     360               END DO 
     361            END DO 
     362            DO jj = 1, jpj                                ! set the associated depth 
     363               DO ji = 1, jpi 
    364364                  h_rnf(ji,jj) = 0._wp 
    365365                  DO jk = 1, nk_rnf(ji,jj)                         
  • trunk/NEMOGCM/NEMO/OPA_SRC/SBC/tide_mod.F90

    r1 r72  
    2121       jpmax_harmo = 19             ! maximum number of harmonic 
    2222 
    23   TYPE tide 
     23  TYPE,PUBLIC:: tide 
    2424     CHARACTER(LEN=4)  :: cname_tide 
    2525     REAL(wp) :: equitide 
  • trunk/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_eiv.F90

    r7 r72  
    168168                  DO ji = fs_2, fs_jpim1   ! vector opt. 
    169169                     z2d(ji,jj) = z2d(ji,jj) + zztmp * u_eiv(ji,jj,jk) & 
    170                        &         * (tsn(ji,jj,jk,jp_tem)+tsn(ji+1,jj,jk,jp_tem)) * e1u(ji,jj) * fse3u(ji,jj,jk)  
     170                       &         * (tsn(ji,jj,jk,jp_tem)+tsn(ji+1,jj,jk,jp_tem)) * e2u(ji,jj) * fse3u(ji,jj,jk)  
    171171                  END DO 
    172172               END DO 
     
    179179                  DO ji = fs_2, fs_jpim1   ! vector opt. 
    180180                     z2d(ji,jj) = z2d(ji,jj) + zztmp * v_eiv(ji,jj,jk) & 
    181                      &           * (tsn(ji,jj,jk,jp_tem)+tsn(ji,jj+1,jk,jp_tem)) * e2v(ji,jj) * fse3v(ji,jj,jk)  
     181                     &           * (tsn(ji,jj,jk,jp_tem)+tsn(ji,jj+1,jk,jp_tem)) * e1v(ji,jj) * fse3v(ji,jj,jk)  
    182182                  END DO 
    183183               END DO 
  • trunk/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_ubs.F90

    r7 r72  
    2424   USE wrk_nemo        ! Memory Allocation 
    2525   USE timing          ! Timing 
     26   USE lib_fortran     ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 
    2627 
    2728   IMPLICIT NONE 
     
    5051      !!      and add it to the general trend of passive tracer equations. 
    5152      !! 
    52       !! ** Method  :   The upstream biased third (UBS) is order scheme based  
    53       !!      on an upstream-biased parabolic interpolation (Shchepetkin and McWilliams 2005) 
     53      !! ** Method  :   The upstream biased 3rd order scheme (UBS) is based on an  
     54      !!      upstream-biased parabolic interpolation (Shchepetkin and McWilliams 2005) 
    5455      !!      It is only used in the horizontal direction. 
    5556      !!      For example the i-component of the advective fluxes are given by : 
    56       !!                !  e1u e3u un ( mi(Tn) - zltu(i  ) )   if un(i) >= 0 
     57      !!                !  e2u e3u un ( mi(Tn) - zltu(i  ) )   if un(i) >= 0 
    5758      !!          zwx = !  or  
    58       !!                !  e1u e3u un ( mi(Tn) - zltu(i+1) )   if un(i) < 0 
     59      !!                !  e2u e3u un ( mi(Tn) - zltu(i+1) )   if un(i) < 0 
    5960      !!      where zltu is the second derivative of the before temperature field: 
    6061      !!          zltu = 1/e3t di[ e2u e3u / e1u di[Tb] ] 
     
    6768      !!      of the scheme, is evaluated using the before velocity (forward in time).  
    6869      !!      Note that UBS is not positive. Do not use it on passive tracers. 
    69       !!      On the vertical, the advection is evaluated using a TVD scheme, as 
    70       !!      the UBS have been found to be too diffusive. 
     70      !!                On the vertical, the advection is evaluated using a TVD scheme, 
     71      !!      as the UBS have been found to be too diffusive. 
    7172      !! 
    7273      !! ** Action : - update (pta) with the now advective tracer trends 
     
    8283      INTEGER                              , INTENT(in   ) ::   kjpt            ! number of tracers 
    8384      REAL(wp), DIMENSION(        jpk     ), INTENT(in   ) ::   p2dt            ! vertical profile of tracer time-step 
    84       REAL(wp), DIMENSION(jpi,jpj,jpk     ), INTENT(in   ) ::   pun, pvn, pwn   ! 3 ocean velocity components 
     85      REAL(wp), DIMENSION(jpi,jpj,jpk     ), INTENT(in   ) ::   pun, pvn, pwn   ! 3 ocean transport components 
    8586      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   ptb, ptn        ! before and now tracer fields 
    8687      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pta             ! tracer trend  
     
    141142            DO jj = 1, jpjm1 
    142143               DO ji = 1, fs_jpim1   ! vector opt. 
    143                   ! upstream transport 
     144                  ! upstream transport (x2) 
    144145                  zfp_ui = pun(ji,jj,jk) + ABS( pun(ji,jj,jk) ) 
    145146                  zfm_ui = pun(ji,jj,jk) - ABS( pun(ji,jj,jk) ) 
    146147                  zfp_vj = pvn(ji,jj,jk) + ABS( pvn(ji,jj,jk) ) 
    147148                  zfm_vj = pvn(ji,jj,jk) - ABS( pvn(ji,jj,jk) ) 
    148                   ! centered scheme 
    149                   zcenut = 0.5 * pun(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji+1,jj  ,jk,jn) ) 
    150                   zcenvt = 0.5 * pvn(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji  ,jj+1,jk,jn) ) 
    151                   ! UBS scheme 
    152                   zwx(ji,jj,jk) =  zcenut - zfp_ui * zltu(ji,jj,jk) - zfm_ui * zltu(ji+1,jj,jk)  
    153                   zwy(ji,jj,jk) =  zcenvt - zfp_vj * zltv(ji,jj,jk) - zfm_vj * zltv(ji,jj+1,jk)  
     149                  ! 2nd order centered advective fluxes (x2) 
     150                  zcenut = pun(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji+1,jj  ,jk,jn) ) 
     151                  zcenvt = pvn(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji  ,jj+1,jk,jn) ) 
     152                  ! UBS advective fluxes 
     153                  zwx(ji,jj,jk) = 0.5 * ( zcenut - zfp_ui * zltu(ji,jj,jk) - zfm_ui * zltu(ji+1,jj,jk) ) 
     154                  zwy(ji,jj,jk) = 0.5 * ( zcenvt - zfp_vj * zltv(ji,jj,jk) - zfm_vj * zltv(ji,jj+1,jk) ) 
    154155               END DO 
    155156            END DO 
     
    198199 
    199200         ! Surface value 
    200          IF( lk_vvl ) THEN   ;   ztw(:,:,1) = 0.e0                      ! variable volume : flux set to zero 
    201          ELSE                ;   ztw(:,:,1) = pwn(:,:,1) * ptb(:,:,1,jn)   ! free constant surface  
     201         IF( lk_vvl ) THEN   ;   ztw(:,:,1) = 0.e0                         ! variable volume : flux set to zero 
     202         ELSE                ;   ztw(:,:,1) = pwn(:,:,1) * ptb(:,:,1,jn)   ! constant volume : non zero flux though z=0  
    202203         ENDIF 
    203204         !  upstream advection with initial mass fluxes & intermediate update 
  • trunk/NEMOGCM/NEMO/OPA_SRC/TRA/traldf_bilapg.F90

    r7 r72  
    248248         ! "Poleward" diffusive heat or salt transport 
    249249         IF( cdtype == 'TRA' .AND. ln_diaptr .AND. ( kaht == 2 ) .AND. ( MOD( kt, nn_fptr ) == 0 ) ) THEN 
    250             IF( jn == jp_tem)   htr_ldf(:) = ptr_vj( zftv(:,:,:) ) 
    251             IF( jn == jp_sal)   str_ldf(:) = ptr_vj( zftv(:,:,:) ) 
     250            ! note sign is reversed to give down-gradient diffusive transports (#1043) 
     251            IF( jn == jp_tem)   htr_ldf(:) = ptr_vj( -zftv(:,:,:) ) 
     252            IF( jn == jp_sal)   str_ldf(:) = ptr_vj( -zftv(:,:,:) ) 
    252253         ENDIF 
    253254 
  • trunk/NEMOGCM/NEMO/OPA_SRC/TRA/traldf_iso.F90

    r7 r72  
    212212         ! "Poleward" diffusive heat or salt transports (T-S case only) 
    213213         IF( cdtype == 'TRA' .AND. ln_diaptr .AND. ( MOD( kt, nn_fptr ) == 0 ) ) THEN 
    214             IF( jn == jp_tem)   htr_ldf(:) = ptr_vj( zftv(:,:,:) ) 
    215             IF( jn == jp_sal)   str_ldf(:) = ptr_vj( zftv(:,:,:) ) 
     214            ! note sign is reversed to give down-gradient diffusive transports (#1043) 
     215            IF( jn == jp_tem)   htr_ldf(:) = ptr_vj( -zftv(:,:,:) ) 
     216            IF( jn == jp_sal)   str_ldf(:) = ptr_vj( -zftv(:,:,:) ) 
    216217         ENDIF 
    217218  
     
    219220         IF( cdtype == 'TRA' .AND. jn == jp_tem  ) THEN 
    220221            z2d(:,:) = 0._wp  
    221             zztmp = rau0 * rcp  
     222            ! note sign is reversed to give down-gradient diffusive transports (#1043) 
     223            zztmp = -1.0_wp * rau0 * rcp  
    222224            DO jk = 1, jpkm1 
    223225               DO jj = 2, jpjm1 
  • trunk/NEMOGCM/NEMO/OPA_SRC/TRA/trasbc.F90

    r43 r72  
    205205      !---------------------------------------- 
    206206      ! 
    207       zfact = 0.5e0 
    208  
    209       ! Effect on (t,s) due to river runoff (dilution effect automatically applied via vertical tracer advection)  
    210       IF( ln_rnf ) THEN   
     207      IF( ln_rnf ) THEN         ! input of heat and salt due to river runoff  
     208         zfact = 0.5_wp 
    211209         DO jj = 2, jpj  
    212210            DO ji = fs_2, fs_jpim1 
    213                zdep = 1. / h_rnf(ji,jj) 
    214                zdep = zfact * zdep   
    215                IF ( rnf(ji,jj) /= 0._wp ) THEN 
     211               IF( rnf(ji,jj) /= 0._wp ) THEN 
     212                  zdep = zfact / h_rnf(ji,jj) 
    216213                  DO jk = 1, nk_rnf(ji,jj) 
    217214                                        tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem)   & 
     
    223220            END DO   
    224221         END DO   
    225       ENDIF   
    226 !!gm  It should be useless 
    227       CALL lbc_lnk( tsa(:,:,:,jp_tem), 'T', 1. )    ;    CALL lbc_lnk( tsa(:,:,:,jp_sal), 'T', 1. ) 
    228  
     222      ENDIF 
     223  
    229224      IF( l_trdtra )   THEN                      ! save the horizontal diffusive trends for further diagnostics 
    230225         ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:) 
  • trunk/NEMOGCM/NEMO/OPA_SRC/TRD/trdmld.F90

    r43 r72  
    760760      ! ------------------------------------------------- 
    761761 
    762       IF( ( lk_trdmld ) .AND. ( MOD( nitend, nn_trd ) /= 0 ) ) THEN 
     762      IF( ( lk_trdmld ) .AND. ( MOD( nitend-nit000+1, nn_trd ) /= 0 ) ) THEN 
    763763         WRITE(numout,cform_err) 
    764764         WRITE(numout,*) '                Your nitend parameter, nitend = ', nitend 
  • trunk/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfevd.F90

    r7 r72  
    1919   USE zdf_oce         ! ocean vertical physics variables 
    2020   USE zdfkpp          ! KPP vertical mixing 
     21   USE zdfgls          ! GLS vertical mixing 
    2122   USE in_out_manager  ! I/O manager 
    2223   USE iom             ! for iom_put 
     
    6768         IF(lwp) WRITE(numout,*) '~~~~~~~ ' 
    6869         IF(lwp) WRITE(numout,*) 
     70         ! 
     71         IF(lwp .AND. lk_zdfgls )   CALL ctl_warn(' No need zdf_evd with GLS closures ') 
     72         ! 
    6973      ENDIF 
    7074 
  • trunk/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfgls.F90

    r43 r72  
    4343   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   mxln    !: now mixing length 
    4444   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   zwall   !: wall function 
     45   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avt_k   ! not enhanced Kz 
     46   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avm_k   ! not enhanced Kz 
     47   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avmu_k  ! not enhanced Kz 
     48   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avmv_k  ! not enhanced Kz 
    4549   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   ustars2 !: Squared surface velocity scale at T-points 
    4650   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   ustarb2 !: Squared bottom  velocity scale at T-points 
     
    118122      !!---------------------------------------------------------------------- 
    119123      ALLOCATE( en(jpi,jpj,jpk),  mxln(jpi,jpj,jpk), zwall(jpi,jpj,jpk) ,     & 
     124         &      avt_k (jpi,jpj,jpk) , avm_k (jpi,jpj,jpk),                    & 
     125         &      avmu_k(jpi,jpj,jpk) , avmv_k(jpi,jpj,jpk),                    & 
    120126         &      ustars2(jpi,jpj), ustarb2(jpi,jpj)                      , STAT= zdf_gls_alloc ) 
    121127         ! 
     
    158164 
    159165      ustars2 = 0._wp   ;   ustarb2 = 0._wp   ;   psi  = 0._wp   ;   zwall_psi = 0._wp 
     166 
     167      IF( kt /= nit000 ) THEN   ! restore before value to compute tke 
     168         avt (:,:,:) = avt_k (:,:,:) 
     169         avm (:,:,:) = avm_k (:,:,:) 
     170         avmu(:,:,:) = avmu_k(:,:,:) 
     171         avmv(:,:,:) = avmv_k(:,:,:)  
     172      ENDIF 
    160173 
    161174      ! Compute surface and bottom friction at T-points 
     
    881894      ENDIF 
    882895      ! 
     896      avt_k (:,:,:) = avt (:,:,:) 
     897      avm_k (:,:,:) = avm (:,:,:) 
     898      avmu_k(:,:,:) = avmu(:,:,:) 
     899      avmv_k(:,:,:) = avmv(:,:,:) 
     900      ! 
    883901      CALL wrk_dealloc( jpi,jpj, zdep, zflxs, zhsro ) 
    884902      CALL wrk_dealloc( jpi,jpj,jpk, eb, mxlb, shear, eps, zwall_psi, z_elem_a, z_elem_b, z_elem_c, psi ) 
     
    12441262         IF(lwp) WRITE(numout,*) '---- gls-rst ----' 
    12451263         CALL iom_rstput( kt, nitrst, numrow, 'en'   , en    ) 
    1246          CALL iom_rstput( kt, nitrst, numrow, 'avt'  , avt   ) 
    1247          CALL iom_rstput( kt, nitrst, numrow, 'avm'  , avm   ) 
    1248          CALL iom_rstput( kt, nitrst, numrow, 'avmu' , avmu  ) 
    1249          CALL iom_rstput( kt, nitrst, numrow, 'avmv' , avmv  ) 
     1264         CALL iom_rstput( kt, nitrst, numrow, 'avt'  , avt_k  ) 
     1265         CALL iom_rstput( kt, nitrst, numrow, 'avm'  , avm_k  ) 
     1266         CALL iom_rstput( kt, nitrst, numrow, 'avmu' , avmu_k ) 
     1267         CALL iom_rstput( kt, nitrst, numrow, 'avmv' , avmv_k ) 
    12501268         CALL iom_rstput( kt, nitrst, numrow, 'mxln' , mxln  ) 
    12511269         ! 
  • trunk/NEMOGCM/NEMO/TOP_SRC/PISCES/p4zsink.F90

    r46 r72  
    1919   USE prtctl_trc      !  print control for debugging 
    2020   USE iom             !  I/O manager 
     21   USE lib_fortran     ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 
    2122 
    2223   IMPLICIT NONE 
     
    614615      INTEGER  ::   ji, jj, jk, jn 
    615616      REAL(wp) ::   zigma,zew,zign, zflx, zstep 
    616       REAL(wp), POINTER, DIMENSION(:,:,:) :: ztraz, zakz, zwsink2  
     617      REAL(wp), POINTER, DIMENSION(:,:,:) :: ztraz, zakz, zwsink2, ztrb  
    617618      !!--------------------------------------------------------------------- 
    618619      ! 
     
    620621      ! 
    621622      ! Allocate temporary workspace 
    622       CALL wrk_alloc( jpi, jpj, jpk, ztraz, zakz, zwsink2 ) 
     623      CALL wrk_alloc( jpi, jpj, jpk, ztraz, zakz, zwsink2, ztrb ) 
    623624 
    624625      zstep = rfact2 / 2. 
     
    626627      ztraz(:,:,:) = 0.e0 
    627628      zakz (:,:,:) = 0.e0 
     629      ztrb (:,:,:) = trn(:,:,:,jp_tra) 
    628630 
    629631      DO jk = 1, jpkm1 
     
    695697            DO ji = 1, jpi 
    696698               zflx = ( psinkflx(ji,jj,jk) - psinkflx(ji,jj,jk+1) ) / fse3t(ji,jj,jk) 
    697                trb(ji,jj,jk,jp_tra) = trb(ji,jj,jk,jp_tra) + 2. * zflx 
     699               ztrb(ji,jj,jk) = ztrb(ji,jj,jk) + 2. * zflx 
    698700            END DO 
    699701         END DO 
    700702      END DO 
    701703 
    702       trn     (:,:,:,jp_tra) = trb(:,:,:,jp_tra) 
     704      trn     (:,:,:,jp_tra) = ztrb(:,:,:) 
    703705      psinkflx(:,:,:)        = 2. * psinkflx(:,:,:) 
    704706      ! 
    705       CALL wrk_dealloc( jpi, jpj, jpk, ztraz, zakz, zwsink2 ) 
     707      CALL wrk_dealloc( jpi, jpj, jpk, ztraz, zakz, zwsink2, ztrb ) 
    706708      ! 
    707709      IF( nn_timing == 1 )  CALL timing_stop('p4z_sink2') 
  • trunk/NEMOGCM/NEMO/TOP_SRC/TRP/trcsbc.F90

    r7 r72  
    5050      !!            tra = tra + emp * trn / e3t   for k=1 
    5151      !!         where emp, the surface freshwater budget (evaporation minus 
    52       !!         precipitation minus runoff) given in kg/m2/s is divided 
     52      !!         precipitation ) given in kg/m2/s is divided 
    5353      !!         by 1035 kg/m3 (density of ocean water) to obtain m/s. 
    5454      !! 
     
    7979      ENDIF 
    8080 
     81      ! Coupling online : river runoff is added to the horizontal divergence (hdivn) in the subroutine sbc_rnf_div  
     82      ! one only consider the concentration/dilution effect due to evaporation minus precipitation + freezing/melting of sea-ice 
    8183 
    82       IF( lk_offline ) THEN          ! emps in dynamical files contains emps - rnf 
    83          zemps(:,:) = emps(:,:)   
    84       ELSE                           ! Concentration dilution effect on tracer due to evaporation, precipitation, and river runoff 
    85          IF( lk_vvl ) THEN                      ! volume variable 
    86             zemps(:,:) = emps(:,:) - emp(:,:)    
    87 !!ch         zemps(:,:) = 0. 
    88          ELSE                                   ! linear free surface 
    89             IF( ln_rnf ) THEN  ;  zemps(:,:) = emps(:,:) - rnf(:,:)   !  E-P-R 
    90             ELSE               ;  zemps(:,:) = emps(:,:) 
    91             ENDIF  
    92          ENDIF  
     84      ! Coupling in offline, hdivn is computed from ocean horizontal velocities only ; the runoff are not included. 
     85      ! emps in dynamical files contains (emps - rnf) 
     86      IF( .NOT. lk_offline .AND. lk_vvl ) THEN  ! online coupling + volume variable  
     87         zemps(:,:) = emps(:,:) - emp(:,:)    
     88      ELSE 
     89         zemps(:,:) = emps(:,:) 
    9390      ENDIF  
    9491 
  • trunk/NEMOGCM/NEMO/TOP_SRC/TRP/trdmld_trc.F90

    r46 r72  
    11741174      ! ------------------------------------------------- 
    11751175 
    1176       IF( ( lk_trdmld_trc ) .AND. ( MOD( nitend, nn_trd_trc ) /= 0 ) ) THEN 
     1176      IF( ( lk_trdmld_trc ) .AND. ( MOD( nitend-nit000+1, nn_trd_trc ) /= 0 ) ) THEN 
    11771177         WRITE(numout,cform_err) 
    11781178         WRITE(numout,*) '                Your nitend parameter, nitend = ', nitend 
Note: See TracChangeset for help on using the changeset viewer.