<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 => Zoltan_Create(MPI_COMM_SELF)
+
+ !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+ !! General Zoltan Parameters
+ !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+ ierr = Zoltan_Set_Param(zz_obj, "ORDER_METHOD", "LOCAL_HSFC")
+
+ !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+ !! 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, &
+ 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, &
+ 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>