Anymatrix Examples

This page provides comprehensive examples of using the anymatrix module.

Basic Matrix Generation

Example 1: Simple Generation

from matrix_toolkit.anymatrix import AnyMatrix

am = AnyMatrix()

# Generate different matrices
beta = am.generate('core/beta', 10)
fourier = am.generate('core/fourier', 8)
hilbert = am.generate('matlab/hilbert', 5)

print(f"Beta shape: {beta.shape}")
print(f"Fourier shape: {fourier.shape}")
print(f"Hilbert shape: {hilbert.shape}")

Output:

Beta shape: (10, 10)
Fourier shape: (8, 8)
Hilbert shape: (5, 5)

Example 2: Matrix with Parameters

# KMS matrix with different parameters
kms_05 = am.generate('gallery/kms', 10, rho=0.5)
kms_09 = am.generate('gallery/kms', 10, rho=0.9)

# Compare condition numbers
import numpy as np
cond_05 = np.linalg.cond(kms_05)
cond_09 = np.linalg.cond(kms_09)

print(f"Condition number (rho=0.5): {cond_05:.2e}")
print(f"Condition number (rho=0.9): {cond_09:.2e}")

Output:

Condition number (rho=0.5): 1.38e+01
Condition number (rho=0.9): 1.26e+03

Property Verification

Example 3: Checking Properties

from matrix_toolkit.anymatrix import MatrixProperties

# Generate a matrix
matrix = am.generate('core/beta', 20)

# Check individual properties
checks = {
    'Symmetric': MatrixProperties.is_symmetric(matrix),
    'Positive Definite': MatrixProperties.is_positive_definite(matrix),
    'Orthogonal': MatrixProperties.is_orthogonal(matrix),
    'Tridiagonal': MatrixProperties.is_tridiagonal(matrix)
}

print("Property Checks:")
for prop, result in checks.items():
    status = "✓" if result else "✗"
    print(f"  {status} {prop}")

Output:

Property Checks:
  ✓ Symmetric
  ✓ Positive Definite
  ✗ Orthogonal
  ✗ Tridiagonal

Example 4: Automated Verification

# Get expected properties
expected = am.properties('core/fourier')
print(f"Expected properties: {expected}")

# Generate and verify
fourier = am.generate('core/fourier', 16)
results = MatrixProperties.check_properties(fourier, expected)

print("\nVerification Results:")
for prop, passed in results.items():
    print(f"  {prop}: {'PASS' if passed else 'FAIL'}")

Output:

Expected properties: ['unitary', 'symmetric']

Verification Results:
  unitary: PASS
  symmetric: PASS

Searching for Matrices

Example 5: Property-Based Search

# Find all symmetric matrices
symmetric = am.search(['symmetric'])
print(f"Found {len(symmetric)} symmetric matrices")
print("First 10:", symmetric[:10])

# Find symmetric AND positive definite
sym_pd = am.search(['symmetric', 'positive definite'])
print(f"\nSymmetric + PD: {len(sym_pd)} matrices")

# Find matrices with ANY of the properties
either = am.search(['symmetric', 'orthogonal'], match_all=False)
print(f"\nSymmetric OR Orthogonal: {len(either)} matrices")

Output:

Found 23 symmetric matrices
First 10: ['core/beta', 'core/wilson', 'core/indef', ...]

Symmetric + PD: 8 matrices

Symmetric OR Orthogonal: 31 matrices

Example 6: Finding Similar Matrices

def find_similar_matrices(matrix_id, am):
    """Find matrices with same properties"""
    # Get properties of reference matrix
    props = am.properties(matrix_id)

    # Search for matrices with same properties
    similar = am.search(props)

    # Remove the reference matrix itself
    similar = [m for m in similar if m != matrix_id]

    return similar

# Find matrices similar to beta
similar = find_similar_matrices('core/beta', am)
print(f"Matrices similar to 'core/beta':")
for matrix_id in similar:
    print(f"  - {matrix_id}")

Output:

Matrices similar to 'core/beta':
  - core/wilson
  - gallery/lehmer
  - gallery/minij
  - matlab/pascal
  - matlab/hilbert

Working with Specific Groups

Example 7: Nilpotent Matrices

