      module communicate
      implicit none
      save 
 
      integer, public :: rcv_p   ! The amount of data received from the parent
      integer, public :: snd_p   ! The amount of data sent to the parent
      integer, public :: snd_c   ! The amount of data sent to ea. child
      integer, public :: rcv_c   ! The amount of data recieved from ea. child
      real(8), public, dimension(10000) :: r8_frm_p !r8_snd_c
      real(8), public, dimension(10000) :: r8_frm_c !r8_rcv_c

! Public Module Procedures
      public :: initp_info,           & 
                get_me,               &
                get_cmytid,           & 
                get_mytid,            &
                get_nproc,            &
                myinit,               & 
                x_cque_c,             & 
                x_rdat_c,             &
                x_rdat_p,             &
                x_rflg_p,             &
                x_rtids_p,            & 
                x_riflg_c,            &
                x_riflg_p,            &
                x_sdat_c,             &
                x_sdat_p,             &
                x_stids,              &
                x_sflg_c,             &
                x_siflg_p,            &
                x_skflg_c,            &
                openf,                &
                x_sp_c,               &
                pk,                   &
                upk,                  &
                destin,               &
                recv                             

     interface pk 
        module procedure pk_i1
        module procedure pk_r1
        module procedure pk_c1
        module procedure pk_in
        module procedure pk_rn
        module procedure pk_cn
     end interface pk

     interface upk
        module procedure upk_i1
        module procedure upk_r1
        module procedure upk_c1
        module procedure upk_in
        module procedure upk_rn
        module procedure upk_cn
     end interface upk

     data pkset /0/

! Private Module Procedure
      private ::  dir, newelm, findelm, add2p_info, &
                  pk_i1, pk_r1, pk_c1, pk_in, pk_rn ,pk_cn,  &
                  upk_i1, upk_r1, upk_c1, upk_in, upk_rn, upk_cn


! Private Global Data
      private
         integer :: pkset  !pk & upks
         integer :: info
         integer :: me
         integer :: mytid
         integer, dimension(0:256):: mptid
         integer :: who
         character(len=10) :: cmytid 

      type, private :: procinfo
           character(len=30) :: exec 
           integer :: status 
           integer :: nproc 
           integer, dimension(0:256) :: tids 
      end type procinfo

      integer,parameter :: elms = 50
      type(procinfo), dimension(elms) :: p_info

      contains

! -------------------------------------------------------------------- !
      subroutine initp_info
! -------------------------------------------------------------------- !
      implicit none
      integer :: i
      p_info(:)%exec="xxx"
      p_info(:)%status = 0
      p_info(:)%nproc  = 0
      do i=0,256
         p_info(:)%tids(i) = 0
      enddo
      end subroutine initp_info

! -------------------------------------------------------------------- !
      function get_nproc(exec)
! -------------------------------------------------------------------- !
      implicit none
      character(len=30), intent(in) :: exec
      integer :: get_nproc
      integer :: i
      i=findelm(exec)
      get_nproc=p_info(i)%nproc
      end function get_nproc
! -------------------------------------------------------------------- !
      function get_me() 
! -------------------------------------------------------------------- !
      implicit none

      integer :: get_me
      get_me = me 

      end function get_me

! -------------------------------------------------------------------- ! 
      function get_mytid()
! -------------------------------------------------------------------- !
      implicit none
      integer get_mytid
      get_mytid = mytid
      end function get_mytid
! -------------------------------------------------------------------- !
      function get_cmytid()
! -------------------------------------------------------------------- !
      implicit none
      character(len=10) :: get_cmytid
      get_cmytid = cmytid
      end function get_cmytid

! -------------------------------------------------------------------- !
      subroutine myinit(exec,parent,envstring) 
! -------------------------------------------------------------------- !
      implicit none
      include 'fpvm3.h'

      character(len=30), intent(in) :: exec
      character(len=70), intent(in) :: envstring

      integer :: i, nnproc
      integer, parameter  :: slen=30
      character(len=30)   :: aexec
      character(len=30)   :: fname
      character(len=30)   :: parent
      character(len=slen) :: mycomp
      character(len= 4)   :: strng = '_std'

      external fbuff
      external getcompname

      call pvmfmytid( mytid )
      call pvmfparent( mptid )

      nnproc=1
      call add2p_info(parent,nnproc,mptid)

      write(cmytid,'(i10)') mytid
      cmytid=trim(adjustl(cmytid))
      fname=trim(adjustl(exec))//strng//cmytid
      call openf( 1, fname, envstring)
      write(1,*) 'hey'
      call getcompname(mycomp, slen)
      mycomp = trim(adjustl(mycomp))
      write(1,*) 'This computer =',mycomp
      write(1,*) 'mptid=',mptid
      call fbuff()
      end subroutine myinit

