Member, IEEE Imène Trigui

Senior Member, IEEE Sofiène Affes

Fellow, IEEE Marco Di Renzo

Senior Member, IEEE Nalin K Jayakody

Coverage Analysis and Scaling Laws in Ultra-Dense Networks

Keywords: Network densification, antenna pattern, stochastic geometry, millimeter wave, antenna height, coverage probability, Fox's H-fading

In this paper, we develop an innovative approach to quantitatively characterize the performance of ultra-dense wireless networks in a plethora of propagation environments.

The proposed framework has the potential of simplifying the cumbersome procedure of analyzing the coverage probability and allowing the unification of single-and multi-antenna networks through compact analytical representations. By harnessing this key feature, we develop a novel statistical machinery to study the scaling laws of wireless networks densification considering general channel power distributions including small-scale fading and shadowing as well as associated beamforming and array gains due to the use of multiple antenna. We further formulate the relationship between network density, antenna height, antenna array seize and carrier frequency showing how the coverage probability can be maintained with ultra-densification. From a system design perspective, we show that, if multiple antenna base stations are deployed at higher frequencies, monotonically increasing the coverage probability by means of ultra-densification is possible, and this without lowering the antenna height. Simulation results substantiate performance trends leveraging network densification and antenna deployment and configuration against path loss models and signal-to-noise plus interference thresholds.

I. INTRODUCTION

Chiefly urged by the unfolding mobile data deluge, a radical design make-over of cellular systems enabled by the so-called network densification and heterogeneity, primarily through the provisioning of small cells, has become an extremely active and promising research topic [START_REF] Andrews | What will 5G be?[END_REF]- [START_REF] Trigui | Unified stochastic geometry modeling and analysis of cellular networks in LOS/NLOS and shadowed fading[END_REF]. While smallcell densification has been recognized as a promising solution to boost capacity and enhance coverage with low cost and power-efficient infrastructure in 5G networks, it also paves the way for reliable and high capacity millimeter wave (mmWave) communication and directional beamforming [START_REF] Andrews | What will 5G be?[END_REF]. Nevertheless, there has been noticeable divergence between the above outlook and conclusions of various studies on the fundamental limits of network densification, according to which the latter may eventually stop, at a certain point, delivering significant capacity gains [START_REF] Andrews | A tractable approach to coverage and rate in cellular networks[END_REF]- [START_REF] Lee | Spectral efficiency scaling laws in dense random wireless networks with multiple receive antennas[END_REF].

In this respect, several valuable contributions leverage stochastic geometry (SG) to investigate ultra-dense networks performance under various pathloss and propagation models [START_REF] Andrews | A tractable approach to coverage and rate in cellular networks[END_REF]- [START_REF] Di | Stochastic geometry modeling of coverage and rate of cellular networks using the Gil-Pelaez inversion theorem[END_REF]. In the single-input single-output (SISO) context, Work supported by the Discovery Grants and the CREATE PERSWADE (www.create-perswade.ca) programs of NSERC and a Discovery Accelerator Supplement (DAS) Award from NSERC. Part of this work has been published in the IEEE WCNC 2020 [START_REF] Trigui | SINR Coverage Analysis of Dense HetNetsOver Fox's H-Fading Channels[END_REF].

conflicting findings based on various choices of path-loss models have identified that the signal-to-noise plus interference (SINR) invariance property, which enables a potentially infinite aggregated data rate resulting from network densification based on the power-law model [START_REF] Andrews | A tractable approach to coverage and rate in cellular networks[END_REF]- [START_REF] Trigui | Unified analysis and optimization of D2D communications in cellular Networks over fading channels[END_REF], vanishes once a more physically feasible path loss model is considered. In the latter case, [START_REF] Liu | Effect of densification on cellular network performance with bounded path-loss model[END_REF]- [START_REF] Ding | Performance impact of LoS and NLoS transmissions in dense cellular networks[END_REF] showed that the coverage probability attains a maximum point before starting to decay when the network becomes denser. Most recently, the authors of [START_REF] Ding | Performance impact of LoS and NLoS transmissions in dense cellular networks[END_REF], [START_REF] Atzeni | Downlink cellular network analysis with LOS/NLOS propagation and elevated base stations[END_REF] and [START_REF] Filo | Stochastic geometry analysis of ultra-dense networks: Impact of antenna height and performance limits[END_REF] have investigated the limits of network densification when the path-loss model includes the antenna height. Besides invalidating the SINR invariance property in this case, these works find that by lowering the antenna height the coverage drop due to ultra densification can be totally offset, thereby improving the network capacity. Motivated by the tractability of the considered system models, most of the previous works assumed the scenario of exponential-based distributions for the channel gains (e.g., integer fading parameter-based power series [START_REF] Di Renzo | Average rate of downlink heterogeneous cellular networks over generalized fading channels: A stochastic geometry approach[END_REF], [START_REF] Gupta | Downlink multi-antenna heterogeneous cellular network with load balancing[END_REF], [START_REF] Khoshkholgh | Coverage analysis of max-SIR cell association in hetNets under nakagami fading[END_REF], and Laguerre polynomial series in [START_REF] Chun | A stochastic geometric analysis of device-to-device communications operating over generalized fading channels[END_REF]) and unbounded power law models, while the few noteworthy studies that incorporate general fading, shadowing and path-loss models often lead to complex mathematical frameworks that fail to explicitly unveil the relationship between network density and system performance [START_REF] Di | Stochastic geometry modeling of coverage and rate of cellular networks using the Gil-Pelaez inversion theorem[END_REF], [START_REF] Trigui | Unified stochastic geometry modeling and analysis of cellular networks in LOS/NLOS and shadowed fading[END_REF], [START_REF] Trigui | Unified analysis and optimization of D2D communications in cellular Networks over fading channels[END_REF]. Moreover, although some works investigated the effect of pathloss singularity [START_REF] Elsawy | Modeling and analysis of cellular networks using stochastic geometry: A tutorial[END_REF], [START_REF] Zhang | Downlink cellular network analysis with multi-slope path loss models[END_REF], [START_REF] Ganti | Asymptotics and approximation of the SIR distribution in general cellular networks[END_REF] or boundedness [START_REF] Li | Success probability and area spectral efficiency in multiuser MIMO HetNets[END_REF], [START_REF] Dhillon | Modeling and analysis of K-tier downlink heterogeneous cellular networks[END_REF], the incorporation of the combined effect of path-loss and generalized fading channel models is usually ignored. This has entailed divergent or even contrasting conclusions on the fundamental limits of network densification [START_REF] Li | Success probability and area spectral efficiency in multiuser MIMO HetNets[END_REF], [START_REF] Ding | Performance impact of LoS and NLoS transmissions in dense cellular networks[END_REF], [START_REF] Atzeni | Downlink cellular network analysis with LOS/NLOS propagation and elevated base stations[END_REF]. More importantly, additional work is necessary to investigate advanced communication and signal processing techniques, e.g., massive multiple-input-multipleoutput (MIMO), coordinated multipoint (CoMP) and mmWave communications [START_REF] Chandrasekhar | Coverage in multi-antenna two-tier networks[END_REF]- [START_REF] George | Ergodic spectral efficiency in MIMO cellular networks[END_REF] that are expected to enhance the channel gain.

Motivated by the above background, our work proposes a unified and comprehensive multiple-parameter Fox's H fading model for general multi-path and/or shadowing distributions, which is proved to enable a tractable analysis of dense networks. The main objective of this paper, in particular, is to introduce a non model-specific channel model that leads to a new unified approach to asses the performance of dense networks. To this end, the proposed framework is based on Fox's H transform theory and the Mellin-Barnes integrals along with SG to investigate the performance limits of network densification under realistic pathloss models and general channel power distributions, including propagation impediments and transmission gains due to the antenna pattern and beamforming, which are particulary relevant in multiantenna settings. Several works [START_REF] Lee | Spectral efficiency scaling laws in dense random wireless networks with multiple receive antennas[END_REF]- [START_REF] Yu | A unified framework for the tractable analysis of multi-antenna wireless networks[END_REF] studied different performance metrics to characterize the performance of multiantenna cellular networks, yet under the assumption of the standard power-law path-loss model, since it leads to tractable analysis. In this paper, by leveraging a novel methodology of analysis that is compatible with a wide class of pathloss models, including the antenna height, we are able to study the achievable performance of multi-antenna networks and understand how scaling the deployment density of the base stations (BSs) helps maintain the per user-coverage in dense networks. The main contributions of this paper are the following:

