0
# Linear Algebra Utilities
1

2
The linalg module provides linear algebra methods and functions essential for quantum computing primitives, including matrix decompositions, quantum-specific operations, and mathematical utilities.
3

4
## Matrix Decompositions
5

6
### KAK Decomposition
7

8
```python { .api }
9
def kak_decomposition(unitary: np.ndarray, *, atol: float = 1e-8) -> 'KakDecomposition':
10
    """KAK decomposition of two-qubit unitary matrices.
11
    
12
    Args:
13
        unitary: 4x4 unitary matrix to decompose
14
        atol: Absolute tolerance for comparisons
15
        
16
    Returns:
17
        KAK decomposition result
18
    """
19

20
class KakDecomposition:
21
    """Result of KAK decomposition: U = (A1 ⊗ A0) exp(i(x XX + y YY + z ZZ)) (A3 ⊗ A2)."""
22
    
23
    def __init__(self, global_phase: complex,
24
                 single_qubit_operations_before: Tuple[np.ndarray, np.ndarray],
25
                 interaction_coefficients: Tuple[float, float, float], 
26
                 single_qubit_operations_after: Tuple[np.ndarray, np.ndarray]) -> None:
27
        """Initialize KAK decomposition."""
28
        
29
    @property
30
    def global_phase(self) -> complex:
31
        """Global phase factor."""
32
        
33
    @property
34
    def single_qubit_operations_before(self) -> Tuple[np.ndarray, np.ndarray]:
35
        """Single-qubit operations before interaction."""
36
        
37
    @property
38
    def interaction_coefficients(self) -> Tuple[float, float, float]:
39
        """Interaction coefficients (x, y, z) for XX, YY, ZZ terms."""
40
        
41
    @property
42
    def single_qubit_operations_after(self) -> Tuple[np.ndarray, np.ndarray]:
43
        """Single-qubit operations after interaction."""
44
        
45
    def __mul__(self, other: 'KakDecomposition') -> 'KakDecomposition':
46
        """Multiply two KAK decompositions."""
47

48
def kak_vector(unitary: np.ndarray, *, atol: float = 1e-8) -> np.ndarray:
49
    """Extract KAK interaction coefficients from two-qubit unitary."""
50

51
def kak_canonicalize_vector(x: float, y: float, z: float, *, atol: float = 1e-8) -> Tuple[float, float, float]:
52
    """Canonicalize KAK vector to preferred coordinate system."""
53
```
54

55
### Single-Qubit Matrix Decomposition
56

57
```python { .api }
58
def axis_angle(matrix: np.ndarray, *, atol: float = 1e-8) -> 'AxisAngleDecomposition':
59
    """Axis-angle decomposition of single-qubit unitary matrix.
60
    
61
    Args:
62
        matrix: 2x2 unitary matrix
63
        atol: Absolute tolerance
64
        
65
    Returns:
66
        Axis-angle decomposition
67
    """
68

69
class AxisAngleDecomposition:
70
    """Axis-angle decomposition of rotation: exp(-i θ/2 n⃗·σ⃗)."""
71
    
72
    def __init__(self, angle: float, axis: np.ndarray, global_phase: complex = 1) -> None:
73
        """Initialize axis-angle decomposition."""
74
        
75
    @property
76
    def angle(self) -> float:
77
        """Rotation angle in radians."""
78
        
79
    @property  
80
    def axis(self) -> np.ndarray:
81
        """Rotation axis as 3D unit vector."""
82
        
83
    @property
84
    def global_phase(self) -> complex:
85
        """Global phase factor."""
86

87
def deconstruct_single_qubit_matrix_into_angles(mat: np.ndarray) -> Tuple[float, float, float]:
88
    """Deconstruct single-qubit matrix into ZYZ Euler angles."""
89
```
90

91
### Matrix Factorization
92

