!!======================================================================

  SUBROUTINE dyn2d_dt(u,v)

!&E---------------------------------------------------------------------
!&E                 ***  ROUTINE dyn2d_dt  ***
!&E
!&E ** Purpose : Estimate of optimal time step
!&E
!&E ** Description : CFL criterion
!&E
!&E ** Called by : step2dxy.F90
!&E
!&E ** External calls :
!&E
!&E ** Used ij-arrays : u,v
!&E
!&E ** Modified variables : dt
!&E
!&E ** Reference :
!&E
!&E ** History :
!&E       !  1987-02  (J.-C. Salomon, P. Lazure)  Original code
!&E       !  2007-09  (V. Garnier) homogeneneity of dx dy...
!&E       !  2009-02  (V. Garnier T. Odaka) MPI2D
!&E       !  2009-05  (F. Dumas)  AGRIF
!&E
!&E---------------------------------------------------------------------





























































!! * Modules used
   USE parameters
   USE comvars2d
   USE toolmpi, ONLY : findmax_mpi_3,MASTER

   IMPLICIT NONE

!! * Arguments
   REAL(KIND=rsh), DIMENSION(liminm2:limaxp2,ljminm2:ljmaxp2), INTENT( in ) :: u,v

!! * Local declarations
   INTEGER                           :: i,j,icrt,jcrt,IERR_MPI
   REAL(KIND=rsh)                    :: gamm,gamp
   REAL(KIND=rsh)                    :: vitu,vitv,vit2
   REAL(KIND=rsh)                    :: crtmax,rap,facor
   REAL(KIND=rlg)                    :: dtmaxp,dtint

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

   gamm=-1.0_rsh
   gamp=-0.5_rsh
   crtmax=0.0_rsh

! Estimate of the maximum velocity and its relative CFL

   DO j=ljminp2,ljmaxm1
     DO i=liminp2,limaxm1
       vitu=(u(i,j)+u(i-1,j))*0.5_rsh
       vitv=(v(i,j)+v(i,j-1))*0.5_rsh
       vit2=(vitu*vitu+vitv*vitv)*dt**2/(MIN(dy(i,j),dx(i,j))**2)
       IF (vit2 >= crtmax) THEN
         crtmax=vit2
         icrt = i
         jcrt = j
       END IF
     END DO
   END DO

   crtmax=SQRT(crtmax)

   CALL FINDMAX_MPI_3(crtmax,icrt,jcrt)
      
! Estimate of the new time step

   rap=(crtmax-cflcrt)/(0.1_rsh*cflcrt)

   IF (rap >= 0.0_rsh) THEN
     facor = MAX(1.0_rsh+rap*gamm,0.8_rsh)
   ELSE
     facor = MIN(1.0_rsh+rap*gamp,1.02_rsh)
   END IF

   crtmax=crtmax/dt*MAX(MIN(dt*facor,dtmax),dtmin)

   dtint=MAX(MIN(dt*facor,dtmax),dtmin)
   !IF (.NOT.Agrif_Root()) THEN
     !IF ( Agrif_NbStepint().EQ.(Agrif_Rhot()-1)) THEN
       ! dtmaxp = dt*facor
       !CALL Agrif_Set_Parent(dt, MIN(Agrif_parent(dtmax),dtmaxp*Agrif_Rhot()))
     !ENDIF
   !ELSE
     dt=dtint
   !ENDIF

   IF (MASTER) THEN 
   WRITE(iscreen,10) dt,crtmax,icrt,jcrt
   ENDIF

   IF (dt < dt_minallow) THEN
     WRITE(ierrorlog,*) 'The time step is too small ( < ',dt_minallow
     WRITE(ierrorlog,*) 'The simulation is stopped'
     CALL MPI_FINALIZE(IERR_MPI)
     STOP
   END IF

   halfdt=0.5_rlg*dt

10 FORMAT ('dt=',f6.1,'   Nb Courant max=',f5.2,' en (',2i4,')')

  END SUBROUTINE dyn2d_dt
