mirror of
https://gitlab.com/scemama/qp_plugins_scemama.git
synced 2024-08-30 15:53:46 +02:00
OpenMP in random svd
This commit is contained in:
parent
16d5b14f36
commit
9ee6182443
@ -13,8 +13,7 @@ subroutine run
|
|||||||
include 'constants.include.F'
|
include 'constants.include.F'
|
||||||
double precision, allocatable :: U(:,:), V(:,:), D(:), A(:,:)
|
double precision, allocatable :: U(:,:), V(:,:), D(:), A(:,:)
|
||||||
integer :: i, j, k, l, q, r, m, n, iter
|
integer :: i, j, k, l, q, r, m, n, iter
|
||||||
double precision,allocatable :: Z(:,:), P(:,:), Yt(:,:), UYt(:,:)
|
double precision,allocatable :: Z(:,:), P(:,:), Yt(:,:), UYt(:,:), r1(:,:)
|
||||||
double precision :: r1,r2
|
|
||||||
|
|
||||||
m = n_det_alpha_unique
|
m = n_det_alpha_unique
|
||||||
n = n_det_beta_unique
|
n = n_det_beta_unique
|
||||||
@ -25,52 +24,75 @@ subroutine run
|
|||||||
|
|
||||||
! Z(m,r) = A(m,n).P(n,r)
|
! Z(m,r) = A(m,n).P(n,r)
|
||||||
Z(:,:) = 0.d0
|
Z(:,:) = 0.d0
|
||||||
|
allocate(P(n,r))
|
||||||
|
|
||||||
|
!$OMP PARALLEL DEFAULT(SHARED) PRIVATE(i,j,k,l,r1)
|
||||||
|
allocate(r1(N_det,2))
|
||||||
|
!$OMP DO
|
||||||
do l=1,r
|
do l=1,r
|
||||||
|
call random_number(r1)
|
||||||
|
r1(:,1) = dsqrt(-2.d0*dlog(r1(:,1)))
|
||||||
|
r1(:,1) = r1(:,1) * dcos(dtwo_pi*r1(:,2))
|
||||||
do k=1,N_det
|
do k=1,N_det
|
||||||
i = psi_bilinear_matrix_rows(k)
|
i = psi_bilinear_matrix_rows(k)
|
||||||
j = psi_bilinear_matrix_columns(k)
|
j = psi_bilinear_matrix_columns(k)
|
||||||
call random_number(r1)
|
Z(i,l) = Z(i,l) + psi_bilinear_matrix_values(k,1) * r1(k,1)
|
||||||
call random_number(r2)
|
|
||||||
r1 = dsqrt(-2.d0*dlog(r1))
|
|
||||||
r2 = dtwo_pi*r2
|
|
||||||
Z(i,l) = Z(i,l) + psi_bilinear_matrix_values(k,1) * r1*dcos(r2)
|
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
!$OMP END DO
|
||||||
|
deallocate(r1)
|
||||||
|
!$OMP END PARALLEL
|
||||||
|
|
||||||
! Power iterations
|
! Power iterations
|
||||||
allocate(P(n,r))
|
|
||||||
do iter=1,20
|
do iter=1,20
|
||||||
|
print *, 'Power iteration ', iter, '/', 20
|
||||||
|
|
||||||
! P(n,r) = At(n,m).Z(m,r)
|
! P(n,r) = At(n,m).Z(m,r)
|
||||||
P(:,:) = 0.d0
|
!$OMP PARALLEL DEFAULT(SHARED) PRIVATE(i,j,k,l)
|
||||||
|
!$OMP DO
|
||||||
do l=1,r
|
do l=1,r
|
||||||
|
P(:,l) = 0.d0
|
||||||
do k=1,N_det
|
do k=1,N_det
|
||||||
i = psi_bilinear_matrix_rows(k)
|
i = psi_bilinear_matrix_rows(k)
|
||||||
j = psi_bilinear_matrix_columns(k)
|
j = psi_bilinear_matrix_columns(k)
|
||||||
P(j,l) = P(j,l) + psi_bilinear_matrix_values(k,1) * Z(i,l)
|
P(j,l) = P(j,l) + psi_bilinear_matrix_values(k,1) * Z(i,l)
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
Z(:,:) = 0.d0
|
!$OMP END DO
|
||||||
|
|
||||||
|
!$OMP BARRIER
|
||||||
|
|
||||||
|
!$OMP DO
|
||||||
do l=1,r
|
do l=1,r
|
||||||
|
Z(:,l) = 0.d0
|
||||||
do k=1,N_det
|
do k=1,N_det
|
||||||
i = psi_bilinear_matrix_rows(k)
|
i = psi_bilinear_matrix_rows(k)
|
||||||
j = psi_bilinear_matrix_columns(k)
|
j = psi_bilinear_matrix_columns(k)
|
||||||
Z(i,l) = Z(i,l) + psi_bilinear_matrix_values(k,1) * P(j,l)
|
Z(i,l) = Z(i,l) + psi_bilinear_matrix_values(k,1) * P(j,l)
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
!$OMP END DO
|
||||||
|
|
||||||
|
!$OMP END PARALLEL
|
||||||
|
|
||||||
! Compute QR
|
! Compute QR
|
||||||
call ortho_qr(Z,size(Z,1),m,r)
|
call ortho_qr(Z,size(Z,1),m,r)
|
||||||
enddo
|
enddo
|
||||||
|
|
||||||
! Y(r,n) = Zt(r,m).A(m,n)
|
! Y(r,n) = Zt(r,m).A(m,n)
|
||||||
allocate(Yt(n,r))
|
allocate(Yt(n,r))
|
||||||
Yt(:,:) = 0.d0
|
!$OMP PARALLEL DEFAULT(SHARED) PRIVATE(i,j,k,l)
|
||||||
|
!$OMP DO
|
||||||
do l=1,r
|
do l=1,r
|
||||||
|
Yt(:,l) = 0.d0
|
||||||
do k=1,N_det
|
do k=1,N_det
|
||||||
i = psi_bilinear_matrix_rows(k)
|
i = psi_bilinear_matrix_rows(k)
|
||||||
j = psi_bilinear_matrix_columns(k)
|
j = psi_bilinear_matrix_columns(k)
|
||||||
Yt(j,l) = Yt(j,l) + Z(i,l) * psi_bilinear_matrix_values(k,1)
|
Yt(j,l) = Yt(j,l) + Z(i,l) * psi_bilinear_matrix_values(k,1)
|
||||||
enddo
|
enddo
|
||||||
enddo
|
enddo
|
||||||
|
!$OMP END DO
|
||||||
|
!$OMP END PARALLEL
|
||||||
|
|
||||||
allocate(D(r),V(n,r), UYt(r,r))
|
allocate(D(r),V(n,r), UYt(r,r))
|
||||||
call svd(Yt,size(Yt,1),V,size(V,1),D,UYt,size(UYt,1),n,r)
|
call svd(Yt,size(Yt,1),V,size(V,1),D,UYt,size(UYt,1),n,r)
|
||||||
|
Loading…
Reference in New Issue
Block a user