93
```python { .api }
94
def kron_factor_4x4_to_2x2s(matrix: np.ndarray, *, atol: float = 1e-8) -> Optional[Tuple[np.ndarray, np.ndarray]]:
95
    """Factor 4x4 matrix into Kronecker product of 2x2 matrices.
96
    
97
    Args:
98
        matrix: 4x4 matrix to factor
99
        atol: Tolerance for factorization
100
        
101
    Returns:
102
        Tuple of 2x2 matrices (A, B) such that matrix ≈ A ⊗ B, or None if not factorizable
103
    """
104

105
def unitary_eig(matrix: np.ndarray, *, atol: float = 1e-8) -> Tuple[np.ndarray, np.ndarray]:
106
    """Eigendecomposition of unitary matrix.
107
    
108
    Args:
109
        matrix: Unitary matrix
110
        atol: Tolerance for unitarity check
111
        
112
    Returns:
113
        Tuple of (eigenvalues, eigenvectors)
114
    """
115

116
def map_eigenvalues(matrix: np.ndarray, func: Callable[[complex], complex], *, atol: float = 1e-8) -> np.ndarray:
117
    """Apply function to eigenvalues of matrix.
118
    
119
    Args:
120
        matrix: Input matrix
121
        func: Function to apply to each eigenvalue
122
        atol: Tolerance for computations
123
        
124
    Returns:
125
        Matrix with transformed eigenvalues
126
    """
127

128
def so4_to_magic_su2s(so4_matrix: np.ndarray, *, atol: float = 1e-8) -> Tuple[np.ndarray, np.ndarray]:
129
    """Convert SO(4) matrix to pair of SU(2) matrices using magic basis."""
130
```
131

132
## Matrix Operations
133

134
### Basic Operations
135

136
```python { .api }
137
def dot(*values: np.ndarray) -> np.ndarray:
138
    """Enhanced matrix multiplication for multiple matrices.
139
    
140
    Args:
141
        *values: Matrices to multiply in sequence
142
        
143
    Returns:
144
        Product of all matrices
145
    """
146

147
def kron(*matrices: np.ndarray) -> np.ndarray:
148
    """Kronecker product of multiple matrices.
149
    
150
    Args:
151
        *matrices: Matrices to take Kronecker product of
152
        
153
    Returns:
154
        Kronecker product ⊗ᵢ matrices[i]
155
    """
156

157
def kron_with_controls(matrices: Sequence[np.ndarray], controls: Sequence[Sequence[int]]) -> np.ndarray:
158
    """Kronecker product with control structure.
159
    
160
    Args:
161
        matrices: Matrices for each subsystem
162
        controls: Control qubit specifications
163
        
164
    Returns:
165
        Controlled Kronecker product
166
    """
167

168
def block_diag(*blocks: np.ndarray) -> np.ndarray:
169
    """Block diagonal matrix from given blocks.
170
    
171
    Args:
172
        *blocks: Diagonal blocks
173
        
174
    Returns:
175
        Block diagonal matrix
176
    """
177
```
178

179
### Quantum-Specific Operations
180

181
```python { .api }
182
def apply_matrix_to_slices(target: np.ndarray,
183
                          matrix: np.ndarray, 
184
                          slices: Sequence[Union[slice, int, Ellipsis]],
185
                          *, out: Optional[np.ndarray] = None) -> np.ndarray:
186
    """Apply matrix to specified slices of target tensor.
187
    
188
    Args:
189
        target: Target tensor to modify
190
        matrix: Matrix to apply
191
        slices: Slices specifying which parts to operate on
192
        out: Output buffer (optional)
193
        
194
    Returns:
195
        Modified tensor
196
    """
197

198
def targeted_left_multiply(left_matrix: np.ndarray,
199
                          right_target: np.ndarray,
200
                          target_axes: Sequence[int],
201
                          *, out: Optional[np.ndarray] = None) -> np.ndarray:
202
    """Left-multiply target tensor by matrix on specified axes."""
203

204
def targeted_conjugate_about(tensor: np.ndarray,
205
                           matrix: np.ndarray,
206
                           indices: Sequence[int],
207
                           *, out: Optional[np.ndarray] = None) -> np.ndarray:
208
    """Conjugate tensor by matrix: M† @ tensor @ M on specified indices."""
209

210
def partial_trace(tensor: np.ndarray,
211
                 keep_indices: Sequence[int],
212
                 qid_shape: Tuple[int, ...]) -> np.ndarray:
213
    """Compute partial trace of tensor, keeping specified subsystems.
214
    
215
    Args:
216
        tensor: Input tensor (state vector or density matrix)
217
        keep_indices: Indices of subsystems to keep
218
        qid_shape: Shape of each subsystem
219
        
220
    Returns:
221
        Partial trace over discarded subsystems
222
    """
223

224
def partial_trace_of_state_vector_as_mixture(state_vector: np.ndarray,
225
                                           keep_indices: Sequence[int],
226
                                           qid_shape: Tuple[int, ...]) -> Tuple[Tuple[float, np.ndarray], ...]:
227
    """Partial trace of state vector as mixture representation."""
228

229
def sub_state_vector(state_vector: np.ndarray,
230
                    keep_indices: Sequence[int],
231
                    qid_shape: Tuple[int, ...],
232
                    *, default: Optional[np.ndarray] = None,
233
                    atol: float = 1e-8) -> Optional[np.ndarray]:
234
    """Extract sub-state vector for specified subsystems."""
235
```
236

