write_pops_det Function

public function write_pops_det(iunit, iunit_2, det, j, gdata) result(bWritten)

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: iunit
integer, intent(in) :: iunit_2
integer(kind=n_int), intent(inout) :: det(0:NIfTot)
integer :: j
real(kind=dp), intent(in) :: gdata(:,:)

Return Value logical


Contents

Source Code


Source Code

    function write_pops_det(iunit, iunit_2, det, j, gdata) result(bWritten)

        ! Output a particle to a popsfile in format acceptable for popsfile v4

        integer, intent(in) :: iunit, iunit_2
        integer(n_int), intent(inout) :: det(0:NIfTot)
        real(dp), intent(in) :: gdata(:, :)
        real(dp) :: real_sgn(lenof_sign), detenergy
        HElement_t(dp) :: hf_helemt, hf_helemt_trans
        integer :: flg, j, k, ex_level, nopen, nI(nel), ex(2, nel)
        logical :: bWritten, is_init, is_init_tmp
        integer :: gdata_size
        character(12) :: format_string
        character(*), parameter :: this_routine = "write_pops_det"
        type(ExcitationInformation_t) :: excitInfo

        bWritten = .false.

        call extract_sign(det, real_sgn)
        !real_sgn = real_sgn * 1000
        ! We don't want to bother outputting empty particles, or those
        ! with a weight which is lower than specified as the cutoff
        if (sum(abs(real_sgn)) > binarypops_min_weight) then
            if (mod(j, iPopsPartEvery) == 0) then
                bWritten = .true.
            end if
        end if

        if (bWritten) then

            gdata_size = size(gdata, dim=1)

            ! Write output in the desired format. If INT64_, we are
            ! including the flag information with the signs in storage in
            ! memory --> need to extract these before outputting them.
            flg = extract_flags(det)
            if (tBinPops) then
                ! All write statements MUST be on the same line, or we end
                ! up with multiple records.
                ! TODO: For POPSFILE V5 --> stream output.
                write (iunit) det(0:NIfD), real_sgn, int(flg, n_int), gdata(1:gdata_size, j)
            else
                do k = 0, nifd
                    write (iunit, '(i24)', advance='no') det(k)
                end do
                do k = 1, lenof_sign
                    write (iunit, '(f30.8)', advance='no') real_sgn(k)
                end do
                write (iunit, '(i24)', advance='no') flg
                do k = 1, gdata_size
                    write (iunit, '(f30.8)', advance='no') gdata(k, j)
                end do
                write (iunit, *)
            end if

            if (tPrintInitiators) then
                is_init = .false.
                do k = 1, inum_runs
                    ! Testing with the TestInititator routine to prevent code
                    ! duplication
                    is_init_tmp = test_flag(det, get_initiator_flag_by_run(k))
                    is_init = is_init .or. TestInitiator(det, j, is_init_tmp, k)
                enddo
                if (is_init) then
                    call decode_bit_det(nI, det)
                    nopen = count_open_orbs(det)

                    hf_helemt = 0.0_dp
                    hf_helemt_trans = 0.0_dp

                    if (tGUGA) then
                        ASSERT(.not. t_non_hermitian_2_body)
                        call calc_guga_matrix_element(&
                                det, CSF_Info_t(det), iLutRef(:, 1), CSF_Info_t(iLutRef(:, 1)), &
                                excitInfo, hf_helemt, .true.)
                        ex_level = excitInfo%excitLvl
                        if (ex_level == -1) ex_level = 0
                    else
                        ex_level = FindBitExcitLevel(ilutRef(:, 1), det, nel, .true.)
                        if (ex_level <= 2 .or. (ex_level == 3 .and. t_3_body_excits)) then
                            if (tHPHF) then
                                hf_helemt = hphf_off_diag_helement(ProjEDet(:, 1), &
                                                                   nI, iLutRef(:, 1), det)

                                if (t_non_hermitian_2_body) then
                                    hf_helemt_trans = hphf_off_diag_helement(nI, &
                                                                             ProjEDet(:, 1), det, iLutRef(:, 1))

                                end if
                            else
                                if (t_lattice_model) then
                                    hf_helemt = get_helement_lattice(ProjEDet(:, 1), &
                                                                     nI, ex_level)
                                    if (t_non_hermitian_2_body) then
                                        hf_helemt_trans = get_helement_lattice(nI, &
                                                                               ProjEDet(:, 1), ex_level)
                                    end if
                                else
                                    hf_helemt = get_helement(ProjEDet(:, 1), nI, &
                                                             ex_level, iLutRef(:, 1), det)
                                    if (t_non_hermitian_2_body) then
                                        hf_helemt_trans = get_helement(nI, ProjEDet(:, 1), &
                                                                       ex_level, det, iLutRef(:, 1))
                                    end if
                                end if
                            end if
                        end if
                    end if

                    ! for singles and doubles also include the excitation matrix
                    ! for GUGA this does not make really sense..
                    ex = 0
                    if (ex_level <= 2) then
                        if (tGUGA) then
                            call getExcitation_guga(ProjEDet(:, 1), nI, ex)
                        else
                            call get_bit_excitmat(ilutRef(:, 1), det, ex, ex_level)
                        end if
                    end if

                    if (tHPHF) then
                        detenergy = hphf_diag_helement(nI, det)
                    else
                        detenergy = get_helement(nI, nI, 0)
                    endif

                    write (format_string, '(a,i0,a)') '(', lenof_sign, 'f20.10,a2)'

                    write (iunit_2, format_string, advance='no') real_sgn, ''

                    ! If energy-scaled walkers are used, also print the scaled number of
                    ! walkers
                    if (tEScaleWalkers) then
                        write (iunit_2, '(f20.10,a20)', advance='no') &
                            real_sgn(1) / scaleFunction(get_diagonal_matel(nI, det) - Hii), ''
                    endif

                    call writebitdet(iunit_2, det, .false.)
                    if (t_non_hermitian_2_body) then
                        write (iunit_2, '(i5,i5,3f20.10,4i5)') &
                            ex_level, nopen, detenergy, hf_helemt, &
                            hf_helemt_trans, ex(1, 1), ex(1, 2), ex(2, 1), ex(2, 2)
                    else
                        write (iunit_2, '(i5,i5,2f20.10,4i5)') &
                            ex_level, nopen, detenergy, hf_helemt, ex(1, 1), ex(1, 2), ex(2, 1), ex(2, 2)
                    end if
                endif
            end if
        end if

    end function write_pops_det