Actual source code: cycliccuda.cu
slepc-3.13.4 2020-09-02
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2020, 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: */
13: #include <slepc/private/svdimpl.h>
14: #include "../src/svd/impls/cyclic/cyclic.h"
16: PetscErrorCode MatMult_Cyclic_CUDA(Mat B,Vec x,Vec y)
17: {
19: SVD svd;
20: SVD_CYCLIC *cyclic;
21: PetscScalar *d_px,*d_py;
22: PetscInt m;
25: MatShellGetContext(B,(void**)&svd);
26: cyclic = (SVD_CYCLIC*)svd->data;
27: SVDMatGetLocalSize(svd,&m,NULL);
28: VecCUDAGetArrayRead(x,(const PetscScalar**)&d_px);
29: VecCUDAGetArray(y,&d_py);
30: VecCUDAPlaceArray(cyclic->x1,d_px);
31: VecCUDAPlaceArray(cyclic->x2,d_px+m);
32: VecCUDAPlaceArray(cyclic->y1,d_py);
33: VecCUDAPlaceArray(cyclic->y2,d_py+m);
34: SVDMatMult(svd,PETSC_FALSE,cyclic->x2,cyclic->y1);
35: SVDMatMult(svd,PETSC_TRUE,cyclic->x1,cyclic->y2);
36: VecCUDAResetArray(cyclic->x1);
37: VecCUDAResetArray(cyclic->x2);
38: VecCUDAResetArray(cyclic->y1);
39: VecCUDAResetArray(cyclic->y2);
40: VecCUDARestoreArrayRead(x,(const PetscScalar**)&d_px);
41: VecCUDARestoreArray(y,&d_py);
42: return(0);
43: }
45: PetscErrorCode MatCreateVecs_Cyclic_CUDA(Mat B,Vec *right,Vec *left)
46: {
48: SVD svd;
49: SVD_CYCLIC *cyclic;
50: PetscInt M,N,m,n;
51: PetscMPIInt size;
54: MatShellGetContext(B,(void**)&svd);
55: cyclic = (SVD_CYCLIC*)svd->data;
56: SVDMatGetSize(svd,&M,&N);
57: SVDMatGetLocalSize(svd,&m,&n);
58: MPI_Comm_size(PetscObjectComm((PetscObject)cyclic->mat),&size);
59: if (right) {
60: VecCreate(PetscObjectComm((PetscObject)cyclic->mat),right);
61: VecSetSizes(*right,m+n,M+N);
62: if (size>1) {
63: VecSetType(*right,VECMPICUDA);
64: } else {
65: VecSetType(*right,VECSEQCUDA);
66: }
67: }
68: if (left) {
69: VecCreate(PetscObjectComm((PetscObject)cyclic->mat),left);
70: VecSetSizes(*left,m+n,M+N);
71: if (size>1) {
72: VecSetType(*left,VECMPICUDA);
73: } else {
74: VecSetType(*left,VECSEQCUDA);
75: }
76: }
77: return(0);
78: }