Kinetic Handover Modeling in Nonstationary Cellular Networks with Fluid Spatiotemporal Stochastic Geometry

1

Kinetic Handover Modeling in Nonstationary Cellular Networks with Fluid Spatiotemporal Stochastic Geometry

arXiv:2610.34211v1 [cs.NI] 25 Oct 2026

Wen-Yu Dong, Member, IEEE, Zheng Lin, Jinfeng Xu, Yuxin Zhang, Zihan Fang Qi Bi, Fellow, IEEE, and Sheng Chen, Life Fellow, IEEE

Abstract—Inferring directional handover traffic from aggregate mobility observations remains challenging because population density evolution does not fully describe the movement of individual users across cell boundaries. This paper develops a kinetic framework for geometry-induced serving-region crossings in nonstationary cellular networks, connecting boundarynormal velocity statistics with fluid spatiotemporal stochastic geometry (F-STSG) for density-based transport reconstruction. The framework distinguishes signed net population transfer from gross bidirectional crossings and reveals two information limits: density evolution does not uniquely determine physical pairwise transport, while first-order transport alone cannot resolve counter-moving crossings. We establish exact recovery conditions for potential transport and develop coherent and Gaussian mobility closures for estimating directional crossing rates. Sharp moment-constrained bounds further characterize the remaining uncertainty under bounded normal velocities without assuming a specific velocity distribution. An error bound links directional estimation accuracy to transport reconstruction and unresolved mobility information. The formulation also recovers the classical stationary Poisson–Voronoi handover law. Numerical experiments examine the information limits, reconstruction sensitivity, and closure accuracy. On 47 held-out vehicular simulation windows, F-STSG achieves lower net-flow reconstruction errors than three baselines, while auxiliary velocity statistics improve directional estimates under matched Gaussian mobility. Index Terms—Fluid spatiotemporal stochastic geometry, handover modeling, mobility analysis, nonstationary cellular networks, stochastic geometry.

I. I NTRODUCTION Mobility management is a core function of cellular networks because changes in serving-cell association affect signaling overhead, service continuity, and radio resource utilization [1, 2]. Network densification, heterogeneous cell sizes, and highfrequency operation further increase the frequency and complexity of mobility-related association changes [3, 4]. HanW.-Y. Dong and Q. Bi are with Future Technology Research Center, China Telecom Research Institute, Beijing 102209, China (E-mails: dongwy@chinatelecom.cn; qibi@chinatelecom.cn) Z. Lin is with the Interdisciplinary Centre for Security, Reliability and Trust (SnT), University of Luxembourg, L-1855 Luxembourg, Luxembourg (E-mail: zhenglin@ieee.org). J. Xu is with the Department of Electrical and Computer Engineering, The University of British Columbia, Vancouver, Canada (E-mail: jinfeng.xu@ubc.ca). Y. Zhang is with the Institute of Space Internet, Fudan University, Shanghai 200433, China (E-mail: yxzhang24@m.fudan.edu.cn). Z. Fang is with the Hong Kong JC STEM Lab of Smart City and Department of Computer Science, City University of Hong Kong, Kowloon, Hong Kong, China (E-mail: zihanfang3-c@my.cityu.edu.hk). S. Chen is with School of Electronics and Computer Science, University of Southampton, Southampton SO17 1BJ, U.K. (E-mail: sqc@ecs.soton.ac.uk).

dover (HO) traffic can consequently vary substantially across neighboring cells and over time, particularly in the presence of collective mobility associated with commuting corridors, transportation hubs, and large public events. Characterizing when, where, and in which direction such association changes arise is therefore relevant to mobility-aware network dimensioning and control. At the trajectory-prediction level, diffusion-based methods have also been developed to forecast individual agent motion [5, 6]. Such methods address user-level motion prediction, whereas the present work characterizes directional cellboundary crossings from aggregate population observations. The modeling difficulty arises from the information lost between individual boundary crossings and aggregate mobility observations. At the microscopic level, a geometryinduced serving-region transition occurs when an individual user equipment (UE) crosses a boundary between neighboring serving regions. At the macroscopic level, however, a population density field retains neither individual trajectories nor the complete local velocity distribution. This leads to two distinct information losses. First, the evolution of the population density constrains the divergence of the physical population flux but does not, in general, uniquely determine the physical transfer between individual cell pairs. Second, even if the physical first-order flux is known, counter-moving users can cancel in the signed flux while still contributing to the total number of boundary crossings. A populationlevel HO model must therefore distinguish signed net transfer from gross crossing traffic and, separately, distinguish the physical mobility flux from a transport field reconstructed from aggregate observations. This paper develops a kinetic characterization of mobilityinduced serving-region crossings. For each neighboring cell pair, the directional crossing rate is defined from the local population density and the positive part of the boundarynormal UE velocity. This representation separates the signed net transfer between two cells from their gross bidirectional crossing traffic. The former is determined by the boundarynormal component of the physical first-order population flux, whereas the latter generally depends on additional local velocity information. The resulting formulation makes explicit which HO quantities are determined by first-order population transport and which remain unresolved when only aggregate motion is available. Fluid spatiotemporal stochastic geometry (F-STSG) [7] is then used to connect this characterization with density-based observation. From successive spatial density snapshots, F-

STSG constructs a canonical first-order transport field that is consistent with the observed density evolution through a minimum-energy formulation. This reconstructed field is not assumed to coincide with the physical flux in general. Rather, the paper identifies conditions under which exact recovery is possible and quantifies how transport-reconstruction error propagates to directional crossing estimates. When counterflow effects are relevant, additional local velocity statistics are introduced through a dispersion-aware closure. The baseline model concerns geometry-induced serving-region crossings. Radio measurements, hysteresis, time-to-trigger, target selection, and other protocol mechanisms governing executed HOs are outside the crossing model itself.

[33]. More recent countrywide measurements have characterized the spatial and temporal distributions of horizontal and vertical HOs and examined their relation to population distribution and network conditions [34, 38]. Such measurements provide direct evidence of how association changes appear in operational networks and complement analytical models based on idealized mobility assumptions. In these studies, however, the HO event is already part of the observation. They do not determine how directional cell-pair crossing rates can be obtained from an evolving population field before individual HO events are observed, nor which components of the resulting crossing traffic are identifiable from aggregate population information. F-STSG [7] starts from a different state description, representing a nonstationary network by a time-varying spatial density and a macroscopic transport field linked through a continuity equation. This field representation has subsequently been applied to flux-aware UAV infrastructure provisioning [35], vorticity-based routing in ultra-dense networks [36], and fluid-dynamic interference modeling for low-Earth-orbit (LEO) mega-constellations [37]. These studies describe network dynamics at the field level but do not formulate directional mobility-induced crossings between neighboring serving regions. Density evolution alone does not specify the physical population transfer between individual cell pairs, while firstorder population transport does not, in general, determine the total crossing traffic generated by counter-moving users. The present work addresses these two distinctions by connecting microscopic boundary crossings with macroscopic population transport and by establishing when directional HO traffic can be reconstructed from successive density observations.

A. Related Work Analytical studies of mobility in cellular networks commonly use stochastic geometry to relate user motion, association geometry, and HO statistics. For hexagonal and Poisson–Voronoi deployments, classical results established the dependence of the mean HO rate on UE speed and base-station (BS) density and derived corresponding cell residence and sojourn-time statistics [8]. Subsequent extensions considered biased multi-tier association, general UE trajectories, network densification, and control-plane/user-plane separation [9– 12]. Related analyses addressed movement-aware cooperation and mobility in ultra-dense and radio-frequency/visible-lightcommunication (RF/VLC) networks [13–15], with a broader overview of mobility-aware stochastic geometry given in [3]. The methodology has also been applied beyond terrestrial cellular networks, including satellite–aerial, cooperative satellite–unmanned-aerial-vehicle (UAV), and heterogeneous non-terrestrial systems [16–22]. Despite the different network geometries considered, mobility in these models is generally specified through individual UE trajectories or prescribed mobility distributions rather than through a population field evolving jointly in space and time. Beyond association geometry, HO models have incorporated radio propagation and mobility-management mechanisms. Multi-connectivity, blockage, beam management, and directional transmission have been considered in millimeterwave, terahertz, and vehicular networks [23–26], while more detailed association models have been used to characterize HO trigger regions, failures, and ping-pong behavior [27]. At the protocol and control level, previous studies have optimized hysteresis and time-to-trigger parameters, incorporated signal prediction and target-cell selection, and applied reinforcement learning and open radio access network (O-RAN) methods to mobility management [28–32]. These models describe how radio measurements, association rules, and protocol parameters affect HO decisions and execution. Their primary state variables, however, remain user-level measurements, trajectories, or serving-cell decisions rather than the aggregate transport of a mobile population. Large-scale measurement studies approach the problem from observed HO records rather than assumed mobility models. Cellular HO traces have been used to infer urban mobility patterns and relate inter-cell transitions to population activity

