B-Spline NN/SA - Structural Vibrational Analysis

bspline_sa_v3.py

"""


bspline_sa_v3.py / self-aware b-spline nn for structural vibrational analytics



___________________________________________________________________________________________




             Input Trajectory Window / B, T, SEQ_LEN
                                │
                                ▼
                  ┌───────────────────────────┐
                  │  batch_hilbert_transform  │  ◄── Fast GPU FFT Block
                  └─────────────┬─────────────┘
                                │
               ┌────────────────┴────────────────┐
               ▼ (Real Component)                ▼ (Quadrature Phase Component)
         [Raw Trajectory]                  [Hilbert Companion]
               │                                 │
               └────────────────┬────────────────┘
                                ▼
                   Stacked Input (B, 2, SEQ_LEN)
                                │
                                ▼
                  ┌───────────────────────────┐
                  │   ControlPointExtractor   │  ◄── Upgraded 2-Channel CNN
                  └─────────────┬─────────────┘
                                │
                                ▼
                     Control Points (C) (B, CP_K, D)
                                │
         ┌──────────────────────┼──────────────────────┐
         ▼                      ▼                      ▼
   Truncated SVD       Orthonormal Alignment    Spectral STFT Features
 (U_r, S_r matrices)   (V_r matrix + QR mask)     (Frequency Energy)
         │                      │                      │
         ▼                      ▼                      ▼
   [SVD Vector]         [Rotation Vector]      [Spectral Vector]
    (128-dim)               (128-dim)              (128-dim)
         │                      │                      │
         └──────────────────────┼──────────────────────┘
                                ▼
                ┌───────────────────────────────┐
                │       Gating And Fusion       │ ◄── Dynamic Softmax Weighting
                └───────────────┬───────────────┘
                                │
                                ▼
                     Fused Window Embedding (128-dim)
                                │
                                ▼
                ┌───────────────────────────────┐
                │       SmallTransformer        │ ◄── Temporal Sequence Modeling
                └───────────────┬───────────────┘
                                │
                                ▼
                       Final Latent State
                                │
                ┌───────────────┴───────────────┐
                ▼                               ▼
       ┌─────────────────┐             ┌─────────────────┐
       │   Action Head   │             │ Prediction Head │ ◄── Self-Awareness Target
       └─────────────────┘             └─────────────────┘



___________________________________________________________________________________________




While processing the data through this pipeline, the following parallel loss calculations
run dynamically during training to ensure structural integrity and convergence:


                              Active Training Loss Loop
                                          │
       ┌──────────────────┬───────────────┼───────────────┬──────────────────┐
       ▼                  ▼               ▼               ▼                  ▼
┌──────────────┐   ┌──────────────┐┌──────────────┐┌──────────────┐   ┌──────────────┐
│ CrossEntropy │   │   MSE Loss   ││  Orthogonal  ││ SVD Entropy  │   │  Gate Consistency
│ (Action Head)│   │ (Self-Aware) ││  Reg (V_r)   ││  Reg (S_r)   │   │  Variance Loss
└──────┬───────┘   └──────┬───────┘└──────┬───────┘└──────┬───────┘   └──────┬───────┘
       │                  │               │               │                  │
       └──────────────────┴───────────────┼───────────────┴──────────────────┘
                                          ▼
                               [Total Aggregated Loss] ──► Scaler Backprop



___________________________________________________________________________________________





Active Optimization & Self-Awareness Targets:



>>Low-Rank Fallback Autoencoder. If the truncated_svd calculation hits a matrix singularity
  on the GPU during a specific step, the system catches the error and passes the control
  points through the LowRankAutoencoder reconstruction layout instead.
  
>>Orthonormal QR Head Constraints. The continuous coordinate mapping framework checks
  Q^T Q=I at every time-step.

>>SVD Entropy Control. The singular values continue to be compressed into a tight, high-energy
  profile, forcing clean representation manifolds.
  
>>Gating Variance Smoothing. The model continues to penalize erratic frame-by-frame expert
  switching via the sequence variance penalty.
  
>>Auxiliary Predictive Target. The model challenges its internal Transformer states to
  actively forecast what the control points of its last trajectory step look like.
  In this version, it evaluates this self-prediction against a target extracted from
  the dual-channel phase-coupled information.

"""

