Skip to content

qlroms.utils.metrics

Lightweight diagnostic metrics used to compare ROM trajectories against the FOM.

dominant_frequency(signal, dt)

Frequency of the highest peak of the Hann-windowed periodogram of a 1D signal (mean removed).

The resolution is 1 / (len(signal) * dt); use a long record (many periods) for a sharp estimate.

Source code in qlroms/utils/metrics.py 
39
40
41
42
43
44
45
46
47
def dominant_frequency(signal, dt: float) -> float:
    """Frequency of the highest peak of the Hann-windowed periodogram of a 1D signal (mean removed).

    The resolution is 1 / (len(signal) * dt); use a long record (many periods) for a sharp estimate.
    """
    x = np.asarray(signal, dtype=float)
    x = (x - x.mean()) * np.hanning(x.size)
    power = np.abs(np.fft.rfft(x)) ** 2
    return float(np.fft.rfftfreq(x.size, dt)[1:][power[1:].argmax()])

relative_error_series(X_ref, X_approx, eps=1e-12, *, Mw=None)

Per-snapshot relative error of (N, T) column-trajectory tensors (torch), in the M metric of Mw: ||x - x_hat||_M / ||x||_M (plain 2-norms for None / uniform weights).

Source code in qlroms/utils/metrics.py 
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
def relative_error_series(X_ref, X_approx, eps: float = 1e-12, *, Mw=None) -> np.ndarray:
    """Per-snapshot relative error of (N, T) column-trajectory tensors (torch), in the
    M metric of Mw: ||x - x_hat||_M / ||x||_M (plain 2-norms for None / uniform weights)."""
    import torch

    # cached snapshots are usually CPU while model output sits on the model's device
    X_approx = torch.as_tensor(X_approx, dtype=X_ref.dtype, device=X_ref.device)
    if X_approx.shape != X_ref.shape:
        hint = (" -- looks time-major; pass .T" if X_approx.shape == X_ref.shape[::-1] else "")
        raise ValueError(f"Both trajectories must be (N, T) with COLUMNS = time: got reference "
                         f"{tuple(X_ref.shape)} and approximation {tuple(X_approx.shape)}{hint}.")
    kind = weight_kind(Mw)
    denom = col_norm(Mw, X_ref, kind).clamp_min(eps)
    return (col_norm(Mw, X_ref - X_approx, kind) / denom).cpu().numpy()

one_step_errors_series(X, model)

Open-loop one-step-ahead error: at every snapshot i, step the model ONE step forward from the TRUE state x_i (not from its own previous prediction) and compare to the TRUE next snapshot x_{i+1}. Unlike a free multi-step run, errors here cannot accumulate across steps -- this isolates the one-step fit quality itself (and, for a quantized-local model, shows up as error spikes right at cluster-boundary snapshots).

Source code in qlroms/utils/metrics.py 
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
def one_step_errors_series(X, model) -> np.ndarray:
    """Open-loop one-step-ahead error: at every snapshot i, step the model ONE
    step forward from the TRUE state x_i (not from its own previous prediction) and
    compare to the TRUE next snapshot x_{i+1}. Unlike a free multi-step run,
    errors here cannot accumulate across steps -- this isolates the one-step fit
    quality itself (and, for a quantized-local model, shows up as error spikes right
    at cluster-boundary snapshots)."""
    import torch

    X = X.to(model.device, model.rdtype)
    Nt = X.shape[1]
    Mw = getattr(model, "Mw", None)
    kind = weight_kind(Mw)                 # classified once, not per snapshot
    errs = np.zeros(Nt - 1)
    for i in range(Nt - 1):
        apod_next = model.step(model.project_state(X[:, i:i + 1]))
        x_pred = model.recover_state(apod_next)
        x_true = X[:, i + 1:i + 2]
        errs[i] = float(fro_norm(Mw, x_true - x_pred, kind)
                        / fro_norm(Mw, x_true, kind).clamp_min(1e-12))
    return errs

reconstruction_error_series(X, model)

Pure POD projection+reconstruction relative error per snapshot -- NO dynamics/ time-stepping: for the nearest cluster k, a = Phi_k^T Mw (x - c_k), x_hat = c_k + Phi_k a. Isolates how well the fitted basis alone represents each snapshot, independent of the model's fitted dynamics (contrast one_step_errors_series).

Source code in qlroms/utils/metrics.py 
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
def reconstruction_error_series(X, model) -> np.ndarray:
    """Pure POD projection+reconstruction relative error per snapshot -- NO dynamics/
    time-stepping: for the nearest cluster k, a = Phi_k^T Mw (x - c_k), x_hat = c_k +
    Phi_k a. Isolates how well the fitted basis alone represents each snapshot,
    independent of the model's fitted dynamics (contrast one_step_errors_series)."""
    import torch

    X = X.to(model.device, model.rdtype)
    Nt = X.shape[1]
    if model.K == 1:
        cids = torch.zeros(Nt, dtype=torch.long, device=model.device)
    else:
        cids_np = nearest_labels(model.tree, X.cpu().numpy().T)
        cids = torch.as_tensor(cids_np, dtype=torch.long, device=model.device).reshape(-1)
    X_hat = torch.empty_like(X)
    for k in cids.unique():
        mask = cids == int(k)
        rom_k = model._get_cluster_rom(int(k))
        X_hat[:, mask] = rom_k.recover_state(rom_k.project_state(X[:, mask]))
    return relative_error_series(X, X_hat, Mw=getattr(model, "Mw", None))

prediction_horizon(rel_err, dt, threshold=0.5)

First time (in physical units) at which rel_err exceeds threshold, matching the ql-ROM papers' T_ph definition (normalized-error threshold, default tau=0.5).

Source code in qlroms/utils/metrics.py 
159
160
161
162
163
164
def prediction_horizon(rel_err: np.ndarray, dt: float, threshold: float = 0.5) -> float:
    """First time (in physical units) at which rel_err exceeds `threshold`, matching
    the ql-ROM papers' T_ph definition (normalized-error threshold, default tau=0.5)."""
    over = np.nonzero(np.asarray(rel_err) > threshold)[0]
    n = int(over[0]) if over.size else len(rel_err)
    return n * dt