B. Contributions and Organization The main contributions are summarized as follows. • We develop a kinetic formulation of directional servingregion crossings based on boundary-normal UE velocity moments. It distinguishes signed net population transfer from gross bidirectional crossings, preserves cell-level population balance, and recovers the classical stationary Poisson–Voronoi crossing-rate law as a special case. • We identify two information limits: density evolution does not uniquely determine physical pairwise transport, and first-order transport does not generally determine gross crossings. We establish exact recovery conditions for potential transport and derive sharp distributionfree bounds on directional crossing rates given the first two normal-velocity moments and a bounded-speed constraint. Coherent and Gaussian mobility closures provide conditional estimates under additional mobility information. • We connect the kinetic formulation with densitysnapshot-based F-STSG reconstruction and establish an error bound accounting for reconstruction error, physicalflux ambiguity, auxiliary-statistic error, and closure mismatch. Controlled experiments examine these effects, while held-out vehicular simulations compare net-flow reconstruction against three baselines and evaluate the role of auxiliary velocity statistics. 2

TABLE I C OMPARISON OF REPRESENTATIVE APPROACHES WITH RESPECT TO THE MOBILITY AND HO MODELING ASPECTS CONSIDERED IN THIS PAPER . Directed cell-pair mobility

Time-varying population field

Net/gross distinction

Density-to-HO reconstruction

Stochastic-geometry mobility and HO modeling [8–15] Radio- and beam-aware HO modeling [23–25, 27] HO control, prediction, and learning [28–31] Empirical HO measurement and mobility analysis [33, 34] F-STSG framework [7] Field-based F-STSG applications [35–37]

△ △ ✓ ✓ – –

– – – △ ✓ △

– – – – – –

– – – – – –

This work

✓

✓

✓

✓

Representative approaches

Note: HO: handover; F-STSG: fluid spatiotemporal stochastic geometry. ✓: explicitly treated; △: partially or indirectly treated; –: not explicitly addressed.

The operators [·]+ , E[·], ∇, and ∇· denote the positive part, expectation, gradient, and divergence, respectively. A hat denotes a reconstructed or estimated quantity, and the corresponding unhatted symbol denotes the reference quantity unless stated otherwise. The principal symbols used throughout the paper are summarized in Table II, while symbols used only in individual derivations are defined locally when introduced.

Table I compares the proposed approach with the closely related studies in terms of mobility and HO modeling aspects. The remainder of this paper is organized as follows. Section II introduces the system model and directional crossing definitions. Section III establishes the kinetic relations, information limits, recovery conditions, and mobility closures. Section IV develops density-based reconstruction and directional crossing estimation. Section V presents numerical results, and Section VI concludes the paper.

II. S YSTEM M ODEL AND P ROBLEM F ORMULATION

TABLE II K EY NOTATION . Symbol

Definition

Ω, Vi Γij , nij λ(x, t) Ni (t) J(x, t) b J(x, t) Uij Hij (t) Fij (t) Gij (t) mij , s2ij U GL ij , Gij λk , qk bk J bij,k , H b ij,k F

Service region and serving region of BS i Common cell boundary and unit normal from i to j Spatial UE population density UE population in cell i Physical first-order population flux Canonical flux reconstructed from exact density evolution UE velocity component normal to Γij Directional crossing rate from cell i to j Signed net flow, Hij − Hji Gross crossing rate, Hij + Hji Conditional mean and variance of normal velocity Lower and upper bounds on gross crossing rate Density snapshot and corrected temporal derivative Flux reconstructed from density snapshots Estimated signed net flow and directional crossing rate

We consider mobility-induced association changes in a cellular network from two complementary descriptions. The physical population flux characterizes the first-order transport generated by UE motion, whereas F-STSG reconstructs a canonical transport field from aggregate density evolution. Directional crossings of neighboring serving regions are then defined from the local boundary-normal velocity distribution. Fig. 1 summarizes the resulting model: physical UE motion generates serving-region crossings, while successive density observations are used to reconstruct first-order transport and, under an appropriate local mobility closure, estimate the corresponding crossing rates. A. Cellular Geometry and Population Transport Consider a bounded planar service region Ω ⊂ R2 containing M fixed BSs at distinct locations {bi }i∈I , where I = {1, . . . , M }. To isolate mobility-induced association changes from radio fluctuations and protocol filtering, the

C. Notation Scalars are denoted by italic letters, vectors by bold lowercase letters, and sets or spatial regions by calligraphic letters.

Fig. 1. System model for mobility-induced serving-region crossings and density-snapshot-based reconstruction. Physical population transport generates directional crossings between neighboring serving regions, whereas aggregate density observations are used to reconstruct first-order transport and estimate the corresponding crossing rates under a local mobility closure.

3

baseline model adopts nearest-BS association. The serving region of BS i is Vi = {x ∈ Ω : ∥x − bi ∥2 ≤ ∥x − bj ∥2 , ∀j ∈ I} .

individual cell boundary.

(1)

B. Boundary-Crossing Model

The regions {Vi }i∈I form a partition of Ω up to their boundaries. Two cells are defined as neighbors when their common boundary has nonzero length. For neighboring cells i and j, let Γij = ∂Vi ∩ ∂Vj ∩ Ω denote the common boundary and let nij denote the unit normal directed from Vi to Vj . Hence, Γij = Γji and nji = −nij . The neighbor set of cell i is denoted by Ni . Let λ(x, t) ≥ 0 denote the spatial UE population density at location x and time t. The UE population in cell i, Ni (t), the local conditional mean velocity, v(x, t), and the corresponding physical population flux, J(x, t), are given respectively: Z Ni (t) = λ(x, t) dx, (2)

Under nearest-BS association, a continuous UE trajectory changes serving region when it crosses the common boundary of two neighboring cells. We refer to this event as a geometryinduced serving-region crossing and denote its directional rate from cell i to cell j by Hij (t). This rate characterizes mobility-induced association crossings rather than protocollevel HO execution, which may additionally depend on radio measurements and mobility-management parameters. For neighboring cells i and j, define the boundary-normal component of the local UE velocity as Uij = UT nij .

A positive value of Uij represents motion from Vi toward Vj , whereas a negative value represents motion in the reverse direction. Assuming continuous trajectories and an integrable first normal-velocity moment along Γij , the directional crossing rate is Z Hij (t) = λ(x, t)E [[Uij ]+ | x, t] dℓ, [a]+ = max{a, 0}, (7)

Vi

v(x, t) = E [U | x, t] ,

J(x, t) = λ(x, t)v(x, t),

(6)

(3)

