# Spectral Energy Density in an Axisymmetric Galaxy as Predicted by an Analytical Model for the Maxwell Field

Authors: Mayeul Arminjon
Journal: Advances in Astronomy
 (2021)

---

## Abstract

An analytical model for the Maxwell radiation field in an axisymmetric galaxy, proposed previously, is first checked for its predictions of the spatial variation of the spectral energy distributions (SEDs) in our Galaxy. First, the model is summarized. It is now shown how to compute the SED with this model. Then the model is adjusted by asking that the SED predicted at our local position in the Galaxy coincides with the available observations. Finally, the first predictions of the model for the spatial variation of the SED in the Galaxy are compared with those of a radiation transfer model. We find that the two predictions do not differ too much. This indicates that, in a future work, it should be possible with the present model to check if the “interaction energy” predicted by an alternative, scalar theory of gravitation, contributes to the dark matter.


---

## Body

## 1. Introduction

The interstellar radiation field (ISRF) is a very important physical characteristic of a galaxy, and there is indeed a rich astrophysical literature that studies the modeling of the ISRF (see, e.g., [1–6] and references therein). This subject is of importance not only in itself but also because the radiation field produced by the stars and other bright objects has a strong interaction with the cosmic rays, which are another important component of the interstellar medium. The necessity of a modeling is obvious, since we can directly measure the radiation field only in the Solar System, in fact in the neighborhood of the Earth [7–9]. The ISRF is significantly affected by its interaction with interstellar dust (e.g., [10, 11]). In order to account for this interaction, equations for radiative transfer are solved either by using Monte Carlo methods, which simulate in a probabilistic framework the propagation of photons and their absorption and scattering by dust, in a discretized spatial domain (e.g., [1, 4]), or by using “ray-tracing” algorithms that calculate the variation, due to absorption and scattering by dust, in the radiation’s specific intensity along a finite set of directions (e.g., [3, 4, 6]). Thus, until now, the electromagnetic (EM) field, made of coupled electric and magnetic fields, has not been considered in the studies of the ISRF in galaxies, even less an EM field that would be a solution of the Maxwell equations—as should actually be the case for the real ISRF, however.Obtaining for the ISRF such an exact solution of the Maxwell equations is somewhat demanding from both the mathematical and the numerical points of view and has been the main purpose of a recent work [12]. A decisive motivation for that recent work was to prepare for checking an exotic prediction of an alternative, scalar theory of gravitation, namely, that the total energy(-momentum-stress) tensor, which provides the source of the gravitational field, is not in general the sum of the energy tensor of the EM field and the energy tensor of the material medium producing it. Instead, each time that an EM field is present together with a variable gravitational field (which is the general case), there should appear an additional energy tensor: the “interaction” energy tensor, Tinter [13]. This speculative energy tensor Tinter, being thus present virtually everywhere according to that theory and having clearly an exotic character, could contribute to the dark matter. To calculate Tinter, and to see if its distribution could be similar to dark matter halos, it is necessary to be able to calculate the EM field and its first-order derivatives [13], whence the work [12]. In the latter work, a model has been formulated. It is based on an analytical representation of the general solution of the source-free Maxwell equations in the axisymmetric case [14]. Moreover, the model has been numerically implemented and has passed a numerical validation test based on an exact solution [12].The aim of the present paper is to make a first check of the predictions of that “Maxwell model of the interstellar radiation field” for the spectral energy distribution (SED) and its spatial variation in a galaxy. Such SEDs are a very important output of the mainstream models for the analysis of the ISRF and, as was recalled at the beginning of this paper, those models depend essentially on the radiation transfer (mainly through dust). Note that, in contrast, the present model [12] does not directly take into account any radiation transfer. Section 2 summarizes that model. Section 3 details the calculation of the SED in that model. Section 4 presents the results obtained with that model for a simulation of the Milky Way. To obtain these results, the model is adjusted on the locally observed SED. In this way, it does indirectly take into account the radiation transfer because the observed SED is indeed affected by the latter.


## 2. The Maxwell Model of the Interstellar Radiation Field

In this section, we give a new presentation of the model that was built in Ref. [12] and that uses results of Ref. [14] and references therein. This presentation uses an easy extension of the result [14]. Such an extension was used implicitly in Ref. [12] and is explicitly formulated in Appendix A.


### 2.1. General Assumptions for the Maxwell Field

