<p><b>mmwolf@sandia.gov</b> 2009-12-21 10:11:51 -0700 (Mon, 21 Dec 2009)</p><p>My first pass at providing an interface to Zoltan functionality,<br>
in particular the ordering functionality.  Currently solves very simple<br>
test HSFC ordering function.  Still need to pass correct objects to order.<br>
Changes to Makefiles to optionally compile module_zoltan_interface.F<br>
and link Zoltan if ZOLTAN_HOME environmental variable is set.<br>
</p><hr noshade><pre><font color="gray">Modified: trunk/swmodel/Makefile
===================================================================
--- trunk/swmodel/Makefile        2009-12-19 00:11:18 UTC (rev 93)
+++ trunk/swmodel/Makefile        2009-12-21 17:11:51 UTC (rev 94)
@@ -57,14 +57,14 @@
    CFLAGS += -DHAVE_ZOLTAN
    CPPFLAGS += -DHAVE_ZOLTAN
 
-   ifdef ZOLTAN_INCLUDE
-      FCINCLUDES += -I$(ZOLTAN_INCLUDE)
+   ifdef ZOLTAN_INC_PATH
+      FCINCLUDES += -I$(ZOLTAN_INC_PATH)
    else
       FCINCLUDES += -I$(ZOLTAN_HOME)/include
    endif
 
-   ifdef ZOLTAN_LIBPATH
-      LIBS += -L$(ZOLTAN_LIBPATH) -lzoltan
+   ifdef ZOLTAN_LIB_PATH
+      LIBS += -L$(ZOLTAN_LIB_PATH) -lzoltan
    else
       LIBS += -L$(ZOLTAN_HOME)/lib -lzoltan
    endif

Modified: trunk/swmodel/src/Makefile
===================================================================
--- trunk/swmodel/src/Makefile        2009-12-19 00:11:18 UTC (rev 93)
+++ trunk/swmodel/src/Makefile        2009-12-21 17:11:51 UTC (rev 94)
@@ -1,5 +1,10 @@
 .SUFFIXES: .F .c .o
 
+ifdef ZOLTAN_HOME
+   ZOLTANOBJ = module_zoltan_interface.o
+endif
+
+
 OBJS = swmodel.o \
        module_timer.o \
        module_configure.o \
@@ -14,6 +19,7 @@
        module_io_input.o \
        module_io_output.o \
        module_time_integration.o \
+       $(ZOLTANOBJ) \
        streams.o
 
 all: swmodel

