From 3998949cda69d0724e7c588144fad6b964a1ef9f Mon Sep 17 00:00:00 2001 From: Abishek Gopal Date: Tue, 8 Apr 2025 16:53:52 -0600 Subject: [PATCH 1/7] Working with upto 24km global mesh --- Makefile | 5 + src/core_atmosphere/Registry.xml | 2 +- src/framework/mpas_block_decomp.F | 209 +++++++++++++++++++++-------- src/framework/mpas_bootstrapping.F | 18 ++- 4 files changed, 176 insertions(+), 58 deletions(-) diff --git a/Makefile b/Makefile index c34823ab5a..7606b47032 100644 --- a/Makefile +++ b/Makefile @@ -759,6 +759,11 @@ endif LIBS += $(NCLIB) endif +export SCOTCH_ROOT=/glade/derecho/scratch/agopal/scotch/build + +FCINCLUDES += -I$(SCOTCH_ROOT)/src/include + +LIBS += -L$(SCOTCH_ROOT)/lib -lscotch -lscotcherr ifneq "$(PNETCDF)" "" ifneq ($(wildcard $(PNETCDF)/lib/libpnetcdf.*), ) diff --git a/src/core_atmosphere/Registry.xml b/src/core_atmosphere/Registry.xml index bdf47ec7b4..ade3916a16 100644 --- a/src/core_atmosphere/Registry.xml +++ b/src/core_atmosphere/Registry.xml @@ -319,7 +319,7 @@ - diff --git a/src/framework/mpas_block_decomp.F b/src/framework/mpas_block_decomp.F index 4f3d197d5d..3079aed4f1 100644 --- a/src/framework/mpas_block_decomp.F +++ b/src/framework/mpas_block_decomp.F @@ -25,6 +25,7 @@ module mpas_block_decomp use mpas_derived_types use mpas_io_units use mpas_log + include "scotchf.h" type graph integer :: nVerticesTotal @@ -49,7 +50,7 @@ module mpas_block_decomp ! !----------------------------------------------------------------------- subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, local_cell_list, block_id, block_start, & - block_count, numBlocks, explicitProcDecomp, blockFilePrefix, procFilePrefix)!{{{ + block_count, cellsOnCellFull, nEdgesOnCellFull, numBlocks, explicitProcDecomp, blockFilePrefix, procFilePrefix)!{{{ implicit none @@ -59,6 +60,8 @@ subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, l integer, dimension(:), pointer :: block_id !< Output: list of global block id's this processor owns integer, dimension(:), pointer :: block_start !< Output: offset in local_cell_list for this blocks list of cells integer, dimension(:), pointer :: block_count !< Output: number of cells in blocks + type (field2dInteger), intent(in) :: cellsOnCellFull + type (field1dInteger), intent(in) :: nEdgesOnCellFull integer, intent(in) :: numBlocks !< Input: Number of blocks (from config_num_blocks) logical, intent(in) :: explicitProcDecomp !< Input: Logical flag controlling if blocks are explicitly assigned to processors @@ -71,12 +74,19 @@ subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, l integer, dimension(:), allocatable :: local_block_list integer, dimension(:,:), allocatable :: sorted_local_cell_list - integer :: i, global_block_id, local_block_id, owning_proc, iunit, istatus - integer :: blocks_per_proc, err + integer, dimension(:), allocatable :: global_block_id_arr, owning_proc_arr + integer :: i, global_block_id, local_block_id, owning_proc, iunit, istatus, vind, j, k + integer :: blocks_per_proc, err, ierr integer, dimension(:), pointer :: local_nvertices - character (len=StrKIND) :: filename + character (len=StrKIND) :: filename, msg logical :: no_blocks + integer :: nTotalEdgesGraph = 0 + integer, dimension(:), allocatable :: edgetab, verttab, parttab + + doubleprecision :: stradat (SCOTCH_STRATDIM) + doubleprecision :: SCOTCHGRAPH (SCOTCH_GRAPHDIM) + integer :: n_size no_blocks = .false. @@ -95,71 +105,158 @@ subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, l allocate(local_nvertices(dminfo % nprocs)) allocate(global_start(dminfo % nprocs)) allocate(global_list(partial_global_graph_info % nVerticesTotal)) + allocate(global_block_id_arr(partial_global_graph_info % nVerticesTotal)) + allocate(owning_proc_arr(partial_global_graph_info % nVerticesTotal)) if (dminfo % my_proc_id == IO_NODE) then + if ( trim(blockFilePrefix) == '' ) then + call mpas_log_write('blockFilePrefix is not set: Using LibScotch for graph partitioning') + do vind=1,partial_global_graph_info % nVerticesTotal + nTotalEdgesGraph = nTotalEdgesGraph + nEdgesOnCellFull% array(vind) + ! do j=1,nEdgesOnCellFull% array(vind) + ! call mpas_log_write('vind=$i j=$i adj= $i', intArgs=(/vind,j,cellsOnCellFull%array(j,vind)/) ) + ! end do + end do + + allocate(edgetab(nTotalEdgesGraph)) + allocate(verttab(partial_global_graph_info % nVerticesTotal + 1)) + !allocate(parttab(partial_global_graph_info % nVerticesTotal)) + + !do vind=1,partial_global_graph_info % nVerticesTotal + !call mpas_log_write('proc=$i vind= $i part= $i', intArgs=(/dminfo % my_proc_id, vind,parttab(vind)/)) + !call mpas_log_write('vind=$i j=$i adj= $i', intArgs=(/vind,j,cellsOnCellFull%array(j,vind)/) ) + !end do + + k = 1 + do vind=1,partial_global_graph_info % nVerticesTotal + verttab(vind) = k + !call mpas_log_write('vind=$i verttab= $i', intArgs=(/vind,verttab(vind)/) ) + do j=1,nEdgesOnCellFull% array(vind) + edgetab(k) = cellsOnCellFull%array(j,vind) + !call mpas_log_write('k=$i edgetab= $i', intArgs=(/k,edgetab(k)/) ) + k = k + 1 + end do + end do + verttab(partial_global_graph_info % nVerticesTotal+1) = nTotalEdgesGraph + 1 + + !call mpas_log_write('nvertices =$i nTotalEdgesGraph= $i', intArgs=(/partial_global_graph_info % nVerticesTotal, nTotalEdgesGraph/)) + + CALL scotchfstratinit (stradat (1), IERR) + CALL SCOTCHFGRAPHINIT (SCOTCHGRAPH (1), IERR) + IF (IERR .NE. 0) THEN + call mpas_log_write('Cannot initialize Scotch Graph', MPAS_LOG_CRIT) + ENDIF + + CALL SCOTCHFGRAPHBUILD (SCOTCHGRAPH (1), 1, partial_global_graph_info % nVerticesTotal, & + verttab (1), verttab (2), & + verttab (1), verttab (1), & + nTotalEdgesGraph, edgetab (1), edgetab (1), IERR) + IF (IERR .NE. 0) THEN + call mpas_log_write('Cannot build Scotch Graph', MPAS_LOG_CRIT) + ENDIF + + CALL SCOTCHFGRAPHCHECK (SCOTCHGRAPH (1), IERR) + IF (IERR .NE. 0) THEN + call mpas_log_write('Cannot check Scotch Graph', MPAS_LOG_CRIT) + ENDIF + + !CALL scotchfgraphpart( SCOTCHGRAPH (1), dminfo % total_blocks, stradat (1) ,parttab(1), IERR) + CALL scotchfgraphpart( SCOTCHGRAPH (1), dminfo % total_blocks, stradat (1) ,global_block_id_arr(1), IERR) + + call scotchfgraphexit (SCOTCHGRAPH (1)) + call scotchfstratexit (stradat (1)) + + local_nvertices(:) = 0 + do i=1,partial_global_graph_info % nVerticesTotal + !global_block_id_arr(i) = parttab(i) + call mpas_get_owning_proc(dminfo, global_block_id_arr(i), owning_proc) + owning_proc_arr(i) = owning_proc + !call mpas_log_write('vert: $i, global_block_id: $i, owning_proc: $i ', intArgs=(/i,global_block_id_arr(i),owning_proc/) ) + local_nvertices(owning_proc+1) = local_nvertices(owning_proc+1) + 1 + end do + + else + + if (dminfo % total_blocks < 10) then + write(filename,'(a,i1)') trim(blockFilePrefix), dminfo % total_blocks + else if (dminfo % total_blocks < 100) then + write(filename,'(a,i2)') trim(blockFilePrefix), dminfo % total_blocks + else if (dminfo % total_blocks < 1000) then + write(filename,'(a,i3)') trim(blockFilePrefix), dminfo % total_blocks + else if (dminfo % total_blocks < 10000) then + write(filename,'(a,i4)') trim(blockFilePrefix), dminfo % total_blocks + else if (dminfo % total_blocks < 100000) then + write(filename,'(a,i5)') trim(blockFilePrefix), dminfo % total_blocks + else if (dminfo % total_blocks < 1000000) then + write(filename,'(a,i6)') trim(blockFilePrefix), dminfo % total_blocks + else if (dminfo % total_blocks < 10000000) then + write(filename,'(a,i7)') trim(blockFilePrefix), dminfo % total_blocks + else if (dminfo % total_blocks < 100000000) then + write(filename,'(a,i8)') trim(blockFilePrefix), dminfo % total_blocks + end if - if (dminfo % total_blocks < 10) then - write(filename,'(a,i1)') trim(blockFilePrefix), dminfo % total_blocks - else if (dminfo % total_blocks < 100) then - write(filename,'(a,i2)') trim(blockFilePrefix), dminfo % total_blocks - else if (dminfo % total_blocks < 1000) then - write(filename,'(a,i3)') trim(blockFilePrefix), dminfo % total_blocks - else if (dminfo % total_blocks < 10000) then - write(filename,'(a,i4)') trim(blockFilePrefix), dminfo % total_blocks - else if (dminfo % total_blocks < 100000) then - write(filename,'(a,i5)') trim(blockFilePrefix), dminfo % total_blocks - else if (dminfo % total_blocks < 1000000) then - write(filename,'(a,i6)') trim(blockFilePrefix), dminfo % total_blocks - else if (dminfo % total_blocks < 10000000) then - write(filename,'(a,i7)') trim(blockFilePrefix), dminfo % total_blocks - else if (dminfo % total_blocks < 100000000) then - write(filename,'(a,i8)') trim(blockFilePrefix), dminfo % total_blocks - end if - - call mpas_new_unit(iunit) - open(unit=iunit, file=trim(filename), form='formatted', status='old', iostat=istatus) - - if (istatus /= 0) then - call mpas_log_write('Could not open block decomposition file for $i blocks.', MPAS_LOG_ERR, intArgs=(/dminfo % total_blocks/) ) - call mpas_log_write('Filename: '//trim(filename), MPAS_LOG_CRIT) - end if - - local_nvertices(:) = 0 - do i=1,partial_global_graph_info % nVerticesTotal - read(unit=iunit, fmt=*, iostat=err) global_block_id - - if ( err .ne. 0 ) then - call mpas_log_write('Decomoposition file: ' // trim(filename) // ' contains less than $i cells', & - MPAS_LOG_CRIT, intArgs=(/partial_global_graph_info % nVerticesTotal/) ) - end if - call mpas_get_owning_proc(dminfo, global_block_id, owning_proc) - local_nvertices(owning_proc+1) = local_nvertices(owning_proc+1) + 1 - end do + call mpas_new_unit(iunit) + open(unit=iunit, file=trim(filename), form='formatted', status='old', iostat=istatus) + + if (istatus /= 0) then + call mpas_log_write('Could not open block decomposition file for $i blocks.', MPAS_LOG_ERR, intArgs=(/dminfo % total_blocks/) ) + call mpas_log_write('Filename: '//trim(filename), MPAS_LOG_CRIT) + end if - read(unit=iunit, fmt=*, iostat=err) + call mpas_log_write('Using block decomposition file: '//trim(filename)) + + call mpas_log_write('First read pass ') + local_nvertices(:) = 0 + do i=1,partial_global_graph_info % nVerticesTotal + read(unit=iunit, fmt=*, iostat=err) global_block_id + global_block_id_arr(i) = global_block_id + + if ( err .ne. 0 ) then + call mpas_log_write('Decomoposition file: ' // trim(filename) // ' contains less than $i cells', & + MPAS_LOG_CRIT, intArgs=(/partial_global_graph_info % nVerticesTotal/) ) + end if + call mpas_get_owning_proc(dminfo, global_block_id, owning_proc) + owning_proc_arr(i) = owning_proc + !call mpas_log_write('vert: $i, global_block_id: $i, owning_proc: $i ', intArgs=(/i,global_block_id,owning_proc/) ) + local_nvertices(owning_proc+1) = local_nvertices(owning_proc+1) + 1 + !call mpas_log_write('owning_proc+1: $i local_nvertices= $i', intArgs=(/owning_proc+1, local_nvertices(owning_proc+1)/)) + end do + + read(unit=iunit, fmt=*, iostat=err) + + if ( err == 0 ) then + call mpas_log_write('Decomposition file: ' // trim(filename) // ' contains more than $i cells', & + MPAS_LOG_CRIT, intArgs=(/partial_global_graph_info % nVerticesTotal/) ) + end if - if ( err == 0 ) then - call mpas_log_write('Decomposition file: ' // trim(filename) // ' contains more than $i cells', & - MPAS_LOG_CRIT, intArgs=(/partial_global_graph_info % nVerticesTotal/) ) + close(unit=iunit) + call mpas_release_unit(iunit) end if global_start(1) = 1 do i=2,dminfo % nprocs global_start(i) = global_start(i-1) + local_nvertices(i-1) + !call mpas_log_write('i: $i global_start= $i', intArgs=(/i, global_start(i)/)) end do - rewind(unit=iunit) + !rewind(unit=iunit) + call mpas_log_write('Second read pass ') do i=1,partial_global_graph_info % nVerticesTotal - read(unit=iunit, fmt=*, iostat=err) global_block_id - call mpas_get_owning_proc(dminfo, global_block_id, owning_proc) + !read(unit=iunit, fmt=*, iostat=err) global_block_id + !call mpas_get_owning_proc(dminfo, global_block_id, owning_proc) + global_block_id = global_block_id_arr(i) + owning_proc = owning_proc_arr(i) global_list(global_start(owning_proc+1)) = i + !call mpas_log_write('owning_proc+1: $i global_start= $i global_list=$i', intArgs=(/owning_proc+1, global_start(owning_proc+1), global_list(global_start(owning_proc+1))/)) global_start(owning_proc+1) = global_start(owning_proc+1) + 1 + !call mpas_log_write('global_start(owning_proc+1): $i', intArgs=(/global_start(owning_proc+1)/)) end do global_start(1) = 0 do i=2,dminfo % nprocs global_start(i) = global_start(i-1) + local_nvertices(i-1) + !call mpas_log_write('i: $i global_start= $i', intArgs=(/i, global_start(i)/)) end do call mpas_dmpar_bcast_ints(dminfo, dminfo % nprocs, local_nvertices) @@ -173,28 +270,34 @@ subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, l global_start(1) = 1 do i=2,dminfo % nprocs global_start(i) = global_start(i-1) + local_nvertices(i-1) + !call mpas_log_write('i: $i global_start= $i', intArgs=(/i, global_start(i)/)) end do - rewind(unit=iunit) + !rewind(unit=iunit) + call mpas_log_write('Third read pass ') do i=1,partial_global_graph_info % nVerticesTotal - read(unit=iunit, fmt=*) global_block_id - call mpas_get_owning_proc(dminfo, global_block_id, owning_proc) + !read(unit=iunit, fmt=*) global_block_id + !call mpas_get_owning_proc(dminfo, global_block_id, owning_proc) + global_block_id = global_block_id_arr(i) + owning_proc = owning_proc_arr(i) global_list(global_start(owning_proc+1)) = global_block_id + !call mpas_log_write('owning_proc+1: $i global_start= $i global_list=$i', intArgs=(/owning_proc+1, global_start(owning_proc+1), global_list(global_start(owning_proc+1))/)) global_start(owning_proc+1) = global_start(owning_proc+1) + 1 + !call mpas_log_write('global_start(owning_proc+1): $i', intArgs=(/global_start(owning_proc+1)/)) end do ! Recompute global start after second read of global_block_list global_start(1) = 0 do i=2,dminfo % nprocs global_start(i) = global_start(i-1) + local_nvertices(i-1) + !call mpas_log_write('i: $i global_start= $i', intArgs=(/i, global_start(i)/)) end do call mpas_dmpar_scatter_ints(dminfo, dminfo % nprocs, local_nvertices(dminfo % my_proc_id + 1), & global_start, local_nvertices, global_list, local_block_list) - close(unit=iunit) - call mpas_release_unit(iunit) + else diff --git a/src/framework/mpas_bootstrapping.F b/src/framework/mpas_bootstrapping.F index 4241255e2a..f444032e3b 100644 --- a/src/framework/mpas_bootstrapping.F +++ b/src/framework/mpas_bootstrapping.F @@ -106,8 +106,8 @@ subroutine mpas_bootstrap_framework_phase1(domain, mesh_filename, mesh_iotype, p type (field1dInteger), pointer :: indexToCellIDField type (field1dInteger), pointer :: indexToEdgeIDField type (field1dInteger), pointer :: indexToVertexIDField - type (field1dInteger), pointer :: nEdgesOnCellField - type (field2dInteger), pointer :: cellsOnCellField + type (field1dInteger), pointer :: nEdgesOnCellField, nEdgesOnCellFull + type (field2dInteger), pointer :: cellsOnCellField, cellsOnCellFull type (field2dInteger), pointer :: edgesOnCellField type (field2dInteger), pointer :: verticesOnCellField type (field2dInteger), pointer :: cellsOnEdgeField @@ -144,7 +144,7 @@ subroutine mpas_bootstrap_framework_phase1(domain, mesh_filename, mesh_iotype, p integer, pointer :: config_num_halos, config_number_of_blocks logical, pointer :: config_explicit_proc_decomp character (len=StrKIND), pointer :: config_block_decomp_file_prefix, config_proc_decomp_file_prefix - integer :: nHalos + integer :: nHalos, j, vind call mpas_pool_get_config(domain % configs, 'config_num_halos', config_num_halos) @@ -197,6 +197,16 @@ subroutine mpas_bootstrap_framework_phase1(domain, mesh_filename, mesh_iotype, p ! which cells/edges/vertices are owned by each block, and which are ghost ! + allocate(cellsOnCellFull) + allocate(cellsOnCellFull % array(maxEdges,nCells)) + call MPAS_io_inq_var(inputHandle, 'cellsOnCell', ierr=ierr) + call mpas_io_get_var(inputHandle, 'cellsOnCell', cellsOnCellFull % array, ierr) + + allocate(nEdgesOnCellFull) + allocate(nEdgesOnCellFull % array(nCells)) + call MPAS_io_inq_var(inputHandle, 'nEdgesOnCell', ierr=ierr) + call mpas_io_get_var(inputHandle, 'nEdgesOnCell', nEdgesOnCellFull % array, ierr) + call mpas_io_setup_cell_block_fields(inputHandle, nreadCells, readCellStart, readingBlock, maxEdges, & indexTocellIDField, xCellField, yCellField, zCellField, nEdgesOnCellField, & cellsOnCellField, edgesOnCellField, verticesOnCellField, nHalos) @@ -235,7 +245,7 @@ subroutine mpas_bootstrap_framework_phase1(domain, mesh_filename, mesh_iotype, p ! file, but in principle it could call some online, distributed mesh partitioning library. ! call mpas_block_decomp_cells_for_proc(domain % dminfo, partial_global_graph_info, local_cell_list, block_id, block_start, & - block_count, config_number_of_blocks, config_explicit_proc_decomp, & + block_count, cellsOnCellFull, nEdgesOnCellFull, config_number_of_blocks, config_explicit_proc_decomp, & config_block_decomp_file_prefix, config_proc_decomp_file_prefix) deallocate(partial_global_graph_info % vertexID) From 4e5564b1b5cd045016dc56fb4f9f6d3baf6664d0 Mon Sep 17 00:00:00 2001 From: Abishek Gopal Date: Wed, 9 Apr 2025 10:08:41 -0600 Subject: [PATCH 2/7] Now checking for zero connectivities. Works with limited area grid --- src/framework/mpas_block_decomp.F | 42 ++++++++++++++++++++----------- 1 file changed, 28 insertions(+), 14 deletions(-) diff --git a/src/framework/mpas_block_decomp.F b/src/framework/mpas_block_decomp.F index 3079aed4f1..81da41292b 100644 --- a/src/framework/mpas_block_decomp.F +++ b/src/framework/mpas_block_decomp.F @@ -75,7 +75,7 @@ subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, l integer, dimension(:,:), allocatable :: sorted_local_cell_list integer, dimension(:), allocatable :: global_block_id_arr, owning_proc_arr - integer :: i, global_block_id, local_block_id, owning_proc, iunit, istatus, vind, j, k + integer :: i, global_block_id, local_block_id, owning_proc, iunit, istatus, j, k integer :: blocks_per_proc, err, ierr integer, dimension(:), pointer :: local_nvertices character (len=StrKIND) :: filename, msg @@ -111,28 +111,36 @@ subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, l if (dminfo % my_proc_id == IO_NODE) then if ( trim(blockFilePrefix) == '' ) then call mpas_log_write('blockFilePrefix is not set: Using LibScotch for graph partitioning') - do vind=1,partial_global_graph_info % nVerticesTotal - nTotalEdgesGraph = nTotalEdgesGraph + nEdgesOnCellFull% array(vind) - ! do j=1,nEdgesOnCellFull% array(vind) - ! call mpas_log_write('vind=$i j=$i adj= $i', intArgs=(/vind,j,cellsOnCellFull%array(j,vind)/) ) + do i=1,partial_global_graph_info % nVerticesTotal + do j=1,nEdgesOnCellFull% array(i) + + if (cellsOnCellFull%array(j,i) == 0) cycle + nTotalEdgesGraph = nTotalEdgesGraph + 1 + !nTotalEdgesGraph = nTotalEdgesGraph + nEdgesOnCellFull% array(i) + ! do j=1,nEdgesOnCellFull% array(i) + ! call mpas_log_write('i=$i j=$i adj= $i', intArgs=(/i,j,cellsOnCellFull%array(j,i)/) ) ! end do + end do end do - + call mpas_log_write('nTotalEdgesGraph is $i', intArgs=(/nTotalEdgesGraph/)) allocate(edgetab(nTotalEdgesGraph)) allocate(verttab(partial_global_graph_info % nVerticesTotal + 1)) !allocate(parttab(partial_global_graph_info % nVerticesTotal)) - !do vind=1,partial_global_graph_info % nVerticesTotal - !call mpas_log_write('proc=$i vind= $i part= $i', intArgs=(/dminfo % my_proc_id, vind,parttab(vind)/)) - !call mpas_log_write('vind=$i j=$i adj= $i', intArgs=(/vind,j,cellsOnCellFull%array(j,vind)/) ) + !do i=1,partial_global_graph_info % nVerticesTotal + !call mpas_log_write('proc=$i i= $i part= $i', intArgs=(/dminfo % my_proc_id, i,parttab(i)/)) + !call mpas_log_write('i=$i j=$i adj= $i', intArgs=(/i,j,cellsOnCellFull%array(j,i)/) ) !end do k = 1 - do vind=1,partial_global_graph_info % nVerticesTotal - verttab(vind) = k - !call mpas_log_write('vind=$i verttab= $i', intArgs=(/vind,verttab(vind)/) ) - do j=1,nEdgesOnCellFull% array(vind) - edgetab(k) = cellsOnCellFull%array(j,vind) + do i=1,partial_global_graph_info % nVerticesTotal + verttab(i) = k + !call mpas_log_write('i=$i verttab= $i', intArgs=(/i,verttab(i)/) ) + do j=1,nEdgesOnCellFull% array(i) + + if (cellsOnCellFull%array(j,i) == 0) cycle + + edgetab(k) = cellsOnCellFull%array(j,i) !call mpas_log_write('k=$i edgetab= $i', intArgs=(/k,edgetab(k)/) ) k = k + 1 end do @@ -154,6 +162,12 @@ subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, l IF (IERR .NE. 0) THEN call mpas_log_write('Cannot build Scotch Graph', MPAS_LOG_CRIT) ENDIF + + ! CALL SCOTCHFGRAPHSAVE (SCOTCHGRAPH (1), 1, IERR) + ! IF (IERR .NE. 0) THEN + ! PRINT *, 'ERROR : MAIN : Invalid graph output' + ! STOP + ! ENDIF CALL SCOTCHFGRAPHCHECK (SCOTCHGRAPH (1), IERR) IF (IERR .NE. 0) THEN From 0d626dbce2c6399cb93577d0eb7d7fe273675845 Mon Sep 17 00:00:00 2001 From: Abishek Gopal Date: Tue, 24 Jun 2025 09:35:48 -0600 Subject: [PATCH 3/7] Modifying the logic for when Scotch partitioning is active --- src/framework/mpas_block_decomp.F | 134 ++++++++++++++++++------------ 1 file changed, 81 insertions(+), 53 deletions(-) diff --git a/src/framework/mpas_block_decomp.F b/src/framework/mpas_block_decomp.F index 81da41292b..e26230f2d9 100644 --- a/src/framework/mpas_block_decomp.F +++ b/src/framework/mpas_block_decomp.F @@ -25,7 +25,10 @@ module mpas_block_decomp use mpas_derived_types use mpas_io_units use mpas_log + + #ifdef MPAS_SCOTCH include "scotchf.h" + #endif type graph integer :: nVerticesTotal @@ -75,7 +78,7 @@ subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, l integer, dimension(:,:), allocatable :: sorted_local_cell_list integer, dimension(:), allocatable :: global_block_id_arr, owning_proc_arr - integer :: i, global_block_id, local_block_id, owning_proc, iunit, istatus, j, k + integer :: i, global_block_id, local_block_id, owning_proc, iunit, ounit, istatus, ostatus, j, k integer :: blocks_per_proc, err, ierr integer, dimension(:), pointer :: local_nvertices character (len=StrKIND) :: filename, msg @@ -84,9 +87,11 @@ subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, l integer :: nTotalEdgesGraph = 0 integer, dimension(:), allocatable :: edgetab, verttab, parttab - doubleprecision :: stradat (SCOTCH_STRATDIM) - doubleprecision :: SCOTCHGRAPH (SCOTCH_GRAPHDIM) - integer :: n_size + logical :: useScotch +#ifdef MPAS_SCOTCH + doubleprecision :: stradat (scotch_stratdim) + doubleprecision :: scotchgraph (scotch_graphdim) +#endif no_blocks = .false. @@ -109,8 +114,50 @@ subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, l allocate(owning_proc_arr(partial_global_graph_info % nVerticesTotal)) if (dminfo % my_proc_id == IO_NODE) then + useScotch = .false. if ( trim(blockFilePrefix) == '' ) then - call mpas_log_write('blockFilePrefix is not set: Using LibScotch for graph partitioning') + call mpas_log_write('Namelist option config_block_decomp_file_prefix is set to \'\' ', MPAS_LOG_ERR) +#ifdef MPAS_SCOTCH + useScotch = .true. +#else + call mpas_log_write('Either build MPAS with the Scotch library or provide a valid file prefix for config_block_decomp_file_prefix', MPAS_LOG_CRIT) +#endif + else + if (dminfo % total_blocks < 10) then + write(filename,'(a,i1)') trim(blockFilePrefix), dminfo % total_blocks + else if (dminfo % total_blocks < 100) then + write(filename,'(a,i2)') trim(blockFilePrefix), dminfo % total_blocks + else if (dminfo % total_blocks < 1000) then + write(filename,'(a,i3)') trim(blockFilePrefix), dminfo % total_blocks + else if (dminfo % total_blocks < 10000) then + write(filename,'(a,i4)') trim(blockFilePrefix), dminfo % total_blocks + else if (dminfo % total_blocks < 100000) then + write(filename,'(a,i5)') trim(blockFilePrefix), dminfo % total_blocks + else if (dminfo % total_blocks < 1000000) then + write(filename,'(a,i6)') trim(blockFilePrefix), dminfo % total_blocks + else if (dminfo % total_blocks < 10000000) then + write(filename,'(a,i7)') trim(blockFilePrefix), dminfo % total_blocks + else if (dminfo % total_blocks < 100000000) then + write(filename,'(a,i8)') trim(blockFilePrefix), dminfo % total_blocks + end if + + call mpas_new_unit(iunit) + open(unit=iunit, file=trim(filename), form='formatted', status='old', iostat=istatus) + + if (istatus /= 0) then + call mpas_log_write('Could not open block decomposition file for $i blocks.', MPAS_LOG_WARN, intArgs=(/dminfo % total_blocks/) ) + call mpas_log_write('Filename: '//trim(filename),MPAS_LOG_WARN) +#ifdef MPAS_SCOTCH + useScotch = .true. +#else + call mpas_log_write('Either build MPAS with Scotch library or provide a valid file prefix for config_block_decomp_file_prefix', MPAS_LOG_CRIT) +#endif + end if + end if + + if (useScotch) then +#ifdef MPAS_SCOTCH + call mpas_log_write('Using LibScotch for graph partitioning') do i=1,partial_global_graph_info % nVerticesTotal do j=1,nEdgesOnCellFull% array(i) @@ -149,33 +196,28 @@ subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, l !call mpas_log_write('nvertices =$i nTotalEdgesGraph= $i', intArgs=(/partial_global_graph_info % nVerticesTotal, nTotalEdgesGraph/)) - CALL scotchfstratinit (stradat (1), IERR) - CALL SCOTCHFGRAPHINIT (SCOTCHGRAPH (1), IERR) - IF (IERR .NE. 0) THEN - call mpas_log_write('Cannot initialize Scotch Graph', MPAS_LOG_CRIT) - ENDIF + call scotchfstratinit (stradat (1), ierr) + call scotchfgraphinit (scotchgraph (1), ierr) + + if (ierr .ne. 0) then + call mpas_log_write('Cannot initialize Scotch Graph', MPAS_LOG_CRIT) + endif - CALL SCOTCHFGRAPHBUILD (SCOTCHGRAPH (1), 1, partial_global_graph_info % nVerticesTotal, & + CALL scotchfgraphbuild (scotchgraph (1), 1, partial_global_graph_info % nVerticesTotal, & verttab (1), verttab (2), & verttab (1), verttab (1), & - nTotalEdgesGraph, edgetab (1), edgetab (1), IERR) - IF (IERR .NE. 0) THEN - call mpas_log_write('Cannot build Scotch Graph', MPAS_LOG_CRIT) - ENDIF - - ! CALL SCOTCHFGRAPHSAVE (SCOTCHGRAPH (1), 1, IERR) - ! IF (IERR .NE. 0) THEN - ! PRINT *, 'ERROR : MAIN : Invalid graph output' - ! STOP - ! ENDIF + nTotalEdgesGraph, edgetab (1), edgetab (1), ierr) + if (ierr .ne. 0) then + call mpas_log_write('Cannot build Scotch Graph', MPAS_LOG_CRIT) + endif + + ! Only needed during development/debugging. + !call scotchfgraphcheck (scotchgraph (1), ierr) + !if (ierr .ne. 0) then + ! call mpas_log_write('Cannot check Scotch Graph', MPAS_LOG_CRIT) + !endif - CALL SCOTCHFGRAPHCHECK (SCOTCHGRAPH (1), IERR) - IF (IERR .NE. 0) THEN - call mpas_log_write('Cannot check Scotch Graph', MPAS_LOG_CRIT) - ENDIF - - !CALL scotchfgraphpart( SCOTCHGRAPH (1), dminfo % total_blocks, stradat (1) ,parttab(1), IERR) - CALL scotchfgraphpart( SCOTCHGRAPH (1), dminfo % total_blocks, stradat (1) ,global_block_id_arr(1), IERR) + call scotchfgraphpart( scotchgraph (1), dminfo % total_blocks, stradat (1) ,global_block_id_arr(1), IERR) call scotchfgraphexit (SCOTCHGRAPH (1)) call scotchfstratexit (stradat (1)) @@ -188,34 +230,20 @@ subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, l !call mpas_log_write('vert: $i, global_block_id: $i, owning_proc: $i ', intArgs=(/i,global_block_id_arr(i),owning_proc/) ) local_nvertices(owning_proc+1) = local_nvertices(owning_proc+1) + 1 end do - - else - if (dminfo % total_blocks < 10) then - write(filename,'(a,i1)') trim(blockFilePrefix), dminfo % total_blocks - else if (dminfo % total_blocks < 100) then - write(filename,'(a,i2)') trim(blockFilePrefix), dminfo % total_blocks - else if (dminfo % total_blocks < 1000) then - write(filename,'(a,i3)') trim(blockFilePrefix), dminfo % total_blocks - else if (dminfo % total_blocks < 10000) then - write(filename,'(a,i4)') trim(blockFilePrefix), dminfo % total_blocks - else if (dminfo % total_blocks < 100000) then - write(filename,'(a,i5)') trim(blockFilePrefix), dminfo % total_blocks - else if (dminfo % total_blocks < 1000000) then - write(filename,'(a,i6)') trim(blockFilePrefix), dminfo % total_blocks - else if (dminfo % total_blocks < 10000000) then - write(filename,'(a,i7)') trim(blockFilePrefix), dminfo % total_blocks - else if (dminfo % total_blocks < 100000000) then - write(filename,'(a,i8)') trim(blockFilePrefix), dminfo % total_blocks - end if - - call mpas_new_unit(iunit) - open(unit=iunit, file=trim(filename), form='formatted', status='old', iostat=istatus) - if (istatus /= 0) then - call mpas_log_write('Could not open block decomposition file for $i blocks.', MPAS_LOG_ERR, intArgs=(/dminfo % total_blocks/) ) - call mpas_log_write('Filename: '//trim(filename), MPAS_LOG_CRIT) + call mpas_log_write('Writing out Scotch Graph paritioning data to '//trim(filename)) + call mpas_new_unit(ounit) + open(unit=ounit, file=trim(filename), form='formatted', status='new', action="write", iostat=ostatus) + do i=1,partial_global_graph_info % nVerticesTotal + write(unit=ounit, fmt='(i0)', iostat=err) global_block_id_arr(i) + !write(unit=ounit, fmt='(*(i0,1x))', iostat=err) global_block_id_arr(i) + end do + close(unit=ounit) + call mpas_release_unit(ounit) end if + #endif + else call mpas_log_write('Using block decomposition file: '//trim(filename)) From e19d74d7c382104e0936c916a6a5323e273ab00d Mon Sep 17 00:00:00 2001 From: Abishek Gopal Date: Tue, 24 Jun 2025 09:36:43 -0600 Subject: [PATCH 4/7] WIP: Makefile changes to enable SCOTCH --- Makefile | 14 +++++++++----- 1 file changed, 9 insertions(+), 5 deletions(-) diff --git a/Makefile b/Makefile index 7606b47032..6fe7bbcfbe 100644 --- a/Makefile +++ b/Makefile @@ -759,11 +759,14 @@ endif LIBS += $(NCLIB) endif -export SCOTCH_ROOT=/glade/derecho/scratch/agopal/scotch/build - -FCINCLUDES += -I$(SCOTCH_ROOT)/src/include - -LIBS += -L$(SCOTCH_ROOT)/lib -lscotch -lscotcherr +ifneq "$(SCOTCH)" "" + override CPPFLAGS += "-DMPAS_SCOTCH" + FCINCLUDES += -I$(SCOTCH)/src/include + LIBS += -L$(SCOTCH)/lib -lscotch -lscotcherr + SCOTCH_MESSAGE = "MPAS has been linked with the Scotch Graph Paritioning library." +else + SCOTCH_MESSAGE = "MPAS was NOT linked with the Scotch Graph Paritioning library." +endif ifneq "$(PNETCDF)" "" ifneq ($(wildcard $(PNETCDF)/lib/libpnetcdf.*), ) @@ -1513,6 +1516,7 @@ mpas_main: $(MAIN_DEPS) @echo $(OPENMP_OFFLOAD_MESSAGE) @echo $(OPENACC_MESSAGE) @echo $(MUSICA_MESSAGE) + @echo $(SCOTCH_MESSAGE) @echo $(SHAREDLIB_MESSAGE) ifeq "$(AUTOCLEAN)" "true" @echo $(AUTOCLEAN_MESSAGE) From e75ea0654d1e18126431340485fbe83b40902a80 Mon Sep 17 00:00:00 2001 From: Abishek Gopal Date: Tue, 24 Jun 2025 16:26:27 -0600 Subject: [PATCH 5/7] Adding compile-time link test for the Scotch library --- Makefile | 47 +++++++++++++++++++++++++++++++++++++++++------ 1 file changed, 41 insertions(+), 6 deletions(-) diff --git a/Makefile b/Makefile index 6fe7bbcfbe..b7d1b9f3dd 100644 --- a/Makefile +++ b/Makefile @@ -760,12 +760,13 @@ endif endif ifneq "$(SCOTCH)" "" - override CPPFLAGS += "-DMPAS_SCOTCH" - FCINCLUDES += -I$(SCOTCH)/src/include - LIBS += -L$(SCOTCH)/lib -lscotch -lscotcherr - SCOTCH_MESSAGE = "MPAS has been linked with the Scotch Graph Paritioning library." -else - SCOTCH_MESSAGE = "MPAS was NOT linked with the Scotch Graph Paritioning library." + SCOTCH_FCINCLUDES += -I$(SCOTCH)/src/include + SCOTCH_LIBS += -L$(SCOTCH)/lib -lscotch -lscotcherr + SCOTCH_FFLAGS = -DMPAS_SCOTCH + + FCINCLUDES += $(SCOTCH_FCINCLUDES) + LIBS += $(SCOTCH_LIBS) + override CPPFLAGS += $(SCOTCH_FFLAGS) endif ifneq "$(PNETCDF)" "" @@ -1423,6 +1424,33 @@ musica_fortran_test: $(eval MUSICA_FORTRAN_VERSION := $(shell pkg-config --modversion musica-fortran)) $(if $(findstring 1,$(MUSICA_FORTRAN_TEST)), $(info Built a simple test program with MUSICA-Fortran version $(MUSICA_FORTRAN_VERSION)), ) +scotch_fortran_test: + @# + @# Create a Fortran test program that will link against the SCOTCH library + @# + $(info Checking for a working MUSICA-Fortran library...) + $(eval SCOTCH_FORTRAN_TEST := $(shell $\ + printf "program test_scotch_fortran\n$\ + & include \"scotchf.h\"\n$\ + & doubleprecision :: scotchgraph (scotch_graphdim)\n$\ + & integer :: ierr\n$\ + & ierr = 0\n$\ + & call scotchfgraphinit(scotchgraph (1), ierr)\n$\ + & call scotchfgraphexit(scotchgraph(1))\n$\ + end program test_scotch_fortran\n" | sed 's/&/ /' > test_scotch_fortran.f90; $\ + $\ + $(FC) $(SCOTCH_FCINCLUDES) $(SCOTCH_FFLAGS) test_scotch_fortran.f90 -o test_scotch_fortran.x $(SCOTCH_LIBS) > /dev/null 2>&1; $\ + scotch_fortran_status=$$?; $\ + rm -f test_scotch_fortran.f90 test_scotch_fortran.x; $\ + if [ $$scotch_fortran_status -eq 0 ]; then $\ + printf "1"; $\ + else $\ + printf "0"; $\ + fi $\ + )) + $(if $(findstring 0,$(SCOTCH_FORTRAN_TEST)), $(error Could not build a simple test program with Scotch)) + $(if $(findstring 1,$(SCOTCH_FORTRAN_TEST)), $(info Built a simple test program with Scotch )) + pnetcdf_test: @# @# Create test C programs that look for PNetCDF header file and some symbols in it @@ -1479,6 +1507,13 @@ else MUSICA_MESSAGE = "MPAS was not linked with the MUSICA-Fortran library." endif +ifneq "$(SCOTCH_FFLAGS)" "" +MAIN_DEPS += scotch_fortran_test +SCOTCH_MESSAGE = "MPAS has been linked with the Scotch graph partitioning library." +else +SCOTCH_MESSAGE = "MPAS was NOT linked with the Scotch graph partitioning library." +endif + mpas_main: $(MAIN_DEPS) cd src; $(MAKE) FC="$(FC)" \ CC="$(CC)" \ From 987ba44b82ef64365abc6839771c74953699e2df Mon Sep 17 00:00:00 2001 From: Abishek Gopal Date: Wed, 2 Jul 2025 16:51:47 -0600 Subject: [PATCH 6/7] few fixes for GNU and Intel builds --- src/framework/mpas_block_decomp.F | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/src/framework/mpas_block_decomp.F b/src/framework/mpas_block_decomp.F index e26230f2d9..67c285d14c 100644 --- a/src/framework/mpas_block_decomp.F +++ b/src/framework/mpas_block_decomp.F @@ -26,9 +26,9 @@ module mpas_block_decomp use mpas_io_units use mpas_log - #ifdef MPAS_SCOTCH - include "scotchf.h" - #endif +#ifdef MPAS_SCOTCH +#include "scotchf.h" +#endif type graph integer :: nVerticesTotal @@ -116,7 +116,7 @@ subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, l if (dminfo % my_proc_id == IO_NODE) then useScotch = .false. if ( trim(blockFilePrefix) == '' ) then - call mpas_log_write('Namelist option config_block_decomp_file_prefix is set to \'\' ', MPAS_LOG_ERR) + call mpas_log_write("Namelist option config_block_decomp_file_prefix is set to ''", MPAS_LOG_ERR) #ifdef MPAS_SCOTCH useScotch = .true. #else @@ -242,7 +242,7 @@ subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, l close(unit=ounit) call mpas_release_unit(ounit) end if - #endif +#endif else call mpas_log_write('Using block decomposition file: '//trim(filename)) From 44d3421d775cd585c7f73abae928adb69a87e7df Mon Sep 17 00:00:00 2001 From: Abishek Gopal Date: Wed, 2 Jul 2025 20:05:12 -0600 Subject: [PATCH 7/7] adding timers --- src/framework/mpas_block_decomp.F | 30 +++++++++++++++++++----------- src/framework/mpas_bootstrapping.F | 5 +++++ 2 files changed, 24 insertions(+), 11 deletions(-) diff --git a/src/framework/mpas_block_decomp.F b/src/framework/mpas_block_decomp.F index 67c285d14c..573220a9fb 100644 --- a/src/framework/mpas_block_decomp.F +++ b/src/framework/mpas_block_decomp.F @@ -55,6 +55,8 @@ module mpas_block_decomp subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, local_cell_list, block_id, block_start, & block_count, cellsOnCellFull, nEdgesOnCellFull, numBlocks, explicitProcDecomp, blockFilePrefix, procFilePrefix)!{{{ + use mpas_timer, only : mpas_timer_start, mpas_timer_stop + implicit none type (dm_info), intent(inout) :: dminfo !< Input: domain information @@ -95,6 +97,8 @@ subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, l no_blocks = .false. + call mpas_timer_start('mpas_block_decomp_cells_for_proc') + if (numBlocks == 0) then dminfo % total_blocks = dminfo % nProcs else @@ -157,6 +161,7 @@ subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, l if (useScotch) then #ifdef MPAS_SCOTCH + call mpas_timer_start('scotch_graph_partitioning') call mpas_log_write('Using LibScotch for graph partitioning') do i=1,partial_global_graph_info % nVerticesTotal do j=1,nEdgesOnCellFull% array(i) @@ -231,17 +236,18 @@ subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, l local_nvertices(owning_proc+1) = local_nvertices(owning_proc+1) + 1 end do - if (istatus /= 0) then - call mpas_log_write('Writing out Scotch Graph paritioning data to '//trim(filename)) - call mpas_new_unit(ounit) - open(unit=ounit, file=trim(filename), form='formatted', status='new', action="write", iostat=ostatus) - do i=1,partial_global_graph_info % nVerticesTotal - write(unit=ounit, fmt='(i0)', iostat=err) global_block_id_arr(i) - !write(unit=ounit, fmt='(*(i0,1x))', iostat=err) global_block_id_arr(i) - end do - close(unit=ounit) - call mpas_release_unit(ounit) - end if + ! if (istatus /= 0) then + ! call mpas_log_write('Writing out Scotch Graph paritioning data to '//trim(filename)) + ! call mpas_new_unit(ounit) + ! open(unit=ounit, file=trim(filename), form='formatted', status='new', action="write", iostat=ostatus) + ! do i=1,partial_global_graph_info % nVerticesTotal + ! write(unit=ounit, fmt='(i0)', iostat=err) global_block_id_arr(i) + ! !write(unit=ounit, fmt='(*(i0,1x))', iostat=err) global_block_id_arr(i) + ! end do + ! close(unit=ounit) + ! call mpas_release_unit(ounit) + ! end if + call mpas_timer_stop('scotch_graph_partitioning') #endif else @@ -430,6 +436,8 @@ subroutine mpas_block_decomp_cells_for_proc(dminfo, partial_global_graph_info, l end if end if + call mpas_timer_stop('mpas_block_decomp_cells_for_proc') + end subroutine mpas_block_decomp_cells_for_proc!}}} !*********************************************************************** diff --git a/src/framework/mpas_bootstrapping.F b/src/framework/mpas_bootstrapping.F index f444032e3b..909b0386ec 100644 --- a/src/framework/mpas_bootstrapping.F +++ b/src/framework/mpas_bootstrapping.F @@ -81,6 +81,7 @@ subroutine mpas_bootstrap_framework_phase1(domain, mesh_filename, mesh_iotype, p #ifdef MPAS_PIO_SUPPORT use pio, only : file_desc_t #endif + use mpas_timer, only : mpas_timer_start, mpas_timer_stop implicit none @@ -155,6 +156,7 @@ subroutine mpas_bootstrap_framework_phase1(domain, mesh_filename, mesh_iotype, p nHalos = config_num_halos + call mpas_timer_start('bootstrap_framework_phase1') inputHandle = MPAS_io_open(trim(mesh_filename), MPAS_IO_READ, mesh_iotype, domain % ioContext, & pio_file_desc=pio_file_desc, ierr=ierr) @@ -440,6 +442,9 @@ subroutine mpas_bootstrap_framework_phase1(domain, mesh_filename, mesh_iotype, p deallocate(block_count) deallocate(readingBlock) + + call mpas_timer_stop('bootstrap_framework_phase1') + end subroutine mpas_bootstrap_framework_phase1 !}}}