!------------------------------------------------------------------------------|
!------------------------------------------------------------------------------|
!                                                                              |
!              ||===\\                                                         | 
!              ||    \\                                                        |
!              ||     ||   //==\\   ||  ||   //==||  ||/==\\                   |
!              ||     ||  ||    ||  ||  ||  ||   ||  ||    ||                  |
!              ||    //   ||    ||  ||  ||  ||   ||  ||                        |
!              ||===//     \\==//    \\==\\  \\==\\  ||                        |
!                                                                              |
!------------------------------------------------------------------------------|
!------------------------------------------------------------------------------|
!                                                                              |
!              DEFINE_BC    Apr. 2007                                          |
!                                                                              |
!------------------------------------------------------------------------------|
!------------------------------------------------------------------------------|

subroutine define_bc (params,osolve,vo)

!------------------------------------------------------------------------------|
!(((((((((((((((( Purpose of the routine  ))))))))))))))))))))))))))))))))))))))
!------------------------------------------------------------------------------|
! subroutine to implement the boundary conditions.
! Here the user should define if the dof is fixed by setting
! kfix((inode-1)*3+idof) to 1 (node inode and dof idof)
! and setting the corresponding value of u, v, or w to the  set fixed value
! Same operation for kfixt and temp but for a single dof
! vo is the structure that contains the void information


!------------------------------------------------------------------------------|
!((((((((((((((((  declaration of the subroutine arguments  ))))))))))))))))))))
!------------------------------------------------------------------------------|

use definitions

implicit none

type (parameters) params
type (octreesolve) osolve
type (void) vo

