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 10377 for NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src – NEMO

Ignore:
Timestamp:
2018-12-10T08:45:39+01:00 (5 years ago)
Author:
smasson
Message:

dev_r10164_HPC09_ESIWACE_PREP_MERGE: merge with trunk@10376, see #2133

Location:
NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src
Files:
8 edited
1 copied

Legend:

Unmodified
Added
Removed

  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/IOM/iom_nf90.F90

    r10358 r10377  
    691691      INTEGER, DIMENSION(4) :: idimid               ! dimensions id 
    692692      CHARACTER(LEN=256)    :: clinfo               ! info character 
    693       CHARACTER(LEN= 12), DIMENSION(4) :: cltmp     ! temporary character 
     693      CHARACTER(LEN= 12), DIMENSION(5) :: cltmp     ! temporary character 
    694694      INTEGER               :: if90id               ! nf90 file identifier 
    695695      INTEGER               :: idmy                 ! dummy variable 
    … …  
    716716         ENDIF 
    717717         ! define the dimension variables if it is not already done 
    718          IF(iom_file(kiomid)%nlev == jpk ) THEN 
    719           cltmp = (/ 'nav_lon     ', 'nav_lat     ', 'nav_lev     ', 'time_counter' /) 
    720          ELSE 
    721           cltmp = (/ 'nav_lon     ', 'nav_lat     ', 'numcat      ', 'time_counter' /) 
    722          ENDIF 
     718         cltmp = (/ 'nav_lon', 'nav_lat', 'nav_lev', 'time_counter', 'numcat' /) 
    723719         CALL iom_nf90_check(NF90_DEF_VAR( if90id, TRIM(cltmp(1)), NF90_FLOAT , (/ 1, 2 /), iom_file(kiomid)%nvid(1) ), clinfo) 
    724720         CALL iom_nf90_check(NF90_DEF_VAR( if90id, TRIM(cltmp(2)), NF90_FLOAT , (/ 1, 2 /), iom_file(kiomid)%nvid(2) ), clinfo) 
    … …  
    728724         iom_file(kiomid)%nvars       = 4 
    729725         iom_file(kiomid)%luld(1:4)   = (/ .FALSE., .FALSE., .FALSE., .TRUE. /) 
    730          iom_file(kiomid)%cn_var(1:4) = cltmp 
    731          iom_file(kiomid)%ndims(1:4)  = (/ 2, 2, 1, 1 /)   
     726         iom_file(kiomid)%cn_var(1:4) = cltmp(1:4) 
     727         iom_file(kiomid)%ndims(1:4)  = (/ 2, 2, 1, 1 /) 
     728         IF( NF90_INQ_DIMID( if90id, 'numcat', idmy ) == nf90_noerr ) THEN   ! add a 5th variable corresponding to the 5th dimension 
     729            CALL iom_nf90_check(NF90_DEF_VAR( if90id, TRIM(cltmp(5)), NF90_FLOAT , (/ 5 /), iom_file(kiomid)%nvid(5) ), clinfo) 
     730            iom_file(kiomid)%nvars     = 5 
     731            iom_file(kiomid)%luld(5)   = .FALSE. 
     732            iom_file(kiomid)%cn_var(5) = cltmp(5) 
     733            iom_file(kiomid)%ndims(5)  = 1 
     734         ENDIF 
    732735         ! trick: defined to 0 to say that dimension variables are defined but not yet written 
    733736         iom_file(kiomid)%dimsz(1, 1)  = 0    
    … …  
    841844               CALL iom_nf90_check( NF90_INQ_VARID( if90id, 'nav_lat'     , idmy )         , clinfo ) 
    842845               CALL iom_nf90_check( NF90_PUT_VAR  ( if90id, idmy, gphit(ix1:ix2, iy1:iy2) ), clinfo ) 
    843                IF(iom_file(kiomid)%nlev == jpk ) THEN  
    844                   !NEMO 
    845                   CALL iom_nf90_check( NF90_INQ_VARID( if90id, 'nav_lev'     , idmy ), clinfo ) 
    846                   CALL iom_nf90_check( NF90_PUT_VAR  ( if90id, idmy, gdept_1d       ), clinfo ) 
    847                ELSE 
    848                   CALL iom_nf90_check( NF90_INQ_VARID( if90id, 'numcat'     , idmy ), clinfo) 
     846               CALL iom_nf90_check( NF90_INQ_VARID( if90id, 'nav_lev'     , idmy ), clinfo ) 
     847               CALL iom_nf90_check( NF90_PUT_VAR  ( if90id, idmy, gdept_1d       ), clinfo ) 
     848               IF( NF90_INQ_VARID( if90id, 'numcat', idmy ) == nf90_noerr ) THEN 
    849849                  CALL iom_nf90_check( NF90_PUT_VAR  ( if90id, idmy, (/ (idlv, idlv = 1,iom_file(kiomid)%nlev) /)), clinfo ) 
    850850               ENDIF 

  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/SBC/geo2ocean.F90

    r10297 r10377  
    1010   !!            3.7  !  11-2015  (G. Madec)  remove the unused repere and repcmo routines 
    1111   !!---------------------------------------------------------------------- 
    12 #if defined key_agrif 
    13 !clem: these lines do not seem necessary anymore 
    14 !!DIR$ OPTIMIZE (-O 1)  ! cray formulation 
    15 # if defined __INTEL_COMPILER 
    16 !acc: still breaks on at least one Ivybridge cluster with ifort 17.0.4 without this directive 
    17 !DIR$ OPTIMIZE:1        ! intel formulation 
    18 # endif 
    19 #endif 
    2012   !!---------------------------------------------------------------------- 
    2113   !!   rot_rep       : Rotate the Repere: geographic grid <==> stretched coordinates grid 
    … …  
    8173         IF(lwp) WRITE(numout,*) ' ~~~~~~~~    ' 
    8274         ! 
    83          CALL angle       ! initialization of the transformation 
     75         CALL angle( glamt, gphit, glamu, gphiu, glamv, gphiv, glamf, gphif )       ! initialization of the transformation 
    8476         lmust_init = .FALSE. 
    8577      ENDIF 
    … …  
    126118 
    127119 
    128    SUBROUTINE angle 
     120   SUBROUTINE angle( plamt, pphit, plamu, pphiu, plamv, pphiv, plamf, pphif ) 
    129121      !!---------------------------------------------------------------------- 
    130122      !!                  ***  ROUTINE angle  *** 
    … …  
    138130      !! ** Action  : - gsint, gcost, gsinu, gcosu, gsinv, gcosv, gsinf, gcosf 
    139131      !!---------------------------------------------------------------------- 
     132      ! WARNING: for an unexplained reason, we need to pass all glam, gphi arrays as input parameters in 
     133      !          order to get AGRIF working with -03 compilation option 
     134      REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) :: plamt, pphit, plamu, pphiu, plamv, pphiv, plamf, pphif   
     135      ! 
    140136      INTEGER  ::   ji, jj   ! dummy loop indices 
    141137      INTEGER  ::   ierr     ! local integer 
    … …  
    167163         DO ji = fs_2, jpi   ! vector opt. 
    168164            !                   
    169             zlam = glamt(ji,jj)     ! north pole direction & modulous (at t-point) 
    170             zphi = gphit(ji,jj) 
     165            zlam = plamt(ji,jj)     ! north pole direction & modulous (at t-point) 
     166            zphi = pphit(ji,jj) 
    171167            zxnpt = 0. - 2. * COS( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. ) 
    172168            zynpt = 0. - 2. * SIN( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. ) 
    173169            znnpt = zxnpt*zxnpt + zynpt*zynpt 
    174170            ! 
    175             zlam = glamu(ji,jj)     ! north pole direction & modulous (at u-point) 
    176             zphi = gphiu(ji,jj) 
     171            zlam = plamu(ji,jj)     ! north pole direction & modulous (at u-point) 
     172            zphi = pphiu(ji,jj) 
    177173            zxnpu = 0. - 2. * COS( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. ) 
    178174            zynpu = 0. - 2. * SIN( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. ) 
    179175            znnpu = zxnpu*zxnpu + zynpu*zynpu 
    180176            ! 
    181             zlam = glamv(ji,jj)     ! north pole direction & modulous (at v-point) 
    182             zphi = gphiv(ji,jj) 
     177            zlam = plamv(ji,jj)     ! north pole direction & modulous (at v-point) 
     178            zphi = pphiv(ji,jj) 
    183179            zxnpv = 0. - 2. * COS( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. ) 
    184180            zynpv = 0. - 2. * SIN( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. ) 
    185181            znnpv = zxnpv*zxnpv + zynpv*zynpv 
    186182            ! 
    187             zlam = glamf(ji,jj)     ! north pole direction & modulous (at f-point) 
    188             zphi = gphif(ji,jj) 
     183            zlam = plamf(ji,jj)     ! north pole direction & modulous (at f-point) 
     184            zphi = pphif(ji,jj) 
    189185            zxnpf = 0. - 2. * COS( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. ) 
    190186            zynpf = 0. - 2. * SIN( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. ) 
    191187            znnpf = zxnpf*zxnpf + zynpf*zynpf 
    192188            ! 
    193             zlam = glamv(ji,jj  )   ! j-direction: v-point segment direction (around t-point) 
    194             zphi = gphiv(ji,jj  ) 
    195             zlan = glamv(ji,jj-1) 
    196             zphh = gphiv(ji,jj-1) 
     189            zlam = plamv(ji,jj  )   ! j-direction: v-point segment direction (around t-point) 
     190            zphi = pphiv(ji,jj  ) 
     191            zlan = plamv(ji,jj-1) 
     192            zphh = pphiv(ji,jj-1) 
    197193            zxvvt =  2. * COS( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. )   & 
    198194               &  -  2. * COS( rad*zlan ) * TAN( rpi/4. - rad*zphh/2. ) 
    … …  
    202198            znvvt = MAX( znvvt, 1.e-14 ) 
    203199            ! 
    204             zlam = glamf(ji,jj  )   ! j-direction: f-point segment direction (around u-point) 
    205             zphi = gphif(ji,jj  ) 
    206             zlan = glamf(ji,jj-1) 
    207             zphh = gphif(ji,jj-1) 
     200            zlam = plamf(ji,jj  )   ! j-direction: f-point segment direction (around u-point) 
     201            zphi = pphif(ji,jj  ) 
     202            zlan = plamf(ji,jj-1) 
     203            zphh = pphif(ji,jj-1) 
    208204            zxffu =  2. * COS( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. )   & 
    209205               &  -  2. * COS( rad*zlan ) * TAN( rpi/4. - rad*zphh/2. ) 
    … …  
    213209            znffu = MAX( znffu, 1.e-14 ) 
    214210            ! 
    215             zlam = glamf(ji  ,jj)   ! i-direction: f-point segment direction (around v-point) 
    216             zphi = gphif(ji  ,jj) 
    217             zlan = glamf(ji-1,jj) 
    218             zphh = gphif(ji-1,jj) 
     211            zlam = plamf(ji  ,jj)   ! i-direction: f-point segment direction (around v-point) 
     212            zphi = pphif(ji  ,jj) 
     213            zlan = plamf(ji-1,jj) 
     214            zphh = pphif(ji-1,jj) 
    219215            zxffv =  2. * COS( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. )   & 
    220216               &  -  2. * COS( rad*zlan ) * TAN( rpi/4. - rad*zphh/2. ) 
    … …  
    224220            znffv = MAX( znffv, 1.e-14 ) 
    225221            ! 
    226             zlam = glamu(ji,jj+1)   ! j-direction: u-point segment direction (around f-point) 
    227             zphi = gphiu(ji,jj+1) 
    228             zlan = glamu(ji,jj  ) 
    229             zphh = gphiu(ji,jj  ) 
     222            zlam = plamu(ji,jj+1)   ! j-direction: u-point segment direction (around f-point) 
     223            zphi = pphiu(ji,jj+1) 
     224            zlan = plamu(ji,jj  ) 
     225            zphh = pphiu(ji,jj  ) 
    230226            zxuuf =  2. * COS( rad*zlam ) * TAN( rpi/4. - rad*zphi/2. )   & 
    231227               &  -  2. * COS( rad*zlan ) * TAN( rpi/4. - rad*zphh/2. ) 
    … …  
    257253      DO jj = 2, jpjm1 
    258254         DO ji = fs_2, jpi   ! vector opt. 
    259             IF( MOD( ABS( glamv(ji,jj) - glamv(ji,jj-1) ), 360. ) < 1.e-8 ) THEN 
     255            IF( MOD( ABS( plamv(ji,jj) - plamv(ji,jj-1) ), 360. ) < 1.e-8 ) THEN 
    260256               gsint(ji,jj) = 0. 
    261257               gcost(ji,jj) = 1. 
    262258            ENDIF 
    263             IF( MOD( ABS( glamf(ji,jj) - glamf(ji,jj-1) ), 360. ) < 1.e-8 ) THEN 
     259            IF( MOD( ABS( plamf(ji,jj) - plamf(ji,jj-1) ), 360. ) < 1.e-8 ) THEN 
    264260               gsinu(ji,jj) = 0. 
    265261               gcosu(ji,jj) = 1. 
    266262            ENDIF 
    267             IF(      ABS( gphif(ji,jj) - gphif(ji-1,jj) )         < 1.e-8 ) THEN 
     263            IF(      ABS( pphif(ji,jj) - pphif(ji-1,jj) )         < 1.e-8 ) THEN 
    268264               gsinv(ji,jj) = 0. 
    269265               gcosv(ji,jj) = 1. 
    270266            ENDIF 
    271             IF( MOD( ABS( glamu(ji,jj) - glamu(ji,jj+1) ), 360. ) < 1.e-8 ) THEN 
     267            IF( MOD( ABS( plamu(ji,jj) - plamu(ji,jj+1) ), 360. ) < 1.e-8 ) THEN 
    272268               gsinf(ji,jj) = 0. 
    273269               gcosf(ji,jj) = 1. 
    … …  
    457453         IF(lwp) WRITE(numout,*) ' obs_rot : geographic <--> stretched' 
    458454         IF(lwp) WRITE(numout,*) ' ~~~~~~~   coordinate transformation' 
    459          CALL angle       ! initialization of the transformation 
     455         CALL angle( glamt, gphit, glamu, gphiu, glamv, gphiv, glamf, gphif )       ! initialization of the transformation 
    460456         lmust_init = .FALSE. 
    461457      ENDIF 

  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p4zmeso.F90

    r10368 r10377  
    252252         DEALLOCATE( zw3d ) 
    253253      ENDIF 
     254      ! 
     255      IF (ln_ligand)  DEALLOCATE( zz2ligprod ) 
    254256      ! 
    255257      IF(ln_ctl)   THEN  ! print mean trends (used for debugging) 

  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p4zmicro.F90

    r10368 r10377  
    205205      ENDIF 
    206206      ! 
     207      IF (ln_ligand)  DEALLOCATE( zzligprod ) 
     208      ! 
    207209      IF(ln_ctl) THEN      ! print mean trends (used for debugging) 
    208210         WRITE(charout, FMT="('micro')") 

  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p4zsink.F90

    r10368 r10377  
    1616   USE trc             !  passive tracers common variables  
    1717   USE sms_pisces      !  PISCES Source Minus Sink variables 
     18   USE trcsink         !  General routine to compute sedimentation 
    1819   USE prtctl_trc      !  print control for debugging 
    1920   USE iom             !  I/O manager 
    … …  
    5960      !!--------------------------------------------------------------------- 
    6061      INTEGER, INTENT(in) :: kt, knt 
    61       INTEGER  ::   ji, jj, jk, jit 
    62       INTEGER  ::   iiter1, iiter2 
    63       REAL(wp) ::   zagg1, zagg2, zagg3, zagg4 
    64       REAL(wp) ::   zagg , zaggfe, zaggdoc, zaggdoc2, zaggdoc3 
    65       REAL(wp) ::   zfact, zwsmax, zmax 
     62      INTEGER  ::   ji, jj, jk 
    6663      CHARACTER (len=25) :: charout 
     64      REAL(wp) :: zmax, zfact 
    6765      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d 
    6866      REAL(wp), ALLOCATABLE, DIMENSION(:,:  ) :: zw2d 
    … …  
    7068      ! 
    7169      IF( ln_timing )   CALL timing_start('p4z_sink') 
    72  
    7370 
    7471      ! Initialization of some global variables 
    … …  
    9794 
    9895      ! 
    99       ! OA This is (I hope) a temporary solution for the problem that may  
    100       ! OA arise in specific situation where the CFL criterion is broken  
    101       ! OA for vertical sedimentation of particles. To avoid this, a time 
    102       ! OA splitting algorithm has been coded. A specific maximum 
    103       ! OA iteration number is provided and may be specified in the namelist  
    104       ! OA This is to avoid very large iteration number when explicit free 
    105       ! OA surface is used (for instance). When niter?max is set to 1,  
    106       ! OA this computation is skipped. The crude old threshold method is  
    107       ! OA then applied. This also happens when niter exceeds nitermax. 
    108       IF( MAX( niter1max, niter2max ) == 1 ) THEN 
    109         iiter1 = 1 
    110         iiter2 = 1 
    111       ELSE 
    112         iiter1 = 1 
    113         iiter2 = 1 
    114         DO jk = 1, jpkm1 
    115           DO jj = 1, jpj 
    116              DO ji = 1, jpi 
    117                 IF( tmask(ji,jj,jk) == 1) THEN 
    118                    zwsmax =  0.5 * e3t_n(ji,jj,jk) / xstep 
    119                    iiter1 =  MAX( iiter1, INT( wsbio3(ji,jj,jk) / zwsmax ) ) 
    120                    iiter2 =  MAX( iiter2, INT( wsbio4(ji,jj,jk) / zwsmax ) ) 
    121                 ENDIF 
    122              END DO 
    123           END DO 
    124         END DO 
    125         IF( lk_mpp ) THEN 
    126            CALL mpp_max( 'p4zsink', iiter1 ) 
    127            CALL mpp_max( 'p4zsink', iiter2 ) 
    128         ENDIF 
    129         iiter1 = MIN( iiter1, niter1max ) 
    130         iiter2 = MIN( iiter2, niter2max ) 
    131       ENDIF 
    132  
    133       DO jk = 1,jpkm1 
    134          DO jj = 1, jpj 
    135             DO ji = 1, jpi 
    136                IF( tmask(ji,jj,jk) == 1 ) THEN 
    137                  zwsmax = 0.5 * e3t_n(ji,jj,jk) / xstep 
    138                  wsbio3(ji,jj,jk) = MIN( wsbio3(ji,jj,jk), zwsmax * REAL( iiter1, wp ) ) 
    139                  wsbio4(ji,jj,jk) = MIN( wsbio4(ji,jj,jk), zwsmax * REAL( iiter2, wp ) ) 
    140                ENDIF 
    141             END DO 
    142          END DO 
    143       END DO 
    144  
    14596      !  Initializa to zero all the sinking arrays  
    14697      !   ----------------------------------------- 
    … …  
    154105      !   Compute the sedimentation term using p4zsink2 for all the sinking particles 
    155106      !   ----------------------------------------------------- 
    156       DO jit = 1, iiter1 
    157         CALL p4z_sink2( wsbio3, sinking , jppoc, iiter1 ) 
    158         CALL p4z_sink2( wsbio3, sinkfer , jpsfe, iiter1 ) 
    159       END DO 
    160  
    161       DO jit = 1, iiter2 
    162         CALL p4z_sink2( wsbio4, sinking2, jpgoc, iiter2 ) 
    163         CALL p4z_sink2( wsbio4, sinkfer2, jpbfe, iiter2 ) 
    164         CALL p4z_sink2( wsbio4, sinksil , jpgsi, iiter2 ) 
    165         CALL p4z_sink2( wsbio4, sinkcal , jpcal, iiter2 ) 
    166       END DO 
     107      CALL trc_sink( kt, wsbio3, sinking , jppoc, rfact2 ) 
     108      CALL trc_sink( kt, wsbio3, sinkfer , jpsfe, rfact2 ) 
     109      CALL trc_sink( kt, wsbio4, sinking2, jpgoc, rfact2 ) 
     110      CALL trc_sink( kt, wsbio4, sinkfer2, jpbfe, rfact2 ) 
     111      CALL trc_sink( kt, wsbio4, sinksil , jpgsi, rfact2 ) 
     112      CALL trc_sink( kt, wsbio4, sinkcal , jpcal, rfact2 ) 
    167113 
    168114      IF( ln_p5z ) THEN 
    … …  
    174120         !   Compute the sedimentation term using p4zsink2 for all the sinking particles 
    175121         !   ----------------------------------------------------- 
    176          DO jit = 1, iiter1 
    177            CALL p4z_sink2( wsbio3, sinkingn , jppon, iiter1 ) 
    178            CALL p4z_sink2( wsbio3, sinkingp , jppop, iiter1 ) 
    179          END DO 
    180  
    181          DO jit = 1, iiter2 
    182            CALL p4z_sink2( wsbio4, sinking2n, jpgon, iiter2 ) 
    183            CALL p4z_sink2( wsbio4, sinking2p, jpgop, iiter2 ) 
    184          END DO 
     122         CALL trc_sink( kt, wsbio3, sinkingn , jppon, rfact2 ) 
     123         CALL trc_sink( kt, wsbio3, sinkingp , jppop, rfact2 ) 
     124         CALL trc_sink( kt, wsbio4, sinking2n, jpgon, rfact2 ) 
     125         CALL trc_sink( kt, wsbio4, sinking2p, jpgop, rfact2 ) 
    185126      ENDIF 
    186127 
    187128      IF( ln_ligand ) THEN 
    188129         wsfep (:,:,:) = wfep 
    189          DO jk = 1,jpkm1 
    190             DO jj = 1, jpj 
    191                DO ji = 1, jpi 
    192                   IF( tmask(ji,jj,jk) == 1 ) THEN 
    193                     zwsmax = 0.5 * e3t_n(ji,jj,jk) / xstep 
    194                     wsfep(ji,jj,jk) = MIN( wsfep(ji,jj,jk), zwsmax * REAL( iiter1, wp ) ) 
    195                   ENDIF 
    196                END DO 
    197             END DO 
    198          END DO 
    199130         ! 
    200131         sinkfep(:,:,:) = 0.e0 
    201          DO jit = 1, iiter1 
    202            CALL p4z_sink2( wsfep, sinkfep , jpfep, iiter1 ) 
    203          END DO 
     132         CALL trc_sink( kt, wsfep, sinkfep , jpfep, rfact2 ) 
    204133      ENDIF 
    205134 
    … …  
    281210   END SUBROUTINE p4z_sink_init 
    282211 
    283  
    284    SUBROUTINE p4z_sink2( pwsink, psinkflx, jp_tra, kiter ) 
    285       !!--------------------------------------------------------------------- 
    286       !!                     ***  ROUTINE p4z_sink2  *** 
    287       !! 
    288       !! ** Purpose :   Compute the sedimentation terms for the various sinking 
    289       !!     particles. The scheme used to compute the trends is based 
    290       !!     on MUSCL. 
    291       !! 
    292       !! ** Method  : - this ROUTINE compute not exactly the advection but the 
    293       !!      transport term, i.e.  div(u*tra). 
    294       !!--------------------------------------------------------------------- 
    295       INTEGER , INTENT(in   )                         ::   jp_tra    ! tracer index index       
    296       INTEGER , INTENT(in   )                         ::   kiter     ! number of iterations for time-splitting  
    297       REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj,jpk) ::   pwsink    ! sinking speed 
    298       REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   psinkflx  ! sinking fluxe 
    299       ! 
    300       INTEGER  ::   ji, jj, jk, jn 
    301       REAL(wp) ::   zigma,zew,zign, zflx, zstep 
    302       REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztraz, zakz, zwsink2, ztrb  
    303       !!--------------------------------------------------------------------- 
    304       ! 
    305       IF( ln_timing )   CALL timing_start('p4z_sink2') 
    306       ! 
    307       zstep = rfact2 / REAL( kiter, wp ) / 2. 
    308  
    309       ztraz(:,:,:) = 0.e0 
    310       zakz (:,:,:) = 0.e0 
    311       ztrb (:,:,:) = trb(:,:,:,jp_tra) 
    312  
    313       DO jk = 1, jpkm1 
    314          zwsink2(:,:,jk+1) = -pwsink(:,:,jk) / rday * tmask(:,:,jk+1)  
    315       END DO 
    316       zwsink2(:,:,1) = 0.e0 
    317  
    318  
    319       ! Vertical advective flux 
    320       DO jn = 1, 2 
    321          !  first guess of the slopes interior values 
    322          DO jk = 2, jpkm1 
    323             ztraz(:,:,jk) = ( trb(:,:,jk-1,jp_tra) - trb(:,:,jk,jp_tra) ) * tmask(:,:,jk) 
    324          END DO 
    325          ztraz(:,:,1  ) = 0.0 
    326          ztraz(:,:,jpk) = 0.0 
    327  
    328          ! slopes 
    329          DO jk = 2, jpkm1 
    330             DO jj = 1,jpj 
    331                DO ji = 1, jpi 
    332                   zign = 0.25 + SIGN( 0.25, ztraz(ji,jj,jk) * ztraz(ji,jj,jk+1) ) 
    333                   zakz(ji,jj,jk) = ( ztraz(ji,jj,jk) + ztraz(ji,jj,jk+1) ) * zign 
    334                END DO 
    335             END DO 
    336          END DO 
    337           
    338          ! Slopes limitation 
    339          DO jk = 2, jpkm1 
    340             DO jj = 1, jpj 
    341                DO ji = 1, jpi 
    342                   zakz(ji,jj,jk) = SIGN( 1., zakz(ji,jj,jk) ) *        & 
    343                      &             MIN( ABS( zakz(ji,jj,jk) ), 2. * ABS(ztraz(ji,jj,jk+1)), 2. * ABS(ztraz(ji,jj,jk) ) ) 
    344                END DO 
    345             END DO 
    346          END DO 
    347           
    348          ! vertical advective flux 
    349          DO jk = 1, jpkm1 
    350             DO jj = 1, jpj       
    351                DO ji = 1, jpi     
    352                   zigma = zwsink2(ji,jj,jk+1) * zstep / e3w_n(ji,jj,jk+1) 
    353                   zew   = zwsink2(ji,jj,jk+1) 
    354                   psinkflx(ji,jj,jk+1) = -zew * ( trb(ji,jj,jk,jp_tra) - 0.5 * ( 1 + zigma ) * zakz(ji,jj,jk) ) * zstep 
    355                END DO 
    356             END DO 
    357          END DO 
    358          ! 
    359          ! Boundary conditions 
    360          psinkflx(:,:,1  ) = 0.e0 
    361          psinkflx(:,:,jpk) = 0.e0 
    362           
    363          DO jk=1,jpkm1 
    364             DO jj = 1,jpj 
    365                DO ji = 1, jpi 
    366                   zflx = ( psinkflx(ji,jj,jk) - psinkflx(ji,jj,jk+1) ) / e3t_n(ji,jj,jk) 
    367                   trb(ji,jj,jk,jp_tra) = trb(ji,jj,jk,jp_tra) + zflx 
    368                END DO 
    369             END DO 
    370          END DO 
    371  
    372       ENDDO 
    373  
    374       DO jk = 1,jpkm1 
    375          DO jj = 1,jpj 
    376             DO ji = 1, jpi 
    377                zflx = ( psinkflx(ji,jj,jk) - psinkflx(ji,jj,jk+1) ) / e3t_n(ji,jj,jk) 
    378                ztrb(ji,jj,jk) = ztrb(ji,jj,jk) + 2. * zflx 
    379             END DO 
    380          END DO 
    381       END DO 
    382  
    383       trb(:,:,:,jp_tra) = ztrb(:,:,:) 
    384       psinkflx(:,:,:)   = 2. * psinkflx(:,:,:) 
    385       ! 
    386       IF( ln_timing )  CALL timing_stop('p4z_sink2') 
    387       ! 
    388    END SUBROUTINE p4z_sink2 
    389  
    390  
    391212   INTEGER FUNCTION p4z_sink_alloc() 
    392213      !!---------------------------------------------------------------------- 

  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p4zsms.F90

    r10345 r10377  
    7474            CALL p4z_che                              ! initialize the chemical constants 
    7575            CALL ahini_for_at(hi)   !  set PH at kt=nit000 
    76             t_oce_co2_flx_cum = 0._wp 
    7776        ELSE 
    7877            CALL p4z_rst( nittrc000, 'READ' )  !* read or initialize all required fields 
    … …  
    189188      !! 
    190189      NAMELIST/nampisbio/ nrdttrc, wsbio, xkmort, ferat3, wsbio2, wsbio2max, wsbio2scale,    & 
    191          &                   niter1max, niter2max, wfep, ldocp, ldocz, lthet,  & 
     190         &                   wfep, ldocp, ldocz, lthet,  & 
    192191         &                   no3rat3, po4rat3 
    193192         ! 
    … …  
    223222         WRITE(numout,*) '      Big particles maximum sinking speed       wsbio2max   =', wsbio2max 
    224223         WRITE(numout,*) '      Big particles sinking speed length scale  wsbio2scale =', wsbio2scale 
    225          WRITE(numout,*) '      Maximum number of iterations for POC      niter1max   =', niter1max 
    226          WRITE(numout,*) '      Maximum number of iterations for GOC      niter2max   =', niter2max 
    227224         IF( ln_ligand ) THEN 
    228225            WRITE(numout,*) '      FeP sinking speed                              wfep   =', wfep 

  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/sms_pisces.F90

    r10368 r10377  
    3636 
    3737   !!*  Biological parameters  
    38    INTEGER  ::   niter1max, niter2max   !: Maximum number of iterations for sinking 
    3938   REAL(wp) ::   rno3              !: ??? 
    4039   REAL(wp) ::   o2ut              !: ??? 

  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/trcini.F90

    r10372 r10377  
    205205      USE trcldf , ONLY:  trc_ldf_ini 
    206206      USE trcrad , ONLY:  trc_rad_ini 
     207      USE trcsink, ONLY:  trc_sink_ini 
    207208      ! 
    208209      INTEGER :: ierr 
    … …  
    214215                       !                          ! vertical diffusion: always implicit time stepping scheme 
    215216                       CALL  trc_rad_ini          ! positivity of passive tracers  
     217                       CALL  trc_sink_ini         ! Vertical sedimentation of particles 
    216218      ! 
    217219   END SUBROUTINE trc_ini_trp 
Note: See TracChangeset for help on using the changeset viewer.