probability_calculator

Dual exact/sampling probability calculator for TE-QSCI excitation analysis.

This module implements both exact and sampling-based probability calculation methods for analyzing excitation probabilities in time-evolved quantum states.

Classes:

Name	Description
ProbabilityCalculator	Dual exact/sampling probability calculator for quantum states.
H6FigureOneAnalyzer	Specialized analyzer for reproducing Figure 1 with H6 molecule.

ProbabilityCalculator

ProbabilityCalculator(method='auto', verbose=True)

Dual exact/sampling probability calculator for quantum states.

Parameters:

Name	Type	Description	Default
method	str	"exact", "sampling", or "auto" for automatic selection	'auto'
verbose	bool	Whether to print progress information	True

Methods:

Name	Description
select_states_at_t1	Select top R states at t=1 and classify by excitation order.
calculate_probabilities	Calculate probabilities for selected states using chosen method.
calculate_grouped_probabilities	Calculate grouped probabilities for Figure 1 reproduction.

Attributes:

Name	Type	Description
method
verbose
sampler
exact_qubit_threshold
exact_state_threshold

Source code in src/probability_calculator.py

def __init__(self, method: str = "auto", verbose: bool = True):
    """Initialize probability calculator.

    Args:
        method: "exact", "sampling", or "auto" for automatic selection
        verbose: Whether to print progress information
    """
    self.method = method
    self.verbose = verbose
    self.sampler = None

    # Performance thresholds for auto-selection
    self.exact_qubit_threshold = 14  # Use exact for ≤14 qubits
    self.exact_state_threshold = 16384  # Use exact for ≤16K states

method instance-attribute

method = method

verbose instance-attribute

verbose = verbose

sampler instance-attribute

sampler = None

exact_qubit_threshold instance-attribute

exact_qubit_threshold = 14

exact_state_threshold instance-attribute

exact_state_threshold = 16384

select_states_at_t1

select_states_at_t1(
    evolved_state_t1, n_electrons, n_qubits, R=850
)

Select top R states at t=1 and classify by excitation order.

Always uses exact calculation for state selection regardless of method.

Parameters:

Name	Type	Description	Default
evolved_state_t1	GeneralCircuitQuantumState	Quantum state evolved to t=1	required
n_electrons	int	Number of electrons in the system	required
n_qubits	int	Number of qubits in the system	required
R	int	Number of top states to select	850

Returns:

Type	Description
Tuple[ndarray, Dict[int, List[int]]]	Tuple of (selected_state_indices, state_orders_dict)

Source code in src/probability_calculator.py

def select_states_at_t1(
    self, 
    evolved_state_t1: GeneralCircuitQuantumState,
    n_electrons: int,
    n_qubits: int,
    R: int = 850
) -> Tuple[np.ndarray, Dict[int, List[int]]]:
    """Select top R states at t=1 and classify by excitation order.

    Always uses exact calculation for state selection regardless of method.

    Args:
        evolved_state_t1: Quantum state evolved to t=1
        n_electrons: Number of electrons in the system
        n_qubits: Number of qubits in the system
        R: Number of top states to select

    Returns:
        Tuple of (selected_state_indices, state_orders_dict)
    """
    if self.verbose:
        print(f"Selecting top R={R} states at t=1 using exact calculation...")

    # Get exact state vector
    state_vector = evaluate_state_to_vector(evolved_state_t1)
    if hasattr(state_vector, 'vector'):
        probabilities = np.abs(state_vector.vector)**2
    else:
        probabilities = np.abs(state_vector)**2

    # Filter to only valid electron configurations
    valid_states = []
    valid_probabilities = []

    for i in range(len(probabilities)):
        if bin(i).count('1') == n_electrons:
            valid_states.append(i)
            valid_probabilities.append(probabilities[i])

    valid_states = np.array(valid_states)
    valid_probabilities = np.array(valid_probabilities)

    if self.verbose:
        print(f"✓ Found {len(valid_states)} valid electron configurations")

    # Select top R states from valid configurations
    if len(valid_states) <= R:
        selected_indices = valid_states
        if self.verbose:
            print(f"✓ Using all {len(valid_states)} valid states (fewer than R={R})")
    else:
        top_indices = np.argsort(valid_probabilities)[-R:]
        selected_indices = valid_states[top_indices]
        if self.verbose:
            selected_probs = valid_probabilities[top_indices]
            print(f"✓ Selected top {R} states")
            print(f"✓ Probability range: {selected_probs.min():.6f} to {selected_probs.max():.6f}")

    # Classify selected states by excitation order
    hf_state_idx = (1 << n_electrons) - 1  # |11111...000...⟩
    selected_state_orders = {}

    for state_idx in selected_indices:
        order = self._classify_excitation_order(state_idx, hf_state_idx, n_electrons)
        if order >= 0:
            if order not in selected_state_orders:
                selected_state_orders[order] = []
            selected_state_orders[order].append(state_idx)

    if self.verbose:
        print("✓ Selected state distribution by excitation order:")
        for order in sorted(selected_state_orders.keys()):
            count = len(selected_state_orders[order])
            print(f"   Order {order}: {count} states")

    return selected_indices, selected_state_orders

