Pivot-safe discrete LU initialization

Integrate this at initialization of an invertible linear layer, coupling layer, reversible MLP block, or normalizing-flow transformation whose weight is represented as an LU factorization. Use matrix dimension n small enough for an exact startup check; for larger blocks, apply the procedure independently to block-diagonal or grouped channel matrices. Pseudocode: (1) choose a finite support S, for example S={-1,0,+1} or S={-a,+a}, with a<1, and a sampling distribution p over S; (2) repeat up to R times: sample A[i,j] independently from p, compute unpivoted LU(A) in float64, record each pivot U[k,k], and reject if any pivot is exactly zero or |U[k,k]| is below a numerical threshold epsilon_pivot; (3) compute g=max(max_abs(L), max_abs(U)/max_abs(A)); reject if g>g_max, such as n^(5/2+delta) or a tighter engineering limit; (4) accept A and store L and U using the network's usual parameterization, often storing log|diag(U)| separately and optionally renormalizing A by a scalar; (5) train the layer normally, but log the same growth metric every K optimizer steps. At each step, the zero-pivot condition is computed by the LU routine, while growth is estimated from the current factors; the probability lower bound and n^(5/2+delta) scaling come from the paper, but the constants c_{xi,delta}, epsilon_pivot, and useful g_max are empirical. Start with a 4-layer RealNVP or reversible MLP on CIFAR-10, using 32- or 64-dimensional channel blocks. Compare Gaussian LU initialization, standard random orthogonal initialization, and discrete-resampled LU at equal parameter count. Measure failed initializations, log-determinant variance, forward/inverse reconstruction error, gradient norm spikes, and validation likelihood. The expected signal is fewer unstable or unusable initializations, lower early-training gradient amplification, and reduced inverse reconstruction error at equal training compute.