"""Restricted metriplectic product of substrate pulse and EPI relaxation."""

import networkx as nx
import numpy as np
import pytest

from tnfr.physics import verify_metriplectic_product
from tnfr.physics.structural_diffusion import compute_diffusion_energy


def _state(graph, epi, frequency):
    for node, value, nu_f in zip(graph, epi, frequency):
        graph.nodes[node].update(
            EPI=float(value),
            nu_f=float(nu_f),
            theta=0.3 * node,
            delta_nfr=(-1.0) ** node * 0.2,
        )
    return graph


def test_product_preserves_hamiltonian_and_dissipates_dirichlet_energy():
    graph = _state(nx.path_graph(4), [2.0, -1.0, 3.0, 0.5], [0.4, 0.8, 1.2, 2.0])

    result = verify_metriplectic_product(graph)

    assert result.is_decoupled_metriplectic_bridge
    assert not result.couples_sectors
    assert result.state_dimension == 5 * len(graph)
    assert result.hamiltonian_derivative == pytest.approx(0.0, abs=1e-12)
    assert result.dirichlet_derivative < 0.0
    assert result.poisson_antisymmetry_residual < 1e-12
    assert result.minimum_dissipative_eigenvalue >= 0.0
    assert not result.stored_pressure_matches_epi_channel


def test_product_reports_when_stored_pressure_is_the_epi_channel():
    graph = _state(nx.path_graph(2), [0.0, 1.0], [0.5, 2.0])
    graph.nodes[0]["delta_nfr"] = 1.0
    graph.nodes[1]["delta_nfr"] = -1.0

    result = verify_metriplectic_product(graph)

    assert result.stored_pressure_matches_epi_channel
    assert result.stored_pressure_consistency_residual < 1e-12
    assert result.is_decoupled_metriplectic_bridge


def test_product_aligns_state_with_canonical_phase_space_node_order():
    graph = nx.Graph()
    graph.add_nodes_from([2, 1, 0])
    graph.add_edges_from([(2, 1), (1, 0)])
    for node in graph:
        graph.nodes[node].update(
            EPI=float(node),
            nu_f=1.0 + 0.25 * node,
            theta=0.3 * node,
            delta_nfr=(-1.0) ** node * 0.2,
        )

    result = verify_metriplectic_product(graph)

    assert result.nodes == (0, 1, 2)
    assert result.is_decoupled_metriplectic_bridge


def test_metriplectic_degeneracy_conditions_are_exact_block_identities():
    graph = _state(nx.cycle_graph(5), np.arange(5.0), np.linspace(0.5, 1.5, 5))

    result = verify_metriplectic_product(graph)

    assert result.poisson_dissipative_degeneracy == pytest.approx(0.0)
    assert result.metric_hamiltonian_degeneracy == pytest.approx(0.0)
    assert result.substrate_velocity_residual < 1e-12
    assert result.epi_velocity_residual < 1e-12


def test_uniform_epi_is_dissipative_equilibrium_while_pulse_can_continue():
    graph = _state(nx.path_graph(3), [7.0] * 3, [0.5, 1.0, 2.0])

    result = verify_metriplectic_product(graph)

    assert result.dirichlet_functional == pytest.approx(0.0)
    assert result.dirichlet_derivative == pytest.approx(0.0)
    assert np.linalg.norm(result.state_velocity[: 4 * len(graph)]) > 0.0
    np.testing.assert_allclose(result.state_velocity[4 * len(graph) :], 0.0)


def test_product_rejects_nonpositive_capacity_and_directed_transport():
    zero = _state(nx.path_graph(3), [0.0, 1.0, 2.0], [1.0, 0.0, 1.0])
    with pytest.raises(ValueError, match="positive finite capacity"):
        verify_metriplectic_product(zero)

    directed = _state(nx.DiGraph([(0, 1), (1, 0), (1, 2)]), [0.0, 1.0, 2.0], [1.0] * 3)
    with pytest.raises(ValueError, match="symmetric adjacency"):
        verify_metriplectic_product(directed)


@pytest.mark.parametrize("graph", [nx.Graph(), nx.empty_graph(1)])
def test_product_rejects_empty_or_singleton_support_explicitly(graph):
    with pytest.raises(ValueError, match="at least two supported nodes"):
        verify_metriplectic_product(graph)