calculate_probabilities

calculate_probabilities(
    evolved_state,
    selected_states,
    n_qubits,
    method_override=None,
)

Calculate probabilities for selected states using chosen method.

Parameters:

Name	Type	Description	Default
evolved_state	GeneralCircuitQuantumState	Time-evolved quantum state	required
selected_states	ndarray	Array of selected state indices	required
n_qubits	int	Number of qubits	required
method_override	Optional[str]	Override the default method selection	None

Returns:

Type	Description
ndarray	Array of probabilities for selected states

Source code in src/probability_calculator.py

def calculate_probabilities(
    self,
    evolved_state: GeneralCircuitQuantumState,
    selected_states: np.ndarray,
    n_qubits: int,
    method_override: Optional[str] = None
) -> np.ndarray:
    """Calculate probabilities for selected states using chosen method.

    Args:
        evolved_state: Time-evolved quantum state
        selected_states: Array of selected state indices
        n_qubits: Number of qubits
        method_override: Override the default method selection

    Returns:
        Array of probabilities for selected states
    """
    # Determine method to use
    if method_override:
        use_method = method_override
    elif self.method == "auto":
        use_method = self._auto_select_method(n_qubits, len(selected_states))
    else:
        use_method = self.method

    if use_method == "exact":
        return self._exact_probabilities(evolved_state, selected_states)
    elif use_method == "sampling":
        return self._sampling_probabilities(evolved_state, selected_states, n_qubits)
    else:
        raise ValueError(f"Unknown method: {use_method}")

calculate_grouped_probabilities

calculate_grouped_probabilities(
    probabilities, selected_state_orders, selected_states
)

Calculate grouped probabilities for Figure 1 reproduction.

Parameters:

Name	Type	Description	Default
probabilities	ndarray	Probabilities for selected states	required
selected_state_orders	Dict[int, List[int]]	Dictionary mapping excitation order to state indices	required
selected_states	ndarray	Array of selected state indices	required

Returns:

Type	Description
Tuple[float, float, float]	Tuple of (prob_1_2, prob_3_4, prob_5_6) for grouped excitation orders

Source code in src/probability_calculator.py

def calculate_grouped_probabilities(
    self,
    probabilities: np.ndarray,
    selected_state_orders: Dict[int, List[int]],
    selected_states: np.ndarray
) -> Tuple[float, float, float]:
    """Calculate grouped probabilities for Figure 1 reproduction.

    Args:
        probabilities: Probabilities for selected states
        selected_state_orders: Dictionary mapping excitation order to state indices
        selected_states: Array of selected state indices

    Returns:
        Tuple of (prob_1_2, prob_3_4, prob_5_6) for grouped excitation orders
    """
    # Create mapping from state index to probability
    state_to_prob = {selected_states[i]: probabilities[i] for i in range(len(selected_states))}

    # Group 1: One/Two-electron excitations (orders 1,2)
    prob_1_2 = 0.0
    for order in [1, 2]:
        if order in selected_state_orders:
            prob_1_2 += sum(state_to_prob.get(state_idx, 0.0) for state_idx in selected_state_orders[order])

    # Group 2: Three/Four-electron excitations (orders 3,4)
    prob_3_4 = 0.0
    for order in [3, 4]:
        if order in selected_state_orders:
            prob_3_4 += sum(state_to_prob.get(state_idx, 0.0) for state_idx in selected_state_orders[order])

    # Group 3: Five/Six-electron excitations (orders 5,6)
    prob_5_6 = 0.0
    for order in [5, 6]:
        if order in selected_state_orders:
            prob_5_6 += sum(state_to_prob.get(state_idx, 0.0) for state_idx in selected_state_orders[order])

    return prob_1_2, prob_3_4, prob_5_6

H6FigureOneAnalyzer

H6FigureOneAnalyzer(
    hamiltonian,
    hartree_fock_state,
    active_space,
    method="auto",
)

Specialized analyzer for reproducing Figure 1 with H6 molecule.

Parameters:

Name	Type	Description	Default
hamiltonian		H6 Hamiltonian operator	required
hartree_fock_state		Hartree-Fock initial state	required
active_space		Active space information	required
method	str	Probability calculation method	'auto'

Methods:

Name	Description
run_figure_one_analysis	Run complete Figure 1 analysis for H6.

Attributes:

Name	Type	Description
hamiltonian
hartree_fock_state
active_space
calculator
n_electrons
n_qubits

Source code in src/probability_calculator.py