import numpy as np

nilpotent_matrices = [
    'core/nilpot_triang',
    'core/nilpot_tridiag',
    'core/creation'
]

print("Nilpotency Tests:")
print("=" * 50)

for matrix_id in nilpotent_matrices:
    A = am.generate(matrix_id, 6)

    # Find nilpotency index (smallest k where A^k = 0)
    for k in range(1, 8):
        Ak = np.linalg.matrix_power(A, k)
        if np.allclose(Ak, 0, atol=1e-10):
            print(f"{matrix_id:25} A^{k} = 0")
            break

Output:

Nilpotency Tests:
==================================================
core/nilpot_triang        A^6 = 0
core/nilpot_tridiag       A^6 = 0
core/creation             A^6 = 0

Example 8: Stochastic Matrices

stochastic = [
    'core/stoch_cesaro',
    'core/stoch_perfect',
    'core/symmstoch'
]

print("Stochastic Matrix Properties:")
print("=" * 60)

for matrix_id in stochastic:
    S = am.generate(matrix_id, 8)

    # Check row sums
    row_sums = S.sum(axis=1)

    # Check column sums
    col_sums = S.sum(axis=0)

    # Determine type
    row_stoch = np.allclose(row_sums, 1.0)
    col_stoch = np.allclose(col_sums, 1.0)

    if row_stoch and col_stoch:
        matrix_type = "Doubly Stochastic"
    elif row_stoch:
        matrix_type = "Row Stochastic"
    else:
        matrix_type = "Unknown"

    print(f"{matrix_id:25} {matrix_type}")

Output:

Stochastic Matrix Properties:
============================================================
core/stoch_cesaro         Row Stochastic
core/stoch_perfect        Doubly Stochastic
core/symmstoch            Doubly Stochastic

Regularization Problems

Example 9: Solving Ill-Posed Problems

import numpy as np
from scipy.linalg import lstsq
import matplotlib.pyplot as plt

# Get Shaw test problem
A, b, x_true = am.generate('regtools/shaw', 100)

# Solve with least squares (no regularization)
x_lstsq = lstsq(A, b)[0]

# Solve with Tikhonov regularization
lambda_reg = 0.01
n = A.shape[1]
A_reg = np.vstack([A, lambda_reg * np.eye(n)])
b_reg = np.concatenate([b, np.zeros(n)])
x_tikhonov = lstsq(A_reg, b_reg)[0]

# Compare solutions
print(f"True solution norm: {np.linalg.norm(x_true):.4f}")
print(f"LSTSQ error: {np.linalg.norm(x_lstsq - x_true):.4f}")
print(f"Tikhonov error: {np.linalg.norm(x_tikhonov - x_true):.4f}")

# Plot
plt.figure(figsize=(12, 4))

plt.subplot(131)
plt.plot(x_true, 'k-', label='True')
plt.title('True Solution')
plt.legend()

plt.subplot(132)
plt.plot(x_lstsq, 'r-', label='LSTSQ')
plt.plot(x_true, 'k--', alpha=0.5, label='True')
plt.title('Least Squares')
plt.legend()

plt.subplot(133)
plt.plot(x_tikhonov, 'b-', label='Tikhonov')
plt.plot(x_true, 'k--', alpha=0.5, label='True')
plt.title('Tikhonov Regularization')
plt.legend()

plt.tight_layout()
plt.savefig('regularization_comparison.png')

Example 10: Discrete Ill-Posed Problems

# Phillips problem
A, b, x = am.generate('regtools/phillips', 80)

# Analyze singular values
U, s, Vt = np.linalg.svd(A)

print("Singular Value Decay:")
print(f"  σ_1 / σ_n = {s[0] / s[-1]:.2e}")
print(f"  Condition number: {np.linalg.cond(A):.2e}")

# Plot singular values
plt.figure(figsize=(10, 4))

plt.subplot(121)
plt.semilogy(s, 'o-')
plt.xlabel('Index')
plt.ylabel('Singular Value')
plt.title('Singular Value Decay')
plt.grid(True)

plt.subplot(122)
plt.plot(x, 'k-', linewidth=2)
plt.xlabel('Index')
plt.ylabel('Value')
plt.title('True Solution')
plt.grid(True)

