finalize_double_occ_and_spin_diff Subroutine

    subroutine finalize_double_occ_and_spin_diff()
        ! routine to communicate all the data from all nodes and also
        ! print out the information i guess..
        character(*), parameter :: this_routine = "finalize_double_occ_and_spin_diff"

        real(dp), allocatable :: all_double_occ_vec(:), all_spin_diff(:), &
                                 all_sum_double_occ_vec(:), all_sum_spin_diff(:)

        allocate(all_double_occ_vec(nBasis / 2))
        allocate(all_spin_diff(nBasis / 2))
        allocate(all_sum_double_occ_vec(nBasis / 2))
        allocate(all_sum_spin_diff(nBasis / 2))

        all_double_occ_vec = 0.0_dp
        all_spin_diff = 0.0_dp
        all_sum_double_occ_vec = 0.0_dp
        all_sum_spin_diff = 0.0_dp

        call MPIAllreduce(spin_diff, MPI_SUM, all_spin_diff)
        call MPIAllreduce(double_occ_vec, MPI_SUM, all_double_occ_vec)
        call MPIAllreduce(sum_double_occ_vec, MPI_SUM, all_sum_double_occ_vec)
        call MPIAllreduce(sum_spin_diff, MPI_SUM, all_sum_spin_diff)

        if (iProcIndex == Root) then
            ! first ouput the double occupancy and spin difference for each
            ! spatial orbital
            print *, "double occupancy for each orbital: "
            print *, sum_double_occ_vec / (sum_norm_psi_squared * real(StepsSft, dp))

            print *, "spin difference for each orbital: "
            print *, sum_spin_diff / (sum_norm_psi_squared * real(StepsSft, dp))

            ! and also calculate the summed value over all orbitals to
            ! compare it with the other method and two_rdms

            print *, "total double double occupancy: "
            print *, 2.0_dp * sum(sum_double_occ_vec) / &
                (sum_norm_psi_squared * real(nBasis, dp) * real(StepsSft, dp))

            print *, "total spin difference (should have to do with Ms!)"
            print *, 2.0_dp * sum(sum_spin_diff) / &
                (sum_norm_psi_squared * real(nbasis, dp) * real(StepsSft, dp))

        end if

        call deallocate_spin_measurements()

        deallocate(all_double_occ_vec)
        deallocate(all_spin_diff)
        deallocate(all_sum_double_occ_vec)

    end subroutine finalize_double_occ_and_spin_diff