  SUBROUTINE cas_instbtjet

!&E---------------------------------------------------------------------
!&E                 ***  ROUTINE cas_instbtjet  ***
!&E
!&E ** Purpose : test case for an instable barotropic current
!&E
!&E ** Description :
!&E
!&E ** Called by : cas_init
!&E
!&E ** External calls :
!&E
!&E ** Used ij-arrays :
!&E
!&E ** Modified variables :
!&E
!&E ** Reference :
!&E
!&E ** History :
!&E       !  2007-03  (T. Duhaut)  Original code
!&E       !  2008-12  (V. Garnier) shift of sigma values : -1 < sig < 0
!&E       !  2011-11  (V. Garnier) Spatial extension of ssh(liminm2:limaxp2,ljminm2:ljmaxp2)
!&E
!&E---------------------------------------------------------------------





























































!! * Modules used
   USE parameters,ONLY : pi,rsh,rlg,riolg,kmax,jmax,imax,jmin,imin, &
                         limin,ljmin,limaxp1,ljmaxp1,ljmaxp2,limax, &
                         ljmax,ljmaxm1,liminp1,limaxp2,liminm2,ljminm2,ljminp2
   USE comvars2d, ONLY : fisac,gwac,fu,fv,dtmin,dtmax,dt,halfdt,suffix,tobs,  &
                         l_corio_cst,icon,l_modele2d,radearth,qmax,    &
                         dfiac,dgac,radearth,l_obc_cycl_x,dy,ig,id
   USE comvarp2d, ONLY : ssh,u,v
   USE comvars3d, ONLY : rhoref,saliref,kzinit,nzinit,l_stateeq_lin,qmaxz
   USE comvarp3d, ONLY : uz,vz
   USE obccombine, ONLY : l_obc_east,l_obc_west,l_obc_cycl
   USE output,    ONLY : pasor,l_out_sal,l_out_temp
   USE toolmpi,ONLY : ex_i_rsh,ex_j_rsh,ex_peri_x_rsh,ex_peri_y_rsh
   USE casinit,   ONLY : cas_init_param,cas_init_grid,cas_init_bathy, &
                         cas_init_io,cas_init_obc,cas_init_meteo,     &
                         cas_init_visc,cas_init_turb,cas_init_tspar,  &
                         cas_init_progvar,cas_init_misc,cas_init_river

   IMPLICIT NONE

!! * Arguments

!! * Local declarations
   CHARACTER(LEN=19) :: date_start,date_end,date_start3d,date_startagrif
   INTEGER           :: i,j,k,taille_germe,optionjet
   INTEGER, ALLOCATABLE, DIMENSION(:) :: germe
   REAL(KIND=rsh)    :: radi,tir
   REAL(KIND=rlg)    :: f,yo,y,uo,et0
   REAL(KIND=rlg)    :: tool_datosec
   INTEGER           :: ni_wall_e,ni_wall_w,nj_wall_n,nj_wall_s
   REAL(KIND=rsh)    :: dxval,dyval,z0bottom,z0surface,depthcst

!!----------------------------------------------------------------------
!! * Executable part


! general parametrization (each parameter can be reassigned)
! -----------------------

! dynamics
   CALL cas_init_param

! over-assignement
!l_corio_cst=.FALSE.  ! beta-plane
   l_modele2d=.TRUE.

! tracers : among others, saliref=saliref_lin=35.5_rsh, tetaref=tetaref_lin=10.0_rsh
   CALL cas_init_tspar

! over-assignement
   icon=0
   qmax=12.0_rsh
   qmaxz=10.0_rsh
   rhoref=1000.0_rsh
   saliref=31.0_rsh
   l_stateeq_lin=.FALSE.


! Bathymetry
! ----------

! assignement : input variables for cas_init_bathy routine
   depthcst = 200.0_rsh
   nj_wall_s = 2  ! land over nj_wall_s cells (from jmin row) at the southern part of the domain
   nj_wall_n = 1  ! land over nj_wall_n cells at the northern part of the domain

! flat bottom, plus walls
   CALL cas_init_bathy(depthcst,nj_wall_n=nj_wall_n,nj_wall_s=nj_wall_s)

! over-assignement


! grid parametrization
! --------------------

! assignement : input variables for cas_init_grid routine
   dxval=10000.0_rsh
   dyval=10000.0_rsh
   radi=57.29578_rsh
   dfiac=dyval*radi/radearth
   dgac=dyval*radi/radearth/COS(45.0_rsh/radi)
   fisac=45.0_riolg-REAL(jmax,riolg)/2.0_riolg*dfiac
   gwac=20.0_riolg

   CALL cas_init_grid(dxval,dyval)

! over-assignement
   DO j = ljmin,ljmax
     DO i = limin,limax
       fu(i,j)=0.0001458_rsh*SIN((30.0_rsh+(float(j-(jmax-1)/2))*dy(i,j)/radearth)/radi)
       fv(i,j)=0.0001458_rsh*SIN((30.0_rsh+(float(j-(jmax-1)/2)-0.5_rsh)*dy(i,j)/radearth)/radi)
     ENDDO
   ENDDO
   dfiac=dyval*radi/radearth
   dgac=dyval*radi/radearth/COS(45.0_rsh/radi)


! IO parametrization
! ------------------

! assignement : input variables for cas_init_io routine
   date_start=   '01/01/1900 00:00:00'
   date_end=   '01/05/1900 00:00:00'
   date_start3d=  '01/01/1900 00:00:00'
   date_startagrif='01/01/1900 00:00:00'
   dtmin=100.0_rlg
   dtmax=800.0_rlg
   dt=800.0_rlg
   halfdt=dt*0.5_rlg

   suffix='JetBt'

   CALL cas_init_io(date_start,date_end,date_start3d,date_startagrif,suffix)

