SALSA Analysis Modules
|
Functions for computing the JPL multicolour structure of a matrix. More...
#include <stdlib.h>
#include "anamod.h"
#include "anamodsalsamodules.h"
#include "anamatrix.h"
#include "petscerror.h"
#include "petscmat.h"
#include "petscis.h"
#include "src/mat/impls/aij/mpi/mpiaij.h"
#include "src/mat/impls/aij/seq/aij.h"
Go to the source code of this file.
Functions | |
static PetscErrorCode | LookAtUncolouredVar (int compare, int this_var, Mat A, PetscScalar *clr_array, PetscScalar *ran_array, PetscReal this_rand, int *neighb, int *nneigh, PetscTruth *colour_now) |
static PetscErrorCode | JonesPlassmannColouring (Mat A, PetscTruth *flg) |
static PetscErrorCode | NColours (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | ColourSizes (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | ColourOffsets (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | Colours (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
PetscErrorCode | RegisterJPLModules (void) |
Functions for computing the JPL multicolour structure of a matrix.
Activate this module with
PetscErrorCode RegisterJPLModules();
Compute these elements with
ComputeQuantity("jpl","element",A,(void*)&res,&flg);
Available elements are:
This routine works in parallel; the stored result is distributed over all participating processors.
@article{JoPl:heuristic, author = {M.T. Jones and P.E. Plassmann}, title = {A parallel graph coloring heuristic}, journal = SJSSC, volume = {14}, year = {1993}, keywords = {incomplete factorization, multicolouring, matrix graph} } @inproceedings{JoPl:parallelunstructured, author = {M.T. Jones and P.E. Plassmann}, title = {Parallel solution of unstructed, sparse systems of linear equations}, booktitle = {Proceedings of the Sixth SIAM conference on Parallel Processing for Scientific Computing}, editor = {R.F. Sincovec and D.E. Keyes and M.R. Leuze and L.R. Petzold and D.A. Reed}, publisher = {SIAM}, address = {Philadelphia}, pages = {471--475}, year = {1993} \section Disclaimer Unabashed use of British spelling of `colour' needs no disclaimer. }
Definition in file jpl.c.
static PetscErrorCode ColourOffsets | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscTruth * | flg | ||
) | [static] |
Compute an array that has the offsets of the locations of the colours, on this processor. If computing the colouring is done in parallel, then the multicolouring is stored in parallel: it is not gathered onto any processor.
Definition at line 451 of file jpl.c.
References GetDataID(), HASTOEXIST, id, AnalysisItem::ii, and JonesPlassmannColouring().
Referenced by RegisterJPLModules().
{ Mat A = (Mat)prob; int llv=0,*v = NULL; PetscTruth has; int id,id2; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("jpl","n-colours",&id2,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id2,llv,*flg); CHKERRQ(ierr); if (*flg) { if (lv) *lv = llv; ierr = GetDataID("jpl","colour-offsets",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetIntstar ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); } if (*flg) rv->ii = v; else { ierr = JonesPlassmannColouring(A,flg); CHKERRQ(ierr); if (*flg) { ierr = PetscObjectComposedDataGetIntstar ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->ii = v; ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id2,llv,*flg); CHKERRQ(ierr); if (*flg) if (lv) *lv = llv; } } PetscFunctionReturn(0); }
static PetscErrorCode Colours | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscTruth * | flg | ||
) | [static] |
Compute an array that has the colour sets on this processor. In order to know where a set starts or ends, the result of ColourOffsets() is needed.
If computing the colouring is done in parallel, then the multicolouring is stored in parallel: it is not gathered onto any processor.
Definition at line 497 of file jpl.c.
References GetDataID(), HASTOEXIST, id, AnalysisItem::ii, and JonesPlassmannColouring().
Referenced by RegisterJPLModules().
{ Mat A = (Mat)prob; int N; int *v = NULL; PetscTruth has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = MatGetLocalSize(A,&N,PETSC_NULL); CHKERRQ(ierr); if (lv) *lv = N; ierr = GetDataID("jpl","colours",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetIntstar ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->ii = v; else { ierr = JonesPlassmannColouring(A,flg); CHKERRQ(ierr); if (*flg) { ierr = PetscObjectComposedDataGetIntstar ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->ii = v; } } PetscFunctionReturn(0); }
static PetscErrorCode ColourSizes | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscTruth * | flg | ||
) | [static] |
Compute an array that has the number of points of each colour, on this processor. If computing the colouring is done in parallel, then the multicolouring is stored in parallel: it is not gathered onto any processor.
Definition at line 408 of file jpl.c.
References GetDataID(), HASTOEXIST, id, AnalysisItem::ii, and JonesPlassmannColouring().
Referenced by RegisterJPLModules().
{ Mat A = (Mat)prob; int llv=0,*v = NULL; PetscTruth has; int id,id2; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("jpl","n-colours",&id2,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id2,llv,*flg); CHKERRQ(ierr); if (*flg) { if (lv) *lv = llv; ierr = GetDataID("jpl","colour-set-sizes",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetIntstar ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); } if (*flg) rv->ii = v; else { ierr = JonesPlassmannColouring(A,flg); CHKERRQ(ierr); if (*flg) { ierr = PetscObjectComposedDataGetIntstar ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->ii = v; ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id2,llv,*flg); CHKERRQ(ierr); if (*flg) if (lv) *lv = llv; } } PetscFunctionReturn(0); }
static PetscErrorCode JonesPlassmannColouring | ( | Mat | A, |
PetscTruth * | flg | ||
) | [static] |
Compute a multicolour ordering of a matrix, in parallel.
This is the main computational routine.
Definition at line 129 of file jpl.c.
References GetDataID(), id, LookAtUncolouredVar(), and MatrixComputeQuantity().
Referenced by ColourOffsets(), Colours(), ColourSizes(), and NColours().
{ MPI_Comm comm; const MatType type; Mat dia,off; Vec rand,ran_bord, clrs,clrs_bkp,clr_bord; VecScatter sctr; PetscScalar *ran_array,*bran_array,*clr_array,*clr_bkp_array,*bclr_array; int *neighb,ncoloured,total_coloured,pass,max_clr, global_maxclr,local_size,total_size, id; PetscTruth parallel,sequential; PetscErrorCode ierr; PetscFunctionBegin; *flg = PETSC_FALSE; /* * What is the matrix type? We only handle AIJ (nonblocked) types */ ierr = PetscObjectGetComm((PetscObject)A,&comm); CHKERRQ(ierr); ierr = MatGetType(A,&type); CHKERRQ(ierr); ierr = PetscStrcmp(type,MATSEQAIJ,&sequential); CHKERRQ(ierr); ierr = PetscStrcmp(type,MATMPIAIJ,¶llel); CHKERRQ(ierr); if (!(sequential || parallel)) { *flg = PETSC_FALSE; PetscFunctionReturn(0); } /* * Set up the matrices and vectors */ if (parallel) { Mat_MPIAIJ *aij = (Mat_MPIAIJ *) A->data; dia=aij->A; off=aij->B; clr_bord = aij->lvec; sctr = aij->Mvctx; ierr = VecDuplicate(aij->lvec,&ran_bord); CHKERRQ(ierr); } else { dia = A; off = NULL; clr_bord = NULL; sctr = NULL; } ierr = MatGetSize(A,&total_size,PETSC_NULL); CHKERRQ(ierr); ierr = MatGetLocalSize(A,&local_size,PETSC_NULL); CHKERRQ(ierr); /* A temporary array for storing neighbour information */ { int l1,l2,len; PetscTruth f; ierr = MatrixComputeQuantity (dia,"structure","max-nnzeros-per-row",(AnalysisItem*)&l1,PETSC_NULL, &f); CHKERRQ(ierr); if (parallel) { ierr = MatrixComputeQuantity (off,"structure","max-nnzeros-per-row",(AnalysisItem*)&l2,PETSC_NULL, &f); CHKERRQ(ierr); if (f) len=l1+l2; else len = local_size; } else len = l1; ierr = PetscMalloc((len+1)*sizeof(int),&neighb); CHKERRQ(ierr); } /* * Create the vectors for random and colours */ { if (parallel) { ierr = VecCreateMPI(comm,local_size,PETSC_DECIDE,&rand); CHKERRQ(ierr); } else { ierr = VecCreateSeq(MPI_COMM_SELF,local_size,&rand); CHKERRQ(ierr); } { PetscRandom rctx; ierr = PetscRandomCreate(MPI_COMM_SELF,&rctx); CHKERRQ(ierr); ierr = PetscRandomSetType(rctx,PETSCRAND48); CHKERRQ(ierr); ierr = VecSetRandom(rand,rctx); CHKERRQ(ierr); ierr = PetscRandomDestroy(rctx); CHKERRQ(ierr); } ierr = VecDuplicate(rand,&clrs); CHKERRQ(ierr); ierr = VecDuplicate(rand,&clrs_bkp); CHKERRQ(ierr); ierr = VecSet(clrs,0.); CHKERRQ(ierr); if (parallel) { ierr = VecScatterBegin (sctr,rand,ran_bord,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr); ierr = VecScatterEnd (sctr,rand,ran_bord,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr); ierr = VecGetArray(ran_bord,&bran_array); CHKERRQ(ierr); ierr = VecGetArray(clr_bord,&bclr_array); CHKERRQ(ierr); } ierr = VecGetArray(rand,&ran_array); CHKERRQ(ierr); ierr = VecGetArray(clrs,&clr_array); CHKERRQ(ierr); ierr = VecGetArray(clrs_bkp,&clr_bkp_array); CHKERRQ(ierr); } /* * Here is the main loop. * Iterate until everything is coloured. */ ncoloured = 0; max_clr = 0; for (pass=0; ; pass++) { PetscTruth busy = PETSC_FALSE,gbusy; int var; /*printf(".. pass %d :",pass);*/ if (max_clr>pass) SETERRQ2(1,"This can not happen: maxclr=%d in pass %d",max_clr,pass); ierr = PetscMemcpy (clr_bkp_array,clr_array,local_size*sizeof(PetscScalar)); CHKERRQ(ierr); /* * Loop over all variables in every pass */ for (var=0; var<local_size; var++) { int nn = 0; PetscTruth now; /* Skip already coloured variables */ if (clr_array[var]>0) continue; /* All others, look at all connections */ ierr = LookAtUncolouredVar (var,var,dia,clr_bkp_array,ran_array,ran_array[var], neighb,&nn,&now); CHKERRQ(ierr); if (now && parallel) { if (parallel) { ierr = LookAtUncolouredVar (-1,var,off,bclr_array,bran_array,ran_array[var], neighb,&nn,&now); CHKERRQ(ierr); } } /* If this variable has the highest value among uncoloured neighbours, * colour it now. */ if (now) { int clr=1,i; attempt: for (i=0; i<nn; i++) if (neighb[i]==clr) {clr++; goto attempt;} clr_array[var] = (PetscScalar) clr; if (clr>max_clr) max_clr = clr; ncoloured++; busy = PETSC_TRUE; /*printf(" %d->%d",var,clr);*/ } else { /*printf(" %d.",var);*/ } } /*printf("\n");*/ /* Exit if all points are coloured */ MPI_Allreduce ((void*)&ncoloured,(void*)&total_coloured,1,MPI_INT,MPI_SUM,comm); if (total_coloured==total_size) break; /* If nothing was coloured in this pass, we're stuck */ MPI_Allreduce ((void*)&busy,(void*)&gbusy,1,MPI_INT,MPI_MAX,comm); if (!gbusy) SETERRQ(1,"We are stuck"); if (parallel) { ierr = VecScatterBegin (sctr,clrs,clr_bord,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr); ierr = VecScatterEnd (sctr,clrs,clr_bord,INSERT_VALUES,SCATTER_FORWARD); CHKERRQ(ierr); } } /* * Store the number of colours */ MPI_Allreduce ((void*)&max_clr,(void*)&global_maxclr,1,MPI_INT,MPI_MAX,comm); ierr = GetDataID("jpl","n-colours",&id,PETSC_NULL); CHKERRQ(ierr); ierr = PetscObjectComposedDataSetInt ((PetscObject)A,id,global_maxclr); CHKERRQ(ierr); /* * Count how many points of each colour */ { int *colour_count,*colour_offsets,*ptr,*colours,clr,i; ierr = PetscMalloc ((global_maxclr+1)*sizeof(int),&colour_count); CHKERRQ(ierr); colour_count[0] = global_maxclr; for (clr=0; clr<global_maxclr; clr++) colour_count[clr+1] = 0; for (i=0; i<local_size; i++) colour_count[(int)(clr_array[i])]++; /* colours are from 1 */ /* printf("coloured:"); for (i=0; i<local_size; i++) printf(" %e",clr_array[i]); printf("\n"); printf("counts:"); for (i=0; i<colour_count[0]; i++) printf(" %d",colour_count[i+1]); printf("\n"); */ ierr = GetDataID("jpl","colour-set-sizes",&id,PETSC_NULL); CHKERRQ(ierr); ierr = PetscObjectComposedDataSetIntstar ((PetscObject)A,id,colour_count); CHKERRQ(ierr); ierr = PetscMalloc ((global_maxclr+1)*sizeof(int),&colour_offsets); CHKERRQ(ierr); colour_offsets[0] = global_maxclr; colour_offsets[1] = 0; for (clr=1; clr<global_maxclr; clr++) colour_offsets[1+clr] = colour_offsets[clr]+colour_count[clr]; /* printf("offsets:"); for (i=0; i<colour_offsets[0]; i++) printf(" %d",colour_offsets[i+1]); printf("\n"); */ ierr = GetDataID("jpl","colour-offsets",&id,PETSC_NULL); CHKERRQ(ierr); ierr = PetscObjectComposedDataSetIntstar ((PetscObject)A,id,colour_offsets); CHKERRQ(ierr); ierr = PetscMalloc((local_size+1)*sizeof(int),&colours); CHKERRQ(ierr); colours[0] = local_size; ierr = PetscMalloc(global_maxclr*sizeof(int),&ptr); CHKERRQ(ierr); PetscMemcpy(ptr,colour_offsets+1,global_maxclr*sizeof(int)); for (i=0; i<local_size; i++) { int c=(int)(clr_array[i]); colours[1+ptr[c-1]++] = i; } ierr = PetscFree(ptr); CHKERRQ(ierr); /* printf("colours:"); for (i=0; i<colours[0]; i++) printf(" %d",colours[i+1]); printf("\n"); */ ierr = GetDataID("jpl","colours",&id,PETSC_NULL); CHKERRQ(ierr); ierr = PetscObjectComposedDataSetIntstar ((PetscObject)A,id,colours); CHKERRQ(ierr); } /* Clean up */ ierr = PetscFree(neighb); CHKERRQ(ierr); ierr = VecRestoreArray(rand,&ran_array); CHKERRQ(ierr); ierr = VecDestroy(rand); CHKERRQ(ierr); ierr = VecRestoreArray(clrs,&clr_array); CHKERRQ(ierr); ierr = VecDestroy(clrs); CHKERRQ(ierr); ierr = VecRestoreArray(clrs_bkp,&clr_bkp_array); CHKERRQ(ierr); ierr = VecDestroy(clrs_bkp); CHKERRQ(ierr); if (parallel) { ierr = VecRestoreArray(ran_bord,&bran_array); CHKERRQ(ierr); ierr = VecDestroy(ran_bord); CHKERRQ(ierr); ierr = VecRestoreArray(clr_bord,&bclr_array); CHKERRQ(ierr); /* clr_bord is the aij->lvec array; not to be destroyed */ } *flg = PETSC_TRUE; PetscFunctionReturn(0); }
static PetscErrorCode LookAtUncolouredVar | ( | int | compare, |
int | this_var, | ||
Mat | A, | ||
PetscScalar * | clr_array, | ||
PetscScalar * | ran_array, | ||
PetscReal | this_rand, | ||
int * | neighb, | ||
int * | nneigh, | ||
PetscTruth * | colour_now | ||
) | [static] |
Investigate the local and remote connections of a variable to see whether it needs to be coloured, and if so, with what colour.
Definition at line 77 of file jpl.c.
Referenced by JonesPlassmannColouring().
{ int nnz,j; const int *idx; PetscErrorCode ierr; PetscFunctionBegin; *colour_now = PETSC_TRUE; ierr = MatGetRow(A,this_var,&nnz,&idx,PETSC_NULL); CHKERRQ(ierr); for (j=0; j<nnz; j++) { /* loop over all points connected to 'this_var' */ int other,clr; other=idx[j]; clr = (int)(clr_array[other]); /* no need to look at ourselves */ if (other==compare) continue; /* this really shouldn't happen */ /* if (ran_array[other]==this_rand) SETERRQ3(1,"Use a better random: %e @ %d,%d",this_rand,compare,other); */ /* * If the other variable is uncoloured and has a higher random value, * we can not colour the current variable. */ if ( (clr==0 && ran_array[other]>this_rand) || (ran_array[other]==this_rand && other>this_var) ) { *colour_now = PETSC_FALSE; goto exit; } /* * Otherwise, we take note of the colour of the other variable * for the computation of the colour of this variable later */ neighb[(*nneigh)++] = clr; } exit: ierr = MatRestoreRow(A,this_var,&nnz,&idx,PETSC_NULL); CHKERRQ(ierr); PetscFunctionReturn(0); }
static PetscErrorCode NColours | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscTruth * | flg | ||
) | [static] |
Compute the global number of colours. Note that any given processor does not need to have all the colours.
Definition at line 377 of file jpl.c.
References GetDataID(), HASTOEXIST, AnalysisItem::i, id, and JonesPlassmannColouring().
Referenced by RegisterJPLModules().
{ Mat A = (Mat)prob; int v = 0; PetscTruth has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("jpl","n-colours",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->i = v; else { ierr = JonesPlassmannColouring(A,flg); CHKERRQ(ierr); if (*flg) { ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->i = v; } } PetscFunctionReturn(0); }
PetscErrorCode RegisterJPLModules | ( | void | ) |
Definition at line 523 of file jpl.c.
References ANALYSISINTARRAY, ANALYSISINTEGER, ColourOffsets(), Colours(), ColourSizes(), NColours(), and RegisterModule().
Referenced by AnaModRegisterSalsaModules(), and AnaModRegisterStandardModules().
{ PetscErrorCode ierr; PetscFunctionBegin; ierr = RegisterModule ("jpl","n-colours",ANALYSISINTEGER,&NColours); CHKERRQ(ierr); ierr = RegisterModule ("jpl","colour-set-sizes",ANALYSISINTARRAY,&ColourSizes); CHKERRQ(ierr); ierr = RegisterModule ("jpl","colour-offsets",ANALYSISINTARRAY,&ColourOffsets); CHKERRQ(ierr); ierr = RegisterModule ("jpl","colours",ANALYSISINTARRAY,&Colours); CHKERRQ(ierr); PetscFunctionReturn(0); }