Changeset 8489


Ignore:
Timestamp:
2017-09-01T16:55:53+02:00 (3 years ago)
Author:
jpalmier
Message:

JPALM — gmed ticket #346 : improve MEDUSA conservation — import BBL bug fix from NEMO ticket #1932 to GO6 branch

Location:
branches/NERC/dev_r5518_GO6_MEDUSA_conserv/NEMOGCM/NEMO
Files:
9 edited

Legend:

Unmodified
Added
Removed
  • branches/NERC/dev_r5518_GO6_MEDUSA_conserv/NEMOGCM/NEMO/OPA_SRC/TRA/trabbl.F90

    r7771 r8489  
    549549      zmbk(:,:) = REAL( mbkv_d(:,:), wp )   ;   CALL lbc_lnk(zmbk,'V',1.)   ;   mbkv_d(:,:) = MAX( INT( zmbk(:,:) ), 1 ) 
    550550 
    551                                         !* sign of grad(H) at u- and v-points 
    552       mgrhu(jpi,:) = 0   ;   mgrhu(:,jpj) = 0   ;   mgrhv(jpi,:) = 0   ;   mgrhv(:,jpj) = 0 
     551      !! AXY (16/08/17): remove the following per George and Andrew bug-hunt 
     552      !!                                   !* sign of grad(H) at u- and v-points 
     553      !! mgrhu(jpi,:) = 0   ;   mgrhu(:,jpj) = 0   ;   mgrhv(jpi,:) = 0   ;   mgrhv(:,jpj) = 0 
     554      !! DO jj = 1, jpjm1 
     555      !!    DO ji = 1, jpim1 
     556      !!       mgrhu(ji,jj) = INT(  SIGN( 1.e0, gdept_0(ji+1,jj,mbkt(ji+1,jj)) - gdept_0(ji,jj,mbkt(ji,jj)) )  ) 
     557      !!       mgrhv(ji,jj) = INT(  SIGN( 1.e0, gdept_0(ji,jj+1,mbkt(ji,jj+1)) - gdept_0(ji,jj,mbkt(ji,jj)) )  ) 
     558      !!    END DO 
     559      !! END DO 
     560 
     561      !! AXY (16/08/17): add the following replacement per George and Andrew bug-hunt 
     562                                        !* sign of grad(H) at u- and  v-points; zero if grad(H) = 0 
     563      mgrhu(:,:) = 0   ;   mgrhv(:,:) = 0 
    553564      DO jj = 1, jpjm1 
    554565         DO ji = 1, jpim1 
    555             mgrhu(ji,jj) = INT(  SIGN( 1.e0, gdept_0(ji+1,jj,mbkt(ji+1,jj)) - gdept_0(ji,jj,mbkt(ji,jj)) )  ) 
    556             mgrhv(ji,jj) = INT(  SIGN( 1.e0, gdept_0(ji,jj+1,mbkt(ji,jj+1)) - gdept_0(ji,jj,mbkt(ji,jj)) )  ) 
     566            IF( gdept_0(ji+1,jj,mbkt(ji+1,jj)) - gdept_0(ji,jj,mbkt(ji,jj)) /= 0._wp ) THEN 
     567               mgrhu(ji,jj) = INT(  SIGN( 1.e0, & 
     568               gdept_0(ji+1,jj,mbkt(ji+1,jj)) - gdept_0(ji,jj,mbkt(ji,jj)) )  ) 
     569            ENDIF 
     570            !      
     571            IF( gdept_0(ji,jj+1,mbkt(ji,jj+1)) - gdept_0(ji,jj,mbkt(ji,jj)) /= 0._wp ) THEN 
     572               mgrhv(ji,jj) = INT(  SIGN( 1.e0, & 
     573               gdept_0(ji,jj+1,mbkt(ji,jj+1)) - gdept_0(ji,jj,mbkt(ji,jj)) )  ) 
     574            ENDIF 
    557575         END DO 
    558576      END DO 
  • branches/NERC/dev_r5518_GO6_MEDUSA_conserv/NEMOGCM/NEMO/TOP_SRC/MEDUSA/air_sea.F90

    r8442 r8489  
    66   !! History : 
    77   !!   -   ! 2017-04 (M. Stringer)        Code taken from trcbio_medusa.F90 
     8   !!   -   ! 2017-08 (A. Yool)            Add air-sea flux kill switch 
    89   !!---------------------------------------------------------------------- 
    910#if defined key_medusa 
     
    277278         ENDDO 
    278279      ENDDO 
     280#  endif 
     281 
     282#  if defined key_axy_killco2flux 
     283      !! AXY (18/08/17): single kill switch on air-sea CO2 flux for budget checking 
     284      f_co2flux(:,:) = 0. 
    279285#  endif 
    280286 
  • branches/NERC/dev_r5518_GO6_MEDUSA_conserv/NEMOGCM/NEMO/TOP_SRC/MEDUSA/bio_medusa_fin.F90

    r8442 r8489  
    66   !! History : 
    77   !!   -   ! 2017-04 (M. Stringer)        Code taken from trcbio_medusa.F90 
     8   !!   -   ! 2017-08 (A. Yool)            Amend bethic reservoir updating 
    89   !!---------------------------------------------------------------------- 
    910#if defined key_medusa 
     
    6263      REAL(wp) :: fq0,fq1,fq2,fq3 
    6364 
     65# if defined key_roam                      
     66      !!---------------------------------------------------------------------- 
     67      !! AXY (09/08/17): fix benthic submodel 
    6468      !!---------------------------------------------------------------------- 
    6569      !! Process benthic in/out fluxes 
    6670      !! These can be handled outside of the 3D calculations since the 
    67       !! benthic pools (and fluxes) are 2D in nature; this code is 
    68       !! (shamelessly) borrowed from corresponding code in the LOBSTER 
    69       !! model 
     71      !! benthic pools (and fluxes) are 2D in nature; this code was 
     72      !! developed with help from George Nurser (NOC); it cannot be run 
     73      !! in a configuration with variable time-stepping with depth 
    7074      !!---------------------------------------------------------------------- 
    7175      !! 
    72       !! IF(lwp) WRITE(numout,*) 'AXY: rdt = ', rdt 
     76      !! time-step calculation 
    7377      if (jorgben.eq.1) then 
    74          za_sed_n(:,:)  = zn_sed_n(:,:)  +                                 &  
    75                           ( f_sbenin_n(:,:)  + f_fbenin_n(:,:)  -          & 
    76                             f_benout_n(:,:)  ) * (rdt / 86400.) 
     78         za_sed_n(:,:)  = zb_sed_n(:,:)  + ((2. * (rdt / 86400.)) * & 
     79                          ( f_sbenin_n(:,:)  + f_fbenin_n(:,:)  - f_benout_n(:,:)  )) 
     80         za_sed_fe(:,:) = zb_sed_fe(:,:) + ((2. * (rdt / 86400.)) * & 
     81                          ( f_sbenin_fe(:,:) + f_fbenin_fe(:,:) - f_benout_fe(:,:) )) 
     82         za_sed_c(:,:)  = zb_sed_c(:,:)  + ((2. * (rdt / 86400.)) * & 
     83                          ( f_sbenin_c(:,:)  + f_fbenin_c(:,:)  - f_benout_c(:,:)  )) 
     84      endif 
     85      if (jinorgben.eq.1) then 
     86         za_sed_si(:,:) = zb_sed_si(:,:) + ((2. * (rdt / 86400.)) * & 
     87                          ( f_fbenin_si(:,:) - f_benout_si(:,:) )) 
     88         za_sed_ca(:,:) = zb_sed_ca(:,:) + ((2. * (rdt / 86400.)) * & 
     89                          ( f_fbenin_ca(:,:) - f_benout_ca(:,:) )) 
     90      endif 
     91      !! 
     92      !! time-level calculation 
     93      if (jorgben.eq.1) then 
     94         zb_sed_n(:,:)  = zn_sed_n(:,:)  + (atfp * & 
     95                          ( za_sed_n(:,:)  - (2. * zn_sed_n(:,:))  + zb_sed_n(:,:)  )) 
    7796         zn_sed_n(:,:)  = za_sed_n(:,:) 
    78          !! 
    79          za_sed_fe(:,:) = zn_sed_fe(:,:) +                                 & 
    80                           ( f_sbenin_fe(:,:) + f_fbenin_fe(:,:) -          & 
    81                             f_benout_fe(:,:) ) * (rdt / 86400.) 
     97         zb_sed_fe(:,:) = zn_sed_fe(:,:) + (atfp * & 
     98                          ( za_sed_fe(:,:) - (2. * zn_sed_fe(:,:)) + zb_sed_fe(:,:) )) 
    8299         zn_sed_fe(:,:) = za_sed_fe(:,:) 
    83          !! 
    84          za_sed_c(:,:)  = zn_sed_c(:,:)  +                                 & 
    85                           ( f_sbenin_c(:,:)  + f_fbenin_c(:,:)  -          & 
    86                             f_benout_c(:,:)  ) * (rdt / 86400.) 
     100         zb_sed_c(:,:)  = zn_sed_c(:,:)  + (atfp * & 
     101                          ( za_sed_c(:,:)  - (2. * zn_sed_c(:,:))  + zb_sed_c(:,:)  )) 
    87102         zn_sed_c(:,:)  = za_sed_c(:,:) 
    88103      endif 
    89104      if (jinorgben.eq.1) then 
    90          za_sed_si(:,:) = zn_sed_si(:,:) +                                 &  
    91                           ( f_fbenin_si(:,:) - f_benout_si(:,:) ) *        & 
    92                           (rdt / 86400.) 
     105         zb_sed_si(:,:) = zn_sed_si(:,:) + (atfp * & 
     106                          ( za_sed_si(:,:) - (2. * zn_sed_si(:,:)) + zb_sed_si(:,:) )) 
    93107         zn_sed_si(:,:) = za_sed_si(:,:) 
    94          !! 
    95          za_sed_ca(:,:) = zn_sed_ca(:,:) +                                 & 
    96                           ( f_fbenin_ca(:,:) - f_benout_ca(:,:) ) *        & 
    97                           (rdt / 86400.) 
     108         zb_sed_ca(:,:) = zn_sed_ca(:,:) + (atfp * & 
     109                          ( za_sed_ca(:,:) - (2. * zn_sed_ca(:,:)) + zb_sed_ca(:,:) )) 
    98110         zn_sed_ca(:,:) = za_sed_ca(:,:) 
    99111      endif 
    100       !! 
    101       if (ibenthic.eq.2) then 
    102          !! The code below (in this if ... then ... endif loop) is 
    103          !! effectively commented out because it does not work as  
    104          !! anticipated; it can be deleted at a later date 
    105          if (jorgben.eq.1) then 
    106             za_sed_n(:,:)  = ( f_sbenin_n(:,:)  + f_fbenin_n(:,:)  -       & 
    107                                f_benout_n(:,:)  ) * rdt 
    108             za_sed_fe(:,:) = ( f_sbenin_fe(:,:) + f_fbenin_fe(:,:) -       & 
    109                                f_benout_fe(:,:) ) * rdt 
    110             za_sed_c(:,:)  = ( f_sbenin_c(:,:)  + f_fbenin_c(:,:)  -       & 
    111                                f_benout_c(:,:)  ) * rdt 
    112          endif 
    113          if (jinorgben.eq.1) then 
    114             za_sed_si(:,:) = ( f_fbenin_si(:,:) - f_benout_si(:,:) ) * rdt 
    115             za_sed_ca(:,:) = ( f_fbenin_ca(:,:) - f_benout_ca(:,:) ) * rdt 
    116          endif 
    117          !! 
    118          !! Leap-frog scheme - only in explicit case, otherwise the  
    119          !! time stepping is already being done in trczdf 
    120          !! IF( l_trczdf_exp .AND. (ln_trcadv_cen2 .OR. ln_trcadv_tvd) ) THEN 
    121          !!    zfact = 2. * rdttra(jk) * FLOAT( ndttrc ) 
    122          !!    IF( neuler == 0 .AND. kt == nittrc000 )  zfact = rdttra(jk) *  
    123          !!                                             FLOAT(ndttrc) 
    124          !!    if (jorgben.eq.1) then 
    125          !!       za_sed_n(:,:)  = zb_sed_n(:,:)  + ( zfact * za_sed_n(:,:)  ) 
    126          !!      za_sed_fe(:,:) = zb_sed_fe(:,:) + ( zfact * za_sed_fe(:,:) ) 
    127          !!       za_sed_c(:,:)  = zb_sed_c(:,:)  + ( zfact * za_sed_c(:,:)  ) 
    128          !!    endif 
    129          !!    if (jinorgben.eq.1) then 
    130          !!       za_sed_si(:,:) = zb_sed_si(:,:) + ( zfact * za_sed_si(:,:) ) 
    131          !!       za_sed_ca(:,:) = zb_sed_ca(:,:) + ( zfact * za_sed_ca(:,:) ) 
    132          !!    endif 
    133          !! ENDIF 
    134          !!  
    135          !! Time filter and swap of arrays 
    136          IF( ln_trcadv_cen2 .OR. ln_trcadv_tvd  ) THEN ! centred or tvd scheme 
    137             IF( neuler == 0 .AND. kt == nittrc000 ) THEN 
    138                if (jorgben.eq.1) then 
    139                   zb_sed_n(:,:)  = zn_sed_n(:,:) 
    140                   zn_sed_n(:,:)  = za_sed_n(:,:) 
    141                   za_sed_n(:,:)  = 0.0 
    142                   !! 
    143                   zb_sed_fe(:,:) = zn_sed_fe(:,:) 
    144                   zn_sed_fe(:,:) = za_sed_fe(:,:) 
    145                   za_sed_fe(:,:) = 0.0 
    146                   !! 
    147                   zb_sed_c(:,:)  = zn_sed_c(:,:) 
    148                   zn_sed_c(:,:)  = za_sed_c(:,:) 
    149                   za_sed_c(:,:)  = 0.0 
    150                endif 
    151                if (jinorgben.eq.1) then 
    152                   zb_sed_si(:,:) = zn_sed_si(:,:) 
    153                   zn_sed_si(:,:) = za_sed_si(:,:) 
    154                   za_sed_si(:,:) = 0.0 
    155                   !! 
    156                   zb_sed_ca(:,:) = zn_sed_ca(:,:) 
    157                   zn_sed_ca(:,:) = za_sed_ca(:,:) 
    158                   za_sed_ca(:,:) = 0.0 
    159                endif 
    160             ELSE 
    161                if (jorgben.eq.1) then 
    162                   zb_sed_n(:,:)  = (atfp  *                                 & 
    163                                     ( zb_sed_n(:,:)  + za_sed_n(:,:)  )) +  & 
    164                                       (atfp1 * zn_sed_n(:,:) ) 
    165                   zn_sed_n(:,:)  = za_sed_n(:,:) 
    166                   za_sed_n(:,:)  = 0.0 
    167                   !! 
    168                   zb_sed_fe(:,:) = (atfp  *                                 & 
    169                                     ( zb_sed_fe(:,:) + za_sed_fe(:,:) )) +  & 
    170                                       (atfp1 * zn_sed_fe(:,:)) 
    171                   zn_sed_fe(:,:) = za_sed_fe(:,:) 
    172                   za_sed_fe(:,:) = 0.0 
    173                   !! 
    174                   zb_sed_c(:,:)  = (atfp  *                                 & 
    175                                     ( zb_sed_c(:,:)  + za_sed_c(:,:)  )) +  & 
    176                                       (atfp1 * zn_sed_c(:,:) ) 
    177                   zn_sed_c(:,:)  = za_sed_c(:,:) 
    178                   za_sed_c(:,:)  = 0.0 
    179                endif 
    180                if (jinorgben.eq.1) then 
    181                   zb_sed_si(:,:) = (atfp  *                                 & 
    182                                     ( zb_sed_si(:,:) + za_sed_si(:,:) )) +  & 
    183                                       (atfp1 * zn_sed_si(:,:)) 
    184                   zn_sed_si(:,:) = za_sed_si(:,:) 
    185                   za_sed_si(:,:) = 0.0 
    186                   !! 
    187                   zb_sed_ca(:,:) = (atfp  *                                 & 
    188                                     ( zb_sed_ca(:,:) + za_sed_ca(:,:) )) +  & 
    189                                       (atfp1 * zn_sed_ca(:,:)) 
    190                   zn_sed_ca(:,:) = za_sed_ca(:,:) 
    191                   za_sed_ca(:,:) = 0.0 
    192                endif 
    193             ENDIF 
    194          ELSE                   !  case of smolar scheme or muscl 
    195             if (jorgben.eq.1) then 
    196                zb_sed_n(:,:)  = za_sed_n(:,:) 
    197                zn_sed_n(:,:)  = za_sed_n(:,:) 
    198                za_sed_n(:,:)  = 0.0 
    199                !! 
    200                zb_sed_fe(:,:) = za_sed_fe(:,:) 
    201                zn_sed_fe(:,:) = za_sed_fe(:,:) 
    202                za_sed_fe(:,:) = 0.0 
    203                !! 
    204                zb_sed_c(:,:)  = za_sed_c(:,:) 
    205                zn_sed_c(:,:)  = za_sed_c(:,:) 
    206                za_sed_c(:,:)  = 0.0 
    207             endif 
    208             if (jinorgben.eq.1) then 
    209                zb_sed_si(:,:) = za_sed_si(:,:) 
    210                zn_sed_si(:,:) = za_sed_si(:,:) 
    211                za_sed_si(:,:) = 0.0 
    212                !! 
    213                zb_sed_ca(:,:) = za_sed_ca(:,:) 
    214                zn_sed_ca(:,:) = za_sed_ca(:,:) 
    215                za_sed_ca(:,:) = 0.0 
    216             endif 
    217          ENDIF 
    218       endif 
     112# endif       
    219113 
    220114#  if defined key_debug_medusa 
     
    253147                  fq1 = f_sbenin_n(ji,jj) + f_fbenin_n(ji,jj) 
    254148                  fq2 = fq0 + fq1 
    255                   IF (lwp) write (numout,'(a,2i3,a,3f15.10)')               & 
    256                      'AXY N   cons: (i,j)=',ji,jj,', (flx,ben,err)=',      & 
    257                      fq0,fq1,fq2 
     149                  fq3 = f_benout_n(ji,jj) 
     150                  if (lwp) write (numout,'a,2i3,a,4f15,5)')                   & 
     151                     'AXY N   cons: (i,j)=',ji,jj,', (flx,ben,err,out)=',      & 
     152                     fq0,fq1,fq2,fq3 
    258153               ENDIF 
    259154            ENDDO 
     
    266161                  fq1 = f_fbenin_si(ji,jj) 
    267162                  fq2 = fq0 + fq1 
    268                   IF (lwp) write (numout,'(a,2i3,a,3f15.10)')               & 
    269                      'AXY Si  cons: (i,j)=',ji,jj,', (flx,ben,err)=',     & 
    270                      fq0,fq1,fq2 
     163                  fq3 = f_benout_si(ji,jj) 
     164                  if (lwp) write (numout,'a,2i3,a,4f15,5)')                   & 
     165                     'AXY Si  cons: (i,j)=',ji,jj,', (flx,ben,err,out)=',     & 
     166                     fq0,fq1,fq2,fq3 
    271167               ENDIF 
    272168            ENDDO 
     
    278174                  fq0 = fflx_c(ji,jj) 
    279175                  fq1 = f_sbenin_c(ji,jj) + f_fbenin_c(ji,jj) + f_fbenin_ca(ji,jj) 
    280                   fq2 = f_co2flux(ji,jj) * e3t_n(ji,jj,1) 
     176                  fq2 = f_co2flux(ji,jj) * fse3t(ji,jj,1) 
    281177                  fq3 = fq0 + fq1 
    282                   IF (lwp) write (numout,'(a,2i3,a,4f15.10)')               & 
    283                     'AXY C   cons: (i,j)=',ji,jj,', (flx,ben,asf,err)=',  & 
    284                     fq0,fq1,fq2,fq3 
    285                ENDIF 
    286             ENDDO 
    287          ENDDO    
    288          !! alkalinity 
    289          DO jj = 2,jpjm1 
    290             DO ji = 2,jpim1 
    291                if (tmask(ji,jj,1) == 1) then 
    292                   fq0 = fflx_a(ji,jj) 
    293                   fq1 = 2.0 * f_fbenin_ca(ji,jj) 
    294                   fq2 = fq0 + fq1 
    295                   IF (lwp) write (numout,'(a,2i3,a,3f15.10)')               & 
    296                      'AXY alk cons: (i,j)=',ji,jj,', (flx,ben,err)=',     & 
    297                      fq0,fq1,fq2 
     178                  fq4 = f_benout_c(ji,jj) + f_benout_ca(ji,jj) 
     179                  if (lwp) write (numout,'a,2i3,a,5f15,5)')                   & 
     180                     'AXY C   cons: (i,j)=',ji,jj,', (flx,ben,asf,err,out)=', & 
     181                     fq0,fq1,fq2,fq3,fq4 
     182                ENDIF 
     183             ENDDO 
     184          ENDDO    
     185          !! alkalinity 
     186          DO jj = 2,jpjm1 
     187             DO ji = 2,jpim1 
     188                if (tmask(ji,jj,1) == 1) then 
     189                   fq0 = fflx_a(ji,jj) 
     190                   fq1 = 2.0 * f_fbenin_ca(ji,jj) 
     191                   fq2 = fq0 + fq1 
     192                   fq3 = 2.0 * f_benout_ca(ji,jj) 
     193                   if (lwp) write (numout,'a,2i3,a,4f15,5)')                   & 
     194                      'AXY alk cons: (i,j)=',ji,jj,', (flx,ben,err,out)=',     & 
     195                      fq0,fq1,fq2,fq3 
    298196               ENDIF 
    299197            ENDDO 
     
    358256         ENDIF 
    359257         IF ( med_diag%OCAL_LVL%dgsave ) THEN 
    360             CALL iom_put( "OCAL_LVL"  , fccd ) 
     258       !! AXY (16/08/17): hijack! 
     259            CALL iom_put( "OCAL_LVL"  , fchl_ml ) 
     260            !! CALL iom_put( "OCAL_LVL"  , fccd ) 
    361261         ENDIF 
    362262         IF ( med_diag%PN_JLIM%dgsave ) THEN 
  • branches/NERC/dev_r5518_GO6_MEDUSA_conserv/NEMOGCM/NEMO/TOP_SRC/MEDUSA/bio_medusa_init.F90

    r8442 r8489  
    66   !! History : 
    77   !!   -   ! 2017-04 (M. Stringer)        Code taken from trcbio_medusa.F90 
     8   !!   -   ! 2017-08 (A. Yool)            Add slow-sinking detrius variables 
    89   !!---------------------------------------------------------------------- 
    910#if defined key_medusa 
     
    161162      fprn_ml(:,:) = 0.0        !! mixed layer PP diagnostics 
    162163      fprd_ml(:,:) = 0.0        !! mixed layer PP diagnostics 
     164      !! AXY (16/08/17) 
     165      fchl_ml(:,:) = 0.0   !! mixed layer chlorophyll diagnostics 
    163166      !!  
    164167      fslownflux(:,:) = 0.0  
    165168      fslowcflux(:,:) = 0.0  
     169      !! 
     170      !! AXY (08/08/17): zero slow detritus fluxes 
     171      fslowsink(:,:)  = 0.0 
     172# if defined key_roam 
     173      fslowsinkc(:,:) = 0.0 
     174# endif       
    166175      !! 
    167176      !! allocate and initiate 2D diag 
  • branches/NERC/dev_r5518_GO6_MEDUSA_conserv/NEMOGCM/NEMO/TOP_SRC/MEDUSA/bio_medusa_mod.F90

    r8442 r8489  
    77   !! History : 
    88   !!   -   ! 2017-04 (M. Stringer)        Code taken from trcbio_medusa.F90 
     9   !!   -   ! 2017-08 (A. Yool)            Slow detritus, ML-avg chl variables 
    910   !!---------------------------------------------------------------------- 
    1011#if defined key_medusa 
     
    5455   !! AXY (01/03/10): add in mixed layer PP diagnostics 
    5556   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fprn_ml,fprd_ml 
     57   !! AXY (16/08/17): add in mixed layer chlorophyll diagnostic 
     58   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fchl_ml 
    5659   !! 
    5760   !! nutrient limiting factors 
     
    9497   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fregenfastc 
    9598# endif 
    96  
     99   !! 
     100   !! AXY (08/08/17): sinking of detritus moved here 
     101   REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::    fslowsink, fslowgain, fslowloss 
     102# if defined key_roam 
     103   REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::    fslowsinkc, fslowgainc, fslowlossc 
     104# endif 
     105   !! 
    97106   !! Particle flux 
    98107   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fdep1 
     
    287296               fjlim_pn(jpi,jpj),fjlim_pd(jpi,jpj),                   & 
    288297               fun_T(jpi,jpj),fun_Q10(jpi,jpj),                       & 
    289                fprn_ml(jpi,jpj),fprd_ml(jpi,jpj),                     & 
     298               fprn_ml(jpi,jpj),fprd_ml(jpi,jpj),fchl_ml(jpi,jpj),    & 
    290299               fnln(jpi,jpj),ffln2(jpi,jpj),                          & 
    291300               fnld(jpi,jpj),ffld(jpi,jpj),fsld(jpi,jpj),             & 
     
    316325               fregenfastc(jpi,jpj),                                  & 
    317326# endif 
     327          fslowsink(jpi,jpj),fslowgain(jpi,jpj),                 & 
     328               fslowloss(jpi,jpj),                                    & 
     329# if defined key_roam 
     330          fslowsinkc(jpi,jpj),fslowgainc(jpi,jpj),               & 
     331               fslowlossc(jpi,jpj),                                   & 
     332# endif 
    318333               fdep1(jpi,jpj),                                        & 
    319334               ftempn(jpi,jpj),ftempsi(jpi,jpj),ftempfe(jpi,jpj),     & 
  • branches/NERC/dev_r5518_GO6_MEDUSA_conserv/NEMOGCM/NEMO/TOP_SRC/MEDUSA/bio_medusa_update.F90

    r8442 r8489  
    66   !! History : 
    77   !!   -   ! 2017-04 (M. Stringer)        Code taken from trcbio_medusa.F90 
     8   !!   -   ! 2017-08 (A. Yool)            Amend slow-detritus bug 
     9   !!   -   ! 2017-08 (A. Yool)            Reformatting for clarity 
    810   !!---------------------------------------------------------------------- 
    911#if defined key_medusa 
     
    6062                                   fsil_cons, fsil_prod, fsdiss,             & 
    6163                                   ftempca, fthetad, fthetan,                & 
     64                                   fslowsink, fslowgain, fslowloss,          & ! AXY (22/08/17) 
     65                                   f_sbenin_n, f_sbenin_c,                   & 
    6266# if defined key_roam 
     67                                   fslowsinkc, fslowgainc, fslowlossc,       & ! AXY (22/08/17) 
    6368                                   fcar_cons, fcar_prod, fcomm_resp,         & 
    6469                                   fddc, fflx_a, fflx_c, fflx_o2, zoxy,      & 
     
    181186      ENDDO 
    182187 
    183       DO jj = 2,jpjm1 
    184          DO ji = 2,jpim1 
    185             if (tmask(ji,jj,jk) == 1) then 
    186                !! 
    187                !!---------------------------------------------------------- 
    188                !! detritus 
    189                btra(ji,jj,jpdet_lc) = b0 *                                      & 
    190                                           ! mort. losses  
    191                                    (fdpn(ji,jj) + ((1.0 - xfdfrac1) *        & 
    192                                                    fdpd(ji,jj)) +            & 
    193                                     fdzmi(ji,jj) +                           & 
    194                                     ((1.0 - xfdfrac2) * fdzme(ji,jj)) +      & 
    195                                           ! assim. inefficiency 
    196                                     ((1.0 - xbetan) * (finmi(ji,jj) +        & 
    197                                                        finme(ji,jj))) -      & 
    198                                           ! grazing and remin. 
    199                                     fgmid(ji,jj) - fgmed(ji,jj) -            & 
    200                                     fdd(ji,jj) +                             & 
    201                                           ! seafloor fast->slow 
    202                                     ffast2slown(ji,jj)) 
    203                !! 
     188      !!---------------------------------------------------------- 
     189      !! detritus 
     190      DO jj = 2,jpjm1 
     191         DO ji = 2,jpim1 
     192            if (tmask(ji,jj,jk) == 1) then 
     193               !! 
     194               btra(ji,jj,jpdet_lc) = b0 * (                           & 
     195                   fdpn(ji,jj)                                         & ! mort. losses  
     196                 + ((1.0 - xfdfrac1) * fdpd(ji,jj))                    & ! mort. losses  
     197                 + fdzmi(ji,jj)                                        & ! mort. losses 
     198                 + ((1.0 - xfdfrac2) * fdzme(ji,jj))                   & ! mort. losses 
     199                 + ((1.0 - xbetan) * (finmi(ji,jj) + finme(ji,jj)))    & ! assim. inefficiency 
     200                 - fgmid(ji,jj) - fgmed(ji,jj)                         & ! grazing 
     201                 - fdd(ji,jj)                                          & ! remin. 
     202                 + fslowgain(ji,jj) - fslowloss(ji,jj)                 & ! slow-sinking 
     203                 - (f_sbenin_n(ji,jj) / fse3t(ji,jj,jk))               & ! slow-sinking loss to seafloor 
     204                 + ffast2slown(ji,jj) )                                  ! seafloor fast->slow  
    204205            ENDIF 
    205206         ENDDO 
     
    305306                                          ffetop(ji,jj) + ffebot(ji,jj) -    & 
    306307                                          ffescav(ji,jj) ) 
     308            ENDIF 
     309         ENDDO 
     310      ENDDO 
     311 
    307312# if defined key_roam 
    308                !! 
    309                !!---------------------------------------------------------- 
    310                !! AXY (26/11/08): implicit detrital carbon change 
    311                btra(ji,jj,jpdtc_lc) = b0 * (                                    & 
    312                                             ! mort. losses 
    313                                          (xthetapn * fdpn(ji,jj)) +          & 
    314                                          ((1.0 - xfdfrac1) *                 & 
    315                                           (xthetapd * fdpd(ji,jj))) +        & 
    316                                          (xthetazmi * fdzmi(ji,jj)) +        & 
    317                                          ((1.0 - xfdfrac2) *                 & 
    318                                           (xthetazme * fdzme(ji,jj))) +      & 
    319                                              ! assim. inefficiency 
    320                                          ((1.0 - xbetac) *                   & 
    321                                           (ficmi(ji,jj) + ficme(ji,jj))) -   & 
    322                                              ! grazing and remin. 
    323                                          fgmidc(ji,jj) - fgmedc(ji,jj) -     & 
    324                                          fddc(ji,jj) +                       & 
    325                                              ! seafloor fast->slow 
    326                                          ffast2slowc(ji,jj) ) 
     313      !!---------------------------------------------------------- 
     314      !! AXY (26/11/08): implicit detrital carbon change 
     315      DO jj = 2,jpjm1 
     316         DO ji = 2,jpim1 
     317            if (tmask(ji,jj,jk) == 1) then  
     318               !! 
     319               btra(ji,jj,jpdtc_lc) = b0 * (                           & 
     320                   (xthetapn * fdpn(ji,jj))                            & ! mort. losses  
     321                 + ((1.0 - xfdfrac1) * (xthetapd * fdpd(ji,jj)))       & ! mort. losses  
     322                 + (xthetazmi * fdzmi(ji,jj))                          & ! mort. losses  
     323                 + ((1.0 - xfdfrac2) * (xthetazme * fdzme(ji,jj)))     & ! mort. losses  
     324                 + ((1.0 - xbetac) * (ficmi(ji,jj) + ficme(ji,jj)))    & ! assim. inefficiency 
     325                 - fgmidc(ji,jj) - fgmedc(ji,jj)                       & ! grazing 
     326                 - fddc(ji,jj)                                         & ! remin. 
     327                 + fslowgainc(ji,jj) - fslowlossc(ji,jj)               & ! slow-sinking 
     328                 - (f_sbenin_c(ji,jj) / fse3t(ji,jj,jk))               & ! slow-sinking loss to seafloor 
     329                 + ffast2slowc(ji,jj) )                                  ! seafloor fast->slow 
    327330            ENDIF 
    328331         ENDDO 
     
    575578                                                           f_o2flux(ji,jj)) 
    576579               endif 
     580            ENDIF 
     581         ENDDO 
     582      ENDDO 
    577583# endif 
    578             ENDIF 
    579          ENDDO 
    580       ENDDO 
    581584 
    582585# if defined key_debug_medusa 
  • branches/NERC/dev_r5518_GO6_MEDUSA_conserv/NEMOGCM/NEMO/TOP_SRC/MEDUSA/detritus.F90

    r8441 r8489  
    66   !! History : 
    77   !!   -   ! 2017-04 (M. Stringer)        Code taken from trcbio_medusa.F90 
     8   !!   -   ! 2017-08 (A. Yool)            Revise slow-sinking of detritus 
    89   !!---------------------------------------------------------------------- 
    910#if defined key_medusa 
     
    3536                                        f_sbenin_n, fdd,                   & 
    3637                                        idf, idfval,                       &    
    37 # if defined key_roam 
     38                                        fslowsink,                         & 
     39                                        fslowgain, fslowloss,              & 
     40# if defined key_roam 
     41                                        fslowsinkc,                        & 
     42                                        fslowgainc, fslowlossc,            & 
    3843                                        fddc,                              & 
    3944# endif 
    4045                                        fun_T, fun_Q10, zdet, zdtc 
    4146      USE detritus_fast_sink_mod, ONLY: detritus_fast_sink 
    42       USE dom_oce,                ONLY: mbathy, tmask 
     47      USE dom_oce,                ONLY: mbathy, e3t_0, e3t_n, gphit, tmask 
    4348      USE in_out_manager,         ONLY: lwp, numout 
    4449      USE par_oce,                ONLY: jpim1, jpjm1 
    4550      USE sms_medusa,             ONLY: jmd, jorgben, jsfd, vsed,          & 
    4651                                        xrfn, xmd, xmdc, xthetad 
     52 
     53   !!* Substitution 
     54#  include "domzgr_substitute.h90" 
    4755 
    4856      !! Level 
     
    123131         DO ji = 2,jpim1 
    124132            if (tmask(ji,jj,jk) == 1) then 
     133               !!---------------------------------------------------------------------- 
     134               !! Detritus sinking (AXY, 08/08/18) 
     135          !! Replaces slow-sinking done in trcsed_medusa.F90 
     136               !! 
     137               !! Uses the fslowsink variable to carry slow-sinking detritus from one 
     138               !! grid level to the next, variable fslowgain to "add" detritus sinking 
     139               !! from above and variable fslowloss to "subtract" detritus sinking out 
     140               !! to below; these variables appear in the differential equations of 
     141               !! detrital nitrogen and carbon below 
     142               !!---------------------------------------------------------------------- 
     143               !! 
     144               fslowgain(ji,jj)  = fslowsink(ji,jj) / fse3t(ji,jj,jk)                  ! = mmol N / m3 / d 
     145               if (jk.lt.mbathy(ji,jj)) then 
     146                  fslowloss(ji,jj)  = (zdet(ji,jj) * vsed * 86400.) / fse3t(ji,jj,jk)  ! = mmol N / m3 / d 
     147               else 
     148                  fslowloss(ji,jj)  = 0.                                               ! = mmol N / m3 / d 
     149               endif 
     150               fslowsink(ji,jj) = fslowloss(ji,jj) * fse3t(ji,jj,jk)                   ! = mmol N / m2 / d 
     151               !! 
     152#  if defined key_roam 
     153               fslowgainc(ji,jj) = fslowsinkc(ji,jj) / fse3t(ji,jj,jk)                 ! = mmol C / m3 / d 
     154               if (jk.lt.mbathy(ji,jj)) then 
     155                  fslowlossc(ji,jj) = (zdtc(ji,jj) * vsed * 86400.) / fse3t(ji,jj,jk)  ! = mmol C / m3 / d 
     156               else 
     157                  fslowlossc(ji,jj) = 0.                                               ! = mmol C / m3 / d 
     158               endif 
     159               fslowsinkc(ji,jj) = fslowlossc(ji,jj) * fse3t(ji,jj,jk)                 ! = mmol C / m2 / d 
     160#  endif 
     161            ENDIF 
     162         ENDDO 
     163      ENDDO 
     164 
     165      DO jj = 2,jpjm1 
     166         DO ji = 2,jpim1 
     167            if (tmask(ji,jj,jk) == 1) then 
    125168               !!--------------------------------------------------------- 
    126169               !! Detritus addition to benthos 
  • branches/NERC/dev_r5518_GO6_MEDUSA_conserv/NEMOGCM/NEMO/TOP_SRC/MEDUSA/phytoplankton.F90

    r8441 r8489  
    66   !! History : 
    77   !!   -   ! 2017-04 (M. Stringer)        Code taken from trcbio_medusa.F90 
     8   !!   -   ! 2017-08 (A. Yool)            Mean mixed layer chlorophyll 
    89   !!---------------------------------------------------------------------- 
    910#if defined key_medusa 
     
    4243                                   zchd, zchn, zdet, zdin, zdtc,         & 
    4344                                   zfer, zpds, zphd, zphn, zsil,         & 
    44                                    zzme, zzmi 
     45                                   zzme, zzmi, fchl_ml 
    4546      USE dom_oce,           ONLY: e3t_0, e3t_n, gdepw_0, gdepw_n, tmask 
    4647      USE in_out_manager,    ONLY: lwp, numout 
     
    373374               fprd_ml(ji,jj) = fprd_ml(ji,jj) + (fprd(ji,jj) * zphd(ji,jj) * & 
    374375                                                  fse3t(ji,jj,jk) * fq0) 
     376          !! AXY (16/08/17) 
     377          fchl_ml(ji,jj) = fchl_ml(ji,jj) + ((zchn(ji,jj) + zchd(ji,jj)) * & 
     378                                             (fse3t(ji,jj,jk) * fq0) / hmld(ji,jj)) 
    375379            ENDIF 
    376380         ENDDO 
  • branches/NERC/dev_r5518_GO6_MEDUSA_conserv/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcsms_medusa.F90

    r8074 r8489  
    88   !!              -   !  2008-11  (A. Yool) continuing adaptation for MEDUSA 
    99   !!              -   !  2010-03  (A. Yool) updated for branch inclusion 
     10   !!              -   !  2017-08  (A. Yool) amend for slow detritus bug 
    1011   !!---------------------------------------------------------------------- 
    1112#if defined key_medusa 
     
    8889#  endif 
    8990       
    90       CALL trc_sed_medusa( kt ) ! sedimentation model 
    91 #  if defined key_debug_medusa 
    92       IF(lwp) WRITE(numout,*) ' MEDUSA done trc_sed_medusa' 
    93       CALL flush(numout) 
    94 #  endif 
     91!! AXY (08/08/2017): remove call to buggy subroutine (now handled by detritus.F90) 
     92!!       CALL trc_sed_medusa( kt ) ! sedimentation model 
     93!! #  if defined key_debug_medusa 
     94!!       IF(lwp) WRITE(numout,*) ' MEDUSA done trc_sed_medusa' 
     95!!       CALL flush(numout) 
     96!! #  endif 
    9597# endif 
    9698 
Note: See TracChangeset for help on using the changeset viewer.