<p><b>duda</b> 2009-11-16 16:17:45 -0700 (Mon, 16 Nov 2009)</p><p>Add code to perform halo exchanges for 1d real fields.<br>
<br>
M    module_dmpar.F<br>
</p><hr noshade><pre><font color="gray">Modified: trunk/swmodel/src/module_dmpar.F
===================================================================
--- trunk/swmodel/src/module_dmpar.F        2009-11-16 23:05:46 UTC (rev 71)
+++ trunk/swmodel/src/module_dmpar.F        2009-11-16 23:17:45 UTC (rev 72)
@@ -1228,6 +1228,66 @@
    end subroutine unpackRecvBuf3dReal
 
 
+   subroutine dmpar_exch_halo_field1dReal(dminfo, array, dim1, sendList, recvList)
+
+      implicit none
+
+      type (dm_info), intent(in) :: dminfo
+      integer, intent(in) :: dim1
+      real (kind=RKIND), dimension(dim1), intent(inout) :: array
+      type (exchange_list), pointer :: sendList, recvList
+
+      type (exchange_list), pointer :: sendListPtr, recvListPtr
+      integer :: lastPackedIdx, lastUnpackedIdx, nPacked, nUnpacked
+      integer :: mpi_ierr
+
+#ifdef _MPI
+
+      recvListPtr =&gt; recvList
+      do while (associated(recvListPtr))
+         if (recvListPtr % procID /= dminfo % my_proc_id) then
+            allocate(recvListPtr % rbuffer(recvListPtr % nlist))
+            call MPI_Irecv(recvListPtr % rbuffer, recvListPtr % nlist, MPI_REALKIND, &amp;
+                           recvListPtr % procID, recvListPtr % procID, dminfo % comm, recvListPtr % reqID, mpi_ierr)
+         end if
+         recvListPtr =&gt; recvListPtr % next
+      end do
+      
+      sendListPtr =&gt; sendList
+      do while (associated(sendListPtr))
+         if (sendListPtr % procID /= dminfo % my_proc_id) then
+            allocate(sendListPtr % rbuffer(sendListPtr % nlist))
+            call packSendBuf1dReal(dim1, array, sendListPtr, 1, sendListPtr % nlist, sendListPtr % rbuffer, nPacked, lastPackedIdx)
+            call MPI_Isend(sendListPtr % rbuffer, sendListPtr % nlist, MPI_REALKIND, &amp;
+                           sendListPtr % procID, dminfo % my_proc_id, dminfo % comm, sendListPtr % reqID, mpi_ierr)
+         end if
+         sendListPtr =&gt; sendListPtr % next
+      end do
+
+      recvListPtr =&gt; recvList
+      do while (associated(recvListPtr))
+         if (recvListPtr % procID /= dminfo % my_proc_id) then
+            call MPI_Wait(recvListPtr % reqID, MPI_STATUS_IGNORE, mpi_ierr)
+            call unpackRecvBuf1dReal(dim1, array, recvListPtr, 1, recvListPtr % nlist, recvListPtr % rbuffer, nUnpacked, lastUnpackedIdx)
+            deallocate(recvListPtr % rbuffer)
+         end if
+         recvListPtr =&gt; recvListPtr % next
+      end do
+      
+      sendListPtr =&gt; sendList
+      do while (associated(sendListPtr))
+         if (sendListPtr % procID /= dminfo % my_proc_id) then
+            call MPI_Wait(sendListPtr % reqID, MPI_STATUS_IGNORE, mpi_ierr)
+            deallocate(sendListPtr % rbuffer)
+         end if
+         sendListPtr =&gt; sendListPtr % next
+      end do
+
+#endif
+
+   end subroutine dmpar_exch_halo_field1dReal
+
+
    subroutine dmpar_exch_halo_field2dReal(dminfo, array, dim1, dim2, sendList, recvList)
 
       implicit none

</font>
</pre>