plt.tight_layout()
plt.savefig('phillips_analysis.png')

MATLAB Gallery Matrices

Example 11: Magic Squares

# Generate magic squares of different sizes
sizes = [3, 4, 5, 6]

fig, axes = plt.subplots(2, 2, figsize=(10, 10))
axes = axes.flatten()

for idx, n in enumerate(sizes):
    M = am.generate('matlab/magic', n)

    # Verify magic property
    magic_sum = n * (n**2 + 1) / 2

    row_sums = M.sum(axis=1)
    col_sums = M.sum(axis=0)
    diag1 = np.trace(M)
    diag2 = np.trace(np.fliplr(M))

    all_equal = (
        np.allclose(row_sums, magic_sum) and
        np.allclose(col_sums, magic_sum) and
        np.isclose(diag1, magic_sum) and
        np.isclose(diag2, magic_sum)
    )

    # Plot
    ax = axes[idx]
    im = ax.imshow(M, cmap='viridis')
    ax.set_title(f'{n}×{n} Magic Square\n(sum={magic_sum:.0f})')

    # Add numbers
    for i in range(n):
        for j in range(n):
            ax.text(j, i, f'{M[i,j]:.0f}',
                   ha='center', va='center',
                   color='white' if M[i,j] < magic_sum else 'black')

    plt.colorbar(im, ax=ax)

plt.tight_layout()
plt.savefig('magic_squares.png')

Example 12: Hilbert Matrix Ill-Conditioning

# Study ill-conditioning of Hilbert matrices
sizes = range(2, 13)
cond_numbers = []

for n in sizes:
    H = am.generate('matlab/hilbert', n)
    cond = np.linalg.cond(H)
    cond_numbers.append(cond)
    print(f"n={n:2d}: cond(H) = {cond:.2e}")

# Plot
plt.figure(figsize=(10, 6))
plt.semilogy(sizes, cond_numbers, 'o-', linewidth=2, markersize=8)
plt.xlabel('Matrix Size')
plt.ylabel('Condition Number')
plt.title('Ill-Conditioning of Hilbert Matrices')
plt.grid(True, which='both', alpha=0.3)
plt.xticks(sizes)
plt.savefig('hilbert_conditioning.png')

Advanced Examples

Example 13: Batch Property Testing

def test_matrix_properties(matrix_id, sizes, am):
    """Test properties across different sizes"""
    results = {}

    for size in sizes:
        try:
            # Generate matrix
            matrix = am.generate(matrix_id, size)

            # Get expected properties
            expected_props = am.properties(matrix_id)

            # Verify properties
            prop_results = MatrixProperties.check_properties(
                matrix, expected_props
            )

            # Store results
            results[size] = {
                'success': True,
                'all_passed': all(prop_results.values()),
                'properties': prop_results
            }

        except Exception as e:
            results[size] = {
                'success': False,
                'error': str(e)
            }

    return results

# Test beta matrix at different sizes
sizes = [10, 20, 50, 100, 200]
results = test_matrix_properties('core/beta', sizes, am)

print("Property Testing Results for 'core/beta':")
print("=" * 60)

for size, result in results.items():
    if result['success']:
        status = "✓" if result['all_passed'] else "✗"
        print(f"Size {size:3d}: {status} All properties verified")
    else:
        print(f"Size {size:3d}: ERROR - {result['error']}")

Example 14: Performance Comparison

import time

def benchmark_matrix_generation(matrix_id, sizes, am):
    """Benchmark matrix generation times"""
    times = []

    for size in sizes:
        start = time.time()
        matrix = am.generate(matrix_id, size)
        elapsed = time.time() - start
        times.append(elapsed)

    return times

# Benchmark different matrices
matrices_to_test = [
    'core/beta',
    'core/fourier',
    'matlab/hilbert',
    'regtools/deriv2'
]

sizes = [100, 200, 500, 1000, 2000]

plt.figure(figsize=(12, 6))

for matrix_id in matrices_to_test:
    times = benchmark_matrix_generation(matrix_id, sizes, am)
    label = matrix_id.split('/')[-1]
    plt.loglog(sizes, times, 'o-', label=label, linewidth=2)

