else if (iusclb(kzone).eq.jbord1) then


! Calculation of the criterion
!===============================================================================
  uxn = ettp(ip,jup)*xnn
  vyn = ettp(ip,jvp)*ynn
  wzn = ettp(ip,jwp)*znn

!    incident velocity calculation
      vi =  abs(ettp(ip,jup) * surfbo(1,kface)                       &
           + ettp(ip,jvp) * surfbo(2,kface)                       &
           + ettp(ip,jwp) * surfbo(3,kface) )                      &
           / surfbn(kface)
     
     v=sqrt(ettp(ip,jup)**2+ettp(ip,jvp)**2+ettp(ip,jwp)**2)



!  energ = 0.5d0*ettp(ip,jmp)*(uxn**2+vyn**2+wzn**2) ! ok au 3 débits et avec différentes tailles mais ne vuet rien dire physiquement
!  energt   = 7d-8*ettp(ip,jdp)

!  energ = 0.5d0*ettp(ip,jmp)*(v*v) ! energie cinétique totale et non la seule composante normale : ne fonctionne plus si on baisse le débit
!  energt   = 6d-7*ettp(ip,jdp)

!  energ = 0.5d0*ettp(ip,jmp)*(uxn+vyn+wzn)**2 ! fonctionne avec 5d-8 en changeant le débit
!  energt   = 5d-8*ettp(ip,jdp)

!     energ = 0.5d0*ettp(ip,jmp)*(vi*vi) ! ne fonctionne plus si on baisse le débit
!     energt   = 5.9d-8*ettp(ip,jdp)


     energ = 0.5d0*ettp(ip,jmp)*(vi*vi) !
     energt   =6.d-14*v/vi ! prendre 5.d-14 si l'efficacité est basée sur la chicane seule, et 5.5d-14 si chicane + plenum



!     energ = 0.5d0*ettp(ip,jmp)*(vi**2) ! 
!     nu=9.02d-5
!     re=vi*ettp(ip,jdp)/nu ! Reynolds number
!     energt   = 1.6d-12*re

!B/ Sommerfeld parameter analysis



!    Reynolds number calculation
!    dynamic viscosity taken as constant : µ=roh*nu, rho = ruslag(1,1,iropt), nu=9.02e-5m²/s
!    nu=9.02d-5
!    re=vi*ettp(ip,jdp)/nu ! Reynolds number
!    we=((ruslag(1,1,iropt))*ettp(ip,jdp)*(vi**2))/(32.d-3) ! weber number
!    oh=sqrt(we)/re ! Ohnesorge number
!    k=oh*(re**(5/4)) ! Sommerfeld parameter correlation avec 1.2 ok pour 120Pa, ne fonctionne plus si on baisse le débit
!    k=we*(oh**(-2/5)) ! Wazel number calculation ok avec 1.3 pour les 3 débits a 5microns

!    re=v*ettp(ip,jdp)/nu ! Reynolds number
!    we=((ruslag(1,1,iropt))*ettp(ip,jdp)*(v**2))/(32.d-3) ! weber number
!    oh=sqrt(we)/re ! Ohnesorge number
!    k=we*(oh**(-2/5)) ! ok avec k=15.5 et les v, mais difficile pour les débit moindre


!     energ = 0.5d0*ettp(ip,jmp)*v**2
!     energt   = 2.095d-13*we ! ne fonctionne pas très bien car c du tout ou rien, veleur comprise entre 2.09 et 2.1


! conditions
!===============================================================================
!
!  if (k .ge. 0.5) then
! the figure 1.2 is adapted from Horter works on particle rebounds/deposition/splashing in ONERA

  if ( energ .ge. energt )then

!  if (we .ge. 1.8) then


! The particle deposits:

    isuivi = 0
!    itepa(ip,jisor) = -itepa(ip,jisor)
    itepa(ip,jisor) = 0

    ettp(ip,jxp) = xk
    ettp(ip,jyp) = yk
    ettp(ip,jzp) = zk

    vitpar(ip,1) = 0.d0
    vitpar(ip,2) = 0.d0
    vitpar(ip,3) = 0.d0

    vitflu(ip,1) = 0.d0
    vitflu(ip,2) = 0.d0
    vitflu(ip,3) = 0.d0

  else

! The particle does not deposit:
! It 'rebounds' on the energy barrier:

  isuivi = 1
  itepa(ip,jisor) = ifabor(kface)

!  The mass flux is not calculated
!
    depch = 0

!-->Modification of the starting point

  ettpa(ip,jxp) = xk
  ettpa(ip,jyp) = yk
  ettpa(ip,jzp) = zk

  if (iensi1.eq.1) then
     nfin = 0
     call enslag                                                 &
          !==========
          ( idebia, idebra  ,                                        &
          nbpmax , nvp    , nvp1   , nvep   , nivep  ,             &
          nfin   , ip     ,                                        &
          itepa  ,                                                 &
          ettpa  , tepa   , ra)
  endif

  !-->Modification of the arrival point
  !   (the absolute value is intended to avoid the negative scalar products
  !  that may occur due to computer round-off error

  aa = 2.d0 * abs( (xq-xk)*xnn + (yq-yk)*ynn + (zq-zk)*znn )

  ettp(ip,jxp) = xq - aa*xnn
  ettp(ip,jyp) = yq - aa*ynn
  ettp(ip,jzp) = zq - aa*znn

  !--> Modification of the particle velocity at the arrival point

!-->Modification of the particle velocity at the impaction point
!   (like an elastic rebound)

    aa = abs(( vitpar(ip,1)*xnn                                   &
              +vitpar(ip,2)*ynn                                   &
              +vitpar(ip,3)*znn) )*2.d0

    vitpar(ip,1) = vitpar(ip,1) - aa*xnn
    vitpar(ip,2) = vitpar(ip,2) - aa*ynn
    vitpar(ip,3) = vitpar(ip,3) - aa*znn

  !--> Modification of the velocity of the flow seen at the arrival point

    aa = abs( (vitflu(ip,1)*xnn                                   &
             + vitflu(ip,2)*ynn                                   &
             + vitflu(ip,3)*znn) ) * 2.d0

  vitflu(ip,1) = vitflu(ip,1) - aa*xnn
  vitflu(ip,2) = vitflu(ip,2) - aa*ynn
  vitflu(ip,3) = vitflu(ip,3) - aa*znn


  endif