Spherically Restricted Random Hyperbolic Diffusion.

This paper investigates solutions of hyperbolic diffusion equations in R 3 with random initial conditions. The solutions are given as spatial-temporal random fields. Their restrictions to the unit sphere S 2 are studied. All assumptions are formulated in terms of the angular power spectrum or the spectral measure of the random initial conditions. Approximations to the exact solutions are given. Upper bounds for the mean-square convergence rates of the approximation fields are obtained. The smoothness properties of the exact solution and its approximation are also investigated. It is demonstrated that the Hölder-type continuity of the solution depends on the decay of the angular power spectrum. Conditions on the spectral measure of initial conditions that guarantee short- or long-range dependence of the solutions are given. Numerical studies are presented to verify the theoretical findings.

1. Introduction

Numerous environmental, biological and astrophysical applications require the modelling of changes in data on the unit sphere or in the space [1,2,3,4,5,6]. One of the conventional tools for such modelling is stochastic partial differential equations (SPDEs)—see, for example [1,5,7,8], and the references therein. Random fields that are solutions of such SPDEs often exhibit dynamics dependent on initial conditions. Investigating the properties of these random fields is important for theoretical insight and practical applications.

SPDEs on surfaces and Riemannian manifolds have found numerous applications to problems in cosmology, physics, biology, fluid dynamics and pattern formation on surfaces, just to mention a few. See [7,9,10,11] and the references therein.

Random fields on a sphere, one of the simplest two-dimensional manifolds, have been used as a standard model in the astrophysical and cosmological literature in the last several decades [3,4,8,12]. NASA and ESA space missions [4] obtained very detailed measurements of Cosmic Microwave Background radiation (CMB), which are interpreted as a realisation of a spherical random field superimposed on an underlying signal of large-scale acoustic waves in plasma near the time of recombination. The theory of the standard inflation scenario uses a Gaussian model for the density fluctuation of this field [3,4,6]. Several new cosmological models were proposed using non-Gaussian assumptions and employed sophisticated statistical tests to justify possible departures from Gaussianity. The understanding of changes in CMB temperature fluctuations is important to predict future cosmological evolution and accurately reconstruct past states of the Universe. It also can help in the estimation and statistical inference of physical parameters obtained from the CMB data. SPDEs on the sphere can be used to model changes in CMB temperature fluctuations; see [7,8].

The pronounced spectral peaks at very large wavelengths in CMB temperature data are evidence of acoustic waves that were seeded by earlier superluminal inflation, with remnant coherent waves remaining until last scattering of photons and recombination of atoms around 340,000 years after the big bang [6,13]. In the plasma Universe there was chaotic mixing, but it is problematic to represent the underlying particle kinematics as standard Brownian motion, which leads to the standard diffusion equation for particle densities. Under standard diffusion, density disturbances have unbounded propagation speeds, which is unacceptable in relativistic cosmological contexts wherein remnant structures that are coherent over space-like domains have not been smeared away by diffusion. Reimberg [14] has directly modelled a sequence of photon–electron collisions backwards in time from the last scattering, with random changes of direction, and with the same distance travelled over equal time steps—unlike in usual random flight theory. Giona has developed a Feynman–Kac stochastic dynamics by which a particle undergoes a succession of collisions with speed-limited jumps. Consequently, the diffusion coefficients decrease as a power of the Lorentz–Fitzgerald contraction factor of an inertial observer. This has the interesting consequence that different observers may disagree on whether a process is deterministic or stochastic. Non-trivial continuous relativistic Markov processes on position space are simply not possible [15]. A simpler alternative phenomenological model is effected by replacing the standard diffusion equation by the simplest hyperbolic diffusion equation that has a variable but bounded speed of propagation. Ali and Zhang [16] recast the hyperbolic diffusion equation as a Lorentz-invariant Liouville conservation equation in one time and four space dimensions, before restricting to be Ali and Zhang then retain the second law of thermodynamics but only as a reaction-diffusion equation in 1+4 dimensions. Section 8 shows how information entropy may decrease by a small amount during the propagation of a point source by hyperbolic diffusion, whereas the overall increase is much larger.

