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 5682 for branches/2015/dev_r5072_UKMO2_OBS_simplification/NEMOGCM/TOOLS/NESTING/src/agrif_partial_steps.f90 – NEMO

Ignore:
Timestamp:
2015-08-12T17:46:45+02:00 (9 years ago)
Author:
mattmartin
Message:

OBS simplification changes committed to branch after running SETTE tests to make sure we get the same results as the trunk for ORCA2_LIM_OBS.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2015/dev_r5072_UKMO2_OBS_simplification/NEMOGCM/TOOLS/NESTING/src/agrif_partial_steps.f90

    r2143 r5682  
    3636    !        
    3737    TYPE(Coordinates) :: Grid                      
    38     REAL*8 :: za1,za0,zsur,zacr,zkth,zdepth,zdepwp,zmin,zmax,zdiff,ze3tp,ze3wp 
     38    REAL*8 :: za2,za1,za0,zsur,zacr,zkth,zacr2,zkth2,zdepth,zdepwp,zmin,zmax,zdiff,ze3tp,ze3wp 
    3939    INTEGER :: i,j,jk,jj,ji,jpj,jpi,ik,ii,ipt,jpt 
    4040    INTEGER, DIMENSION(1) :: k 
     
    7676       za0  = pa0 
    7777       za1  = pa1 
     78       za2  = pa2 
    7879       ! 
    7980    ELSE 
     
    8889    ENDIF 
    8990 
    90     zacr = ppacr 
    91     zkth = ppkth        
    92     ! 
    93     DO i = 1,N 
    94        !  
    95        gdepw(i) = (zsur+za0*i+za1*zacr*LOG(COSH((i-zkth)/zacr)))  
    96        gdept(i) = (zsur+za0*(i+0.5)+za1*zacr*LOG(COSH(((i+0.5)-zkth)/zacr))) 
    97        e3w(i)   = (za0 + za1 * TANH((i-zkth)/zacr)) 
    98        e3t(i)   = (za0 + za1 * TANH(((i+0.5)-zkth)/zacr)) 
    99        ! 
    100     END DO 
    101     ! 
    102  
     91    zacr  = ppacr 
     92    zkth  = ppkth      
     93    zacr2 = ppacr2 
     94    zkth2 = ppkth2   
     95    ! 
     96    IF( ppkth == 0. ) THEN            !  uniform vertical grid  
     97         za1 = pphmax / FLOAT(N-1)  
     98         DO i = 1, N 
     99            gdepw(i) = ( i - 1   ) * za1 
     100            gdept(i) = ( i - 0.5 ) * za1 
     101            e3w  (i) =  za1 
     102            e3t  (i) =  za1 
     103         END DO 
     104    ELSE                            ! Madec & Imbard 1996 function 
     105       IF( .NOT. ldbletanh ) THEN 
     106          DO i = 1,N 
     107             !  
     108             gdepw(i) = (zsur+za0*i+za1*zacr*LOG(COSH((i-zkth)/zacr)))  
     109             gdept(i) = (zsur+za0*(i+0.5)+za1*zacr*LOG(COSH(((i+0.5)-zkth)/zacr))) 
     110             e3w(i)   = (za0 + za1 * TANH((i-zkth)/zacr)) 
     111             e3t(i)   = (za0 + za1 * TANH(((i+0.5)-zkth)/zacr)) 
     112             ! 
     113          END DO 
     114       ELSE 
     115            DO i = 1,N 
     116               ! Double tanh function 
     117               gdepw(i) = ( zsur + za0*i  + za1 * zacr * LOG ( COSH( (i-zkth ) / zacr  ) )               & 
     118                  &                       + za2 * zacr2* LOG ( COSH( (i-zkth2) / zacr2 ) )  ) 
     119               gdept(i) = ( zsur + za0*(i+0.5) + za1 * zacr * LOG ( COSH( ((i+0.5)-zkth ) / zacr  ) )    & 
     120                  &                            + za2 * zacr2* LOG ( COSH( ((i+0.5)-zkth2) / zacr2 ) )  ) 
     121               e3w  (i) =          za0         + za1        * TANH(       (i-zkth ) / zacr  )            & 
     122                  &                            + za2        * TANH(       (i-zkth2) / zacr2 ) 
     123               e3t  (i) =          za0         + za1        * TANH(       ((i+0.5)-zkth ) / zacr  )      & 
     124                  &                            + za2        * TANH(       ((i+0.5)-zkth2) / zacr2 ) 
     125            END DO 
     126       ENDIF 
     127    ENDIF 
    103128    gdepw(1) = 0.0   
    104129    ! 
     
    106131    ! 
    107132    zmax = gdepw(N) + e3t(N) 
    108     zmin = gdepw(4) 
     133    IF( rn_hmin < 0. ) THEN  ;   i = - INT( rn_hmin )                                  ! from a nb of level 
     134    ELSE                     ;   i = MINLOC( gdepw, mask = gdepw > rn_hmin, dim = 1 )  ! from a depth 
     135    ENDIF 
     136    zmin = gdepw(i+1) 
    109137    ! 
    110138    ! Initialize bathy_level to the maximum ocean level available 
     
    235263    INTEGER :: i,j,ji,ij,ii,jj,jpt,ipt 
    236264    REAL,DIMENSION(N) :: gdepw,e3t 
    237     REAL :: za0,za1,zsur,zacr,zkth,zmin,zmax,zdepth 
     265    REAL :: za0,za1,za2,zsur,zacr,zacr2,zkth,zkth2,zmin,zmax,zdepth 
    238266    INTEGER :: kbathy,jk,diff 
    239267    INTEGER :: bornex,borney,bornex2,borney2 
    240     !        
     268    ! 
    241269    IF ( ( pa0 == 0 .OR. pa1 == 0 .OR. psur == 0 ) & 
    242270         .AND. ppdzmin.NE.0 .AND. pphmax.NE.0 ) THEN  
    243271       !     
     272       WRITE(*,*) 'psur,pa0,pa1 computed' 
    244273       za1=( ppdzmin - pphmax / (N-1) )          & 
    245274            / ( TANH((1-ppkth)/ppacr) - ppacr/(N-1) & 
     
    253282         pa0.NE.0 .AND. pa1.NE.0 ) THEN 
    254283       !        
     284       WRITE(*,*) 'psur,pa0,pa1 given by namelist' 
    255285       zsur = psur 
    256286       za0  = pa0 
    257287       za1  = pa1 
     288       za2  = pa2 
    258289       ! 
    259290    ELSE 
     
    263294       WRITE(*,*) 'please check values of variables' 
    264295       WRITE(*,*) 'in namelist vertical_grid section' 
    265        WRITE(*,*) ' '       
    266        !        
    267     ENDIF 
    268     !        
    269     zacr = ppacr 
    270     zkth = ppkth        
    271     ! 
    272     DO i = 1,N 
    273        !  
    274        gdepw(i) = (zsur+za0*i+za1*zacr*LOG(COSH((i-zkth)/zacr))) 
    275        e3t(i)   = (za0 + za1 * TANH(((i+0.5)-zkth)/zacr))        
    276     END DO 
    277     ! 
     296       WRITE(*,*) ' '   
     297       STOP     
     298       !        
     299    ENDIF 
     300 
     301    zacr  = ppacr 
     302    zkth  = ppkth      
     303    zacr2 = ppacr2 
     304    zkth2 = ppkth2   
     305    ! 
     306    IF( ppkth == 0. ) THEN            !  uniform vertical grid  
     307         za1 = pphmax / FLOAT(N-1)  
     308         DO i = 1, N 
     309            gdepw(i) = ( i - 1   ) * za1 
     310            e3t  (i) =  za1 
     311         END DO 
     312    ELSE                            ! Madec & Imbard 1996 function 
     313       IF( .NOT. ldbletanh ) THEN 
     314          DO i = 1,N 
     315             !  
     316             gdepw(i) = (zsur+za0*i+za1*zacr*LOG(COSH((i-zkth)/zacr)))  
     317             e3t(i)   = (za0 + za1 * TANH(((i+0.5)-zkth)/zacr)) 
     318             ! 
     319          END DO 
     320       ELSE 
     321            DO i = 1,N 
     322               ! Double tanh function 
     323               gdepw(i) = ( zsur + za0*i  + za1 * zacr * LOG ( COSH( (i-zkth ) / zacr  ) )               & 
     324                  &                       + za2 * zacr2* LOG ( COSH( (i-zkth2) / zacr2 ) )  ) 
     325               e3t  (i) =          za0         + za1        * TANH(       ((i+0.5)-zkth ) / zacr  )      & 
     326                  &                            + za2        * TANH(       ((i+0.5)-zkth2) / zacr2 ) 
     327            END DO 
     328       ENDIF 
     329    ENDIF 
    278330    gdepw(1) = 0.0 
    279     ! 
    280331    ! 
    281332    diff = 0       
     
    344395    ! 
    345396    zmax = gdepw(N) + e3t(N) 
    346     zmin = gdepw(4)  
     397    IF( rn_hmin < 0. ) THEN  ;   i = - INT( rn_hmin )                                  ! from a nb of level 
     398    ELSE                     ;   i = MINLOC( gdepw, mask = gdepw > rn_hmin, dim = 1 )  ! from a depth 
     399    ENDIF 
     400    zmin = gdepw(i+1) 
    347401    ! 
    348402    ! check that interpolated value stays at the same level          
     
    646700    REAL*8, DIMENSION(:,:,:) :: fse3u,fse3t,fse3v 
    647701    !                                   
    648     REAL*8 :: za1,za0,zsur,zacr,zkth,zdepth,zdepwp,zmin,zmax,zdiff,ze3tp,ze3wp 
     702    REAL*8 :: za2,za1,za0,zsur,zacr,zkth,zacr2,zkth2,zdepth,zdepwp,zmin,zmax,zdiff,ze3tp,ze3wp 
    649703    INTEGER :: i,j,jk,jj,ji,jpj,jpi,ik,ii,ipt,jpt,jpk 
    650704    INTEGER, DIMENSION(1) :: k 
     
    660714    ALLOCATE(gdepw(jpk),e3t(jpk)) 
    661715    ALLOCATE(gdepw_ps(jpi,jpj,jpk))                   
    662     !        
     716    !   
    663717    IF ( ( pa0 == 0 .OR. pa1 == 0 .OR. psur == 0 ) & 
    664718         .AND. ppdzmin.NE.0 .AND. pphmax.NE.0 ) THEN  
     
    668722            / ( TANH((1-ppkth)/ppacr) - ppacr/(jpk-1) & 
    669723            *  (  LOG( COSH( (jpk - ppkth) / ppacr) )      & 
    670             - LOG( COSH( ( 1  - ppkth) / ppacr) )  )  )  
    671        ! 
     724            - LOG( COSH( ( 1  - ppkth) / ppacr) )  )  ) 
     725 
    672726       za0  = ppdzmin - za1 * TANH( (1-ppkth) / ppacr ) 
    673727       zsur = - za0 - za1 * ppacr * LOG( COSH( (1-ppkth) / ppacr )  ) 
     
    676730         pa0.NE.0 .AND. pa1.NE.0 ) THEN 
    677731       !        
     732       WRITE(*,*) 'psur,pa0,pa1 given by namelist' 
    678733       zsur = psur 
    679734       za0  = pa0 
    680        za1  = pa1    
    681        !        
    682     ENDIF 
    683  
    684     zacr = ppacr 
    685     zkth = ppkth        
     735       za1  = pa1 
     736       za2  = pa2 
     737       ! 
     738    ELSE 
     739       !        
     740       WRITE(*,*) 'ERROR ***** bad vertical grid parameters ...'  
     741       WRITE(*,*) ' ' 
     742       WRITE(*,*) 'please check values of variables' 
     743       WRITE(*,*) 'in namelist vertical_grid section' 
     744       WRITE(*,*) ' '   
     745       STOP     
     746       !        
     747    ENDIF 
     748 
     749    zacr  = ppacr 
     750    zkth  = ppkth      
     751    zacr2 = ppacr2 
     752    zkth2 = ppkth2   
     753    ! 
     754    IF( ppkth == 0. ) THEN            !  uniform vertical grid  
     755         za1 = pphmax / FLOAT(jpk-1)  
     756         DO i = 1, jpk 
     757            gdepw(i) = ( i - 1   ) * za1 
     758            e3t  (i) =  za1 
     759         END DO 
     760    ELSE                            ! Madec & Imbard 1996 function 
     761       IF( .NOT. ldbletanh ) THEN 
     762          DO i = 1,jpk 
     763             !  
     764             gdepw(i) = (zsur+za0*i+za1*zacr*LOG(COSH((i-zkth)/zacr)))  
     765             e3t(i)   = (za0 + za1 * TANH(((i+0.5)-zkth)/zacr)) 
     766             ! 
     767          END DO 
     768       ELSE 
     769            DO i = 1,jpk 
     770               ! Double tanh function 
     771               gdepw(i) = ( zsur + za0*i  + za1 * zacr * LOG ( COSH( (i-zkth ) / zacr  ) )               & 
     772                  &                       + za2 * zacr2* LOG ( COSH( (i-zkth2) / zacr2 ) )  ) 
     773               e3t  (i) =          za0         + za1        * TANH(       ((i+0.5)-zkth ) / zacr  )      & 
     774                  &                            + za2        * TANH(       ((i+0.5)-zkth2) / zacr2 ) 
     775            END DO 
     776       ENDIF 
     777    ENDIF          
    686778    !          
    687779    !                 
    688780    DO i = 1,jpk 
    689        !  
    690        gdepw(i) = (zsur+za0*i+za1*zacr*LOG(COSH((i-zkth)/zacr)))  
    691        e3t(i)   = (za0 + za1 * TANH(((i+0.5)-zkth)/zacr)) 
    692781       !        
    693782       fse3t(:,:,i) = e3t(i) 
     
    700789    ! 
    701790    zmax = gdepw(jpk) + e3t(jpk) 
    702     zmin = gdepw(4) 
     791    IF( rn_hmin < 0. ) THEN  ;   i = - INT( rn_hmin )                                  ! from a nb of level 
     792    ELSE                     ;   i = MINLOC( gdepw, mask = gdepw > rn_hmin, dim = 1 )  ! from a depth 
     793    ENDIF 
     794    zmin = gdepw(i+1) 
    703795    ! 
    704796    DO jj = 1, jpj 
Note: See TracChangeset for help on using the changeset viewer.