LCPQ/qp2
9
1
Fork 0
mirror of https://github.com/QuantumPackage/qp2.git synced 2024-12-22 11:33:29 +01:00
Code Releases Activity

OpenMP in Hessian

Browse Source
This commit is contained in:
Anthony Scemama 2019-07-06 01:04:50 +02:00
parent a744bc30d4
commit 970fd8837a
1 changed files with 60 additions and 34 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

84
src/casscf/hessian.irp.f
View File

@ -189,7 +189,7 @@ BEGIN_PROVIDER [real*8, hessmat2, (nMonoEx,nMonoEx)]
! !
END_DOC END_DOC
implicit none implicit none
integer :: i,j,t,u,a,b,indx,jndx,bstart,ustart integer :: i,j,t,u,a,b,indx,jndx,bstart,ustart,indx_shift
real*8 :: hessmat_itju real*8 :: hessmat_itju
real*8 :: hessmat_itja real*8 :: hessmat_itja
@ -203,9 +203,14 @@ BEGIN_PROVIDER [real*8, hessmat2, (nMonoEx,nMonoEx)]
write(6,*) ' nMonoEx = ',nMonoEx write(6,*) ' nMonoEx = ',nMonoEx
endif endif
indx=1 !$OMP PARALLEL DEFAULT(NONE) &
!$OMP SHARED(hessmat2,n_core_inact_orb,n_act_orb,n_virt_orb,nMonoEx) &
!$OMP PRIVATE(i,indx,jndx,j,ustart,t,u,a,bstart,indx_shift)
!$OMP DO
do i=1,n_core_inact_orb do i=1,n_core_inact_orb
do t=1,n_act_orb do t=1,n_act_orb
indx = t + (i-1)*n_act_orb
jndx=indx jndx=indx
do j=i,n_core_inact_orb do j=i,n_core_inact_orb
if (i.eq.j) then if (i.eq.j) then
@ -214,31 +219,31 @@ BEGIN_PROVIDER [real*8, hessmat2, (nMonoEx,nMonoEx)]
ustart=1 ustart=1
end if end if
do u=ustart,n_act_orb do u=ustart,n_act_orb
hessmat2(indx,jndx)=hessmat_itju(i,t,j,u) hessmat2(jndx,indx)=hessmat_itju(i,t,j,u)
hessmat2(jndx,indx)=hessmat2(indx,jndx)
jndx+=1 jndx+=1
end do end do
end do end do
do j=1,n_core_inact_orb do j=1,n_core_inact_orb
do a=1,n_virt_orb do a=1,n_virt_orb
hessmat2(indx,jndx)=hessmat_itja(i,t,j,a) hessmat2(jndx,indx)=hessmat_itja(i,t,j,a)
hessmat2(jndx,indx)=hessmat2(indx,jndx)
jndx+=1 jndx+=1
end do end do
end do end do
do u=1,n_act_orb do u=1,n_act_orb
do a=1,n_virt_orb do a=1,n_virt_orb
hessmat2(indx,jndx)=hessmat_itua(i,t,u,a) hessmat2(jndx,indx)=hessmat_itua(i,t,u,a)
hessmat2(jndx,indx)=hessmat2(indx,jndx)
jndx+=1 jndx+=1
end do end do
end do end do
indx+=1
end do end do
end do end do
!$OMP END DO NOWAIT
do i=1,n_core_inact_orb indx_shift = n_core_inact_orb*n_act_orb
!$OMP DO
do a=1,n_virt_orb do a=1,n_virt_orb
do i=1,n_core_inact_orb
indx = a + (i-1)*n_virt_orb + indx_shift
jndx=indx jndx=indx
do j=i,n_core_inact_orb do j=i,n_core_inact_orb
if (i.eq.j) then if (i.eq.j) then
@ -247,24 +252,25 @@ BEGIN_PROVIDER [real*8, hessmat2, (nMonoEx,nMonoEx)]
bstart=1 bstart=1
end if end if
do b=bstart,n_virt_orb do b=bstart,n_virt_orb
hessmat2(indx,jndx)=hessmat_iajb(i,a,j,b) hessmat2(jndx,indx)=hessmat_iajb(i,a,j,b)
hessmat2(jndx,indx)=hessmat2(indx,jndx)
jndx+=1 jndx+=1
end do end do
end do end do
do t=1,n_act_orb do t=1,n_act_orb
do b=1,n_virt_orb do b=1,n_virt_orb
hessmat2(indx,jndx)=hessmat_iatb(i,a,t,b) hessmat2(jndx,indx)=hessmat_iatb(i,a,t,b)
hessmat2(jndx,indx)=hessmat2(indx,jndx)
jndx+=1 jndx+=1
end do end do
end do end do
indx+=1
end do end do
end do end do
!$OMP END DO NOWAIT
do t=1,n_act_orb indx_shift += n_core_inact_orb*n_virt_orb
!$OMP DO
do a=1,n_virt_orb do a=1,n_virt_orb
do t=1,n_act_orb
indx = a + (t-1)*n_virt_orb + indx_shift
jndx=indx jndx=indx
do u=t,n_act_orb do u=t,n_act_orb
if (t.eq.u) then if (t.eq.u) then
@ -273,14 +279,22 @@ BEGIN_PROVIDER [real*8, hessmat2, (nMonoEx,nMonoEx)]
bstart=1 bstart=1
end if end if
do b=bstart,n_virt_orb do b=bstart,n_virt_orb
hessmat2(indx,jndx)=hessmat_taub(t,a,u,b) hessmat2(jndx,indx)=hessmat_taub(t,a,u,b)
hessmat2(jndx,indx)=hessmat2(indx,jndx)
jndx+=1 jndx+=1
end do end do
end do end do
indx+=1
end do end do
end do end do
!$OMP END DO
!$OMP END PARALLEL
do jndx=1,nMonoEx
do indx=1,jndx-1
hessmat2(indx,jndx) = hessmat2(jndx,indx)
enddo
enddo
END_PROVIDER END_PROVIDER
@ -524,8 +538,8 @@ real*8 function hessmat_taub(t,a,u,b)
tt=list_act(t) tt=list_act(t)
aa=list_virt(a) aa=list_virt(a)
if (t.eq.u) then if (t == u) then
if (a.eq.b) then if (a == b) then
! ta/ta ! ta/ta
t1=occnum(tt)*Fipq(aa,aa) t1=occnum(tt)*Fipq(aa,aa)
t2=0.D0 t2=0.D0
@ -581,8 +595,8 @@ real*8 function hessmat_taub(t,a,u,b)
if (a.eq.b) then if (a.eq.b) then
term-=0.5D0*(occnum(tt)*Fipq(uu,tt)+occnum(uu)*Fipq(tt,uu)) term-=0.5D0*(occnum(tt)*Fipq(uu,tt)+occnum(uu)*Fipq(tt,uu))
do v=1,n_act_orb do v=1,n_act_orb
do x=1,n_act_orb
do y=1,n_act_orb do y=1,n_act_orb
do x=1,n_act_orb
term-=P0tuvx_no(t,v,x,y)*bielecCI_no(x,y,v,uu) term-=P0tuvx_no(t,v,x,y)*bielecCI_no(x,y,v,uu)
term-=P0tuvx_no(u,v,x,y)*bielecCI_no(x,y,v,tt) term-=P0tuvx_no(u,v,x,y)*bielecCI_no(x,y,v,tt)
end do end do
@ -602,29 +616,41 @@ BEGIN_PROVIDER [real*8, hessdiag, (nMonoEx)]
! the diagonal of the Hessian, needed for the Davidson procedure ! the diagonal of the Hessian, needed for the Davidson procedure
END_DOC END_DOC
implicit none implicit none
integer :: i,t,a,indx integer :: i,t,a,indx,indx_shift
real*8 :: hessmat_itju,hessmat_iajb,hessmat_taub real*8 :: hessmat_itju,hessmat_iajb,hessmat_taub
indx=0 !$OMP PARALLEL DEFAULT(NONE) &
!$OMP SHARED(hessdiag,n_core_inact_orb,n_act_orb,n_virt_orb,nMonoEx) &
!$OMP PRIVATE(i,indx,t,a,indx_shift)
!$OMP DO
do i=1,n_core_inact_orb do i=1,n_core_inact_orb
do t=1,n_act_orb do t=1,n_act_orb
indx+=1 indx = t + (i-1)*n_act_orb
hessdiag(indx)=hessmat_itju(i,t,i,t) hessdiag(indx)=hessmat_itju(i,t,i,t)
end do end do
end do end do
!$OMP END DO NOWAIT
do i=1,n_core_inact_orb indx_shift = n_core_inact_orb*n_act_orb
!$OMP DO
do a=1,n_virt_orb do a=1,n_virt_orb
indx+=1 do i=1,n_core_inact_orb
indx = a + (i-1)*n_virt_orb + indx_shift
hessdiag(indx)=hessmat_iajb(i,a,i,a) hessdiag(indx)=hessmat_iajb(i,a,i,a)
end do end do
end do end do
!$OMP END DO NOWAIT
do t=1,n_act_orb indx_shift += n_core_inact_orb*n_virt_orb
!$OMP DO
do a=1,n_virt_orb do a=1,n_virt_orb
indx+=1 do t=1,n_act_orb
indx = a + (t-1)*n_virt_orb + indx_shift
hessdiag(indx)=hessmat_taub(t,a,t,a) hessdiag(indx)=hessmat_taub(t,a,t,a)
end do end do
end do end do
!$OMP END DO
!$OMP END PARALLEL
END_PROVIDER END_PROVIDER