!******************** start file bmzini.for ; group bmzini *************
!..........................................
!
subroutine bmzini(ierr)
  !
  use nbi_com
  use nbi_idhdr
  use nbi_mpi_com
  use nbi_halo_mod_def
  !
  !
  implicit none
  integer, intent(out) :: ierr  ! completion code, 0=ok
  integer :: isb
  !
  !  beam code version of old transp routine
  !  get 2d mc grid data from xplasma...
  !  dmc 7 dec 2001
  !
  !--------------------------------------------------
  integer :: iznbmr,iznbmri
  integer :: iudsym,n_split
  !
  integer :: izones,jznbmr,iabort,i,ith,iwarn,ier
  !
  integer, parameter :: mjth = 40
  real(fp) :: z_split
  integer :: i_gen_req,i_gen,i_nseed

  !
  !
  !--------------------------------------------------
  !

  ierr = 0
  !
  call mcgrid_getsym(id_mcgrid,iudsym)
  if(iudsym.eq.1) then
     nlsym2b=.true.
  else
     nlsym2b=.false.
  end if

  !
  call mcgrid_getnumr(id_mcgrid,iznbmr,iznbmri)
  nxskpb=nzones/iznbmri
  nxtbzn=iznbmr-iznbmri

  !
  !  test...
  !
  izones = (nzones/nxskpb) * nxskpb
  if(izones.ne.nzones) then
     write(nonlin,9901) nzones,nxskpb
     ierr=1
     return
  end if
9901 format(' ?bmzini -- nzones=',i3,' and nxskpb=',i3/ &
       '   nxskpb must divide nzones without remainder for transp'/ &
       '   beam code indexing to work properly.')
  !
  dxitk = 1.0_fp/iznbmri
  nznbmr = iznbmr
  nznbmri = iznbmri
  jznbmr = nznbmr+lcm1
  !
  !  count the zones; double if not updown symmetric.
  !


  if(allocated(nthzsm)) deallocate(nthzsm)
  allocate(nthzsm(jznbmr))


  nthzsm(1:lcm1) = 0
  call mcgrid_getnumzns(id_mcgrid,nthzsm(lcentr:jznbmr))
  !
  !  convert to sum
  !
  do i=lcentr,jznbmr
     nthzsm(i)=nthzsm(i-1)+nthzsm(i)
  end do
  !
  nfbzns=nthzsm(jznbmr)
  nfbznsi=nthzsm(jznbmr-nxtbzn)
  nthzsm0=nthzsm(lcentr)



#ifdef __MPI
  !
  !--------------------------------------------------------
  ! *** allocating "track" data for a new time step on all cpu for a sublist
  ! of track data i.e. ~ ndep0/numprocd,
  ! at this point, if it's step>1 then cpu =0 allocated for all track list
  ! and cpu>0 is/are deallocated (see depall_dpt_buf_gather)
  !

  if(ndep0>numprocd) then
     mdep0 = ndep0/numprocd
     if(mod(ndep0,numprocd).ne.0) then
        mdep0=mdep0+1       ! round up
     end if
  else
     mdep0 =1
  end if
  ndep0_use=mdep0*numprocd

#else
  !--------------------------------------------------------
  ! ***if serial run  allocating "track" data for a first time step on  cpu=0 for all list
  ! of track data
  !
  mdep0=max(1,ndep0)
  ndep0_use=max(1,ndep0)
