1 2 3 /* 4 Defines the basic matrix operations for the AIJ (compressed row) 5 matrix storage format. 6 */ 7 8 #include <petscconf.h> 9 #include <../src/mat/impls/aij/seq/aij.h> /*I "petscmat.h" I*/ 10 #include <petscbt.h> 11 #include <../src/vec/vec/impls/dvecimpl.h> 12 #include <petsc/private/vecimpl.h> 13 14 #include <../src/mat/impls/aij/seq/seqviennacl/viennaclmatimpl.h> 15 16 17 #include <algorithm> 18 #include <vector> 19 #include <string> 20 21 #include "viennacl/linalg/prod.hpp" 22 23 PETSC_INTERN PetscErrorCode MatConvert_SeqAIJ_SeqAIJViennaCL(Mat A, MatType type, MatReuse reuse, Mat *newmat); 24 PETSC_INTERN PetscErrorCode MatGetFactor_seqaij_petsc(Mat,MatFactorType,Mat*); 25 26 27 PetscErrorCode MatViennaCLCopyToGPU(Mat A) 28 { 29 30 Mat_SeqAIJViennaCL *viennaclstruct = (Mat_SeqAIJViennaCL*)A->spptr; 31 Mat_SeqAIJ *a = (Mat_SeqAIJ*)A->data; 32 PetscErrorCode ierr; 33 34 PetscFunctionBegin; 35 if (A->rmap->n > 0 && A->cmap->n > 0 && a->nz) { //some OpenCL SDKs have issues with buffers of size 0 36 if (A->valid_GPU_matrix == PETSC_OFFLOAD_UNALLOCATED || A->valid_GPU_matrix == PETSC_OFFLOAD_CPU) { 37 ierr = PetscLogEventBegin(MAT_ViennaCLCopyToGPU,A,0,0,0);CHKERRQ(ierr); 38 39 try { 40 if (a->compressedrow.use) { 41 if (!viennaclstruct->compressed_mat) viennaclstruct->compressed_mat = new ViennaCLCompressedAIJMatrix(); 42 43 // Since PetscInt is different from cl_uint, we have to convert: 44 viennacl::backend::mem_handle dummy; 45 46 viennacl::backend::typesafe_host_array<unsigned int> row_buffer; row_buffer.raw_resize(dummy, a->compressedrow.nrows+1); 47 for (PetscInt i=0; i<=a->compressedrow.nrows; ++i) 48 row_buffer.set(i, (a->compressedrow.i)[i]); 49 50 viennacl::backend::typesafe_host_array<unsigned int> row_indices; row_indices.raw_resize(dummy, a->compressedrow.nrows); 51 for (PetscInt i=0; i<a->compressedrow.nrows; ++i) 52 row_indices.set(i, (a->compressedrow.rindex)[i]); 53 54 viennacl::backend::typesafe_host_array<unsigned int> col_buffer; col_buffer.raw_resize(dummy, a->nz); 55 for (PetscInt i=0; i<a->nz; ++i) 56 col_buffer.set(i, (a->j)[i]); 57 58 viennaclstruct->compressed_mat->set(row_buffer.get(), row_indices.get(), col_buffer.get(), a->a, A->rmap->n, A->cmap->n, a->compressedrow.nrows, a->nz); 59 } else { 60 if (!viennaclstruct->mat) viennaclstruct->mat = new ViennaCLAIJMatrix(); 61 62 // Since PetscInt is in general different from cl_uint, we have to convert: 63 viennacl::backend::mem_handle dummy; 64 65 viennacl::backend::typesafe_host_array<unsigned int> row_buffer; row_buffer.raw_resize(dummy, A->rmap->n+1); 66 for (PetscInt i=0; i<=A->rmap->n; ++i) 67 row_buffer.set(i, (a->i)[i]); 68 69 viennacl::backend::typesafe_host_array<unsigned int> col_buffer; col_buffer.raw_resize(dummy, a->nz); 70 for (PetscInt i=0; i<a->nz; ++i) 71 col_buffer.set(i, (a->j)[i]); 72 73 viennaclstruct->mat->set(row_buffer.get(), col_buffer.get(), a->a, A->rmap->n, A->cmap->n, a->nz); 74 } 75 ViennaCLWaitForGPU(); 76 } catch(std::exception const & ex) { 77 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"ViennaCL error: %s", ex.what()); 78 } 79 80 // Create temporary vector for v += A*x: 81 if (viennaclstruct->tempvec) { 82 if (viennaclstruct->tempvec->size() != static_cast<std::size_t>(A->rmap->n)) { 83 delete (ViennaCLVector*)viennaclstruct->tempvec; 84 viennaclstruct->tempvec = new ViennaCLVector(A->rmap->n); 85 } else { 86 viennaclstruct->tempvec->clear(); 87 } 88 } else { 89 viennaclstruct->tempvec = new ViennaCLVector(A->rmap->n); 90 } 91 92 A->valid_GPU_matrix = PETSC_OFFLOAD_BOTH; 93 94 ierr = PetscLogEventEnd(MAT_ViennaCLCopyToGPU,A,0,0,0);CHKERRQ(ierr); 95 } 96 } 97 PetscFunctionReturn(0); 98 } 99 100 PetscErrorCode MatViennaCLCopyFromGPU(Mat A, const ViennaCLAIJMatrix *Agpu) 101 { 102 Mat_SeqAIJ *a = (Mat_SeqAIJ*)A->data; 103 PetscInt m = A->rmap->n; 104 PetscErrorCode ierr; 105 106 107 PetscFunctionBegin; 108 if (A->valid_GPU_matrix == PETSC_OFFLOAD_UNALLOCATED) { 109 try { 110 if (a->compressedrow.use) SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "ViennaCL: Cannot handle row compression for GPU matrices"); 111 else { 112 113 if ((PetscInt)Agpu->size1() != m) SETERRQ2(PETSC_COMM_WORLD, PETSC_ERR_ARG_SIZ, "GPU matrix has %d rows, should be %d", Agpu->size1(), m); 114 a->nz = Agpu->nnz(); 115 a->maxnz = a->nz; /* Since we allocate exactly the right amount */ 116 A->preallocated = PETSC_TRUE; 117 if (a->singlemalloc) { 118 if (a->a) {ierr = PetscFree3(a->a,a->j,a->i);CHKERRQ(ierr);} 119 } else { 120 if (a->i) {ierr = PetscFree(a->i);CHKERRQ(ierr);} 121 if (a->j) {ierr = PetscFree(a->j);CHKERRQ(ierr);} 122 if (a->a) {ierr = PetscFree(a->a);CHKERRQ(ierr);} 123 } 124 ierr = PetscMalloc3(a->nz,&a->a,a->nz,&a->j,m+1,&a->i);CHKERRQ(ierr); 125 ierr = PetscLogObjectMemory((PetscObject)A, a->nz*(sizeof(PetscScalar)+sizeof(PetscInt))+(m+1)*sizeof(PetscInt));CHKERRQ(ierr); 126 127 a->singlemalloc = PETSC_TRUE; 128 129 /* Setup row lengths */ 130 if (a->imax) {ierr = PetscFree2(a->imax,a->ilen);CHKERRQ(ierr);} 131 ierr = PetscMalloc2(m,&a->imax,m,&a->ilen);CHKERRQ(ierr); 132 ierr = PetscLogObjectMemory((PetscObject)A, 2*m*sizeof(PetscInt));CHKERRQ(ierr); 133 134 /* Copy data back from GPU */ 135 viennacl::backend::typesafe_host_array<unsigned int> row_buffer; row_buffer.raw_resize(Agpu->handle1(), Agpu->size1() + 1); 136 137 // copy row array 138 viennacl::backend::memory_read(Agpu->handle1(), 0, row_buffer.raw_size(), row_buffer.get()); 139 (a->i)[0] = row_buffer[0]; 140 for (PetscInt i = 0; i < (PetscInt)Agpu->size1(); ++i) { 141 (a->i)[i+1] = row_buffer[i+1]; 142 a->imax[i] = a->ilen[i] = a->i[i+1] - a->i[i]; //Set imax[] and ilen[] arrays at the same time as i[] for better cache reuse 143 } 144 145 // copy column indices 146 viennacl::backend::typesafe_host_array<unsigned int> col_buffer; col_buffer.raw_resize(Agpu->handle2(), Agpu->nnz()); 147 viennacl::backend::memory_read(Agpu->handle2(), 0, col_buffer.raw_size(), col_buffer.get()); 148 for (PetscInt i=0; i < (PetscInt)Agpu->nnz(); ++i) 149 (a->j)[i] = col_buffer[i]; 150 151 // copy nonzero entries directly to destination (no conversion required) 152 viennacl::backend::memory_read(Agpu->handle(), 0, sizeof(PetscScalar)*Agpu->nnz(), a->a); 153 154 ViennaCLWaitForGPU(); 155 /* TODO: Once a->diag is moved out of MatAssemblyEnd(), invalidate it here. */ 156 } 157 } catch(std::exception const & ex) { 158 SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_LIB, "ViennaCL error: %s", ex.what()); 159 } 160 161 /* This assembly prevents resetting the flag to PETSC_OFFLOAD_CPU and recopying */ 162 ierr = MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 163 ierr = MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 164 165 A->valid_GPU_matrix = PETSC_OFFLOAD_BOTH; 166 } else SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "ViennaCL error: Only valid for unallocated GPU matrices"); 167 PetscFunctionReturn(0); 168 } 169 170 PetscErrorCode MatMult_SeqAIJViennaCL(Mat A,Vec xx,Vec yy) 171 { 172 Mat_SeqAIJ *a = (Mat_SeqAIJ*)A->data; 173 PetscErrorCode ierr; 174 Mat_SeqAIJViennaCL *viennaclstruct = (Mat_SeqAIJViennaCL*)A->spptr; 175 const ViennaCLVector *xgpu=NULL; 176 ViennaCLVector *ygpu=NULL; 177 178 PetscFunctionBegin; 179 if (A->rmap->n > 0 && A->cmap->n > 0 && a->nz) { 180 ierr = VecViennaCLGetArrayRead(xx,&xgpu);CHKERRQ(ierr); 181 ierr = VecViennaCLGetArrayWrite(yy,&ygpu);CHKERRQ(ierr); 182 try { 183 if (a->compressedrow.use) { 184 *ygpu = viennacl::linalg::prod(*viennaclstruct->compressed_mat, *xgpu); 185 } else { 186 *ygpu = viennacl::linalg::prod(*viennaclstruct->mat,*xgpu); 187 } 188 ViennaCLWaitForGPU(); 189 } catch (std::exception const & ex) { 190 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"ViennaCL error: %s", ex.what()); 191 } 192 ierr = VecViennaCLRestoreArrayRead(xx,&xgpu);CHKERRQ(ierr); 193 ierr = VecViennaCLRestoreArrayWrite(yy,&ygpu);CHKERRQ(ierr); 194 ierr = PetscLogFlops(2.0*a->nz - a->nonzerorowcnt);CHKERRQ(ierr); 195 } else { 196 ierr = VecSet(yy,0);CHKERRQ(ierr); 197 } 198 PetscFunctionReturn(0); 199 } 200 201 PetscErrorCode MatMultAdd_SeqAIJViennaCL(Mat A,Vec xx,Vec yy,Vec zz) 202 { 203 Mat_SeqAIJ *a = (Mat_SeqAIJ*)A->data; 204 PetscErrorCode ierr; 205 Mat_SeqAIJViennaCL *viennaclstruct = (Mat_SeqAIJViennaCL*)A->spptr; 206 const ViennaCLVector *xgpu=NULL,*ygpu=NULL; 207 ViennaCLVector *zgpu=NULL; 208 209 PetscFunctionBegin; 210 if (A->rmap->n > 0 && A->cmap->n > 0 && a->nz) { 211 try { 212 ierr = VecViennaCLGetArrayRead(xx,&xgpu);CHKERRQ(ierr); 213 ierr = VecViennaCLGetArrayRead(yy,&ygpu);CHKERRQ(ierr); 214 ierr = VecViennaCLGetArrayWrite(zz,&zgpu);CHKERRQ(ierr); 215 216 if (a->compressedrow.use) { 217 ViennaCLVector temp = viennacl::linalg::prod(*viennaclstruct->compressed_mat, *xgpu); 218 *zgpu = *ygpu + temp; 219 ViennaCLWaitForGPU(); 220 } else { 221 if (zz == xx || zz == yy) { //temporary required 222 ViennaCLVector temp = viennacl::linalg::prod(*viennaclstruct->mat, *xgpu); 223 *zgpu = *ygpu; 224 *zgpu += temp; 225 ViennaCLWaitForGPU(); 226 } else { 227 *viennaclstruct->tempvec = viennacl::linalg::prod(*viennaclstruct->mat, *xgpu); 228 *zgpu = *ygpu + *viennaclstruct->tempvec; 229 ViennaCLWaitForGPU(); 230 } 231 } 232 233 ierr = VecViennaCLRestoreArrayRead(xx,&xgpu);CHKERRQ(ierr); 234 ierr = VecViennaCLRestoreArrayRead(yy,&ygpu);CHKERRQ(ierr); 235 ierr = VecViennaCLRestoreArrayWrite(zz,&zgpu);CHKERRQ(ierr); 236 237 } catch(std::exception const & ex) { 238 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"ViennaCL error: %s", ex.what()); 239 } 240 ierr = PetscLogFlops(2.0*a->nz);CHKERRQ(ierr); 241 } else { 242 ierr = VecCopy(yy,zz);CHKERRQ(ierr); 243 } 244 PetscFunctionReturn(0); 245 } 246 247 PetscErrorCode MatAssemblyEnd_SeqAIJViennaCL(Mat A,MatAssemblyType mode) 248 { 249 PetscErrorCode ierr; 250 251 PetscFunctionBegin; 252 ierr = MatAssemblyEnd_SeqAIJ(A,mode);CHKERRQ(ierr); 253 if (!A->pinnedtocpu) { 254 ierr = MatViennaCLCopyToGPU(A);CHKERRQ(ierr); 255 } 256 PetscFunctionReturn(0); 257 } 258 259 /* --------------------------------------------------------------------------------*/ 260 /*@ 261 MatCreateSeqAIJViennaCL - Creates a sparse matrix in AIJ (compressed row) format 262 (the default parallel PETSc format). This matrix will ultimately be pushed down 263 to GPUs and use the ViennaCL library for calculations. For good matrix 264 assembly performance the user should preallocate the matrix storage by setting 265 the parameter nz (or the array nnz). By setting these parameters accurately, 266 performance during matrix assembly can be increased substantially. 267 268 269 Collective 270 271 Input Parameters: 272 + comm - MPI communicator, set to PETSC_COMM_SELF 273 . m - number of rows 274 . n - number of columns 275 . nz - number of nonzeros per row (same for all rows) 276 - nnz - array containing the number of nonzeros in the various rows 277 (possibly different for each row) or NULL 278 279 Output Parameter: 280 . A - the matrix 281 282 It is recommended that one use the MatCreate(), MatSetType() and/or MatSetFromOptions(), 283 MatXXXXSetPreallocation() paradigm instead of this routine directly. 284 [MatXXXXSetPreallocation() is, for example, MatSeqAIJSetPreallocation] 285 286 Notes: 287 If nnz is given then nz is ignored 288 289 The AIJ format (also called the Yale sparse matrix format or 290 compressed row storage), is fully compatible with standard Fortran 77 291 storage. That is, the stored row and column indices can begin at 292 either one (as in Fortran) or zero. See the users' manual for details. 293 294 Specify the preallocated storage with either nz or nnz (not both). 295 Set nz=PETSC_DEFAULT and nnz=NULL for PETSc to control dynamic memory 296 allocation. For large problems you MUST preallocate memory or you 297 will get TERRIBLE performance, see the users' manual chapter on matrices. 298 299 Level: intermediate 300 301 .seealso: MatCreate(), MatCreateAIJ(), MatCreateAIJCUSPARSE(), MatSetValues(), MatSeqAIJSetColumnIndices(), MatCreateSeqAIJWithArrays(), MatCreateAIJ() 302 303 @*/ 304 PetscErrorCode MatCreateSeqAIJViennaCL(MPI_Comm comm,PetscInt m,PetscInt n,PetscInt nz,const PetscInt nnz[],Mat *A) 305 { 306 PetscErrorCode ierr; 307 308 PetscFunctionBegin; 309 ierr = MatCreate(comm,A);CHKERRQ(ierr); 310 ierr = MatSetSizes(*A,m,n,m,n);CHKERRQ(ierr); 311 ierr = MatSetType(*A,MATSEQAIJVIENNACL);CHKERRQ(ierr); 312 ierr = MatSeqAIJSetPreallocation_SeqAIJ(*A,nz,(PetscInt*)nnz);CHKERRQ(ierr); 313 PetscFunctionReturn(0); 314 } 315 316 317 PetscErrorCode MatDestroy_SeqAIJViennaCL(Mat A) 318 { 319 PetscErrorCode ierr; 320 Mat_SeqAIJViennaCL *viennaclcontainer = (Mat_SeqAIJViennaCL*)A->spptr; 321 322 PetscFunctionBegin; 323 try { 324 if (viennaclcontainer) { 325 delete viennaclcontainer->tempvec; 326 delete viennaclcontainer->mat; 327 delete viennaclcontainer->compressed_mat; 328 delete viennaclcontainer; 329 } 330 A->valid_GPU_matrix = PETSC_OFFLOAD_UNALLOCATED; 331 } catch(std::exception const & ex) { 332 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"ViennaCL error: %s", ex.what()); 333 } 334 335 ierr = PetscObjectComposeFunction((PetscObject)A,"MatConvert_seqaij_seqaijviennacl_C",NULL);CHKERRQ(ierr); 336 337 /* this next line is because MatDestroy tries to PetscFree spptr if it is not zero, and PetscFree only works if the memory was allocated with PetscNew or PetscMalloc, which don't call the constructor */ 338 A->spptr = 0; 339 ierr = MatDestroy_SeqAIJ(A);CHKERRQ(ierr); 340 PetscFunctionReturn(0); 341 } 342 343 344 PETSC_EXTERN PetscErrorCode MatCreate_SeqAIJViennaCL(Mat B) 345 { 346 PetscErrorCode ierr; 347 348 PetscFunctionBegin; 349 ierr = MatCreate_SeqAIJ(B);CHKERRQ(ierr); 350 ierr = MatConvert_SeqAIJ_SeqAIJViennaCL(B,MATSEQAIJVIENNACL,MAT_INPLACE_MATRIX,&B); 351 PetscFunctionReturn(0); 352 } 353 354 static PetscErrorCode MatPinToCPU_SeqAIJViennaCL(Mat,PetscBool); 355 static PetscErrorCode MatDuplicate_SeqAIJViennaCL(Mat A,MatDuplicateOption cpvalues,Mat *B) 356 { 357 PetscErrorCode ierr; 358 Mat C; 359 360 PetscFunctionBegin; 361 ierr = MatDuplicate_SeqAIJ(A,cpvalues,B);CHKERRQ(ierr); 362 C = *B; 363 364 ierr = MatPinToCPU_SeqAIJViennaCL(A,PETSC_FALSE);CHKERRQ(ierr); 365 C->ops->pintocpu = MatPinToCPU_SeqAIJViennaCL; 366 367 C->spptr = new Mat_SeqAIJViennaCL(); 368 ((Mat_SeqAIJViennaCL*)C->spptr)->tempvec = NULL; 369 ((Mat_SeqAIJViennaCL*)C->spptr)->mat = NULL; 370 ((Mat_SeqAIJViennaCL*)C->spptr)->compressed_mat = NULL; 371 372 ierr = PetscObjectChangeTypeName((PetscObject)C,MATSEQAIJVIENNACL);CHKERRQ(ierr); 373 374 C->valid_GPU_matrix = PETSC_OFFLOAD_UNALLOCATED; 375 376 /* If the source matrix is already assembled, copy the destination matrix to the GPU */ 377 if (C->assembled) { 378 ierr = MatViennaCLCopyToGPU(C);CHKERRQ(ierr); 379 } 380 381 PetscFunctionReturn(0); 382 } 383 384 static PetscErrorCode MatPinToCPU_SeqAIJViennaCL(Mat A,PetscBool flg) 385 { 386 PetscFunctionBegin; 387 A->pinnedtocpu = flg; 388 if (flg) { 389 A->ops->mult = MatMult_SeqAIJ; 390 A->ops->multadd = MatMultAdd_SeqAIJ; 391 A->ops->assemblyend = MatAssemblyEnd_SeqAIJ; 392 A->ops->destroy = MatDestroy_SeqAIJ; 393 A->ops->duplicate = MatDuplicate_SeqAIJ; 394 } else { 395 A->ops->mult = MatMult_SeqAIJViennaCL; 396 A->ops->multadd = MatMultAdd_SeqAIJViennaCL; 397 A->ops->assemblyend = MatAssemblyEnd_SeqAIJViennaCL; 398 A->ops->destroy = MatDestroy_SeqAIJViennaCL; 399 A->ops->duplicate = MatDuplicate_SeqAIJViennaCL; 400 } 401 PetscFunctionReturn(0); 402 } 403 404 PETSC_INTERN PetscErrorCode MatConvert_SeqAIJ_SeqAIJViennaCL(Mat A,MatType type,MatReuse reuse,Mat *newmat) 405 { 406 PetscErrorCode ierr; 407 Mat B; 408 Mat_SeqAIJ *aij; 409 410 PetscFunctionBegin; 411 412 if (reuse == MAT_REUSE_MATRIX) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"MAT_REUSE_MATRIX is not supported. Consider using MAT_INPLACE_MATRIX instead"); 413 414 if (reuse == MAT_INITIAL_MATRIX) { 415 ierr = MatDuplicate(A,MAT_COPY_VALUES,newmat);CHKERRQ(ierr); 416 } 417 418 B = *newmat; 419 420 aij = (Mat_SeqAIJ*)B->data; 421 aij->inode.use = PETSC_FALSE; 422 423 B->spptr = new Mat_SeqAIJViennaCL(); 424 425 ((Mat_SeqAIJViennaCL*)B->spptr)->tempvec = NULL; 426 ((Mat_SeqAIJViennaCL*)B->spptr)->mat = NULL; 427 ((Mat_SeqAIJViennaCL*)B->spptr)->compressed_mat = NULL; 428 429 ierr = MatPinToCPU_SeqAIJViennaCL(A,PETSC_FALSE);CHKERRQ(ierr); 430 A->ops->pintocpu = MatPinToCPU_SeqAIJViennaCL; 431 432 ierr = PetscObjectChangeTypeName((PetscObject)B,MATSEQAIJVIENNACL);CHKERRQ(ierr); 433 ierr = PetscFree(B->defaultvectype);CHKERRQ(ierr); 434 ierr = PetscStrallocpy(VECVIENNACL,&B->defaultvectype);CHKERRQ(ierr); 435 436 ierr = PetscObjectComposeFunction((PetscObject)A,"MatConvert_seqaij_seqaijviennacl_C",MatConvert_SeqAIJ_SeqAIJViennaCL);CHKERRQ(ierr); 437 438 B->valid_GPU_matrix = PETSC_OFFLOAD_UNALLOCATED; 439 440 /* If the source matrix is already assembled, copy the destination matrix to the GPU */ 441 if (B->assembled) { 442 ierr = MatViennaCLCopyToGPU(B);CHKERRQ(ierr); 443 } 444 445 PetscFunctionReturn(0); 446 } 447 448 449 /*MC 450 MATSEQAIJVIENNACL - MATAIJVIENNACL = "aijviennacl" = "seqaijviennacl" - A matrix type to be used for sparse matrices. 451 452 A matrix type type whose data resides on GPUs. These matrices are in CSR format by 453 default. All matrix calculations are performed using the ViennaCL library. 454 455 Options Database Keys: 456 + -mat_type aijviennacl - sets the matrix type to "seqaijviennacl" during a call to MatSetFromOptions() 457 . -mat_viennacl_storage_format csr - sets the storage format of matrices for MatMult during a call to MatSetFromOptions(). 458 - -mat_viennacl_mult_storage_format csr - sets the storage format of matrices for MatMult during a call to MatSetFromOptions(). 459 460 Level: beginner 461 462 .seealso: MatCreateSeqAIJViennaCL(), MATAIJVIENNACL, MatCreateAIJViennaCL() 463 M*/ 464 465 466 PETSC_EXTERN PetscErrorCode MatSolverTypeRegister_ViennaCL(void) 467 { 468 PetscErrorCode ierr; 469 470 PetscFunctionBegin; 471 ierr = MatSolverTypeRegister(MATSOLVERPETSC, MATSEQAIJVIENNACL, MAT_FACTOR_LU,MatGetFactor_seqaij_petsc);CHKERRQ(ierr); 472 ierr = MatSolverTypeRegister(MATSOLVERPETSC, MATSEQAIJVIENNACL, MAT_FACTOR_CHOLESKY,MatGetFactor_seqaij_petsc);CHKERRQ(ierr); 473 ierr = MatSolverTypeRegister(MATSOLVERPETSC, MATSEQAIJVIENNACL, MAT_FACTOR_ILU,MatGetFactor_seqaij_petsc);CHKERRQ(ierr); 474 ierr = MatSolverTypeRegister(MATSOLVERPETSC, MATSEQAIJVIENNACL, MAT_FACTOR_ICC,MatGetFactor_seqaij_petsc);CHKERRQ(ierr); 475 PetscFunctionReturn(0); 476 } 477