Commit 403f611f authored by Marie bernadette Lepetit's avatar Marie bernadette Lepetit
Browse files

dens dans le cas semble marcher mais doit etre verifiee plus avant

parent 122fb6e1
......@@ -218,6 +218,22 @@ Program RelaxSE
flush(f_output)
endif
if (lprop) then
write(f_output,*)
call wrt_info(prog_info,g_info,o_info,v_info,det_info,d,sym_info,f_info)
write(f_output,*)">>> x_info written on file"
flush(f_info)
flush(f_output)
write(f_output,*)
call deter_init(det,det_info%ndet)
call fill_detd(det, d)
call wrt_bdet(det,det_info%ndet,f_bdet)
write(f_output,*)">>> Determinants written on file"
flush(f_output)
flush(f_bdet)
call deter_free(det)
endif
!Initialisation of the bare h matrix
! contains ntot * ntot elements
! their indices run from ngel+1 to ngel+ntot
......@@ -243,24 +259,6 @@ Program RelaxSE
#ifdef VAR_MPI
call MPI_BARRIER(MPI_COMM_WORLD,ierr)
#endif
! write info for prop
if (lprop) then
write(f_output,*)
call wrt_info(prog_info,g_info,o_info,v_info,det_info,int_info,d,sym_info,f_info)
write(f_output,*)">>> x_info written on file"
flush(f_info)
flush(f_output)
write(f_output,*)
call deter_init(det,det_info%ndet)
call fill_detd(det, d)
call wrt_bdet(det,det_info%ndet,f_bdet)
write(f_output,*)">>> Determinants written on file"
flush(f_output)
flush(f_bdet)
call deter_free(det)
endif
!Build the fock matrix
!Add the 2e- part to the bare h_pq
......
......@@ -41,6 +41,7 @@ Module files
!!$ f_det : tous les determinants en clair
!!$ f_bdet : tous les determinants en binaire
!!$ f_info : variables globales x_info
!!$ f_dens : matrice densite
!!$ Autres
!!$ f_hcore : Hcore (a voir si vraiment nécessaire)
......@@ -52,7 +53,7 @@ Module files
!!$ -----------------------------------------------------------
integer, parameter :: f_input=1, f_output=7, f_ref0=8, f_det=9, f_bdet=39
integer, parameter :: f_fock=10
integer, parameter :: f_info=24
integer, parameter :: f_info=24, f_dens=25
integer, parameter :: f_gen=14, f_restart = 15, f_gen0 = 16
integer, parameter :: f_tone=31, f_tint=32, f_mat=33, f_bmat=34, f_mat2=35
......
......@@ -28,7 +28,7 @@
!!----
subroutine init_prop(g_info, prog_info, o_info, v_info, det_info, int_info, sym_info, prop_info)
subroutine init_prop(g_info, prog_info, o_info, v_info, det_info, sym_info, prop_info)
!!$ Initialisation des variables
!!$ -------- Donness globales ---------------------------------
use dimensions
......@@ -44,7 +44,6 @@ subroutine init_prop(g_info, prog_info, o_info, v_info, det_info, int_info, sym_
type(o_infotype), intent(inout) :: o_info
type(v_infotype), intent(inout) :: v_info
type(det_infotype), intent(inout) :: det_info
type(int_infotype), intent(inout) :: int_info
type(sym_infotype), intent(inout) :: sym_info
type(prop_infotype), intent(inout) :: prop_info
......@@ -97,41 +96,6 @@ subroutine init_prop(g_info, prog_info, o_info, v_info, det_info, int_info, sym_
det_info%nmonoref0 = 0
det_info%ndiref0 = 0
! int_info
int_info%n_1int = 0
int_info%n_oooo = 0
int_info%n_aaaa = 0
int_info%n_aaao = 0
!int_info%n_aoaa = 0
int_info%n_aaoo = 0
int_info%n_aoao = 0
int_info%n_aooo = 0
int_info%n_vooo = 0
int_info%n_vaoo = 0
int_info%n_voao = 0
int_info%n_vaao = 0
int_info%n_voaa = 0
int_info%n_vaaa = 0
int_info%n_vvoo = 0
int_info%n_vovo = 0
int_info%n_vvao = 0
int_info%n_vavo = 0
!int_info%n_vova = 0
int_info%n_vvaa = 0
int_info%n_vava = 0
int_info%n_vvvo = 0
int_info%n_vvva = 0
int_info%n_vvvv = 0
int_info%n_gint = 0
int_info%n_2int = 0
int_info%nintkind = 0
int_info%CASS_nintkind = 21
allocate(int_info%CASS_intkind(int_info%CASS_nintkind))
int_info%CASS_intkind(:) &
= (/ 'fock', 'aaaa', 'aaao', 'vaaa', 'aaoo', 'vaao', 'vvaa','vaoo', &
'vvao','vvoo', 'vava','vvvo','vvva','vvvv','oooo','vovo','aooo','vavo',&
'vooo','aoao','voao'/)
!sym_info
sym_info%iChTb(:,:) = 0
sym_info%iIrTb(:,:) = 0
......@@ -186,6 +150,7 @@ subroutine init_prop(g_info, prog_info, o_info, v_info, det_info, int_info, sym_
read(iunit) int_info%CASS_intkind(1:int_info%CASS_nintkind)
end subroutine lire_intinfo
......
......@@ -56,7 +56,6 @@ Program proprietes
type(o_infotype) :: o_info
type(v_infotype) :: v_info
type(det_infotype) :: det_info
type(int_infotype) :: int_info
type(sym_infotype) :: sym_info
......@@ -102,7 +101,7 @@ Program proprietes
!!$-----
!!$----- Initialisations
!!$-----
call init_prop(g_info, prog_info, o_info, v_info, det_info, int_info, sym_info, prop_info)
call init_prop(g_info, prog_info, o_info, v_info, det_info, sym_info, prop_info)
#ifdef VAR_OMP
!$OMP PARALLEL SHARED(nb_thread)
......@@ -126,12 +125,11 @@ Program proprietes
read(f_info) o_info
read(f_info) v_info
read(f_info) det_info
call lire_intinfo(f_info, int_info)
read(f_info) nblock
allocate(shtblkdet(nblock), nblkdet(nblock), deter_index(nblock))
shtblkdet(:) = 0
nblkdet(:) = 0
deter_index = 99
shtblkdet(:) = 0
nblkdet(:) = 0
deter_index(:) = 99
read(f_info) shtblkdet(1:nblock), nblkdet(1:nblock), deter_index(1:nblock)
read(f_info) sym_info
......@@ -232,23 +230,27 @@ Program proprietes
allocate(psi_S(nvec))
psi_S(:) = 0.d0
call s2(psi,psi_S,ndet,nvec,nelact, det, det_info, d, o_info, f_output)
nq=nvec/10
nr=nvec - nq*pas
do i =1,nq
write(f_output,'(" S ")',advance='no')
do iwr = 1, max(2*nact +1-6,0)
write(f_output,'(a1)',advance='no') A1
end do
write(f_output,9012) psi_S((i-1)*pas+1:i*pas)
end do
write(f_output,'(" S ")',advance='no')
do iwr = 1, max(2*nact +1-6,0)
write(f_output,'(a1)',advance='no') A1
end do
write(f_output,9012) psi_S(nq*pas+1:nvec)
write(f_output,*)
deallocate(psi_S)
9012 format(10(F18.8,2x))
nq=nvec/10
nr=nvec - nq*pas
do i =1,nq
write(f_output,'(" S ")',advance='no')
do iwr = 1, max(2*nact +1-6,0)
write(f_output,'(a1)',advance='no') A1
end do
write(f_output,9012) psi_S((i-1)*pas+1:i*pas)
end do
write(f_output,'(" S ")',advance='no')
do iwr = 1, max(2*nact +1-6,0)
write(f_output,'(a1)',advance='no') A1
end do
write(f_output,9012) psi_S(nq*pas+1:nvec)
write(f_output,*)
deallocate(psi_S)
9012 format(10(F18.8,2x))
case("wf")
write(f_output,*)
write(f_output,*) " >>> Total WFs"
call wrt_WF(ndet,nvec,1,nvec,det,psi,o_info,f_output)
case("pref1")
! Projection sur les ref1
write(f_output,*)
......@@ -262,10 +264,9 @@ Program proprietes
case("dens")
! Calcul de la matrice densite
write(f_output,*)
write(f_output,*) " >>> One particule density matrix"
write(f_output,*) " Not yet implemented"
write(f_output,*) ">>> One particule density matrix"
#ifdef VAR_DEV
call dens(ndet,nvec, psi,det, d,rspin, o_info,g_info, rho_tot,rho_spin)
call dens(ndet,nvec, psi,det, d,rspin, o_info,g_info,prog_info)
#endif
end select
end do
......
......@@ -48,7 +48,8 @@ module utils_wrt
wrtdet, &
wrtact, &
wrtactalpha, &
wrtmat, wrtmat_accurate, wrtmatE, wrt_bdet, wrt_info
wrtmat, wrtmat_accurate, wrtmatE, wrt_bdet, wrt_info, &
wrt_WF
! interface wrt_bdet
! module procedure wrt_bdet
......@@ -153,7 +154,7 @@ contains
!> @param
!> @param
!$============================================================
subroutine wrt_info(prog_info,g_info,o_info,v_info,det_info,int_info,d,sym_info,iunit)
subroutine wrt_info(prog_info,g_info,o_info,v_info,det_info,d,sym_info,iunit)
! -------- Donnes locales -----------------------------------
Implicit none
type(prog_infotype), intent(in) :: prog_info
......@@ -161,7 +162,6 @@ contains
type(o_infotype), intent(in) :: o_info
type(v_infotype), intent(in) :: v_info
type(det_infotype), intent(in) :: det_info
type(int_infotype), intent(in) :: int_info
type(deter_dblocklist),intent(in):: d
type(sym_infotype), intent(in) :: sym_info
integer, intent(in) :: iunit
......@@ -182,55 +182,11 @@ contains
write(iunit) o_info
write(iunit) v_info
write(iunit) det_info
call wrt_intinfo(iunit,int_info)
write(iunit) nblock
write(iunit) shtblkdet(1:nblock), nblkdet(1:nblock), deter_index(1:nblock)
write(iunit) sym_info
end subroutine wrt_info
!$============================================================
!> @brief Write int_info types in binary file
!> @author MBL
!> @date juin 2021
!> @param
!> @param
!$============================================================
subroutine wrt_intinfo(iunit,int_info)
Implicit none
type(int_infotype), intent(in) :: int_info
integer, intent(in) :: iunit
write(iunit) int_info%n_1int, &
int_info%n_oooo , &
int_info%n_aaaa , &
int_info%n_aaao , &
int_info%n_aaoo , &
int_info%n_aoao , &
int_info%n_aooo , &
int_info%n_vooo , &
int_info%n_vaoo , &
int_info%n_voao , &
int_info%n_vaao , &
int_info%n_voaa , &
int_info%n_vaaa , &
int_info%n_vvoo , &
int_info%n_vovo , &
int_info%n_vvao , &
int_info%n_vavo , &
int_info%n_vvaa , &
int_info%n_vava , &
int_info%n_vvvo , &
int_info%n_vvva , &
int_info%n_vvvv , &
int_info%n_gint , &
int_info%n_2int , &
int_info%nintkind , &
int_info%CASS_nintkind
write(iunit) int_info%CASS_intkind(1:int_info%CASS_nintkind)
end subroutine wrt_intinfo
!$============================================================
!> @brief Write a determinant using the global variables
!> @author
......@@ -358,6 +314,7 @@ contains
return
end if
end if
9001 format(17x," -> ",I5,a1)
9002 format(17x," -> ",I5,a1," , ",i5,a1)
......@@ -495,6 +452,191 @@ contains
9993 format (5x,10(E15.6,1x))
end Subroutine wrtmatE
!!$============================================================
!> @brief Write WF determinants in formated file
!> @author MBL
!> @date July 2021
!!$============================================================
subroutine wrt_WF(ndet,nvec,first_vec,last_vec,det,vect,o_info,iunit)
implicit none
Integer, parameter :: pas = 5
Integer (KIND=kd_int), intent(in) :: ndet,nvec,first_vec,last_vec,iunit
type(deter), dimension(ndet), intent(in) :: det
real(kd_dble), dimension(ndet,nvec), intent(in) :: vect
type(o_infotype), intent(in) :: o_info
Integer (KIND=kd_int) :: idet,j,ivec,npas, kpas, shtrou, shpart, shtu, shtd, shpu, shpd, ibit
Character*1, dimension(2) :: spin
Integer(KIND=kd_int), dimension(2) :: sh
! -------- Code ---------------------------------------------
npas = (last_vec - first_vec + 1)/pas
shtrou = 2*o_info%ngel + o_info%nocc + o_info%nligo
shtu = o_info%ngel
shtd = o_info%ngel + o_info%nocc + o_info%nligo
shpart = 2*(o_info%ngel + o_info%nocc + o_info%nligo + o_info%nact) + o_info%nligv + o_info%nvirt
shpu = o_info%ngel + o_info%nocc + o_info%nligo + o_info%nact
shpd = o_info%ngel + o_info%nocc + o_info%nligo + o_info%nact + o_info%nligv + o_info%nvirt
! if (npas.eq.0) goto 10
! do kpas = 1,npas
! do i = 1, ndet
! call wrtdet_local(det(i), iunit, o_info%ngel, o_info%nocc, &
! o_info%nligo, o_info%nact, o_info%nligv, &
! o_info%nvirt, o_info%ndel)
! write(iunit,9990) (vect(i,ivec), ivec = first_vec + (kpas-1)*pas+1, first_vec + kpas*pas)
! end do
! write(iunit,*)
! end do
! 10 if (mod(last_vec - first_vec + 1,pas) .eq. 0) return
write(iunit,*) "***** numeros particules faux"
if (npas.eq.0) goto 10
kp : do kpas = 1,npas
do idet = 1, ndet
! actives
do ibit = 0,o_info%nact-1
write(iunit,'(i1)',advance='no') ibits(det(idet)%detact, ibit, 1)
end do
write(iunit,'(" ")',advance='no')
do ibit = o_info%nact, 2*o_info%nact-1
write(iunit,'(i1)',advance='no') ibits(det(idet)%detact, ibit, 1)
end do
write(iunit,'(" ")',advance='no')
! trous
if (det(idet)%dettr(1).eq.0) then ! 0 trou
write(iunit,'(11x," -> ")',advance='no')
else if (det(idet)%dettr(2).eq.0) then ! 1 trou
if (det(idet)%dettr(1).le.shtrou) then
spin(1)="u"
sh(1) = shtu
else
spin(1)="d"
sh(1) = shtd
end if
write(iunit,'(I4,a1,6x," -> ")',advance='no') det(idet)%dettr(1)-sh(1), spin(1)
else ! 2 trou
if (det(idet)%dettr(1).le.shtrou) then
spin(:)="u"
sh(:) = shtu
else if (det(idet)%dettr(2).le.shtrou) then
spin(1)="d"
sh(1) = shtd
spin(2)="u"
sh(2) = shtu
else
spin(:)="d"
sh(:) = shtd
end if
write(iunit,'(I4,a1,I4,a1," -> ")',advance='no') (det(idet)%dettr(j)-sh(j), spin(j), j=2,1,-1)
end if
! particules
if (det(idet)%detprt(1).eq.0) then ! 0 particules
write(iunit,'(11x,2x)',advance='no')
else if (det(idet)%detprt(2).eq.0) then ! 1 particules
if (det(idet)%detprt(1).le.shpart) then
spin(1)="u"
sh(1) = shpu
else
spin(1)="d"
sh(1) = shpd
end if
write(iunit,'(I4,a1,6x,2x)',advance='no') det(idet)%detprt(1)-sh(1), spin(1)
else ! 2 particules
if (det(idet)%detprt(2).le.shpart) then
spin(:)="u"
sh(:) = shpu
else if (det(idet)%detprt(1).le.shpart) then
spin(1)="u"
sh(1) = shpu
spin(2)="d"
sh(2) = shpd
else
spin(:)="d"
sh(:) = shpd
end if
write(iunit,'(I4,a1,I4,a1,2x)',advance='no') (det(idet)%detprt(j)-sh(j), spin(j), j=1,2)
end if
! vecteur
write(iunit,9990) (vect(idet,ivec), ivec = first_vec + (kpas-1)*pas, first_vec-1 + kpas*pas)
end do
write(iunit,*)
end do kp
10 if (mod(last_vec - first_vec + 1,pas) .eq. 0) return
do idet = 1, ndet
! actives
do ibit = 0,o_info%nact-1
write(iunit,'(i1)',advance='no') ibits(det(idet)%detact, ibit, 1)
end do
write(iunit,'(" ")',advance='no')
do ibit = o_info%nact, 2*o_info%nact-1
write(iunit,'(i1)',advance='no') ibits(det(idet)%detact, ibit, 1)
end do
write(iunit,'(" ")',advance='no')
! trous
if (det(idet)%dettr(1).eq.0) then ! 0 trou
write(iunit,'(11x," -> ")',advance='no')
else if (det(idet)%dettr(2).eq.0) then ! 1 trou
if (det(idet)%dettr(1).le.shtrou) then
spin(1)="u"
sh(1) = shtu
else
spin(1)="d"
sh(1) = shtd
end if
write(iunit,'(I4,a1,6x," -> ")',advance='no') det(idet)%dettr(1)-sh(1), spin(1)
else ! 2 trou
if (det(idet)%dettr(1).le.shtrou) then
spin(:)="u"
sh(:) = shtu
else if (det(idet)%dettr(2).le.shtrou) then
spin(1)="d"
sh(1) = shtd
spin(2)="u"
sh(2) = shtu
else
spin(:)="d"
sh(:) = shtd
end if
write(iunit,'(I4,a1,I4,a1," -> ")',advance='no') (det(idet)%dettr(j)-sh(j), spin(j), j=2,1,-1)
end if
! particules
if (det(idet)%detprt(1).eq.0) then ! 0 particules
write(iunit,'(11x,2x)',advance='no')
else if (det(idet)%detprt(2).eq.0) then ! 1 particules
if (det(idet)%detprt(1).le.shpart) then
spin(1)="u"
sh(1) = shpu
else
spin(1)="d"
sh(1) = shpd
end if
write(iunit,'(I4,a1,6x,2x)',advance='no') det(idet)%detprt(1)-sh(1), spin(1)
else ! 2 particules
if (det(idet)%detprt(2).le.shpart) then
spin(:)="u"
sh(:) = shpu
else if (det(idet)%detprt(1).le.shpart) then
spin(1)="u"
sh(1) = shpu
spin(2)="d"
sh(2) = shpd
else
spin(:)="d"
sh(:) = shpd
end if
write(iunit,'(I4,a1,I4,a1,2x)',advance='no') (det(idet)%detprt(j)-sh(j), spin(j), j=1,2)
end if
! vecteur
write(iunit,9990) (vect(idet,ivec), ivec = first_vec + npas*pas, last_vec )
end do
write(iunit,*)
9990 format(5(2x,F16.12))
end subroutine wrt_WF
!!$============================================================
!!$============================================================
End Module utils_wrt
!!$ Local Variables:
......
......@@ -40,3 +40,7 @@
tol_conv = 1.d-16,
&end
&propinp
nprop = 1,
whichprop="dens"
&end
......@@ -39,3 +39,7 @@
tol_conv = 1.d-16,
&end
&propinp
nprop = 1,
whichprop="dens"
&end
......@@ -40,3 +40,7 @@
tol_conv = 1.d-16,
&end
&propinp
nprop = 2,
whichprop="WF", "dens"
&end
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment