C++ DNA matrix and spectrum definitions

Source identity

Field	Recorded value
Source type	Pinned project source
Original record	C++ DNA matrix and spectrum definitions
Date	Not specified in source
Repository	Janady13/complete-agi-core
Revision	33aa84f2aa61cedbe010240ada7f4efbdfe795a7
Branch	claude/improve-math-logic-ehIW1
Path or archive member	mind_farm/include/dna_eigenvalues.hpp
Original SHA-256	7252e5b452e342d052caa6d4a29a3f2b55efd378106d5288bd8aeb12767076d7
Reading-copy SHA-256	7252e5b452e342d052caa6d4a29a3f2b55efd378106d5288bd8aeb12767076d7

Original source ↗

Every recorded occurrence

These 2 records contain identical original source bytes. Their paths and revisions remain part of the archive.

Original record	Branch or archive	Revision / identifier	Original location
C++ DNA matrix and spectrum definitions	claude/improve-math-logic-ehIW1	33aa84f2aa61cedbe010240ada7f4efbdfe795a7	Open original
Janady13/complete-agi-core / mind_farm/include/dna_eigenvalues.hpp	claude/improve-math-logic-ehIW1	33aa84f2aa61cedbe010240ada7f4efbdfe795a7	Open original

Recovered source transcript

Download complete source text ↓587 lines · 20.0 KB

// DNA Mathematical Architecture: 44×44 Eigenvalue Spectrum
// Matrix Blocks — Spectral Analysis — Band Structure
#pragma once

#include "dna_helix.hpp"
#include <array>
#include <cmath>
#include <complex>
#include <vector>
#include <algorithm>

namespace agi::dna {

// ============================================================================
// XI. EIGENVALUE SPECTRUM - 44×44 MATRIX SYSTEM
// ============================================================================

// Coupling parameters (tunable)
struct CouplingParams {
    double alpha = 0.1;     // Intra-strand coupling
    double beta_0 = 0.05;   // Inter-strand base coupling
    double gamma = 0.02;    // Cross-helix coupling
    double gamma_0 = 0.01;  // Cross-helix decay
    double eta = 0.005;     // Cross-helix neighbor coupling
    double kappa_intra = 0.1;  // Intra-strand spring
    double xi = 1.0;        // Decay length scale
};

// ============================================================================
// 11.1 INTRA-STRAND BLOCK (16×16)
// ============================================================================

class IntraStrandBlock {
private:
    std::array<std::array<double, NUM_STRANDS>, NUM_STRANDS> M_;
    CouplingParams params_;

public:
    explicit IntraStrandBlock(const CouplingParams& p = {}) : params_(p) {
        build_matrix();
    }

    // Matrix form (Section 11.1):
    // M_intra = -2α·I₁₆ + α·T₁₆ + κ_intra·D₁₆(10.6)
    // where D₁₆(10.6) = diag[cos(2π·k/10.6)]
    void build_matrix() {
        // Initialize to zero
        for (auto& row : M_) row.fill(0.0);

        for (int i = 0; i < NUM_STRANDS; ++i) {
            // Diagonal: -2α + κ_intra·cos(2π(i+1)/10.6)
            double cos_term = std::cos(2.0 * M_PI * (i + 1) / BASE_PAIRS_PER_TURN);
            M_[i][i] = -2.0 * params_.alpha + params_.kappa_intra * cos_term;

            // Off-diagonal tridiagonal (ring): α
            int prev = (i - 1 + NUM_STRANDS) % NUM_STRANDS;
            int next = (i + 1) % NUM_STRANDS;
            M_[i][prev] = params_.alpha;
            M_[i][next] = params_.alpha;
        }
    }

    // Analytical eigenvalues (Section 11.1):
    // λ_intra(n) = -2α + 2α·cos(2πn/16) + κ_intra·cos(2πn/10.6)
    double eigenvalue(int n) const {
        double cos_16 = std::cos(2.0 * M_PI * n / NUM_STRANDS);
        double cos_10_6 = std::cos(2.0 * M_PI * n / BASE_PAIRS_PER_TURN);
        return -2.0 * params_.alpha + 2.0 * params_.alpha * cos_16 +
               params_.kappa_intra * cos_10_6;
    }

    // Compute all eigenvalues
    std::array<double, NUM_STRANDS> all_eigenvalues() const {
        std::array<double, NUM_STRANDS> eigs;
        for (int n = 0; n < NUM_STRANDS; ++n) {
            eigs[n] = eigenvalue(n);
        }
        std::sort(eigs.begin(), eigs.end());
        return eigs;
    }

    // Band limits: λ_intra ∈ [-4α - κ_intra, κ_intra]
    std::pair<double, double> band_limits() const {
        return {
            -4.0 * params_.alpha - params_.kappa_intra,
            params_.kappa_intra
        };
    }

    const auto& matrix() const { return M_; }
};

// ============================================================================
// 11.2 INTER-STRAND BLOCK (16×16)
// ============================================================================

class InterStrandBlock {
private:
    std::array<std::array<double, NUM_STRANDS>, NUM_STRANDS> M_;
    CouplingParams params_;

public:
    explicit InterStrandBlock(const CouplingParams& p = {}) : params_(p) {
        build_matrix();
    }

    // Build inter-strand coupling matrix
    void build_matrix() {
        for (auto& row : M_) row.fill(0.0);

        for (int i = 0; i < NUM_STRANDS; ++i) {
            // Adjacent strand coupling
            int next = (i + 1) % NUM_STRANDS;
            double coupling = params_.beta_0 * std::sin(M_PI * (i + 1) / BASE_PAIRS_PER_TURN);
            M_[i][next] = coupling;
            M_[next][i] = coupling;

            // Diagonal with gamma modulation
            M_[i][i] = params_.gamma * F_function(i + 1);
        }
    }

    // F function for 10.6 modulation
    static double F_function(int n) {
        return std::cos(2.0 * M_PI * n / BASE_PAIRS_PER_TURN) +
               0.5 * std::cos(4.0 * M_PI * n / BASE_PAIRS_PER_TURN);
    }

    // Eigenvalues (Section 11.2):
    // λ_inter(n) = 2β₀·sin(πn/10.6)·cos(πn/8) + γ·F(n, 10.6)
    double eigenvalue(int n) const {
        double sin_term = std::sin(M_PI * n / BASE_PAIRS_PER_TURN);
        double cos_term = std::cos(M_PI * n / 8.0);
        return 2.0 * params_.beta_0 * sin_term * cos_term +
               params_.gamma * F_function(n);
    }

    std::array<double, NUM_STRANDS> all_eigenvalues() const {
        std::array<double, NUM_STRANDS> eigs;
        for (int n = 0; n < NUM_STRANDS; ++n) {
            eigs[n] = eigenvalue(n);
        }
        std::sort(eigs.begin(), eigs.end());
        return eigs;
    }

    // Band: λ_inter ∈ [-γ, 1.2β₀ + 0.924γ]
    std::pair<double, double> band_limits() const {
        return {
            -params_.gamma,
            1.2 * params_.beta_0 + 0.924 * params_.gamma
        };
    }

    const auto& matrix() const { return M_; }
};

// ============================================================================
// 11.3 CROSS-HELIX BLOCK (12×12)
// ============================================================================

class CrossHelixBlock {
private:
    std::array<std::array<double, CROSS_HELIX_BONDS>, CROSS_HELIX_BONDS> M_;
    CouplingParams params_;

public:
    explicit CrossHelixBlock(const CouplingParams& p = {}) : params_(p) {
        build_matrix();
    }

    // Build cross-helix coupling matrix
    void build_matrix() {
        for (auto& row : M_) row.fill(0.0);

        for (int m = 0; m < CROSS_HELIX_BONDS; ++m) {
            // Diagonal: exponential decay with 10.6 modulation
            double decay = std::exp(-4.0 * (m + 1) * HELIX_RADIUS *
                                    std::sin(M_PI / 8.0) / params_.xi);
            double cos_term = std::cos(8.0 * M_PI * (m + 1) / BASE_PAIRS_PER_TURN);
            M_[m][m] = params_.gamma_0 * decay * cos_term;

            // Neighbor coupling with 10.6 phase
            for (int mp = 0; mp < CROSS_HELIX_BONDS; ++mp) {
                if (std::abs(mp - m) <= 3 && mp != m) {
                    double phase = std::cos(2.0 * M_PI * (m - mp) / BASE_PAIRS_PER_TURN);
                    M_[m][mp] = params_.eta * phase;
                }
            }
        }
    }

