1 /* 2 GAMG geometric-algebric multigrid PC - Mark Adams 2011 3 */ 4 #include <petsc/private/matimpl.h> 5 #include <../src/ksp/pc/impls/gamg/gamg.h> /*I "petscpc.h" I*/ 6 #include <../src/ksp/pc/impls/bjacobi/bjacobi.h> /* Hack to access same_local_solves */ 7 8 #if defined PETSC_GAMG_USE_LOG 9 PetscLogEvent petsc_gamg_setup_events[NUM_SET]; 10 #endif 11 12 #if defined PETSC_USE_LOG 13 PetscLogEvent PC_GAMGGraph_AGG; 14 PetscLogEvent PC_GAMGGraph_GEO; 15 PetscLogEvent PC_GAMGCoarsen_AGG; 16 PetscLogEvent PC_GAMGCoarsen_GEO; 17 PetscLogEvent PC_GAMGProlongator_AGG; 18 PetscLogEvent PC_GAMGProlongator_GEO; 19 PetscLogEvent PC_GAMGOptProlongator_AGG; 20 #endif 21 22 #define GAMG_MAXLEVELS 30 23 24 /* #define GAMG_STAGES */ 25 #if (defined PETSC_GAMG_USE_LOG && defined GAMG_STAGES) 26 static PetscLogStage gamg_stages[GAMG_MAXLEVELS]; 27 #endif 28 29 static PetscFunctionList GAMGList = 0; 30 static PetscBool PCGAMGPackageInitialized; 31 32 /* ----------------------------------------------------------------------------- */ 33 #undef __FUNCT__ 34 #define __FUNCT__ "PCReset_GAMG" 35 PetscErrorCode PCReset_GAMG(PC pc) 36 { 37 PetscErrorCode ierr; 38 PC_MG *mg = (PC_MG*)pc->data; 39 PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx; 40 41 PetscFunctionBegin; 42 if (pc_gamg->data) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_PLIB,"This should not happen, cleaned up in SetUp\n"); 43 pc_gamg->data_sz = 0; 44 ierr = PetscFree(pc_gamg->orig_data);CHKERRQ(ierr); 45 PetscFunctionReturn(0); 46 } 47 48 /* -------------------------------------------------------------------------- */ 49 /* 50 PCGAMGCreateLevel_GAMG: create coarse op with RAP. repartition and/or reduce number 51 of active processors. 52 53 Input Parameter: 54 . pc - parameters + side effect: coarse data in 'pc_gamg->data' and 55 'pc_gamg->data_sz' are changed via repartitioning/reduction. 56 . Amat_fine - matrix on this fine (k) level 57 . cr_bs - coarse block size 58 In/Output Parameter: 59 . a_P_inout - prolongation operator to the next level (k-->k-1) 60 . a_nactive_proc - number of active procs 61 Output Parameter: 62 . a_Amat_crs - coarse matrix that is created (k-1) 63 */ 64 65 #undef __FUNCT__ 66 #define __FUNCT__ "PCGAMGCreateLevel_GAMG" 67 static PetscErrorCode PCGAMGCreateLevel_GAMG(PC pc,Mat Amat_fine,PetscInt cr_bs,Mat *a_P_inout,Mat *a_Amat_crs,PetscMPIInt *a_nactive_proc,IS * Pcolumnperm) 68 { 69 PetscErrorCode ierr; 70 PC_MG *mg = (PC_MG*)pc->data; 71 PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx; 72 Mat Cmat,Pold=*a_P_inout; 73 MPI_Comm comm; 74 PetscMPIInt rank,size,new_size,nactive=*a_nactive_proc; 75 PetscInt ncrs_eq,ncrs,f_bs; 76 77 PetscFunctionBegin; 78 ierr = PetscObjectGetComm((PetscObject)Amat_fine,&comm);CHKERRQ(ierr); 79 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 80 ierr = MPI_Comm_size(comm, &size);CHKERRQ(ierr); 81 ierr = MatGetBlockSize(Amat_fine, &f_bs);CHKERRQ(ierr); 82 ierr = MatPtAP(Amat_fine, Pold, MAT_INITIAL_MATRIX, 2.0, &Cmat);CHKERRQ(ierr); 83 84 /* set 'ncrs' (nodes), 'ncrs_eq' (equations)*/ 85 ierr = MatGetLocalSize(Cmat, &ncrs_eq, NULL);CHKERRQ(ierr); 86 if (pc_gamg->data_cell_rows>0) { 87 ncrs = pc_gamg->data_sz/pc_gamg->data_cell_cols/pc_gamg->data_cell_rows; 88 } else { 89 PetscInt bs; 90 ierr = MatGetBlockSize(Cmat, &bs);CHKERRQ(ierr); 91 ncrs = ncrs_eq/bs; 92 } 93 94 /* get number of PEs to make active 'new_size', reduce, can be any integer 1-P */ 95 { 96 PetscInt ncrs_eq_glob; 97 ierr = MatGetSize(Cmat, &ncrs_eq_glob, NULL);CHKERRQ(ierr); 98 new_size = (PetscMPIInt)((float)ncrs_eq_glob/(float)pc_gamg->min_eq_proc + 0.5); /* hardwire min. number of eq/proc */ 99 if (!new_size) new_size = 1; /* not likely, posible? */ 100 else if (new_size >= nactive) new_size = nactive; /* no change, rare */ 101 } 102 103 if (Pcolumnperm) *Pcolumnperm = NULL; 104 105 if (!pc_gamg->repart && new_size==nactive) *a_Amat_crs = Cmat; /* output - no repartitioning or reduction - could bail here */ 106 else { 107 PetscInt *counts,*newproc_idx,ii,jj,kk,strideNew,*tidx,ncrs_new,ncrs_eq_new,nloc_old; 108 IS is_eq_newproc,is_eq_num,is_eq_num_prim,new_eq_indices; 109 110 nloc_old = ncrs_eq/cr_bs; 111 if (ncrs_eq % cr_bs) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"ncrs_eq %D not divisible by cr_bs %D",ncrs_eq,cr_bs); 112 #if defined PETSC_GAMG_USE_LOG 113 ierr = PetscLogEventBegin(petsc_gamg_setup_events[SET12],0,0,0,0);CHKERRQ(ierr); 114 #endif 115 /* make 'is_eq_newproc' */ 116 ierr = PetscMalloc1(size, &counts);CHKERRQ(ierr); 117 if (pc_gamg->repart) { 118 /* Repartition Cmat_{k} and move colums of P^{k}_{k-1} and coordinates of primal part accordingly */ 119 Mat adj; 120 121 ierr = PetscInfo3(pc,"Repartition: size (active): %D --> %D, neq = %D\n",*a_nactive_proc,new_size,ncrs_eq);CHKERRQ(ierr); 122 123 /* get 'adj' */ 124 if (cr_bs == 1) { 125 ierr = MatConvert(Cmat, MATMPIADJ, MAT_INITIAL_MATRIX, &adj);CHKERRQ(ierr); 126 } else { 127 /* make a scalar matrix to partition (no Stokes here) */ 128 Mat tMat; 129 PetscInt Istart_crs,Iend_crs,ncols,jj,Ii; 130 const PetscScalar *vals; 131 const PetscInt *idx; 132 PetscInt *d_nnz, *o_nnz, M, N; 133 static PetscInt llev = 0; 134 MatType mtype; 135 136 ierr = PetscMalloc2(ncrs, &d_nnz,ncrs, &o_nnz);CHKERRQ(ierr); 137 ierr = MatGetOwnershipRange(Cmat, &Istart_crs, &Iend_crs);CHKERRQ(ierr); 138 ierr = MatGetSize(Cmat, &M, &N);CHKERRQ(ierr); 139 for (Ii = Istart_crs, jj = 0; Ii < Iend_crs; Ii += cr_bs, jj++) { 140 ierr = MatGetRow(Cmat,Ii,&ncols,0,0);CHKERRQ(ierr); 141 d_nnz[jj] = ncols/cr_bs; 142 o_nnz[jj] = ncols/cr_bs; 143 ierr = MatRestoreRow(Cmat,Ii,&ncols,0,0);CHKERRQ(ierr); 144 if (d_nnz[jj] > ncrs) d_nnz[jj] = ncrs; 145 if (o_nnz[jj] > (M/cr_bs-ncrs)) o_nnz[jj] = M/cr_bs-ncrs; 146 } 147 148 ierr = MatGetType(Amat_fine,&mtype);CHKERRQ(ierr); 149 ierr = MatCreate(comm, &tMat);CHKERRQ(ierr); 150 ierr = MatSetSizes(tMat, ncrs, ncrs,PETSC_DETERMINE, PETSC_DETERMINE);CHKERRQ(ierr); 151 ierr = MatSetType(tMat,mtype);CHKERRQ(ierr); 152 ierr = MatSeqAIJSetPreallocation(tMat,0,d_nnz);CHKERRQ(ierr); 153 ierr = MatMPIAIJSetPreallocation(tMat,0,d_nnz,0,o_nnz);CHKERRQ(ierr); 154 ierr = PetscFree2(d_nnz,o_nnz);CHKERRQ(ierr); 155 156 for (ii = Istart_crs; ii < Iend_crs; ii++) { 157 PetscInt dest_row = ii/cr_bs; 158 ierr = MatGetRow(Cmat,ii,&ncols,&idx,&vals);CHKERRQ(ierr); 159 for (jj = 0; jj < ncols; jj++) { 160 PetscInt dest_col = idx[jj]/cr_bs; 161 PetscScalar v = 1.0; 162 ierr = MatSetValues(tMat,1,&dest_row,1,&dest_col,&v,ADD_VALUES);CHKERRQ(ierr); 163 } 164 ierr = MatRestoreRow(Cmat,ii,&ncols,&idx,&vals);CHKERRQ(ierr); 165 } 166 ierr = MatAssemblyBegin(tMat,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 167 ierr = MatAssemblyEnd(tMat,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 168 169 if (llev++ == -1) { 170 PetscViewer viewer; char fname[32]; 171 ierr = PetscSNPrintf(fname,sizeof(fname),"part_mat_%D.mat",llev);CHKERRQ(ierr); 172 PetscViewerBinaryOpen(comm,fname,FILE_MODE_WRITE,&viewer); 173 ierr = MatView(tMat, viewer);CHKERRQ(ierr); 174 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 175 } 176 ierr = MatConvert(tMat, MATMPIADJ, MAT_INITIAL_MATRIX, &adj);CHKERRQ(ierr); 177 ierr = MatDestroy(&tMat);CHKERRQ(ierr); 178 } /* create 'adj' */ 179 180 { /* partition: get newproc_idx */ 181 char prefix[256]; 182 const char *pcpre; 183 const PetscInt *is_idx; 184 MatPartitioning mpart; 185 IS proc_is; 186 PetscInt targetPE; 187 188 ierr = MatPartitioningCreate(comm, &mpart);CHKERRQ(ierr); 189 ierr = MatPartitioningSetAdjacency(mpart, adj);CHKERRQ(ierr); 190 ierr = PCGetOptionsPrefix(pc, &pcpre);CHKERRQ(ierr); 191 ierr = PetscSNPrintf(prefix,sizeof(prefix),"%spc_gamg_",pcpre ? pcpre : "");CHKERRQ(ierr); 192 ierr = PetscObjectSetOptionsPrefix((PetscObject)mpart,prefix);CHKERRQ(ierr); 193 ierr = MatPartitioningSetFromOptions(mpart);CHKERRQ(ierr); 194 ierr = MatPartitioningSetNParts(mpart, new_size);CHKERRQ(ierr); 195 ierr = MatPartitioningApply(mpart, &proc_is);CHKERRQ(ierr); 196 ierr = MatPartitioningDestroy(&mpart);CHKERRQ(ierr); 197 198 /* collect IS info */ 199 ierr = PetscMalloc1(ncrs_eq, &newproc_idx);CHKERRQ(ierr); 200 ierr = ISGetIndices(proc_is, &is_idx);CHKERRQ(ierr); 201 targetPE = 1; /* bring to "front" of machine */ 202 /*targetPE = size/new_size;*/ /* spread partitioning across machine */ 203 for (kk = jj = 0 ; kk < nloc_old ; kk++) { 204 for (ii = 0 ; ii < cr_bs ; ii++, jj++) { 205 newproc_idx[jj] = is_idx[kk] * targetPE; /* distribution */ 206 } 207 } 208 ierr = ISRestoreIndices(proc_is, &is_idx);CHKERRQ(ierr); 209 ierr = ISDestroy(&proc_is);CHKERRQ(ierr); 210 } 211 ierr = MatDestroy(&adj);CHKERRQ(ierr); 212 213 ierr = ISCreateGeneral(comm, ncrs_eq, newproc_idx, PETSC_COPY_VALUES, &is_eq_newproc);CHKERRQ(ierr); 214 ierr = PetscFree(newproc_idx);CHKERRQ(ierr); 215 } else { /* simple aggreagtion of parts -- 'is_eq_newproc' */ 216 PetscInt rfactor,targetPE; 217 218 /* find factor */ 219 if (new_size == 1) rfactor = size; /* easy */ 220 else { 221 PetscReal best_fact = 0.; 222 jj = -1; 223 for (kk = 1 ; kk <= size ; kk++) { 224 if (!(size%kk)) { /* a candidate */ 225 PetscReal nactpe = (PetscReal)size/(PetscReal)kk, fact = nactpe/(PetscReal)new_size; 226 if (fact > 1.0) fact = 1./fact; /* keep fact < 1 */ 227 if (fact > best_fact) { 228 best_fact = fact; jj = kk; 229 } 230 } 231 } 232 if (jj != -1) rfactor = jj; 233 else rfactor = 1; /* does this happen .. a prime */ 234 } 235 new_size = size/rfactor; 236 237 if (new_size==nactive) { 238 *a_Amat_crs = Cmat; /* output - no repartitioning or reduction, bail out because nested here */ 239 ierr = PetscFree(counts);CHKERRQ(ierr); 240 ierr = PetscInfo2(pc,"Aggregate processors noop: new_size=%D, neq(loc)=%D\n",new_size,ncrs_eq);CHKERRQ(ierr); 241 #if defined PETSC_GAMG_USE_LOG 242 ierr = PetscLogEventEnd(petsc_gamg_setup_events[SET12],0,0,0,0);CHKERRQ(ierr); 243 #endif 244 PetscFunctionReturn(0); 245 } 246 247 ierr = PetscInfo1(pc,"Number of equations (loc) %D with simple aggregation\n",ncrs_eq);CHKERRQ(ierr); 248 targetPE = rank/rfactor; 249 ierr = ISCreateStride(comm, ncrs_eq, targetPE, 0, &is_eq_newproc);CHKERRQ(ierr); 250 } /* end simple 'is_eq_newproc' */ 251 252 /* 253 Create an index set from the is_eq_newproc index set to indicate the mapping TO 254 */ 255 ierr = ISPartitioningToNumbering(is_eq_newproc, &is_eq_num);CHKERRQ(ierr); 256 is_eq_num_prim = is_eq_num; 257 /* 258 Determine how many equations/vertices are assigned to each processor 259 */ 260 ierr = ISPartitioningCount(is_eq_newproc, size, counts);CHKERRQ(ierr); 261 ncrs_eq_new = counts[rank]; 262 ierr = ISDestroy(&is_eq_newproc);CHKERRQ(ierr); 263 ncrs_new = ncrs_eq_new/cr_bs; /* eqs */ 264 265 ierr = PetscFree(counts);CHKERRQ(ierr); 266 #if defined PETSC_GAMG_USE_LOG 267 ierr = PetscLogEventEnd(petsc_gamg_setup_events[SET12],0,0,0,0);CHKERRQ(ierr); 268 #endif 269 /* data movement scope -- this could be moved to subclasses so that we don't try to cram all auxilary data into some complex abstracted thing */ 270 { 271 Vec src_crd, dest_crd; 272 const PetscInt *idx,ndata_rows=pc_gamg->data_cell_rows,ndata_cols=pc_gamg->data_cell_cols,node_data_sz=ndata_rows*ndata_cols; 273 VecScatter vecscat; 274 PetscScalar *array; 275 IS isscat; 276 277 /* move data (for primal equations only) */ 278 /* Create a vector to contain the newly ordered element information */ 279 ierr = VecCreate(comm, &dest_crd);CHKERRQ(ierr); 280 ierr = VecSetSizes(dest_crd, node_data_sz*ncrs_new, PETSC_DECIDE);CHKERRQ(ierr); 281 ierr = VecSetType(dest_crd,VECSTANDARD);CHKERRQ(ierr); /* this is needed! */ 282 /* 283 There are 'ndata_rows*ndata_cols' data items per node, (one can think of the vectors of having 284 a block size of ...). Note, ISs are expanded into equation space by 'cr_bs'. 285 */ 286 ierr = PetscMalloc1(ncrs*node_data_sz, &tidx);CHKERRQ(ierr); 287 ierr = ISGetIndices(is_eq_num_prim, &idx);CHKERRQ(ierr); 288 for (ii=0,jj=0; ii<ncrs; ii++) { 289 PetscInt id = idx[ii*cr_bs]/cr_bs; /* get node back */ 290 for (kk=0; kk<node_data_sz; kk++, jj++) tidx[jj] = id*node_data_sz + kk; 291 } 292 ierr = ISRestoreIndices(is_eq_num_prim, &idx);CHKERRQ(ierr); 293 ierr = ISCreateGeneral(comm, node_data_sz*ncrs, tidx, PETSC_COPY_VALUES, &isscat);CHKERRQ(ierr); 294 ierr = PetscFree(tidx);CHKERRQ(ierr); 295 /* 296 Create a vector to contain the original vertex information for each element 297 */ 298 ierr = VecCreateSeq(PETSC_COMM_SELF, node_data_sz*ncrs, &src_crd);CHKERRQ(ierr); 299 for (jj=0; jj<ndata_cols; jj++) { 300 const PetscInt stride0=ncrs*pc_gamg->data_cell_rows; 301 for (ii=0; ii<ncrs; ii++) { 302 for (kk=0; kk<ndata_rows; kk++) { 303 PetscInt ix = ii*ndata_rows + kk + jj*stride0, jx = ii*node_data_sz + kk*ndata_cols + jj; 304 PetscScalar tt = (PetscScalar)pc_gamg->data[ix]; 305 ierr = VecSetValues(src_crd, 1, &jx, &tt, INSERT_VALUES);CHKERRQ(ierr); 306 } 307 } 308 } 309 ierr = VecAssemblyBegin(src_crd);CHKERRQ(ierr); 310 ierr = VecAssemblyEnd(src_crd);CHKERRQ(ierr); 311 /* 312 Scatter the element vertex information (still in the original vertex ordering) 313 to the correct processor 314 */ 315 ierr = VecScatterCreate(src_crd, NULL, dest_crd, isscat, &vecscat);CHKERRQ(ierr); 316 ierr = ISDestroy(&isscat);CHKERRQ(ierr); 317 ierr = VecScatterBegin(vecscat,src_crd,dest_crd,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 318 ierr = VecScatterEnd(vecscat,src_crd,dest_crd,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 319 ierr = VecScatterDestroy(&vecscat);CHKERRQ(ierr); 320 ierr = VecDestroy(&src_crd);CHKERRQ(ierr); 321 /* 322 Put the element vertex data into a new allocation of the gdata->ele 323 */ 324 ierr = PetscFree(pc_gamg->data);CHKERRQ(ierr); 325 ierr = PetscMalloc1(node_data_sz*ncrs_new, &pc_gamg->data);CHKERRQ(ierr); 326 327 pc_gamg->data_sz = node_data_sz*ncrs_new; 328 strideNew = ncrs_new*ndata_rows; 329 330 ierr = VecGetArray(dest_crd, &array);CHKERRQ(ierr); 331 for (jj=0; jj<ndata_cols; jj++) { 332 for (ii=0; ii<ncrs_new; ii++) { 333 for (kk=0; kk<ndata_rows; kk++) { 334 PetscInt ix = ii*ndata_rows + kk + jj*strideNew, jx = ii*node_data_sz + kk*ndata_cols + jj; 335 pc_gamg->data[ix] = PetscRealPart(array[jx]); 336 } 337 } 338 } 339 ierr = VecRestoreArray(dest_crd, &array);CHKERRQ(ierr); 340 ierr = VecDestroy(&dest_crd);CHKERRQ(ierr); 341 } 342 /* move A and P (columns) with new layout */ 343 #if defined PETSC_GAMG_USE_LOG 344 ierr = PetscLogEventBegin(petsc_gamg_setup_events[SET13],0,0,0,0);CHKERRQ(ierr); 345 #endif 346 347 /* 348 Invert for MatGetSubMatrix 349 */ 350 ierr = ISInvertPermutation(is_eq_num, ncrs_eq_new, &new_eq_indices);CHKERRQ(ierr); 351 ierr = ISSort(new_eq_indices);CHKERRQ(ierr); /* is this needed? */ 352 ierr = ISSetBlockSize(new_eq_indices, cr_bs);CHKERRQ(ierr); 353 if (is_eq_num != is_eq_num_prim) { 354 ierr = ISDestroy(&is_eq_num_prim);CHKERRQ(ierr); /* could be same as 'is_eq_num' */ 355 } 356 if (Pcolumnperm) { 357 ierr = PetscObjectReference((PetscObject)new_eq_indices);CHKERRQ(ierr); 358 *Pcolumnperm = new_eq_indices; 359 } 360 ierr = ISDestroy(&is_eq_num);CHKERRQ(ierr); 361 #if defined PETSC_GAMG_USE_LOG 362 ierr = PetscLogEventEnd(petsc_gamg_setup_events[SET13],0,0,0,0);CHKERRQ(ierr); 363 ierr = PetscLogEventBegin(petsc_gamg_setup_events[SET14],0,0,0,0);CHKERRQ(ierr); 364 #endif 365 /* 'a_Amat_crs' output */ 366 { 367 Mat mat; 368 ierr = MatGetSubMatrix(Cmat, new_eq_indices, new_eq_indices, MAT_INITIAL_MATRIX, &mat);CHKERRQ(ierr); 369 *a_Amat_crs = mat; 370 } 371 ierr = MatDestroy(&Cmat);CHKERRQ(ierr); 372 373 #if defined PETSC_GAMG_USE_LOG 374 ierr = PetscLogEventEnd(petsc_gamg_setup_events[SET14],0,0,0,0);CHKERRQ(ierr); 375 #endif 376 /* prolongator */ 377 { 378 IS findices; 379 PetscInt Istart,Iend; 380 Mat Pnew; 381 382 ierr = MatGetOwnershipRange(Pold, &Istart, &Iend);CHKERRQ(ierr); 383 #if defined PETSC_GAMG_USE_LOG 384 ierr = PetscLogEventBegin(petsc_gamg_setup_events[SET15],0,0,0,0);CHKERRQ(ierr); 385 #endif 386 ierr = ISCreateStride(comm,Iend-Istart,Istart,1,&findices);CHKERRQ(ierr); 387 ierr = ISSetBlockSize(findices,f_bs);CHKERRQ(ierr); 388 ierr = MatGetSubMatrix(Pold, findices, new_eq_indices, MAT_INITIAL_MATRIX, &Pnew);CHKERRQ(ierr); 389 ierr = ISDestroy(&findices);CHKERRQ(ierr); 390 391 #if defined PETSC_GAMG_USE_LOG 392 ierr = PetscLogEventEnd(petsc_gamg_setup_events[SET15],0,0,0,0);CHKERRQ(ierr); 393 #endif 394 ierr = MatDestroy(a_P_inout);CHKERRQ(ierr); 395 396 /* output - repartitioned */ 397 *a_P_inout = Pnew; 398 } 399 ierr = ISDestroy(&new_eq_indices);CHKERRQ(ierr); 400 401 *a_nactive_proc = new_size; /* output */ 402 } 403 PetscFunctionReturn(0); 404 } 405 406 /* -------------------------------------------------------------------------- */ 407 /* 408 PCSetUp_GAMG - Prepares for the use of the GAMG preconditioner 409 by setting data structures and options. 410 411 Input Parameter: 412 . pc - the preconditioner context 413 414 */ 415 #undef __FUNCT__ 416 #define __FUNCT__ "PCSetUp_GAMG" 417 PetscErrorCode PCSetUp_GAMG(PC pc) 418 { 419 PetscErrorCode ierr; 420 PC_MG *mg = (PC_MG*)pc->data; 421 PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx; 422 Mat Pmat = pc->pmat; 423 PetscInt fine_level,level,level1,bs,M,qq,lidx,nASMBlocksArr[GAMG_MAXLEVELS]; 424 MPI_Comm comm; 425 PetscMPIInt rank,size,nactivepe; 426 Mat Aarr[GAMG_MAXLEVELS],Parr[GAMG_MAXLEVELS]; 427 IS *ASMLocalIDsArr[GAMG_MAXLEVELS]; 428 PetscLogDouble nnz0=0.,nnztot=0.; 429 MatInfo info; 430 431 PetscFunctionBegin; 432 ierr = PetscObjectGetComm((PetscObject)pc,&comm);CHKERRQ(ierr); 433 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 434 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 435 436 if (pc_gamg->setup_count++ > 0) { 437 if ((PetscBool)(!pc_gamg->reuse_prol)) { 438 /* reset everything */ 439 ierr = PCReset_MG(pc);CHKERRQ(ierr); 440 pc->setupcalled = 0; 441 } else { 442 PC_MG_Levels **mglevels = mg->levels; 443 /* just do Galerkin grids */ 444 Mat B,dA,dB; 445 446 if (!pc->setupcalled) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PCSetUp() has not been called yet"); 447 if (pc_gamg->Nlevels > 1) { 448 /* currently only handle case where mat and pmat are the same on coarser levels */ 449 ierr = KSPGetOperators(mglevels[pc_gamg->Nlevels-1]->smoothd,&dA,&dB);CHKERRQ(ierr); 450 /* (re)set to get dirty flag */ 451 ierr = KSPSetOperators(mglevels[pc_gamg->Nlevels-1]->smoothd,dA,dB);CHKERRQ(ierr); 452 453 for (level=pc_gamg->Nlevels-2; level>=0; level--) { 454 /* the first time through the matrix structure has changed from repartitioning */ 455 if (pc_gamg->setup_count==2) { 456 ierr = MatPtAP(dB,mglevels[level+1]->interpolate,MAT_INITIAL_MATRIX,1.0,&B);CHKERRQ(ierr); 457 ierr = MatDestroy(&mglevels[level]->A);CHKERRQ(ierr); 458 459 mglevels[level]->A = B; 460 } else { 461 ierr = KSPGetOperators(mglevels[level]->smoothd,NULL,&B);CHKERRQ(ierr); 462 ierr = MatPtAP(dB,mglevels[level+1]->interpolate,MAT_REUSE_MATRIX,1.0,&B);CHKERRQ(ierr); 463 } 464 ierr = KSPSetOperators(mglevels[level]->smoothd,B,B);CHKERRQ(ierr); 465 dB = B; 466 } 467 } 468 469 ierr = PCSetUp_MG(pc);CHKERRQ(ierr); 470 PetscFunctionReturn(0); 471 } 472 } 473 474 if (!pc_gamg->data) { 475 if (pc_gamg->orig_data) { 476 ierr = MatGetBlockSize(Pmat, &bs);CHKERRQ(ierr); 477 ierr = MatGetLocalSize(Pmat, &qq, NULL);CHKERRQ(ierr); 478 479 pc_gamg->data_sz = (qq/bs)*pc_gamg->orig_data_cell_rows*pc_gamg->orig_data_cell_cols; 480 pc_gamg->data_cell_rows = pc_gamg->orig_data_cell_rows; 481 pc_gamg->data_cell_cols = pc_gamg->orig_data_cell_cols; 482 483 ierr = PetscMalloc1(pc_gamg->data_sz, &pc_gamg->data);CHKERRQ(ierr); 484 for (qq=0; qq<pc_gamg->data_sz; qq++) pc_gamg->data[qq] = pc_gamg->orig_data[qq]; 485 } else { 486 if (!pc_gamg->ops->createdefaultdata) SETERRQ(comm,PETSC_ERR_PLIB,"'createdefaultdata' not set(?) need to support NULL data"); 487 ierr = pc_gamg->ops->createdefaultdata(pc,Pmat);CHKERRQ(ierr); 488 } 489 } 490 491 /* cache original data for reuse */ 492 if (!pc_gamg->orig_data && (PetscBool)(!pc_gamg->reuse_prol)) { 493 ierr = PetscMalloc1(pc_gamg->data_sz, &pc_gamg->orig_data);CHKERRQ(ierr); 494 for (qq=0; qq<pc_gamg->data_sz; qq++) pc_gamg->orig_data[qq] = pc_gamg->data[qq]; 495 pc_gamg->orig_data_cell_rows = pc_gamg->data_cell_rows; 496 pc_gamg->orig_data_cell_cols = pc_gamg->data_cell_cols; 497 } 498 499 /* get basic dims */ 500 ierr = MatGetBlockSize(Pmat, &bs);CHKERRQ(ierr); 501 ierr = MatGetSize(Pmat, &M, &qq);CHKERRQ(ierr); 502 503 ierr = MatGetInfo(Pmat,MAT_GLOBAL_SUM,&info);CHKERRQ(ierr); /* global reduction */ 504 nnz0 = info.nz_used; 505 nnztot = info.nz_used; 506 ierr = PetscInfo6(pc,"level %d) N=%D, n data rows=%d, n data cols=%d, nnz/row (ave)=%d, np=%d\n",0,M,pc_gamg->data_cell_rows,pc_gamg->data_cell_cols,(int)(nnz0/(PetscReal)M+0.5),size);CHKERRQ(ierr); 507 508 /* Get A_i and R_i */ 509 for (level=0, Aarr[0]=Pmat, nactivepe = size; level < (pc_gamg->Nlevels-1) && (!level || M>pc_gamg->coarse_eq_limit); level++) { 510 pc_gamg->current_level = level; 511 level1 = level + 1; 512 #if defined PETSC_GAMG_USE_LOG 513 ierr = PetscLogEventBegin(petsc_gamg_setup_events[SET1],0,0,0,0);CHKERRQ(ierr); 514 #if (defined GAMG_STAGES) 515 ierr = PetscLogStagePush(gamg_stages[level]);CHKERRQ(ierr); 516 #endif 517 #endif 518 { /* construct prolongator */ 519 Mat Gmat; 520 PetscCoarsenData *agg_lists; 521 Mat Prol11; 522 523 ierr = pc_gamg->ops->graph(pc,Aarr[level], &Gmat);CHKERRQ(ierr); 524 ierr = pc_gamg->ops->coarsen(pc, &Gmat, &agg_lists);CHKERRQ(ierr); 525 ierr = pc_gamg->ops->prolongator(pc,Aarr[level],Gmat,agg_lists,&Prol11);CHKERRQ(ierr); 526 527 /* could have failed to create new level */ 528 if (Prol11) { 529 /* get new block size of coarse matrices */ 530 ierr = MatGetBlockSizes(Prol11, NULL, &bs);CHKERRQ(ierr); 531 532 if (pc_gamg->ops->optprolongator) { 533 /* smooth */ 534 ierr = pc_gamg->ops->optprolongator(pc, Aarr[level], &Prol11);CHKERRQ(ierr); 535 } 536 537 Parr[level1] = Prol11; 538 } else Parr[level1] = NULL; 539 540 if (pc_gamg->use_aggs_in_asm) { 541 PetscInt bs; 542 ierr = MatGetBlockSizes(Prol11, &bs, NULL);CHKERRQ(ierr); 543 ierr = PetscCDGetASMBlocks(agg_lists, bs, Gmat, &nASMBlocksArr[level], &ASMLocalIDsArr[level]);CHKERRQ(ierr); 544 } 545 546 ierr = MatDestroy(&Gmat);CHKERRQ(ierr); 547 ierr = PetscCDDestroy(agg_lists);CHKERRQ(ierr); 548 } /* construct prolongator scope */ 549 #if defined PETSC_GAMG_USE_LOG 550 ierr = PetscLogEventEnd(petsc_gamg_setup_events[SET1],0,0,0,0);CHKERRQ(ierr); 551 #endif 552 if (!level) Aarr[0] = Pmat; /* use Pmat for finest level setup */ 553 if (!Parr[level1]) { 554 ierr = PetscInfo1(pc,"Stop gridding, level %D\n",level);CHKERRQ(ierr); 555 #if defined PETSC_GAMG_USE_LOG && defined GAMG_STAGES 556 ierr = PetscLogStagePop();CHKERRQ(ierr); 557 #endif 558 break; 559 } 560 #if defined PETSC_GAMG_USE_LOG 561 ierr = PetscLogEventBegin(petsc_gamg_setup_events[SET2],0,0,0,0);CHKERRQ(ierr); 562 #endif 563 564 ierr = pc_gamg->ops->createlevel(pc, Aarr[level], bs,&Parr[level1], &Aarr[level1], &nactivepe, NULL);CHKERRQ(ierr); 565 566 #if defined PETSC_GAMG_USE_LOG 567 ierr = PetscLogEventEnd(petsc_gamg_setup_events[SET2],0,0,0,0);CHKERRQ(ierr); 568 #endif 569 ierr = MatGetSize(Aarr[level1], &M, &qq);CHKERRQ(ierr); 570 ierr = MatGetInfo(Aarr[level1], MAT_GLOBAL_SUM, &info);CHKERRQ(ierr); 571 nnztot += info.nz_used; 572 ierr = PetscInfo5(pc,"%d) N=%D, n data cols=%d, nnz/row (ave)=%d, %d active pes\n",level1,M,pc_gamg->data_cell_cols,(int)(info.nz_used/(PetscReal)M),nactivepe);CHKERRQ(ierr); 573 574 #if (defined PETSC_GAMG_USE_LOG && defined GAMG_STAGES) 575 ierr = PetscLogStagePop();CHKERRQ(ierr); 576 #endif 577 /* stop if one node or one proc -- could pull back for singular problems */ 578 if ( (pc_gamg->data_cell_cols && M/pc_gamg->data_cell_cols < 2) || (!pc_gamg->data_cell_cols && M/bs < 2) ) { 579 ierr = PetscInfo2(pc,"HARD stop of coarsening on level %D. Grid too small: %D block nodes\n",level,M/bs);CHKERRQ(ierr); 580 level++; 581 break; 582 } 583 } /* levels */ 584 ierr = PetscFree(pc_gamg->data);CHKERRQ(ierr); 585 586 ierr = PetscInfo2(pc,"%D levels, grid complexity = %g\n",level+1,nnztot/nnz0);CHKERRQ(ierr); 587 pc_gamg->Nlevels = level + 1; 588 fine_level = level; 589 ierr = PCMGSetLevels(pc,pc_gamg->Nlevels,NULL);CHKERRQ(ierr); 590 591 if (pc_gamg->Nlevels > 1) { /* don't setup MG if one level */ 592 /* set default smoothers & set operators */ 593 for (lidx = 1, level = pc_gamg->Nlevels-2; lidx <= fine_level; lidx++, level--) { 594 KSP smoother; 595 PC subpc; 596 597 ierr = PCMGGetSmoother(pc, lidx, &smoother);CHKERRQ(ierr); 598 ierr = KSPGetPC(smoother, &subpc);CHKERRQ(ierr); 599 600 ierr = KSPSetNormType(smoother, KSP_NORM_NONE);CHKERRQ(ierr); 601 /* set ops */ 602 ierr = KSPSetOperators(smoother, Aarr[level], Aarr[level]);CHKERRQ(ierr); 603 ierr = PCMGSetInterpolation(pc, lidx, Parr[level+1]);CHKERRQ(ierr); 604 605 /* set defaults */ 606 ierr = KSPSetType(smoother, KSPCHEBYSHEV);CHKERRQ(ierr); 607 608 /* set blocks for ASM smoother that uses the 'aggregates' */ 609 if (pc_gamg->use_aggs_in_asm) { 610 PetscInt sz; 611 IS *iss; 612 613 sz = nASMBlocksArr[level]; 614 iss = ASMLocalIDsArr[level]; 615 ierr = PCSetType(subpc, PCASM);CHKERRQ(ierr); 616 ierr = PCASMSetOverlap(subpc, 0);CHKERRQ(ierr); 617 ierr = PCASMSetType(subpc,PC_ASM_BASIC);CHKERRQ(ierr); 618 if (!sz) { 619 IS is; 620 ierr = ISCreateGeneral(PETSC_COMM_SELF, 0, NULL, PETSC_COPY_VALUES, &is);CHKERRQ(ierr); 621 ierr = PCASMSetLocalSubdomains(subpc, 1, NULL, &is);CHKERRQ(ierr); 622 ierr = ISDestroy(&is);CHKERRQ(ierr); 623 } else { 624 PetscInt kk; 625 ierr = PCASMSetLocalSubdomains(subpc, sz, NULL, iss);CHKERRQ(ierr); 626 for (kk=0; kk<sz; kk++) { 627 ierr = ISDestroy(&iss[kk]);CHKERRQ(ierr); 628 } 629 ierr = PetscFree(iss);CHKERRQ(ierr); 630 } 631 ASMLocalIDsArr[level] = NULL; 632 nASMBlocksArr[level] = 0; 633 } else { 634 ierr = PCSetType(subpc, PCSOR);CHKERRQ(ierr); 635 } 636 } 637 { 638 /* coarse grid */ 639 KSP smoother,*k2; PC subpc,pc2; PetscInt ii,first; 640 Mat Lmat = Aarr[(level=pc_gamg->Nlevels-1)]; lidx = 0; 641 ierr = PCMGGetSmoother(pc, lidx, &smoother);CHKERRQ(ierr); 642 ierr = KSPSetOperators(smoother, Lmat, Lmat);CHKERRQ(ierr); 643 ierr = KSPSetNormType(smoother, KSP_NORM_NONE);CHKERRQ(ierr); 644 ierr = KSPGetPC(smoother, &subpc);CHKERRQ(ierr); 645 ierr = PCSetType(subpc, PCBJACOBI);CHKERRQ(ierr); 646 ierr = PCSetUp(subpc);CHKERRQ(ierr); 647 ierr = PCBJacobiGetSubKSP(subpc,&ii,&first,&k2);CHKERRQ(ierr); 648 if (ii != 1) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_PLIB,"ii %D is not one",ii); 649 ierr = KSPGetPC(k2[0],&pc2);CHKERRQ(ierr); 650 ierr = PCSetType(pc2, PCLU);CHKERRQ(ierr); 651 ierr = PCFactorSetShiftType(pc2,MAT_SHIFT_INBLOCKS);CHKERRQ(ierr); 652 ierr = KSPSetTolerances(k2[0],PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,1);CHKERRQ(ierr); 653 ierr = KSPSetType(k2[0], KSPPREONLY);CHKERRQ(ierr); 654 /* This flag gets reset by PCBJacobiGetSubKSP(), but our BJacobi really does the same algorithm everywhere (and in 655 * fact, all but one process will have zero dofs), so we reset the flag to avoid having PCView_BJacobi attempt to 656 * view every subdomain as though they were different. */ 657 ((PC_BJacobi*)subpc->data)->same_local_solves = PETSC_TRUE; 658 } 659 660 /* should be called in PCSetFromOptions_GAMG(), but cannot be called prior to PCMGSetLevels() */ 661 ierr = PetscObjectOptionsBegin((PetscObject)pc);CHKERRQ(ierr); 662 ierr = PCSetFromOptions_MG(PetscOptionsObject,pc);CHKERRQ(ierr); 663 ierr = PetscOptionsEnd();CHKERRQ(ierr); 664 if (!mg->galerkin) SETERRQ(comm,PETSC_ERR_USER,"PCGAMG must use Galerkin for coarse operators."); 665 if (mg->galerkin == 1) mg->galerkin = 2; 666 667 /* clean up */ 668 for (level=1; level<pc_gamg->Nlevels; level++) { 669 ierr = MatDestroy(&Parr[level]);CHKERRQ(ierr); 670 ierr = MatDestroy(&Aarr[level]);CHKERRQ(ierr); 671 } 672 ierr = PCSetUp_MG(pc);CHKERRQ(ierr); 673 } else { 674 KSP smoother; 675 ierr = PetscInfo(pc,"One level solver used (system is seen as DD). Using default solver.\n");CHKERRQ(ierr); 676 ierr = PCMGGetSmoother(pc, 0, &smoother);CHKERRQ(ierr); 677 ierr = KSPSetOperators(smoother, Aarr[0], Aarr[0]);CHKERRQ(ierr); 678 ierr = KSPSetType(smoother, KSPPREONLY);CHKERRQ(ierr); 679 ierr = PCSetUp_MG(pc);CHKERRQ(ierr); 680 } 681 PetscFunctionReturn(0); 682 } 683 684 /* ------------------------------------------------------------------------- */ 685 /* 686 PCDestroy_GAMG - Destroys the private context for the GAMG preconditioner 687 that was created with PCCreate_GAMG(). 688 689 Input Parameter: 690 . pc - the preconditioner context 691 692 Application Interface Routine: PCDestroy() 693 */ 694 #undef __FUNCT__ 695 #define __FUNCT__ "PCDestroy_GAMG" 696 PetscErrorCode PCDestroy_GAMG(PC pc) 697 { 698 PetscErrorCode ierr; 699 PC_MG *mg = (PC_MG*)pc->data; 700 PC_GAMG *pc_gamg= (PC_GAMG*)mg->innerctx; 701 702 PetscFunctionBegin; 703 ierr = PCReset_GAMG(pc);CHKERRQ(ierr); 704 if (pc_gamg->ops->destroy) { 705 ierr = (*pc_gamg->ops->destroy)(pc);CHKERRQ(ierr); 706 } 707 ierr = PetscRandomDestroy(&pc_gamg->random);CHKERRQ(ierr); 708 ierr = PetscFree(pc_gamg->ops);CHKERRQ(ierr); 709 ierr = PetscFree(pc_gamg->gamg_type_name);CHKERRQ(ierr); 710 ierr = PetscFree(pc_gamg);CHKERRQ(ierr); 711 ierr = PCDestroy_MG(pc);CHKERRQ(ierr); 712 PetscFunctionReturn(0); 713 } 714 715 #undef __FUNCT__ 716 #define __FUNCT__ "PCGAMGSetProcEqLim" 717 /*@ 718 PCGAMGSetProcEqLim - Set number of equations to aim for per process on the coarse grids via processor reduction. 719 720 Logically Collective on PC 721 722 Input Parameters: 723 + pc - the preconditioner context 724 - n - the number of equations 725 726 727 Options Database Key: 728 . -pc_gamg_process_eq_limit <limit> 729 730 Notes: GAMG will reduce the number of MPI processes used directly on the coarse grids so that there are around <limit> equations on each process 731 that has degrees of freedom 732 733 Level: intermediate 734 735 Concepts: Unstructured multigrid preconditioner 736 737 .seealso: PCGAMGSetCoarseEqLim() 738 @*/ 739 PetscErrorCode PCGAMGSetProcEqLim(PC pc, PetscInt n) 740 { 741 PetscErrorCode ierr; 742 743 PetscFunctionBegin; 744 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 745 ierr = PetscTryMethod(pc,"PCGAMGSetProcEqLim_C",(PC,PetscInt),(pc,n));CHKERRQ(ierr); 746 PetscFunctionReturn(0); 747 } 748 749 #undef __FUNCT__ 750 #define __FUNCT__ "PCGAMGSetProcEqLim_GAMG" 751 static PetscErrorCode PCGAMGSetProcEqLim_GAMG(PC pc, PetscInt n) 752 { 753 PC_MG *mg = (PC_MG*)pc->data; 754 PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx; 755 756 PetscFunctionBegin; 757 if (n>0) pc_gamg->min_eq_proc = n; 758 PetscFunctionReturn(0); 759 } 760 761 #undef __FUNCT__ 762 #define __FUNCT__ "PCGAMGSetCoarseEqLim" 763 /*@ 764 PCGAMGSetCoarseEqLim - Set maximum number of equations on coarsest grid. 765 766 Collective on PC 767 768 Input Parameters: 769 + pc - the preconditioner context 770 - n - maximum number of equations to aim for 771 772 Options Database Key: 773 . -pc_gamg_coarse_eq_limit <limit> 774 775 Level: intermediate 776 777 Concepts: Unstructured multigrid preconditioner 778 779 .seealso: PCGAMGSetProcEqLim() 780 @*/ 781 PetscErrorCode PCGAMGSetCoarseEqLim(PC pc, PetscInt n) 782 { 783 PetscErrorCode ierr; 784 785 PetscFunctionBegin; 786 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 787 ierr = PetscTryMethod(pc,"PCGAMGSetCoarseEqLim_C",(PC,PetscInt),(pc,n));CHKERRQ(ierr); 788 PetscFunctionReturn(0); 789 } 790 791 #undef __FUNCT__ 792 #define __FUNCT__ "PCGAMGSetCoarseEqLim_GAMG" 793 static PetscErrorCode PCGAMGSetCoarseEqLim_GAMG(PC pc, PetscInt n) 794 { 795 PC_MG *mg = (PC_MG*)pc->data; 796 PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx; 797 798 PetscFunctionBegin; 799 if (n>0) pc_gamg->coarse_eq_limit = n; 800 PetscFunctionReturn(0); 801 } 802 803 #undef __FUNCT__ 804 #define __FUNCT__ "PCGAMGSetRepartition" 805 /*@ 806 PCGAMGSetRepartition - Repartition the degrees of freedom across the processors on the coarser grids 807 808 Collective on PC 809 810 Input Parameters: 811 + pc - the preconditioner context 812 - n - PETSC_TRUE or PETSC_FALSE 813 814 Options Database Key: 815 . -pc_gamg_repartition <true,false> 816 817 Notes: this will generally improve the loading balancing of the work on each level 818 819 Level: intermediate 820 821 Concepts: Unstructured multigrid preconditioner 822 823 .seealso: () 824 @*/ 825 PetscErrorCode PCGAMGSetRepartition(PC pc, PetscBool n) 826 { 827 PetscErrorCode ierr; 828 829 PetscFunctionBegin; 830 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 831 ierr = PetscTryMethod(pc,"PCGAMGSetRepartition_C",(PC,PetscBool),(pc,n));CHKERRQ(ierr); 832 PetscFunctionReturn(0); 833 } 834 835 #undef __FUNCT__ 836 #define __FUNCT__ "PCGAMGSetRepartition_GAMG" 837 static PetscErrorCode PCGAMGSetRepartition_GAMG(PC pc, PetscBool n) 838 { 839 PC_MG *mg = (PC_MG*)pc->data; 840 PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx; 841 842 PetscFunctionBegin; 843 pc_gamg->repart = n; 844 PetscFunctionReturn(0); 845 } 846 847 #undef __FUNCT__ 848 #define __FUNCT__ "PCGAMGSetReuseInterpolation" 849 /*@ 850 PCGAMGSetReuseInterpolation - Reuse prolongation when rebuilding algebraic multigrid preconditioner 851 852 Collective on PC 853 854 Input Parameters: 855 + pc - the preconditioner context 856 - n - PETSC_TRUE or PETSC_FALSE 857 858 Options Database Key: 859 . -pc_gamg_reuse_interpolation <true,false> 860 861 Level: intermediate 862 863 Notes: this may negatively affect the convergence rate of the method on new matrices if the matrix entries change a great deal, but allows 864 rebuilding the preconditioner quicker. 865 866 Concepts: Unstructured multigrid preconditioner 867 868 .seealso: () 869 @*/ 870 PetscErrorCode PCGAMGSetReuseInterpolation(PC pc, PetscBool n) 871 { 872 PetscErrorCode ierr; 873 874 PetscFunctionBegin; 875 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 876 ierr = PetscTryMethod(pc,"PCGAMGSetReuseInterpolation_C",(PC,PetscBool),(pc,n));CHKERRQ(ierr); 877 PetscFunctionReturn(0); 878 } 879 880 #undef __FUNCT__ 881 #define __FUNCT__ "PCGAMGSetReuseInterpolation_GAMG" 882 static PetscErrorCode PCGAMGSetReuseInterpolation_GAMG(PC pc, PetscBool n) 883 { 884 PC_MG *mg = (PC_MG*)pc->data; 885 PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx; 886 887 PetscFunctionBegin; 888 pc_gamg->reuse_prol = n; 889 PetscFunctionReturn(0); 890 } 891 892 #undef __FUNCT__ 893 #define __FUNCT__ "PCGAMGASMSetUseAggs" 894 /*@ 895 PCGAMGASMSetUseAggs - Have the PCGAMG smoother on each level use the aggregates defined by the coarsening process as the subdomains for the additive Schwarz preconditioner. 896 897 Collective on PC 898 899 Input Parameters: 900 + pc - the preconditioner context 901 - flg - PETSC_TRUE to use aggregates, PETSC_FALSE to not 902 903 Options Database Key: 904 . -pc_gamg_asm_use_agg 905 906 Level: intermediate 907 908 Concepts: Unstructured multigrid preconditioner 909 910 .seealso: () 911 @*/ 912 PetscErrorCode PCGAMGASMSetUseAggs(PC pc, PetscBool flg) 913 { 914 PetscErrorCode ierr; 915 916 PetscFunctionBegin; 917 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 918 ierr = PetscTryMethod(pc,"PCGAMGASMSetUseAggs_C",(PC,PetscBool),(pc,flg));CHKERRQ(ierr); 919 PetscFunctionReturn(0); 920 } 921 922 #undef __FUNCT__ 923 #define __FUNCT__ "PCGAMGASMSetUseAggs_GAMG" 924 static PetscErrorCode PCGAMGASMSetUseAggs_GAMG(PC pc, PetscBool flg) 925 { 926 PC_MG *mg = (PC_MG*)pc->data; 927 PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx; 928 929 PetscFunctionBegin; 930 pc_gamg->use_aggs_in_asm = flg; 931 PetscFunctionReturn(0); 932 } 933 934 #undef __FUNCT__ 935 #define __FUNCT__ "PCGAMGSetNlevels" 936 /*@ 937 PCGAMGSetNlevels - Sets the maximum number of levels PCGAMG will use 938 939 Not collective on PC 940 941 Input Parameters: 942 + pc - the preconditioner 943 - n - the maximum number of levels to use 944 945 Options Database Key: 946 . -pc_mg_levels 947 948 Level: intermediate 949 950 Concepts: Unstructured multigrid preconditioner 951 952 .seealso: () 953 @*/ 954 PetscErrorCode PCGAMGSetNlevels(PC pc, PetscInt n) 955 { 956 PetscErrorCode ierr; 957 958 PetscFunctionBegin; 959 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 960 ierr = PetscTryMethod(pc,"PCGAMGSetNlevels_C",(PC,PetscInt),(pc,n));CHKERRQ(ierr); 961 PetscFunctionReturn(0); 962 } 963 964 #undef __FUNCT__ 965 #define __FUNCT__ "PCGAMGSetNlevels_GAMG" 966 static PetscErrorCode PCGAMGSetNlevels_GAMG(PC pc, PetscInt n) 967 { 968 PC_MG *mg = (PC_MG*)pc->data; 969 PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx; 970 971 PetscFunctionBegin; 972 pc_gamg->Nlevels = n; 973 PetscFunctionReturn(0); 974 } 975 976 #undef __FUNCT__ 977 #define __FUNCT__ "PCGAMGSetThreshold" 978 /*@ 979 PCGAMGSetThreshold - Relative threshold to use for dropping edges in aggregation graph 980 981 Not collective on PC 982 983 Input Parameters: 984 + pc - the preconditioner context 985 - threshold - the threshold value, 0.0 means keep all nonzero entries in the graph; negative means keep even zero entries in the graph 986 987 Options Database Key: 988 . -pc_gamg_threshold <threshold> 989 990 Notes: Before aggregating the graph GAMG will remove small values from the graph thus reducing the coupling in the graph and a different 991 (perhaps better) coarser set of points. 992 993 Level: intermediate 994 995 Concepts: Unstructured multigrid preconditioner 996 997 .seealso: () 998 @*/ 999 PetscErrorCode PCGAMGSetThreshold(PC pc, PetscReal n) 1000 { 1001 PetscErrorCode ierr; 1002 1003 PetscFunctionBegin; 1004 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 1005 ierr = PetscTryMethod(pc,"PCGAMGSetThreshold_C",(PC,PetscReal),(pc,n));CHKERRQ(ierr); 1006 PetscFunctionReturn(0); 1007 } 1008 1009 #undef __FUNCT__ 1010 #define __FUNCT__ "PCGAMGSetThreshold_GAMG" 1011 static PetscErrorCode PCGAMGSetThreshold_GAMG(PC pc, PetscReal n) 1012 { 1013 PC_MG *mg = (PC_MG*)pc->data; 1014 PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx; 1015 1016 PetscFunctionBegin; 1017 pc_gamg->threshold = n; 1018 PetscFunctionReturn(0); 1019 } 1020 1021 #undef __FUNCT__ 1022 #define __FUNCT__ "PCGAMGSetType" 1023 /*@ 1024 PCGAMGSetType - Set solution method 1025 1026 Collective on PC 1027 1028 Input Parameters: 1029 + pc - the preconditioner context 1030 - type - PCGAMGAGG, PCGAMGGEO, or PCGAMGCLASSICAL 1031 1032 Options Database Key: 1033 . -pc_gamg_type <agg,geo,classical> - type of algebraic multigrid to apply 1034 1035 Level: intermediate 1036 1037 Concepts: Unstructured multigrid preconditioner 1038 1039 .seealso: PCGAMGGetType(), PCGAMG, PCGAMGType 1040 @*/ 1041 PetscErrorCode PCGAMGSetType(PC pc, PCGAMGType type) 1042 { 1043 PetscErrorCode ierr; 1044 1045 PetscFunctionBegin; 1046 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 1047 ierr = PetscTryMethod(pc,"PCGAMGSetType_C",(PC,PCGAMGType),(pc,type));CHKERRQ(ierr); 1048 PetscFunctionReturn(0); 1049 } 1050 1051 #undef __FUNCT__ 1052 #define __FUNCT__ "PCGAMGGetType" 1053 /*@ 1054 PCGAMGGetType - Get solution method 1055 1056 Collective on PC 1057 1058 Input Parameter: 1059 . pc - the preconditioner context 1060 1061 Output Parameter: 1062 . type - the type of algorithm used 1063 1064 Level: intermediate 1065 1066 Concepts: Unstructured multigrid preconditioner 1067 1068 .seealso: PCGAMGSetType(), PCGAMGType 1069 @*/ 1070 PetscErrorCode PCGAMGGetType(PC pc, PCGAMGType *type) 1071 { 1072 PetscErrorCode ierr; 1073 1074 PetscFunctionBegin; 1075 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 1076 ierr = PetscUseMethod(pc,"PCGAMGGetType_C",(PC,PCGAMGType*),(pc,type));CHKERRQ(ierr); 1077 PetscFunctionReturn(0); 1078 } 1079 1080 #undef __FUNCT__ 1081 #define __FUNCT__ "PCGAMGGetType_GAMG" 1082 static PetscErrorCode PCGAMGGetType_GAMG(PC pc, PCGAMGType *type) 1083 { 1084 PC_MG *mg = (PC_MG*)pc->data; 1085 PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx; 1086 1087 PetscFunctionBegin; 1088 *type = pc_gamg->type; 1089 PetscFunctionReturn(0); 1090 } 1091 1092 #undef __FUNCT__ 1093 #define __FUNCT__ "PCGAMGSetType_GAMG" 1094 static PetscErrorCode PCGAMGSetType_GAMG(PC pc, PCGAMGType type) 1095 { 1096 PetscErrorCode ierr,(*r)(PC); 1097 PC_MG *mg = (PC_MG*)pc->data; 1098 PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx; 1099 1100 PetscFunctionBegin; 1101 pc_gamg->type = type; 1102 ierr = PetscFunctionListFind(GAMGList,type,&r);CHKERRQ(ierr); 1103 if (!r) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown GAMG type %s given",type); 1104 if (pc_gamg->ops->destroy) { 1105 ierr = (*pc_gamg->ops->destroy)(pc);CHKERRQ(ierr); 1106 ierr = PetscMemzero(pc_gamg->ops,sizeof(struct _PCGAMGOps));CHKERRQ(ierr); 1107 pc_gamg->ops->createlevel = PCGAMGCreateLevel_GAMG; 1108 /* cleaning up common data in pc_gamg - this should disapear someday */ 1109 pc_gamg->data_cell_cols = 0; 1110 pc_gamg->data_cell_rows = 0; 1111 pc_gamg->orig_data_cell_cols = 0; 1112 pc_gamg->orig_data_cell_rows = 0; 1113 ierr = PetscFree(pc_gamg->data);CHKERRQ(ierr); 1114 pc_gamg->data_sz = 0; 1115 } 1116 ierr = PetscFree(pc_gamg->gamg_type_name);CHKERRQ(ierr); 1117 ierr = PetscStrallocpy(type,&pc_gamg->gamg_type_name);CHKERRQ(ierr); 1118 ierr = (*r)(pc);CHKERRQ(ierr); 1119 PetscFunctionReturn(0); 1120 } 1121 1122 #undef __FUNCT__ 1123 #define __FUNCT__ "PCView_GAMG" 1124 static PetscErrorCode PCView_GAMG(PC pc,PetscViewer viewer) 1125 { 1126 PetscErrorCode ierr; 1127 PC_MG *mg = (PC_MG*)pc->data; 1128 PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx; 1129 1130 PetscFunctionBegin; 1131 ierr = PetscViewerASCIIPrintf(viewer," GAMG specific options\n");CHKERRQ(ierr); 1132 ierr = PetscViewerASCIIPrintf(viewer," Threshold for dropping small values from graph %g\n",(double)pc_gamg->threshold);CHKERRQ(ierr); 1133 if (pc_gamg->use_aggs_in_asm) { 1134 ierr = PetscViewerASCIIPrintf(viewer," Using aggregates from coarsening process to define subdomains for PCASM\n");CHKERRQ(ierr); 1135 } 1136 if (pc_gamg->ops->view) { 1137 ierr = (*pc_gamg->ops->view)(pc,viewer);CHKERRQ(ierr); 1138 } 1139 PetscFunctionReturn(0); 1140 } 1141 1142 #undef __FUNCT__ 1143 #define __FUNCT__ "PCSetFromOptions_GAMG" 1144 PetscErrorCode PCSetFromOptions_GAMG(PetscOptionItems *PetscOptionsObject,PC pc) 1145 { 1146 PetscErrorCode ierr; 1147 PC_MG *mg = (PC_MG*)pc->data; 1148 PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx; 1149 PetscBool flag; 1150 MPI_Comm comm; 1151 char prefix[256]; 1152 const char *pcpre; 1153 1154 PetscFunctionBegin; 1155 ierr = PetscObjectGetComm((PetscObject)pc,&comm);CHKERRQ(ierr); 1156 ierr = PetscOptionsHead(PetscOptionsObject,"GAMG options");CHKERRQ(ierr); 1157 { 1158 char tname[256]; 1159 ierr = PetscOptionsFList("-pc_gamg_type","Type of AMG method","PCGAMGSetType",GAMGList, pc_gamg->gamg_type_name, tname, sizeof(tname), &flag);CHKERRQ(ierr); 1160 if (flag) { 1161 ierr = PCGAMGSetType(pc,tname);CHKERRQ(ierr); 1162 } 1163 ierr = PetscOptionsBool("-pc_gamg_repartition","Repartion coarse grids","PCGAMGSetRepartition",pc_gamg->repart,&pc_gamg->repart,NULL);CHKERRQ(ierr); 1164 ierr = PetscOptionsBool("-pc_gamg_reuse_interpolation","Reuse prolongation operator","PCGAMGReuseInterpolation",pc_gamg->reuse_prol,&pc_gamg->reuse_prol,NULL);CHKERRQ(ierr); 1165 ierr = PetscOptionsBool("-pc_gamg_asm_use_agg","Use aggregation agragates for ASM smoother","PCGAMGASMSetUseAggs",pc_gamg->use_aggs_in_asm,&pc_gamg->use_aggs_in_asm,NULL);CHKERRQ(ierr); 1166 ierr = PetscOptionsInt("-pc_gamg_process_eq_limit","Limit (goal) on number of equations per process on coarse grids","PCGAMGSetProcEqLim",pc_gamg->min_eq_proc,&pc_gamg->min_eq_proc,NULL);CHKERRQ(ierr); 1167 ierr = PetscOptionsInt("-pc_gamg_coarse_eq_limit","Limit on number of equations for the coarse grid","PCGAMGSetCoarseEqLim",pc_gamg->coarse_eq_limit,&pc_gamg->coarse_eq_limit,NULL);CHKERRQ(ierr); 1168 ierr = PetscOptionsReal("-pc_gamg_threshold","Relative threshold to use for dropping edges in aggregation graph","PCGAMGSetThreshold",pc_gamg->threshold,&pc_gamg->threshold,&flag);CHKERRQ(ierr); 1169 ierr = PetscOptionsInt("-pc_mg_levels","Set number of MG levels","PCGAMGSetNlevels",pc_gamg->Nlevels,&pc_gamg->Nlevels,NULL);CHKERRQ(ierr); 1170 1171 /* set options for subtype */ 1172 if (pc_gamg->ops->setfromoptions) {ierr = (*pc_gamg->ops->setfromoptions)(PetscOptionsObject,pc);CHKERRQ(ierr);} 1173 } 1174 ierr = PCGetOptionsPrefix(pc, &pcpre);CHKERRQ(ierr); 1175 ierr = PetscSNPrintf(prefix,sizeof(prefix),"%spc_gamg_",pcpre ? pcpre : "");CHKERRQ(ierr); 1176 ierr = PetscObjectSetOptionsPrefix((PetscObject)pc_gamg->random,prefix);CHKERRQ(ierr); 1177 ierr = PetscRandomSetFromOptions(pc_gamg->random);CHKERRQ(ierr); 1178 ierr = PetscOptionsTail();CHKERRQ(ierr); 1179 PetscFunctionReturn(0); 1180 } 1181 1182 /* -------------------------------------------------------------------------- */ 1183 /*MC 1184 PCGAMG - Geometric algebraic multigrid (AMG) preconditioner 1185 1186 Options Database Keys: 1187 + -pc_gamg_type <type> - one of agg, geo, or classical 1188 . -pc_gamg_repartition <true,default=false> - repartition the degrees of freedom accross the coarse grids as they are determined 1189 . -pc_gamg_reuse_interpolation <true,default=false> - when rebuilding the algebraic multigrid preconditioner reuse the previously computed interpolations 1190 . -pc_gamg_asm_use_agg <true,default=false> - use the aggregates from the coasening process to defined the subdomains on each level for the PCASM smoother 1191 . -pc_gamg_process_eq_limit <limit, default=50> - GAMG will reduce the number of MPI processes used directly on the coarse grids so that there are around <limit> 1192 equations on each process that has degrees of freedom 1193 . -pc_gamg_coarse_eq_limit <limit, default=50> - Set maximum number of equations on coarsest grid to aim for. 1194 - -pc_gamg_threshold <thresh,default=0> - Before aggregating the graph GAMG will remove small values from the graph thus reducing the coupling in the graph and a different 1195 1196 Options Database Keys for default Aggregation: 1197 + -pc_gamg_agg_nsmooths <nsmooth, default=1> - number of smoothing steps to use with smooth aggregation 1198 . -pc_gamg_sym_graph <true,default=false> - symmetrize the graph before computing the aggregation 1199 - -pc_gamg_square_graph <n,default=1> - number of levels to square the graph before aggregating it 1200 1201 Multigrid options(inherited): 1202 + -pc_mg_cycles <v>: v or w (PCMGSetCycleType()) 1203 . -pc_mg_smoothup <1>: Number of post-smoothing steps (PCMGSetNumberSmoothUp) 1204 . -pc_mg_smoothdown <1>: Number of pre-smoothing steps (PCMGSetNumberSmoothDown) 1205 . -pc_mg_type <multiplicative>: (one of) additive multiplicative full kascade 1206 - -pc_mg_levels <levels> - Number of levels of multigrid to use. 1207 1208 1209 Notes: In order to obtain good performance for PCGAMG for vector valued problems you must 1210 $ Call MatSetBlockSize() to indicate the number of degrees of freedom per grid point 1211 $ Call MatSetNearNullSpace() (or PCSetCoordinates() if solving the equations of elasticity) to indicate the near null space of the operator 1212 $ See the Users Manual Chapter 4 for more details 1213 1214 Level: intermediate 1215 1216 Concepts: algebraic multigrid 1217 1218 .seealso: PCCreate(), PCSetType(), MatSetBlockSize(), PCMGType, PCSetCoordinates(), MatSetNearNullSpace(), PCGAMGSetType(), PCGAMGAGG, PCGAMGGEO, PCGAMGCLASSICAL, PCGAMGSetProcEqLim(), 1219 PCGAMGSetCoarseEqLim(), PCGAMGSetRepartition(), PCGAMGRegister(), PCGAMGSetReuseInterpolation(), PCGAMGASMSetUseAggs(), PCGAMGSetNlevels(), PCGAMGSetThreshold(), PCGAMGGetType(), PCGAMGSetReuseInterpolation() 1220 M*/ 1221 1222 #undef __FUNCT__ 1223 #define __FUNCT__ "PCCreate_GAMG" 1224 PETSC_EXTERN PetscErrorCode PCCreate_GAMG(PC pc) 1225 { 1226 PetscErrorCode ierr; 1227 PC_GAMG *pc_gamg; 1228 PC_MG *mg; 1229 1230 PetscFunctionBegin; 1231 /* register AMG type */ 1232 ierr = PCGAMGInitializePackage();CHKERRQ(ierr); 1233 1234 /* PCGAMG is an inherited class of PCMG. Initialize pc as PCMG */ 1235 ierr = PCSetType(pc, PCMG);CHKERRQ(ierr); 1236 ierr = PetscObjectChangeTypeName((PetscObject)pc, PCGAMG);CHKERRQ(ierr); 1237 1238 /* create a supporting struct and attach it to pc */ 1239 ierr = PetscNewLog(pc,&pc_gamg);CHKERRQ(ierr); 1240 mg = (PC_MG*)pc->data; 1241 mg->galerkin = 2; /* Use Galerkin, but it is computed externally from PCMG by GAMG code */ 1242 mg->innerctx = pc_gamg; 1243 1244 ierr = PetscNewLog(pc,&pc_gamg->ops);CHKERRQ(ierr); 1245 1246 pc_gamg->setup_count = 0; 1247 /* these should be in subctx but repartitioning needs simple arrays */ 1248 pc_gamg->data_sz = 0; 1249 pc_gamg->data = 0; 1250 1251 /* overwrite the pointers of PCMG by the functions of base class PCGAMG */ 1252 pc->ops->setfromoptions = PCSetFromOptions_GAMG; 1253 pc->ops->setup = PCSetUp_GAMG; 1254 pc->ops->reset = PCReset_GAMG; 1255 pc->ops->destroy = PCDestroy_GAMG; 1256 mg->view = PCView_GAMG; 1257 1258 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCGAMGSetProcEqLim_C",PCGAMGSetProcEqLim_GAMG);CHKERRQ(ierr); 1259 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCGAMGSetCoarseEqLim_C",PCGAMGSetCoarseEqLim_GAMG);CHKERRQ(ierr); 1260 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCGAMGSetRepartition_C",PCGAMGSetRepartition_GAMG);CHKERRQ(ierr); 1261 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCGAMGSetReuseInterpolation_C",PCGAMGSetReuseInterpolation_GAMG);CHKERRQ(ierr); 1262 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCGAMGASMSetUseAggs_C",PCGAMGASMSetUseAggs_GAMG);CHKERRQ(ierr); 1263 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCGAMGSetThreshold_C",PCGAMGSetThreshold_GAMG);CHKERRQ(ierr); 1264 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCGAMGSetType_C",PCGAMGSetType_GAMG);CHKERRQ(ierr); 1265 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCGAMGGetType_C",PCGAMGGetType_GAMG);CHKERRQ(ierr); 1266 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCGAMGSetNlevels_C",PCGAMGSetNlevels_GAMG);CHKERRQ(ierr); 1267 pc_gamg->repart = PETSC_FALSE; 1268 pc_gamg->reuse_prol = PETSC_FALSE; 1269 pc_gamg->use_aggs_in_asm = PETSC_FALSE; 1270 pc_gamg->min_eq_proc = 50; 1271 pc_gamg->coarse_eq_limit = 50; 1272 pc_gamg->threshold = 0.; 1273 pc_gamg->Nlevels = GAMG_MAXLEVELS; 1274 pc_gamg->current_level = 0; /* don't need to init really */ 1275 pc_gamg->ops->createlevel = PCGAMGCreateLevel_GAMG; 1276 1277 ierr = PetscRandomCreate(PetscObjectComm((PetscObject)pc),&pc_gamg->random);CHKERRQ(ierr); 1278 1279 /* PCSetUp_GAMG assumes that the type has been set, so set it to the default now */ 1280 ierr = PCGAMGSetType(pc,PCGAMGAGG);CHKERRQ(ierr); 1281 PetscFunctionReturn(0); 1282 } 1283 1284 #undef __FUNCT__ 1285 #define __FUNCT__ "PCGAMGInitializePackage" 1286 /*@C 1287 PCGAMGInitializePackage - This function initializes everything in the PCGAMG package. It is called 1288 from PetscDLLibraryRegister() when using dynamic libraries, and on the first call to PCCreate_GAMG() 1289 when using static libraries. 1290 1291 Level: developer 1292 1293 .keywords: PC, PCGAMG, initialize, package 1294 .seealso: PetscInitialize() 1295 @*/ 1296 PetscErrorCode PCGAMGInitializePackage(void) 1297 { 1298 PetscErrorCode ierr; 1299 1300 PetscFunctionBegin; 1301 if (PCGAMGPackageInitialized) PetscFunctionReturn(0); 1302 PCGAMGPackageInitialized = PETSC_TRUE; 1303 ierr = PetscFunctionListAdd(&GAMGList,PCGAMGGEO,PCCreateGAMG_GEO);CHKERRQ(ierr); 1304 ierr = PetscFunctionListAdd(&GAMGList,PCGAMGAGG,PCCreateGAMG_AGG);CHKERRQ(ierr); 1305 ierr = PetscFunctionListAdd(&GAMGList,PCGAMGCLASSICAL,PCCreateGAMG_Classical);CHKERRQ(ierr); 1306 ierr = PetscRegisterFinalize(PCGAMGFinalizePackage);CHKERRQ(ierr); 1307 1308 /* general events */ 1309 ierr = PetscLogEventRegister("PCGAMGGraph_AGG", 0, &PC_GAMGGraph_AGG);CHKERRQ(ierr); 1310 ierr = PetscLogEventRegister("PCGAMGGraph_GEO", PC_CLASSID, &PC_GAMGGraph_GEO);CHKERRQ(ierr); 1311 ierr = PetscLogEventRegister("PCGAMGCoarse_AGG", PC_CLASSID, &PC_GAMGCoarsen_AGG);CHKERRQ(ierr); 1312 ierr = PetscLogEventRegister("PCGAMGCoarse_GEO", PC_CLASSID, &PC_GAMGCoarsen_GEO);CHKERRQ(ierr); 1313 ierr = PetscLogEventRegister("PCGAMGProl_AGG", PC_CLASSID, &PC_GAMGProlongator_AGG);CHKERRQ(ierr); 1314 ierr = PetscLogEventRegister("PCGAMGProl_GEO", PC_CLASSID, &PC_GAMGProlongator_GEO);CHKERRQ(ierr); 1315 ierr = PetscLogEventRegister("PCGAMGPOpt_AGG", PC_CLASSID, &PC_GAMGOptProlongator_AGG);CHKERRQ(ierr); 1316 1317 #if defined PETSC_GAMG_USE_LOG 1318 ierr = PetscLogEventRegister("GAMG: createProl", PC_CLASSID, &petsc_gamg_setup_events[SET1]);CHKERRQ(ierr); 1319 ierr = PetscLogEventRegister(" Graph", PC_CLASSID, &petsc_gamg_setup_events[GRAPH]);CHKERRQ(ierr); 1320 /* PetscLogEventRegister(" G.Mat", PC_CLASSID, &petsc_gamg_setup_events[GRAPH_MAT]); */ 1321 /* PetscLogEventRegister(" G.Filter", PC_CLASSID, &petsc_gamg_setup_events[GRAPH_FILTER]); */ 1322 /* PetscLogEventRegister(" G.Square", PC_CLASSID, &petsc_gamg_setup_events[GRAPH_SQR]); */ 1323 ierr = PetscLogEventRegister(" MIS/Agg", PC_CLASSID, &petsc_gamg_setup_events[SET4]);CHKERRQ(ierr); 1324 ierr = PetscLogEventRegister(" geo: growSupp", PC_CLASSID, &petsc_gamg_setup_events[SET5]);CHKERRQ(ierr); 1325 ierr = PetscLogEventRegister(" geo: triangle", PC_CLASSID, &petsc_gamg_setup_events[SET6]);CHKERRQ(ierr); 1326 ierr = PetscLogEventRegister(" search-set", PC_CLASSID, &petsc_gamg_setup_events[FIND_V]);CHKERRQ(ierr); 1327 ierr = PetscLogEventRegister(" SA: col data", PC_CLASSID, &petsc_gamg_setup_events[SET7]);CHKERRQ(ierr); 1328 ierr = PetscLogEventRegister(" SA: frmProl0", PC_CLASSID, &petsc_gamg_setup_events[SET8]);CHKERRQ(ierr); 1329 ierr = PetscLogEventRegister(" SA: smooth", PC_CLASSID, &petsc_gamg_setup_events[SET9]);CHKERRQ(ierr); 1330 ierr = PetscLogEventRegister("GAMG: partLevel", PC_CLASSID, &petsc_gamg_setup_events[SET2]);CHKERRQ(ierr); 1331 ierr = PetscLogEventRegister(" repartition", PC_CLASSID, &petsc_gamg_setup_events[SET12]);CHKERRQ(ierr); 1332 ierr = PetscLogEventRegister(" Invert-Sort", PC_CLASSID, &petsc_gamg_setup_events[SET13]);CHKERRQ(ierr); 1333 ierr = PetscLogEventRegister(" Move A", PC_CLASSID, &petsc_gamg_setup_events[SET14]);CHKERRQ(ierr); 1334 ierr = PetscLogEventRegister(" Move P", PC_CLASSID, &petsc_gamg_setup_events[SET15]);CHKERRQ(ierr); 1335 1336 /* PetscLogEventRegister(" PL move data", PC_CLASSID, &petsc_gamg_setup_events[SET13]); */ 1337 /* PetscLogEventRegister("GAMG: fix", PC_CLASSID, &petsc_gamg_setup_events[SET10]); */ 1338 /* PetscLogEventRegister("GAMG: set levels", PC_CLASSID, &petsc_gamg_setup_events[SET11]); */ 1339 /* create timer stages */ 1340 #if defined GAMG_STAGES 1341 { 1342 char str[32]; 1343 PetscInt lidx; 1344 sprintf(str,"MG Level %d (finest)",0); 1345 ierr = PetscLogStageRegister(str, &gamg_stages[0]);CHKERRQ(ierr); 1346 for (lidx=1; lidx<9; lidx++) { 1347 sprintf(str,"MG Level %d",lidx); 1348 ierr = PetscLogStageRegister(str, &gamg_stages[lidx]);CHKERRQ(ierr); 1349 } 1350 } 1351 #endif 1352 #endif 1353 PetscFunctionReturn(0); 1354 } 1355 1356 #undef __FUNCT__ 1357 #define __FUNCT__ "PCGAMGFinalizePackage" 1358 /*@C 1359 PCGAMGFinalizePackage - This function frees everything from the PCGAMG package. It is 1360 called from PetscFinalize() automatically. 1361 1362 Level: developer 1363 1364 .keywords: Petsc, destroy, package 1365 .seealso: PetscFinalize() 1366 @*/ 1367 PetscErrorCode PCGAMGFinalizePackage(void) 1368 { 1369 PetscErrorCode ierr; 1370 1371 PetscFunctionBegin; 1372 PCGAMGPackageInitialized = PETSC_FALSE; 1373 ierr = PetscFunctionListDestroy(&GAMGList);CHKERRQ(ierr); 1374 PetscFunctionReturn(0); 1375 } 1376 1377 #undef __FUNCT__ 1378 #define __FUNCT__ "PCGAMGRegister" 1379 /*@C 1380 PCGAMGRegister - Register a PCGAMG implementation. 1381 1382 Input Parameters: 1383 + type - string that will be used as the name of the GAMG type. 1384 - create - function for creating the gamg context. 1385 1386 Level: advanced 1387 1388 .seealso: PCGAMGType, PCGAMG, PCGAMGSetType() 1389 @*/ 1390 PetscErrorCode PCGAMGRegister(PCGAMGType type, PetscErrorCode (*create)(PC)) 1391 { 1392 PetscErrorCode ierr; 1393 1394 PetscFunctionBegin; 1395 ierr = PCGAMGInitializePackage();CHKERRQ(ierr); 1396 ierr = PetscFunctionListAdd(&GAMGList,type,create);CHKERRQ(ierr); 1397 PetscFunctionReturn(0); 1398 } 1399 1400