1 #include <petscconf.h> 2 #include <../src/mat/impls/aij/mpi/mpiaij.h> /*I "petscmat.h" I*/ 3 #include <../src/mat/impls/aij/seq/seqviennacl/viennaclmatimpl.h> 4 5 PetscErrorCode MatMPIAIJSetPreallocation_MPIAIJViennaCL(Mat B,PetscInt d_nz,const PetscInt d_nnz[],PetscInt o_nz,const PetscInt o_nnz[]) 6 { 7 Mat_MPIAIJ *b = (Mat_MPIAIJ*)B->data; 8 PetscErrorCode ierr; 9 10 PetscFunctionBegin; 11 ierr = PetscLayoutSetUp(B->rmap);CHKERRQ(ierr); 12 ierr = PetscLayoutSetUp(B->cmap);CHKERRQ(ierr); 13 if (!B->preallocated) { 14 /* Explicitly create the two MATSEQAIJVIENNACL matrices. */ 15 ierr = MatCreate(PETSC_COMM_SELF,&b->A);CHKERRQ(ierr); 16 ierr = MatSetSizes(b->A,B->rmap->n,B->cmap->n,B->rmap->n,B->cmap->n);CHKERRQ(ierr); 17 ierr = MatSetType(b->A,MATSEQAIJVIENNACL);CHKERRQ(ierr); 18 ierr = PetscLogObjectParent((PetscObject)B,(PetscObject)b->A);CHKERRQ(ierr); 19 ierr = MatCreate(PETSC_COMM_SELF,&b->B);CHKERRQ(ierr); 20 ierr = MatSetSizes(b->B,B->rmap->n,B->cmap->N,B->rmap->n,B->cmap->N);CHKERRQ(ierr); 21 ierr = MatSetType(b->B,MATSEQAIJVIENNACL);CHKERRQ(ierr); 22 ierr = PetscLogObjectParent((PetscObject)B,(PetscObject)b->B);CHKERRQ(ierr); 23 } 24 ierr = MatSeqAIJSetPreallocation(b->A,d_nz,d_nnz);CHKERRQ(ierr); 25 ierr = MatSeqAIJSetPreallocation(b->B,o_nz,o_nnz);CHKERRQ(ierr); 26 B->preallocated = PETSC_TRUE; 27 PetscFunctionReturn(0); 28 } 29 30 PetscErrorCode MatAssemblyEnd_MPIAIJViennaCL(Mat A,MatAssemblyType mode) 31 { 32 Mat_MPIAIJ *b = (Mat_MPIAIJ*)A->data; 33 PetscErrorCode ierr; 34 PetscBool v; 35 36 PetscFunctionBegin; 37 ierr = MatAssemblyEnd_MPIAIJ(A,mode);CHKERRQ(ierr); 38 ierr = PetscObjectTypeCompare((PetscObject)b->lvec,VECSEQVIENNACL,&v);CHKERRQ(ierr); 39 if (!v) { 40 PetscInt m; 41 ierr = VecGetSize(b->lvec,&m);CHKERRQ(ierr); 42 ierr = VecDestroy(&b->lvec);CHKERRQ(ierr); 43 ierr = VecCreateSeqViennaCL(PETSC_COMM_SELF,m,&b->lvec);CHKERRQ(ierr); 44 } 45 PetscFunctionReturn(0); 46 } 47 48 PetscErrorCode MatDestroy_MPIAIJViennaCL(Mat A) 49 { 50 PetscErrorCode ierr; 51 52 PetscFunctionBegin; 53 ierr = MatDestroy_MPIAIJ(A);CHKERRQ(ierr); 54 PetscFunctionReturn(0); 55 } 56 57 PETSC_EXTERN PetscErrorCode MatCreate_MPIAIJViennaCL(Mat A) 58 { 59 PetscErrorCode ierr; 60 61 PetscFunctionBegin; 62 ierr = MatCreate_MPIAIJ(A);CHKERRQ(ierr); 63 ierr = PetscFree(A->defaultvectype);CHKERRQ(ierr); 64 ierr = PetscStrallocpy(VECVIENNACL,&A->defaultvectype);CHKERRQ(ierr); 65 ierr = PetscObjectComposeFunction((PetscObject)A,"MatMPIAIJSetPreallocation_C",MatMPIAIJSetPreallocation_MPIAIJViennaCL);CHKERRQ(ierr); 66 A->ops->assemblyend = MatAssemblyEnd_MPIAIJViennaCL; 67 ierr = PetscObjectChangeTypeName((PetscObject)A,MATMPIAIJVIENNACL);CHKERRQ(ierr); 68 PetscFunctionReturn(0); 69 } 70 71 72 /*@ 73 MatCreateAIJViennaCL - Creates a sparse matrix in AIJ (compressed row) format 74 (the default parallel PETSc format). This matrix will ultimately be pushed down 75 to GPUs and use the ViennaCL library for calculations. For good matrix 76 assembly performance the user should preallocate the matrix storage by setting 77 the parameter nz (or the array nnz). By setting these parameters accurately, 78 performance during matrix assembly can be increased substantially. 79 80 81 Collective 82 83 Input Parameters: 84 + comm - MPI communicator, set to PETSC_COMM_SELF 85 . m - number of rows 86 . n - number of columns 87 . nz - number of nonzeros per row (same for all rows) 88 - nnz - array containing the number of nonzeros in the various rows 89 (possibly different for each row) or NULL 90 91 Output Parameter: 92 . A - the matrix 93 94 It is recommended that one use the MatCreate(), MatSetType() and/or MatSetFromOptions(), 95 MatXXXXSetPreallocation() paradigm instead of this routine directly. 96 [MatXXXXSetPreallocation() is, for example, MatSeqAIJSetPreallocation] 97 98 Notes: 99 If nnz is given then nz is ignored 100 101 The AIJ format (also called the Yale sparse matrix format or 102 compressed row storage), is fully compatible with standard Fortran 77 103 storage. That is, the stored row and column indices can begin at 104 either one (as in Fortran) or zero. See the users' manual for details. 105 106 Specify the preallocated storage with either nz or nnz (not both). 107 Set nz=PETSC_DEFAULT and nnz=NULL for PETSc to control dynamic memory 108 allocation. For large problems you MUST preallocate memory or you 109 will get TERRIBLE performance, see the users' manual chapter on matrices. 110 111 Level: intermediate 112 113 .seealso: MatCreate(), MatCreateAIJ(), MatCreateAIJCUSPARSE(), MatSetValues(), MatSeqAIJSetColumnIndices(), MatCreateSeqAIJWithArrays(), MatCreateAIJ(), MATMPIAIJVIENNACL, MATAIJVIENNACL 114 @*/ 115 PetscErrorCode MatCreateAIJViennaCL(MPI_Comm comm,PetscInt m,PetscInt n,PetscInt M,PetscInt N,PetscInt d_nz,const PetscInt d_nnz[],PetscInt o_nz,const PetscInt o_nnz[],Mat *A) 116 { 117 PetscErrorCode ierr; 118 PetscMPIInt size; 119 120 PetscFunctionBegin; 121 ierr = MatCreate(comm,A);CHKERRQ(ierr); 122 ierr = MatSetSizes(*A,m,n,M,N);CHKERRQ(ierr); 123 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 124 if (size > 1) { 125 ierr = MatSetType(*A,MATMPIAIJVIENNACL);CHKERRQ(ierr); 126 ierr = MatMPIAIJSetPreallocation(*A,d_nz,d_nnz,o_nz,o_nnz);CHKERRQ(ierr); 127 } else { 128 ierr = MatSetType(*A,MATSEQAIJVIENNACL);CHKERRQ(ierr); 129 ierr = MatSeqAIJSetPreallocation(*A,d_nz,d_nnz);CHKERRQ(ierr); 130 } 131 PetscFunctionReturn(0); 132 } 133 134 /*MC 135 MATAIJVIENNACL - MATMPIAIJVIENNACL= "aijviennacl" = "mpiaijviennacl" - A matrix type to be used for sparse matrices. 136 137 A matrix type (CSR format) whose data resides on GPUs. 138 All matrix calculations are performed using the ViennaCL library. 139 140 This matrix type is identical to MATSEQAIJVIENNACL when constructed with a single process communicator, 141 and MATMPIAIJVIENNACL otherwise. As a result, for single process communicators, 142 MatSeqAIJSetPreallocation is supported, and similarly MatMPIAIJSetPreallocation is supported 143 for communicators controlling multiple processes. It is recommended that you call both of 144 the above preallocation routines for simplicity. 145 146 Options Database Keys: 147 + -mat_type mpiaijviennacl - sets the matrix type to "mpiaijviennacl" during a call to MatSetFromOptions() 148 149 Level: beginner 150 151 .seealso: MatCreateAIJViennaCL(), MATSEQAIJVIENNACL, MatCreateSeqAIJVIENNACL() 152 M*/ 153