The so-called Cattaneo hyperbolic diffusion equation [17,18] has been used to explain outcomes of heat conduction experiments in liquid in the super-fluid state [19,20], and solid and solid at very low temperatures. In these materials, and in some nanotubes and other graphite structures [21], heat energy propagates as a “second sound” wave mediated by phonons, at a propagation speed of around one-tenth the normal speed of sound. Since the experiments on graphite have been conducted on nano-scale samples, we expect that waves of second sound could likewise be detected in the spherical surface of a or larger fullerene ball, as formulated in our previous paper on hyperbolic diffusion on a spherical surface [8].

SPDEs on have been extensively studied. However, SPDEs on manifolds have only recently attracted a great deal of attention [7,8,22,23]. The results in these papers demonstrate that the continuity properties of solutions and the convergence rates of approximations to solutions are determined by decay rates of the angular power spectrum of initial random conditions. This article continues studies of solutions of SPDEs on the sphere. However, in contrast to the above publications that directly model spherical random fields using Laplace or Laplace–Beltrami operators on the sphere, we employ another approach. Namely, we consider the restriction of the stochastic hyperbolic diffusion in to the unit sphere. Compared to the available literature this approach is more consistent with real CMB observations that exist in space but are measured only on . From a mathematical point of view, additional investigations are required to show that solutions of known models on the sphere admit physically meaningful extensions to that are consistent with observations. By its construction, our model directly provides this consistency. The proposed model may find new applications for the next generation of CMB experiments, CMB-, which will be collecting observations. A very detailed discussion of SPDEs on manifolds and their physical and mathematical justification for CMB problems can be found in [8]. The hyperbolic diffusion equation prohibits the superluminal propagation of density disturbances, which is an unwanted feature of pure diffusion models over super-galactic distances. In addition, the linear hyperbolic diffusion equation, expressed in terms of co-moving material space coordinates and conformal time coordinate, is a good approximation to the field equation of a scalar field minimally coupled to an expanding Robertson–Walker space-time. However, speed-limited diffusion raises some interesting questions about the dynamics of Shannon entropy. For physical concentrations governed by linear or nonlinear heat diffusion equations of parabolic type, Shannon entropy is fully analogous to thermodynamic entropy and it increases monotonically [24]. It will be explained that at low wave numbers, the hyperbolic diffusion equation behaves as a dissipative diffusion equation. However, above some cut-off wave number it behaves as a bi-directional wave equation which has increasing entropy when twin pulses separate but has decreasing entropy when pulses approach each other constructively.

This paper is organised as follows. Section 2 presents definitions and results about spatial-temporal random fields in . It also introduces hyperbolic diffusion equations with random initial conditions and their solutions. Section 3 investigates the spatial-temporal hyperbolic diffusion field on the unit sphere. The Hölder-type continuity of the exact solution of the spatial-temporal hyperbolic diffusion field on the sphere is investigated in Section 4. In Section 5 we study the dependence structures of the spherical hyperbolic diffusion random fields. Section 6 obtains the mean-square convergence rate to the diffusion field in terms of the angular power spectrum. Section 7 provides some numerical results. Finally, Shannon entropy behaviour is discussed in Section 8, followed by some conclusions.

In the following sections we will use the symbol C to denote constants that are not important for our exposition. The same symbol may be used for different constants appearing in the same proof.

2. Spatial Random Hyperbolic Diffusion

This section reviews the basic theory of random fields in and introduces a hyperbolic diffusion with random initial conditions. Then, the solution of the diffusion equation is derived and analysed.

We consider the hyperbolic diffusion equation where is a random field satisfying the random initial conditions: where is the Laplacian in and the random field , defined on a suitable complete probability space , is assumed to be measurable, mean-square continuous, wide-sense homogeneous and isotropic with zero mean and the covariance function .

The covariance function has the following representation: for some bounded, non-negative measures on and on , such that See [25], pp. 1–5 and [26], pp. 10–15 for more details.

Then, there exists a complex-valued orthogonally scattered random measure such that for every , the field itself has the spectral representation Let , , , , , be complex spherical harmonics defined by the relation where are the associated Legendre polynomials with indices l and m. For spherical harmonics it holds that where the symbol * denotes the complex conjugation and is the Kronecker delta function. The addition formula for spherical harmonics gives The Bessel function of the first kind of order is defined by where is the Gamma function.

It admits the following representation by the Poisson integral, see (10.9.4) in [27]: By the addition theorem for Bessel functions(e.g., [26], p. 14), where is a family of random measures on , such that The stochastic integrals in (3) and (4) are viewed as integrals with the structural measures F and G, correspondingly.