• We introduce a unified analytical framework for analyzing heterogeneous networks under general Fox's H distributed channel models and both unbounded and bounded path-loss models. Closed-form expressions for the coverage probability and corresponding scaling laws allow us to confirm that the path-loss model plays a significant role in determining the network performance, a result corroborated by other recent works on ultra-dense networks [START_REF] Andrews | A tractable approach to coverage and rate in cellular networks[END_REF]- [START_REF] Yu | A unified framework for the tractable analysis of multi-antenna wireless networks[END_REF].

• By exploiting the proposed Fox's H-based channel power representation, we obtain an analytical framework for the coverage probability in multi-antenna networks that is shown to preserve the tractability of the single-antenna case. The asymptotic performance limits of multi-antenna networks are derived in closed-form showing that there is potential for improving the scaling laws of the coverage by increasing the number of BS antennas. • Harnessing the tractability of the developed analytical model, the impact of network densification is investigated by considering advanced transmission techniques, such as MIMO and directional beamforming, and by considering the effect of high transmission frequencies (e.g., mmWave). We show that maintaining the maximum coverage is possible by deploying multiple antenna at the BS and by operating at higher frequency bands, and without lowering the BS height. The obtained scaling laws provide valuable system design guidelines for optimizing general networks deployment. The rest of the paper is organized as follows. In Section II, we present the system model and the modeling assumptions. In Section III, we introduce our approach to obtain exact closedform expressions and scaling laws for the coverage probability. Section IV is focused on multi-antenna BSs and single-antenna users, i.e., MISO networks. Applications of the obtained coverage expressions in different wireless communication scenarios are detailed in order to leverage the full potential of network densification. Numerical and simulation results are illustrated in Section V. Finally, Section IV concludes the paper.

II. SYSTEM MODEL

We consider the downlink transmission of a T -tier heterogeneous wireless network. We focus on the performance analysis of a typical user equipment (UE) which is assumed, without loss of generality, to be located at the origin and to be served by the k-th tier. Hence, its SINR is given by

SINR k = L(r k )g x k T j=1 ri∈Φj \r k P i L(r i )g xi + σ 2 k , (1)
where the following notation is used:

• L(r i ) is the large-scale channel gain between the typical UE and the BS at distance r i , where L(r) = r -α for an unbounded path-loss and L(r) = (1+r) -α for a bounded one.

• Pi = Pi P k is the power of the i-th BS normalized by the power of the BS with index k serving the typical UE.

• σ 2 k is the normalized noise power defined as

σ 2 k = σ 2 P k
• g x k is the channel power gain for the desired signal from the associated transmitter located at x k . Different channel distributions and MIMO techniques lead to different distributions for g x k . In this paper, a general type of distribution is assumed for g x k , as in Assumption 1.

Assumption 1: The channel power gain g x k for the typical UE has a Fox's H distribution, i.e., g x k ∼ H u,v p,q (x; P k ), with the parameter sequence

P k = (κ k , c k , a k , b k , A k , B k
) and probability density function (pdf) [START_REF] Mathai | The H-function: Theory and Applications[END_REF]

f gx k (x) = κ k H u,v p,q c k x (a k , A k ) p (b k , B k ) q , x ≥ 0. (2)
The main advantage of the H-function representation for statistical distributions is that any algebraic combination involving products, quotients, or powers of any number of independent positive continuous random variables can be written as an H-function distribution. Indeed, the Fox's H function distribution captures composite effects of multipath fading and shadowing, subsuming a wide variety of important or generalized fading distributions adopted in wireless communications such as α-µ1 , N -Nakagami-m, (generalized) K-fading, and Weibull/gamma fading, and the Fisher-Snedecor F-S F ( [START_REF] Yoo | The Fisher-Snedecor F distribution: A simple and accurate composite fading model[END_REF] and [START_REF] Yilmaz | A novel unified expression for the capacity and bit error probability of wireless communication systems over generalized fading channels[END_REF] and references therein). • α is the path loss exponent.

• g xi is the interferer's power gain from the interfering transmitter located at x i . In the proposed framework, we assume that g xi , i ∈ {1, . . . T } are non-negative random variables that are independent and identically distributed according to (2).

III. UNIFIED ANALYTICAL FRAMEWORK

In this section, we derive the complementary cumulative distribution function (ccdf) of the SINR, also called the coverage probability, in single-antenna networks. The obtained framework is utilized to analyze multi-antenna networks in general settings in the next section.

C B = T k=1 λ k ∞ 0 e -j∈T πλj P δ j δξ(Ψ1-Ψ2) H v,u q,p+1 ξ, P k B ξ 2 T j=1 πλ j P δ j δξ(Ψ 1 + Ψ 2 ) H 1,1 1,1 T j=1 λ j P δ j (1+δξΨ 1 ) T j=1 λ j P δ j δξ(2Ψ 1 +Ψ 2 ) , P δ dξ, (6)
where P δ = (1, 1, -1, 0, 2, 1), P k B = P k U and Ψ x = H v+1,u+2 q+2,p+3 ξ, P x,I B , x ∈ {1, 2} with

P 1,j,I B = κ j c 2 j , 1 c j , (1 -b j -2B j , 0, δ), (0, 1 -a j -2A j , -1, δ -1), (B j , 1, 1), (1, A j , 1, 1) (7)
and

P 2,j,I B = κ j c 2 j , 1 c j , 1-b j -2B j , 0, δ 2 , 0, 1-a j -2A j , -1, δ 2 -1 , (B j , 1, 1), (1, A j , 1, 1) (8)

A. Coverage Analysis in Closest-BS-Association-Based Cellular Networks

Proposition 1: When the locations of the BSs are modeled as a Poisson point process (PPP) [START_REF] Di Renzo | A mathematical framework to the computation of the error probability of downlink MIMO cellular networks by using stochastic geometry[END_REF] and the nearest-BS association is adopted, the SINR coverage probability at the typical UE for an unbounded path-loss model and the SINR thresholds β k , k ∈ {1, . . . , T }, is given by

C U = πδ T k=1 λ k P k σ 2 k δ ∞ 0 1 ξ 2+δ H v,u q,p+1 ξ, P k U H 1,1 1,1   P k ξσ 2 k δ T j=1 πλ j P δ j 1+δξH v+1,u+2 q+2,p+3 ξ, P Ij U ,P δ   dξ, (3)
where δ =2 α , P δ = (1, 1, 1 -δ, 0, δ, 1), with

P k U = κ k β k , 1 c k β k , 1-b k , (1-a k , 1), B k , (A k , 1) , (4)
and

P Ij U = κ j c 2 j , 1 c j , (1
-b j -2B j , 0, δ), (0, 1-a j -2A j , -1, δ-1),

(B j , 1, 1), (1, A j , 1, 1) . (5)
Proof: See Appendix A. The main assumptions in Proposition 1 are the Fox's H distributed signal and interference channel power gains and the standard power-law unbounded path loss model. The unbounded power-law path loss is known to be inaccurate for short distances, due to the singularity at the origin, which affects the scaling laws of the coverage probability [START_REF] Elsawy | Modeling and analysis of cellular networks using stochastic geometry: A tutorial[END_REF]. Next, a more physically feasible path-loss model is considered.

Proposition 2: When a bounded path-loss model is adopted, the coverage probability of cellular networks based on the nearest-BS association strategy is given in ( 6)- [START_REF] Gupta | Downlink multi-antenna heterogeneous cellular network with load balancing[END_REF], that are shown at the top of this page.

Proof: See Appendix B. Remark 1: For arbitrary distributions for the channel gain, the coverage expressions in (3) and ( 6) are independent of the n-th derivative of the Laplace transform of the aggregate interference, n ∈ [0, ∞), while accurately reflecting the behavior of multi-tiers networks in all operating regimes without the need of applying approximations or upper bounds. Compared with the coverage approximations in [START_REF] Chun | A stochastic geometric analysis of device-to-device communications operating over generalized fading channels[END_REF], [START_REF] Gupta | Downlink multi-antenna heterogeneous cellular network with load balancing[END_REF], and [START_REF] Khoshkholgh | Coverage analysis of max-SIR cell association in hetNets under nakagami fading[END_REF] and expressions in [START_REF] Di | Stochastic geometry modeling of coverage and rate of cellular networks using the Gil-Pelaez inversion theorem[END_REF], [START_REF] Di Renzo | Average rate of downlink heterogeneous cellular networks over generalized fading channels: A stochastic geometry approach[END_REF], the proposed approach yields a more compact analytical result for the coverage probability, where only an integration of Fox's-H functions is needed thanks to the novel handling of fading distributions. Table II lists some commonly-used channel fading distributions and the corresponding expression for C. It is worthy to note that the proposed framework can be extended to other network models, for example, where the transmitters are spatially distributed according to other point processes [START_REF] Di Renzo | System-level analysis/optimization of cellular networks with simultaneous wireless information and power transfer: Stochastic geometry modeling[END_REF], [START_REF] Di Renzo | Inhomogeneous double thinning Modeling and analysis of cellular networks by using inhomogeneous Poisson point processes[END_REF], notably including non-Poisson models [START_REF] Ganti | Asymptotics and approximation of the SIR distribution in general cellular networks[END_REF], or under multi-slope path loss models [START_REF] Nguyen | Performance limits of network densification[END_REF], [START_REF] Zhang | Downlink cellular network analysis with multi-slope path loss models[END_REF]. Hence, the results of this paper allow an exact and tractable approximation for the coverage probability of any stationary and ergodic point process [START_REF] Ganti | Asymptotics and approximation of the SIR distribution in general cellular networks[END_REF], [START_REF] Di Renzo | System-level analysis/optimization of cellular networks with simultaneous wireless information and power transfer: Stochastic geometry modeling[END_REF]. In the next section, we show the usefulness of the proposed approach for obtaining insightful design guidelines for multi-antenna and mmWave networks.

B. Coverage Analysis in Strongest-BS-Association-Based Cellular Networks

The strongest-BS association rule, according to which the serving BS is the one that provides the maximum signal-tointerference (SIR) 2 , can be particularly advantageous for application to scenarios in which the closest-BS association strategy may provide poor performance due to severe blockage. Also, the strongest-BS association criterion may yield performance bounds for other, more practical, cell association strategies.

Proposition 3: When the strongest-BS association is adopted, the SIR coverage probability of the typical UE, given the SIR thresholds β k , k ∈ {1, . . . , T }, is given by

C = 2π T k=1 κ k λ k c k ∞ 0 r k Υ(r k )d r k , = π C(δ) T k=1 λ k β -δ k Λ k j∈T λ j P δ j Λ j , (9)

TABLE I COVERAGE PROBABILITY OF SOME WELL-KNOWN FADING CHANNEL MODELS BASED ON THE CLOSEST-BS STRATEGY

Instantaneous Fading Distribution

Coverage Probability C U Gamma Fading

fg(z) = m Γ(m) H 1,0 0,1 mz - (m -1, 1) C U = πδ T k=1 λ k m k β K Γ(m k ) P k σ 2 k δ ∞ 0 H 0,1 1,1 ξ m k β k (2 -m k , 1) (1, 1) ξ 2+δ H 1,1 1,1      P k σ 2 k δ ξ δ T j=1 πλ j P δ j 1+ δξ m j Γ(m j ) H 1,3 3,3 ξ m j
(-m j , 1), (0, 1), (δ, 1) (0, 1), (-1, 1), (δ -1, 1)

,P δ      dξ Generalized Gamma fg(z) = µ Γ(m) H 1,0 0,1 µz - (m -1 η , 1 η )
where µ =

Γ(m+ 1 η ) Γ(m) . C U = πδ T k=1 λ k µ k β K Γ(m k ) P k σ 2 k δ ∞ 0 H 0,1 1,1 ξ µ k β k (1 + 1 η k -m k , 1 η k ) (1, 1) ξ 2+δ H 1,1 1,1      P k σ 2 k δ ξ δ T j=1 πλ j P δ j 1+ δξ µ j Γ(m j ) H 1,3 3,3 ξ µ j (1 -1 η j -m j , 1), (0, 1), (δ, 1) (0, 1), (-1, 1), (δ -1, 1) ,P δ      dξ Power of Nakagami-n (Rice) fg(z) = ∞ k=0 Ψ k m k Γ(m k )Ω k H 1,0 0,1 m k Ω k z - (m k -1, 1) with m k = k + 1 and Ω k = k+1
1+K R , where KR is the Rician factor.

C U = πδ T k=1 lim K k -→∞ K k t=0 Ψtmt Γ(mt)Ωtλ k β k P k σ 2 k δ ∞ 0 H 0,1 1,1 ξΩt mtβ k (2 -mt, 1) (1, 1) ξ 2+δ H 1,1 1,1      P k σ 2 k δ ξ δ T j=1 πλ j P δ j 1+ ∞ t=0 ΨtΩtδξ Γ(mt)mt H 1,3 3,3
ξΩt mt (-mt, 1), (0, 1), (δ, 1) (0, 1), (-1, 1), (δ -1, 1)

,P δ      dξ , where Ψ k = K k R e -K R /Γ(k + 1).
Lognormal Fading

fg(z) = N n=0 wn ωn H 0,0 0,0 z ωn - -
where ωn = 10 √ 2σun+µ , while un and wn represrent the weight factors and the zeros of the N -order Hermite polynomial [START_REF] Di Renzo | Average rate of downlink heterogeneous cellular networks over generalized fading channels: A stochastic geometry approach[END_REF]Table 25.10].

C U = πδ T k=1 N k t=0 wt ωtλ k β k P k σ 2 k δ ∞ 0 H 0,0 0,1 ξωt β k - (1, 1) ξ 2+δ H 1,1 1,1      P k σ 2 k δ ξ δ T j=1 πλ j P δ j   1+ N j t=0 wtδξ ω 2 t H 1,2
2,3 ξωt (0, 1), (δ, 1) (0, 1), (-1, 1), (δ -1, 1)

  ,P δ      dξ Fisher-Snedecor Fading fg(z) = m msΓ(ms)Γ(m) H 1,1 1,1 mz ms (-ms, 1) (m -1, 1) C U = πδ T k=1 λ k m k β K ms k Γ(m k )Γ(ms k ) P k σ 2 k δ ∞ 0 H 1,1 1,2 ξms k m k β k (2 -m k , 1) (1 + ms k , 1), (1, 1) ξ 2+δ H 1,1 1,1      P k σ 2 k δ ξ δ T j=1 πλ j P δ j 1+ δξ m j Γ(m j ) H 1,3 3,3 ξ m j (-m j , 1), (0, 1), (δ, 1) (0, 1), (-1, 1), (δ -1, 1) ,P δ      dξ with Υ(r k ) = H u+1,v p+1,q+1   T j=1 πr 2 k λ j Γ(1-δ)Λ j P -δ j (c k β k ) -δ (a k +A k , δA k ), (1, δ) (0, 1), (b+B k , δB k )   , (10)
where C(δ) = π 2 δ csc(πδ) and

Λ j = κ j c δ+1 j u t=1 Γ (b jt + (1 + δ)B jt ) p t=u+1 Γ (1 -b jt -(1 + δ)B jt ) × v k=1 Γ (1 -a j k -(1 + δ)A j k ) p k=v+1 Γ (a j k + (1 + δ)A j k ) . ( 11
)
Proof: The proof follows from Appendix C along with the fact that

E r k [Υ(r k )] = 2πλ k ∞ 0 r k H u+1,v p+1,q+1 T j=1 πr 2 k λ j Γ(1-δ)Λ j P -δ j (c k β k ) -δ (a k +A k , δA k ), (1, δ) (0, 1), (b+B k , δB k ) d r k . (12)
Then, applying the transformation H m,n p,q x

(a i , kA j ) p (b i , kB j ) q = 1 k H m,n p,q x 1 k (a i , A j ) p (b i , B j ) q
, k > 0 and the Mellin transform in [START_REF] Mathai | The H-function: Theory and Applications[END_REF], we obatin

E r k [Υ(r k )] = Γ(1 -δ) -1 Γ(1 + δ) c k β -δ k Λ k κ k j∈T λ j P δ j Λ j . (13)
Finally, plugging ( 13) into ( 14) yields the desired result after some manipulations. Remark 2: As shown in [START_REF] Khoshkholgh | Coverage analysis of max-SIR cell association in hetNets under nakagami fading[END_REF], the main task in deriving the coverage probability in cellular networks under the strongest-BS cell association criterion is to calculate Λ. In Table II, we show the coverage probability for the strongest-BS association criterion when various special cases of the Fox's H-function distribution are considered. Notably, ( 9) is instrumental in evaluating the impact of the number of tiers or their relative densities, transmit powers, and target SIR over generalized fading scenarios. This result complements existing valuable coverage studies of cellular networks over generalized fading [14, Proposition 1], [9, Corollary 1].