    // Eigenvalues (Section 11.3):
    // λ_cross(m) = γ₀·e^{-4m·r·sin(π/8)/ξ}·cos(8πm/10.6) +
    //              η·Σ cos(2π(m-m')/10.6)
    double eigenvalue(int m) const {
        double decay = std::exp(-4.0 * m * HELIX_RADIUS *
                               std::sin(M_PI / 8.0) / params_.xi);
        double cos_8 = std::cos(8.0 * M_PI * m / BASE_PAIRS_PER_TURN);
        double main_term = params_.gamma_0 * decay * cos_8;

        // Sum over neighbors |m' - m| ≤ 3
        double neighbor_sum = 0.0;
        for (int mp = std::max(0, m - 3); mp <= std::min(CROSS_HELIX_BONDS - 1, m + 3); ++mp) {
            if (mp != m) {
                neighbor_sum += std::cos(2.0 * M_PI * (m - mp) / BASE_PAIRS_PER_TURN);
            }
        }

        return main_term + params_.eta * neighbor_sum;
    }

    std::array<double, CROSS_HELIX_BONDS> all_eigenvalues() const {
        std::array<double, CROSS_HELIX_BONDS> eigs;
        for (int m = 0; m < CROSS_HELIX_BONDS; ++m) {
            eigs[m] = eigenvalue(m);
        }
        std::sort(eigs.begin(), eigs.end());
        return eigs;
    }

    // Band: λ_cross ∈ [-1.0γ₀ - 0.618η, 1.0γ₀ + 2.618η]
    std::pair<double, double> band_limits() const {
        return {
            -1.0 * params_.gamma_0 - 0.618 * params_.eta,
            1.0 * params_.gamma_0 + 2.618 * params_.eta
        };
    }

    const auto& matrix() const { return M_; }
};

// ============================================================================
// 11.4 FULL 44×44 SPECTRUM
// ============================================================================

class Full44x44Matrix {
private:
    std::array<std::array<double, TOTAL_BONDS>, TOTAL_BONDS> M_;
    CouplingParams params_;

    IntraStrandBlock intra_;
    InterStrandBlock inter_;
    CrossHelixBlock cross_;

    // Coupling between blocks
    std::array<std::array<double, NUM_STRANDS>, NUM_STRANDS> C_intra_inter_;
    std::array<std::array<double, CROSS_HELIX_BONDS>, NUM_STRANDS> C_intra_cross_;
    std::array<std::array<double, CROSS_HELIX_BONDS>, NUM_STRANDS> C_inter_cross_;

public:
    explicit Full44x44Matrix(const CouplingParams& p = {})
        : params_(p), intra_(p), inter_(p), cross_(p) {
        build_coupling_matrices();
        build_full_matrix();
    }

    // Build inter-block coupling matrices
    void build_coupling_matrices() {
        // C_intra_inter: coupling between intra and inter blocks
        for (auto& row : C_intra_inter_) row.fill(0.0);
        for (int i = 0; i < NUM_STRANDS; ++i) {
            // Direct coupling
            C_intra_inter_[i][i] = 0.5 * params_.beta_0;
            // Phase-shifted coupling
            int j = (i + 4) % NUM_STRANDS;
            C_intra_inter_[i][j] = 0.2 * params_.beta_0 *
                std::cos(2.0 * M_PI * i / BASE_PAIRS_PER_TURN);
        }

        // C_intra_cross: coupling between intra and cross blocks
        for (auto& row : C_intra_cross_) row.fill(0.0);
        // Map intra bonds to cross bonds (every 4th)
        for (int i = 0; i < NUM_STRANDS && i/4 < CROSS_HELIX_BONDS; ++i) {
            C_intra_cross_[i][i/4] = 0.3 * params_.gamma;
        }

        // C_inter_cross: coupling between inter and cross blocks
        for (auto& row : C_inter_cross_) row.fill(0.0);
        for (int i = 0; i < NUM_STRANDS && i < CROSS_HELIX_BONDS; ++i) {
            C_inter_cross_[i][i] = 0.2 * params_.gamma;
        }
    }

