Skip to content

Commit b62ef9a

Browse files
brunowubrunowu
committed
restruct lanczos for CUDA-SEQ
1 parent d453ba9 commit b62ef9a

2 files changed

Lines changed: 137 additions & 21 deletions

File tree

‎ChASE-MPI/impl/chase_mpidla_cuda_seq.hpp‎

Lines changed: 119 additions & 21 deletions
Original file line numberDiff line numberDiff line change
@@ -28,6 +28,52 @@
2828
* multiGPUs (both block-block and block-cyclic distributions)
2929
* @{
3030
*/
31+
//! generate `n` random float numbers in normal distribution on each GPU device.
32+
//!
33+
//! @param[in] seed the seed of random number generator
34+
//! @param[in] states the states of the sequence of random number generator
35+
//! @param[in,out] v a pointer to the device memory to store the random
36+
//! generated numbers
37+
//! @param[in] stream_ an asynchronous CUDA stream which allows to run this
38+
//! function asynchronously
39+
void chase_rand_normal(unsigned long long seed, curandStatePhilox4_32_10_t* states, float* v,
40+
int n, cudaStream_t stream_);
41+
//! generate `n` random double numbers in normal distribution on each GPU
42+
//! device.
43+
//!
44+
//! @param[in] seed the seed of random number generator
45+
//! @param[in] states the states of the sequence of random number generator
46+
//! @param[in,out] v a pointer to the device memory to store the random
47+
//! generated numbers
48+
//! @param[in] stream_ an asynchronous CUDA stream which allows to run this
49+
//! function asynchronously
50+
void chase_rand_normal(unsigned long long seed, curandStatePhilox4_32_10_t* states, double* v,
51+
int n, cudaStream_t stream_);
52+
//! generate `n` random complex float numbers in normal distribution on each GPU
53+
//! device. The real part and the imaginary part of each individual random
54+
//! number are the same.
55+
//!
56+
//! @param[in] seed the seed of random number generator
57+
//! @param[in] states the states of the sequence of random number generator
58+
//! @param[in,out] v a pointer to the device memory to store the random
59+
//! generated numbers
60+
//! @param[in] stream_ an asynchronous CUDA stream which allows to run this
61+
//! function asynchronously
62+
void chase_rand_normal(unsigned long long seed, curandStatePhilox4_32_10_t* states,
63+
std::complex<float>* v, int n, cudaStream_t stream_);
64+
//! generate `n` random complex double numbers in normal distribution on each
65+
//! GPU device. The real part and the imaginary part of each individual random
66+
//! number are the same.
67+
//!
68+
//! @param[in] seed the seed of random number generator
69+
//! @param[in] states the states of the sequence of random number generator
70+
//! @param[in,out] v a pointer to the device memory to store the random
71+
//! generated numbers
72+
//! @param[in] stream_ an asynchronous CUDA stream which allows to run this
73+
//! function asynchronously
74+
void chase_rand_normal(unsigned long long seed, curandStatePhilox4_32_10_t* states,
75+
std::complex<double>* v, int n, cudaStream_t stream_);
76+
3177

3278
//! shift the diagonal of a `nxn` square matrix `A` in float real data type on a
3379
//! single GPU.
@@ -96,7 +142,6 @@ class ChaseMpiDLACudaSeq : public ChaseMpiDLAInterface<T>
96142
ldh_(matrices.get_ldh())
97143
{
98144
cuda_exec(cudaSetDevice(0));
99-
100145
cuda_exec(cudaMalloc((void**)&(d_V1_), N_ * (nev_ + nex_) * sizeof(T)));
101146
cuda_exec(cudaMalloc((void**)&(d_V2_), N_ * (nev_ + nex_) * sizeof(T)));
102147
cuda_exec(cudaMalloc((void**)&(d_H_), N_ * N_ * sizeof(T)));
@@ -195,41 +240,28 @@ class ChaseMpiDLACudaSeq : public ChaseMpiDLAInterface<T>
195240
}
196241
void initVecs() override
197242
{
198-
cuda_exec(cudaMemcpy(d_V1_, V1_, (nev_ + nex_) * N_ * sizeof(T),
199-
cudaMemcpyHostToDevice));
200-
cuda_exec(cudaMemcpy(d_V2_, d_V1_, (nev_ + nex_) * N_ * sizeof(T),
243+
cuda_exec(cudaMemcpy(d_V2_, d_V1_, (nev_ + nex_) * N_ * sizeof(T),
201244
cudaMemcpyDeviceToDevice));
202-
//cuda_exec(
203-
// cudaMemcpy(d_H_, H_, N_ * N_ * sizeof(T), cudaMemcpyHostToDevice));
204245
cublasSetMatrix(N_, N_, sizeof(T), H_, ldh_, d_H_, N_);
205246
}
206247
void initRndVecs() override
207248
{
208-
209-
std::mt19937 gen(1337.0);
210-
std::normal_distribution<> d;
211-
for (auto j = 0; j < (nev_ + nex_); j++)
212-
{
213-
for (auto i = 0; i < N_; i++)
214-
{
215-
V1_[i + j * N_] = getRandomT<T>([&]() { return d(gen); });
216-
}
217-
}
249+
unsigned long long seed = 24141;
250+
chase_rand_normal(seed, states_, d_V1_, N_ * (nev_ + nex_),
251+
(cudaStream_t)0);
218252
}
219253

220-
// host->device: v1 on host, v2 on device
221254
void V2C(T* v1, std::size_t off1, T* v2, std::size_t off2,
222255
std::size_t block) override
223256
{
224257
cuda_exec(cudaMemcpy(v2 + off2 * N_, v1 + off1 * N_,
225-
block * N_ * sizeof(T), cudaMemcpyHostToDevice));
258+
block * N_ * sizeof(T), cudaMemcpyDeviceToDevice));
226259
}
227-
// device->host: v1 on device, v2 on host
228260
void C2V(T* v1, std::size_t off1, T* v2, std::size_t off2,
229261
std::size_t block) override
230262
{
231263
cuda_exec(cudaMemcpy(v2 + off2 * N_, v1 + off1 * N_,
232-
block * N_ * sizeof(T), cudaMemcpyDeviceToHost));
264+
block * N_ * sizeof(T), cudaMemcpyDeviceToDevice));
233265
}
234266