We make two assumptions regarding the form of the Maxwell radiation fieldE,B (where E is the electric field and B is the magnetic field):(i)
We assume that it is anaxially symmetric field. This is a relevant first approximation for many galaxies, even though it neglects some aspects like the arm structure of a spiral galaxy. Axial symmetry is often assumed in the literature on modeling the ISRF in galaxies and/or its interaction with dust (e.g., [2, 3, 6, 10]).(ii)
We search for asource-free solution of the Maxwell equations. This is because we seek to describe the EM field at the galactic scale, not the local EM field in the neighborhood of the stars that are the primary sources of that field.


### 2.2. The Explicit Representation for Axisymmetric Free Maxwell Fields

This is based on the following recentTheorem [14]: any time-harmonic axisymmetric solution E,B of the source-free Maxwell equations has a unique decomposition(1)E=E1+E2′,B=B1+B2′,where E1,B1 has the form (in cylindrical coordinates ρ,ϕ,z corresponding to the symmetry axis Oz):(2)Bϕ=−∂Ψ∂ρ,Eϕ=0,(3)Eρ=ic2ω∂2Ψ∂ρ∂z,Bρ=0,(4)Ez=ic2ω∂2Ψ∂z2+iωΨ,Bz=0,and where E2′,B2′ is deduced from a solution E2,B2 having just that same form (2)–(4), by the EM duality—i.e.,(5)E′=cB,B′=−Ec.Here,c is the velocity of light. Note that the dual solution (5) has the complementary zero and nonzero components as compared with the direct solution (2)–(4), i.e., Bϕ′=Eρ′=Ez′=0, while Eϕ′,Bρ′,Bz′≠0 in general. Both solutions E1,B1 and E2′,B2′ are time-harmonic with the same angular frequency ω as has E,B. In equations (2)–(4), the potential Ψt,ρ,z=e−iωtΨ^ρ,z is a time-harmonic axisymmetric solution of the standard scalar wave equation, which also has the same frequency ω as has E,B; in general, the two solutions E1,B1 and E2,B2 correspond to two different scalar potentials Ψ and Ψ′ [14]. The solution (2)–(4) actually derives from a vector potential having just its axial component Az=Ψ nonzero, i.e., A=Ψez [14, 15]. The Theorem in [14] states the existence of the two potentials Ψ and Ψ′, not their uniqueness. In Appendix A, we extend this result to EM fields having a finite spectrum, which is easy and allows us to give a more rigorous presentation of the model in Section 2.3.If a time-harmonic, axisymmetric solution of the scalar wave equation istotally propagating, it has the following explicit form [15, 16]:(6)ψωSt,ρ,z=e−iωt∫−K+KJ0ρK2−k2eikzSkdk,where ω is the angular frequency, K≔ω/c, J0 is the first-kind Bessel function of order 0, and S is some (generally complex) function of k∈−K,+K. “Totally propagating” means, in short, that Ψ has no “evanescent” component [15]; a detailed explanation is given in Ref. [12]. In turn, we shall define simply that a time-harmonic, axisymmetric EM field F=E,B is totally propagating if the two potentials Ψ and Ψ′, whose existence is guaranteed by the Theorem recalled above, can themselves be chosen totally propagating, i.e., having the form (6).


### 2.3. Statement of the Model