237
### Kronecker Products for Quantum States
238

239
```python { .api }
240
def state_vector_kronecker_product(*factors: np.ndarray) -> np.ndarray:
241
    """Kronecker product of state vectors.
242
    
243
    Args:
244
        *factors: State vectors to tensor together
245
        
246
    Returns:
247
        Tensor product state vector
248
    """
249

250
def density_matrix_kronecker_product(*factors: np.ndarray) -> np.ndarray:
251
    """Kronecker product of density matrices.
252
    
253
    Args:
254
        *factors: Density matrices to tensor together
255
        
256
    Returns:
257
        Tensor product density matrix
258
    """
259
```
260

261
## Matrix Properties and Predicates
262

263
### Unitary and Orthogonal Matrices
264

265
```python { .api }
266
def is_unitary(matrix: np.ndarray, *, atol: float = 1e-8) -> bool:
267
    """Check if matrix is unitary (U† U = I).
268
    
269
    Args:
270
        matrix: Matrix to test
271
        atol: Absolute tolerance
272
        
273
    Returns:
274
        True if matrix is unitary
275
    """
276

277
def is_orthogonal(matrix: np.ndarray, *, atol: float = 1e-8) -> bool:
278
    """Check if matrix is orthogonal (O^T O = I)."""
279

280
def is_special_unitary(matrix: np.ndarray, *, atol: float = 1e-8) -> bool:
281
    """Check if matrix is special unitary (unitary with determinant 1)."""
282

283
def is_special_orthogonal(matrix: np.ndarray, *, atol: float = 1e-8) -> bool:
284
    """Check if matrix is special orthogonal (orthogonal with determinant 1)."""
285
```
286

287
### Hermitian and Normal Matrices
288

289
```python { .api }
290
def is_hermitian(matrix: np.ndarray, *, atol: float = 1e-8) -> bool:
291
    """Check if matrix is Hermitian (M† = M).
292
    
293
    Args:
294
        matrix: Matrix to test
295
        atol: Absolute tolerance
296
        
297
    Returns:
298
        True if matrix is Hermitian
299
    """
300

301
def is_normal(matrix: np.ndarray, *, atol: float = 1e-8) -> bool:
302
    """Check if matrix is normal (M M† = M† M)."""
303

304
def is_diagonal(matrix: np.ndarray, *, atol: float = 1e-8) -> bool:
305
    """Check if matrix is diagonal."""
306
```
307

308
### Quantum Channel Properties
309

310
```python { .api }
311
def is_cptp(choi: np.ndarray, *, atol: float = 1e-8) -> bool:
312
    """Check if Choi matrix represents completely positive trace-preserving map.
313
    
314
    Args:
315
        choi: Choi matrix representation
316
        atol: Tolerance for tests
317
        
318
    Returns:
319
        True if map is CPTP
320
    """
321
```
322

323
### Matrix Comparison
324

325
```python { .api }
326
def allclose_up_to_global_phase(a: np.ndarray, b: np.ndarray, *, 
327
                               rtol: float = 1e-5, atol: float = 1e-8) -> bool:
328
    """Check if matrices are equal up to global phase.
329
    
330
    Args:
331
        a, b: Matrices to compare
332
        rtol: Relative tolerance
333
        atol: Absolute tolerance
334
        
335
    Returns:
336
        True if matrices are equal up to global phase
337
    """
338

339
def matrix_commutes(m1: np.ndarray, m2: np.ndarray, *, atol: float = 1e-8) -> bool:
340
    """Check if two matrices commute (AB = BA).
341
    
342
    Args:
343
        m1, m2: Matrices to test
344
        atol: Absolute tolerance
345
        
346
    Returns:
347
        True if matrices commute
348
    """
349
```
350