@pytest.mark.parametrize("tolerance", [True, "1e-10", 1e-10j, None])
def test_product_rejects_nonreal_tolerance(tolerance):
    graph = _state(nx.path_graph(2), [0.0, 1.0], [1.0, 1.0])

    with pytest.raises(ValueError, match="finite and positive"):
        verify_metriplectic_product(graph, tolerance=tolerance)


@pytest.mark.parametrize("attribute", ["EPI", "nu_f", "theta", "delta_nfr"])
def test_product_rejects_boolean_nodal_channels_before_extraction(
    attribute, monkeypatch
):
    graph = _state(nx.path_graph(2), [0.0, 1.0], [1.0, 1.0])
    graph.nodes[1][attribute] = True

    def unexpected_extraction(_graph):
        raise AssertionError("invalid nodal channels reached field extraction")

    monkeypatch.setattr(
        "tnfr.physics.metriplectic.extract_phase_space_point", unexpected_extraction
    )

    with pytest.raises(ValueError, match="finite real"):
        verify_metriplectic_product(graph)


def test_product_rejects_nonrepresentable_finite_intermediates():
    graph = _state(nx.path_graph(2), [1e300, -1e300], [1.0, 1.0])
    graph.edges[0, 1]["weight"] = 1e100

    with pytest.raises(ValueError, match="floating-point range"):
        verify_metriplectic_product(graph)


def test_positive_mobility_underflow_cannot_certify_zero_nodal_velocity():
    graph = _state(nx.path_graph(2), [1.0, 0.0], [1e-200, 1e-200])
    graph.edges[0, 1].update(weight=1e200, length=1.0)
    graph.nodes[0]["delta_nfr"] = -1.0
    graph.nodes[1]["delta_nfr"] = 1.0
    # x'=(-1e-200,1e-200) and E_D'=-2 are representable, but the product
    # tensor's positive mobility 1e-400 cannot be represented in binary64.
    with pytest.raises(ValueError, match="mobility.*floating-point range"):
        verify_metriplectic_product(graph)


@pytest.mark.parametrize("epi_power,capacity_power", [(-600, 0), (-300, -900)])
def test_nonzero_diffusion_energy_or_rate_cannot_be_certified_as_zero(
    epi_power, capacity_power
):
    graph = _state(nx.path_graph(2), [2.0**epi_power, 0.0], [2.0**capacity_power] * 2)
    # E=x^2/2 and x'=(-nu*x,nu*x): the first case loses positive energy;
    # the second loses a nonzero nodal rate. Both require explicit refusal.
    with pytest.raises(ValueError, match="floating-point range"):
        verify_metriplectic_product(graph)


@pytest.mark.parametrize("offset", [0.0, 1e12])
def test_product_inherits_edge_balance_with_explicit_node_alignment(offset):
    graph = nx.Graph()
    graph.add_nodes_from([2, 0, 1])
    graph.add_weighted_edges_from([(2, 0, 2.0), (0, 1, 0.5), (0, 0, 3.0)])
    _state(graph, [offset + 1.5, offset - 0.25, offset + 0.75], [0.25, 1.0, 2.0])
    expected = compute_diffusion_energy(graph)
    result = verify_metriplectic_product(graph)

    assert expected.nodes == [2, 0, 1]
    assert result.nodes == (0, 1, 2)
    order = [expected.nodes.index(node) for node in result.nodes]
    np.testing.assert_array_equal(result.state_velocity[12:], expected.epi_rate[order])
    np.testing.assert_array_equal(
        np.diag(result.dissipative_tensor)[12:], expected.mobility[order]
    )
    # Independent edge energy: 1*(7/4)^2 + (1/4)*1^2 = 53/16.
    assert result.dirichlet_functional == expected.energy == 53 / 16
    assert result.dirichlet_derivative == expected.energy_rate
    assert result.dirichlet_derivative == pytest.approx(-1915 / 352)
    np.testing.assert_allclose(result.state_velocity[12:], [8 / 11, -2.0, -7 / 16])
    assert result.is_decoupled_metriplectic_bridge