where U denotes the velocity of a randomly selected UE conditioned on location and time. The conditional velocity distribution is assumed to have a finite first moment. Thus, J retains the first-order population transport induced by the conditional mean velocity but not the complete distribution of individual UE velocities. The baseline model considers population redistribution within Ω without physical arrivals or departures through its outer boundary. Assuming sufficient regularity for the differential and boundary operations, the density and physical flux satisfy ∂t λ(x, t) + ∇· J(x, t) = 0, x ∈ Ω, (4) J(x, t) · nΩ (x) = 0, x ∈ ∂Ω. where nΩ denotes the outward unit normal on ∂Ω. The noflux boundary condition implies conservation of the total UE population in Ω. Open boundaries or physical arrivals and departures can instead be represented through boundary fluxes or source terms but are not included in the baseline model. Equation (4) determines the divergence of the physical flux but not the complete field J. F-STSG [39] selects R a canonical flux by minimizing the weighted kinetic energy Ω ∥j∥2 /λ dx over admissible flux fields j satisfying the continuity equation and no-flux boundary condition. Under the regularity conditions of the weighted elliptic problem, the potential ϕ(x, t) b t) satisfy and reconstructed flux J(x, ∂λ(x, t) , x ∈ Ω, ∇· (λ(x, t)∇ϕ(x, t)) = ∂t (5) λ(x, t)∇ϕ(x, t) · nΩ (x) = 0, x ∈ ∂Ω, b t) = −λ(x, t)∇ϕ(x, t). J(x,

Γij

where dℓ denotes the line element along Γij . Reversing the boundary normal gives the corresponding rate Hji . Equation (7) has units of crossings per unit time and depends only on the velocity component normal to the serving-region boundary. Tangential motion does not contribute to the crossing rate. Intersections shared by three or more Voronoi cells have zero line measure and do not contribute to the boundary integral. The problem considered in this paper is to determine which directional crossing quantities can be characterized or estimated from aggregate population observations without reconstructing individual UE trajectories. The physical firstorder flux J contains only the conditional mean velocity in (3), whereas (7) depends on a one-sided moment of the boundary-normal velocity distribution. Moreover, density evolution does not, in general, uniquely determine the physical pairwise transport represented by J. The next section therefore separates the information contained in density evolution, firstorder population transport, and unresolved local mobility, and identifies the conditions under which directional crossing rates can be recovered. III. K INETIC H ANDOVER C HARACTERIZATION The directional crossing rate in (7) depends on both first-order population transport and the local distribution of boundary-normal UE velocities. This section first relates directional crossings to the physical population flux and then identifies two distinct information limitations. The first concerns the recovery of physical pairwise transport from density evolution, whereas the second concerns the recovery of gross crossing traffic from first-order transport. Conditions under which these ambiguities can be resolved are then established, followed by a stationary isotropic limit that recovers the classical Poisson– Voronoi HO result.

For the closed system, the Neumann compatibility condition follows from conservation of the total population. The potential is determined up to an additive constant, which does not afb Throughout the paper, J denotes the physical first-order fect J. b denotes the F-STSG reconstruction population flux, whereas J inferred from density evolution. The two quantities need not coincide in general, even though both are consistent with the same population balance. In particular, density evolution alone does not determine the physical transport across each 4

Then J + C produces the same density evolution and outer boundary flux as J, while its flow across Γij changes by Z C(x, t) · nij dℓ. (16)

A. Net and Gross Crossing Flows The two directional crossing rates associated with a common cell boundary contain different information about population transport. Define the total outgoing and incoming crossing rates of cell i respectively as X X Hiout (t) ≜ Hij (t), Hiin (t) ≜ Hji (t). (8) j∈Ni

Γij

b satisfy the same continuity equation Proof. Both J and J and outer no-flux condition. Integrating their divergences over Vi and applying the divergence theorem gives (14). A field satisfying (15) changes neither the divergence of the transport field nor its normal component on the outer boundary. Hence, provided that J + C remains admissible, it produces the same density evolution. Its contribution to the flow across an individual common boundary is exactly (16).

j∈Ni

Proposition 1. For any neighboring cells i and j satisfying the assumptions of (7), define the signed net flow Fij (t) and gross crossing rate Gij (t) by Z Fij (t) ≜ Hij (t)−Hji (t) = J(x, t) · nij dℓ, ZΓij (9) Gij (t) ≜ Hij (t)+Hji (t) = λ(x, t)E [|Uij | | x, t] dℓ.

Proposition 2 separates cell population balance from physical pairwise transport. Density evolution fixes the former, but divergence-free transport components can alter the latter without changing the observed density evolution. The pairwise transport ambiguity is resolved for the following potential-transport class.

Γij

Then Z |J(x, t) · nij | dℓ ≥ |Fij (t)|,

Gij (t) ≥

(10)

Γij

and the cell population satisfies X dNi (t) = Hiin (t) − Hiout (t) = − Fij (t). dt

(11)

Corollary 1. Assume that the physical population flux admits the representation

j∈Ni

Proof. Because nji = −nij , we have Uji = −Uij . The identities [a]+ − [−a]+ = a and [a]+ + [−a]+ = |a|, together with J = λE[U | x, t], give (9). The pointwise inequality |E[Uij | x, t]| ≤ E[|Uij | | x, t] gives the first inequality in (10), and the triangle inequality gives the second. Integrating (4) over Vi and applying the divergence theorem, with the outer no-flux condition where ∂Vi intersects ∂Ω, yields (11).

J(x, t) = −λ(x, t)∇ψ(x, t),

(17)

satisfies the same outer no-flux condition as (5), and assume that the corresponding weighted Neumann problem admits a solution that is unique up to an additive constant. Then b t) = J(x, t). J(x, (18) Consequently, Fbij (t) = Fij (t) for every neighboring cell pair.

The quantity Fij is the signed net population transfer from cell i to cell j, whereas Gij is the gross crossing rate across their common boundary. By construction, Fji = −Fij and Gji = Gij , and Gij − Fij Gij + Fij , Hji = . (12) Hij = 2 2 Equality in the first inequality of (10) holds when the conditional normal velocity has one sign almost surely at almost every point of Γij . Equality in the second additionally requires the mean normal flux to have a common sign along the boundary. When only the density evolution is available, the physical flux J is not directly observed. Let the reconstructed pairwise flow associated with (5) be Z b t) · nij dℓ. Fbij (t) ≜ J(x, (13)

Proof. Substituting (17) into (4) gives ∂λ(x, t) ∇· [λ(x, t)∇ψ(x, t)] = . ∂t The physical no-flux condition also gives λ∇ψ·nΩ = 0 on ∂Ω. Hence, ψ and ϕ satisfy the same weighted Neumann problem. Uniqueness up to an additive constant gives ∇ψ = ∇ϕ, and b therefore J = J. The preceding result identifies a class for which the F-STSG reconstruction coincides with the physical first-order transport. Outside this class, agreement of cell population balance does not imply agreement of individual pairwise flows. A second information limitation remains even when the physical first-order population flux is known.

Γij

Proposition 3. Suppose that Γij contains a subset of positive line measure on which λ(x, t) > 0. Then the local density λ and physical first-order population flux J do not, in general, uniquely determine the gross crossing rate Gij .

Proposition 2. For a given density evolution and the same outer no-flux condition, X X dNi (t) Fbij (t) = Fij (t) = − . (14) dt j∈Ni

j∈Ni

Proof. Consider a measurable subset of Γij with positive line measure and positive density. On this subset, let two conditional normal velocity distributions be U1 ∼ Unif[−a, a] and U2 ∼ Unif[−b, b], where a, b > 0 and a ̸= b, while keeping the tangential conditional mean and all mobility statistics outside the subset identical. Both normal distributions

The individual physical pairwise flows are not uniquely determined by the density evolution in general. In particular, let C(x, t) be a sufficiently regular field such that J + C remains an admissible population flux and ∇· C(x, t) = 0,

C(x, t) · nΩ (x) = 0,

x ∈ ∂Ω. (15) 5

Let φ and Φ denote the probability density function and cumulative distribution function of a standard normal random variable. For sij > 0, define mij . (23) zij ≜ sij The one-sided Gaussian moments give Z Hij (t) = λ (sij φ(zij ) + mij Φ(zij )) dℓ, Γij Z (24) Gij (t) = λ (2sij φ(zij ) + mij (2Φ(zij ) − 1)) dℓ.

have zero mean, and they produce the same conditional mean velocity and hence the same physical first-order flux. However, b a E|U2 | = . E|U1 | = , 2 2 Their contributions to the integral defining Gij are therefore different on a set of positive line measure. Thus, the same λ and J can correspond to different gross crossing rates. Propositions 2 and 3 concern different information losses. The first occurs when physical pairwise transport is inferred from density evolution. The second occurs when gross crossing traffic is inferred from first-order population transport.

Γij

All quantities inside the integrals in (24) are evaluated at (x, t). The reverse rate Hji follows by replacing mij with −mij . Hence, Hij (t) − Hji (t) = Fij (t), (25) and the closure preserves the prescribed first-order population transport. As sij → 0, (24) approaches the coherent result in Corollary 2. When mij = 0 and sij > 0, the net normal flux is zero while the gross crossing rate remains positive. The covariance Σ cannot be inferred from density snapshots alone. It must therefore be supplied by an auxiliary mobility statistic, a trajectory subsample, or a prescribed residualmobility model. The Gaussian closure augments the firstorder transport with this independently supplied second-order information. As an alternative to the Gaussian closure, prescribed velocity moments and a bounded-speed assumption yield distribution-free crossing-rate intervals. The following result applies classical moment-bound principles [40] to directional crossing rates.

B. Closure for Unresolved Normal Mobility We next consider the second information limitation while treating the physical first-order population flux as given. Corollary 2. If Uij is almost surely nonnegative or almost surely nonpositive at almost every point of Γij , then Z Gij (t) = |J(x, t) · nij | dℓ. (19) Γij

If the mean normal flux also has a common sign along Γij , then Gij (t) = |Fij (t)|. (20) Proof. Under the stated condition, E [|Uij | | x, t] = |E [Uij | x, t]| at almost every boundary point. Substitution into (9) gives (19). If the mean normal flux has a common sign along Γij , the integral of its absolute value equals the absolute value of its integral, which gives (20).

Proposition 4. Suppose that |Uij | ≤ vmax almost surely on Γij , where vmax > 0, and the conditional moments mij and s2ij in (21) are given. Let ηij = m2ij + s2ij and define   Z ηij L dℓ, λ max |mij |, Gij = vmax ZΓij (26) √ GU λ ηij dℓ. ij =

Corollary 2 does not require identical UE velocities. It requires only that the conditional normal-velocity distribution at almost every boundary point be supported on one side of zero. When this condition is violated, an additional description of local velocity variation is required. For counter-moving populations, assume that the conditional UE velocity has finite second moments wherever λ(x, t) > 0. Define the conditional velocity covariance, normal mean, and normal variance respectively as h i T Σ(x, t) ≜ E (U − v) (U − v) | x, t , mij (x, t) ≜ v(x, t) · nij =

J(x, t) · nij , λ(x, t)

Γij

Then GLij ≤ Gij ≤ GU ij and Fij + GLij Fij + GU ij ≤ Hij ≤ . (27) 2 2 Both endpoints are attainable within the class of boundarylocal conditional velocity distributions satisfying the prescribed moments and support. Thus, the interval is sharp within this class.

(21)

s2ij (x, t) ≜ nT ij Σ(x, t)nij . These definitions apply wherever λ(x, t) > 0, and zero-density points do not contribute to the boundary integrals. We adopt the Gaussian moment closure
Uij | x, t ∼ N mij (x, t), s2ij (x, t) . (22)

Proof. Write m = mij and η = ηij at a fixed boundary point. The inequalities √ |m| ≤ E|Uij |, η ≤ vmax E|Uij |, E|Uij | ≤ η give (26). Equation (12) then yields (27). For η > 0, the upper √ √ bound is attained by Uij ∈ {− η, η} with probabilities √ √ (1 − m/ η)/2 and (1 + m/ η)/2. For 0 < η ≤ vmax |m|, the lower bound is attained by Uij = sgn(m)η/|m| with probability m2 /η and zero otherwise. For η ≥ vmax |m|, it is attained by Uij ∈ {−vmax , 0, vmax } with probabilities 2 2 2 (η−vmax m)/(2vmax ), 1−η/vmax , and (η+vmax m)/(2vmax ),

This closure is fully specified by the prescribed first two moments and yields closed-form one-sided velocity moments. It is also the maximum-entropy distribution on R for fixed mean and variance. The closure is imposed only on the boundary-normal velocity and does not require the complete two-dimensional UE velocity distribution to be Gaussian. 6

respectively. When η = 0, both bounds vanish. Applying these constructions along the boundary establishes attainability of the integrated endpoints.

therefore require an additional spatial interpolation or density model before continuous reconstruction can be applied. The effect of density observation errors is examined in Section V-D.

L The interval width is (GU ij −Gij )/2. It vanishes when sij = 2 0 or ηij = vmax almost everywhere on the boundary. The bounds do not require a prescribed velocity distribution.

A. Population Flux Reconstruction Suppose that density snapshots are available at uniformly spaced times tk = t0 +k∆t, where k = 0, . . . , K and ∆t > 0. Let λk (x) denote the density snapshot used for reconstruction at time tk . Under ideal observation, λk (x) = λ(x, tk ), and perturbed snapshots are considered separately in the numerical evaluation. For each interior snapshot k = 1, . . . , K − 1, the temporal derivative is approximated by the centered difference λk+1 (x) − λk−1 (x) . (29) qek (x) = 2∆t For sufficiently smooth density fields, (29) is second-order accurate. In the midpoint-based numerical experiments, we instead use (λk+1 −λk )/∆t with the midpoint density (λk+1 + λk )/2, retaining second-order temporal accuracy at the interval midpoint. For the closed system in Section II, the weighted Neumann problem requires the temporal derivative to have zero spatial integral. This condition follows from constant total population, while small numerical inconsistencies are corrected in the corrected temporal derivative qk given by Z Z 1 qk (x) = qek (x)− qek (y) dy with qk (x) dx = 0, (30) |Ω| Ω Ω where |Ω| denotes the area of the service region. This correction is intended only for small compatibility errors. A physical change in the total UE population must instead be represented through an open boundary or a source term and should not be removed by (30). Using λk and qk , the F-STSG reconstruction at time tk is obtained from ∇· (λk (x)∇ϕk (x)) = qk (x), x ∈ Ω,

C. Stationary Isotropic Limit The stationary isotropic case provides a consistency check with the classical stochastic-geometry result and does not use the Gaussian closure. Consider an infinite plane in which BSs form a homogeneous Poisson point process with intensity µBS . Let the UE density be constant, denoted by λu , and let each UE move with constant speed v0 in an independently and uniformly distributed direction. Nearest-BS association then induces a stationary Poisson–Voronoi tessellation. For an arbitrary cell boundary, let Θ denote the angle between the UE velocity and a selected boundary normal. Isotropy gives Θ ∼ Unif[0, 2π). Let Le denote the mean Voronoi edge length per unit area when each geometric edge is counted once. For a homogeneous Poisson–Voronoi √ tessellation, Le = 2 µBS . The mean absolute normal speed and the resulting crossing rate per UE are therefore given respectively by 2v0 , E [|v0 cos Θ|] = π √ (28) 4v0 µBS λu E[|v0 cos Θ|]Le . hUE = = λu π Equation (28) is the classical mean HO rate under isotropic mobility and nearest-BS association in a homogeneous Poisson cellular network [8]. This special case also illustrates the information absent from the first-order population field. Isotropy gives E[U] = 0 and hence J = 0, even though E|Uij | > 0. Because the density is stationary, the minimum-energy F-STSG reconstruction also b = 0. The nonzero crossing rate in (28) therefore arises gives J from velocity variation that is not represented by the first-order population flux.

λk (x)∇ϕk (x) · nΩ (x) = 0, x ∈ ∂Ω, b Jk (x) = −λk (x)∇ϕk (x).

(31)

The potential is determined up to an additive constant, which bk . We assume that λk is sufficiently regular does not affect J and bounded away from zero on the reconstruction domain, as required by the weighted elliptic formulation. In numerical implementation, zero-density regions can be excluded from the active reconstruction support or handled through the regularization adopted by the F-STSG solver. After spatial discretization, (31) yields a sparse linear system with one potential unknown per active grid point. For the local finite-volume discretization used here, the number of nonzero matrix entries grows linearly with the number of active grid points. This sparsity property does not imply linear solution time, which also depends on matrix conditioning and iterative convergence. Reconstructions at different snapshot times are independent and can be processed in parallel. The reconstructed signed flow across a common cell boundary is Z

IV. H ANDOVER E STIMATION FROM D ENSITY S NAPSHOTS The preceding analysis separates HO estimation from aggregate observations into two steps. Density evolution is first used to reconstruct an admissible first-order population transport field, after which directional crossing rates are obtained under a specified local mobility closure. This section implements these steps for discrete density snapshots. The F-STSG reconstruction selects the minimum-energy transport field consistent with the observed population evolution. Proposition 2 and Corollary 1 characterize when this field can represent the physical pairwise transport, while the local mobility closure determines its conversion to directional crossing rates. A density snapshot is a spatial estimate of the UE population over the service region at a given observation time. Such a field may be constructed from gridded positioning estimates, spatially aggregated measurement reports, or other locationaware population measurements [41]. Cell-level UE counts provide only population integrals over the serving regions and

bk (x) · nij dℓ. J

Fbij,k = Γij

7

(32)

velocity covariance. For locations where λk (x) > 0, define bk (x) · nij J b , sb2ij,k (x) = nT m b ij,k (x) = ij Σk (x)nij . (37) λk (x) Points with λk = 0 make no contribution to the crossing integrals. For sbij,k > 0, let m b ij,k zbij,k = . (38) sbij,k Using the Gaussian moment closure in (22), the dispersionaware directional estimate is Z

Because (31) satisfies the population balance associated with qk , the reconstructed cell flows obey Z X b Fij,k = − qk (x) dx. (33) Vi

j∈Ni

Equation (33) provides both a population-consistency relation and a numerical check on the boundary integration. It ensures that the reconstructed pairwise flows reproduce the cell population change implied by the density snapshots, but it does not imply Fbij,k = Fij (tk ) for every neighboring pair. Such equality depends on the identifiability of the physical transport from density evolution. In particular, Corollary 1 identifies a potential-transport class for which the continuous F-STSG reconstruction coincides with the physical first-order flux when the exact density derivative is available.

b disp = H ij,k

λk (b sij,k φ(b zij,k ) + m b ij,k Φ(b zij,k )) dℓ.

(39)

All quantities inside (39) are evaluated at x and time tk . The covariance field is additional mobility information and cannot be inferred from density snapshots alone. It may be obtained from an auxiliary mobility measurement, a trajectory subsample, or a prescribed residual-mobility model. The dispersion-aware estimate also admits an error bound that separates transport reconstruction from the information absent in density evolution and first-order mobility. Define the Gaussian positive-part moment ( aφ(f /a) + f Φ(f /a), a > 0, K(f, a) = (40) [f ]+ , a = 0.

B. Directional Handover Estimation The reconstructed field in (31) determines signed first-order transport consistent with the observed density evolution. The remaining step is to map this transport to directional mobilityinduced crossing rates. As established in Section III-B, this mapping depends on the unresolved local normal-velocity distribution and is separate from the transport reconstruction. When local normal motion is coherent, the reconstructed directional rate is Z h i coh bk (x) · nij b ij,k H = J dℓ. (34) Γij

Γij

b · nij , fbk = J bk · nij , At time tk , let f = J · nij , f c = J b a = λsij , and b ak = λk sbij,k . Here J is the canonical field bk is obtained from the exact density evolution in (5), whereas J reconstructed from snapshots. Write rij = λE[[Uij ]+ | x, tk ]. All local quantities are evaluated on Γij at tk . Then b disp − Hij (tk ) ≤ Erec + Eid + Estat + Ecl , (41) H

+

The reverse rate is obtained by replacing nij with −nij , and hence coh coh b ij,k b ji,k H −H = Fbij,k . (35)

ij,k

The effect of first-order transport reconstruction error on the directional estimate is bounded directly.

where the boundary-integrated terms are Z Erec = |fbk − f c | dℓ, Γij Z Eid = |f c − f | dℓ, Γij Z 1 |b ak − a| dℓ, Estat = √ 2π Γij Z Ecl = |K(f, a) − rij | dℓ.

Proposition 5. Suppose that the physical normal mobility is coherent on Γij at time tk . Then Z   coh bk (x)−J(x, tk ) · nij dℓ. (36) b ij,k J H − Hij (tk ) ≤ Γij

Consequently, convergence of the reconstructed boundarynormal transport in L1 (Γij ) implies convergence of the corresponding coherent directional crossing estimate.

(42)

Γij

√ Indeed, |∂f K| ≤ 1 and |∂a K| ≤ 1/ 2π for a > 0, and the same Lipschitz bounds extend to a = 0 by continuity. Equation (41) follows by the triangle inequality. The terms quantify bounds associated with snapshot reconstruction, canonical-tophysical flux mismatch, auxiliary statistics, and distributional closure, respectively, and they are not additive components of the realized estimation error. For each snapshot, the estimated neighboring-cell rates form an M × M geometric crossing-rate matrix h i bk = H b ij,k . H (43)

Proof. Under coherent normal mobility, Corollary 2 gives Z Hij (tk ) = [J(x, tk ) · nij ]+ dℓ. Γij

The positive-part mapping is nonexpansive, and |[a]+ −[b]+ | ≤ |a − b|. Applying this inequality pointwise on Γij and integrating gives (36). Proposition 5 bounds coherent directional estimation error by the difference between reconstructed and physical boundary-normal fluxes. This difference includes both snapshot reconstruction error and, when present, canonical-tophysical transport mismatch. When local velocity dispersion cannot be neglected, let b k (x) denote an independently supplied estimate of the local Σ

b ii,k = 0 and entries corresponding to nonneighboring Here, H b ij,k is obtained from either cell pairs are zero. The entry H the coherent closure in (34) or the dispersion-aware closure in 8

(39), according to the adopted local mobility model. The overall computational chain is bk , {λk−1 , λk , λk+1 } −→ qk −→ J n o (44) bk −→ Fbij,k , H b ij,k (with closure) . J The first part of (44) reconstructs an F-STSG transport field consistent with the observed density evolution, whereas the final mapping determines directional crossing rates under the adopted local mobility closure. Under coherent normal mobility, no additional velocity statistic is required, and Proposition 5 relates the directional estimation error directly to the boundary-normal transport error. Under unresolved counterflow, the dispersion-aware estimate additionally requires the b k and local normal-velocity variance. Agreement between H the physical mobility-induced crossing-rate matrix therefore depends separately on transport reconstruction and on the adequacy of the local mobility closure.

Fig. 2. Mean gross crossing rate per UE under nonstationary coherent mobility, comparing microscopic trajectory counts, physical boundary integration, and F-STSG reconstruction.

boundary-normal mobility condition in Corollary 2. The simulation region is a 2 km×2 km square containing 100 fixed BSs and 10,000 UEs. The initial population consists of a uniform background and a Gaussian hotspot with mixture weight 0.65 and standard deviation 180 m. A horizontal potential flow with maximum speed 15 m/s drives the population evolution under a no-flux outer boundary condition. Microscopic trajectories are sampled every 0.05 s over 30 independent UE realizations. The F-STSG reconstruction uses noise-free analytical density snapshots at 1 s intervals on a grid with 20 m spacing. The physical reference is obtained directly from the prescribed mobility field, whereas the F-STSG estimate uses density snapshots without access to the physical velocity field. Fig. 2 compares the three gross crossing-rate estimates. The rate first rises from 0.064 to 0.073 s−1 and then falls to 0.061 s−1 as the UE distribution moves across fixed serving regions. The physical and F-STSG curves nearly coincide, whereas trajectory counts fluctuate with finite sampling. Relative to the physical boundary reference, trajectories yield ϵH = ϵF = 9.60 × 10−3 , versus 7.32 × 10−4 for F-STSG, and the cell-balance residual is 1.19 × 10−2 . The smaller FSTSG error reflects analytical density inputs rather than a general advantage over trajectory observations. The results are consistent with the potential-transport recovery condition and Proposition 5.

V. N UMERICAL R ESULTS This section evaluates the kinetic crossing model, its information limits, and density-snapshot reconstruction. Controlled mobility fields test the conditions in Section III, while an independent vehicular simulation provides external validation. Unless otherwise stated, association follows the nearest BS. In synthetic experiments, density, flux, and crossing rates are normalized by the total UE population NUE . For a directed neighboring-cell quantity Qij , where Q denotes H or F , we use the normalized pairwise error P est ref (i,j)∈Ed |Qij − Qij | P , Q ∈ {H, F }, (45) ϵQ = ref (i,j)∈Ed |Qij | where Ed = {(i, j) : j ∈ Ni } is the set of ordered neighboring pairs. In the controlled experiments, the reference quantities are calculated from the prescribed physical mobility field. The b − J∥L2 (Ω) /∥J∥L2 (Ω) , evaluated at flux-field error is ϵJ = ∥J the corresponding reconstruction time. We also evaluate celllevel population conservation using the normalized balance residual ϵbal , defined as the accumulated absolute mismatch between the cell-population change and the net outgoing flow, normalized by the corresponding absolute magnitudes. We distinguish ϵphys bal , computed from the physical flux, from ϵF−STSG , computed from the reconstructed flux. Both use bal the same cell-population changes but their respective boundary fluxes, and measure conservation consistency rather than pairwise reconstruction accuracy. For external evaluation, Qij denotes window-integrated counts. The window-level net-flow (w) error is denoted as ϵF , and the paired baseline difference is: (w)

∆ϵF

(w)

(w)

= ϵF,baseline − ϵF,F-STSG ,

B. Information Limits of First-Order Population Fields We next examine two information limits of first-order population fields: the nonuniqueness of physical pairwise transport given density evolution and the inability of first-order transport to determine gross crossing rates. To illustrate the first limit, we consider a disk of radius 1000 m with 80 BSs and 100,000 uniformly distributed UEs undergoing rigid rotation at 0.01 s−1 . Over 60 s, the divergence-free transport preserves the density while producing inter-cell crossings. Trajectory-derived pairwise flows agree with physical boundary integration to within a normalized error of 2.11 × 10−2 . In Fig. 3(a), the physical pairwise flows remain nonzero, whereas the canonical F-STSG flux vanishes under stationary density, illustrating Proposition 2. We further examine how unobserved divergence-free transport affects pairwise reconstruction. Starting from the potential flow in Section V-A, we add a divergence-free

(46)

(w) (w) where ϵF,baseline and ϵF,F-STSG are the normalized pairwise

net-flow errors of the baseline and F-STSG, respectively, evaluated on the same directed cell pairs in window w. Mean paired differences across windows are distinct from the pooled normalized error over all directed pairs. A. Nonstationary Coherent Mobility We first consider a potential-transport field satisfying the exact-recovery condition in Corollary 1 and the coherent 9

component generated by the stream function χ(x, y) = sin2 (πx/L) sin2 (πy/L), where L = 2 km. Its L2 norm is matched to that of the potential flux before scaling by α. Because this component has zero divergence and is tangent to the outer boundary, varying α changes the physical pairwise flows without changing the density evolution or the F-STSG reconstruction. Fig. 3(b) shows that increasing α from zero to one raises the normalized net-flow and directional crossing errors from approximately 7.32 × 10−4 to 0.854 and 0.817, respectively. In contrast, the physical and F-STSG cell-balance residuals remain near 1.2 × 10−2 , mainly due to differences

between grid-based population integration and Voronoi boundary flux integration. Thus, both fields can satisfy nearly the same population balance while having substantially different pairwise flows, as characterized in Proposition 2. For the second information limit, Fig. 4(a) compares the Gaussian moments in (24) with Monte Carlo estimates using 106 independent samples per point, sij = 5 m/s, and zij = mij /sij ∈ [−3, 3]. The analytical and sampled positive-part and absolute-velocity moments agree closely. Fig. 4(b) considers velocities of +10 and −10 m/s, with reverse-moving fraction ρ ∈ [0, 0.5]. The exact normalized gross rate remains one, whereas the coherent estimate falls to zero as opposite flows cancel. At ρ = 0.5, Gaussian and Laplace closures matched to the first two velocity moments yield 0.798 and 0.707, respectively. For this fixed-speed case, 2 ηij = vmax with vmax = 10 m/s, and therefore Proposition 4 recovers the exact gross crossing rate. C. Stationary Isotropic Limit

Fig. 3. Effect of divergence-free transport on pairwise flow reconstruction. (a) Physical and F-STSG pairwise net-flow magnitudes under stationary density. (b) Net-flow and directional crossing errors versus the divergence-free weight α, compared with the physical and reconstructed cell-balance residuals. The balance residuals measure conservation consistency, not pairwise reconstruction accuracy. Flows in (a) are normalized per UE.

We examine the stationary Poisson–Voronoi limit derived in Section III-C. BS locations follow a homogeneous Poisson point process, and UEs move at constant speed in independent, uniformly distributed directions. Statistics are collected in a 2 km × 2 km central region with a 1500 m guard margin. Each operating point uses 5000 UEs and 30 independent BS realizations. The trajectory Monte Carlo method directly counts intersections between UE trajectories and Voronoi boundaries, whereas the kinetic integral evaluates boundarycrossing rates using the realized Voronoi geometry and the

Fig. 4. Normal-velocity closure and counterflow. (a) Analytical Gaussian positive-part and absolute-velocity moments compared with Monte Carlo estimates. (b) Normalized gross crossing rates under bimodal counterflow, comparing the exact result with coherent, Gaussian, and Laplace closures. The Gaussian and Laplace closures are matched to the mean and variance of the bimodal distribution.

Fig. 5. Stationary Poisson–Voronoi crossing-rate limit. (a) Mean crossing rate per UE versus speed at µBS = 50 BS/km2 , comparing the theoretical result with trajectory Monte Carlo estimates. (b) Absolute relative errors of trajectory Monte Carlo and kinetic boundary integration with respect to the theoretical rate across BS densities at v0 = 10 m/s.

10

Fig. 6. Sensitivity of density-based reconstruction to observation resolution and accuracy. (a) Normalized flux and directional crossing errors versus snapshot interval ∆t at h = 20 m. (b) Reconstruction errors versus spatial grid spacing h at ∆t = 0.5 s. (c) Reconstruction errors versus relative root-mean-square (RMS) density error at ∆t = 1 s.

prescribed velocity distribution. Both methods are evaluated against the theoretical result in (28). Fig. 5(a) confirms the predicted linear dependence on UE speed. Across the BS densities in Fig. 5(b), the maximum trajectory-to-theory deviation is 1.02%, while the largest trajectory-to-kinetic-integral difference is 0.184% of the theoretical rate. Both results are consistent with the square-root BS-density scaling in (28).

TABLE III C ONTROLLED DIRECTIONAL - ERROR DIAGNOSIS UNDER NATIVE 1 S EVALUATION . S TAGE DISCREPANCIES ARE NOT ADDITIVE CONTRIBUTIONS TO THE TOTAL ERROR . Mobility Shape Ident. Recon. Total Potential coherent 0 0 0.305 0.305 Potential Gaussian 0 0 0.222 0.225 Nonpotential (α = 0.5) 0 0.633 0.195 0.688 Bimodal (ρ = 0.5) 0.202 0 0.510 0.551

kinetic reference, physical-mean Gaussian closure, canonicalmean Gaussian closure, and snapshot-based reconstruction. The final estimate includes a closure amplitude based on observed density and a variance estimated from 20% of an independent auxiliary pool. Table III gives normalized adjacent-stage L1 discrepancies and total error, using the same kinetic-reference denominator within each scenario. Unlike the boundary-integral bounds in (42), these are differences between successive estimates. The nonpotential case shows a canonical-to-physical discrepancy of 0.633, whereas balanced bimodal motion gives a Gaussian shape mismatch of 0.202. Under noisy potential transport, the reconstruction discrepancies are 0.305 for coherent motion and 0.222 for Gaussian motion. The total also includes density-amplitude and variance-estimation effects. Unlike Fig. 6(c), which evaluates errors after temporal aggregation, Table III accumulates absolute errors over native 1 s intervals. For the same 2% noise realizations used in Fig. 6(c), these two evaluation orders yield 0.0555 and 0.3021, respectively.

D. Snapshot Resolution and Observation Error We next evaluate the sensitivity of density-based reconstruction to temporal resolution, spatial discretization, and observation errors under the coherent nonstationary model. The temporal experiment uses analytical density snapshots with intervals ∆t ranging from 0.25 to 8 s at a fixed grid spacing of h = 20 m. The spatial experiment varies h ∈ {10, 20, 40} m at ∆t = 0.5 s. To assess observation sensitivity, spatially correlated perturbations with a correlation length of 100 m are added to the analytical density snapshots, with 30 independent realizations at each noise level. In Fig. 6(a), errors change little below ∆t = 1–2 s: for 0.25 ≤ ∆t ≤ 2 s, ϵJ < 1.18 × 10−3 and ϵH < 8.12 × 10−4 . At ∆t = 8 s they increase to 9.33 × 10−3 and 4.25 × 10−3 , respectively, indicating diminishing returns from finer temporal sampling at fixed spatial resolution. Fig. 6(b) shows that reducing h from 40 to 10 m decreases ϵJ from 4.64 × 10−3 to 2.61 × 10−4 and ϵH from 2.87 × 10−3 to 1.86 × 10−4 . For h = 40, 20, and 10 m (2,500, 10,000, and 40,000 unknowns), mean elliptic-solve times are 8.37, 38.75, and 457.52 ms per snapshot. Timings use a single-threaded Jacobi-preconditioned conjugate-gradient solver on an AMD Ryzen 5 9600X and exclude density acquisition and boundary integration. Fig. 6(c) shows stronger sensitivity to density-observation error: at 1%, 2%, and 10% relative RMS density error, the respective (ϵJ , ϵH ) values are (0.120, 0.0175), (0.243, 0.0555), and (1.219, 0.825). The growth reflects error amplification in temporal differentiation, even when the physical transport is of potential form. We next separate error sources using 60 native 1 s intervals on the same 100-BS geometry, 12 independent observation seeds, h = 20 m, and nominal 2% spatially correlated density noise. The staged estimates use, in turn, the physical

E. External Validation and Baseline Comparison We evaluate F-STSG using the public FUMD-AI Kraków vehicular dataset [42], generated with Simulation of Urban MObility (SUMO) and OMNeT++/Simu5G. The dataset contains vehicle trajectories, serving-cell records, and nine fixed BS locations. We construct nearest-BS Voronoi regions and use nonoverlapping 60 s windows with fixed vehicle cohorts, obtaining 99 valid windows and 7523 geometric crossings. The reference crossings are extracted from consecutive 1 s trajectory segments. Density snapshots are estimated from vehicle positions without using serving-cell records or crossing labels. Reconstruction uses grid spacing h ≃ dBS /30 and Gaussian density bandwidth dBS /10, where dBS is the median nearest11

TABLE IV E FFECTS OF TEMPORAL AGGREGATION AND KINETIC CLOSURE ON DIRECTIONAL CROSSING ERRORS . (a) Temporal aggregation Window (s) Oracle EH (95% CI) 1 0.068 [0.058, 0.079] 5 0.235 [0.213, 0.256] 15 0.415 [0.386, 0.441] 58 0.584 [0.556, 0.608]

Graph-MinNorm

Graph-Weighted

(b) Controlled kinetic closure: normalized error EH Mobility Kinetic P-G F-C F-G Low 0.0153 0.0153 0.0189 0.0155 Gaussian 0.0246 0.0246 0.5332 0.0250 Bimodal 0.0232 0.0559 0.5462 0.0574 Correlated 0.0256 0.0256 0.5426 0.0765

Continuity-OF

-0.2

0

0.2

0.4

0.6

P-G: Gaussian closure with physical mean transport; F-C: F-STSG with coherent closure; F-G: F-STSG with Gaussian closure and correct auxiliary variance. Kinetic denotes the physical kinetic reference evaluated against finite-sample trajectory counts.

0.6

high-crossing-edge recall from 0.418 to 0.539 and reduces directional error from 0.772 to 0.611. Thus, aggregation order affects directional estimates even when the reconstructed net flow is unchanged. Across all 99 valid windows, raw Simu5G serving-cell transitions occur 2.567 times as often as geometric crossings, whereas dataset-defined actionable events have a count ratio of 0.928. Geometric crossings thus characterize mobility-related boundary activity rather than executed protocol handovers.

0.4 Graph-MinNorm Graph-Weighted Continuity-OF F-STSG-Common F-STSG-Native * Nominal Top-10%: K=3 (actual 12.5%)

0.2

0 2

3*

4

5

Fig. 7. Baseline comparison on 47 held-out windows. (a) Window-level paired net-flow error differences between each baseline and F-STSG. Dots indicate individual windows, while diamonds and horizontal intervals show mean differences and 95% confidence intervals. Positive differences favor FSTSG. (b) Mean recall of high-crossing directed edges versus the selection budget K among 24 eligible edges. The two F-STSG variants differ in their directional closure and aggregation order.

F. Temporal Aggregation and Kinetic Closure Table IV(a) examines temporal aggregation using trajectoryderived net crossings as oracle inputs on 47 held-out windows. Evaluation uses the common 58 s interior interval of each 60 s window, with 95% confidence intervals. The coherent-closure error increases from 0.068 at 1 s to 0.584 at 58 s as opposite directional crossings cancel over longer intervals. Table IV(b) compares kinetic closures under controlled nonstationary mobility. P-G uses the physical mean transport, whereas F-G uses the prescribed true velocity variance, and these auxiliary quantities serve as oracle information for evaluating the closures. Under matched Gaussian mobility, FG reduces the directional error from 0.5332 (F-C) to 0.0250, close to the physical kinetic reference. The errors increase to 0.0574 under bimodal mobility and 0.0765 under correlated mobility, reflecting distributional mismatch and, in the latter case, canonical-to-physical transport differences.

BS distance. The iterative solver converges on 86 of the 99 valid windows, giving net-flow and directional errors of 0.592 and 0.606, respectively, with an estimated-to-reference crossing-count ratio of 49.9%. An equivalent direct solver resolves all 99 windows without materially changing the original results. The baseline comparison uses 47 held-out windows with complete outputs from F-STSG and the three baselines. We compare F-STSG with three baselines. Graph-MinNorm reconstructs pairwise flows from cell-population balance, while Graph-Weighted additionally accounts for Voronoi boundary lengths. Continuity-OF (continuity-constrained optical flow) estimates a spatial transport field using the population continuity constraint and smoothness regularization. All methods are evaluated on the same held-out windows and directed cell pairs. Unlike the graph-based baselines, F-STSG and Continuity-OF use spatial density information. Fig. 7(a) shows positive mean paired net-flow error differences of 0.224, 0.231, and 0.123 against Graph-MinNorm, Graph-Weighted, and Continuity-OF, respectively, and all paired 95% confidence intervals lie above zero. The pooled normalized net-flow error on the same held-out subset is nevertheless 0.585, distinct from the window-level paired differences and indicating substantial remaining absolute error. Fig. 7(b) compares two aggregation orders using the same reconstructed flux: F-STSG-Common takes the positive part after pairwise aggregation, whereas F-STSG-Native integrates the local positive part first. At K = 3, Native improves

VI. C ONCLUSION This paper has developed a kinetic formulation of geometric serving-region crossings in nonstationary cellular networks. By linking directional crossing rates to boundary-normal velocity moments, the formulation distinguishes signed net transfer from gross bidirectional crossings. We have established two information limits concerning physical transport recovery from density evolution and gross crossing recovery from firstorder transport. Potential transport permits exact canonical flux recovery, while coherent and Gaussian closures characterize directional rates under specified mobility conditions. Sharp moment-constrained bounds further quantify crossingrate uncertainty without prescribing a velocity distribution. The stationary Poisson–Voronoi crossing law follows as a special case. 12

Controlled experiments have examined transport ambiguity, snapshot reconstruction error, and velocity-closure mismatch. On held-out vehicular simulations, F-STSG achieves lower net-flow error than the considered baselines, although substantial absolute error remains. Auxiliary velocity statistics improve directional estimates under matched Gaussian mobility. Future work will extend the framework to radio-dependent association and handover mechanisms.

IEEE Trans. Wireless Commun., vol. 20, no. 3, pp. 1816– 1830, 2021. [13] W. Sun, L. Wang, J. Liu, N. Kato, and Y. Zhang, “Movement Aware CoMP Handover in Heterogeneous Ultra-Dense Networks,” IEEE Trans. Commun., vol. 69, no. 1, pp. 340–352, 2021. [14] R. Arshad and L. Lampe, “Stochastic Geometry Analysis of User Mobility in RF/VLC Hybrid Networks,” IEEE Trans. Wireless Commun., vol. 20, no. 11, pp. 7404– 7419, 2021. [15] ——, “Stochastic Geometry Analysis of Sojourn Time in RF/VLC Hybrid Networks,” IEEE Trans. Wireless Commun., vol. 21, no. 12, pp. 11 237–11 251, 2022. [16] Y. Tian, G. Pan, H. ElSawy, and M.-S. Alouini, “SatelliteAerial Communications with Multi-Aircraft Interference,” IEEE Trans. Wireless Commun., vol. 22, no. 10, pp. 7008–7024, 2023. [17] Y. Tian, G. Pan, M. A. Kishk, and M.-S. Alouini, “Stochastic Analysis of Cooperative Satellite-UAV Communications,” IEEE Trans. Wireless Commun., vol. 21, no. 6, pp. 3570–3586, 2022. [18] W.-Y. Dong, S. Yang, P. Zhang, and S. Chen, “Stochastic Geometry Based Modeling and Analysis of Uplink Cooperative Satellite-Aerial-Terrestrial Networks for Nomadic Communications with Weak Satellite Coverage,” IEEE J. Sel. Areas Commun., vol. 42, no. 12, pp. 3428–3444, 2024. [19] W.-Y. Dong, S. Yang, and S. Chen, “Uplink Performance Analysis of Heterogeneous Non-Terrestrial Networks in Harsh Environments: A Novel Stochastic Geometry Model,” IEEE Trans. Commun., vol. 73, no. 8, pp. 6734– 6747, 2025. [20] F. Teng, S. Mou, C. Yang, J. Yuan, Q. Zhang, Z. Lin, and W. Ni, “Communication-Efficient Satellite-Ground Federated Learning over Intermittent LEO Satellite Networks,” IEEE Trans. Mobile Comput., 2026. [21] Q. Zhang, X. Li, L. Shi, Z. Lin, J. Wen, A. Zhou, and S. Wang, “Resource-Aware Parameter-Efficient Model Adaptation for Onboard High-Dimensional Data,” arXiv preprint arXiv:2609.33687, 2026. [22] Z. Lin, Z. Chen, Z. Fang, X. Chen, X. Wang, and Y. Gao, “FedSN: A Federated Learning Framework over Heterogeneous LEO Satellite Networks,” IEEE Trans. Mobile Comput., vol. 24, no. 3, pp. 1293–1307, 2024. [23] M. F. Özkoç, A. Koutsaftis, R. Kumar, P. Liu, and S. S. Panwar, “The Impact of Multi-Connectivity and Handover Constraints on Millimeter Wave and Terahertz Cellular Networks,” IEEE J. Sel. Areas Commun., vol. 39, no. 6, pp. 1833–1853, 2021. [24] S. S. Kalamkar, F. Baccelli, F. M. Abinader, Jr., A. S. Marcano Fani, and L. G. Uzeda Garcia, “Beam Management in 5G: A Stochastic Geometry Analysis,” IEEE Trans. Wireless Commun., vol. 21, no. 4, pp. 2275–2290, 2022. [25] S. Aghashahi, S. Aghashahi, Z. Zeinalpour-Yazdi, A. Tadaion, and A. Asadi, “Stochastic Modeling of Beam Management in mmWave Vehicular Networks,” IEEE Trans. Mobile Comput., vol. 22, no. 6, pp. 3665–3676,

R EFERENCES [1] I. Shayea et al., “Key Challenges, Drivers and Solutions for Mobility Management in 5G Networks: A Survey,” IEEE Access, vol. 8, pp. 172 534–172 552, 2020. [2] E. Gures, I. Shayea, A. Alhammadi, M. Ergen, and H. Mohamad, “A Comprehensive Survey on Mobility Management in 5G Heterogeneous Networks: Architectures, Challenges and Solutions,” IEEE Access, vol. 8, pp. 195 883–195 913, 2020. [3] H. Tabassum, M. Salehi, and E. Hossain, “Fundamentals of Mobility-Aware Performance Characterization of Cellular Networks: A Tutorial,” IEEE Commun. Surveys Tuts., vol. 21, no. 3, pp. 2288–2308, 2019. [4] Z. Lin, L. Wang, B. Tan, and X. Li, “Spatial-Spectral Terahertz Networks,” IEEE Trans. Wireless Commun., 2021. [5] C. Yang, H. Pan, J. Wang, and Y. Hong, “TrajDiff: Trajectory Prediction with Diffusion Probabilistic Models,” IEEE Trans. Image Process., vol. 34, pp. 8257–8270, 2025. [6] S. Zhang, J. Yang, Z. Song, Z. Lin, S. Prajapat, Z. Xia, H. Ghayvat, H. Ding, L. Y. Por, and A. K. Das, “Large Language Model Enhanced Differentiable Trajectory Planning for IoT-Enabled Autonomous Driving,” IEEE Internet Things J., 2026. [7] W.-Y. Dong, W. Jiang, S. Zhao, Q. Bi, and S. Chen, “Fluid-Spatiotemporal Stochastic Geometry: Information Flow in Non-Stationary Fields,” arXiv preprint arXiv:2607.00616, 2026. [8] X. Lin, R. K. Ganti, P. J. Fleming, and J. G. Andrews, “Towards Understanding the Fundamentals of Mobility in Cellular Networks,” IEEE Trans. Wireless Commun., vol. 12, no. 4, pp. 1686–1698, 2013. [9] S. Sadr and R. S. Adve, “Handoff Rate and Coverage Analysis in Multi-Tier Heterogeneous Networks,” IEEE Trans. Wireless Commun., vol. 14, no. 5, pp. 2626–2638, 2015. [10] W. Bao and B. Liang, “Stochastic Geometric Analysis of User Mobility in Heterogeneous Wireless Networks,” IEEE J. Sel. Areas Commun., vol. 33, no. 10, pp. 2212– 2225, 2015. [11] H. Ibrahim, H. ElSawy, U. T. Nguyen, and M.-S. Alouini, “Mobility-Aware Modeling and Analysis of Dense Cellular Networks with C-Plane/U-Plane Split Architecture,” IEEE Trans. Commun., vol. 64, no. 11, pp. 4879–4894, 2016. [12] M. Salehi and E. Hossain, “Stochastic Geometry Analysis of Sojourn Time in Multi-Tier Cellular Networks,” 13

2023. [26] Z. Lin, L. Wang, J. Ding, Y. Xu, and B. Tan, “Tracking and Transmission Design in Terahertz V2I Networks,” IEEE Trans. Wireless Commun., vol. 22, no. 6, pp. 3586– 3598, 2022. [27] H. Wei and H. Zhang, “Handover Optimization for IRSAided Networks with Experimental Verification,” IEEE Trans. Wireless Commun., vol. 24, no. 11, pp. 9217– 9232, 2025. [28] W. Huang et al., “Self-Adapting Handover Parameters Optimization for SDN-Enabled UDN,” IEEE Trans. Wireless Commun., vol. 21, no. 8, pp. 6434–6447, 2022. [29] R. Karmakar, G. Kaddoum, and S. Chattopadhyay, “Mobility Management in 5G and Beyond: A Novel Smart Handover with Adaptive Time-to-Trigger and Hysteresis Margin,” IEEE Trans. Mobile Comput., vol. 22, no. 10, pp. 5995–6010, 2023. [30] A. Prado, F. Stöckeler, F. Mehmeti, P. Krämer, and W. Kellerer, “Enabling Proportionally-Fair Mobility Management with Reinforcement Learning in 5G Networks,” IEEE J. Sel. Areas Commun., vol. 41, no. 6, pp. 1845–1858, 2023. [31] M. Kalntis, G. Iosifidis, J. Suárez-Varela, A. Lutu, and F. A. Kuipers, “Meta-Learning-Based Handover Management in NextG O-RAN,” IEEE J. Sel. Areas Commun., vol. 44, pp. 3251–3267, 2026. [32] T. Duan, Z. Zhang, Z. Lin, S. Guo, X. Guan, G. Wu, Z. Fang, H. Meng, X. Du, J.-Z. Zhou et al., “LLMDriven Stationarity-Aware Expert Demonstrations for Multi-Agent Reinforcement Learning in Mobile Systems,” IEEE Trans. Mobile Comput., 2026. [33] M. G. Demissie, G. H. de A. Correia, and C. Bento, “Exploring Cellular Network Handover Information for Urban Mobility Analysis,” J. Transp. Geogr., vol. 31, pp. 164–170, 2013. [34] M. Kalntis et al., “Through the Telco Lens: A Countrywide Empirical Study of Cellular Handovers,” in Proc. IMC, 2024, pp. 51–67. [35] W.-Y. Dong, S. Zhao, R.-S. Han, Q. Bi, and S. Chen, “Digital Tides: A Fluid-Dynamic Framework for FluxAware Infrastructure Provisioning in UAV Logistics Networks,” IEEE Trans. Mobile Comput., vol. 25, no. 10, pp. 16 705–16 720, 2026. [36] W.-Y. Dong et al., “Vorticity Dissipation Based Routing: A Fluid-Kinetic Framework for Loop-Free Transport in Ultra-Dense Networks,” IEEE Trans. Mobile Comput., vol. 25, no. 11, pp. 22 370–22 387, 2026. [37] ——, “Fluid-Dynamic Interference Modeling for LEO Mega-Constellations: A Spatiotemporal Kinetic Field Approach,” arXiv preprint arXiv:2608.20651, 2026. [38] X. Wang, Y. Du, B. Gu, Z. Lin, X. Li, Y. Gong, W. Gong, and L. Kong, “Exploring the Initial Performance of NBIoT NTN over GEO: Measurement and Analysis,” arXiv preprint arXiv:2609.06280, 2026. [39] W.-Y. Dong, W. Jiang, S. Zhao, Q. Bi, and S. Chen, “Fluid-Spatiotemporal Stochastic Geometry for Information Flow in Non-Stationary Fields: Kinematics and Field Reconstruction,” arXiv preprint arXiv:2610.23277, 2026.

[40] D. Bertsimas and I. Popescu, “Optimal Inequalities in Probability Theory: A Convex Optimization Approach,” SIAM J. Optim., vol. 15, no. 3, pp. 780–804, 2005. [41] H. Yuan, Z. Chen, Z. Lin, J. Peng, Y. Zhong, X. Hu, S. Xue, W. Li, and Y. Gao, “Constructing 4D Radio Map in LEO Satellite Networks with Limited Samples,” in Proc. INFOCOM, 2025, pp. 1–10. [42] C. Bernad, S. Filiposka, and K. Gilly, “FUMD-AI Kraków Vehicular Handover Dataset,” Hugging Face Hub, 2026, https://huggingface.co/datasets/FUMD-AI/ krakow-handover-dataset.

14

Related documents

• Kinetic Handover Modeling in Nonstationary Cellular Networks with Fluid Spatiotemporal Stochastic Geometry#1457034
• A Stochastic Geometry Framework for Performance Analysis of RIS-Assisted OFDM Cellular Networks#1459697
• Heterogeneity and Aggregate Fluctuations#21691
• Data Free Training of Diffusion Models for Lattice Field Theories#1344352
• Data Free Training of Diffusion Models for Lattice Field Theories#1319003