Proceedings of Machine Learning Research vol :1–35, 2026
39th Annual Conference on Learning Theory
Deep Q-Learning on Hölder Spaces Qian Qi
QIQIAN @ PKU . EDU . CN
Peking University
Editors: Steve Hanneke and Tor Lattimore
arXiv:2606.16846v1 [cs.LG] 15 Jun 2026
Abstract We study the operator-theoretic core of Q-learning in continuous-time stochastic control with continuous states and actions. In value-based reinforcement learning, each Q-learning or DQN update is built from a Bellman optimality target; our analysis isolates this target in a diffusion setting and studies its regularity and approximation complexity. Under uniform ellipticity and Hölder-regular coefficients, we show that a Bellman update maps bounded inputs into an anisotropic regularity class, smoothing the state variable while leaving only Lipschitz dependence on the action variable. This yields a compact family of Bellman iterates and motivates a tensor-product DeepONet architecture adapted to the mixed regularity of the problem. We then derive explicit approximation and resource bounds, together with a stiffness–complexity trade-off as the time step δ → 0. The resulting theory makes a direct contribution to Q-learning theory at the level of Bellman target regularity and approximation in continuous stochastic control. At the same time, we do not claim a full convergence theorem for practical sampled Q-learning with exploration, replay, and stochastic gradient updates. Keywords: Q-Learning Theory, Bellman Value Iteration, Neural Operators, Parabolic Smoothing, Stochastic Optimal Control, Approximation Theory.
1. Introduction This paper studies the operator-theoretic core of Q-learning in continuous-time stochastic control with continuous states and actions. In value-based reinforcement learning, the central object is the Bellman optimality update for the action-value function; every Q-learning or DQN step is, at heart, an attempt to approximate this update from data. Our focus is on this Bellman target itself in a continuous-time diffusion model, where the coefficients are known and the operator can be analyzed directly. This does not yet cover all ingredients of practical Q-learning, but it does address what we view as one of its core theoretical questions: what regularity and approximation complexity does the Bellman target class possess? The relationship with Q-learning can be stated concretely. In discrete time, a typical valuebased target has the form r + β maxa′ Q(x′ , a′ ), namely a stochastic approximation to a Bellman optimality operator. Our continuous-time discretization replaces the single sampled transition by its conditional expectation over a short diffusion step of size δ, thereby isolating the operator-level part of the theory from the separate issues of exploration, replay, and stochastic optimization. This is why we view the paper as a contribution to Q-learning theory, while also being explicit that it operates at the Bellman-operator level. Classical analyses of value-based methods with function approximation often assume that the transition and reward coefficients are globally Lipschitz in state and action. While convenient, this hypothesis misses a central feature of stochastic control: non-degenerate diffusion can regularize *. Accepted for presentation at the Conference on Learning Theory (COLT) 2026. © 2026 Q. Qi.
QI
the Bellman update through parabolic smoothing. In PDE language, the Bellman step is tied to a controlled parabolic equation, so rough inputs may become more regular after one application of the operator (Krylov, 1987). Motivated by this observation, we study the time-discretized Bellman optimality operator Tδ under uniform ellipticity and Hölder-regular coefficients. The main message is that the Bellman targets underlying value-based learning in this regime lie in an anisotropic regularity class substantially smaller than a generic Lipschitz class. This creates an approximation-theoretic route toward neural architectures tailored to the regularity of continuous-control Bellman iterates, and hence toward a sharper theory of what Q-learning-style methods are trying to approximate in this setting. 1.1. Main Contributions Our contributions are fourfold. First, we prove that under uniform ellipticity the Bellman operator maps bounded inputs into an anisotropic Hölder class: smooth in the state variable and only Lipschitz in the action variable. This identifies a nontrivial regularity structure for the Bellman targets that arise in continuous-time Q-learning problems. Second, we show that this regularity statement is already a core theoretical contribution to Qlearning: it describes the target class that any value-function approximator must represent in this diffusion regime. The operator-theoretic viewpoint is therefore not peripheral to Q-learning, but rather isolates its mathematical backbone. Third, we use this regularity to derive approximation guarantees for a tensor-product DeepONet architecture that approximates Bellman iterates. The architecture is chosen because it matches the anisotropic regularity revealed by the Bellman analysis; the regularity statement itself is architectureindependent, while the theorem provides one constructive approximation route. Fourth, we quantify the stiffness–complexity trade-off induced by the time step δ. As δ → 0, the Bellman operator approaches the identity and the smoothing effect weakens, which increases the approximation complexity. This explains why continuous-time limits can become numerically stiff even when the diffusion is regularizing. Overall, the paper aims to make a precise claim: in continuous stochastic control with uniformly elliptic diffusion, the Bellman targets at the heart of Q-learning enjoy exploitable anisotropic regularity, and this has concrete consequences for approximation theory and architecture design. What we do not claim is a full convergence theorem for sampled Q-learning with exploration, replay, target-network heuristics, and stochastic gradient optimization. 1.2. Literature Review Our work sits at the intersection of approximation theory for Deep Learning, stochastic analysis of PDEs, and reinforcement learning. We briefly survey the relevant literature in these domains. The convergence of Q-learning with function approximation has been a central topic since the early work of Tsitsiklis and Van Roy (1996). The empirical success of DQN was popularized by Mnih et al. (2015), while rigorous analyses of value-based methods with function approximation often impose strong smoothness or contraction assumptions; see, for example, Fan et al. (2020). Our work addresses a complementary but, in our view, foundational theoretical question: what regularity and approximation complexity do the Bellman targets themselves possess in a continuous-time stochastic-control regime? Because Q-learning and DQN are built around repeated approximation 2
D EEP Q-L EARNING ON H OLDER S PACES
of Bellman optimality targets, understanding this target class is already a direct contribution to Qlearning theory, even before one adds finite-sample and optimization analysis. The observations of Tallec et al. (2019) about stiffness as δ → 0, of Qi (2025) about regularity collapse near decision boundaries, and of Qi (2026) about viscosity-consistency constraints for deep Q-network discretizations help motivate our complexity discussion. At the same time, to avoid overstatement, we emphasize that our theorem is an operator-level approximation result for Bellman iterates in a diffusion setting, not a full convergence guarantee for practical sampled Q-learning training dynamics. The correspondence between the Bellman equation and parabolic PDEs is classical, established by Krylov (1980) and Fleming and Soner (2006). While the existence of viscosity solutions for HJB equations is well-known (Crandall et al., 1992), the application of interior regularity estimates to analyze the sample complexity of RL is relatively new. Most existing finite-sample bounds for RL (Azar et al., 2017; Jin et al., 2018) operate in the tabular or linear setting. In contrast, our work explicitly utilizes the anisotropic nature of the Schauder estimates (C 2,α in state, C 0,1 in action) induced by the diffusion term Tr(σσ ⊤ D2 V ), a concept deeply rooted in the theory of regularization by noise. Universal Approximation Theorems (UATs) for neural networks have evolved from qualitative density results (Cybenko, 1989; Hornik, 1991) to quantitative complexity bounds. Yarotsky (2017) established the relationship between network depth/width and the Sobolev regularity of the target function, showing that W ∼ ϵ−d/s . However, applying these isotropic bounds to Q∗ yields pessimistic estimates due to the joint dimension dstate + daction . Our analysis circumvents this by adopting the anisotropic viewpoint of Mhaskar and Poggio (2016) and Poggio et al. (2017), decoupling the complexity of the smooth state manifold from the non-smooth policy landscape. The architecture proposed in this paper draws upon the emerging field of Operator Learning. The Deep Operator Network (DeepONet), introduced by Lu et al. (2021), generalizes the universal approximation theorem to non-linear operators on Banach spaces. Similarly, the Fourier Neural Operator (FNO) (Li et al., 2021) has shown success in solving PDEs by learning the solution operator directly. While these methods are primarily applied to forward problems in fluid dynamics and solid mechanics, we adapt the framework to the backward Bellman flow. Our Tensor-Product Trunk modification addresses the specific mixed-regularity requirements of the Control problem, which are not typically present in physical conservation laws.
2. Preliminaries and Problem Formulation 2.1. The Stochastic Control Problem Let (Ω, F , F = {Fs }s≥0 , P) be a complete filtered probability space satisfying the usual conditions, supporting a standard d-dimensional Brownian motion W = (Ws )s≥0 . We consider a finite-horizon Markov Decision Process over the time interval t ∈ [0, T ]. The state process Xs ∈ Rd is controlled by an action process a = (as )s∈[t,T ] taking values in a compact metric space A ⊂ Rda , where da denotes the action dimension. We denote the set of admissible control policies by U [t, T ], defined as the space of all F-progressively measurable processes taking values in A such that the associated SDE has a unique strong solution. 3
QI
For a fixed initial time t ∈ [0, T ) and state x ∈ Rd , the dynamics of the system are governed by the forward Itô Stochastic Differential Equation (SDE): ( dXst,x,a = b(s, Xst,x,a , as ) ds + σ(s, Xst,x,a , as ) dWs , s ∈ [t, T ], (1) Xtt,x,a = x. The objective of the agent is to maximize the expected cumulative reward defined by the cost functional: Z T
J(t, x, a) = E t
r(s, Xst,x,a , as ) ds + g(XTt,x,a ) Ft .
(2)
We assume the drift b : [0, T ] × Rd × A → Rd , diffusion σ : [0, T ] × Rd × A → Rd×d , running reward r : [0, T ] × Rd × A → R, and terminal reward g : Rd → R satisfy standard measurability and boundedness conditions (detailed in Assumption 1). The optimal value function V ∗ (t, x) is defined via the maximization over admissible controls: V ∗ (t, x) =
sup J(t, x, a).
(3)
a∈U [t,T ]
It is well-known from classical stochastic control theory (Krylov, 1980) that, under suitable regularity conditions, V ∗ is the unique viscosity solution to the backward Hamilton-Jacobi-Bellman (HJB) equation: ( −∂t V ∗ (t, x) − supa∈A {La V ∗ (t, x) + r(t, x, a)} = 0, (t, x) ∈ [0, T ) × Rd , (4) V ∗ (T, x) = g(x), where La is the second-order linear differential operator (generator) associated with a constant action a: 1 (5) La ϕ(x) := Tr σ(t, x, a)σ(t, x, a)⊤ Dx2 ϕ(x) + b(t, x, a)⊤ ∇x ϕ(x). 2 In the present idealized Bellman-iteration model, we do not analyze data-driven Q-learning or exploration directly. Instead, we approximate the action-value function Q∗ (t, x, a) associated with the semi-discrete Bellman operator obtained by holding an action constant for a time step δ and acting optimally thereafter. 2.2. Assumptions: Ellipticity and Regularity To invoke the Schauder interior regularity estimates, we require specific structural assumptions on the MDP. A critical distinction in our analysis is the separation between regularity in the statetime domain (t, x), where diffusion generates smoothing, and the action space A, where no such smoothing occurs. We first define the parabolic metric on [0, T ] × Rd : dP ((t, x), (s, y)) := |t − s|1/2 + ∥x − y∥.
(6)
A function f : [0, T ] × Rd → R is said to be parabolic Hölder continuous with exponent α ∈ (0, 1), denoted f ∈ C α/2,α ([0, T ] × Rd ), if: ∥f ∥C α/2,α := sup |f (t, x)| + (t,x)
|f (t, x) − f (s, y)| < ∞. α (t,x)̸=(s,y) dP ((t, x), (s, y)) sup
This scaling reflects the natural heat-kernel behavior where time scales as the square of space. 4
(7)
D EEP Q-L EARNING ON H OLDER S PACES
Assumption 1 (Global Coefficients) The drift b, diffusion σ, and reward r satisfy the following conditions on [0, T ] × Rd × A: (i) Global Boundedness and Hölder Continuity: The functions b, σ, r are bounded and uniformly Hölder continuous with exponent α ∈ (0, 1) in the parabolic metric on the entire domain DT = [0, T ] × Rd . They are also globally Lipschitz in the state variable, uniformly over time and actions, so that the controlled SDE and the flow-sensitivity estimates used below are well posed. Specifically, there exists a constant K such that for all a ∈ A: ∥b(·, ·, a)∥C α/2,α (DT ) + ∥σ(·, ·, a)∥C α/2,α (DT ) + ∥r(·, ·, a)∥C α/2,α (DT ) ≤ K.
(8)
The terminal payoff g is bounded and globally Lipschitz in x. (ii) Uniform Ellipticity on Rd : There exists a constant ν > 0 such that for all (t, x, a) ∈ [0, T ] × Rd × A and ξ ∈ Rd : ξ ⊤ (σσ ⊤ )(t, x, a)ξ ≥ ν∥ξ∥2 . (9) (iii) Lipschitz Continuity in Actions: The coefficients are globally Lipschitz continuous with respect to the action parameter a, uniformly in (t, x). Remark 2 (Domain Definition) By assuming the coefficients are well-behaved on Rd , we treat the Bellman update as a Cauchy problem. This avoids the technical complications of boundary conditions (e.g., reflecting or absorbing boundaries) on a compact domain Ω, where interior estimates would diverge as dist(x, ∂Ω) → 0. The compactness required for the approximation theorem will be recovered via the restriction of the global solution to compact reachable sets. Assumption 3 (Uniform Ellipticity) The diffusion is uniformly non-degenerate; this repeats the ellipticity condition separately because it is the assumption invoked directly by the PDE estimates. Let Σ(t, x, a) := σ(t, x, a)σ(t, x, a)⊤ . There exists a constant ν > 0 such that for all (t, x, a) ∈ [0, T ] × Rd × A and all ξ ∈ Rd : ξ ⊤ Σ(t, x, a)ξ ≥ ν∥ξ∥2 .
(10)
Remark 4 (Necessity of Assumptions) Assumption 1(i) ensures that the coefficients possess sufficient regularity to support classical solutions to the parabolic PDE. If the coefficients were merely bounded measurable, solutions would only exist in Sobolev spaces Wp1,2 (Krylov, 1980), which is insufficient for the pointwise gradient estimates required for our operator network bounds. Assumption 3 prevents the diffusions from becoming singular, ensuring the existence of a smooth Gaussian transition density (heat kernel) which drives the regularization mechanism. Example 1 (A controlled diffusion satisfying the assumptions) A concrete example is a bounded stochastic tracking problem with compact action set A ⊂ Rda : dXs = tanh(AXs ) + Bas + β(s) ds + Σ0 dWs , where A and B are fixed matrices, β is bounded and smooth, and Σ0 Σ⊤ 0 ≥ νId for some ν > 0. Let the running reward and terminal payoff be bounded Hölder-continuous functions, for instance r(t, x, a) = r̄ t, tanh(x), a , g(x) = ḡ(tanh(x)), 5
QI
with r̄ and ḡ smooth and bounded and with A compact. Then b, σ, and r satisfy Assumption 1, Assumption 3 holds by construction, and the model retains the key features relevant to our setting: continuous states, continuous actions, and genuinely stochastic transitions. The example should be read as a bounded-coefficient proxy for continuous-control problems in which diffusion regularizes the Bellman update. 2.3. Anisotropic Function Spaces Standard isotropic function spaces are ill-suited for the analysis of the Bellman operator in stochastic control. The diffusion term σσ ⊤ provides parabolic smoothing in the state-time domain (t, x), but the dependence on the action a is merely parametric. The optimal policy, and consequently the value function, may exhibit non-smooth behavior with respect to a (e.g., kinks at points where the optimal action switches). To rigorously capture this mixed regularity, we introduce anisotropic Hölder spaces. Let DT = [0, T ]×Rd . We denote by C 1,2,α (DT ) the classical parabolic Hölder space consisting of functions u(t, x) such that ∂t u, ∇x u, and Dx2 u are continuous and bounded, and the highestorder derivatives satisfy the α-Hölder condition with respect to the parabolic metric dP . The norm is defined as: ∥u∥C 1,2,α (DT ) := ∥u∥∞ + ∥∂t u∥C α/2,α + ∥∇x u∥C (1+α)/2,1+α + ∥Dx2 u∥C α/2,α .
(11)
Here, the Hölder seminorms are defined with respect to dP ((t, x), (s, y)) = |t−s|1/2 +∥x−y∥. This definition ensures that if u ∈ C 1,2,α , it satisfies the requisite regularity for the classical Schauder estimates applied to ∂t u + La u = f . Definition 5 (Anisotropic Space Hα,Lip ) A function Q : [0, T ] × Rd × A → R belongs to the anisotropic space Hα,Lip if: (i) Smoothing Domain: For every fixed action a ∈ A, the map (t, x) 7→ Q(t, x, a) belongs to the parabolic space C 1,2,α (DT ). (ii) Parametric Domain: For every fixed (t, x) ∈ DT , the map a 7→ Q(t, x, a) is Lipschitz continuous (i.e., belongs to C 0,1 (A)). The associated norm is: ∥Q∥Hα,Lip := sup ∥Q(·, ·, a)∥C 1,2,α (DT ) + a∈A
|Q(t, x, a) − Q(t, x, a′ )| . ∥a − a′ ∥A (t,x)∈DT a̸=a′ sup
sup
(12)
This construction is critical for our main result. While the Bellman operator lifts the regularity of the value function in state space (from L∞ to C 2,α ), it preserves (but does not improve) the regularity with respect to actions. Since the coefficients b, σ, r are assumed to be Lipschitz in a (Assumption 1), the solution Q naturally resides in Hα,Lip . 2.4. The Bellman Semi-Group The central object of our analysis is the evolution of the action-value function under the optimal control dynamics. Unlike the infinite-horizon setting which relies on contraction properties of a static operator, the finite-horizon setting is naturally described by a backward parabolic evolution. 6
D EEP Q-L EARNING ON H OLDER S PACES
Let Lat denote the infinitesimal generator associated with a fixed action a ∈ A at time t, acting on functions ϕ ∈ C 2 (Rd ): 1 Lat ϕ(x) := Tr σ(t, x, a)σ(t, x, a)⊤ Dx2 ϕ(x) + b(t, x, a)⊤ ∇x ϕ(x). (13) 2 We define the discrete-time backward Bellman operator Tδ which maps a terminal value function approximation at time t + δ to a value function at time t. Let Φ : Rd × A → R be a bounded measurable function representing the Q-function at the future time step t + δ. The update consists of two stages: (i) Non-linear Maximization (Optimization Step): We first compute the optimal value at the future state: VΦ (y) := sup Φ(y, a′ ). (14) a′ ∈A
Note that even if Φ is smooth, VΦ is typically only Lipschitz continuous due to the supremum operation. (ii) Linear Evolution (Smoothing Step): We solve the backward Kolmogorov equation (FeynmanKac formula) over the interval [t, t + δ). For a fixed current action a, let u(·, ·; a) be the solution to the linear parabolic PDE: ( ∂s u(s, x) + Las u(s, x) + r(s, x, a) = 0, s ∈ [t, t + δ), (15) u(t + δ, x) = VΦ (x). The Bellman Operator Tt,δ is defined as the evaluation of this solution at the initial time s = t: (Tt,δ Φ)(x, a) := u(t, x; a).
(16)
Remark 6 (Probabilistic Interpretation) Under the uniform ellipticity condition (Assumption 3), there exists a smooth transition density function (fundamental solution) Γa (s, x; τ, y) associated with the operator ∂s + Las . The operator admits the integral representation: Z Z t+δ Z (Tt,δ Φ)(x, a) = Γa (t, x; t + δ, y)VΦ (y) dy + Γa (t, x; s, y)r(s, y, a) dy ds. (17) Rd
t
Rd
The first term represents the diffusion of the terminal value (convolution with a heat kernel), while the second term accumulates the running reward. It is the smoothing property of the kernel Γa , specifically its Gaussian decay and singularity behavior, that allows us to recover C 2,α regularity for u(t, ·) even when the terminal data VΦ is non-differentiable.
3. Regularity of the Bellman Operator In this section, we establish the smoothing properties of the Bellman operator. Unlike standard RL theory which relies on the contraction property to preserve Lipschitz continuity, we utilize the regularizing effect of the diffusion term σσ ⊤ . We analyze the single-step backward operator Tδ defined in (16). This operator can be decomposed into two distinct stages: 7
QI
(i) Non-linear Maximization: V (x) = supa′ ∈A Φ(x, a′ ). This step generally destroys differentiability, even if Φ is smooth. (ii) Linear Smoothing (Convolution): The evolution of V (x) backwards in time via the parabolic PDE. This step restores differentiability with respect to x. 3.1. Interior Schauder Estimates We analyze the regularity of the single-step Bellman update by casting it as a solution to a linear parabolic Cauchy problem. A critical improvement over standard analysis arises from the observation that while the maximization step destroys differentiability (preventing C 1 propagation), it preserves Lipschitz continuity. We leverage this to derive tighter smoothing estimates than those available for generic bounded measurable initial data. Let a ∈ A be a fixed action. Let Φ be a bounded Q-function approximation at the next time step, uniformly Lipschitz in the state variable as in Theorem 7. The update involves first computing the optimal value at the future state: VΦ (x) := sup Φ(x, a′ ).
(18)
a′ ∈A
We emphasize that due to the supremum operator, VΦ typically lacks differentiability (it is generally only Lipschitz continuous even if Φ is smooth). This roughness prevents the existence of classical solutions at τ = 0. However, the diffusion term induces immediate smoothing for τ > 0. Let v(τ, x) := (Tt+δ−τ,τ Φ)(x, a). Then v satisfies the forward parabolic PDE with rough initial data: ( ∂τ v(τ, x) = La v(τ, x) + r(t + δ − τ, x, a) in (0, δ] × Rd , (19) v(0, x) = VΦ (x). Here, the operator La acts on the spatial variable x. Under Assumption 3 (Uniform Ellipticity) and Assumption 1 (Hölder continuous coefficients), the operator ∂τ −La is strictly parabolic. We invoke interior Schauder estimates to bound the solution at time τ = δ. Theorem 7 (Refined Interior Parabolic Smoothing) Grant Assumptions 1 and 3. Let the input function Φ be uniformly Lipschitz continuous in the state variable, i.e., supa′ ∈A ∥Φ(·, a′ )∥C 0,1 (Rd ) ≤ LΦ . Then, for any fixed δ > 0, the function x 7→ (Tδ Φ)(t, x, a) belongs to the Hölder space C 2,α (Rd ). Moreover, there exists a constant C > 0 (depending only on d, ν, and coefficient bounds K) such that the following tighter smoothing estimate holds: 1+α ∥(Tδ Φ)(t, ·, a)∥C 2,α (Rd ) ≤ C ∥r∥C α/2,α + δ − 2 sup ∥Φ(·, a′ )∥C 0,1 (Rd ) . (20) a′ ∈A
1+α
2+α
Critically, the singularity as δ → 0 scales as O(δ − 2 ), which is strictly milder than the O(δ − 2 ) scaling associated with discontinuous initial data. Proof See Appendix B.1.
8
D EEP Q-L EARNING ON H OLDER S PACES
Remark 8 (Improved Stiffness Scaling) The refined estimate in Theorem 7 has significant impli1+α cations for the complexity analysis. The radius of the reachable function set scales as Rδ ∼ δ − 2 . Accordingly, throughout the remainder of the paper we use the stiffness scaling factor γ(α) =
3+α . 2(2 + α)
(21)
This is the quantity appearing in Theorem 15 and in the entropy bounds of the appendix. For α ≈ 1, it reduces the exponent of the complexity explosion by approximately 17%, reflecting the benefit of propagating Lipschitz continuity through the max-operator. Remark 9 (Regularity Bootstrap) Crucially, Theorem 7 guarantees that even if the input Φ is initialized as a rough approximation (e.g., from a neural network with random weights), the output after a single Bellman step lies in C 2,α (Rd ). This healing property ensures that the sequence of value iterates remains within a compact subset of the function space (for fixed δ), satisfying the prerequisites for the Operator Approximation Theorem. 3.2. Domain Truncation and Concentration of Measure A central technical challenge in applying global PDE estimates to Deep Reinforcement Learning is the bounded domain of the training data. While the theoretical diffusion process defined in (1) is supported on the entire space Rd , practical neural networks are trained on compact sets. Applying the Bellman operator restricted to a compact domain introduces a boundary pollution error, as the true value Q(x, a) depends on the probability mass that diffuses outside the domain. To resolve the inconsistency between the global Cauchy problem (analyzed in Section 3.1) and the compact function support required for complexity analysis, we invoke the concentration of measure phenomenon associated with uniform ellipticity. We define a computational domain strictly larger than the target inference domain and prove that the error induced by truncating the integral outside this √ inflated domain decays exponentially with the margin width relative to the diffusion horizon δ. Definition 10 (Nested Domains) Let S ⊂ Rd be the compact Target Domain where we require an accurate approximation of Q∗ (e.g., the region of interest for the policy). For a margin parameter ρ > 0, we define the Computational Domain Sρ as the closed ρ-neighborhood of S: d Sρ := y ∈ R : inf ∥y − x∥2 ≤ ρ . (22) x∈S
We assume the sensor set Zm used for training covers the inflated domain Sρ , while the error analysis is restricted to the interior S. We define the Truncated Bellman Operator Tδρ , which integrates the target value only over trajectories that remain within the computational domain. Let τρ := inf{s > t : Xs ∈ / Sρ } be the first exit time. The truncated operator is defined implicitly by evaluating the conditional expectation on the event {τρ > t + δ}. The following lemma quantifies the probability of escaping the computational domain within one time-step δ. 9
QI
Lemma 11 (Sub-Gaussian Exit Probability) Under Assumption 3 (Uniform Ellipticity) and the boundedness of the drift b and diffusion σ, the process Xs exhibits sub-Gaussian concentration. Let Mb := sup(t,x,a) ∥b(t, x, a)∥2 and Λ := sup(t,x,a) ∥σ(t, x, a)∥F . For any x ∈ S, fixed δ > 0, and margin ρ > Mb δ, the probability that the diffusion exits Sρ satisfies: (ρ − Mb δ)2 P (τρ ≤ t + δ | Xt = x) ≤ 2d · exp − . (23) 2dΛ2 δ Proof Consider the integral form of the dynamics: Z t+δ Z t+δ Xt+δ − x = b(s, Xs , a) ds + σ(s, Xs , a) dWs . t
t
By the triangle inequality, the exit event implies a large deviation of the running martingale supremum: ( ) Z t+s
{τρ ≤ t + δ} ⊆
σ(u, Xu , a) dWu ≥ ρ − Mb δ
sup 0≤s≤δ
.
t
R t+s Let Ms = t σ(u, Xu , a) dWu . This is a continuous vector-valued local martingale with coordinate quadratic variations bounded by Λ2 δ. If sups ∥Ms ∥2 ≥ λ, then some coordinate satisfies √ (k) sups |Ms | ≥ λ/ d. Applying the maximal Bernstein inequality for continuous martingales to each coordinate (see, e.g., Revuz and Yor (2013), Ch. IV) and taking a union bound gives ! λ2 P sup ∥Ms ∥2 ≥ λ ≤ 2d exp − . 2dΛ2 δ 0≤s≤δ Setting λ = ρ − Mb δ yields the claimed bound. Theorem 12 (Truncation Error and Interior Validity) Let Q ∈ L∞ (Rd × A) be bounded, and set Mval := ∥VΦ ∥∞ + δ∥r∥∞ for the one-step payoff bound. Let Tδ be the true global operator and Tδρ be the operator truncated to Sρ . To ensure the boundary pollution error satisfies ∥Tδ Q − Tδρ Q∥L∞ (S) ≤ εtrunc , it suffices to choose the computational margin ρ scaling as: s 2dMval ρ(δ, εtrunc ) ≥ Mb δ + Λ 2dδ log . (24) εtrunc Critically, with this margin, the Interior Schauder estimates on S are preserved. The global smoothing estimates from Theorem 7 apply to Tδ Q|S with negligible perturbation from the boundary ∂Sρ . Proof See Appendix B.2. Remark 13 (Validity of Interior Estimates) Standard interior Schauder estimates on a domain S depend on the distance to the boundary √ of the domain where the PDE is solved, scaling roughly as −k dist(S, ∂Sρ ) . By choosing ρ ∼ δ, we ensure that the boundary of the computational domain is located several standard deviations away from the target domain. Consequently, the transition density Γ(t, x; t + δ, y) is effectively zero for x ∈ S and y ∈ ∂Sρ . This justifies using the global regularity constant K (from Theorem 7) for functions restricted to S, as the boundary layer effects are suppressed by the Gaussian tail decay. 10
D EEP Q-L EARNING ON H OLDER S PACES
3.3. Compactness of Q-Iterates We now consider the sequence of Bellman iterates defined by backward induction on a time grid tk = T − kδ for k = 0, . . . , N . Let Q(0) (x, a) = g(x). The approximate value functions are generated by: Q(k+1) = Tδ Q(k) . (25) A naive analysis might suggest that repeated applications of the Bellman operator degrade regularity due to the non-smooth maximization step supa′ Q(k) (x, a′ ). However, the parabolic smoothing established in Theorem 7 ensures that regularity is not only preserved but reset at each step. Corollary 14 (Anisotropic Compactness of Iterates) Let S ⊂ Rd be the compact state domain. Under Assumptions 1 and 3, for a fixed time-step δ > 0 and finite horizon T = N δ: (k) approximates (i) Uniform Spatial Regularity: The sequence of iterates {Q(k) }N k=0 (where Q Q∗ (tk , ·, ·)) is uniformly bounded in the anisotropic Hölder space on the state domain. Specifically, there exists Cδ such that for all k:
sup ∥Q(k) (·, a)∥C 2,α (S) + sup ∥Q(k) (x, ·)∥C 0,1 (A) ≤ Cδ .
a∈A
(26)
x∈S
(ii) Sequential Compactness: The finite family Q∗δ = {Q(k) (·, ·)}0≤k≤N lies in a pre-compact subset of the Banach space C 1 (S) × C 0 (A). Equivalently, any sequence drawn from the corresponding uniformly bounded family admits a limit point in the topology of uniform convergence of values and state-gradients. Proof See Appendix B.3.
3.4. Implication for Neural Approximation The compactness result (Corollary 14) provides the rigorous foundation for approximating Q∗ using deep neural networks. Standard universal approximation theorems often assume the target function lies in a generic Sobolev space W k,p , leading to pessimistic bounds on the required network size due to the Curse of Dimensionality. Our PDE analysis offers a more precise characterization via the concept of Metric Entropy. M ETRIC E NTROPY AND C OVERING N UMBERS Let Q∗δ ⊂ Hα,Lip be the set of reachable Q-functions for a fixed time-step δ. The complexity of learning this class is governed by its Kolmogorov ϵ-entropy, log N (ϵ, Q∗δ ), where N is the covering number (the minimal number of balls of radius ϵ required to cover the set). For the Hölder ball of radius R in dimension d with smoothness s, the entropy scales as: d/s R . (27) log N (ϵ, BallR (C s )) ≍ ϵ In our context, the effective smoothness in the state space is s = 2 + α. This implies that the diffusion term effectively reduces the complexity of the function class compared to the standard Lipschitz assumption (s = 1), where the exponent would be d. 11
QI
T HE C OMPLEXITY-S TABILITY T RADE - OFF However, we must address the dependency on the time discretization δ. From the smoothing esti1+α mate (20), the radius of the function class scales as R ≈ O(δ − 2 ). Consequently, the network complexity required to achieve an error ϵ scales roughly as:
d
1
Complexity(δ, ϵ) ∼ δ
1+α 2
2+α
.
(28)
ϵ
This reveals a fundamental trade-off: • Smoothing Benefit: For a fixed δ, the parabolic nature of the Bellman operator constrains the target functions to a smooth manifold, making them efficiently learnable. • Stiffness Cost: As δ → 0 (approximating the continuous-time control limit), the Bellman operator approaches the Identity, and the smoothing effect acts over an infinitesimally small time. The functions become stiff (large second derivatives), inflating the constant R and requiring larger networks (increased width/depth) to resolve the sharp curvatures near decision boundaries. A NISOTROPIC I NDUCTIVE B IAS The distinction between state and action regularity dictates a specific inductive bias for the network architecture. A standard Multi-Layer Perceptron (MLP) treating (x, a) as a concatenated vector represents an isotropic prior. Corollary 14 suggests that the optimal architecture should be anisotropic: • The state subnet should utilize smooth activation functions (e.g., tanh, SiLU) to match the Cx2,α regularity. • The action subnet should utilize piecewise-linear activations (e.g., ReLU) to model the potentially non-smooth (Ca0,1 ) dependence and the sharp transitions induced by the max operator. This motivates the specific design of the Operator-Network proposed in Appendix A.
4. Universal Approximation Analysis In this section, we establish the convergence rate of the proposed Neural Operator scheme. A rigorous analysis must address the anisotropic nature of the value function regularity established in Section 3: the Bellman iterates exhibit parabolic smoothing in the state space S but remain merely Lipschitz continuous in the action space A. Standard isotropic approximations (e.g., standard MLPs) suffer from the curse of dimensionality over the combined dimension d + da . We prove that our Tensor-Product DeepONet architecture (Section A.2) successfully decouples these complexities, allowing the parabolic smoothing to reduce the effective dimension of the state space without being bottlenecked by the action space geometry. 4.1. Main Theorem: Decoupled Complexity Analysis We now state the main approximation result. We seek a global uniform error bound sup0≤k≤N ∥Q̂(k) − Q∗(k) ∥∞ ≤ ε for approximation of the Bellman-value-iteration sequence. The following theorem quantifies the explicit computational cost (network size and sensor density) required to achieve this 12
D EEP Q-L EARNING ON H OLDER S PACES
precision for the direct learning of the Bellman map. We interpret this as a theorem about the complexity of the Q-learning target class in the present diffusion setting: it analyzes the Bellman objects that value-based methods aim to approximate, while stopping short of a full sampling-andoptimization analysis for stochastic Q-learning training. Theorem 15 (Global Approximation of Bellman Iterates and Q-Learning Target Complexity) Grant Assumptions 1, 3, 23, 24, and 25. Let T be the finite horizon and δ the time-step. Let α ∈ (0, 1) be the Hölder exponent. Define the stiffness scaling factor γ(α) derived from the heat kernel singularity as: 3+α γ(α) := . (29) 2(2 + α) For any target tolerance ε > 0, there exists a Tensor-Product DeepONet with width W and a sensor set of size m such that the global error satisfies sup0≤k≤N ∥Q̂(k) − Q∗(k) ∥∞ ≤ ε. The required resources scale as follows, decoupling the state diffusion from the action geometry: d 3+α d (30) W ∼ O ε− 2+α δ −γ(α)d + (εδ)−da , m ∼ O ε− 2 δ − 4 d · (εδ)−da . This bound assumes the Branch Network approximates a monotone (positivity-preserving) reconstruction operator ensuring strict L∞ stability, and that the target is the direct map Q(k+1) = Tδ Q(k) . Proof See Appendix C.1.
5. Conclusion We have established an approximation-theoretic framework for Bellman value iteration in continuous stochastic control by grounding the analysis in the uniform ellipticity of the environment noise rather than restrictive global Lipschitz assumptions. By exploiting interior Schauder estimates, we proved that the Bellman operator maps bounded iterates into the anisotropic Hölder space Cx2,α (S) × Ca0,1 (A). Our proposed Tensor-Product DeepONet architecture explicitly leverages this mixed regularity, decoupling the smooth state manifold from the non-smooth policy landscape to mitigate the curse of dimensionality. Crucially, we quantified the stiffness–complexity trade-off, proving that the required network capacity scales as O(δ −γ(α)d ) due to the singularity of the heat kernel as δ → 0. We further validated the restriction to compact domains via the exponential decay of diffusion exit probabilities. We view these results as a core theoretical contribution to Q-learning in the continuous-time diffusion regime: they characterize the regularity and approximation complexity of the Bellman targets that value-based methods are built to learn. What the paper does not provide is a full convergence theorem for practical Deep Q-Learning with exploration, replay, and stochastic optimization. Rather, it isolates the Bellman side of the theory and shows that, in this regime, those targets enjoy exploitable anisotropic regularity and admit quantifiable approximation guarantees.
Acknowledgments We thank a bunch of people and funding agency. 13
QI
References Mohammad Gheshlaghi Azar, Ian Osband, and Rémi Munos. Minimax regret bounds for reinforcement learning. In International Conference on Machine Learning (ICML), pages 263–272. PMLR, 2017. Michael G Crandall, Hitoshi Ishii, and Pierre-Louis Lions. User’s guide to viscosity solutions of second order partial differential equations. Bulletin of the American Mathematical Society, 27 (1):1–67, 1992. George Cybenko. Approximation by superpositions of a sigmoidal function. Mathematics of control, signals and systems, 2(4):303–314, 1989. Jianqing Fan, Zhaoran Wang, Yuchen Xie, and Zhuoran Yang. A theoretical analysis of deep qlearning. In Learning for dynamics and control, pages 486–489. PMLR, 2020. Wendell H Fleming and Halil M Soner. Controlled Markov processes and viscosity solutions, volume 25. Springer Science & Business Media, 2006. Kurt Hornik. Approximation capabilities of multilayer feedforward networks. Neural networks, 4 (2):251–257, 1991. Chi Jin, Zeyuan Allen-Zhu, Sebastien Bubeck, and Michael I Jordan. Is q-learning provably efficient? In Advances in Neural Information Processing Systems (NeurIPS), volume 31, 2018. Nicolai V Krylov. Controlled diffusion processes. Springer, 1980. Nikolai Vladimirovich Krylov. Nonlinear elliptic and parabolic equations of the second order. Springer, 1987. Zongyi Li, Nikola Kovachki, Kamyar Azizzadenesheli, Burigede Liu, Kaushik Bhattacharya, Andrew Stuart, and Anima Anandkumar. Fourier neural operator for parametric partial differential equations. In International Conference on Learning Representations (ICLR), 2021. Lu Lu, Pengzhan Jin, Guofei Pang, Zhongqiang Zhang, and George Em Karniadakis. Learning nonlinear operators via deeponet based on the universal approximation theorem of operators. Nature machine intelligence, 3(3):218–229, 2021. Alessandra Lunardi. Analytic semigroups and optimal regularity in parabolic problems. Springer Science & Business Media, 2012. Hrushikesh N Mhaskar and Tomaso Poggio. Deep vs. shallow networks: An approximation theory perspective. Analysis and Applications, 14(06):829–848, 2016. Volodymyr Mnih, Koray Kavukcuoglu, David Silver, Andrei A Rusu, Joel Veness, Marc G Bellemare, Alex Graves, Martin Riedmiller, Andreas K Fidjeland, Georg Ostrovski, et al. Human-level control through deep reinforcement learning. Nature, 518(7540):529–533, 2015. Tomaso Poggio, Hrushikesh Mhaskar, Lorenzo Rosasco, Brando Miranda, and Qianli Liao. Why and when can deep-but not shallow-networks avoid the curse of dimensionality: A review. International Journal of Automation and Computing, 14(5):503–519, 2017. 14
D EEP Q-L EARNING ON H OLDER S PACES
Qian Qi. Universal approximation theorem of deep q-networks. In Proceedings of the 42nd International Conference on Machine Learning, volume 267 of Proceedings of Machine Learning Research, pages 50024–50046. PMLR, 2025. Qian Qi. Viscosity convergence analysis for deep q-networks. Revision submitted to the Journal of Machine Learning Research, 2026. Daniel Revuz and Marc Yor. Continuous Martingales and Brownian Motion, volume 293 of Grundlehren der mathematischen Wissenschaften. Springer Science & Business Media, 3rd edition, 2013. ISBN 978-3-642-08400-3. Reprint of the 1999 Edition. Corentin Tallec, Léonard Blier, and Yann Ollivier. Making deep q-learning methods robust to time discretization. In International Conference on Machine Learning (ICML), pages 6096–6104. PMLR, 2019. Matus Telgarsky. Neural networks and rational functions. In Proceedings of the 34th International Conference on Machine Learning, pages 3387–3393. PMLR, 2017. John Tsitsiklis and Benjamin Van Roy. Analysis of temporal-diffference learning with function approximation. Advances in neural information processing systems, 9, 1996. Holger Wendland. Scattered data approximation, volume 17. Cambridge university press, 2004. Dmitry Yarotsky. Error bounds for approximations with deep relu networks. Neural networks, 94: 103–114, 2017.
15
QI
Appendix A. Operator-Based Network Architecture To leverage the smoothing properties established in Section 3, we propose a network architecture designed to approximate the operator map Tδ : C(D) → Cx2,α (D). Consequently, we adopt a meshfree approach based on the Deep Operator Network (DeepONet) framework (Lu et al., 2021), which generalizes the universal approximation theorem to non-linear operators on Banach spaces. A.1. The DeepONet Framework We seek a neural operator Nθ that approximates the Bellman transition map Tδ : Hα,Lip → Hα,Lip . Following the Deep Operator Network (DeepONet) framework (Lu et al., 2021), we decompose the architecture into a Branch Net (which encodes the input function space) and a Trunk Net (which encodes the output domain). A.1.1. S ENSOR D ISCRETIZATION AND E NCODING A fundamental challenge in operator learning is that the input to the network is a function, an infinite-dimensional object. In practice, we access the function via a finite set of observations. Let {zj }m j=1 ⊂ S × A be a fixed set of m sensor locations. We define the discretization operator Pm : C(S × A) → Rm as: Pm (Q) := [Q(z1 ), Q(z2 ), . . . , Q(zm )]⊤ .
(31)
In the context of RL, these sensors correspond to the states stored in a Replay Buffer or a fixed quadrature grid used for fitting the value iteration. A.1.2. N ETWORK D EFINITION The Deep Bellman Operator Network is defined as the map Nθ : Rm × (S × A) → R. Let Qin be the value function from the previous iteration. The approximate next iterate Q̂out (x, a) is given by the inner product of the Branch and Trunk outputs: Q̂out (x, a) ≈ Nθ (Pm (Qin ))(x, a) =
p X k=1
bk (Qin (z1 ), . . . , Qin (zm )) · | {z } Branch Net β:Rm →Rp
τk (x, a) | {z }
+b0 .
Trunk Net τ :S×A→Rp
(32) Here, θ represents the trainable weights of the neural networks β and τ . • The Branch Net β learns to extract a latent embedding (coefficients) of the current value landscape from the sensor readings. • The Trunk Net τ learns a continuous basis for the image space of the Bellman operator. Remark 16 (The Role of Sensor Density) The standard Universal Approximation Theorem for Operators assumes the input is the entire function Qin . By using Pm (Qin ), we introduce an encoding error (or aliasing error). However, Corollary 14 guarantees that Qin belongs to a compact subset of the Hölder space Hα,Lip . Let h = supy∈S×A minj ∥y − zj ∥ be the fill-distance of the sensors. Due to the Hölder continuity of the inputs, the reconstruction error is bounded by Lhα . Thus, as the number of sensors m → ∞ (and h → 0), the discrete representation Pm (Qin ) uniquely characterizes Qin up to arbitrary precision, making the operator approximation well-posed. 16
D EEP Q-L EARNING ON H OLDER S PACES
A.2. Anisotropic Trunk Architecture The structural results of Corollary 14 indicate that the target functions exhibit mixed regularity: they are locally C 2,α in the state variable x (due to parabolic smoothing) but only Lipschitz continuous in the action variable a (due to the parametric nature of the operator and the supremum envelope). A standard isotropic MLP (e.g., fully connected layers on concatenated inputs [x, a]) fails to exploit this structure efficiently. Furthermore, a full tensor-product basis (via the Kronecker product of feature vectors) would result in a feature dimension scaling quadratically with the width. To remedy this, we propose a Separable Rank-p Trunk Architecture that explicitly forces a low-rank decomposition of the state-action coupling. Definition 17 (Separable Rank-p Trunk) The Trunk Network τ (x, a) : S × A → Rp is constructed as the Hadamard product (element-wise multiplication) of two independent sub-networks. This enforces that the learned basis functions are separable: τ (x, a) = Nstate (x; θx ) ⊙ Naction (a; θa ).
(33)
Here, ⊙ denotes the element-wise product. Consequently, the k-th component of the trunk output is given by τk (x, a) = ϕk (x)ψk (a), where: (i) The State Basis ϕ(·) = Nstate (·) : Rd → Rp employs smooth activation functions σs ∈ ∞ C (R) (e.g., tanh, SiLU) to match the C 2,α (S) regularity of the heat kernel eigenfunctions. (ii) The Action Basis ψ(·) = Naction (·) : Rda → Rp employs piecewise-linear activation functions σr (e.g., ReLU) to capture the non-differentiable transitions of the optimal policy. T HEORETICAL J USTIFICATION The Hadamard formulation in (33) combined with the DeepONet inner product (32) implies that the network approximates the target function Q̂ via a Canonical Polyadic (CP) decomposition of rank p: p X Q̂(x, a) = βk (Qin ) · ϕk (x) · ψk (a). (34) k=1
This is structurally distinct from a full Tensor Product basis (which would imply terms of the form ϕi (x)ψj (a) for all pairs). While a full tensor product avoids any structural constraints on the stateaction coupling, it suffers from an interaction parameter complexity of O(p2 ). To justify the restriction to a diagonal (rank-p) structure, we must carefully distinguish between the spatial smoothing of the PDE and the parametric dependence on the control variable. For a fixed action a, the dominant term in the backward Bellman update is the convolution with the fundamental solution Γa (t, x, y). Because the elliptic operator is compact on bounded domains, the spatial update admits a rapidly decaying spectral expansion involving the eigenfunctions ek (x; a) of the generator La . Critically, the true spatial eigenfunctions ek (x; a) depend on a because the diffusion tensor Σ(x, a) and drift b(x, a) vary with the action. Therefore, the exact global value function is not mathematically separable. However, because the dependence of the environment coefficients on a 17
QI
is uniformly Lipschitz (Assumption 1), the manifold of spatial solutions M = {Q(·, a) : a ∈ A} exhibits low ε-rank. By employing a CP decomposition, the network learns a finite, universal set of spatial basis functions {ϕk (x)}pk=1 that span the dominant eigenspaces of La across all a ∈ A. Concurrently, the action network {ψk (a)}pk=1 learns the non-smooth, action-dependent mixing coefficients required to linearly combine these spatial modes. This formulation enforces a principled anisotropic inductive bias: it allocates the smooth network capacity (e.g., tanh activations) exclusively to the spatial factors ϕk ∈ C 2,α (S) to emulate parabolic smoothing, and allocates the non-smooth capacity (e.g., ReLU activations) exclusively to the action factors ψk ∈ C 0,1 (A) to capture the non-differentiable kinks of the policy landscape. This decoupling reduces the interaction complexity to O(p) while efficiently parameterizing the mixed-regularity hypothesis space Hα,Lip defined in Corollary 14. A.3. Complexity Analysis: Dimensionality and Discretization Standard grid-based value iteration suffers from the classic Curse of Dimensionality: to halve the approximation error, the computational cost scales as 2d . Our mesh-free Neural Operator approach mitigates this scaling through the regularizing effect of diffusion, but it is subject to trade-offs regarding the time-discretization δ and sensor density m. A.3.1. R EDUCTION OF M ETRIC E NTROPY The fundamental difficulty of learning the Q-function is governed by the size of the reachable function class Q∗ . As established in Corollary 14, the uniform ellipticity of the SDE ensures that Q∗ is a bounded subset of the anisotropic space Cx2,α ∩ Ca0,1 . The Kolmogorov ϵ-entropy (log-covering number) of the Hölder ball of radius R in dimension d with smoothness s scales as log N ≍ (R/ϵ)d/s . In our setting, the effective smoothness in the state space is s = 2 + α. Thus, the complexity scales as: d ! R 2+α Complexity(ϵ) ∼ exp C . ϵ
(35)
Compared to standard Lipschitz continuity assumptions (s = 1, exponent d), the diffusion term reduces the effective dimensionality exponent by a factor of roughly 3 (since 2 + α ≈ 3). While the dependence on d remains exponential in the worst case (covering the entire uniform norm), the volume of the function space is significantly compressed. A.3.2. S ENSOR A LIASING AND S TATISTICAL R ATES The DeepONet replaces the spatial grid with a set of m scattered sensors {zj }m j=1 . The total error decomposes into the operator approximation error and the input encoding (aliasing) error: ∥N − T ∥ ≤ ∥N (Pm Q) − T (RPm Q)∥ + ∥T (RPm Q − Q)∥, | {z } | {z } Neural approx. error
Sensor encoding error
where Pm is the evaluation map and R is an optimal reconstruction operator. 18
(36)
D EEP Q-L EARNING ON H OLDER S PACES
• Worst-Case Error: To guarantee uniformly bounded error over the compact set, the sensor fill-distance h must satisfy Lhα < ϵ. This requires m ∼ ϵ−d/α sensors, which still suffers from the curse of dimensionality. • Generalization Error: However, in Deep RL, we minimize the L2µ error over an exploration distribution µ. By standard statistical p learning theory (e.g., Rademacher complexity bounds), the generalization error scales as O(1/ Nsamples ), independent of the dimension d. The Neural Operator efficiently learns the map on the support of µ, effectively breaking the curse for the visited states. A.3.3. T HE C OMPLEXITY-S TABILITY T RADE - OFF (δ → 0) 1+α
From the Schauder estimate (20), the radius of the function class scales as Rδ ∼ δ − 2 . Substituting this into the entropy bound, the required network capacity W behaves as:
d
1
W (δ, ϵ) ∼ δ
1+α 2
2+α
.
(37)
ϵ
This reveals a fundamental Stability-Complexity Trade-off: • Smoothing Regime (Fixed δ): For a discrete MDP with fixed step δ, the operator is strongly smoothing, and the Q-functions lie on a low-complexity manifold. The network size remains manageable. • Stiff Regime (δ → 0): As we approximate the continuous-time limit, the operator approaches the Identity, and the smoothing vanishes. The target functions become stiff (large second derivatives), and the required network capacity explodes. This suggests that value-based discretizations built from Bellman iteration are best conditioned when δ is chosen large enough to allow diffusion to act, rather than trying to approximate the SDE with infinitesimal steps directly. A.4. Recursive Residual Formulation For problems with long horizons (large N = T /δ), we construct the solution by recursively applying the operator. Direct composition of neural operators can lead to optimization difficulties akin to the vanishing gradient problem. However, the parabolic nature of the Bellman equation suggests a natural residual structure. Since the value function evolves continuously in time, for small δ, the linear diffusion part of Tδ is close to the Identity; the nonlinear maximization term also requires the current action to be neargreedy for the normalized residual to remain O(1). We propose a Residual Operator architecture that learns the temporal difference rather than the absolute value. A.4.1. D ISCRETE E ULER S CHEME IN F UNCTION S PACE Let Q̂(k) denote the approximate Q-function at time-step k (corresponding to physical time tk = T − kδ). We formulate the update as: Q̂(k+1) (x, a) = Q̂(k) (x, a) + δ · Nres Q̂(k) (x, a), (38) 19
QI
with the initialization Q̂(0) (x, a) = g(x). Here, Nres is a DeepONet trained to approximate the normalized displacement operator: Rδ (Q) :=
Tδ Q − Q . δ
(39)
From the PDE characterization in Eq. (15), the normalized residual decomposes into a local generator contribution plus a maximization-gap contribution: Rδ (Q)(x, a) = Lat Q(x, a) + r(t, x, a) +
1 (VQ (x) − Q(x, a)) + o(1), δ
(40)
where VQ (x) := supa′ ∈A Q(x, a′ ), with the expansion interpreted locally for time-dependent coefficients through the parabolic evolution family. Thus the residual converges to the Hamiltonian only on the zero-gap or near-greedy regime where VQ (x) − Q(x, a) = o(δ). The residual operator network should therefore be read as learning the Hamiltonian on this active/near-active set, not as removing the Bellman maximization gap for arbitrary off-policy actions. A.4.2. A DVANTAGES OF THE R ESIDUAL A PPROACH (i) Well-Conditioned Learning: In the stiff regime where δ is small, the diffusion-evolution part of Q 7→ Tδ Q is dominated by the Identity. Learning the full map directly forces the network to expend capacity approximating this nearly identity component. The residual formulation subtracts it and focuses capacity on the dynamics, provided the maximization gap is separately controlled on the active or near-greedy action set. (ii) Flow Consistency: This architecture effectively implements a Forward Euler discretization of the abstract Ordinary Differential Equation (ODE) in the Banach space Hα,Lip : d Q(τ ) = F(Q(τ )), dτ
Q(0) = g,
(41)
where τ is time-to-go and F is the non-linear Bellman generator. This ensures that the inductive bias of the network aligns with the physical flow of information. (iii) Iterative Stability: By sharing the weights of Nres across time-steps (Recurrent Operator Learning), we enforce that the laws of physics (the environment dynamics) remain constant, while the value function evolves. This dramatically reduces the number of parameters compared to training separate networks for each time step.
Appendix B. Detailed Proofs and Mathematical Derivations In this appendix, we provide rigorous proofs for the theorems and lemmas presented in the main text. We explicitly address the technical concerns raised regarding the spectral separability of the diffusion operator and the dependence of Lipschitz constants on the time horizon. 20
D EEP Q-L EARNING ON H OLDER S PACES
B.1. Proof of Theorem 7 (Refined Interior Parabolic Smoothing) Proof We fix the current action a and introduce the time-to-go variable τ = t + δ − s. Let v(τ, x) := u(t + δ − τ, x). Then v solves the forward parabolic equation on (0, δ] × Rd : ∂τ v = La v + r(t + δ − τ, ·, a),
v(0, ·) = VΦ .
Since the coefficients may depend on time, the solution is represented by the parabolic evolution family U a (τ, s) associated with the non-autonomous operators Lat+δ−τ : Z δ v(δ) = U a (δ, 0)VΦ + U a (δ, s)r(t + δ − s, ·, a) ds. (42) 0
(i) Regularity inheritance of the Maximization (Initial Data). First, we establish the regularity of the initial condition v(0) = VΦ . Since Φ(·, a′ ) is uniformly Lipschitz with constant LΦ , for any x, y ∈ Rd : |VΦ (x) − VΦ (y)| = sup Φ(x, a′ ) − sup Φ(y, a′ ) a′
a′
′
≤ sup |Φ(x, a ) − Φ(y, a′ )| ≤ LΦ ∥x − y∥. a′
Thus, VΦ ∈ C 0,1 (Rd ) with semi-norm [VΦ ]C 0,1 ≤ LΦ . (ii) Smoothing of Lipschitz Data (The Homogeneous Term). We estimate ∥U a (δ, 0)VΦ ∥C 2,α . The parabolic evolution family satisfies the same local smoothing scale as the frozen-coefficient heat kernel under the uniform Hölder and ellipticity bounds: it maps C k (Rd ) to C m (Rd ) for m > k with norm scaling τ −(m−k)/2 . We decompose the C 2,α norm into the integer and fractional parts. The dominant singularity arises from the highest order derivative. We regard VΦ as an element of the interpolation space roughly corresponding to C 1 . To be precise, since VΦ ∈ C 0,1 (Rd ), it resides in the Besov space 1 B∞,∞ . Using the standard smoothing property of the heat kernel (see Lunardi (2012, Prop 2.2)), the cost of lifting regularity from C 0,1 (regularity index 1) to C 2,α (regularity index 2+α) is determined by the difference in indices: Diff = (2 + α) − 1 = 1 + α. (43) Consequently, the bound scales as τ −Diff/2 : 1+α
∥U a (δ, 0)VΦ ∥C 2,α (Rd ) ≤ C1 δ − 2 ∥VΦ ∥C 0,1 (Rd ) .
(44) 2+α
Note that if we had only assumed VΦ ∈ L∞ (Rd ) (regularity 0), the cost would have been δ − 2 . The Lipschitz assumption improves the exponent by 1/2. (iii) Regularity of the Source Potential (The Inhomogeneous Term). The second term involves the running reward r. By Assumption 1, r ∈ C α/2,α . Standard Schauder estimates for the inhomogeneous problem guarantee that the integral operator preserves the regularity of the source without singular scaling in time (modulo the bounded horizon). There exists C2 independent of small δ such that: Z δ
U a (δ, s)r(t + δ − s, ·, a) ds
0
C 2,α
21
≤ C2 ∥r∥C α/2,α .
(45)
QI
(iv) Synthesis. Combining (44) and (45), we obtain: 1+α ∥(Tδ Φ)(t, ·, a)∥C 2,α ≤ C δ − 2 ∥VΦ ∥C 0,1 + ∥r∥C α/2,α .
(46)
1+α
For δ < 1, the term δ − 2 dominates.
B.2. Proof of Theorem 12 (Truncation Error) Proof Let S be the target domain and Sρ = {y : d(y, S) ≤ ρ} be the computational domain. Let u(t, x) be the solution to the global Cauchy problem and uρ (t, x) be the solution with the process killed upon exiting Sρ . By the Feynman-Kac formula: Z t+δ u(t, x) − uρ (t, x) = Ex r(s, Xs )ds + VΦ (Xt+δ ) 1{τρ ≤t+δ} . (47) t
Let Mval = ∥VΦ ∥∞ + δ∥r∥∞ . |u(t, x) − uρ (t, x)| ≤ Mval · Px (τρ ≤ t + δ).
(48)
Since x ∈ S, the distance to the boundary is at Rleast ρ. We R apply the concentration inequality for Itô diffusions. The displacement is Xt+δ − x = bds + σdWs . Z Px (sup ∥Xt+s − x∥ ≥ ρ) ≤ P(sup ∥ σdW ∥ ≥ ρ − ∥b∥∞ δ). (49) s≤δ
s≤δ
for continuous local martingales Ms = RAssuming ρ > ∥b∥∞ δ, we apply the Bernstein 2inequality σdW . The quadratic variation is ⟨M ⟩δ ≤ ∥σ∥F δ =: Λ2 δ. λ2 . (50) P(sup ∥Ms ∥ ≥ λ) ≤ 2d exp − 2dΛ2 δ s≤δ Setting λ = ρ − ∥b∥∞ δ and equating the probability to εtrunc /Mval , we solve for ρ: s εtrunc (ρ − ∥b∥∞ δ)2 2dMval = exp − =⇒ ρ ≥ ∥b∥∞ δ + Λ 2dδ log . 2dMval 2dΛ2 δ εtrunc
(51)
B.3. Proof of Corollary 14 (Anisotropic Compactness) Proof We proceed by induction on the iteration index k, analyzing the spatial and parametric regularity separately. Let δ > 0 be fixed. While the PDE estimates in Theorem 7 are derived globally on Rd to ensure well-posedness without boundary effects, the compactness result requires restricting the function sequence to a compact subdomain. Let K ⊂ S be a compact set with Lipschitz boundary (the region of interest). 22
D EEP Q-L EARNING ON H OLDER S PACES
(i) Uniform L∞ Boundedness. We first establish that the sequence is uniformly bounded in the supremum norm. Since the reward r and terminal condition g are bounded, and the discount factor implies a finite effective horizon (or finite T ), there exists a constant M independent of k such that: sup ∥Q(k) ∥L∞ (Rd ×A) ≤ M.
0≤k≤N
Let V (k) (x) := supa′ ∈A Q(k) (x, a′ ). It follows immediately that ∥V (k) ∥L∞ (Rd ) ≤ M . (ii) Uniform Spatial Lipschitz Continuity. Before invoking higher-order interior estimates, we must guarantee that the spatial Lipschitz constant of the sequence does not explode. By Assumption 1(i), the coefficients b, σ are globally Lipschitz in x. Using standard estimates for the stochastic flow of diffeomorphisms, the spatial gradients of the expected reward and terminal value grow by at most a factor of eCb δ per step, where Cb = ∥∇x b∥∞ . The maximization step V (k) (x) = supa Q(k) (x, a) preserves the spatial Lipschitz constant. Over the finite horizon T = N δ, the spatial Lipschitz constant accumulates as: sup ∥Q(k) (·, a)∥C 0,1 (Rd ) ≤ eCb T ∥g∥C 0,1 (Rd ) + T eCb T ∥r∥C 0,1 =: LS .
(52)
a∈A
Crucially, this uniform Lipschitz bound LS depends only on the finite horizon T and the coefficient bounds, independent of δ. (iii) Spatial Regularity (The Healing Step). Consider the update Q(k+1) (·, a) = Tδ Q(k) (·, a). As defined in (15), this requires solving the parabolic Cauchy problem with terminal condition V (k) . Applying the refined global Schauder estimate from Theorem 7 (Eq. 20) which leverages the uniform C 0,1 bound established above: 1+α ∥Q(k+1) (·, a)∥C 2,α (Rd ) ≤ C ∥r(·, a)∥C α/2,α + δ − 2 ∥V (k) ∥C 0,1 (Rd ) . (53) Since ∥V (k) ∥C 0,1 ≤ LS for all k, we successfully propagate the tighter singularity scaling without recursive explosion, obtaining a uniform spatial bound: 1+α
sup sup ∥Q(k) (·, a)∥C 2,α (Rd ) ≤ KS (δ) := C ′ · δ − 2 .
(54)
0≤k≤N a∈A
Note that while KS (δ) diverges as δ → 0, it remains finite and constant for the fixed time-step analysis of this Corollary. (iv) Parametric Regularity in Action Space. We now address the Lipschitz continuity with respect to a. The value Q(k+1) (t, x, a) is the solution at time t of the linear PDE: ∂τ u = La u + r(·, a),
u(0) = V (k) .
Notice that the action a appears only as a parameter in the coefficients (b, σ) and the source term r. The initial data V (k) (x) depends on the previous iteration’s optimization and is independent of the current action parameter a. By Assumption 1(iii), the coefficients are globally Lipschitz in a. Using standard estimates for the sensitivity of parabolic PDEs to parameters, the gradient with respect to the parameter satisfies: ∥∇a Q(k+1) ∥∞ ≤ CT ∥∇a r∥∞ + ∥∇a L∥∞ ∥∇x Q(k+1) ∥∞ . (55) 23
QI
From step (iii), ∥∇x Q(k+1) ∥∞ ≤ ∥Q(k+1) ∥C 2,α is uniformly bounded by KS (δ). Consequently, the Lipschitz constant with respect to actions is uniformly bounded: sup sup ∥Q(k) (x, ·)∥C 0,1 (A) ≤ LA (δ, T ).
(56)
0≤k≤N x∈Rd
(v) Topological Compactness. Combining (iii) and (iv), the sequence Q∗δ = {Q(k) }0≤k≤N resides in the bounded set: B = Q : sup ∥Q(·, a)∥C 2,α (Rd ) ≤ KS and sup ∥Q(x, ·)∥C 0,1 (A) ≤ LA . a
x
To establish pre-compactness, we restrict our attention to the compact state domain K ⊂ Rd and the compact action space A. • By the compact embedding of Hölder spaces, the inclusion C 2,α (K) ,→ C 1 (K) is compact. • By the Arzelà-Ascoli theorem, the embedding C 0,1 (A) ,→ C 0 (A) is compact. Since Q∗δ is uniformly bounded in the stronger product norm Hα,Lip , any sequence drawn from this uniformly bounded family admits a limit point in the topology of C 1 (K) × C 0 (A), proving the stated compactness.
Appendix C. Appendix for Universal Approximation Analysis C.1. Proof of Theorem 15 Proof We derive the complexity bounds by establishing the local consistency tolerance required to satisfy the global stability constraint, and then determining the network capacity and sensor density necessary to achieve this tolerance given the anisotropic regularity of the reachable set Q∗δ . (i) Global Error Stability and Local Budget. We first determine the maximum permissible singlestep error ηmax . Let Ek := ∥Q̂(k) − Q∗(k) ∥∞ . From Lemma 29, the error propagates according to: Ek+1 ≤ eρstab δ · 1 · Ek + ηk , where ηk = εnet + εenc comprises the network approximation error and the sensor encoding error, and we explicitly leverage the stable Lebesgue constant Λm = 1. The constant ρstab ≥ 0 is any one-step sup-norm stability exponent for the Bellman propagation and reconstruction pipeline; in the ideal Markov expectation part alone one may take the non-expansive case ρstab = 0. Iterating over N = T /δ steps with E0 = 0: EN ≤ ηmax
N −1 X k=0
24
ekρstab δ .
(57)
D EEP Q-L EARNING ON H OLDER S PACES
To ensure the final error EN ≤ ε, the local error budget must scale linearly with the time-step: ( ρstab /(eρstab T − 1), ρstab > 0, ηmax ≤ CT · δ · ε, CT := (58) 1/T, ρstab = 0, We allocate this budget equally among the error sources: εnet ≤ 13 CT δε, and the interpolation errors for state and action spaces. (ii) Anisotropic Resolution Requirements. We determine the required sensor fill-distances hS and hA . Action Space (Lipschitz Regime): The dependence on actions is strictly Lipschitz with constant LA . The encoding error constraint is LA hA ≤ 31 CT δε. Thus: hA ∼ O(δε).
(59)
State Space (Refined Smoothing Regime): Using the Refined Smoothing Estimate (Theorem 7), 1+α the regularity constant of the output scales as KS (δ) ≍ δ − 2 . To ensure L∞ stability, the Branch Net models a monotone reconstruction operator bounded by O(h2 ) consistency. The reconstruction error scales as KS (δ)h2S . The constraint is: 1+α 1 Cdiff δ − 2 h2S ≤ CT δ ε. 3
Solving for hS : hS ≲
δε 1+α
δ− 2
1/2
3+α 1+α 1/2 = ε1/2 · δ 4 . = ε · δ 1+ 2
(60)
(iii) Resource Complexity Analysis. We map the resolutions to Sensor Count m and Network Width W . Sensor Count (m): For the quasi-uniform cover, m ≈ (1/h)d . 3+α
d
3+α
1/2 δ 4 )−d = ε− 2 δ − 4 d . • State sensors: mS ∼ h−d S ∼ (ε a ∼ (δε)−da . • Action sensors: mA ∼ h−d A
Network Width (W ) and Branch Net Capacity: By Assumption 24, W scales with the metric entropy of the target class relative to the tolerance τnet ∼ δε. Note that while the data complexity (sensor count) is capped by the O(h2 ) interpolation stability limit, the network’s universal approximation capability still explores the full C 2,α smoothness of the function manifold. 1+α
The covering radius is Rδ ≈ KS (δ) ∼ δ − 2 . The effective entropy ratio for the state component is: 1+α Rδ δ− 2 1 3+α ≍ = δ− 2 . (61) τnet δε ε The required width for the state trunk WS follows the entropy scaling (R/τ )d/s with s = 2 + α: d d 3+α 2+α WS ∼ ε−1 δ − 2 = ε− 2+α δ −γ(α)d , 25
(62)
QI
where the stiffness scaling factor is explicitly identified as: γ(α) :=
3+α . 2(2 + α)
(63)
The Action Space Curse: For the action trunk, the Lipschitz complexity yields WA ∼ (LA /δε)da . Combining the components via the decoupled DeepONet architecture yields the final resource bounds.
C.2. Approximation Setup and Assumptions To derive explicit convergence rates, we must characterize the geometric complexity of the function class generated by the Bellman recursion. We rely on the anisotropic regularity established in Corollary 14 and formalized via the notion of Metric Entropy. C.2.1. A NISOTROPIC R EACHABLE S ETS Let S ⊂ Rd be a compact domain satisfying the cone condition, and let A ⊂ Rda be a compact metric space. We define the set of reachable action-value functions at time-step δ. Definition 18 (Reachable Set Q∗δ ) Let Q∗δ denote the image of the unit ball of L∞ (S × A) under the Bellman operator Tδ . Based on the interior Schauder estimate (20), any Q ∈ Q∗δ satisfies the following disjoint regularity bounds: (i) State Regularity (Smoothing): For every fixed a ∈ A, the spatial section Q(·, a) belongs to the Hölder space C 2,α (S). The semi-norm is controlled by the diffusion time-scale: 1+α
sup ∥Q(·, a)∥C 2,α (S) ≤ KS (δ) := Cdiff δ − 2 ,
(64)
a∈A
where Cdiff depends on the uniform ellipticity constant ν and dimension d. (ii) Action Regularity (Lipschitz): For every fixed x ∈ S, the parametric section Q(x, ·) belongs to C 0,1 (A). Unlike the local coefficients, the Lipschitz constant of the value function depends on the global stability of the SDE flow over the horizon T . By standard estimates for controlled diffusion processes (Krylov, 1980, Corollary 2.10), there exists a constant LA such that: sup
∥Q(t, x, ·)∥C 0,1 (A) ≤ LA (T, KLip ),
(65)
(t,x)∈DT
where LA (T, KLip ) depends exponentially on the horizon T and polynomially on the global Lipschitz constants KLip of the coefficients b, σ, r. Remark 19 (Stability with respect to Control Parameters) We caution that LA cannot be equated simply to sup ∥∇a coeff∥. Since the action a modifies the drift and diffusion, a perturbation in a propagates through the trajectory Xs . By Gronwall’s inequality, the sensitivity of the value function scales as O(eCT ). However, since T is fixed and finite, and the coefficients satisfy Assumption 1(iii), LA remains a dimension-independent constant distinct from the time-step regularization parameter δ. 26
D EEP Q-L EARNING ON H OLDER S PACES
C.2.2. A NISOTROPIC S ENSOR E NCODING AND S TABILITY A central challenge in applying Operator Learning to Reinforcement Learning is the discretization of the input function Qin . While the standard theory often assumes access to the full function, in practice we observe Qin via a finite set of sensors Zm . As noted in recent stability analyses of scattered data interpolation (Wendland, 2004), standard i.i.d. random sampling leads to the formation of clusters (Poisson clumping), causing the separation distance qZ to decay much faster than the fill distance hZ . This results in a diverging mesh ratio ρ = hZ /qZ → ∞, which causes the Lebesgue constant of the reconstruction operator to explode, destabilizing the iterative learning scheme. To guarantee the stability required for high-order Schauder estimates, we must replace generic random sampling with a controlled Quasi-Uniform covering of the domain. Definition 20 (Geometric Properties of Scattered Sensors) Let Ω ⊂ Rd be a compact domain. Let Zm = {zi }m i=1 ⊂ Ω be a set of m sensor locations. (i) The Fill-Distance (covering radius) is defined as: hZm ,Ω := sup min ∥x − zi ∥2 . x∈Ω zi ∈Zm
(66)
(ii) The Separation Radius is defined as: qZm :=
1 min ∥zi − zj ∥2 . 2 i̸=j
(67)
(iii) The Mesh Ratio is defined as ρZm := hZm ,Ω /qZm . A sequence of sensor sets is called Quasi-Uniform if there exists a constant cqu ≥ 1 such that ρZm ≤ cqu for all m. We define the sampling strategy via a greedy packing mechanism (Farthest Point Sampling) acting separately on the state and action spaces to preserve the tensor structure required by the DeepONet. Lemma 21 (Deterministic Quasi-Uniform Covering) Let Ω be a compact domain in Rd with Lipschitz boundary. There exists a deterministic construction of sensor sets Zm (via Farthest Point Sampling) such that for all m: clow m−1/d ≤ qZm ≤ hZm ,Ω ≤ chigh m−1/d .
(68)
Consequently, the mesh ratio is uniformly bounded: ρZm ≤ chigh /clow := cqu . Proof We construct the set Zm inductively. Let z1 ∈ Ω be arbitrary. For k = 1, . . . , m − 1, select: zk+1 = argmax min ∥x − zj ∥2 . x∈Ω
1≤j≤k
By definition, the separation distance at step k + 1 is exactly the fill distance of the set Zk . Standard sphere-packing arguments in Rd imply that the covering radius of an optimal packing scales as m−1/d . Specifically, since Ω is bounded, the union of disjoint balls B(zi , qZm ) must fit within Ω (modulo boundary effects), implying m · Vol(B1 )qZd m ≤ Vol(Ω), which gives the lower bound on q. Conversely, the covering property gives the upper bound on h. The ratio remains bounded by a constant depending only on the dimension d and the geometry of Ω.
27
QI
Remark 22 (Rejection of Random Sampling) We explicitly reject the use of Replay Buffers populated by i.i.d. exploration for the definition of the sensor set Zm . For i.i.d. uniform samples, qZm ∼ m−2/d (or worse) almost surely, while hZm ∼ m−1/d . This implies ρZm → ∞, which would cause the error coefficients in Lemma 29 to diverge. The set Zm in our architecture represents a fixed quadrature cores set, distinct from the transient training data. With the geometry of the sensors fixed, we rely on the stability of the reconstruction operator. Assumption 23 (Monotone Reconstruction and the Godunov Barrier) Let Z = ZS × ZA be the tensor product of two quasi-uniform sets constructed via Lemma 21. To prevent unconditional instability as δ → 0, the reconstruction operator must not amplify uniform errors. By Godunov’s theorem equivalents for scattered data, any linear scheme with consistency order > 2 must have a T /δ Lebesgue constant Λm > 1, which would cause the error to explode as Λm → ∞. Therefore, we assume the Branch Network β approximates a Monotone Averager (e.g., simplicial interpolation or a positive partition of unity) satisfying: (i) Unconditional L∞ Stability: Because R is a positive operator (convex combination of sensor values), its Lebesgue constant is exactly 1: Λm := sup ∥R(v)∥∞ = 1.
(69)
∥v∥∞ ≤1
(ii) Consistency (Capped at Second Order): While the target functions belong to C 2,α (S), positive operators can only reproduce polynomials up to degree 1. Consequently, the spatial interpolation error is bounded by O(h2 ) rather than O(h2+α ). For f ∈ C 2,α (S) and g ∈ C 0,1 (A): ∥f − RS (PZS f )∥∞ ≤ Capprox h2ZS ∥f ∥C 2,α ,
(70)
∥g − RA (PZA g)∥∞ ≤ Capprox hZA ∥g∥C 0,1 .
(71)
C.2.3. C APACITY D ECOUPLING VIA S EPARATION R ANK A standard Multi-Layer Perceptron (MLP) approximates functions in the joint Sobolev space W k,p (S× A), resulting in a width requirement scaling exponentially with d+da . However, the proposed SplitTrunk DeepONet (Eq. 33) approximates Q(x, a) via a finite sum of separable functions: Q̂(x, a) =
W X
ψk (x)ϕk (a).
(72)
k=1
The efficiency of this architecture relies on the concept of ε-Separation Rank. Since Tδ involves convolution with the Gaussian heat kernel (which has rapidly decaying singular values), the target functions Q ∈ Q∗δ admit low-rank approximations. This allows us to decompose the metric entropy of the hypothesis space. Assumption 24 (Decoupled Network Capacity) Let N (ε, F) denote the covering number of a function class F at scale ε. We assume the DeepONet with width W can approximate any Q ∈ Q∗δ such that the required width scales with the sum of the component entropies rather than the product: W (εnet ) ∼ O log N εnet , BKS (δ) (C 2,α (S)) + log N εnet , BLA (C 0,1 (A)) . (73) 28
D EEP Q-L EARNING ON H OLDER S PACES
Invoking the standard entropy bounds for Hölder balls, log N (ε, BR (C s )) ≍ (R/ε)d/s , we obtain the rigorous width bound: " d # KS (δ) 2+α LA da W (εnet , δ) ≤ Cnet + . (74) εnet εnet Assumption 25 (Regularized Hypothesis Space) The neural network training algorithm is constrained to the compact hypothesis class defined by the reachable set Q∗δ . Specifically, we assume the optimization includes a projection step or explicit regularization (e.g., spectral normalization of the weights) such that for every iteration k, the learned function Q̂(k) satisfies the anisotropic regularity bounds of the target class: Q̂(k) ∈ BKS (δ) (C 2,α (S)) ∩ BLA (C 0,1 (A)).
(75)
Consequently, the approximation produced by the network is admissible as an input for the sensor reconstruction operator R in the subsequent Bellman update. Remark 26 Assumption 24 formalizes the benefit of the specific inductive bias. The term (KS /ε)d/(2+α) represents the cost of learning the smooth manifold in state space, while (LA /ε)da represents the cost of resolving the non-smooth control policies. By decoupling them, we avoid the curse of dimensionality associated with the joint dimension d + da in the exponent. C.2.4. C ONSISTENCY AND S TABILITY D ECOMPOSITION A significant technical challenge in the analysis of Neural Operators for PDEs is the regularity mismatch between the target function and the hypothesis class. While Corollary 14 guarantees that the true Bellman iterate Q∗(k) belongs to the smooth anisotropic space C 2,α (S), the neural approximation Q̂(k) , typically a composition of piecewise-linear activations (e.g., ReLU), does not possess classical higher-order derivatives. Consequently, applying interpolation error bounds of the form ∥RPm f − f ∥ ≤ Chs ∥f ∥C s directly to f = Q̂(k) is ill-posed, as the norm ∥Q̂(k) ∥C s may be infinite or diverge with network depth. To resolve this, we employ a Consistency-Stability analysis. We decompose the total error into a consistency error (depending on the smoothness of the target Q∗ ) and a stability term (depending on the Lebesgue constant of the reconstruction operator). To ensure the stability of high-order reconstruction on scattered data, we must strengthen the sampling assumption from generic i.i.d. distributions (which may produce arbitrarily small separation distances) to quasi-uniform configurations. Assumption 27 (Quasi-Uniform Sensor Distribution) Let hZ,Ω be the fill-distance and qZ := 1 2 mini̸=j ∥zi −zj ∥ be the separation radius of the sensor set Z. We assume the sensor configuration is quasi-uniform, meaning there exists a constant cqu ≥ 1 such that: hZ,Ω ≤ cqu · qZ .
(76)
This prevents the clustering of sensors that leads to the ill-conditioning of the Gram matrix in Moving Least Squares (MLS) reconstruction. Such configurations can be generated via Poisson Disk Sampling or low-discrepancy sequences. 29
QI
Definition 28 (Lebesgue Stability) Let R : Rm → L∞ (S × A) be the linear reconstruction operator mapping sensor values to functions. The Lebesgue constant Λm is the operator norm induced by the uniform metric: Λm := sup ∥R(v)∥∞ . (77) v∈Rm ,∥v∥∞ ≤1
Under Assumption 27, for Moving Least Squares with appropriate weight functions, Λm is uniformly bounded (or grows at most logarithmically with m) (Wendland, 2004). We denote this bound by ΛM LS . C.3. Error Recurrence Analysis We analyze the propagation of error through the iterative approximate Value Iteration scheme. The bound requires a careful decomposition into consistency (approximation power) and stability (error amplification). Let Q∗(k) denote the true solution at time step tk , and let Q̂(k) denote the Neural Operator approximation. We define the uniform error at step k as Ek := ∥Q̂(k) − Q∗(k) ∥∞ . Lemma 29 (Stable Error Recurrence) Grant Assumptions 3, 23, and 24. Let LT ≤ eρstab δ be a one-step L∞ stability bound for the Bellman propagation composed with the reconstruction pipeline. In the ideal Markov expectation and maximization operator without reconstruction amplification, one may take ρstab = 0 and LT = 1. The error at the next iteration satisfies the recurrence: Ek+1 ≤ LT Ek + εnet + LT Capprox KS (δ) h2S + LA hA . (78) | {z } Local Consistency Error ηloc
Proof We decompose the error Ek+1 = ∥Nθ (Pm Q̂(k) ) − Tδ Q∗(k) ∥∞ using the triangle inequality: Ek+1 ≤ ∥Nθ (Pm Q̂(k) ) − Tδ (RPm Q̂(k) )∥∞ | {z } (I) Network Approx.
+ ∥Tδ (RPm Q̂(k) ) − Tδ (RPm Q∗(k) )∥∞ | {z } (II) Stability
+ ∥Tδ (RPm Q |
∗(k)
) − Tδ (Q∗(k) )∥∞ . {z }
(III) Consistency
(i) Term (I): By Assumption 24, the network fits the composed operator within tolerance εnet . (ii) Term (II): This term captures the amplification of previous errors. (II) ≤ LT ∥RPm Q̂(k) − RPm Q∗(k) ∥∞ = LT ∥R Pm (Q̂(k) − Q∗(k) ) ∥∞ ≤ LT · Λm · ∥Pm (Q̂(k) − Q∗(k) )∥∞ ≤ LT · 1 · Ek . Crucially, because we restricted R to a Monotone Averager (Assumption 23), Λm = 1 exactly. This prevents the unconditional instability that plagues high-order collocation methods. 30
D EEP Q-L EARNING ON H OLDER S PACES
(iii) Term (III): This term represents the interpolation error on the true function. By Corollary 14, Q∗(k) lies in the anisotropic Hölder space. Using the capped second-order bound: (III) ≤ LT ∥RPm Q∗(k) − Q∗(k) ∥∞ ≤ LT Capprox KS (δ) h2S + LA hA . Summing these yields the result. Remark 30 (Global Stability and Sample Curation) For the global error EN to remain bounded over the horizon T = N δ, the amplification factor (LT Λm )N must not explode. Since LT ≤ eρstab δ , the physical propagation contributes at most eρstab T over the horizon, and equals 1 in the non-expansive Markov case. However, the Lebesgue constant Λm acts as a multiplier at every step. If Λm > 1, the error grows as (Λm )T /δ . This highlights a strict requirement for the sensor set: we must employ Sample Curation (or regularization) to ensure Λm ≈ 1. The Quasi-Uniform construction (Lemma 21) is a necessary condition for keeping Λm controlled, whereas random sampling guarantees instability (Λm ≫ 1) for large m. C.4. The Residual Limit and Conservation of Difficulty The stiffness scaling in Theorem 15 seemingly contradicts the intuition that the continuous-time limit (the HJB PDE) is a well-posed problem involving smooth generators. This apparent paradox is resolved by analyzing the Residual Operator formulation proposed in Appendix A.4. Recall the normalized residual operator acting on Q ∈ Hα,Lip : Rδ [Q](x, a) :=
Tδ Q(x, a) − Q(x, a) . δ
(79)
In the residual learning framework, the network Nres targets Rδ [Q] rather than Tδ Q. W ELL -C ONDITIONED TARGET VS . I LL -P OSED I NPUTS For time-dependent coefficients, the local parabolic evolution family gives the expansion U a (t + δ, t)Q − Q → Lat Q + r(t, ·, a) δ
(80)
for Q in the generator domain. For the Bellman optimality operator, this Hamiltonian limit is valid only after controlling the maximization gap VQ (x) − Q(x, a). In particular, on a zero-gap or neargreedy set satisfying VQ (x)−Q(x, a) = o(δ), the target function for the residual network converges to H(Q) = Lat Q + r, which is O(1) with respect to δ. Off that set, the normalized gap contributes (VQ − Q)/δ and can dominate the residual. However, the complexity does not vanish; it is merely shifted from width (W ) to sensor resolution (m). To approximate the generator La Q(x) = 12 Tr(ΣD2 Q(x)) + b⊤ ∇Q(x), the Branch Net must implicitly compute the Hessian D2 Q from the discrete sensor inputs Pm (Q). Proposition 31 (Conservation of Difficulty in the Residual Limit) Let Nres be the branch-trunk network approximating the residual operator Rδ = 1δ (Tδ − I) on a zero-gap or near-greedy subset where VQ (x) − Q(x, a) = o(δ). Let hS be the fill-distance (covering radius) of the quasi-uniform sensor set Zm ⊂ S, and let ϵenc be the precision of the input encoding. 31
QI
For the network to consistently approximate the infinitesimal generator La as δ → 0 (i.e., to distinguish the dynamics from the identity map), the sensor resolution must satisfy the Diffusion Nyquist Condition: √ hS ≲ νδ. (81) Consequently, the required number of sensors m scales as: m ∼ O (νδ)−d/2 .
(82)
Violating this condition results in an O(1) aliasing error where the diffusion operator is numerically indistinguishable from the Identity on the grid Zm . Proof We analyze the error in identifying the operator La from the discrete samples of the mapped function. We decompose the total approximation error into Temporal Truncation, Spatial Aliasing, and Conditioning errors. (i) Temporal Truncation (The Semigroup Expansion). Let U a (t+δ, t) be the local parabolic evolution family generated by the time-dependent operator Lat . On the near-greedy set where the maximization gap is o(δ), the Bellman update is locally U a (t + δ, t)Q + δr + o(δ). The residual is therefore: U a (t + δ, t) − I Rδ Q = Q + r + o(1). (83) δ A local Taylor expansion of the parabolic evolution family for Q ∈ Dom((Lat )2 ) gives U a (t + δ, t)Q = Q + δLat Q + O(δ 2 ). Thus, the target consistency error is: ∥Rδ Q − (Lat Q + r)∥∞ ≤ C∥Lat Lat Q∥∞ · δ + o(1).
(84)
This term vanishes as δ → 0, suggesting the problem becomes easier. This is the source of the Residual intuition. (ii) Spatial Aliasing (The Spectral Barrier). However, the network does not see Q; it sees Pm Q. It must infer the action of T (δ) from these samples. Consider the simplified diffusion case where La = ν2 ∆. In the Fourier domain, the operator T (δ) acts as a multiplier on frequency ξ ∈ Rd : ν 2 b T\ (δ)Q(ξ) = e− 2 δ∥ξ∥ Q(ξ).
(85)
This is a low-pass filter. The signal of the physics (the change in values) corresponds to the attenuν 2 ation 1 − e− 2 δ∥ξ∥ . The sensor grid Zm with spacing hS supports spatial √ frequencies only up to the Nyquist limit π ξmax ≈ hS . If the grid is too coarse such that hS ≫ νδ, then for all visible frequencies ∥ξ∥ ≤ ξmax : ν νδπ 2 δ∥ξ∥2 ≤ ≪ 1. (86) 2 2h2S ν
2
Consequently, e− 2 δ∥ξ∥ ≈ 1 − O(δ/h2S ). The relative change in the function value on the grid is negligible compared to the observation noise or approximation error. The discrete operator acts as the Identity: Pm Tδ R ≈ Pm R. 32
D EEP Q-L EARNING ON H OLDER S PACES
To resolve the curvature (the Laplacian), the grid must capture the transition where the spectrum rolls off. This requires the Nyquist frequency to exceed the Gaussian width: √ π 1 (87) ≳√ =⇒ hS ≲ π νδ. hS νδ (iii) Conditioning (Numerical Differentiation). Even if the grid satisfies the Nyquist condition, we face the ill-conditioning of inverse problems. The Branch Net implicitly computes the discrete Laplacian ∆h . For inputs with encoding error ϵenc (from the previous step’s network approximation), the error in the discrete Laplacian scales as: ϵenc Error(∆h ) ∼ C1 h2S + . (88) C2 2 | {z } hS {z } | Discretization Noise Amplification
1/4
Balancing these terms requires hS ∼ ϵenc . However, √ combining with condition (ii), we require ϵenc /h2S ≪ 1 to get a stable signal. Substituting hS ∼ δ, the noise amplification factor becomes ϵenc /δ. This implies that to learn the generator with time-step δ, the previous value function must be approximated with precision ϵenc ≪ δ. Combining (ii) and the definition of the sensor count m ≈ Vol(S)h−d S , we obtain: m ≳ Vol(S)(νδ)−d/2 .
(89)
This proves that the complexity reduction in W (width) achieved by the residual formulation is perfectly conserved by the explosion in m (sample complexity). Remark 32 (Implication for Replay Buffers) This proposition rigorously explains the empirical instability of continuous-time RL (small δ) when used with fixed-size replay √ buffers. As δ → 0, a fixed buffer becomes effectively sparse relative to the diffusion length δ. The gradients of the Bellman loss becomes dominated by the high-frequency aliasing noise rather than the physical drift/diffusion signal, leading to value function collapse unless the buffer density m increases dynamically with 1/δ. Thus, Theorem 15 correctly captures the fundamental hardness of the problem. On the zero-gap or near-greedy subset, the Residual architecture allows the weights to converge to a stable generatorlike limit, but the input data requirements (m) still explode as δ −d/2 to avoid aliasing the diffusion process. The residual framework is therefore a method for conditioning optimization, not for circumventing the information-theoretic lower bounds of the parabolic operator.
Appendix D. Discussion Our main result, Theorem 15, establishes a rigorous upper bound on the complexity of approximating Bellman value-iteration iterates, and therefore on the complexity of the idealized Q-learning targets in the present diffusion regime, by identifying the stiffness scaling factor γ(α) which penalizes the network capacity as the time-step δ → 0. In this concluding section, we refine the physical interpretation of these bounds, addressing the limit behavior of the residual operator, the implications of semi-convexity, and the practical realizability of the proposed architecture. 33
QI
D.1. The Residual Limit and Generator Consistency While the direct analysis of Tδ reveals a stiffness explosion, modern Deep RL typically employs Residual Learning. We define the normalized residual operator Rδ acting on a function Q ∈ Hα,Lip as: Tδ Q(x, a) − Q(x, a) . (90) Rδ [Q](x, a) := δ When the DeepONet is trained via the skip-connection ansatz Q̂out = Qin + δ · Nres (Qin ), the network Nres effectively learns Rδ [Q]. Unlike the direct diffusion-evolution operator, which converges to the Identity (a trivial but uninformative limit), the residual operator admits a non-trivial strong limit only after the Bellman maximization gap is controlled. For any Q(·, a) in the domain of the local generator and any sequence of active or near-greedy actions satisfying VQ (x) − Q(x, a) = o(δ): lim Rδ [Q](x, a) = Lat Q(x, a) + r(t, x, a) . {z } |
δ→0
(91)
Hamiltonian
For arbitrary off-policy actions, the additional normalized gap (VQ − Q)/δ may diverge. Thus the residual formulation should be interpreted as a generator-consistency statement on the active/nearactive set, with the remaining difficulty appearing in data resolution as discussed in Section D.3. D.2. Beyond Lipschitz: Semi-Convexity and Regularity In Section 3, we derived smoothing estimates based on the assumption that the value function VΦ (x) = supa Φ(x, a) is merely Lipschitz continuous. While rigorous, this bound is pessimistic. Remark 33 (Semi-Convexity of the Value Function) If the input Q-function Φ(·, a) is C 2 smooth with uniformly bounded Hessian, the upper envelope VΦ is semi-convex. That is, there exists a constant λ ≥ 0 such that x 7→ VΦ (x)+ λ2 ∥x∥2 is convex. Semi-convex functions possess a distributional lower bound on their Hessian: D2 VΦ ≥ −λI. Since the parabolic operator ∂t − Tr(σσ ⊤ D2 ) is order-preserving, the solution to the Bellman equation preserves this semi-convexity. Consequently, the singularities in the second derivative are one-sided (related to the downward kinks of the maxoperator). Standard viscosity solution theory suggests that parabolic smoothing acts more efficiently 3+α on semi-convex data than on general Lipschitz data. While our derived exponent γ(α) = 2(2+α) relies on the conservative Lipschitz embedding, exploiting semi-convexity could theoretically reduce this further, potentially approaching the optimal heat kernel scaling of δ −1/2 . D.3. Conservation of Difficulty: The Sensor Resolution Condition To approximate the generator target La Q ≈ 12 Tr(ΣD2 Q) + b⊤ ∇Q, the Branch Net must implicitly compute second-order derivatives of the input function from the discrete sensor readings {Q(zj )}. This is an inverse problem. Let hS be the sensor fill-distance. By standard numerical analysis, the error in estimating the Laplacian behaves as O(h2S ) + O(η/h2S ), where η is the observation noise (or function approximation error). √ Crucially, the physical diffusion process operates on a spatial scale of νδ. If the sensor spacing √ hS is coarser than the diffusion length scale (hS ≫ δ), the discrete observations cannot resolve the local curvature driving the dynamics. 34
D EEP Q-L EARNING ON H OLDER S PACES
Remark 34 (Necessary Resolution Condition) For the Neural Operator to consistently approximate the generator Lat in the residual framework as δ → 0 on the active or near-greedy set, the sensor fill-distance must satisfy the stricter condition: √ (92) hS = o( δ). Failure to satisfy this leads to aliasing of the high-frequency modes of the Q-function, resulting in a biased estimate of the drift and diffusion terms. This explains the empirical instability of continuoustime RL algorithms when applied with fixed, coarse replay buffers. D.4. Realizability of the Branch Net Our analysis assumes the Branch Net β can approximate the anisotropic Moving Least Squares (MLS) operator (Assumption 23). The MLS shape functions are rational functions of the input coordinates (ratios of polynomials involving the kernel weights). Recent results in approximation theory (Telgarsky, 2017; Yarotsky, 2017) establish that deep ReLU networks can approximate rational functions with exponential convergence w.r.t. depth. Specifically, to approximate a rational function with error ε, a depth of O(log(1/ε)) suffices. Thus, the assumption that standard MLPs can serve as the Branch Net encoder for the proposed architecture is well-founded, provided the network is sufficiently deep to resolve the algebraic structure of the reconstruction operator. D.5. Limitations: The Persistence of the Action Curse Finally, we must clarify the scope of the Curse of Dimensionality mitigation. Theorem 15 demonstrates that the parabolic smoothing reduces the effective dimension of the state space S from d to roughly d/3. However, the Bellman operator provides no smoothing in the action space A. The derived complexity bound: d W ∼ O δ −γd ε− 2+α + (δε)−da , (93) retains an exponential dependence on the action dimension da . For high-dimensional control problems (da ≫ 1), the second term dominates. Our anisotropic DeepONet architecture prevents the state dimension from contaminating the action complexity, but it cannot remove the intrinsic hardness of optimizing non-concave functions over high-dimensional action spaces. This remains a fundamental limitation of value-based methods in continuous control, suggesting that for large da , policy-gradient methods (which avoid explicit maximization over A) may be structurally necessary.
35