c-----------------------------------------------------------------------
c Partial MPI library based on the Program-to-Program Communications
c ToolBox in the Macintosh OS.  No local buffering of messages is
c implemented, so that all messages must be received in the order sent,
c and receives with wildcard sources are not supported.
c the following subroutines are implemented:
c MPI_INIT, MPI_FINALIZE, MPI_SEND, MPI_RECV, MPI_ISEND, MPI_IRECV
c MPI_TEST, MPI_WAIT, MPI_SENDRECV, MPI_SSEND, MPI_ISSEND, MPI_WAITALL
c MPI_WAITANY, MPI_GET_COUNT, MPI_INITIALIZED, MPI_COMM_SIZE
c MPI_COMM_RANK, MPI_BCAST, MPI_BARRIER, MPI_REDUCE, MPI_SCAN
c MPI_ALLREDUCE, MPI_GATHER, MPI_ALLGATHER, MPI_SCATTER, MPI_ALLTOALL
c MPI_GATHERV, MPI_ALLGATHERV, MPI_SCATTERV, MPI_ALLTOALLV
c MPI_REDUCE_SCATTER, MPI_ABORT, MPI_WTIME, MPI_WTICK, MPI_TYPE_EXTENT
c The PPC Toolbox is described in Inside Macintosh: Interapplication
c Communication [Addison-Wesley, Reading, MA, 1993], chapter 11.
c The Message Passing Interface (MPI) is described in the reference,
c M. Snir, S. Otto, S. Huss-Lederman, D. Walker, and J. Dongarra,
c MPI: The Complete Reference [MIT Press, Cambridge, MA,1996]. 
c Fortran unit 2 is used throughout for error messages, and
c Fortran units 3 and 4 are used in MPI_INIT.
c written by viktor k. decyk, ucla
c copyright 1998-1999, regents of the university of california.
c all rights reserved.
c no warranty for proper operation of this software is given or implied.
c software or information may be copied, distributed, and used at own
c risk; it may not be distributed without this notice included verbatim
c with each file.
c update: november 9, 2000
      block data
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
c srnum = array of session reference numbers for each participating node
      common /mpiparms/ nproc, idproc, srnum, prnum
c declare common block for non-blocking messages
      integer MAXM, curreq, monitor
      parameter(MAXM=16,MAXMS=5*MAXM)
      dimension curreq(5,MAXM)
      integer*2 readrec(35,MAXM)
c curreq = request record for transmission parameters
c readrec = PPCReadPB Record
c monitor = (0,1,2) = (suppress,display,display & log) monitor messages
      common /mpisendrec/ curreq, readrec, monitor
c common block for message window
c cpptr = pointer to window structure
c nsp = amount of space between boxes
c nbx = size of box
c nds = number of message sizes monitored
      integer cpptr
      integer*2 crect(4), nsp, nbx, nds
      common /winmess/ cpptr, crect, nsp, nbx, nds
      save /mpiparms/, /mpisendrec/, /winmess/
      data nproc /0/
      data srnum /MAXS*0,0/
      data curreq /MAXMS*0/
      data monitor /1/
      data cpptr, nsp, nbx, nds /0,8,16,24/
      end
c-----------------------------------------------------------------------
      subroutine MPI_INIT(ierror)
c initialize the MPI execution environment
c ierror = error indicator
c input: none, output: ierror
      implicit none
      integer ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1)
      integer*2 prnum
c nproc = number of real or virtual processors obtained
c idproc = processor id
c srnum = array of session reference numbers for each participating node
c prnum = port reference number
      common /mpiparms/ nproc, idproc, srnum, prnum
c declare common block for non-blocking messages
      integer MAXM, curreq, monitor
      parameter(MAXM=16)
      dimension curreq(5,MAXM)
      integer*2 readrec(35,MAXM)
c monitor = (0,1,2) = (suppress,display,display & log) monitor messages
      common /mpisendrec/ curreq, readrec, monitor
c function declarations
      integer*2 Gestalt, PPCInit, PPCOpenSync
      integer*2 PPCInformAsync, PPCStartSync
      integer*2 PPCAcceptSync, PPCRejectSync
      integer*2 ioresult
      logical checkesc
      external Gestalt, PPCInit, PPCOpenSync
      external PPCInformAsync, PPCStartSync
      external PPCAcceptSync, PPCRejectSync
      external ioresult, checkesc
c MPI constants
      integer MPI_STATUS_SIZE, MPI_INTEGER
      parameter(MPI_STATUS_SIZE=5,MPI_INTEGER=18)
c local data
      integer*2 intw(2), ierr
      integer*2 noerr, notdone
      parameter(noerr=0,notdone=1)
      integer*2 openrec(28), ppcport(36), locname(52)
      integer*2 remoteport(36), remotename(52)
      integer*2 informrec(32), startrec(33), rejectrec(29)
      integer longw
      integer response, nerr, nofile, nv, stat(MPI_STATUS_SIZE), i
      character*4 creator, ftype
      character*32 portname
      character*64 location
c longw and intw are used to convert types between integer*2 and integer
      equivalence(longw,intw)
c creator and type correspond to executables produced by Absoft
      data creator,ftype /'MRWE','APPL'/
c
c check status of PPC ToolBox
c
c get information about the operating environment
      ierr = Gestalt(VAL4('ppc '),response)
      if (ierr.ne.noerr) then
         write (2,*) 'PPC not supported on this computer'
         ierror = ierr
         return
      endif
c check if gestaltPPCSupportsRealTime is set
      nerr = response/4096
      if ((nerr-2*(nerr/2)).ne.1) then
         write (2,*) 'Initializing PPC'
c initialize PPC Toolbox
         ierr = PPCInit()
         if (ierr.ne.noerr) then
            write (2,*) 'PPCInit Failed, ierr = ', ierr
            ierror = ierr
            return
         else
c update information about the operating environment
            ierr = Gestalt(VAL4('ppc '),response)
         endif
      endif
      nerr = response/2
c check if gestaltPPCSupportsOutgoing is set
      if ((nerr-2*(nerr/2)).ne.1) then
         write (2,*) 'AppleTalk not enabled in Chooser'
         ierror = 3
         return
      endif
c check if gestaltPPCSupportsIncoming is set
      if ((response-2*nerr).ne.1) then
         write (2,*) 'Program Linking not enabled in Sharing Setup CP'
         ierror = 4
         return
      endif
c
c everyone opens a port
c
c open file containing portname (and possibly participating nodes)
c first line in nodelist file on all nodes contains common portname
c if the file is missing or empty, a default name of ppc_link is used
      open(unit=3,file='nodelist',form='formatted',status='old',iostat=n
     1ofile)
      if (nofile.eq.0) then
         read (3,'(a32)',iostat=nerr) portname
         if ((nerr.ne.0).or.(portname.eq.' ')) then
            portname = 'ppc_link'
         endif
      else
         portname = 'ppc_link'
      endif
c check if port name is already being used
      prnum = 0
      call clrport(portname,prnum)
c set PPCPortRec
c script code = smRoman
      ppcport(1) = 0
c set port name string to portname
      call ctopascl(ppcport(2),portname)
c set selector to Creator and Type
      ppcport(19) = 1
c set creator
      ppcport(20) = 256*ichar(creator(1:1)) + ichar(creator(2:2))
      ppcport(21) = 256*ichar(creator(3:3)) + ichar(creator(4:4))
c set type
      ppcport(22) = 256*ichar(ftype(1:1)) + ichar(ftype(2:2))
      ppcport(23) = 256*ichar(ftype(3:3)) + ichar(ftype(4:4))
c set LocationNameRec
c set selector to NBPTypeLocation
      locname(1) = 2
c set object string to portname
      call ctopascl(locname(2),portname)
c set PPCOpenPBRec structure
c set ioCompletion code to NIL
      openrec(7) = 0
      openrec(8) = 0
c service type = ppcServiceRealTime and resFlag = 0
      openrec(23) = 256
c set pointer to ppcPortRec
      longw = loc(ppcport)
      openrec(24) = intw(1)
      openrec(25) = intw(2)
c set pointer to locationNameRec
      longw = loc(locname)
      openrec(26) = intw(1)
      openrec(27) = intw(2)
c make port network visible
      openrec(28) = 256
c open a PPC port synchronously
      ierr = PPCOpenSync(openrec)
      prnum = openrec(20)
      if (ierr.ne.noerr) then
         write (2,*) 'PPCOpen Failed, ierr = ', ierr
         if (ierr.eq.(-910)) then
            write (2,*) 'Another port is already open with this name'
         endif
         ierror = ierr
         return
c write port reference number to file
      else
         call clrport(portname,prnum)
      endif
c debug
      if (monitor.eq.2) write (2,*) 'local portname=', portname
c
c determine if node is master (idproc=0) or slave (idproc>0).
c on the master node, the second and subsequent lines of nodelist file
c contain names of the nodes participating, in the format hostname@zone.
c if this list of nodes is missing, then the node is a slave.
c every node also makes a connection to itself.
c
      nproc = 0
      if (nofile.eq.0) read (3,'(a64)',iostat=nerr) location
c must be slave
      if ((nofile.ne.0).or.(nerr.ne.0).or.(location.eq.' ')) then
         idproc = 1
c must be master
      else
         idproc = 0
      endif
c 
c * * * begin main iteration loop * * *
c
c prepare to accept connection
c
c set PPCInformPBRec
c set ioCompletion code to NIL
   10 informrec(7) = 0
      informrec(8) = 0
c set portRefNum
      informrec(20) = prnum
c set autoAccept to false
      informrec(23) = 0
c for connections to oneself, set autoAccept to true
      if (idproc.eq.nproc) informrec(23) = 1
c set pointer to ppcPortRec
      longw = loc(remoteport)
      informrec(24) = intw(1)
      informrec(25) = intw(2)
c set locationNameRec pointer to NIL
      informrec(26) = 0
      informrec(27) = 0
c set userName pointer to NIL
      informrec(28) = 0
      informrec(29) = 0
c receive session requests asynchronously
      ierr = PPCInformAsync(informrec)
      if (ierr.ne.noerr) then
         write (2,*) 'PPCInform failed, ierr=', ierr
         ierror = ierr
         call MPI_FINALIZE(nerr)
         return
      endif
c for connections to oneself, set selector to NoLocation
c and jump to PPCStart
      if (idproc.eq.nproc) then
         remotename(1) = 0
         go to 70
      endif
c wait for connection
   20 if (ioresult(informrec).eq.notdone) then
         if (checkesc(60)) then
            call MPI_FINALIZE(nerr)
            stop
         else
            go to 20
         endif
      else
         ierr = ioresult(informrec)
         if (ierr.ne.noerr) then
            write (2,*) 'PPCInform failed, ierr=', ierr
            ierror = ierr
            call MPI_FINALIZE(nerr)
            return
         endif
      endif
c extract processor id on first connection
      if (nproc.eq.0) idproc = informrec(31)
c check if remote portname agrees with local portname
      nv = ppcport(2)/512 + 1
      nerr = 0
      do 30 i = 1, nv
      nerr = nerr + abs(ppcport(i+1) - remoteport(i+1))
   30 continue
c set PPCRejectPBRec (also used for PPCAcceptPBRec)
c set ioCompletion code to NIL
      rejectrec(7) = 0
      rejectrec(8) = 0
c set session reference number
      rejectrec(21) = informrec(21)
      rejectrec(22) = informrec(22)
c accept if portnames agree
      if (nerr.eq.0) then
c indicate willingness to accept incoming session request
         ierr = PPCAcceptSync(rejectrec)
         if (ierr.ne.noerr) then
            write (2,*) 'PPCAccept error, idproc, ierr = ', idproc, ierr
            ierror = ierr
            call MPI_FINALIZE(nerr)
            return
         endif
c reject if portnames disagree
      else
c set rejectInfo code to 1
         rejectrec(28) = 0
         rejectrec(29) = 1
c reject a session request
         ierr = PPCRejectSync(rejectrec)
         if (ierr.ne.noerr) then
            write (2,*) 'PPCReject error, idproc, ierr = ', idproc, ierr
            ierror = ierr
         else
            write (2,*) 'Session rejected, idproc = ', idproc
            ierror = 5
         endif
c make string length visible
         remoteport(2) = remoteport(2) + 256*ichar('0')
         write (2,'(17a2)') 'remoteport = ', (remoteport(1+i),i=1,nv)
         call MPI_FINALIZE(nerr)
         return
      endif
      nproc = nproc + 1
c check for processor number overflow
      if (nproc.gt.MAXS) then
         write (2,*) 'processor number overflow, nproc = ', nproc
         ierror = 6
         call MPI_FINALIZE(nerr)
         return
      endif
c extract and store session reference number
      intw(1) = informrec(21)
      intw(2) = informrec(22)
      srnum(nproc) = longw
c debug
      if (monitor.eq.2) then
         write (2,*) 'connection accepted with idproc=', nproc-1
      endif
c accept more connections
      if (idproc.ge.nproc) go to 10
c
c master prepares to start connection
c
   40 nv = len(location)
c find zone delimiter '@'
      i = 0
   50 i = i + 1
      if (i.gt.nv) go to 60
      if (location(i:i).eq.'@') go to 60
      go to 50
   60 nv = i - 1
c set LocationNameRec
c set selector to NBPLocation
      remotename(1) = 1
c set object string to remote computer name
      call ctopascl(remotename(2),location(1:nv))
c set type string to 'PPCToolBox'
      call ctopascl(remotename(19),'PPCToolBox')
c set zone string to remote zone name
      nv = nv + 1
      i = len(location)
      if (nv.lt.i) then
         call ctopascl(remotename(36),location(nv+1:i))
      else
         remotename(36) = 0
      endif
c set PPCStartPBRec
c set ioCompletion code to NIL
   70 startrec(7) = 0
      startrec(8) = 0
c set portRefNum
      startrec(20) = prnum
c service type = ppcServiceRealTime and resFlag = 0
      startrec(23) = 256
c set pointer to ppcPortRec
      longw = loc(ppcport)
      startrec(24) = intw(1)
      startrec(25) = intw(2)
c set pointer to locationNameRec
      longw = loc(remotename)
      startrec(26) = intw(1)
      startrec(27) = intw(2)
c set rejectInfo to Null
      startrec(28) = 0
      startrec(29) = 0
c set userData to idproc
      startrec(30) = 0
      startrec(31) = nproc
c set userRefNum to guest
      startrec(32) = 0
      startrec(33) = 0
c initiate a PPC session
      ierr = PPCStartSync(startrec)
      if (ierr.ne.noerr) then
         write (2,*) 'PPCStart failed, ierr=', ierr
         write (2,*) 'Trying to start ',location
         write (2,*) 'rejectInfo=', startrec(29)
         if (startrec(29).eq.1) then
            write (2,*) 'Port name not recognized'
         endif
         if (ierr.eq.(-906)) then
            write (2,*) 'Port does not exist at destination'
         elseif (ierr.eq.(-912)) then
            write (2,*) 'Destination rejected the session request'
         elseif (ierr.eq.(-915)) then
            write (2,*) 'Unable to contact application'
         elseif (ierr.eq.(-926)) then
            write (2,*) 'Target application had no inform pending'
            write (2,*) 'nodelist may contain name of master node'
         elseif (ierr.eq.(-932)) then
            write (2,*) 'Destination port requires authentication'
            write (2,*) 'Guest access probably not enabled'
         endif
         ierror = ierr
         call MPI_FINALIZE(nerr)
         return
      endif
c for connections to oneself, sends use the inform session reference
c number and receives use the start session reference number
      if (idproc.eq.nproc) then
c extract and store session reference number for sends to oneself
         intw(1) = informrec(21)
         intw(2) = informrec(22)
         srnum(MAXS+1) = longw
      endif
      nproc = nproc + 1
c check for processor number overflow
      if (nproc.gt.MAXS) then
         write (2,*) 'processor number overflow, nproc = ', nproc
         ierror = 6
         call MPI_FINALIZE(nerr)
         return
      endif
c extract session reference number
      intw(1) = startrec(21)
      intw(2) = startrec(22)
      srnum(nproc) = longw
c debug
      if (monitor.eq.2) then
         write (2,*) 'connection confirmed with idproc=', nproc-1
      endif