Let us consider the initial conditions of the form: where is the Dirac delta function.

Let be the fundamental solution (or the Green’s function) of the initial-value problem (1) and (6), and let be its Fourier transform.

The following theorem derives the Fourier transform Contrary to many models in the literature, for the initial-value problem (1) and (2) it can be explicitly written in terms of elementary functions. Later, this result will be used to obtain the solution and its covariance function.

The Fourier transform (7) of the initial-value problem (1) and (6) is given by the formula where

The Fourier transform (7) is the solution of the initial-value problem The characteristic equation for the ordinary differential equation in (11) is with the roots Therefore, the general solution of the ordinary differential equation in (11) has the form where are some functions that do not depend on t and are given by (12). From the initial conditions in (11) we obtain the system of equations to find these functions which results in Thus, by (12) and (13) the solution of the initial-value problem (11) is The theorem is proved. □

The function (8) is radial. That is, there exists a function

[17]. At higher wave numbers, it can easily be seen from the one-dimensional solutions that the Fourier components may be viewed as travelling waves but with exponentially decaying amplitude. At high wave numbers, the governing PDE may be regarded as a damped wave equation.

Let us denote such that

It holds that and

It follows from (13)–(15) that for it holds Note that for and because and if . Thus,

As , one obtains the upper bound (19) from the representation (17) for . □

The following theorem provides the solution of the initial-value problem and its covariance function in terms of the Fourier transform As the explicit expression of in terms of elementary functions is given in Theorem 1, it can be used to obtain an explicit expression for the solution and then easily investigate various properties of

The solution (1) and (2) can be written as the convolution The covariance function of the spatio-temporal random field (20) is

Notice that where is given by (8), assuming that the random initial condition has the spectral measure F, such that Under the condition (22), the stochastic integral (20) exists in the -sense.

By Lemma 1 the function can be bounded by a constant which depends only on t. Noting that we obtain (22). The representation (21) immediately follows from (20) and the orthogonality of . □

3. Spherical Random Hyperbolic Diffusion

In this section we investigate a restriction of the spatial-temporal hyperbolic diffusion field from Section 2 to the unit sphere.

Consider the sphere in the three-dimensional Euclidean space with the Lebesgue measure

A spatio-temporal spherical random field defined on a probability space is a stochastic function

We consider a real-valued spatio-temporal spherical random field T with zero mean and finite second-order moments and being continuous in the mean-square sense (see, e.g., Marinucci and Peccati [3] for definitions and other details). Under these conditions, the zero-mean random field T can be expanded in the mean-square sense as the Laplace series [25]: where the functions represent the spherical harmonics and the coefficients are given by the formula We assume that the field T is isotropic (in the weak sense), that is, , and for every , the group of rotations in . This is equivalent to the condition that the covariance function depends only on the angular distance between the points and on for every .

The field is isotropic if and only if Hence,

The functional series is called the angular time-dependent power spectrum of the isotropic random field .

We can define a covariance function between two locations with the angular distance at times t and by where is the l-th Legendre polynomial.

If is a zero-mean isotropic Gaussian field, then the coefficients are complex-valued Gaussian stochastic processes with

By Remark 1, the random field given by (20) is homogeneous and isotropic in , and hence its covariance function (21) can be written in the form: where and are spherical coordinates of and respectively.

Using the Karhunen theorem we obtain the following spectral representation of the random field: where the random measures are given in (5).

Similarly to the condition (22) the isotropic measure satisfies the following condition if the field has a finite variance:

Subclasses of covariance functions of the isotropic fields on the sphere can be obtained from covariance functions of homogeneous isotropic random fields in Euclidean space, since a restriction of the homogeneous and isotropic random field to the sphere yields an isotropic spherical field (e.g., [25], p. 76).

Consider two locations and on the unit sphere with the angle between them. Then, the Euclidean distance between these two points is , the inner product is , which gives a direct correspondence between the covariance function in the Euclidean space and the covariance function on the sphere for every fixed . Thus, the restriction of the homogeneous and isotropic hyperbolic diffusion field (25) to is an isotropic spherical random field for every fixed . We will call it the spherical hyperbolic diffusion isotropic random field .

Its covariance function is of the form: By the addition theorem for Bessel functions, the random field has the following spectral representation: where and the random measure satisfies (5).

