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.

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., [16] 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 [79]. 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, [13]. This speculative energy tensor , being thus present virtually everywhere according to that theory and having clearly an exotic character, could contribute to the dark matter. To calculate , 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 field (where is the electric field and is the magnetic field):(i)We assume that it is an axially 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 a source-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 recent Theorem [14]: any time-harmonic axisymmetric solution of the source-free Maxwell equations has a unique decompositionwhere has the form (in cylindrical coordinates corresponding to the symmetry axis ):and where is deduced from a solution having just that same form (2)–(4), by the EM duality—i.e.,

Here, 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., , while in general. Both solutions and are time-harmonic with the same angular frequency as has . In equations (2)–(4), the potential is a time-harmonic axisymmetric solution of the standard scalar wave equation, which also has the same frequency as has ; in general, the two solutions and correspond to two different scalar potentials and [14]. The solution (2)–(4) actually derives from a vector potential having just its axial component nonzero, i.e., [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 is totally propagating, it has the following explicit form [15, 16]:where is the angular frequency, , is the first-kind Bessel function of order 0, and is some (generally complex) function of . “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 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 of triplets was used as in Ref. [12], thus 9216 points. The distribution of and 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:

Thus, the sum of the potentials emitted by the different stars at the different frequencies is

Here , where is the spatial position at which the function is evaluated and is the spatial position of the th star; the numbers , with , are the weights affected to the different frequencies (in Reference [12], these weights were noted , which now would clash with the functions 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 (with restricted to the plane , due to the axial symmetry), making a regular spatiotemporal grid . We thus want to determine complex functions , such thatwhere 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), , 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 —that is,(cf. equation (A.3) in the Appendix).(ii)Or, the dual (5) of that field—that is,(iii)Or still, the sum of these two fields, that is(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 of form (6), we discretize the integration intervals and the functions , and we use the “Simpson composite rule” [12]. We set , with some (arbitrary) reference frequency, and

With this discretization, the unknown functions in equation (9) are characterized by the complex numbers . Thus, the equality in the sense of the least squares, schematized by equation (9), actually becomes [12]with

The ’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):where and

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 the remainder) an exact solution of the source-free Maxwell equations (see §3.3 in Reference [12]).

The arguments of the Bessel functions and 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 . For a galaxy, these are numbers of the order of ; 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, iswhere . 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 . Let be some components of the EM field at some given spatial position . We have, with some complex coefficients, :whence

Unless , the time average is zero. Unless , the time average is zero. But does not occur because , and occurs only when . Hence, (26) gives indeed the additivity for the time average of the square of a component admitting an expansion (25):

For the EM field (10), its components admit expansions (25) with coefficients , , and , respectively, which are directly read on equations (20)–(22):where and