! over-assignement
   pasor(1)=24.0_rlg*3600.0_rlg
   l_out_sal(1)=.FALSE.
   l_out_temp(1)=.FALSE.
   tobs=24.0_rlg*3600.0_rlg


! open boundaries
! ---------------

! no open boundary
   CALL cas_init_obc


! over-assignement
! geometry (where are the open boundaries ?)
   l_obc_west=.TRUE.
   l_obc_east=.TRUE.
! used scheme to manage obc
   l_obc_cycl = .TRUE.
   l_obc_cycl_x=.TRUE.

! viscosity
!-----------

! none
   CALL cas_init_visc


! turbulence
! ----------

   z0bottom   =1.0e-13
   z0surface  =0.0_rsh
   kzinit=0.0_rsh
   nzinit=0.0_rsh

! Gaspard parametrization
   CALL cas_init_turb(kzinit,nzinit,z0surface,z0bottom)

! over-assignement


! meteorology
! -----------

! no forcing
   CALL cas_init_meteo


! River and outflows
! -------------------

! no river
   CALL cas_init_river


! initialization of the pronostic variables
! -----------------------------------------

! set to 0,tetaref,saliref
   CALL cas_init_progvar

! over-assignement

! Sea surface height and barotropic velocities
!---------------------------------------------


   optionjet=2

   IF (optionjet == 1 ) THEN

     uo=1.0_rlg
     yo=100000.0_rlg
!!  profile en cos**2 pour Uo
     DO j=ljminp2,ljmaxm1
       y=float(j-(jmax-1)/2)*dyval
       IF (ABS(y) > yo) THEN
         IF ( y > 0.0_rlg ) THEN
           ssh(limin:limax,j)=-fu(limin:limax,j)*yo*uo/(2.0_rlg*9.81_rlg)
         ELSE
           ssh(limin:limax,j)=fu(limin:limax,j)*yo*uo/(2.0_rlg*9.81_rlg)
         END IF
       ELSE
         ssh(limin:limax,j)=-fu(limin:limax,j)*yo*uo/(pi*9.81_rlg)*(pi*y*0.5_rlg/yo+0.5_rlg*SIN(pi*y/yo))
         u(limin:limax,j)=uo*(COS(pi*y*0.5_rlg/yo))**2
       END IF
     END DO

   ELSEIF (optionjet == 2 ) THEN
     yo=50000.0_rlg
     et0=0.8_rlg
!! profile gausien zonal pour le surface libre
     DO j=ljminp2,ljmaxm1
       y=float(j-(jmax-1)/2)*dyval
       ssh(limin:limax,j)=et0*EXP(-0.5_rlg*(y/yo)**2)
       u(limin:limax,j)=(9.81_rlg*et0*y/fu(limin:limax,j)/yo/yo)*EXP(-0.5_rlg*(y/yo)**2)
     END DO

   ELSEIF (optionjet == 3 ) THEN
     yo=30000.0_rlg
     et0=0.5_rlg
!! profile en tanh pour la surface libre
     DO j=ljminp2,ljmaxm1
       y=float(j-(jmax-1)/2)*dyval
       ssh(limin:limax,j)=-et0*TANH(y/yo)
       u(limin:limax,j)=(9.81_rlg*et0)/(fu(limin:limax,j)*yo)*(1.0_rlg-TANH(y/yo)**2)
     END DO

   END IF

   CALL ex_i_rsh(-2,2,1,liminm2,limaxp2,ljminm2,ljmaxp2,ssh(liminm2:limaxp2,ljminm2:ljmaxp2))
   CALL ex_j_rsh(-2,2,1,liminm2,limaxp2,ljminm2,ljmaxp2,ssh(liminm2:limaxp2,ljminm2:ljmaxp2))
   CALL ex_i_rsh(-2,2,1,liminm2,limaxp2,ljminm2,ljmaxp2,u(liminm2:limaxp2,ljminm2:ljmaxp2))
   CALL ex_j_rsh(-2,2,1,liminm2,limaxp2,ljminm2,ljmaxp2,u(liminm2:limaxp2,ljminm2:ljmaxp2))

! add a small perturbation to enhance
! the instability
!------------------------------------
   CALL RANDOM_SEED(size=taille_germe)
   ALLOCATE(germe(taille_germe))
   germe(:)=1
   CALL RANDOM_SEED(put=germe(:))
   DO j=jmin,jmax
     DO i=imin,imax
       CALL RANDOM_NUMBER(tir)
       IF ( (liminm2 <= i) .AND. (i <= limaxp2) .AND. (ljminm2 <= j) .AND.  (j <= ljmaxp2) ) THEN
       u(i,j)=u(i,j)*(1.0_rsh+0.00001_rsh*(0.5_rsh-tir))
       v(i,j)=v(i,j)*(1.0_rsh+0.00001_rsh*(0.5_rsh-tir))
       ENDIF
     END DO
   END DO



! initialization of vertical scale factor (part relative to dc*sds)
!------------------------------------------------------------------

! Domain of calculation
! Percent of wet surface on each cell
   CALL cas_init_misc


! Periodic conditions use ig,id,jb,jh which are defined in cas_init_misc

   IF (l_obc_cycl_x) THEN

     CALL ex_peri_x_rsh(2,1,liminm2,limaxp2,ljminm2,ljmaxp2,ssh(liminm2:limaxp2,ljminm2:ljmaxp2))
     CALL ex_peri_x_rsh(2,1,liminm2,limaxp2,ljminm2,ljmaxp2,u(liminm2:limaxp2,ljminm2:ljmaxp2))

   END IF

  END SUBROUTINE cas_instbtjet
