scemama/qp_plugins_scemama
1
0
Fork 0
mirror of https://gitlab.com/scemama/qp_plugins_scemama.git synced 2024-12-22 04:13:40 +01:00
Code Issues Releases Wiki Activity

H_oo on GPU

Browse Source
This commit is contained in:
Anthony Scemama 2023-08-03 19:46:41 +02:00
parent 2ea84d4ffc
commit ecd9531ee6
6 changed files with 163 additions and 19 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

7
devel/ccsd_gpu/ccsd_space_orb_sub.irp.f
View File

@ -120,17 +120,18 @@ subroutine run_ccsd_space_orb
gpu_data = gpu_init(nO, nV, cholesky_mo_num, &
cc_space_v_oo_chol, cc_space_v_ov_chol, cc_space_v_vo_chol, cc_space_v_vv_chol, &
cc_space_v_oooo, cc_space_v_vooo, cc_space_v_oovv, cc_space_v_vvoo, &
cc_space_v_oovo, cc_space_v_ovvo, cc_space_v_ovov, cc_space_v_ovoo, cc_space_f_vo)
cc_space_v_oovo, cc_space_v_ovvo, cc_space_v_ovov, cc_space_v_ovoo, &
cc_space_f_oo, cc_space_f_vo, cc_space_f_vv)
do while (not_converged)
! Residue
if (do_ao_cholesky) then
call compute_H_oo_chol(nO,nV,tau_x,H_oo)
call compute_H_vv_chol(nO,nV,tau_x,H_vv)
call compute_H_vo_chol(nO,nV,t1,H_vo)
call gpu_upload(gpu_data, nO, nV, t1, t2, tau, tau_x, H_vv);
call compute_H_oo_chol_gpu(gpu_data,nO,nV,0,H_oo)
call gpu_upload(gpu_data, nO, nV, t1, t2, tau, tau_x, H_oo, H_vv);
call compute_r1_space_chol(nO,nV,t1,t2,tau,H_oo,H_vv,H_vo,r1,max_r1)
call compute_r2_space_chol_gpu(gpu_data, nO, nV, t1, r2, max_r2)
else

38
devel/ccsd_gpu/ccsd_space_orb_sub_chol.irp.f
View File

@ -439,3 +439,41 @@ subroutine compute_H_vo_chol(nO,nV,t1,H_vo)
end
subroutine compute_H_oo_chol2(nO,nV,tau_x,H_oo)
implicit none
integer, intent(in) :: nO,nV
double precision, intent(in) :: tau_x(nO, nO, nV, nV)
double precision, intent(out) :: H_oo(nO, nO)
integer :: a,b,i,j,u,k
double precision, allocatable :: tau_kau(:,:,:), tmp_vov(:,:,:), tmp_ovv(:,:,:)
allocate(tau_kau(cholesky_mo_num,nV,nO))
allocate(tmp_vov(nV,nO,nV) )
allocate(tmp_ovv(nO,nV,nV) )
do u = 1, nO
call dcopy(nO*nV*nV, tau_x(u,1,1,1), nO, tmp_ovv, 1)
print *, u
print *, tmp_ovv
do b=1,nV
do j=1,nO
do a=1,nV
tmp_vov(a,j,b) = tmp_ovv(j,a,b)
enddo
enddo
enddo
call dgemm('N','T',cholesky_mo_num,nV,nO*nV,1.d0, &
cc_space_v_ov_chol, cholesky_mo_num, tmp_vov, nV, &
0.d0, tau_kau(1,1,u), cholesky_mo_num)
enddo
deallocate(tmp_vov)
call dcopy(nO*nO, cc_space_f_oo, 1, H_oo, 1);
call dgemm('T', 'N', nO, nO, cholesky_mo_num*nV, 1.d0, &
tau_kau, cholesky_mo_num*nV, cc_space_v_vo_chol, cholesky_mo_num*nV, &
1.d0, H_oo, nO)
end

98
devel/ccsd_gpu/gpu.c
View File