From (24) and (27), we thus get for field (10) a (time-averaged) energy density with a discrete spectrum:where the coefficients are given by equations (28)–(30), , and . 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 will be imposed to be the ones measured at our local position 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 given by equation (8) depends on the numbers , with . Each of them represents the relative weight of the corresponding frequency —though in some average sense because is common to all “stars” at the points . There is no exact way to assign values to these relative weights , which, together with the set of the positions of the stars, determine the functions on the r.h.s. of equation (9)—or more exactly, determine the discretized values in equations (14) and (18). In turn, the values , for a given value of the frequency index , determine the intensity of the -component of the scalar wave obtained after the fitting (the r.h.s. of (18)) and thus determine the intensity of the -component of the EM field (20)–(22). However, since we start only from relative weights , this determination can be only up to a multiplying factor, and since the very values of the ’s are assigned somewhat arbitrarily, that determination of the ’s is really up to a factor depending on , (this would be even clearer if we would separate the least-squares adjustment (18) into the different frequencies , which is actually easy and even would allow us to eliminate the common time dependence 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 by imposing that the values given by equation (32), calculated at our local position in the Galaxy, be equal to the ones measured. We take  kpc and  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.

Then, at any other spatial position (in the axially symmetric case considered), the model gives a prediction for the values of the SED corresponding to the frequencies (or the wavelengths ) considered in our calculation. The spatiotemporal grid used was (thus, 1050 events ), with varying from 0 to with , varying from 0 to  kpc, and from to  kpc. Other grids have also been tested, giving very similar results. Moreover, unless otherwise mentioned, the discretization number was (this is the number of subintervals in which any of the integration intervals is subdivided, see equations (13) and (14)). Figures 25 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).

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 ( kpc, Figures 4 and 5); and merely at short wavelengths, for the two positions that are farther from the galactic centre ( 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 (, or )—except for the position ( kpc and  kpc), that is, relatively speaking, very close to our local position ( kpc and  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 ( 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 (or rather ) 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 (for the radiation transfer model of Popescu et al. [6]). However, the and 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 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 for the present model and to understand the reason for it.

A surprising prediction of our model is the one for the maximum values of the SED in our model of the Galaxy—the SED being calculated at each point on a regular 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 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 axis, thus on the central axis of the galaxy (actually, the values are fairly independent of the value of in the domain investigated, i.e., ). However, if one refines the frequency mesh, going from to , and then to , one observes that the very high values of are obtained only at shorter and shorter wavelengths (Figures 7 and 8). A similar phenomenon is observed when one refines the discretization parameter : for a given frequency mesh, increasing reduces the domain of the high values of towards smaller values of , or, equivalently, decreasing expands that domain towards larger values of (Figures 911). The likely common reason for the two effects is that, when one refines either the frequency mesh ( increased) or the discretization of the integration intervals ( increased), then the number of the solved-for parameters increases. In most of the settings investigated in the present work, this number, , is actually larger than the number 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 ) and thus . We have for the grid used in this work, except for the two grids used for the first two curves in Figure 6, which give, respectively, and (we have 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 parameters for reason of insufficient numerical accuracy, so that their effective number was 150 to 200, thus significantly smaller than ). Higher values of the ratio correspond to a more overfitted situation, in which one expects the predictions to be less accurate. Now, on each among Figures 711, we saw that the domain of the low values of increases when increases, the number being fixed at for these five figures—thus we saw that the domain of the low values of increases when increases, i.e., when more overfitting occurs. Therefore, it is the low values of which are likely to be a numerical artefact, not the extremely high values. Note that this effect can be seen also in Figure 6: and , hence , are the same for the three curves. So the third and finest grid leads to a lower than do the other two grids and —and accordingly the domain of the high values of 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 from 76 to 23, the SEDs considered in Figures 25 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 : 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.

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  kpc and  kpc, coincides with the observations collected through different spatial missions [79]. 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.

Appendix

A. Extension of the Explicit Representation to a Finite Spectrum

In order to present this extension of the result [14] (the latter being recalled in Section 2.2), it is necessary to write this result in terms of mappings. Consider the vector space (respectively, ) made of those solutions of the scalar wave equation (respectively, made of those source-free Maxwell fields ) that are axisymmetric and time-harmonic with frequency (as the domain of definition of the fields in and , we take the whole spacetime. More exactly, since we assume axial symmetry, we assume that they are defined for the whole range of values of the spacetime variables , thus in the domain : (, )). We define two mappings and from into : given in , is the EM field given by equations (2)–(4), and is the EM field given by equations (2)–(4), followed by equation (5). Obviously, and are linear. We can see that belongs to , the subspace of with , whereas belongs to , the subspace of with . The uniqueness of the decomposition (1) is just that of the EM fields and (if they exist), such that the starting EM field is . This results immediately from the complementarity of the nonzero components, noted in Section 2.2, and inherent in the definition of and : the nonzero components of and must be just the corresponding components of , hence are known from the data of . Henceforth, for any EM field , we will note (respectively, ) the field whose components are zero (respectively, the field whose components are zero), the other components being those of . In contrast with the uniqueness, the existence of the decomposition (1) is a rather strong result, which amounts to state that the mapping from to defined byis surjective. This also amounts to state that, given any EM field , there exist and in , such that and .

The general axisymmetric solutions of either the scalar wave equation or the source-free Maxwell equations are obtained by summing time-harmonic axisymmetric solutions of the respective equations. We restrict ourselves to such solutions having a finite set of frequencies, say , for simplicity [12, 14]. The axisymmetric solutions of the wave equation with this set of frequencies form the vector space . The axisymmetric solutions of the source-free Maxwell equations with this set of frequencies form the vector space . The sums are direct because the decomposition in a finite sum of time-harmonic fields is of course unique. Therefore, the linear mappings and extend immediately to ones from to or, respectively, , by linearity, and so does the theorem. That is, for any EM field in the direct sum , the result valid for the time-harmonic case ensures that there exist scalar potentials: , such that

Hence, for any EM field , there are indeed two scalar potentials and in the direct sum , such thatand thus, finally,

In the totally propagating case, we have and , for some functions , with .

Data Availability

Previously reported SED data were used to support this study and are publicly available as scientific papers, references [2, 6], which are cited at relevant places within the text.

Conflicts of Interest

The author declares no conflicts of interest.