#endif

  !--------------------------------------------------------
  ! ***allocating "beam_in_box" data for a new time step on all cpu for a sublist
  ! of track data i.e. ~ ndepbox/numprocd,
  if(nbbox.gt.0) then
     if(ndepbox>numprocd) then
        mdepbox = ndepbox/numprocd
        if(mod(ndepbox,numprocd).ne.0) then
           mdepbox=mdepbox+1 ! round up
        end if
     else
        mdepbox =1
     end if
     ndepbox_use=mdepbox*numprocd
  end if

  if(mdepbox.gt.0) then
     call boxn0_setup(iwarn) ! beam-in-box calculation setup:  ignore error
     ! all cpus know
     ! set monte carlo statistics control for beam-in-box
     ! 'ndepbox' divisible by 'numproc'
     if(levmod_halo.eq.2) then
        ! nsplt # of spliting in all -- all generations are included
        ! # of 'geo' splitings for a given generation
        ! # of 'kin' splitings for a given generation
        if(nsplt.eq.0.or.nsplt_geo*nsplt_kin.le.0) then
           write(nonlin,*) ' ?bmzini:levmod_halo=2 but nsplt=0'
           write(nonlin,*) '          or nsplt_geo*nsplt_kin=0'
           ierr=1
           return
        end if
        if(nsplt.lt.0) then
           call track_halo_alloc(mdepbox) !'0' generation only
        else
           !     approximate # of halo generation
           z_split=1._fp*nsplt/(nsplt_geo*nsplt_kin)**nsplt_gen
           i_gen=log(z_split)/log(1.4_fp)+1 !empirical~1/cxrate=1.0/0.6
           i_gen=max(i_gen,1)
           i_gen_req=min(i_gen,nsplt_gen,8)
           z_split=1._fp
           !     approximate # of splittings
           if(nsplt_gen.eq.0) then
              do i=1,8      !0 through 8 th generation are included
                 z_split=z_split+(0.6_fp*nsplt_geo*nsplt_kin)**i  !cxrate=0.6
              end do
           else
              do i=1,i_gen_req !0 through 8 th generation are included
                 z_split=z_split+(nsplt_geo*nsplt_kin)**i
              end do
              if(i_gen_req.lt.i_gen) then
                 z_split=z_split+(nsplt_geo*nsplt_kin)** &
                      i_gen_req*(i_gen-i_gen_req)
              end if
           end if
           n_split=max(1.0,z_split)
           n_split=min(n_split,abs(nsplt))

           write(nonlin,*) '  3d halo:n_split=',n_split

           call track_halo_alloc(mdepbox*n_split) !call for 3d halo model
        end if
     end if
  end if

  !
  !--------------------------------------------------------
  !  *** if this is the first call ***
  !
  i_nseed=nseed+int(tbm1*1.0e+5_fp)

  call nbi_random_init(i_nseed,0,1)  ! init c. karney rng

  if(.not.nbi_init_flag) then

     ! set some array dimensions
     ! mdep0 = 0
     ! track data arays will be allocated only if nltrk_dep0 = .true.

     mimxbz=nfbzns
     mimxbzf=nfbznsi
     !
     !
     !

     call nbi_alloc2_init   ! allocate rest of nbi_com ...!


     !  nbi_alloc2_init sets nbi_init_flag

     delt = 0.0_fp
     dt_adj0 = 1.0_fp
     dt_adj0s = 1.0_fp

     if(ierr.ne.0) return

  else if(nbi_idhdr_skip) then

     mimxbz=nfbzns
     mimxbzf=nfbznsi

  end if
  !
  !--------------------------------------------------------
  !
  if(nfbzns.ne.mimxbz) then
     write(nonlin,*) '?bmzini -- nfbzns/mimxbz mismatch: ', &
          nfbzns,mimxbe
     ierr=1
  end if
  if(nfbznsi.ne.mimxbzf) then
     write(nonlin,*) '?bmzini -- nfbznsi/mimxbzf mismatch: ', &
          nfbznsi, mimxbzf
     ierr=1
  end if
  if(ierr.eq.1) return
  !
  xbmbnd = nznbmr * dxitk
  xminbm = xbmbnd
  !
  nznblo=int(xminbm/dxizc(lcentr) + 0.5_fp)

  if(allocated(xiblo)) deallocate(xiblo)
  allocate(xiblo(mtbl))

  xiblo(1:lcm1)=0.0_fp
  do i=1,mtbl-lcm1
     xiblo(lcm1+i)=(i-1)*dxizc(lcentr)
  end do
  xiblo(lep1)=xi(lep1,1)
  !
  !  a poloidal angle grid -- added dmc 6 june 1994
  !
  if(nlsym2b) then
     thbdy0=0.0_fp
     thbdy1=xpi
     nthzons=mjth/2
  else
     thbdy0=-xpi
     thbdy1=xpi
     nthzons=mjth
  end if
  !
  if(allocated(thzons)) deallocate(thzons)
  if(allocated(gxtabl)) deallocate(gxtabl)
  allocate(thzons(nthzons),gxtabl(nthzons,mtbl))
  thzons=0.0_fp
  !
  do ith=1,nthzons
     thzons(ith)=thbdy0 + &
          (thbdy1-thbdy0)*(real(ith)-0.5_fp)/real(nthzons)
  end do
  !
  xbjaco=0.0_fp  ! clear xbjaco
  !
  return
end subroutine bmzini
!******************** end file bmzini.for ; group bmzini ***************
subroutine nbi_shrink
  use nbi_com

  implicit none
  nznbmr=nzones/nxskpb
  nznblo=nzones
  xbmbnd=1.0_fp

end subroutine nbi_shrink
!------------------------------------------
subroutine nbi_rbloat
  use nbi_com

  implicit none
  integer :: iznbmri

  call mcgrid_getnumr(id_mcgrid,nznbmr,iznbmri)
  xbmbnd=xminbm
  nznblo=int(xminbm/dxizc(lcentr) + 0.5_fp)

end subroutine nbi_rbloat
!------------------------------------------

subroutine nbi_check_beam_species
  use iso_c_binding, only: fp => c_double
  use nbi_com

  implicit none
  integer :: isb

  mibs2use = 0
  nltrk_dep0 = nltrk_dep0.and.(ndep0.gt.0)
  if(nltrk_dep0) then
     mibs_use=0
     mtrk0=0
     do isb=1,nsbeam
        if(.not.nlfprod(isb).and.xninjs(isb).gt.0.0_fp) then !do not allow gather track data
           !if beam is 'off'
           mibs_use=mibs_use+1
           mibs2use(isb) = mibs_use
        end if
     end do
  end if
  nltrk_dep0 = nltrk_dep0 .and.(mibs_use.gt.0)
  if(nltrk_dep0) then
     call nbi_alloc_track
     !the rest of trak data  will be allocated in dep.for (call nbi_alloc_track_zone)

  else
     call nbi_dalloc_track
     call nbi_dalloc_track_zone
  end if
 
  return
end subroutine nbi_check_beam_species

! 11sep2004 fgtok