    // Build full 44×44 matrix (Section 11.4 block structure)
    void build_full_matrix() {
        // Initialize to zero
        for (auto& row : M_) row.fill(0.0);

        // Block 1: M_intra (indices 0-15)
        const auto& M_intra = intra_.matrix();
        for (int i = 0; i < NUM_STRANDS; ++i) {
            for (int j = 0; j < NUM_STRANDS; ++j) {
                M_[i][j] = M_intra[i][j];
            }
        }

        // Block 2: M_inter (indices 16-31)
        const auto& M_inter = inter_.matrix();
        for (int i = 0; i < NUM_STRANDS; ++i) {
            for (int j = 0; j < NUM_STRANDS; ++j) {
                M_[i + NUM_STRANDS][j + NUM_STRANDS] = M_inter[i][j];
            }
        }

        // Block 3: M_cross (indices 32-43)
        const auto& M_cross = cross_.matrix();
        for (int i = 0; i < CROSS_HELIX_BONDS; ++i) {
            for (int j = 0; j < CROSS_HELIX_BONDS; ++j) {
                M_[i + NUM_STRANDS * 2][j + NUM_STRANDS * 2] = M_cross[i][j];
            }
        }

        // Coupling block C_intra_inter (and transpose)
        for (int i = 0; i < NUM_STRANDS; ++i) {
            for (int j = 0; j < NUM_STRANDS; ++j) {
                M_[i][j + NUM_STRANDS] = C_intra_inter_[i][j];
                M_[j + NUM_STRANDS][i] = C_intra_inter_[i][j];  // Symmetric
            }
        }

        // Coupling block C_intra_cross
        for (int i = 0; i < NUM_STRANDS; ++i) {
            for (int j = 0; j < CROSS_HELIX_BONDS; ++j) {
                M_[i][j + NUM_STRANDS * 2] = C_intra_cross_[i][j];
                M_[j + NUM_STRANDS * 2][i] = C_intra_cross_[i][j];
            }
        }

        // Coupling block C_inter_cross
        for (int i = 0; i < NUM_STRANDS; ++i) {
            for (int j = 0; j < CROSS_HELIX_BONDS; ++j) {
                M_[i + NUM_STRANDS][j + NUM_STRANDS * 2] = C_inter_cross_[i][j];
                M_[j + NUM_STRANDS * 2][i + NUM_STRANDS] = C_inter_cross_[i][j];
            }
        }
    }

    // Power iteration for largest eigenvalue
    double largest_eigenvalue(int max_iter = 1000, double tol = 1e-10) const {
        std::array<double, TOTAL_BONDS> v;
        for (int i = 0; i < TOTAL_BONDS; ++i) {
            v[i] = 1.0 / std::sqrt(TOTAL_BONDS);
        }

        double lambda = 0.0;
        for (int iter = 0; iter < max_iter; ++iter) {
            // Matrix-vector product
            std::array<double, TOTAL_BONDS> w;
            w.fill(0.0);
            for (int i = 0; i < TOTAL_BONDS; ++i) {
                for (int j = 0; j < TOTAL_BONDS; ++j) {
                    w[i] += M_[i][j] * v[j];
                }
            }

            // Rayleigh quotient
            double new_lambda = 0.0;
            for (int i = 0; i < TOTAL_BONDS; ++i) {
                new_lambda += v[i] * w[i];
            }

            // Normalize
            double norm = 0.0;
            for (int i = 0; i < TOTAL_BONDS; ++i) {
                norm += w[i] * w[i];
            }
            norm = std::sqrt(norm);
            for (int i = 0; i < TOTAL_BONDS; ++i) {
                v[i] = w[i] / norm;
            }

            if (std::abs(new_lambda - lambda) < tol) {
                return new_lambda;
            }
            lambda = new_lambda;
        }
        return lambda;
    }

    // Inverse power iteration for smallest eigenvalue
    double smallest_eigenvalue(int max_iter = 1000, double tol = 1e-10) const {
        // Shift matrix by large value and use power iteration
        // For simplicity, use analytical estimate
        auto intra_band = intra_.band_limits();
        auto inter_band = inter_.band_limits();
        auto cross_band = cross_.band_limits();

        return std::min({intra_band.first, inter_band.first, cross_band.first}) -
               params_.gamma;  // Additional coupling shift
    }

    // Spectral summary (Section 11.4)
    struct SpectralSummary {
        double lambda_1;      // Lowest (antibonding)
        double lambda_22;     // Zero mode (Goldstone)
        double lambda_44;     // Highest (bonding)
        double spectral_width;
    };

