Changeset 10127 for NEMO/trunk/src/TOP

Timestamp:

2018-09-13T17:27:54+02:00 (6 years ago)

Author:

cetlod

Message:

Bugfix to allow the use of on-the-fly interpolation for PISCES ancillary data, see ticket #2003

Location:

NEMO/trunk/src/TOP/PISCES/P4Z

Files:

• p4zsbc.F90 (modified) (7 diffs)
• p4zsed.F90 (modified) (7 diffs)

NEMO/trunk/src/TOP/PISCES/P4Z/p4zsbc.F90

@@ -75 +75 @@
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   hydrofe          !: Hydrothermal vent supply of iron
 
-   REAL(wp), PUBLIC ::   rivalkinput, rivdicinput, nitdepinput, sumdepsi
-   REAL(wp), PUBLIC ::   rivdininput, rivdipinput, rivdsiinput
+   REAL(wp), PUBLIC :: sedsilfrac, sedcalfrac
+   REAL(wp), PUBLIC :: rivalkinput, rivdicinput
+   REAL(wp), PUBLIC :: rivdininput, rivdipinput, rivdsiinput
 
    !! * Substitutions
@@ -110 +111 @@
          IF( kt == nit000 .OR. ( kt /= nit000 .AND. ntimes_dust > 1 ) ) THEN
             CALL fld_read( kt, 1, sf_dust )
-            IF( nn_ice_tr == -1 .AND. .NOT.ln_ironice ) THEN   ;   dust(:,:) = sf_dust(1)%fnow(:,:,1)
-            ELSE                                               ;   dust(:,:) = sf_dust(1)%fnow(:,:,1) * ( 1.-fr_i(:,:) )
+            IF( nn_ice_tr == -1 .AND. .NOT.ln_ironice ) THEN   ;   dust(:,:) = MAX( rtrn, sf_dust(1)%fnow(:,:,1) )
+            ELSE                                               ;   dust(:,:) = MAX( rtrn, sf_dust(1)%fnow(:,:,1) * ( 1.-fr_i(:,:) ) )
             ENDIF
          ENDIF
@@ -174 +175 @@
              zcoef = rno3 * 14E6 * ryyss
              CALL fld_read( kt, 1, sf_ndepo )
-             nitdep(:,:) = sf_ndepo(1)%fnow(:,:,1) / zcoef / e3t_n(:,:,1) 
+             nitdep(:,:) = MAX( rtrn, sf_ndepo(1)%fnow(:,:,1) / zcoef / e3t_n(:,:,1) )
          ENDIF
          IF( .NOT.ln_linssh ) THEN
            zcoef = rno3 * 14E6 * ryyss
-           nitdep(:,:) = sf_ndepo(1)%fnow(:,:,1) / zcoef / e3t_n(:,:,1) 
+           nitdep(:,:) = MAX( rtrn, sf_ndepo(1)%fnow(:,:,1) / zcoef / e3t_n(:,:,1) )
          ENDIF
       ENDIF
@@ -210 +211 @@
       REAL(wp), DIMENSION(nbtimes) :: zsteps                 ! times records
       REAL(wp), DIMENSION(:), ALLOCATABLE :: rivinput
-      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zdust, zndepo, zriver, zcmask
+      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zriver, zcmask
       !
       CHARACTER(len=100) ::  cn_dir          ! Root directory for location of ssr files
@@ -309 +310 @@
             CALL iom_open (  TRIM( sn_dust%clname ) , numdust )
             CALL iom_gettime( numdust, zsteps, kntime=ntimes_dust)  ! get number of record in file
-            ALLOCATE( zdust(jpi,jpj,ntimes_dust) )
-            DO jm = 1, ntimes_dust
-               CALL iom_get( numdust, jpdom_data, TRIM( sn_dust%clvar ), zdust(:,:,jm), jm )
-            END DO
-            CALL iom_close( numdust )
-            ztimes_dust = 1._wp / REAL(ntimes_dust, wp) 
-            sumdepsi = 0.e0
-            DO jm = 1, ntimes_dust
-               sumdepsi = sumdepsi + glob_sum( zdust(:,:,jm) * e1e2t(:,:) * tmask(:,:,1) * ztimes_dust )
-            END DO
-            sumdepsi = sumdepsi / ( nyear_len(1) * rday ) * 12. * 8.8 * 0.075 * mfrac / 28.1 
-            DEALLOCATE( zdust)
-         ENDIF
-      ELSE
-         sumdepsi  = 0._wp
+         END IF
       END IF
 