import os
import time
import math
import numpy as np
from collections import defaultdict

import torch
import torch.nn as nn
import torch.nn.functional as F

# Configurations
device = torch.device("cuda" if torch.cuda.is_available() else "cpu")
PRINT_INTERVAL = 50
SEQ_LEN = 128
CP_K = 8
D = 3
SVD_RANK = 3
BATCH = 32
ACTION_DIM = 4
TRANSFORMER_DIM = 128
TRANSFORMER_LAYERS = 2
TRANSFORMER_HEADS = 4
CHECKPOINT_DIR = "./channel_state_v3"
os.makedirs(CHECKPOINT_DIR, exist_ok=True)

# Training runtime flags
USE_AUTOENCODER_FALLBACK = True
USE_MIXED_PRECISION = True

# Hilbert Transform Feature Engineering
def batch_hilbert_transform(x):
    """
    Computes the 1D Hilbert Transform along the last dimension for a batch of signals
    Utilizes a fast FFT-based analytic signal approximation framework
    Input x shape: (B, SeqLen) or (B, T, SeqLen)
    Output shape: ...as input, representing the quadrature component
    """
    input_shape = x.shape
    # Flatten leading dimensions to handle both (B, SeqLen) and (B, T, SeqLen)
    x_flat = x.view(-1, input_shape[-1])
    N = x_flat.shape[-1]
    
    # Forward real FFT
    Xf = torch.fft.fft(x_flat, n=N, dim=-1)
    
    # Create the Hilbert step function multiplier
    h = torch.zeros(N, dtype=Xf.dtype, device=Xf.device)
    if N % 2 == 0:
        h[0] = 1
        h[N // 2] = 1
        h[1 : N // 2] = 2
    else:
        h[0] = 1
        h[1 : (N + 1) // 2] = 2
        
    # Apply multiplier across the batch sequence
    Xf = Xf * h.unsqueeze(0)
    
    # Inverse FFT yields the full complex Analytic Signal: Analytic = x + i*hilbert(x)
    analytic = torch.fft.ifft(Xf, dim=-1)
    
    # ***Extract the imaginary part as the true quadrature signal***
    quadrature = analytic.imag
    return quadrature.view(*input_shape)

# Utilities
def truncated_svd(C, r):
    Cc = C - C.mean(dim=1, keepdim=True)
    U, S, Vh = torch.linalg.svd(Cc, full_matrices=False)
    U_r = U[:, :, :r]
    S_r = S[:, :r]
    V_r = Vh[:, :r, :]
    return U_r, S_r, V_r

def sv_entropy_loss(S, eps=1e-8):
    norm = S / (S.sum(dim=-1, keepdim=True) + eps)
    ent = - (norm * torch.log(norm + eps)).sum(dim=-1)
    return -ent.mean()

def orthogonality_loss(Q):
    B = Q.shape[0]
    Id = torch.eye(Q.shape[-1], device=Q.device).unsqueeze(0).expand(B, -1, -1)
    return F.mse_loss(Q.transpose(-2,-1) @ Q, Id)

def spectral_features_stft(C, n_fft=16, hop=4, topk=6):
    B, k, d = C.shape
    frames = []
    for start in range(0, max(1, k - n_fft + 1), hop):
        seg = C[:, start:start+n_fft, :]
        Cf = torch.fft.rfft(seg, dim=1)
        mag = torch.abs(Cf)
        phase = torch.angle(Cf)
        frames.append(torch.cat([mag.real, phase.real], dim=-1).reshape(B, -1))
    if len(frames) == 0:
        Cf = torch.fft.rfft(C, dim=1)
        mag = torch.abs(Cf)
        phase = torch.angle(Cf)
        feat = torch.cat([mag, phase], dim=1).reshape(B, -1)
    else:
        feat = torch.cat(frames, dim=-1)
    feat = feat.float()
    feat = feat / (feat.norm(dim=-1, keepdim=True) + 1e-9)
    target_dim = topk * d * 2
    if feat.shape[-1] >= target_dim:
        feat = feat[:, :target_dim]
    else:
        pad = target_dim - feat.shape[-1]
        feat = F.pad(feat, (0, pad))
    pal_loss = torch.tensor(0.0, device=C.device)
    return feat, pal_loss

def save_channel_state(channel_id, state_dict):
    path = os.path.join(CHECKPOINT_DIR, f"channel_{channel_id}.npz")
    np_dict = {}
    for k, v in state_dict.items():
        if isinstance(v, torch.Tensor):
            np_dict[k] = v.detach().cpu().numpy()
        else:
            np_dict[k] = v
    np.savez(path, **np_dict)

def load_channel_state(channel_id):
    path = os.path.join(CHECKPOINT_DIR, f"channel_{channel_id}.npz")
    if not os.path.exists(path):
        return None
    data = np.load(path)
    out = {}
    for k in data.files:
        out[k] = torch.tensor(data[k], device=device)
    return out

# Modules
class ControlPointExtractor(nn.Module):
    def __init__(self, seq_len, k, d, base_channels=64):
        super().__init__()
        self.k = k
        self.d = d
        # Upgraded input channels from 1 to 2 to handle [Raw Trajectory, Hilbert Phase]
        self.net = nn.Sequential(
            nn.Conv1d(2, base_channels, kernel_size=5, padding=2),
            nn.ReLU(),
            nn.Conv1d(base_channels, base_channels*2, kernel_size=5, padding=2),
            nn.ReLU(),
            nn.AdaptiveAvgPool1d(k),
            nn.Flatten(),
            nn.Linear(base_channels*2 * k, k * d))
    def forward(self, x_two_channel):
        # Expected shapes: B, 2, SEQ_LEN
        out = self.net(x_two_channel)
        return out.view(x_two_channel.size(0), self.k, self.d)

class LowRankAutoencoder(nn.Module):
    def __init__(self, k, d, r):
        super().__init__()
        self.encoder = nn.Sequential(
            nn.Flatten(),
            nn.Linear(k*d, max(128, r*(d+k))),
            nn.ReLU(),
            nn.Linear(max(128, r*(d+k)), r * d))
        self.decoder = nn.Sequential(
            nn.Linear(r * d, max(128, r*(d+k))),
            nn.ReLU(),
            nn.Linear(max(128, r*(d+k)), k*d),
            nn.Unflatten(1, (k, d)))
            
    def forward(self, C):
        z = self.encoder(C)
        rec = self.decoder(z)
        return rec, z

def qr_orthonormal_from_matrix(M):
    Q, R = torch.linalg.qr(M)
    diag = torch.sign(torch.diagonal(R, dim1=-2, dim2=-1))
    Q = Q * diag.unsqueeze(-2)
    return Q

class ConditionalOrthonormalHead(nn.Module):
    def __init__(self, r, d, hidden=128):
        super().__init__()
        self.r = r
        self.d = d
        self.mlp = nn.Sequential(
            nn.Linear(r * d, hidden),
            nn.ReLU(),
            nn.Linear(hidden, d * d))
            
    def forward(self, V_r):
        B = V_r.shape[0]
        vflat = V_r.reshape(B, -1)
        M = self.mlp(vflat)
        M = M.view(B, self.d, self.d)
        Q = qr_orthonormal_from_matrix(M)
        V_rot = (Q.unsqueeze(1) @ V_r.permute(0,2,1)).permute(0,2,1)
        return V_rot.reshape(B, -1), Q

class SVDHead(nn.Module):
    def __init__(self, k, r):
        super().__init__()
        self.k = k
        self.r = r
    def forward(self, U_r, S_r):
        return torch.cat([U_r.reshape(U_r.shape[0], -1), S_r], dim=-1)

class SmallTransformer(nn.Module):
    def __init__(self, dim=TRANSFORMER_DIM, layers=2, heads=4, ff=256, max_len=32):
        super().__init__()
        self.dim = dim
        encoder_layer = nn.TransformerEncoderLayer(d_model=dim, nhead=heads, dim_feedforward=ff, activation='relu', batch_first=True)
        self.encoder = nn.TransformerEncoder(encoder_layer, num_layers=layers)
        self.pos_emb = nn.Parameter(torch.randn(max_len, dim) * 0.01)
    def forward(self, x):
        B, T, _ = x.shape
        pe = self.pos_emb[:T, :].unsqueeze(0).expand(B, -1, -1)
        return self.encoder(x + pe)

class GatingAndFusion(nn.Module):
    def __init__(self, in_dims, hidden=256, out_dim=TRANSFORMER_DIM):
        super().__init__()
        self.n = len(in_dims)
        total = sum(in_dims)
        self.gate_mlp = nn.Sequential(
            nn.Linear(total, hidden),
            nn.ReLU(),
            nn.Linear(hidden, self.n))
        self.projs = nn.ModuleList([nn.Linear(in_dims[i], out_dim) for i in range(self.n)])
        
    def forward(self, embeddings):
        concat = torch.cat(embeddings, dim=-1)
        logits = self.gate_mlp(concat)
        weights = F.softmax(logits, dim=-1)
        outs = torch.stack([self.projs[i](embeddings[i]) for i in range(self.n)], dim=1)
        fused = (weights.unsqueeze(-1) * outs).sum(dim=1)
        return fused, weights

class MetaController(nn.Module):
    def __init__(self, embed_dim=TRANSFORMER_DIM, transformer_layers=TRANSFORMER_LAYERS):
        super().__init__()
        self.transformer = SmallTransformer(dim=embed_dim, layers=transformer_layers, heads=TRANSFORMER_HEADS, max_len=64)
        self.fc = nn.Sequential(nn.Linear(embed_dim, 256), nn.ReLU())
        self.action_head = nn.Linear(256, ACTION_DIM)
        self.pred_head = nn.Linear(256, CP_K * D)
    def forward(self, path_emb):
        tr_out = self.transformer(path_emb)
        last = tr_out[:, -1, :]
        h = self.fc(last)
        logits = self.action_head(h)
        pred_cp = self.pred_head(h).view(-1, CP_K, D)
        return logits, pred_cp, tr_out

class FullModelV3(nn.Module):
    def __init__(self):
        super().__init__()
        self.extractor = ControlPointExtractor(SEQ_LEN, CP_K, D)
        self.svd_head = SVDHead(CP_K, SVD_RANK)
        self.conditional_head = ConditionalOrthonormalHead(SVD_RANK, D)
        self.autoencoder = LowRankAutoencoder(CP_K, D, SVD_RANK) if USE_AUTOENCODER_FALLBACK else None
        self.spectral_proj = nn.Linear(6 * D * 2, 64)
        svd_dim = CP_K * SVD_RANK + SVD_RANK
        learned_dim = SVD_RANK * D
        spec_dim = 64
        self.proj_svd = nn.Linear(svd_dim, 128)
        self.proj_learned = nn.Linear(learned_dim, 128)
        self.proj_spec = nn.Linear(spec_dim, 128)
        self.gate = GatingAndFusion([128, 128, 128], hidden=256, out_dim=TRANSFORMER_DIM)
        self.meta = MetaController(embed_dim=TRANSFORMER_DIM, transformer_layers=TRANSFORMER_LAYERS)
        self.orth_weight = 1e-2
        self.sv_ent_weight = 1e-2
        self.gate_consistency_weight = 1e-2
        self.ema_decay = 0.98
        self.register_buffer("_dummy", torch.tensor(0.0))

    def forward(self, x_path):
        B, T, L = x_path.shape
        # Generate the orthogonal quadrature signal for the entire time array block
        x_quadrature = batch_hilbert_transform(x_path)
        per_t_emb = []
        pal_losses = []
        orth_losses = []
        sv_ent_losses = []
        gate_weights_all = []
        for t in range(T):
            # Isolate raw path data frame and its phase companion
            x_raw_t = x_path[:, t, :]   # (B, SEQ_LEN)
            x_quad_t = x_quadrature[:, t, :] # (B, SEQ_LEN)
            # Stack along a new structural channel dimension: B, 2, SEQ_LEN
            x_two_channel = torch.stack([x_raw_t, x_quad_t], dim=1)
            # Extract geometric spline coordinates from the phase-coupled inputs
            C = self.extractor(x_two_channel)
            # Stable SVD Decomposition
            try:
                U_r, S_r, V_r = truncated_svd(C, SVD_RANK)
            except Exception:
                if self.autoencoder is not None:
                    rec, z = self.autoencoder(C)
                    U_r, S_r, V_r = truncated_svd(rec, SVD_RANK)
                else:
                    C_pert = C + 1e-6 * torch.randn_like(C)
                    U_r, S_r, V_r = truncated_svd(C_pert, SVD_RANK)
            # Geometric Loss Calculations
            sv_ent_losses.append(sv_entropy_loss(S_r) * self.sv_ent_weight)
            e_svd = self.svd_head(U_r, S_r)
            e_svd = F.relu(self.proj_svd(e_svd))
            e_learned, Q = self.conditional_head(V_r)
            e_learned = F.relu(self.proj_learned(e_learned))
            orth_losses.append(orthogonality_loss(Q) * self.orth_weight)
            spec_feat, pal_loss = spectral_features_stft(C, n_fft=16, hop=4, topk=6)
            pal_losses.append(pal_loss)
            spec_p = F.relu(self.spectral_proj(spec_feat))
            spec_p = F.relu(self.proj_spec(spec_p))
            fused, weights = self.gate([e_svd, e_learned, spec_p])
            per_t_emb.append(fused)
            gate_weights_all.append(weights)
        path_emb = torch.stack(per_t_emb, dim=1)
        logits, pred_cp, tr_out = self.meta(path_emb)
        # Loss aggregations
        orth_loss = torch.stack(orth_losses).mean() if len(orth_losses) else torch.tensor(0.0, device=device)
        sv_ent_loss = torch.stack(sv_ent_losses).mean() if len(sv_ent_losses) else torch.tensor(0.0, device=device)
        pal_loss_mean = torch.stack(pal_losses).mean() if len(pal_losses) else torch.tensor(0.0, device=device)
        gate_weights = torch.stack(gate_weights_all, dim=1).mean(dim=1)
        gate_var = torch.var(torch.stack(gate_weights_all, dim=1), dim=1).mean()
        gate_consistency = gate_var * self.gate_consistency_weight
        reg_loss = orth_loss + sv_ent_loss + pal_loss_mean * 1e-2 + gate_consistency
        return logits, pred_cp, reg_loss, gate_weights.mean(dim=0).detach().cpu().numpy()

    def save_channel_summary(self, channel_id, C):
        U_r, S_r, V_r = truncated_svd(C, SVD_RANK)
        state = load_channel_state(channel_id)
        summary = {
            "U_mean": U_r.mean(dim=0).detach().cpu().numpy(),
            "S_mean": S_r.mean(dim=0).detach().cpu().numpy(),
            "V_mean": V_r.mean(dim=0).detach().cpu().numpy(),
            "ts": time.time()}
        if state is not None:
            prev_U = state.get("U_mean").cpu().numpy()
            prev_S = state.get("S_mean").cpu().numpy()
            prev_V = state.get("V_mean").cpu().numpy()
            summary["U_mean"] = (self.ema_decay * prev_U + (1-self.ema_decay) * summary["U_mean"])
            summary["S_mean"] = (self.ema_decay * prev_S + (1-self.ema_decay) * summary["S_mean"])
            summary["V_mean"] = (self.ema_decay * prev_V + (1-self.ema_decay) * summary["V_mean"])
        save_channel_state(channel_id, summary)

# Synthetic dataset for path windows
def make_synthetic_path_batch(batch, T=6, seq_len=SEQ_LEN):
    batch_windows = []
    actions = []
    for i in range(batch):
        total_len = seq_len + T - 1
        freq = np.random.uniform(0.5, 3.0)
        phase = np.random.uniform(0, 2*np.pi)
        amp = np.random.uniform(0.5, 1.5)
        t = np.linspace(0, 1, total_len)
        x = amp * np.cos(2*np.pi*freq*t + phase)
        y = amp * np.sin(2*np.pi*freq*t + phase)
        z = 0.2 * np.sin(4*np.pi*freq*t + phase*0.5)
        traj = np.stack([x,y,z], axis=1)
        proj = np.random.randn(3); proj /= np.linalg.norm(proj)
        scalar = traj @ proj
        scalar += 0.02 * np.random.randn(*scalar.shape)
        windows = []
        for s in range(T):
            win = scalar[s:s+seq_len]
            windows.append(win.astype(np.float32))
        batch_windows.append(np.stack(windows, axis=0))
        actions.append(int(min(3, int(freq // 0.75))))
    x = torch.tensor(np.stack(batch_windows, axis=0), device=device)
    actions = torch.tensor(actions, device=device)
    return x, actions

def train(model, steps=1000, batch=BATCH, T=6, lr=1e-3):
    model.train()
    opt = torch.optim.Adam(model.parameters(), lr=lr)
    scaler = torch.cuda.amp.GradScaler(enabled=USE_MIXED_PRECISION and device.type=='cuda')
    for step in range(steps):
        x, actions = make_synthetic_path_batch(batch, T=T)
        with torch.cuda.amp.autocast(enabled=USE_MIXED_PRECISION and device.type=='cuda'):
            logits, pred_cp, reg_loss, gate_avg = model(x)
            action_loss = F.cross_entropy(logits, actions)
            # > Internal self-awareness alignment target verification
            # > Create dummy channel vector to simulate target extraction mapping
            dummy_target_in = torch.stack([x[:, -1, :], batch_hilbert_transform(x[:, -1, :])], dim=1)
            target_cp = model.extractor(dummy_target_in).detach()
            pred_loss = F.mse_loss(pred_cp, target_cp)
            loss = action_loss + pred_loss + reg_loss
        opt.zero_grad()
        scaler.scale(loss).backward()
        scaler.unscale_(opt)
        torch.nn.utils.clip_grad_norm_(model.parameters(), 5.0)
        scaler.step(opt)
        scaler.update()
        if step % PRINT_INTERVAL == 0:
            print(f"step {step:5d} loss {loss.item():.4f} act {action_loss.item():.4f} pred {pred_loss.item():.4f} reg {reg_loss.item():.6f} gate_avg {gate_avg}")
    return model

def demo_save_load(model):
    x, _ = make_synthetic_path_batch(1, T=1)
    quad = batch_hilbert_transform(x[:,0,:])
    ch2 = torch.stack([x[:,0,:], quad], dim=1)
    C = model.extractor(ch2)
    model.save_channel_summary("demo_v3", C)
    st = load_channel_state("demo_v3")
    print("Loaded demo state keys:", list(st.keys()))

if __name__ == "__main__":
    print("Device identified:", device)
    model = FullModelV3().to(device)
    trained = train(model, steps=400, batch=BATCH, T=6, lr=1e-3)
    demo_save_load(trained)

The model above can be placed in authorization matrices with linear chaining across any llm nn, however qbit driven representation preferred.