Thus, the angular spectrum of the isotropic spherical random field is given by the formula

Therefore, we obtained the following result.

Consider the random initial-value problem (1) and (2), in which (5).

Then, the restriction of the spatio-temporal hyperbolic-diffusion random field (25) to the sphere The field and its covariance functions are given by (27) and (26), respectively.

This result investigates the restriction of the spatio-temporal hyperbolic-diffusion random field to the sphere It shows how the angular power spectrum of the restriction depends on the Fourier transform Hence, one can explicitly compute coefficients and study the contributions of different spherical harmonics to the spatial-temporal field

Notice that . The angular power spectrum of , , will be denoted by , For spherical random fields with finite variances, it holds

If (30) holds true, then

By Theorem 3 Now, combining (31) and Lemma 1, one gets as , , and . This completes the proof. □

It follows from Lemma 2 and the estimate (24) is finite if the initial condition

4. Smoothness of Solutions

Numerous problems in mathematical physics and geosciences require studying the regularity properties of solutions of differential equations. Smoothness, boundedness of derivatives or Hölder continuity conditions are often used to describe and investigate local changes and growth rates of solutions. Knowing regularity properties is also essential for an adequate approximation of SPDE solutions. In those cases where solutions are given by infinite series, it is a rather difficult mathematical problem, as the tail terms of such series can accumulate.

In this section, we investigate the Hölder-type continuity of the solution given by (27) on the sphere. Estimations of the closeness of values at spherical points and are obtained. We demonstrate how they depend on the decay of the angular power spectrum and provide some specifications in terms of the spectral measure

First, we obtain the continuity of the solution with respect to the geodesic distance on the sphere. To prove this we use the approach from Corollary 5 in [8].

Let (1)– (2) and the random initial condition

Then, for

where γ is the angle between directions

If the measure

(a) It follows from (23), (24), (27) and (32) that

Applying the property of Legendre polynomials, [23], p. 16, one obtains the statement (a) of the theorem.

(b) It follows from the proof of (32) that in the case of it holds The remaining steps are similar to the proof in (a). □

When the geodesic distance vanishes (i.e., ), it is easy to see that and, therefore, as well.

The next two results provide conditions on the field’s spectrum that guarantee the Hölder-type regularity of

If the measure

even for the case of (33).

Indeed, by (26) we get

For it holds , when , and therefore

Hence, and (34) follows from Lemma 1. □

The next result gives sufficient conditions to guarantee (33).

Suppose that (33) holds true.

By the Poisson integral representation of the Bessel function it follows that as . □

5. Short and Long Memory

Investigating statistical dependence between measurements at two points with increasing time or spatial distance between them is an important issue for practical temporal or spatial predictions. The spatial domain of the considered random fields is restricted to the sphere with the geodesic distance Note that this distance is bounded to the interval , but time t can unboundedly increase and takes values in Hence, this section investigates only temporal statistical dependencies, namely, slow or fast decays of covariance functions in time. The corresponding cases represent long or short memory scenarios.

In this section we use the representation (26) of covariance functions to investigate the structure of dependences of over time. We demonstrate that conditional on the spectral isotropic measure of the initial random condition , the random field can exhibit short- or long-range dependence.

The random field will be called short-range dependent if for all and . If the integral in (35) is divergent, the field is called long-range dependent.

Results that link behaviours of covariance functions at infinity and spectral measures at the origin are called Abelian–Tauberian theorems. A very detailed overview of such results for random fields can be found in [28].

First we investigate the case of in (26) (i.e., the behaviour of ).

For

It follows from (16), (17) and (26) that

Using the upper bound from (19) we get Hence, to study the integrability of the covariance function one has to investigate the integral As for , , it is equivalent to studying the integral or, by (16) and for , to investigating the finiteness of the integral Noting that on , , we obtain that (37) is finite if and only if the following integral converges: The last integral is finite only if , which completes the proof. □

Now we extend Theorem 7 to the case of arbitrary and from .

The random field

Note that by (26) the integrators in and differ only by a multiplier .

Thus,

It follows from the estimates (19), (36) and the inequality that

Now, note that for the interval can be split into a finite number of subintervals where and denotes the integer part of a. The ratio has the same sign on each of these subintervals. Therefore, similar to the proof of Theorem 7 we obtain the sufficient and necessary condition for the integrability of Note that by this condition is equivalent to the one in Theorem 7. This completes the proof. □

