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

Changeset 811


Ignore:
Timestamp:
2008-02-07T17:00:12+01:00 (16 years ago)
Author:
ctlod
Message:

dev_001_SBC: merge with the trunk last changesets: #780, 782, 783, 784, 785, 788, 789, 793, 794

Location:
branches/dev_001_SBC/NEMO
Files:
9 deleted
51 edited
1 copied

Legend:

Unmodified
Added
Removed

  • branches/dev_001_SBC/NEMO/LIM_SRC/limrst.F90

    r717 r811  
    5454      !!---------------------------------------------------------------------- 
    5555      ! 
    56       IF( kt == nit000 )   lrst_ice = .FALSE. 
    57        
    58       IF( kt == nitrst - 2*nn_fsbc + 1 .OR.  nitend - nit000 + 1 <= nn_fsbc ) THEN 
    59          ! beware if model runs less than nn_fsbc + 1 time step 
    60          ! beware of the format used to write kt (default is i8.8, that should be large enough) 
    61          IF( nitrst > 1.0e9 ) THEN    
    62             WRITE(clkt,*) nitrst 
    63          ELSE 
    64             WRITE(clkt,'(i8.8)') nitrst 
     56      IF( kt == nit000 )   lrst_ice = .FALSE.   ! default definition 
     57       
     58      ! to get better performances with NetCDF format: 
     59      ! we open and define the ice restart file one ice time step before writing the data (-> at nitrst - 2*nn_fsbc + 1) 
     60      ! except if we write ice restart files every ice time step or if an ice restart file was writen at nitend - 2*nn_fsbc + 1 
     61      IF( kt == nitrst - 2*nn_fsbc + 1 .OR. nstock == nn_fsbc .OR. ( kt == nitend - nn_fsbc + 1 .AND. .NOT. lrst_ice ) ) THEN 
     62         ! beware of the format used to write kt (default is i8.8, that should be large enough...) 
     63         IF( nitrst > 99999999 ) THEN   ;   WRITE(clkt, *       ) nitrst 
     64         ELSE                           ;   WRITE(clkt, '(i8.8)') nitrst 
    6565         ENDIF 
    6666         ! create the file 
    67          IF(lwp) WRITE(numout,*) 
    6867         clname = TRIM(cexper)//"_"//TRIM(ADJUSTL(clkt))//"_restart_ice" 
    69          IF(lwp) WRITE(numout,*) '             open ice restart.output NetCDF file: '//clname 
     68         IF(lwp) THEN 
     69            WRITE(numout,*) 
     70            SELECT CASE ( jprstlib ) 
     71            CASE ( jprstdimg )   ;   WRITE(numout,*) '             open ice restart binary file: '//clname 
     72            CASE DEFAULT         ;   WRITE(numout,*) '             open ice restart NetCDF file: '//clname 
     73            END SELECT 
     74            IF( kt == nitrst - 2*nn_fsbc + 1 ) THEN    
     75               WRITE(numout,*)         '             kt = nitrst - 2*nn_fsbc + 1 = ', kt,' date= ', ndastp 
     76            ELSE   ;   WRITE(numout,*) '             kt = '                         , kt,' date= ', ndastp 
     77            ENDIF 
     78         ENDIF 
     79 
    7080         CALL iom_open( clname, numriw, ldwrt = .TRUE., kiolib = jprstlib ) 
    7181         lrst_ice = .TRUE. 
    … …  
    8696      !!---------------------------------------------------------------------- 
    8797 
    88       iter = kt + nn_fsbc - 1 
     98      iter = kt + nn_fsbc - 1   ! ice restarts are written at kt == nitrst - nn_fsbc + 1 
    8999 
    90100      IF( iter == nitrst ) THEN 
    91101         IF(lwp) WRITE(numout,*) 
    92          IF(lwp) WRITE(numout,*) 'lim_rst_write : write ice restart.output NetCDF file  kt =', kt 
     102         IF(lwp) WRITE(numout,*) 'lim_rst_write : write ice restart file  kt =', kt 
    93103         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~' 
    94104      ENDIF 

  • branches/dev_001_SBC/NEMO/NST_SRC/agrif_opa_interp.F90

    r699 r811  
    88   USE dom_oce       
    99   USE sol_oce 
     10   USE agrif_oce 
    1011 
    1112   IMPLICIT NONE 
    … …  
    1617   CONTAINS 
    1718    
    18    SUBROUTINE Agrif_tra( kt ) 
     19   SUBROUTINE Agrif_tra 
    1920      !!--------------------------------------------- 
    2021      !!   *** ROUTINE Agrif_Tra *** 
    … …  
    2324#  include "vectopt_loop_substitute.h90" 
    2425       
    25       INTEGER, INTENT(in) :: kt 
    26  
    2726      INTEGER :: ji,jj,jk 
    2827      REAL(wp) :: zrhox 
    … …  
    178177 
    179178      Agrif_SpecialValue=0. 
    180       Agrif_UseSpecialValue = .TRUE. 
     179      Agrif_UseSpecialValue = ln_spc_dyn 
     180 
    181181      zua = 0. 
    182182      zva = 0. 
    … …  
    187187 
    188188      Agrif_SpecialValue=0. 
    189       Agrif_UseSpecialValue = .TRUE. 
     189      Agrif_UseSpecialValue = ln_spc_dyn 
    190190      CALL Agrif_Bc_variable(zua2d,e1u,calledweight=1.,procname=interpu2d) 
    191191      CALL Agrif_Bc_variable(zva2d,e2v,calledweight=1.,procname=interpv2d) 

  • branches/dev_001_SBC/NEMO/NST_SRC/agrif_opa_sponge.F90

    r699 r811  
    1010   USE dom_oce 
    1111   USE in_out_manager 
     12   USE agrif_oce 
    1213 
    1314   IMPLICIT NONE 
    … …  
    1617   PUBLIC Agrif_Sponge_Tra, Agrif_Sponge_Dyn, interptn, interpsn, interpun, interpvn 
    1718 
    18    !! * Namelist (namagrif) 
    19    REAL(wp) :: visc_tra = rdt 
    20    REAL(wp) :: visc_dyn = rdt 
    2119 
    2220   CONTAINS 
    2321 
    24    SUBROUTINE Agrif_Sponge_Tra( kt ) 
     22   SUBROUTINE Agrif_Sponge_Tra 
    2523      !!--------------------------------------------- 
    2624      !!   *** ROUTINE Agrif_Sponge_Tra *** 
    … …  
    2826#include "domzgr_substitute.h90" 
    2927 
    30       INTEGER, INTENT(in) :: kt 
    31  
    3228      INTEGER :: ji,jj,jk 
    33       REAL(wp), DIMENSION(jpi,jpj,jpk) :: umasktemp,vmasktemp 
    3429      INTEGER :: spongearea 
    3530      REAL(wp) :: timecoeff 
    3631      REAL(wp) :: zta, zsa, zabe1, zabe2, zbtr 
    3732      REAL(wp), DIMENSION(jpi,jpj) :: localviscsponge 
    38       REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztab, tbdiff, sbdiff 
     33      REAL(wp), DIMENSION(jpi,jpj,jpk) :: tbdiff, sbdiff 
    3934      REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztu ,ztv, zsu ,zsv 
     35      REAL(wp), DIMENSION(jpi,jpj,jpk) ::  ztab 
    4036 
    4137#if defined SPONGE 
    42  
    43       IF( kt == nit000  )   CALL agrif_sponge_init 
    4438 
    4539      timecoeff = REAL(Agrif_NbStepint(),wp)/Agrif_rhot() 
    … …  
    6155      sbdiff(:,:,:) = sb(:,:,:) - ztab(:,:,:) 
    6256 
     57 
    6358      spongearea = 2 + 2 * Agrif_irhox() 
    6459 
    6560      localviscsponge = 0. 
    66       umasktemp = 0. 
    67       vmasktemp = 0. 
     61       
     62      IF (.NOT. spongedoneT) THEN 
     63         spe1ur(:,:) = 0. 
     64         spe2vr(:,:) = 0. 
    6865 
    6966      IF ((nbondi == -1).OR.(nbondi == 2)) THEN 
    … …  
    7168            localviscsponge(ji,:) = visc_tra * (spongearea-ji)/real(spongearea-2) 
    7269         ENDDO 
    73          DO jk = 1, jpkm1 
    74             umasktemp(2:spongearea-1,:,jk) = umask(2:spongearea-1,:,jk) & 
    75                * 0.5 * (localviscsponge(2:spongearea-1,:) + localviscsponge(3:spongearea,:)) 
    76          ENDDO 
    77          DO jk = 1, jpkm1 
    78             vmasktemp(2:spongearea,1:jpjm1,jk) = vmask(2:spongearea,1:jpjm1,jk) & 
    79                * 0.5 * (localviscsponge(2:spongearea,1:jpjm1) + localviscsponge(2:spongearea,2:jpj)) 
    80          ENDDO 
     70     
     71    spe1ur(2:spongearea-1,:)=0.5 * (localviscsponge(2:spongearea-1,:) + localviscsponge(3:spongearea,:)) & 
     72          * e2u(2:spongearea-1,:) / e1u(2:spongearea-1,:) 
     73 
     74         spe2vr(2:spongearea,1:jpjm1) = 0.5 * (localviscsponge(2:spongearea,1:jpjm1) + & 
     75             localviscsponge(2:spongearea,2:jpj)) & 
     76           * e1v(2:spongearea,1:jpjm1) / e2v(2:spongearea,1:jpjm1) 
    8177      ENDIF 
    8278 
    … …  
    8581            localviscsponge(ji,:) = visc_tra * (ji - (nlci-spongearea+1))/real(spongearea-2) 
    8682         ENDDO 
    87          DO jk = 1, jpkm1 
    88             umasktemp(nlci-spongearea + 1:nlci-2,:,jk) = umask(nlci-spongearea + 1:nlci-2,:,jk) & 
    89                * 0.5 * (localviscsponge(nlci-spongearea + 1:nlci-2,:) + localviscsponge(nlci-spongearea + 2:nlci-1,:)) 
    90          ENDDO 
    91          DO jk = 1, jpkm1 
    92             vmasktemp(nlci-spongearea + 1:nlci-1,1:jpjm1,jk) = vmask(nlci-spongearea + 1:nlci-1,1:jpjm1,jk) & 
    93                * 0.5 * (localviscsponge(nlci-spongearea + 1:nlci-1,1:jpjm1) + localviscsponge(nlci-spongearea + 1:nlci-1,2:jpj)) 
    94          ENDDO 
     83     
     84    spe1ur(nlci-spongearea + 1:nlci-2,:)=0.5 * (localviscsponge(nlci-spongearea + 1:nlci-2,:) + & 
     85           localviscsponge(nlci-spongearea + 2:nlci-1,:)) & 
     86          * e2u(nlci-spongearea + 1:nlci-2,:) / e1u(nlci-spongearea + 1:nlci-2,:) 
     87 
     88         spe2vr(nlci-spongearea + 1:nlci-1,1:jpjm1) = 0.5 * (localviscsponge(nlci-spongearea + 1:nlci-1,1:jpjm1) & 
     89              + localviscsponge(nlci-spongearea + 1:nlci-1,2:jpj)) & 
     90           * e1v(nlci-spongearea + 1:nlci-1,1:jpjm1) / e2v(nlci-spongearea + 1:nlci-1,1:jpjm1) 
    9591      ENDIF 
    9692 
    … …  
    10096            localviscsponge(:,jj) = visc_tra * (spongearea-jj)/real(spongearea-2) 
    10197         ENDDO 
    102          DO jk = 1, jpkm1 
    103             vmasktemp(:,2:spongearea-1,jk) = vmask(:,2:spongearea-1,jk) & 
    104                * 0.5 * (localviscsponge(:,2:spongearea-1) + localviscsponge(:,3:spongearea)) 
    105          ENDDO 
    106          DO jk = 1, jpkm1 
    107             umasktemp(1:jpim1,2:spongearea,jk) = umask(1:jpim1,2:spongearea,jk) & 
    108                * 0.5 * (localviscsponge(1:jpim1,2:spongearea) + localviscsponge(2:jpi,2:spongearea)) 
    109          ENDDO 
     98     
     99    spe1ur(1:jpim1,2:spongearea)=0.5 * (localviscsponge(1:jpim1,2:spongearea) + & 
     100           localviscsponge(2:jpi,2:spongearea)) & 
     101          * e2u(1:jpim1,2:spongearea) / e1u(1:jpim1,2:spongearea) 
     102 
     103         spe2vr(:,2:spongearea-1) = 0.5 * (localviscsponge(:,2:spongearea-1) + & 
     104             localviscsponge(:,3:spongearea)) & 
     105           * e1v(:,2:spongearea-1) / e2v(:,2:spongearea-1) 
    110106      ENDIF 
    111107 
    … …  
    114110            localviscsponge(:,jj) = visc_tra * (jj - (nlcj-spongearea+1))/real(spongearea-2) 
    115111         ENDDO 
    116          DO jk = 1, jpkm1 
    117             vmasktemp(:,nlcj-spongearea + 1:nlcj-2,jk) = vmask(:,nlcj-spongearea + 1:nlcj-2,jk) & 
    118                * 0.5 * (localviscsponge(:,nlcj-spongearea + 1:nlcj-2) + localviscsponge(:,nlcj-spongearea + 2:nlcj-1)) 
    119          ENDDO 
    120          DO jk = 1, jpkm1 
    121             umasktemp(1:jpim1,nlcj-spongearea + 1:nlcj-1,jk) = umask(1:jpim1,nlcj-spongearea + 1:nlcj-1,jk) & 
    122                * 0.5 * (localviscsponge(1:jpim1,nlcj-spongearea + 1:nlcj-1) + localviscsponge(2:jpi,nlcj-spongearea + 1:nlcj-1)) 
    123          ENDDO 
    124       ENDIF 
    125  
    126       IF (.NOT. spongedoneT) THEN 
    127          zspe1ur(:,:) = e2u(:,:) / e1u(:,:) 
    128          zspe2vr(:,:) = e1v(:,:) / e2v(:,:) 
    129          zspbtr2(:,:) = 1. / ( e1t(:,:) * e2t(:,:)) 
     112     
     113    spe1ur(1:jpim1,nlcj-spongearea + 1:nlcj-1)=0.5 * (localviscsponge(1:jpim1,nlcj-spongearea + 1:nlcj-1) + & 
     114            localviscsponge(2:jpi,nlcj-spongearea + 1:nlcj-1)) & 
     115          * e2u(1:jpim1,nlcj-spongearea + 1:nlcj-1) / e1u(1:jpim1,nlcj-spongearea + 1:nlcj-1) 
     116 
     117         spe2vr(:,nlcj-spongearea + 1:nlcj-2) = 0.5 * (localviscsponge(:,nlcj-spongearea + 1:nlcj-2) + & 
     118            localviscsponge(:,nlcj-spongearea + 2:nlcj-1)) & 
     119           * e1v(:,nlcj-spongearea + 1:nlcj-2) / e2v(:,nlcj-spongearea + 1:nlcj-2) 
     120      ENDIF 
     121       
     122         spbtr2(:,:) = 1. / ( e1t(:,:) * e2t(:,:)) 
    130123 
    131124         spongedoneT = .TRUE. 
    … …  
    136129            DO ji = 1, jpim1 
    137130#if defined key_zco 
    138                zabe1 = umasktemp(ji,jj,jk) * zspe1ur(ji,jj) 
    139                zabe2 = vmasktemp(ji,jj,jk) * zspe2vr(ji,jj) 
    140 #else 
    141                zabe1 = umasktemp(ji,jj,jk) * zspe1ur(ji,jj) * fse3u(ji,jj,jk) 
    142                zabe2 = vmasktemp(ji,jj,jk) * zspe2vr(ji,jj) * fse3v(ji,jj,jk) 
     131               zabe1 = umask(ji,jj,jk) * spe1ur(ji,jj) 
     132               zabe2 = vmask(ji,jj,jk) * spe2vr(ji,jj) 
     133#else 
     134               zabe1 = umask(ji,jj,jk) * spe1ur(ji,jj) * fse3u(ji,jj,jk) 
     135               zabe2 = vmask(ji,jj,jk) * spe2vr(ji,jj) * fse3v(ji,jj,jk) 
    143136#endif 
    144137               ztu(ji,jj,jk) = zabe1 * ( tbdiff(ji+1,jj  ,jk) - tbdiff(ji,jj,jk) ) 
    … …  
    152145            DO ji = 2,jpim1 
    153146#if defined key_zco 
    154                zbtr = zspbtr2(ji,jj) 
    155 #else 
    156                zbtr = zspbtr2(ji,jj) / fse3t(ji,jj,jk) 
     147               zbtr = spbtr2(ji,jj) 
     148#else 
     149               zbtr = spbtr2(ji,jj) / fse3t(ji,jj,jk) 
    157150#endif 
    158151               ! horizontal diffusive trends 
    … …  
    173166   END SUBROUTINE Agrif_Sponge_Tra 
    174167 
    175    SUBROUTINE Agrif_Sponge_dyn( kt ) 
     168   SUBROUTINE Agrif_Sponge_dyn 
    176169      !!--------------------------------------------- 
    177170      !!   *** ROUTINE Agrif_Sponge_dyn *** 
    178171      !!--------------------------------------------- 
    179172#include "domzgr_substitute.h90" 
    180  
    181       INTEGER,INTENT(in) :: kt 
    182173 
    183174      INTEGER :: ji,jj,jk 
    184175      INTEGER :: spongearea 
    185176      REAL(wp) :: timecoeff 
    186       REAL(wp) :: ze2u, ze1v, zua, zva 
     177      REAL(wp) :: ze2u, ze1v, zua, zva, zbtr 
    187178      REAL(wp), DIMENSION(jpi,jpj) :: localviscsponge 
    188179      REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztab, ubdiff, vbdiff,rotdiff,hdivdiff 
    189       REAL(wp), DIMENSION(jpi,jpj,jpk) :: umasktemp,vmasktemp 
    190180 
    191181#if defined SPONGE 
    … …  
    194184 
    195185      Agrif_SpecialValue=0. 
    196       Agrif_UseSpecialValue = .TRUE. 
     186      Agrif_UseSpecialValue = ln_spc_dyn 
    197187      ztab = 0.e0 
    198188      CALL Agrif_Bc_Variable(ztab, ua,calledweight=timecoeff,procname=interpun) 
    199189      Agrif_UseSpecialValue = .FALSE. 
    200190 
    201       ubdiff(:,:,:) = ub(:,:,:) - ztab(:,:,:) 
     191      ubdiff(:,:,:) = (ub(:,:,:) - ztab(:,:,:))*umask(:,:,:) 
    202192 
    203193      ztab = 0.e0 
    204194      Agrif_SpecialValue=0. 
    205       Agrif_UseSpecialValue = .TRUE. 
     195      Agrif_UseSpecialValue = ln_spc_dyn 
    206196      CALL Agrif_Bc_Variable(ztab, va,calledweight=timecoeff,procname=interpvn) 
    207197      Agrif_UseSpecialValue = .FALSE. 
    208198 
    209       vbdiff(:,:,:) = vb(:,:,:) - ztab(:,:,:) 
     199      vbdiff(:,:,:) = (vb(:,:,:) - ztab(:,:,:))*vmask(:,:,:) 
    210200 
    211201      spongearea = 2 + 2 * Agrif_irhox() 
    212202 
    213203      localviscsponge = 0. 
    214       umasktemp = 0. 
    215       vmasktemp = 0. 
     204 
     205      IF (.NOT. spongedoneU) THEN 
     206         spe1ur2(:,:) = 0. 
     207         spe2vr2(:,:) = 0. 
    216208 
    217209      IF ((nbondi == -1).OR.(nbondi == 2)) THEN 
    … …  
    219211            localviscsponge(ji,:) = visc_dyn * (spongearea-ji)/real(spongearea-2) 
    220212         ENDDO 
    221          DO jk = 1, jpkm1 
    222             umasktemp(2:spongearea-1,:,jk) = umask(2:spongearea-1,:,jk) & 
    223                * 0.5 * (localviscsponge(2:spongearea-1,:) + localviscsponge(3:spongearea,:)) 
    224          ENDDO 
    225          DO jk = 1, jpkm1 
    226             vmasktemp(2:spongearea,1:jpjm1,jk) = vmask(2:spongearea,1:jpjm1,jk) & 
    227                * 0.5 * (localviscsponge(2:spongearea,1:jpjm1) + localviscsponge(2:spongearea,2:jpj)) 
    228          ENDDO 
     213     
     214    spe1ur2(2:spongearea-1,:)=0.5 * (localviscsponge(2:spongearea-1,:) + localviscsponge(3:spongearea,:)) 
     215 
     216         spe2vr2(2:spongearea,1:jpjm1) = 0.5 * (localviscsponge(2:spongearea,1:jpjm1) + & 
     217             localviscsponge(2:spongearea,2:jpj)) 
    229218      ENDIF 
    230219 
    … …  
    233222            localviscsponge(ji,:) = visc_dyn * (ji - (nlci-spongearea+1))/real(spongearea-2) 
    234223         ENDDO 
    235          DO jk = 1, jpkm1 
    236             umasktemp(nlci-spongearea + 1:nlci-2,:,jk) = umask(nlci-spongearea + 1:nlci-2,:,jk) & 
    237                * 0.5 * (localviscsponge(nlci-spongearea + 1:nlci-2,:) + localviscsponge(nlci-spongearea + 2:nlci-1,:)) 
    238          ENDDO 
    239          DO jk = 1, jpkm1 
    240             vmasktemp(nlci-spongearea + 1:nlci-1,1:jpjm1,jk) = vmask(nlci-spongearea + 1:nlci-1,1:jpjm1,jk) & 
    241                * 0.5 * (localviscsponge(nlci-spongearea + 1:nlci-1,1:jpjm1) + localviscsponge(nlci-spongearea + 1:nlci-1,2:jpj)) 
    242          ENDDO 
    243       ENDIF 
     224     
     225    spe1ur2(nlci-spongearea + 1:nlci-2,:)=0.5 * (localviscsponge(nlci-spongearea + 1:nlci-2,:) + & 
     226           localviscsponge(nlci-spongearea + 2:nlci-1,:)) 
     227 
     228         spe2vr2(nlci-spongearea + 1:nlci-1,1:jpjm1) = 0.5 * (localviscsponge(nlci-spongearea + 1:nlci-1,1:jpjm1) & 
     229              + localviscsponge(nlci-spongearea + 1:nlci-1,2:jpj)) 
     230      ENDIF 
     231 
    244232 
    245233      IF ((nbondj == -1).OR.(nbondj == 2)) THEN 
    … …  
    247235            localviscsponge(:,jj) = visc_dyn * (spongearea-jj)/real(spongearea-2) 
    248236         ENDDO 
    249          DO jk = 1, jpkm1 
    250             vmasktemp(:,2:spongearea-1,jk) = vmask(:,2:spongearea-1,jk) & 
    251                * 0.5 * (localviscsponge(:,2:spongearea-1) + localviscsponge(:,3:spongearea)) 
    252          ENDDO 
    253          DO jk = 1, jpkm1 
    254             umasktemp(1:jpim1,2:spongearea,jk) = umask(1:jpim1,2:spongearea,jk) & 
    255                * 0.5 * (localviscsponge(1:jpim1,2:spongearea) + localviscsponge(2:jpi,2:spongearea)) 
    256          ENDDO 
     237     
     238    spe1ur2(1:jpim1,2:spongearea)=0.5 * (localviscsponge(1:jpim1,2:spongearea) + & 
     239           localviscsponge(2:jpi,2:spongearea)) 
     240 
     241         spe2vr2(:,2:spongearea-1) = 0.5 * (localviscsponge(:,2:spongearea-1) + & 
     242             localviscsponge(:,3:spongearea)) 
    257243      ENDIF 
    258244 
    … …  
    261247            localviscsponge(:,jj) = visc_dyn * (jj - (nlcj-spongearea+1))/real(spongearea-2) 
    262248         ENDDO 
    263          DO jk = 1, jpkm1 
    264             vmasktemp(:,nlcj-spongearea + 1:nlcj-2,jk) = vmask(:,nlcj-spongearea + 1:nlcj-2,jk) & 
    265                * 0.5 * (localviscsponge(:,nlcj-spongearea + 1:nlcj-2) + localviscsponge(:,nlcj-spongearea + 2:nlcj-1)) 
    266          ENDDO 
    267          DO jk = 1, jpkm1 
    268             umasktemp(1:jpim1,nlcj-spongearea + 1:nlcj-1,jk) = umask(1:jpim1,nlcj-spongearea + 1:nlcj-1,jk) & 
    269                * 0.5 * (localviscsponge(1:jpim1,nlcj-spongearea + 1:nlcj-1) + localviscsponge(2:jpi,nlcj-spongearea + 1:nlcj-1)) 
    270          ENDDO 
    271       ENDIF 
    272  
    273       ubdiff = ubdiff * umasktemp 
    274       vbdiff = vbdiff * vmasktemp 
    275  
     249     
     250    spe1ur2(1:jpim1,nlcj-spongearea + 1:nlcj-1)=0.5 * (localviscsponge(1:jpim1,nlcj-spongearea + 1:nlcj-1) + & 
     251            localviscsponge(2:jpi,nlcj-spongearea + 1:nlcj-1)) 
     252 
     253         spe2vr2(:,nlcj-spongearea + 1:nlcj-2) = 0.5 * (localviscsponge(:,nlcj-spongearea + 1:nlcj-2) + & 
     254            localviscsponge(:,nlcj-spongearea + 2:nlcj-1)) 
     255      ENDIF 
     256 
     257         spongedoneU = .TRUE. 
     258     
     259    spbtr3(:,:) = 1./( e1f(:,:) * e2f(:,:)) 
     260      ENDIF 
     261       
     262      IF (.NOT. spongedoneT) THEN 
     263        spbtr2(:,:) = 1. / ( e1t(:,:) * e2t(:,:))       
     264      ENDIF 
     265       
     266      DO jk=1,jpkm1 
     267      ubdiff(:,:,jk) = ubdiff(:,:,jk) * spe1ur2(:,:) 
     268      vbdiff(:,:,jk) = vbdiff(:,:,jk) * spe2vr2(:,:) 
     269      ENDDO 
     270       
    276271      hdivdiff = 0. 
    277272      rotdiff = 0. 
    … …  
    286281            DO ji = 2, jpim1   ! vector opt. 
    287282#if defined key_zco 
     283               zbtr = spbtr2(ji,jj) 
    288284               hdivdiff(ji,jj,jk) = (  e2u(ji,jj) * ubdiff(ji,jj,jk) & 
    289285                  - e2u(ji-1,jj  ) * ubdiff(ji-1,jj  ,jk)      & 
    290286                  &               + e1v(ji,jj) * vbdiff(ji,jj,jk) - & 
    291                   &              e1v(ji  ,jj-1) * vbdiff(ji  ,jj-1,jk)  )   & 
    292                   &            / ( e1t(ji,jj) * e2t(ji,jj) ) 
    293 #else 
     287                  &              e1v(ji  ,jj-1) * vbdiff(ji  ,jj-1,jk)  ) * zbtr 
     288#else 
     289               zbtr = spbtr2(ji,jj) / fse3t(ji,jj,jk) 
    294290               hdivdiff(ji,jj,jk) =   & 
    295291                  (  e2u(ji,jj)*fse3u(ji,jj,jk) * &  
    … …  
    298294                  + e1v(ji,jj)*fse3v(ji,jj,jk) * & 
    299295                  vbdiff(ji,jj,jk) - e1v(ji  ,jj-1)* & 
    300                   fse3v(ji  ,jj-1,jk)  * vbdiff(ji  ,jj-1,jk)  )    & 
    301                   / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 
     296                  fse3v(ji  ,jj-1,jk)  * vbdiff(ji  ,jj-1,jk)  ) * zbtr 
    302297#endif 
    303298            END DO 
    … …  
    306301         DO jj = 1, jpjm1 
    307302            DO ji = 1, jpim1   ! vector opt. 
     303#if defined key_zco       
     304               zbtr = spbtr3(ji,jj) 
    308305               rotdiff(ji,jj,jk) = (  e2v(ji+1,jj  ) * vbdiff(ji+1,jj  ,jk) - e2v(ji,jj) * vbdiff(ji,jj,jk)    & 
    309306                  &              - e1u(ji  ,jj+1) * ubdiff(ji  ,jj+1,jk) + e1u(ji,jj) * ubdiff(ji,jj,jk)  ) & 
    310                   &           * fmask(ji,jj,jk) / ( e1f(ji,jj) * e2f(ji,jj) ) 
     307                  &           * fmask(ji,jj,jk) * zbtr 
     308#else 
     309               zbtr = spbtr3(ji,jj) * fse3f(ji,jj,jk) 
     310               rotdiff(ji,jj,jk) = (  e2v(ji+1,jj  ) * vbdiff(ji+1,jj  ,jk) - e2v(ji,jj) * vbdiff(ji,jj,jk)    & 
     311                  &              - e1u(ji  ,jj+1) * ubdiff(ji  ,jj+1,jk) + e1u(ji,jj) * ubdiff(ji,jj,jk)  ) & 
     312                  &           * fmask(ji,jj,jk) * zbtr 
     313#endif         
    311314            END DO 
    312315         END DO 
    … …  
    333336                  ze1v                  ) / e2v(ji,jj) 
    334337#else 
    335                ze2u = rotdiff (ji,jj,jk)*fse3f(ji,jj,jk) 
     338               ze2u = rotdiff (ji,jj,jk) 
    336339               ze1v = hdivdiff(ji,jj,jk) 
    337340               ! horizontal diffusive trends 
    338                zua = - ( ze2u - rotdiff (ji,jj-1,jk)* & 
    339                   fse3f(ji,jj-1,jk) ) / ( e2u(ji,jj) * fse3u(ji,jj,jk) )   & 
     341               zua = - ( ze2u - rotdiff (ji,jj-1,jk)) / ( e2u(ji,jj) * fse3u(ji,jj,jk) )   & 
    340342                  + ( hdivdiff(ji+1,jj,jk) - ze1v      & 
    341343                  ) / e1u(ji,jj) 
    342344 
    343                zva = + ( ze2u - rotdiff (ji-1,jj,jk)* & 
    344                   fse3f(ji-1,jj,jk) ) / ( e1v(ji,jj) * fse3v(ji,jj,jk) )   & 
     345               zva = + ( ze2u - rotdiff (ji-1,jj,jk)) / ( e1v(ji,jj) * fse3v(ji,jj,jk) )   & 
    345346                  + ( hdivdiff(ji,jj+1,jk) - ze1v    & 
    346347                  ) / e2v(ji,jj) 
    … …  
    360361   END SUBROUTINE Agrif_Sponge_dyn 
    361362 
    362    SUBROUTINE agrif_sponge_init 
    363       !!--------------------------------------------- 
    364       !!   *** ROUTINE agrif_sponge_init *** 
    365       !!--------------------------------------------- 
    366       NAMELIST/namagrif/ visc_tra, visc_dyn 
    367       REWIND ( numnam ) 
    368       READ   ( numnam, namagrif ) 
    369  
    370       IF(lwp) THEN 
    371          WRITE(numout,*) 
    372          WRITE(numout,*) 'agrif_sponge_init : agrif sponge parameters' 
    373          WRITE(numout,*) '~~~~~~~~~~~~' 
    374          WRITE(numout,*) '          Namelist namagrif : set sponge parameters' 
    375          WRITE(numout,*) '             sponge coefficient for tracers =  ', visc_tra 
    376          WRITE(numout,*) '             sponge coefficient for dynamics = ', visc_dyn 
    377       ENDIF 
    378  
    379    END SUBROUTINE agrif_sponge_init 
    380  
    381363   SUBROUTINE interptn(tabres,i1,i2,j1,j2,k1,k2) 
    382364      !!--------------------------------------------- 
    … …  
    405387   END SUBROUTINE interpsn 
    406388 
    407  
    408389   SUBROUTINE interpun(tabres,i1,i2,j1,j2,k1,k2) 
    409390      !!--------------------------------------------- 

  • branches/dev_001_SBC/NEMO/NST_SRC/agrif_opa_update.F90

    r699 r811  
    99   USE oce 
    1010   USE dom_oce 
     11   USE agrif_oce 
    1112 
    1213   IMPLICIT NONE 
    … …  
    1516   PUBLIC Agrif_Update_Tra, Agrif_Update_Dyn 
    1617 
    17    INTEGER, PARAMETER :: nbclineupdate = 3 
    1818   INTEGER :: nbcline 
    1919 
    … …  
    7171      nbcline = nbcline + 1 
    7272 
    73       Agrif_UseSpecialValueInUpdate = .TRUE. 
     73      Agrif_UseSpecialValueInUpdate = ln_spc_dyn 
    7474      Agrif_SpecialValueFineGrid = 0. 
    7575      CALL Agrif_Update_Variable(ztab2d,sshn,procname = updateSSH) 

  • branches/dev_001_SBC/NEMO/NST_SRC/agrif_top_interp.F90

    r699 r811  
    88   USE dom_oce       
    99   USE sol_oce 
     10   USE agrif_oce 
    1011   USE trcstp 
    1112   USE sms 

  • branches/dev_001_SBC/NEMO/NST_SRC/agrif_top_update.F90

    r699 r811  
    1010   USE oce 
    1111   USE dom_oce 
     12   USE agrif_oce 
    1213   USE trcstp 
    1314   USE sms 
    … …  
    1819   PUBLIC Agrif_Update_Trc 
    1920 
    20    INTEGER, PARAMETER :: nbclineupdate = 3 
    2121   INTEGER :: nbcline 
    2222 

  • branches/dev_001_SBC/NEMO/NST_SRC/agrif_user.F90

    r699 r811  
    6868      USE ice_oce 
    6969#endif 
    70 #if defined key_agrif 
    7170      USE agrif_opa_update 
    7271      USE agrif_opa_interp 
    … …  
    7473      USE agrif_top_update 
    7574      USE agrif_top_interp 
    76 #endif 
    7775 
    7876      IMPLICIT NONE 
    … …  
    9290 
    9391      Call opa_init  ! Initializations of each fine grid 
     92      Call agrif_opa_init 
    9493 
    9594      ! Specific fine grid Initializations 
    … …  
    314313   End SUBROUTINE Agrif_detect 
    315314 
     315   SUBROUTINE agrif_opa_init 
     316      !!--------------------------------------------- 
     317      !!   *** ROUTINE agrif_init *** 
     318      !!--------------------------------------------- 
     319      USE agrif_oce  
     320      USE in_out_manager 
     321 
     322      IMPLICIT NONE 
     323 
     324      NAMELIST/namagrif/ nbclineupdate, visc_tra, visc_dyn, ln_spc_dyn 
     325 
     326      REWIND ( numnam ) 
     327      READ   ( numnam, namagrif ) 
     328      IF(lwp) THEN 
     329         WRITE(numout,*) 
     330         WRITE(numout,*) 'agrif_opa_init : agrif parameters' 
     331         WRITE(numout,*) '~~~~~~~~~~~~' 
     332         WRITE(numout,*) '          Namelist namagrif : set agrif parameters' 
     333         WRITE(numout,*) '             baroclinic update frequency          =  ', nbclineupdate 
     334         WRITE(numout,*) '             sponge coefficient for tracers       =  ', visc_tra 
     335         WRITE(numout,*) '             sponge coefficient for dynamics      =  ', visc_dyn 
     336         WRITE(numout,*) '             use special values for dynamics      =  ', ln_spc_dyn 
     337         WRITE(numout,*)  
     338      ENDIF 
     339 
     340    END SUBROUTINE agrif_opa_init 
    316341#if defined key_mpp_mpi 
    317  
    318342   SUBROUTINE Agrif_InvLoc(indloc,nprocloc,i,indglob) 
    319343      !!------------------------------------------ 
    … …  
    338362 
    339363   END SUBROUTINE Agrif_InvLoc 
    340  
    341 #endif 
    342  
     364#endif 
    343365#else 
    344366   SUBROUTINE Subcalledbyagrif 

  • branches/dev_001_SBC/NEMO/OPA_SRC/DOM/dom_oce.F90

    r699 r811  
    216216      !                        ! parameterize exchanges through straits 
    217217 
    218 #if defined key_agrif 
    219    !!---------------------------------------------------------------------- 
    220    !! agrif sponge layer 
    221    !!---------------------------------------------------------------------- 
    222    LOGICAL, PUBLIC :: spongedoneT = .FALSE.  !: ??? 
    223    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   & 
    224       zspe1ur, zspe2vr ,zspbtr2    !: ??? 
    225    !!---------------------------------------------------------------------- 
    226 #endif 
    227  
    228218END MODULE dom_oce 

  • branches/dev_001_SBC/NEMO/OPA_SRC/DOM/domain.F90

    r717 r811  
    179179      ndastp = ndate0                ! Assign initial date to current date 
    180180 
    181 ! ... Control of output frequency 
    182       IF ( nstock == 0 ) THEN 
    183           IF(lwp)WRITE(numout,cform_war) 
    184           IF(lwp)WRITE(numout,*) '           nstock = ', nstock, ' it is forced to ', nitend 
    185           nstock = nitend 
    186           nwarn = nwarn + 1 
     181      ! ... Control of output frequency 
     182      IF ( nstock == 0 .OR. nstock > nitend - nit000 + 1 ) THEN 
     183         WRITE(ctmp1,*) '           nstock = ', nstock, ' it is forced to ', nitend - nit000 + 1 
     184         CALL ctl_warn( ctmp1 ) 
     185         nstock = nitend - nit000 + 1 
    187186      ENDIF 
    188187      IF ( nwrite == 0 ) THEN 
    189           IF(lwp)WRITE(numout,cform_war) 
    190           IF(lwp)WRITE(numout,*) '           nwrite = ', nwrite, ' it is forced to ', nitend 
    191           nwrite = nitend 
    192           nwarn = nwarn + 1 
     188         WRITE(ctmp1,*) '           nwrite = ', nwrite, ' it is forced to ', nitend 
     189         CALL ctl_warn( ctmp1 ) 
     190         nwrite = nitend 
    193191      ENDIF 
    194192 

  • branches/dev_001_SBC/NEMO/OPA_SRC/DYN/divcur.F90

    r699 r811  
    122122 
    123123#if defined key_obc 
     124#if defined key_agrif 
     125         IF (Agrif_Root() ) THEN 
     126#endif 
    124127         ! open boundaries (div must be zero behind the open boundary) 
    125128         !  mpp remark: The zeroing of hdivn can probably be extended to 1->jpi/jpj for the correct row/column 
    … …  
    128131         IF( lp_obc_north )   hdivn(nin0  :nin1  ,njn0p1:njn1p1,jk) = 0.e0      ! north 
    129132         IF( lp_obc_south )   hdivn(nis0  :nis1  ,njs0  :njs1  ,jk) = 0.e0      ! south 
     133#if defined key_agrif 
     134         ENDIF 
     135#endif 
    130136#endif          
     137#if defined key_agrif 
     138         if ( .NOT. AGRIF_Root() ) then 
     139           IF ((nbondi ==  1).OR.(nbondi == 2)) hdivn(nlci-1 , :     ,jk) = 0.e0      ! east 
     140           IF ((nbondi == -1).OR.(nbondi == 2)) hdivn(2      , :     ,jk) = 0.e0      ! west 
     141           IF ((nbondj ==  1).OR.(nbondj == 2)) hdivn(:      ,nlcj-1 ,jk) = 0.e0      ! north 
     142           IF ((nbondj == -1).OR.(nbondj == 2)) hdivn(:      ,2      ,jk) = 0.e0      ! south 
     143         endif 
     144#endif     
    131145 
    132146         !                                             ! -------- 
    … …  
    317331 
    318332#if defined key_obc 
     333#if defined key_agrif 
     334         IF ( Agrif_Root() ) THEN 
     335#endif 
    319336         ! open boundaries (div must be zero behind the open boundary) 
    320337         !  mpp remark: The zeroing of hdivn can probably be extended to 1->jpi/jpj for the correct row/column 
    … …  
    323340         IF( lp_obc_north )   hdivn(nin0  :nin1  ,njn0p1:njn1p1,jk) = 0.e0      ! north 
    324341         IF( lp_obc_south )   hdivn(nis0  :nis1  ,njs0  :njs1  ,jk) = 0.e0      ! south 
     342#if defined key_agrif 
     343         ENDIF 
     344#endif 
    325345#endif          
     346#if defined key_agrif 
     347         if ( .NOT. AGRIF_Root() ) then 
     348            IF ((nbondi ==  1).OR.(nbondi == 2)) hdivn(nlci-1 , :     ,jk) = 0.e0      ! east 
     349            IF ((nbondi == -1).OR.(nbondi == 2)) hdivn(2      , :     ,jk) = 0.e0      ! west 
     350            IF ((nbondj ==  1).OR.(nbondj == 2)) hdivn(:      ,nlcj-1 ,jk) = 0.e0      ! north 
     351            IF ((nbondj == -1).OR.(nbondj == 2)) hdivn(:      ,2      ,jk) = 0.e0      ! south 
     352         endif 
     353#endif     
     354 
    326355         !                                             ! -------- 
    327356         ! relative vorticity                          !   rot  

  • branches/dev_001_SBC/NEMO/OPA_SRC/DYN/dynhpg.F90

    r699 r811  
    1818   !!   dyn_hpg      : update the momentum trend with the now horizontal 
    1919   !!                  gradient of the hydrostatic pressure 
    20    !!                  default case : k-j-i loops (vector opt. available) 
    2120   !!       hpg_ctl  : initialisation and control of options 
    2221   !!       hpg_zco  : z-coordinate scheme 
    … …  
    3029   USE oce             ! ocean dynamics and tracers 
    3130   USE dom_oce         ! ocean space and time domain 
    32    USE dynhpg_jki      ! 
    3331   USE phycst          ! physical constants 
    3432   USE in_out_manager  ! I/O manager 
    … …  
    4240 
    4341   PUBLIC   dyn_hpg    ! routine called by step module 
    44  
    45 #if defined key_mpp_omp 
    46    !!---------------------------------------------------------------------- 
    47    !!   'key_mpp_omp' :                                 j-k-i loop (j-slab) 
    48    !!---------------------------------------------------------------------- 
    49    LOGICAL, PUBLIC, PARAMETER ::   lk_dynhpg_jki = .TRUE.    !: OpenMP hpg flag 
    50    LOGICAL, PUBLIC, PARAMETER ::   lk_dynhpg     = .FALSE.   !: vector hpg flag 
    51 #else 
    52    !!---------------------------------------------------------------------- 
    53    !!   default case :                             k-j-i loop (vector opt.) 
    54    !!----------------------------------------------------------------------    
    55    LOGICAL, PUBLIC, PARAMETER ::   lk_dynhpg_jki = .FALSE.   !: OpenMP hpg flag 
    56    LOGICAL, PUBLIC, PARAMETER ::   lk_dynhpg     = .TRUE.    !: vector hpg flag 
    57 #endif 
    5842 
    5943   !!* Namelist nam_dynhpg : Choice of horizontal pressure gradient computation 
    … …  
    11195      CASE (  5 )   ;   CALL hpg_djc    ( kt )      ! s-coordinate (Density Jacobian with Cubic polynomial) 
    11296      CASE (  6 )   ;   CALL hpg_rot    ( kt )      ! s-coordinate (ROTated axes scheme) 
    113       CASE ( 10 )   ;   CALL hpg_zco_jki( kt )      ! z-coordinate (k-j-i) 
    114       CASE ( 11 )   ;   CALL hpg_zps_jki( kt )      ! z-coordinate plus partial steps (interpolation) (k-j-i) 
    115       CASE ( 12 )   ;   CALL hpg_sco_jki( kt )      ! s-coordinate (standard jacobian formulation) (k-j-i) 
    11697      END SELECT 
    11798 
    … …  
    186167      IF ( ioptio /= 1 )   CALL ctl_stop( ' NO or several hydrostatic pressure gradient options used' ) 
    187168 
    188       IF( lk_dynhpg_jki ) THEN 
    189          nhpg = nhpg + 10 
    190          IF(lwp) WRITE(numout,*) 
    191          IF(lwp) WRITE(numout,*) '          Autotasking or OPENMP: use j-k-i loops (i.e. _jki routines)' 
    192       ENDIF 
    193169      ! 
    194170   END SUBROUTINE hpg_ctl 

  • branches/dev_001_SBC/NEMO/OPA_SRC/DYN/dynspg.F90

    r699 r811  
    2020   USE dynspg_flt     ! surface pressure gradient     (dyn_spg_flt routine) 
    2121   USE dynspg_rl      ! surface pressure gradient     (dyn_spg_rl  routine) 
    22    USE dynspg_exp_jki ! surface pressure gradient (dyn_spg_exp_jki routine) 
    23    USE dynspg_ts_jki  ! surface pressure gradient (dyn_spg_ts_jki  routine) 
    24    USE dynspg_flt_jki ! surface pressure gradient (dyn_spg_flt_jki routine) 
    2522   USE trdmod         ! ocean dynamics trends 
    2623   USE trdmod_oce     ! ocean variables trends 
    … …  
    6865 
    6966      SELECT CASE ( nspg )                       ! compute surf. pressure gradient trend and add it to the general trend 
    70       !                                                      ! k-j-i loops 
     67      !                                                      
    7168      CASE (  0 )   ;   CALL dyn_spg_exp    ( kt )              ! explicit 
    7269      CASE (  1 )   ;   CALL dyn_spg_ts     ( kt )              ! time-splitting 
    7370      CASE (  2 )   ;   CALL dyn_spg_flt    ( kt, kindic )      ! filtered 
    7471      CASE (  3 )   ;   CALL dyn_spg_rl     ( kt, kindic )      ! rigid lid 
    75       !                                                      ! j-k-i loops 
    76       CASE ( 10 )   ;   CALL dyn_spg_exp_jki( kt )              ! explicit with j-k-i loop 
    77       CASE ( 11 )   ;   CALL dyn_spg_ts_jki ( kt )              ! time-splitting with j-k-i loop 
    78       CASE ( 12 )   ;   CALL dyn_spg_flt_jki( kt, kindic )      ! filtered with j-k-i loop 
    79       ! 
     72      !                                                     
    8073      CASE ( -1 )                                       ! esopa: test all possibility with control print 
    8174         ;              CALL dyn_spg_exp    ( kt ) 
    … …  
    8780         ;              CALL dyn_spg_flt  ( kt, kindic ) 
    8881         ;              CALL prt_ctl( tab3d_1=ua, clinfo1=' spg2 - Ua: ', mask1=umask, & 
    89             &                         tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' ) 
    90          ;              CALL dyn_spg_exp_jki( kt ) 
    91          ;              CALL prt_ctl( tab3d_1=ua, clinfo1=' spg10- Ua: ', mask1=umask, & 
    92             &                         tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' ) 
    93          ;              CALL dyn_spg_ts_jki ( kt ) 
    94          ;              CALL prt_ctl( tab3d_1=ua, clinfo1=' spg12- Ua: ', mask1=umask, & 
    95             &                         tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' ) 
    96          ;              CALL dyn_spg_flt_jki( kt, kindic ) 
    97          ;              CALL prt_ctl( tab3d_1=ua, clinfo1=' spg13- Ua: ', mask1=umask, & 
    9882            &                         tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' ) 
    9983      END SELECT 
    … …  
    159143      IF( lk_dynspg_flt)   nspg =  2 
    160144      IF( lk_dynspg_rl )   nspg =  3 
    161       IF( lk_jki       )   nspg =  nspg + 10 
    162145      IF( nspg == 13   )   nspg =  3 
    163146 
    … …  
    171154         IF( nspg ==  2 )   WRITE(numout,*) '     filtered free surface' 
    172155         IF( nspg ==  3 )   WRITE(numout,*) '     rigid-lid' 
    173          IF( nspg == 10 )   WRITE(numout,*) '     explicit free surface with j-k-i loop' 
    174          IF( nspg == 11 )   WRITE(numout,*) '     time splitting free surface with j-k-i loop' 
    175          IF( nspg == 12 )   WRITE(numout,*) '     filtered free surface with j-k-i loop' 
    176156      ENDIF 
    177157 

  • branches/dev_001_SBC/NEMO/OPA_SRC/DYN/dynspg_exp.F90

    r762 r811  
    3333   !! * Accessibility 
    3434   PUBLIC dyn_spg_exp  ! routine called by step.F90 
    35    PUBLIC exp_rst      ! routine called j-k-i subroutine 
     35   PUBLIC exp_rst      ! routine called by istate.F90 
    3636 
    3737   !! * Substitutions 

  • branches/dev_001_SBC/NEMO/OPA_SRC/DYN/dynspg_flt.F90

    r762 r811  
    3434   USE solsor          ! Successive Over-relaxation solver 
    3535   USE solfet          ! FETI solver 
    36    USE solsor_e        ! Successive Over-relaxation solver with MPP optimization 
    3736   USE obcdyn          ! ocean open boundary condition (obc_dyn routines) 
    3837   USE obcvol          ! ocean open boundary condition (obc_vol routines) 
    … …  
    5150 
    5251   PUBLIC dyn_spg_flt  ! routine called by step.F90 
    53    PUBLIC flt_rst      ! routine called by j-k-i subroutine 
     52   PUBLIC flt_rst      ! routine called by istate.F90 
    5453 
    5554   !! * Substitutions 
    … …  
    315314      END DO 
    316315      ! applied the lateral boundary conditions 
    317       IF( nsolv == 4 )  CALL lbc_lnk_e( gcb, c_solver_pt, 1. )    
     316      IF( nsolv == 2 .AND. MAX( jpr2di, jpr2dj ) > 0 ) CALL lbc_lnk_e( gcb, c_solver_pt, 1. )    
    318317 
    319318#if defined key_agrif 
    … …  
    360359         ELSEIF( nsolv == 3 ) THEN     ! FETI solver 
    361360            CALL sol_fet( kindic ) 
    362          ELSEIF( nsolv == 4 ) THEN     ! successive-over-relaxation with extra outer halo 
    363             CALL sol_sor_e( kindic ) 
    364361         ELSE                          ! e r r o r in nsolv namelist parameter 
    365362            WRITE(ctmp1,*) ' ~~~~~~~~~~~                not = ', nsolv 
    366             CALL ctl_stop( ' dyn_spg_flt : e r r o r, nsolv = 1, 2 ,3 or 4', ctmp1 ) 
     363            CALL ctl_stop( ' dyn_spg_flt : e r r o r, nsolv = 1, 2 or 3', ctmp1 ) 
    367364         ENDIF 
    368365      ENDIF 

  • branches/dev_001_SBC/NEMO/OPA_SRC/DYN/dynspg_rl.F90

    r762 r811  
    3636   USE solsor          ! Successive Over-relaxation solver 
    3737   USE solfet          ! FETI solver 
    38    USE solsor_e        ! Successive Over-relaxation solver with MPP optimization 
    3938   USE solisl          ! ??? 
    4039   USE obc_oce         ! Lateral open boundary condition 
    … …  
    177176      END DO 
    178177      ! applied the lateral boundary conditions 
    179       IF( nsolv == 4)  CALL lbc_lnk_e( gcb, c_solver_pt, 1. )    
     178      IF( nsolv == 2 .AND. MAX( jpr2di, jpr2dj ) > 0 ) CALL lbc_lnk_e( gcb, c_solver_pt, 1. )    
    180179 
    181180      ! Relative precision (computation on one processor) 
    … …  
    204203         CASE( 3 )                     ! FETI solver 
    205204            CALL sol_fet( kindic ) 
    206          CASE( 4 )                     ! successive-over-relaxation with extra outer halo 
    207             CALL sol_sor_e( kindic ) 
    208205         CASE DEFAULT                  ! e r r o r in nsolv namelist parameter 
    209206            WRITE(ctmp1,*) ' ~~~~~~~~~~                not = ', nsolv 
    210             CALL ctl_stop( ' dyn_spg_rl : e r r o r, nsolv = 1, 2 ,3 or 4', ctmp1 ) 
     207            CALL ctl_stop( ' dyn_spg_rl : e r r o r, nsolv = 1, 2 or 3', ctmp1 ) 
    211208         END SELECT 
    212209      ENDIF 

  • branches/dev_001_SBC/NEMO/OPA_SRC/DYN/dynspg_ts.F90

    r762 r811  
    77   !!             " "  !  05-11  (V. Garnier, G. Madec)  optimization 
    88   !!             9.0  !  06-08  (S. Masson)  distributed restart using iom 
     9   !!             9.0  !  08-01  (R. Benshila)  change averaging method 
    910   !!--------------------------------------------------------------------- 
    1011#if defined key_dynspg_ts   ||   defined key_esopa 
    … …  
    4041 
    4142   PUBLIC dyn_spg_ts  ! routine called by step.F90 
    42    PUBLIC ts_rst      ! routine called by j-k-i subroutine 
     43   PUBLIC ts_rst      ! routine called by istate.F90 
    4344 
    4445   REAL(wp), DIMENSION(jpi,jpj) ::  ftnw, ftne,   &  ! triad of coriolis parameter 
    … …  
    9293      !! * Local declarations 
    9394      INTEGER  ::  ji, jj, jk, jit             ! dummy loop indices 
    94       INTEGER  ::  icycle                      ! temporary scalar 
     95      INTEGER  ::  icycle, ibaro               ! temporary scalar 
    9596      REAL(wp) ::                           & 
    9697         zraur, zcoef, z2dt_e, z2dt_b, zfac25,   &  ! temporary scalars 
    … …  
    295296      !---------------- 
    296297      ! Number of iteration of the barotropic loop 
    297       icycle = FLOOR( z2dt_b / rdtbt ) 
     298      ibaro = FLOOR( rdt / rdtbt ) 
     299      icycle = 3 /2 * ibaro  
    298300 
    299301      ! variables for the barotropic equations 
    … …  
    304306      zun_e  (:,:) = un_b  (:,:) 
    305307      zvn_e  (:,:) = vn_b  (:,:) 
    306       zssha_b(:,:) = sshn  (:,:)       ! time averaged variables over all sub-timesteps 
    307       zua_b  (:,:) = un_b  (:,:)    
    308       zva_b  (:,:) = vn_b  (:,:) 
     308      zssha_b(:,:) = 0.e0 
     309      zua_b  (:,:) = 0.e0 
     310      zva_b  (:,:) = 0.e0 
    309311      IF( lk_vvl ) THEN 
    310312         zsshun_e(:,:)    = sshu (:,:)     ! (barotropic) sea surface height (now) at u-point 
    … …  
    464466         ! temporal sum 
    465467         !------------- 
    466          zssha_b(:,:) = zssha_b(:,:) + ssha_e(:,:) 
    467          zua_b  (:,:) = zua_b  (:,:) + ua_e  (:,:) 
    468          zva_b  (:,:) = zva_b  (:,:) + va_e  (:,:)  
     468         IF( jit >= ibaro/2 ) THEN 
     469            zssha_b(:,:) = zssha_b(:,:) + ssha_e(:,:) 
     470            zua_b  (:,:) = zua_b  (:,:) + ua_e  (:,:) 
     471            zva_b  (:,:) = zva_b  (:,:) + va_e  (:,:)  
     472         ENDIF 
    469473 
    470474         ! Time filter and swap of dynamics arrays 
    … …  
    530534 
    531535      ! Time average of after barotropic variables 
    532       zcoef =  1.e0 / (  FLOAT( icycle +1 ) ) 
     536      zcoef =  1.e0 / ( ibaro + 1 ) 
    533537      zssha_b(:,:) = zcoef * zssha_b(:,:)  
    534538      zua_b  (:,:) = zcoef *  zua_b (:,:)  

  • branches/dev_001_SBC/NEMO/OPA_SRC/DYN/dynvor.F90

    r699 r811  
    226226 
    227227!CDIR PARALLEL DO PRIVATE( zwx, zwy, zwz ) 
    228 !$OMP PARALLEL DO PRIVATE( zwx, zwy, zwz ) 
    229228      !                                                ! =============== 
    230229      DO jk = 1, jpkm1                                 ! Horizontal slab 
    … …  
    333332 
    334333!CDIR PARALLEL DO PRIVATE( zwx, zwy, zwz, zww ) 
    335 !$OMP PARALLEL DO PRIVATE( zwx, zwy, zwz, zww ) 
    336334      !                                                ! =============== 
    337335      DO jk = 1, jpkm1                                 ! Horizontal slab 
    … …  
    444442 
    445443!CDIR PARALLEL DO PRIVATE( zwx, zwy, zwz ) 
    446 !$OMP PARALLEL DO PRIVATE( zwx, zwy, zwz ) 
    447444      !                                                ! =============== 
    448445      DO jk = 1, jpkm1                                 ! Horizontal slab 
    … …  
    567564       
    568565!CDIR PARALLEL DO PRIVATE( zwx, zwy, zwz, ztnw, ztne, ztsw, ztse ) 
    569 !$OMP PARALLEL DO PRIVATE( zwx, zwy, zwz, ztnw, ztne, ztsw, ztse ) 
    570566      !                                                ! =============== 
    571567      DO jk = 1, jpkm1                                 ! Horizontal slab 

  • branches/dev_001_SBC/NEMO/OPA_SRC/DYN/dynzad.F90

    r708 r811  
    3838CONTAINS 
    3939 
    40 #if defined key_mpp_omp 
    41    !!---------------------------------------------------------------------- 
    42    !!   'key_mpp_omp'        OpenMP / NEC autotasking: j-k-i loops (j-slab) 
    43    !!---------------------------------------------------------------------- 
    44  
    45    SUBROUTINE dyn_zad( kt ) 
    46       !!---------------------------------------------------------------------- 
    47       !!                  ***  ROUTINE dynzad  *** 
    48       !! 
    49       !! ** Purpose :   Compute the now vertical momentum advection trend and  
    50       !!      add it to the general trend of momentum equation. 
    51       !! 
    52       !! ** Method  :   Use j-slab (j-k-i loops) for OpenMP / NEC autotasking 
    53       !!      The now vertical advection of momentum is given by: 
    54       !!         w dz(u) = ua + 1/(e1u*e2u*e3u) mk+1[ mi(e1t*e2t*wn) dk(un) ] 
    55       !!         w dz(v) = va + 1/(e1v*e2v*e3v) mk+1[ mj(e1t*e2t*wn) dk(vn) ] 
    56       !!      Add this trend to the general trend (ua,va): 
    57       !!         (ua,va) = (ua,va) + w dz(u,v) 
    58       !! 
    59       !! ** Action  : - Update (ua,va) with the vert. momentum advection trends 
    60       !!              - Save the trends in (ztrdu,ztrdv) ('key_trddyn') 
    61       !!---------------------------------------------------------------------- 
    62       USE oce, ONLY:   zwuw => ta   ! use ta as 3D workspace 
    63       USE oce, ONLY:   zwvw => sa   ! use sa as 3D workspace 
    64       !! 
    65       INTEGER, INTENT(in) ::   kt   ! ocean time-step inedx 
    66       !! 
    67       INTEGER  ::   ji, jj, jk      ! dummy loop indices 
    68       REAL(wp) ::   zvn, zua, zva   ! temporary scalars 
    69       REAL(wp), DIMENSION(jpi)         ::   zww            ! 1D workspace 
    70       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   ztrdu, ztrdv   ! 3D workspace 
    71       !!---------------------------------------------------------------------- 
    72        
    73       IF( kt == nit000 ) THEN 
    74          IF(lwp) WRITE(numout,*) 
    75          IF(lwp) WRITE(numout,*) 'dyn_zad : arakawa advection scheme' 
    76          IF(lwp) WRITE(numout,*) '~~~~~~~   Auto-tasking case, j-slab, no vector opt.' 
    77       ENDIF 
    78  
    79       IF( l_trddyn )   THEN         ! Save ua and va trends 
    80          ztrdu(:,:,:) = ua(:,:,:)  
    81          ztrdv(:,:,:) = va(:,:,:)  
    82       ENDIF 
    83  
    84       !                                                ! =============== 
    85       DO jj = 2, jpjm1                                 !  Vertical slab 
    86          !                                             ! =============== 
    87          DO jk = 2, jpkm1         ! Vertical momentum advection at uw and vw-pts 
    88             DO ji = 2, jpi              ! vertical fluxes  
    89                zww(ji) = 0.25 * e1t(ji,jj) * e2t(ji,jj) * wn(ji,jj,jk) 
    90             END DO 
    91             DO ji = 2, jpim1            ! vertical momentum advection at w-point 
    92                zvn = 0.25 * e1t(ji,jj+1) * e2t(ji,jj+1) * wn(ji,jj+1,jk) 
    93                zwuw(ji,jj,jk) = ( zww(ji+1) + zww(ji) ) * ( un(ji,jj,jk-1)-un(ji,jj,jk) ) 
    94                zwvw(ji,jj,jk) = ( zvn       + zww(ji) ) * ( vn(ji,jj,jk-1)-vn(ji,jj,jk) ) 
    95             END DO   
    96          END DO    
    97          DO ji = 2, jpim1               ! Surface and bottom values set to zero 
    98             zwuw(ji,jj, 1 ) = 0.e0 
    99             zwvw(ji,jj, 1 ) = 0.e0 
    100             zwuw(ji,jj,jpk) = 0.e0 
    101             zwvw(ji,jj,jpk) = 0.e0 
    102          END DO   
    103          ! 
    104          DO jk = 1, jpkm1         ! Vertical momentum advection at u- and v-points 
    105             DO ji = 2, jpim1 
    106                !                        ! vertical momentum advective trends 
    107                zua = - ( zwuw(ji,jj,jk) + zwuw(ji,jj,jk+1) ) / ( e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) ) 
    108                zva = - ( zwvw(ji,jj,jk) + zwvw(ji,jj,jk+1) ) / ( e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) ) 
    109                !                        ! add the trends to the general momentum trends 
    110                ua(ji,jj,jk) = ua(ji,jj,jk) + zua 
    111                va(ji,jj,jk) = va(ji,jj,jk) + zva 
    112             END DO   
    113          END DO   
    114          !                                             ! =============== 
    115       END DO                                           !   End of slab 
    116       !                                                ! =============== 
    117       ! 
    118       IF( l_trddyn ) THEN         ! save the vertical advection trends for diagnostic 
    119          ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:) 
    120          ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:) 
    121          CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_zad, 'DYN', kt ) 
    122       ENDIF 
    123       !                           ! Control print 
    124       IF(ln_ctl)   CALL prt_ctl( tab3d_1=ua, clinfo1=' zad  - Ua: ', mask1=umask,   & 
    125          &                       tab3d_2=va, clinfo2=' Va: ', mask2=vmask, clinfo3='dyn' ) 
    126       ! 
    127    END SUBROUTINE dyn_zad 
    128  
    129 #else 
    130    !!---------------------------------------------------------------------- 
    131    !!   Default option                             k-j-i loop (vector opt.) 
    132    !!---------------------------------------------------------------------- 
    133  
    13440   SUBROUTINE dyn_zad ( kt ) 
    13541      !!---------------------------------------------------------------------- 
    … …  
    16268         IF(lwp)WRITE(numout,*) 
    16369         IF(lwp)WRITE(numout,*) 'dyn_zad : arakawa advection scheme' 
    164          IF(lwp)WRITE(numout,*) '~~~~~~~   vector optimization k-j-i loop' 
    16570      ENDIF 
    16671 
    … …  
    215120      ! 
    216121   END SUBROUTINE dyn_zad 
    217 #endif 
    218122 
    219123   !!====================================================================== 

  • branches/dev_001_SBC/NEMO/OPA_SRC/DYN/dynzdf.F90

    r699 r811  
    1717   USE dynzdf_exp      ! vertical diffusion: explicit (dyn_zdf_exp     routine) 
    1818   USE dynzdf_imp      ! vertical diffusion: implicit (dyn_zdf_imp     routine) 
    19    USE dynzdf_imp_jki  ! vertical diffusion  implicit (dyn_zdf_imp_jki routine) 
    2019 
    2120   USE ldfdyn_oce      ! ocean dynamics: lateral physics 
    … …  
    7473      ! 
    7574      CASE ( 0 )   ;   CALL dyn_zdf_exp    ( kt, r2dt )      ! explicit scheme 
    76       CASE ( 1 )   ;   CALL dyn_zdf_imp    ( kt, r2dt )      ! implicit scheme (k-j-i loop) 
    77       CASE ( 2 )   ;   CALL dyn_zdf_imp_jki( kt, r2dt )      ! implicit scheme (j-k-i loop) 
     75      CASE ( 1 )   ;   CALL dyn_zdf_imp    ( kt, r2dt )      ! implicit scheme 
    7876      ! 
    7977      CASE ( -1 )                                      ! esopa: test all possibility with control print 
    … …  
    8381                       CALL dyn_zdf_imp    ( kt, r2dt ) 
    8482                       CALL prt_ctl( tab3d_1=ua, clinfo1=' zdf1 - Ua: ', mask1=umask,               & 
    85             &                        tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' ) 
    86                        CALL dyn_zdf_imp_jki( kt, r2dt ) 
    87                        CALL prt_ctl( tab3d_1=ua, clinfo1=' zdf2 - Ua: ', mask1=umask,               & 
    8883            &                        tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' ) 
    8984      END SELECT 
    … …  
    109104      !! ** Method  :   implicit (euler backward) scheme (default) 
    110105      !!                explicit (time-splitting) scheme if ln_zdfexp=T 
    111       !!                OpenMP / NEC autotasking: use j-k-i loops 
    112106      !!---------------------------------------------------------------------- 
    113107      USE zdftke 
    … …  
    125119      IF( ln_dynldf_hor .AND. ln_sco )   nzdf = 1   ! horizontal lateral physics in s-coordinate 
    126120 
    127       ! OpenMP / NEC autotasking 
    128 #if defined key_mpp_omp 
    129       IF( nzdf == 1 )   nzdf = 2                    ! j-k-i loop 
    130 #endif 
    131  
    132121      IF( lk_esopa )    nzdf = -1                   ! Esopa key: All schemes used 
    133122 
    … …  
    139128         IF( nzdf ==  0 )   WRITE(numout,*) '              Explicit time-splitting scheme' 
    140129         IF( nzdf ==  1 )   WRITE(numout,*) '              Implicit (euler backward) scheme' 
    141          IF( nzdf ==  2 )   WRITE(numout,*) '              Implicit (euler backward) scheme with j-k-i loops' 
    142130      ENDIF 
    143131      ! 

  • branches/dev_001_SBC/NEMO/OPA_SRC/DYN/wzvmod.F90

    r708 r811  
    3535 
    3636CONTAINS 
    37  
    38 #if defined key_mpp_omp 
    39    !!---------------------------------------------------------------------- 
    40    !!   'key_mpp_omp'                                   j-k-i loop (j-slab) 
    41    !!---------------------------------------------------------------------- 
    42  
    43    SUBROUTINE wzv( kt ) 
    44       !!---------------------------------------------------------------------- 
    45       !!                    ***  ROUTINE wzv  *** 
    46       !!                      
    47       !! ** Purpose :   Compute the now vertical velocity after the array swap 
    48       !! 
    49       !! ** Method  :   Using the incompressibility hypothesis, the vertical 
    50       !!     velocity is computed by integrating the horizontal divergence  
    51       !!     from the bottom to the surface. 
    52       !!     The boundary conditions are w=0 at the bottom (no flux) and, 
    53       !!     in rigid-lid case, w=0 at the sea surface. 
    54       !! 
    55       !! ** action  :    wn array : the now vertical velocity 
    56       !!---------------------------------------------------------------------- 
    57       !! * Arguments 
    58       INTEGER, INTENT( in ) ::   kt      ! ocean time-step index 
    59  
    60       !! * Local declarations 
    61       INTEGER ::   jj, jk      ! dummy loop indices 
    62       !!---------------------------------------------------------------------- 
    63  
    64       IF( kt == nit000 ) THEN 
    65          IF(lwp) WRITE(numout,*) 
    66          IF(lwp) WRITE(numout,*) 'wzv     : vertical velocity from continuity eq.' 
    67          IF(lwp) WRITE(numout,*) '~~~~~~~                     j-k-i loops' 
    68  
    69          ! bottom boundary condition: w=0 (set once for all) 
    70          wn(:,:,jpk) = 0.e0 
    71       ENDIF 
    72  
    73       !                                                ! =============== 
    74       DO jj = 1, jpj                                   !  Vertical slab 
    75          !                                             ! =============== 
    76          ! Computation from the bottom 
    77          DO jk = jpkm1, 1, -1 
    78             wn(:,jj,jk) = wn(:,jj,jk+1) - fse3t(:,jj,jk) * hdivn(:,jj,jk) 
    79          END DO 
    80          !                                             ! =============== 
    81       END DO                                           !   End of slab 
    82       !                                                ! =============== 
    83  
    84       IF(ln_ctl)   CALL prt_ctl(tab3d_1=wn, clinfo1=' w**2 -   : ', mask1=wn) 
    85  
    86    END SUBROUTINE wzv 
    87  
    88 #else 
    89    !!---------------------------------------------------------------------- 
    90    !!   Default option                                           k-j-i loop 
    91    !!---------------------------------------------------------------------- 
    9237 
    9338   SUBROUTINE wzv( kt ) 
    … …  
    188133 
    189134   END SUBROUTINE wzv 
    190 #endif 
    191135 
    192136   !!====================================================================== 

  • branches/dev_001_SBC/NEMO/OPA_SRC/LDF/ldfeiv.F90

    r717 r811  
    1010   !!---------------------------------------------------------------------- 
    1111   !!   ldf_eiv      : compute the eddy induced velocity coefficients 
    12    !!                  Same results but not same routine if 'key_mpp_omp' 
    13    !!                  is defined or not 
    1412   !!---------------------------------------------------------------------- 
    1513   !! * Modules used 
    … …  
    4139 
    4240CONTAINS 
    43  
    44 # if defined key_mpp_omp 
    45    !!---------------------------------------------------------------------- 
    46    !!   'key_mpp_omp' :                  OpenMP /  NEC autotasking (j-slab) 
    47    !!---------------------------------------------------------------------- 
    48  
    49    SUBROUTINE ldf_eiv( kt ) 
    50       !!---------------------------------------------------------------------- 
    51       !!                  ***  ROUTINE ldf_eiv  *** 
    52       !! 
    53       !! ** Purpose :   Compute the eddy induced velocity coefficient from the 
    54       !!      growth rate of baroclinic instability. 
    55       !! 
    56       !! ** Method : 
    57       !! 
    58       !! ** Action :   uslp(),   : i- and j-slopes of neutral surfaces 
    59       !!               vslp()      at u- and v-points, resp. 
    60       !!               wslpi(),  : i- and j-slopes of neutral surfaces 
    61       !!               wslpj()     at w-points.  
    62       !! 
    63       !! History : 
    64       !!   8.1  !  99-03  (G. Madec, A. Jouzeau)  Original code 
    65       !!   8.5  !  02-06  (G. Madec)  Free form, F90 
    66       !!---------------------------------------------------------------------- 
    67       !! * Arguments 
    68       INTEGER, INTENT( in ) ::   kt     ! ocean time-step inedx 
    69        
    70       !! * Local declarations 
    71       INTEGER ::   ji, jj, jk           ! dummy loop indices 
    72       REAL(wp) ::   & 
    73          zfw, ze3w, zn2, zf20,       &  ! temporary scalars 
    74          zaht, zaht_min 
    75       REAL(wp), DIMENSION(jpi,jpj) ::   & 
    76          zn, zah, zhw, zross            ! workspace 
    77       !!---------------------------------------------------------------------- 
    78  
    79       IF( kt == nit000 ) THEN 
    80          IF(lwp) WRITE(numout,*) 
    81          IF(lwp) WRITE(numout,*) 'ldf_eiv : eddy induced velocity coefficients' 
    82          IF(lwp) WRITE(numout,*) '~~~~~~~   NEC autotasking / OpenMP : j-slab' 
    83       ENDIF 
    84        
    85       !                                                ! =============== 
    86       DO jj = 2, jpjm1                                 !  Vertical slab 
    87          !                                             ! =============== 
    88           
    89          ! 0. Local initialization 
    90          ! ----------------------- 
    91          zn   (:,jj) = 0.e0 
    92          zhw  (:,jj) = 5.e0 
    93          zah  (:,jj) = 0.e0 
    94          zross(:,jj) = 0.e0 
    95           
    96          ! 1. Compute lateral diffusive coefficient  
    97          ! ---------------------------------------- 
    98  
    99 !CDIR NOVERRCHK  
    100          DO jk = 1, jpk 
    101 !CDIR NOVERRCHK  
    102             DO ji = 2, jpim1 
    103                ! Take the max of N^2 and zero then take the vertical sum  
    104                ! of the square root of the resulting N^2 ( required to compute  
    105                ! internal Rossby radius Ro = .5 * sum_jpk(N) / f  
    106                zn2 = MAX( rn2(ji,jj,jk), 0.e0 ) 
    107                ze3w = fse3w(ji,jj,jk) * tmask(ji,jj,jk) 
    108                zn(ji,jj) = zn(ji,jj) + SQRT( zn2 ) * fse3w(ji,jj,jk) 
    109                ! Compute elements required for the inverse time scale of baroclinic 
    110                ! eddies using the isopycnal slopes calculated in ldfslp.F :  
    111                ! T^-1 = sqrt(m_jpk(N^2*(r1^2+r2^2)*e3w)) 
    112                zah(ji,jj) = zah(ji,jj) + zn2   & 
    113                               * ( wslpi(ji,jj,jk) * wslpi(ji,jj,jk)    & 
    114                                 + wslpj(ji,jj,jk) * wslpj(ji,jj,jk) )   & 
    115                               * ze3w 
    116                zhw(ji,jj) = zhw(ji,jj) + ze3w 
    117             END DO  
    118          END DO  
    119   
    120 !CDIR NOVERRCHK  
    121          DO ji = 2, jpim1 
    122             zfw = MAX( ABS( 2. * omega * SIN( rad * gphit(ji,jj) ) ) , 1.e-10 ) 
    123             ! Rossby radius at w-point taken < 40km and  > 2km 
    124             zross(ji,jj) = MAX( MIN( .4 * zn(ji,jj) / zfw, 40.e3 ), 2.e3 ) 
    125             ! Compute aeiw by multiplying Ro^2 and T^-1 
    126             aeiw(ji,jj) = zross(ji,jj) * zross(ji,jj) * SQRT( zah(ji,jj) / zhw(ji,jj) ) * tmask(ji,jj,1) 
    127             IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN   ! ORCA R02 
    128                ! Take the minimum between aeiw and aeiv0 for depth levels 
    129                ! lower than 20 (21 in w- point) 
    130                IF( mbathy(ji,jj) <= 21. ) aeiw(ji,jj) = MIN( aeiw(ji,jj), 1000. ) 
    131             ENDIF 
    132          END DO 
    133  
    134          ! Decrease the coefficient in the tropics (20N-20S)  
    135          zf20 = 2. * omega * sin( rad * 20. ) 
    136          DO ji = 2, jpim1 
    137             aeiw(ji,jj) = MIN( 1., ABS( ff(ji,jj) / zf20 ) ) * aeiw(ji,jj) 
    138          END DO 
    139    
    140          ! ORCA R05: Take the minimum between aeiw and aeiv0 
    141          IF( cp_cfg == "orca" .AND. jp_cfg == 05 ) THEN   ! ORCA R05 
    142             DO ji = 2, jpim1 
    143                aeiw(ji,jj) = MIN( aeiw(ji,jj), aeiv0 ) 
    144             END DO 
    145          ENDIF 
    146          !                                             ! =============== 
    147       END DO                                           !   End of slab 
    148       !                                                ! =============== 
    149  
    150       !,,,,,,,,,,,,,,,,,,,,,,,,,,,,,synchro,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,, 
    151  
    152       ! lateral boundary condition on aeiw  
    153       CALL lbc_lnk( aeiw, 'W', 1. ) 
    154  
    155       ! Average the diffusive coefficient at u- v- points  
    156       DO jj = 2, jpjm1 
    157          DO ji = fs_2, fs_jpim1   ! vector opt. 
    158             aeiu(ji,jj) = .5 * (aeiw(ji,jj) + aeiw(ji+1,jj  )) 
    159             aeiv(ji,jj) = .5 * (aeiw(ji,jj) + aeiw(ji  ,jj+1)) 
    160          END DO  
    161       END DO  
    162       !,,,,,,,,,,,,,,,,,,,,,,,,,,,,,synchro,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,, 
    163  
    164       ! lateral boundary condition on aeiu, aeiv  
    165       CALL lbc_lnk( aeiu, 'U', 1. ) 
    166       CALL lbc_lnk( aeiv, 'V', 1. ) 
    167  
    168       IF(ln_ctl)   THEN 
    169          CALL prt_ctl(tab2d_1=aeiu, clinfo1=' eiv  - u: ', ovlap=1) 
    170          CALL prt_ctl(tab2d_1=aeiv, clinfo1=' eiv  - v: ', ovlap=1) 
    171       ENDIF 
    172        
    173       ! ORCA R05: add a space variation on aht (=aeiv except at the equator and river mouth) 
    174       IF( cp_cfg == "orca" .AND. jp_cfg == 05 ) THEN 
    175          zf20     = 2. * omega * SIN( rad * 20. ) 
    176          zaht_min = 100.                              ! minimum value for aht 
    177          DO jj = 1, jpj 
    178             DO ji = 1, jpi 
    179                zaht      = ( 1. -  MIN( 1., ABS( ff(ji,jj) / zf20 ) ) ) * ( aht0 - zaht_min )  & 
    180                   &      + aht0 * rnfmsk(ji,jj)                          ! enhanced near river mouths 
    181                ahtu(ji,jj) = MAX( MAX( zaht_min, aeiu(ji,jj) ) + zaht, aht0 ) 
    182                ahtv(ji,jj) = MAX( MAX( zaht_min, aeiv(ji,jj) ) + zaht, aht0 ) 
    183                ahtw(ji,jj) = MAX( MAX( zaht_min, aeiw(ji,jj) ) + zaht, aht0 ) 
    184             END DO 
    185          END DO 
    186          IF(ln_ctl) THEN 
    187             CALL prt_ctl(tab2d_1=ahtu, clinfo1=' aht  - u: ', ovlap=1) 
    188             CALL prt_ctl(tab2d_1=ahtv, clinfo1=' aht  - v: ', ovlap=1) 
    189             CALL prt_ctl(tab2d_1=ahtw, clinfo1=' aht  - w: ', ovlap=1) 
    190          ENDIF 
    191       ENDIF 
    192  
    193       IF( aeiv0 == 0.e0 ) THEN 
    194          aeiu(:,:) = 0.e0 
    195          aeiv(:,:) = 0.e0 
    196          aeiw(:,:) = 0.e0 
    197       ENDIF 
    198  
    199    END SUBROUTINE ldf_eiv 
    200  
    201 # else 
    202    !!---------------------------------------------------------------------- 
    203    !!   Default key                                             k-j-i loops 
    204    !!---------------------------------------------------------------------- 
    20541 
    20642   SUBROUTINE ldf_eiv( kt ) 
    … …  
    25288 
    25389      DO jk = 1, jpk 
    254 #  if defined key_vectopt_loop  &&  ! defined key_mpp_omp 
     90#  if defined key_vectopt_loop   
    25591!CDIR NOVERRCHK  
    25692         DO ji = 1, jpij   ! vector opt. 
    … …  
    374210   END SUBROUTINE ldf_eiv 
    375211 
    376 # endif 
    377  
    378212#else 
    379213   !!---------------------------------------------------------------------- 

  • branches/dev_001_SBC/NEMO/OPA_SRC/LDF/ldfslp.F90

    r699 r811  
    141141 
    142142      IF( ln_zps ) THEN      ! partial steps correction at the bottom ocean level (zps_hde routine) 
    143 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     143# if defined key_vectopt_loop   
    144144         jj = 1 
    145145         DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    … …  
    153153               zgru(ji,jj,iku) = gru(ji,jj)  
    154154               zgrv(ji,jj,ikv) = grv(ji,jj)                
    155 # if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     155# if ! defined key_vectopt_loop  
    156156            END DO 
    157157# endif 
    … …  
    492492      ! mask for mixed layer 
    493493      DO jk = 1, jpk 
    494 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     494# if defined key_vectopt_loop 
    495495         jj = 1 
    496496         DO ji = 1, jpij   ! vector opt. (forced unrolling) 
    … …  
    506506                  omlmask(ji,jj,jk) = 0.e0 
    507507               ENDIF 
    508 # if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     508# if ! defined key_vectopt_loop 
    509509            END DO 
    510510# endif 
    … …  
    524524      zwy(:,jpj) = 0.e0 
    525525      zwy(jpi,:) = 0.e0 
    526 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     526# if defined key_vectopt_loop 
    527527      jj = 1 
    528528      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    … …  
    537537               &             * ( pn2(ji,jj,ik) + pn2(ji,jj,ik+1) )   & 
    538538               &             / MAX( tmask(ji,jj,ik) + tmask (ji,jj,ik+1), 1. ) 
    539 # if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     539# if ! defined key_vectopt_loop 
    540540         END DO 
    541541# endif 
    … …  
    545545 
    546546      ! Slope at u points 
    547 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     547# if defined key_vectopt_loop 
    548548      jj = 1 
    549549      DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling) 
    … …  
    562562            ! uslpml 
    563563            uslpml (ji,jj) = zau / ( zbu - zeps ) * umask (ji,jj,ik) 
    564 # if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     564# if ! defined key_vectopt_loop 
    565565         END DO 
    566566# endif 
    … …  
    574574      zwy ( :, jpj) = 0.e0 
    575575      zwy ( jpi, :) = 0.e0 
    576 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     576# if defined key_vectopt_loop 
    577577      jj = 1 
    578578      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    … …  
    586586               &             * ( pn2(ji,jj,ik) + pn2(ji,jj,ik+1) )   & 
    587587               &             / MAX( tmask(ji,jj,ik) + tmask (ji,jj,ik+1), 1. ) 
    588 # if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     588# if ! defined key_vectopt_loop 
    589589         END DO 
    590590# endif 
    … …  
    595595 
    596596      ! Slope at v points 
    597 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     597# if defined key_vectopt_loop 
    598598      jj = 1 
    599599      DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling) 
    … …  
    612612            ! vslpml 
    613613            vslpml (ji,jj) = zav / ( zbv - zeps ) * vmask (ji,jj,ik) 
    614 # if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     614# if ! defined key_vectopt_loop 
    615615         END DO 
    616616# endif 
    … …  
    626626      ! Local vertical density gradient evaluated from N^2 
    627627      ! zwy = d/dz(prd)= - mk ( prd ) / grav * pn2 -- at w point 
    628 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     628# if defined key_vectopt_loop 
    629629      jj = 1 
    630630      DO ji = 1, jpij   ! vector opt. (forced unrolling) 
    … …  
    638638            zwy (ji,jj) = zm05g * pn2 (ji,jj,ik) *     & 
    639639               &             ( prd (ji,jj,ik) + prd (ji,jj,ikm1) + 2. ) 
    640 # if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     640# if ! defined key_vectopt_loop 
    641641         END DO 
    642642# endif 
    … …  
    644644 
    645645      ! Slope at w point 
    646 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     646# if defined key_vectopt_loop 
    647647      jj = 1 
    648648      DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling) 
    … …  
    674674            wslpiml (ji,jj) = zai / ( zbi - zeps) * tmask (ji,jj,ik) 
    675675            wslpjml (ji,jj) = zaj / ( zbj - zeps) * tmask (ji,jj,ik) 
    676 # if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     676# if ! defined key_vectopt_loop 
    677677         END DO 
    678678# endif 

  • branches/dev_001_SBC/NEMO/OPA_SRC/SOL/sol_oce.F90

    r699 r811  
    2424   !! ---------------------------------- 
    2525   INTEGER , PUBLIC ::      & !!: namsol   elliptic solver / island / free surface 
    26       nsolv    =    1 ,     &  !: = 1/2/3/4 type of elliptic solver 
     26      nsolv    =    1 ,     &  !: = 1/2/3 type of elliptic solver 
    2727      nsol_arp =    0 ,     &  !: = 0/1 absolute/relative precision convergence test 
    2828      nmin     =  300 ,     &  !: minimum of iterations for the SOR solver 

  • branches/dev_001_SBC/NEMO/OPA_SRC/SOL/solmat.F90

    r699 r811  
    291291         gcp(:,:,3) = gcp(:,:,3) * gcdprc(:,:) 
    292292         gcp(:,:,4) = gcp(:,:,4) * gcdprc(:,:) 
    293          IF( ( nsolv == 2 ) .OR. ( nsolv == 4 ) )  gccd(:,:) = sor * gcp(:,:,2) 
    294  
    295          IF( nsolv == 4 ) THEN 
     293         IF( nsolv == 2 )  gccd(:,:) = sor * gcp(:,:,2) 
     294 
     295         IF( nsolv == 2 .AND. MAX( jpr2di, jpr2dj ) > 0) THEN 
    296296            CALL lbc_lnk_e( gcp   (:,:,1), c_solver_pt, 1. )   ! lateral boundary conditions 
    297297            CALL lbc_lnk_e( gcp   (:,:,2), c_solver_pt, 1. )   ! lateral boundary conditions 

  • branches/dev_001_SBC/NEMO/OPA_SRC/SOL/solpcg.F90

    r699 r811  
    5252      !!      where q is the preconditioning matrix = diagonal matrix of the 
    5353      !!                                              diagonal elements of a 
    54       !!      Initialization: 
     54      !!      Initialization  : 
    5555      !!         x(o) = gcx 
    5656      !!         r(o) = d(o) = pgcb - pa.x(o) 
    5757      !!         rr(o)= < r(o) , r(o) >_q 
    58       !!      Iteration n   : 
    59       !!         z(n)   = pa.d(n) 
    60       !!         alp(n) = rr(n) / < z(n) , d(n) >_q 
     58      !!      Iteration 1     : 
     59      !!         standard PCG algorithm 
     60      !!      Iteration n > 1 : 
     61      !!         s(n)   = pa.r(n) 
     62      !!         gam(n) = < r(n) , r(n) >_q 
     63      !!         del(n) = < r(n) , s(n) >_q 
     64      !!         bet(n) = gam(n) / gam(n-1) 
     65      !!         d(n)   = r(n) + bet(n) d(n-1) 
     66      !!         z(n)   = s(n) + bet(n) z(n-1)  
     67      !!         sig(n) = del(n) - bet(n)*bet(n)*sig(n-1)  
     68      !!         alp(n) = gam(n) / sig(n)  
    6169      !!         x(n+1) = x(n) + alp(n) d(n) 
    6270      !!         r(n+1) = r(n) - alp(n) z(n) 
    63       !!         rr(n+1)= < r(n+1) , r(n+1) >_q 
    64       !!         bet(n) = rr(n+1) / rr(n) 
    65       !!         r(n+1) = r(n+1) + bet(n+1) d(n) 
    6671      !!      Convergence test : 
    6772      !!         rr(n+1) / < gcb , gcb >_q   =< epsr 
    … …  
    7378      !! References : 
    7479      !!      Madec et al. 1988, Ocean Modelling, issue 78, 1-6. 
     80      !!      D Azevedo et al. 1993, Computer Science Technical Report, Tennessee U. 
    7581      !! 
    7682      !! History : 
    … …  
    8288      !!        !  96-11  (A. Weaver)  correction to preconditioning 
    8389      !!   8.5  !  02-08  (G. Madec)  F90: Free form 
     90      !!        !  08-01  (R. Benshila) mpp optimization 
    8491      !!---------------------------------------------------------------------- 
    8592      !! * Arguments 
    … …  
    9198      !! * Local declarations 
    9299      INTEGER ::   ji, jj, jn                ! dummy loop indices 
    93       REAL(wp) ::   zgcad                    ! temporary scalars 
     100      REAL(wp) ::  zgcad                     ! temporary scalars 
     101      REAL(wp), DIMENSION(2) :: zsum 
     102      REAL(wp), DIMENSION(jpi,jpj) :: zgcr 
    94103      !!---------------------------------------------------------------------- 
    95104 
     105      ! Initialization of the algorithm with standard PCG 
     106      ! ------------------------------------------------- 
     107 
     108      CALL lbc_lnk( gcx, c_solver_pt, 1. )   ! lateral boundary condition 
     109 
     110      ! gcr   = gcb-a.gcx 
     111      ! gcdes = gcr 
     112 
     113      DO jj = 2, jpjm1 
     114         DO ji = fs_2, fs_jpim1   ! vector opt. 
     115            zgcad = bmask(ji,jj) * ( gcb(ji,jj  ) -                gcx(ji  ,jj  )   & 
     116               &                                  - gcp(ji,jj,1) * gcx(ji  ,jj-1)   & 
     117               &                                  - gcp(ji,jj,2) * gcx(ji-1,jj  )   & 
     118               &                                  - gcp(ji,jj,3) * gcx(ji+1,jj  )   & 
     119               &                                  - gcp(ji,jj,4) * gcx(ji  ,jj+1)   ) 
     120            gcr  (ji,jj) = zgcad 
     121            gcdes(ji,jj) = zgcad 
     122         END DO 
     123      END DO 
     124 
     125      ! rnorme = (gcr,gcr) 
     126      rnorme = SUM(  gcr(:,:) * gcdmat(:,:) * gcr(:,:)  ) 
     127      IF( lk_mpp )   CALL mpp_sum( rnorme )   ! sum over the global domain 
     128 
     129      CALL lbc_lnk( gcdes, c_solver_pt, 1. )   ! lateral boundary condition 
     130 
     131      ! gccd = matrix . gcdes 
     132      DO jj = 2, jpjm1 
     133         DO ji = fs_2, fs_jpim1   ! vector opt. 
     134            gccd(ji,jj) = bmask(ji,jj)*( gcdes(ji,jj)   & 
     135               &        +gcp(ji,jj,1)*gcdes(ji,jj-1)+gcp(ji,jj,2)*gcdes(ji-1,jj)   & 
     136               &        +gcp(ji,jj,4)*gcdes(ji,jj+1)+gcp(ji,jj,3)*gcdes(ji+1,jj)   ) 
     137         END DO 
     138      END DO  
     139 
     140      ! alph = (gcr,gcr)/(gcdes,gccd) 
     141      radd = SUM(  gcdes(:,:) * gcdmat(:,:) * gccd(:,:)  ) 
     142      IF( lk_mpp )   CALL mpp_sum( radd )   ! sum over the global domain 
     143      alph = rnorme /radd 
     144 
     145      ! gcx = gcx + alph * gcdes 
     146      ! gcr = gcr - alph * gccd 
     147      DO jj = 2, jpjm1 
     148         DO ji = fs_2, fs_jpim1   ! vector opt. 
     149            gcx(ji,jj) = bmask(ji,jj) * ( gcx(ji,jj) + alph * gcdes(ji,jj) ) 
     150            gcr(ji,jj) = bmask(ji,jj) * ( gcr(ji,jj) - alph * gccd (ji,jj) ) 
     151         END DO 
     152      END DO 
     153 
     154      ! Algorithm wtih Eijkhout rearrangement 
     155      ! ------------------------------------- 
     156         
    96157      !                                                !================ 
    97158      DO jn = 1, nmax                                  ! Iterative loop 
    98159         !                                             !================ 
    99160 
    100          IF( jn == 1 ) THEN           ! Initialization of the algorithm 
    101  
    102             CALL lbc_lnk( gcx, c_solver_pt, 1. )   ! lateral boundary condition 
    103  
    104             ! gcr   = gcb-a.gcx 
    105             ! gcdes = gsr 
    106             DO jj = 2, jpjm1 
    107                DO ji = fs_2, fs_jpim1   ! vector opt. 
    108                   zgcad = bmask(ji,jj) * ( gcb(ji,jj  ) -                gcx(ji  ,jj  )   & 
    109                      &                                  - gcp(ji,jj,1) * gcx(ji  ,jj-1)   & 
    110                      &                                  - gcp(ji,jj,2) * gcx(ji-1,jj  )   & 
    111                      &                                  - gcp(ji,jj,3) * gcx(ji+1,jj  )   & 
    112                      &                                  - gcp(ji,jj,4) * gcx(ji  ,jj+1)   ) 
    113                   gcr  (ji,jj) = zgcad 
    114                   gcdes(ji,jj) = zgcad 
    115                END DO 
    116             END DO 
    117              
    118             rnorme = SUM(  gcr(:,:) * gcdmat(:,:) * gcr(:,:)  ) 
    119             IF( lk_mpp )   CALL mpp_sum( rnorme )   ! sum over the global domain 
    120             rr = rnorme 
    121  
    122          ENDIF 
    123          
    124          !                             ! Algorithm 
    125          
    126          CALL lbc_lnk( gcdes, c_solver_pt, 1. )   ! lateral boundary condition 
    127          
    128          ! ... gccd = matrix . gcdes 
     161         CALL lbc_lnk( gcr, c_solver_pt, 1. )   ! lateral boundary condition 
     162 
     163         ! zgcr = matrix . gcr 
    129164         DO jj = 2, jpjm1 
    130165            DO ji = fs_2, fs_jpim1   ! vector opt. 
    131                gccd(ji,jj) = bmask(ji,jj)*( gcdes(ji,jj)   & 
    132                   &        +gcp(ji,jj,1)*gcdes(ji,jj-1)+gcp(ji,jj,2)*gcdes(ji-1,jj)   & 
    133                   &        +gcp(ji,jj,4)*gcdes(ji,jj+1)+gcp(ji,jj,3)*gcdes(ji+1,jj)   ) 
     166               zgcr(ji,jj) = bmask(ji,jj)*( gcr(ji,jj)   & 
     167                  &        +gcp(ji,jj,1)*gcr(ji,jj-1)+gcp(ji,jj,2)*gcr(ji-1,jj)   & 
     168                  &        +gcp(ji,jj,4)*gcr(ji,jj+1)+gcp(ji,jj,3)*gcr(ji+1,jj)   ) 
    134169            END DO 
    135170         END DO 
    136171  
    137          ! alph = (gcr,gcr)/(gcdes,gccd) 
    138          radd = SUM(  gcdes(:,:) * gcdmat(:,:) * gccd(:,:)  ) 
    139          IF( lk_mpp )   CALL mpp_sum( radd )   ! sum over the global domain 
    140          alph = rr / radd 
    141           
    142          ! gcx = gcx + alph * gcdes 
    143          ! gcr = gcr - alph * gccd 
    144          DO jj = 2, jpjm1 
    145             DO ji = fs_2, fs_jpim1   ! vector opt. 
    146                gcx(ji,jj) = bmask(ji,jj) * ( gcx(ji,jj) + alph * gcdes(ji,jj) ) 
    147                gcr(ji,jj) = bmask(ji,jj) * ( gcr(ji,jj) - alph * gccd (ji,jj) ) 
    148             END DO 
    149          END DO 
    150          
    151172         ! rnorme = (gcr,gcr) 
    152          rnorme = SUM(  gcr(:,:) * gcdmat(:,:) * gcr(:,:)  ) 
    153          IF( lk_mpp )   CALL  mpp_sum( rnorme )   ! sum over the global domain 
    154          
     173         rr = rnorme 
     174         zsum(1) = SUM(  gcr(:,:) * gcdmat(:,:) * gcr(:,:)  ) 
     175 
     176         ! zgcad = (zgcr,gcr)  
     177         zsum(2) = SUM( gcr(2:jpim1,2:jpjm1) * gcdmat(2:jpim1,2:jpjm1) * zgcr(2:jpim1,2:jpjm1) ) 
     178 
     179         IF( lk_mpp )   CALL mpp_sum( zsum, 2 )   ! sum over the global domain 
     180         rnorme = zsum(1)   
     181         zgcad  = zsum(2) 
     182 
    155183         ! test of convergence 
    156184         IF( rnorme < epsr .OR. jn == nmax ) THEN 
    … …  
    159187            ncut = 999 
    160188         ENDIF 
    161          
     189 
    162190         ! beta = (rk+1,rk+1)/(rk,rk) 
    163191         beta = rnorme / rr 
    164          rr   = rnorme 
    165  
     192         radd = zgcad - beta*beta*radd 
     193         alph = rnorme / radd 
     194 
     195         ! gcx = gcx + alph * gcdes 
     196         ! gcr = gcr - alph * gccd 
     197         DO jj = 2, jpjm1 
     198            DO ji = fs_2, fs_jpim1   ! vector opt. 
     199               gcdes(ji,jj) = gcr (ji,jj) + beta * gcdes(ji,jj)  
     200               gccd (ji,jj) = zgcr(ji,jj) + beta * gccd (ji,jj)  
     201               gcx  (ji,jj) = gcx (ji,jj) + alph * gcdes(ji,jj)  
     202               gcr  (ji,jj) = gcr (ji,jj) - alph * gccd (ji,jj)  
     203            END DO 
     204         END DO 
     205         
    166206         ! indicator of non-convergence or explosion 
    167207         IF( jn == nmax .OR. SQRT(epsr)/eps > 1.e+20 ) kindic = -2 
    168208         IF( ncut == 999 ) GOTO 999 
    169209 
    170          ! gcdes = gcr + beta * gcdes 
    171          DO jj = 2, jpjm1 
    172             DO ji = fs_2, fs_jpim1   ! vector opt. 
    173                gcdes(ji,jj) = bmask(ji,jj)*( gcr(ji,jj) + beta * gcdes(ji,jj) ) 
    174             END DO 
    175          END DO 
    176          
    177210         !                                             !================ 
    178211      END DO                                           !    End Loop 

  • branches/dev_001_SBC/NEMO/OPA_SRC/SOL/solsor.F90

    r699 r811  
    4343      !!     as well as islands) while in the latter the boundary condition 
    4444      !!     specification is not required. 
    45       !! 
     45      !!       This routine provides a MPI optimization to the existing solsor 
     46      !!     by reducing the number of call to lbc. 
     47      !!  
    4648      !! ** Method  :   Successive-over-relaxation method using the red-black  
    4749      !!      technique. The former technique used was not compatible with  
    4850      !!      the north-fold boundary condition used in orca configurations. 
    49       !! 
     51      !!      Compared to the classical sol_sor, this routine provides a  
     52      !!      mpp optimization by reducing the number of calls to lnc_lnk 
     53      !!      The solution is computed on a larger area and the boudary 
     54      !!      conditions only when the inside domain is reached. 
     55      !!  
    5056      !! References : 
    5157      !!      Madec et al. 1988, Ocean Modelling, issue 78, 1-6. 
     58      !!      Beare and Stevens 1997 Ann. Geophysicae 15, 1369-1377 
    5259      !! 
    5360      !! History : 
    … …  
    5865      !!        !  96-11  (A. Weaver)  correction to preconditioning 
    5966      !!   9.0  !  03-04  (C. Deltel, G. Madec)  Red-Black SOR in free form 
     67      !!   9.0  !  05-09  (R. Benshila, G. Madec)  MPI optimization 
    6068      !!---------------------------------------------------------------------- 
    6169      !! * Arguments 
    … …  
    6674      !! * Local declarations 
    6775      INTEGER  ::   ji, jj, jn               ! dummy loop indices 
    68       INTEGER  ::   ishift 
     76      INTEGER  ::   ishift, icount 
    6977      REAL(wp) ::   ztmp, zres, zres2 
    7078 
    7179      INTEGER  ::   ijmppodd, ijmppeven 
     80      INTEGER  ::   ijpr2d 
    7281      !!---------------------------------------------------------------------- 
    7382       
    74       ijmppeven = MOD(nimpp+njmpp  ,2) 
    75       ijmppodd  = MOD(nimpp+njmpp+1,2) 
     83      ijmppeven = MOD(nimpp+njmpp+jpr2di+jpr2dj,2) 
     84      ijmppodd  = MOD(nimpp+njmpp+jpr2di+jpr2dj+1,2) 
     85      ijpr2d = MAX(jpr2di,jpr2dj) 
     86      icount = 0 
    7687      !                                                       ! ============== 
    7788      DO jn = 1, nmax                                         ! Iterative loop  
    7889         !                                                    ! ============== 
    7990 
    80          CALL lbc_lnk( gcx, c_solver_pt, 1. )   ! applied the lateral boundary conditions 
    81           
     91         ! applied the lateral boundary conditions 
     92         IF( MOD(icount,ijpr2d+1) == 0 ) CALL lbc_lnk_e( gcx, c_solver_pt, 1. )    
     93         
    8294         ! Residus 
    8395         ! ------- 
    8496 
    8597         ! Guess black update 
    86          DO jj = 2, jpjm1 
    87             ishift = MOD( jj-ijmppodd, 2 ) 
    88             DO ji = 2+ishift, jpim1, 2 
     98         DO jj = 2-jpr2dj, nlcj-1+jpr2dj 
     99            ishift = MOD( jj-ijmppodd-jpr2dj, 2 ) 
     100            DO ji = 2-jpr2di+ishift, nlci-1+jpr2di, 2 
    89101               ztmp =                  gcb(ji  ,jj  )   & 
    90102                  &   - gcp(ji,jj,1) * gcx(ji  ,jj-1)   & 
    … …  
    99111            END DO 
    100112         END DO 
    101  
    102          CALL lbc_lnk( gcx, c_solver_pt, 1. )   ! applied the lateral boubary conditions 
     113         icount = icount + 1  
     114  
     115         ! applied the lateral boundary conditions 
     116         IF( MOD(icount,ijpr2d+1) == 0 ) CALL lbc_lnk_e( gcx, c_solver_pt, 1. )   
    103117 
    104118         ! Guess red update 
    105          DO jj = 2, jpjm1 
    106             ishift = MOD( jj-ijmppeven, 2 ) 
    107             DO ji = 2+ishift, jpim1, 2 
     119         DO jj = 2-jpr2dj, nlcj-1+jpr2dj 
     120            ishift = MOD( jj-ijmppeven-jpr2dj, 2 ) 
     121            DO ji = 2-jpr2di+ishift, nlci-1+jpr2di, 2 
    108122               ztmp =                  gcb(ji  ,jj  )   & 
    109123                  &   - gcp(ji,jj,1) * gcx(ji  ,jj-1)   & 
    … …  
    118132            END DO 
    119133         END DO 
     134         icount = icount + 1 
    120135 
    121136         ! test of convergence 
    … …  
    124139            SELECT CASE ( nsol_arp ) 
    125140            CASE ( 0 )                 ! absolute precision (maximum value of the residual) 
    126                zres2 = MAXVAL( gcr(2:jpim1,2:jpjm1) ) 
     141               zres2 = MAXVAL( gcr(2:nlci-1,2:nlcj-1) ) 
    127142               IF( lk_mpp )   CALL mpp_max( zres2 )   ! max over the global domain 
    128143               ! test of convergence 
    … …  
    133148               ENDIF 
    134149            CASE ( 1 )                 ! relative precision 
    135                rnorme = SUM( gcr(2:jpim1,2:jpjm1) ) 
     150               rnorme = SUM( gcr(2:nlci-1,2:nlcj-1) ) 
    136151               IF( lk_mpp )   CALL mpp_sum( rnorme )   ! sum over the global domain 
    137152               ! test of convergence 
    … …  
    163178      !  ------------- 
    164179 
    165       CALL lbc_lnk( gcx, c_solver_pt, 1. )    ! boundary conditions 
     180      CALL lbc_lnk_e( gcx, c_solver_pt, 1. )    ! boundary conditions 
    166181 
    167182       

  • branches/dev_001_SBC/NEMO/OPA_SRC/SOL/solver.F90

    r699 r811  
    152152 
    153153      CASE ( 2 )                ! successive-over-relaxation solver 
    154          IF(lwp) WRITE(numout,*) '          a successive-over-relaxation solver is used' 
    155          IF( jpr2di /= 0 .AND. jpr2dj /= 0 ) & 
    156              CALL ctl_stop( ' jpr2di and jpr2dj should be equal to zero' ) 
     154         IF(lwp) WRITE(numout,*) '          a successive-over-relaxation solver with extra outer halo is used' 
     155         IF(lwp) WRITE(numout,*) '          with jpr2di =', jpr2di, ' and  jpr2dj =', jpr2dj 
     156         IF( .NOT. lk_mpp .AND. jpr2di /= 0 .AND. jpr2dj /= 0 ) THEN 
     157             CALL ctl_stop( ' jpr2di and jpr2dj are not equal to zero',   & 
     158             &              ' In this case this algorithm should be used only with the key_mpp_... option' ) 
     159         ELSE 
     160            IF( ( ( jperio == 1 .OR. jperio == 4 .OR. jperio == 6 ) .OR. ( jpni /= 1 ) ) & 
     161              &  .AND. ( jpr2di /= jpr2dj ) ) CALL ctl_stop( '          jpr2di should be equal to jpr2dj' ) 
     162         ENDIF 
    157163 
    158164      CASE ( 3 )                ! FETI solver 
    … …  
    169175         ENDIF 
    170176          
    171       CASE ( 4 )                ! successive-over-relaxation solver with extra outer halo 
    172          IF(lwp) WRITE(numout,*) '          a successive-over-relaxation solver with extra outer halo is used' 
    173          IF(lwp) WRITE(numout,*) '          with jpr2di =', jpr2di, ' and  jpr2dj =', jpr2dj 
    174          IF( .NOT. lk_mpp .AND. jpr2di /= 0 .AND. jpr2dj /= 0 ) THEN 
    175              CALL ctl_stop( ' jpr2di and jpr2dj are not equal to zero',   & 
    176              &              ' In this case this algorithm should be used only with the key_mpp_... option' ) 
    177          ELSE 
    178             IF( ( ( jperio == 1 .OR. jperio == 4 .OR. jperio == 6 ) .OR. ( jpni /= 1 ) ) & 
    179               &  .AND. ( jpr2di /= jpr2dj ) ) CALL ctl_stop( '          jpr2di should be equal to jpr2dj' ) 
    180          ENDIF 
    181  
    182177      CASE DEFAULT 
    183178         WRITE(ctmp1,*) '          bad flag value for nsolv = ', nsolv 
    … …  
    186181      END SELECT 
    187182 
    188       IF( nbit_cmp == 1 .AND. nsolv /= 2 ) THEN 
    189          CALL ctl_stop( ' Reproductibility tests (nbit_cmp=1) require the SOR solver: nsolv = 2' ) 
     183      IF( nbit_cmp == 1 ) THEN 
     184         IF( nsolv /= 2 ) THEN 
     185            CALL ctl_stop( ' Reproductibility tests (nbit_cmp=1) require the SOR solver: nsolv = 2' ) 
     186         ELSE IF( MAX( jpr2di, jpr2dj ) > 0 ) THEN 
     187            CALL ctl_stop( ' Reproductibility tests (nbit_cmp=1) require jpr2di = jpr2dj = 0' ) 
     188         END IF  
    190189      END IF 
    191190 

  • branches/dev_001_SBC/NEMO/OPA_SRC/TRA/traadv.F90

    r699 r811  
    1414   USE dom_oce         ! ocean space and time domain 
    1515   USE traadv_cen2     ! 2nd order centered scheme (tra_adv_cen2   routine) 
    16    USE traadv_cen2_jki ! 2nd order centered scheme (tra_adv_cen2   routine) 
    1716   USE traadv_tvd      ! TVD      scheme           (tra_adv_tvd    routine) 
    1817   USE traadv_muscl    ! MUSCL    scheme           (tra_adv_muscl  routine) 
    … …  
    8887 
    8988      SELECT CASE ( nadv )                           ! compute advection trend and add it to general trend 
    90       CASE ( 0 )   ;   CALL tra_adv_cen2    ( kt, zun, zvn, zwn )    ! 2nd order centered scheme k-j-i loops 
    91       CASE ( 1 )   ;   CALL tra_adv_cen2_jki( kt, zun, zvn, zwn )    ! 2nd order centered scheme 
     89      CASE ( 1 )   ;   CALL tra_adv_cen2    ( kt, zun, zvn, zwn )    ! 2nd order centered scheme 
    9290      CASE ( 2 )   ;   CALL tra_adv_tvd     ( kt, zun, zvn, zwn )    ! TVD      scheme 
    9391      CASE ( 3 )   ;   CALL tra_adv_muscl   ( kt, zun, zvn, zwn )    ! MUSCL    scheme 
    … …  
    9997                       CALL tra_adv_cen2    ( kt, zun, zvn, zwn ) 
    10098                       CALL prt_ctl( tab3d_1=ta, clinfo1=' adv0 - Ta: ', mask1=tmask,               & 
    101             &                        tab3d_2=sa, clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' ) 
    102                        CALL tra_adv_cen2_jki( kt, zun, zvn, zwn ) 
    103                        CALL prt_ctl( tab3d_1=ta, clinfo1=' adv1 - Ta: ', mask1=tmask,               & 
    10499            &                        tab3d_2=sa, clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' ) 
    105100                       CALL tra_adv_tvd     ( kt, zun, zvn, zwn ) 
    … …  
    171166 
    172167      !                              ! Set nadv 
    173       IF( ln_traadv_cen2   )   nadv =  0 
    174 #if defined key_mpp_omp 
    175168      IF( ln_traadv_cen2   )   nadv =  1 
    176 #endif 
    177169      IF( ln_traadv_tvd    )   nadv =  2 
    178170      IF( ln_traadv_muscl  )   nadv =  3 
    … …  
    184176      IF(lwp) THEN                   ! Print the choice 
    185177         WRITE(numout,*) 
    186          IF( nadv ==  0 )   WRITE(numout,*) '         2nd order scheme is used' 
    187          IF( nadv ==  1 )   WRITE(numout,*) '         2nd order scheme is usedi, k-j-i case' 
     178         IF( nadv ==  1 )   WRITE(numout,*) '         2nd order scheme is used' 
    188179         IF( nadv ==  2 )   WRITE(numout,*) '         TVD       scheme is used' 
    189180         IF( nadv ==  3 )   WRITE(numout,*) '         MUSCL     scheme is used' 
    190181         IF( nadv ==  4 )   WRITE(numout,*) '         MUSCL2    scheme is used' 
    191182         IF( nadv ==  5 )   WRITE(numout,*) '         UBS       scheme is used' 
    192          IF( nadv ==  6 )   WRITE(numout,*) '         PPM       scheme is used' 
    193          IF( nadv ==  7 )   WRITE(numout,*) '         QUICKEST  scheme is used' 
     183         IF( nadv ==  6 )   WRITE(numout,*) '         QUICKEST  scheme is used' 
    194184         IF( nadv == -1 )   WRITE(numout,*) '         esopa test: use all advection scheme' 
    195185      ENDIF 

  • branches/dev_001_SBC/NEMO/OPA_SRC/TRA/traadv_cen2.F90

    r717 r811  
    1515   !!   tra_adv_cen2 : update the tracer trend with the horizontal and 
    1616   !!                  vertical advection trends using a seconder order 
    17    !!                  centered scheme. (k-j-i loops) 
    1817   !!   ups_orca_set : allow mixed upstream/centered scheme in specific 
    1918   !!                  area (set for orca 2 and 4 only) 

  • branches/dev_001_SBC/NEMO/OPA_SRC/TRA/traadv_qck.F90

    r708 r811  
    5050 
    5151CONTAINS 
    52  
    53 #if ! defined key_mpp_omp 
    54    !!---------------------------------------------------------------------- 
    55    !!   Default option :             quickest advection scheme (k-j-i loop) 
    56    !!---------------------------------------------------------------------- 
    5752 
    5853   SUBROUTINE tra_adv_qck( kt, pun, pvn, pwn ) 
    … …  
    579574   END FUNCTION 
    580575 
    581 #else 
    582    !!---------------------------------------------------------------------- 
    583    !!   'key_mpp_omp' :      quickest advection (k- and j-slabs) 
    584    !!---------------------------------------------------------------------- 
    585    SUBROUTINE tra_adv_qck( kt, pun, pvn, pwn  ) 
    586       !!---------------------------------------------------------------------- 
    587       !! * Arguments 
    588       INTEGER, INTENT( in ) ::   kt             ! ocean time-step index 
    589       REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) ::  pun   ! effective ocean velocity, u_component 
    590       REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) ::  pvn   ! effective ocean velocity, v_component 
    591       REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) ::  pwn   ! effective ocean velocity, w_component 
    592       !!----------------------------------------------------------------------   
    593          IF(lwp) WRITE(numout,*) 
    594          IF(lwp) WRITE(numout,*) 'tra_adv_ qck:3st order quickest advection scheme'  
    595          IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~   Vector optimization case' 
    596          IF(lwp) WRITE(numout,*) 'WITH AUTOTASKING =>this routine doesn t exist for the moment' 
    597     IF(lwp) WRITE(numout,*) ' EMPTY ROUTINE!!!!!!'       
    598  
    599    END SUBROUTINE tra_adv_qck 
    600  
    601 #endif 
    602  
    603576   !!====================================================================== 
    604577END MODULE traadv_qck 

  • branches/dev_001_SBC/NEMO/OPA_SRC/TRA/traadv_tvd.F90

    r699 r811  
    7575      INTEGER  ::   ji, jj, jk              ! dummy loop indices 
    7676      REAL(wp) ::                        &  ! temporary scalar 
    77          ztai, ztaj, ztak,               &  !    "         "    
    78          zsai, zsaj, zsak,               &  !    "         "    
     77         ztat, zsat,                     &  !    "         "    
    7978         z_hdivn_x, z_hdivn_y, z_hdivn 
    8079      REAL(wp) ::   & 
    … …  
    158157      ! total advective trend 
    159158      DO jk = 1, jpkm1 
     159         z2dtt = z2 * rdttra(jk) 
    160160         DO jj = 2, jpjm1 
    161161            DO ji = fs_2, fs_jpim1   ! vector opt. 
    162162               zbtr = 1./ ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 
    163                ! i- j- horizontal & k- vertical advective trends 
    164                ztai = - ( ztu(ji,jj,jk) - ztu(ji-1,jj  ,jk  ) ) * zbtr 
    165                ztaj = - ( ztv(ji,jj,jk) - ztv(ji  ,jj-1,jk  ) ) * zbtr 
    166                ztak = - ( ztw(ji,jj,jk) - ztw(ji  ,jj  ,jk+1) ) * zbtr 
    167                zsai = - ( zsu(ji,jj,jk) - zsu(ji-1,jj  ,jk  ) ) * zbtr 
    168                zsaj = - ( zsv(ji,jj,jk) - zsv(ji  ,jj-1,jk  ) ) * zbtr 
    169                zsak = - ( zsw(ji,jj,jk) - zsw(ji  ,jj  ,jk+1) ) * zbtr 
    170163               ! total intermediate advective trends 
    171                zti(ji,jj,jk) = ztai + ztaj + ztak 
    172                zsi(ji,jj,jk) = zsai + zsaj + zsak 
     164               ztat = - ( ztu(ji,jj,jk) - ztu(ji-1,jj  ,jk  )   & 
     165                  &     + ztv(ji,jj,jk) - ztv(ji  ,jj-1,jk  )   & 
     166                  &     + ztw(ji,jj,jk) - ztw(ji  ,jj  ,jk+1) ) * zbtr 
     167               zsat = - ( zsu(ji,jj,jk) - zsu(ji-1,jj  ,jk  )   &  
     168                  &     + zsv(ji,jj,jk) - zsv(ji  ,jj-1,jk  )   & 
     169                  &     + zsw(ji,jj,jk) - zsw(ji  ,jj  ,jk+1) ) * zbtr 
     170               ! update and guess with monotonic sheme 
     171               ta(ji,jj,jk) =  ta(ji,jj,jk) + ztat 
     172               sa(ji,jj,jk) =  sa(ji,jj,jk) + zsat 
     173               zti (ji,jj,jk) = ( tb(ji,jj,jk) + z2dtt * ztat ) * tmask(ji,jj,jk) 
     174               zsi (ji,jj,jk) = ( sb(ji,jj,jk) + z2dtt * zsat ) * tmask(ji,jj,jk) 
    173175            END DO 
    174176         END DO 
    … …  
    222224         ! 
    223225      ENDIF 
    224  
    225       ! update and guess with monotonic sheme 
    226       DO jk = 1, jpkm1 
    227          z2dtt = z2 * rdttra(jk) 
    228          DO jj = 2, jpjm1 
    229             DO ji = fs_2, fs_jpim1   ! vector opt. 
    230                ta(ji,jj,jk) =  ta(ji,jj,jk) + zti(ji,jj,jk) 
    231                sa(ji,jj,jk) =  sa(ji,jj,jk) + zsi(ji,jj,jk) 
    232                zti (ji,jj,jk) = ( tb(ji,jj,jk) + z2dtt * zti(ji,jj,jk) ) * tmask(ji,jj,jk) 
    233                zsi (ji,jj,jk) = ( sb(ji,jj,jk) + z2dtt * zsi(ji,jj,jk) ) * tmask(ji,jj,jk) 
    234             END DO 
    235          END DO 
    236       END DO 
    237226 
    238227      ! Lateral boundary conditions on zti, zsi   (unchanged sign) 
    … …  
    290279            DO ji = fs_2, fs_jpim1   ! vector opt.   
    291280               zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 
    292                ! i- j- horizontal & k- vertical advective trends 
    293                ztai = - ( ztu(ji,jj,jk) - ztu(ji-1,jj  ,jk  )) * zbtr 
    294                ztaj = - ( ztv(ji,jj,jk) - ztv(ji  ,jj-1,jk  )) * zbtr 
    295                ztak = - ( ztw(ji,jj,jk) - ztw(ji  ,jj  ,jk+1)) * zbtr 
    296                zsai = - ( zsu(ji,jj,jk) - zsu(ji-1,jj  ,jk  )) * zbtr 
    297                zsaj = - ( zsv(ji,jj,jk) - zsv(ji  ,jj-1,jk  )) * zbtr 
    298                zsak = - ( zsw(ji,jj,jk) - zsw(ji  ,jj  ,jk+1)) * zbtr 
    299  
     281               ! total advective trends 
     282               ztat = - ( ztu(ji,jj,jk) - ztu(ji-1,jj  ,jk  )   & 
     283                  &     + ztv(ji,jj,jk) - ztv(ji  ,jj-1,jk  )   & 
     284                  &     + ztw(ji,jj,jk) - ztw(ji  ,jj  ,jk+1) ) * zbtr 
     285               zsat = - ( zsu(ji,jj,jk) - zsu(ji-1,jj  ,jk  )   & 
     286                  &     + zsv(ji,jj,jk) - zsv(ji  ,jj-1,jk  )   & 
     287                  &     + zsw(ji,jj,jk) - zsw(ji  ,jj  ,jk+1) ) * zbtr 
    300288               ! add them to the general tracer trends 
    301                ta(ji,jj,jk) = ta(ji,jj,jk) + ztai + ztaj + ztak 
    302                sa(ji,jj,jk) = sa(ji,jj,jk) + zsai + zsaj + zsak 
     289               ta(ji,jj,jk) = ta(ji,jj,jk) + ztat 
     290               sa(ji,jj,jk) = sa(ji,jj,jk) + zsat 
    303291            END DO 
    304292         END DO 
    … …  
    396384      !!       in-space based differencing for fluid 
    397385      !!---------------------------------------------------------------------- 
    398       REAL(wp), INTENT( in ) ::   prdt                               ! ??? 
     386      REAL(wp), INTENT( in ) ::   prdt   ! ??? 
    399387      REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT( inout ) ::   & 
    400388         pbef,                            & ! before field 
    … …  
    407395      INTEGER ::   ikm1 
    408396      REAL(wp), DIMENSION (jpi,jpj,jpk) ::   zbetup, zbetdo 
     397      REAL(wp), DIMENSION (jpi,jpj,jpk) ::   zbup, zbdo 
    409398      REAL(wp) ::   zpos, zneg, zbt, za, zb, zc, zbig, zrtrn, z2dtt 
     399      REAL(wp) ::   zau, zbu, zcu, zav, zbv, zcv 
     400      REAL(wp) ::   zup, zdo 
    410401      !!---------------------------------------------------------------------- 
    411402 
    412403      zbig = 1.e+40 
    413404      zrtrn = 1.e-15 
    414       zbetup(:,:,:) = 0.e0   ;   zbetdo(:,:,:) = 0.e0   
     405      zbetup(:,:,jpk) = 0.e0   ;   zbetdo(:,:,jpk) = 0.e0 
     406 
    415407 
    416408      ! Search local extrema 
    417409      ! -------------------- 
    418       ! large negative value (-zbig) inside land 
    419       pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) - zbig * ( 1.e0 - tmask(:,:,:) ) 
    420       paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) - zbig * ( 1.e0 - tmask(:,:,:) ) 
    421       ! search maximum in neighbourhood 
     410      ! max/min of pbef & paft with large negative/positive value (-/+zbig) inside land 
     411      zbup = MAX( pbef * tmask - zbig * ( 1.e0 - tmask ),   & 
     412         &        paft * tmask - zbig * ( 1.e0 - tmask )  ) 
     413      zbdo = MIN( pbef * tmask + zbig * ( 1.e0 - tmask ),   & 
     414         &        paft * tmask + zbig * ( 1.e0 - tmask )  ) 
     415 
    422416      DO jk = 1, jpkm1 
    423417         ikm1 = MAX(jk-1,1) 
    424          DO jj = 2, jpjm1 
    425             DO ji = fs_2, fs_jpim1   ! vector opt. 
    426                zbetup(ji,jj,jk) = MAX(  pbef(ji  ,jj  ,jk  ), paft(ji  ,jj  ,jk  ),   & 
    427                   &                     pbef(ji-1,jj  ,jk  ), pbef(ji+1,jj  ,jk  ),   & 
    428                   &                     paft(ji-1,jj  ,jk  ), paft(ji+1,jj  ,jk  ),   & 
    429                   &                     pbef(ji  ,jj-1,jk  ), pbef(ji  ,jj+1,jk  ),   & 
    430                   &                     paft(ji  ,jj-1,jk  ), paft(ji  ,jj+1,jk  ),   & 
    431                   &                     pbef(ji  ,jj  ,ikm1), pbef(ji  ,jj  ,jk+1),   & 
    432                   &                     paft(ji  ,jj  ,ikm1), paft(ji  ,jj  ,jk+1)  ) 
    433             END DO 
    434          END DO 
    435       END DO 
    436       ! large positive value (+zbig) inside land 
    437       pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) + zbig * ( 1.e0 - tmask(:,:,:) ) 
    438       paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) + zbig * ( 1.e0 - tmask(:,:,:) ) 
    439       ! search minimum in neighbourhood 
    440       DO jk = 1, jpkm1 
    441          ikm1 = MAX(jk-1,1) 
    442          DO jj = 2, jpjm1 
    443             DO ji = fs_2, fs_jpim1   ! vector opt. 
    444                zbetdo(ji,jj,jk) = MIN(  pbef(ji  ,jj  ,jk  ), paft(ji  ,jj  ,jk  ),   & 
    445                   &                     pbef(ji-1,jj  ,jk  ), pbef(ji+1,jj  ,jk  ),   & 
    446                   &                     paft(ji-1,jj  ,jk  ), paft(ji+1,jj  ,jk  ),   & 
    447                   &                     pbef(ji  ,jj-1,jk  ), pbef(ji  ,jj+1,jk  ),   & 
    448                   &                     paft(ji  ,jj-1,jk  ), paft(ji  ,jj+1,jk  ),   & 
    449                   &                     pbef(ji  ,jj  ,ikm1), pbef(ji  ,jj  ,jk+1),   & 
    450                   &                     paft(ji  ,jj  ,ikm1), paft(ji  ,jj  ,jk+1)  ) 
    451             END DO 
    452          END DO 
    453       END DO 
    454  
    455       ! restore masked values to zero 
    456       pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) 
    457       paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) 
    458   
    459  
    460       ! 2. Positive and negative part of fluxes and beta terms 
    461       ! ------------------------------------------------------ 
    462  
    463       DO jk = 1, jpkm1 
    464418         z2dtt = prdt * rdttra(jk) 
    465419         DO jj = 2, jpjm1 
    466420            DO ji = fs_2, fs_jpim1   ! vector opt. 
    467                ! positive & negative part of the flux 
     421 
     422               ! search maximum in neighbourhood 
     423               zup = MAX(  zbup(ji  ,jj  ,jk  ),   & 
     424                  &        zbup(ji-1,jj  ,jk  ), zbup(ji+1,jj  ,jk  ),   & 
     425                  &        zbup(ji  ,jj-1,jk  ), zbup(ji  ,jj+1,jk  ),   & 
     426                  &        zbup(ji  ,jj  ,ikm1), zbup(ji  ,jj  ,jk+1)  ) 
     427 
     428               ! search minimum in neighbourhood 
     429               zdo = MIN(  zbdo(ji  ,jj  ,jk  ),   & 
     430                  &        zbdo(ji-1,jj  ,jk  ), zbdo(ji+1,jj  ,jk  ),   & 
     431                  &        zbdo(ji  ,jj-1,jk  ), zbdo(ji  ,jj+1,jk  ),   & 
     432                  &        zbdo(ji  ,jj  ,ikm1), zbdo(ji  ,jj  ,jk+1)  ) 
     433 
     434               ! positive part of the flux 
    468435               zpos = MAX( 0., paa(ji-1,jj  ,jk  ) ) - MIN( 0., paa(ji  ,jj  ,jk  ) )   & 
    469436                  & + MAX( 0., pbb(ji  ,jj-1,jk  ) ) - MIN( 0., pbb(ji  ,jj  ,jk  ) )   & 
    470437                  & + MAX( 0., pcc(ji  ,jj  ,jk+1) ) - MIN( 0., pcc(ji  ,jj  ,jk  ) ) 
     438 
     439               ! negative part of the flux 
    471440               zneg = MAX( 0., paa(ji  ,jj  ,jk  ) ) - MIN( 0., paa(ji-1,jj  ,jk  ) )   & 
    472441                  & + MAX( 0., pbb(ji  ,jj  ,jk  ) ) - MIN( 0., pbb(ji  ,jj-1,jk  ) )   & 
    473442                  & + MAX( 0., pcc(ji  ,jj  ,jk  ) ) - MIN( 0., pcc(ji  ,jj  ,jk+1) ) 
     443 
    474444               ! up & down beta terms 
    475445               zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) / z2dtt 
    476                zbetup(ji,jj,jk) = ( zbetup(ji,jj,jk) - paft(ji,jj,jk) ) / (zpos+zrtrn) * zbt 
    477                zbetdo(ji,jj,jk) = ( paft(ji,jj,jk) - zbetdo(ji,jj,jk) ) / (zneg+zrtrn) * zbt 
     446               zbetup(ji,jj,jk) = ( zup            - paft(ji,jj,jk) ) / ( zpos + zrtrn ) * zbt 
     447               zbetdo(ji,jj,jk) = ( paft(ji,jj,jk) - zdo            ) / ( zneg + zrtrn ) * zbt 
    478448            END DO 
    479449         END DO 
    … …  
    490460         DO jj = 2, jpjm1 
    491461            DO ji = fs_2, fs_jpim1   ! vector opt. 
    492                za = MIN( 1.e0, zbetdo(ji,jj,jk), zbetup(ji+1,jj,jk) ) 
    493                zb = MIN( 1.e0, zbetup(ji,jj,jk), zbetdo(ji+1,jj,jk) ) 
    494                zc = 0.5 * ( 1.e0 + SIGN( 1.e0, paa(ji,jj,jk) ) ) 
    495                paa(ji,jj,jk) = paa(ji,jj,jk) * ( zc * za + ( 1.e0 - zc) * zb ) 
    496  
    497                za = MIN( 1.e0, zbetdo(ji,jj,jk), zbetup(ji,jj+1,jk) ) 
    498                zb = MIN( 1.e0, zbetup(ji,jj,jk), zbetdo(ji,jj+1,jk) ) 
    499                zc = 0.5 * ( 1.e0 + SIGN( 1.e0, pbb(ji,jj,jk) ) ) 
    500                pbb(ji,jj,jk) = pbb(ji,jj,jk) * ( zc * za + ( 1.e0 - zc) * zb ) 
    501             END DO 
    502          END DO 
    503       END DO 
    504  
     462               zau = MIN( 1.e0, zbetdo(ji,jj,jk), zbetup(ji+1,jj,jk) ) 
     463               zbu = MIN( 1.e0, zbetup(ji,jj,jk), zbetdo(ji+1,jj,jk) ) 
     464               zcu =       ( 0.5  + SIGN( 0.5 , paa(ji,jj,jk) ) ) 
     465               paa(ji,jj,jk) = paa(ji,jj,jk) * ( zcu * zau + ( 1.e0 - zcu) * zbu ) 
     466 
     467               zav = MIN( 1.e0, zbetdo(ji,jj,jk), zbetup(ji,jj+1,jk) ) 
     468               zbv = MIN( 1.e0, zbetup(ji,jj,jk), zbetdo(ji,jj+1,jk) ) 
     469               zcv =       ( 0.5  + SIGN( 0.5 , pbb(ji,jj,jk) ) ) 
     470               pbb(ji,jj,jk) = pbb(ji,jj,jk) * ( zcv * zav + ( 1.e0 - zcv) * zbv ) 
    505471 
    506472      ! monotonic flux in the k direction, i.e. pcc 
    507473      ! ------------------------------------------- 
    508       DO jk = 2, jpkm1 
    509          DO jj = 2, jpjm1 
    510             DO ji = fs_2, fs_jpim1   ! vector opt. 
    511  
    512                za = MIN( 1., zbetdo(ji,jj,jk), zbetup(ji,jj,jk-1) ) 
    513                zb = MIN( 1., zbetup(ji,jj,jk), zbetdo(ji,jj,jk-1) ) 
    514                zc = 0.5 * ( 1.e0 + SIGN( 1.e0, pcc(ji,jj,jk) ) ) 
    515                pcc(ji,jj,jk) = pcc(ji,jj,jk) * ( zc * za + ( 1.e0 - zc) * zb ) 
     474               za = MIN( 1., zbetdo(ji,jj,jk+1), zbetup(ji,jj,jk) ) 
     475               zb = MIN( 1., zbetup(ji,jj,jk+1), zbetdo(ji,jj,jk) ) 
     476               zc =       ( 0.5  + SIGN( 0.5 , pcc(ji,jj,jk+1) ) ) 
     477               pcc(ji,jj,jk+1) = pcc(ji,jj,jk+1) * ( zc * za + ( 1.e0 - zc) * zb ) 
    516478            END DO 
    517479         END DO 
    … …  
    521483      CALL lbc_lnk( paa, 'U', -1. )      ! changed sign 
    522484      CALL lbc_lnk( pbb, 'V', -1. )      ! changed sign 
    523       CALL lbc_lnk( pcc, 'W',  1. )      ! NO changed sign 
    524485      ! 
    525486   END SUBROUTINE nonosc 

  • branches/dev_001_SBC/NEMO/OPA_SRC/TRA/trabbc.F90

    r699 r811  
    7575      INTEGER, INTENT( in ) ::   kt   ! ocean time-step index 
    7676      !! 
    77 #if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     77#if defined key_vectopt_loop 
    7878      INTEGER ::   ji       ! dummy loop indices 
    7979#else 
    … …  
    9595      ! 
    9696      CASE ( 1:2 )                !  geothermal heat flux 
    97 #if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     97#if defined key_vectopt_loop 
    9898         DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling) 
    9999            zqgh_trd = ro0cpr * qgh_trd0(ji,1) / fse3t(ji,1,nbotlevt(ji,1) ) 

  • branches/dev_001_SBC/NEMO/OPA_SRC/TRA/trabbl.F90

    r699 r811  
    154154      ! ----------------------------------------------------------------- 
    155155      ! mbathy= number of w-level, minimum value=1 (cf dommsk.F) 
    156 #  if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     156#  if defined key_vectopt_loop 
    157157      jj = 1 
    158158      DO ji = 1, jpij   ! vector opt. (forced unrolling) 
    … …  
    167167            zsbb(ji,jj) = sb(ji,jj,ik) * tmask(ji,jj,1) 
    168168            zdep(ji,jj) = fsdept(ji,jj,ik)                ! depth of the ocean bottom T-level 
    169 #  if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     169#  if ! defined key_vectopt_loop 
    170170         END DO 
    171171#  endif 
    … …  
    173173 
    174174      IF( ln_zps ) THEN      ! partial steps correction  
    175 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     175# if defined key_vectopt_loop 
    176176         jj = 1 
    177177         DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    … …  
    188188               zahu(ji,jj) = atrbbl * e2u(ji,jj) * ze3u / e1u(ji,jj) * umask(ji,jj,1) 
    189189               zahv(ji,jj) = atrbbl * e1v(ji,jj) * ze3v / e2v(ji,jj) * vmask(ji,jj,1) 
    190 # if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     190# if ! defined key_vectopt_loop 
    191191            END DO 
    192192# endif 
    193193         END DO 
    194194      ELSE                    ! z-coordinate - full steps or s-coordinate 
    195 #   if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     195#   if defined key_vectopt_loop 
    196196         jj = 1 
    197197         DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    … …  
    204204               zahu(ji,jj) = atrbbl * e2u(ji,jj) * fse3u(ji,jj,iku) / e1u(ji,jj) * umask(ji,jj,1) 
    205205               zahv(ji,jj) = atrbbl * e1v(ji,jj) * fse3v(ji,jj,ikv) / e2v(ji,jj) * vmask(ji,jj,1) 
    206 #   if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     206#   if ! defined key_vectopt_loop 
    207207            END DO 
    208208#   endif 
    … …  
    219219      CASE ( 0 )                 ! Jackett and McDougall (1994) formulation 
    220220 
    221 #  if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     221#  if defined key_vectopt_loop 
    222222      jj = 1 
    223223      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    … …  
    238238            zsign = SIGN( 0.5, - zgdrho * ( zdep(ji+1,jj) - zdep(ji,jj) ) ) 
    239239            zki(ji,jj) = ( 0.5 - zsign ) * zahu(ji,jj) 
    240 #  if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
    241          END DO 
    242 #  endif 
    243       END DO 
    244  
    245 #  if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     240#  if ! defined key_vectopt_loop 
     241         END DO 
     242#  endif 
     243      END DO 
     244 
     245#  if defined key_vectopt_loop 
    246246      jj = 1 
    247247      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    … …  
    262262            zsign = sign( 0.5, -zgdrho * ( zdep(ji,jj+1) - zdep(ji,jj) ) ) 
    263263            zkj(ji,jj) = ( 0.5 - zsign ) * zahv(ji,jj) 
    264 #  if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     264#  if ! defined key_vectopt_loop 
    265265         END DO 
    266266#  endif 
    … …  
    269269      CASE ( 1 )               ! Linear formulation function of temperature only 
    270270                               !  
    271 #  if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     271#  if defined key_vectopt_loop 
    272272      jj = 1 
    273273      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    … …  
    281281            zsign = SIGN( 0.5, - zgdrho * ( zdep(ji+1,jj) - zdep(ji,jj) ) ) 
    282282            zki(ji,jj) = ( 0.5 - zsign ) * zahu(ji,jj) 
    283 #  if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
    284          END DO 
    285 #  endif 
    286       END DO 
    287  
    288 #  if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     283#  if ! defined key_vectopt_loop 
     284         END DO 
     285#  endif 
     286      END DO 
     287 
     288#  if defined key_vectopt_loop 
    289289      jj = 1 
    290290      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    … …  
    298298            zsign = sign( 0.5, -zgdrho * ( zdep(ji,jj+1) - zdep(ji,jj) ) ) 
    299299            zkj(ji,jj) = ( 0.5 - zsign ) * zahv(ji,jj) 
    300 #  if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     300#  if ! defined key_vectopt_loop 
    301301         END DO 
    302302#  endif 
    … …  
    305305      CASE ( 2 )               ! Linear formulation function of temperature and salinity 
    306306 
    307 #  if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     307#  if defined key_vectopt_loop 
    308308      jj = 1 
    309309      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    … …  
    318318            zsign = SIGN( 0.5, - zgdrho * ( zdep(ji+1,jj) - zdep(ji,jj) ) ) 
    319319            zki(ji,jj) = ( 0.5 - zsign ) * zahu(ji,jj) 
    320 #  if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
    321          END DO 
    322 #  endif 
    323       END DO 
    324  
    325 #  if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     320#  if ! defined key_vectopt_loop 
     321         END DO 
     322#  endif 
     323      END DO 
     324 
     325#  if defined key_vectopt_loop 
    326326      jj = 1 
    327327      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    … …  
    336336            zsign = sign( 0.5, -zgdrho * ( zdep(ji,jj+1) - zdep(ji,jj) ) ) 
    337337            zkj(ji,jj) = ( 0.5 - zsign ) * zahv(ji,jj) 
    338 #  if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     338#  if ! defined key_vectopt_loop 
    339339         END DO 
    340340#  endif 
    … …  
    352352 
    353353      ! first derivative (gradient) 
    354 #  if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     354#  if defined key_vectopt_loop 
    355355      jj = 1 
    356356      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    … …  
    364364            zky(ji,jj) = zkj(ji,jj) * ( ztbb(ji,jj+1) - ztbb(ji,jj) ) 
    365365            zkw(ji,jj) = zkj(ji,jj) * ( zsbb(ji,jj+1) - zsbb(ji,jj) ) 
    366 #  if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     366#  if ! defined key_vectopt_loop 
    367367         END DO 
    368368#  endif 
    … …  
    402402 
    403403      ! second derivative (divergence) and add to the general tracer trend 
    404 #  if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     404#  if defined key_vectopt_loop 
    405405      jj = 1 
    406406      DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling) 
    … …  
    417417            ta(ji,jj,ik) = ta(ji,jj,ik) + zta 
    418418            sa(ji,jj,ik) = sa(ji,jj,ik) + zsa 
    419 #  if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     419#  if ! defined key_vectopt_loop 
    420420         END DO 
    421421#  endif 

  • branches/dev_001_SBC/NEMO/OPA_SRC/TRA/trabbl_adv.h90

    r699 r811  
    9999      ! mbathy= number of w-level, minimum value=1 (cf dommsk.F) 
    100100 
    101 #if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     101#if defined key_vectopt_loop 
    102102      jj = 1 
    103103      DO ji = 1, jpij   ! vector opt. (forced unrolling) 
    … …  
    115115            zunb(ji,jj) = un(ji,jj,mbku(ji,jj))  
    116116            zvnb(ji,jj) = vn(ji,jj,mbkv(ji,jj))  
    117 #if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     117#if ! defined key_vectopt_loop 
    118118         END DO 
    119119#endif 
    … …  
    229229      IF( ln_zps ) THEN     ! partial steps correction    
    230230       
    231 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     231# if defined key_vectopt_loop 
    232232         jj = 1 
    233233         DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    … …  
    249249                  v_bbl(ji,jj,ikv) = zalphay(ji,jj) * vn(ji,jj,ikv) * ze3v / fse3v(ji,jj,ikv)        
    250250               ENDIF 
    251 # if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     251# if ! defined key_vectopt_loop 
    252252            END DO 
    253253# endif 
    … …  
    259259      ELSE       ! if not partial step loop over the whole domain no lbc call 
    260260 
    261 #if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     261#if defined key_vectopt_loop 
    262262         jj = 1 
    263263         DO ji = 1, jpij   ! vector opt. (forced unrolling) 
    … …  
    272272                  v_bbl(ji,jj,ikv) = zalphay(ji,jj) * vn(ji,jj,ikv)        
    273273               ENDIF 
    274 #if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     274#if ! defined key_vectopt_loop 
    275275            END DO 
    276276#endif 
    … …  
    284284      ! ... Second order centered tracer flux at u and v-points 
    285285 
    286 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     286# if defined key_vectopt_loop 
    287287      jj = 1 
    288288      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    … …  
    309309            zwz(ji,jj) = ( ( zfvj + ABS( zfvj ) ) * zsbb(ji  ,jj  )   & 
    310310               &          +( zfvj - ABS( zfvj ) ) * zsbb(ji  ,jj+1) ) * 0.5 
    311 #if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     311#if ! defined key_vectopt_loop 
    312312         END DO 
    313313#endif 
    314314        END DO 
    315 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     315# if defined key_vectopt_loop 
    316316      jj = 1 
    317317      DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling) 
    … …  
    331331            ta(ji,jj,ik) = ta(ji,jj,ik) + zta 
    332332            sa(ji,jj,ik) = sa(ji,jj,ik) + zsa 
    333 #if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     333#if ! defined key_vectopt_loop 
    334334         END DO 
    335335#endif 
    … …  
    365365      IF( ln_zps ) THEN 
    366366      
    367 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     367# if defined key_vectopt_loop 
    368368         jj = 1 
    369369         DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    … …  
    383383               zwu(ji,jj) = zalphax(ji,jj) * e2u(ji,jj) * ze3u   
    384384               zwv(ji,jj) = zalphay(ji,jj) * e1v(ji,jj) * ze3v 
    385 #if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     385#if ! defined key_vectopt_loop 
    386386            END DO 
    387387#endif 
    … …  
    390390      ELSE 
    391391 
    392 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     392# if defined key_vectopt_loop 
    393393         jj = 1 
    394394         DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    … …  
    401401               zwu(ji,jj) = zalphax(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,iku)  
    402402               zwv(ji,jj) = zalphay(ji,jj) * e1v(ji,jj) * fse3v(ji,jj,ikv)  
    403 #if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     403#if ! defined key_vectopt_loop 
    404404            END DO 
    405405#endif 
    … …  
    409409  
    410410 
    411 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     411# if defined key_vectopt_loop 
    412412      jj = 1 
    413413      DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling) 
    … …  
    425425               &   ) / zbt 
    426426 
    427 # if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     427# if ! defined key_vectopt_loop 
    428428         END DO 
    429429# endif 

  • branches/dev_001_SBC/NEMO/OPA_SRC/TRA/traldf_iso.F90

    r699 r811  
    1818   !!                  and with the vertical part of 
    1919   !!                  the isopycnal or geopotential s-coord. operator  
    20    !!                  vector optimization, use k-j-i loops. 
    2120   !!---------------------------------------------------------------------- 
    2221   USE oce             ! ocean dynamics and active tracers 
    … …  
    151150       
    152151!CDIR PARALLEL DO PRIVATE( zdk1t, zdk1s, zftu, zfsu )  
    153 !$OMP PARALLEL DO PRIVATE( zdk1t, zdk1s, zftu, zfsu ) 
    154152      !                                                ! =============== 
    155153      DO jk = 1, jpkm1                                 ! Horizontal slab 

  • branches/dev_001_SBC/NEMO/OPA_SRC/TRA/tranxt.F90

    r708 r811  
    168168      !                                                ! =============== 
    169169      ! Update tracers on open boundaries. 
    170       CALL Agrif_tra( kt ) 
     170      CALL Agrif_tra 
    171171      !                                                ! =============== 
    172172      DO jk = 1, jpkm1                                 ! Horizontal slab 

  • branches/dev_001_SBC/NEMO/OPA_SRC/TRA/trazdf.F90

    r708 r811  
    1919   USE trazdf_exp      ! vertical diffusion: explicit (tra_zdf_exp     routine) 
    2020   USE trazdf_imp      ! vertical diffusion: implicit (tra_zdf_imp     routine) 
    21    USE trazdf_imp_jki  ! vertical diffusion  implicit (tra_zdf_imp_jki routine) 
    2221 
    2322   USE ldftra_oce      ! ocean active tracers: lateral physics 
    … …  
    9089         CALL prt_ctl( tab3d_1=ta, clinfo1=' zdf1 - Ta: ', mask1=tmask,               & 
    9190            &          tab3d_2=sa, clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' ) 
    92          CALL tra_zdf_imp_jki( kt, r2dt ) 
    93          CALL prt_ctl( tab3d_1=ta, clinfo1=' zdf2 - Ta: ', mask1=tmask,               & 
    94             &          tab3d_2=sa, clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' ) 
    9591 
    9692      CASE ( 0 )                                       ! explicit scheme 
    … …  
    10298         ENDIF 
    10399 
    104       CASE ( 1 )                                       ! implicit scheme (k-j-i loop) 
     100      CASE ( 1 )                                       ! implicit scheme  
    105101         CALL tra_zdf_imp    ( kt, r2dt ) 
    106          IF( l_trdtra )   THEN                         ! save the vertical diffusive trends for further diagnostics 
    107             DO jk = 1, jpkm1 
    108                ztrdt(:,:,jk) = ( ( ta(:,:,jk) - tb(:,:,jk) ) / r2dt(jk) ) - ztrdt(:,:,jk) 
    109                ztrds(:,:,jk) = ( ( sa(:,:,jk) - sb(:,:,jk) ) / r2dt(jk) ) - ztrds(:,:,jk) 
    110             END DO 
    111             CALL trd_mod( ztrdt, ztrds, jptra_trd_zdf, 'TRA', kt ) 
    112          ENDIF 
    113  
    114       CASE ( 2 )                                       ! implicit scheme (j-k-i loop) 
    115          CALL tra_zdf_imp_jki( kt, r2dt ) 
    116102         IF( l_trdtra )   THEN                         ! save the vertical diffusive trends for further diagnostics 
    117103            DO jk = 1, jpkm1 
    … …  
    137123      !! ** Purpose :   Choose the vertical mixing scheme 
    138124      !! 
    139       !! ** Method  :   Set nzdf from ln_zdfexp and 'key_mpp_omp'. 
     125      !! ** Method  :   Set nzdf from ln_zdfexp 
    140126      !!      nzdf = 0   explicit (time-splitting) scheme (ln_zdfexp=T) 
    141127      !!           = 1   implicit (euler backward) scheme (ln_zdfexp=F) 
    142       !!           = 2   implicit (euler backward) scheme with j-k-i loops 
    143       !!                 (ln_zdfexp=T and 'key_mpp_omp') 
    144128      !!      NB: rotation of lateral mixing operator or TKE or KPP scheme, 
    145129      !!      the implicit scheme is required. 
    … …  
    169153      ENDIF 
    170154 
    171       ! NEC autotasking / OpenMP 
    172 #if defined key_mpp_omp 
    173       IF( nzdf == 1 )   nzdf = 2                       ! j-k-i loop 
    174 #endif 
    175  
    176155      ! Test: esopa 
    177156      IF( lk_esopa )    nzdf = -1                      ! All schemes used 
    … …  
    184163         IF( nzdf ==  0 )   WRITE(numout,*) '              Explicit time-splitting scheme' 
    185164         IF( nzdf ==  1 )   WRITE(numout,*) '              Implicit (euler backward) scheme' 
    186          IF( nzdf ==  2 )   WRITE(numout,*) '              Implicit (euler backward) scheme with j-k-i loops' 
    187165      ENDIF 
    188166 

  • branches/dev_001_SBC/NEMO/OPA_SRC/TRA/trazdf_exp.F90

    r699 r811  
    8585      IF( kt == nit000 ) THEN 
    8686         IF(lwp) WRITE(numout,*) 
    87          IF(lwp) WRITE(numout,*) 'tra_zdf_exp : explicit vertical mixing (j-k-i loops)' 
     87         IF(lwp) WRITE(numout,*) 'tra_zdf_exp : explicit vertical mixing' 
    8888         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 
    8989      ENDIF 

  • branches/dev_001_SBC/NEMO/OPA_SRC/TRA/trazdf_imp.F90

    r699 r811  
    1717   !!   tra_zdf_imp : Update the tracer trend with the diagonal vertical   
    1818   !!                 part of the mixing tensor. 
    19    !!                 Vector optimization, use k-j-i loops. 
    2019   !!---------------------------------------------------------------------- 
    2120   !! * Modules used 
    … …  
    108107      IF( kt == nit000 ) THEN 
    109108         IF(lwp)WRITE(numout,*) 
    110          IF(lwp)WRITE(numout,*) 'tra_zdf_imp : implicit vertical mixing (k-j-i loops)' 
     109         IF(lwp)WRITE(numout,*) 'tra_zdf_imp : implicit vertical mixing' 
    111110         IF(lwp)WRITE(numout,*) '~~~~~~~~~~~ ' 
    112111         zavi = 0.e0      ! avoid warning at compilation phase when lk_ldfslp=F 

  • branches/dev_001_SBC/NEMO/OPA_SRC/TRA/zpshde.F90

    r699 r811  
    130130 
    131131      ! Interpolation of T and S at the last ocean level 
    132 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     132# if defined key_vectopt_loop 
    133133         jj = 1 
    134134         DO ji = 1, jpij-jpi   ! vector opt. (forced unrolled) 
    … …  
    193193               pgsv(ji,jj) = vmask(ji,jj,1) * ( psal(ji,jj+1,ikv) - zsj(ji,jj) ) 
    194194            ENDIF 
    195 # if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     195# if ! defined key_vectopt_loop 
    196196         END DO 
    197197# endif 
    … …  
    205205 
    206206      ! Gradient of density at the last level  
    207 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     207# if defined key_vectopt_loop 
    208208         jj = 1 
    209209         DO ji = 1, jpij-jpi   ! vector opt. (forced unrolled) 
    … …  
    226226               pgrv(ji,jj) = vmask(ji,jj,1) * ( prd(ji,jj+1,ikv) - zrj(ji,jj) ) 
    227227            ENDIF 
    228 # if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     228# if ! defined key_vectopt_loop 
    229229         END DO 
    230230# endif 

  • branches/dev_001_SBC/NEMO/OPA_SRC/TRD/trdmld_rst.F90

    r699 r811  
    4646      !!-------------------------------------------------------------------------------- 
    4747 
    48       IF( kt == nitrst-1 ) THEN 
    49          IF( nitrst > 1.0e9 ) THEN    
    50             WRITE(clkt,*) nitrst 
    51          ELSE 
    52             WRITE(clkt,'(i8.8)') nitrst 
     48      ! to get better performances with NetCDF format: 
     49      ! we open and define the ocean restart_mld file one time step before writing the data (-> at nitrst - 1) 
     50      ! except if we write ocean restart_mld files every time step or if an ocean restart_mld file was writen at nitend - 1 
     51      IF( kt == nitrst - 1 .OR. nstock == 1 .OR. ( kt == nitend .AND. MOD( nitend - 1, nstock ) == 0 ) ) THEN 
     52         ! beware of the format used to write kt (default is i8.8, that should be large enough...) 
     53         IF( nitrst > 999999999 ) THEN   ;   WRITE(clkt, *       ) nitrst 
     54         ELSE                            ;   WRITE(clkt, '(i8.8)') nitrst 
    5355         ENDIF 
     56         ! create the file 
    5457         clname = TRIM(cexper)//"_"//TRIM(ADJUSTL(clkt))//"_restart_mld" 
    55          IF(lwp) WRITE(numout,*) '             open ocean restart_mld NetCDF file: '//clname 
     58         IF(lwp) THEN 
     59            WRITE(numout,*) 
     60            SELECT CASE ( jprstlib ) 
     61            CASE ( jprstdimg )   ;   WRITE(numout,*) '             open ocean restart_mld binary file: '//clname 
     62            CASE DEFAULT         ;   WRITE(numout,*) '             open ocean restart_mld NetCDF file: '//clname 
     63            END SELECT 
     64            IF( kt == nitrst - 1 ) THEN   ;   WRITE(numout,*) '             kt = nitrst - 1 = ', kt,' date= ', ndastp 
     65            ELSE                          ;   WRITE(numout,*) '             kt = '             , kt,' date= ', ndastp 
     66            ENDIF 
     67         ENDIF 
     68 
    5669         CALL iom_open( clname, nummldw, ldwrt = .TRUE., kiolib = jprstlib ) 
    5770      ENDIF 
    … …  
    5972      IF( kt == nitrst .AND. lwp ) THEN 
    6073         WRITE(numout,*) 
    61          WRITE(numout,*) 'trdmld_rst: output for ML diags. restart, with trd_mld_rst_write routine' 
     74         WRITE(numout,*) 'trdmld_rst: output for ML diags. restart, with trd_mld_rst_write routine kt =', kt 
    6275         WRITE(numout,*) '~~~~~~~~~~' 
    6376         WRITE(numout,*) 
    … …  
    8194         CALL iom_rstput( kt, nitrst, nummldw, 'tml_sumb'        , tml_sumb        ) 
    8295         DO jk = 1, jpltrd 
    83             IF( jk < 10 )   THEN 
    84                WRITE(charout,FMT="('tmltrd_csum_ub_', I1)") jk 
    85             ELSE 
    86                WRITE(charout,FMT="('tmltrd_csum_ub_', I2)") jk 
     96            IF( jk < 10 ) THEN   ;   WRITE(charout,FMT="('tmltrd_csum_ub_', I1)") jk 
     97            ELSE                 ;   WRITE(charout,FMT="('tmltrd_csum_ub_', I2)") jk 
    8798            ENDIF 
    8899            CALL iom_rstput( kt, nitrst, nummldw, charout,  tmltrd_csum_ub(:,:,jk) ) 
    … …  
    94105         CALL iom_rstput( kt, nitrst, nummldw, 'sml_sumb'        , sml_sumb        ) 
    95106         DO jk = 1, jpltrd 
    96             IF( jk < 10 )   THEN 
    97                WRITE(charout,FMT="('smltrd_csum_ub_', I1)") jk 
    98             ELSE 
    99                WRITE(charout,FMT="('smltrd_csum_ub_', I2)") jk 
     107            IF( jk < 10 ) THEN   ;   WRITE(charout,FMT="('smltrd_csum_ub_', I1)") jk 
     108            ELSE                 ;   WRITE(charout,FMT="('smltrd_csum_ub_', I2)") jk 
    100109            ENDIF 
    101110            CALL iom_rstput( kt, nitrst, nummldw, charout , smltrd_csum_ub(:,:,jk) ) 

  • branches/dev_001_SBC/NEMO/OPA_SRC/ZDF/zdfbfr.F90

    r699 r811  
    8080 
    8181      CASE( 0 )                 ! no-slip boundary condition 
    82 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     82# if defined key_vectopt_loop 
    8383         jj = 1 
    8484         DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling) 
    … …  
    9393               avmu(ji,jj,ikbu) = 2. * avmu(ji,jj,ikbum1) 
    9494               avmv(ji,jj,ikbv) = 2. * avmv(ji,jj,ikbvm1) 
    95 # if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     95# if ! defined key_vectopt_loop 
    9696            END DO 
    9797# endif 
    … …  
    9999 
    100100      CASE( 1 )                 ! linear botton friction 
    101 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     101# if defined key_vectopt_loop 
    102102         jj = 1 
    103103         DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling) 
    … …  
    110110               avmu(ji,jj,ikbu) = bfri1 * fse3uw(ji,jj,ikbu) 
    111111               avmv(ji,jj,ikbv) = bfri1 * fse3vw(ji,jj,ikbv) 
    112 # if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     112# if ! defined key_vectopt_loop 
    113113            END DO 
    114114# endif 
    … …  
    116116 
    117117      CASE( 2 )                 ! quadratic botton friction 
    118 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     118# if defined key_vectopt_loop 
    119119         jj = 1 
    120120!CDIR NOVERRCHK 
    … …  
    142142               avmu(ji,jj,ikbu) = bfri2 * zecu * fse3uw(ji,jj,ikbu) 
    143143               avmv(ji,jj,ikbv) = bfri2 * zecv * fse3vw(ji,jj,ikbv) 
    144 # if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     144# if ! defined key_vectopt_loop 
    145145            END DO 
    146146# endif 
    … …  
    148148 
    149149      CASE( 3 )                 ! free-slip boundary condition 
    150 # if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     150# if defined key_vectopt_loop 
    151151         jj = 1 
    152152         DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling) 
    … …  
    159159               avmu(ji,jj,ikbu) = 0.e0 
    160160               avmv(ji,jj,ikbv) = 0.e0 
    161 # if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     161# if ! defined key_vectopt_loop 
    162162            END DO 
    163163# endif 

  • branches/dev_001_SBC/NEMO/OPA_SRC/ZDF/zdfevd.F90

    r699 r811  
    7979         DO jk = 1, jpkm1                                 ! Horizontal slab 
    8080            !                                             ! =============== 
    81 #   if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     81#   if defined key_vectopt_loop 
    8282!!!         WHERE( rn2(:,:,jk) <= -1.e-12 ) avt(:,:,jk) = tmask(:,:,jk) * avevd   ! agissant sur T SEUL! 
    8383            jj = 1                     ! big loop forced 
    … …  
    129129            !                                             ! =============== 
    130130!!!         WHERE( rn2(:,:,jk) <= -1.e-12 ) avt(:,:,jk) = tmask(:,:,jk) * avevd   ! agissant sur T SEUL!  
    131 #   if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     131#   if defined key_vectopt_loop 
    132132            jj = 1                     ! big loop forced 
    133133            DO ji = 1, jpij    

  • branches/dev_001_SBC/NEMO/OPA_SRC/ZDF/zdfmxl.F90

    r699 r811  
    4444 
    4545CONTAINS 
    46  
    47 # if defined key_mpp_omp 
    48    !!---------------------------------------------------------------------- 
    49    !!   'key_mpp_omp'                                   j-k-i loop (j-slab) 
    50    !!---------------------------------------------------------------------- 
    51  
    52    SUBROUTINE zdf_mxl( kt ) 
    53       !!---------------------------------------------------------------------- 
    54       !!                    ***  ROUTINE zdfmxl  *** 
    55       !!                    
    56       !! ** Purpose :   Compute the turbocline depth and the mixed layer depth 
    57       !!      with a density criteria. 
    58       !! 
    59       !! ** Method  :   The turbocline depth is the depth at which the vertical  
    60       !!      eddy diffusivity coefficient (resulting from the vertical physics 
    61       !!      alone, not the isopycnal part, see trazdf.F) fall below a given 
    62       !!      value defined locally (avt_c here taken equal to 5 cm/s2) 
    63       !! 
    64       !! ** Action  : 
    65       !! 
    66       !!---------------------------------------------------------------------- 
    67       !! * Arguments 
    68       INTEGER, INTENT( in ) ::   kt         ! ocean time-step index 
    69  
    70       !! * Local declarations 
    71       INTEGER ::   ji, jj, jk     ! dummy loop indices 
    72       INTEGER ::   ik             ! temporary integer 
    73       INTEGER, DIMENSION(jpi,jpj) ::   & 
    74          imld                     ! temporary workspace 
    75       !!---------------------------------------------------------------------- 
    76  
    77       IF( kt == nit000 ) THEN 
    78          IF(lwp) WRITE(numout,*) 
    79          IF(lwp) WRITE(numout,*) 'zdf_mxl : mixed layer depth, j-k-i loops' 
    80          IF(lwp) WRITE(numout,*) '~~~~~~~' 
    81       ENDIF 
    82  
    83       !                                                ! =============== 
    84       DO jj = 1, jpj                                   !  Vertical slab 
    85          !                                             ! =============== 
    86  
    87          ! 1. Turbocline depth 
    88          ! ------------------- 
    89          ! last w-level at which avt<avt_c (starting from the bottom jk=jpk) 
    90          ! (since avt(.,.,jpk)=0, we have jpk=< imld =< 2 ) 
    91          DO jk = jpk, 2, -1 
    92             DO ji = 1, jpi 
    93                IF( avt(ji,jj,jk) < avt_c ) imld(ji,jj) = jk  
    94             END DO 
    95          END DO 
    96  
    97          ! Turbocline depth and sub-turbocline temperature 
    98          DO ji = 1, jpi 
    99             ik = imld(ji,jj) 
    100             hmld (ji,jj) = fsdepw(ji,jj,ik) * tmask(ji,jj,1) 
    101          END DO 
    102  
    103 !!gm idea 
    104 !!    
    105 !!gm     DO jk = jpk, 2, -1 
    106 !!gm        DO ji = 1, jpi 
    107 !!gm           IF( avt(ji,jj,jk) < avt_c ) hmld(ji,jj) = fsdepw(ji,jj,jk) * tmask(ji,jj,1) 
    108 !!gm        END DO 
    109 !!gm     END DO 
    110 !!gm 
    111  
    112          ! 2. Mixed layer depth 
    113          ! -------------------- 
    114          ! Initialization to the number of w ocean point mbathy 
    115          nmln(:,jj) = mbathy(:,jj) 
    116  
    117          ! Last w-level at which rhop>=rho surf+rho_c (starting from jpk-1) 
    118          ! (rhop defined at t-point, thus jk-1 for w-level just above) 
    119          DO jk = jpkm1, 2, -1 
    120             DO ji = 1, jpi 
    121                IF( rhop(ji,jj,jk) > rhop(ji,jj,1) + rho_c )   nmln(ji,jj) = jk 
    122             END DO 
    123          END DO 
    124  
    125          ! Mixed layer depth 
    126          DO ji = 1, jpi 
    127             ik = nmln(ji,jj) 
    128             hmlp (ji,jj) = fsdepw(ji,jj,ik) * tmask(ji,jj,1) 
    129             hmlpt(ji,jj) = fsdept(ji,jj,ik-1) 
    130          END DO 
    131          !                                             ! =============== 
    132       END DO                                           !   End of slab 
    133       !                                                ! =============== 
    134  
    135       IF(ln_ctl)   CALL prt_ctl( tab2d_1=REAL(nmln,wp), clinfo1=' nmln : ', tab2d_2=hmld, clinfo2=' hmld : ', ovlap=1 ) 
    136  
    137    END SUBROUTINE zdf_mxl 
    138  
    139 # else 
    140    !!---------------------------------------------------------------------- 
    141    !!   Default option :                                         k-j-i loop 
    142    !!---------------------------------------------------------------------- 
    14346 
    14447   SUBROUTINE zdf_mxl( kt ) 
    … …  
    237140 
    238141   END SUBROUTINE zdf_mxl 
    239 #endif 
    240142 
    241143   !!====================================================================== 

  • branches/dev_001_SBC/NEMO/OPA_SRC/ZDF/zdftke.F90

    r762 r811  
    4343 
    4444   PUBLIC   zdf_tke        ! routine called in step module 
    45    PUBLIC   zdf_tke_init   ! routine also called in zdftke_jki module 
    46    PUBLIC   tke_rst        ! routine also called in zdftke_jki module 
    4745 
    4846   LOGICAL , PUBLIC, PARAMETER              ::   lk_zdftke = .TRUE.  !: TKE vertical mixing flag 

  • branches/dev_001_SBC/NEMO/OPA_SRC/eosbn2.F90

    r703 r811  
    444444         DO jj = 1, jpjm1 
    445445!CDIR NOVERRCHK 
    446 #if defined key_mpp_omp 
    447             DO ji = 1, jpim1 
    448 #else 
    449446            DO ji = 1, fs_jpim1   ! vector opt. 
    450 #endif 
    451447               zws(ji,jj) = SQRT( ABS( psal(ji,jj) ) ) 
    452448            END DO 
    … …  
    455451         DO jj = 1, jpjm1                                 ! Horizontal slab 
    456452            !                                             ! =============== 
    457 #if defined key_mpp_omp 
    458             DO ji = 1, jpim1 
    459 #else 
    460453            DO ji = 1, fs_jpim1   ! vector opt. 
    461 #endif 
    462454               zmask = tmask(ji,jj,1)      ! land/sea bottom mask = surf. mask 
    463455 
    … …  
    508500         DO jj = 1, jpjm1                                 ! Horizontal slab 
    509501            !                                             ! =============== 
    510 #if defined key_mpp_omp 
    511             DO ji = 1, jpim1 
    512 #else 
    513502            DO ji = 1, fs_jpim1   ! vector opt. 
    514 #endif 
    515503               prd(ji,jj) = ( 0.0285 - ralpha * ptem(ji,jj) ) * tmask(ji,jj,1) 
    516504            END DO 
    … …  
    524512         DO jj = 1, jpjm1                                 ! Horizontal slab 
    525513            !                                             ! =============== 
    526 #if defined key_mpp_omp 
    527             DO ji = 1, jpim1 
    528 #else 
    529514            DO ji = 1, fs_jpim1   ! vector opt. 
    530 #endif 
    531515               prd(ji,jj) = ( rbeta * psal(ji,jj) - ralpha * ptem(ji,jj) ) * tmask(ji,jj,1)  
    532516            END DO 

  • branches/dev_001_SBC/NEMO/OPA_SRC/par_oce.F90

    r699 r811  
    219219#endif 
    220220 
    221 #if defined key_mpp_omp 
    222    LOGICAL, PUBLIC, PARAMETER ::   lk_jki = .TRUE.   !: j-k-i loop flag 
    223 #else 
    224    LOGICAL, PUBLIC, PARAMETER ::   lk_jki = .FALSE.  !: k-j-i loop flag 
    225 #endif 
    226  
    227221   !!====================================================================== 
    228222END MODULE par_oce 

  • branches/dev_001_SBC/NEMO/OPA_SRC/restart.F90

    r762 r811  
    6363      !!---------------------------------------------------------------------- 
    6464      ! 
    65       IF( kt == nit000 ) THEN   ! default initialization, to do: should be read in the namelist... 
    66          nitrst = nitend        ! to do: should be read in the namelist in a cleaver way... 
    67          lrst_oce = .FALSE. 
    68       ENDIF 
    69  
    70       IF    ( kt == nitrst-1 .AND. lrst_oce         ) THEN 
    71          CALL ctl_stop( 'rst_opn: we cannot create an ocean restart at every time step',    & 
    72             &           'if the run has more than one time step!!!' ) 
    73          numrow = 0 
    74       ELSEIF( kt == nitrst-1 .OR.  nitend == nit000 ) THEN   ! beware if model runs only one time step 
    75          ! beware of the format used to write kt (default is i8.8, that should be large enough) 
    76          IF( nitrst > 1.0e9 ) THEN    
    77             WRITE(clkt,*) nitrst 
    78          ELSE 
    79             WRITE(clkt,'(i8.8)') nitrst 
     65      IF( kt == nit000 ) THEN   ! default definitions 
     66         lrst_oce = .FALSE.    
     67         nitrst = nitend 
     68      ENDIF 
     69      IF( MOD( kt - 1, nstock ) == 0 ) THEN    
     70         ! we use kt - 1 and not kt - nit000 to keep the same periodicity from the beginning of the experiment 
     71         nitrst = kt + nstock - 1                  ! define the next value of nitrst for restart writing 
     72         IF( nitrst > nitend )   nitrst = nitend   ! make sure we write a restart at the end of the run 
     73      ENDIF 
     74 
     75      ! to get better performances with NetCDF format: 
     76      ! we open and define the ocean restart file one time step before writing the data (-> at nitrst - 1) 
     77      ! except if we write ocean restart files every time step or if an ocean restart file was writen at nitend - 1 
     78      IF( kt == nitrst - 1 .OR. nstock == 1 .OR. ( kt == nitend .AND. .NOT. lrst_oce ) ) THEN 
     79         ! beware of the format used to write kt (default is i8.8, that should be large enough...) 
     80         IF( nitrst > 999999999 ) THEN   ;   WRITE(clkt, *       ) nitrst 
     81         ELSE                            ;   WRITE(clkt, '(i8.8)') nitrst 
    8082         ENDIF 
    8183         ! create the file 
    … …  
    8486            WRITE(numout,*) 
    8587            SELECT CASE ( jprstlib ) 
    86             CASE ( jpnf90 ) 
    87                WRITE(numout,*) '             open ocean restart.output NetCDF file: '//clname 
    88             CASE ( jprstdimg ) 
    89                WRITE(numout,*) '             open ocean restart.output binary file: '//clname 
     88            CASE ( jprstdimg )   ;   WRITE(numout,*) '             open ocean restart binary file: '//clname 
     89            CASE DEFAULT         ;   WRITE(numout,*) '             open ocean restart NetCDF file: '//clname 
    9090            END SELECT 
    91             IF( kt == nitrst-1 ) THEN 
    92                WRITE(numout,*) '             kt = nitrst - 1 = ', kt,' date= ', ndastp 
    93             ELSE 
    94                WRITE(numout,*) '             kt = ', kt,' date= ', ndastp 
     91            IF( kt == nitrst - 1 ) THEN   ;   WRITE(numout,*) '             kt = nitrst - 1 = ', kt,' date= ', ndastp 
     92            ELSE                          ;   WRITE(numout,*) '             kt = '             , kt,' date= ', ndastp 
    9593            ENDIF 
    9694         ENDIF 
    … …  
    114112      INTEGER, INTENT(in) ::   kt   ! ocean time-step 
    115113      !!---------------------------------------------------------------------- 
     114 
     115      IF( kt == nitrst ) THEN 
     116         IF(lwp) WRITE(numout,*) 
     117         IF(lwp) WRITE(numout,*) 'rst_write : write oce restart file  kt =', kt 
     118         IF(lwp) WRITE(numout,*) '~~~~~~~'          
     119      ENDIF 
    116120 
    117121      ! calendar control 

  • branches/dev_001_SBC/NEMO/OPA_SRC/step.F90

    r709 r811  
    8484   USE zdfbfr          ! bottom friction                  (zdf_bfr routine) 
    8585   USE zdftke          ! TKE vertical mixing              (zdf_tke routine) 
    86    USE zdftke_jki      ! TKE vertical mixing              (zdf_tke routine) 
    8786   USE zdfkpp          ! KPP vertical mixing              (zdf_kpp routine) 
    8887   USE zdfddm          ! double diffusion mixing          (zdf_ddm routine) 
    … …  
    210209      !                                                     ! Vertical eddy viscosity and diffusivity coefficients 
    211210      IF( lk_zdfric )   CALL zdf_ric( kstp )                       ! Richardson number dependent Kz 
    212 #if defined key_mpp_omp 
    213       IF( lk_zdftke )   CALL zdf_tke_jki( kstp )                   ! TKE closure scheme for Kz - j-k-i loops 
    214 #else 
     211 
    215212      IF( lk_zdftke )   CALL zdf_tke( kstp )                       ! TKE closure scheme for Kz 
    216 #endif 
     213 
    217214      IF( lk_zdfkpp )   CALL zdf_kpp( kstp )                       ! KPP closure scheme for Kz 
    218215 
    … …  
    275272                             CALL tra_ldf    ( kstp )       ! lateral mixing 
    276273#if defined key_agrif 
    277       IF(.NOT. Agrif_Root()) CALL Agrif_Sponge_tra( kstp )  ! tracers sponge 
     274      IF(.NOT. Agrif_Root()) CALL Agrif_Sponge_tra          ! tracers sponge 
    278275#endif 
    279276                             CALL tra_zdf    ( kstp )       ! vertical mixing 
    … …  
    306303                               CALL dyn_ldf( kstp )           ! lateral mixing 
    307304#if defined key_agrif 
    308       IF(.NOT. Agrif_Root())   CALL Agrif_Sponge_dyn( kstp )  ! momemtum sponge 
     305      IF(.NOT. Agrif_Root())   CALL Agrif_Sponge_dyn          ! momemtum sponge 
    309306#endif 
    310307                               CALL dyn_hpg( kstp )           ! horizontal gradient of Hydrostatic pressure 

  • branches/dev_001_SBC/NEMO/TOP_SRC/initrc.F90

    r760 r811  
    2525   USE trcini 
    2626   USE prtctl_trc      ! Print control passive tracers (prt_ctl_trc_init routine) 
     27   USE lib_mpp         ! distributed memory computing 
    2728    
    2829   IMPLICIT NONE 

  • branches/dev_001_SBC/NEMO/TOP_SRC/trcrst.F90

    r760 r811  
    335335      ENDDO 
    336336#endif 
    337       trb(:,:,:,:) = trn(:,:,:,:) 
     337!CT comment the line below which doesn't ensure  
     338!CT restartability of the GYRE_LOBSTER configuration 
     339!CT      trb(:,:,:,:) = trn(:,:,:,:) 
    338340 
    339341      CALL iom_close( numrtr ) 
Note: See TracChangeset for help on using the changeset viewer.