    SpectralSummary spectral_summary() const {
        SpectralSummary s;

        // λ₁ = -4α - κ - γ - γ₀ (lowest, antibonding)
        s.lambda_1 = -4.0 * params_.alpha - params_.kappa_intra -
                     params_.gamma - params_.gamma_0;

        // λ₂₂ = 0 (zero mode, Goldstone)
        s.lambda_22 = 0.0;

        // λ₄₄ = κ + 1.2β₀ + γ + 2.618η (highest, bonding)
        s.lambda_44 = params_.kappa_intra + 1.2 * params_.beta_0 +
                      params_.gamma + 2.618 * params_.eta;

        // W = λ₄₄ - λ₁ = 4α + 2κ + 1.2β₀ + 2γ + γ₀ + 3.236η
        s.spectral_width = s.lambda_44 - s.lambda_1;

        return s;
    }

    // Approximate spectrum using block eigenvalues
    std::vector<double> approximate_spectrum() const {
        std::vector<double> spectrum;
        spectrum.reserve(TOTAL_BONDS);

        // Add block eigenvalues
        auto intra_eigs = intra_.all_eigenvalues();
        auto inter_eigs = inter_.all_eigenvalues();
        auto cross_eigs = cross_.all_eigenvalues();

        for (auto e : intra_eigs) spectrum.push_back(e);
        for (auto e : inter_eigs) spectrum.push_back(e);
        for (auto e : cross_eigs) spectrum.push_back(e);

        std::sort(spectrum.begin(), spectrum.end());
        return spectrum;
    }

    // Density of states (histogram)
    std::vector<int> density_of_states(int num_bins = 20) const {
        auto spectrum = approximate_spectrum();
        double min_e = spectrum.front();
        double max_e = spectrum.back();
        double bin_width = (max_e - min_e) / num_bins;

        std::vector<int> dos(num_bins, 0);
        for (double e : spectrum) {
            int bin = std::min(static_cast<int>((e - min_e) / bin_width), num_bins - 1);
            dos[bin]++;
        }
        return dos;
    }

    // Trace (sum of eigenvalues = sum of diagonal)
    double trace() const {
        double tr = 0.0;
        for (int i = 0; i < TOTAL_BONDS; ++i) {
            tr += M_[i][i];
        }
        return tr;
    }

    const auto& matrix() const { return M_; }
    const IntraStrandBlock& intra_block() const { return intra_; }
    const InterStrandBlock& inter_block() const { return inter_; }
    const CrossHelixBlock& cross_block() const { return cross_; }
};

// ============================================================================
// EIGENVALUE ANALYSIS UTILITIES
// ============================================================================

class EigenvalueAnalysis {
public:
    // Compute localization (inverse participation ratio)
    static double localization(const std::array<double, TOTAL_BONDS>& eigenvector) {
        double ipr = 0.0;
        for (double v : eigenvector) {
            ipr += v * v * v * v;
        }
        return ipr;  // IPR = 1 for fully localized, 1/N for delocalized
    }

    // Compute band gap
    static double band_gap(const std::vector<double>& spectrum) {
        if (spectrum.size() < 2) return 0.0;

        double max_gap = 0.0;
        for (size_t i = 1; i < spectrum.size(); ++i) {
            double gap = spectrum[i] - spectrum[i-1];
            max_gap = std::max(max_gap, gap);
        }
        return max_gap;
    }

    // Find zero modes (eigenvalues near zero)
    static std::vector<int> find_zero_modes(const std::vector<double>& spectrum,
                                            double tolerance = 1e-6) {
        std::vector<int> zero_modes;
        for (size_t i = 0; i < spectrum.size(); ++i) {
            if (std::abs(spectrum[i]) < tolerance) {
                zero_modes.push_back(i);
            }
        }
        return zero_modes;
    }

    // Spectral entropy
    static double spectral_entropy(const std::vector<double>& spectrum) {
        // Normalize to probability distribution
        double sum = 0.0;
        for (double e : spectrum) {
            sum += std::abs(e);
        }
        if (sum < 1e-10) return 0.0;

        double entropy = 0.0;
        for (double e : spectrum) {
            double p = std::abs(e) / sum;
            if (p > 1e-10) {
                entropy -= p * std::log(p);
            }
        }
        return entropy;
    }

