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