  SUBROUTINE cas_1DVconso

!&E---------------------------------------------------------------------
!&E                 ***  ROUTINE cas_1DVconso  ***
!&E
!&E ** Purpose : definition of a test case for 1DV sediment consolidation
!&E              (containing sandy and muddy sediments)
!&E
!&E ** Description :
!&E                  1DV test case to simulate sedimentation and consolidation processes of mud-sand
!&E                  mixtures, as sediment settling column experiments. It takes into account segregation
!&E                  processes and simulate hindered settling, permeability and effective stress regimes,
!&E                  as presented and validated in Grasso et al. (Ocean Dynamics 2014).
!&E                  This test case is based on:
!&E                  - 1DV hydrodynamics test case "cas_turbkatophil"
!&E                  - 2DV sediment test case "cas_tidalflat2DV"
!&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       !  2011-03 (F. Cayocca) new test case issued from test caschansed
!&E       !  2011-11  (V. Garnier) Spatial extension of ssh(liminm2:limaxp2,ljminm2:ljmaxp2)
!&E       !  2014-05  (F. Grasso) new test case for sediment consolidation
!&E
!&E---------------------------------------------------------------------





























































!! * Modules used
   USE parameters, ONLY : rsh,rlg,riolg,valmanq,kmax,imin,imax,lchain,&
                          ljminm1,ljmaxp2,liminm1,limaxp2
   USE comvarp3d, ONLY : ect

   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
   USE output,    ONLY : pasor,l_out_sal,l_out_temp,fileresu,tdebsor,tfinsor, &
                         l_out_buoy,l_out_wind,l_out_zaxis_w,l_out_turb,   &
                         iextmin,iextmax,jextmin,jextmax
   USE toolgeom,  ONLY : f_zu
   USE comvarp2d, ONLY : temp
   USE comvars3d, ONLY : tetaref,saliref,nzinit,kzinit,qmaxz,l_stateeq_lin,  &
                         turb_0eq_option,turb_nbeq,turb_2eq_option,bgdiff,OptionSF
   USE comvars2d, ONLY : fisac,gwac,dtmin,dtmax,dt,suffix,fvisc,kx,ky,qmax,halfdt, &
                         icon,fu,fv,h0,tobs,hminfrot,iscreenlog, &
                         h0_g,hx_g,hy_g,h0fond,hm,vis2d,vish_xe,vish_phi
   USE comdiag,   ONLY : l_diag,file_diag,tdebdiag,tfindiag
   USE comturb,   ONLY : eps

   USE comsubstance
   USE casinit,   ONLY : cas_init_subs,cas_init_subs_disspart,cas_init_subs_budget


   USE parameters, ONLY : ksdmax
   USE comMUSTANG
   USE casinit,   ONLY : cas_init_sedim


   IMPLICIT NONE

!! * Local declarations
   CHARACTER(LEN=19) :: date_start,date_end,date_start3d,date_startagrif
   CHARACTER(LEN=19) :: date_startout,date_endout
   INTEGER           :: i,j,k,isubs,iv,ig
   INTEGER           :: ni_wall_e,ni_wall_w,nj_wall_n,nj_wall_s
   INTEGER           :: IERR_MPI
   REAL(KIND=rsh)    :: dxval,dyval,dzsmaxinit
   REAL(KIND=rsh)    :: z0bottom,z0surface,depthcst,zr
   REAL(KIND=rlg)    :: tool_datosec
   LOGICAL           :: l_filesubs,l_casinit_hsed


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

   !PRINT*, 'enter CAS_CONSO routine'

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

! Dynamics
! ----------
   CALL cas_init_param

! Tracers
! ----------
   CALL cas_init_tspar
   kx=1.0_rsh
   ky=1.0_rsh
   qmax=8.0_rsh
   qmaxz=10.0_rsh
   saliref=0.0_rsh
   l_stateeq_lin=.FALSE.


! Bathymetry
! ----------
! assignement : input variables for cas_init_bathy routine
   depthcst = 1.0_rsh
   ni_wall_e = 1  ! land over ni_wall_e cells at the eastern part of the domain
   ni_wall_w = 2  ! land over ni_wall_w cells (from imin column) at the western part of the domain
   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


   hm=0.5_rsh

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

! Grid parametrization
! --------------------
! assignement : input variables for cas_init_grid routine

! parametres de forme
!---------------------
   dxval = 1.0_rsh
   dyval = 1.0_rsh
   fisac=1.0_riolg
   gwac=1.0_riolg



   CALL cas_init_grid(dxval,dyval)
   fu(:,:)=0.0_rsh
   fv(:,:)=0.0_rsh