We model a disk galaxy as a disk-like set of point-like “stars” (the corresponding points are obtained by pseudorandom generation of their cylindrical coordinates. Here, the same set of16×16×36 triplets ρ,z,ϕ was used as in Ref. [12], thus 9216 points. The distribution of ρ and z is approximately valid for the star distribution in the Milky Way, and the distribution of ϕ is uniform, thus ensuring approximate axial symmetry). Each of them is assumed to produce an EM field with a finite spectrum: that field derives, in the way defined in the previous section, from potentials that obey the scalar wave equation, more precisely ones with spherical symmetry:(7)ψωjt,x=eiKjr−ωjtKjr,Kj≔ωjc,r≔x.Thus, the sum of the potentials emitted by the different stars at the different frequencies is(8)Ψxiωjwjt,x=∑i=1imax∑j=1NωwjeiKjri−ωjtKjri,Hereri≔x−xi, where x is the spatial position at which the function is evaluated and xi is the spatial position of the ith star; the numbers wj>0, with ∑jwj=1, are the weights affected to the different frequencies (in Reference [12], these weights were noted Sj′, which now would clash with the functions Sj′ mentioned at the end of the present study). Since the disk-like distribution of the “stars” is built (approximately) axisymmetric [12], function (8) is also (approximately) axisymmetric. Moreover, this is a solution of the scalar wave equation, and it is clearly a totally propagating one, as is each of the individual spherical scalar waves in (8). It should thus be possible to put function (8)—at least approximately—in the form of a sum, over the frequencies, of functions having the form (6). This is done by the least-squares method applied on a finite set of events t,x (with x=xρ,ϕ=0,z restricted to the plane ϕ=0, due to the axial symmetry), making a regular spatiotemporal grid tl,ρm,zp. We thus want to determine complex functions Sj=Sjk−Kj≤k≤+Kj, such that(9)Ψxiωjwjt,ρ,z≅∑jψωjSjt,ρ,z≔ΨωjSjt,ρ,z,where the sign ≅ indicates that the equality is in the sense of the least squares (on the relevant spatiotemporal grid, just mentioned). To the unique function on the r.h.s. of (9), Ψ=∑jψωjSj, we may then associate the following:(i)
Either the EM field that is obtained by summing the time-harmonic contributions given by equations (2)–(4), with Az≔ψωjSj—that is,(10)F1≔Z1Ψ.(cf. equation (A.3) in the Appendix).(ii)
Or, the dual (5) of that field—that is,(11)F2≔Z2Ψ.(iii)
Or still, the sum of these two fields, that is(12)F=F1+F2=Z1Ψ+Z2Ψ≔ZΨ,Ψ.(cf. equation (A.4) in the Appendix).


### 2.4. Discretization

To calculate the integrals in equation (6) and in equations (2)–(4) applied with functions ψωjSj of form (6), we discretize the integration intervals −Kj,+Kj and the functions Sj, and we use the “Simpson 3/8 composite rule” [12]. We set K0≔ω0/c, with ω0 some (arbitrary) reference frequency, and(13)kn≔−K0+nδ0,n=0,…,N,δ0≔2K0N,(14)Snj≔Sjωjω0kn,n=0,…,N,j=1,…,Nω,(15)an≔38δ0n=0orn=N,(16)an≔2×38δ0,modn,3=0andn≠0andn≠N,(17)an≔3×38δ0,otherwise.With this discretization, the unknown functionsSjj=1,…,Nω in equation (9) are characterized by the complex numbers Snjn=0,…,N;j=1,…,Nω. Thus, the equality in the sense of the least squares, schematized by equation (9), actually becomes [12](18)Ψxiωjwjt,ρ,z≅∑n=0N∑j=1Nωfnjt,ρ,zSnj,with(19)fnjt,ρ,z=ωjω0anJ0ρωjω0K02−kn2expiωjω0knz−ωjt.TheSnj’s are the solved-for parameters in the least-squares problem (18).This gives us for the EM field (10) [12], as also for the field (12):(20)Bϕt,ρ,z=∑n=0N∑j=1NωRnJ1ρωjω0RnℛeFnjt,z+O1N4,(21)Eρt,ρ,z=∑n=0N∑j=1Nωc2ω0knRnJ1ρωjω0RnℛeFnjt,z+O1N4,(22)Ezt,ρ,z=∑n=0N∑j=1Nωc2ω0kn2−ω0J0ρωjω0RnℐmFnjt,z+O1N4,where Rn=K02−kn2 and(23)Fnjt,z=ωjω02anexpiωjω0knz−ωjtSnj.Even though these formulae follow from a discretization and from using an integration rule that is only approximate for the starting integral, they provide (without theO1/N4 remainder) an exact solution of the source-free Maxwell equations (see §3.3 in Reference [12]).The arguments of the Bessel functionsJ0 and J1 in equations (19) and (20)–(22), as well as the spatial part of the argument of the complex exponential in equations (19) and (23), are of the order of ρ/λ or z/λ. For a galaxy, these are numbers of the order of 1025; therefore, the numerical implementation of the model needs to use a precision better than quadruple precision—which leads to long computation times [12].


## 2.1. General Assumptions for the Maxwell Field

We make two assumptions regarding the form of the Maxwell radiation fieldE,B (where E is the electric field and B is the magnetic field):(i)
We assume that it is anaxially symmetric field. This is a relevant first approximation for many galaxies, even though it neglects some aspects like the arm structure of a spiral galaxy. Axial symmetry is often assumed in the literature on modeling the ISRF in galaxies and/or its interaction with dust (e.g., [2, 3, 6, 10]).(ii)
We search for asource-free solution of the Maxwell equations. This is because we seek to describe the EM field at the galactic scale, not the local EM field in the neighborhood of the stars that are the primary sources of that field.