C. Coverage Analysis in Ad Hoc Networks

Ad hoc networks with short range transmission are, from an architecture perspective, similar to device-to-device (D2D) communication networks where Internet of Things (IoT) devices communicate directly over the regular cellular spectrum but without using the BSs. In ad hoc networks, the communication distance r k between the typical receiver and its associated transmitter in the k-th tier is assumed to be fixed and independent of the set of interfering transmitters and their densities.

Proposition 4: The coverage probability of ad hoc networks over the Fox's H fading channel is given by

C = T k=1 κ k c k Υ(r k ). ( 14
)
where Υ(r k ) is given in [START_REF] Trigui | Unified stochastic geometry modeling and analysis of cellular networks in LOS/NLOS and shadowed fading[END_REF]. We note that the coverage probability in ad hoc networks involves finite summation of Fox's H functions which can be efficiently evaluated [START_REF] Trigui | Unified analysis and optimization of D2D communications in cellular Networks over fading channels[END_REF]. Overall, the obtained analytical expressions are easier to compute than existing results [START_REF] Andrews | A tractable approach to coverage and rate in cellular networks[END_REF]- [START_REF] Di Renzo | Average rate of downlink heterogeneous cellular networks over generalized fading channels: A stochastic geometry approach[END_REF], [START_REF] Ding | Performance impact of LoS and NLoS transmissions in dense cellular networks[END_REF]- [START_REF] Lee | Spectral efficiency scaling laws in dense random wireless networks with multiple receive antennas[END_REF] that contain multiple nested integrals.

D. The Impact of Network Densification

In this section, we exploit the derived analytical framework to analyze the coverage scaling laws for single-antenna multitiers cellular and ad hoc networks. Assuming λ k = λ → ∞, k = 1, . . . , T , the coverage scaling laws are given in the following text.

1) Coverage Scaling Law in Cellular Networks: The coverage probability of single antenna-cellular networks with an unbounded path-loss model is invariant to the BS density λ. Specifically, we have

C U ,∞ = T k=1 ∞ 0 H v,u q,p+1 ξ, P k U dξ ξ 2 T j=1 P δ j 1+δξH v+1,u+2 q+2,p+3 ξ, P I U , (15)
which is obtained by letting λ → ∞ in Proposition 1 and resorting to the asymptotic expansion of the Fox's H function

H 1,1 1,1 (x; P δ ) ≈ x→∞ 1 δ
x -1 along with applying [37, Eq. (1.5.9)]. We note that (15) generalizes the SINR invariance property that has been revealed in some specific settings, e.g., [START_REF] Li | Success probability and area spectral efficiency in multiuser MIMO HetNets[END_REF], [START_REF] Andrews | A tractable approach to coverage and rate in cellular networks[END_REF], [START_REF] Di | Stochastic geometry modeling of coverage and rate of cellular networks using the Gil-Pelaez inversion theorem[END_REF], and [START_REF] Dhillon | Modeling and analysis of K-tier downlink heterogeneous cellular networks[END_REF]. Contrary to what the standard unbounded path-loss model predicts, the coverage probability under the bounded path-loss model scales with e -λ and approaches zero with increasing λ for general values of δ. This is readily shown in the following asymptotic coverage expression obtained by letting λ → ∞ in Proposition 2, as3

C B,∞ = T k=1 ∞ 0 e -λ T j=1 π P δ j δξ(Ψ1-Ψ2) × H v,u q,p+1 ξ, P k B H 1,1 1,1 T j=1 P δ j (1+δξΨ1) j∈T P δ j δξ(2Ψ1+Ψ2) , P δ ξ 2 T j=1 π P δ j δξ(Ψ 1 + Ψ 2 )
dξ.( 16)

Due to the complexity of the bounded model, its impact was only understood through approximations in [START_REF] Khoshkholgh | Coverage analysis of max-SIR cell association in hetNets under nakagami fading[END_REF] and [START_REF] Liu | Effect of densification on cellular network performance with bounded path-loss model[END_REF] and for fading scenarios with integer parameters. Thanks to our proposed unified approach, the impact of ultra densification can be scrutinized in the most comprehensive setting of multitier networks under the Fox's H fading channel.

2) Coverage Scaling Law in Ad Hoc Networks: In ad hoc networks, to the best of our knowledge, there exists no works that quantified the effect of densification over generalized fading channels. By exploiting the proposed analytical framework, the coverage scaling law in ad hoc networks is revealed in this paper. First, its is pertinent to remark that g k ∼ H-{(q, 0, p, q), P} can be assumed in the majority of fading distributions as shown in Table I. In this case, applying the asymptotic expansion of the Fox's H function [START_REF] Kilbas | H-Transforms: Theory and Applications[END_REF]Eq. (1.7.14)] H q,0 p,q (x) ∼ x

ν+ 1 2 ∆ exp -∆ x ρ 1/∆ to (14), we obtain C ≈ λ→∞ T k=1 κ k c k (λA) ν k + 1 2 ∆ k exp -∆ k λ A ρ k 1/∆ k , (17)
where [START_REF] Kilbas | H-Transforms: Theory and Applications[END_REF]Eq. (1.1.9)], and [37, Eq. (1.1.10)], respectively. In the special case of Gamma fading, i.e.,

A = T j=1 πr 2 k Γ(1-δ)Λj P -δ j (c k β k ) -δ , ∆ k = 1 + δ q j=1 B j k - p j=1 A j k -1 , ρ k = δ δ p j=1 (δA j k ) -δAj k q j=1 (δB j k ) -δBj k , and ν k = q j=1 b j k - p j=1 a j k + q j=1 B j k - p j=1 A j k + p-q 2 -1 are constants defined in [37, Eq. (1.1.8)],
g x k ∼ Gamma(m k , 1) ∼ H-{(1, 0, 0, 1), P}, it can be shown that ∆ k = 1, ρ k = 1, and ν k = m k -3 2 , which results in C ≈ λ→∞ T k=1 (λA) m k -1 Γ(m k ) exp (-λA) . (18)
It turns out that the coverage probability of ad hoc networks in arbitrary Nakagami-m fading (i.e., g x k ∼ Gamma(m k , 1),

Gamma Fading

C U = π C(δ) T k=1 λ k β -δ k Γ(m k +δ) Γ(m k )m δ k T j=1 λ j P δ j Γ(m j +δ) Γ(m j )m δ j . Generalized Gamma C U = π C(δ) T k=1 λ k β -δ k Γ(µ k ) δ-1 Γ µ k + 1 α k δ Γ µ k + δ α k T j=1 λ j P δ j Γ(µ j ) δ-1 Γ µ j + 1 α j δ Γ µ j + δ α j . Power of Nakagami-n (Rice) C U = π C(δ) T k=1 λ k β -δ k e -K R k ∞ t=0 K t R k Γ(mt+δ) Γ(t+1)Γ(mt) Ωt mt δ T j=1 λ j e -K R j P δ j ∞ t=0 K t R j Γ(mt+δ) Γ(t+1)Γ(mt) Ωt mt δ . Lognormal Fading C U = π C(δ) T k=1 λ k β -δ k N k n=0 wn10 δ( √ 2σ k un+µ k ) T j=1 λ j P δ j N j n=0 wn10 δ( √ 2σ j un+µ j )
.

Fisher-Snedecor Fading of the transmitter density λ. When T = 1, i.e., in singletier networks, the coverage probability is a product of an exponential function and a power function of order m -1.

C U = π C(δ) T k=1 λ k β -δ k m δ s k Γ(m k +δ)Γ(ms k -δ) Γ(ms k )Γ(m k )m δ k T j=1 λ j P δ j m δ s j Γ(m j +δ)Γ(ms j -δ) Γ(ms j )Γ(m j )m δ j . k = 1, . . . , T ) is
In the special case when T = m = 1, i.e., in single-tier ad hoc networks over Rayleigh fading channel, the coverage probability reduces to an exponential function.

IV. MULTI-ANTENNA VS. SINGLE-ANTENNA NETWORKS

A. Coverage Analysis