@ -12,7 +12,6 @@ void gpu_upload(gpu_data* data,
double* t2,
double* tau,
double* tau_x,
double* H_oo,
double* H_vv)
{
int lda;
@ -21,7 +20,6 @@ void gpu_upload(gpu_data* data,
int ngpus = 1;
if (MULTIGPU == 1) cudaGetDeviceCount(&ngpus);
#pragma omp parallel num_threads(ngpus)
{
int igpu = omp_get_thread_num();
@ -43,10 +41,6 @@ void gpu_upload(gpu_data* data,
lda = nO*nO;
cublasSetMatrix(nO*nO, nV*nV, sizeof(double), t2, lda, d_t2, lda);
double* d_H_oo = data[igpu].H_oo;
lda = nO;
cublasSetMatrix(nO, nO, sizeof(double), H_oo, lda, d_H_oo, lda);
double* d_H_vv = data[igpu].H_vv;
lda = nV;
cublasSetMatrix(nV, nV, sizeof(double), H_vv, lda, d_H_vv, lda);
@ -75,8 +69,8 @@ void compute_r2_space_chol_gpu(gpu_data* data, int nO, int nV, double* t1, doubl
int igpu = omp_get_thread_num();
cudaSetDevice(igpu);
cublasHandle_t handle;
cublasHandle_t handle;
cublasCreate(&handle);
double* d_r2;
@ -97,7 +91,9 @@ void compute_r2_space_chol_gpu(gpu_data* data, int nO, int nV, double* t1, doubl
double* d_cc_space_v_ovvo = data[igpu].cc_space_v_ovvo;
double* d_cc_space_v_ovov = data[igpu].cc_space_v_ovov;
double* d_cc_space_v_ovoo = data[igpu].cc_space_v_ovoo;
double* d_cc_space_f_oo = data[igpu].cc_space_f_oo;
double* d_cc_space_f_vo = data[igpu].cc_space_f_vo;
double* d_cc_space_f_vv = data[igpu].cc_space_f_vv;
double* d_tau = data[igpu].tau;
double* d_t1 = data[igpu].t1;
double* d_t2 = data[igpu].t2;
@ -1353,3 +1349,91 @@ void compute_r2_space_chol_gpu(gpu_data* data, int nO, int nV, double* t1, doubl
}
void compute_h_oo_chol_gpu(gpu_data* data, int nO, int nV, int igpu, double* H_oo)
{
int ngpus = 1;
if (MULTIGPU == 1) cudaGetDeviceCount(&ngpus);
const int cholesky_mo_num = data[igpu].cholesky_mo_num;
cudaSetDevice(igpu);
int m,n,k, lda, ldb, ldc;
double alpha, beta;
double* A;
double* B;
double* C;
cudaStream_t stream[nV];
cublasHandle_t handle;
cublasCreate(&handle);
double* d_H_oo = data[igpu].H_oo;
double* d_tau_x = data[igpu].tau_x;
double* d_cc_space_f_oo = data[igpu].cc_space_f_oo;
double* d_cc_space_v_vo_chol = data[igpu].cc_space_v_vo_chol;
double* d_cc_space_v_ov_chol = data[igpu].cc_space_v_ov_chol;
double* d_tau_kau;
cudaMalloc((void **)&d_tau_kau, cholesky_mo_num*nV*nO * sizeof(double));
double* d_tmp_ovv;
cudaMalloc((void **)&d_tmp_ovv, nO*nV*nV * sizeof(double));
double* d_tmp_vov;
cudaMalloc((void **)&d_tmp_vov, nV*nO*nV * sizeof(double));
alpha = 1.0;
beta = 0.0;
for (int u=0 ; u<nO ; ++u) {
cublasDcopy(handle, nO*nV*nV, &(d_tau_x[u]), nO, d_tmp_ovv, 1);
for (int i=0 ; i<nV ; ++i) {
cudaStreamCreate(&(stream[i]));
}
for (int b=0 ; b<nV ; ++b) {
cublasSetStream(handle, stream[b]);
A = &(d_tmp_ovv[nO*nV*b]); lda = nO;
B = &(d_tmp_ovv[nO*nV*b]); ldb = nO;
C = &(d_tmp_vov[nV*nO*b]); ldc = nV;
cublasDgeam(handle, CUBLAS_OP_T, CUBLAS_OP_T, nV, nO, &alpha, A, lda, &beta, B, ldb, C, ldc);
}
for (int i=0 ; i<nV ; ++i) {
cudaStreamDestroy(stream[i]);
}
cublasSetStream(handle, NULL);
alpha = 1.0;
beta = 0.0;
m=cholesky_mo_num; n=nV; k=nO*nV;
A=d_cc_space_v_ov_chol; lda=cholesky_mo_num;
B=d_tmp_vov; ldb=nV;
C=&(d_tau_kau[cholesky_mo_num*nV*u]); ldc=cholesky_mo_num;
cublasDgemm(handle, CUBLAS_OP_N, CUBLAS_OP_T, m, n, k, &alpha, A, lda, B, ldb, &beta, C, ldc);
}
cudaFree(d_tmp_vov);
cudaFree(d_tmp_ovv);
cublasDcopy(handle, nO*nO, d_cc_space_f_oo, 1, d_H_oo, 1);
alpha = 1.0;
beta = 1.0;
m=nO; n=nO; k=cholesky_mo_num*nV;
A=d_tau_kau; lda=cholesky_mo_num*nV;
B=d_cc_space_v_vo_chol; ldb=cholesky_mo_num*nV;
C=d_H_oo; ldc=nO;
cublasDgemm(handle, CUBLAS_OP_T, CUBLAS_OP_N, m, n, k, &alpha, A, lda, B, ldb, &beta, C, ldc);
cudaFree(d_tau_kau);
// double* H_oo = malloc(nO*nO*sizeof(double));
cublasGetMatrix(nO, nO, sizeof(double), d_H_oo, nO, H_oo, nO);
for (int i=0 ; i<ngpus ; ++i) {
if (i != igpu) {
double* d_H_oo = data[i].H_oo;
cudaSetDevice(i);
cublasSetMatrix(nO, nO, sizeof(double), H_oo, nO, d_H_oo, nO);
}
}
// free(H_oo);
cublasDestroy(handle);
}

2
devel/ccsd_gpu/gpu.h
View File

@ -11,7 +11,9 @@ typedef struct {
double* cc_space_v_ovvo;
double* cc_space_v_ovov;
double* cc_space_v_ovoo;
double* cc_space_f_oo;
double* cc_space_f_vo;
double* cc_space_f_vv;
double* tau;
double* tau_x;
double* t1;

13
devel/ccsd_gpu/gpu_init.c
View File

@ -14,7 +14,8 @@ gpu_data* gpu_init(
double* cc_space_v_oovv, double* cc_space_v_vvoo,
double* cc_space_v_oovo, double* cc_space_v_ovvo,
double* cc_space_v_ovov, double* cc_space_v_ovoo,
double* cc_space_f_vo)
double* cc_space_f_oo, double* cc_space_f_vo,
double* cc_space_f_vv)
{
int ngpus = 1;
cudaGetDeviceCount(&ngpus);
@ -86,10 +87,18 @@ gpu_data* gpu_init(
cudaMalloc((void **)&d_cc_space_v_ovoo, nO*nV*nO*nO * sizeof(double));
cublasSetMatrix(lda, nO*nO, sizeof(double), cc_space_v_ovoo, lda, d_cc_space_v_ovoo, lda);
double* d_cc_space_f_oo;
cudaMalloc((void**)&d_cc_space_f_oo, nO*nO*sizeof(double));
cublasSetMatrix(nO, nO, sizeof(double), cc_space_f_oo, nO, d_cc_space_f_oo, nO);
double* d_cc_space_f_vo;
cudaMalloc((void**)&d_cc_space_f_vo, nV*nO*sizeof(double));
cublasSetMatrix(nV, nO, sizeof(double), cc_space_f_vo, nV, d_cc_space_f_vo, nV);
double* d_cc_space_f_vv;
cudaMalloc((void**)&d_cc_space_f_vv, nV*nV*sizeof(double));
cublasSetMatrix(nV, nV, sizeof(double), cc_space_f_vv, nV, d_cc_space_f_vv, nV);
double* d_tau;
lda = nO * nO;
cudaMalloc((void **)&d_tau, lda * nV * nV * sizeof(double));
@ -126,7 +135,9 @@ gpu_data* gpu_init(
data[igpu].cc_space_v_ovvo = d_cc_space_v_ovvo;
data[igpu].cc_space_v_ovov = d_cc_space_v_ovov;
data[igpu].cc_space_v_ovoo = d_cc_space_v_ovoo;
data[igpu].cc_space_f_oo = d_cc_space_f_oo;
data[igpu].cc_space_f_vo = d_cc_space_f_vo;
data[igpu].cc_space_f_vv = d_cc_space_f_vv;
data[igpu].tau = d_tau;
data[igpu].tau_x = d_tau_x;
data[igpu].t1 = d_t1;

14
devel/ccsd_gpu/gpu_module.f90
View File

@ -8,7 +8,7 @@ module gpu_module
cc_space_v_oo_chol, cc_space_v_ov_chol, cc_space_v_vo_chol, cc_space_v_vv_chol, &
cc_space_v_oooo, cc_space_v_vooo, cc_space_v_oovv, cc_space_v_vvoo, &
cc_space_v_oovo, cc_space_v_ovvo, cc_space_v_ovov, cc_space_v_ovoo, &
cc_space_f_vo) bind(C)
cc_space_f_oo, cc_space_f_vo, cc_space_f_vv) bind(C)
import c_int, c_double, c_ptr
integer(c_int), intent(in), value :: nO, nV, cholesky_mo_num
real(c_double), intent(in) :: cc_space_v_oo_chol(cholesky_mo_num,nO,nO)
@ -23,10 +23,12 @@ module gpu_module
real(c_double), intent(in) :: cc_space_v_ovvo(nO,nV,nV,nO)
real(c_double), intent(in) :: cc_space_v_ovov(nO,nV,nO,nV)
real(c_double), intent(in) :: cc_space_v_ovoo(nO,nV,nO,nO)
real(c_double), intent(in) :: cc_space_f_oo(nO,nO)
real(c_double), intent(in) :: cc_space_f_vo(nV,nO)
real(c_double), intent(in) :: cc_space_f_vv(nV,nV)
end function
subroutine gpu_upload(gpu_data, nO, nV, t1, t2, tau, tau_x, H_oo, H_vv) bind(C)
subroutine gpu_upload(gpu_data, nO, nV, t1, t2, tau, tau_x, H_vv) bind(C)
import c_int, c_double, c_ptr
type(c_ptr), value :: gpu_data
integer(c_int), intent(in), value :: nO, nV
@ -34,10 +36,16 @@ module gpu_module
real(c_double), intent(in) :: t2(nO,nO,nV,nV)
real(c_double), intent(in) :: tau(nO,nO,nV,nV)
real(c_double), intent(in) :: tau_x(nO,nO,nV,nV)
real(c_double), intent(in) :: H_oo(nO,nO)
real(c_double), intent(in) :: H_vv(nV,nV)
end subroutine
subroutine compute_H_oo_chol_gpu(gpu_data, nO, nV, igpu, H_oo) bind(C)
import c_int, c_double, c_ptr
type(c_ptr), value :: gpu_data
integer(c_int), intent(in), value :: nO, nV, igpu
real(c_double), intent(out) :: H_oo(nO,nO)
end subroutine
subroutine compute_r2_space_chol_gpu(gpu_data, nO, nV, t1, r2, max_r2) bind(C)
import c_int, c_double, c_ptr
type(c_ptr), value :: gpu_data