updates indexing in scotch partitioning for newer gcc versions (>= 13.x)

src/decompose_mesh/partition_scotch.F90

@@ -173,7 +173,7 @@ subroutine partition_scotch()
       stop 'Scotch ERROR : MAIN : Invalid check'
     endif
 
-    call scotchfgraphpart (scotchgraph(1), nparts, scotchstrat(1),part(1),ier)
+    call scotchfgraphpart (scotchgraph(1), nparts, scotchstrat(1), part(1),ier)
     if (ier /= 0) then
       stop 'Scotch ERROR : MAIN : Cannot part graph'
     endif
@@ -271,11 +271,15 @@ subroutine lts_partition_scotch()
     enddo
     if (inum /= nspec_p) stop 'Error number of p-vertices not equals to number of p-elements'
 
+    ! note: newer gcc compilers (version >= 13.x) seem to work only if the adjaceny arrays start indexing from 0.
+    !       that is, we use scotchfgraphbuild(..) with a baseval == 0 and use index values from 0 to nspec_p-1 in
+    !       xadj_tmp() and adjncy_tmp() arrays.
+
     ! builds up local arrays only for p-elements
     xadj_tmp(:) = 0
     adjncy_tmp(:) = -1
     elmnts_load_tmp(:) = 0
-    nb_edges_tmp = 1
+    nb_edges_tmp = 0
     do inum = 1,nspec_p
       ispec = ispec_global(inum)
 
@@ -284,11 +288,9 @@ subroutine lts_partition_scotch()
 
       ! number of neighbors
       ! note: xadj() values start from 0; adjncy() values start from 0 to nspec-1
-      !       we shift these by +1 since arrays in this subroutine are defined between 1 to ..
       !
       !       xadj(i)  -> adjncy( . ) : values in xadj point to the first index in adjncy() for vertex i
       !                   adjncy() holds all indices of neighbor vertices(=elements)
-      ! note: for xadj_tmp and adjncy_tmp we start indexing at 1
       xadj_tmp(inum) = nb_edges_tmp
       do j = xadj(ispec),xadj(ispec+1)-1
         ! gets neighbor index
@@ -299,14 +301,20 @@ subroutine lts_partition_scotch()
         if (ispec_p_refine(k) == p) then
           ! we store the local index between 1 and nspec_p
           i = ispec_local( k )
-          adjncy_tmp(nb_edges_tmp) = i
+          if (i < 1 .or. i > nspec_p) then
+            print *,'Error: invalid local index ',i,'ispec k',k,'p',ispec_p_refine(k),'ilevel/p',ilevel,p
+            stop 'Error local index'
+          endif
+          adjncy_tmp(nb_edges_tmp+1) = i - 1  ! shift index to start at 0
           nb_edges_tmp = nb_edges_tmp + 1
         endif
       enddo
     enddo
     ! last entry for xadj for a contiguous range of indices ("compact edge array")
     xadj_tmp(nspec_p+1) = nb_edges_tmp
-    nb_edges_tmp = nb_edges_tmp - 1
+
+    ! since we start counting at 0, the total number is nb_edges_tmp not nb_edges_tmp-1
+    !nb_edges_tmp = nb_edges_tmp - 1
 
     ! debug
     !print *,'xadj:',xadj(ispec_global(1):ispec_global(1)+1),xadj(nspec-1:nspec+1)
@@ -324,8 +332,8 @@ subroutine lts_partition_scotch()
     !  print *,j-xadj_tmp(1)+1,j,adjncy_tmp(j),ispec_global(adjncy_tmp(j))
     !enddo
 
-    ! checks ranges
-    if (minval(adjncy_tmp(1:nb_edges_tmp)) < 1 .or. maxval(adjncy_tmp(1:nb_edges_tmp)) > nspec_p) then
+    ! checks ranges (index values in adjncy_tmp() start at 0)
+    if (minval(adjncy_tmp(1:nb_edges_tmp)) < 0 .or. maxval(adjncy_tmp(1:nb_edges_tmp)) > nspec_p - 1) then
       print *,'Error adjncy bounds invalid', minval(adjncy_tmp(1:nb_edges_tmp)),maxval(adjncy_tmp(1:nb_edges_tmp))
       stop 'Error adjncy bounds invalid'
     endif
@@ -349,7 +357,7 @@ subroutine lts_partition_scotch()
     ! Setting vendtab to refer to one cell after verttab yields the same result,
     ! as it is the exact semantics of a compact vertex array.
 
-    call scotchfgraphbuild (scotchgraph(1), 1, nspec_p, &
+    call scotchfgraphbuild (scotchgraph(1), 0, nspec_p, &
                             xadj_tmp(1), xadj_tmp(1), &
                             elmnts_load_tmp(1), xadj_tmp(1), &
                             nb_edges_tmp, adjncy_tmp(1), &
@@ -372,19 +380,22 @@ subroutine lts_partition_scotch()
         'using subset maximum',P_LEVEL_PARTIAL_SUBSET_MINIMUM
 
       ! partitions among subset nparts_partial processes
-      call scotchfgraphpart (scotchgraph(1), nparts_partial, scotchstrat(1),part_tmp(1),ier)
+      call scotchfgraphpart (scotchgraph(1), nparts_partial, scotchstrat(1), part_tmp(1), ier)
       if (ier /= 0) stop 'Scotch ERROR : MAIN : Cannot part graph'
 
     else
+      ! parition over all nparts processes
+      ! user output
+      print *,'    partitioning: nparts =',nparts
 
       ! partitions among all nparts processes
-      call scotchfgraphpart (scotchgraph(1), nparts, scotchstrat(1),part_tmp(1),ier)
+      call scotchfgraphpart (scotchgraph(1), nparts, scotchstrat(1), part_tmp(1), ier)
       if (ier /= 0) stop 'Scotch ERROR : MAIN : Cannot part graph'
 
     endif
 
     ! frees scotch graph for subsequent calls of scotch again
-    call scotchfgraphfree (scotchgraph(1),ier)
+    call scotchfgraphfree (scotchgraph(1), ier)
     if (ier /= 0) stop 'Scotch ERROR : MAIN : Cannot free graph'
 
     ! stitch partitioning together