!------------------------------------------------------------------------------|
!(((((((((((((((( declaration of the subroutine internal variables )))))))))))))
!------------------------------------------------------------------------------|

integer i,iproc,nproc,ierr
double precision eps

!------------------------------------------------------------------------------|
!------------------------------------------------------------------------------|

INCLUDE 'mpif.h'

call mpi_comm_size (mpi_comm_world,nproc,ierr)
call mpi_comm_rank (mpi_comm_world,iproc,ierr)

eps=1.d-10
osolve%kfix=0
osolve%kfixt=0
      
select case(trim(params%infile))

case('input.txt')
 
   do i=1,osolve%nnode
      if (osolve%x(i).lt.eps) then
         osolve%kfix((i-1)*3+1)=1 ; osolve%u(i)=0.d0
         !osolve%kfix((i-1)*3+2)=1 ; osolve%v(i)=0.d0
         !osolve%kfix((i-1)*3+3)=1 ; osolve%w(i)=0.d0
      endif
      if (osolve%x(i).gt.1.d0-eps) then
         osolve%kfix((i-1)*3+1)=1 ; osolve%u(i)=0.d0
         !osolve%kfix((i-1)*3+2)=1 ; osolve%v(i)=0.d0
         !osolve%kfix((i-1)*3+3)=1 ; osolve%w(i)=0.d0
      endif
      if (osolve%y(i).lt.eps) then
         !osolve%kfix((i-1)*3+1)=1 ; osolve%u(i)=0.d0
         osolve%kfix((i-1)*3+2)=1 ; osolve%v(i)=0.d0
         if (osolve%x(i).lt.0.5001d0) osolve%v(i)=1.0d0
         osolve%kfix((i-1)*3+3)=1 ; osolve%w(i)=0.d0
      endif
      if (osolve%y(i).gt.1.d0-eps) then
         !osolve%kfix((i-1)*3+1)=1 ; osolve%u(i)=0.d0
         osolve%kfix((i-1)*3+2)=1 ; osolve%v(i)=0.d0
         !osolve%kfix((i-1)*3+3)=1 ; osolve%w(i)=0.d0
      endif
      if (osolve%z(i).lt.eps) then
         osolve%kfix((i-1)*3+1)=1 ; osolve%u(i)=0.d0
         osolve%kfix((i-1)*3+2)=1 ; osolve%v(i)=0.d0
         osolve%kfix((i-1)*3+3)=1 ; osolve%w(i)=0.d0
         if (osolve%x(i).lt.0.5001d0.and.osolve%y(i).lt.0.5001d0) osolve%v(i)=1.0d0
         osolve%kfixt(i)=1
         osolve%temp(i)=1.d0
      endif
      if (osolve%z(i).gt.1.d0-eps) then
         !osolve%kfix((i-1)*3+1)=1 ; osolve%u(i)=0.d0
         !osolve%kfix((i-1)*3+2)=1 ; osolve%v(i)=0.d0
!         osolve%kfix((i-1)*3+3)=1 ; osolve%w(i)=sin(osolve%y(i)*2.*3.141592654)*1.d-3
         osolve%kfixt(i)=1
         osolve%temp(i)=0.d0
      endif
      if (.not.vo%influid(i)) then
         osolve%kfixt(i)=1
         osolve%temp(i)=0.d0
      endif
   enddo

case ('input.jgr')
   call define_bc_jgr (osolve%nnode,osolve%kfix,osolve%kfixt,osolve%x,osolve%y,osolve%z,osolve%u,osolve%v,osolve%w,osolve%temp,vo) 

case ('input.sphere','input.cube')
   call define_bc_sphere (osolve%nnode,osolve%kfix,osolve%kfixt,osolve%x,osolve%y,osolve%z,osolve%u,osolve%v,osolve%w,osolve%temp,vo) 

case ('input.3Dpunch')
   call define_bc_3Dpunch (osolve%nnode,osolve%kfix,osolve%kfixt,osolve%x,osolve%y,osolve%z,osolve%u,osolve%v,osolve%w,osolve%temp,vo) 

case ('input.2Dpunch')
   call define_bc_2Dpunch (osolve%nnode,osolve%kfix,osolve%kfixt,osolve%x,osolve%y,osolve%z,osolve%u,osolve%v,osolve%w,osolve%temp,vo) 

case ('input.2Dpunch_vert')
   call define_bc_2Dpunch_vert (osolve%nnode,osolve%kfix,osolve%kfixt,osolve%x,osolve%y,osolve%z,osolve%u,osolve%v,osolve%w,osolve%temp,vo) 

case ('input.folding')
   call define_bc_folding (osolve%nnode,osolve%kfix,osolve%kfixt,osolve%x,osolve%y,osolve%z,osolve%u,osolve%v,osolve%w,osolve%temp,vo) 

case ('input.pipo')
   call define_bc_pipo (osolve%nnode,osolve%kfix,osolve%kfixt,osolve%x,osolve%y,osolve%z,osolve%u,osolve%v,osolve%w,osolve%temp,vo) 

case ('input.ritske')
   call define_bc_ritske (osolve%nnode,osolve%kfix,osolve%kfixt,osolve%x,osolve%y,osolve%z,osolve%u,osolve%v,osolve%w,osolve%temp,vo) 

case ('input.ritske_isurf')
   call define_bc_ritske_isurf (osolve%nnode,osolve%kfix,osolve%kfixt,osolve%x,osolve%y,osolve%z,osolve%u,osolve%v,osolve%w,osolve%temp,vo) 

case ('input.subduction')
   call define_bc_subduction (osolve%nnode,osolve%kfix,osolve%kfixt,osolve%x,osolve%y,osolve%z,osolve%u,osolve%v,osolve%w,osolve%temp,vo) 

case ('input.riedel')
   call define_bc_riedel (osolve%nnode,osolve%kfix,osolve%kfixt,osolve%x,osolve%y,osolve%z,osolve%u,osolve%v,osolve%w,osolve%temp,vo) 

case ('input.parallipipede')
   call define_bc_parallipipede (osolve%nnode,osolve%kfix,osolve%kfixt,osolve%x,osolve%y,osolve%z,osolve%u,osolve%v,osolve%w,osolve%temp,vo) 

case ('input.model1','input.model2')
   call define_bc_model1 (osolve%nnode,osolve%kfix,osolve%kfixt,osolve%x,osolve%y,osolve%z,osolve%u,osolve%v,osolve%w,osolve%temp,vo) 

case default
   if (iproc.eq.0) print *,params%infile
   call stop_run (' pb with input file in define_bc.f90')

end select


if (params%debug==2) call output_bc (osolve)

return
end
!------------------------------------------------------------------------------|
!------------------------------------------------------------------------------|