In multi-antenna networks, the analysis of the coverage probability is more difficult due to more complicated signal and interference distributions. However, we emphasize that, for several MIMO techniques, the associated post-processing signal power gain can include Gamma-type fading [START_REF] Lee | Spectral efficiency scaling laws in dense random wireless networks with multiple receive antennas[END_REF], [START_REF] Hunter | Transmission capacity of ad hoc networks with spatial diversity[END_REF], [START_REF] Di Renzo | Stochastic geometry modeling and systemlevel analysis of uplink heterogeneous cellular networks with multiantenna base stations[END_REF], [START_REF] Yu | A unified framework for the tractable analysis of multi-antenna wireless networks[END_REF] with g x ∼ Gamma(M, θ) where M is typically related to the number of antennas (e.g. M = N t , θ = 1 for maximum-ratio-transmission (MRT) and M = N t , θ = 1/N t for millimeter wave analog beamforming, where N t is the number of antennas at the transmitter) [START_REF] Yu | A unified framework for the tractable analysis of multi-antenna wireless networks[END_REF]Table II]. Hence, assuming that the signal power is gamma distributed in multi-antenna networks and recognizing that f gx (x) =

θ -1 Γ(M) H 1,0 0,1 x θ - (M -1, 1)
, then the Fox's H-based modeling of the coverage presented in Section III can be generalized to muti-anetnna networks analysis.

1) Multi-Antenna Cellular Networks:

We consider multiple-input-single-output (MISO) networks using MRT where the BSs in the k-th tier are equipped with N t k antennas. We assume that the channel power gain g x k for the desired signal is gamma distributed such that g x k ∼ Gamma(N t k , 1) [START_REF] Yu | A unified framework for the tractable analysis of multi-antenna wireless networks[END_REF]. As far as the interference distribution is concerned, we assume that g xi are identically distributed according to an arbitrary Fox's H distribution. Hence, Proposition 1 can be generalized to obtain the coverage probability in multi-antenna cellular networks with arbitrary interference, as

C = πδ T k=1 P k σ 2 k δ λ k β K Γ(N t k ) ∞ 0 η(ξ) ξ 2+δ H 0,1 1,1 ξ β k (2-N t k , 1) (1, 1) dξ, (19)
where

η(ξ) = H 1,1 1,1    T j=1 πλ j P δ j P k ξσ 2 k -δ 1+δξH v+1,u+2 q+2,p+3 ξ, P Ij U ,P δ    .
(20) Compared with existing approaches in [START_REF] Lee | Spectral efficiency scaling laws in dense random wireless networks with multiple receive antennas[END_REF]- [START_REF] George | Ergodic spectral efficiency in MIMO cellular networks[END_REF], which requires the calculation of N t k -1 derivatives of η(ξ) when g x k is gamma distributed as Gamma(N t k , 1), the framework in [START_REF] Filo | Stochastic geometry analysis of ultra-dense networks: Impact of antenna height and performance limits[END_REF] and [START_REF] Lee | Spectral efficiency scaling laws in dense random wireless networks with multiple receive antennas[END_REF] adds no computational complexity and thus preserves the tractability of single-antenna settings. We note that assuming a Gamma distribution for the interferers' power gain, i.e. g xj ∼ Gamma(χ j , φ j ), j ∈ {1, . . . , T }, is commonly encountered in multi-antenna networks [START_REF] Yu | A unified framework for the tractable analysis of multi-antenna wireless networks[END_REF], [START_REF] Yu | Coverage analysis for millimeter wave networks: The impact of directional antenna arrays[END_REF], [START_REF] Deng | Millimeter-wave deviceto-device networks with heterogeneous antenna arrays[END_REF]. In this case, we only need to modify the parameters of η(ξ) by replacing in [START_REF] Filo | Stochastic geometry analysis of ultra-dense networks: Impact of antenna height and performance limits[END_REF] the following equation

P I,j U = φ j Γ(χ j )
, φ j , (-χ j , 0, δ), (0, 1, , δ-1),

(1, 1, 1), (1, 1, 1, 1) . ( 21)

2) Multi-Antenna Ad Hoc Networks: The coverage probability of ad hoc networks for different multi-antenna transmission strategies for which g x k ∼ Gamma(M k , θ k ), k = 1, . . . , T is directly obtained from [START_REF] Khoshkholgh | Coverage analysis of max-SIR cell association in hetNets under nakagami fading[END_REF] as

C = T k=1 1 Γ(M k ) ×H 2,0 1,2   T j=1 πr 2 k λ j Γ(1 -δ)Λ j β δ k P -δ j (θ k ) -δ (1, δ) (0, 1), (M k , δ)   ,( 22
)
where Λ j accounts for Fox's H identically distributed interferences and is given in [START_REF] Trigui | Unified analysis and optimization of D2D communications in cellular Networks over fading channels[END_REF]. The coverage probability scaling law of multi-antenna ad hoc networks using MRT with N t k antenna at the k-th tier BS is obtained from applying ( 18) to [START_REF] Zhang | Downlink cellular network analysis with multi-slope path loss models[END_REF] as

C ≈ λ→∞ T k=1 (λA) Nt k -1 Γ(N t k ) exp (-λA) -→ λ-→∞ 0 (23)
where

A = T j=1 πr 2 k Γ(1-δ)Λj P -δ j β -δ k
. This last result reveals that, although the SIR increases in multi-antenna ad hoc networks, it will continue to drop to zero as the transmitter density increases.

B. The Impact of the Antenna Size

In this subsection, we consider multiple-input-single-output single-tier networks (i.e., T = 1) in which the BSs are equipped with N t antennas. Next we exploit the expressions and tools of the previous sections to derive the scaling laws for different multi-antenna networks including ad hoc, cellular, mmWave and networks with elevated BSs.

1) Antenna Scaling in Ad Hoc Networks: For the multiantenna case, the coverage expression in [START_REF] Zhang | Downlink cellular network analysis with multi-slope path loss models[END_REF] can be used to find the asymptotic scaling laws summarized as follows.

Proposition 5: Consider a multiple-input-single-output ad hoc network with N t transmit antennas such that lim λ→∞ Nt λ

1 δ = γ,
where γ ∈ [0, ∞], then the asymptotic coverage probability has the following scaling law

lim λ→∞ C =        0, γ = 0; H 1,0 1,1 T γ δ (1, δ) (0, 1) , γ ∈ R * + ; 1, γ = ∞. (24)
where γ = 0, ∈ R * + , ∞ stands for asymptotically sublinear, linear and super-linear scaling of N t and T = πr 2 Γ(1δ)β δ θ δ Λ.

Proof: Resorting to the Mellin-Barnes representation of the Fox's H-function in [START_REF] Zhang | Downlink cellular network analysis with multi-slope path loss models[END_REF], it follows that

