make_1e_rdm_hermitian Subroutine

    subroutine make_1e_rdm_hermitian(matrix, norm_1rdm)

        ! Simply average the 1-RDM(i,j) and 1-RDM(j,i) elements which should
        ! be equal in a perfect world.

        use RotateOrbsData, only: SymLabelListInv_rot, NoOrbs

        real(dp), intent(inout) :: matrix(:, :)
        real(dp), intent(in) :: norm_1rdm

        real(dp) :: max_error_herm, sum_error_herm
        integer :: i, j
        real(dp) :: temp

        max_error_herm = 0.0_dp
        sum_error_herm = 0.0_dp

        associate(ind => SymLabelListInv_rot)
            do i = 1, NoOrbs
                do j = i, NoOrbs
                    if ((abs((matrix(ind(i), ind(j)) * norm_1rdm) - (matrix(ind(j), ind(i)) * norm_1rdm))) > max_error_herm) then
                        max_error_herm = abs(matrix(ind(i), ind(j)) * norm_1rdm - matrix(ind(j), ind(i)) * norm_1rdm)
                    end if

                    sum_error_herm = sum_error_herm + abs(matrix(ind(i), ind(j)) * norm_1rdm - matrix(ind(j), ind(i)) * norm_1rdm)

                    temp = (matrix(ind(i), ind(j)) + matrix(ind(j), ind(i))) / 2.0_dp
                    matrix(ind(i), ind(j)) = temp
                    matrix(ind(j), ind(i)) = temp
                end do
            end do
        end associate

        ! Output the hermiticity errors.
        write(stdout, '(1X,"MAX ABS ERROR IN 1-RDM HERMITICITY",F20.13)') max_error_herm
        write(stdout, '(1X,"MAX SUM ERROR IN 1-RDM HERMITICITY",F20.13)') sum_error_herm

    end subroutine make_1e_rdm_hermitian