def __init__(self, hamiltonian, hartree_fock_state, active_space, method: str = "auto"):
    """Initialize H6 Figure 1 analyzer.

    Args:
        hamiltonian: H6 Hamiltonian operator
        hartree_fock_state: Hartree-Fock initial state
        active_space: Active space information
        method: Probability calculation method
    """
    self.hamiltonian = hamiltonian
    self.hartree_fock_state = hartree_fock_state
    self.active_space = active_space
    self.calculator = ProbabilityCalculator(method=method)

    # H6 system parameters
    self.n_electrons = active_space.n_active_ele if active_space else 6
    self.n_qubits = 2 * active_space.n_active_orb if active_space else 12

hamiltonian instance-attribute

hamiltonian = hamiltonian

hartree_fock_state instance-attribute

hartree_fock_state = hartree_fock_state

active_space instance-attribute

active_space = active_space

calculator instance-attribute

calculator = ProbabilityCalculator(method=method)

n_electrons instance-attribute

n_electrons = n_active_ele if active_space else 6

n_qubits instance-attribute

n_qubits = 2 * n_active_orb if active_space else 12

run_figure_one_analysis

run_figure_one_analysis(R=850, time_points=None)

Run complete Figure 1 analysis for H6.

Parameters:

Name	Type	Description	Default
R	int	Number of states to select (paper uses R=850)	850
time_points	Optional[ndarray]	Custom time points (default: focus on small-t regime)	None

Returns:

Type	Description
Dict	Dictionary containing all analysis results

Source code in src/probability_calculator.py

def run_figure_one_analysis(self, R: int = 850, time_points: Optional[np.ndarray] = None) -> Dict:
    """Run complete Figure 1 analysis for H6.

    Args:
        R: Number of states to select (paper uses R=850)
        time_points: Custom time points (default: focus on small-t regime)

    Returns:
        Dictionary containing all analysis results
    """
    if time_points is None:
        # Focus on small-t regime for scaling analysis, then broader range
        small_times = np.logspace(-2, 0, 20)  # 0.01 to 1.0
        large_times = np.linspace(1.2, 3.0, 10)  # 1.2 to 3.0
        time_points = np.concatenate([small_times, large_times])
        time_points = np.sort(time_points)

    print(f"=== H6 FIGURE 1 ANALYSIS ===")
    print(f"System: {self.n_electrons} electrons, {self.n_qubits} qubits")
    print(f"States to select: R = {R}")
    print(f"Time points: {len(time_points)} (range {time_points.min():.3f} to {time_points.max():.3f})")
    print(f"Method: {self.calculator.method}")

    # Step 1: State selection at t=1
    print(f"\n=== STEP 1: STATE SELECTION AT t=1 ===")
    start_time = time.time()

    evolved_state_t1 = self._evolve_exact(1.0)
    selected_states, selected_state_orders = self.calculator.select_states_at_t1(
        evolved_state_t1, self.n_electrons, self.n_qubits, R
    )

    selection_time = time.time() - start_time
    print(f"✓ State selection completed in {selection_time:.2f} seconds")

    # Step 2: Time evolution analysis
    print(f"\n=== STEP 2: TIME EVOLUTION ANALYSIS ===")
    start_time = time.time()

    results_1_2 = []
    results_3_4 = []
    results_5_6 = []

    for i, t in enumerate(time_points):
        if i % max(1, len(time_points) // 10) == 0:
            print(f"  Progress: {i+1}/{len(time_points)} (t = {t:.4f})")

        # Evolve state to time t
        evolved_state = self._evolve_exact(t)

        # Calculate probabilities for selected states
        probabilities = self.calculator.calculate_probabilities(
            evolved_state, selected_states, self.n_qubits
        )

        # Group by excitation orders
        prob_1_2, prob_3_4, prob_5_6 = self.calculator.calculate_grouped_probabilities(
            probabilities, selected_state_orders, selected_states
        )

        results_1_2.append(prob_1_2)
        results_3_4.append(prob_3_4)
        results_5_6.append(prob_5_6)

    evolution_time = time.time() - start_time
    print(f"✓ Time evolution analysis completed in {evolution_time:.2f} seconds")

    # Package results
    results = {
        'times': time_points,
        'prob_1_2': np.array(results_1_2),
        'prob_3_4': np.array(results_3_4),
        'prob_5_6': np.array(results_5_6),
        'selected_states': selected_states,
        'selected_state_orders': selected_state_orders,
        'system_info': {
            'n_electrons': self.n_electrons,
            'n_qubits': self.n_qubits,
            'R': R,
            'method': self.calculator.method
        },
        'timing': {
            'selection_time': selection_time,
            'evolution_time': evolution_time,
            'total_time': selection_time + evolution_time
        }
    }

    print(f"\n=== ANALYSIS COMPLETE ===")
    print(f"Total time: {results['timing']['total_time']:.2f} seconds")

    return results