<p><b>duda</b> 2009-08-17 18:10:26 -0600 (Mon, 17 Aug 2009)</p><p>When reading namelist and block decomposition files, have only one process<br>
read the files and distribute the information to all other processes, rather <br>
than allow all processes to simultaneously open and read files.<br>
<br>
M module_dmpar.F<br>
M swmodel.F<br>
M module_configure.F<br>
M module_block_decomp.F<br>
M Makefile<br>
</p><hr noshade><pre><font color="gray">Modified: trunk/swmodel/Makefile
===================================================================
--- trunk/swmodel/Makefile        2009-08-15 04:57:41 UTC (rev 30)
+++ trunk/swmodel/Makefile        2009-08-18 00:10:26 UTC (rev 31)
@@ -50,6 +50,8 @@
swmodel.o: module_configure.o module_dmpar.o module_grid_types.o module_test_cases.o module_io_input.o module_sw_solver.o module_timer.o
+module_configure.o: module_dmpar.o
+
module_grid_types.o: module_dmpar.o
module_dmpar.o: module_sort.o streams.o
Modified: trunk/swmodel/module_block_decomp.F
===================================================================
--- trunk/swmodel/module_block_decomp.F        2009-08-15 04:57:41 UTC (rev 30)
+++ trunk/swmodel/module_block_decomp.F        2009-08-18 00:10:26 UTC (rev 31)
@@ -23,53 +23,84 @@
type (dm_info), intent(in) :: dminfo
type (graph), intent(in) :: partial_global_graph_info
integer, dimension(:), pointer :: local_cell_list
+ integer, dimension(:), pointer :: global_cell_list
+ integer, dimension(:), pointer :: global_start
- integer :: i, j, owner, iunit, istatus, local_nvertices
+ integer :: i, j, owner, iunit, istatus
+ integer, dimension(:), pointer :: local_nvertices
character (len=256) :: filename
if (dminfo % nprocs > 1) then
- iunit = 50 + dminfo % my_proc_id
- if (dminfo % nprocs < 10) then
- write(filename,'(a,i1)') 'graph.info.part.', dminfo % nprocs
- else if (dminfo % nprocs < 100) then
- write(filename,'(a,i2)') 'graph.info.part.', dminfo % nprocs
- else if (dminfo % nprocs < 1000) then
- write(filename,'(a,i3)') 'graph.info.part.', dminfo % nprocs
- else if (dminfo % nprocs < 10000) then
- write(filename,'(a,i4)') 'graph.info.part.', dminfo % nprocs
- else if (dminfo % nprocs < 100000) then
- write(filename,'(a,i5)') 'graph.info.part.', dminfo % nprocs
+ allocate(local_nvertices(dminfo % nprocs))
+ allocate(global_start(dminfo % nprocs))
+
+ if (dminfo % my_proc_id == IO_NODE) then
+
+ iunit = 50 + dminfo % my_proc_id
+ if (dminfo % nprocs < 10) then
+ write(filename,'(a,i1)') 'graph.info.part.', dminfo % nprocs
+ else if (dminfo % nprocs < 100) then
+ write(filename,'(a,i2)') 'graph.info.part.', dminfo % nprocs
+ else if (dminfo % nprocs < 1000) then
+ write(filename,'(a,i3)') 'graph.info.part.', dminfo % nprocs
+ else if (dminfo % nprocs < 10000) then
+ write(filename,'(a,i4)') 'graph.info.part.', dminfo % nprocs
+ else if (dminfo % nprocs < 100000) then
+ write(filename,'(a,i5)') 'graph.info.part.', dminfo % nprocs
+ end if
+
+ open(unit=iunit, file=trim(filename), form='formatted', status='old', iostat=istatus)
+
+ if (istatus /= 0) then
+ write(0,*) 'Could not open block decomposition file for ',dminfo % nprocs,' tasks.'
+ write(0,*) 'Filename: ',trim(filename)
+ call dmpar_abort(dminfo)
+ end if
+
+ local_nvertices(:) = 0
+ do i=1,partial_global_graph_info % nVerticesTotal
+ read(unit=iunit, fmt=*) owner
+ local_nvertices(owner+1) = local_nvertices(owner+1) + 1
+ end do
+
+ allocate(global_cell_list(partial_global_graph_info % nVerticesTotal))
+
+ global_start(1) = 1
+ do i=2,dminfo % nprocs
+ global_start(i) = global_start(i-1) + local_nvertices(i-1)
+ end do
+
+ rewind(unit=iunit)
+
+ do i=1,partial_global_graph_info % nVerticesTotal
+ read(unit=iunit, fmt=*) owner
+ global_cell_list(global_start(owner+1)) = i
+ global_start(owner+1) = global_start(owner+1) + 1
+ end do
+
+ global_start(1) = 0
+ do i=2,dminfo % nprocs
+ global_start(i) = global_start(i-1) + local_nvertices(i-1)
+ end do
+
+ close(unit=iunit)
+
+ call dmpar_bcast_ints(dminfo, dminfo % nprocs, local_nvertices)
+ allocate(local_cell_list(local_nvertices(dminfo % my_proc_id + 1)))
+
+ call dmpar_scatter_ints(dminfo, dminfo % nprocs, local_nvertices(dminfo % my_proc_id + 1), &
+ global_start, local_nvertices, global_cell_list, local_cell_list)
+
+ else
+
+ call dmpar_bcast_ints(dminfo, dminfo % nprocs, local_nvertices)
+ allocate(local_cell_list(local_nvertices(dminfo % my_proc_id + 1)))
+
+ call dmpar_scatter_ints(dminfo, dminfo % nprocs, local_nvertices(dminfo % my_proc_id + 1), &
+ global_start, local_nvertices, global_cell_list, local_cell_list)
+
end if
-
- open(unit=iunit, file=trim(filename), form='formatted', status='old', iostat=istatus)
-
- if (istatus /= 0) then
- write(0,*) 'Could not open block decomposition file for ',dminfo % nprocs,' tasks.'
- write(0,*) 'Filename: ',trim(filename)
- call dmpar_abort(dminfo)
- end if
-
- local_nvertices = 0
- do i=1,partial_global_graph_info % nVerticesTotal
- read(unit=iunit, fmt=*) owner
- if (owner == dminfo % my_proc_id) local_nvertices = local_nvertices + 1
- end do
-
- allocate(local_cell_list(local_nvertices))
-
- rewind(unit=iunit)
-
- j = 0
- do i=1,partial_global_graph_info % nVerticesTotal
- read(unit=iunit, fmt=*) owner
- if (owner == dminfo % my_proc_id) then
- j = j + 1
- local_cell_list(j) = i
- end if
- end do
-
- close(unit=iunit)
else
allocate(local_cell_list(partial_global_graph_info % nVerticesTotal))
local_cell_list(:) = dminfo % my_proc_id
Modified: trunk/swmodel/module_configure.F
===================================================================
--- trunk/swmodel/module_configure.F        2009-08-15 04:57:41 UTC (rev 30)
+++ trunk/swmodel/module_configure.F        2009-08-18 00:10:26 UTC (rev 31)
@@ -1,5 +1,7 @@
module configure
+ use dmpar
+
integer :: config_test_case
character (len=32) :: config_time_integration
real (kind=RKIND) :: config_dt
@@ -10,10 +12,12 @@
contains
- subroutine read_namelist()
+ subroutine read_namelist(dminfo)
implicit none
+ type (dm_info), intent(in) :: dminfo
+
integer :: funit
namelist /sw_model/ config_test_case, &
@@ -25,15 +29,22 @@
funit = 21
! Set default values for namelist options
- config_test_case = 1
+ config_test_case = 5
config_time_integration = "RK4"
- config_dt = 0.0001
- config_ntimesteps = 100
- config_output_interval = 100
+ config_dt = 172.8
+ config_ntimesteps = 7500
+ config_output_interval = 500
- open(funit,file='namelist.input',status='old',form='formatted')
- read(funit,sw_model)
- close(funit)
+ if (dminfo % my_proc_id == IO_NODE) then
+ open(funit,file='namelist.input',status='old',form='formatted')
+ read(funit,sw_model)
+ close(funit)
+ end if
+ call dmpar_bcast_int(dminfo, config_test_case)
+ call dmpar_bcast_char(dminfo, config_time_integration)
+ call dmpar_bcast_real(dminfo, config_dt)
+ call dmpar_bcast_int(dminfo, config_ntimesteps)
+ call dmpar_bcast_int(dminfo, config_output_interval)
end subroutine read_namelist
Modified: trunk/swmodel/module_dmpar.F
===================================================================
--- trunk/swmodel/module_dmpar.F        2009-08-15 04:57:41 UTC (rev 30)
+++ trunk/swmodel/module_dmpar.F        2009-08-18 00:10:26 UTC (rev 31)
@@ -107,6 +107,88 @@
end subroutine dmpar_abort
+ subroutine dmpar_bcast_int(dminfo, i)
+
+ implicit none
+
+ type (dm_info), intent(in) :: dminfo
+ integer, intent(inout) :: i
+
+#ifdef _MPI
+ integer :: mpi_ierr
+
+ call MPI_Bcast(i, 1, MPI_INTEGER, IO_NODE, dminfo % comm, mpi_ierr)
+#endif
+
+ end subroutine dmpar_bcast_int
+
+
+ subroutine dmpar_bcast_ints(dminfo, n, iarray)
+
+ implicit none
+
+ type (dm_info), intent(in) :: dminfo
+ integer, intent(in) :: n
+ integer, dimension(n), intent(inout) :: iarray
+
+#ifdef _MPI
+ integer :: mpi_ierr
+
+ call MPI_Bcast(iarray, n, MPI_INTEGER, IO_NODE, dminfo % comm, mpi_ierr)
+#endif
+
+ end subroutine dmpar_bcast_ints
+
+
+ subroutine dmpar_bcast_real(dminfo, r)
+
+ implicit none
+
+ type (dm_info), intent(in) :: dminfo
+ real (kind=RKIND), intent(inout) :: r
+
+#ifdef _MPI
+ integer :: mpi_ierr
+
+ call MPI_Bcast(r, 1, MPI_REALKIND, IO_NODE, dminfo % comm, mpi_ierr)
+#endif
+
+ end subroutine dmpar_bcast_real
+
+
+ subroutine dmpar_bcast_reals(dminfo, n, rarray)
+
+ implicit none
+
+ type (dm_info), intent(in) :: dminfo
+ integer, intent(in) :: n
+ real (kind=RKIND), dimension(n), intent(inout) :: rarray
+
+#ifdef _MPI
+ integer :: mpi_ierr
+
+ call MPI_Bcast(rarray, n, MPI_REALKIND, IO_NODE, dminfo % comm, mpi_ierr)
+#endif
+
+ end subroutine dmpar_bcast_reals
+
+
+ subroutine dmpar_bcast_char(dminfo, c)
+
+ implicit none
+
+ type (dm_info), intent(in) :: dminfo
+ character (len=*), intent(inout) :: c
+
+#ifdef _MPI
+ integer :: mpi_ierr
+
+ call MPI_Bcast(c, len(c), MPI_CHARACTER, IO_NODE, dminfo % comm, mpi_ierr)
+#endif
+
+ end subroutine dmpar_bcast_char
+
+
subroutine dmpar_sum_int(dminfo, i, isum)
implicit none
@@ -126,6 +208,25 @@
end subroutine dmpar_sum_int
+ subroutine dmpar_scatter_ints(dminfo, nprocs, noutlist, displs, counts, inlist, outlist)
+
+ implicit none
+
+ type (dm_info), intent(in) :: dminfo
+ integer, intent(in) :: nprocs, noutlist
+ integer, dimension(nprocs), intent(in) :: displs, counts
+ integer, dimension(:), pointer :: inlist
+ integer, dimension(noutlist), intent(inout) :: outlist
+
+#ifdef _MPI
+ integer :: mpi_ierr
+
+ call MPI_Scatterv(inlist, counts, displs, MPI_INTEGER, outlist, noutlist, MPI_INTEGER, IO_NODE, dminfo % comm, mpi_ierr)
+#endif
+
+ end subroutine dmpar_scatter_ints
+
+
subroutine dmpar_get_index_range(dminfo, &
global_start, global_end, &
local_start, local_end)
Modified: trunk/swmodel/swmodel.F
===================================================================
--- trunk/swmodel/swmodel.F        2009-08-15 04:57:41 UTC (rev 30)
+++ trunk/swmodel/swmodel.F        2009-08-18 00:10:26 UTC (rev 31)
@@ -13,12 +13,11 @@
type (dm_info), pointer :: dminfo
type (domain_type), pointer :: domain
-
- call read_namelist()
-
allocate(dminfo)
call dmpar_init(dminfo)
+ call read_namelist(dminfo)
+
call timer_start("total time")
call timer_start("initialize")
</font>
</pre>