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 3527 – NEMO

Changeset 3527


Ignore:
Timestamp:
2012-11-05T16:23:05+01:00 (11 years ago)
Author:
poddo
Message:

include the up-stream in muscl.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2012/dev_r3521_INGV7_muscl/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_muscl.F90

    r3294 r3527  
    88   !!   NEMO     1.0  !  2002-06  (G. Madec)  F90: Free form and module 
    99   !!            3.2  !  2010-05  (C. Ethe, G. Madec)  merge TRC-TRA + switch from velocity to transport 
     10   !!            3.4  !  2012-06  (P. Oddo) include the upstream where needed 
    1011   !!---------------------------------------------------------------------- 
    1112 
     
    1819   USE trdmod_oce      ! tracers trends  
    1920   USE trdtra      ! tracers trends  
     21   USE eosbn2          ! equation of state  
    2022   USE in_out_manager  ! I/O manager 
    2123   USE dynspg_oce      ! choice/control of key cpp for surface pressure gradient 
    2224   USE trabbl          ! tracers: bottom boundary layer 
     25   USE sbcrnf          ! river runoffs 
    2326   USE lib_mpp         ! distribued memory computing 
    2427   USE lbclnk          ! ocean lateral boundary condition (or mpp link)  
     
    2730   USE wrk_nemo        ! Memory Allocation 
    2831   USE timing          ! Timing 
     32   USE eosbn2          ! equation of state 
     33   USE sbcrnf          ! river runoffs 
    2934 
    3035   IMPLICIT NONE 
     
    3540   LOGICAL  :: l_trd       ! flag to compute trends 
    3641 
     42   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: upsmsk !: mixed upstream/centered scheme near some straits 
     43   !                                                             !  and in closed seas (orca 2 and 4 configurations) 
    3744   !! * Substitutions 
    3845#  include "domzgr_substitute.h90" 
     
    7481      ! 
    7582      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices 
     83      INTEGER  ::   ierr 
    7684      REAL(wp) ::   zu, z0u, zzwx, zw         ! local scalars 
    7785      REAL(wp) ::   zv, z0v, zzwy, z0w        !   -      - 
    7886      REAL(wp) ::   ztra, zbtr, zdt, zalpha   !   -      - 
    7987      REAL(wp), POINTER, DIMENSION(:,:,:) :: zslpx, zslpy 
     88      INTEGER  ::   ierr 
     89      REAL(wp) ::   zice                             ! temporary scalars 
     90      REAL(wp), POINTER, DIMENSION(:,:  ) :: ztfreez 
     91      REAL(wp), POINTER, DIMENSION(:,:,:) :: zind 
    8092      !!---------------------------------------------------------------------- 
    8193      ! 
     
    8395      ! 
    8496      CALL wrk_alloc( jpi, jpj, jpk, zslpx, zslpy ) 
     97      CALL wrk_alloc( jpi, jpj, ztfreez ) 
     98      CALL wrk_alloc( jpi, jpj, jpk, zind ) 
    8599      ! 
    86100 
     
    89103         IF(lwp) WRITE(numout,*) 'tra_adv : MUSCL advection scheme on ', cdtype 
    90104         IF(lwp) WRITE(numout,*) '~~~~~~~' 
     105         IF(lwp) WRITE(numout,*) 
     106         ! 
     107         ! 
     108         IF (.not. ALLOCATED(upsmsk))THEN 
     109             ALLOCATE( upsmsk(jpi,jpj), STAT=ierr ) 
     110             IF( ierr /= 0 )   CALL ctl_stop('STOP', 'tra_adv_muscl: unable to allocate array') 
     111         ENDIF 
     112         ! 
     113         upsmsk(:,:) = 0._wp                             ! not upstream by default 
    91114         ! 
    92115         l_trd = .FALSE. 
     
    94117      ENDIF 
    95118 
     119      ! 
     120      ! Upstream / centered scheme indicator 
     121      ! ------------------------------------ 
     122      ztfreez(:,:) = tfreez( tsn(:,:,1,jp_sal) ) 
     123      DO jk = 1, jpk 
     124         DO jj = 1, jpj 
     125            DO ji = 1, jpi 
     126               !                                        ! below ice covered area (if tn < "freezing"+0.1 ) 
     127               IF( tsn(ji,jj,jk,jp_tem) <= ztfreez(ji,jj) + 0.1_wp ) THEN   ;   zice = 1.e0 
     128               ELSE                                                         ;   zice = 0.e0 
     129               ENDIF 
     130               zind(ji,jj,jk) = MAX (   & 
     131                  rnfmsk(ji,jj) * rnfmsk_z(jk),      &  ! near runoff mouths (& closed sea outflows) 
     132                  upsmsk(ji,jj)               ,      &  ! some of some straits 
     133                  zice                               &  ! below ice covered area (if tn < "freezing"+0.1 ) 
     134                  &                  ) * tmask(ji,jj,jk) 
     135               zind(ji,jj,jk) = 1 - zind(ji,jj,jk) 
     136            END DO 
     137         END DO 
     138      END DO 
    96139      !                                                     ! =========== 
    97140      DO jn = 1, kjpt                                       ! tracer loop 
     
    148191                  zalpha = 0.5 - z0u 
    149192                  zu  = z0u - 0.5 * pun(ji,jj,jk) * zdt / ( e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) ) 
    150                   zzwx = ptb(ji+1,jj,jk,jn) + zu * zslpx(ji+1,jj,jk) 
    151                   zzwy = ptb(ji  ,jj,jk,jn) + zu * zslpx(ji  ,jj,jk) 
     193                  zzwx = ptb(ji+1,jj,jk,jn) + zind(ji,jj,jk) * (zu * zslpx(ji+1,jj,jk)) 
     194                  zzwy = ptb(ji  ,jj,jk,jn) + zind(ji,jj,jk) * (zu * zslpx(ji  ,jj,jk)) 
    152195                  zwx(ji,jj,jk) = pun(ji,jj,jk) * ( zalpha * zzwx + (1.-zalpha) * zzwy ) 
    153196                  ! 
     
    155198                  zalpha = 0.5 - z0v 
    156199                  zv  = z0v - 0.5 * pvn(ji,jj,jk) * zdt / ( e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) ) 
    157                   zzwx = ptb(ji,jj+1,jk,jn) + zv * zslpy(ji,jj+1,jk) 
    158                   zzwy = ptb(ji,jj  ,jk,jn) + zv * zslpy(ji,jj  ,jk)  
     200                  zzwx = ptb(ji,jj+1,jk,jn) + zind(ji,jj,jk) * (zv * zslpy(ji,jj+1,jk)) 
     201                  zzwy = ptb(ji,jj  ,jk,jn) + zind(ji,jj,jk) * (zv * zslpy(ji,jj  ,jk)) 
    159202                  zwy(ji,jj,jk) = pvn(ji,jj,jk) * ( zalpha * zzwx + (1.-zalpha) * zzwy ) 
    160203               END DO 
     
    230273                  zalpha = 0.5 + z0w 
    231274                  zw  = z0w - 0.5 * pwn(ji,jj,jk+1) * zdt * zbtr  
    232                   zzwx = ptb(ji,jj,jk+1,jn) + zw * zslpx(ji,jj,jk+1) 
    233                   zzwy = ptb(ji,jj,jk  ,jn) + zw * zslpx(ji,jj,jk  ) 
     275                  zzwx = ptb(ji,jj,jk+1,jn) + zind(ji,jj,jk) * (zw * zslpx(ji,jj,jk+1)) 
     276                  zzwy = ptb(ji,jj,jk  ,jn) + zind(ji,jj,jk) * (zw * zslpx(ji,jj,jk  )) 
    234277                  zwx(ji,jj,jk+1) = pwn(ji,jj,jk+1) * ( zalpha * zzwx + (1.-zalpha) * zzwy ) 
    235278               END DO  
     
    255298      ! 
    256299      CALL wrk_dealloc( jpi, jpj, jpk, zslpx, zslpy ) 
     300      CALL wrk_dealloc( jpi, jpj, ztfreez ) 
     301      CALL wrk_dealloc( jpi, jpj, jpk, zind ) 
    257302      ! 
    258303      IF( nn_timing == 1 )  CALL timing_stop('tra_adv_muscl') 
Note: See TracChangeset for help on using the changeset viewer.