! IO parametrization
! ------------------
! assignement : input variables for cas_init_io routine
   date_start=     '01/01/1900 00:00:00'
   date_end=       '01/02/1900 00:00:00'
   date_start3d=   '01/01/1900 00:00:00'
   date_startagrif='01/01/1900 00:00:00'
   dtmin=30.0_rlg
   dtmax=50.0_rlg
   dt   =30.0_rlg
   halfdt=dt*0.5_rlg
   suffix='conso'


   date_startsub='01/01/1900 00:00:00'
   tdebsubs=tool_datosec(date_startsub)

! consolidation only if  l_consolid=.true.  or l_diffused or l_bioturb
   date_start_sed=date_start
   tstart_sed=tool_datosec(date_start_sed)
   tstart_sed=MAX(tstart_sed,tdebsubs)
   t_sed=tstart_sed
   dt_sed=500.0_rlg



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

! over-assignement
   iextmin(1)=2
   iextmax(1)=2
   jextmin(1)=2
   jextmax(1)=2
   pasor(1)=600.0_rlg
   tobs=2.0_rlg*3600.0_rlg

   date_startout='01/01/1900 00:00:00'
   date_endout=date_end
   tdebsor(1)=tool_datosec(date_startout)
   tfinsor(1)=tool_datosec(date_endout)
!pasor(1)=3.0_rsh*86400.0_rlg
   l_diag=.FALSE.

! Open boundaries
! ---------------
   CALL cas_init_obc

! Viscosity  (none by default)
!-----------
   CALL cas_init_visc
   fvisc=0.1_rsh                                          ! =0 in casturb
   vis2d(:,:) = fvisc*0.01_rsh*dyval**1.15_rsh            ! =0 in casturb
   vish_xe(:,:,:) = fvisc*0.01_rsh*dyval**1.15_rsh        ! =0 in casturb
   vish_phi(:,:,:)= fvisc*0.01_rsh*dyval**1.15_rsh        ! =0 in casturb

! Turbulence
! ----------
   z0bottom=0.0001
   z0surface  =0.0001_rsh
   nzinit=0.001_rsh
   kzinit=0.00001_rsh

! Gaspard parametrization
   CALL cas_init_turb(kzinit,nzinit,z0surface,z0bottom)
   turb_nbeq=2
   turb_0eq_option=1
   turb_2eq_option=2
   ect(:,:,:)=0.000001_rsh
   eps(:,:,:)=0.000001_rsh
   bgdiff=0.000001_rsh


! Meteorology (none by default)
! -----------
! no forcing
   CALL cas_init_meteo

! Initialization of the pronostic variables
! -----------------------------------------
! set to 0,tetaref,saliref
   CALL cas_init_progvar

! over-assignement


! Miscellaneous
! -------------

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

! Substances
!
!--------------

! Initialization of the number of variables depending on their type
! (gravel,sand,mud,constitutive,non-constitutive,dissolved,fixed,driving or intermediate)
!----------------------------------------------------------------------------------------
   nv_grav=0
   nv_sand=1
   nv_mud=1
   nv_cp=0
   nv_ncp=0
   nv_sorb=0
   nv_dis=0
   nv_fix=0
   nv_dri=0
   nv_int=0
   nvp=2
   nvpc=2
   l_filesubs=.TRUE.
   IF (l_filesubs) THEN  ! we use "sub_read_var" instead of "cas_init_subs"
     name_filesubs = '../inputs/variable_Cancale_2cl.dat'
     CALL cas_init_subs_budget(date_start)
   ELSE
! allocation et initialisation des tableaux cv_wat,f90,t90...
     CALL cas_init_subs(date_start)   ! nvp must be defined before

   ENDIF
   l_bilan=.TRUE.
   dtout_budget=1.0_rlg

   IF(.NOT.l_filesubs) THEN
! caracteristics of each variable (settling velocity)
! ---------------------------------------------------
! typart=0: non-const. part., typart=1: const. part.
! itypv==1: GRAV, itypv==2: SAND, itypv==3: MUDS/PART, itypv==4: NCP,itypv==5: SORB,
! itypv==6: DISS, itypv==7: FIXE, itypv==8: DRIV, itypv==9: INTE,

!-------------------------------- beginning of substance definition
     isubs=1
     name_var(isubs)='Light Sand'
     unit_var(isubs)='kg.m-3'
     typart(isubs)=1
     itypv(isubs)=2
     ws_free_opt(2)=0
     ws_free_min(2)=.0001_rsh
     ws_free_max(2)=.0001_rsh
     ws_free_para(1,2)=0.0_rsh
     ws_free_para(2,2)=0.0_rsh
     ws_free_para(3,2)=0.0_rsh
     ws_free_para(4,2)=0.0_rsh
     ws_hind_opt(2)=0
     ws_hind_para(1,2)=.0_rsh
     ws_hind_para(2,2)=.0_rsh 
     tocd(isubs)=1.0_rsh
     cini_wat(isubs)=0.0_rsh
     cv_wat(isubs,:,:,:) = cini_wat(isubs)

     isubs=2
     name_var(isubs)='Mud'
     unit_var(isubs)='kg.m-3'

     typart(isubs)=1
     itypv(isubs)=3
     ws_free_opt(2)=0
     ws_free_min(2)=0003_rsh
     ws_free_max(2)=.006_rsh
     ws_free_para(1,2)=0.0_rsh
     ws_free_para(2,2)=0.0_rsh
     ws_free_para(3,2)=0.0_rsh
     ws_free_para(4,2)=0.0_rsh
     ws_hind_opt(2)=0
     ws_hind_para(1,2)=.0_rsh
     ws_hind_para(2,2)=.0_rsh 
     tocd(isubs)=1.0_rsh
     cini_wat(isubs)=0.0_rsh
     cv_wat(isubs,:,:,:) = cini_wat(isubs)