@@ -421 +408 @@
             CALL iom_open ( TRIM( sn_ndepo%clname ), numdepo )
             CALL iom_gettime( numdepo, zsteps, kntime=ntimes_ndep)
-            ALLOCATE( zndepo(jpi,jpj,ntimes_ndep) )
-            DO jm = 1, ntimes_ndep
-               CALL iom_get( numdepo, jpdom_data, TRIM( sn_ndepo%clvar ), zndepo(:,:,jm), jm )
-            END DO
-            CALL iom_close( numdepo )
-            ztimes_ndep = 1._wp / REAL(ntimes_ndep, wp) 
-            nitdepinput = 0._wp
-            DO jm = 1, ntimes_ndep
-              nitdepinput = nitdepinput + glob_sum( zndepo(:,:,jm) * e1e2t(:,:) * tmask(:,:,1) * ztimes_ndep )
-            ENDDO
-            nitdepinput = nitdepinput / rno3 / 14E6 
-            DEALLOCATE( zndepo)
-         ENDIF
-      ELSE
-         nitdepinput = 0._wp
+         ENDIF
       ENDIF
 
@@ -527 +500 @@
          WRITE(numout,*) '    DIC Supply : ', rivdicinput*1E3*12./1E12     ,' TgC/yr'
          WRITE(numout,*) 
-         WRITE(numout,*) '    Total input of elements from atmospheric supply'
-         WRITE(numout,*) '    ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
-         WRITE(numout,*) '    N Supply   : ', nitdepinput*rno3*1E3/1E12*14.,' TgN/yr'
-         WRITE(numout,*) 
-      ENDIF
+      ENDIF
+      !
+      sedsilfrac = 0.03     ! percentage of silica loss in the sediments
+      sedcalfrac = 0.6      ! percentage of calcite loss in the sediments
       !
    END SUBROUTINE p4z_sbc_init

NEMO/trunk/src/TOP/PISCES/P4Z/p4zsed.F90

@@ -49 +49 @@
       !
       INTEGER, INTENT(in) ::   kt, knt ! ocean time step
-      INTEGER  ::   ji, jj, jk, ikt
-      REAL(wp) ::   zsumsedsi, zsumsedpo4, zsumsedcal
-      REAL(wp) ::   zrivalk, zrivsil, zrivno3
+      INTEGER  ::  ji, jj, jk, ikt
+      REAL(wp) ::  zrivalk, zrivsil, zrivno3
       REAL(wp) ::  zwflux, zfminus, zfplus
       REAL(wp) ::  zlim, zfact, zfactcal
@@ -62 +61 @@
       !
       CHARACTER (len=25) :: charout
-      REAL(wp), DIMENSION(jpi,jpj    ) :: zwork1, zwork2, zwork3
-      REAL(wp), DIMENSION(jpi,jpj    ) :: zdenit2d, zbureff
+      REAL(wp), DIMENSION(jpi,jpj    ) :: zdenit2d, zbureff, zwork
       REAL(wp), DIMENSION(jpi,jpj    ) :: zwsbio3, zwsbio4, zwscal
       REAL(wp), DIMENSION(jpi,jpj    ) :: zsedcal, zsedsi, zsedc
@@ -82 +80 @@
       zdenit2d(:,:) = 0.e0
       zbureff (:,:) = 0.e0
-      zwork1  (:,:) = 0.e0
-      zwork2  (:,:) = 0.e0
-      zwork3  (:,:) = 0.e0
+      zwork   (:,:) = 0.e0
       zsedsi  (:,:) = 0.e0
       zsedcal (:,:) = 0.e0
@@ -255 +251 @@
                    &     + trb(ji,jj,ikt,jppoc) * zwsbio3(ji,jj) ) * 1E6
                  zbureff(ji,jj) = 0.013 + 0.53 * zflx**2 / ( 7.0 + zflx )**2
