New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 7646 for trunk/NEMOGCM/NEMO/OPA_SRC/DIA – NEMO

Ignore:
Timestamp:
2017-02-06T10:25:03+01:00 (7 years ago)
Author:
timgraham
Message:

Merge of dev_merge_2016 into trunk. UPDATE TO ARCHFILES NEEDED for XIOS2.
LIM_SRC_s/limrhg.F90 to follow in next commit due to change of kind (I'm unable to do it in this commit).
Merged using the following steps:

1) svn merge --reintegrate ​svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/trunk .
2) Resolve minor conflicts in sette.sh and namelist_cfg for ORCA2LIM3 (due to a change in trunk after branch was created)
3) svn commit
4) svn switch ​svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/trunk
5) svn merge ​svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/branches/2016/dev_merge_2016 .
6) At this stage I checked out a clean copy of the branch to compare against what is about to be committed to the trunk.
6) svn commit #Commit code to the trunk

In this commit I have also reverted a change to Fcheck_archfile.sh which was causing problems on the Paris machine.

Location:
trunk/NEMOGCM/NEMO/OPA_SRC/DIA
Files:
1 deleted
7 edited

Legend:

Unmodified
Added
Removed

  • trunk/NEMOGCM/NEMO/OPA_SRC/DIA/diaar5.F90

    r6665 r7646  
    66   !! History :  3.2  !  2009-11  (S. Masson)  Original code 
    77   !!            3.3  !  2010-10  (C. Ethe, G. Madec) reorganisation of initialisation phase + merge TRC-TRA 
    8    !!---------------------------------------------------------------------- 
    9 #if defined key_diaar5 
    10    !!---------------------------------------------------------------------- 
    11    !!   'key_diaar5'  :                           activate ar5 diagnotics 
    128   !!---------------------------------------------------------------------- 
    139   !!   dia_ar5       : AR5 diagnostics 
    … …  
    2420   USE phycst         ! physical constant 
    2521   USE in_out_manager  ! I/O manager 
     22   USE zdfddm 
     23   USE zdf_oce 
    2624 
    2725   IMPLICIT NONE 
    … …  
    2927 
    3028   PUBLIC   dia_ar5        ! routine called in step.F90 module 
    31    PUBLIC   dia_ar5_init   ! routine called in opa.F90 module 
    3229   PUBLIC   dia_ar5_alloc  ! routine called in nemogcm.F90 module 
    33  
    34    LOGICAL, PUBLIC, PARAMETER :: lk_diaar5 = .TRUE.   ! coupled flag 
     30   PUBLIC   dia_ar5_hst    ! heat/salt transport 
    3531 
    3632   REAL(wp)                         ::   vol0         ! ocean volume (interior domain) 
    … …  
    3935   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:  ) ::   thick0       ! ocean thickness (interior domain) 
    4036   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   sn0          ! initial salinity 
     37 
     38   LOGICAL  :: l_ar5 
    4139       
     40   !! * Substitutions 
     41#  include "zdfddm_substitute.h90" 
     42#  include "vectopt_loop_substitute.h90" 
    4243   !!---------------------------------------------------------------------- 
    4344   !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
    … …  
    7374      INTEGER  ::   ji, jj, jk                      ! dummy loop arguments 
    7475      REAL(wp) ::   zvolssh, zvol, zssh_steric, zztmp, zarho, ztemp, zsal, zmass 
     76      REAL(wp) ::   zaw, zbw, zrw 
    7577      ! 
    7678      REAL(wp), POINTER, DIMENSION(:,:)     :: zarea_ssh , zbotpres       ! 2D workspace  
     79      REAL(wp), POINTER, DIMENSION(:,:)     :: zpe                         ! 2D workspace  
    7780      REAL(wp), POINTER, DIMENSION(:,:,:)   :: zrhd , zrhop               ! 3D workspace 
    7881      REAL(wp), POINTER, DIMENSION(:,:,:,:) :: ztsn                       ! 4D workspace 
    … …  
    8083      IF( nn_timing == 1 )   CALL timing_start('dia_ar5') 
    8184  
    82       CALL wrk_alloc( jpi , jpj              , zarea_ssh , zbotpres ) 
    83       CALL wrk_alloc( jpi , jpj , jpk        , zrhd      , zrhop    ) 
    84       CALL wrk_alloc( jpi , jpj , jpk , jpts , ztsn                 ) 
    85  
    86       zarea_ssh(:,:) = area(:,:) * sshn(:,:) 
    87  
    88       !                                         ! total volume of liquid seawater 
    89       zvolssh = SUM( zarea_ssh(:,:) )  
    90       IF( lk_mpp )   CALL mpp_sum( zvolssh ) 
    91       zvol = vol0 + zvolssh 
     85      IF( kt == nit000 )     CALL dia_ar5_init 
     86 
     87      IF( l_ar5 ) THEN 
     88         CALL wrk_alloc( jpi , jpj              , zarea_ssh , zbotpres ) 
     89         CALL wrk_alloc( jpi , jpj , jpk        , zrhd      , zrhop    ) 
     90         CALL wrk_alloc( jpi , jpj , jpk , jpts , ztsn                 ) 
     91         zarea_ssh(:,:) = area(:,:) * sshn(:,:) 
     92      ENDIF 
     93      ! 
     94      IF( iom_use( 'voltot' ) .OR. iom_use( 'sshtot' )  .OR. iom_use( 'sshdyn' )  ) THEN     
     95         !                                         ! total volume of liquid seawater 
     96         zvolssh = SUM( zarea_ssh(:,:) )  
     97         IF( lk_mpp )   CALL mpp_sum( zvolssh ) 
     98         zvol = vol0 + zvolssh 
    9299       
    93       CALL iom_put( 'voltot', zvol               ) 
    94       CALL iom_put( 'sshtot', zvolssh / area_tot ) 
    95  
    96       !                      
    97       ztsn(:,:,:,jp_tem) = tsn(:,:,:,jp_tem)                    ! thermosteric ssh 
    98       ztsn(:,:,:,jp_sal) = sn0(:,:,:) 
    99       CALL eos( ztsn, zrhd, gdept_n(:,:,:) )                       ! now in situ density using initial salinity 
    100       ! 
    101       zbotpres(:,:) = 0._wp                        ! no atmospheric surface pressure, levitating sea-ice 
    102       DO jk = 1, jpkm1 
    103          zbotpres(:,:) = zbotpres(:,:) + e3t_n(:,:,jk) * zrhd(:,:,jk) 
    104       END DO 
    105       IF( ln_linssh ) THEN 
    106          IF( ln_isfcav ) THEN 
    107             DO ji=1,jpi 
    108                DO jj=1,jpj 
    109                   zbotpres(ji,jj) = zbotpres(ji,jj) + sshn(ji,jj) * zrhd(ji,jj,mikt(ji,jj)) + riceload(ji,jj) 
    110                END DO 
    111             END DO 
    112          ELSE 
    113             zbotpres(:,:) = zbotpres(:,:) + sshn(:,:) * zrhd(:,:,1) 
    114          END IF 
     100         CALL iom_put( 'voltot', zvol               ) 
     101         CALL iom_put( 'sshtot', zvolssh / area_tot ) 
     102         CALL iom_put( 'sshdyn', sshn(:,:) - (zvolssh / area_tot) ) 
     103         ! 
     104      ENDIF 
     105 
     106      IF( iom_use( 'botpres' ) .OR. iom_use( 'sshthster' )  .OR. iom_use( 'sshsteric' )  ) THEN     
     107         !                      
     108         ztsn(:,:,:,jp_tem) = tsn(:,:,:,jp_tem)                    ! thermosteric ssh 
     109         ztsn(:,:,:,jp_sal) = sn0(:,:,:) 
     110         CALL eos( ztsn, zrhd, gdept_n(:,:,:) )                       ! now in situ density using initial salinity 
     111         ! 
     112         zbotpres(:,:) = 0._wp                        ! no atmospheric surface pressure, levitating sea-ice 
     113         DO jk = 1, jpkm1 
     114            zbotpres(:,:) = zbotpres(:,:) + e3t_n(:,:,jk) * zrhd(:,:,jk) 
     115         END DO 
     116         IF( ln_linssh ) THEN 
     117            IF( ln_isfcav ) THEN 
     118               DO ji = 1, jpi 
     119                  DO jj = 1, jpj 
     120                     zbotpres(ji,jj) = zbotpres(ji,jj) + sshn(ji,jj) * zrhd(ji,jj,mikt(ji,jj)) + riceload(ji,jj) 
     121                  END DO 
     122               END DO 
     123            ELSE 
     124               zbotpres(:,:) = zbotpres(:,:) + sshn(:,:) * zrhd(:,:,1) 
     125            END IF 
    115126!!gm 
    116127!!gm   riceload should be added in both ln_linssh=T or F, no? 
    117128!!gm 
    118       END IF 
    119       !                                          
    120       zarho = SUM( area(:,:) * zbotpres(:,:) )  
    121       IF( lk_mpp )   CALL mpp_sum( zarho ) 
    122       zssh_steric = - zarho / area_tot 
    123       CALL iom_put( 'sshthster', zssh_steric ) 
     129         END IF 
     130         !                                          
     131         zarho = SUM( area(:,:) * zbotpres(:,:) )  
     132         IF( lk_mpp )   CALL mpp_sum( zarho ) 
     133         zssh_steric = - zarho / area_tot 
     134         CALL iom_put( 'sshthster', zssh_steric ) 
    124135       
    125       !                                         ! steric sea surface height 
    126       CALL eos( tsn, zrhd, zrhop, gdept_n(:,:,:) )                 ! now in situ and potential density 
    127       zrhop(:,:,jpk) = 0._wp 
    128       CALL iom_put( 'rhop', zrhop ) 
    129       ! 
    130       zbotpres(:,:) = 0._wp                        ! no atmospheric surface pressure, levitating sea-ice 
     136         !                                         ! steric sea surface height 
     137         CALL eos( tsn, zrhd, zrhop, gdept_n(:,:,:) )                 ! now in situ and potential density 
     138         zrhop(:,:,jpk) = 0._wp 
     139         CALL iom_put( 'rhop', zrhop ) 
     140         ! 
     141         zbotpres(:,:) = 0._wp                        ! no atmospheric surface pressure, levitating sea-ice 
     142         DO jk = 1, jpkm1 
     143            zbotpres(:,:) = zbotpres(:,:) + e3t_n(:,:,jk) * zrhd(:,:,jk) 
     144         END DO 
     145         IF( ln_linssh ) THEN 
     146            IF ( ln_isfcav ) THEN 
     147               DO ji = 1,jpi 
     148                  DO jj = 1,jpj 
     149                     zbotpres(ji,jj) = zbotpres(ji,jj) + sshn(ji,jj) * zrhd(ji,jj,mikt(ji,jj)) + riceload(ji,jj) 
     150                  END DO 
     151               END DO 
     152            ELSE 
     153               zbotpres(:,:) = zbotpres(:,:) + sshn(:,:) * zrhd(:,:,1) 
     154            END IF 
     155         END IF 
     156         !     
     157         zarho = SUM( area(:,:) * zbotpres(:,:) )  
     158         IF( lk_mpp )   CALL mpp_sum( zarho ) 
     159         zssh_steric = - zarho / area_tot 
     160         CALL iom_put( 'sshsteric', zssh_steric ) 
     161       
     162         !                                         ! ocean bottom pressure 
     163         zztmp = rau0 * grav * 1.e-4_wp               ! recover pressure from pressure anomaly and cover to dbar = 1.e4 Pa 
     164         zbotpres(:,:) = zztmp * ( zbotpres(:,:) + sshn(:,:) + thick0(:,:) ) 
     165         CALL iom_put( 'botpres', zbotpres ) 
     166         ! 
     167      ENDIF 
     168 
     169      IF( iom_use( 'masstot' ) .OR. iom_use( 'temptot' )  .OR. iom_use( 'saltot' )  ) THEN     
     170         !                                         ! Mean density anomalie, temperature and salinity 
     171         ztemp = 0._wp 
     172         zsal  = 0._wp 
     173         DO jk = 1, jpkm1 
     174            DO jj = 1, jpj 
     175               DO ji = 1, jpi 
     176                  zztmp = area(ji,jj) * e3t_n(ji,jj,jk) 
     177                  ztemp = ztemp + zztmp * tsn(ji,jj,jk,jp_tem) 
     178                  zsal  = zsal  + zztmp * tsn(ji,jj,jk,jp_sal) 
     179               END DO 
     180            END DO 
     181         END DO 
     182         IF( ln_linssh ) THEN 
     183            IF( ln_isfcav ) THEN 
     184               DO ji = 1, jpi 
     185                  DO jj = 1, jpj 
     186                     ztemp = ztemp + zarea_ssh(ji,jj) * tsn(ji,jj,mikt(ji,jj),jp_tem)  
     187                     zsal  = zsal  + zarea_ssh(ji,jj) * tsn(ji,jj,mikt(ji,jj),jp_sal)  
     188                  END DO 
     189               END DO 
     190            ELSE 
     191               ztemp = ztemp + SUM( zarea_ssh(:,:) * tsn(:,:,1,jp_tem) ) 
     192               zsal  = zsal  + SUM( zarea_ssh(:,:) * tsn(:,:,1,jp_sal) ) 
     193            END IF 
     194         ENDIF 
     195         IF( lk_mpp ) THEN   
     196            CALL mpp_sum( ztemp ) 
     197            CALL mpp_sum( zsal  ) 
     198         END IF 
     199         ! 
     200         zmass = rau0 * ( zarho + zvol )                 ! total mass of liquid seawater 
     201         ztemp = ztemp / zvol                            ! potential temperature in liquid seawater 
     202         zsal  = zsal  / zvol                            ! Salinity of liquid seawater 
     203         ! 
     204         CALL iom_put( 'masstot', zmass ) 
     205         CALL iom_put( 'temptot', ztemp ) 
     206         CALL iom_put( 'saltot' , zsal  ) 
     207         ! 
     208      ENDIF 
     209 
     210      IF( iom_use( 'tnpeo' )) THEN     
     211      ! Work done against stratification by vertical mixing 
     212      ! Exclude points where rn2 is negative as convection kicks in here and 
     213      ! work is not being done against stratification 
     214          CALL wrk_alloc( jpi, jpj, zpe ) 
     215          zpe(:,:) = 0._wp 
     216          IF( lk_zdfddm ) THEN 
     217             DO ji=1,jpi 
     218                DO jj=1,jpj 
     219                   DO jk=1,jpk 
     220                      zrw =   ( gdepw_n(ji,jj,jk  ) - gdept_n(ji,jj,jk) )   & 
     221                         &  / ( gdept_n(ji,jj,jk-1) - gdept_n(ji,jj,jk) ) 
     222                      ! 
     223                      zaw = rab_n(ji,jj,jk,jp_tem) * (1. - zrw) + rab_n(ji,jj,jk-1,jp_tem)* zrw 
     224                      zbw = rab_n(ji,jj,jk,jp_sal) * (1. - zrw) + rab_n(ji,jj,jk-1,jp_sal)* zrw 
     225                      ! 
     226                      zpe(ji, jj) = zpe(ji, jj) - MIN(0._wp, rn2(ji,jj,jk)) * & 
     227                           &       grav * (avt(ji,jj,jk) * zaw * (tsn(ji,jj,jk-1,jp_tem) - tsn(ji,jj,jk,jp_tem) )  & 
     228                           &       - fsavs(ji,jj,jk) * zbw * (tsn(ji,jj,jk-1,jp_sal) - tsn(ji,jj,jk,jp_sal) ) ) 
     229 
     230                   ENDDO 
     231                ENDDO 
     232             ENDDO 
     233          ELSE 
     234             DO ji = 1, jpi 
     235                DO jj = 1, jpj 
     236                   DO jk = 1, jpk 
     237                       zpe(ji,jj) = zpe(ji,jj) + avt(ji, jj, jk) * MIN(0._wp,rn2(ji, jj, jk)) * rau0 * e3w_n(ji, jj, jk) 
     238                   ENDDO 
     239                ENDDO 
     240             ENDDO 
     241          ENDIF 
     242          CALL lbc_lnk( zpe, 'T', 1._wp)          
     243          CALL iom_put( 'tnpeo', zpe ) 
     244          CALL wrk_dealloc( jpi, jpj, zpe ) 
     245      ENDIF 
     246      ! 
     247      IF( l_ar5 ) THEN 
     248        CALL wrk_dealloc( jpi , jpj              , zarea_ssh , zbotpres ) 
     249        CALL wrk_dealloc( jpi , jpj , jpk        , zrhd      , zrhop    ) 
     250        CALL wrk_dealloc( jpi , jpj , jpk , jpts , ztsn                 ) 
     251      ENDIF 
     252      ! 
     253      IF( nn_timing == 1 )   CALL timing_stop('dia_ar5') 
     254      ! 
     255   END SUBROUTINE dia_ar5 
     256 
     257   SUBROUTINE dia_ar5_hst( ktra, cptr, pua, pva )  
     258      !!---------------------------------------------------------------------- 
     259      !!                    ***  ROUTINE dia_ar5_htr *** 
     260      !!---------------------------------------------------------------------- 
     261      !! Wrapper for heat transport calculations 
     262      !! Called from all advection and/or diffusion routines 
     263      !!---------------------------------------------------------------------- 
     264      INTEGER                         , INTENT(in )  :: ktra  ! tracer index 
     265      CHARACTER(len=3)                , INTENT(in)   :: cptr  ! transport type  'adv'/'ldf' 
     266      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in)   :: pua   ! 3D input array of advection/diffusion 
     267      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in)   :: pva   ! 3D input array of advection/diffusion 
     268      ! 
     269      INTEGER    ::  ji, jj, jk 
     270      REAL(wp), POINTER, DIMENSION(:,:)  :: z2d 
     271 
     272     
     273 
     274      CALL wrk_alloc( jpi, jpj, z2d ) 
     275      z2d(:,:) = pua(:,:,1)  
    131276      DO jk = 1, jpkm1 
    132          zbotpres(:,:) = zbotpres(:,:) + e3t_n(:,:,jk) * zrhd(:,:,jk) 
    133       END DO 
    134       IF( ln_linssh ) THEN 
    135          IF ( ln_isfcav ) THEN 
    136             DO ji=1,jpi 
    137                DO jj=1,jpj 
    138                   zbotpres(ji,jj) = zbotpres(ji,jj) + sshn(ji,jj) * zrhd(ji,jj,mikt(ji,jj)) + riceload(ji,jj) 
    139                END DO 
     277         DO jj = 2, jpjm1 
     278            DO ji = fs_2, fs_jpim1   ! vector opt. 
     279               z2d(ji,jj) = z2d(ji,jj) + pua(ji,jj,jk)  
    140280            END DO 
    141          ELSE 
    142             zbotpres(:,:) = zbotpres(:,:) + sshn(:,:) * zrhd(:,:,1) 
    143          END IF 
    144       END IF 
    145       !     
    146       zarho = SUM( area(:,:) * zbotpres(:,:) )  
    147       IF( lk_mpp )   CALL mpp_sum( zarho ) 
    148       zssh_steric = - zarho / area_tot 
    149       CALL iom_put( 'sshsteric', zssh_steric ) 
    150        
    151       !                                         ! ocean bottom pressure 
    152       zztmp = rau0 * grav * 1.e-4_wp               ! recover pressure from pressure anomaly and cover to dbar = 1.e4 Pa 
    153       zbotpres(:,:) = zztmp * ( zbotpres(:,:) + sshn(:,:) + thick0(:,:) ) 
    154       CALL iom_put( 'botpres', zbotpres ) 
    155  
    156       !                                         ! Mean density anomalie, temperature and salinity 
    157       ztemp = 0._wp 
    158       zsal  = 0._wp 
    159       DO jk = 1, jpkm1 
    160          DO jj = 1, jpj 
    161             DO ji = 1, jpi 
    162                zztmp = area(ji,jj) * e3t_n(ji,jj,jk) 
    163                ztemp = ztemp + zztmp * tsn(ji,jj,jk,jp_tem) 
    164                zsal  = zsal  + zztmp * tsn(ji,jj,jk,jp_sal) 
    165             END DO 
    166          END DO 
    167       END DO 
    168       IF( ln_linssh ) THEN 
    169          IF( ln_isfcav ) THEN 
    170             DO ji=1,jpi 
    171                DO jj=1,jpj 
    172                   ztemp = ztemp + zarea_ssh(ji,jj) * tsn(ji,jj,mikt(ji,jj),jp_tem)  
    173                   zsal  = zsal  + zarea_ssh(ji,jj) * tsn(ji,jj,mikt(ji,jj),jp_sal)  
    174                END DO 
    175             END DO 
    176          ELSE 
    177             ztemp = ztemp + SUM( zarea_ssh(:,:) * tsn(:,:,1,jp_tem) ) 
    178             zsal  = zsal  + SUM( zarea_ssh(:,:) * tsn(:,:,1,jp_sal) ) 
    179          END IF 
    180       ENDIF 
    181       IF( lk_mpp ) THEN   
    182          CALL mpp_sum( ztemp ) 
    183          CALL mpp_sum( zsal  ) 
    184       END IF 
    185       ! 
    186       zmass = rau0 * ( zarho + zvol )                 ! total mass of liquid seawater 
    187       ztemp = ztemp / zvol                            ! potential temperature in liquid seawater 
    188       zsal  = zsal  / zvol                            ! Salinity of liquid seawater 
    189       ! 
    190       CALL iom_put( 'masstot', zmass ) 
    191       CALL iom_put( 'temptot', ztemp ) 
    192       CALL iom_put( 'saltot' , zsal  ) 
    193       ! 
    194       CALL wrk_dealloc( jpi , jpj              , zarea_ssh , zbotpres ) 
    195       CALL wrk_dealloc( jpi , jpj , jpk        , zrhd      , zrhop    ) 
    196       CALL wrk_dealloc( jpi , jpj , jpk , jpts , ztsn                 ) 
    197       ! 
    198       IF( nn_timing == 1 )   CALL timing_stop('dia_ar5') 
    199       ! 
    200    END SUBROUTINE dia_ar5 
     281         END DO 
     282       END DO 
     283       CALL lbc_lnk( z2d, 'U', -1. ) 
     284       IF( cptr == 'adv' ) THEN 
     285          IF( ktra == jp_tem ) CALL iom_put( "uadv_heattr" , rau0_rcp * z2d )  ! advective heat transport in i-direction 
     286          IF( ktra == jp_sal ) CALL iom_put( "uadv_salttr" , rau0     * z2d )  ! advective salt transport in i-direction 
     287       ENDIF 
     288       IF( cptr == 'ldf' ) THEN 
     289          IF( ktra == jp_tem ) CALL iom_put( "udiff_heattr" , rau0_rcp * z2d ) ! diffusive heat transport in i-direction 
     290          IF( ktra == jp_sal ) CALL iom_put( "udiff_salttr" , rau0     * z2d ) ! diffusive salt transport in i-direction 
     291       ENDIF 
     292       ! 
     293       z2d(:,:) = pva(:,:,1)  
     294       DO jk = 1, jpkm1 
     295          DO jj = 2, jpjm1 
     296             DO ji = fs_2, fs_jpim1   ! vector opt. 
     297                z2d(ji,jj) = z2d(ji,jj) + pva(ji,jj,jk)  
     298             END DO 
     299          END DO 
     300       END DO 
     301       CALL lbc_lnk( z2d, 'V', -1. ) 
     302       IF( cptr == 'adv' ) THEN 
     303          IF( ktra == jp_tem ) CALL iom_put( "vadv_heattr" , rau0_rcp * z2d )  ! advective heat transport in j-direction 
     304          IF( ktra == jp_sal ) CALL iom_put( "vadv_salttr" , rau0     * z2d )  ! advective salt transport in j-direction 
     305       ENDIF 
     306       IF( cptr == 'ldf' ) THEN 
     307          IF( ktra == jp_tem ) CALL iom_put( "vdiff_heattr" , rau0_rcp * z2d ) ! diffusive heat transport in j-direction 
     308          IF( ktra == jp_sal ) CALL iom_put( "vdiff_salttr" , rau0     * z2d ) ! diffusive salt transport in j-direction 
     309       ENDIF 
     310           
     311       CALL wrk_dealloc( jpi, jpj, z2d ) 
     312 
     313   END SUBROUTINE dia_ar5_hst 
    201314 
    202315 
    … …  
    217330      IF( nn_timing == 1 )   CALL timing_start('dia_ar5_init') 
    218331      ! 
    219       CALL wrk_alloc( jpi , jpj , jpk, jpts, zsaldta ) 
    220       !                                      ! allocate dia_ar5 arrays 
    221       IF( dia_ar5_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'dia_ar5_init : unable to allocate arrays' ) 
    222  
    223       area(:,:) = e1e2t(:,:) * tmask_i(:,:) 
    224  
    225       area_tot = SUM( area(:,:) )   ;   IF( lk_mpp )   CALL mpp_sum( area_tot ) 
    226  
    227       vol0        = 0._wp 
    228       thick0(:,:) = 0._wp 
    229       DO jk = 1, jpkm1 
    230          vol0        = vol0        + SUM( area (:,:) * tmask(:,:,jk) * e3t_0(:,:,jk) ) 
    231          thick0(:,:) = thick0(:,:) +    tmask_i(:,:) * tmask(:,:,jk) * e3t_0(:,:,jk) 
    232       END DO 
    233       IF( lk_mpp )   CALL mpp_sum( vol0 ) 
    234  
    235  
    236       CALL iom_open ( 'sali_ref_clim_monthly', inum ) 
    237       CALL iom_get  ( inum, jpdom_data, 'vosaline' , zsaldta(:,:,:,1), 1  ) 
    238       CALL iom_get  ( inum, jpdom_data, 'vosaline' , zsaldta(:,:,:,2), 12 ) 
    239       CALL iom_close( inum ) 
    240  
    241       sn0(:,:,:) = 0.5_wp * ( zsaldta(:,:,:,1) + zsaldta(:,:,:,2) )         
    242       sn0(:,:,:) = sn0(:,:,:) * tmask(:,:,:) 
    243       IF( ln_zps ) THEN               ! z-coord. partial steps 
    244          DO jj = 1, jpj               ! interpolation of salinity at the last ocean level (i.e. the partial step) 
    245             DO ji = 1, jpi 
    246                ik = mbkt(ji,jj) 
    247                IF( ik > 1 ) THEN 
    248                   zztmp = ( gdept_1d(ik) - gdept_0(ji,jj,ik) ) / ( gdept_1d(ik) - gdept_1d(ik-1) ) 
    249                   sn0(ji,jj,ik) = ( 1._wp - zztmp ) * sn0(ji,jj,ik) + zztmp * sn0(ji,jj,ik-1) 
    250                ENDIF 
     332      l_ar5 = .FALSE. 
     333      IF(   iom_use( 'voltot'  ) .OR. iom_use( 'sshtot'    )  .OR. iom_use( 'sshdyn' )  .OR.  &  
     334         &  iom_use( 'masstot' ) .OR. iom_use( 'temptot'   )  .OR. iom_use( 'saltot' ) .OR.  &     
     335         &  iom_use( 'botpres' ) .OR. iom_use( 'sshthster' )  .OR. iom_use( 'sshsteric' )  ) L_ar5 = .TRUE. 
     336   
     337      IF( l_ar5 ) THEN 
     338         ! 
     339         CALL wrk_alloc( jpi , jpj , jpk, jpts, zsaldta ) 
     340         !                                      ! allocate dia_ar5 arrays 
     341         IF( dia_ar5_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'dia_ar5_init : unable to allocate arrays' ) 
     342 
     343         area(:,:) = e1e2t(:,:) * tmask_i(:,:) 
     344 
     345         area_tot = SUM( area(:,:) )   ;   IF( lk_mpp )   CALL mpp_sum( area_tot ) 
     346 
     347         vol0        = 0._wp 
     348         thick0(:,:) = 0._wp 
     349         DO jk = 1, jpkm1 
     350            vol0        = vol0        + SUM( area (:,:) * tmask(:,:,jk) * e3t_0(:,:,jk) ) 
     351            thick0(:,:) = thick0(:,:) +    tmask_i(:,:) * tmask(:,:,jk) * e3t_0(:,:,jk) 
     352         END DO 
     353         IF( lk_mpp )   CALL mpp_sum( vol0 ) 
     354 
     355         CALL iom_open ( 'sali_ref_clim_monthly', inum ) 
     356         CALL iom_get  ( inum, jpdom_data, 'vosaline' , zsaldta(:,:,:,1), 1  ) 
     357         CALL iom_get  ( inum, jpdom_data, 'vosaline' , zsaldta(:,:,:,2), 12 ) 
     358         CALL iom_close( inum ) 
     359 
     360         sn0(:,:,:) = 0.5_wp * ( zsaldta(:,:,:,1) + zsaldta(:,:,:,2) )         
     361         sn0(:,:,:) = sn0(:,:,:) * tmask(:,:,:) 
     362         IF( ln_zps ) THEN               ! z-coord. partial steps 
     363            DO jj = 1, jpj               ! interpolation of salinity at the last ocean level (i.e. the partial step) 
     364               DO ji = 1, jpi 
     365                  ik = mbkt(ji,jj) 
     366                  IF( ik > 1 ) THEN 
     367                     zztmp = ( gdept_1d(ik) - gdept_0(ji,jj,ik) ) / ( gdept_1d(ik) - gdept_1d(ik-1) ) 
     368                     sn0(ji,jj,ik) = ( 1._wp - zztmp ) * sn0(ji,jj,ik) + zztmp * sn0(ji,jj,ik-1) 
     369                  ENDIF 
     370               END DO 
    251371            END DO 
    252          END DO 
    253       ENDIF 
    254       ! 
    255       CALL wrk_dealloc( jpi , jpj , jpk, jpts, zsaldta ) 
     372         ENDIF 
     373         ! 
     374         CALL wrk_dealloc( jpi , jpj , jpk, jpts, zsaldta ) 
     375         ! 
     376      ENDIF 
    256377      ! 
    257378      IF( nn_timing == 1 )   CALL timing_stop('dia_ar5_init') 
    258379      ! 
    259380   END SUBROUTINE dia_ar5_init 
    260  
    261 #else 
    262    !!---------------------------------------------------------------------- 
    263    !!   Default option :                                         NO diaar5 
    264    !!---------------------------------------------------------------------- 
    265    LOGICAL, PUBLIC, PARAMETER :: lk_diaar5 = .FALSE.   ! coupled flag 
    266 CONTAINS 
    267    SUBROUTINE dia_ar5_init    ! Dummy routine 
    268    END SUBROUTINE dia_ar5_init 
    269    SUBROUTINE dia_ar5( kt )   ! Empty routine 
    270       INTEGER ::   kt 
    271       WRITE(*,*) 'dia_ar5: You should not have seen this print! error?', kt 
    272    END SUBROUTINE dia_ar5 
    273 #endif 
    274381 
    275382   !!====================================================================== 

  • trunk/NEMOGCM/NEMO/OPA_SRC/DIA/diadct.F90

    r6981 r7646  
    392392        ENDIF                
    393393 
    394         IF( iptglo .NE. 0 )THEN 
     394        IF( iptglo /= 0 )THEN 
    395395              
    396396           !read points'coordinates and directions  
    … …  
    399399           directemp(:) = 0                  !value of directions of each points 
    400400           DO jpt=1,iptglo 
    401               READ(numdct_in)i1,i2 
     401              READ(numdct_in) i1, i2 
    402402              coordtemp(jpt)%I = i1  
    403403              coordtemp(jpt)%J = i2 
    404404           ENDDO 
    405            READ(numdct_in)directemp(1:iptglo) 
     405           READ(numdct_in) directemp(1:iptglo) 
    406406     
    407407           !debug 
    … …  
    416416           !Now each proc selects only points that are in its domain: 
    417417           !-------------------------------------------------------- 
    418            iptloc = 0                    !initialize number of points selected 
    419            DO jpt=1,iptglo               !loop on listpoint read in the file 
    420                      
     418           iptloc = 0                    ! initialize number of points selected 
     419           DO jpt = 1, iptglo            ! loop on listpoint read in the file 
     420              !       
    421421              iiglo=coordtemp(jpt)%I          ! global coordinates of the point 
    422422              ijglo=coordtemp(jpt)%J          !  "  
    423423 
    424               IF( iiglo==jpidta .AND. nimpp==1 ) iiglo = 2 
    425  
    426               iiloc=iiglo-jpizoom+1-nimpp+1   ! local coordinates of the point 
    427               ijloc=ijglo-jpjzoom+1-njmpp+1   !  " 
     424              IF( iiglo==jpiglo .AND. nimpp==1 )   iiglo = 2         !!gm BUG: Hard coded periodicity ! 
     425 
     426              iiloc=iiglo-nimpp+1   ! local coordinates of the point 
     427              ijloc=ijglo-njmpp+1   !  " 
    428428 
    429429              !verify if the point is on the local domain:(1,nlei)*(1,nlej) 
    430               IF( iiloc .GE. 1 .AND. iiloc .LE. nlei .AND. & 
    431                   ijloc .GE. 1 .AND. ijloc .LE. nlej       )THEN 
     430              IF( iiloc >= 1 .AND. iiloc <= nlei .AND. & 
     431                  ijloc >= 1 .AND. ijloc <= nlej       )THEN 
    432432                 iptloc = iptloc + 1                                                 ! count local points 
    433433                 secs(jsec)%listPoint(iptloc) = POINT_SECTION(mi0(iiglo),mj0(ijglo)) ! store local coordinates 
    434434                 secs(jsec)%direction(iptloc) = directemp(jpt)                       ! store local direction 
    435435              ENDIF 
    436  
    437            ENDDO 
     436              ! 
     437           END DO 
    438438      
    439439           secs(jsec)%nb_point=iptloc !store number of section's points 
    … …  
    444444              WRITE(numout,*)"      List of points selected by the proc:" 
    445445              DO jpt = 1,iptloc 
    446                  iiglo = secs(jsec)%listPoint(jpt)%I + jpizoom - 1 + nimpp - 1 
    447                  ijglo = secs(jsec)%listPoint(jpt)%J + jpjzoom - 1 + njmpp - 1 
     446                 iiglo = secs(jsec)%listPoint(jpt)%I + nimpp - 1 
     447                 ijglo = secs(jsec)%listPoint(jpt)%J + njmpp - 1 
    448448                 WRITE(numout,*)'         # I J : ',iiglo,ijglo 
    449449              ENDDO 
    … …  
    452452              IF(jsec==nn_secdebug .AND. secs(jsec)%nb_point .NE. 0)THEN 
    453453              DO jpt = 1,iptloc 
    454                  iiglo = secs(jsec)%listPoint(jpt)%I + jpizoom - 1 + nimpp - 1 
    455                  ijglo = secs(jsec)%listPoint(jpt)%J + jpjzoom - 1 + njmpp - 1 
     454                 iiglo = secs(jsec)%listPoint(jpt)%I + nimpp - 1 
     455                 ijglo = secs(jsec)%listPoint(jpt)%J + njmpp - 1 
    456456              ENDDO 
    457457              ENDIF 
    … …  
    468468           IF(jsec==nn_secdebug .AND. secs(jsec)%nb_point .NE. 0)THEN 
    469469              DO jpt = 1,secs(jsec)%nb_point 
    470                  iiglo = secs(jsec)%listPoint(jpt)%I + jpizoom - 1 + nimpp - 1 
    471                  ijglo = secs(jsec)%listPoint(jpt)%J + jpjzoom - 1 + njmpp - 1 
     470                 iiglo = secs(jsec)%listPoint(jpt)%I + nimpp - 1 
     471                 ijglo = secs(jsec)%listPoint(jpt)%J + njmpp - 1 
    472472              ENDDO 
    473473           ENDIF 
    … …  
    479479              iptloc = secs(jsec)%nb_point 
    480480              DO jpt = 1,iptloc 
    481                  iiglo = secs(jsec)%listPoint(jpt)%I + jpizoom - 1 + nimpp - 1 
    482                  ijglo = secs(jsec)%listPoint(jpt)%J + jpjzoom - 1 + njmpp - 1 
     481                 iiglo = secs(jsec)%listPoint(jpt)%I + nimpp - 1 
     482                 ijglo = secs(jsec)%listPoint(jpt)%J + njmpp - 1 
    483483                 WRITE(numout,*)'         # I J : ',iiglo,ijglo 
    484484                 CALL FLUSH(numout) 

  • trunk/NEMOGCM/NEMO/OPA_SRC/DIA/diaharm.F90

    r6140 r7646  
    66   !! History :  3.1  !  2007  (O. Le Galloudec, J. Chanut)  Original code 
    77   !!---------------------------------------------------------------------- 
    8 #if defined key_diaharm && defined key_tide 
     8#if defined key_diaharm 
    99   !!---------------------------------------------------------------------- 
    1010   !!   'key_diaharm' 
    11    !!   'key_tide' 
    1211   !!---------------------------------------------------------------------- 
    1312   USE oce             ! ocean dynamics and tracers variables 
    … …  
    1615   USE daymod 
    1716   USE tide_mod 
     17   USE sbctide         ! Tidal forcing or not 
    1818   ! 
    1919   USE in_out_manager  ! I/O units 
    … …  
    8282         WRITE(numout,*) '~~~~~~~ ' 
    8383      ENDIF 
     84      ! 
     85      IF( .NOT. ln_tide )   CALL ctl_stop( 'dia_harm_init : ln_tide must be true for harmonic analysis') 
    8486      ! 
    8587      CALL tide_init_Wave 

  • trunk/NEMOGCM/NEMO/OPA_SRC/DIA/diahsb.F90

    r6140 r7646  
    2323   USE trabbc          ! bottom boundary condition  
    2424   USE trabbc          ! bottom boundary condition 
    25    USE bdy_par         ! (for lk_bdy) 
    2625   USE restart         ! ocean restart 
     26   USE bdy_oce   , ONLY: ln_bdy 
    2727   ! 
    2828   USE iom             ! I/O manager 
    … …  
    3838   PUBLIC   dia_hsb        ! routine called by step.F90 
    3939   PUBLIC   dia_hsb_init   ! routine called by nemogcm.F90 
    40    PUBLIC   dia_hsb_rst    ! routine called by step.F90 
    4140 
    4241   LOGICAL, PUBLIC ::   ln_diahsb   !: check the heat and salt budgets 
    … …  
    8685      !!--------------------------------------------------------------------------- 
    8786      IF( nn_timing == 1 )   CALL timing_start('dia_hsb')       
     87      ! 
    8888      CALL wrk_alloc( jpi,jpj,   z2d0, z2d1 ) 
    8989      ! 
    … …  
    171171      ENDDO 
    172172 
    173       ! Substract forcing from heat content, salt content and volume variations 
     173      ! ------------------------ ! 
     174      ! 3 -  Drifts              ! 
     175      ! ------------------------ ! 
    174176      zdiff_v1 = zdiff_v1 - frc_v 
    175177      IF( .NOT.ln_linssh )   zdiff_v2 = zdiff_v2 - frc_v 
    … …  
    184186 
    185187      ! ----------------------- ! 
    186       ! 3 - Diagnostics writing ! 
     188      ! 4 - Diagnostics writing ! 
    187189      ! ----------------------- ! 
    188190      zvol_tot = 0._wp                    ! total ocean volume (calculated with scale factors) 
    … …  
    197199!!gm end 
    198200 
    199       IF( ln_linssh ) THEN 
    200         CALL iom_put( 'bgtemper' , zdiff_hc1 / zvol_tot)              ! Heat content variation (C)  
    201         CALL iom_put( 'bgsaline' , zdiff_sc1 / zvol_tot)              ! Salt content variation (psu) 
    202         CALL iom_put( 'bgheatco' , zdiff_hc1 * 1.e-20 * rau0 * rcp )  ! Heat content variation (1.e20 J)  
    203         CALL iom_put( 'bgsaltco' , zdiff_sc1 * 1.e-9   )              ! Salt content variation (psu*km3) 
    204         CALL iom_put( 'bgvolssh' , zdiff_v1  * 1.e-9   )              ! volume ssh variation (km3)   
    205         CALL iom_put( 'bgfrcvol' , frc_v    * 1.e-9    )              ! vol - surface forcing (km3)  
    206         CALL iom_put( 'bgfrctem' , frc_t / zvol_tot    )              ! hc  - surface forcing (C)  
    207         CALL iom_put( 'bgfrcsal' , frc_s / zvol_tot    )              ! sc  - surface forcing (psu)  
     201      CALL iom_put(   'bgfrcvol' , frc_v    * 1.e-9    )              ! vol - surface forcing (km3)  
     202      CALL iom_put(   'bgfrctem' , frc_t    * rau0 * rcp * 1.e-20 )   ! hc  - surface forcing (1.e20 J)  
     203      CALL iom_put(   'bgfrchfx' , frc_t    * rau0 * rcp /  &         ! hc  - surface forcing (W/m2)  
     204         &                       ( surf_tot * kt * rdt )        ) 
     205      CALL iom_put(   'bgfrcsal' , frc_s    * 1.e-9    )              ! sc  - surface forcing (psu*km3)  
     206 
     207      IF( .NOT. ln_linssh ) THEN 
     208        CALL iom_put( 'bgtemper' , zdiff_hc / zvol_tot )              ! Temperature drift     (C)  
     209        CALL iom_put( 'bgsaline' , zdiff_sc / zvol_tot )              ! Salinity    drift     (pss) 
     210        CALL iom_put( 'bgheatco' , zdiff_hc * 1.e-20 * rau0 * rcp )   ! Heat content drift    (1.e20 J)  
     211        CALL iom_put( 'bgheatfx' , zdiff_hc * rau0 * rcp /  &         ! Heat flux drift       (W/m2)  
     212           &                       ( surf_tot * kt * rdt )        ) 
     213        CALL iom_put( 'bgsaltco' , zdiff_sc * 1.e-9    )              ! Salt content drift    (psu*km3) 
     214        CALL iom_put( 'bgvolssh' , zdiff_v1 * 1.e-9    )              ! volume ssh drift      (km3)   
     215        CALL iom_put( 'bgvole3t' , zdiff_v2 * 1.e-9    )              ! volume e3t drift      (km3)   
     216      ELSE 
     217        CALL iom_put( 'bgtemper' , zdiff_hc1 / zvol_tot)              ! Heat content drift    (C)  
     218        CALL iom_put( 'bgsaline' , zdiff_sc1 / zvol_tot)              ! Salt content drift    (pss) 
     219        CALL iom_put( 'bgheatco' , zdiff_hc1 * 1.e-20 * rau0 * rcp )  ! Heat content drift    (1.e20 J)  
     220        CALL iom_put( 'bgheatfx' , zdiff_hc1 * rau0 * rcp /  &        ! Heat flux drift       (W/m2)  
     221           &                       ( surf_tot * kt * rdt )         ) 
     222        CALL iom_put( 'bgsaltco' , zdiff_sc1 * 1.e-9    )             ! Salt content drift    (psu*km3) 
     223        CALL iom_put( 'bgvolssh' , zdiff_v1 * 1.e-9    )              ! volume ssh drift      (km3)   
    208224        CALL iom_put( 'bgmistem' , zerr_hc1 / zvol_tot )              ! hc  - error due to free surface (C) 
    209225        CALL iom_put( 'bgmissal' , zerr_sc1 / zvol_tot )              ! sc  - error due to free surface (psu) 
    210       ELSE 
    211         CALL iom_put( 'bgtemper' , zdiff_hc / zvol_tot )              ! Temperature variation (C)  
    212         CALL iom_put( 'bgsaline' , zdiff_sc / zvol_tot )              ! Salinity    variation (psu) 
    213         CALL iom_put( 'bgheatco' , zdiff_hc * 1.e-20 * rau0 * rcp )   ! Heat content variation (1.e20 J)  
    214         CALL iom_put( 'bgsaltco' , zdiff_sc * 1.e-9    )              ! Salt content variation (psu*km3) 
    215         CALL iom_put( 'bgvolssh' , zdiff_v1 * 1.e-9    )              ! volume ssh variation (km3)   
    216         CALL iom_put( 'bgvole3t' , zdiff_v2 * 1.e-9    )              ! volume e3t variation (km3)   
    217         CALL iom_put( 'bgfrcvol' , frc_v    * 1.e-9    )              ! vol - surface forcing (km3)  
    218         CALL iom_put( 'bgfrctem' , frc_t / zvol_tot    )              ! hc  - surface forcing (C)  
    219         CALL iom_put( 'bgfrcsal' , frc_s / zvol_tot    )              ! sc  - surface forcing (psu)  
    220226      ENDIF 
    221227      ! 
    … …  
    231237   SUBROUTINE dia_hsb_rst( kt, cdrw ) 
    232238     !!--------------------------------------------------------------------- 
    233      !!                   ***  ROUTINE limdia_rst  *** 
     239     !!                   ***  ROUTINE dia_hsb_rst  *** 
    234240     !!                      
    235241     !! ** Purpose :   Read or write DIA file in restart file 
    … …  
    241247     ! 
    242248     INTEGER ::   ji, jj, jk   ! dummy loop indices 
    243      INTEGER ::   id1          ! local integers 
    244249     !!---------------------------------------------------------------------- 
    245250     ! 
    246251     IF( TRIM(cdrw) == 'READ' ) THEN        ! Read/initialise  
    247252        IF( ln_rstart ) THEN                   !* Read the restart file 
    248            !id1 = iom_varid( numror, 'frc_vol'  , ldstop = .FALSE. ) 
    249253           ! 
    250254           IF(lwp) WRITE(numout,*) '~~~~~~~' 
    … …  
    259263           ENDIF 
    260264           CALL iom_get( numror, jpdom_autoglo, 'surf_ini', surf_ini ) ! ice sheet coupling 
    261            CALL iom_get( numror, jpdom_autoglo, 'ssh_ini', ssh_ini ) 
    262            CALL iom_get( numror, jpdom_autoglo, 'e3t_ini', e3t_ini ) 
    263            CALL iom_get( numror, jpdom_autoglo, 'hc_loc_ini', hc_loc_ini ) 
    264            CALL iom_get( numror, jpdom_autoglo, 'sc_loc_ini', sc_loc_ini ) 
     265           CALL iom_get( numror, jpdom_autoglo, 'ssh_ini', ssh_ini(:,:) ) 
     266           CALL iom_get( numror, jpdom_autoglo, 'e3t_ini', e3t_ini(:,:,:) ) 
     267           CALL iom_get( numror, jpdom_autoglo, 'hc_loc_ini', hc_loc_ini(:,:,:) ) 
     268           CALL iom_get( numror, jpdom_autoglo, 'sc_loc_ini', sc_loc_ini(:,:,:) ) 
    265269           IF( ln_linssh ) THEN 
    266               CALL iom_get( numror, jpdom_autoglo, 'ssh_hc_loc_ini', ssh_hc_loc_ini ) 
    267               CALL iom_get( numror, jpdom_autoglo, 'ssh_sc_loc_ini', ssh_sc_loc_ini ) 
     270              CALL iom_get( numror, jpdom_autoglo, 'ssh_hc_loc_ini', ssh_hc_loc_ini(:,:) ) 
     271              CALL iom_get( numror, jpdom_autoglo, 'ssh_sc_loc_ini', ssh_sc_loc_ini(:,:) ) 
    268272           ENDIF 
    269273       ELSE 
    … …  
    313317        ENDIF 
    314318        CALL iom_rstput( kt, nitrst, numrow, 'surf_ini', surf_ini )      ! ice sheet coupling 
    315         CALL iom_rstput( kt, nitrst, numrow, 'ssh_ini', ssh_ini ) 
    316         CALL iom_rstput( kt, nitrst, numrow, 'e3t_ini', e3t_ini ) 
    317         CALL iom_rstput( kt, nitrst, numrow, 'hc_loc_ini', hc_loc_ini ) 
    318         CALL iom_rstput( kt, nitrst, numrow, 'sc_loc_ini', sc_loc_ini ) 
     319        CALL iom_rstput( kt, nitrst, numrow, 'ssh_ini', ssh_ini(:,:) ) 
     320        CALL iom_rstput( kt, nitrst, numrow, 'e3t_ini', e3t_ini(:,:,:) ) 
     321        CALL iom_rstput( kt, nitrst, numrow, 'hc_loc_ini', hc_loc_ini(:,:,:) ) 
     322        CALL iom_rstput( kt, nitrst, numrow, 'sc_loc_ini', sc_loc_ini(:,:,:) ) 
    319323        IF( ln_linssh ) THEN 
    320            CALL iom_rstput( kt, nitrst, numrow, 'ssh_hc_loc_ini', ssh_hc_loc_ini ) 
    321            CALL iom_rstput( kt, nitrst, numrow, 'ssh_sc_loc_ini', ssh_sc_loc_ini ) 
     324           CALL iom_rstput( kt, nitrst, numrow, 'ssh_hc_loc_ini', ssh_hc_loc_ini(:,:) ) 
     325           CALL iom_rstput( kt, nitrst, numrow, 'ssh_sc_loc_ini', ssh_sc_loc_ini(:,:) ) 
    322326        ENDIF 
    323327        ! 
    … …  
    339343      !!             - Compute coefficients for conversion 
    340344      !!--------------------------------------------------------------------------- 
    341       INTEGER ::   jk       ! dummy loop indice 
    342345      INTEGER ::   ierror   ! local integer 
    343346      INTEGER ::   ios 
    344       ! 
     347      !! 
    345348      NAMELIST/namhsb/ ln_diahsb 
    346349      !!---------------------------------------------------------------------- 
    347  
    348       IF(lwp) THEN 
    349          WRITE(numout,*) 
    350          WRITE(numout,*) 'dia_hsb_init : check the heat and salt budgets' 
    351          WRITE(numout,*) '~~~~~~~~ ' 
    352       ENDIF 
    353  
     350      ! 
    354351      REWIND( numnam_ref )              ! Namelist namhsb in reference namelist 
    355352      READ  ( numnam_ref, namhsb, IOSTAT = ios, ERR = 901) 
    … …  
    361358      IF(lwm) WRITE ( numond, namhsb ) 
    362359 
    363       ! 
    364       IF(lwp) THEN                   ! Control print 
     360      IF(lwp) THEN 
    365361         WRITE(numout,*) 
    366          WRITE(numout,*) 'dia_hsb_init : check the heat and salt budgets' 
    367          WRITE(numout,*) '~~~~~~~~~~~~' 
    368          WRITE(numout,*) '   Namelist namhsb : set hsb parameters' 
    369          WRITE(numout,*) '      Switch for hsb diagnostic (T) or not (F)  ln_diahsb  = ', ln_diahsb 
    370          WRITE(numout,*) 
    371       ENDIF 
    372  
     362         WRITE(numout,*) 'dia_hsb_init' 
     363         WRITE(numout,*) '~~~~~~~~ ' 
     364         WRITE(numout,*) '  check the heat and salt budgets (T) or not (F)       ln_diahsb = ', ln_diahsb 
     365      ENDIF 
     366      ! 
    373367      IF( .NOT. ln_diahsb )   RETURN 
    374          !      IF( .NOT. lk_mpp_rep ) & 
    375          !        CALL ctl_stop (' Your global mpp_sum if performed in single precision - 64 bits -', & 
    376          !             &         ' whereas the global sum to be precise must be done in double precision ',& 
    377          !             &         ' please add key_mpp_rep') 
    378368 
    379369      ! ------------------- ! 
    … …  
    383373         &      e3t_ini(jpi,jpj,jpk), surf(jpi,jpj),  ssh_ini(jpi,jpj), STAT=ierror  ) 
    384374      IF( ierror > 0 ) THEN 
    385          CALL ctl_stop( 'dia_hsb: unable to allocate hc_loc_ini' )   ;   RETURN 
     375         CALL ctl_stop( 'dia_hsb_init: unable to allocate hc_loc_ini' )   ;   RETURN 
    386376      ENDIF 
    387377 
    388378      IF( ln_linssh )   ALLOCATE( ssh_hc_loc_ini(jpi,jpj), ssh_sc_loc_ini(jpi,jpj),STAT=ierror ) 
    389379      IF( ierror > 0 ) THEN 
    390          CALL ctl_stop( 'dia_hsb: unable to allocate hc_loc_ini' )   ;   RETURN 
     380         CALL ctl_stop( 'dia_hsb: unable to allocate ssh_hc_loc_ini' )   ;   RETURN 
    391381      ENDIF 
    392382 
    … …  
    394384      ! 2 - Time independant variables and file opening ! 
    395385      ! ----------------------------------------------- ! 
    396       IF(lwp) WRITE(numout,*) "dia_hsb: heat salt volume budgets activated" 
    397       IF(lwp) WRITE(numout,*) '~~~~~~~' 
    398386      surf(:,:) = e1t(:,:) * e2t(:,:) * tmask_i(:,:)      ! masked surface grid cell area 
    399       surf_tot  = glob_sum( surf(:,:) )                                       ! total ocean surface area 
    400  
    401       IF( lk_bdy ) CALL ctl_warn( 'dia_hsb does not take open boundary fluxes into account' )          
     387      surf_tot  = glob_sum( surf(:,:) )                   ! total ocean surface area 
     388 
     389      IF( ln_bdy ) CALL ctl_warn( 'dia_hsb_init: heat/salt budget does not consider open boundary fluxes' )          
    402390      ! 
    403391      ! ---------------------------------- ! 

  • trunk/NEMOGCM/NEMO/OPA_SRC/DIA/diaptr.F90

    r6140 r7646  
    99   !!            3.3  ! 2010-10  (G. Madec)  dynamical allocation 
    1010   !!            3.6  ! 2014-12  (C. Ethe) use of IOM 
     11   !!            3.6  ! 2016-06  (T. Graham) Addition of diagnostics for CMIP6 
    1112   !!---------------------------------------------------------------------- 
    1213 
    … …  
    3839   PUBLIC   dia_ptr_init   ! call in step module 
    3940   PUBLIC   dia_ptr        ! call in step module 
     41   PUBLIC   dia_ptr_hst    ! called from tra_ldf/tra_adv routines 
    4042 
    4143   !                                  !!** namelist  namptr  ** 
    42    REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:) ::   htr_adv, htr_ldf   !: Heat TRansports (adv, diff, overturn.) 
    43    REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:) ::   str_adv, str_ldf   !: Salt TRansports (adv, diff, overturn.) 
    44     
     44   REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) ::   htr_adv, htr_ldf, htr_eiv   !: Heat TRansports (adv, diff, Bolus.) 
     45   REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) ::   str_adv, str_ldf, str_eiv   !: Salt TRansports (adv, diff, Bolus.) 
     46   REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) ::   htr_ove, str_ove   !: heat Salt TRansports ( overturn.) 
     47   REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) ::   htr_btr, str_btr   !: heat Salt TRansports ( barotropic ) 
    4548 
    4649   LOGICAL, PUBLIC ::   ln_diaptr   !  Poleward transport flag (T) or not (F) 
    4750   LOGICAL, PUBLIC ::   ln_subbas   !  Atlantic/Pacific/Indian basins calculation 
    48    INTEGER        ::   nptr        ! = 1 (l_subbas=F) or = 5 (glo, atl, pac, ind, ipc) (l_subbas=T)  
     51   INTEGER, PUBLIC ::   nptr        ! = 1 (l_subbas=F) or = 5 (glo, atl, pac, ind, ipc) (l_subbas=T)  
    4952 
    5053   REAL(wp) ::   rc_sv    = 1.e-6_wp   ! conversion from m3/s to Sverdrup 
    … …  
    7578      ! 
    7679      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices 
    77       REAL(wp) ::   zv, zsfc               ! local scalar 
     80      REAL(wp) ::   zsfc,zvfc               ! local scalar 
    7881      REAL(wp), DIMENSION(jpi,jpj)     ::  z2d   ! 2D workspace 
    7982      REAL(wp), DIMENSION(jpi,jpj,jpk) ::  z3d   ! 3D workspace 
    8083      REAL(wp), DIMENSION(jpi,jpj,jpk) ::  zmask   ! 3D workspace 
    8184      REAL(wp), DIMENSION(jpi,jpj,jpk,jpts) ::  zts   ! 3D workspace 
    82       CHARACTER( len = 10 )  :: cl1 
     85      REAL(wp), DIMENSION(jpj)     ::  vsum   ! 1D workspace 
     86      REAL(wp), DIMENSION(jpj,jpts)     ::  tssum   ! 1D workspace 
     87  
     88      ! 
     89      !overturning calculation 
     90      REAL(wp), DIMENSION(jpj,jpk,nptr) ::   sjk  , r1_sjk ! i-mean i-k-surface and its inverse 
     91      REAL(wp), DIMENSION(jpj,jpk,nptr) ::   v_msf, sn_jk  , tn_jk ! i-mean T and S, j-Stream-Function 
     92      REAL(wp), DIMENSION(jpi,jpj,jpk) ::  zvn   ! 3D workspace 
     93 
     94 
     95      CHARACTER( len = 12 )  :: cl1 
    8396      !!---------------------------------------------------------------------- 
    8497      ! 
    … …  
    109122            END DO 
    110123         ENDIF 
     124         IF( iom_use("sopstove") .OR. iom_use("sophtove") .OR. iom_use("sopstbtr") .OR. iom_use("sophtbtr") ) THEN 
     125            ! define fields multiplied by scalar 
     126            zmask(:,:,:) = 0._wp 
     127            zts(:,:,:,:) = 0._wp 
     128            zvn(:,:,:) = 0._wp 
     129            DO jk = 1, jpkm1 
     130               DO jj = 1, jpjm1 
     131                  DO ji = 1, jpi 
     132                     zvfc = e1v(ji,jj) * e3v_n(ji,jj,jk) 
     133                     zmask(ji,jj,jk)      = vmask(ji,jj,jk)      * zvfc 
     134                     zts(ji,jj,jk,jp_tem) = (tsn(ji,jj,jk,jp_tem)+tsn(ji,jj+1,jk,jp_tem)) * 0.5 * zvfc  !Tracers averaged onto V grid 
     135                     zts(ji,jj,jk,jp_sal) = (tsn(ji,jj,jk,jp_sal)+tsn(ji,jj+1,jk,jp_sal)) * 0.5 * zvfc 
     136                     zvn(ji,jj,jk)        = vn(ji,jj,jk)         * zvfc 
     137                  ENDDO 
     138               ENDDO 
     139             ENDDO 
     140         ENDIF 
     141         IF( iom_use("sopstove") .OR. iom_use("sophtove") ) THEN 
     142             sjk(:,:,1) = ptr_sjk( zmask(:,:,:), btmsk(:,:,1) ) 
     143             r1_sjk(:,:,1) = 0._wp 
     144             WHERE( sjk(:,:,1) /= 0._wp )   r1_sjk(:,:,1) = 1._wp / sjk(:,:,1) 
     145 
     146             ! i-mean T and S, j-Stream-Function, global 
     147             tn_jk(:,:,1) = ptr_sjk( zts(:,:,:,jp_tem) ) * r1_sjk(:,:,1) 
     148             sn_jk(:,:,1) = ptr_sjk( zts(:,:,:,jp_sal) ) * r1_sjk(:,:,1) 
     149             v_msf(:,:,1) = ptr_sjk( zvn(:,:,:) ) 
     150 
     151             htr_ove(:,1) = SUM( v_msf(:,:,1)*tn_jk(:,:,1) ,2 ) 
     152             str_ove(:,1) = SUM( v_msf(:,:,1)*sn_jk(:,:,1) ,2 ) 
     153 
     154             z2d(1,:) = htr_ove(:,1) * rc_pwatt        !  (conversion in PW) 
     155             DO ji = 1, jpi 
     156               z2d(ji,:) = z2d(1,:) 
     157             ENDDO 
     158             cl1 = 'sophtove' 
     159             CALL iom_put( TRIM(cl1), z2d ) 
     160             z2d(1,:) = str_ove(:,1) * rc_ggram        !  (conversion in Gg) 
     161             DO ji = 1, jpi 
     162               z2d(ji,:) = z2d(1,:) 
     163             ENDDO 
     164             cl1 = 'sopstove' 
     165             CALL iom_put( TRIM(cl1), z2d ) 
     166             IF( ln_subbas ) THEN 
     167                DO jn = 2, nptr 
     168                    sjk(:,:,jn) = ptr_sjk( zmask(:,:,:), btmsk(:,:,jn) ) 
     169                    r1_sjk(:,:,jn) = 0._wp 
     170                    WHERE( sjk(:,:,jn) /= 0._wp )   r1_sjk(:,:,jn) = 1._wp / sjk(:,:,jn) 
     171 
     172                    ! i-mean T and S, j-Stream-Function, basin 
     173                    tn_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) * r1_sjk(:,:,jn) 
     174                    sn_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) * r1_sjk(:,:,jn) 
     175                    v_msf(:,:,jn) = ptr_sjk( zvn(:,:,:), btmsk(:,:,jn) )  
     176                    htr_ove(:,jn) = SUM( v_msf(:,:,jn)*tn_jk(:,:,jn) ,2 ) 
     177                    str_ove(:,jn) = SUM( v_msf(:,:,jn)*sn_jk(:,:,jn) ,2 ) 
     178 
     179                    z2d(1,:) = htr_ove(:,jn) * rc_pwatt !  (conversion in PW) 
     180                    DO ji = 1, jpi 
     181                        z2d(ji,:) = z2d(1,:) 
     182                    ENDDO 
     183                    cl1 = TRIM('sophtove_'//clsubb(jn)) 
     184                    CALL iom_put( cl1, z2d ) 
     185                    z2d(1,:) = str_ove(:,jn) * rc_ggram        ! (conversion in Gg) 
     186                    DO ji = 1, jpi 
     187                        z2d(ji,:) = z2d(1,:) 
     188                    ENDDO 
     189                    cl1 = TRIM('sopstove_'//clsubb(jn)) 
     190                    CALL iom_put( cl1, z2d ) 
     191                END DO 
     192             ENDIF 
     193         ENDIF 
     194         IF( iom_use("sopstbtr") .OR. iom_use("sophtbtr") ) THEN 
     195         ! Calculate barotropic heat and salt transport here  
     196             sjk(:,1,1) = ptr_sj( zmask(:,:,:), btmsk(:,:,1) ) 
     197             r1_sjk(:,1,1) = 0._wp 
     198             WHERE( sjk(:,1,1) /= 0._wp )   r1_sjk(:,1,1) = 1._wp / sjk(:,1,1) 
     199             
     200            vsum = ptr_sj( zvn(:,:,:), btmsk(:,:,1)) 
     201            tssum(:,jp_tem) = ptr_sj( zts(:,:,:,jp_tem), btmsk(:,:,1) ) 
     202            tssum(:,jp_sal) = ptr_sj( zts(:,:,:,jp_sal), btmsk(:,:,1) ) 
     203            htr_btr(:,1) = vsum * tssum(:,jp_tem) * r1_sjk(:,1,1) 
     204            str_btr(:,1) = vsum * tssum(:,jp_sal) * r1_sjk(:,1,1) 
     205            z2d(1,:) = htr_btr(:,1) * rc_pwatt        !  (conversion in PW) 
     206            DO ji = 2, jpi 
     207               z2d(ji,:) = z2d(1,:) 
     208            ENDDO 
     209            cl1 = 'sophtbtr' 
     210            CALL iom_put( TRIM(cl1), z2d ) 
     211            z2d(1,:) = str_btr(:,1) * rc_ggram        !  (conversion in Gg) 
     212            DO ji = 2, jpi 
     213              z2d(ji,:) = z2d(1,:) 
     214            ENDDO 
     215            cl1 = 'sopstbtr' 
     216            CALL iom_put( TRIM(cl1), z2d ) 
     217            IF( ln_subbas ) THEN 
     218                DO jn = 2, nptr 
     219                    sjk(:,1,jn) = ptr_sj( zmask(:,:,:), btmsk(:,:,jn) ) 
     220                    r1_sjk(:,1,jn) = 0._wp 
     221                    WHERE( sjk(:,1,jn) /= 0._wp )   r1_sjk(:,1,jn) = 1._wp / sjk(:,1,jn) 
     222                    vsum = ptr_sj( zvn(:,:,:), btmsk(:,:,jn)) 
     223                    tssum(:,jp_tem) = ptr_sj( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) 
     224                    tssum(:,jp_sal) = ptr_sj( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) 
     225                    htr_btr(:,jn) = vsum * tssum(:,jp_tem) * r1_sjk(:,1,jn) 
     226                    str_btr(:,jn) = vsum * tssum(:,jp_sal) * r1_sjk(:,1,jn) 
     227                    z2d(1,:) = htr_btr(:,jn) * rc_pwatt !  (conversion in PW) 
     228                    DO ji = 1, jpi 
     229                        z2d(ji,:) = z2d(1,:) 
     230                    ENDDO 
     231                    cl1 = TRIM('sophtbtr_'//clsubb(jn)) 
     232                    CALL iom_put( cl1, z2d ) 
     233                    z2d(1,:) = str_btr(:,jn) * rc_ggram        ! (conversion in Gg) 
     234                    DO ji = 1, jpi 
     235                        z2d(ji,:) = z2d(1,:) 
     236                    ENDDO 
     237                    cl1 = TRIM('sopstbtr_'//clsubb(jn)) 
     238                    CALL iom_put( cl1, z2d ) 
     239               ENDDO 
     240            ENDIF !ln_subbas 
     241         ENDIF !iom_use("sopstbtr....) 
    111242         ! 
    112243      ELSE 
    … …  
    148279         !                                ! Advective and diffusive heat and salt transport 
    149280         IF( iom_use("sophtadv") .OR. iom_use("sopstadv") ) THEN    
    150             z2d(1,:) = htr_adv(:) * rc_pwatt        !  (conversion in PW) 
     281            z2d(1,:) = htr_adv(:,1) * rc_pwatt        !  (conversion in PW) 
    151282            DO ji = 1, jpi 
    152283               z2d(ji,:) = z2d(1,:) 
    … …  
    154285            cl1 = 'sophtadv'                  
    155286            CALL iom_put( TRIM(cl1), z2d ) 
    156             z2d(1,:) = str_adv(:) * rc_ggram        ! (conversion in Gg) 
     287            z2d(1,:) = str_adv(:,1) * rc_ggram        ! (conversion in Gg) 
    157288            DO ji = 1, jpi 
    158289               z2d(ji,:) = z2d(1,:) 
    … …  
    160291            cl1 = 'sopstadv' 
    161292            CALL iom_put( TRIM(cl1), z2d ) 
     293            IF( ln_subbas ) THEN 
     294              DO jn=2,nptr 
     295               z2d(1,:) = htr_adv(:,jn) * rc_pwatt        !  (conversion in PW) 
     296               DO ji = 1, jpi 
     297                 z2d(ji,:) = z2d(1,:) 
     298               ENDDO 
     299               cl1 = TRIM('sophtadv_'//clsubb(jn))                  
     300               CALL iom_put( cl1, z2d ) 
     301               z2d(1,:) = str_adv(:,jn) * rc_ggram        ! (conversion in Gg) 
     302               DO ji = 1, jpi 
     303                  z2d(ji,:) = z2d(1,:) 
     304               ENDDO 
     305               cl1 = TRIM('sopstadv_'//clsubb(jn))                  
     306               CALL iom_put( cl1, z2d )               
     307              ENDDO 
     308            ENDIF 
    162309         ENDIF 
    163310         ! 
    164311         IF( iom_use("sophtldf") .OR. iom_use("sopstldf") ) THEN    
    165             z2d(1,:) = htr_ldf(:) * rc_pwatt        !  (conversion in PW)  
     312            z2d(1,:) = htr_ldf(:,1) * rc_pwatt        !  (conversion in PW)  
    166313            DO ji = 1, jpi 
    167314               z2d(ji,:) = z2d(1,:) 
    … …  
    169316            cl1 = 'sophtldf' 
    170317            CALL iom_put( TRIM(cl1), z2d ) 
    171             z2d(1,:) = str_ldf(:) * rc_ggram        !  (conversion in Gg) 
     318            z2d(1,:) = str_ldf(:,1) * rc_ggram        !  (conversion in Gg) 
    172319            DO ji = 1, jpi 
    173320               z2d(ji,:) = z2d(1,:) 
    … …  
    175322            cl1 = 'sopstldf' 
    176323            CALL iom_put( TRIM(cl1), z2d ) 
     324            IF( ln_subbas ) THEN 
     325              DO jn=2,nptr 
     326               z2d(1,:) = htr_ldf(:,jn) * rc_pwatt        !  (conversion in PW) 
     327               DO ji = 1, jpi 
     328                 z2d(ji,:) = z2d(1,:) 
     329               ENDDO 
     330               cl1 = TRIM('sophtldf_'//clsubb(jn))                  
     331               CALL iom_put( cl1, z2d ) 
     332               z2d(1,:) = str_ldf(:,jn) * rc_ggram        ! (conversion in Gg) 
     333               DO ji = 1, jpi 
     334                  z2d(ji,:) = z2d(1,:) 
     335               ENDDO 
     336               cl1 = TRIM('sopstldf_'//clsubb(jn))                  
     337               CALL iom_put( cl1, z2d )               
     338              ENDDO 
     339            ENDIF 
     340         ENDIF 
     341 
     342         IF( iom_use("sophteiv") .OR. iom_use("sopsteiv") ) THEN  
     343            z2d(1,:) = htr_eiv(:,1) * rc_pwatt        !  (conversion in PW)  
     344            DO ji = 1, jpi 
     345               z2d(ji,:) = z2d(1,:) 
     346            ENDDO 
     347            cl1 = 'sophteiv' 
     348            CALL iom_put( TRIM(cl1), z2d ) 
     349            z2d(1,:) = str_eiv(:,1) * rc_ggram        !  (conversion in Gg) 
     350            DO ji = 1, jpi 
     351               z2d(ji,:) = z2d(1,:) 
     352            ENDDO 
     353            cl1 = 'sopsteiv' 
     354            CALL iom_put( TRIM(cl1), z2d ) 
     355            IF( ln_subbas ) THEN 
     356               DO jn=2,nptr 
     357                  z2d(1,:) = htr_eiv(:,jn) * rc_pwatt        !  (conversion in PW) 
     358                  DO ji = 1, jpi 
     359                     z2d(ji,:) = z2d(1,:) 
     360                  ENDDO 
     361                  cl1 = TRIM('sophteiv_'//clsubb(jn))                  
     362                  CALL iom_put( cl1, z2d ) 
     363                  z2d(1,:) = str_eiv(:,jn) * rc_ggram        ! (conversion in Gg) 
     364                  DO ji = 1, jpi 
     365                     z2d(ji,:) = z2d(1,:) 
     366                  ENDDO 
     367                  cl1 = TRIM('sopsteiv_'//clsubb(jn))  
     368                  CALL iom_put( cl1, z2d )               
     369               ENDDO 
     370            ENDIF 
    177371         ENDIF 
    178372         ! 
    … …  
    254448         ! Initialise arrays to zero because diatpr is called before they are first calculated 
    255449         ! Note that this means diagnostics will not be exactly correct when model run is restarted. 
    256          htr_adv(:) = 0._wp  ;  str_adv(:) =  0._wp   
    257          htr_ldf(:) = 0._wp  ;  str_ldf(:) =  0._wp  
     450         htr_adv(:,:) = 0._wp  ;  str_adv(:,:) =  0._wp  
     451         htr_ldf(:,:) = 0._wp  ;  str_ldf(:,:) =  0._wp  
     452         htr_eiv(:,:) = 0._wp  ;  str_eiv(:,:) =  0._wp  
     453         htr_ove(:,:) = 0._wp  ;   str_ove(:,:) =  0._wp 
     454         htr_btr(:,:) = 0._wp  ;   str_btr(:,:) =  0._wp 
    258455         ! 
    259456      ENDIF  
    … …  
    261458   END SUBROUTINE dia_ptr_init 
    262459 
     460   SUBROUTINE dia_ptr_hst( ktra, cptr, pva )  
     461      !!---------------------------------------------------------------------- 
     462      !!                    ***  ROUTINE dia_ptr_hst *** 
     463      !!---------------------------------------------------------------------- 
     464      !! Wrapper for heat and salt transport calculations to calculate them for each basin 
     465      !! Called from all advection and/or diffusion routines 
     466      !!---------------------------------------------------------------------- 
     467      INTEGER                         , INTENT(in )  :: ktra  ! tracer index 
     468      CHARACTER(len=3)                , INTENT(in)   :: cptr  ! transport type  'adv'/'ldf'/'eiv' 
     469      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in)   :: pva   ! 3D input array of advection/diffusion 
     470      INTEGER                                        :: jn    ! 
     471 
     472      IF( cptr == 'adv' ) THEN 
     473         IF( ktra == jp_tem )  htr_adv(:,1) = ptr_sj( pva(:,:,:) ) 
     474         IF( ktra == jp_sal )  str_adv(:,1) = ptr_sj( pva(:,:,:) ) 
     475      ENDIF 
     476      IF( cptr == 'ldf' ) THEN 
     477         IF( ktra == jp_tem )  htr_ldf(:,1) = ptr_sj( pva(:,:,:) ) 
     478         IF( ktra == jp_sal )  str_ldf(:,1) = ptr_sj( pva(:,:,:) ) 
     479      ENDIF 
     480      IF( cptr == 'eiv' ) THEN 
     481         IF( ktra == jp_tem )  htr_eiv(:,1) = ptr_sj( pva(:,:,:) ) 
     482         IF( ktra == jp_sal )  str_eiv(:,1) = ptr_sj( pva(:,:,:) ) 
     483      ENDIF 
     484      ! 
     485      IF( ln_subbas ) THEN 
     486         ! 
     487         IF( cptr == 'adv' ) THEN 
     488             IF( ktra == jp_tem ) THEN  
     489                DO jn = 2, nptr 
     490                   htr_adv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
     491                END DO 
     492             ENDIF 
     493             IF( ktra == jp_sal ) THEN  
     494                DO jn = 2, nptr 
     495                   str_adv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
     496                END DO 
     497             ENDIF 
     498         ENDIF 
     499         IF( cptr == 'ldf' ) THEN 
     500             IF( ktra == jp_tem ) THEN  
     501                DO jn = 2, nptr 
     502                    htr_ldf(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
     503                 END DO 
     504             ENDIF 
     505             IF( ktra == jp_sal ) THEN  
     506                DO jn = 2, nptr 
     507                   str_ldf(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
     508                END DO 
     509             ENDIF 
     510         ENDIF 
     511         IF( cptr == 'eiv' ) THEN 
     512             IF( ktra == jp_tem ) THEN  
     513                DO jn = 2, nptr 
     514                    htr_eiv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
     515                 END DO 
     516             ENDIF 
     517             IF( ktra == jp_sal ) THEN  
     518                DO jn = 2, nptr 
     519                   str_eiv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 
     520                END DO 
     521             ENDIF 
     522         ENDIF 
     523         ! 
     524      ENDIF 
     525   END SUBROUTINE dia_ptr_hst 
     526 
    263527 
    264528   FUNCTION dia_ptr_alloc() 
    … …  
    271535      ierr(:) = 0 
    272536      ! 
    273       ALLOCATE( btmsk(jpi,jpj,nptr) ,           & 
    274          &      htr_adv(jpj) , str_adv(jpj) ,   & 
    275          &      htr_ldf(jpj) , str_ldf(jpj) , STAT=ierr(1)  ) 
     537      ALLOCATE( btmsk(jpi,jpj,nptr) ,              & 
     538         &      htr_adv(jpj,nptr) , str_adv(jpj,nptr) ,   & 
     539         &      htr_eiv(jpj,nptr) , str_eiv(jpj,nptr) ,   & 
     540         &      htr_ove(jpj,nptr) , str_ove(jpj,nptr) ,   & 
     541         &      htr_btr(jpj,nptr) , str_btr(jpj,nptr) ,   & 
     542         &      htr_ldf(jpj,nptr) , str_ldf(jpj,nptr) , STAT=ierr(1)  ) 
    276543         ! 
    277544      ALLOCATE( p_fval1d(jpj), p_fval2d(jpj,jpk), Stat=ierr(2)) 

  • trunk/NEMOGCM/NEMO/OPA_SRC/DIA/diatmb.F90

    r6140 r7646  
    44   !! Harmonic analysis of tidal constituents  
    55   !!====================================================================== 
    6    !! History :  3.6  !  2014  (E O'Dea)  Original code 
     6   !! History :  3.6  !  08-2014  (E O'Dea)  Original code 
     7   !!            3.7  !  05-2016  (G. Madec)  use mbkt, mikt to account for ocean cavities 
    78   !!---------------------------------------------------------------------- 
    89   USE oce             ! ocean dynamics and tracers variables 
    910   USE dom_oce         ! ocean space and time domain 
     11   ! 
    1012   USE in_out_manager  ! I/O units 
    1113   USE iom             ! I/0 library 
    … …  
    3133      !!                  ***  ROUTINE dia_tmb_init  *** 
    3234      !!      
    33       !! ** Purpose: Initialization of tmb namelist  
     35      !! ** Purpose :  Initialization of tmb namelist  
    3436      !!         
    35       !! ** Method : Read namelist 
    36       !!   History 
    37       !!   3.6  !  08-14  (E. O'Dea) Routine to initialize dia_tmb 
     37      !! ** Method  :   Read namelist 
    3838      !!--------------------------------------------------------------------------- 
    39       !! 
    4039      INTEGER ::   ios                 ! Local integer output status for namelist read 
    4140      ! 
    … …  
    5958         WRITE(numout,*) 'Switch for TMB diagnostics (T) or not (F)  ln_diatmb  = ', ln_diatmb 
    6059      ENDIF 
    61  
     60      ! 
    6261   END SUBROUTINE dia_tmb_init 
    6362 
    64    SUBROUTINE dia_calctmb( pinfield,pouttmb ) 
     63 
     64   SUBROUTINE dia_calctmb( pfield, ptmb ) 
    6565      !!--------------------------------------------------------------------- 
    6666      !!                  ***  ROUTINE dia_tmb  *** 
    … …  
    6868      !! ** Purpose :    Find the Top, Mid and Bottom fields of water Column 
    6969      !! 
    70       !! ** Method  :    
    71       !!      use mbathy to find surface, mid and bottom of model levels 
     70      !! ** Method  :    use mbkt, mikt to find surface, mid and bottom of  
     71      !!              model levels due to potential existence of ocean cavities 
    7272      !! 
    73       !! History : 
    74       !!   3.6  !  08-14  (E. O'Dea) Routine based on dia_wri_foam 
    7573      !!---------------------------------------------------------------------- 
    76       !! * Modules used 
    77  
    78       ! Routine to map 3d field to top, middle, bottom 
    79       IMPLICIT NONE 
    80  
    81  
    82       ! Routine arguments 
    83       REAL(wp), DIMENSION(jpi, jpj, jpk), INTENT(IN   ) :: pinfield    ! Input 3d field and mask 
    84       REAL(wp), DIMENSION(jpi, jpj, 3  ), INTENT(  OUT) :: pouttmb     ! Output top, middle, bottom 
    85  
    86  
    87  
    88       ! Local variables 
    89       INTEGER :: ji,jj,jk  ! Dummy loop indices 
    90  
    91       ! Local Real 
    92       REAL(wp)                         ::   zmdi  !  set masked values 
    93  
    94       zmdi=1.e+20 !missing data indicator for masking 
    95  
    96       ! Calculate top 
    97       pouttmb(:,:,1) = pinfield(:,:,1)*tmask(:,:,1)  + zmdi*(1.0-tmask(:,:,1)) 
    98  
    99       ! Calculate middle 
    100       DO jj = 1,jpj 
    101          DO ji = 1,jpi 
    102             jk              = max(1,mbathy(ji,jj)/2) 
    103             pouttmb(ji,jj,2) = pinfield(ji,jj,jk)*tmask(ji,jj,jk)  + zmdi*(1.0-tmask(ji,jj,jk)) 
     74      REAL(wp), DIMENSION(jpi, jpj, jpk), INTENT(in   ) :: pfield   ! Input 3d field and mask 
     75      REAL(wp), DIMENSION(jpi, jpj,  3 ), INTENT(  out) :: ptmb     ! top, middle, bottom extracted from pfield 
     76      ! 
     77      INTEGER  ::   ji, jj  ! Dummy loop indices 
     78      INTEGER  ::   itop, imid, ibot  ! local integers 
     79      REAL(wp) ::   zmdi = 1.e+20_wp  ! land value 
     80      !!--------------------------------------------------------------------- 
     81      ! 
     82      DO jj = 1, jpj 
     83         DO ji = 1, jpi 
     84            itop = mikt(ji,jj)                        ! top    ocean  
     85            ibot = mbkt(ji,jj)                        ! bottom ocean  
     86            imid =  itop + ( ibot - itop + 1 ) / 2    ! middle ocean           
     87            !                     
     88            ptmb(ji,jj,1) = pfield(ji,jj,itop)*tmask(ji,jj,itop)  + zmdi*( 1._wp-tmask(ji,jj,itop) ) 
     89            ptmb(ji,jj,2) = pfield(ji,jj,imid)*tmask(ji,jj,imid)  + zmdi*( 1._wp-tmask(ji,jj,imid) ) 
     90            ptmb(ji,jj,3) = pfield(ji,jj,ibot)*tmask(ji,jj,ibot)  + zmdi*( 1._wp-tmask(ji,jj,ibot) ) 
    10491         END DO 
    10592      END DO 
    106  
    107       ! Calculate bottom 
    108       DO jj = 1,jpj 
    109          DO ji = 1,jpi 
    110             jk              = max(1,mbathy(ji,jj) - 1) 
    111             pouttmb(ji,jj,3) = pinfield(ji,jj,jk)*tmask(ji,jj,jk)  + zmdi*(1.0-tmask(ji,jj,jk)) 
    112          END DO 
    113       END DO 
    114  
     93      ! 
    11594   END SUBROUTINE dia_calctmb 
    116  
    11795 
    11896 
    … …  
    122100      !! ** Purpose :   Write diagnostics for Top, Mid and Bottom of water Column 
    123101      !! 
    124       !! ** Method  :    
    125       !!      use mbathy to find surface, mid and bottom of model levels 
     102      !! ** Method  :  use mikt,mbkt to find surface, mid and bottom of model levels 
    126103      !!      calls calctmb to retrieve TMB values before sending to iom_put 
    127104      !! 
    128       !! History : 
    129       !!   3.6  !  08-14  (E. O'Dea)  
    130       !!          
    131105      !!-------------------------------------------------------------------- 
    132       REAL(wp), POINTER, DIMENSION(:,:,:) :: zwtmb    ! temporary workspace  
    133       REAL(wp)                         ::   zmdi      ! set masked values 
    134  
    135       zmdi=1.e+20 !missing data indicator for maskin 
    136  
     106      REAL(wp) ::   zmdi =1.e+20     ! land value 
     107      REAL(wp), POINTER, DIMENSION(:,:,:) :: zwtmb    ! workspace  
     108      !!-------------------------------------------------------------------- 
     109      ! 
    137110      IF (ln_diatmb) THEN 
    138          CALL wrk_alloc( jpi , jpj, 3 , zwtmb ) 
     111         CALL wrk_alloc( jpi,jpj,3  , zwtmb ) 
    139112         CALL dia_calctmb(  tsn(:,:,:,jp_tem),zwtmb ) 
    140113         !ssh already output but here we output it masked 
    141          CALL iom_put( "sshnmasked" , sshn(:,:)*tmask(:,:,1) + zmdi*(1.0 - tmask(:,:,1)) )   ! tmb Temperature 
     114         CALL iom_put( "sshnmasked" , sshn(:,:)*tmask(:,:,1) + zmdi*(1.0 - tmask(:,:,1)) ) 
    142115         CALL iom_put( "top_temp" , zwtmb(:,:,1) )    ! tmb Temperature 
    143116         CALL iom_put( "mid_temp" , zwtmb(:,:,2) )    ! tmb Temperature 
    … …  
    161134         CALL iom_put( "bot_v" , zwtmb(:,:,3) )    ! tmb  V Velocity 
    162135!Called in  dynspg_ts.F90       CALL iom_put( "baro_v" , vn_b )    ! Barotropic  V Velocity 
     136         CALL wrk_dealloc( jpi,jpj,3   , zwtmb ) 
    163137      ELSE 
    164138         CALL ctl_warn('dia_tmb: tmb diagnostic is set to false you should not have seen this') 
    165139      ENDIF 
    166  
     140      ! 
    167141   END SUBROUTINE dia_tmb 
    168142   !!====================================================================== 

  • trunk/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90

    r6387 r7646  
    4141   USE zdfddm          ! vertical  physics: double diffusion 
    4242   USE diahth          ! thermocline diagnostics 
     43   USE wet_dry         ! wetting and drying 
     44   USE sbcwave         ! wave parameters 
    4345   ! 
    4446   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    … …  
    153155 
    154156      CALL iom_put( "ssh" , sshn )                 ! sea surface height 
     157      IF( iom_use("wetdep") )   &                  ! wet depth 
     158         CALL iom_put( "wetdep" , ht_wd(:,:) + sshn(:,:) ) 
    155159       
    156160      CALL iom_put( "toce", tsn(:,:,:,jp_tem) )    ! 3D temperature 
    … …  
    302306      CALL iom_put( "hdiv", hdivn )                  ! Horizontal divergence 
    303307      ! 
    304       IF( iom_use("u_masstr") .OR. iom_use("u_heattr") .OR. iom_use("u_salttr") ) THEN 
     308      IF( iom_use("u_masstr") .OR. iom_use("u_masstr_vint") .OR. iom_use("u_heattr") .OR. iom_use("u_salttr") ) THEN 
    305309         z3d(:,:,jpk) = 0.e0 
     310         z2d(:,:) = 0.e0 
    306311         DO jk = 1, jpkm1 
    307312            z3d(:,:,jk) = rau0 * un(:,:,jk) * e2u(:,:) * e3u_n(:,:,jk) * umask(:,:,jk) 
     313            z2d(:,:) = z2d(:,:) + z3d(:,:,jk) 
    308314         END DO 
    309315         CALL iom_put( "u_masstr", z3d )                  ! mass transport in i-direction 
     316         CALL iom_put( "u_masstr_vint", z2d )             ! mass transport in i-direction vertical sum 
    310317      ENDIF 
    311318       
    … …  
    370377         CALL iom_put( "v_salttr", 0.5 * z2d )            !  heat transport in j-direction 
    371378      ENDIF 
     379 
     380      ! Vertical integral of temperature 
     381      IF( iom_use("tosmint") ) THEN 
     382         z2d(:,:)=0._wp 
     383         DO jk = 1, jpkm1 
     384            DO jj = 2, jpjm1 
     385               DO ji = fs_2, fs_jpim1   ! vector opt. 
     386                  z2d(ji,jj) = z2d(ji,jj) + rau0 * e3t_n(ji,jj,jk) *  tsn(ji,jj,jk,jp_tem) 
     387               END DO 
     388            END DO 
     389         END DO 
     390         CALL lbc_lnk( z2d, 'T', -1. ) 
     391         CALL iom_put( "tosmint", z2d )  
     392      ENDIF 
     393 
     394      ! Vertical integral of salinity 
     395      IF( iom_use("somint") ) THEN 
     396         z2d(:,:)=0._wp 
     397         DO jk = 1, jpkm1 
     398            DO jj = 2, jpjm1 
     399               DO ji = fs_2, fs_jpim1   ! vector opt. 
     400                  z2d(ji,jj) = z2d(ji,jj) + rau0 * e3t_n(ji,jj,jk) * tsn(ji,jj,jk,jp_sal) 
     401               END DO 
     402            END DO 
     403         END DO 
     404         CALL lbc_lnk( z2d, 'T', -1. ) 
     405         CALL iom_put( "somint", z2d )  
     406      ENDIF 
     407 
     408      CALL iom_put( "bn2", rn2 )  !Brunt-Vaisala buoyancy frequency (N^2) 
    372409      ! 
    373410      CALL wrk_dealloc( jpi , jpj      , z2d ) 
    … …  
    666703         CALL histdef( nid_T, "so28chgt", "Depth of 28C isotherm"              , "m"      ,   & ! hd28 
    667704            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
    668          CALL histdef( nid_T, "sohtc300", "Heat content 300 m"                 , "W"      ,   & ! htc3 
     705         CALL histdef( nid_T, "sohtc300", "Heat content 300 m"                 , "J/m2"   ,   & ! htc3 
    669706            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
    670707#endif 
    … …  
    682719         CALL histdef( nid_U, "vozocrtx", "Zonal Current"                      , "m/s"    ,   &  ! un 
    683720            &          jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout ) 
     721         IF( ln_wave .AND. ln_sdw) THEN 
     722            CALL histdef( nid_U, "sdzocrtx", "Stokes Drift Zonal Current"         , "m/s"    ,   &  ! usd 
     723               &          jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout ) 
     724         ENDIF 
    684725         !                                                                                      !!! nid_U : 2D 
    685726         CALL histdef( nid_U, "sozotaux", "Wind Stress along i-axis"           , "N/m2"   ,   &  ! utau 
    … …  
    691732         CALL histdef( nid_V, "vomecrty", "Meridional Current"                 , "m/s"    ,   &  ! vn 
    692733            &          jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout ) 
     734         IF( ln_wave .AND. ln_sdw) THEN 
     735            CALL histdef( nid_V, "sdmecrty", "Stokes Drift Meridional Current"    , "m/s"    ,   &  ! vsd 
     736               &          jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout ) 
     737         ENDIF 
    693738         !                                                                                      !!! nid_V : 2D 
    694739         CALL histdef( nid_V, "sometauy", "Wind Stress along j-axis"           , "N/m2"   ,   &  ! vtau 
    … …  
    707752         IF( lk_zdfddm ) THEN 
    708753            CALL histdef( nid_W,"voddmavs","Salt Vertical Eddy Diffusivity"    , "m2/s"   ,   &  ! avs 
     754               &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout ) 
     755         ENDIF 
     756          
     757         IF( ln_wave .AND. ln_sdw) THEN 
     758            CALL histdef( nid_W, "sdvecrtz", "Stokes Drift Vertical Current"   , "m/s"    ,   &  ! wsd 
    709759               &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout ) 
    710760         ENDIF 
    … …  
    829879      ENDIF 
    830880 
     881      IF( ln_wave .AND. ln_sdw ) THEN 
     882         CALL histwrite( nid_U, "sdzocrtx", it, usd           , ndim_U , ndex_U )    ! i-StokesDrift-current 
     883         CALL histwrite( nid_V, "sdmecrty", it, vsd           , ndim_V , ndex_V )    ! j-StokesDrift-current 
     884         CALL histwrite( nid_W, "sdvecrtz", it, wsd           , ndim_T , ndex_T )    ! StokesDrift vert. current 
     885      ENDIF 
     886 
    831887      ! 3. Close all files 
    832888      ! --------------------------------------- 
    … …  
    912968         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )  
    913969         ! 
    914       CALL histdef( id_i, "ahtu"    , "u-eddy diffusivity"    , "m2/s"    ,   &   ! zonal current 
    915          &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
    916       CALL histdef( id_i, "ahtv"    , "v-eddy diffusivity"    , "m2/s"    ,   &   ! meridonal current 
    917          &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )  
    918       CALL histdef( id_i, "ahmt"    , "t-eddy viscosity"      , "m2/s"    ,   &   ! zonal current 
    919          &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
    920       CALL histdef( id_i, "ahmf"    , "f-eddy viscosity"      , "m2/s"    ,   &   ! meridonal current 
    921          &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )  
     970      IF( ALLOCATED(ahtu) ) THEN 
     971         CALL histdef( id_i, "ahtu"    , "u-eddy diffusivity"    , "m2/s"    ,   &   ! zonal current 
     972            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
     973         CALL histdef( id_i, "ahtv"    , "v-eddy diffusivity"    , "m2/s"    ,   &   ! meridonal current 
     974            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )  
     975      ENDIF 
     976      IF( ALLOCATED(ahmt) ) THEN  
     977         CALL histdef( id_i, "ahmt"    , "t-eddy viscosity"      , "m2/s"    ,   &   ! zonal current 
     978            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
     979         CALL histdef( id_i, "ahmf"    , "f-eddy viscosity"      , "m2/s"    ,   &   ! meridonal current 
     980            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )  
     981      ENDIF 
    922982         ! 
    923983      CALL histdef( id_i, "sowaflup", "Net Upward Water Flux" , "Kg/m2/S",   &   ! net freshwater  
    … …  
    939999            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
    9401000      ENDIF 
     1001      ! 
     1002      IF( ln_wave .AND. ln_sdw ) THEN 
     1003         CALL histdef( id_i, "sdzocrtx", "Stokes Drift Zonal"    , "m/s"    , &   ! StokesDrift zonal current 
     1004            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
     1005         CALL histdef( id_i, "sdmecrty", "Stokes Drift Merid"    , "m/s"    , &   ! StokesDrift meridonal current 
     1006            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
     1007         CALL histdef( id_i, "sdvecrtz", "Stokes Drift Vert"     , "m/s"    , &   ! StokesDrift vertical current 
     1008            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
     1009      ENDIF 
    9411010 
    9421011#if defined key_lim2 
    … …  
    9631032      CALL histwrite( id_i, "vovecrtz", kt, wn               , jpi*jpj*jpk, idex )    ! now k-velocity 
    9641033      ! 
    965       CALL histwrite( id_i, "ahtu"    , kt, ahtu             , jpi*jpj*jpk, idex )    ! aht at u-point 
    966       CALL histwrite( id_i, "ahtv"    , kt, ahtv             , jpi*jpj*jpk, idex )    !  -  at v-point 
    967       CALL histwrite( id_i, "ahmt"    , kt, ahmt             , jpi*jpj*jpk, idex )    ! ahm at t-point 
    968       CALL histwrite( id_i, "ahmf"    , kt, ahmf             , jpi*jpj*jpk, idex )    !  -  at f-point 
    969       ! 
    970       CALL histwrite( id_i, "sowaflup", kt, emp-rnf          , jpi*jpj    , idex )    ! freshwater budget 
     1034      IF( ALLOCATED(ahtu) ) THEN 
     1035         CALL histwrite( id_i, "ahtu"    , kt, ahtu             , jpi*jpj*jpk, idex )    ! aht at u-point 
     1036         CALL histwrite( id_i, "ahtv"    , kt, ahtv             , jpi*jpj*jpk, idex )    !  -  at v-point 
     1037      ENDIF 
     1038      IF( ALLOCATED(ahmt) ) THEN 
     1039         CALL histwrite( id_i, "ahmt"    , kt, ahmt             , jpi*jpj*jpk, idex )    ! ahm at t-point 
     1040         CALL histwrite( id_i, "ahmf"    , kt, ahmf             , jpi*jpj*jpk, idex )    !  -  at f-point 
     1041      ENDIF 
     1042      ! 
     1043      CALL histwrite( id_i, "sowaflup", kt, emp - rnf        , jpi*jpj    , idex )    ! freshwater budget 
    9711044      CALL histwrite( id_i, "sohefldo", kt, qsr + qns        , jpi*jpj    , idex )    ! total heat flux 
    9721045      CALL histwrite( id_i, "soshfldo", kt, qsr              , jpi*jpj    , idex )    ! solar heat flux 
    … …  
    9791052         CALL histwrite( id_i, "vovvle3t", kt, e3t_n (:,:,:) , jpi*jpj*jpk, idex )!  T-cell thickness   
    9801053      END IF  
     1054  
     1055      IF( ln_wave .AND. ln_sdw ) THEN 
     1056         CALL histwrite( id_i, "sdzocrtx", kt, usd           , jpi*jpj*jpk, idex)     ! now StokesDrift i-velocity 
     1057         CALL histwrite( id_i, "sdmecrty", kt, vsd           , jpi*jpj*jpk, idex)     ! now StokesDrift j-velocity 
     1058         CALL histwrite( id_i, "sdvecrtz", kt, wsd           , jpi*jpj*jpk, idex)     ! now StokesDrift k-velocity 
     1059      ENDIF 
     1060 
    9811061      ! 3. Close the file 
    9821062      ! ----------------- 
Note: See TracChangeset for help on using the changeset viewer.