!-------------------------------- end of substance definition

   ENDIF


! si il existe des variables particulaires
   IF (nvp>0) THEN

      CALL cas_init_subs_disspart(l_filesubs)

! over-assignement
!  DO j=ljminm1,ljmaxp2
!    DO i=liminm1,limaxp2
!      IF (h0(i,j) > -valmanq) THEN
!        DO k=1,kmax
!          zr=-f_zu(h0(i,j),0.0_rsh,k,i,j)
!          cv_wat(1,k,i,j)= 0.1 - 0.00001_rsh/(9.81_rsh*0.0016)*zr
!          cv_wat(2,k,i,j)= 50. - 0.001_rsh/(9.81_rsh*0.0016)*zr
!        cv_wat(:,k,i,j) = cini_wat(:)
!       END DO

!     ELSE
!        cv_wat(:,:,i,j) = valmanq
!     END IF
!   END DO
!  END DO
    
! obligation de definir les ros et diam des particules apres l appel a cas_init_subs_disspart dans lasuelle ces variable sont allouees
! vase ne se sert pas du diametre mais il ne doit pas etre nul et il se sert de ros

   IF(.NOT.l_filesubs) THEN
     ros(1)=2600.0_rsh
     diam_sed(1)=0.000100_rsh
     cini_sed(1)=0.5       ! initial fraction [0-1] of total sediment concentration (cseduni)
     ros(2)=2600.0_rsh
     diam_sed(2)=0.000020_rsh
     cini_sed(2)=0.5       ! initial fraction [0-1] of total sediment concentration (cseduni)
   ENDIF



   ENDIF

! fin substance

!   VARIABLES DIAGNOSTIQUES
! ---------------------------------------
    
   
! River and outflows
! -------------------



! Initialization of sedimentary parameters
! ----------------------------------------
! pre-assignement
   l_consolid = .TRUE.
   hseduni = 1.0_rsh
   cseduni = 106.0_rsh
   ksmiuni = 1
   ksmauni =200
   dzsmaxinit=5.0e-3
   l_casinit_hsed=.FALSE.

! Initialization of sediment layers
! ---------------------------------
   CALL cas_init_sedim(l_filesubs,dzsmaxinit=dzsmaxinit,l_casinit_hsed=l_casinit_hsed)

! over-assignement
   dzsmin = 1.0e-5
   cfreshmud = 100.0_rsh

   
   IF(l_diffused .or. l_consolid .OR. l_bioturb)l_dyn_insed = .TRUE.

! Initialization of sediment output
! ---------------------------------
!  **********choice output in sediment *******************************
!  sediment output over all layers                     --- option 1, default ---
!  choice_nivsed_out=1
!
!  sediment output over a fixed number of layers from the surface ----option 2---
!  choice_nivsed_out=2
!   nk_nivsed_out=3
!   ep_nivsed_out(:)=0.0_rsh
!
!  sediment output : max thickness (in m)                       ---- option 3---
!  choice_nivsed_out=3
!   thickness max MUST be > dzsmaxuni
!   nk_nivsed_out=0
!   epmax_nivsed_out=0.030_rsh
!
!  sediment output : 5 layers max - thicknesses in m            ---- option 4---
!   donnera une couche supplmentaire avec resultats integres jusqu au fond
!   choice_nivsed_out=4
!   nk_nivsed_out=4
!   ep_nivsed_out(1)=0.001_rsh
!   ep_nivsed_out(2)=0.005_rsh
!   ep_nivsed_out(3)=0.010_rsh
!   ep_nivsed_out(4)=0.050_rsh



! fin sedim


!--------------------------------------------------
   PRINT*, 'Number of substances :',nv_tot
   PRINT*, 'cseduni (kg.m-3) =',cseduni
   DO iv=1,nv_tot
     PRINT*, iv,': "',TRIM(name_var(iv)),'", cini_sed (fraction) =',cini_sed(iv)
     PRINT*, iv,': "',TRIM(name_var(iv)),'", cini_wat (',TRIM(unit_var(iv)),') =',cini_wat(iv)
   END DO
!--------------------------------------------------


   !PRINT*, 'exit CAS_1DVCONSO routine'

  END SUBROUTINE cas_1DVconso