! -------------------------------------------------------------------- !
      function newelm()
! -------------------------------------------------------------------- !
! Purpose: Determines the index of the first unused element of p_info
! -------------------------------------------------------------------- !
      implicit none
      integer :: newelm
      integer :: i
      integer :: nsize
      external fbuff

      nsize = size(p_info)
      do i=1,nsize
        if(p_info(i)%status.eq.0) then 
           newelm=i        
           return
        endif 
      enddo    
      write(1,*) 'ERROR newelm: no more new elements'
      call fbuff()

      end function newelm

! -------------------------------------------------------------------- !
      function findelm(exec)
! -------------------------------------------------------------------- !
! Purpose: Determines the index of the array p_info corresponding to
!          exec) 
! -------------------------------------------------------------------- !
      implicit none
      integer :: findelm
      integer :: i
      integer :: nsize
      character(len=30) :: exec
      external fbuff

      nsize = size(p_info)
      do i=1,nsize
         if(p_info(i)%exec .eq. exec) then 
            findelm = i
            return
          endif
      enddo
      write(1,*) 'ERROR findelm: p_info index corr. to exec not found'
      call fbuff()

      end function findelm

! -------------------------------------------------------------------- !
      subroutine add2p_info(exec,nproc,tids)
! -------------------------------------------------------------------- !
      implicit none
      character(len=30), intent(in) :: exec 
      integer, intent(in) :: nproc
      integer, intent(in), dimension(0:256) :: tids 
      integer :: i 
      integer,save :: ifst
      data ifst /1/

      if(ifst==1) then 
        call initp_info
        ifst=0
      endif

      i=newelm()
      p_info(i)%status = 1
      p_info(i)%exec   = exec 
      p_info(i)%nproc  = nproc 
      p_info(i)%tids   = tids
    
      end subroutine add2p_info

! ----------------------------------------------------------------- !
      subroutine x_cque_c(msgtype,exec,ipsz)
! ----------------------------------------------------------------- !
! Purpose :  This subroutine acts as a shell containing message 
!            passing code to communicated with the subroutines 
!            found in the spawned slaves
! Glossary:
! ipsz       -- Total number of evaluations
! snd_c      -- Number of real numbers sent to ea. child process 
! rcv_c      -- Number of real numbers received from ea. child process 
! paraminlin -- Array storing data to be sent to ea. child process
! fitarrya   -- Array storing data returned from ea. child process
! ----------------------------------------------------------------- !

      implicit none
      include 'fpvm3.h'
      
      integer,intent(in) :: msgtype
      character(len=30),intent(in) :: exec
      integer,intent(in) :: ipsz

      integer :: i,ip,j,k,kkk,icy,jcy
      integer :: istat
      integer :: job
      integer :: jobindex
      integer, dimension(100)   :: jobline
!     integer, dimension(10000) :: jobs
     
      ip = findelm(exec)
      istat = 1
      jobline(:) = -1
      jobindex = 0

! ---- Send a job to each slave process

      do i=0,(p_info(ip)%nproc-1)
         jobindex=jobindex+1
         j=snd_c*jobindex-(snd_c-1)
         call pvmfinitsend( PVMDEFAULT, info )
         call pvmfpack( INTEGER4, istat, 1, 1, info )
         call pvmfpack( INTEGER4, snd_c, 1, 1, info )
         call pvmfpack( REAL8, r8_frm_p(j), snd_c, 1, info )
         call pvmfsend( p_info(ip)%tids(i), msgtype, info )
         jobline(jobindex)=i
      enddo

      if(jobindex.eq.ipsz) goto 20   !#jobs = # slaves

! ---- Retrieve results and submit new job until queue is empty 
      do while(jobindex .lt. ipsz)
         call pvmfrecv( -1, msgtype, info )
         call pvmfunpack( INTEGER4, who, 1, 1, info )
         do k=1,ipsz+1
            if(jobline(k).eq.who)then
               job=k
               jobline(job)=-2 ! -2-> complete 
               goto 10
            endif
            if(k.eq.ipsz+1) then
               print*,'ERROR TRAP1: x_cdat_c k=',k 
               stop
            endif 
         enddo
10       continue
         call pvmfunpack( INTEGER4, rcv_c, 1, 1, info )
         j=rcv_c*job-(rcv_c-1)
         call pvmfunpack( REAL8, r8_frm_c(j), rcv_c, 1, info )