c pass current location to next node
      if (nproc.gt.(idproc+2)) then
         call MPI_SEND(location,16,MPI_INTEGER,idproc+1,1,0,nerr)
      endif
c read location of next node from file
      if (idproc.eq.0) then
         if ((nproc.ge.2).or.(location.eq.'self')) then
            if (nofile.ne.0) go to 80
            read (3,'(a64)',iostat=nerr) location
c end of file
            if ((nerr.ne.0).or.(location.eq.' ')) go to 80
         endif
c receive location of next node from another processor
      else
         call MPI_RECV(location,16,MPI_INTEGER,idproc-1,1,0,stat,nerr)
c end of file marker received
         if (stat(4).eq.0) go to 80
      endif
c start another connection
      go to 40
c 
c * * * end main iteration loop * * *
c
c all expected nodes activated
   80 nv = nproc - 1
c debug
      if (monitor.eq.2) then
         write (2,*) 'all nodes activated: idproc, nproc=', idproc,nproc
      endif
c send null record to next processor
      if (idproc.lt.nv) then
         call MPI_SEND(location,0,MPI_INTEGER,idproc+1,1,0,nerr)
      endif
      if (nofile.eq.0) close(unit=3)
c check number of processors
      if (idproc.eq.nv) then
         do 90 i = 1, nv
         call MPI_SEND(nproc,1,MPI_INTEGER,nv-i,2,0,nerr)
   90    continue
      else
         call MPI_RECV(response,1,MPI_INTEGER,nv,2,0,stat,nerr)
c local processor does not agree with last processor on total number
         if (response.ne.nproc) then
            write (2,*) 'processor number error, local/remote nproc = '
     1, nproc, response
            ierror = 7
            call MPI_FINALIZE(nerr)
            return
         endif
      endif
c create window for showing MPI message status
      if (monitor.gt.0) then
         call messwin(nproc)
         call checkesc(1)
      endif
c set errror code to success
      ierror = 0
      return
      end
c-----------------------------------------------------------------------
      subroutine clrport(portname,prnum)
c this subroutine checks to see if requested port name has already been
c used, and if so, it closes the old port.  this is necessary because
c applications which terminate abnormally can leave open ports.
c a file is created whose name is the port name and whose only contents
c is the port reference number (needed to close the port)
c if prnum = 0, the file is read (if it exists), and the port closed.
c otherwise, the current port reference number is written to the file
      implicit none
      character*(*) portname
      integer*2 prnum
c function declarations
      integer*2 PPCCloseSync
      external PPCCloseSync
c local data
      integer*2 closerec(20), ierr
      integer lc, nofile, nerr
      lc = len(portname)
      if (lc.lt.1) return