6. Approximations to Solutions

The results in the previous sections were based on the series representation of the random field In applications and numerical studies, only a finite number of series terms are available. Hence, one has to investigate the behaviours of finite cumulative sums. This section provides an analysis of truncated series expansions of the solution field and shows the role of the decay rate of the angular power spectrum. These results can be used to determine the number of terms for a given accuracy of approximate solutions.

This section introduces and studies approximate solutions of the initial value problem (1) and (2). A mean-square convergence rate to the diffusion field in terms of the angular power spectrum is obtained. Then, several specifications in terms of the measure are discussed.

We define the approximation of the truncation degree to the solution given by (27) as The next result provides the convergence rate of to when

Let (1) and (2) and

Note that by properties of we get for all , , and .

Then, by (23) and (27) it follows that

Using the addition formula for spherical harmonics one gets

Finally, by (32)  □

For the SPDE model studied in [8] it was shown that its solution has an exponential decay in t and the corresponding approximation error can be bounded as see (36) in [8].

The following result shows that the considered model is more complex. In the general case of an arbitrary measure it is impossible to get a bound similar to (39) even for a sufficiently large L.

For any fixed (39).

Indeed, let us consider some .

Then, if .

Let the measure be concentrated on the interval . By (16), if then

Hence, by (38) and Theorem 3 for any , there exists , and the measure such that for the corresponding and it holds  □

However, it is possible to obtain a rate of convergence that is exponential in t if the measure has a bounded support.

Let (1) and (2) and its approximation

As is an increasing function on (0, ∞) it follows from (16) that for Notice that for and it holds that .

Then, using the definition of in (17) we get for

Hence, if it follows from Theorem 3 that Applying this bound to (38) we obtain the statement of the theorem. □

The next result follows from (38) and the upper bound (19) for .

If

The rates of convergence in Theorems 9, 11 and Corollary 1 are sharp. Indeed, for

The angular power spectrum , of the initial random field is determined by the measure . The following results provide some insight into the behaviour of in terms of the spectral measure .

Let the angular power spectrum of

Then, it holds that

If

If the measure

(a) It follows from the representation that By von Lommel’s formula (see (2.60) in [29]), where and , we obtain

Now, (40) follows by substituting the last expression in (43).

(b) Using the inequality from [30] we obtain that for which after the substitution in (43) gives (41).

(c) By the Poisson integral formula and the identity one obtains

If , is the support of the measure , then it follows from (43)–(45) that which completes the proof. □

7. Numerical Studies

This section presents numerical studies of the solution , its angular spectrum, and covariance functions over time. We also provide some numerical analysis of approximation errors.