!        print*,'rcv_c',rcv_c,' who=', who,' job=',job  !TESTT
!        write(6,900) rcv_c, who,job,r8_frm_c(j+1),r8_frm_c(j)
900      format('rcv_c',i6,' who=',i6,' job=',i6,' r8=',f14.3,2x,e14.7)
!        print*,'jcy=',jcy
!        jcy=snd_c*job-(snd_c-1)
!        do icy=jcy,jcy+snd_c-1
!           write(6,*) icy, r8_frm_p(icy)
!        enddo
!        do kkk=0,rcv_c-1
!            print*,'r8_frm_c(',j+kkk,')=',r8_frm_c(j+kkk)
!        enddo
!        print*,' '

         jobindex = jobindex +1
         j=snd_c*jobindex -(snd_c-1)
         call pvmfinitsend( PVMDEFAULT, info )
         call pvmfpack( INTEGER4, istat, 1, 1, info )
         call pvmfpack( INTEGER4, snd_c, 1, 1, info )
         call pvmfpack( REAL8, r8_frm_p(j), snd_c, 1, info )
         call pvmfsend( p_info(ip)%tids(who), msgtype, info ) 
         jobline(jobindex)=who
      enddo

20    continue

! ----- Retrieve all remaining results from slaves

      do i=0,(p_info(ip)%nproc-1)
         call pvmfrecv( -1, msgtype, info )
         call pvmfunpack( INTEGER4, who, 1, 1, info )
         do k=1,ipsz+1
            if(jobline(k).eq.who)then
               job=k
               jobline(job)=-2 ! -2-> complete
               goto 30
            endif
            if(k.eq.ipsz+1) then
               print*,'ERROR TRAP2: x_cdat_c who=',who
               stop
            endif
         enddo
30       continue
         call pvmfunpack( INTEGER4, rcv_c, 1, 1, info )
         j=rcv_c*job-(rcv_c-1)
         call pvmfunpack( REAL8, r8_frm_c(j), rcv_c, 1, info )

!        write(6,901) rcv_c, who,job,r8_frm_c(j+1),r8_frm_c(j)
901      format('rcv_c',i6,' who=',i6,' job=',i6,' r8=',f14.3,2x,e14.7)

!        print*,'rcv_c',rcv_c,' who=', who,' job=',job !TESTT
!        do kkk=0,rcv_c-1
!           print*,'r8_frm_c(',j+kkk,')=',r8_frm_c(j+kkk)
!        enddo
!        print*,' '

      enddo

      end subroutine x_cque_c

! -------------------------------------------------------------
      subroutine x_rdat_c(msgtype,exec)
! -------------------------------------------------------------

      implicit none
      include 'fpvm3.h'

      integer,intent(in) :: msgtype
      character(len=30), intent(in) :: exec
      integer :: i,ip,j,k
      external fbuff

      ip = findelm(exec)

      write(1,*) 'ip=',ip
      write(1,*) 'p_info(ip)%nproc=',p_info(ip)%nproc
      do i = 1,p_info(ip)%nproc
         call pvmfrecv( -1, msgtype, info )
         call pvmfunpack( integer4, who, 1, 1, info )
         call pvmfunpack( integer4, rcv_c, 1, 1, info )
         k= 1+rcv_c*(who)
         write(1,*) 'i=',i,'rcv from =',who,' rcv_c=',rcv_c,' info=',info
         call pvmfunpack( real8, r8_frm_c(k), rcv_c, 1, info )
      enddo

      end subroutine x_rdat_c

! --------------------------------------------------------------  !
     subroutine x_rdat_p(msgtype,istat)
! -------------------------------------------------------------- !
! Purpose: This subroutine receives data from ct that will be 
!          be distributed to ms for atomic kinetic cal.
! -------------------------------------------------------------- !
! NOTE: It is the parents responsibility to state the size of the
!       array section (rcv_p)
! ______________________________________________________________ !
     implicit none

     integer,intent(in)  :: msgtype
     integer,intent(out) :: istat

     integer :: i

     include 'fpvm3.h'
     external fbuff

     call pvmfrecv( mptid, msgtype, info )
     call pvmfunpack( integer4, istat, 1, 1, info )
!    write(1,*) 'istat=',istat
     if(istat.eq.0) return 
     call pvmfunpack( integer4, rcv_p, 1, 1, info )
     call pvmfunpack( real8, r8_frm_p, rcv_p, 1, info )

!     write(1,*) 'rcv_p=',rcv_p
!     do i=1,rcv_p
!        write(1,*) 'r8_frm_p(',i,')=',r8_frm_p(i)
!        call fbuff ()
!     enddo     

     end subroutine x_rdat_p

