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