C (a) = 1 2πj C Γ(N t + δs)Γ(s) Γ(N t )Γ(1 + δs) (T λ) -s ds (b) ≃ λ→∞ 1 2πj C T λ N δ t -s Γ(1 + δs) ds, (25)
where (a) follows from using [32, Eq. (2. Hence, based on Proposition 5, we evince that scaling the number of antennas linearly with the density does not prevent the SINR from dropping to zero for high BSs densities (as δ -1 = α/2 with α > 2) thereby hindering the SINR invariance property. Interestingly, when the number of antennas scales super-linearly with the BSs density, the coverage approaches a finite constant which is desirable since it guarantees a certain quality of service (QoS) for the users in the dense regime.

2) Antenna Scaling in Cellular Networks: Before delving into the analysis, it is important to recall that, in the single antenna case, the coverage probability under a practical bounded path-loss model drops to zero as λ → ∞ (see Section II.C). In the multi-antenna case, the asymptotic coverage scaling laws are summarized in the following proposition.

lim λ→∞ C =        0, ζ = 0; H 1,0 1,1 2π β ηζ (1, 1) (0, 1) , ζ ∈ R * + ; 1, ζ = ∞, (26)
where η = 2 -3α + α 2 and ζ = 0, ∈ R * + , ∞ stands for asymptotically sub-linear, linear and super-linear scaling of N t .

Proof: Due to the intricacy of L I in (47), we resort to an analytically tractable tight lower bound. Under the bounded path-loss model, the coverage probability in (46) involves the interference Laplace transform L I (s(1 + r) α ) = exp(2πλΘ(s(1 + r) α ), where we have

Θ(s) = E g ∞ r 1 -exp -sg(1 + t) -α ) tdt (a) ≤ sE[g] ∞ r t(1 + t) -α dt (b) ≈ s η (1 + r) 1-α (1 -r + αr) , (27)
where the inequality in (a) follows from the fact that 1e -x ≤ x, ∀x ≥ 0 and (b) holds since g has a unit mean. Note that when r becomes smaller, the inequality in [START_REF] Di Renzo | Stochastic geometry modeling and systemlevel analysis of uplink heterogeneous cellular networks with multiantenna base stations[END_REF] becomes tighter. This is typically the case in ultra-dense networks, where the closest distance to the origin tends to be infinitesimally small. Accordingly, by relabeling r ← r λ , we obtain

Θ s 1 + r λ α ≈ λ→∞ s η . (28)
Hence, the coverage probability can be obtained by merging [START_REF] Di Renzo | Stochastic geometry modeling and systemlevel analysis of uplink heterogeneous cellular networks with multiantenna base stations[END_REF] and (46) as

C ≈ β Γ(N t ) ∞ 0 e -2π λ η ξ ξ 2 H 0,1 1,1 ξ β (2 -N t , 1) (1, 1) dξ, (a) ≈ 1 Γ(N t ) H 2,0 1,2 2π λβ η (1, 1) (0, 1), (N t , 1) , (29)
where (a) follows from applying

∞ 0 f r (r)dr = 1 and [32, Eq. (2.3)]. As N t (λ) → λ→∞ ∞, we obtain C (b) ≃ λ→∞ H 1,0 1,1 2π λβ N t η (1, 1) (0, 1) , ( 30
)
where (b) follows along the same lines of [START_REF] Hunter | Transmission capacity of ad hoc networks with spatial diversity[END_REF]. The proof follows by resorting to the asymptotic expansions of the Fox' H function in [START_REF] Cheng | Coverage analysis for millimeter wave cellular networks with imperfect beam alignment[END_REF] when ζ = Nt λ is near zero [37, Eq. (1.7.14)] and infinity [START_REF] Kilbas | H-Transforms: Theory and Applications[END_REF]Eq. (1.8.7)].

The obtained result in [START_REF] Yu | Coverage analysis for millimeter wave networks: The impact of directional antenna arrays[END_REF] allows us to conclude that monotonically increasing the per-user coverage performance by means of ultra-densification is theoretically possible through the deployment of multi-antenna BSs. Specifically, [START_REF] Yu | Coverage analysis for millimeter wave networks: The impact of directional antenna arrays[END_REF] unveils that scaling linearly the number of antennas with the BS density constitutes a solution for the coverage drop in traditional dense networks.

C. The Impact of Antenna Gain in mmWave Networks

In multiple-input-single-output mmWave networks, the channel gain for the signal g x follows a gamma distribution g x ∼ Gamma(N t , 1

Nt ), where N t is the number of antennas at the BS [START_REF] Yu | A unified framework for the tractable analysis of multi-antenna wireless networks[END_REF]- [START_REF] Yu | Coverage analysis for millimeter wave networks: The impact of directional antenna arrays[END_REF]. As for the interference received at the typical user, the total channel gain is the product of an arbitrary unit mean small-scale fading gain g [START_REF] Di Renzo | Stochastic geometry modeling and systemlevel analysis of uplink heterogeneous cellular networks with multiantenna base stations[END_REF], [START_REF] Yu | Coverage analysis for millimeter wave networks: The impact of directional antenna arrays[END_REF] and the directional antenna array gain G( d λt θ x ), where d and λ t are the antenna spacing and wavelength, respectively, and θ x is a uniformly distributed random variable over [-1, 1]. An example of antenna pattern based on the cosine function is given by [START_REF] Yu | Coverage analysis for millimeter wave networks: The impact of directional antenna arrays[END_REF], [START_REF] Deng | Millimeter-wave deviceto-device networks with heterogeneous antenna arrays[END_REF]

G(x) = cos 2 πNt 2 x , |x| ≤ 1 Nt ; 0, otherwise. (31)
In dense mmWave network deployments, it is reasonable to assume that the link between any serving BS and the user is in line-of-sight (LOS). Mathematically, the probability of being in a LOS propagation can be formulated as p(r) = e -τ r , where τ is the blockage parameter determined by the density and average size of the spatial blockage [START_REF] Atzeni | Downlink cellular network analysis with LOS/NLOS propagation and elevated base stations[END_REF], [START_REF] Yu | A unified framework for the tractable analysis of multi-antenna wireless networks[END_REF]. Accordingly, [START_REF] Yu | A unified framework for the tractable analysis of multi-antenna wireless networks[END_REF] can be derived, based on the cosine antenna pattern and the blockage model, as

Θ(s) ≤ λ t s πdN t ∞ r t e -τ t (1 + t) αL dt+ ∞ (1+r) α L α N -1 t(1-e -τ t ) (1 + t) αN dt × π 0 cos 2 x 2 dx (a) = λ→∞ λ t e τ s 2dN t (P + J (τ )) , (32)
where (a) follows along the same lines of ( 28) with α L (α N ) is the path-loss exponent of the LOS (NLOS) link, P = αN -αL-1

(1-αL)(αN -2) , J (τ ) = (αL-1+τ )Eα L -1 (τ ) αL-1 - (αN -1+τ )Eα N -1 (τ ) αN -1
, with E ν (•) denoting the Exponential Integral function [START_REF] Gradshteyn | Table of Integrals, Series and Products[END_REF]. Then, using [START_REF] Mathai | The H-function: Theory and Applications[END_REF] and following the same steps as in [START_REF] Cheng | Coverage analysis for millimeter wave cellular networks with imperfect beam alignment[END_REF], we obtain

C ≈ 1 Γ(N t ) H 2,0 1,2
πλλ t e τ β d (P +J (τ )) (1, 1) (0, 1), (N t , 1)

.

Hence, the coverage scaling laws in mmWave networks are given in the following proposition. Proposition 7: In mmWave networks in which lim λ→∞ λ t = 0, and lim

lim λ→∞ C =        0, ρ = 0; H 1,0 1,1 π βe τ (P+J (τ )) dρ (1, 1) (0, 1) , ρ ∈ R * + ; 1, ρ = ∞. (34)
Proof: The proof is similar to those of Propositions 5 and 6. The obtained result unveils that the scaling laws derived for mmWave cellular networks are similar to those obtained for legacy frequency bands (see Proposition 6). Specifically, maintaining a linear scaling between the density of BSs and the number of antennas is sufficient to prevent the SINR from dropping to zero and to guarantee a certain QOS to the UE. In addition, the optimal coverage can be achieved by linearly scaling the number of antennas and the mmWave carrier frequency, which reduces both cost and power consumption. This result provides evidence that moving toward higher frequency bands may be an attractive solution for high capacity ultradense networks.

Achieving optimal coverage rely on determining the optimal scaling factor below which further densification becomes destructive or cost-ineffective. This operating point will depend on properties of the channel power distribution and pathloss and is of cardinal importance for the successful deployment of ultra-dense networks.

Corollary 1 (Optimal Scaling Factor in Dense mmWave Networks): Capitalizing on Proposition 7, the optimal scaling factor ρ that prevents the outage drop in dense mmWave networks is given by

ρ = N t f c λ (a) = πβe τ (P + J (τ )) d , (35)
where f c is the mmWave carrier frequency and (a) follows from recognizing that H 1,0 1,1 x

(1, 1) (0, 1) = U (1 -x), where U (x) = 1, x ≥ 0; 0, otherwise. stands for the Heaviside function [START_REF] Gradshteyn | Table of Integrals, Series and Products[END_REF]. In particular, [START_REF] Di Renzo | System-level analysis/optimization of cellular networks with simultaneous wireless information and power transfer: Stochastic geometry modeling[END_REF] unveils that under a full-blockage scenario (i.e., τ → ∞), a super-linear scaling of N t f c is required to offset the coverage drop. However, in the noblockage regime (i.e., τ → 0), only a linear scaling is needed. Using this framework, enhanced antenna models can be considered to investigate the impact of beam alignment errors on the coverage probability of mmWave dense networks [START_REF] Cheng | Coverage analysis for millimeter wave cellular networks with imperfect beam alignment[END_REF].

D. The Impact of Antenna Height in 3D Networks

The vast majority of spatial models for cellular networks are usually 2D and ignore the impact of the BS height. Recent papers have, however, tackled this issue and have highlighted the importance of taking this parameter into account to appropriately estimate the network performance [START_REF] Ding | Performance impact of LoS and NLoS transmissions in dense cellular networks[END_REF], [START_REF] Atzeni | Downlink cellular network analysis with LOS/NLOS propagation and elevated base stations[END_REF], [START_REF] Filo | Stochastic geometry analysis of ultra-dense networks: Impact of antenna height and performance limits[END_REF]. In 3D cellular networks, the distance between a BS and the typical UE can be expressed as √ h 2 + r 2 , where h is the absolute antenna height difference between the serving BS and the typical UE. Adapting the coverage probability expression in (46) to the 3D context results in an interference distribution whose Laplace transform is of the form

L I (s(r 2 + h 2 ) α 2 ) = exp(2πλΘ(s(r 2 + h 2 ) α 2 )) where Θ(s) ≤ sE[g] ∞ r t(h 2 + t 2 ) -α 2 dt ( 36
) (a) = s α -2 h 2 + r 2 1-α/2 .
By employing the change of variable x ← λr 2 , we obtain

L I s r λ + h 2 α 2 ≈ λ→∞ e 2πλsh 2
α-2 . Hence, the coverage probability in 3D multi-antenna cellular networks can be formulated similar to [START_REF] Cheng | Coverage analysis for millimeter wave cellular networks with imperfect beam alignment[END_REF] and [START_REF] Yu | Coverage analysis for millimeter wave networks: The impact of directional antenna arrays[END_REF] as

C ≃ λ→∞ H 1,0 1,1 2π λh 2 β N t (α -2) (1, 1) (0, 1) . ( 37
)
The obtained analytical expression for the coverage probability unveils the impact of the antenna height coupled with other design parameters. Recent works [START_REF] Ding | Performance impact of LoS and NLoS transmissions in dense cellular networks[END_REF]- [START_REF] Filo | Stochastic geometry analysis of ultra-dense networks: Impact of antenna height and performance limits[END_REF] proposed to maintain the SINR invariance of the coverage probability by lowering the height of the BSs. Based on (37), we evince that the SINR invariance of the coverage probability in 3D networks can be maintained by enforcing a super-linear scaling with the number of antennas. Corollary 2 (Optimal Scaling Factor in Dense 3D Networks): The optimal scaling factor for the successful deployment of dense 3D networks is

N t λ = 2πβh 2 α -2 , (38)
which exploits [START_REF] Kilbas | H-Transforms: Theory and Applications[END_REF] and follows along the same lines of Corollary 1. In particular, the last result shows that the coverage probability monotonically decreases as the BS density increases, if

lim λ→∞ N t λ = ζ ∈ R * + , and if h > (α-2)ζ
2π . Interestingly, it is possible to counteract this decay by tuning the antenna number according to BS density in order to maintain the per-user coverage performance.

V. NUMERICAL RESULTS

In this section, we substantiate our theoretical coverage expressions and scaling laws using system level simulations. Unless otherwise stated, the noise power is set to σ 2 = -70 dBm and the path loss is given by L(r) = r -α for powerlaw unbounded model and L(r) = (1 + r) -α for physically feasible bounded model, with α = 3.

The performance comparisons between strongest-BS-and closest-BS-association-based two-tier (i.e., T = 2) cellular networks with unbounded power-law path-loss model are illustrated in Fig. 1. Overall, the strongest-BS strategy provides significant performance gain over the closest-BS strategy especially in low density range. Furthermore, depending on the target SINR thresholds, the effect of increasing densification is beneficial while, in some cases, tends to be negligible. Indeed, since using the power law model, the coverage saturates to a non-zero finite constant in the limit of λ 1 , λ 2 → ∞.

Fig. 2 plots the scaling of the coverage probability with BS densities for both bounded and unbounded path-loss models. Analytical and experimental curves are in full agreement. It shows that the unbounded model (i.e., r -α ) guarantees a certain QoS or coverage for the users in the dense regime by preventing the SINR form dropping to zero. However, this SINR-invariance property is unattainable because the unbounded model is physically impracticable and unrealistic. The figure also highlights the diminishing gains achieved with the more realistic bounded, as anticipated by Eq. ( 16). In this case, new densification strategies are required to prevent the SINR from dropping to zero and avoid the densification plateau. This will be discussed later in Fig. 6. Fig. 3 shows the scaling of the SIR coverage probability of ad hoc networks against the transmitter density for various common fading distributions stemming from the general Fox's H fading model. In particular, we corroborate the result of Eq. [START_REF] Ding | Performance impact of LoS and NLoS transmissions in dense cellular networks[END_REF] stipulating that increasing the transmitter density degrades the coverage probability in ad hoc networks, and that the coverage probability is a product of an exponential function and a polynomial function of order T (m -1) of the transmitter density. Moreover, the multi-path fading model has a less noticeable impact on the coverage performance than the path-loss model (cf. Fig. 2) and the number of tiers.

Fig. 4 shows the SIR outage probability of cellular networks for an unbounded path-loss model versus the antenna size when assuming that the interferers' power gain follows a Gamma distribution, i.e., g ∼ Gamma(χ, φ). Fig. 4 demonstrates that increasing the antenna size keeps improving the coverage probability, less significantly, however, as the number of antennas grows large. Fig. 5 illustrates the SIR coverage probability of a two-tiers cellular network over Rician fading with closed-BS association obtained from Proposition 3 for different Rician power factors. We observe a substantial increase of the coverage probability only in the non-asymptotic regime, i.e., K 1 = K 2 . Moreover, we observe that the two extreme regimes of pure fading channel with (K → 0) and pure LOS propagation (K → ∞) achieve worse coverage performance. Fig. 6 shows the scaling of the coverage probability in ad hoc networks against the transmitter density for different scaling rates of the number of antennas; super-linear, linear, sub-linear, and constant (i.e, single antenna). We notice that the coverage decreases with the density for the single antenna case, as anticipated in [START_REF] Ganti | Asymptotics and approximation of the SIR distribution in general cellular networks[END_REF], and also when the number of antennas is scaled sublinearly or linearly with the density, as predicted by Proposition 5. We observe that a super-linear scaling of the number of antennas with the BS density is required to prevent the SIR from dropping to zero, and thereby restore the SIR invariance property. The impact of the BS height on the coverage probability is illustrated in Fig. 7. As predicted in Section IV.D, we note that a linear scaling of the number of antennas is required to maintain a non-zero SINR for low value of h. When the BS height increases, the coverage probability decreases due to the increase of the path-loss and the linear scaling becomes insufficient.

VI. CONCLUSION

By leveraging the properties of Fox's H random variables, we developed a unifying framework to characterize the coverage probability of heterogeneous and muti-antenna networks under both the closest-BS and the strongest-BS cell association strategies. We studied the impact of BS densification on the coverage performance both under bounded and unbounded path loss models. By direct inspection of the obtained analytical framework, we have been able to derive exact closedform formulations and scaling laws of the coverage probability for two typical network models, i.e., heterogeneous and multiantenna cellular and ad hoc networks, while incorporating generalized fading distributions. The obtained results encompass insightful relationships between the BS density and the relative antenna array size, gain and height, showing how the coverage can be maintained whilst increasing the network density. The insights provided in this work are of cardinal importance for optimally deploying general ultra-dense networks.

VII. APPENDIX A: PROOF OF PROPOSITION 1

With the closest-BS association strategy, the coverage probability is given by

C T k=1 θ k P (SINR k > β k ) , (39)
where θ k denotes the association probability and is expressed as

θ k = λ k j∈T λj P δ j
, and P j = Pj P k . Using [10, Theorem 1] and [START_REF] Trigui | Unified analysis and optimization of D2D communications in cellular Networks over fading channels[END_REF]Eq. (39)] and applying the Fox's H-transform in [32, Eq. (2.