351
## Diagonalization and Eigenvalue Problems
352

353
### Real Symmetric Matrices
354

355
```python { .api }
356
def diagonalize_real_symmetric_matrix(matrix: np.ndarray, *,
357
                                     atol: float = 1e-8) -> Tuple[np.ndarray, np.ndarray]:
358
    """Diagonalize real symmetric matrix.
359
    
360
    Args:
361
        matrix: Real symmetric matrix
362
        atol: Tolerance for symmetry check
363
        
364
    Returns:
365
        Tuple of (eigenvalues, eigenvectors)
366
    """
367

368
def diagonalize_real_symmetric_and_sorted_diagonal_matrices(symmetric_matrix: np.ndarray,
369
                                                          diagonal_matrix: np.ndarray, *,
370
                                                          atol: float = 1e-8) -> Tuple[np.ndarray, np.ndarray]:
371
    """Simultaneously diagonalize symmetric and diagonal matrices with sorted eigenvalues."""
372
```
373

374
### Bidiagonalization
375

376
```python { .api }
377
def bidiagonalize_real_matrix_pair_with_symmetric_products(mat1: np.ndarray, mat2: np.ndarray, *,
378
                                                         atol: float = 1e-8) -> Tuple[np.ndarray, np.ndarray, np.ndarray, np.ndarray]:
379
    """Bidiagonalize pair of real matrices with symmetric products."""
380

381
def bidiagonalize_unitary_with_special_orthogonals(mat: np.ndarray, *,
382
                                                  atol: float = 1e-8) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
383
    """Bidiagonalize unitary matrix using special orthogonal matrices."""
384
```
385

386
## Operator Space Operations
387

388
### Pauli Basis Operations
389

390
```python { .api }
391
PAULI_BASIS: Tuple[np.ndarray, np.ndarray, np.ndarray, np.ndarray]
392
"""Standard Pauli basis matrices [I, X, Y, Z]."""
393

394
def expand_matrix_in_orthogonal_basis(matrix: np.ndarray,
395
                                    basis: Sequence[np.ndarray]) -> np.ndarray:
396
    """Expand matrix in orthogonal basis.
397
    
398
    Args:
399
        matrix: Matrix to expand
400
        basis: Orthogonal basis matrices
401
        
402
    Returns:
403
        Coefficients in the basis expansion
404
    """
405

406
def matrix_from_basis_coefficients(coefficients: Sequence[complex],
407
                                 basis: Sequence[np.ndarray]) -> np.ndarray:
408
    """Construct matrix from basis coefficients.
409
    
410
    Args:
411
        coefficients: Expansion coefficients
412
        basis: Basis matrices
413
        
414
    Returns:
415
        Matrix = Σᵢ coefficients[i] * basis[i]
416
    """
417

418
def kron_bases(*bases: Sequence[np.ndarray]) -> List[np.ndarray]:
419
    """Kronecker product of multiple bases.
420
    
421
    Args:
422
        *bases: Bases to take Kronecker product of
423
        
424
    Returns:
425
        List of all Kronecker products
426
    """
427

428
def hilbert_schmidt_inner_product(mat1: np.ndarray, mat2: np.ndarray) -> complex:
429
    """Hilbert-Schmidt inner product: ⟨A, B⟩ = tr(A† B).
430
    
431
    Args:
432
        mat1, mat2: Matrices for inner product
433
        
434
    Returns:
435
        Inner product value
436
    """
437

438
def pow_pauli_combination(pauli_map: Dict[str, float], exponent: float) -> Dict[str, complex]:
439
    """Raise Pauli combination to a power."""
440
```
441

442
## Tolerance and Comparison Functions
443

444
### Near-Zero Checks
445

446
```python { .api }
447
def all_near_zero(a: np.ndarray, *, atol: float = 1e-8) -> bool:
448
    """Check if all elements are near zero.
449
    
450
    Args:
451
        a: Array to check
452
        atol: Absolute tolerance
453
        
454
    Returns:
455
        True if all elements have magnitude < atol
456
    """
457

458
def all_near_zero_mod(a: np.ndarray, period: float, *, atol: float = 1e-8) -> bool:
459
    """Check if all elements are near zero modulo period."""
460
```
461

462
### Utility Functions
463