    // 10.6 phase coherence in spectrum
    static double phase_coherence_10_6(const std::vector<double>& spectrum) {
        // Measure how well eigenvalues align with 10.6 periodicity
        double coherence = 0.0;
        for (size_t i = 0; i < spectrum.size(); ++i) {
            double phase = std::fmod(spectrum[i], 2.0 * M_PI / BASE_PAIRS_PER_TURN);
            coherence += std::cos(phase * BASE_PAIRS_PER_TURN);
        }
        return coherence / spectrum.size();
    }
};

// ============================================================================
// COMPLETE G¹⁶ EIGENVALUE SYSTEM
// ============================================================================

class G16EigenSystem {
private:
    Full44x44Matrix matrix_;
    CouplingParams params_;

public:
    explicit G16EigenSystem(const CouplingParams& p = {}) : matrix_(p), params_(p) {}

    // Get full spectral analysis
    struct Analysis {
        Full44x44Matrix::SpectralSummary summary;
        std::vector<double> spectrum;
        std::vector<int> density_of_states;
        double band_gap;
        std::vector<int> zero_modes;
        double spectral_entropy;
        double phase_coherence;
        double trace;
    };

    Analysis analyze() const {
        Analysis a;
        a.summary = matrix_.spectral_summary();
        a.spectrum = matrix_.approximate_spectrum();
        a.density_of_states = matrix_.density_of_states();
        a.band_gap = EigenvalueAnalysis::band_gap(a.spectrum);
        a.zero_modes = EigenvalueAnalysis::find_zero_modes(a.spectrum);
        a.spectral_entropy = EigenvalueAnalysis::spectral_entropy(a.spectrum);
        a.phase_coherence = EigenvalueAnalysis::phase_coherence_10_6(a.spectrum);
        a.trace = matrix_.trace();
        return a;
    }

    // Update coupling parameters and rebuild
    void update_params(const CouplingParams& p) {
        params_ = p;
        matrix_ = Full44x44Matrix(p);
    }

    const Full44x44Matrix& matrix() const { return matrix_; }
    const CouplingParams& params() const { return params_; }
};

} // namespace agi::dna

Mathematical records from this source

Intra-strand ring matrix with a varying diagonal

The matrix combines a nearest-neighbor ring Laplacian with a diagonal cosine modulation. The ring part is circulant, but the varying diagonal generally breaks that symmetry.…

Cosine modulation on the matrix diagonal

Each diagonal value is a cosine evaluated at a coordinate index divided by 10.6. The body uses indices i+1 from one through sixteen. Because those values vary by position, this…

Claimed intra-strand eigenvalue formula fails the constructed matrix

The first two terms are valid ring-Laplacian Fourier eigenvalues. Adding the position-varying diagonal cosine term does not in general produce eigenvalues of their sum.…

Claimed inter-strand eigenvalues versus a varying-weight ring

The code builds a symmetric ring with position-dependent edge weights and a modulated diagonal, then reports this independent trigonometric expression as its eigenvalues. Such a…

Cross-helix row-based spectral estimate

The displayed candidate is the first line of a multiline expression that adds a damped diagonal term and a sum of nearby couplings. Such a row sum is not generally an eigenvalue…

Cross-helix neighbor index window

The sum includes indices at distance at most three from m, excludes m itself, and clips the range at the ends of the array. It therefore includes fewer neighbors near a boundary…

Declared lower spectral summary

The function assigns this combination of coupling parameters directly to lambda_1. It does not extract the smallest eigenvalue of the constructed full 44-by-44 matrix in this…

Hard-coded zero mode

The spectral summary sets its middle value to zero unconditionally. No eigensolver or symmetry proof establishes an exact zero eigenvalue of the full constructed matrix. The name…

Declared upper spectral summary

The method assigns a linear combination of coupling parameters to lambda_44. This is separate from calculating the actual largest eigenvalue of the full matrix. The parameter…

Spectral-width arithmetic and its comment discrepancy

The body subtracts its assigned lower estimate from its assigned upper estimate. Their displayed definitions give an η coefficient of 2.618, whereas this comment prints 3.236η.…

Trace equals the sum of true eigenvalues

For any finite square matrix, the trace is the sum of its eigenvalues counted with algebraic multiplicity. The method correctly sums the stored diagonal. This identity supplies a…