## 2.2. The Explicit Representation for Axisymmetric Free Maxwell Fields

This is based on the following recentTheorem [14]: any time-harmonic axisymmetric solution E,B of the source-free Maxwell equations has a unique decomposition(1)E=E1+E2′,B=B1+B2′,where E1,B1 has the form (in cylindrical coordinates ρ,ϕ,z corresponding to the symmetry axis Oz):(2)Bϕ=−∂Ψ∂ρ,Eϕ=0,(3)Eρ=ic2ω∂2Ψ∂ρ∂z,Bρ=0,(4)Ez=ic2ω∂2Ψ∂z2+iωΨ,Bz=0,and where E2′,B2′ is deduced from a solution E2,B2 having just that same form (2)–(4), by the EM duality—i.e.,(5)E′=cB,B′=−Ec.Here,c is the velocity of light. Note that the dual solution (5) has the complementary zero and nonzero components as compared with the direct solution (2)–(4), i.e., Bϕ′=Eρ′=Ez′=0, while Eϕ′,Bρ′,Bz′≠0 in general. Both solutions E1,B1 and E2′,B2′ are time-harmonic with the same angular frequency ω as has E,B. In equations (2)–(4), the potential Ψt,ρ,z=e−iωtΨ^ρ,z is a time-harmonic axisymmetric solution of the standard scalar wave equation, which also has the same frequency ω as has E,B; in general, the two solutions E1,B1 and E2,B2 correspond to two different scalar potentials Ψ and Ψ′ [14]. The solution (2)–(4) actually derives from a vector potential having just its axial component Az=Ψ nonzero, i.e., A=Ψez [14, 15]. The Theorem in [14] states the existence of the two potentials Ψ and Ψ′, not their uniqueness. In Appendix A, we extend this result to EM fields having a finite spectrum, which is easy and allows us to give a more rigorous presentation of the model in Section 2.3.If a time-harmonic, axisymmetric solution of the scalar wave equation istotally propagating, it has the following explicit form [15, 16]:(6)ψωSt,ρ,z=e−iωt∫−K+KJ0ρK2−k2eikzSkdk,where ω is the angular frequency, K≔ω/c, J0 is the first-kind Bessel function of order 0, and S is some (generally complex) function of k∈−K,+K. “Totally propagating” means, in short, that Ψ has no “evanescent” component [15]; a detailed explanation is given in Ref. [12]. In turn, we shall define simply that a time-harmonic, axisymmetric EM field F=E,B is totally propagating if the two potentials Ψ and Ψ′, whose existence is guaranteed by the Theorem recalled above, can themselves be chosen totally propagating, i.e., having the form (6).


## 2.3. Statement of the Model

We model a disk galaxy as a disk-like set of point-like “stars” (the corresponding points are obtained by pseudorandom generation of their cylindrical coordinates. Here, the same set of16×16×36 triplets ρ,z,ϕ was used as in Ref. [12], thus 9216 points. The distribution of ρ and z is approximately valid for the star distribution in the Milky Way, and the distribution of ϕ is uniform, thus ensuring approximate axial symmetry). Each of them is assumed to produce an EM field with a finite spectrum: that field derives, in the way defined in the previous section, from potentials that obey the scalar wave equation, more precisely ones with spherical symmetry:(7)ψωjt,x=eiKjr−ωjtKjr,Kj≔ωjc,r≔x.Thus, the sum of the potentials emitted by the different stars at the different frequencies is(8)Ψxiωjwjt,x=∑i=1imax∑j=1NωwjeiKjri−ωjtKjri,Hereri≔x−xi, where x is the spatial position at which the function is evaluated and xi is the spatial position of the ith star; the numbers wj>0, with ∑jwj=1, are the weights affected to the different frequencies (in Reference [12], these weights were noted Sj′, which now would clash with the functions Sj′ mentioned at the end of the present study). Since the disk-like distribution of the “stars” is built (approximately) axisymmetric [12], function (8) is also (approximately) axisymmetric. Moreover, this is a solution of the scalar wave equation, and it is clearly a totally propagating one, as is each of the individual spherical scalar waves in (8). It should thus be possible to put function (8)—at least approximately—in the form of a sum, over the frequencies, of functions having the form (6). This is done by the least-squares method applied on a finite set of events t,x (with x=xρ,ϕ=0,z restricted to the plane ϕ=0, due to the axial symmetry), making a regular spatiotemporal grid tl,ρm,zp. We thus want to determine complex functions Sj=Sjk−Kj≤k≤+Kj, such that(9)Ψxiωjwjt,ρ,z≅∑jψωjSjt,ρ,z≔ΨωjSjt,ρ,z,where the sign ≅ indicates that the equality is in the sense of the least squares (on the relevant spatiotemporal grid, just mentioned). To the unique function on the r.h.s. of (9), Ψ=∑jψωjSj, we may then associate the following:(i)
Either the EM field that is obtained by summing the time-harmonic contributions given by equations (2)–(4), with Az≔ψωjSj—that is,(10)F1≔Z1Ψ.(cf. equation (A.3) in the Appendix).(ii)
Or, the dual (5) of that field—that is,(11)F2≔Z2Ψ.(iii)
Or still, the sum of these two fields, that is(12)F=F1+F2=Z1Ψ+Z2Ψ≔ZΨ,Ψ.(cf. equation (A.4) in the Appendix).