! --------------------------------------------------------------  !
     subroutine x_rflg_p(msgtype,flg)
! --------------------------------------------------------------  !
     implicit none
     include 'fpvm3.h'

     integer,intent(in)  :: msgtype
     integer,intent(out) :: flg

     call pvmfrecv( mptid, msgtype, info )
     call pvmfunpack( INTEGER4, flg, 1, 1, info )

     end subroutine x_rflg_p

! ------------------------------------------------------------ !
      subroutine x_rtids_p(msgtype,exec)
! ------------------------------------------------------------ !
! Purpose: Receive initial data from parent
! ------------------------------------------------------------ !

      implicit none
      include 'fpvm3.h'

      integer, intent(in) :: msgtype
      character(len=30), intent(in)  :: exec 

      integer :: i
      integer :: nproc
      integer, dimension(0:256) :: thslyrstids
      external fbuff

      call pvmfparent( mptid )
      call pvmfmytid( mytid )

      write(1,*) 'mptid=',mptid
      call fbuff()
      call pvmfrecv( mptid, msgtype, info )
      call pvmfunpack( integer4, nproc, 1, 1, info )
      call pvmfunpack( integer4, thslyrstids, nproc, 1, info )
      do i=0,nproc-1
         if(mytid.eq.thslyrstids(i)) me = i 
      enddo
      write(1,*) 'me=',me
      call fbuff()
      call add2p_info(exec,nproc,thslyrstids)

      end subroutine x_rtids_p

! -----------------------------------------------------------  !
     subroutine x_riflg_c(msgtype,exec)
! ------------------------------------------------------------ !
! Purpose:
! --  wait for notice of initialization from nodes

      implicit none
      include 'fpvm3.h'

      integer,intent(in) :: msgtype
      character(len=30), intent(in) :: exec

      integer :: i,ip
      integer :: igo

      ip = findelm(exec)

      igo = 1
      do i=1,p_info(ip)%nproc
         call pvmfinitsend( PVMDEFAULT, info )
         call pvmfpack( INTEGER4, igo, 1, 1, info )
         call pvmfsend( p_info(ip)%tids(i-1), msgtype, info )

         call pvmfrecv( -1, msgtype, info )
         call pvmfunpack( INTEGER4, who, 1, 1, info )
         write(1,*) 'who=',who
      enddo
      end subroutine x_riflg_c

! --------------------------------------------------------------  !
     subroutine x_riflg_p(msgtype)
! -------------------------------------------------------------- !
! Purpose: This subroutine receives initialization flag from parent 
! -------------------------------------------------------------- !
     implicit none
     include 'fpvm3.h'

     integer,intent(in) :: msgtype

     integer :: igo

     call pvmfrecv( mptid, msgtype, info )
     call pvmfunpack( integer4, igo, 1, 1, info )
     end subroutine x_riflg_p

! -------------------------------------------------------------
      subroutine x_sdat_c(msgtype,exec,istat)
! -------------------------------------------------------------
    
      implicit none
      include 'fpvm3.h'

      integer,intent(in) :: msgtype

      integer :: i,ip,j,k
      character(len=30), intent(in) :: exec
      integer :: istat
      external fbuff

      ip = findelm(exec)
       
      do i= 1,p_info(ip)%nproc
         j= 1+snd_c*(i-1)
!        write(1,*) 'i=',i,' j=',j
         call fbuff()
         call pvmfinitsend( PVMDEFAULT, info )
         call pvmfpack( INTEGER4, istat, 1, 1, info )
         if(istat.ne.0) then
           call pvmfpack( integer4, snd_c, 1, 1, info )
           call pvmfpack( real8, r8_frm_p(j), snd_c, 1, info )
!          write(1,*) p_info(ip)%tids(i-1),' r8(',j,')=',r8_frm_p(j),' r8(',j+1,')=',r8_frm_p(j+1)
         endif
         call pvmfsend( p_info(ip)%tids(i-1), msgtype, info )
      enddo

      end subroutine x_sdat_c

! ----------------------------------------------------- !
      subroutine x_sdat_p(msgtype)
! ----------------------------------------------------- !
! Purpose: This program sends data from the current
!          executable (x) to the parent (p).
! ----------------------------------------------------- !

      implicit none
      include 'fpvm3.h'

      integer,intent(in) :: msgtype
      call pvmfinitsend( PVMDEFAULT, info )
      call pvmfpack( integer4, me, 1, 1, info )
      call pvmfpack( integer4, snd_p, 1, 1, info )
      call pvmfpack( REAL8, r8_frm_c, snd_p, 1, info )
      call pvmfsend( mptid, msgtype, info )

      end subroutine x_sdat_p