-                ENDIF
-              END DO
-           END DO 
-
-           ! Loss of biogenic silicon, Caco3 organic carbon in the sediments. 
-           ! First, the total loss is computed.
-           ! The factor for calcite comes from the alkalinity effect
-           ! -------------------------------------------------------------
-           DO jj = 1, jpj
-              DO ji = 1, jpi
-                 IF( tmask(ji,jj,1) == 1 ) THEN
-                    ikt = mbkt(ji,jj) 
-                    zwork1(ji,jj) = trb(ji,jj,ikt,jpgsi) * zwsbio4(ji,jj)
-                    zwork2(ji,jj) = trb(ji,jj,ikt,jpgoc) * zwsbio4(ji,jj) + trb(ji,jj,ikt,jppoc) * zwsbio3(ji,jj) 
-                    ! For calcite, burial efficiency is made a function of saturation
-                    zfactcal      = MIN( excess(ji,jj,ikt), 0.2 )
-                    zfactcal      = MIN( 1., 1.3 * ( 0.2 - zfactcal ) / ( 0.4 - zfactcal ) )
-                    zwork3(ji,jj) = trb(ji,jj,ikt,jpcal) * zwscal(ji,jj) * 2.e0 * zfactcal
-                ENDIF
-            END DO
-         END DO
-         zsumsedsi  = glob_sum( zwork1(:,:) * e1e2t(:,:) ) * r1_rday
-         zsumsedpo4 = glob_sum( zwork2(:,:) * e1e2t(:,:) ) * r1_rday
-         zsumsedcal = glob_sum( zwork3(:,:) * e1e2t(:,:) ) * r1_rday
+              ENDIF
+            END DO
+         END DO 
          !
       ENDIF
@@ -285 +260 @@
       ! Thus, the amount of silica lost in the sediments equal the supply at the surface (dust+rivers)
       ! ------------------------------------------------------
-      IF( .NOT.lk_sed )  zrivsil =  1._wp - ( sumdepsi + rivdsiinput * r1_ryyss ) / ( zsumsedsi + rtrn )
+      IF( .NOT.lk_sed )  zrivsil = 1._wp - sedsilfrac
 
       DO jj = 1, jpj
@@ -312 +287 @@
                zfactcal = MIN( excess(ji,jj,ikt), 0.2 )
                zfactcal = MIN( 1., 1.3 * ( 0.2 - zfactcal ) / ( 0.4 - zfactcal ) )
-               zrivalk  =  1._wp - ( rivalkinput * r1_ryyss ) * zfactcal / ( zsumsedcal + rtrn )
+               zrivalk  = sedcalfrac * zfactcal
                tra(ji,jj,ikt,jptal) =  tra(ji,jj,ikt,jptal) + zcaloss * zrivalk * 2.0
                tra(ji,jj,ikt,jpdic) =  tra(ji,jj,ikt,jpdic) + zcaloss * zrivalk
@@ -492 +467 @@
             IF( iom_use("Nfix"   ) ) CALL iom_put( "Nfix", nitrpot(:,:,:) * nitrfix * rno3 * zfact * tmask(:,:,:) )  ! nitrogen fixation 
             IF( iom_use("INTNFIX") ) THEN   ! nitrogen fixation rate in ocean ( vertically integrated )
-               zwork1(:,:) = 0.
+               zwork(:,:) = 0.
                DO jk = 1, jpkm1
-                 zwork1(:,:) = zwork1(:,:) + nitrpot(:,:,jk) * nitrfix * rno3 * zfact * e3t_n(:,:,jk) * tmask(:,:,jk)
+                 zwork(:,:) = zwork(:,:) + nitrpot(:,:,jk) * nitrfix * rno3 * zfact * e3t_n(:,:,jk) * tmask(:,:,jk)
                ENDDO
-               CALL iom_put( "INTNFIX" , zwork1 ) 
+               CALL iom_put( "INTNFIX" , zwork ) 
             ENDIF
             IF( iom_use("SedCal" ) ) CALL iom_put( "SedCal", zsedcal(:,:) * zfact )