All numerical computations and simulations in this paper were performed using the software R version 3.6.1 and Python version 3.7.5. The results were derived using the HEALPix representation of spherical data (see [31] and http://healpix.sourceforge.net). In particular, the R package [32,33] was used for computations and visualisations of the obtained results. The Python package was used for fast spherical harmonics generation of spherical maps from Laplace series coefficients. The R and Python code used for numerical examples in Section 7 are freely available in the folder “Research materials” from the website https://sites.google.com/site/olenkoandriy/.

It is important to clarify that the numerical analysis in this paper is rather different from the one in [8] and requires more advanced approximation approaches. Namely, the stochastic model in [8] yielded the representation of the Laplace series coefficients for some functions and which can be explicitly computed in terms of elementary functions. However, for the model (1) and (2) there is no such simple functional relation that links and . As a result, there are no explicit elementary functional relations between , and , respectively. To compute spectral and covariance functions of at time one has to use formulae (26), (28) and (29). These integral representations are given in terms of the spectral measure and stochastic measures of the initial random condition field .

By (1.2.5) in [26], it follows from that which can be used to compute (26), (29) and simulate for computations in (28). However, obtaining a reliable approximation of the integral in (46) and stochastic measures requires the estimation of the empirical covariance function on a dense grid. Moreover, for data observed on bounded subsets of covariance functions can be estimated only for distances that do not exceed their diameters. Thus, it is important to verify that empirical covariance functions are sufficiently quickly decaying to be assumed negligible for distances greater than these diameters. We postpone the solution of these technical problems and analysis of real data to future publications.

In the following examples, we study the properties of solutions and their approximations using simulated data. The examples were constructed to demonstrate that the model is sufficiently powerful to imitate behaviours of the empirical CMB covariance function and oscillating angular spectrum (see [8,34]). The actual CMB covariance function and angular spectrum are shown in Figure 1a,b. Note that the estimated angular CMB spectrum shown in Figure 1b was obtained by a piecewise fitting of several physical models and interpolation techniques for different intervals of the spectrum [6,13]. Some actual spectrum estimates deviate substantially from the fitted curve in Figure 1b [34]. Therefore, in predicting CMB and spectrum changes over time, small details can be ignored and one needs to focus on a general pattern. Thus, the following examples with the analysis of simulated data which spectrum is analogous to the real one can offer important insights on the future evolution of CMB and its spectral properties.

Figure 1

CMB: Cosmic Microwave Background radiation.

In the examples, the case of a discrete measure is considered (i.e., the support of is a finite set ). We employ real-valued stochastic measures that are concentrated on this set and satisfy the condition We assume that the random field is centred Gaussian. Hence, we can choose that are independent for different l, m and i.

In these settings, formulae (26), (28) and (29) take the following discrete forms: which are convenient for simulations.

This approach can also be used to approximate absolutely continuous spectral measures by considering a sufficiently large I, small and ,

This example illustrates changes over time of the covariance function

To understand the effect of the parameters c and D on the covariance function, we also produced

In this example we use a discrete spectrum concentrated on the two intervals

To produce the plots and computations in this paper we used the first 100 coefficients

For the following numerical studies we used simulated data from two windows shown in for the empirical distributions

The q-statistic, see [

8. Entropy and Hyperbolic Diffusion

This section discusses the evolution of Shannon entropy for hyperbolic diffusion. Theoretical analysis and several numerical examples are presented. To simplify the exposition and plots, only the case of and various non-random initial conditions are studied.

For diffusive transport that arises from the random motion of particles, the mass distribution may indeed be regarded as a probability distribution, after which the Shannon entropy may be calculated. For a simple thermodynamic system governed purely by linear or nonlinear heat conduction, there is a close analogy between thermodynamic entropy and Shannon entropy (e.g., [24,36]). When the transport mechanism is modified to hyperbolic diffusion, the behaviour of entropy requires more scrutiny. In order to illustrate this, consider one-dimensional solutions on , subject to Neumann boundary conditions This may represent transport in the x-direction through a linear conduit of cross-sectional area A, with the variation of density in each cross section being effectively zero. It will be seen that the total mass M is constant. Therefore, the scaled density has constant unit integral on , from which physically relevant non-negative solutions may be regarded as distributions. By choosing length scale and time scale , it may be assumed that the coefficients in the hyperbolic diffusion equation are normalised to .

Let , and . Then, subject to boundary conditions and initial conditions

Defining Shannon entropy density to be , the hyperbolic diffusion equation for implies The case of unbounded speed of propagation is obtained by taking the limit , which results in a positive entropy production rate . This is familiar from the theory of heat conduction, for which the entropy production rate is , where T is absolute temperature and is the Lewis number, which is the order-1 ratio of thermal diffusivity to mass diffusivity.

For uni-directional waves of velocity , the entropy production rate is zero. For bi-directional waves, the total Shannon entropy is constant when opposite-travelling waves are not superposing, increasing when opposite-travelling superposing waves are separating, and decreasing when they are superposing and approaching. However, non-constant travelling wave solutions of the hyperbolic diffusion equation must have speed less than c and they must have an amplitude that decreases with time. For the remainder of this section, the asterisk superscripts will be conveniently omitted.

Some solutions of the hyperbolic diffusion equation may be of dissipative diffusive type, while others may be dissipative bi-directional waves. In order to illustrate this, by the completeness of the Fourier transform, the general even solution by the separation of variables is where , , and .

The first summation covers modes that are purely dissipative in character, just as for the linear heat diffusion equation. However, in this case, the dissipative modes exist only when . The second summation covers standing wave modes with decaying amplitude. These may be regarded as a superposition of a decaying left-travelling wave and a decaying right-travelling wave. Note that the dissipative mode with logarithmic decay rate decays more slowly than all other modes.

The above solution is mass-conserving with mean value and constant mass integral by normalisation. For a single decaying standing wave mode of a hyperbolic diffusion equation distribution, for some value of t, Then, the total Shannon entropy is At times , the distribution is uniform, which is the state of maximum entropy . Overall, the total entropy oscillates as it approaches the limiting equilibrium state. However the negative excursions of entropy may be quite small since the amplitude of oscillation decreases exponentially.

Figure 10 plots the total entropy for a wave with single harmonic, calculated by trapezoidal integration with 400 intervals, versus time.

Figure 10

Total entropy for standing wave with single harmonic. , wave number

It would be helpful to have a point-source solution for the hyperbolic diffusion equation. As far as we are aware, there is no known simple expression for the point source evolution but it has the standard uniform Fourier spectrum that evolves according to (50). It is plotted in Figure 11 after truncating the Fourier series at 100 terms. As in the d’Alembert wave equation, two separating travelling delta waves emerge but now the amplitudes of the truncated spikes are decreasing and there is an additional central symmetric hump due to the purely diffusive terms. The leading edges of the spikes are travelling at maximum speed c. In two and three dimensions there would be similar solutions with a single travelling cylindrical or spherical shock wave surrounding a central hump.

Figure 11

Evolving spike solution for .

It is instructive also to view the motion of an initial rectangular disturbance of finite amplitude. This is approximated in Figure 12 and Figure 13 by a Fourier series of 200 terms. The truncated Fourier series is an exact solution, but due to the truncation and the boundary conditions, the solution is negative at some values of the domain, so that Shannon entropy cannot be calculated. However, the solution is indicative of the behaviour of a non-negative solution with initial rectangle. As in the bidirectional wave equation, the symmetric solution consists of two superposed rectangles that increase entropy as they begin to separate by travelling in opposite directions. After they have separated, their amplitude decreases, which leads to further entropy increase. The height of the leading edge decreases more rapidly than the trailing edge, so each rectangle evolves to a trapezoid. The leading edge—which is the boundary of the disturbance—continues to move at maximum speed c. Between the trapezoids, there is a central hump that eventually dominates, and resembles a diffusive Gaussian, increasing entropy further. With this kind of peaked initial condition, there is no indication of any significant period of entropy decrease.

Figure 12

Evolving symmetric rectangle: emergent bidirectional wave.

Figure 13

Evolving rectangle: dominant diffusive hump at large t, with leading edge of remnant rectangle demarcating the extent of the disturbance.

9. Future Research Problems

This paper investigated evolutions of random fields determined by hyperbolic diffusion equations with random initial conditions. Spherical random fields were modelled as restrictions of 3D solution fields to the sphere. Compared to the previous publications, it resulted in a more realistic physical model. However, the solution field for the new model cannot be represented by using Laplace series coefficients of the initial condition directly. The more complicated representation involves spectral measures of initial random conditions.

Detailed studies of the solutions and their approximations were presented. In particular, regularity properties and temporal dependencies of solutions were investigated. Approximations to the SPDE solutions were proposed, and the upper bound analysis of approximation errors was provided. It was demonstrated that the magnitude of approximation errors is determined by the angular power spectrum and decreases at the rate of the cumulative tail sums

The numerical studies investigated the dependence of solution fields on parameters of the SPDE model and provided some insight into the evolution of Shannon entropy for hyperbolic diffusion.

Some important problems and extensions for future research are:

Investigating the sharpness of the obtained upper bounds on approximation errors (see [8]);

Developing statistical estimators of the equation parameters and studying their asymptotic properties;

Extending the methodology to tangent spherical vector fields (see [37]);

Developing numerical methods for the obtained representations to deal with spectra of initial conditions;

Extending the analysis and numerical studies in Section 8 to other scenarios;

In line with the theme of this Special Issue, in future we intend to study the effect of nonlinear diffusivity in the equation

For example, if q is the electron density in a plasma, is typically decreasing [38].

Table 1

Statistics for windows 1 and 2.

Time	0	0.05	10
Mean for window 1	1.353×10^−5	−5.62×10^−6	3.501×10^−7
Mean for window 2	7.083×10^−6	−1.132×10^−5	−5.166×10^−8
IQR for window 1	2.877×10^−4	1.307×10^−4	6.78×10^−6
IQR for window 2	3.252×10^−4	1.452×10^−4	7.11×10^−6
Entropy for window 1	2.193	2.116	2.369
Entropy for window 2	2.302	2.221	2.387
q-statistics	1.986×10^−4	7.272×10^−4	1.5×10^−3