1 #include "petscconf.h" 2 #include <../src/mat/impls/aij/mpi/mpiaij.h> /*I "petscmat.h" I*/ 3 4 #undef __FUNCT__ 5 #define __FUNCT__ "MatMPIAIJSetPreallocation_MPIAIJViennaCL" 6 PetscErrorCode MatMPIAIJSetPreallocation_MPIAIJViennaCL(Mat B,PetscInt d_nz,const PetscInt d_nnz[],PetscInt o_nz,const PetscInt o_nnz[]) 7 { 8 Mat_MPIAIJ *b = (Mat_MPIAIJ*)B->data; 9 PetscErrorCode ierr; 10 PetscInt i; 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(B,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(B,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 #undef __FUNCT__ 33 #define __FUNCT__ "MatGetVecs_MPIAIJViennaCL" 34 PetscErrorCode MatGetVecs_MPIAIJViennaCL(Mat mat,Vec *right,Vec *left) 35 { 36 PetscErrorCode ierr; 37 38 PetscFunctionBegin; 39 if (right) { 40 ierr = VecCreate(PetscObjectComm((PetscObject)mat),right);CHKERRQ(ierr); 41 ierr = VecSetSizes(*right,mat->cmap->n,PETSC_DETERMINE);CHKERRQ(ierr); 42 ierr = VecSetBlockSize(*right,mat->rmap->bs);CHKERRQ(ierr); 43 ierr = VecSetType(*right,VECVIENNACL);CHKERRQ(ierr); 44 ierr = VecSetLayout(*right,mat->cmap);CHKERRQ(ierr); 45 } 46 if (left) { 47 ierr = VecCreate(PetscObjectComm((PetscObject)mat),left);CHKERRQ(ierr); 48 ierr = VecSetSizes(*left,mat->rmap->n,PETSC_DETERMINE);CHKERRQ(ierr); 49 ierr = VecSetBlockSize(*left,mat->rmap->bs);CHKERRQ(ierr); 50 ierr = VecSetType(*left,VECVIENNACL);CHKERRQ(ierr); 51 ierr = VecSetLayout(*left,mat->rmap);CHKERRQ(ierr); 52 } 53 PetscFunctionReturn(0); 54 } 55 56 57 #undef __FUNCT__ 58 #define __FUNCT__ "MatDestroy_MPIAIJViennaCL" 59 PetscErrorCode MatDestroy_MPIAIJViennaCL(Mat A) 60 { 61 PetscErrorCode ierr; 62 63 PetscFunctionBegin; 64 ierr = MatDestroy_MPIAIJ(A);CHKERRQ(ierr); 65 PetscFunctionReturn(0); 66 } 67 68 #undef __FUNCT__ 69 #define __FUNCT__ "MatCreate_MPIAIJViennaCL" 70 PETSC_EXTERN PetscErrorCode MatCreate_MPIAIJViennaCL(Mat A) 71 { 72 PetscErrorCode ierr; 73 74 PetscFunctionBegin; 75 ierr = MatCreate_MPIAIJ(A);CHKERRQ(ierr); 76 ierr = PetscObjectComposeFunction((PetscObject)A,"MatMPIAIJSetPreallocation_C",MatMPIAIJSetPreallocation_MPIAIJViennaCL);CHKERRQ(ierr); 77 A->ops->getvecs = MatGetVecs_MPIAIJViennaCL; 78 79 ierr = PetscObjectChangeTypeName((PetscObject)A,MATMPIAIJVIENNACL);CHKERRQ(ierr); 80 PetscFunctionReturn(0); 81 } 82 83 84 /*@ 85 MatCreateAIJViennaCL - Creates a sparse matrix in AIJ (compressed row) format 86 (the default parallel PETSc format). This matrix will ultimately pushed down 87 to GPUs and use the ViennaCL library for calculations. For good matrix 88 assembly performance the user should preallocate the matrix storage by setting 89 the parameter nz (or the array nnz). By setting these parameters accurately, 90 performance during matrix assembly can be increased substantially. 91 92 93 Collective on MPI_Comm 94 95 Input Parameters: 96 + comm - MPI communicator, set to PETSC_COMM_SELF 97 . m - number of rows 98 . n - number of columns 99 . nz - number of nonzeros per row (same for all rows) 100 - nnz - array containing the number of nonzeros in the various rows 101 (possibly different for each row) or NULL 102 103 Output Parameter: 104 . A - the matrix 105 106 It is recommended that one use the MatCreate(), MatSetType() and/or MatSetFromOptions(), 107 MatXXXXSetPreallocation() paradigm instead of this routine directly. 108 [MatXXXXSetPreallocation() is, for example, MatSeqAIJSetPreallocation] 109 110 Notes: 111 If nnz is given then nz is ignored 112 113 The AIJ format (also called the Yale sparse matrix format or 114 compressed row storage), is fully compatible with standard Fortran 77 115 storage. That is, the stored row and column indices can begin at 116 either one (as in Fortran) or zero. See the users' manual for details. 117 118 Specify the preallocated storage with either nz or nnz (not both). 119 Set nz=PETSC_DEFAULT and nnz=NULL for PETSc to control dynamic memory 120 allocation. For large problems you MUST preallocate memory or you 121 will get TERRIBLE performance, see the users' manual chapter on matrices. 122 123 By default, this format uses inodes (identical nodes) when possible, to 124 improve numerical efficiency of matrix-vector products and solves. We 125 search for consecutive rows with the same nonzero structure, thereby 126 reusing matrix information to achieve increased efficiency. 127 128 Level: intermediate 129 130 .seealso: MatCreate(), MatCreateAIJ(), MatCreateAIJCUSP(), MatSetValues(), MatSeqAIJSetColumnIndices(), MatCreateSeqAIJWithArrays(), MatCreateAIJ(), MATMPIAIJVIENNACL, MATAIJVIENNACL 131 @*/ 132 #undef __FUNCT__ 133 #define __FUNCT__ "MatCreateAIJViennaCL" 134 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) 135 { 136 PetscErrorCode ierr; 137 PetscMPIInt size; 138 139 PetscFunctionBegin; 140 ierr = MatCreate(comm,A);CHKERRQ(ierr); 141 ierr = MatSetSizes(*A,m,n,M,N);CHKERRQ(ierr); 142 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 143 if (size > 1) { 144 ierr = MatSetType(*A,MATMPIAIJVIENNACL);CHKERRQ(ierr); 145 ierr = MatMPIAIJSetPreallocation(*A,d_nz,d_nnz,o_nz,o_nnz);CHKERRQ(ierr); 146 } else { 147 ierr = MatSetType(*A,MATSEQAIJVIENNACL);CHKERRQ(ierr); 148 ierr = MatSeqAIJSetPreallocation(*A,d_nz,d_nnz);CHKERRQ(ierr); 149 } 150 PetscFunctionReturn(0); 151 } 152 153 /*M 154 MATAIJVIENNACL - MATMPIAIJVIENNACL= "aijviennacl" = "mpiaijviennacl" - A matrix type to be used for sparse matrices. 155 156 A matrix type (CSR format) whose data resides on GPUs. 157 All matrix calculations are performed using the ViennaCL library. 158 159 This matrix type is identical to MATSEQAIJVIENNACL when constructed with a single process communicator, 160 and MATMPIAIJVIENNACL otherwise. As a result, for single process communicators, 161 MatSeqAIJSetPreallocation is supported, and similarly MatMPIAIJSetPreallocation is supported 162 for communicators controlling multiple processes. It is recommended that you call both of 163 the above preallocation routines for simplicity. 164 165 Options Database Keys: 166 + -mat_type mpiaijviennacl - sets the matrix type to "mpiaijviennacl" during a call to MatSetFromOptions() 167 168 Level: beginner 169 170 .seealso: MatCreateAIJViennaCL(), MATSEQAIJVIENNACL, MatCreateSeqAIJVIENNACL() 171 M*/ 172 173