SUBROUTINE INITFROMFCID(NEL, NBASISMAX, LEN, LMS, TBIN)
logical, intent(in) :: tbin
integer, intent(out) :: nBasisMax(5, *), LEN, LMS
integer, intent(inout) :: NEL
integer SYMLZ(1000), ST, III
integer OCC(nIrreps), CLOSED(nIrreps), FROZEN(nIrreps)
integer(int64) :: ORBSYM(1000)
INTEGER NORB, NELEC, MS2, ISYM, i, SYML(1000), iunit, iuhf
LOGICAL exists
logical :: uhf, trel, tDetectSym
character(*), parameter :: this_routine = 'INITFROMFCID'
real(dp), allocatable :: FOCK(:)
CHARACTER(len=3) :: fmat
NAMELIST /FCI/ NORB, NELEC, MS2, ORBSYM, OCC, CLOSED, FROZEN, &
ISYM, IUHF, UHF, TREL, SYML, SYMLZ, PROPBITLEN, NPROP, ST, III, &
FOCK
UHF = .FALSE.
fmat = 'NO'
PROPBITLEN = 0
NPROP = 0
IUHF = 0
TREL = .false.
SYMLZ(:) = 0
OCC = -1
CLOSED = -1
FROZEN = -1
allocate(FOCK(1000), source=0.0_dp)
! [W.D. 15.5.2017:]
! with the new relativistic calculations, withoug a ms value in the
! FCIDUMP, we have to set some more defaults..
MS2 = 0
IF (iProcIndex == 0) THEN
iunit = get_free_unit()
IF (TBIN) THEN
INQUIRE (FILE='FCISYM', EXIST=exists)
IF (.not. exists) THEN
CALL Stop_All(this_routine, 'FCISYM file does not exist')
end if
INQUIRE (FILE=FCIDUMP_name, EXIST=exists, FORMATTED=fmat)
IF (.not. exists) THEN
CALL Stop_All(this_routine, 'FCIDUMP file does not exist')
end if
open(iunit, FILE='FCISYM', STATUS='OLD', FORM='FORMATTED')
read(iunit, FCI)
ELSE
INQUIRE (FILE=FCIDUMP_name, EXIST=exists, UNFORMATTED=fmat)
IF (.not. exists) THEN
CALL Stop_All(this_routine, 'FCIDUMP file does not exist')
end if
open(iunit, FILE=FCIDUMP_name, STATUS='OLD', FORM='FORMATTED')
read(iunit, FCI)
end if
close(iunit)
end if
! Currently, FOCK array has been introduced to prevent crashing calculations
! It is deallocated here, as it is not used in the code anymore.
! If you want to use FOCK, don't forget to add MPIBCast(FOCK).
deallocate(FOCK)
call reorder_sym_labels(ORBSYM, SYML, SYMLZ)
!Now broadcast these values to the other processors
CALL MPIBCast(NORB, 1)
CALL MPIBCast(NELEC, 1)
CALL MPIBCast(MS2, 1)
CALL MPIBCast(ORBSYM, 1000)
CALL MPIBCast(SYML, 1000)
CALL MPIBCast(SYMLZ, 1000)
CALL MPIBCast(ISYM, 1)
CALL MPIBCast(UHF)
CALL MPIBCast(tRel)
CALL MPIBCast(PROPBITLEN, 1)
CALL MPIBCast(NPROP, 3)
call MPIBCast(OCC, 8)
call MPIBCast(CLOSED, nIrreps)
call MPIBCast(FROZEN, nIrreps)
if (UHF .and. .not. (tUHF .or. tROHF)) then
! unfortunately, the `UHF` keyword in the FCIDUMP namelist indicates
! spin-orbital-resolved integrals, not necessarily UHF
call stop_all(this_routine, 'UHF in FCIDUMP but neither uhf nor rohf in input.')
end if
! If PropBitLen has been set then assume we're not using an Abelian
! symmetry group which has two cycle generators (ie the group has
! complex representations).
TwoCycleSymGens = PropBitLen == 0
tReltvy = tRel
nOccOrbs = OCC
nClosedOrbs = CLOSED
IF (.not. TwoCycleSymGens .and. ((NPROP(1) + NPROP(2) + NPROP(3)) > 3)) THEN
!We are using abelian k-point symmetry. Turn it on.
tKPntSym = .true.
write(stdout, "(A,I5,A)") "Using abelian k-point symmetry - ", NPROP(1) * NPROP(2) * NPROP(3), " kpoints found."
ELSE
tKPntSym = .false.
end if
#ifndef CMPLX_
if (PropBitLen /= 0) then
!We have compiled the code in real, but we are looking at a complex FCIDUMP, potentially
!even with multiple k-points. However, these orbitals must be real, and ensure that the
!integrals are read in as real.
write(stdout, "(A)") "Neci compiled in real mode, but kpoints detected."
write(stdout, "(A)") "This will be ok if kpoints are all at Gamma or BZ boundary and are correctly rotated."
tRotatedOrbsReal = .true.
else if (tRel) then
write(stdout, "(A)") "Relativistic integrals, but neci compiled in real mode."
end if
#endif
if (.not. tMolpro) then
IF (tROHF .and. (.not. UHF)) THEN
CALL Stop_All(this_routine, "ROHF specified, but FCIDUMP is not in a high-spin format.")
end if
end if
tDetectSym = .true.
DO i = 1, NORB
IF (ORBSYM(i) == 0 .and. TwoCycleSymGens) THEN
write(stdout, *) "** WARNING **"
write(stdout, *) "** Unconverged symmetry of orbitals **"
write(stdout, *) "** Turning point group symmetry off for rest of run **"
if (.not. (tFixLz .or. tKPntSym .or. tUEG .or. tHub)) then
write(stdout, *) "** No symmetry at all will be used in excitation generators **"
tNoSymGenRandExcits = .true.
end if
lNoSymmetry = .true.
EXIT
end if
end do
if (tDetectSym) then
!Check that there is more than one irrep
tDetectSym = .false.
do i = 1, Norb
if (OrbSym(i) /= 1) then
tDetectSym = .true.
exit
end if
end do
end if
if (.not. tDetectSym) then
write(stdout, *) "Only one irrep found. Turning off symmetry for rest of calculation."
if (.not. tFixLz .or. tKPntSym) then
tNoSymGenRandExcits = .true.
lNoSymmetry = .true.
end if
end if
IF (NELEC /= NEL) THEN
if( NEL == NEL_UNINITIALIZED ) then
write(stdout,*) "No number of electrons given, using NEL in FCIDUMP"
NEL = NELEC
else
write(stdout, *) &
& '*** WARNING: NEL in FCIDUMP differs from input file ***'
write(stdout, *) ' NUMBER OF ELECTRONS : ', NEL
endif
end if
! NEL=NELEC
IF (LMS /= MS2) THEN
write(stdout, *) &
& '*** WARNING: LMS in FCIDUMP differs from input file ***'
LMS = MS2
write(stdout, *) ' BASIS MS : ', LMS
end if
if (tMolpro) then
!Molpros FCIDUMPs always indicate the number of spatial orbitals.
!Need to x2 to get spin orbitals
LEN = 2 * NORB
else
IF (UHF .or. tRel) then
LEN = NORB
ELSE
LEN = 2 * NORB
end if
end if
NBASISMAX(1:5, 1:3) = 0
! NBASISMAX(1,1)=1
! NBASISMAX(1,2)=NORB
NBASISMAX(1, 1) = 0
NBASISMAX(1, 2) = 0
NBASISMAX(1, 3) = 2
NBASISMAX(4, 1) = -1
NBASISMAX(4, 2) = 1
!the indicator of UHF (which becomes ISpinSkip or ISS later
if (tMolpro) then
if (tUHF) then
NBASISMAX(2, 3) = 1
else
NBASISMAX(2, 3) = 2
end if
else
tStoreSpinOrbs = tStoreSpinOrbs .or. tRel
IF ((UHF .and. (.not. tROHF)) .or. tRel) then
NBASISMAX(2, 3) = 1
else
NBASISMAX(2, 3) = 2
end if
end if
! Show that there's no momentum conservation
NBASISMAX(3, 3) = 1
RETURN
END SUBROUTINE INITFROMFCID