text:
module gsw_mod_netcdf

use netcdf

use gsw_mod_kinds

integer, private :: ncid

public :: ncdf_open
public :: ncdf_close
public :: ncdf_get_dim
public :: ncdf_get_var
public :: ncdf_get_var_att

private :: ncdf_handle_err

contains

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

    subroutine ncdf_open (file_name)
    implicit none

    character (*), intent(in) :: file_name

    integer :: istat

    character (*), parameter :: fname = 'ncdf_open'

    istat = nf90_open(file_name, nf90_nowrite, ncid)
    if (istat /= nf90_noerr) call ncdf_handle_err(istat,fname, &
                                                  'nf90_open',file_name)

    return
    end subroutine ncdf_open

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

    subroutine ncdf_close ()
    implicit none

    integer :: istat

    character (*), parameter :: fname = 'ncdf_close'

    istat = nf90_close(ncid)
    if (istat /= nf90_noerr) call ncdf_handle_err(istat,fname,'nf90_close')

    return
    end subroutine ncdf_close

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

    integer function ncdf_get_dim (dim_name)
    implicit none

    character (*), intent(in) :: dim_name

    integer :: istat, dimid

    character (*), parameter :: fname = 'ncdf_get_dim'

    istat = nf90_inq_dimid(ncid, dim_name, dimid)
    if (istat /= nf90_noerr) call ncdf_handle_err(istat,fname, &
                                                  'nf90_inq_dimid',dim_name)

    istat = nf90_inquire_dimension(ncid, dimid, len=ncdf_get_dim)
    if (istat /= nf90_noerr) call ncdf_handle_err(istat,fname, &
                                             'nf90_inquire_dimension',dim_name)
    return
    end function ncdf_get_dim

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

    subroutine ncdf_get_var (var_name, var0, var1, var2, var3)
    implicit none

    character (*), intent(in) :: var_name
    real (r8), intent(out), optional :: var0, var1(:), var2(:,:), var3(:,:,:)

    integer :: istat, varid

    character (*), parameter :: fname = 'ncdf_get_var'

    istat = nf90_inq_varid(ncid, var_name, varid)
    if (istat /= nf90_noerr) call ncdf_handle_err(istat,fname, &
                                                  'nf90_inq_varid',var_name)

    if (present(var0)) then

       istat = nf90_get_var(ncid, varid, var0)
       if (istat /= nf90_noerr) call ncdf_handle_err(istat,fname, &
                                                     'nf90_get_var',var_name)

    else if (present(var1)) then

       istat = nf90_get_var(ncid, varid, var1)
       if (istat /= nf90_noerr) call ncdf_handle_err(istat,fname, &
                                                     'nf90_get_var',var_name)

    else if (present(var2)) then

       istat = nf90_get_var(ncid, varid, var2)
       if (istat /= nf90_noerr) call ncdf_handle_err(istat,fname, &
                                                     'nf90_get_var',var_name)

    else if (present(var3)) then

       istat = nf90_get_var(ncid, varid, var3)
       if (istat /= nf90_noerr) call ncdf_handle_err(istat,fname, &
                                                     'nf90_get_var',var_name)

    end if

    return
    end subroutine ncdf_get_var

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

    subroutine ncdf_get_var_att (var_name, var, att_name, att)
    implicit none

    character (*), intent(in) :: var_name, att_name
    real (r8), intent(out) :: var(:,:), att

    integer :: istat, varid

    character (*), parameter :: fname = 'ncdf_get_var_att'

    istat = nf90_inq_varid(ncid, var_name, varid)
    if (istat /= nf90_noerr) call ncdf_handle_err(istat,fname, &
                                                  'nf90_inq_varid',var_name)

    istat = nf90_get_var(ncid, varid, var)
    if (istat /= nf90_noerr) call ncdf_handle_err(istat,fname, &
                                                  'nf90_get_var',var_name)

    istat = nf90_get_att(ncid, varid, att_name, att)
    if (istat /= nf90_noerr) call ncdf_handle_err(istat,fname, &
                                                  'nf90_get_att',att_name)

    return
    end subroutine ncdf_get_var_att

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

    subroutine ncdf_handle_err (status, fname, nf90_name, var_name)
    implicit none

    integer, intent(in) :: status
    character (*), intent(in) :: fname, nf90_name
    character (*), intent(in), optional :: var_name

    if (status == nf90_noerr) return

    if (present(var_name)) then
        print*, '**Error return from ', nf90_name, ' (', var_name, ') in ', &
	        fname, ' ...'
    else
        print*, '**Error return from ', nf90_name, ' in ', fname, ' ...'
    end if
    print '(3x,a," (status=",i0,")")', trim(nf90_strerror(status)), status

    stop
    end subroutine ncdf_handle_err

end module gsw_mod_netcdf

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

max_stars_repo_name: minsukji/ci-debug

max_issues_repo_name: minsukji/ci-debug

max_forks_repo_name: minsukji/ci-debug


text:
C*************************************************************
C     MINIMIZATION OF THE NONLINEAR FUNCTIONAL FUNC(U)
C     IN GIVEN DIRECTION  ZZ
C*************************************************************
      Subroutine minim(
C group(F)
     &       nU, U, ZZ, Lamda1, Lamda2, fMin, U1,
     &       icnt, rMove, flagResult,
C group (ANI)
     &       XYP, IPE, IEE, HesP, hStar, status,
     &       lE, iEs, XYPs, IPEs, detGs, HesPs, qEs,
     &       nPw, nEw, XYPw, HesPw, IPEw,
     &       MetricFunction, flagAnalytic,
     &       iSE, rSE, iP1, iFNCs, calCrv,
     &       L1Et, L2Et, tE,
     &       nL2t, nStept, nEt, nCrvFnc, LFnc, ILt,
C group (Tangle)
     &       flagTM, nBad)
C ===========================================================
      include 'makS.fd'
      include 'status.fd'
C ===========================================================
C group (F)
      real  U(*), ZZ(*), U1(*)
      real  Lamda1, Lamda2, fMin, rMove
      Logical flagResult

C group (ANI)
      Integer IPE(3, *), IEE(3, *)
      Integer iEs(*), IPEs(3, *), IPEw(3, *), iSE(*), status
      real  XYP(2, *),  HesP(3, *), hStar

      real  XYPs(2, *), HesPs(*),   detGs, qEs(*)
      real  XYPw(2, *), HesPw(3, *), rSE(*)
      real  tE(*)

      EXTERNAL calCrv

      Integer L1Et(2, *), L2Et(*)
      Integer nL2t(*), nStept(4, *), nEt(*)
      Integer LFnc(*), ILt(*)

      Logical flagAnalytic

      Integer  MetricFunction
      EXTERNAL MetricFunction

C group (Tangle)
      Logical flagTM
      Integer nBad

C group (Local variables)
      real  NLnFnc
      real  Lamda, dLamda, step, x, f0, f1
      real  aStep

      real  qEt(MaxS), XYPt(2), HesPt(3), detGt
      real  XYPu(2)

      Integer iOs(MaxS)
      Logical flagOrient

C ===========================================================
      flagResult = .FALSE.

      x = 0D0
      Do i = 1, nU
         x = x + ZZ(i) ** 2
      End do
      x = sqrt(x)

      Do i = 1, nU
         ZZ(i) = ZZ(i) / x
      End do

      x=Lamda1
      Lamda1 = min(x, Lamda2)
      Lamda2 = max(x, Lamda2)

      icnt = 0
      step = 5D-3
      aStep = 0D0
      dLamda = Lamda2 - Lamda1

      f0 = fMin

      Call copySQ(lE, qEs, XYPs, HesPs, detGs,
     &                qEt, XYPt, HesPt, detGt)

 10   Lamda = step  * dLamda
      If(Lamda.GT.(1D0 - aStep - step) * dLamda) Then
         f1 = f0 + 1D0
      Else
         Do i = 1, nU
            U1(i) = U(i) + Lamda * ZZ(i)
         End do

         icnt = icnt + 1
         f1 = NLnFnc(nU, U1,
C group (ANI)
     &        XYP, IPE, IEE, HesP, hStar, status,
     &        lE, iEs, XYPs, IPEs, HesPs, detGs, qEs,
     &        nPw, nEw, XYPw, HesPw, IPEw,
     &        MetricFunction, flagAnalytic,
     &        iSE, rSE, iP1, iFNCs, calCrv,
     &        L1Et, L2Et, tE,
     &        nL2t, nStept, nEt, nCrvFnc, LFnc, ILt)

c  ...  checking for inverted elements
         If(flagTM) Then
            Do i = 1, 2
               XYPu(i) = XYP(i, iP1)
               XYP(i, iP1) = XYPs(i, 1)
            End do

            Do n = 1, lE
               If(qEs(n).GT.0D0) Then
                  Call updQb(n, lE, iEs, XYP, IPEs, qEs)
               End if
            End do

            Do i = 1, 2
               XYP(i, iP1) = XYPu(i)
            End do

            mBad = 0
            Do n = 1, lE
               If(qEs(n).LE.0D0) Then
                  mBad = mBad + 1
                  f1 = max(f1, 1D0 - qEs(n))
               End if
            End do
            If(mBad.GT.nBad) f1 = f0 + 1D0

c  ...  checking orientation of triangles
         Else
            Call calSO(XYP, IPE, lE, iEs, iOs)
            Call chkSO(iP1, iP1, XYPs, XYP, IPE, lE, iEs, iOs,
     &                 flagOrient)
            If(.NOT.flagOrient) f1 = f0 + 1D0
         End if
      End if

      If(f1.LT.f0) Then
         f0 = f1
         Do i = 1, nU
            U(i) = U1(i)
         End do

         flagResult = .TRUE.

         Call copySQ(lE, qEs, XYPs, HesPs, detGs,
     &                   qEt, XYPt, HesPt, detGt)

         aStep = aStep + step
         step = 4D0 * step
      Else
         step = step / 4D0
         If(step.LT.1D-5) Then
            fMin = f0
            rMove = aStep * dLamda
            goto 1000
         End if
      End if
      goto 10

 1000 Call copySQ(lE, qEt, XYPt, HesPt, detGt,
     &                qEs, XYPs, HesPs, detGs)

 9000 Return
      End

ext: f