! ----------------------------------------------------------------- !
      subroutine x_stids(msgtype,whosinfo,where,itm_info,itm_where)
! ----------------------------------------------------------------- !
! Purpose:  Send initialization data to child.
! ----------------------------------------------------------------- !
! NOTE: Default send all "whosinfo" to target process "where".
!       The info send is all the tids aquired by x_sp_c up to the time
!       of this sub. call. Any further info added by calls 
!       to x_sp_c will not be included. It is not advisable to call
!       x_sp_c several times to spawn several copies of the same 
!       executable
!
!     msgtype     - Message tag
!     whosinfo    - The name of the processes whos info will be sent
!     where       - The name of the processes that will recv the info 
!     itm_info    - Optional which singular proc info will be sent 
!                   (Index)
!     itm_where   - Optional which singular proc will recv the info
!                   (Index)
! ----------------------------------------------------------------- !

      implicit none
      include 'fpvm3.h'

      integer, intent(in) :: msgtype
      character(len=30), intent(in) :: whosinfo
      character(len=30), intent(in) :: where 
      integer, intent(in), optional :: itm_info
      integer, intent(in), optional :: itm_where

      integer :: i,j
      logical :: entry_chk,chk1,chk2
      external fbuff

      i=findelm(whosinfo)
      j=findelm(where) 

      chk1 = present(itm_info)
      chk2 = present(itm_where)
      if (chk1 .neqv. chk2) then
         write(1,*) 'ERROR x_stds: Both optional entries NOT present'
         call fbuff
         stop
      endif     
      entry_chk = present(itm_info)
      select case(entry_chk)
      case(.false.)
        call pvmfinitsend( PVMDEFAULT, info )
        call pvmfpack( INTEGER4, p_info(i)%nproc, 1, 1, info )
        call pvmfpack( INTEGER4, p_info(i)%tids, p_info(i)%nproc, 1, info )
!       print*,'p_info(',i,')%nproc=',p_info(i)%nproc 
!       print*,'p_info(',i,')%tids=',p_info(i)%tids
        call pvmfmcast( p_info(j)%nproc, p_info(j)%tids, msgtype, info )
!       print*,'aft mcast 1'
      case(.true.)
        call pvmfinitsend( PVMDEFAULT, info )
        call pvmfpack( INTEGER4, 1, 1, 1, info )
        call pvmfpack( INTEGER4, p_info(i)%tids(itm_info), 1, 1, info )
        call pvmfsend( p_info(j)%tids(itm_where), msgtype, info )
!       print*,'aft mcast2'
      end select 

      end subroutine x_stids

! -------------------------------------------------------------  !
     subroutine x_siflg_p(msgtype)
! -------------------------------------------------------------- !
! Purpose: This subroutine sends initialization flag back to
!          parent. 
! -------------------------------------------------------------- !
     implicit none
     include 'fpvm3.h'

     integer,intent(in) :: msgtype

     call pvmfinitsend( pvmdefault, info )
     call pvmfpack( integer4, me, 1, 1, info )
     call pvmfsend( mptid, msgtype, info )

     end subroutine x_siflg_p

! ------------------------------------------------------------- !
      subroutine x_sflg_c(msgtype,exec,flg)
! ------------------------------------------------------------- !
! Purpose: Sends flag to child process
! ------------------------------------------------------------- !

      implicit none
      include 'fpvm3.h'

      integer,intent(in) :: msgtype
      character(len=30), intent(in) :: exec
      integer, intent(in) :: flg
      integer :: i,ip
      ip=findelm(exec)
      do i=0,(p_info(ip)%nproc-1)
         call pvmfinitsend( PVMDEFAULT, info )
         call pvmfpack( INTEGER4, flg, 1, 1, info )
         call pvmfsend( p_info(ip)%tids(i), msgtype, info )
      enddo
      end subroutine x_sflg_c
! ------------------------------------------------------------- ! 
      subroutine x_skflg_c(msgtype,exec) 
! ------------------------------------------------------------- ! 
! Purpose: Sends kill flag to child process
! ------------------------------------------------------------- !

      implicit none
      include 'fpvm3.h'

      integer,intent(in) :: msgtype
      character(len=30), intent(in) :: exec

      integer :: i,ip 
      integer :: istat

      ip=findelm(exec)
      