3)], the coverage probability under unbounded path-loss model, denoted as C U , is given by

C U (r k ) = ∞ 0 1 √ ξ L -1 1 √ s H u,v p,q {f (t); P} (sξ); s; β k e -σ 2 k ξ r α k P k j∈T L Ij ξ r α k P k dξ, (40)
where

f (t) = √ tJ 1 2 √ stξ , H u,v p,q {f (t); P} (s) is the Mellin transform [32, Eq. (2.3)], J 1 (x) = H 1,0 0,2 x 2 4 ; (1, 1, 1 2 , -1 2 , 1, 1)
is the Bessel function of the first kind [38, Eq. (8.402)], and L -1 is the inverse Laplace transform. Moreover L Ij , in (46), is the Laplace transform of the aggregate interference from the j-th tier evaluated as

L Ij (s) = exp(2πλ j Θ(s)), (41)
where

Θ | Hj =y (s) (a) = ∞ P j P k 2 α r k 1 -exp -syP j r -α rdr (b) = syP j α ∞ P k P j r -α k x -2 α e -syx 1 F 1 (1, 2, syx) dx (c) = syP j 2 F 2 1, -2 α + 1; 2; -2 α + 2; -syP k r -α k P k Pj 2 α -1 r (α-2) α 2 α -1 , (42)
where (a) follows from the probability generating functional [START_REF] Dhillon | Modeling and analysis of K-tier downlink heterogeneous cellular networks[END_REF], [START_REF] Haenggi | Stochastic Geometry for Wireless Networks[END_REF], while relabeling x as r -α k and (1e -x )/x = e -x 1 F 1 (1, 2; x) is applied in (b), and (c) follows from applying

