!------------------------------------------------------------------------ ! FILE: pt2pt-combo.f ! ! PROGRAM DESCRIPTION: ! Based on least-squares-pt2pt.f, this program combines all the ! solutions to problems 2, 3, and 4 of Exercise 4 in the ! "MPI Point to Point Communication" lab. ! ! USAGE: Run with number of processes that equals to the power of 2. ! Tested with 2, 4, and 8 processes. ! ! AUTHOR: Dora Abdullah (MPI version, 11/96) !------------------------------------------------------------------------ program linear implicit none include 'mpif.h' real(8), allocatable :: x(:),y(:) real(8) mySUMx, mySUMy, mySUMxy, mySUMxx, SUMx, SUMy, SUMxy, & SUMxx, SUMres, res, slope, y_intercept, y_estimate real value integer i,j,k,n,myid,numprocs,maxprocs parameter (maxprocs = 10) integer ierr,request,status(MPI_STATUS_SIZE), & irequest(4),istatus_array(MPI_STATUS_SIZE,4), & nrequest(maxprocs), & nstatus_array(MPI_STATUS_SIZE,maxprocs), & naverage,nremain,mypoints,ishift, & igroup,npower,isendto,irecvfrom open (unit=10, file='xydata') call MPI_INIT(ierr) call MPI_COMM_RANK(MPI_COMM_WORLD, myid, ierr) call MPI_COMM_SIZE(MPI_COMM_WORLD, numprocs, ierr) ! ---------------------------------------------------------- ! Step 1: Process 0 reads data and sends the value of n ! ---------------------------------------------------------- if (myid.eq.0) then write(6,*) write(6,*) 'Number of processes used: ', numprocs write(6,*) '-------------------------------------' write(6,*) 'The x coordinates on worker processes:' ! -------------------------------------------------------- ! This call is used to achieve a consistent output format call new_sleep (%val(3)) ! -------------------------------------------------------- read(10,*) n allocate (x(n),y(n)) do i=1,n read(10,*) x(i),y(i) end do write(6,*) do i = 1,numprocs-1 call MPI_ISEND(n,1,MPI_INTEGER,i,10,MPI_COMM_WORLD, & nrequest(i),ierr) ! ---------------------------------------------------- ! We don't need to wait for this ISEND before the end ! of the program, because we don't need to write to n ! ---------------------------------------------------- end do else call MPI_IRECV(n,1,MPI_INTEGER,0,10,MPI_COMM_WORLD, & request,ierr) call MPI_WAIT(request,status,ierr) allocate (x(n),y(n)) end if ! ---------------------------------------------------------- naverage = n/numprocs nremain = mod(n,numprocs) ! ---------------------------------------------------------- ! Step 2: Process 0 sends subsets of x and y ! ---------------------------------------------------------- if (myid.eq.0) then do i = 1,numprocs-1 if (i.lt.nremain) then ishift = i*(naverage+1) + 1 mypoints = naverage + 1 else ishift = i*(naverage+1) + 1 + (nremain-i) mypoints = naverage end if call MPI_ISEND(ishift,1,MPI_INTEGER,i,1, & MPI_COMM_WORLD, irequest(1), ierr) call MPI_ISEND(mypoints,1,MPI_INTEGER,i,2, & MPI_COMM_WORLD, irequest(2), ierr) call MPI_ISEND(x(ishift),mypoints,MPI_DOUBLE_PRECISION,i,3, & MPI_COMM_WORLD, irequest(3), ierr) call MPI_ISEND(y(ishift),mypoints,MPI_DOUBLE_PRECISION,i,4, & MPI_COMM_WORLD, irequest(4), ierr) call MPI_WAITALL(4, irequest, istatus_array, ierr) end do else ! ---------------the other processes receive---------------- call MPI_IRECV(ishift,1,MPI_INTEGER,0,1, & MPI_COMM_WORLD, irequest(1), ierr) call MPI_IRECV(mypoints,1,MPI_INTEGER,0,2, & MPI_COMM_WORLD, irequest(2), ierr) call MPI_WAITALL(2, irequest, istatus_array, ierr) call MPI_IRECV(x(ishift),mypoints,MPI_DOUBLE_PRECISION,0,3, & MPI_COMM_WORLD, irequest(3), ierr) call MPI_IRECV(y(ishift),mypoints,MPI_DOUBLE_PRECISION,0,4, & MPI_COMM_WORLD, irequest(4), ierr) call MPI_WAITALL(2, irequest(3), istatus_array, ierr) ! ---------------------------------------------------------- end if if (myid.ne.0) then write(6,'(A,i2,15f5.2)') 'id',myid, (x(j),j=1,n) end if ! ---------------------------------------------------------- ! Step 3: Each process calculates its partial sum ! ---------------------------------------------------------- mySUMx = 0; mySUMy = 0; mySUMxy = 0; mySUMxx = 0 if (myid.eq.0) then ishift = 1 if (nremain.eq.0) then mypoints = naverage else mypoints = naverage + 1 end if end if j = 0 do while (j.lt.mypoints) mySUMx = mySUMx + x(ishift+j) mySUMy = mySUMy + y(ishift+j) mySUMxy = mySUMxy + x(ishift+j)*y(ishift+j) mySUMxx = mySUMxx + x(ishift+j)*x(ishift+j) j = j + 1 end do ! ---------------------------------------------------------- ! Step 4: Process 0 receives partial sums from the others ! ---------------------------------------------------------- value = numprocs npower = alog(value) / alog(2.) igroup = numprocs/2 SUMx = mySUMx; SUMy = mySUMy SUMxy = mySUMxy; SUMxx = mySUMxx do i = 1,npower isendto = myid-igroup irecvfrom = myid+igroup if (myid.ge.igroup.and.myid.lt.2*igroup) then call MPI_ISEND(SUMx,1,MPI_DOUBLE_PRECISION, & isendto, 5, MPI_COMM_WORLD, irequest(1), ierr) call MPI_ISEND(SUMy,1,MPI_DOUBLE_PRECISION, & isendto, 6, MPI_COMM_WORLD, irequest(2), ierr) call MPI_ISEND(SUMxy,1,MPI_DOUBLE_PRECISION, & isendto, 7, MPI_COMM_WORLD, irequest(3), ierr) call MPI_ISEND(SUMxx,1,MPI_DOUBLE_PRECISION, & isendto, 8, MPI_COMM_WORLD, irequest(4), ierr) call MPI_WAITALL(4, irequest, istatus_array, ierr) else if (myid.lt.igroup) then call MPI_IRECV(mySUMx,1,MPI_DOUBLE_PRECISION, & irecvfrom, 5, MPI_COMM_WORLD, irequest(1), ierr) call MPI_IRECV(mySUMy,1,MPI_DOUBLE_PRECISION, & irecvfrom, 6, MPI_COMM_WORLD, irequest(2), ierr) call MPI_IRECV(mySUMxy,1,MPI_DOUBLE_PRECISION, & irecvfrom, 7, MPI_COMM_WORLD, irequest(3), ierr) call MPI_IRECV(mySUMxx,1,MPI_DOUBLE_PRECISION, & irecvfrom, 8, MPI_COMM_WORLD, irequest(4), ierr) call MPI_WAITALL(4, irequest, istatus_array, ierr) SUMx = SUMx + mySUMx SUMy = SUMy + mySUMy SUMxy = SUMxy + mySUMxy SUMxx = SUMxx + mySUMxx end if igroup = igroup/2 end do ! ---------------------------------------------------------- ! Step 5: Process 0 does the final steps ! ---------------------------------------------------------- if (myid.eq.0) then slope = ( SUMx*SUMy - n*SUMxy ) / ( SUMx*SUMx - n*SUMxx ) y_intercept = ( SUMy - slope*SUMx ) / n call new_sleep (%val(3)) write(6,*) write(6,*) 'The linear equation that best fits the given data:' write(6,11) slope, y_intercept 11 format (10x,'y = ',f6.2,'x + ',f6.2) write(6,*) '--------------------------------------------------' write(6,*) ' Original (x,y) Estimated y Residual' write(6,*) '--------------------------------------------------' SUMres = 0 do i = 1,n y_estimate = slope*x(i) + y_intercept res = y(i) - y_estimate SUMres = SUMres + res*res write(6,12) x(i), y(i), y_estimate, res end do write(6,*) '--------------------------------------------------' write(6,13) SUMres 12 format (3x,'(',f6.2,',',f6.2,')',6x,f6.2,8x,f6.2) 13 format (1x,'Residual sum = ',f6.2/) ! Complete any outstanding communication ! -------------------------------------- call MPI_WAITALL(numprocs-1,nrequest,nstatus_array,ierr) endif ! ---------------------------------------------------------- call MPI_FINALIZE(ierr) end