!     Sending kill command to slaves (istat=0)
      istat=0
      do i=0,(p_info(ip)%nproc-1)
         call pvmfinitsend( PVMDEFAULT, info )
         call pvmfpack( INTEGER4, istat, 1, 1, info )
         call pvmfsend( p_info(ip)%tids(i), msgtype, info )
         print*,'Kill notice from proc =',i,' has been sent!'
      enddo

      end subroutine x_skflg_c

! -------------------------------------------------------- !
      subroutine x_sp_c(exec,nproc,strt,shft)
! -------------------------------------------------------- !
! Purpose:  This program spawns child process  
! -------------------------------------------------------- !

      implicit none
      include 'fpvm3.h'  

      integer, parameter :: s_len = 30
      character(len=s_len), intent(in) :: exec 
      integer,intent(in)               :: nproc
      integer :: strt   ! Starting computer (0=This computer)
      integer :: shft   ! Shift from starting computer (0=Zero shift)
      integer :: i,clen,len,lmin
      integer :: nhost
      integer :: narch 
      integer :: numtot, numt
      integer, dimension(256) :: dtid
      integer, dimension(256) :: speed
      integer, dimension(0:256) :: tids
      character(len=s_len) mycomp
      character(len=s_len) comp
      character(len=s_len) arch(256)
      character(len=s_len) host(256)
      character(len=s_len) prmhost(256)
      external getcompname
      call getcompname(mycomp,s_len)
      len=len_trim(mycomp)
      write(1,*) 'mycomp=',mycomp,' len=',len

      i=1
      numtot=0
      do
        call pvmfconfig(nhost,narch,dtid(i),host(i),arch(i),speed(i),info)
        prmhost(i)=host(i)
        write(1,*) 'nhost=',nhost,' host(',i,')=',host(i)
        call fbuff()
        i=i+1
        if(i.gt.nhost) exit 
      enddo
     
      do i=1,nhost
        comp=host(i)
        clen=len_trim(comp) 
        lmin=min(clen,len)
          if(mycomp(1:lmin).eq.comp(1:lmin)) then
            write(1,*) 'I am host(',i,')=',host(i)
            call fbuff()
            prmhost=cshift(host,i-1,1)   ! prmhost(1)=This host
            exit
          endif
      enddo
 
      prmhost=cshift(prmhost,strt,1)
      do i=1,nproc
        write(1,*) 'prmhost(1)= ',prmhost(1)
        call fbuff()
        call pvmfspawn( exec, PVMHOST, prmhost(1), 1, tids(i-1), numt )
        call cmpshft(prmhost,shft,nhost)
        numtot=numtot+numt
      enddo

! --- Check for problems
      if( numtot .lt. nproc ) then
         write(1,*) 'trouble spawning ',exec
         write(1,*) ' Check tids for error code'
         call fbuff()
         call shutdown(nproc,tids)
      endif
      
      call add2p_info(exec,numtot,tids)
      
      end subroutine x_sp_c

! ------------------------------------------------- !
      subroutine shutdown(nproc,tids)
! ------------------------------------------------- !
! Purpose: Kill all tasks I spawned and then myself
! ------------------------------------------------- !

      integer :: i,nproc
      integer, dimension(0:256) :: tids
      do i=0, nproc
         call pvmfkill( tids(i), info )
      enddo
      call pvmfexit( info )
      stop
      return

      end subroutine shutdown

! ------------------------------------------------- !
      subroutine cmpshft(array,ishft,ncomp)
! ------------------------------------------------- !
! Purpose: Shifts list of computers for spawning 
!          new processes
! ------------------------------------------------- !
      implicit none

      integer, parameter :: s_len = 30
      character(len=s_len), intent(inout), dimension(256) :: array
      integer, intent(in) :: ishft
      integer, intent(in) :: ncomp

      integer :: i,j
      character(len=s_len) :: aux

      if(ishft .gt.0) then
        do j=1,ishft
           aux=array(1)
           do i=1,(ncomp-1)
              array(i) = array(i+1)
           enddo
           array(ncomp) = aux
        enddo
      elseif(ishft.lt.0)then
        write(1,*) 'ERROR: x_sp_c ishft < 0'
        stop
      endif

      end subroutine cmpshft

! -------------------------------------------------------------- ! 
    subroutine openf(unitnum,ctemp,envstring)