Added: trunk/swmodel/src/module_zoltan_interface.F
===================================================================
--- trunk/swmodel/src/module_zoltan_interface.F                                (rev 0)
+++ trunk/swmodel/src/module_zoltan_interface.F        2009-12-21 17:11:51 UTC (rev 94)
@@ -0,0 +1,183 @@
+module zoltan_interface
+   use zoltan
+   use mpi
+
+   implicit none
+
+   contains
+
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+   !! Perhaps not necessary, but implemented in case it helps
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+   subroutine zoltanStart()  
+
+      integer(Zoltan_INT) :: error
+      real(Zoltan_FLOAT) :: version
+
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+      !! Body of subroutine
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+      error = Zoltan_Initialize(version)
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+      
+   end subroutine
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+   subroutine zoltanOrderLocHSFC_Cells()
+
+      implicit none
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+      !! local variables
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+      type(Zoltan_Struct), pointer :: zz_obj
+      integer(ZOLTAN_INT) :: ierr
+
+      integer :: numGidEntries, myrank, locNumObj, i
+      integer(ZOLTAN_INT), allocatable :: global_ids(:)
+      integer(ZOLTAN_INT), allocatable :: permGIDs(:)
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+      !! Body of subroutine
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+      nullify(zz_obj)
+      zz_obj =&gt; Zoltan_Create(MPI_COMM_SELF)
+
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+      !! General Zoltan Parameters
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+      ierr = Zoltan_Set_Param(zz_obj, &quot;ORDER_METHOD&quot;, &quot;LOCAL_HSFC&quot;)
+
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+      !! register query functions
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+      ierr = Zoltan_Set_Fn(zz_obj, ZOLTAN_NUM_OBJ_FN_TYPE,zqfNumCells)
+      ierr = Zoltan_Set_Fn(zz_obj, ZOLTAN_OBJ_LIST_FN_TYPE,zqfGetCells)
+      ierr = Zoltan_Set_Fn(zz_obj, ZOLTAN_NUM_GEOM_FN_TYPE,zqfGeomDim)
+      ierr =  Zoltan_Set_Fn(zz_obj, ZOLTAN_GEOM_FN_TYPE, zqfGetCellGeom)
+
+      locNumObj=16
+      numGidEntries=1
+
+      allocate(global_ids(locNumObj))
+      allocate(permGIDs(locNumObj))
+
+      call MPI_Comm_rank(MPI_COMM_WORLD, myrank, ierr)
+
+      do i=1,locNumObj
+        global_ids(i) = i
+      end do
+
+      ierr = Zoltan_Order(zz_obj, numGidEntries, locNumObj, global_ids, permGIDs);
+
+
+      do i=1,locNumObj
+         write(*,*) global_ids(i), permGIDs(i)
+      end do
+
+      deallocate(global_ids)
+      deallocate(permGIDs)
+
+      call Zoltan_Destroy(zz_obj)
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!                                                                                                                                                 
+
+   end subroutine zoltanOrderLocHSFC_Cells
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+   !! zoltan query function:
+   !!    Returns number of cells
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+   integer function zqfNumCells(data, ierr)
+
+      ! Local declarations
+      integer(ZOLTAN_INT), intent(in) :: data(*)
+      integer(ZOLTAN_INT), intent(out) :: ierr
+
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+      zqfNumCells = 16
+      ierr = ZOLTAN_OK
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+   end function zqfNumCells
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+   !! zoltan query function: 
+   !!    Returns lists of Cell IDs
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+   subroutine zqfGetCells (data, num_gid_entries, num_lid_entries, global_ids, &amp;
+                           local_ids, wgt_dim, obj_wgts, ierr)
+     use mpi
+
+     integer(ZOLTAN_INT), intent(in) :: data(*)
+     integer(ZOLTAN_INT), intent(in) :: num_gid_entries 
+     integer(ZOLTAN_INT), intent(in) ::  num_lid_entries
+     integer(ZOLTAN_INT), intent(out) :: global_ids(*)
+     integer(ZOLTAN_INT), intent(out) :: local_ids(*)
+     integer(ZOLTAN_INT), intent(in) :: wgt_dim 
+     real(ZOLTAN_FLOAT), intent(out) :: obj_wgts(*)
+     integer(ZOLTAN_INT), intent(out) :: ierr
+
+     ! local declarations
+     integer :: i
+
+     !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+     !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+     do i= 1, 16
+       global_ids(i) = i
+       local_ids(i) = i
+     end do
+
+     ierr = ZOLTAN_OK
+     !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+   end subroutine zqfGetCells
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+   !! Zoltan Query Function:
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+   integer function zqfGeomDim(data, ierr)
+      !use zoltan
+      implicit none
+      integer(ZOLTAN_INT), intent(in) :: data(*)
+      integer(ZOLTAN_INT) :: ierr
+
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+      zqfGeomDim = 2
+      ierr = ZOLTAN_OK
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+   end function zqfGeomDim
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+   !! Zoltan Query Function:
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+   subroutine zqfGetCellGeom(data, num_gid_entries, num_lid_entries, global_id, &amp;
+                             local_id, geom_vec, ierr)
+      !use zoltan
+      implicit none
+
+      integer(ZOLTAN_INT), intent(in) :: data(*)
+      integer(ZOLTAN_INT), intent(in) :: num_gid_entries 
+      integer(ZOLTAN_INT), intent(in) ::  num_lid_entries
+      integer(ZOLTAN_INT), intent(in) :: global_id
+      integer(ZOLTAN_INT), intent(in) :: local_id
+      real(ZOLTAN_DOUBLE), intent(out) :: geom_vec(*)
+      integer(ZOLTAN_INT), intent(out) :: ierr
+
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+      geom_vec(1) =   4*( (global_id-1)/32 ) + MODULO(global_id-1,4) + 1;
+      geom_vec(2) =  MODULO(global_id-1,32)   /4   +   1;
+
+      ierr = ZOLTAN_OK
+      !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+   end subroutine zqfGetCellGeom
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+end module zoltan_interface

</font>
</pre>