subroutine Standalone_Errors() use sym_mod, only: getsym USE MolproPlugin, only: MolproPluginResult real(dp) :: mean_ProjE_re, mean_ProjE_im, mean_Shift real(dp) :: ProjE_Err_re, ProjE_Err_im, Shift_Err logical :: tNoProjEValue, tNoShiftValue integer :: iroot, isymh TYPE(BasisFn) RefSym HElement_t(dp) :: h_tmp real(dp) :: Hii, BestEnergy, EnergyDiff !Automatic error analysis call error_analysis(tSinglePartPhase(1), iBlockingIter(1), mean_ProjE_re, ProjE_Err_re, & mean_ProjE_im, ProjE_Err_im, mean_Shift, Shift_Err, tNoProjEValue, tNoShiftValue, & equilshift=iBlockEquilShift, equilproje=iBlockEquilProjE) call MPIBCast(ProjectionE) call MPIBCast(mean_ProjE_re) call MPIBCast(ProjE_Err_re) call MPIBCast(mean_ProjE_im) call MPIBCast(ProjE_Err_im) call MPIBCast(mean_Shift) call MPIBCast(Shift_Err) call MPIBCast(tNoProjEValue) call MPIBCast(tNoShiftValue) h_tmp = get_helement(FDet, FDet, 0) Hii = real(h_tmp, dp) iroot = 1 CALL GetSym(FDet, NEl, G1, NBasisMax, RefSym) isymh = int(RefSym%Sym%S) + 1 write(stdout, 10101) iroot, isymh 10101 format(//'RESULTS FOR STATE', i2, '.', i1/'====================='/) write(stdout, '('' Current reference energy'',T52,F19.12)') Hii if (tNoProjEValue) then write(stdout, '('' No projected energy value could be obtained'',T52)') else write(stdout, '('' Projected correlation energy'',T52,F19.12)') mean_ProjE_re write(stdout, '('' Estimated error in Projected correlation energy'',T52,F19.12)') ProjE_Err_re end if if (lenof_sign == 2 .and. inum_runs == 1) then write(stdout, '('' Projected imaginary energy'',T52,F19.12)') mean_ProjE_im write(stdout, '('' Estimated error in Projected imaginary energy'',T52,F19.12)') ProjE_Err_im end if if (tNoShiftValue) then write(stdout, '('' No shift energy value could be obtained'',T52)') else write(stdout, '('' Shift correlation energy'',T52,F19.12)') mean_Shift write(stdout, '('' Estimated error in shift correlation energy'',T52,F19.12)') shift_err end if !Do shift and projected energy agree? write(stdout, "(A)") if (tNoProjEValue .and. tNoShiftValue) return EnergyDiff = abs(mean_Shift - mean_ProjE_re) if (EnergyDiff <= sqrt(shift_err**2 + ProjE_Err_re**2)) then write(stdout, "(A,F15.8)") " Projected and shift energy estimates agree " & //"within errorbars: EDiff = ", EnergyDiff else if (EnergyDiff <= sqrt((max(shift_err, ProjE_Err_re) * 2)**2 + min(shift_err, ProjE_Err_re)**2)) then write(stdout, "(A,F15.8)") " Projected and shift energy estimates agree to within two sigma " & //"of largest error: EDiff = ", EnergyDiff else write(stdout, "(A,F15.8)") " Projected and shift energy estimates do not agree to " & //"within approximate errorbars: EDiff = ", EnergyDiff end if if (ProjE_Err_re < shift_err) then BestEnergy = mean_ProjE_re + Hii else BestEnergy = mean_shift + Hii end if write(stdout, "(A)") write(stdout, "(A,F20.8,A,G15.6)") " Total projected energy ", & mean_ProjE_re + Hii, " +/- ", ProjE_Err_re write(stdout, "(A,F20.8,A,G15.6)") " Total shift energy ", & mean_shift + Hii, " +/- ", shift_err CALL MolproPluginResult('ENERGY', [BestEnergy]) CALL MolproPluginResult('FCIQMC_ERR', [min(ProjE_Err_re, shift_err)]) write(stdout, "(/)") end subroutine Standalone_Errors