plt.xlabel('Matrix Size')
plt.ylabel('Generation Time (seconds)')
plt.title('Matrix Generation Performance')
plt.legend()
plt.grid(True, which='both', alpha=0.3)
plt.savefig('generation_performance.png')

Example 15: Custom Matrix Collection

class CustomMatrixCollection:
    """Create a custom collection of matrices for specific application"""

    def __init__(self, am):
        self.am = am
        self.collection = {}

    def add_matrix(self, name, matrix_id, size, **params):
        """Add matrix to collection"""
        matrix = self.am.generate(matrix_id, size, **params)
        self.collection[name] = {
            'matrix': matrix,
            'matrix_id': matrix_id,
            'size': size,
            'params': params,
            'properties': self.am.properties(matrix_id)
        }

    def verify_all(self):
        """Verify properties of all matrices"""
        results = {}

        for name, info in self.collection.items():
            props = info['properties']
            matrix = info['matrix']

            verification = MatrixProperties.check_properties(matrix, props)
            results[name] = all(verification.values())

        return results

    def summary(self):
        """Print collection summary"""
        print("Custom Matrix Collection")
        print("=" * 60)

        for name, info in self.collection.items():
            print(f"\n{name}:")
            print(f"  Matrix ID: {info['matrix_id']}")
            print(f"  Size: {info['size']}")
            print(f"  Properties: {', '.join(info['properties'])}")

# Create custom collection
collection = CustomMatrixCollection(am)

# Add matrices for a specific application
collection.add_matrix('test_spd', 'core/beta', 100)
collection.add_matrix('test_ortho', 'core/fourier', 64)
collection.add_matrix('test_sparse', 'regtools/deriv2', 200)
collection.add_matrix('test_stoch', 'core/symmstoch', 50)

# Verify and summarize
collection.summary()

verification_results = collection.verify_all()
print("\n\nVerification Results:")
for name, passed in verification_results.items():
    status = "✓" if passed else "✗"
    print(f"  {status} {name}")

Integration Examples

Example 16: With NumPy/SciPy

import numpy as np
import scipy.linalg as la

# Generate test matrix
A = am.generate('core/beta', 100)

# NumPy operations
eigenvalues = np.linalg.eigvals(A)
print(f"Eigenvalue range: [{eigenvalues.min():.2e}, {eigenvalues.max():.2e}]")

# SciPy operations
U, s, Vt = la.svd(A)
print(f"Condition number: {s[0]/s[-1]:.2e}")

# Matrix functions
sqrtA = la.sqrtm(A)
logA = la.logm(A)

# Verify
print(f"||sqrtA^2 - A|| = {np.linalg.norm(sqrtA @ sqrtA - A):.2e}")

Example 17: Machine Learning Application

from sklearn.decomposition import PCA
from sklearn.preprocessing import StandardScaler

# Generate multiple test matrices
matrices = []
labels = []

for matrix_id in ['core/beta', 'gallery/lehmer', 'matlab/hilbert']:
    for size in [10, 20, 30]:
        M = am.generate(matrix_id, size)

        # Extract features
        features = [
            np.linalg.norm(M, 'fro'),  # Frobenius norm
            np.trace(M),                # Trace
            np.linalg.det(M),           # Determinant
            np.linalg.cond(M),          # Condition number
        ]

        matrices.append(features)
        labels.append(matrix_id)

# Analyze with PCA
scaler = StandardScaler()
X_scaled = scaler.fit_transform(matrices)

pca = PCA(n_components=2)
X_pca = pca.fit_transform(X_scaled)

# Plot
plt.figure(figsize=(10, 6))

for label in set(labels):
    mask = [l == label for l in labels]
    plt.scatter(X_pca[mask, 0], X_pca[mask, 1],
               label=label, s=100, alpha=0.6)

plt.xlabel('PC1')
plt.ylabel('PC2')
plt.title('PCA of Matrix Features')
plt.legend()
plt.grid(True, alpha=0.3)
plt.savefig('matrix_pca.png')

See Also

• Working with Anymatrix - Anymatrix user guide

• Anymatrix Module - Complete API reference

• basic_examples - Basic SuiteSparse examples

• advanced_examples - Advanced usage patterns