## 2.4. Discretization

To calculate the integrals in equation (6) and in equations (2)–(4) applied with functions ψωjSj of form (6), we discretize the integration intervals −Kj,+Kj and the functions Sj, and we use the “Simpson 3/8 composite rule” [12]. We set K0≔ω0/c, with ω0 some (arbitrary) reference frequency, and(13)kn≔−K0+nδ0,n=0,…,N,δ0≔2K0N,(14)Snj≔Sjωjω0kn,n=0,…,N,j=1,…,Nω,(15)an≔38δ0n=0orn=N,(16)an≔2×38δ0,modn,3=0andn≠0andn≠N,(17)an≔3×38δ0,otherwise.With this discretization, the unknown functionsSjj=1,…,Nω in equation (9) are characterized by the complex numbers Snjn=0,…,N;j=1,…,Nω. Thus, the equality in the sense of the least squares, schematized by equation (9), actually becomes [12](18)Ψxiωjwjt,ρ,z≅∑n=0N∑j=1Nωfnjt,ρ,zSnj,with(19)fnjt,ρ,z=ωjω0anJ0ρωjω0K02−kn2expiωjω0knz−ωjt.TheSnj’s are the solved-for parameters in the least-squares problem (18).This gives us for the EM field (10) [12], as also for the field (12):(20)Bϕt,ρ,z=∑n=0N∑j=1NωRnJ1ρωjω0RnℛeFnjt,z+O1N4,(21)Eρt,ρ,z=∑n=0N∑j=1Nωc2ω0knRnJ1ρωjω0RnℛeFnjt,z+O1N4,(22)Ezt,ρ,z=∑n=0N∑j=1Nωc2ω0kn2−ω0J0ρωjω0RnℐmFnjt,z+O1N4,where Rn=K02−kn2 and(23)Fnjt,z=ωjω02anexpiωjω0knz−ωjtSnj.Even though these formulae follow from a discretization and from using an integration rule that is only approximate for the starting integral, they provide (without theO1/N4 remainder) an exact solution of the source-free Maxwell equations (see §3.3 in Reference [12]).The arguments of the Bessel functionsJ0 and J1 in equations (19) and (20)–(22), as well as the spatial part of the argument of the complex exponential in equations (19) and (23), are of the order of ρ/λ or z/λ. For a galaxy, these are numbers of the order of 1025; therefore, the numerical implementation of the model needs to use a precision better than quadruple precision—which leads to long computation times [12].


## 3. Calculation of the Spectral Energy Density

