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 5956 for branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DOM/istate.F90 – NEMO

Ignore:
Timestamp:
2015-11-30T20:55:41+01:00 (8 years ago)
Author:
mathiot
Message:

ISF : merged trunk (5936) into branch

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DOM/istate.F90

    r5621 r5956  
    2929   USE daymod          ! calendar 
    3030   USE eosbn2          ! eq. of state, Brunt Vaisala frequency (eos     routine) 
    31    USE ldftra_oce      ! ocean active tracers: lateral physics 
     31   USE ldftra          ! lateral physics: ocean active tracers 
    3232   USE zdf_oce         ! ocean vertical physics 
    3333   USE phycst          ! physical constants 
    3434   USE dtatsd          ! data temperature and salinity   (dta_tsd routine) 
    3535   USE dtauvd          ! data: U & V current             (dta_uvd routine) 
    36    USE zpshde          ! partial step: hor. derivative (zps_hde routine) 
    37    USE eosbn2          ! equation of state            (eos bn2 routine) 
    3836   USE domvvl          ! varying vertical mesh 
    39    USE dynspg_oce      ! pressure gradient schemes 
    40    USE dynspg_flt      ! filtered free surface 
    41    USE sol_oce         ! ocean solver variables 
    4237   ! 
    4338   USE in_out_manager  ! I/O manager 
     
    7671      ! 
    7772 
    78       IF(lwp) WRITE(numout,*) ' ' 
     73      IF(lwp) WRITE(numout,*) 
    7974      IF(lwp) WRITE(numout,*) 'istate_ini : Initialization of the dynamics and tracers' 
    8075      IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
    8176 
    82       CALL dta_tsd_init                       ! Initialisation of T & S input data 
    83       IF( lk_c1d ) CALL dta_uvd_init          ! Initialization of U & V input data 
     77                     CALL dta_tsd_init        ! Initialisation of T & S input data 
     78      IF( lk_c1d )   CALL dta_uvd_init        ! Initialization of U & V input data 
    8479 
    8580      rhd  (:,:,:  ) = 0._wp   ;   rhop (:,:,:  ) = 0._wp      ! set one for all to 0 at level jpk 
     
    10398         ub   (:,:,:) = 0._wp   ;   un   (:,:,:) = 0._wp 
    10499         vb   (:,:,:) = 0._wp   ;   vn   (:,:,:) = 0._wp   
    105          rotb (:,:,:) = 0._wp   ;   rotn (:,:,:) = 0._wp 
    106          hdivb(:,:,:) = 0._wp   ;   hdivn(:,:,:) = 0._wp 
     100                                    hdivn(:,:,:) = 0._wp 
    107101         ! 
    108102         IF( cp_cfg == 'eel' ) THEN 
     
    119113            ENDIF 
    120114            IF ( ln_uvd_init .AND. lk_c1d ) THEN ! read 3D U and V data at nit000 
    121                CALL wrk_alloc( jpi, jpj, jpk, 2, zuvd ) 
     115               CALL wrk_alloc( jpi,jpj,jpk,2,  zuvd ) 
    122116               CALL dta_uvd( nit000, zuvd ) 
    123117               ub(:,:,:) = zuvd(:,:,:,1) ;  un(:,:,:) = ub(:,:,:) 
    124118               vb(:,:,:) = zuvd(:,:,:,2) ;  vn(:,:,:) = vb(:,:,:) 
    125                CALL wrk_dealloc( jpi, jpj, jpk, 2, zuvd ) 
     119               CALL wrk_dealloc( jpi,jpj,jpk,2,  zuvd ) 
    126120            ENDIF 
    127121         ENDIF 
    128          ! 
    129          CALL eos( tsb, rhd, rhop, gdept_0(:,:,:) )        ! before potential and in situ densities 
    130 #if ! defined key_c1d 
    131          IF( ln_zps .AND. .NOT. ln_isfcav)                                 & 
    132             &            CALL zps_hde    ( nit000, jpts, tsb, gtsu, gtsv,  &    ! Partial steps: before horizontal gradient 
    133             &                                            rhd, gru , grv    )  ! of t, s, rd at the last ocean level 
    134          IF( ln_zps .AND.       ln_isfcav)                                 & 
    135             &            CALL zps_hde_isf( nit000, jpts, tsb, gtsu, gtsv, gtui, gtvi,    &    ! Partial steps for top cell (ISF) 
    136             &                                            rhd, gru , grv , grui, grvi     )    ! of t, s, rd at the last ocean level 
    137 #endif 
    138122         !    
    139123         ! - ML - sshn could be modified by istate_eel, so that initialization of fse3t_b is done here 
     
    141125            DO jk = 1, jpk 
    142126               fse3t_b(:,:,jk) = fse3t_n(:,:,jk) 
    143             ENDDO 
     127            END DO 
    144128         ENDIF 
    145129         !  
    146       ENDIF 
    147       ! 
    148       IF( lk_agrif ) THEN                  ! read free surface arrays in restart file 
    149          IF( ln_rstart ) THEN 
    150             IF( lk_dynspg_flt )  THEN      ! read or initialize the following fields 
    151                !                           ! gcx, gcxb for agrif_opa_init 
    152                IF( sol_oce_alloc()  > 0 )   CALL ctl_stop('agrif sol_oce_alloc: allocation of arrays failed') 
    153                CALL flt_rst( nit000, 'READ' ) 
    154             ENDIF 
    155          ENDIF                             ! explicit case not coded yet with AGRIF 
    156130      ENDIF 
    157131      ! 
     
    162136      ! 
    163137      ! 
    164       un_b(:,:) = 0._wp ; vn_b(:,:) = 0._wp 
    165       ub_b(:,:) = 0._wp ; vb_b(:,:) = 0._wp 
     138      un_b(:,:) = 0._wp   ;  vn_b(:,:) = 0._wp 
     139      ub_b(:,:) = 0._wp   ;  vb_b(:,:) = 0._wp 
    166140      ! 
    167141      DO jk = 1, jpkm1 
     
    201175      !! References :  Philander ??? 
    202176      !!---------------------------------------------------------------------- 
    203       INTEGER  :: ji, jj, jk 
    204       REAL(wp) ::   zsal = 35.50 
     177      INTEGER  ::   ji, jj, jk 
     178      REAL(wp) ::   zsal = 35.50_wp 
    205179      !!---------------------------------------------------------------------- 
    206180      ! 
     
    232206      !!                and relative vorticity fields 
    233207      !!---------------------------------------------------------------------- 
    234       USE divcur     ! hor. divergence & rel. vorticity      (div_cur routine) 
     208      USE divhor     ! hor. divergence      (div_hor routine) 
    235209      USE iom 
    236210      ! 
     
    281255            tsb(:,:,:,jp_sal) = tsn(:,:,:,jp_sal) 
    282256            ! 
    283             ! set the dynamics: U,V, hdiv, rot (and ssh if necessary) 
     257            ! set the dynamics: U,V, hdiv (and ssh if necessary) 
    284258            ! ---------------- 
    285259            ! Start EEL5 configuration with barotropic geostrophic velocities  
     
    317291            ENDIF 
    318292            ! 
    319             CALL div_cur( nit000 )                  ! horizontal divergence and relative vorticity (curl) 
     293!!gm  Check  here call to div_hor should not be necessary 
     294!!gm         div_hor call runoffs  not sure they are defined at that level 
     295            CALL div_hor( nit000 )                  ! horizontal divergence and relative vorticity (curl) 
    320296            ! N.B. the vertical velocity will be computed from the horizontal divergence field 
    321297            ! in istate by a call to wzv routine 
     
    370346      !! 
    371347      !! ** Method  : - set temprature field 
    372       !!              - set salinity field 
     348      !!              - set salinity   field 
    373349      !!---------------------------------------------------------------------- 
    374350      INTEGER :: ji, jj, jk  ! dummy loop indices 
     
    456432      !!                 p=integral [ rau*g dz ] 
    457433      !!---------------------------------------------------------------------- 
    458       USE dynspg          ! surface pressure gradient             (dyn_spg routine) 
    459       USE divcur          ! hor. divergence & rel. vorticity      (div_cur routine) 
     434      USE divhor          ! hor. divergence                       (div_hor routine) 
    460435      USE lbclnk          ! ocean lateral boundary condition (or mpp link) 
    461436      ! 
    462437      INTEGER ::   ji, jj, jk        ! dummy loop indices 
    463       INTEGER ::   indic             ! ??? 
    464438      REAL(wp) ::   zmsv, zphv, zmsu, zphu, zalfg     ! temporary scalars 
    465439      REAL(wp), POINTER, DIMENSION(:,:,:) :: zprn 
    466440      !!---------------------------------------------------------------------- 
    467441      ! 
    468       CALL wrk_alloc( jpi, jpj, jpk, zprn) 
     442      CALL wrk_alloc( jpi,jpj,jpk,  zprn) 
    469443      ! 
    470444      IF(lwp) WRITE(numout,*)  
     
    528502      vb(:,:,:) = vn(:,:,:) 
    529503       
    530       ! WARNING !!!!! 
    531       ! after initializing u and v, we need to calculate the initial streamfunction bsf. 
    532       ! Otherwise, only the trend will be computed and the model will blow up (inconsistency). 
    533       ! to do that, we call dyn_spg with a special trick: 
    534       ! we fill ua and va with the velocities divided by dt, and the streamfunction will be brought to the 
    535       ! right value assuming the velocities have been set up in one time step. 
    536       ! we then set bsfd to zero (first guess for next step is d(psi)/dt = 0.) 
    537       !  sets up s false trend to calculate the barotropic streamfunction. 
    538  
    539       ua(:,:,:) = ub(:,:,:) / rdt 
    540       va(:,:,:) = vb(:,:,:) / rdt 
    541  
    542       ! calls dyn_spg. we assume euler time step, starting from rest. 
    543       indic = 0 
    544       CALL dyn_spg( nit000, indic )       ! surface pressure gradient 
    545  
    546       ! the new velocity is ua*rdt 
    547  
    548       CALL lbc_lnk( ua, 'U', -1. ) 
    549       CALL lbc_lnk( va, 'V', -1. ) 
    550  
    551       ub(:,:,:) = ua(:,:,:) * rdt 
    552       vb(:,:,:) = va(:,:,:) * rdt 
    553       ua(:,:,:) = 0.e0 
    554       va(:,:,:) = 0.e0 
    555       un(:,:,:) = ub(:,:,:) 
    556       vn(:,:,:) = vb(:,:,:) 
    557         
    558       ! Compute the divergence and curl 
    559  
    560       CALL div_cur( nit000 )            ! now horizontal divergence and curl 
    561  
    562       hdivb(:,:,:) = hdivn(:,:,:)       ! set the before to the now value 
    563       rotb (:,:,:) = rotn (:,:,:)       ! set the before to the now value 
    564       ! 
    565       CALL wrk_dealloc( jpi, jpj, jpk, zprn) 
     504      ! 
     505!!gm  Check  here call to div_hor should not be necessary 
     506!!gm         div_hor call runoffs  not sure they are defined at that level 
     507      CALL div_hor( nit000 )            ! now horizontal divergence 
     508      ! 
     509      CALL wrk_dealloc( jpi,jpj,jpk,   zprn) 
    566510      ! 
    567511   END SUBROUTINE istate_uvg 
Note: See TracChangeset for help on using the changeset viewer.