module MpiExchangeGwfModule

  use KindModule, only: DP, I4B
  use ConstantsModule, only: LENMODELNAME, LENMEMPATH, LENVARNAME, LENFTYPE,    &
                             LENPACKAGENAME
  use ArrayHandlersModule, only: ifind
  use MpiExchangeGenModule, only: serialrun
  use MpiWrapper, only: mpiwrpallgather, mpiwrpbarrier,  mpiwrptypefree
  use MpiWrapperMemory, only: ColMemoryType, mpiwrpcolstruct, mpiwrpallgathervcol
  use MpiExchangeModule, only: MpiWorld, mpi_to_colmem
  use SimModule, only: store_error, store_error_unit, ustop ! @@@@@ debug
  !
  use MawModule, only: mawftype => ftype
  use SfrModule, only: sfrftype => ftype
  use LakModule, only: lakftype => ftype
  use UzfModule, only: uzfftype => ftype
  !
  implicit none

  integer, parameter :: imaw = 1
  integer, parameter :: isfr = 2
  integer, parameter :: ilak = 3
  integer, parameter :: iuzf = 4
  integer, parameter :: mxpack = iuzf

  character(LENFTYPE), dimension(2,mxpack), save :: packftype

  private

  ! -- Public functions
  public :: mpi_gwfhalo_world
  !public :: mpi_set_gwfhalo_world
  public :: mpi_set_gwfhalo_world_var_int
  public :: mpi_set_gwfhalo_world_var_int1d
  public :: mpi_set_gwfhalo_world_dis
  public :: mpi_set_gwfhalo_world_mvr

  !interface mpi_set_gwfhalo_world_var
  !  module procedure mpi_set_gwfhalo_world_var_int,                             &
  !                   mpi_set_gwfhalo_world_var_int1d
  !end interface mpi_set_gwfhalo_world_var

  contains

  subroutine mpi_gwfhalo_world(iopt)
! ******************************************************************************
! This subroutine gathers DIS information for all models.
! ******************************************************************************
!
!    SPECIFICATIONS:
! ------------------------------------------------------------------------------
    ! -- modules
    use MpiExchangeColModule, only: ciopt, n_recv, cmt_recv
    use MemoryManagerModule, only: mem_setptr, mem_get_ptr, mem_setval
    use MemoryHelperModule, only: split_mem_path_ff
    use MemoryTypeModule, only: MemoryType
    use BaseModelModule, only: BaseModelType, GetBaseModelFromList
    use ListsModule, only: basemodellist
    use NumericalModelModule, only: NumericalModelType
    use BndModule, only: BndType, GetBndFromList
    use GwfModule, only: GwfModelType, cunit, niunit

    ! -- dummy
    integer, intent(in) :: iopt
    ! -- local
    character(len=1) :: cdum
    character(len=4) :: dis_type
    integer :: im, ipo, ip, is, iact
    class(BaseModelType), pointer :: mb
    class(NumericalModelType), pointer :: mp
    class(GwfModelType), pointer :: gp
    class(BndType), pointer :: packobj
    type(MemoryType), pointer :: mt
    !
    integer, parameter :: cntypes = 6
    integer :: i, n, n_send, newtype
    integer, dimension(1) :: iwrk
    integer, dimension(:), allocatable :: recvcounts, displs
    type(ColMemoryType), dimension(:), allocatable :: cmt_send
    character(len=LENFTYPE) :: ft, cft
    integer :: rank