The volume energy density of an EM field, in SI units, is(24)U=δWδV=ε02E2+c2B2,where ε0=1/4π×9×109. In the relevant case with a frequency spectrum, the fields corresponding to the different frequencies do add; hence, the quadratic expression (24) is not additive. The time average of such a quadratic expression is yet known to be additive in the case of a Fourier expansion [17]. Let us check if this remains true in the present case of an arbitrary finite frequency spectrum ωjj=1,…,Nω. Let Ft be some components of the EM field at some given spatial position x. We have, with some complex coefficients, Cj=Cjx:(25)Ft=ℛe∑j=1NωCje−iωjt,whence(26)4F2=∑j=1NωCje−iωjt+Cj∗eiωjt∑k=1NωCke−iωkt+Ck∗eiωkt.Unlessωj+ωk=0, the time average e−iωj+ωkt¯ is zero. Unless ωj−ωk=0, the time average e−iωj−ωkt¯ is zero. But ωj+ωk=0 does not occur because ωj>0, and ωj−ωk=0 occurs only when j=k. Hence, (26) gives indeed the additivity for the time average of the square of a component admitting an expansion (25):(27)F2¯=14∑j=1NωCjCj∗+Cj∗Cj=12∑j=1NωCj2.For the EM field (10), its components Bϕ,Eρ,andEz admit expansions (25) with coefficients Cj1, Cj2, and Cj3, respectively, which are directly read on equations (20)–(22):(28)Cj1ρ,z=∑n=0NRnJ1ρωjω0RnGnj,(29)Cj2ρ,z=∑n=0Nc2ω0knRnJ1ρωjω0RnGnj,(30)Cj3ρ,z=∑n=0Nc2ω0kn2−ω0J0ρωjω0Rn−iGnj,where Rn=K02−kn2 and(31)Gnjz=expiωjtFnjt,z=ωjω02anexpiωjω0knzSnj.From (24) and (27), we thus get for field (10) a (time-averaged) energy density with a discrete spectrum:(32)U¯ρ,z≔δWδV¯=∑j=1Nωuj,uj≔14∑q=13αqCjq2,where the coefficients Cjqq=1,2,3 are given by equations (28)–(30), α1=ε0c2, and α2=α3=ε0. As to field (11), it is the dual (5) of field (10). Therefore, the energy density (24) for field (11) is given by the same equation (32). Finally, the energy density for field (12) is hence just the double of this. Now, in the next section, the model will be adjusted so as to exactly describe the measured local EM spectrum, i.e., the values ujxloc will be imposed to be the ones measured at our local position xloc in the Galaxy. In this context, the precise choice of the model, i.e., equations (10) or (11) or (12), is therefore totally neutral.


## 4. Results and Comparison with a Radiation Transfer Model