464
```python { .api }
465
def match_global_phase(a: np.ndarray, b: np.ndarray) -> np.ndarray:
466
    """Multiply b by phase to match global phase of a.
467
    
468
    Args:
469
        a: Reference matrix
470
        b: Matrix to adjust
471
        
472
    Returns:
473
        Phase-adjusted version of b
474
    """
475

476
def slice_for_qubits_equal_to(target_qubit_axes: Sequence[int],
477
                             little_endian_qureg_value: int = 0,
478
                             qid_shape: Optional[Sequence[int]] = None) -> Tuple[Union[slice, int], ...]:
479
    """Create slice for computational basis state."""
480

481
def phase_delta(mat: np.ndarray, mat2: np.ndarray) -> float:
482
    """Phase difference between matrices."""
483

484
def reflection_matrix_pow(reflection_matrix: np.ndarray, exponent: float) -> np.ndarray:
485
    """Raise reflection matrix to a power."""
486

487
def to_special(matrix: np.ndarray) -> np.ndarray:
488
    """Convert to special unitary/orthogonal (determinant 1)."""
489

490
def extract_right_diag(left: np.ndarray, middle: np.ndarray, right: np.ndarray) -> np.ndarray:
491
    """Extract right diagonal matrix from decomposition."""
492
```
493

494
## Utility and Validation Functions
495

496
```python { .api }
497
def num_cnots_required(mat: np.ndarray) -> int:
498
    """Minimum number of CNOTs required to implement two-qubit unitary.
499
    
500
    Args:
501
        mat: 4x4 unitary matrix
502
        
503
    Returns:
504
        Minimum CNOT count (0, 1, 2, or 3)
505
    """
506

507
def can_numpy_support_shape(shape: Sequence[int]) -> bool:
508
    """Check if numpy can support given tensor shape."""
509

510
def transpose_flattened_array(array: np.ndarray,
511
                            shape: Tuple[int, ...],
512
                            permutation: Sequence[int]) -> np.ndarray:
513
    """Transpose flattened array according to permutation."""
514

515
def scatter_plot_normalized_kak_interaction_coefficients(mat: np.ndarray) -> None:
516
    """Create scatter plot of normalized KAK interaction coefficients."""
517

518
CONTROL_TAG: str
519
"""Tag for control operations in linear algebra contexts."""
520
```
521

522
## Usage Examples
523

524
### Matrix Decompositions
525

526
```python
527
import cirq
528
import numpy as np
529

530
# KAK decomposition of two-qubit gate
531
print("=== KAK Decomposition ===")
532
# Create a two-qubit unitary (CNOT gate)
533
cnot_matrix = np.array([[1, 0, 0, 0],
534
                       [0, 1, 0, 0], 
535
                       [0, 0, 0, 1],
536
                       [0, 0, 1, 0]], dtype=complex)
537

538
kak = cirq.kak_decomposition(cnot_matrix)
539
print(f"Global phase: {kak.global_phase}")
540
print(f"Interaction coefficients: {kak.interaction_coefficients}")
541
print(f"Single-qubit ops before shape: {[m.shape for m in kak.single_qubit_operations_before]}")
542

543
# Verify decomposition
544
reconstructed = kak.global_phase * np.kron(kak.single_qubit_operations_before[1], 
545
                                          kak.single_qubit_operations_before[0])
546
# Apply interaction term...
547
print(f"Decomposition valid: {np.allclose(cnot_matrix, reconstructed, atol=1e-10)}")
548

549
# Axis-angle decomposition for single qubit
550
print("\n=== Axis-Angle Decomposition ===")
551
h_matrix = np.array([[1, 1], [1, -1]]) / np.sqrt(2)  # Hadamard
552
axis_angle = cirq.axis_angle(h_matrix)
553
print(f"Rotation angle: {axis_angle.angle:.4f} radians")
554
print(f"Rotation axis: {axis_angle.axis}")
555
```
556

557
### Matrix Properties and Validation
558

