Recurrent neural networks (RNNs) offer linear-time scaling with sequence length while requiring only constant memory, yet they struggle to capture long-range dependencies due to vanishing gradients and limited receptive fields. To address these limitations, we introduce a second-order recurrent model in which the standard neuron-to-neuron communication is replaced by a spatially evolving field governed by (discretized) partial differential equations. Drawing inspiration from the role of cortical waves in brain computation, this mechanism allows structured spatiotemporal patterns to serve as an implicit, high-capacity memory. We show that the resulting model is equivalent to a structured infinite-order RNN in which the current state depends explicitly on its entire history of past states, yielding an effectively unbounded receptive field with a fixed number of parameters. We further derive constructive conditions to ensure marginal stability, constraining the gradient spectrum on the unit circle and thereby eliminating vanishing and exploding gradients. Empirically, the proposed architecture outperforms other recurrent models on long-horizon benchmarks while using substantially fewer parameters, demonstrating that spatial dynamics can effectively bridge the gap between efficient inference and long-term memory.
Figures & tables
Figure 1: Schematic of the model architectures. ( a ) SpatialRNN: The neuron states ht (represented by white nodes) are coupled through a shared spatiotemporal communication medium ψt (color coded), where information propagates locally over space and time. Each node writes to and reads from the medium, inducing wave-like interactions that encode past activity. ( b ) HORNN: Nodes correspond to neuron states ht , and interactions are explicitly represented by edges connecting the current state to a finite set of past states ht−k , for k=1,…,d . For d=1 , the model reduces to a standard RNN. Block diagrams of the corresponding discretized computations are shown on the right of each panel.
Figure 2: Gradient analysis of SpatialRNNs and HORNNs. ( a ) Eigenvalue spectrum of the gradient dynamics for standard RNNs (top) and SpatialRNNs (bottom). Red and green markers indicate the eigenvalues computed for D=I and D∈D , respectively. The spectrum of standard RNNs is sensitive to variations in D , whereas SpatialRNNs satisfying Theorem 3.3 allow for arbitrary eigenvalue placement with invariant magnitude, including on the unit circle. ( b ) Gradient norm versus time horizon. ( c ) Gradient norm (at lag m=20 ) versus the number of parameters (with 8 neurons). The results in (b,c) are averaged over 10 trajectories with uniformly sampled inputs.
Figure 3: Performance of SpatialRNN and HORNNs on the copy task problem, shown as a function of sequence length T .
Model (# units)
Acc. (%)
Params
unitary RNN (512) ( Arjovsky et al., 2016 )
91.4 (max)
9k
LSTM (256) ( Helfrich et al., 2018 )
92.9 (max)
270k
GRU (256) ( Chang et al., 2017 )
94.1 (max)
200k
NWM (256) ( Keller and Welling, 2023 )
94.8±1.1
50k
coRNN (256) ( Keller and Welling, 2023 )
95.0±2.4
134k
SpatialRNN (88)
95.5±0.2
13k
Table 1: psMNIST performance. Mean ± std over 3 runs when available; otherwise, maximum accuracy. Only single-layer architectures are accounted for.
Model (# units)
Acc. (%)
Params
LSTM (128) ( Kag et al., 2020 )
11.6
116k
GRU (128) ( Keller et al., 2024 )
43.8
88k
iRNN (529) ( Keller et al., 2024 )
51.3
336k
wRNN (256) ( Keller et al., 2024 )
55.0
435k
SpatialRNN (88)
57.2
21k
LRNN (128) ( Keller et al., 2024 )
57.4
46k
Table 2: npCIFAR10 performance for T=1000 . Maximum test accuracy over 3 trials.
Figure 4: SpatialRNN validation accuracy on the npCIFAR10 task for varying T .
Appendix figures & tables8 assets
Supplementary material from the paper’s appendix.
Appendix
Figure 5: Evolution of the parameters {αi}i=1k during training on the copy task with sequence length T=128 . The dashed lines show the learned values of αi for a SpatialRNN architecture with block sizes [8, 8, 8, 8], while the solid line indicates training accuracy. Performance remains at chance-level until α4→1 , after which learning rapidly progresses.
Figure 6: Gradient norm ∥∂hT−t∂L∥ as a function of the time horizon, shown for the psMNIST classification task. Matrix norms are set such that ∣∣Bi∣∣2<0.4,∣∣Wψ,i∣∣2<0.4 and ∣∣Wi∣∣2<1.8 .
Figure 7: Frobenius norm of the exact state-to-state Jacobians Gmτ=∂hτ+m/∂hτ , propagated via Eqs. ( 13 )–( 14 ) along 100 psMNIST validation sequences ( τ=0 , T=784 ), for a trained SpatialRNN (left) and the same architecture at initialization (right). Thin lines show individual sequences, while the black line and the shaded band show the geometric mean and ±1 standard deviation (in log scale).
Model
Log-linear slope (per step)
Change over 783 steps
CV of ∥Gmτ∥F
Trained SpatialRNN
(7.98±7.75)⋅10−5
≈6%
0.071±0.010
Untrained SpatialRNN
(2.60±1.51)⋅10−4
≈23%
0.072±0.032
Appendix
Table 3: Empirical LTV gradient statistics over 100 psMNIST validation sequences (mean ± standard deviation across sequences). Slopes are fitted for m≥50 , and the coefficient of variation (CV) measures the dispersion of ∥Gmτ∥F along each sequence.
Figure 8: Relative Frobenius error ∥Gmexact−Gmrec∥F/∥Gmexact∥F between the exact Jacobians Gm=∂hT−1/∂hT−1−m , computed from the stored forward states, and those computed from the backward-reconstructed states, as a function of the time lag m , in single (float32) and double (float64) precision. (a) Non-dissipative SpatialRNN (single block, n=40 , α1=1 ): the error remains at machine precision over the entire horizon, with no exponential accumulation. (b) SpatialRNN with dissipative blocks (block sizes {22,6,6,6} , α=[1,0.5,0.5,0.5] ): round-off is amplified exponentially, at a measured rate of 1.42 per backward step, matching the predicted αmin−1/2=1.414 . Vertical lines mark the lags at which the median error first exceeds 10−3 , yielding m=44 in single and m=101 in double precision. Since the reconstruction starts from (ψT−1,ψT−2) and proceeds backward, the state at time t is recovered after T−1−t backward steps, so m is also the depth of the backward recursion. Curves show the median over 20 input sequences of length T=784 and shaded bands the 10 – 90 percentile range; the same effective matrices are used in both precisions, so the two curves differ only in the arithmetic precision of the recursions. Both models are randomly initialized, with the diagonal blocks rescaled so that condition 2 of Theorem 3.3 holds.
Model
Neurons
Parameters
Mean Accuracy(%) ± std.
SpatialRNN [32, 8, 8, 8]
56
6,815
99.04 ± 0.07
SpatialRNN [8, 8, 8, 8]
32
2,915
96.90 ± 4.34
SpatialRNN [4, 4, 4, 4]
16
1,179
69.33 ± 16.54
HORNN (T=[1,5,10,15])
64
17,744
13.39 ± 3.07
HORNN (T=[1,5,10])
64
13,648
14.37 ± 2.72
HORNN (T=[1,5])
64
9,552
11.05 ± 1.20
Appendix
Table 4: Copy task architecture configurations and performance. Mean accuracy is evaluated at the maximum sequence length T=1024 . Average and standard deviation over 10 runs.
Model
Neurons
Parameters
Mean Accuracy
SpatialRNN [128, 8, 8, 8]
152
36,497
96.34±0.09
SpatialRNN [64, 8, 8, 8]
88
12,849
95.49±0.22
SpatialRNN [32, 8, 8, 8]
56
5,633
94.37±0.13
Appendix
Table 5: psMNIST task architecture configurations and performance. Mean accuracy is evaluated over 3 runs, ordered by mean accuracy.
Model
Neurons
Parameters
Max Accuracy
SpatialRNN [128, 8, 8, 8]
152
50,937
57.8
SpatialRNN [64, 8, 8, 8]
88
21,209
57.2
SpatialRNN [32, 8, 8, 8]
56
10,953
54.1
SpatialRNN [16, 8, 8, 8]
40
6,977
50.8
SpatialRNN [8, 8, 8, 8]
32
5,277
45.6
Appendix
Table 6: npCIFAR10 task architecture configurations and performance. Maximum accuracy is evaluated at the sequence length T=1000 over 3 runs.
Training recurrent neural networks (RNNs) requires assigning credit across long sequences of computations. Standard backpropagation through time (BPTT) addresses this problem poorly: it is sequential in time, limiting parallelism, and suffers from vanishing or exploding gradients, making long-range associations difficult to learn. We propose Supervised Memory Training (SMT), a method for training nonlinear RNNs that sidesteps recurrent credit propagation entirely by reducing RNN training to supervised learning on one-step memory transition labels (mt,xt+1)→mt+1. SMT acquires these memory labels by training a Transformer-based encoder on a predictive state objective--retaining only information from the past necessary to predict the future. By decoupling what to remember from how to update memory, SMT enables time-parallel RNN training with a stable O(1) length gradient path between any two tokens--without ever unrolling the RNN. We find that SMT outperforms BPTT when pretraining various RNN architectures on tasks like language modeling and pixel sequence modeling. SMT enables nonlinear RNNs to better capture long-range dependencies and train in parallel, potentially unlocking the scaling of models that build temporal abstractions of past experience.
Recurrent models provide a natural path to long-context modeling, yet models trained with backpropagation through time (BPTT) often fail beyond their training horizon. Classical analyses emphasize gradients that vanish or explode along temporal paths. However, dense per-token losses can still train a shared recurrent rule despite severe decay, showing that decay alone does not determine whether learning fails. We instead study state credit: the signal through which future losses reach earlier recurrent states before contributing to parameter updates. Accordingly, we intervene directly on state credit and propose Credit Stabilization through Time (CST). During backward propagation, CST locally rescales the state-credit signal to stabilize its norm without rotating the component being corrected, while leaving the forward computation unchanged. Because controlled synthetic tasks and real data exhibit different credit dynamics, we specialize CST to each regime. In both settings, CST improves performance beyond the training horizon, with gains observed at up to 128x the training length.
For over a decade, explicit memory architectures like the Neural Turing Machine have remained theoretically appealing yet practically intractable for language modeling due to catastrophic gradient instability during Backpropagation Through Time. In this work, we break this stalemate with \textit{Phasor Memory Network} (PMNet), a novel architecture that structurally resolves memory volatility through \textit{Unitary Phasor Dynamics} and \textit{Hierarchical Learnable Anchors}. Rather than relying on brute-force scaling, we present a mechanistic proof-of-concept in a controlled byte-level setting. By constraining recurrent state updates to phase rotations on a complex unit circle, PMNet preserves gradient norms and inherently prevents divergence without the need for specialized initialization. We empirically demonstrate the active actuation of the memory module through a synthetic Copy-Paste task, where PMNet utilizes an expansive \textit{85-slot hierarchical memory tree} (=∑h=144h−1) to achieve near 100% exact retrieval across temporal distances that completely exceed the local sliding window attention's receptive field. Furthermore, despite being a compact 119M parameter model trained on 18.8B tokens, PMNet matches the zero-shot long-context robustness of a Mamba model that is three times larger. Our ablation studies and gradient analyses confirm that the historical failure of explicit memory was a structural alignment problem, which PMNet effectively overcomes, providing a theoretically grounded foundation for scalable sequence modeling.
Sungwoo Goo, Hwi-yeol Yun, Sangkeun Jung
College of Pharmacy Chungnam National University Daejeon, South Korea · Department of Computer Science & Engineering Chungnam National University Daejeon, South Korea