Actual source code: cycliccuda.cu

slepc-3.11.2 2019-07-30
Report Typos and Errors
  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: }