235267
void preApplication(T* V, std::size_t locked, std::size_t block) override
@@ -442,7 +474,73 @@ class ChaseMpiDLACudaSeq : public ChaseMpiDLAInterface<T>
442474
cudaMemcpyDeviceToDevice));
443475
}
444476
void Lanczos(std::size_t M, int idx, Base<T>* d, Base<T>* e, Base<T> *r_beta) override
445-
{}
477+
{
478+
Base<T> real_beta;
479+
480+
T alpha = T(1.0);
481+
T beta = T(0.0);
482+
483+
cudaMemset(v0_, 0, sizeof(T) * N_);
484+
485+
#ifdef USE_NSIGHT
486+
nvtxRangePushA("C2V");
487+
#endif
488+
if(idx >= 0)
489+
{
490+
this->C2V(d_V2_, idx, v1_, 0, 1);
491+
}else
492+
{
493+
unsigned long long seed = 2342;
494+
chase_rand_normal(seed, states_, v1_, N_, (cudaStream_t)0);
495+
}
496+
497+
#ifdef USE_NSIGHT
498+
nvtxRangePop();
499+
#endif
500+
// ENSURE that v1 has one norm
501+
#ifdef USE_NSIGHT
502+
nvtxRangePushA("Lanczos: loop");
503+
#endif
504+
Base<T> real_alpha = this->nrm2(N_, v1_, 1);
505+
alpha = T(1 / real_alpha);
506+
this->scal(N_, &alpha, v1_, 1);
507+
for (std::size_t k = 0; k < M; k = k + 1)
508+
{
509+
if(idx >= 0){
510+
this->V2C(v1_, 0, d_V1_, k, 1);
511+
}
512+
this->applyVec(v1_, w_);
513+
alpha = this->dot(N_, v1_, 1, w_, 1);
514+
alpha = -alpha;
515+
this->axpy(N_, &alpha, v1_, 1, w_, 1);
516+
alpha = -alpha;
517+
518+
d[k] = std::real(alpha);
519+
520+
if (k == M - 1)
521+
break;
522+
523+
beta = T(-real_beta);
524+
this->axpy(N_, &beta, v0_, 1, w_, 1);
525+
beta = -beta;
526+
527+
real_beta = this->nrm2(N_, w_, 1);
528+
529+
beta = T(1.0 / real_beta);
530+
531+
this->scal(N_, &beta, w_, 1);
532+
533+
e[k] = real_beta;
534+
535+
std::swap(v1_, v0_);
536+
std::swap(v1_, w_);
537+
}
538+
#ifdef USE_NSIGHT
539+
nvtxRangePop();
540+
#endif
541+
*r_beta = real_beta;
542+
543+
}
446544

447545
void B2C(T* B, std::size_t off1, T* C, std::size_t off2, std::size_t block) override
448546
{}

‎cmake/FindSCALAPACK.cmake‎

Lines changed: 18 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -354,6 +354,24 @@ if(BLAS_FOUND
354354
)
355355
endif(NOT SCALAPACK_LIBRARIES)
356356
endif()
357+
358+
#Fujitsu ScaLAPACK
359+
if(BLA_VENDOR STREQUAL "Fujitsu_SSL2BLAMP" OR BLA_VENDOR STREQUAL "All")
360+
message(STATUS "generic ")
361+
if(NOT SCALAPACK_LIBRARIES)
362+
check_scalapack_libraries(
363+
SCALAPACK_LIBRARIES
364+
SCALAPACK
365+
pdgemm
366+
"-Kopenmp"
367+
"fjscalapack" # scalapack lib to look for
368+
"${LAPACK_LIBRARIES};${BLAS_LIBRARIES}" # blas and lapack libs
369+
"${MPI_Fortran_LIBRARIES}" # mpi libs
370+
"" # threads libs
371+
)
372+
endif(NOT SCALAPACK_LIBRARIES)
373+
endif()
374+
357375
# intel scalapack
358376
if(BLA_VENDOR MATCHES "Intel" OR BLA_VENDOR STREQUAL "All")
359377
if(UNIX AND NOT WIN32)

0 commit comments

Comments
 (0)