import numpy as np from scipy.io.wavfile import read # ========================================== # HARDWARE ABSTRACTION LAYER (SOFT NPU) # ========================================== # Simulates the systolic grid, on-chip SRAM, and dataflow waves. # In real hardware, these would be fixed-point operations in silicon gates. class NPUTensor: """Wrapper representing data in Grid SRAM.""" def __init__(self, data, name="tensor"): self.data = np.array(data, dtype=np.float32) # Simulating Fixed-Point self.name = name self.shape = self.data.shape def __repr__(self): return f"NPU-Tensor[{self.name}] {self.shape}" class SoftNPUGrid: """ The Dataflow Runtime. All operations here are 'waves' that pass through the grid without CPU intervention. """ @staticmethod def MATMUL(A, B): """General Matrix Multiply (GEMM) - Systolic Array Pattern.""" return NPUTensor(np.dot(A.data, B.data), name="gemm_out") @staticmethod def BROADCAST_ADD(A, b): """Adds bias vector b to every row of matrix A.""" # In hardware: b is broadcast from the edge of the grid. return NPUTensor(A.data + b.data, name="broadcast_out") @staticmethod def RELU(X): """Element-wise ReLU Activation.""" return NPUTensor(np.maximum(0, X.data), name="relu_out") @staticmethod def RELU_DERIV(X): """Derivative of ReLU for Backpropagation.""" return NPUTensor(np.where(X.data > 0, 1.0, 0.0), name="relu_grad") @staticmethod def SOFTMAX(X): """Softmax with numerical stability.""" # Hardware: Requires exponential function units. exp_x = np.exp(X.data - np.max(X.data, axis=1, keepdims=True)) return NPUTensor(exp_x / np.sum(exp_x, axis=1, keepdims=True), name="softmax_out") @staticmethod def MEAS(X): """ COLLAPSE operation. Reduces a probability distribution to a single scalar value (ArgMax). Equivalent to: final_class = argmax(probs) """ return np.argmax(X.data, axis=1) @staticmethod def TRANSPOSE(X): return NPUTensor(X.data.T, name="transpose") @staticmethod def REDUCE_MEAN_AXIS0(X): """Reduces a matrix to a vector by averaging columns.""" return NPUTensor(np.sum(X.data, axis=0, keepdims=True) / X.data.shape[0], name="reduce_mean") # ========================================== # DATA DIVISION # ========================================== # 1. Load Data (I/O Phase: Main RAM -> Grid Interface) # Note: Loading from WAV files is treated as a high-latency I/O burst. X_train_np = read('../X_train.wav')[1].reshape(-1, 784) y_train_np = (read('../y_train.wav')[1] * 9).astype(int) X_test_np = read('../X_test.wav')[1].reshape(-1, 784) y_test_np = (read('../y_test.wav')[1] * 9).astype(int) yt20 = y_test_np[:1000] X20 = X_test_np[:1000] # 2. Define the Probabilistic / Dataflow Neural Network class MLP_ProbolSC: """ Implements an MLP using PROBOL-SC (Scientific) semantics. Weights and Biases are Grid Registers (SRAM). """ def __init__(self, input_size, hidden_size, output_size, learning_rates): # Initialize Grid Weights (Gaussian initialization as degenerate prob tensors) self.W1 = NPUTensor(np.random.randn(input_size, hidden_size) * 0.01, name="W1") self.b1 = NPUTensor(np.zeros((1, hidden_size)), name="b1") self.W2 = NPUTensor(np.random.randn(hidden_size, output_size) * 0.01, name="W2") self.b2 = NPUTensor(np.zeros((1, output_size)), name="b2") # Learning rates as a schedule self.lr = learning_rates # Pipeline Registers (holds intermediate activation waves) self.z1 = None self.a1 = None self.z2 = None def FORWARD_WAVE(self, X_batch): """ Executes the forward pass as a series of systolic waves. Input X_batch flows from grid edge to output edge. """ # --- Layer 1 --- self.z1 = SoftNPUGrid.MATMUL(X_batch, self.W1) self.z1 = SoftNPUGrid.BROADCAST_ADD(self.z1, self.b1) self.a1 = SoftNPUGrid.RELU(self.z1) # --- Layer 2 --- self.z2 = SoftNPUGrid.MATMUL(self.a1, self.W2) self.z2 = SoftNPUGrid.BROADCAST_ADD(self.z2, self.b2) output = SoftNPUGrid.SOFTMAX(self.z2) return output def BACKWARD_AND_UPDATE_WAVE(self, X_batch, y_true_onehot): """ Executes Gradient Descent using Reverse-Mode Differentiation (Backprop). Updates weights via Grid Feedback Registers. """ m = X_batch.shape[0] # Batch size (for normalization) # --- Gradient Wave: Output Layer --- # Loss = CrossEntropy(y_true, softmax(z2)) # dL/dz2 = softmax(z2) - y_true dz2 = NPUTensor(self.z2.data - y_true_onehot.data, name="grad_z2") # dW2 = (A1.T @ dz2) / m dW2 = SoftNPUGrid.MATMUL(SoftNPUGrid.TRANSPOSE(self.a1), dz2) dW2 = NPUTensor(dW2.data * (1.0/m), name="grad_W2") # db2 = sum(dz2) / m db2 = SoftNPUGrid.REDUCE_MEAN_AXIS0(dz2) # Normalized sum # --- Gradient Wave: Hidden Layer --- # dL/da1 = dz2 @ W2.T da1 = SoftNPUGrid.MATMUL(dz2, SoftNPUGrid.TRANSPOSE(self.W2)) # dL/dz1 = dL/da1 * relu'(z1) # Element-wise multiplication wave dz1_deriv = SoftNPUGrid.RELU_DERIV(self.z1) dz1 = NPUTensor(da1.data * dz1_deriv.data, name="grad_z1") # dW1 = (X.T @ dz1) / m dW1 = SoftNPUGrid.MATMUL(SoftNPUGrid.TRANSPOSE(X_batch), dz1) dW1 = NPUTensor(dW1.data * (1.0/m), name="grad_W1") # db1 = sum(dz1) / m db1 = SoftNPUGrid.REDUCE_MEAN_AXIS0(dz1) # --- Parameter Update (Feedback to SRAM) --- # W = W - lr * dW # Note: In Soft NPU theory, this is a 'DEGENERATE DISTRIBUTION' update (deterministic) self.W2.data -= self.lr[2] * dW2.data self.b2.data -= self.lr[3] * db2.data self.W1.data -= self.lr[0] * dW1.data self.b1.data -= self.lr[1] * db1.data def TRAIN_CYCLE(self, X_np, y_np): """ One full iteration on the NPU: 1. Load Batch (I/O) 2. Forward Wave 3. Collapse/Measure (Accuracy) 4. Backward Wave + Update """ # 1. Load Data into Grid Tensors X_tensor = NPUTensor(X_np, name="X_input") y_onehot_np = np.eye(10)[y_np] y_tensor = NPUTensor(y_onehot_np, name="y_target") # 2. Forward Pass probs = self.FORWARD_WAVE(X_tensor) # 3. Measurement (MEAS) - Check Accuracy # This collapses the probability distribution to a discrete label predictions = SoftNPUGrid.MEAS(probs) accuracy = np.mean(predictions == y_np) # 4. Backward Pass & Optimization self.BACKWARD_AND_UPDATE_WAVE(X_tensor, y_tensor) return accuracy # ========================================== # PROCEDURE DIVISION # ========================================== if __name__ == "__main__": # Initialize Learning Rate Schedule (Simulating Optimizer State) learning_rates = np.random.rand(4) # Instantiate the Neural Grid f = MLP_ProbolSC(input_size=784, hidden_size=100, output_size=10, learning_rates=learning_rates) i = 0 # Main Training Loop # In theory, this is a continuous dataflow stream feeding the grid. while True: # I/O Operation: Fetch a random minibatch from main memory idx = np.random.randint(0, 60000, 1000) X_batch = X_train_np[idx] y_batch = y_train_np[idx] # Execute one Grid Cycle for _ in range(50): idx = np.random.randint(0, 60000, 100) X0 = X_train_np[idx] y0 = y_train_np[idx] current_acc = f.TRAIN_CYCLE(X_batch, y_batch) current_acc = f.TRAIN_CYCLE(X0, y0) y = f.FORWARD_WAVE(X_test_np).data.argmax(1) acc = np.mean(y==y_test_np) # Reporting (Host CPU side) print(f"Iteration {i}: Grid Accuracy = {current_acc:.4f}: Accuracy {acc}") i += 1