printProjEContrib Subroutine

public subroutine printProjEContrib()

Arguments

None

Contents

Source Code


Source Code

    subroutine printProjEContrib()
        implicit none
        character(*), parameter :: refFilename = "kMatProjE"
        integer :: iunit
        integer :: i, j, k, l
        real(dp) :: tmp
        integer(MPIArg) :: ierror
        integer :: kMatSize
        integer :: nBI

        nBI = numBasisIndices(nBasis)
        kMatSize = int(determineKMatSize())

! I/O only done by root
        if (iProcIndex == root) then
            ! do an in-place reduction to save memory
            call MPI_Reduce(MPI_IN_PLACE, kMatProjEContrib, kMatSize, MPI_DOUBLE_PRECISION, &
                            root, MPI_SUM, MPI_COMM_WORLD, ierror)
            ! open the file
            iunit = get_free_unit()
            call open_new_file(iunit, refFilename)

            ! for each possible double excit, print the contribution
            ! TODO: In principle, we only need to loop over i,j occ and k,l virtual
            do i = 1, nBI
                do j = i, nBI
                    do k = 1, nBI
                        do l = k, nBI
                            tmp = kMatProjEContrib(UMatInd(i, j, k, l))
                            if (abs(tmp) > eps) then
                                ! stored is the full accumulated contribution, but we are
                                ! only interested in the average
                                tmp = tmp / sum(all_sum_proje_denominator)
                                write(iunit, *) i, j, k, l, tmp
                            end if
                        end do
                    end do
                end do
            end do

            tmp = sum(kMatProjEContrib) / sum(all_sum_proje_denominator)
            write(stdout, *) "Total transcorrelated 2-Body correlation energy", tmp
        else
            call MPI_Reduce(kMatProjEContrib, kMatProjEContrib, kMatSize, MPI_DOUBLE_PRECISION, &
                            root, MPI_SUM, MPI_COMM_WORLD, ierror)
        end if
    end subroutine printProjEContrib