c omit trailing blanks
   10 if (portname(lc:lc).eq.' ')  then
         lc = lc - 1
         if (lc.gt.0) go to 10
      endif
      if (lc.lt.1) return
      if (prnum.ne.0) then
         open(unit=4,file=portname(1:lc),form='formatted',status='unknow
     1n')
         write (4,*) prnum
         close(unit=4)
      else
         open(unit=4,file=portname(1:lc),form='formatted',status='old',i
     1ostat=nofile)
         if (nofile.eq.0) then
            read (4,*,iostat=nerr) prnum
c set PPCClosePBRec structure
c set ioCompletion to NIL
            closerec(7) = 0
            closerec(8) = 0
c set port reference number
            closerec(20) = prnum
c close a PPC Port synchronously
            ierr = PPCCloseSync(closerec)
            if (ierr.eq.0) write (2,*) 'closed old port:',portname(1:lc)
            close(unit=4)
         endif
      endif
      return
      end
c-----------------------------------------------------------------------
      subroutine ctopascl(ichr,chr)
c this subroutine converts fortran characters into pascal style strings
c stored as integer*2 structures.  trailing blanks are omitted.
c ichr = output integer structure
c chr = input characters
      implicit none
      integer*2 ichr(*)
      character*(*) chr
c local data
c nd = number of words of output written
      integer*2 nd
      integer lc, nds, i, ii
      lc = len(chr)
      nd = 0
      if (lc.lt.1) return
c omit trailing blanks
   10 if (chr(lc:lc).eq.' ')  then
         lc = lc - 1
         if (lc.gt.0) go to 10
      endif
      lc = lc - 1
      if (lc.lt.0) return
      nds = lc/2
c set length and copy first character
      ichr(1) = 256*(lc+1) + ichar(chr(1:1))
c copy remaining characters in blocks of two
      do 20 i = 1, nds
      ii = 2*i
      ichr(i+1) = 256*ichar(chr(ii:ii)) + ichar(chr(ii+1:ii+1))
   20 continue
c set output length
      nd = nds + 1
c copy last character if necessary
      if (lc.gt.(2*nds)) then
         nd = nd + 1
         ichr(nd) = 256*ichar(chr(lc+1:lc+1))
      endif
      return
      end
c-----------------------------------------------------------------------
      function ioresult(pblock)
c this function returns ioResult for asynchronous PPC procedures
c input: pblock
      implicit none
      integer*2 ioresult
      integer*2 pblock(19)
      ioresult = pblock(9)
      return
      end
c-----------------------------------------------------------------------
      logical function checkesc(stk)
c this procedure allows user to abort a procedure by checking for
c escape, Cmd-. or Ctrl-C keystrokes.  Calling EventAvail also
c permits an idle procedure to time-share and checks for Quit Events
c returns true if an escape event occurred.
c stk = maximum number of sleepTicks (sixtieths of a second) that
c application agrees to relinquish the processor if no events are 
c pending for it.
c input: stk
      implicit none
      integer stk
c function declarations
      integer*2 WaitNextEvent, FindWindow
      external WaitNextEvent, FindWindow
c common block for message window
c cpptr = pointer to window structure
c crect = current drag region
      integer cpptr
      integer*2 crect(4), nsp, nbx, nds
      common /winmess/ cpptr, crect, nsp, nbx, nds
c MPI constants
      integer MPI_COMM_WORLD
      parameter(MPI_COMM_WORLD=0)
c local data
      integer*2 mouseDown, keyDown, autoKey, updateEvt, kHighLevelEvent
c myEventMask looks for mouse, keyboard, update, and quit events
      parameter(mouseDown=1,keyDown=3,autoKey=5)
      parameter(updateEvt=6,kHighLevelEvent=23)
      integer*2 myEventMask
      integer*2 event(8), intw(2)
      integer key, longw, nvp
c longw and intw are used to convert types between integer*2 and integer
      equivalence(longw,intw)
      checkesc = .false.
c myEventMask looks for mouse, keyboard, and quit events
      myEventMask = 1086
c if monitor window is open, look for update events also
      if (cpptr.ne.0) myEventMask = myEventMask + 64
c receive next event from event manager
      if (WaitNextEvent(val2(myEventMask),event,val4(stk),val4(0))) then
         if ((event(1).eq.keyDown).or.(event(1).eq.autoKey)) then
c check for escape key
            key = event(3) - 256*(event(3)/256)
            if (key.eq.27) then
               checkesc = .true.
c check for Cmd-.
            elseif (key.eq.46) then
               if ((event(8)/256).ne.(2*(event(8)/512))) then
                  checkesc = .true.
               endif
c check for Ctrl-C
            elseif (key.eq.3) then
               checkesc = .true.
            endif
c check for 'QuitApplication' Apple Event
         elseif (event(1).eq.kHighLevelEvent) then
c check if event(6) = 'qu' and event(7) = 'it'
            if ((event(6).eq.29045).and.(event(7).eq.26996)) then
               write (2,*) 'Quit Application Apple Event received'
               call MPI_ABORT(MPI_COMM_WORLD,998,key)
               stop
            endif
c check for update events
         elseif (event(1).eq.updateEvt) then
c get window pointer
            intw(1) = event(2)
            intw(2) = event(3)
            if (cpptr.eq.longw) then
c signal start of window update
               call BeginUpdate(val4(cpptr))
               call MPI_COMM_SIZE(MPI_COMM_WORLD,nvp,key)
               call messwin(nvp)
c signal end of update after BeginUpdate
               call EndUpdate(val4(cpptr))
            endif
c check for drag window event
         elseif (event(1).eq.mouseDown) then
c longw = position in global coordinates where mouse event took place
            intw(1) = event(6)
            intw(2) = event(7)
c see which window part, including menu bar, is at a point
            if (FindWindow(val4(longw),key).eq.4) then
               if (cpptr.eq.key) then
c track the mouse and move a window
                  call DragWindow(val4(cpptr),val4(longw),crect)
               endif
            endif
         endif
      endif
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_FINALIZE(ierror)
c terminate MPI execution environment
c ierror = error indicator
c output: ierror
      implicit none
      integer ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
c srnum = array of session reference numbers for each participating node
c prnum = port reference number
      common /mpiparms/ nproc, idproc, srnum, prnum
c declare common block for non-blocking messages
      integer MAXM, curreq, monitor
      parameter(MAXM=16)
      dimension curreq(5,MAXM)
      integer*2 readrec(35,MAXM)
c monitor = (0,1,2) = (suppress,display,display & log) monitor messages
      common /mpisendrec/ curreq, readrec, monitor
c function declarations
      integer*2 PPCEndSync, PPCCloseSync
      external PPCEndSync, PPCCloseSync
c local data
      integer*2 intw(2), endrec(22), closerec(20), ierr, noerr
      integer longw, i
      parameter(noerr=0)
c longw and intw are used to convert types between integer*2 and integer
      equivalence(longw,intw)
c set PPCEndPBRec
      do 10 i = 1, nproc
c set ioCompletion code to NIL
      endrec(7) = 0
      endrec(8) = 0
c set session reference number
      longw = srnum(i)
      endrec(21) = intw(1)
      endrec(22) = intw(2)
c end a PPC Session
      ierr = PPCEndSync(endrec)
      if (ierr.ne.noerr) then
         write (2,*) 'PPCEnd failed, i, ierr=', i, ierr
      endif
  10  continue
c set PPCClosePBRec structure
c set ioCompletion to NIL
      closerec(7) = 0
      closerec(8) = 0
c set port reference number
      closerec(20) = prnum
      ierror = 0
c MPI already finalized
      if (nproc.le.0) ierror = 1
c close a PPC Port synchronously
      ierr = PPCCloseSync(closerec)
      if (ierr.ne.noerr) then
         write (2,*) 'PPCClose Failed, ierr, prnum = ', ierr, prnum
         ierror = ierr
      endif
c close window for showing MPI message status
      if (monitor.gt.0) then
         call logmess(0,0,-1,0)
         call delmess()
      endif
c nullify nproc
      nproc = 0
c nullify session reference numbers
      do 20 i = 1, MAXS
      srnum(i) = 0
   20 continue
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_SEND(buf,count,datatype,dest,tag,comm,ierror)
c blocking standard mode send
c buf = initial address of send buffer
c count = number of entries to send
c datatype = datatype of each entry
c dest = rank of destination
c tag = message tag
c comm = communicator (only MPI_COMM_WORLD currently supported)
c ierror = error indicator
c input: buf, count, datatype, dest, tag, comm
c output: ierror
      implicit none
      integer buf(*)
      integer count, datatype, dest, tag, comm, ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
c idproc = processor id
c srnum = array of session reference numbers for each participating node
      common /mpiparms/ nproc, idproc, srnum, prnum
c declare common block for non-blocking messages
      integer MAXM, curreq, monitor
      parameter(MAXM=16)
      dimension curreq(5,MAXM)
      integer*2 readrec(35,MAXM)
c monitor = (0,1,2) = (suppress,display,display & log) monitor messages
      common /mpisendrec/ curreq, readrec, monitor
c function declarations
      integer*2 PPCWriteAsync, ioresult
      logical checkesc
      external PPCWriteAsync
c MPI constants
      integer MPI_COMM_WORLD
      parameter(MPI_COMM_WORLD=0)
      integer MPI_INTEGER, MPI_REAL, MPI_DOUBLE_PRECISION
      parameter(MPI_INTEGER=18,MPI_REAL=19,MPI_DOUBLE_PRECISION=20)
      integer MPI_COMPLEX, MPI_DOUBLE_COMPLEX, MPI_BYTE
      parameter(MPI_COMPLEX=22,MPI_DOUBLE_COMPLEX=23,MPI_BYTE=2)
c local data
      integer*2 intw(2), ierr
      integer*2 noerr, notdone
      parameter(noerr=0,notdone=1)
      integer*2 writerec(35)
      integer longw
      integer slen, rlen, nerr
c longw and intw are used to convert types between integer*2 and integer
      equivalence(longw,intw)
      ierror = 0
c check for error conditions
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c invalid comm
      elseif (comm.gt.MPI_COMM_WORLD) then
         ierror = 2
c invalid count
      elseif (count.lt.0) then
         ierror = 3
c invalid destination
      elseif ((dest.lt.0).or.(dest.ge.nproc)) then
         write (2,*) 'destination = ', dest
         ierror = 4
c invalid tag
      elseif (tag.lt.(-1)) then
         ierror = 6
      endif
c handle errors
      if (ierror.ne.0) then
         call writerrs('MPI_SEND: ',ierror)
         return
      endif
c set PPCWritePBRec
c set ioCompletion code to NIL
      writerec(7) = 0
      writerec(8) = 0
c for connections to oneself, sends use the inform session reference
      if (idproc.eq.dest) then
         longw = srnum(MAXS+1)
      else
c normal reference number
         longw = srnum(dest+1)
      endif
c set session reference number
      writerec(21) = intw(1)
      writerec(22) = intw(2)
c find buffer length
      if ((datatype.eq.MPI_INTEGER).or.(datatype.eq.MPI_REAL)) then
         slen = 4*count
      elseif ((datatype.eq.MPI_DOUBLE_PRECISION).or.(datatype.eq.MPI_COM
     1PLEX)) then
         slen = 8*count
      elseif (datatype.eq.MPI_DOUBLE_COMPLEX) then
         slen = 16*count
      elseif (datatype.eq.MPI_BYTE) then
         slen = count
c invalid datatype
      else
         ierror = 7
         call writerrs('MPI_SEND: ',ierror)
         return
      endif
c set buffer length
      longw = slen
      writerec(23) = intw(1)
      writerec(24) = intw(2)
c set buffer pointer
      longw = loc(buf)
      writerec(27) = intw(1)
      writerec(28) = intw(2)
c set more flag to FALSE
      writerec(29) = 0
c set userData to tag
      longw = tag
      writerec(30) = intw(1)
      writerec(31) = intw(2)
c set blockCreator to comm
      longw = comm
      writerec(32) = intw(1)
      writerec(33) = intw(2)
c set blockType to datatype
      longw = datatype
      writerec(34) = intw(1)
      writerec(35) = intw(2)
c write to an application during a ppc session
      ierr = PPCWriteAsync(writerec)
c check for write errors
      if (ierr.ne.noerr) then
         write (2,*) 'PPCWrite Error, ierr, dest, tag = ', ierr, dest
     1, tag
         if (ierr.eq.(-917)) then
            write (2,*) 'The session has closed'
         endif
         ierror = ierr
         call writerrs('MPI_SEND: ',ierror)
         return
      endif
c log MPI message state change and displays status
      if (monitor.gt.0) call logmess(dest,1,slen,tag)
c wait for connection
   10 if (ioresult(writerec).eq.notdone) then
         if (checkesc(0)) then
            call MPI_ABORT(MPI_COMM_WORLD,999,nerr)
            stop
         else
            go to 10
         endif
      else
         ierr = ioresult(writerec)
      endif
c get actual length
      intw(1) = writerec(25)
      intw(2) = writerec(26)
      rlen = longw
c check for write errors
      if (ierr.ne.noerr) then
         write (2,*) 'PPCWrite Error, ierr, dest, tag = ', ierr, dest
     1, tag
         if (ierr.eq.(-917)) then
            write (2,*) 'The session has closed'
         endif
         ierror = ierr
         call writerrs('MPI_SEND: ',ierror)
      else
         if (rlen.ne.slen) then
            write(2,*) 'Send Length Error, dest, tag, requested/actual l
     1ength = ', dest, tag, slen, rlen
            ierror = 8
            call writerrs('MPI_SEND: ',ierror)
         endif
      endif
c log MPI message state change and display status
      if (monitor.gt.0) then
         if (ierror.eq.0) call logmess(dest,-1,rlen,tag)
      endif
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_RECV(buf,count,datatype,source,tag,comm,status,ierr
     1or)
c blocking receive
c buf = initial address of receive buffer
c count = maximum number of entries to receive
c datatype = datatype of each entry
c source = rank of source
c tag = message tag
c comm = communicator (only MPI_COMM_WORLD currently supported)
c status = return status
c ierror = error indicator
c input: count, datatype, source, tag, comm
c output: buf, status, ierror
      implicit none
      integer buf(*), status(*)
      integer count, datatype, source, tag, comm, ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
c srnum = array of session reference numbers for each participating node
      common /mpiparms/ nproc, idproc, srnum, prnum
c declare common block for non-blocking messages
      integer MAXM, curreq, monitor
      parameter(MAXM=16)
      dimension curreq(5,MAXM)
      integer*2 readwrec(35,MAXM)
c monitor = (0,1,2) = (suppress,display,display & log) monitor messages
      common /mpisendrec/ curreq, readwrec, monitor
c function declarations
      integer*2 PPCReadAsync, ioresult
      logical checkesc
      external PPCReadAsync
c MPI constants
      integer MPI_COMM_WORLD, MPI_ANY_SOURCE
      parameter(MPI_COMM_WORLD=0,MPI_ANY_SOURCE=-1)
      integer MPI_INTEGER, MPI_REAL, MPI_DOUBLE_PRECISION
      parameter(MPI_INTEGER=18,MPI_REAL=19,MPI_DOUBLE_PRECISION=20)
      integer MPI_COMPLEX, MPI_DOUBLE_COMPLEX, MPI_BYTE
      parameter(MPI_COMPLEX=22,MPI_DOUBLE_COMPLEX=23,MPI_BYTE=2)
c local data
      integer*2 intw(2), ierr
      integer*2 noerr, notdone
      parameter(noerr=0,notdone=1)
      integer*2 readrec(35)
      integer longw
      integer slen, rlen, rcomm, rtag, rdatat, nerr
c longw and intw are used to convert types between integer*2 and integer
      equivalence(longw,intw)
      ierror = 0
c check for error conditions
c set status to empty
      status(1) = -1
      status(2) = -1
      status(4) = 0
      status(5) = 0
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c invalid comm
      elseif (comm.gt.MPI_COMM_WORLD) then
         ierror = 2
c invalid count
      elseif (count.lt.0) then
         ierror = 3
c invalid source
      elseif ((source.lt.0).or.(source.ge.nproc)) then
         if (source.eq.MPI_ANY_SOURCE) then
            write (2,*) 'MPI_ANY_SOURCE not supported'
         else
            write (2,*) 'source = ', source
         endif
         ierror = 5
c invalid tag
      elseif (tag.lt.(-1)) then
         ierror = 6
      endif
c handle errors
      if (ierror.ne.0) then
         call writerrs('MPI_RECV: ',ierror)
         return
      endif
c set PPCReadPBRec
c set ioCompletion code to NIL
      readrec(7) = 0
      readrec(8) = 0
c set session reference number
      longw = srnum(source+1)
      readrec(21) = intw(1)
      readrec(22) = intw(2)
c find buffer length
      if ((datatype.eq.MPI_INTEGER).or.(datatype.eq.MPI_REAL)) then
         slen = 4*count
      elseif ((datatype.eq.MPI_DOUBLE_PRECISION).or.(datatype.eq.MPI_COM
     1PLEX)) then
         slen = 8*count
      elseif (datatype.eq.MPI_DOUBLE_COMPLEX) then
         slen = 16*count
      elseif (datatype.eq.MPI_BYTE) then
         slen = count
c invalid datatype
      else
         ierror = 7
         call writerrs('MPI_RECV: ',ierror)
         return
      endif
c set buffer length
      longw = slen
      readrec(23) = intw(1)
      readrec(24) = intw(2)
c set buffer pointer
      longw = loc(buf)
   10 readrec(27) = intw(1)
      readrec(28) = intw(2)
c read incoming data from an application
      ierr = PPCReadAsync(readrec)
c check for read errors
      if (ierr.ne.noerr) then
         write (2,*) 'PPCRead Error, ierr, source, tag = ', ierr, source
     1, tag
         if (ierr.eq.(-917)) then
            write (2,*) 'The session has closed'
         endif
         ierror = ierr
         status(3) = ierror
         call writerrs('MPI_RECV: ',ierror)
         return
      endif
      status(1) = source
c log MPI message state change and displays status
      if (monitor.gt.0) call logmess(source,2,slen,tag)
c wait for connection
   20 if (ioresult(readrec).eq.notdone) then
         if (checkesc(0)) then
            call MPI_ABORT(MPI_COMM_WORLD,999,nerr)
            stop
         else
            go to 20
         endif
      else
         ierr = ioresult(readrec)
      endif
c get received tag from userData
      intw(1) = readrec(30)
      intw(2) = readrec(31)
      rtag = longw
      status(2) = rtag
c get actual length
      intw(1) = readrec(25)
      intw(2) = readrec(26)
      rlen = longw
      status(4) = status(4) + rlen
c get received comm from blockCreator
      intw(1) = readrec(32)
      intw(2) = readrec(33)
      rcomm = longw
c get received datatype from blockType
      intw(1) = readrec(34)
      intw(2) = readrec(35)
      rdatat = longw
      status(5) = rdatat
c check for read errors
      if (ierr.ne.noerr) then
         write (2,*) 'PPCRead Error, ierr, source, tag = ', ierr, source
     1, tag
         if (ierr.eq.(-917)) then
            write (2,*) 'The session has closed'
         endif
         ierror = ierr
c comm error
      elseif (rcomm.ne.comm) then
         write (2,*) 'Read Comm Error, source, tag, expected/received co
     1mm = ', source, tag, comm, rcomm
         ierror = 9
c tag error
      elseif ((tag.ge.0).and.(rtag.ne.tag)) then
         write (2,*) 'Read Tag Error, source, expected/received tag = '
     1, source, tag, rtag
         ierror = 10
c type error
      elseif (rdatat.ne.datatype) then
         write (2,*) 'Read Type Error, source, tag, expected/received ty
     1pe = ', source, tag, datatype, rdatat
         ierror = 11
c incomplete data error
      elseif ((readrec(29)/256).ne.0) then
         write (2,*) 'Incomplete Read, source, tag, requested/actual = '
     1, source, tag, slen, rlen
         slen = slen - rlen
c this case is a workaround for a PPC Toolbox bug
         if (slen.eq.0) then
            readrec(23) = 0
            readrec(24) = 4
            longw = loc(nerr)
            go to 10
c get more data
         elseif (slen.gt.0) then
            longw = slen
            readrec(23) = intw(1)
            readrec(24) = intw(2)
            longw = loc(buf) + rlen
            go to 10
         else
            ierror = 12
         endif
c length error
      else
         if (rlen.gt.slen) then
            write(2,*) 'Read Length Error, source, tag, requested/actual
     1 = ', source, tag, slen, rlen
            ierror = 13
         endif
      endif
c store error code
      status(3) = ierror
c log MPI message state change and display status
      if (monitor.gt.0) then
         if (ierror.eq.0) call logmess(source,-2,rlen,tag)
      endif
c handle read errors
      if (ierror.ne.0) call writerrs('MPI_RECV: ',ierror)
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_ISEND(buf,count,datatype,dest,tag,comm,request,ierr
     1or)
c start a non-blocking send
c buf = initial address of send buffer
c count = number of entries to send
c datatype = datatype of each entry
c dest = rank of destination
c tag = message tag
c comm = communicator (only MPI_COMM_WORLD currently supported)
c request = request handle
c ierror = error indicator
c input: buf, count, datatype, dest, tag, comm
c output: request, ierror
      implicit none
      integer buf(*)
      integer count, datatype, dest, tag, comm, request, ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
c idproc = processor id
c srnum = array of session reference numbers for each participating node
      common /mpiparms/ nproc, idproc, srnum, prnum
c declare common block for non-blocking messages
      integer MAXM, curreq, monitor
      parameter(MAXM=16)
      dimension curreq(5,MAXM)
      integer*2 writerec(35,MAXM)
c curreq = request record for transmission parameters
c writerec = PPCWritePB Record
c monitor = (0,1,2) = (suppress,display,display & log) monitor messages
      common /mpisendrec/ curreq, writerec, monitor
c function declarations
      integer*2 PPCWriteAsync
      external PPCWriteAsync
c MPI constants
      integer MPI_COMM_WORLD
      parameter(MPI_COMM_WORLD=0)
      integer MPI_INTEGER, MPI_REAL, MPI_DOUBLE_PRECISION
      parameter(MPI_INTEGER=18,MPI_REAL=19,MPI_DOUBLE_PRECISION=20)
      integer MPI_COMPLEX, MPI_DOUBLE_COMPLEX, MPI_BYTE
      parameter(MPI_COMPLEX=22,MPI_DOUBLE_COMPLEX=23,MPI_BYTE=2)
c local data
      integer*2 intw(2), ierr
      integer*2 noerr
      parameter(noerr=0)
      integer longw
      integer i, slen
c longw and intw are used to convert types between integer*2 and integer
      equivalence(longw,intw)
      ierror = 0
c find space for record
      i = 0
   10 i = i + 1
      if (i.gt.MAXM) then
         write (2,*) 'too many sends waiting, dest, tag = ', dest, tag
         request = -1
         ierror = 14
         call writerrs('MPI_ISEND: ',ierror)
         return
      elseif (curreq(1,i).ne.0) then
         go to 10
      endif 
c check for error conditions
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c invalid comm
      elseif (comm.gt.MPI_COMM_WORLD) then
         ierror = 2
c invalid count
      elseif (count.lt.0) then
         ierror = 3
c invalid destination
      elseif ((dest.lt.0).or.(dest.ge.nproc)) then
         write (2,*) 'destination = ', dest
         ierror = 4
c invalid tag
      elseif (tag.lt.(-1)) then
         ierror = 6
      endif
c handle errors
      if (ierror.ne.0) then
         call writerrs('MPI_ISEND: ',ierror)
         return
      endif
c set PPCWritePBRec
c set ioCompletion code to NIL
      writerec(7,i) = 0
      writerec(8,i) = 0
c for connections to oneself, sends use the inform session reference
      if (idproc.eq.dest) then
         longw = srnum(MAXS+1)
      else
c normal reference number
         longw = srnum(dest+1)
      endif
c set session reference number
      writerec(21,i) = intw(1)
      writerec(22,i) = intw(2)
c find buffer length
      if ((datatype.eq.MPI_INTEGER).or.(datatype.eq.MPI_REAL)) then
         slen = 4*count
      elseif ((datatype.eq.MPI_DOUBLE_PRECISION).or.(datatype.eq.MPI_COM
     1PLEX)) then
         slen = 8*count
      elseif (datatype.eq.MPI_DOUBLE_COMPLEX) then
         slen = 16*count
      elseif (datatype.eq.MPI_BYTE) then
         slen = count
c invalid datatype
      else
         ierror = 7
         call writerrs('MPI_ISEND: ',ierror)
         return
      endif
c set buffer length
      longw = slen
      writerec(23,i) = intw(1)
      writerec(24,i) = intw(2)
c set buffer pointer
      longw = loc(buf)
      writerec(27,i) = intw(1)
      writerec(28,i) = intw(2)
c set more flag to FALSE
      writerec(29,i) = 0
c set userData to tag
      longw = tag
      writerec(30,i) = intw(1)
      writerec(31,i) = intw(2)
c set blockCreator to comm
      longw = comm
      writerec(32,i) = intw(1)
      writerec(33,i) = intw(2)
c set blockType to datatype
      longw = datatype
      writerec(34,i) = intw(1)
      writerec(35,i) = intw(2)
c write to an application during a ppc session
      ierr = PPCWriteAsync(writerec(1,i))
c check for write errors
      if (ierr.ne.noerr) then
         write (2,*) 'PPCWrite Error, ierr, dest, tag = ', ierr, dest
     1, tag
         if (ierr.eq.(-917)) then
            write (2,*) 'The session has closed'
         endif
         ierror = ierr
         call writerrs('MPI_ISEND: ',ierror)
      endif
c log MPI message state change and displays status
      if (monitor.gt.0) call logmess(dest,1,slen,tag)
c save transmission mode as send
      curreq(1,i) = -1
c save destination/source id
      curreq(2,i) = dest
c save communicator
      curreq(3,i) = comm
c save tag
      curreq(4,i) = tag
c save datatype
      curreq(5,i) = datatype
c assign request handle
      request = i
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_IRECV(buf,count,datatype,source,tag,comm,request,ie
     1rror)
c begin a non-blocking receive
c buf = initial address of receive buffer
c count = maximum number of entries to receive
c datatype = datatype of each entry
c source = rank of source
c tag = message tag
c comm = communicator (only MPI_COMM_WORLD currently supported)
c request = request handle
c ierror = error indicator
c input: count, datatype, source, tag, comm
c output: buf, request, ierror
      implicit none
      integer buf(*)
      integer count, datatype, source, tag, comm, request, ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
c srnum = array of session reference numbers for each participating node
      common /mpiparms/ nproc, idproc, srnum, prnum
c declare common block for non-blocking messages
      integer MAXM, curreq, monitor
      parameter(MAXM=16)
      dimension curreq(5,MAXM)
      integer*2 readrec(35,MAXM)
c curreq = request record for transmission parameters
c readrec = PPCReadPB Record
c monitor = (0,1,2) = (suppress,display,display & log) monitor messages
      common /mpisendrec/ curreq, readrec, monitor
c function declarations
      integer*2 PPCReadAsync
      external PPCReadAsync
c MPI constants
      integer MPI_COMM_WORLD, MPI_ANY_SOURCE
      parameter(MPI_COMM_WORLD=0,MPI_ANY_SOURCE=-1)
      integer MPI_INTEGER, MPI_REAL, MPI_DOUBLE_PRECISION
      parameter(MPI_INTEGER=18,MPI_REAL=19,MPI_DOUBLE_PRECISION=20)
      integer MPI_COMPLEX, MPI_DOUBLE_COMPLEX, MPI_BYTE
      parameter(MPI_COMPLEX=22,MPI_DOUBLE_COMPLEX=23,MPI_BYTE=2)
c local data
      integer*2 intw(2), ierr
      integer*2 noerr
      parameter(noerr=0)
      integer longw
      integer i, slen
c longw and intw are used to convert types between integer*2 and integer
      equivalence(longw,intw)
      ierror = 0
c find space for record
      i = 0
   10 i = i + 1
      if (i.gt.MAXM) then
         write (2,*) 'too many receives waiting, source, tag = ', source
     1, tag
         request = -1
         ierror = 15
         call writerrs('MPI_IRECV: ',ierror)
         return
      elseif (curreq(1,i).ne.0) then
         go to 10
      endif
c check for error conditions
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c invalid comm
      elseif (comm.gt.MPI_COMM_WORLD) then
         ierror = 2
c invalid count
      elseif (count.lt.0) then
         ierror = 3
c invalid source
      elseif ((source.lt.0).or.(source.ge.nproc)) then
         if (source.eq.MPI_ANY_SOURCE) then
            write (2,*) 'MPI_ANY_SOURCE not supported'
         else
            write (2,*) 'source = ', source
         endif
         ierror = 5
c invalid tag
      elseif (tag.lt.(-1)) then
         ierror = 6
      endif
c handle errors
      if (ierror.ne.0) then
         call writerrs('MPI_IRECV: ',ierror)
         return
      endif
c set PPCReadPBRec
c set ioCompletion code to NIL
      readrec(7,i) = 0
      readrec(8,i) = 0
c set session reference number
      longw = srnum(source+1)
      readrec(21,i) = intw(1)
      readrec(22,i) = intw(2)
c find buffer length
      if ((datatype.eq.MPI_INTEGER).or.(datatype.eq.MPI_REAL)) then
         slen = 4*count
      elseif ((datatype.eq.MPI_DOUBLE_PRECISION).or.(datatype.eq.MPI_COM
     1PLEX)) then
         slen = 8*count
      elseif (datatype.eq.MPI_DOUBLE_COMPLEX) then
         slen = 16*count
      elseif (datatype.eq.MPI_BYTE) then
         slen = count
c invalid datatype
      else
         ierror = 7
         call writerrs('MPI_IRECV: ',ierror)
         return
      endif
c set buffer length
      longw = slen
      readrec(23,i) = intw(1)
      readrec(24,i) = intw(2)
c set buffer pointer
      longw = loc(buf)
      readrec(27,i) = intw(1)
      readrec(28,i) = intw(2)
c read incoming data from an application
      ierr = PPCReadAsync(readrec(1,i))
c check for read errors
      if (ierr.ne.noerr) then
         write (2,*) 'PPCRead Error, ierr, source, tag = ', ierr, source
     1, tag
         if (ierr.eq.(-917)) then
            write (2,*) 'The session has closed'
         endif
         ierror = ierr
         call writerrs('MPI_IRECV: ',ierror)
         return
      endif
c log MPI message state change and displays status
      if (monitor.gt.0) call logmess(source,2,slen,tag)
c save transmission mode as receive
      curreq(1,i) = 1
c save destination/source id
      curreq(2,i) = source
c save communicator
      curreq(3,i) = comm
c save tag
      curreq(4,i) = tag
c save datatype
      curreq(5,i) = datatype
c assign request handle
      request = i
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_TEST(request,flag,status,ierror)
c check to see if a nonblocking send or receive operation has completed
c request = request handle
c flag = true if operation completed
c status = status object
c ierror = error indicator
c input: request
c output: request, flag, status, ierror
      implicit none
      integer status(*)
      integer request, ierror
      logical flag
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
      common /mpiparms/ nproc, idproc, srnum, prnum
c declare common block for non-blocking messages
      integer MAXM, curreq, monitor
      parameter(MAXM=16)
      dimension curreq(5,MAXM)
      integer*2 readwrec(35,MAXM)
c curreq = request record for transmission parameters
c readwrec = PPCReadPB or PPCWritePB Record
c monitor = (0,1,2) = (suppress,display,display & log) monitor messages
      common /mpisendrec/ curreq, readwrec, monitor
c function declarations
      integer*2 ioresult, PPCReadAsync
      logical checkesc
      external PPCReadAsync
c MPI constants
      integer MPI_COMM_WORLD
      parameter(MPI_COMM_WORLD=0)
c local data
      integer*2 intw(2), ierr
      integer*2 noerr, notdone
      parameter(noerr=0,notdone=1)
      integer longw
      integer i, dest, source, slen, comm, tag, datatype
      integer rlen, rcomm, rtag, rdatat, nerr
c longw and intw are used to convert types between integer*2 and integer
      equivalence(longw,intw)
      ierror = 0
c check for error conditions
c set status to empty
      status(1) = -1
      status(2) = -1
      status(3) = 0
      status(4) = 0
      status(5) = 0
      i = request
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c null request
      elseif (request.lt.0) then
         flag = .true.
         return
c invalid request handle
      elseif ((i.lt.1).or.(i.gt.MAXM)) then
         ierror = 16
      elseif (curreq(1,i).eq.0) then
         ierror = 16
      endif
c handle errors
      if (ierror.ne.0) then
         status(3) = ierror
         call writerrs('MPI_TEST: ',ierror)
         return
      endif
      flag = .false.
c check for connection
   10 if (ioresult(readwrec(1,i)).eq.notdone) then
         if (checkesc(1)) then
            call MPI_ABORT(MPI_COMM_WORLD,999,nerr)
            stop
c if Incomplete Read bug has occurred, wait for the remainder 
         else
            if (flag) go to 10
            return
         endif
      else
         ierr = ioresult(readwrec(1,i))
         flag = .true.
      endif
c get requested length
      intw(1) = readwrec(23,i)
      intw(2) = readwrec(24,i)
      slen = longw
c get actual length
      intw(1) = readwrec(25,i)
      intw(2) = readwrec(26,i)
      rlen = longw
c read current request record
      dest = curreq(2,i)
      tag = curreq(4,i)
c check for send errors
      if (curreq(1,i).lt.0) then
c check for write errors
         if (ierr.ne.noerr) then
            write (2,*) 'PPCWrite Error, ierr, dest, tag = ', ierr, dest
     1, tag
            if (ierr.eq.(-917)) then
               write (2,*) 'The session has closed'
            endif
            ierror = ierr
         elseif (rlen.ne.slen) then
            write(2,*) 'Send Length Error, dest, tag, requested/actual l
     1ength = ', dest, tag, slen, rlen
            ierror = 8
         endif
c define length and type for MPI_GET_COUNT
         status(4) = rlen
         status(5) = curreq(5,i)
c log MPI message state change and display status
         if (monitor.gt.0) then
            if (ierror.eq.0) call logmess(dest,-1,rlen,tag)
         endif
         go to 30
      endif
c read current request record
      source = curreq(2,i)
      comm = curreq(3,i)
      datatype = curreq(5,i)
      status(1) = source
      status(4) = status(4) + rlen
c get received tag from userData
      intw(1) = readwrec(30,i)
      intw(2) = readwrec(31,i)
      rtag = longw
      status(2) = rtag
c get received comm from blockCreator
      intw(1) = readwrec(32,i)
      intw(2) = readwrec(33,i)
      rcomm = longw
c get received datatype from blockType
      intw(1) = readwrec(34,i)
      intw(2) = readwrec(35,i)
      rdatat = longw
      status(5) = rdatat
c check for read errors
      if (ierr.ne.noerr) then
         write (2,*) 'PPCRead Error, ierr, source, tag = ', ierr, source
     1, tag
         if (ierr.eq.(-917)) then
            write (2,*) 'The session has closed'
         endif
         ierror = ierr
c comm error
      elseif (rcomm.ne.comm) then
         write (2,*) 'Read Comm Error, source, tag, expected/received co
     1mm = ', source, tag, comm, rcomm
         ierror = 9
c tag error
      elseif ((tag.ge.0).and.(rtag.ne.tag)) then
         write (2,*) 'Read Tag Error, source, expected/received tag = '
     1, source, tag, rtag
         ierror = 10
c type error
      elseif (rdatat.ne.datatype) then
         write (2,*) 'Read Type Error, source, tag, expected/received ty
     1pe = ', source, tag, datatype, rdatat
         ierror = 11
c incomplete data error
      elseif ((readwrec(29,i)/256).ne.0) then
         write (2,*) 'Incomplete Read, source, tag, requested/actual = '
     1, source, tag, slen, rlen
         slen = slen - rlen
c this case is a workaround for a PPC Toolbox bug
         if (slen.eq.0) then
            readwrec(23,i) = 0
            readwrec(24,i) = 4
            longw = loc(nerr)
c get more data
         elseif (slen.gt.0) then
            longw = slen
            readwrec(23,i) = intw(1)
            readwrec(24,i) = intw(2)
            intw(1) = readwrec(27,i)
            intw(2) = readwrec(28,i)
            longw = longw + rlen
         else
            ierror = 12
            go to 30
         endif
c set buffer pointer
         readwrec(27,i) = intw(1)
         readwrec(28,i) = intw(2)
c read incoming data from an application
         ierr = PPCReadAsync(readwrec(1,i))
c check for read errors
         if (ierr.ne.noerr) then
            write (2,*) 'PPCRead Error, ierr, source, tag = ', ierr
     1, source, tag
            ierror = ierr
            go to 30
         endif
         go to 10
c length error
      elseif (rlen.gt.slen) then
         write (2,*) 'Read Length Error, source, tag, requested/actual='
     1, source, tag, slen, rlen
         ierror = 13
      endif
c log MPI message state change and display status
      if (monitor.gt.0) then
         if (ierror.eq.0) call logmess(source,-2,rlen,tag)
      endif
c store error code
   30 status(3) = ierror
c nullify transmission mode
      curreq(1,i) = 0
c nullify request handle
      request = -1
c handle read and write errors
      if (ierror.ne.0) call writerrs('MPI_TEST: ',ierror)
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_WAIT(request,status,ierror)
c wait for an MPI send or receive to complete
c request = request handle
c status = status object
c ierror = error indicator
c input: request
c output: request, status, ierror
      implicit none
      integer status(*)
      integer request, ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
      common /mpiparms/ nproc, idproc, srnum, prnum
c declare common block for non-blocking messages
      integer MAXM, curreq, monitor
      parameter(MAXM=16)
      dimension curreq(5,MAXM)
      integer*2 readwrec(35,MAXM)
c curreq = request record for transmission parameters
c readwrec = PPCReadPB or PPCWritePB Record
c monitor = (0,1,2) = (suppress,display,display & log) monitor messages
      common /mpisendrec/ curreq, readwrec, monitor
c function declarations
      integer*2 ioresult, PPCReadAsync
      logical checkesc
      external PPCReadAsync
c MPI constants
      integer MPI_COMM_WORLD
      parameter(MPI_COMM_WORLD=0)
c local data
      integer*2 intw(2), ierr
      integer*2 noerr, notdone
      parameter(noerr=0,notdone=1)
      integer longw
      integer i, dest, source, slen, comm, tag, datatype
      integer rlen, rcomm, rtag, rdatat, nerr
c longw and intw are used to convert types between integer*2 and integer
      equivalence(longw,intw)
      ierror = 0
c check for error conditions
c set status to empty
      status(1) = -1
      status(2) = -1
      status(3) = 0
      status(4) = 0
      status(5) = 0
      i = request
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c null request
      elseif (request.lt.0) then
         return
c invalid request handle
      elseif ((i.lt.1).or.(i.gt.MAXM)) then
         ierror = 16
      elseif (curreq(1,i).eq.0) then
         ierror = 16
      endif
c handle errors
      if (ierror.ne.0) then
         status(3) = ierror
         call writerrs('MPI_WAIT: ',ierror)
         return
      endif
c wait for connection
   10 if (ioresult(readwrec(1,i)).eq.notdone) then
         if (checkesc(1)) then
            call MPI_ABORT(MPI_COMM_WORLD,999,nerr)
            stop
         else
            go to 10
         endif
      else
         ierr = ioresult(readwrec(1,i))
      endif
c get requested length
      intw(1) = readwrec(23,i)
      intw(2) = readwrec(24,i)
      slen = longw
c get actual length
      intw(1) = readwrec(25,i)
      intw(2) = readwrec(26,i)
      rlen = longw
c read current request record
      dest = curreq(2,i)
      tag = curreq(4,i)
c check for send errors
      if (curreq(1,i).lt.0) then
c check for write errors
         if (ierr.ne.noerr) then
            write (2,*) 'PPCWrite Error, ierr, dest, tag = ', ierr, dest
     1, tag
            if (ierr.eq.(-917)) then
               write (2,*) 'The session has closed'
            endif
            ierror = ierr
         elseif (rlen.ne.slen) then
            write(2,*) 'Send Length Error, dest, tag, requested/actual l
     1ength = ', dest, tag, slen, rlen
            ierror = 8
         endif
c define length and type for MPI_GET_COUNT
         status(4) = rlen
         status(5) = curreq(5,i)
c log MPI message state change and display status
         if (monitor.gt.0) then
            if (ierror.eq.0) call logmess(dest,-1,rlen,tag)
         endif
      endif
c read current request record
      source = curreq(2,i)
      comm = curreq(3,i)
      datatype = curreq(5,i)
      status(1) = source
      status(4) = status(4) + rlen
c get received tag from userData
      intw(1) = readwrec(30,i)
      intw(2) = readwrec(31,i)
      rtag = longw
      status(2) = rtag
c get received comm from blockCreator
      intw(1) = readwrec(32,i)
      intw(2) = readwrec(33,i)
      rcomm = longw
c get received datatype from blockType
      intw(1) = readwrec(34,i)
      intw(2) = readwrec(35,i)
      rdatat = longw
      status(5) = rdatat
c check for read errors
      if (ierr.ne.noerr) then
         write (2,*) 'PPCRead Error, ierr, source, tag = ', ierr, source
     1, tag
         if (ierr.eq.(-917)) then
            write (2,*) 'The session has closed'
         endif
         ierror = ierr
c comm error
      elseif (rcomm.ne.comm) then
         write (2,*) 'Read Comm Error, source, tag, expected/received co
     1mm = ', source, tag, comm, rcomm
         ierror = 9
c tag error
      elseif ((tag.ge.0).and.(rtag.ne.tag)) then
         write (2,*) 'Read Tag Error, source, expected/received tag = '
     1, source, tag, rtag
         ierror = 10
c type error
      elseif (rdatat.ne.datatype) then
         write (2,*) 'Read Type Error, source, tag, expected/received ty
     1pe = ', source, tag, datatype, rdatat
         ierror = 11
c incomplete data error
      elseif ((readwrec(29,i)/256).ne.0) then
         write (2,*) 'Incomplete Read, source, tag, requested/actual = '
     1, source, tag, slen, rlen
         slen = slen - rlen
c this case is a workaround for a PPC Toolbox bug
         if (slen.eq.0) then
            readwrec(23,i) = 0
            readwrec(24,i) = 4
            longw = loc(nerr)
c get more data
         elseif (slen.gt.0) then
            longw = slen
            readwrec(23,i) = intw(1)
            readwrec(24,i) = intw(2)
            intw(1) = readwrec(27,i)
            intw(2) = readwrec(28,i)
            longw = longw + rlen
         else
            ierror = 12
            go to 30
         endif
c set buffer pointer
         readwrec(27,i) = intw(1)
         readwrec(28,i) = intw(2)
c read incoming data from an application
         ierr = PPCReadAsync(readwrec(1,i))
c check for read errors
         if (ierr.ne.noerr) then
            write (2,*) 'PPCRead Error, ierr, source, tag = ', ierr
     1, source, tag
            ierror = ierr
            go to 30
         endif
         go to 10
c length error
      elseif (rlen.gt.slen) then
         write (2,*) 'Read Length Error, source, tag, requested/actual='
     1, source, tag, slen, rlen
         ierror = 13
      endif
c log MPI message state change and display status
      if (monitor.gt.0) then
         if (ierror.eq.0) call logmess(source,-2,rlen,tag)
      endif
c store error code
   30 status(3) = ierror
c nullify transmission mode
      curreq(1,i) = 0
c nullify request handle
      request = -1
c handle read and write errors
      if (ierror.ne.0) call writerrs('MPI_WAIT: ',ierror)
      return
      end
c-----------------------------------------------------------------------
      subroutine MPE_WAIT(request,status,ierror)
c wait for an MPI send or receive to complete
c request = request handle
c status = status object
c ierror = error indicator
c input: request
c output: request, status, ierror
      implicit none
      integer status(*)
      integer request, ierror
c local data
      logical flag
   10 call MPI_TEST(request,flag,status,ierror)
      if (.not.flag) go to 10
      return
      end
c-----------------------------------------------------------------------
      subroutine MPE_SEND(buf,count,datatype,dest,tag,comm,ierror)
c blocking standard mode send
c buf = initial address of send buffer
c count = number of entries to send
c datatype = datatype of each entry
c dest = rank of destination
c tag = message tag
c comm = communicator (only MPI_COMM_WORLD currently supported)
c ierror = error indicator
c input: buf, count, datatype, dest, tag, comm
c output: ierror
      implicit none
      integer buf(*)
      integer count, datatype, dest, tag, comm, ierror
c MPI constants
      integer MPI_STATUS_SIZE
      parameter(MPI_STATUS_SIZE=5)
c local data
      integer request, status
      dimension status(MPI_STATUS_SIZE)
      call MPI_ISEND(buf,count,datatype,dest,tag,comm,request,ierror)
      call MPE_WAIT(request,status,ierror)
      return
      end
c-----------------------------------------------------------------------
      subroutine MPE_RECV(buf,count,datatype,source,tag,comm,status,ierr
     1or)
c blocking receive
c buf = initial address of receive buffer
c count = maximum number of entries to receive
c datatype = datatype of each entry
c source = rank of source
c tag = message tag
c comm = communicator (only MPI_COMM_WORLD currently supported)
c status = return status
c ierror = error indicator
c input: count, datatype, source, tag, comm
c output: buf, status, ierror
      implicit none
      integer buf(*), status(*)
      integer count, datatype, source, tag, comm, ierror
c local data
      integer request
      call MPI_IRECV(buf,count,datatype,source,tag,comm,request,ierror)
      call MPE_WAIT(request,status,ierror)
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_SENDRECV (sendbuf,sendcount,sendtype,dest,sendtag,
     1recvbuf,recvcount,recvtype,source,recvtag,comm,status,ierror)
c blocking send and receive operation
c sendbuf = initial address of send buffer
c sendcount = number of entries to send
c sendtype = type of entries in send buffer
c dest = rank of destination
c sendtag = send tag
c recvbuf = initial address of receive buffer
c recvcount = max number of entries to receive
c recvtype = type of entries in receive buffer
c source = rank of source
c recvtag = receive tag
c comm = communicator (only MPI_COMM_WORLD currently supported)
c status = return status
c ierror = error indicator
c input: sendbuf, sendcount, sendtype, dest, sendtag, recvbuf, recvcount
c        recvtype, source, recvtag, comm
c output: recvbuf, status, ierror
      implicit none
      integer sendbuf(*), recvbuf(*), status(*)
      integer sendcount, sendtype, dest, sendtag, recvcount, recvtype
      integer source, recvtag, comm, ierror
c local data
      integer recvreq, sendreq
c post non-blocking receive and send
      call MPI_IRECV(recvbuf,recvcount,recvtype,source,recvtag,comm,recv
     1req,ierror)
      call MPI_ISEND(sendbuf,sendcount,sendtype,dest,sendtag,comm,sendre
     1q,ierror)
c wait for send and receive
      call MPI_WAIT(sendreq,status,ierror)
      call MPI_WAIT(recvreq,status,ierror)
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_SSEND(buf,count,datatype,dest,tag,comm,ierror)
c blocking synchronous mode send
c buf = initial address of send buffer
c count = number of entries to send
c datatype = datatype of each entry
c dest = rank of destination
c tag = message tag
c comm = communicator (only MPI_COMM_WORLD currently supported)
c ierror = error indicator
c input: buf, count, datatype, dest, tag, comm
c output: ierror
      implicit none
      integer buf(*)
      integer count, datatype, dest, tag, comm, ierror
      call MPI_SEND(buf,count,datatype,dest,tag,comm,ierror)
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_ISSEND(buf,count,datatype,dest,tag,comm,request,ier
     1ror)
c start a non-blocking synchronous mode send
c buf = initial address of send buffer
c count = number of entries to send
c datatype = datatype of each entry
c dest = rank of destination
c tag = message tag
c comm = communicator (only MPI_COMM_WORLD currently supported)
c request = request handle
c ierror = error indicator
c input: buf, count, datatype, dest, tag, comm
c output: request, ierror
      implicit none
      integer buf(*)
      integer count, datatype, dest, tag, comm, request, ierror
      call MPI_ISEND(buf,count,datatype,dest,tag,comm,request,ierror)
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_WAITALL (count,array_of_requests,array_of_statuses,
     1ierror)
c wait for a collection of specified MPI sends or receives to complete
c count = list length
c array_of_requests = array of request handles
c array_of_statuses = array of status objects
c ierror = error indicator
c input: count, array_of_requests
c output: array_of_requests, array_of_statuses, ierror
      implicit none
c MPI constants
      integer MPI_STATUS_SIZE, MPI_ERR_IN_STATUS
      parameter(MPI_STATUS_SIZE=5,MPI_ERR_IN_STATUS=67)
      integer array_of_requests(*), array_of_statuses(MPI_STATUS_SIZE,*)
      integer count, ierror
c local data
      integer i, ierr
c invalid count
      if (count.lt.0) then
         write (2,*) 'Invalid list length = ', count
         ierror = 17
         call writerrs('MPI_WAITALL: ',ierror)
         return
      endif
      ierror = 0
      do 10 i = 1, count
      call MPI_WAIT(array_of_requests(i),array_of_statuses(1,i),ierr)
      if (ierr.ne.0) ierror = MPI_ERR_IN_STATUS
   10 continue
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_WAITANY(count,array_of_requests,index,status,ierror
     1)
c wait for any specified MPI send or receive to complete
c count = list length
c array_of_requests = array of request handles
c index = index of request handle that completed
c status = status object
c ierror = error indicator
c input: count, array_of_requests
c output: array_of_requests, index, status, ierror
      implicit none
      integer array_of_requests(*), status(*)
      integer count, index, ierror
c local data
      integer i, k
      logical flag
c invalid count
      if (count.lt.0) then
         write (2,*) 'Invalid list length = ', count
         ierror = 17
         call writerrs('MPI_WAITANY: ',ierror)
         return
      endif
c find number of requests already completed
      k = 0
      do 10 i = 1, count
      if (array_of_requests(i).lt.0) k = k + 1
   10 continue
      if (k.eq.count) then
         index = -1
         ierror = 0
         return
      endif
      i = 1
   20 flag = .false.
      if (array_of_requests(i).ge.0) then
         call MPI_TEST(array_of_requests(i),flag,status,ierror)
      endif
      if (flag) then
         index = i
      else
         i = i + 1
         if (i.gt.count) i = 1
         go to 20
      endif
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_GET_COUNT(status,datatype,count,ierror)
c get the number of "top level" elements
c status = return status of receive operation
c datatype = datatype of each receive buffer entry
c count = number of received entries
c ierror = error indicator
c input: status, datatype
c output: count, ierror
      implicit none
      integer status(*)
      integer datatype, count, ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
      common /mpiparms/ nproc, idproc, srnum, prnum
c MPI constants
      integer MPI_INTEGER, MPI_REAL, MPI_DOUBLE_PRECISION
      parameter(MPI_INTEGER=18,MPI_REAL=19,MPI_DOUBLE_PRECISION=20)
      integer MPI_COMPLEX, MPI_DOUBLE_COMPLEX, MPI_BYTE
      parameter(MPI_COMPLEX=22,MPI_DOUBLE_COMPLEX=23,MPI_BYTE=2)
      integer MPI_UNDEFINED
      parameter(MPI_UNDEFINED=-1)
      ierror = 0
      count = 0
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c mismatched datatype
      elseif (datatype.ne.status(5)) then
         ierror = 18
c calculate count
      elseif ((datatype.eq.MPI_INTEGER).or.(datatype.eq.MPI_REAL)) then
         count = status(4)/4
         if (4*count.ne.status(4)) count = MPI_UNDEFINED
      elseif ((datatype.eq.MPI_DOUBLE_PRECISION).or.(datatype.eq.MPI_COM
     1PLEX)) then
         count = status(4)/8
         if (8*count.ne.status(4)) count = MPI_UNDEFINED
      elseif (datatype.eq.MPI_DOUBLE_COMPLEX) then
         count = status(4)/16
         if (16*count.ne.status(4)) count = MPI_UNDEFINED
      elseif (datatype.eq.MPI_BYTE) then
         count = status(4)
c invalid datatype
      else
         ierror = 7
      endif
c handle errors
      if (ierror.ne.0) call writerrs('MPI_GET_COUNT: ',ierror)
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_INITIALIZED(flag,ierror)
c indicate whether MPI_init has been called
c flag = true if MPI_Init has been called, false otherwise
c ierror = error indicator
c output: flag, ierror
      implicit none
      logical flag
      integer ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
      common /mpiparms/ nproc, idproc, srnum, prnum
      if (nproc.gt.0) then
         flag = .true.
      else
         flag = .false.
      endif
      ierror = 0
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_COMM_SIZE(comm,size,ierror)
c determine the size of the group associated with a communicator
c comm = communicator (this is ignored)
c size = number of processors in the group of comm
c ierror = error indicator
c input: comm
c output: size, ierror
      implicit none
      integer comm, size, ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
      common /mpiparms/ nproc, idproc, srnum, prnum
c MPI constants
      integer MPI_COMM_WORLD
      parameter(MPI_COMM_WORLD=0)
c check for error conditions
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c invalid comm
      elseif (comm.ne.MPI_COMM_WORLD) then
         ierror = 2
c get size
      else
         size = nproc
         ierror = 0
      endif
c handle errors
      if (ierror.ne.0) call writerrs('MPI_COMM_SIZE: ',ierror)
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_COMM_RANK(comm,rank,ierror)
c determine the rank of the calling process in the communicator
c comm = communicator (this is ignored)
c rank = rank of the calling process in group of comm
c ierror = error indicator
c input: comm
c output: rank, ierror
      implicit none
      integer comm, rank, ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c idproc = processor id
      common /mpiparms/ nproc, idproc, srnum, prnum
c MPI constants
      integer MPI_COMM_WORLD
      parameter(MPI_COMM_WORLD=0)
c check for error conditions
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c invalid comm
      elseif (comm.ne.MPI_COMM_WORLD) then
         ierror = 2
c get rank
      else
         rank = idproc
         ierror = 0
      endif
c handle errors
      if (ierror.ne.0) call writerrs('MPI_COMM_RANK: ',ierror)
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_BCAST(buffer,count,datatype,root,comm,ierror)
c broadcast a message from root to all processes in comm
c buffer = starting address of buffer
c count = number of entries in buffer
c datatype = datatype of buffer
c root = rank of broadcast root
c comm = communicator (only MPI_COMM_WORLD currently supported)
c ierror = error indicator
c input: buffer, count, datatype, root, comm
c output: buffer, ierror
      implicit none
      integer buffer(*)
      integer count, datatype, root, comm, ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
c idproc = processor id
      common /mpiparms/ nproc, idproc, srnum, prnum
c MPI constants
      integer MPI_COMM_WORLD, MPI_STATUS_SIZE
      parameter(MPI_COMM_WORLD=0,MPI_STATUS_SIZE=5)
c local data
      integer i, id, status
      dimension status(MPI_STATUS_SIZE)
      ierror = 0
c check for error conditions
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c invalid comm
      elseif (comm.ne.MPI_COMM_WORLD) then
         ierror = 2
c invalid root
      elseif ((root.lt.0).or.(root.ge.nproc)) then
         ierror = 19
      endif
c handle errors
      if (ierror.ne.0) then
         call writerrs('MPI_BCAST: ',ierror)
         return
      endif
c start broadcast
      if (idproc.eq.root) then
         do 10 i = 1, nproc
         id = i - 1
         if (id.ne.root) call MPI_SEND(buffer,count,datatype,id,0,-1,ier
     1ror)
   10    continue
      else
         call MPI_RECV(buffer,count,datatype,root,0,-1,status,ierror)
      endif
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_BARRIER(comm,ierror)
c blocks each process in comm until all processes have called it.
c comm = communicator (only MPI_COMM_WORLD currently supported)
c ierror = error indicator
c input: comm
c output: ierror
      implicit none
      integer comm, ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
c idproc = processor id
      common /mpiparms/ nproc, idproc, srnum, prnum
c MPI constants
      integer MPI_COMM_WORLD, MPI_STATUS_SIZE, MPI_INTEGER
      parameter(MPI_COMM_WORLD=0,MPI_STATUS_SIZE=5,MPI_INTEGER=18)
c local data
      integer ntasks, isync, irync, i, status
      dimension status(MPI_STATUS_SIZE)
      ierror = 0
c check for error conditions
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c invalid comm
      elseif (comm.ne.MPI_COMM_WORLD) then
         ierror = 2
      endif
c handle errors
      if (ierror.ne.0) then
         call writerrs('MPI_BARRIER: ',ierror)
         return
      endif
c begin synchronization
      ntasks = nproc - 1
      isync = -1
      if (idproc.eq.0) then
c processor 0 receives a message from everyone else
         do 10 i = 1, ntasks
         call MPI_RECV(irync,1,MPI_INTEGER,i,0,-2,status,ierror)
         if (irync.ne.isync) write (2,*) 'sync error from proc', i
   10    continue
c then sends an acknowledgment back
         isync = 1
         call MPI_BCAST(isync,1,MPI_INTEGER,0,comm,ierror)
      else
c remaining processors send a message to processor 0
         call MPI_SEND(isync,1,MPI_INTEGER,0,0,-2,ierror)
c then receive an acknowledgement back
         isync = 1
         call MPI_BCAST(irync,1,MPI_INTEGER,0,comm,ierror)
         if (irync.ne.isync) write (2,*) 'rsync error at proc', idproc
      endif
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_REDUCE(sendbuf,recvbuf,count,datatype,op,root,comm,
     1ierror)
c applies a reduction operation to the vector sendbuf over the set of
c processes specified by comm and places the result in recvbuf on root
c sendbuf = address of send buffer
c recvbuf = address of receive buffer
c count = number of elements in send buffer
c datatype = datatype of elements in send buffer
c op = reduce operation (only max, min and sum currently supported)
c root = rank of root process
c comm = communicator (only MPI_COMM_WORLD currently supported)
c ierror = error indicator
c input: sendbuf, count, datatype, op, root, comm
c output: recvbuf, ierror
      implicit none
      integer sendbuf(*), recvbuf(*)
      integer count, datatype, op, root, comm, ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
c idproc = processor id
      common /mpiparms/ nproc, idproc, srnum, prnum
c function declarations
      integer NewPtr
      external NewPtr
c MPI constants
      integer MPI_COMM_WORLD ,MPI_STATUS_SIZE, MPI_SUM
      parameter(MPI_COMM_WORLD=0,MPI_STATUS_SIZE=5,MPI_SUM=2)
      integer MPI_INTEGER, MPI_REAL, MPI_DOUBLE_PRECISION
      parameter(MPI_INTEGER=18,MPI_REAL=19,MPI_DOUBLE_PRECISION=20)
      integer MPI_COMPLEX, MPI_DOUBLE_COMPLEX
      parameter(MPI_COMPLEX=22,MPI_DOUBLE_COMPLEX=23)
c local data
      integer i, j, id, ltmp, loct, nl, lcnt, lsize, tmpbuf, status
      dimension status(MPI_STATUS_SIZE)
      ierror = 0
c check for error conditions
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c invalid comm
      elseif (comm.ne.MPI_COMM_WORLD) then
         ierror = 2
c invalid root
      elseif ((root.lt.0).or.(root.ge.nproc)) then
         ierror = 19
c invalid op
      elseif ((op.lt.0).or.(op.gt.2)) then
         ierror = 20
      endif
c handle errors
      if (ierror.ne.0) then
         call writerrs('MPI_REDUCE: ',ierror)
         return
      endif
c determine size of temporary buffer
      if ((datatype.eq.MPI_INTEGER).or.(datatype.eq.MPI_REAL)) then
         lsize = 4
      elseif (datatype.eq.MPI_DOUBLE_PRECISION) then
         lsize = 8
      elseif ((datatype.eq.MPI_COMPLEX).and.(op.eq.MPI_SUM)) then
         lsize = 8
      elseif ((datatype.eq.MPI_DOUBLE_COMPLEX).and.(op.eq.MPI_SUM)) then
         lsize = 16
c invalid datatype
      else
         ierror = 7
         call writerrs('MPI_REDUCE: ',ierror)
         return
      endif
      if (idproc.eq.root) then
         loct = 0
c initialize by copying from send to receive buffer
         if (datatype.eq.MPI_INTEGER) then
            call iredux(recvbuf,sendbuf,loct,count,-1)
         elseif (datatype.eq.MPI_REAL) then
            call fredux(recvbuf,sendbuf,loct,count,-1)
         elseif (datatype.eq.MPI_DOUBLE_PRECISION) then
            call dredux(recvbuf,sendbuf,loct,count,-1)
         elseif (datatype.eq.MPI_COMPLEX) then
            call fredux(recvbuf,sendbuf,loct,2*count,-1)
         elseif (datatype.eq.MPI_DOUBLE_COMPLEX) then
            call dredux(recvbuf,sendbuf,loct,2*count,-1)
         endif
      else
         loct = 1
      endif
      ltmp = lsize*count
c allocate a nonrelocatable block of memory
      tmpbuf = NewPtr(val4(ltmp))
c memory not available
      if (tmpbuf.eq.0) then
         ierror = 21
         call writerrs('MPI_REDUCE: ',ierror)
         return
      endif
      ltmp = ltmp/lsize
c send messages in groups of ltmp
      nl = (count - 1)/ltmp + 1
      lcnt = ltmp
      lsize = lsize*ltmp/4
      do 20 j = 1, nl
      if (j.eq.nl) lcnt = count - ltmp*(nl - 1)
      if (idproc.eq.root) then
c root receives data from everyone else
         do 10 i = 1, nproc
         id = i - 1
         if (id.ne.root) then
            call MPI_RECV(val4(tmpbuf),lcnt,datatype,id,j,-3,status,ierr
     1or)
c reduce data
            if (datatype.eq.MPI_INTEGER) then
               call iredux(recvbuf,val4(tmpbuf),loct,lcnt,op)
            elseif (datatype.eq.MPI_REAL) then
               call fredux(recvbuf,val4(tmpbuf),loct,lcnt,op)
            elseif (datatype.eq.MPI_DOUBLE_PRECISION) then
               call dredux(recvbuf,val4(tmpbuf),loct,lcnt,op)
            elseif (datatype.eq.MPI_COMPLEX) then
               call fredux(recvbuf,val4(tmpbuf),loct,2*lcnt,op)
            elseif (datatype.eq.MPI_DOUBLE_COMPLEX) then
               call dredux(recvbuf,val4(tmpbuf),loct,2*lcnt,op)
            endif
         endif
   10    continue
         loct = loct + ltmp
      else
c remaining processors send data to root
         call MPI_SEND(sendbuf(loct),lcnt,datatype,root,j,-3,ierror)
         loct = loct + lsize
      endif
   20 continue
c release nonrelocatable memory block
      call DisposePtr(val4(tmpbuf))
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_SCAN(sendbuf,recvbuf,count,datatype,op,comm,ierror)
c performs a parallel prefix reduction on data distributed across a
c group
c sendbuf = address of send buffer
c recvbuf = address of receive buffer
c count = number of elements in send buffer
c datatype = datatype of elements in send buffer
c op = reduce operation (only max, min and sum currently supported)
c comm = communicator (only MPI_COMM_WORLD currently supported)
c ierror = error indicator
c input: sendbuf, count, datatype, op, comm
c output: recvbuf, ierror
      implicit none
      integer sendbuf(*), recvbuf(*)
      integer count, datatype, op, comm, ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
c idproc = processor id
      common /mpiparms/ nproc, idproc, srnum, prnum
c function declarations
      integer NewPtr
      external NewPtr
c MPI constants
      integer MPI_COMM_WORLD ,MPI_STATUS_SIZE, MPI_SUM
      parameter(MPI_COMM_WORLD=0,MPI_STATUS_SIZE=5,MPI_SUM=2)
      integer MPI_INTEGER, MPI_REAL, MPI_DOUBLE_PRECISION
      parameter(MPI_INTEGER=18,MPI_REAL=19,MPI_DOUBLE_PRECISION=20)
      integer MPI_COMPLEX, MPI_DOUBLE_COMPLEX
      parameter(MPI_COMPLEX=22,MPI_DOUBLE_COMPLEX=23)
c local data
      integer status
      integer i, j, root, id, ltmp, loct, nl, lcnt, lsize, tmpbuf
      dimension status(MPI_STATUS_SIZE)
      data root /0/
      ierror = 0
c check for error conditions
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c invalid comm
      elseif (comm.ne.MPI_COMM_WORLD) then
         ierror = 2
c invalid root
      elseif ((root.lt.0).or.(root.ge.nproc)) then
         ierror = 19
c invalid op
      elseif ((op.lt.0).or.(op.gt.2)) then
         ierror = 20
      endif
c handle errors
      if (ierror.ne.0) then
         call writerrs('MPI_SCAN: ',ierror)
         return
      endif
c determine size of temporary buffer
      if ((datatype.eq.MPI_INTEGER).or.(datatype.eq.MPI_REAL)) then
         lsize = 4
      elseif (datatype.eq.MPI_DOUBLE_PRECISION) then
         lsize = 8
      elseif ((datatype.eq.MPI_COMPLEX).and.(op.eq.MPI_SUM)) then
         lsize = 8
      elseif ((datatype.eq.MPI_DOUBLE_COMPLEX).and.(op.eq.MPI_SUM)) then
         lsize = 16
c invalid datatype
      else
         ierror = 7
         call writerrs('MPI_SCAN: ',ierror)
         return
      endif
      if (idproc.eq.root) then
         loct = 0
c initialize by copying from send to receive buffer
         if (datatype.eq.MPI_INTEGER) then
            call iredux(recvbuf,sendbuf,loct,count,-1)
         elseif (datatype.eq.MPI_REAL) then
            call fredux(recvbuf,sendbuf,loct,count,-1)
         elseif (datatype.eq.MPI_DOUBLE_PRECISION) then
            call dredux(recvbuf,sendbuf,loct,count,-1)
         elseif (datatype.eq.MPI_COMPLEX) then
            call fredux(recvbuf,sendbuf,loct,2*count,-1)
         elseif (datatype.eq.MPI_DOUBLE_COMPLEX) then
            call dredux(recvbuf,sendbuf,loct,2*count,-1)
         endif
      else
         loct = 1
      endif
      ltmp = lsize*count
c allocate a nonrelocatable block of memory
      tmpbuf = NewPtr(val4(ltmp))
c memory not available
      if (tmpbuf.eq.0) then
         ierror = 21
         call writerrs('MPI_SCAN: ',ierror)
         return
      endif
      ltmp = ltmp/lsize
c send messages in groups of ltmp
      nl = (count - 1)/ltmp + 1
      lcnt = ltmp
      lsize = lsize*ltmp/4
      do 20 j = 1, nl
      if (j.eq.nl) lcnt = count - ltmp*(nl - 1)
      if (idproc.eq.root) then
c root receives data from everyone else
         do 10 i = 1, nproc
         id = i - 1
         if (id.ne.root) then
            call MPI_RECV(val4(tmpbuf),lcnt,datatype,id,j,-3,status,ierr
     1or)
c reduce data
            if (datatype.eq.MPI_INTEGER) then
               call iredux(recvbuf,val4(tmpbuf),loct,lcnt,op)
            elseif (datatype.eq.MPI_REAL) then
               call fredux(recvbuf,val4(tmpbuf),loct,lcnt,op)
            elseif (datatype.eq.MPI_DOUBLE_PRECISION) then
               call dredux(recvbuf,val4(tmpbuf),loct,lcnt,op)
            elseif (datatype.eq.MPI_COMPLEX) then
               call fredux(recvbuf,val4(tmpbuf),loct,2*lcnt,op)
            elseif (datatype.eq.MPI_DOUBLE_COMPLEX) then
               call dredux(recvbuf,val4(tmpbuf),loct,2*lcnt,op)
            endif
c send partial result data to processor id
            call MPI_SEND(recvbuf(loct+1),lcnt,datatype,id,j+nproc,-3,ie
     1rror)
         endif
   10    continue
         loct = loct + ltmp
      else
c remaining processors send data to root
         call MPI_SEND(sendbuf(loct),lcnt,datatype,root,j,-3,ierror)
c receive partial result data from root
         call MPI_RECV(recvbuf(loct),lcnt,datatype,root,j+nproc,-3,statu
     1s,ierror)
         loct = loct + lsize
      endif
   20 continue
      if (idproc.eq.root) then
c initialize by copying from send to receive buffer
         if (datatype.eq.MPI_INTEGER) then
            call iredux(recvbuf,sendbuf,0,count,-1)
         elseif (datatype.eq.MPI_REAL) then
            call fredux(recvbuf,sendbuf,0,count,-1)
         elseif (datatype.eq.MPI_DOUBLE_PRECISION) then
            call dredux(recvbuf,sendbuf,0,count,-1)
         elseif (datatype.eq.MPI_COMPLEX) then
            call fredux(recvbuf,sendbuf,0,2*count,-1)
         elseif (datatype.eq.MPI_DOUBLE_COMPLEX) then
            call dredux(recvbuf,sendbuf,0,2*count,-1)
         endif
      endif
c release nonrelocatable memory block
      call DisposePtr(val4(tmpbuf))
      return
      end
c-----------------------------------------------------------------------
      subroutine iredux(recvbuf,sendbuf,offset,count,op)
c perform reduction operation for integer types
c recvbuf = address of receive buffer
c sendbuf = address of send buffer
c offset = starting location minus one in receive buffer
c count = number of elements in send buffer
c op = reduce operation (only max, min and sum currently supported)
c input: recvbuf, sendbuf, offset, count, op
c output: recvbuf
      implicit none
      integer offset, count, op
      integer recvbuf(*), sendbuf(*)
c MPI constants
      integer MPI_MAX, MPI_MIN, MPI_SUM
      parameter(MPI_MAX=0,MPI_MIN=1,MPI_SUM=2)
c local data
      integer i
c perform reduction
      if (op.eq.MPI_MAX) then
         do 10 i = 1, count
         recvbuf(i+offset) = max(recvbuf(i+offset),sendbuf(i))
   10    continue
      elseif (op.eq.MPI_MIN) then
         do 20 i = 1, count
         recvbuf(i+offset) = min(recvbuf(i+offset),sendbuf(i))
   20    continue
      elseif (op.eq.MPI_SUM) then
         do 30 i = 1, count
         recvbuf(i+offset) = recvbuf(i+offset) + sendbuf(i)
   30    continue
c copy initial data
      elseif (op.eq.(-1)) then
         do 40 i = 1, count
         recvbuf(i+offset) = sendbuf(i)
   40    continue
      endif
      return
      end
c-----------------------------------------------------------------------
      subroutine fredux(recvbuf,sendbuf,offset,count,op)
c perform reduction operation for real types
c recvbuf = address of receive buffer
c sendbuf = address of send buffer
c offset = starting location minus one in receive buffer
c count = number of elements in send buffer
c op = reduce operation (only max, min and sum currently supported)
c input: recvbuf, sendbuf, offset, count, op
c output: recvbuf
      implicit none
      integer offset, count, op
      real recvbuf(*), sendbuf(*)
c MPI constants
      integer MPI_MAX, MPI_MIN, MPI_SUM
      parameter(MPI_MAX=0,MPI_MIN=1,MPI_SUM=2)
c local data
      integer i
c perform reduction
      if (op.eq.MPI_MAX) then
         do 10 i = 1, count
         recvbuf(i+offset) = max(recvbuf(i+offset),sendbuf(i))
   10    continue
      elseif (op.eq.MPI_MIN) then
         do 20 i = 1, count
         recvbuf(i+offset) = min(recvbuf(i+offset),sendbuf(i))
   20    continue
      elseif (op.eq.MPI_SUM) then
         do 30 i = 1, count
         recvbuf(i+offset) = recvbuf(i+offset) + sendbuf(i)
   30    continue
c copy initial data
      elseif (op.eq.(-1)) then
         do 40 i = 1, count
         recvbuf(i+offset) = sendbuf(i)
   40    continue
      endif
      return
      end
c-----------------------------------------------------------------------
      subroutine dredux(recvbuf,sendbuf,offset,count,op)
c perform reduction operation for double precision types
c recvbuf = address of receive buffer
c sendbuf = address of send buffer
c offset = starting location minus one in receive buffer
c count = number of elements in send buffer
c op = reduce operation (only max, min and sum currently supported)
c input: recvbuf, sendbuf, offset, count, op
c output: recvbuf
      implicit none
      integer offset, count, op
      double precision recvbuf(*), sendbuf(*)
c MPI constants
      integer MPI_MAX, MPI_MIN, MPI_SUM
      parameter(MPI_MAX=0,MPI_MIN=1,MPI_SUM=2)
c local data
      integer i
c perform reduction
      if (op.eq.MPI_MAX) then
         do 10 i = 1, count
         recvbuf(i+offset) = max(recvbuf(i+offset),sendbuf(i))
   10    continue
      elseif (op.eq.MPI_MIN) then
         do 20 i = 1, count
         recvbuf(i+offset) = min(recvbuf(i+offset),sendbuf(i))
   20    continue
      elseif (op.eq.MPI_SUM) then
         do 30 i = 1, count
         recvbuf(i+offset) = recvbuf(i+offset) + sendbuf(i)
   30    continue
c copy initial data
      elseif (op.eq.(-1)) then
         do 40 i = 1, count
         recvbuf(i+offset) = sendbuf(i)
   40    continue
      endif
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_ALLREDUCE(sendbuf,recvbuf,count,datatype,op,comm,ie
     1rror)
c applies a reduction operation to the vector sendbuf over the set of
c processes specified by comm and places result in recvbuf on all nodes
c sendbuf = address of send buffer
c recvbuf = address of receive buffer
c count = number of elements in send buffer
c datatype = datatype of elements in send buffer
c op = reduce operation (only max, min and sum currently supported)
c comm = communicator (only MPI_COMM_WORLD currently supported)
c ierror = error indicator
c input: sendbuf, count, datatype, op, root, comm
c output: recvbuf, ierror
      implicit none
      integer sendbuf(*), recvbuf(*)
      integer count, datatype, op, comm, ierror
c local data
      integer root, ierr
      root = 0
      call MPI_REDUCE(sendbuf,recvbuf,count,datatype,op,root,comm,ierror
     1)
      call MPI_BCAST(recvbuf,count,datatype,root,comm,ierr)
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_GATHER(sendbuf,sendcount,sendtype,recvbuf,recvcount
     1,recvtype,root,comm,ierror)
c collect individual messages from each process in comm at root
c sendbuf = starting address of send buffer
c sendcount = number of elements in send buffer
c sendtype = datatype of send buffer elements
c recvbuf = address of receive buffer
c recvcount = number of elements for any single receive
c recvtype = datatype of recv buffer elements
c root = rank of receiving process
c comm = communicator (only MPI_COMM_WORLD currently supported)
c ierror = error indicator
c input: sendbuf, sendcount, sendtype, recvcount, recvtype, root, comm
c output: recvbuf, ierror
      implicit none
      integer sendbuf(*), recvbuf(*)
      integer sendcount, sendtype, recvcount, recvtype, root, comm
      integer ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
c idproc = processor id
      common /mpiparms/ nproc, idproc, srnum, prnum
c MPI constants
      integer MPI_COMM_WORLD ,MPI_STATUS_SIZE
      parameter(MPI_COMM_WORLD=0,MPI_STATUS_SIZE=5)
      integer MPI_INTEGER, MPI_REAL, MPI_DOUBLE_PRECISION
      parameter(MPI_INTEGER=18,MPI_REAL=19,MPI_DOUBLE_PRECISION=20)
      integer MPI_COMPLEX, MPI_DOUBLE_COMPLEX
      parameter(MPI_COMPLEX=22,MPI_DOUBLE_COMPLEX=23)
c local data
      integer loct, lsize, id, i, j, status
      dimension status(MPI_STATUS_SIZE)
      ierror = 0
c check for error conditions
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c invalid comm
      elseif (comm.ne.MPI_COMM_WORLD) then
         ierror = 2
c invalid root
      elseif ((root.lt.0).or.(root.ge.nproc)) then
         ierror = 19
c invalid count
      elseif (sendcount.lt.0) then
         ierror = 3
      endif
c handle errors
      if (ierror.ne.0) then
         call writerrs('MPI_GATHER: ',ierror)
         return
      endif
c root receives data
      if (idproc.eq.root) then
c invalid count
         if (recvcount.lt.0) ierror = 3
c determine size of data to be sent
         if ((sendtype.eq.MPI_INTEGER).or.(sendtype.eq.MPI_REAL)) then
            loct = sendcount
         elseif ((sendtype.eq.MPI_DOUBLE_PRECISION).or.(sendtype.eq.MPI_
     1COMPLEX)) then
            loct = 2*sendcount
         elseif (sendtype.eq.MPI_DOUBLE_COMPLEX) then
            loct = 4*sendcount
c invalid datatype
         else
            loct = 0
            ierror = 7
         endif
c determine size of data to be received
         if ((recvtype.eq.MPI_INTEGER).or.(recvtype.eq.MPI_REAL)) then
            lsize = recvcount
         elseif ((recvtype.eq.MPI_DOUBLE_PRECISION).or.(recvtype.eq.MPI_
     1COMPLEX)) then
            lsize = 2*recvcount
         elseif (recvtype.eq.MPI_DOUBLE_COMPLEX) then
            lsize = 4*recvcount
c invalid datatype
         else
            lsize = 0
            ierror = 7
         endif
c unequal message length error
         if (loct.ne.lsize) then
            write (2,*) 'Unequal message length, send/receive bytes = ',
     1loct, lsize
            ierror = 22
         endif
c handle count, datatype and length errors
         if (ierror.ne.0) then
            call writerrs('MPI_GATHER: ',ierror)
            return
         endif
         do 20 i = 1, nproc
         id = i - 1
         loct = lsize*id
c root copies its own data directly
         if (id.eq.root) then
            do 10 j = 1, lsize
            recvbuf(j+loct) = sendbuf(j)
   10       continue
c otherwise, root receives data from other processors
         else
            call MPI_RECV(recvbuf(loct+1),recvcount,recvtype,id,1,-4,sta
     1tus,ierror)
         endif
   20    continue
c processors other than root send data to root
      else
         call MPI_SEND(sendbuf,sendcount,sendtype,root,1,-4,ierror)
      endif
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_ALLGATHER(sendbuf,sendcount,sendtype,recvbuf,recvco
     1unt,recvtype,comm,ierror)
c gather individual messages from each process in comm and distribute
c the resulting message to each process.
c sendbuf = starting address of send buffer
c sendcount = number of elements in send buffer
c sendtype = datatype of send buffer elements
c recvbuf = address of receive buffer
c recvcount = number of elements for any process
c recvtype = datatype of receive buffer elements
c comm = communicator (only MPI_COMM_WORLD currently supported)
c ierror = error indicator
c input: sendbuf, sendcount, sendtype, recvcount, recvtype, comm
c output: recvbuf, ierror
      implicit none
      integer sendbuf(*), recvbuf(*)
      integer sendcount, sendtype, recvcount, recvtype, comm, ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
      common /mpiparms/ nproc, idproc, srnum, prnum
c local data
      integer root, ierr
      root = 0
      call MPI_GATHER(sendbuf,sendcount,sendtype,recvbuf,recvcount,recvt
     1ype,root,comm,ierror)
      call MPI_BCAST(recvbuf,nproc*recvcount,recvtype,root,comm,ierr)
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_SCATTER(sendbuf,sendcount,sendtype,recvbuf,recvcoun
     1t,recvtype,root,comm,ierror)
c distribute individual messages from root to each process in comm
c sendbuf = starting address of send buffer
c sendcount = number of elements in send buffer
c sendtype = datatype of send buffer elements
c recvbuf = address of receive buffer
c recvcount = number of elements for any single receive
c recvtype = datatype of recv buffer elements
c root = rank of sending process
c comm = communicator (only MPI_COMM_WORLD currently supported)
c ierror = error indicator
c input: sendbuf, sendcount, sendtype, recvcount, recvtype, root, comm
c output: recvbuf, ierror
      implicit none
      integer sendbuf(*), recvbuf(*)
      integer sendcount, sendtype, recvcount, recvtype, root, comm
      integer ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
c idproc = processor id
      common /mpiparms/ nproc, idproc, srnum, prnum
c MPI constants
      integer MPI_COMM_WORLD ,MPI_STATUS_SIZE
      parameter(MPI_COMM_WORLD=0,MPI_STATUS_SIZE=5)
      integer MPI_INTEGER, MPI_REAL, MPI_DOUBLE_PRECISION
      parameter(MPI_INTEGER=18,MPI_REAL=19,MPI_DOUBLE_PRECISION=20)
      integer MPI_COMPLEX, MPI_DOUBLE_COMPLEX
      parameter(MPI_COMPLEX=22,MPI_DOUBLE_COMPLEX=23)
c local data
      integer lsize, loct, id, i, j, status
      dimension status(MPI_STATUS_SIZE)
      ierror = 0
c check for error conditions
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c invalid comm
      elseif (comm.ne.MPI_COMM_WORLD) then
         ierror = 2
c invalid root
      elseif ((root.lt.0).or.(root.ge.nproc)) then
         ierror = 19
c invalid counts
      elseif (recvcount.lt.0) then
         ierror = 3
      endif
c handle errors
      if (ierror.ne.0) then
         call writerrs('MPI_SCATTER: ',ierror)
         return
      endif
c root sends data
      if (idproc.eq.root) then
c invalid counts
         if (sendcount.lt.0) ierror = 3
c determine size of data to be sent
         if ((sendtype.eq.MPI_INTEGER).or.(sendtype.eq.MPI_REAL)) then
            lsize = sendcount
         elseif ((sendtype.eq.MPI_DOUBLE_PRECISION).or.(sendtype.eq.MPI_
     1COMPLEX)) then
            lsize = 2*sendcount
         elseif (sendtype.eq.MPI_DOUBLE_COMPLEX) then
            lsize = 4*sendcount
c invalid datatype
         else
            lsize = 0
            ierror = 7
         endif
c determine size of data to be received
         if ((recvtype.eq.MPI_INTEGER).or.(recvtype.eq.MPI_REAL)) then
            loct = recvcount
         elseif ((recvtype.eq.MPI_DOUBLE_PRECISION).or.(recvtype.eq.MPI_
     1COMPLEX)) then
            loct = 2*recvcount
         elseif (recvtype.eq.MPI_DOUBLE_COMPLEX) then
            loct = 4*recvcount
c invalid datatype
         else
            loct = 0
            ierror = 7
         endif
c unequal message length error
         if (loct.ne.lsize) then
            write (2,*) 'Unequal message length, send/receive bytes = ',
     1lsize, loct
            ierror = 22
         endif
c handle count, datatype and length errors
         if (ierror.ne.0) then
            call writerrs('MPI_SCATTER: ',ierror)
            return
         endif
         do 20 i = 1, nproc
         id = i - 1
         loct = lsize*id
c root copies its own data directly
         if (id.eq.root) then
            do 10 j = 1, lsize
            recvbuf(j) = sendbuf(j+loct)
   10       continue
c otherwise, root sends data to other processors
         else
            call MPI_SEND(sendbuf(loct+1),sendcount,sendtype,id,1,-5,ier
     1ror)
         endif
   20    continue
c processors other than root receive data from root
      else
         call MPI_RECV(recvbuf,recvcount,recvtype,root,1,-5,status,ierro
     1r)
      endif
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_ALLTOALL(sendbuf,sendcount,sendtype,recvbuf,recvcou
     1nt,recvtype,comm,ierror)
c send a distinct message from each process to every process
c sendbuf = starting address of send buffer
c sendcount = number of elements in send buffer
c sendtype = datatype of send buffer elements
c recvbuf = address of receive buffer
c recvcount = number of elements for any single receive
c recvtype = datatype of recv buffer elements
c comm = communicator (only MPI_COMM_WORLD currently supported)
c ierror = error indicator
c input: sendbuf, sendcount, sendtype, recvcount, recvtype, comm
c output: recvbuf, ierror
      implicit none
      integer sendbuf(*), recvbuf(*)
      integer sendcount, sendtype, recvcount, recvtype, comm
      integer ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
c idproc = processor id
      common /mpiparms/ nproc, idproc, srnum, prnum
c MPI constants
      integer MPI_COMM_WORLD ,MPI_STATUS_SIZE
      parameter(MPI_COMM_WORLD=0,MPI_STATUS_SIZE=5)
      integer MPI_INTEGER, MPI_REAL, MPI_DOUBLE_PRECISION
      parameter(MPI_INTEGER=18,MPI_REAL=19,MPI_DOUBLE_PRECISION=20)
      integer MPI_COMPLEX, MPI_DOUBLE_COMPLEX
      parameter(MPI_COMPLEX=22,MPI_DOUBLE_COMPLEX=23)
c local data
      integer loct, lsize, id, i, j, request, status
      dimension status(MPI_STATUS_SIZE)
      ierror = 0
c check for error conditions
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c invalid comm
      elseif (comm.ne.MPI_COMM_WORLD) then
         ierror = 2
c invalid counts
      elseif ((sendcount.lt.0).or.(recvcount.lt.0)) then
         ierror = 3
      endif
c handle errors
      if (ierror.ne.0) then
         call writerrs('MPI_ALLTOALL: ',ierror)
         return
      endif
c determine size of data to be sent
      if ((sendtype.eq.MPI_INTEGER).or.(sendtype.eq.MPI_REAL)) then
         loct = sendcount
      elseif ((sendtype.eq.MPI_DOUBLE_PRECISION).or.(sendtype.eq.MPI_
     1COMPLEX)) then
         loct = 2*sendcount
      elseif (sendtype.eq.MPI_DOUBLE_COMPLEX) then
         loct = 4*sendcount
c invalid datatype
      else
         loct = 0
         ierror = 7
      endif
c determine size of data to be received
      if ((recvtype.eq.MPI_INTEGER).or.(recvtype.eq.MPI_REAL)) then
         lsize = recvcount
      elseif ((recvtype.eq.MPI_DOUBLE_PRECISION).or.(recvtype.eq.MPI_
     1COMPLEX)) then
         lsize = 2*recvcount
      elseif (recvtype.eq.MPI_DOUBLE_COMPLEX) then
         lsize = 4*recvcount
c invalid datatype
      else
         lsize = 0
         ierror = 7
      endif
c unequal message length error
      if (loct.ne.lsize) then
         write (2,*) 'Unequal message length, send/receive bytes = ',
     1loct, lsize
         ierror = 22
      endif
c handle count, datatype and length errors
      if (ierror.ne.0) then
         call writerrs('MPI_ALLTOALL: ',ierror)
         return
      endif
      do 20 i = 1, nproc
      id = i - idproc - 1
      if (id.lt.0) id = id + nproc
      loct = lsize*id
c each node copies its own data directly
      if (idproc.eq.id) then
         do 10 j = 1, lsize
         recvbuf(j+loct) = sendbuf(j+loct)
   10    continue
c otherwise, each node receives data from other nodes
      else
         call MPI_IRECV(recvbuf(loct+1),recvcount,recvtype,id,i,-6,reque
     1st,ierror)
         call MPI_SEND(sendbuf(loct+1),sendcount,sendtype,id,i,-6,ierror
     1)
         call MPI_WAIT(request,status,ierror)
      endif
   20 continue
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_GATHERV(sendbuf,sendcount,sendtype,recvbuf,recvcoun
     1ts,displs,recvtype,root,comm,ierror)
c collect individual messages from each process in comm at root
c messages can have different sizes and displacements
c sendbuf = starting address of send buffer
c sendcount = number of elements in send buffer
c sendtype = datatype of send buffer elements
c recvbuf = address of receive buffer
c recvcounts = integer array
c displs = integer array of displacements
c recvtype = datatype of recv buffer elements
c root = rank of receiving process
c comm = communicator (only MPI_COMM_WORLD currently supported)
c ierror = error indicator
c input: sendbuf, sendcount, sendtype, recvcounts, displs, recvtype
c input: root, comm
c output: recvbuf, ierror
      implicit none
      integer sendbuf(*), recvbuf(*), recvcounts(*), displs(*)
      integer sendcount, sendtype, recvtype, root, comm
      integer ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
c idproc = processor id
      common /mpiparms/ nproc, idproc, srnum, prnum
c MPI constants
      integer MPI_COMM_WORLD ,MPI_STATUS_SIZE
      parameter(MPI_COMM_WORLD=0,MPI_STATUS_SIZE=5)
      integer MPI_INTEGER, MPI_REAL, MPI_DOUBLE_PRECISION
      parameter(MPI_INTEGER=18,MPI_REAL=19,MPI_DOUBLE_PRECISION=20)
      integer MPI_COMPLEX, MPI_DOUBLE_COMPLEX
      parameter(MPI_COMPLEX=22,MPI_DOUBLE_COMPLEX=23)
c local data
      integer loct, lsize, id, i, j, status
      dimension status(MPI_STATUS_SIZE)
      ierror = 0
c check for error conditions
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c invalid comm
      elseif (comm.ne.MPI_COMM_WORLD) then
         ierror = 2
c invalid root
      elseif ((root.lt.0).or.(root.ge.nproc)) then
         ierror = 19
c invalid count
      elseif (sendcount.lt.0) then
         ierror = 3
      endif
c handle errors
      if (ierror.ne.0) then
         call writerrs('MPI_GATHERV: ',ierror)
         return
      endif
c root receives data
      if (idproc.eq.root) then
c invalid counts
         do 10 i = 1, nproc
         if (recvcounts(i).lt.0) ierror = 3
   10    continue
c determine size of data to be sent
         if ((sendtype.eq.MPI_INTEGER).or.(sendtype.eq.MPI_REAL)) then
            loct = 1
         elseif ((sendtype.eq.MPI_DOUBLE_PRECISION).or.(sendtype.eq.MPI_
     1COMPLEX)) then
            loct = 2
         elseif (sendtype.eq.MPI_DOUBLE_COMPLEX) then
            loct = 4
c invalid datatype
         else
            ierror = 7
         endif
c determine size of data to be received
         if ((recvtype.eq.MPI_INTEGER).or.(recvtype.eq.MPI_REAL)) then
            lsize = 1
         elseif ((recvtype.eq.MPI_DOUBLE_PRECISION).or.(recvtype.eq.MPI_
     1COMPLEX)) then
            lsize = 2
         elseif (recvtype.eq.MPI_DOUBLE_COMPLEX) then
            lsize = 4
c invalid datatype
         else
            ierror = 7
         endif
c unequal message length error
         id = lsize*recvcounts(idproc+1)
         if ((ierror.eq.0).and.(loct*sendcount.ne.id)) then
            write (2,*) 'Unequal self message, send/receive bytes = ',
     1loct*sendcount, id
            ierror = 22
         endif
c handle count and datatype errors
         if (ierror.ne.0) then
            call writerrs('MPI_GATHERV: ',ierror)
            return
         endif
         do 30 i = 1, nproc
         id = i - 1
         loct = lsize*displs(i)
c root copies its own data directly
         if (id.eq.root) then
            do 20 j = 1, lsize*recvcounts(i)
            recvbuf(j+loct) = sendbuf(j)
   20       continue
c otherwise, root receives data from other processors
         else
            call MPI_RECV(recvbuf(loct+1),recvcounts(i),recvtype,id,1,-7
     1,status,ierror)
         endif
   30    continue
c processors other than root send data to root
      else
         call MPI_SEND(sendbuf,sendcount,sendtype,root,1,-7,ierror)
      endif
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_ALLGATHERV(sendbuf,sendcount,sendtype,recvbuf,recvc
     1ounts,displs,recvtype,comm,ierror)
c gather individual messages from each process in comm and distribute
c the resulting message to each process.
c messages can have different sizes and displacements
c sendbuf = starting address of send buffer
c sendcount = number of elements in send buffer
c sendtype = datatype of send buffer elements
c recvbuf = address of receive buffer
c recvcounts = integer array
c displs = integer array of displacements
c recvtype = datatype of receive buffer elements
c comm = communicator (only MPI_COMM_WORLD currently supported)
c ierror = error indicator
c input: sendbuf, sendcount, sendtype, recvcounts, displs, recvtype
c input: comm
c output: recvbuf, ierror
      implicit none
      integer sendbuf(*), recvbuf(*), recvcounts(*), displs(*)
      integer sendcount, sendtype, recvtype, comm, ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
      common /mpiparms/ nproc, idproc, srnum, prnum
c local data
      integer root, i, ierr
      ierror = 0
      do 10 i = 1, nproc
      root = i - 1
      call MPI_GATHERV(sendbuf,sendcount,sendtype,recvbuf,recvcounts,dis
     1pls,recvtype,root,comm,ierr)
      if (ierr.ne.0) ierror = ierr
   10 continue
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_SCATTERV(sendbuf,sendcounts,displs,sendtype,recvbuf
     1,recvcount,recvtype,root,comm,ierror)
c distribute individual messages from root to each process in comm
c messages can have different sizes and displacements
c sendbuf = starting address of send buffer
c sendcounts = integer array
c displs = integer array of displacements
c sendtype = datatype of send buffer elements
c recvbuf = address of receive buffer
c recvcount = number of elements for any single receive
c recvtype = datatype of recv buffer elements
c root = rank of sending process
c comm = communicator (only MPI_COMM_WORLD currently supported)
c ierror = error indicator
c input: sendbuf, sendcounts, displs, sendtype, recvcount, recvtype
c input: root, comm
c output: recvbuf, ierror
      implicit none
      integer sendbuf(*), recvbuf(*), sendcounts(*), displs(*)
      integer sendtype, recvcount, recvtype, root, comm
      integer ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
c idproc = processor id
      common /mpiparms/ nproc, idproc, srnum, prnum
c MPI constants
      integer MPI_COMM_WORLD ,MPI_STATUS_SIZE
      parameter(MPI_COMM_WORLD=0,MPI_STATUS_SIZE=5)
      integer MPI_INTEGER, MPI_REAL, MPI_DOUBLE_PRECISION
      parameter(MPI_INTEGER=18,MPI_REAL=19,MPI_DOUBLE_PRECISION=20)
      integer MPI_COMPLEX, MPI_DOUBLE_COMPLEX
      parameter(MPI_COMPLEX=22,MPI_DOUBLE_COMPLEX=23)
c local data
      integer lsize, loct, id, i, j, status
      dimension status(MPI_STATUS_SIZE)
      ierror = 0
c check for error conditions
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c invalid comm
      elseif (comm.ne.MPI_COMM_WORLD) then
         ierror = 2
c invalid root
      elseif ((root.lt.0).or.(root.ge.nproc)) then
         ierror = 19
c invalid counts
      elseif (recvcount.lt.0) then
         ierror = 3
      endif
c handle errors
      if (ierror.ne.0) then
         call writerrs('MPI_SCATTERV: ',ierror)
         return
      endif
c root sends data
      if (idproc.eq.root) then
c invalid counts
         do 10 i = 1, nproc
         if (sendcounts(i).lt.0) ierror = 3
   10    continue
c determine size of data to be sent
         if ((sendtype.eq.MPI_INTEGER).or.(sendtype.eq.MPI_REAL)) then
            lsize = 1
      elseif ((sendtype.eq.MPI_DOUBLE_PRECISION).or.(sendtype.eq.MPI_
     1COMPLEX)) then
            lsize = 2
         elseif (sendtype.eq.MPI_DOUBLE_COMPLEX) then
            lsize = 4
c invalid datatype
         else
            ierror = 7
         endif
c determine size of data to be received
         if ((recvtype.eq.MPI_INTEGER).or.(recvtype.eq.MPI_REAL)) then
            loct = 1
         elseif ((recvtype.eq.MPI_DOUBLE_PRECISION).or.(recvtype.eq.MPI_
     1COMPLEX)) then
            loct = 2
         elseif (recvtype.eq.MPI_DOUBLE_COMPLEX) then
            loct = 4
c invalid datatype
         else
            ierror = 7
         endif
c unequal message length error
         id = lsize*sendcounts(idproc+1)
         if ((ierror.eq.0).and.(loct*recvcount.ne.id)) then
            write (2,*) 'Unequal self message, send/receive bytes = ',
     1id, loct*recvcount
            ierror = 22
         endif
c handle count and datatype errors
         if (ierror.ne.0) then
            call writerrs('MPI_SCATTERV: ',ierror)
            return
         endif
         do 30 i = 1, nproc
         id = i - 1
         loct = lsize*displs(i)
c root copies its own data directly
         if (id.eq.root) then
            do 20 j = 1, lsize*sendcounts(i)
            recvbuf(j) = sendbuf(j+loct)
   20       continue
c otherwise, root sends data to other processors
         else
            call MPI_SEND(sendbuf(loct+1),sendcounts(i),sendtype,id,1,-8
     1,ierror)
         endif
   30    continue
c processors other than root receive data from root
      else
         call MPI_RECV(recvbuf,recvcount,recvtype,root,1,-8,status,ierro
     1r)
      endif
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_ALLTOALLV(sendbuf,sendcounts,sdispls,sendtype,recvb
     1uf,recvcounts,rdispls,recvtype,comm,ierror)
c send a distinct message from each process to every process
c messages can have different sizes and displacements
c sendbuf = starting address of send buffer
c sendcounts = integer array
c sdispls = integer array of send displacements
c sendtype = datatype of send buffer elements
c recvbuf = address of receive buffer
c recvcounts = integer array
c rdispls = integer array of receive displacements
c recvtype = datatype of recv buffer elements
c comm = communicator (only MPI_COMM_WORLD currently supported)
c ierror = error indicator
c input: sendbuf, sendcount, sendtype, recvcount, recvtype, comm
c output: recvbuf, ierror
      implicit none
      integer sendbuf(*), recvbuf(*)
      integer sendcounts(*), sdispls(*), recvcounts(*), rdispls(*)
      integer sendtype, recvtype, comm
      integer ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
c idproc = processor id
      common /mpiparms/ nproc, idproc, srnum, prnum
c MPI constants
      integer MPI_COMM_WORLD ,MPI_STATUS_SIZE
      parameter(MPI_COMM_WORLD=0,MPI_STATUS_SIZE=5)
      integer MPI_INTEGER, MPI_REAL, MPI_DOUBLE_PRECISION
      parameter(MPI_INTEGER=18,MPI_REAL=19,MPI_DOUBLE_PRECISION=20)
      integer MPI_COMPLEX, MPI_DOUBLE_COMPLEX
      parameter(MPI_COMPLEX=22,MPI_DOUBLE_COMPLEX=23)
c local data
      integer locs, msize, loct, lsize, id, ld, i, j, request, status
      dimension status(MPI_STATUS_SIZE)
      ierror = 0
c check for error conditions
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c invalid comm
      elseif (comm.ne.MPI_COMM_WORLD) then
         ierror = 2
      endif
c invalid counts
      do 10 i = 1, nproc
      if ((sendcounts(i).lt.0).or.(recvcounts(i).lt.0)) ierror = 3
   10 continue
c handle errors
      if (ierror.ne.0) then
         call writerrs('MPI_ALLTOALLV: ',ierror)
         return
      endif
c determine size of data to be sent
      if ((sendtype.eq.MPI_INTEGER).or.(sendtype.eq.MPI_REAL)) then
         msize = 1
      elseif ((sendtype.eq.MPI_DOUBLE_PRECISION).or.(sendtype.eq.MPI_
     1COMPLEX)) then
         msize = 2
      elseif (sendtype.eq.MPI_DOUBLE_COMPLEX) then
         msize = 4
c invalid datatype
      else
         ierror = 7
      endif
c determine size of data to be received
      if ((recvtype.eq.MPI_INTEGER).or.(recvtype.eq.MPI_REAL)) then
         lsize = 1
      elseif ((recvtype.eq.MPI_DOUBLE_PRECISION).or.(recvtype.eq.MPI_
     1COMPLEX)) then
         lsize = 2
      elseif (recvtype.eq.MPI_DOUBLE_COMPLEX) then
         lsize = 4
c invalid datatype
      else
         ierror = 7
      endif
c unequal message length error
      id = msize*sendcounts(idproc+1)
      ld = lsize*recvcounts(idproc+1)
      if ((ierror.eq.0).and.(id.ne.ld)) then
         write (2,*) 'Unequal self message length, send/receive bytes=',
     1id, ld
         ierror = 22
      endif
c handle count and datatype errors
      if (ierror.ne.0) then
         call writerrs('MPI_ALLTOALLV: ',ierror)
         return
      endif
      do 30 i = 1, nproc
      id = i - idproc - 1
      if (id.lt.0) id = id + nproc
      ld = id + 1
      locs = msize*sdispls(ld)
      loct = lsize*rdispls(ld)
c each node copies its own data directly
      if (idproc.eq.id) then
         do 20 j = 1, lsize*recvcounts(ld)
         recvbuf(j+loct) = sendbuf(j+locs)
   20    continue
c otherwise, each node receives data from other nodes
      else
         call MPI_IRECV(recvbuf(loct+1),recvcounts(ld),recvtype,id,i,-9,
     1request,ierror)
         call MPI_SEND(sendbuf(locs+1),sendcounts(ld),sendtype,id,i,-9,i
     1error)
         call MPI_WAIT(request,status,ierror)
      endif
   30 continue
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_REDUCE_SCATTER(sendbuf,recvbuf,recvcounts,datatype,
     1op,comm,ierror)
c applies a reduction operation to the vector sendbuf over the set of
c processes specified by comm and scatters the result according to the
c values in recvcounts
c sendbuf = starting address of send buffer
c recvbuf = starting address of receive buffer
c recvcounts = integer array
c datatype = datatype of elements in input buffer
c op = reduce operation (only max, min and sum currently supported)
c comm = communicator (only MPI_COMM_WORLD currently supported)
c ierror = error indicator
c input: sendbuf, recvcounts, datatype, op, comm
c output: recvbuf, ierror
      implicit none
      integer sendbuf(*), recvbuf(*), recvcounts(*)
      integer datatype, op, comm, ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
c idproc = processor id
      common /mpiparms/ nproc, idproc, srnum, prnum
c function declarations
      integer NewPtr
      external NewPtr
c MPI constants
      integer MPI_SUM
      parameter(MPI_SUM=2)
      integer MPI_INTEGER, MPI_REAL, MPI_DOUBLE_PRECISION
      parameter(MPI_INTEGER=18,MPI_REAL=19,MPI_DOUBLE_PRECISION=20)
      integer MPI_COMPLEX, MPI_DOUBLE_COMPLEX
      parameter(MPI_COMPLEX=22,MPI_DOUBLE_COMPLEX=23)
c local data
      integer root, count, lsize, ltmp, i, displs, tmpbuf
      data root /0/
      ltmp = 4*nproc
c allocate a nonrelocatable block of memory
      displs = NewPtr(val4(ltmp))
c memory not available
      if (displs.eq.0) then
         ierror = 21
         call writerrs('MPI_REDUCE_SCATTER: ',ierror)
         return
      endif
      count = 0
c determines overall count and displacements
      do 10 i = 1, nproc
      long(displs+4*(i-1)) = count
      count = count + recvcounts(i)
   10 continue
c determine size of temporary buffer
      if ((datatype.eq.MPI_INTEGER).or.(datatype.eq.MPI_REAL)) then
         lsize = 4
      elseif (datatype.eq.MPI_DOUBLE_PRECISION) then
         lsize = 8
      elseif ((datatype.eq.MPI_COMPLEX).and.(op.eq.MPI_SUM)) then
         lsize = 8
      elseif ((datatype.eq.MPI_DOUBLE_COMPLEX).and.(op.eq.MPI_SUM)) then
         lsize = 16
      else
         ierror = 7
         call writerrs('MPI_REDUCE_SCATTER: ',ierror)
         return
      endif
      ltmp = lsize*count
c allocate a nonrelocatable block of memory
      tmpbuf = NewPtr(val4(ltmp))
c memory not available
      if (tmpbuf.eq.0) then
         ierror = 21
         call writerrs('MPI_REDUCE_SCATTER: ',ierror)
         return
      endif
      call MPI_REDUCE(sendbuf,val4(tmpbuf),count,datatype,op,root,comm,i
     1error)
      call MPI_SCATTERV(val4(tmpbuf),recvcounts,val4(displs),datatype,re
     1cvbuf,recvcounts(idproc+1),datatype,root,comm,ierror)
c release nonrelocatable memory block
      call DisposePtr(val4(tmpbuf))
c release nonrelocatable memory block
      call DisposePtr(val4(displs))
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_ABORT(comm,errorcode,ierror)
c force all tasks on an MPI environment to terminate
c comm = communicator (only MPI_COMM_WORLD currently supported)
c errorcode = error code to return to invoking environment
c ierror = error indicator
c input: comm, errorcode
c output: ierror
      implicit none
      integer comm, errorcode, ierror
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = number of real or virtual processors obtained
      common /mpiparms/ nproc, idproc, srnum, prnum
c MPI not initialized
      if (nproc.le.0) then
         ierror = 1
c this is just a temporary patch, have not yet notified everyone else
      else
         call MPI_FINALIZE(ierror)
      endif
      return
      end
c-----------------------------------------------------------------------
      function MPI_WTIME()
c return an elapsed time on the calling processor in seconds
      implicit none
      double precision MPI_WTIME, tick
      parameter(tick=1.0d0/60.0d0)
      integer TickCount
      external TickCount
c get current system tick count
      MPI_WTIME = dble(TickCount())*tick
      return
      end
c-----------------------------------------------------------------------
      function MPI_WTICK()
c return the resolution of MPI_WTIME in seconds
      implicit none
      double precision MPI_WTICK, tick
      parameter(tick=1.0d0/60.0d0)
      MPI_WTICK = tick
      return
      end
c-----------------------------------------------------------------------
      subroutine MPI_TYPE_EXTENT(datatype,extent,ierror)
c returns the size of a datatype
c datatype = datatype
c extent = datatype extent
c ierror = error indicator
c input: dataype
c output: extent, ierror
      implicit none
      integer datatype, extent, ierror
c MPI constants
      integer MPI_INTEGER, MPI_REAL, MPI_DOUBLE_PRECISION
      parameter(MPI_INTEGER=18,MPI_REAL=19,MPI_DOUBLE_PRECISION=20)
      integer MPI_COMPLEX, MPI_DOUBLE_COMPLEX, MPI_BYTE
      parameter(MPI_COMPLEX=22,MPI_DOUBLE_COMPLEX=23,MPI_BYTE=2)
      ierror = 0
c find size of datatype
      if ((datatype.eq.MPI_INTEGER).or.(datatype.eq.MPI_REAL)) then
         extent = 4
      elseif ((datatype.eq.MPI_DOUBLE_PRECISION).or.(datatype.eq.MPI_COM
     1PLEX)) then
         extent = 8
      elseif (datatype.eq.MPI_DOUBLE_COMPLEX) then
         extent = 16
      elseif (datatype.eq.MPI_BYTE) then
         extent = 1
c invalid datatype
      else
         ierror = 7
         write (2,*) 'MPI_TYPE_EXTENT: Invalid datatype'
      endif
      return
      end
c-----------------------------------------------------------------------
      subroutine writerrs(source,ierror)
c this subroutine writes out error descriptions from error codes
c source = source subroutine of error message
c ierror = error indicator
c input: source, ierror
      implicit none
      character*(*) source
      integer ierror
c MPI constants
      integer MPI_COMM_WORLD
      parameter(MPI_COMM_WORLD=0)
c local data
      integer ierr
      logical fatal
      fatal = .true.
c check error code and print corresponding message
      if (ierror.eq.1) then
         write (2,*) source, 'MPI not initialized'
      elseif (ierror.eq.2) then
         write (2,*) source, 'Invalid Communicator'
      elseif (ierror.eq.3) then
         write (2,*) source,'Invalid count'
      elseif (ierror.eq.4) then
         write (2,*) source, 'Invalid destination'
      elseif (ierror.eq.5) then
         write (2,*) source, 'Invalid source'
      elseif (ierror.eq.6) then
         write (2,*) source, 'Invalid tag'
      elseif (ierror.eq.7) then
         write (2,*) source, 'Invalid datatype'
      elseif (ierror.eq.12) then
         write (2,*) source, 'Incomplete read'
      elseif (ierror.eq.16) then
         write (2,*) source, 'Invalid request handle'
      elseif (ierror.eq.18) then
         write (2,*) source, 'Mismatched dataype'
      elseif (ierror.eq.19) then
         write (2,*) source, 'Invalid root'
      elseif (ierror.eq.20) then
         write (2,*) source, 'Invalid operation'
      elseif (ierror.eq.21) then
         write (2,*) source, 'Unable to allocate memory'
c unlisted error code
      else
         write (2,*) source, 'Error code = ', ierror
      endif
c abort if error is fatal
      if (fatal) then
         end file 2
         backspace 2
         call MPI_ABORT(MPI_COMM_WORLD,ierror,ierr)
         close(unit=2)
         stop
      endif
      return
      end
c-----------------------------------------------------------------------
      subroutine messwin(nvp)
c this subroutine creates a window for showing MPI message status
c nvp = number of real or virtual processors
c input argument: nvp
      implicit none
      integer nvp
c function declarations
      integer*2 TextWidth
      integer GetMainDevice, NewWindow
      external GetMainDevice, NewWindow, TextWidth
c declare internal mpi common block
      integer nproc, idproc
      integer MAXS, srnum
      parameter(MAXS=16)
      dimension srnum(MAXS+1) 
      integer*2 prnum
c nproc = num