! ------------------------------------------------------------------------------
    if (serialrun) then
      return
    end if
    !
    ! -- initialize package filetype
    packftype(1,imaw) = trim(mawftype)
    packftype(1,isfr) = trim(sfrftype)
    packftype(1,ilak) = trim(lakftype)
    packftype(1,iuzf) = trim(uzfftype)
    !
    do ip = 1, mxpack
      ft = packftype(1,ip)
      n = len_trim(ft)
      do i = 1, niunit
        cft = cunit(i)
        if (ft(1:n) == cft(1:n)) then
          packftype(2,ip) = trim(cunit(i))
        end if
      end do
    end do
    !
    allocate(cmt_send(1))
    do iact = 1, 2
      is = 0
      do im = 1, basemodellist%Count()
        mb => GetBaseModelFromList(basemodellist, im)
        select type (mb)
        class is (NumericalModelType)
          mp => mb
        end select
        select type (mb)
        class is (GwfModelType)
          gp => mb
        end select
        call split_mem_path_ff(mp%dis%memoryPath, cdum, dis_type)
        select case (dis_type)
          case ('DIS')
            if (iopt == 1) then
              call mem_get_ptr('DNDIM', mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
            else
              call mem_get_ptr('DNDIM', mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
              call mem_get_ptr('NLAY', mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
              call mem_get_ptr('NROW', mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
              call mem_get_ptr('NCOL', mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
              call mem_get_ptr('MSHAPE', mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
              call mem_get_ptr('NODESUSER', mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
            end if
          case ('DISV')
            if (iopt == 1) then
              call mem_get_ptr('DNDIM',  mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
            else
              call mem_get_ptr('DNDIM',  mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
              call mem_get_ptr('NLAY',  mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
              call mem_get_ptr('NCPL',  mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
              call mem_get_ptr('NVERT', mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
              call mem_get_ptr('MSHAPE', mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
              call mem_get_ptr('NODESUSER', mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
            end if
          case ('DISU')
            if (iopt == 1) then
              call mem_get_ptr('DNDIM', mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
            else
              call mem_get_ptr('DNDIM', mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
              call mem_get_ptr('NVERT', mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
              call mem_get_ptr('MSHAPE', mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
              call mem_get_ptr('NODESUSER', mp%dis%memoryPath, mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
            end if
        end select
        if (iopt == 2) then
          ! -- Mover package data
          do ipo=1,mp%bndlist%Count()
            packobj => GetBndFromList(mp%bndlist, ipo)
            if (packobj%filtyp == packftype(1,imaw)) then
              call mem_get_ptr('NMAWWELLS', trim(packobj%memoryPath), mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
              call mem_get_ptr('IMOVER', trim(packobj%memoryPath), mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
            end if
            if (packobj%filtyp == packftype(1,isfr)) then
              call mem_get_ptr('MAXBOUND', trim(packobj%memoryPath), mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
              call mem_get_ptr('IMOVER', trim(packobj%memoryPath), mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
            end if
            if (packobj%filtyp == packftype(1,ilak)) then
              call mem_get_ptr('NOUTLETS', trim(packobj%memoryPath), mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
              call mem_get_ptr('NLAKES', trim(packobj%memoryPath), mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
              call mem_get_ptr('IMOVER', trim(packobj%memoryPath), mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
            end if
            if (packobj%filtyp == packftype(1,iuzf)) then
              call mem_get_ptr('MAXBOUND', trim(packobj%memoryPath), mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
              call mem_get_ptr('IMOVER', trim(packobj%memoryPath), mt)
              call mpi_to_colmem(mt, is, cmt_send, iact)
            end if
          end do
          call mem_get_ptr('INNPF', mp%memoryPath, mt)
          call mpi_to_colmem(mt, is, cmt_send, iact)
          call mem_get_ptr('INIC', mp%memoryPath, mt)
          call mpi_to_colmem(mt, is, cmt_send, iact)
        endif
        if (iopt == 3) then
          ! -- NPF related
          call mem_get_ptr('IK22', gp%npf%memoryPath, mt)
          call mpi_to_colmem(mt, is, cmt_send, iact)
          call mem_get_ptr('IK33', gp%npf%memoryPath, mt)
          call mpi_to_colmem(mt, is, cmt_send, iact)
          call mem_get_ptr('IANGLE1', gp%npf%memoryPath, mt)
          call mpi_to_colmem(mt, is, cmt_send, iact)
          call mem_get_ptr('IANGLE2', gp%npf%memoryPath, mt)
          call mpi_to_colmem(mt, is, cmt_send, iact)
          call mem_get_ptr('IANGLE3', gp%npf%memoryPath, mt)
          call mpi_to_colmem(mt, is, cmt_send, iact)
          call mem_get_ptr('IREWET', gp%npf%memoryPath, mt)
          call mpi_to_colmem(mt, is, cmt_send, iact)
        end if
      end do
      if (iact == 1) then
        n_send = is
        if (allocated(cmt_send)) then
          deallocate(cmt_send)
        end if
        allocate(cmt_send(max(n_send,1))) !PAR
      end if
    end do
    !
    ! -- gather sizes
    allocate(recvcounts(MpiWorld%nrproc), displs(MpiWorld%nrproc))
    iwrk(1) = n_send
    call mpiwrpallgather(MpiWorld%comm, iwrk, 1, recvcounts, 1)
    n_recv = sum(recvcounts)
    if (allocated(cmt_recv)) then
      deallocate(cmt_recv)
    end if
    allocate(cmt_recv(n_recv))
    displs(1) = 0
    do i = 2, MpiWorld%nrproc
      displs(i) = displs(i-1) + recvcounts(i-1)
    end do
    ! -- create derived type
    call mpiwrpcolstruct(newtype)
    ! -- gather data
    call mpiwrpallgathervcol(MpiWorld%comm, cmt_send, n_send, newtype, cmt_recv,   &
                             recvcounts, newtype, displs)
    ! -- DEBUG
    if (.false.) then
      do rank = 0, MpiWorld%nrproc-1
        if (MpiWorld%myrank == rank) then
          write(*,*) '=== send myrank', rank
          do i = 1, n_send
            write(*,'(a,1x,a,1x,i0)') trim(cmt_send(i)%name), trim(cmt_send(i)%path), cmt_send(i)%intsclr
          end do
          write(*,*) '=== recv myrank', rank
          do i = 1, n_recv
            write(*,'(a,1x,a,1x,i0)') trim(cmt_recv(i)%name), trim(cmt_recv(i)%path), cmt_recv(i)%intsclr
          end do
        end if
        call MpiWorld%mpi_barrier()
      end do
      call ustop('@@@@')
    end if
    !
    ! -- clean up
    call mpiwrptypefree(newtype)
    deallocate(cmt_send, recvcounts, displs)
    ciopt = iopt
    !
    ! -- return
    return
  end subroutine mpi_gwfhalo_world

  function mpi_set_gwfhalo_world_var_int(tgtname, tgtmempath) result(aint)
! ******************************************************************************

! ******************************************************************************
!
!    SPECIFICATIONS:
! ------------------------------------------------------------------------------
    ! -- modules
    use MpiExchangeColModule, only: n_recv, cmt_recv
    use MemoryTypeModule, only: iintsclr
    ! -- dummy
    character(len=*), intent(in) :: tgtname
    character(len=*), intent(in) :: tgtmempath
    ! -- return
    integer(I4B) :: aint
    ! -- local
    integer(I4B) :: i
    logical :: lok
    character(len=LENVARNAME) :: name        !name of the array
    character(len=LENMEMPATH) :: mempath     !memory path
! ------------------------------------------------------------------------------
    !
    if (serialrun) then
      return
    end if
    !
    ! -- Set other variables
    aint = -1
    lok = .false.
    do i = 1, n_recv
      name = cmt_recv(i)%name
      read(cmt_recv(i)%path,*) mempath
      !write(*,*) '@@@@ name, path: "'//trim(cmt_recv(i)%name)//'" "'//trim(cmt_recv(i)%path)//'"'
      if ((trim(name) == trim(tgtname)) .and. (trim(mempath) == trim(tgtmempath))) then
        lok = .true.
        if (cmt_recv(i)%memitype /= iintsclr) then
          lok = .false.
        else
          aint = cmt_recv(i)%intsclr
        end if
        exit
      end if
    end do
    !
    ! -- check
    if (.not.lok) then
      call store_error('Program error: mpi_set_gwfhalo_world_var_int')
      write(*,*) '@@@ not found: "',trim(tgtname),'" "',trim(tgtmempath),'"'
      call ustop()
    end if
    ! -- return
    return
  end function mpi_set_gwfhalo_world_var_int

  function mpi_set_gwfhalo_world_var_int1d(tgtname, tgtmempath) result(aint1d)
! ******************************************************************************

! ******************************************************************************
!
!    SPECIFICATIONS:
! ------------------------------------------------------------------------------
    ! -- modules
    use MpiExchangeColModule, only: n_recv, cmt_recv
    use MemoryTypeModule, only: iaint1d
    ! -- dummy
    character(len=*), intent(in) :: tgtname
    character(len=*), intent(in) :: tgtmempath
    ! -- return
    integer(I4B), dimension(3) :: aint1d
    ! -- local
    integer(I4B) :: i
    logical :: lok
    character(len=LENVARNAME) :: name        !name of the array
    character(len=LENMEMPATH) :: mempath     !memory path
! ------------------------------------------------------------------------------
    !
    if (serialrun) then
      return
    end if
    !
    ! -- Set other variables
    aint1d = -1
    lok = .false.
    do i = 1, n_recv
      name = cmt_recv(i)%name
      read(cmt_recv(i)%path,*) mempath
      if ((trim(name) == trim(tgtname)) .and. (trim(mempath) == trim(tgtmempath))) then
        lok = .true.
        if (cmt_recv(i)%memitype /= iaint1d) then
          lok = .false.
        else
          aint1d(1) = cmt_recv(i)%aint1d(1)
          aint1d(2) = cmt_recv(i)%aint1d(2)
          aint1d(3) = cmt_recv(i)%aint1d(3)
        end if
        exit
      end if
    end do
    !
    ! -- check
    if (.not.lok) then
      call store_error('Program error: mpi_set_gwfhalo_world_var_int')
      write(*,*) '@@@ not found: "',trim(tgtname),'" "',trim(tgtmempath),'"'
      call ustop()
    end if
    !
    ! -- return
    return
  end function mpi_set_gwfhalo_world_var_int1d
  !
  subroutine mpi_set_gwfhalo_world_dis(mname, m2)
! ******************************************************************************
! ******************************************************************************
!
!    SPECIFICATIONS:
! ------------------------------------------------------------------------------
    ! -- modules
    use GwfDisModule,          only: dis_cr,  GwfDisType !,dis_init_mem TODO
    use GwfDisvModule,         only: disv_cr, GwfDisvType
    use GwfDisuModule,         only: disu_cr, GwfDisuType
    use NumericalModelModule,   only: NumericalModelType
    use BaseDisModule,         only: DisBaseType !PAR
    use MpiExchangeColModule,  only: mpi_get_distype
    use MemoryManagerModule,   only: mem_setval, mem_allocate, mem_setval
    use MemoryHelperModule, only: create_mem_path
    use ConstantsModule, only: IZERO, DZERO, LENMEMPATH
    ! -- dummy
    character(len=*), intent(in) :: mname
    class(NumericalModelType), pointer, intent(out) :: m2
    ! -- local
    integer(I4B) :: in_dum, iout_dum
    integer(I4B) :: ndim, nodesuser
    integer(I4B), dimension(3) :: mshape
    logical ::  ldis, ldisu, ldisv
    character(len=LENMEMPATH) :: disMemPath
    class(DisBaseType), pointer :: bp
    class(GwfDisType), pointer :: disp
    class(GwfDisvType), pointer :: disvp
! ------------------------------------------------------------------------------
    !
    if (serialrun) then
      return
    end if
    !
    disMemPath = 'M2'
    !
    allocate(m2)
    allocate(m2%dis)
    in_dum = -1
    iout_dum = -1
    call mpi_get_distype(mname, ldis, ldisu, ldisv)
    if (ldis) then
      call dis_cr(m2%dis, disMemPath, in_dum, iout_dum)
      disMemPath = create_mem_path(trim(disMemPath), 'DIS')
      bp => m2%dis
      select type (bp)
      class is (GwfDisType)
        disp => bp
      end select
      call mem_allocate(disp%nodereduced, 1, 'NODEREDUCED', disMemPath)
      call mem_allocate(disp%idomain, 1, 1, 1, 'IDOMAIN', disMemPath)
      call mem_allocate(disp%top2d, 1, 1, 'TOP2D', disMemPath)
      call mem_allocate(disp%bot3d, 1, 1, 1, 'BOT3D', disMemPath)
    endif
    if (ldisu) then
      call disu_cr(m2%dis, disMemPath, in_dum, iout_dum)
      disMemPath = create_mem_path(trim(disMemPath), 'DISU')
    endif
    if (ldisv) then
      call disv_cr(m2%dis, disMemPath, in_dum, iout_dum)
      disMemPath = create_mem_path(trim(disMemPath), 'DISV')
      bp => m2%dis
      select type (bp)
      class is (GwfDisvType)
        disvp => bp
      end select
      call mem_allocate(disvp%nodereduced, 1, 'NODEREDUCED', disMemPath)
      call mem_allocate(disvp%idomain, 1, 1, 1, 'IDOMAIN', disMemPath)
      call mem_allocate(disvp%top2d, 1, 1, 'TOP2D', disMemPath)
      call mem_allocate(disvp%bot3d, 1, 1, 1, 'BOT3D', disMemPath)
    endif
    !
    ndim = mpi_set_gwfhalo_world_var_int('DNDIM', mname)
    call mem_setval(ndim, 'NDIM', disMemPath)
    !
    nodesuser = mpi_set_gwfhalo_world_var_int('NODESUSER', mname)
    call mem_setval(nodesuser, 'NODESUSER', disMemPath)
    call mem_setval(nodesuser, 'NODES', disMemPath)
    !
    mshape = mpi_set_gwfhalo_world_var_int1d('MSHAPE', mname)
    call mem_allocate(m2%dis%mshape, ndim, 'MSHAPE', disMemPath)
    call mem_setval(mshape, 'MSHAPE', disMemPath)
    !
    ! -- return
    return
  end subroutine mpi_set_gwfhalo_world_dis

  subroutine mpi_set_gwfhalo_world_mvr()
! ******************************************************************************
! This subroutine sets halo model packages for mover
! ******************************************************************************
!
!    SPECIFICATIONS:
! ------------------------------------------------------------------------------
    ! -- modules
    use MemoryTypeModule, only: iintsclr, iaint1d
    use MpiExchangeColModule, only: n_recv, cmt_recv
    use ArrayHandlersModule, only: ifind
    use ListsModule, only: halomodellist
    use BaseModelModule, only: BaseModelType, GetBaseModelFromListByName
    use NumericalModelModule, only: NumericalModelType
    use BndModule, only: BndType, GetBndFromList
    use GwfModule, only: GwfModelType
    use MemoryManagerModule, only: mem_setval
    use MemoryHelperModule, only: split_mem_path_ff, create_mem_path
    ! -- local
    character(len=LENVARNAME)   :: name         !name of the array
    character(len=LENMEMPATH)   :: mempath      !memory path
    character(len=LENMEMPATH)   :: mempath_halo !memory path halo
    character(len=LENMODELNAME) :: mname        !name of model
    character(len=LENMODELNAME) :: hmname       !name of halo model
    character(len=LENMODELNAME) :: id
    character(len=LENFTYPE)     :: ft
    character(len=LENPACKAGENAME) :: pakname
    integer(I4B) :: i, j, ih, ip, istat, lenid, lenft, ipo, ibcnum
    logical, dimension(mxpack) :: lpack, lcreate
    class(BaseModelType), pointer :: mb
    class(GwfModelType), pointer :: mp
    class(NumericalModelType), pointer :: nmp
    class(BndType), pointer :: packobj
! ------------------------------------------------------------------------------
    !
    if (serialrun) then
      return
    end if
    !
    do i = 1, n_recv
      name   = cmt_recv(i)%name
      mempath = cmt_recv(i)%path
      call split_mem_path_ff(mempath, mname, id)
      lenid = len_trim(id)
      ! -- check if the data is related to packages
      lpack = .false.
      if (lenid > 0) then
        do ip = 1, mxpack
          ft = packftype(1,ip)
          lenft = len_trim(ft)
          if (lenid >= lenft) then
            if (id(1:lenft) == ft(1:lenft)) then
              lpack(ip) = .true.
            end if
          end if
        end do
        if (any(lpack)) then
          pakname = ''
          j = index(mempath,'-',back=.true.)
          read(mempath(j+1:),*,iostat=istat) ibcnum
          if (istat /= 0) then
            ibcnum = 0
            pakname = trim(id)
          end if
        end if
        if (any(lpack)) then
          do ih = 1, MpiWorld%hnmodel
            if (trim(MpiWorld%hmodelm2names(ih)) == trim(mname)) then
              hmname = MpiWorld%hmodelnames(ih)
              mb => GetBaseModelFromListByName(halomodellist, hmname)
              ! -- cast model
              mp => null()
              select type (mb)
              class is (GwfModelType)
                mp => mb
              end select
              nmp => null()
              select type (mb)
              class is (NumericalModelType)
                nmp => mb
              end select
              if((.not.associated(mp)).or.(.not.associated(nmp))) then
                call ustop('Program error 1 mpi_set_gwfhalo_world_mvr')
              end if
              ! -- check if package should be created
              lcreate = .false.
              do ip = 1, mxpack
                if (lpack(ip)) then
                  lcreate(ip) = .true.
                end if
              end do
              do ipo=1,nmp%bndlist%Count()
                packobj => GetBndFromList(nmp%bndlist, ipo)
                do ip = 1, mxpack
                  if (lpack(ip) .and. packobj%filtyp == packftype(1,ip) .and.&
                    packobj%ibcnum == ibcnum) then
                    lcreate(ip) = .false.
                  end if
                end do
              end do
              ! -- create package
              do ip = 1, mxpack
                if (lcreate(ip)) then
                  call mp%package_create(trim(packftype(2,ip)), 0, ibcnum,   &
                                         pakname, 0, 0)
                  packobj => GetBndFromList(nmp%bndlist, nmp%bndlist%Count())
                  call mem_setval(.true., 'P_ISHALO', packobj%memoryPath)
                end if
              end do
            end if
          end do
        end if
      end if
    end do
    !
    ! -- Set other variables
    do i = 1, n_recv
      name   = cmt_recv(i)%name
      mempath = cmt_recv(i)%path
      call split_mem_path_ff(mempath, mname, id)
      lenid = len_trim(id)
      lpack = .false.
      if (lenid > 0) then
        do ip = 1, mxpack
          ft = packftype(1,ip)
          lenft = len_trim(ft)
          if (lenid >= lenft) then
            if (id(1:lenft) == ft(1:lenft)) then
              lpack(ip) = .true.
            end if
          end if
        end do
      end if
      !
      if (any(lpack)) then
        do ih = 1, MpiWorld%hnmodel
          if (trim(MpiWorld%hmodelm2names(ih)) == trim(mname)) then
            hmname = MpiWorld%hmodelnames(ih)
            mempath_halo = create_mem_path(trim(hmname), trim(id))
            select case(cmt_recv(i)%memitype)
              case(iintsclr)
                call mem_setval(cmt_recv(i)%intsclr, name, mempath_halo)
              case(iaint1d)
                call mem_setval(cmt_recv(i)%aint1d, name, mempath_halo)
            end select
          end if
        end do
      end if
    end do
    !
    ! -- return
    return
  end subroutine mpi_set_gwfhalo_world_mvr

end module MpiExchangeGwfModule


text:
      subroutine DSBASI(K, N, T, X1, X2, J1, J2, BASI)
c Copyright (c) 1996 California Institute of Technology, Pasadena, CA.
c ALL RIGHTS RESERVED.
c Based on Government Sponsored Research NAS7-03001.
c>> 1994-10-20 DSBASI Krogh  Changes to use M77CON
c>> 1993-01-12 DSBASI CLL  Using DSBASD in place of _SBAS.
c>> 1992-11-12 C. L. Lawson, JPL Made J1 & J2 inout.  Were out.
c>> 1992-11-02 C. L. Lawson, JPL
c>> 1988-04-01 C. L. Lawson, JPL
C
C   This subroutine computes the integral over [X1,X2] of each of the N
c   members of a Kth order family of B-spline basis functions defined
c   over the knot set T().  The values of these integrals will be
c   returned in BASI(j), j = 1, ..., N.  The indices J1 and J2 will be
c   set to indicate the subset of these that might be nonzero.  Thus
c   BASI(j), j = J1, ..., J2 might be nonzero, whereas the other
c   elements of BASI() will be zero.
c   Orders K as high as 20 are permitted.
c   Permits X1 .le. X2 or X1 gt. X2.  In the former case all returned
c   values will be .ge. 0, while in the latter case all will be .le. 0.
C      The "proper interpolation interval" is from T(K) to T(N+1).
c      For most reliable results, X1 and X2 should lie in this
c      interval, however this subr will give a result even if this
c      is not the case by use of extrapolation.
C      This subr applies a 2, 6, or 10 point Gauss formula on
c      subintervals of [X1,X2] which are formed by the included
c      (distinct) knots.
c     ------------------------------------------------------------------
C                           ARGUMENTS
c
C  K      [in] ORDER OF B-SPLINE, 2 .le. K .le. 20
c               Note that the polynomial degree of the segments of
c               the spline are one less than the order.
C  N      [in] LENGTH OF COEFFICIENT ARRAY
C  T()    [in] KNOT ARRAY OF LENGTH N+K
C  X1,X2  [in] END POINTS OF QUADRATURE INTERVAL.
C               Permit X1 .le. X2 or X1 .gt. X2.
c
c  J1, J2  [inout integers] On entry must contain integer values.
c     Will be passed to DSFIND for lookup.  Any integer values
c     are ok but previous returned values may aid efficiency.
c     On return will satisfy 1 .le. J1 .le. J2 .le. N.
c     Indicates that only the output values BASI(j), j = J1,...,J2
c     might be nonzero.
C  BASI()  [out, floating pt]  Array in which the N integrals will be
c           stored.
C
C     ERROR CONDITIONS: None
c     ------------------------------------------------------------------
C     Based on spline quadrature subroutines called BQUAD or BSQUAD
c     written by D. E. Amos,  Sandia, June, 1979.  Documented in
c     SAND79-1825.
c     Uses the (T, BCOEF, N, K) representation of a B-spline as
c     presented in A PRACTICAL GUIDE TO SPLINES by Carl De Boor,
c     Springer-Verlag, 1978.
c     Current version by C. L. Lawson, JPL, March 1988.
c     ------------------------------------------------------------------
C--D replaces "?": ?SBASI, ?SFIND, ?SBASD
c     ------------------------------------------------------------------
c     GPTS(), GWTS()  Index 1 selects 2-point Gaussian formula.
c                     Indices 2-4 select 6-point Gaussian formula.
c                     Indices 5-9 select 10-point Gaussian formula.
c     ------------------------------------------------------------------
      integer KMAX
      parameter(KMAX = 20)
      integer IK, IL1, IL2, J, J1, J2, JF, K, LEFT, M, MF, MODE, N
      double precision A, AA, B, BASI(N), BB
      double precision BMA, BPA, BVAL1(KMAX), BVAL2(KMAX), C1
      double precision FAC, GPTS(9), GWTS(9), HALF
      double precision T(N+K), TA, TB, X1, X2, ZERO
      parameter(HALF = 0.50D0, ZERO = 0.0D0)
C
*     DATA GPTS            / 5.77350269189626D-01, 2.38619186083197D-01,
*    1 6.61209386466265D-01, 9.32469514203152D-01, 1.48874338981631D-01,
*    2 4.33395394129247D-01, 6.79409568299024D-01, 8.65063366688985D-01,
*    3 9.73906528517172D-01/
C
*     DATA GWTS            / 1.00000000000000D+00, 4.67913934572691D-01,
*    1 3.60761573048139D-01, 1.71324492379170D-01, 2.95524224714753D-01,
*    2 2.69266719309996D-01, 2.19086362515982D-01, 1.49451349150581D-01,
*    3 6.66713443086880D-02/
c
c     Gaussian quadrature abcissae and weights to 40 decimal digits.
      data GPTS / 0.5773502691896257645091487805019574556476d0,
     2            0.2386191860831969086305017216807119354186d0,
     3            0.6612093864662645136613995950199053470064d0,
     4            0.9324695142031520278123015544939946091348d0,
     5            0.1488743389816312108848260011297199846176d0,
     6            0.4333953941292471907992659431657841622001d0,
     7            0.6794095682990244062343273651148735757693d0,
     8            0.8650633666889845107320966884234930485275d0,
     9            0.9739065285171717200779640120844520534283d0/
c
      data GWTS / 1.0d0,
     2            0.4679139345726910473898703439895509948120d0,
     3            0.3607615730481386075698335138377161116612d0,
     4            0.1713244923791703450402961421727328935273d0,
     5            0.2955242247147528701738929946513383294210d0,
     6            0.2692667193099963550912269215694693528586d0,
     7            0.2190863625159820439955349342281631924634d0,
     8            0.1494513491505805931457763396576973323908d0,
     9            0.06667134430868813759356880989333179285686d0/
c     ------------------------------------------------------------------
      if(K .gt. KMAX) then
         call IERM1('DSBASI',1,2,'Require KORDER .le. KMAX.',
     *   'KORDER',K,',')
         call IERV1('KMAX',KMAX,'.')
         return
      endif
C
      AA= min(X1,X2)
      BB= max(X1,X2)
C
      IL1 = J1+K-1
      IL2 = J2
      call DSFIND(T,K, N+1,AA,IL1,MODE)
      call DSFIND(T,K, N+1,BB,IL2,MODE)
c
      if(BB .eq. T(IL2) .and. IL2 .gt. K .and. IL2 .gt. IL1)
     *    IL2 = IL2 -1
      J1 = IL1+1-K
      J2 = IL2
      do 15 J=1,N
         BASI(J) = ZERO
   15 continue
      if(AA .eq. BB) return
c
C                  Selection of Gaussian quadrature formula:
c                  KORDER .le. 4          =>   2 point formula
c                  5 .le. KORDER .le. 12  =>  6 point formula
c                  KORDER .ge. 13         =>  10 point formula
c
      if(K .le. 4) then
         JF=0
         MF=1
      elseif(K .le. 12) then
         JF=1
         MF=3
      else
         JF=4
         MF=5
      endif
c
      do 40 LEFT=IL1,IL2
         TA=T(LEFT)
         TB=T(LEFT+1)
         if(TA .ne. TB) then
            if(LEFT .eq. K) then
               A = AA
            else
               A= max(AA,TA)
            endif
            if(LEFT .eq. N) then
               B = BB
            else
               B= min(BB,TB)
            endif
            BMA= HALF *(B-A)
            BPA= HALF *(B+A)
            DO 30 M=1,MF
               C1=BMA*GPTS(JF+M)
               FAC = GWTS(JF+M)*BMA
               call DSBASD(K, LEFT, T, BPA+C1, 0, BVAL1)
               call DSBASD(K, LEFT, T, BPA-C1, 0, BVAL2)
               J = LEFT-K
               do 20 IK = 1,K
                  J = J+1
                  BASI(J) = BASI(J)+FAC*(BVAL1(IK)+BVAL2(IK))
   20          continue
   30       continue
         endif
   40 continue
      if(X1 .gt. X2) then
         do 50 J=J1, J2
            BASI(J) = -BASI(J)
   50    continue
      endif
      return
      end


text:
c     program DRSRANUS
c>> 1995-10-23 DRSRANUS Krogh  Reduced N to work on PC in D.P.
c>> 1994-10-19 DRSRANUS Krogh  Changes to use M77CON
c>> 1992-03-13 DRSRANUS CLL
c
c     Driver to demonstrate use of SRANUS to generate random numbers
c     from the uniform distribution between A and A+B.
c     Program computes histogram for N numbers
c     ------------------------------------------------------------------
c--S replaces "?": DR?RANUS, ?RANUS, ?STAT1, ?STAT2
c     ------------------------------------------------------------------
      integer N, NCELLS
      parameter(N = 8000, NCELLS = 10+2)
      real             A, B, STATS(5), YTAB(N)
      integer IHIST(NCELLS)
      data A, B / 20.0e0, 10.0e0/
c     ------------------------------------------------------------------
      STATS(1) = 0.0e0
c                          Get array of scaled random numbers.
c
      call SRANUS(YTAB, N, A, B)
c
c                          Compute statistics and histogram.
c
        call SSTAT1(YTAB, N, STATS, IHIST, NCELLS, A, A+B)
c
c                          Print the statistics and histogram.
c
      write(*,'(13x,a//)') 'Uniform random numbers from SRANUS'
      call SSTAT2(STATS, IHIST, NCELLS, A, A+B)
      end