! -------------------------------------------------------------- !
!  Purpose: This is a driver routine used to open files in a 
!           location defined by a given system enviroment 
!           variable (envstring) with a given unit number (unitnum) 
!           with a file name (ctemp)
! 
!  Limitations: 
!           ctemp   <= clen (30)    -File name
!           workdir <= slen (70)    -Working directory 
!           cffoo   <= dlen (100)   -Filename + Working dir.
!
!  ERRORS: Written in unit number 1 (FORCED)
!
! -------------------------------------------------------------- !
    implicit none
    integer, parameter :: slen = 70
    integer, parameter :: clen = 30 
    integer, parameter :: dlen = 100
    integer, parameter :: elen = 100
    integer :: envlen
    integer :: wlen
    integer :: ios
    integer :: flg
    integer :: unitnum
    character(len=clen) :: ctemp
    character(len=slen) :: workdir
    character(len=slen) :: envstring
    character(len=dlen) :: cffoo
    character(len=elen) :: errorstr
    external fbuff

   
    envlen = len_trim(envstring)
    call fbuff()
    call envarc( envstring, envlen, workdir, slen, flg, errorstr, elen )
    if(flg.eq.-1) then
      write(99,*) 'ERROR - sub. openf: incompatible string length ',flg
      call fbuff() 
      stop
    endif
    if(flg.eq.-2) then
      write(99,*) 'ERROR - sub. openf: env. var does not exist ',flg
      call fbuff()
      stop
    endif
    if(flg.eq.-3) then
      write(99,*) 'ERROR - sub. openf: envarc mkdir error  ',errorstr,' ',flg
      call fbuff()
      stop
    endif
    if(flg.eq.-4) then
      write(99,*) 'ERROR - sub. openf:Problem mkdir ',errorstr,' ',flg
      call fbuff()
      stop
    endif
    if(flg.eq.-99) then
      write(1,*) 'ERROR - sub. openf:TRAP -99 '
      call fbuff()
      stop
    endif

    wlen = len_trim(workdir)
    if(workdir(wlen:wlen).ne."/") workdir(wlen+1:wlen+1)="/"
    cffoo = dir(workdir,ctemp)
    close(unitnum)
    open(unit=unitnum,file=cffoo,iostat=ios)
    if(ios /= 0)then
      write(1,*) 'ERROR - sub. openf opening file ',cffoo
      call fbuff()
      stop
    endif
    write(1,*) 'FILE OPENED:',cffoo
    call fbuff()

    end subroutine openf
  
! -------------------------------------------------------------------- !
      function dir(string1,string2)
