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 2528 for trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbcfwb.F90 – NEMO

Ignore:
Timestamp:
2010-12-27T18:33:53+01:00 (13 years ago)
Author:
rblod
Message:

Update NEMOGCM from branch nemo_v3_3_beta

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbcfwb.F90

    r2471 r2528  
    44   !! Ocean fluxes   : domain averaged freshwater budget 
    55   !!====================================================================== 
    6    !! History :  8.2  !  01-02  (E. Durand)  Original code 
    7    !!            8.5  !  02-06  (G. Madec)  F90: Free form and module 
    8    !!            9.0  !  06-08  (G. Madec)  Surface module 
    9    !!            9.2  !  09-07  (C. Talandier) emp mean s spread over erp area  
     6   !! History :  OPA  ! 2001-02  (E. Durand)  Original code 
     7   !!   NEMO     1.0  ! 2002-06  (G. Madec)  F90: Free form and module 
     8   !!            3.0  ! 2006-08  (G. Madec)  Surface module 
     9   !!            3.2  ! 2009-07  (C. Talandier) emp mean s spread over erp area  
    1010   !!---------------------------------------------------------------------- 
    1111 
     
    2323   USE lib_mpp         ! distribued memory computing library 
    2424   USE lbclnk          ! ocean lateral boundary conditions 
     25   USE lib_fortran 
    2526 
    2627   IMPLICIT NONE 
     
    3132   REAL(wp) ::   a_fwb_b            ! annual domain averaged freshwater budget 
    3233   REAL(wp) ::   a_fwb              ! for 2 year before (_b) and before year. 
    33    REAL(wp) ::   empold             ! empold to be suppressed 
     34   REAL(wp) ::   fwfold             ! fwfold to be suppressed 
    3435   REAL(wp) ::   area               ! global mean ocean surface (interior domain) 
    3536 
    36    REAL(wp), DIMENSION(jpi,jpj) ::   e1e2_i    ! area of the interior domain (e1t*e2t*tmask_i) 
     37   REAL(wp), DIMENSION(jpi,jpj) ::   e1e2    ! area of the interior domain (e1t*e2t) 
    3738 
    3839   !! * Substitutions 
     
    4041#  include "vectopt_loop_substitute.h90" 
    4142   !!---------------------------------------------------------------------- 
    42    !!  OPA 9.0 , LOCEAN-IPSL (2006)  
     43   !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
    4344   !! $Id$ 
    44    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     45   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    4546   !!---------------------------------------------------------------------- 
    4647CONTAINS 
     
    6566      INTEGER  ::   inum                  ! temporary logical unit 
    6667      INTEGER  ::   ikty, iyear           !  
    67       REAL(wp) ::   z_emp, z_emp_nsrf, zsum_emp, zsum_erp       ! temporary scalars 
     68      REAL(wp) ::   z_fwf, z_fwf_nsrf, zsum_fwf, zsum_erp       ! temporary scalars 
    6869      REAL(wp) ::   zsurf_neg, zsurf_pos, zsurf_tospread 
    6970      REAL(wp), DIMENSION(jpi,jpj) ::   ztmsk_neg, ztmsk_pos, ztmsk_tospread 
     
    7273      ! 
    7374      IF( kt == nit000 ) THEN 
    74          ! 
    7575         IF(lwp) THEN 
    7676            WRITE(numout,*) 
     
    7979            IF( kn_fwb == 1 )   WRITE(numout,*) '          instantaneously set to zero' 
    8080            IF( kn_fwb == 2 )   WRITE(numout,*) '          adjusted from previous year budget' 
    81             IF( kn_fwb == 3 )   WRITE(numout,*) '          emp set to zero and spread out over erp area' 
     81            IF( kn_fwb == 3 )   WRITE(numout,*) '          fwf set to zero and spread out over erp area' 
    8282         ENDIF 
    8383         ! 
    8484         IF( kn_fwb == 3 .AND. nn_sssr /= 2 )   CALL ctl_stop( 'sbc_fwb: nn_fwb = 3 requires nn_sssr = 2, we stop ' ) 
    8585         ! 
    86          e1e2_i(:,:) = e1t(:,:) * e2t(:,:) * tmask_i(:,:) 
    87          area = SUM( e1e2_i(:,:) ) 
    88          IF( lk_mpp )   CALL  mpp_sum( area    )   ! sum over the global domain 
    89          ! 
     86         e1e2(:,:) = e1t(:,:) * e2t(:,:)  
     87         area = glob_sum( e1e2(:,:) )           ! interior global domain surface 
    9088      ENDIF 
    9189       
     
    9391      SELECT CASE ( kn_fwb ) 
    9492      ! 
    95       CASE ( 0 )                                
    96          WRITE(ctmp1,*)'sbc_fwb : nn_fwb = ', kn_fwb, ' is not yet associated to an option, choose either 1/2' 
    97          CALL ctl_stop( ctmp1 ) 
     93      CASE ( 1 )                             !==  global mean fwf set to zero  ==! 
    9894         ! 
    99           
    100       ! 
    101       CASE ( 1 )                               ! global mean emp set to zero 
    10295         IF( MOD( kt-1, kn_fsbc ) == 0 ) THEN 
    103             z_emp = SUM( e1e2_i(:,:) * emp(:,:) ) / area 
    104             IF( lk_mpp )   CALL  mpp_sum( z_emp    )   ! sum over the global domain 
    105             emp (:,:) = emp (:,:) - z_emp 
    106             emps(:,:) = emps(:,:) - z_emp 
     96            z_fwf = glob_sum( e1e2(:,:) * ( emp(:,:) - rnf(:,:) ) ) / area   ! sum over the global domain 
     97            emp (:,:) = emp (:,:) - z_fwf  
     98            emps(:,:) = emps(:,:) - z_fwf  
    10799         ENDIF 
    108100         ! 
    109       CASE ( 2 )                               ! emp budget adjusted from the previous year 
    110          ! initialisation 
    111          IF( kt == nit000 ) THEN 
    112             ! Read the corrective factor on precipitations (empold) 
     101      CASE ( 2 )                             !==  fwf budget adjusted from the previous year  ==! 
     102         ! 
     103         IF( kt == nit000 ) THEN                   ! initialisation 
     104            !                                         ! Read the corrective factor on precipitations (fwfold) 
    113105            CALL ctl_opn( inum, 'EMPave_old.dat', 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE. ) 
    114106            READ ( inum, "(24X,I8,2ES24.16)" ) iyear, a_fwb_b, a_fwb 
    115107            CLOSE( inum ) 
    116             empold = a_fwb                  ! current year freshwater budget correction 
    117             !                               ! estimate from the previous year budget 
     108            fwfold = a_fwb                            ! current year freshwater budget correction 
     109            !                                         ! estimate from the previous year budget 
    118110            IF(lwp)WRITE(numout,*) 
    119             IF(lwp)WRITE(numout,*)'sbc_fwb : year = ',iyear  , ' freshwater budget correction = ', empold 
     111            IF(lwp)WRITE(numout,*)'sbc_fwb : year = ',iyear  , ' freshwater budget correction = ', fwfold 
    120112            IF(lwp)WRITE(numout,*)'          year = ',iyear-1, ' freshwater budget read       = ', a_fwb 
    121113            IF(lwp)WRITE(numout,*)'          year = ',iyear-2, ' freshwater budget read       = ', a_fwb_b 
    122114         ENDIF    
    123          !  
    124          ! Update empold if new year start 
     115         !                                         ! Update fwfold if new year start 
    125116         ikty = 365 * 86400 / rdttra(1)    !!bug  use of 365 days leap year or 360d year !!!!!!! 
    126117         IF( MOD( kt, ikty ) == 0 ) THEN 
    127118            a_fwb_b = a_fwb 
    128             a_fwb   = SUM( e1e2_i(:,:) * sshn(:,:) ) 
    129             IF( lk_mpp )   CALL  mpp_sum( a_fwb    )   ! sum over the global domain 
     119            a_fwb   = glob_sum( e1e2(:,:) * sshn(:,:) )   ! sum over the global domain 
    130120            a_fwb   = a_fwb * 1.e+3 / ( area * 86400. * 365. )     ! convert in Kg/m3/s = mm/s 
    131121!!gm        !                                                      !!bug 365d year  
    132             empold =  a_fwb                 ! current year freshwater budget correction 
    133             !                               ! estimate from the previous year budget 
     122            fwfold =  a_fwb                           ! current year freshwater budget correction 
     123            !                                         ! estimate from the previous year budget 
    134124         ENDIF 
    135125         !  
    136          ! correct the freshwater fluxes 
    137          IF( MOD( kt-1, kn_fsbc ) == 0 ) THEN 
    138             emp (:,:) = emp (:,:) + empold 
    139             emps(:,:) = emps(:,:) + empold 
     126         IF( MOD( kt-1, kn_fsbc ) == 0 ) THEN      ! correct the freshwater fluxes 
     127            emp (:,:) = emp (:,:) + fwfold 
     128            emps(:,:) = emps(:,:) + fwfold 
    140129         ENDIF 
    141130         ! 
    142          ! save empold value in a file 
    143          IF( kt == nitend .AND. lwp ) THEN 
     131         IF( kt == nitend .AND. lwp ) THEN         ! save fwfold value in a file 
    144132            CALL ctl_opn( inum, 'EMPave.dat', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE., narea ) 
    145133            WRITE( inum, "(24X,I8,2ES24.16)" ) nyear, a_fwb_b, a_fwb 
     
    147135         ENDIF 
    148136         ! 
    149       CASE ( 3 )                               ! global emp set to zero and spread out over erp area 
     137      CASE ( 3 )                             !==  global fwf set to zero and spread out over erp area  ==! 
    150138         ! 
    151139         IF( MOD( kt-1, kn_fsbc ) == 0 ) THEN 
    152             ! Select <0 and >0 area of erp 
    153             ztmsk_pos(:,:) = tmask_i(:,:) 
    154             WHERE( erp < 0.e0 ) ztmsk_pos = 0.e0 
     140            ztmsk_pos(:,:) = tmask_i(:,:)                      ! Select <0 and >0 area of erp 
     141            WHERE( erp < 0._wp )   ztmsk_pos = 0._wp 
    155142            ztmsk_neg(:,:) = tmask_i(:,:) - ztmsk_pos(:,:) 
    156  
    157             ! Area filled by <0 and >0 erp  
    158             zsurf_neg = SUM( e1e2_i(:,:)*ztmsk_neg(:,:) ) 
    159             zsurf_pos = SUM( e1e2_i(:,:)*ztmsk_pos(:,:) ) 
    160          
    161             ! emp global mean  
    162             z_emp = SUM( e1e2_i(:,:) * emp(:,:) ) / area 
    163143            ! 
    164             IF( lk_mpp )   CALL  mpp_sum( z_emp ) 
    165             IF( lk_mpp )   CALL  mpp_sum( zsurf_neg ) 
    166             IF( lk_mpp )   CALL  mpp_sum( zsurf_pos ) 
    167              
    168             IF( z_emp < 0.e0 ) THEN 
    169                 ! to spread out over >0 erp area to increase evaporation damping process 
    170                 zsurf_tospread = zsurf_pos 
     144            zsurf_neg = glob_sum( e1e2(:,:)*ztmsk_neg(:,:) )   ! Area filled by <0 and >0 erp  
     145            zsurf_pos = glob_sum( e1e2(:,:)*ztmsk_pos(:,:) ) 
     146            !                                                  ! fwf global mean  
     147            z_fwf = glob_sum( e1e2(:,:) * ( emp(:,:) - rnf(:,:) ) ) / area 
     148            !             
     149            IF( z_fwf < 0._wp ) THEN         ! spread out over >0 erp area to increase evaporation 
     150                zsurf_tospread      = zsurf_pos 
    171151                ztmsk_tospread(:,:) = ztmsk_pos(:,:) 
    172             ELSE 
    173                 ! to spread out over <0 erp area to increase precipitation damping process 
    174                 zsurf_tospread = zsurf_neg 
     152            ELSE                             ! spread out over <0 erp area to increase precipitation 
     153                zsurf_tospread      = zsurf_neg 
    175154                ztmsk_tospread(:,:) = ztmsk_neg(:,:) 
    176155            ENDIF 
    177  
    178             ! emp global mean over <0 or >0 erp area 
    179             zsum_emp = SUM( e1e2_i(:,:) * z_emp ) 
    180             IF( lk_mpp )   CALL  mpp_sum( zsum_emp ) 
    181             z_emp_nsrf =  zsum_emp / ( zsurf_tospread + rsmall ) 
    182             ! weight to respect erp field 2D structure  
    183             zsum_erp = SUM( ztmsk_tospread(:,:) * erp(:,:) * e1e2_i(:,:) ) 
    184             IF( lk_mpp )   CALL  mpp_sum( zsum_erp ) 
     156            ! 
     157            zsum_fwf   = glob_sum( e1e2(:,:) * z_fwf )         ! fwf global mean over <0 or >0 erp area 
     158!!gm :  zsum_fwf   = z_fwf * area   ???  it is right?  I think so.... 
     159            z_fwf_nsrf =  zsum_fwf / ( zsurf_tospread + rsmall ) 
     160            !                                                  ! weight to respect erp field 2D structure  
     161            zsum_erp = glob_sum( ztmsk_tospread(:,:) * erp(:,:) * e1e2(:,:) ) 
    185162            z_wgt(:,:) = ztmsk_tospread(:,:) * erp(:,:) / ( zsum_erp + rsmall ) 
    186  
    187             ! final correction term to apply 
    188             zerp_cor(:,:) = -1. * z_emp_nsrf * zsurf_tospread * z_wgt(:,:) 
    189  
     163            !                                                  ! final correction term to apply 
     164            zerp_cor(:,:) = -1. * z_fwf_nsrf * zsurf_tospread * z_wgt(:,:) 
     165            ! 
     166!!gm   ===>>>>  lbc_lnk should be useless as all the computation is done over the whole domain ! 
    190167            CALL lbc_lnk( zerp_cor, 'T', 1. ) 
    191  
     168            ! 
    192169            emp (:,:) = emp (:,:) + zerp_cor(:,:) 
    193170            emps(:,:) = emps(:,:) + zerp_cor(:,:) 
    194171            erp (:,:) = erp (:,:) + zerp_cor(:,:) 
    195              
    196             IF( nprint == 1 .AND. lwp ) THEN 
    197                IF( z_emp < 0.e0 ) THEN 
    198                   WRITE(numout,*)'       z_emp < 0' 
    199                   WRITE(numout,*)'       SUM(erp+)        = ', SUM( ztmsk_tospread(:,:)*erp(:,:)*e1e2_i(:,:) )*1.e-3,' m3.s-1' 
     172            ! 
     173            IF( nprint == 1 .AND. lwp ) THEN                   ! control print 
     174               IF( z_fwf < 0._wp ) THEN 
     175                  WRITE(numout,*)'   z_fwf < 0' 
     176                  WRITE(numout,*)'   SUM(erp+)     = ', SUM( ztmsk_tospread(:,:)*erp(:,:)*e1e2(:,:) )*1.e-9,' Sv' 
    200177               ELSE 
    201                    WRITE(numout,*)'      z_emp >= 0' 
    202                    WRITE(numout,*)'      SUM(erp-)        = ', SUM( ztmsk_tospread(:,:)*erp(:,:)*e1e2_i(:,:) )*1.e-3,' m3.s-1' 
     178                  WRITE(numout,*)'   z_fwf >= 0' 
     179                  WRITE(numout,*)'   SUM(erp-)     = ', SUM( ztmsk_tospread(:,:)*erp(:,:)*e1e2(:,:) )*1.e-9,' Sv' 
    203180               ENDIF 
    204                WRITE(numout,*)'      SUM(empG)        = ', SUM( z_emp*e1e2_i(:,:) )*1.e-3,' m3.s-1' 
    205                WRITE(numout,*)'      z_emp            = ', z_emp      ,' mm.s-1' 
    206                WRITE(numout,*)'      z_emp_nsrf       = ', z_emp_nsrf ,' mm.s-1' 
    207                WRITE(numout,*)'      MIN(zerp_cor)    = ', MINVAL(zerp_cor)  
    208                WRITE(numout,*)'      MAX(zerp_cor)    = ', MAXVAL(zerp_cor)  
     181               WRITE(numout,*)'   SUM(empG)     = ', SUM( z_fwf*e1e2(:,:) )*1.e-9,' Sv' 
     182               WRITE(numout,*)'   z_fwf         = ', z_fwf      ,' Kg/m2/s' 
     183               WRITE(numout,*)'   z_fwf_nsrf    = ', z_fwf_nsrf ,' Kg/m2/s' 
     184               WRITE(numout,*)'   MIN(zerp_cor) = ', MINVAL(zerp_cor)  
     185               WRITE(numout,*)'   MAX(zerp_cor) = ', MAXVAL(zerp_cor)  
    209186            ENDIF 
    210             ! 
    211187         ENDIF 
    212188         ! 
    213       CASE DEFAULT    ! you should never be there 
    214          WRITE(ctmp1,*)'sbc_fwb : nn_fwb = ', kn_fwb, ' is not permitted for the FreshWater Budget correction, choose either 1/2' 
    215          CALL ctl_stop( ctmp1 ) 
     189      CASE DEFAULT                           !==  you should never be there  ==! 
     190         CALL ctl_stop( 'sbc_fwb : wrong nn_fwb value for the FreshWater Budget correction, choose either 1, 2 or 3' ) 
    216191         ! 
    217192      END SELECT 
Note: See TracChangeset for help on using the changeset viewer.