NECICalcInit Subroutine

subroutine NECICalcInit(iCacheFlag)
    != Calculation specific initialisation: just a wrapper for the individual
    != initiialisation routines for each part of the calculation.
    != In:
    !=    iCacheFlag: controls the behaviour of the 4-index integral cache.
    !=                Currently relevant only for CPMD and VASP calculations.
    !=                iCacheFlag=0,1: initialise the cache.
    !=                iCacheFlag=2,3: reuse the cache from the previous NECI
    !=                                calculation from within the same CPMD/VASP
    !=                                calculation.
    use constants, only: stdout
    use System, only: SysInit
    use SystemData, only: tRotateOrbs, tFindCINatOrbs, tGUGA, tUEG, &
                          t_ueg_transcorr, t_ueg_dump, tContact, t_mol_3_body
    use Integrals_neci, only: IntInit, IntFreeze, tPostFreezeHF, DumpFCIDUMP
    use IntegralsData, only: tDumpFCIDUMP
    use DetCalc, only: DetCalcInit, DoDetCalc
    use Determinants, only: DetPreFreezeInit, DetInit, DetPreFreezeInit_old
    use Calc, only: CalcInit
    use HFCalc, only: HFDoCalc
    use RotateOrbsMod, only: RotateOrbs
    use replica_data, only: init_replica_arrays
    use gen_coul_ueg_mod, only: GEN_Umat_TC, prep_ueg_dump, GEN_Umat_TC_Contact
    use LMat_mod, only: readLMat, setup_tchint_ints
    use guga_init, only: init_guga
    implicit none
    integer, intent(in) :: iCacheFlag

    ! Initialise the global data arrays, whose size depends on global values
    ! that used to be constant.
    ! These are essentially constant arrays available after the input has
    ! been read
    call init_replica_arrays()
    if (t_mol_3_body) call setup_tchint_ints()

!   Initlialize the system.  Sets up ...
!   Symmetry is a subset of the system
    call SysInit()

!   Initialize the integrals.  This will read in integrals, as well as calculating
!   some relevant integrals if they are calculated
    call IntInit(iCacheFlag)

!   If we need to freeze some integrals, do so.  Previously we used to do an HF
!   calculation here.  Instead we relegate that to the main Calc section, and are
!   required to read in the relevant orbitals if necessary.

!   This will also call SysPostFreezeInit()
    if (tGUGA) then
        call DetPreFreezeInit_old()
    else
        call DetPreFreezeInit()
    end if

    !! we prepare the contribution of the 2 body transcorrelated operator

    If (tUEG .and. t_ueg_transcorr) then
!                  CALL GetUMatSize(nBasis,nEl,UMATINT)
!                  call shared_allocate("umat_TC3", umat_TC3, (/UMatInt/))
!                  !allocate(UMat(UMatInt), stat=ierr)
!                  LogAlloc(ierr, 'UMat_TC3', int(UMatInt),HElement_t_SizeB, tagUMat)
!                  UMat_TC3 = 0.0_dp
!                  write(stdout,*) "Size of UMat_TC3 is: ",UMATINT

        write(stdout, *) 'prepare the convolution part of the 2 body transcorrelated operator'

        If (tContact) then
            call GEN_Umat_TC_contact
        else
            call GEN_Umat_TC
        end if
        write(stdout, *) "The infinite sums for the transcorrelated approach is determined."

        if (t_ueg_dump) call prep_ueg_dump

    !!$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
    end if

    if (.not. tPostFreezeHF) then
        call HFDoCalc()
    end if

    call IntFreeze()
    if(t_mol_3_body) call readLMat()

    if (tPostFreezeHF) call HFDoCalc()

    if (tDumpFCIDUMP) then
        !If wanted, write out a new FCIDUMP file (with the core frozen if necessary)
        call DumpFCIDUMP()
    end if

    call DetInit()

! Deal with the many-electron basis, setting up sym etc.
! If we are doing a determinant-enumerated calc, this will generate lists of determinants.
!    call DetFreeze()

!   Do any determinant-list based calculations if we have them.  (diagonalizations probably)
    call DetCalcInit()
    call DoDetCalc()

!   Do any initializations we need to do for calculations (e.g. ...?)

    call CalcInit()

    ! can i initialize the GUGA stuff here? after freezing? or otherwise
    ! it is incorrectly setup..
    ! try to init guga here..
    if (tGUGA) call init_guga()

    IF (tRotateOrbs .and. (.not. tFindCINatOrbs)) THEN
        CALL RotateOrbs()
    end if

end subroutine NECICalcInit