We first explain the final adjustment of the model. The scalar waveΨxiωjwj given by equation (8) depends on the numbers wj>0, with ∑jwj=1. Each of them represents the relative weight of the corresponding frequency ωj—though in some average sense because wj is common to all “stars” at the points xi. There is no exact way to assign values to these relative weights wj, which, together with the set xi of the positions of the stars, determine the functions Sj on the r.h.s. of equation (9)—or more exactly, determine the discretized values Snj in equations (14) and (18). In turn, the values Snjn=0,…,N, for a given value of the frequency index j, determine the intensity of the ωj-component of the scalar wave obtained after the fitting (the r.h.s. of (18)) and thus determine the intensity of the ωj-component of the EM field (20)–(22). However, since we start only from relative weights wj, this determination can be only up to a multiplying factor, and since the very values of the wj’s are assigned somewhat arbitrarily, that determination of the Snj’s n=0,…,N is really up to a factor depending on j, ξj>0 (this would be even clearer if we would separate the least-squares adjustment (18) into the different frequencies ωj, which is actually easy and even would allow us to eliminate the common time dependence exp−iωjt on both sides, thus eliminating the time variable in the adjustment. This different procedure will be numerically tested in a future work). We determine the values ξj by imposing that the values uj given by equation (32), calculated at our local position xloc in the Galaxy, be equal to the ones measured. We take ρloc=8 kpc and zloc=0.02 kpc (e.g., [18]). We took the local values of the SED from Figure 1 of Reference [2], using the software CurveUnscan to digitalize the curves. Those local values were determined from the Apollo 17 mission [7] and from the COBE mission: the DIRBE experiment [8] and the FIRAS experiment [9]. Figure 1 shows the corresponding SED. For the calculations, in order to save computer time, we extracted 76 values of λ from the complete list of values of λ in the curve; the corresponding points are also shown in Figure 1.Figure 1
SED as measured by space missions. Between the data of Apollo (ending atlog10λ≃−0.8) and DIRBE (beginning at log10λ≃+0.1), the curve has been interpolated.Then, at any other spatial positionxρ,z (in the axially symmetric case considered), the model gives a prediction for the values ujx of the SED corresponding to the frequencies ωj (or the wavelengths λj) considered in our calculation. The spatiotemporal grid used was Nt=5,Nρ=10,andNz=21 (thus, 1050 events t,ρ,z), with t varying from 0 to T0=λ0/c with λ0=10μm, ρ varying from 0 to 10 kpc, and z from −1 to +1 kpc. Other grids have also been tested, giving very similar results. Moreover, unless otherwise mentioned, the discretization number N was N=48 (this is the number of subintervals in which any of the integration intervals −Kj,+Kj is subdivided, see equations (13) and (14)). Figures 2–5 show the comparison of the SEDs, obtained at four different places in our Galaxy, according to either the thus-adjusted present model or the radiation transfer model of Popescu et al. [6] (again, the relevant curves, Figure 9 in Reference [6], have been digitalized using CurveUnscan).Figure 2
SED atρ=8 kpc and z=0 kpc, according to this model or a radiation transfer model.Figure 3
SED atρ=8 kpc and z=1 kpc, according to this model or a radiation transfer model.Figure 4
SED atρ=1 kpc and z=0 kpc, according to this model or a radiation transfer model.Figure 5
SED atρ=1 kpc and z=1 kpc, according to this model or a radiation transfer model.Figure 6
Predicted maximum values of the SED:Nω=23, three different grids.Clearly, the general agreement between these two calculations is far from perfect. In particular, our calculation oscillates markedly as function ofλ: at all wavelengths, for the two positions that are closer to the galactic centre (ρ=1 kpc, Figures 4 and 5); and merely at short wavelengths, for the two positions that are farther from the galactic centre (ρ=8 kpc, Figures 2 and 3). Nevertheless, this agreement is actually surprisingly good if one remembers how different are the two models. The SED predicted by our model is higher than the one predicted by the radiation transfer model in Ref. [6] at the short wavelengths (λ≲0.5μm, or log10λ≲−0.3)—except for the position (ρ=8 kpc and z=0 kpc), that is, relatively speaking, very close to our local position (ρloc=8 kpc and zloc=0.02 kpc), where the adjustment is done. The higher levels found with our model at the short wavelengths might have a relation with the fact that the absorption by dust, not taken into account by our model, is strong for these radiations. However, this relation is not obvious, since our model is adjusted on the measured SED at our local position, and that local SED is of course affected by the dust absorption to which the radiations directed towards here have been subjected along their paths in the Galaxy. Actually, one notes also that, in the long wavelength domain as well as in the very short wavelength domain, the SED predicted by this model is, in average, higher than that predicted by the radiation transfer model, at the two positions which have a higher altitude (z=1 kpc): Figures 3 and 5. Comparing Figures 2 and 3, and comparing Figures 4 and 5, it appears that the present model with its current numerical settings predicts a slower decrease of the energy density as z (or rather z) increases than does the radiation transfer model. Of course the two models are very different, also in their modeling of the Galaxy: a discrete star distribution (for the present model) vs. an analytical, continuous variation of the stellar emissivity as function of ρ and z (for the radiation transfer model of Popescu et al. [6]). However, the ρ and z distribution of the stars in our model has been built to correspond (at least approximately) with that for the Milky Way [12], as has been built the ρ and z dependence of the stellar emissivity function for the radiation transfer model [6]. Of course also, the SEDs of a radiation transfer model (like that of Ref. [6]) are a prediction, too, not a measurement. Nevertheless, a further investigation is necessary to quantify the slower decrease of the SED with increasing z for the present model and to understand the reason for it.A surprising prediction of our model is the one for the maximum valuesujmax of the SED in our model of the Galaxy—the SED being calculated at each point on a regular ρ,z grid, as described above. Figure 6 shows these maximum values for three different spatial grids, with the same set of 23 wavelengths or frequencies. It is seen that the ujmax values are quite similar for the different grids, and these values are extremely high at the short wavelengths, i.e., at high energies of the radiation field. They are found on the ρ=0 axis, thus on the central axis of the galaxy (actually, the values ujρ=0,z are fairly independent of the value of z in the domain investigated, i.e., −1kpc≤z≤1kpc). However, if one refines the frequency mesh, going from Nω=23 to Nω=46, and then to Nω=76, one observes that the very high values of ujmax are obtained only at shorter and shorter wavelengths (Figures 7 and 8). A similar phenomenon is observed when one refines the discretization parameter N: for a given frequency mesh, increasing N reduces the domain of the high values of ujmax towards smaller values of λj, or, equivalently, decreasing N expands that domain towards larger values of λj (Figures 9–11). The likely common reason for the two effects is that, when one refines either the frequency mesh (Nω increased) or the discretization of the integration intervals (N increased), then the number of the solved-for parameters Snjn=0,…,N;j=1,…,Nω increases. In most of the settings investigated in the present work, this number, Npara=N+1×Nω, is actually larger than the number Ndata of “data,” that is, the number of points in the spatiotemporal grid (indeed, at each of those points, one calculates the function on the l.h.s. of equation (18). These are the data that have to be fitted by using the function on the r.h.s., which depends on the solved-for parameters Snj) and thus Ndata=Nt×Nρ×Nz. We have Ndata=1050 for the grid used in this work, except for the two grids used for the first two curves in Figure 6, which give, respectively, Ndata=660 and Ndata=550 (we have Nt=5 also for those two grids). Hence, we are often in a situation of “overfitting” (however, the least-squares solver, the Matlab backslash operator mldivide, puts to zero many of the Snj parameters for reason of insufficient numerical accuracy, so that their effective number was 150 to 200, thus significantly smaller than Ndata). Higher values of the ratio R=Npara/Ndata correspond to a more overfitted situation, in which one expects the predictions to be less accurate. Now, on each among Figures 7–11, we saw that the domain of the low values of ujmax increases when Npara increases, the number Ndata being fixed at Ndata=1050 for these five figures—thus we saw that the domain of the low values of ujmax increases when R increases, i.e., when more overfitting occurs. Therefore, it is the low values of ujmax which are likely to be a numerical artefact, not the extremely high values. Note that this effect can be seen also in Figure 6: N and Nω, hence Npara, are the same for the three curves. So the third and finest grid Nρ=10×Nz=21 leads to a lower R than do the other two grids Nρ=12×Nz=11 and Nρ=10×Nz=11—and accordingly the domain of the high values of ujmax is larger for the third curve. On the other hand, the effect of overfitting seems less important at a distance from the axis (where the values of the SED are not high at all): e.g., if one changes Nω from 76 to 23, the SEDs considered in Figures 2–5 stay similar: they show a smaller number of oscillations due to the smaller number of points, but the amplitude of the oscillations is comparable. Tentatively, we may understand this by finding it natural that the very high numbers are more sensitive to overfitting. In summary, the status of the very high values of the maxima of ujρ,z: physical prediction of the model or numerical artefact, is not really clear at the present stage of this work—but the foregoing discussion, with different examples that all show an effect which can be attributed to the “amount of overfitting,” seems to indicate that those high values might be a true prediction of the model.Figure 7
Predicted maximum values of the SED:Nω=23 vs. Nω=46.Figure 8
Predicted maximum values of the SED:Nω=46 vs. Nω=76.Figure 9
Predicted maximum values of the SED forNω=23: N=48 vs. N=24.Figure 10
Predicted maximum values of the SED forNω=46: N=48 vs. N=24.Figure 11
Predicted maximum values of the SED forNω=46: N=48 vs. N=192.