x β-1 e -cx 1 F 1 (a, b, cx) = (1 -δ) H u,v p,q {g(t); P} (ξ) ,

where g(t) = t 2 F 2 1, 1 -δ; 2; 2 -δ; -ξtr -α k = tH 1,2 2,3 (t; P 1 ), P 1 = (1 -δ, ξ(r 2 k ) -α 2 , (0, δ), (0, -1, δ -1), 1 2 , 1 3 ), and p F q (•) is the generalized hypergeometric function of [START_REF] Gradshteyn | Table of Integrals, Series and Products[END_REF]Eq. (9.14.1)]. In (43), in particular, we first take the expectation over the interferers' locations and then average over the Fox's H distributed channel gains, which is in the reverse order compared to the conventional derivations in [START_REF] Andrews | A tractable approach to coverage and rate in cellular networks[END_REF]- [START_REF] Dhillon | Modeling and analysis of K-tier downlink heterogeneous cellular networks[END_REF]. The reason behind this order swapping is that the Fox's H fading model is more complicated than the conventional exponential model, and therefore averaging over it in a latter step preserves the analytical tractability.

Finally, applying [32, Eq. (1.58)], the Mellin transform [32, Eq. (2.3)], and the inverse Laplace transform of the Fox's Hfunction [START_REF] Mathai | The H-function: Theory and Applications[END_REF]Eq. (2.21)] given by L -1 {x -ρ H u,v p,q (x; P); x; t} = t -ρ-1 H u,v p+1,q 1 t ; P l , (44)

where P l = (κ, c, (a, ρ), b, (A, 1), B), the desired result is obtained after applying the Fox's H reduction formulae in [START_REF] Mathai | The H-function: Theory and Applications[END_REF]Eq. (1.57)]. The coverage probability over Fox's H-fading 4 for a receiver connecting to a k-th tier BS located at x k is given by Recall that the pdf of the link's distance r k is given by f r k (x) = 2πλ k θ k x exp -j∈T πx 2 λ j P δ j [START_REF] Li | Success probability and area spectral efficiency in multiuser MIMO HetNets[END_REF], [START_REF] Haenggi | Stochastic Geometry for Wireless Networks[END_REF]. Then recognizing that exp(-x) = H 1,0 0,1 (x; 1, 1, 0, 1, 1, 1) [32, Eq. (1.125)] in (45), we apply [START_REF] Mathai | The H-function: Theory and Applications[END_REF]Eq. (2.3)] to obtain the average coverage probability in (3) after some manipulations.

C U (r k ) = ∞ 0 1 ξ 2 H v,

VIII. APPENDIX B: PROOF OF PROPOSITION 2

The proof of Proposition 2 relies on the same approach adopted in Appendix A, yielding where rearranging [START_REF] Trigui | Unified analysis and optimization of D2D communications in cellular Networks over fading channels[END_REF]Eq. (39)] after carrying out the change of variable relabeling (1 + x) -α as x, we have L Ij (ξ) = exp -πδλ j ξ (1+r k ) 2-α (1 -δ) H u,v p,q {g 1 (t); P 1 } (ξ) -

C B (r k ) = ∞ 0 1 √ ξ L -
(1 + r k ) 1-α 1 -δ 2 H u,v p,q {g 2 (t); P} (ξ) , (47)
where g 1 (t) = t 2 F 2 (1, 1 -δ; 2; 2 -δ; -ξt(1 + r k ) -α ) and g 2 (t) = t 2 F 2 1, 1 -δ 2 ; 2; 2 -δ 2 ; -ξt(1 + r k ) -α . Finally applying the Mellin transform in [START_REF] Mathai | The H-function: Theory and Applications[END_REF]Eq. (2.3)] and plugging the obtained result into (46), Proposition 2 follows after some manipulations.

IX. APPENDIX C: PROOF OF PROPOSITION 3

Based on [START_REF] Trigui | Unified analysis and optimization of D2D communications in cellular Networks over fading channels[END_REF], the Laplace transform of the aggregate interference from tier j under the max-SINR association strategy is evaluated as L Ij (ξ) = exp -πλ j ξ δ Γ ( 1

where (a) follows from substituting exp(-x) = H 1,0 0,1 (x; 1, 1, -, 0, -, 1) [32, Eq.

(1.125)] and applying the transformation H u,v p,q x (a i , kA j ) p (b i , kB j ) q = 1 k H u,v p,q x 1 k (a i , A j ) p (b i , B j ) q . Finally applying [32, Eq. (2.

3)] yields the strongest-BS based coverage probability as shwon in Proposition 3.

  formulated as the product of an exponential function and a polynomial function of order T (max m k k=1,...,M -1)

  t ) δs . The proof follows by recognizing the Fox's H function definition in [32, Eq. (1.2)], along with its asymptotic expansions near zero [37, Eq. (1.7.14)] and infinity [37, Eq. (1.8.7)].

Proposition 6 :

 6 In multiple-input-single-output cellular networks with N t transmit antennas such that lim λ→∞ Nt λ = ζ, where ζ ∈ [0, ∞], the asymptotic coverage probability under a bounded path-loss model has the following scaling law

  where ζ ∈ [0, ∞], the asymptotic coverage probability has the following scaling law

Fig. 1 .

 1 Fig. 1. Coverage probability vs. the BS density λ for multi-tier cellular networks with T = 2 over arbitrary Nakagami-m fading with m 1 = 1.5, m 2 = 2.5, P 1 = 50 W, and P 2 = 1 W.

Fig. 2 .

 2 Fig. 2. Coverage probability and scaling laws vs. the BS density λ for multitier cellular networks over Nakagami-m fading with T = 2, P 1 = 50 W, and P 2 = 1 W.

Fig. 3 .

 3 Fig. 3. Coverage probability in ad hoc networks vs. the transmitter density λ when β = 0 dB.

Fig. 4 .Fig. 5 .

 45 Fig. 4. Outage probability in MISO cellular networks assuming MRT vs. the number of antennas at the transmitter with λ = 10 -3 and φ = 1.

Fig. 6 .

 6 Fig. 6. Coverage probability of MISO ad hoc networks vs the BS density λ for different scaling of the number of antennas.

Fig. 7 .

 7 Fig. 7. Coverage probability of MISO cellular networks vs the BS height h for different scaling of the antenna number Nt.

x β β 2 F 2

 22 (b -a, β, b, β + 1, -cx). Hence, we obtain L Ij (ξ) = exp(-2πλ j E Hj [Θ(s)]) = exp -πδλ j ξr 2-α k

  δ = 2 α , and the parameter sequencesP k U = κβ k , 1 cβ k , 1-b, (1-a, 1), B, (A, 1) , andP I U = κ c 2 , 1 c , (1-b-2B, 0, δ), (0, 1-a-2A, -1, δ -1), (B, 1, 1), (1, A, 1, 1) .

2 k

 2 -δ) E[H δ ] , where E[H δ ] is the Mellin transform of the Fox's-H function obtained as E[H δ ] = Λ[START_REF] Mathai | The H-function: Theory and Applications[END_REF] Eq. (2.8)]. Then following the analytical steps as in Appendix A, we obtain πλ j P δ j ξ δ Γ(1 -δ)Λ j dξ

TABLE II COVERAGE

 II PROBABILITY OF SOME WELL-KNOWN FADING CHANNEL MODELS BASED ON THE STRONGEST-BS ASSOCIATION

	Instantaneous Fading Distribution	Coverage Probability C U

The α-µ distributions can be attributed to exponential, one-sided Gaussian, Rayleigh, Nakagami-m, Weibull and Gamma fading distributions by assigning specific values for α and µ.

[START_REF] Dhillon | Modeling and analysis of K-tier downlink heterogeneous cellular networks[END_REF] showed that self-interference dominates noise in typical heterogeneous networks under strongest-BS association. Therefore, we ignore noise in the rest of this section.

Using the Mellin-Barnes integral representations of Ψ 1 and Ψ 2[START_REF] Mathai | The H-function: Theory and Applications[END_REF], we can easily show that Ψ 1 -Ψ 2 > 0.

We dropped the index i from Fox's H-distribution {O i , P i } for notation simplicity.