Actual source code: cycliccuda.cu
slepc-3.11.2 2019-07-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2019, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
10: /*
11: SLEPc singular value solver: "cyclic" (CUDA implementation)
12: */
14: #include <slepc/private/svdimpl.h>
15: #include "../cyclicimpl.h"
16: #include <petsccuda.h>
18: PetscErrorCode MatMult_Cyclic_CUDA(Mat B,Vec x,Vec y)
19: {
21: SVD svd;
22: SVD_CYCLIC *cyclic;
23: PetscScalar *d_px,*d_py;
24: PetscInt m;
27: MatShellGetContext(B,(void**)&svd);
28: cyclic = (SVD_CYCLIC*)svd->data;
29: SVDMatGetLocalSize(svd,&m,NULL);
30: VecCUDAGetArrayRead(x,(const PetscScalar**)&d_px);
31: VecCUDAGetArrayReadWrite(y,&d_py);
32: VecCUDAPlaceArray(cyclic->x1,d_px);
33: VecCUDAPlaceArray(cyclic->x2,d_px+m);
34: VecCUDAPlaceArray(cyclic->y1,d_py);
35: VecCUDAPlaceArray(cyclic->y2,d_py+m);
36: SVDMatMult(svd,PETSC_FALSE,cyclic->x2,cyclic->y1);
37: SVDMatMult(svd,PETSC_TRUE,cyclic->x1,cyclic->y2);
38: VecCUDAResetArray(cyclic->x1);
39: VecCUDAResetArray(cyclic->x2);
40: VecCUDAResetArray(cyclic->y1);
41: VecCUDAResetArray(cyclic->y2);
42: VecCUDARestoreArrayRead(x,(const PetscScalar**)&d_px);
43: VecCUDARestoreArrayReadWrite(y,&d_py);
44: return(0);
45: }
47: PetscErrorCode MatCreateVecs_Cyclic_CUDA(Mat B,Vec *right,Vec *left)
48: {
50: SVD svd;
51: SVD_CYCLIC *cyclic;
52: PetscInt M,N,m,n;
53: PetscMPIInt size;
56: MatShellGetContext(B,(void**)&svd);
57: cyclic = (SVD_CYCLIC*)svd->data;
58: SVDMatGetSize(svd,&M,&N);
59: SVDMatGetLocalSize(svd,&m,&n);
60: MPI_Comm_size(PetscObjectComm((PetscObject)cyclic->mat),&size);
61: if (right) {
62: VecCreate(PetscObjectComm((PetscObject)cyclic->mat),right);
63: VecSetSizes(*right,m+n,M+N);
64: if (size>1) {
65: VecSetType(*right,VECMPICUDA);
66: } else {
67: VecSetType(*right,VECSEQCUDA);
68: }
69: }
70: if (left) {
71: VecCreate(PetscObjectComm((PetscObject)cyclic->mat),left);
72: VecSetSizes(*left,m+n,M+N);
73: if (size>1) {
74: VecSetType(*left,VECMPICUDA);
75: } else {
76: VecSetType(*left,VECSEQCUDA);
77: }
78: }
79: return(0);
80: }