! -------------------------------------------------------------------- !
      implicit none
      character(len=100) :: dir
      character(len=*),intent(in)  :: string1,string2
      dir=trim(adjustl(string1))//adjustl(string2)

      end function dir

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

        subroutine destin( msgtype, prcn, ip  )
          implicit none
          include 'fpvm3.h'
          integer, intent(in) ::  msgtype        ! message id
          character(len=30), intent(in) :: prcn  ! proc name
          integer, intent(in), optional :: ip    ! individual proc
          integer :: i
          logical :: pi_chk
          external fbuff 
          i=findelm(prcn)
          pi_chk = present(ip)
          select case(pi_chk)
          case(.false.)
             call pvmfmcast( p_info(i)%nproc, p_info(i)%tids, msgtype, info )
          case(.true.)
             call pvmfsend( p_info(i)%tids(ip), msgtype, info )
          end select
          pkset=0
        end subroutine destin

        subroutine recv( msgtype, prcn, ip )
          implicit none
          integer, intent(in) :: msgtype
          character(len=30), intent(in), optional :: prcn  ! proc name 
          integer, intent(in), optional :: ip              ! indiv. proc indx 
          integer :: i
          logical :: prcn_chk
          logical :: ip_chk
          logical :: and_s , or_s
          external fbuff
          prcn_chk = present(prcn)
          ip_chk   = present(ip)
          and_s = (prcn_chk .and. ip_chk) 
          or_s  = (prcn_chk .or. ip_chk)
          if(or_s .eqv. .false.) then                               !No opt
             call pvmfrecv( -1, msgtype, info )
          elseif((and_s .eqv. .false.).and.(or_s .eqv. .true.)) then !One opt
             i=findelm(prcn) 
             call pvmfrecv( p_info(i)%tids(0), msgtype, info)
          elseif((and_s .eqv. .true.).and.(or_s .eqv. .true.)) then  !Both opt
             i=findelm(prcn)
             call pvmfrecv( p_info(i)%tids(ip), msgtype, info)
          else
             write(1,*) 'ERROR communicate:recv logic falure'
             call fbuff()
          end if
        end subroutine recv

        subroutine pk_i1(i)
          implicit none
          include 'fpvm3.h'
          integer, intent(in)  :: i
          if(pkset.eq.0)then
            call pvmfinitsend( PVMDEFAULT, info ) 
            pkset=1
          endif
          call pvmfpack( INTEGER4, i, 1, 1, info )
        end subroutine pk_i1

        subroutine pk_r1(r)
          implicit none
          include 'fpvm3.h'
          real(8), intent(in) :: r
          if(pkset.eq.0)then
            call pvmfinitsend( PVMDEFAULT, info ) 
            pkset=1
          endif
          call pvmfpack( REAL8, r, 1, 1, info)
        end subroutine pk_r1

        subroutine pk_c1(c)
          implicit none
          include 'fpvm3.h'
          character(len=*), intent(in) :: c
          if(pkset.eq.0)then
            call pvmfinitsend( PVMDEFAULT, info )
            pkset=1
          endif
          call pvmfpack( STRING, c, 1, 1, info)
        end subroutine pk_c1

        subroutine pk_in(i,n,strd)
          implicit none
          include 'fpvm3.h'
          integer, dimension(:), intent(in) :: i
          integer, intent(in) :: n
          integer, intent(in), optional :: strd
          logical :: strd_chk
          if(pkset.eq.0)then
            call pvmfinitsend( PVMDEFAULT, info ) 
            pkset=1
          endif
          strd_chk = present(strd)
          call pvmfpack( INTEGER4, n, 1, 1, info )
          select case(strd_chk)
          case(.true.)
            call pvmfpack( INTEGER4, i, n, strd, info )
          case(.false.)
            call pvmfpack( INTEGER4, i, n, 1, info )
          end select
        end subroutine pk_in 

        subroutine pk_rn(r,n,strd)
          implicit none
          include 'fpvm3.h'
          real(8), intent(in), dimension(:) :: r
          integer, intent(in) :: n
          integer, intent(in), optional :: strd
          logical :: strd_chk
          if(pkset.eq.0)then
            call pvmfinitsend( PVMDEFAULT, info ) 
            pkset=1
          endif
          strd_chk = present(strd)
          call pvmfpack( INTEGER4, n, 1, 1, info )
          select case(strd_chk)
          case(.true.)
          call pvmfpack( REAL8, r, n, strd, info )
          case(.false.)
          call pvmfpack( REAL8, r, n, 1, info )
          end select
        end subroutine pk_rn

        subroutine pk_cn(c,n,strd)
          implicit none
          include 'fpvm3.h'
          character(len=*), intent(in) :: c
          integer, intent(in) :: n
          integer, intent(in), optional :: strd
          logical :: strd_chk
          if(pkset.eq.0)then
            call pvmfinitsend( PVMDEFAULT, info )
            pkset=1
          endif
          strd_chk = present(strd)
          call pvmfpack( INTEGER4, n, 1, 1, info )
          select case(strd_chk)
          case(.true.)
            call pvmfpack( STRING, c, n, strd, info )
          case(.false.)
            call pvmfpack( STRING, c, n, 1, info )
          end select
        end subroutine pk_cn

        subroutine upk_i1(i)
          implicit none
          include 'fpvm3.h'
          integer, intent(out) :: i
          call pvmfunpack( INTEGER4, i, 1, 1, info )          
        end subroutine upk_i1

        subroutine upk_r1(r)
          implicit none
          include 'fpvm3.h'
          real(8), intent(out) :: r 
          call pvmfunpack( REAL8, r, 1, 1, info )
        end subroutine upk_r1

        subroutine upk_c1(c)
          implicit none
          include 'fpvm3.h'
          character(len=*), intent(out) :: c
          call pvmfunpack( STRING, c, 1, 1, info )
        end subroutine upk_c1
 
        subroutine upk_in(i,n)
          implicit none
          include 'fpvm3.h'
          integer, dimension(:), intent(out) :: i
          integer, intent(out) :: n
          call pvmfunpack( INTEGER4, n, 1, 1, info )
          call pvmfunpack( INTEGER4, i, n, 1, info )
        end subroutine upk_in

        subroutine upk_rn(r,n)
          implicit none
          include 'fpvm3.h'
          real(8), dimension(:), intent(out) :: r
          integer, intent(out) :: n
          call pvmfunpack( INTEGER4, n, 1, 1, info )
          call pvmfunpack( REAL8, r, n, 1, info )
        end subroutine upk_rn 

        subroutine upk_cn(c,n)
          implicit none
          include 'fpvm3.h'
          character(len=*), intent(out) :: c
          integer, intent(out) :: n
          call pvmfunpack( INTEGER4, n, 1, 1, info)
          call pvmfunpack( STRING, c, n, 1, info)
       end subroutine upk_cn
 ! ------------------------------------------------------------ ! 

    end module communicate
