∀x' : ‖x' - x‖_∞ ≤ ε → argmax(f(x')) = argmax(f(x))

∀x ∈ S : g(f(x)) ≤ 0

|P(ŷ=1|z=0, y=1) - P(ŷ=1|z=1, y=1)| ≤ δ

def relu_abstract(l, u):
    """区间传播经过 ReLU 抽象转换"""
    return max(0, l), max(0, u)

x̂ = c + Σ α_i · ε_i    where α_i ∈ [-1, 1]

import numpy as np

class Zonotope:
    """Zonotope 抽象域：中心向量 + 生成器矩阵"""
    def __init__(self, center, generators):
        # center: shape (n,)
        # generators: shape (n, m)，每列是一个噪声符号的系数
        self.c = center
        self.G = generators

    def relu_abstract(self):
        """Zonotope 上的 ReLU 抽象转换"""
        n = self.c.shape[0]
        new_G = np.hstack([self.G, np.zeros((n, 1))])
        new_c = self.c.copy()

        for i in range(n):
            l, u = self.bounds_of(i)
            if l >= 0:
                continue  # 肯定为正，直接通过
            elif u <= 0:
                new_c[i] = 0  # 肯定为负，置零
                new_G[i, :] = 0
            else:
                # 引入新的噪声符号处理过零点
                lam = u / (u - l)
                new_c[i] = lam * (new_c[i] + l)
                new_G[i, -1] = -lam * l

        return Zonotope(new_c, new_G)

    def bounds_of(self, dim):
        """计算某一维的精确区间界"""
        abs_sum = np.sum(np.abs(self.G[dim, :]))
        return self.c[dim] - abs_sum, self.c[dim] + abs_sum

P = {x | Ax ≤ b}

from z3 import *

def encode_network(solver, layers, input_vars, output_vars):
    """将前馈神经网络编码为 Z3 约束系统"""
    current_vars = input_vars

    for i, (W, b, activation) in enumerate(layers):
        n_out = W.shape[0]
        next_vars = [Real(f'layer{i}_n{j}') for j in range(n_out)]

        for j in range(n_out):
            # 线性变换: z = Wx + b
            z = Sum([W[j, k] * current_vars[k]
                     for k in range(len(current_vars))]) + b[j]

            if activation == 'relu':
                # ReLU: y = max(0, z)
                solver.add(next_vars[j] == If(z >= 0, z, 0))
            elif activation == 'linear':
                solver.add(next_vars[j] == z)

        current_vars = next_vars

    # 输出变量等于最后一层的激活值
    for j in range(len(output_vars)):
        solver.add(output_vars[j] == current_vars[j])

    return solver

def verify_local_robustness(network, x0, epsilon, target_class):
    """验证给定点的 ε-鲁棒性"""
    solver = Solver()

    # 输入变量
    n_input = x0.shape[0]
    x_vars = [Real(f'x_{i}') for i in range(n_input)]

    # 输入约束：ε-球内
    for i in range(n_input):
        solver.add(x_vars[i] >= x0[i] - epsilon)
        solver.add(x_vars[i] <= x0[i] + epsilon)

    # 输出变量（logits）
    n_output = 10  # 假设 10 分类
    y_vars = [Real(f'y_{i}') for i in range(n_output)]

    # 编码网络
    encode_network(solver, network, x_vars, y_vars)

    # 性质规否：存在某个类 k≠target 使得 y_k > y_target
    target_var = y_vars[target_class]
    violation = Or([y_vars[k] > target_var
                    for k in range(n_output) if k != target_class])
    solver.add(violation)

    # 检查可满足性
    result = solver.check()
    if result == unsat:
        return True   # 鲁棒：不存在对抗样本
    else:
        return False  # 非鲁棒：找到了对抗样本（反例）

# Marabou 验证调用示例（Python API）
from maraboupy import Marabou

def verify_with_marabou(network_path, input_bounds, property):
    """使用 Marabou 验证网络性质"""
    network = Marabou.read_network(network_path)

    # 设置输入边界
    for i, (lo, hi) in enumerate(input_bounds):
        network.set_lower_bound(i, lo)
        network.set_upper_bound(i, hi)

    # 添加性质约束
    output_vars = [network.output_var(i) for i in range(num_classes)]
    property_constraint = encode_property(output_vars, property)

    result, stats = network.solve()

    if result is None:
        return "UNKNOWN"
    elif len(result) == 0:
        return "VERIFIED"  # 证明成立
    else:
        return "COUNTEREXAMPLE", result  # 返回反例

import torch
import torch.nn as nn

class VATLoss(nn.Module):
    """验证感知训练损失：最小化扰动下的输出变化"""
    def __init__(self, epsilon=0.01, xi=1e-6, num_iters=1):
        super().__init__()
        self.epsilon = epsilon
        self.xi = xi
        self.num_iters = num_iters

    def forward(self, model, x):
        # 随机方向初始化
        d = torch.randn_like(x)
        d = self.xi * self._normalize(d)

        # 幂法求最坏扰动方向（近似梯度方向）
        for _ in range(self.num_iters):
            d.requires_grad = True
            pred = model(x)
            pred_pert = model(x + d)
            # 最大化 KL 散度方向
            kl_div = self._kl_divergence(pred, pred_pert)
            kl_div.backward()
            d = d.grad.detach()
            d = self._normalize(d)

        # 在最坏方向施加 ε 扰动
        r_adv = self.epsilon * d
        pred = model(x)
        pred_adv = model(x + r_adv)

        return self._kl_divergence(pred, pred_adv)

    def _normalize(self, v):
        return v / (torch.norm(v, dim=(1,2,3), keepdim=True) + 1e-8)

    def _kl_divergence(self, p, q):
        return torch.sum(p * (torch.log(p + 1e-10) - torch.log(q + 1e-10)),
                         dim=1).mean()

# 集成到训练循环
vat_loss_fn = VATLoss(epsilon=0.01)

for batch in dataloader:
    images, labels = batch
    outputs = model(images)
    loss = criterion(outputs, labels)

    # VAT 作为正则项
    vat_loss = vat_loss_fn(model, images)
    total_loss = loss + 0.1 * vat_loss

    optimizer.zero_grad()
    total_loss.backward()
    optimizer.step()

训练框架 → ONNX 导出 → 验证专用 IR → 验证引擎
                            ↓
                   标准化规格语言 → 性质检查器