559
```python
560
import cirq
561
import numpy as np
562

563
print("=== Matrix Properties ===")
564

565
# Test various matrix properties
566
matrices = {
567
    'Identity': np.eye(2),
568
    'Hadamard': np.array([[1, 1], [1, -1]]) / np.sqrt(2),
569
    'Pauli-X': np.array([[0, 1], [1, 0]]),
570
    'Random': np.random.randn(2, 2) + 1j * np.random.randn(2, 2)
571
}
572

573
for name, matrix in matrices.items():
574
    print(f"\n{name} matrix:")
575
    print(f"  Unitary: {cirq.is_unitary(matrix)}")
576
    print(f"  Hermitian: {cirq.is_hermitian(matrix)}")
577
    print(f"  Normal: {cirq.is_normal(matrix)}")
578
    
579
    if cirq.is_unitary(matrix):
580
        print(f"  Special unitary: {cirq.is_special_unitary(matrix)}")
581

582
# Test quantum channel properties
583
print("\n=== Quantum Channel Properties ===")
584
# Depolarizing channel Choi matrix
585
p = 0.1
586
depol_choi = (1-p) * np.kron([[1, 0], [0, 0]], np.eye(2)) + p/3 * sum(
587
    np.kron(pauli, pauli) for pauli in [[[0, 1], [1, 0]], [[0, -1j], [1j, 0]], [[1, 0], [0, -1]]]
588
)
589
print(f"Depolarizing channel is CPTP: {cirq.is_cptp(depol_choi)}")
590
```
591

592
### Quantum State Operations
593

594
```python
595
import cirq
596
import numpy as np
597

598
print("=== Quantum State Operations ===")
599

600
# Create entangled state
601
bell_state = np.array([1, 0, 0, 1]) / np.sqrt(2)  # (|00⟩ + |11⟩)/√2
602
density_matrix = np.outer(bell_state, bell_state)
603

604
# Partial trace to get reduced density matrix of first qubit  
605
reduced_dm = cirq.partial_trace(density_matrix, keep_indices=[0], qid_shape=(2, 2))
606
print(f"Reduced density matrix of first qubit:")
607
print(reduced_dm)
608

609
# Check if reduced state is pure
610
eigenvals = np.linalg.eigvals(reduced_dm)
611
purity = np.sum(eigenvals**2)
612
print(f"Purity of reduced state: {purity.real:.4f}")
613
print(f"Is pure state: {np.isclose(purity, 1.0)}")
614

615
# Kronecker products
616
state1 = np.array([1, 0])  # |0⟩
617
state2 = np.array([0, 1])  # |1⟩
618
product_state = cirq.state_vector_kronecker_product(state1, state2)
619
print(f"Product state |01⟩: {product_state}")
620
```
621

622
### Linear Algebra Utilities
623

624
```python
625
import cirq
626
import numpy as np
627

628
print("=== Linear Algebra Utilities ===")
629

630
# Matrix operations
631
A = np.array([[1, 2], [3, 4]])
632
B = np.array([[5, 6], [7, 8]])
633
C = np.array([[1, 0], [0, -1]])
634

635
# Enhanced dot product
636
product = cirq.dot(A, B, C)
637
print(f"A @ B @ C shape: {product.shape}")
638

639
# Kronecker product
640
kron_product = cirq.kron(A, B)
641
print(f"A ⊗ B shape: {kron_product.shape}")
642

643
# Try to factor back
644
factors = cirq.kron_factor_4x4_to_2x2s(kron_product)
645
if factors is not None:
646
    A_recovered, B_recovered = factors
647
    print(f"Successfully factored: {np.allclose(A, A_recovered) and np.allclose(B, B_recovered)}")
648

649
# Pauli basis expansion
650
pauli_x = np.array([[0, 1], [1, 0]])
651
coeffs = cirq.expand_matrix_in_orthogonal_basis(pauli_x, cirq.PAULI_BASIS)
652
print(f"Pauli-X coefficients in Pauli basis: {coeffs}")
653

654
# Reconstruct from coefficients
655
reconstructed = cirq.matrix_from_basis_coefficients(coeffs, cirq.PAULI_BASIS)
656
print(f"Reconstruction successful: {np.allclose(pauli_x, reconstructed)}")
657

658
# Near-zero checks
659
small_matrix = np.array([[1e-10, 1e-9], [1e-11, 1e-12]])
660
print(f"All elements near zero (atol=1e-8): {cirq.all_near_zero(small_matrix, atol=1e-8)}")
661
```
662

663
This comprehensive linear algebra toolkit provides the mathematical foundation for quantum computing operations, enabling efficient manipulation of quantum states, gates, and channels.