## 5. Conclusion

In this work, the Maxwell model of the interstellar radiation field [12] has been first checked for its predictions of the spatial variation of the spectral energy distributions (SEDs) in our Galaxy. The model has been adjusted by asking that the SED predicted at our local position in the Galaxy, equation (32) with ρ=ρloc≔8 kpc and z=zloc≔0.02 kpc, coincides with the observations collected through different spatial missions [7–9]. Then the predictions of the model for other positions in the Galaxy (for which, of course, we have no real observation) have been compared with the predictions obtained by using a recent radiation transfer model [6]. The two predictions do not differ too much in magnitude, even though the predictions of the present model oscillate rather strongly as function of the wavelength, especially at short wavelengths and not very far (1 kpc) from the axis of the Galaxy. Also, the SED decreases more slowly with increasing altitude for the present model. The values of the SED reach their maximum on the axis itself of the Galaxy; that maximum is extremely high at short wavelengths. These surprisingly high values occur in a larger or smaller domain of wavelengths, depending on the settings of the calculation. However, this dependence is governed by the “amount of overfitting:” less overfitting increases the range of the high values. This makes it plausible that the high values might be a true prediction of the model. The next work will aim at improving the numerical procedure, in order (i) to try to establish if the high values are indeed a true prediction of the model; (ii) to check if a numerically improved model could lead to oscillations of a smaller amplitude; and (iii) to quantify the slower decrease of the SED with increasing altitude for the present model and to find the reason for it. Nevertheless, the present results (1) show that the prediction, by the Maxwell model of the ISRF, of the spatial variation of the SED in the Galaxy is quite comparable in magnitude to that provided by a radiation transfer model (except perhaps for high energies when the point is very close to the galaxy’s axis) and therefore (2) indicate that an order-of-magnitude estimate of the “interaction energy” predicted [13] by an alternative, scalar theory of gravitation, should be indeed possible with the present model in a future work.



---

*Source: 5524600-2021-08-02.xml*