subroutine initialize_user
  !=======================================================
#define DEBUG_LEVEL 3

  use mpi_parameters      , only: comm,myid
  use dependent_variables
  use grid_info           , only: Grid 
  use counters_flags_etc  , only: user_subs_flag,restart_flag,ierr 
  use dimensional_scales  , only: L_scale,U_scale,s1_scale,s2_scale
  use user_parameters     , only: Km_of_z,Kt_of_z 
  use user_routines       , only: user_ics 

  implicit none

  include '../input/problem_size.h'
  include 'mpif.h'

  integer    :: i,j,k,ifield,num_champ
  real       :: alpha1,alpha2,alpha3,alpha4,Lx,Ly,Lz
  real,dimension(:,:,:,:),allocatable   :: vals


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

  if( user_subs_flag == 'no' ) return

#if DEBUG_LEVEL >= 1
  if( myid==0) then
     write(0,*) ' '
     write(0,*) 'hello world from initialize_user'
  endif
#endif



  !========================================================
  !ENABLE VERTICALLY VARYING EDDY VISCOSITY AND DIFFUSIVITY
  !
  allocate( Km_of_z(locnz), Kt_of_z(locnz) )

  call set_diffusion

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



  if( restart_flag == 'yes' ) return

  ! we can also specify any non-zero initial conditions


  num_champ=3+num_scalars

  allocate ( vals(nx,ny,locnz,num_champ) )

  !! appel de user_ics pour definir la condition initiale
  do ifield=1,num_champ
     call user_ics( Grid(0)%x(:,1)*L_scale , Grid(0)%y(:)*L_scale , Grid(0)%z(1,:)*L_scale, &
                    ifield , nx , ny , locnz, vals(:,:,:,ifield) ) 
  enddo

  !! Orientation de la grille pour U & W, affectation des champs vals(:,:,:,ifield)
  do k=1,locnz                                                                                           
     do i=1,nx
        alpha1=Grid(0)%x_xi(i,k)
        alpha2=Grid(0)%x_zeta(i,k)
        alpha3=Grid(0)%z_xi(i,k)
        alpha4=Grid(0)%z_zeta(i,k)
        if ( i == 10 .and. k == 10 ) print*,'alpha1-4=',alpha1,alpha2,alpha3,alpha4
        do j=1,ny
           U(i,j,k)=1./( alpha1*alpha4-alpha3*alpha2 ) * (  alpha4*vals(i,j,k,1) - alpha2*vals(i,j,k,3) )
           W(i,j,k)=1./( alpha1*alpha4-alpha3*alpha2 ) * ( -alpha3*vals(i,j,k,1) + alpha1*vals(i,j,k,3) )
        enddo
     enddo
  enddo
  V = vals(:,:,:,2)
  s1= vals(:,:,:,4)
  if ( num_scalars > 1 ) s2 = vals(:,:,:,5)

  !! adimensionnement
  U=U/U_scale
  V=V/U_scale
  W=W/U_scale
  s1=s1/s1_scale
  s2=s2/s2_scale

  !! affichage des min/max a l'ecran
  call maxim


  deallocate( vals )

end subroutine initialize_user