That is, at long enough time, the second term, the contribution of the initial value in the solution, Equation (5), is no longer significant. From observing the cosmic microwave background, this background temperature of the universe has been measured at 2.73 K [2]. This information could be used to fine-tune the system. Choosing S_{T} = 1 K and converting to dimensional quantities, then from Equation (4),

$\therefore \frac{{c}_{1}}{{c}_{2}}=\frac{Q}{h}=2.73\text{\hspace{0.17em}}\text{K}$ (7)

As c_{1} and c_{2} are found only as a ratio, we still cannot use Equation (5) to work out the histories as well as the initial temperature. Other assumptions are needed as described in the next model.

2) A zero net energy model

This popular model is based on the argument that although the universe contains billions of energy generation bodies such as our sun, the gravitational forces needed to keep the system together with other energy consuming processes use up equal amount of energy resulting in zero nett energy production. The reason that the universe still has a small temperature of 2.73 K is due to the remnant of its extremely high temperature past.

For a zero net energy universe, the constant, c_{1} = 0, in the field Equation (3) and also in its solution Equation (5). But, since T_{0} the initial value could be chosen arbitrarily, c_{2} can only be found by knowing a particular temperature at a specific time. The widely accepted assumption, based on distances of faraway stars or on the rate of expansion of the universe, is that the universe has an age of 13.8 billion years. There are widely different suggestions about the initial temperature; as an illustration, it is assumed that T_{0} = 10^{14} K. Based on S_{T} = 1 K and S_{t} = 1 year and using Equation (5), c_{2} could be found by solving:

${10}^{14}\times \mathrm{exp}\left(-{c}_{2}\times 13.8\times {10}^{9}\right)=2.73$ (8)

to give c_{2} = 2.273 × 10^{−9}. With c_{1} and c_{2} known, Equation (5) could be used to give the whole temperature history.

3) A composite model

The present observed temperature could be assumed to have come from both the heat source and the remnant of past temperature. In all cases, the system parameters, c_{1} and c_{2} could be fine-tuned so that the universe could start from any given initial value and cool down or heat up to the present temperature after a specific period of time. If T_{0} = 10^{14} K and the age is 13.8 billion years and the presently observed temperature of 2.73 K is coming from a combination of the remnant of T_{0} and heat source Q at different proportions, the system needs to be fine-tuned according to the values given in Table 1, where Percentage refers to the percentage contribution from the remnant of T_{0} and the last two columns are what the system parameters should be fine-tuned to give the present observed background temperature of 2.73 K.

As the time span involved is so large, it is not a surprise to see that c_{2} is varying over such a very narrow range. In fact, it is only from analytical solutions that c_{2} could be so determined.

4) Non-linear Models

Although the linear models could be fine-tuned to fit what have been hypothesized as the temperatures in the universe, namely the background temperature and the age, there is no justification that the system parameters, c_{1}/c_{2} and c_{2} could both be constant over such a huge temperature range. Non-linear models using temperature dependent parameters would be more likely to give better predictions. But it is unlikely that analytical solutions could be found for such nonlinear models. We shall use instead a linearized approximation by subdividing the time span into segments so that over each time interval the system parameters are replaced by constants that represent the average over that particular time interval. To fine-tune such a system, it is necessary that the temperatures at the beginning and at the end of each segment are known or could be estimated.

An example with c_{1} = 0 (that is assuming zero energy source.) and T_{0} = 10^{14} at t = 0 and temperatures at each other period are as given, the values for c_{2} found are shown in Table 2. It should be noted that the chosen time and temperature values in Table 2 are compatible with those often used in Big Bang theory [7].

These results show that the cooling rate would be very high at the beginning and slowed down to a small value only toward the end of the current time history of the universe. To lower c_{2} at the beginning, a plausible model may assume that there is a rapid expansion that consumes so much energy leading to a negative c_{1}. For example, if c_{1}/c_{2} = −10^{11}, then for the first period between 0 and 180 seconds, the value for c_{2} = 1.21 × 10^{6}. However, the above is only one of unlimited number of possibilities. In fact, the system could start with a zero temperature providing that the heat source is positive as in the cases shown in Table 3, for c_{2} = 1.

5) Low temperature models

The problem with high temperature models is that it is difficult to explain how the temperature could be cooled down so rapidly from infinite to say 10^{14} K in just 10^{−42} second. In fact, if the initial temperature is infinite as postulated, there is no way that this temperature could be lowered, no matter how fast and how much energy, that is any amount less than infinite, is removed, by expansion or any other means. Because it is difficult to justify both physically and mathematically, dark matter and dark energy, both unknown to us, have been suggested as the best plausible explanations.

It is much easier to assume a low temperature model. For example, taking the energy source to be the same as in Case (a), that is c_{1}/c_{2} = 2.73, and cooling coefficient in Case (b), that is c_{2} = 2.273 × 10^{−9} the time requires for temperature to

Table 1. Fine-tuning a composite model.

Table 2. Fine-tuning a non-linearized model.

Table 3. Temperature rises from zero at the beginning.

cool down from 1000 K to 3 K is only 114.6 years. If starting at zero temperature, the time taken for the temperature to heat up to 2.70 K is 62.9 years.

It is better still to assume that the starting temperature is 2.73 K, the same as the observed temperature now and Q = 2.73 h, then the universe is always maintained at a constant temperature of 2.73 K all the time.

An often used justification for a high temperature start is based on the second law of thermodynamics. As the temperature is so high, the initial entropy is very low. Then, all spontaneous processes supposed to take place subsequently should bring down the temperature and hence an increase of entropy. But it should be noted that out of all the cases studied above from (a) to (e), with the case of a zero starting temperature as an exception, energy is coming from a higher temperature than the loss term associated with h. The implication has been that entropy is always increasing as requires by the second law.

4. Electromagnetic Waves Transmission

Consider special electromagnetic waves called vortex solitons. Using cylindrical coordinates that is natural for this type of wave transmission, the nonlinear Schrödinger equations [4] in the spatiotemporal domain are:

$\begin{array}{l}i\frac{\partial u}{\partial z}+\frac{1}{2}\left(\frac{{\partial}^{2}u}{\partial {r}^{2}}+\frac{1}{r}\frac{\partial u}{\partial r}+\frac{1}{{r}^{2}}\frac{{\partial}^{2}u}{\partial {\theta}^{2}}\right)+u\ast v-\sigma \left({\left|u\right|}^{2}/4+2{\left|v\right|}^{2}\right)u=0\\ 2i\frac{\partial v}{\partial z}+\frac{1}{2}\left(\frac{{\partial}^{2}v}{\partial {r}^{2}}+\frac{1}{r}\frac{\partial v}{\partial \theta}+\frac{1}{{r}^{2}}\frac{{\partial}^{2}v}{\partial {\theta}^{2}}\right)-qv+{u}^{2}/2-\sigma \left(4{\left|v\right|}^{2}+2{\left|u\right|}^{2}\right)v=0\end{array}$ (11)

where u and v are the amplitudes of the coupled electromagnetic waves, q is a coupling coefficient, σ is the nonlinear coefficient and z the distance travelled.

The above equation could be solved by the Lanczos-Chebyshev pseudospectral method [8], seeking solutions in the form:

$\begin{array}{l}u\left(r,\theta ,z\right)={\displaystyle {\sum}_{i=0,M}{\displaystyle {\sum}_{j=0,N}{u}_{ij}\left(z\right){r}^{i}{\theta}^{j}}}\\ v\left(r,\theta ,z\right)={\displaystyle {\sum}_{i=0,M}{\displaystyle {\sum}_{j=1,N}{v}_{ij}\left(z\right){r}^{i}{\theta}^{j}}}\end{array}$ (12)

For the given computational 2D domain, an affine transformation is carried out to rescale it into new coordinates $\left\{r,\theta \right\}$ to $\left[-1,+1\right]\times \left[-1,+1\right]$. Spatial-discretization of Equation (11) is carried out by using collocations at points corresponding to the roots of Chebyshev polynomials:

$\begin{array}{l}{r}_{i}=-\mathrm{cos}\left[\frac{\left(2i+1\right)\pi}{2\left(M-2\right)}\right],\text{\hspace{0.17em}}\text{\hspace{0.17em}}i=0,\cdots ,M-2\\ {\theta}_{j}=-\mathrm{cos}\left[\frac{\left(2j+1\right)\pi}{2\left(N-2\right)}\right],\text{\hspace{0.17em}}\text{\hspace{0.17em}}j=0,\cdots ,N-2\end{array}$ (13)

This approach would be the most economical choice numerical-wise; more details about the Lanczos-Chebyshev psudospectral method may be found in Reference [9].

Equation (11), together with the boundary conditions, is then transformed to a set of first order ordinary differential equations in matrix form of

$iA{U}_{z}-LU-H\left(U,z\right)=0$ (14)

where $U$ is a vector ${\left({u}_{11},{u}_{12},\cdots ,{u}_{1M},\cdots ,{u}_{N1},{u}_{N2},\cdots ,{u}_{NM},{v}_{11},{v}_{12},\cdots ,{v}_{1M},\cdots ,{v}_{N1},{v}_{N2},\cdots ,{u}_{NM}\right)}^{\prime}$, $A$ and $L$ are linear matrices, $H$ is a nonlinear vector.

Using the unconditional, stable and implicit Crank-Nicholson algorithm, we set

$iA\left({U}^{m+1}-{U}^{m}\right)-\frac{\Delta z}{2}\left[L\left({U}^{m+1}+{U}^{m}\right)+H\left({U}^{m+1},{z}^{m+1}\right)+H\left({U}^{m},{z}^{m}\right)\right]=0$ (15)

where the superscript m refers to the time-step number. Since Equation (14) is nonlinear, we let

${U}^{m+1,0}={U}^{m}$, then ${U}^{m+1,r-1}\to {U}^{m+1,r}$ (16)

where the superscript r refers to the iteration number and the arrow implies that integration is used. It should be noted that the iteration involves only one nonlinear $H$ term and since $A$ and $L$ are linear, only one inversion is required for all the time steps:

${U}^{m+1,r}={\left[iA-\frac{\Delta z}{2}L\right]}^{-1}\left\{\left[iA+\frac{\Delta z}{2}L\right]{U}^{m}+\frac{\Delta z}{2}\left[H\left({U}^{m+1},{z}^{m+1}\right)+H\left({U}^{m},{z}^{m}\right)\right]\right\}$ (17)

The system is a nonlinear initial value problem. In theory any set of arbitrarily chosen initial values could be used as non-equilibrium components present initially would disappear through the marching process. But some unwanted components may grow instead. Therefore, numerically, to ensure success, it is necessary to select as good as possible a set of initial values and also to use a small enough step size. In our case, we took advantage that an approximate solution could be obtained more easily. If we consider that

$u\left(r,\theta ,z\right)=\mathrm{exp}\left(iS\theta \right)f\left(r,z\right),\text{\hspace{0.17em}}\text{\hspace{0.17em}}v=\mathrm{exp}\left(2iS\theta \right)g\left(r,z\right)$ (18)

where $S=0,1,2,\cdots $ is the integer vorticity. Substituting Equation (18) into Equation (11) gives,

$\begin{array}{l}i\frac{\partial f}{\partial z}+\frac{1}{2}\left(\frac{{\partial}^{2}f}{\partial {r}^{2}}+\frac{1}{r}\frac{\partial f}{\partial r}-\frac{{S}^{2}f}{{r}^{2}}\right)+f\ast g-\sigma \left(\left|f\right|/4+2{\left|g\right|}^{2}\right)f=0\\ 2i\frac{\partial g}{\partial z}+\frac{1}{2}\left(\frac{{\partial}^{2}g}{\partial {r}^{2}}-\frac{2iSg}{r}+\frac{4{S}^{2}g}{{r}^{2}}\right)-qg+{g}^{2}/2-\sigma \left(4{\left|g\right|}^{2}+2{\left|f\right|}^{2}\right)g=0\end{array}$ (19)

We used the solutions for S = 0 as our initial input to start our iterative process.

For a nonlinear problem like the one at hand, it has to keep in mind that stationary solutions only exist over a certain combinations and ranges of the system parameters. It is also possible that multiple solutions could exist.

We have obtained a set of stationary vortex solutions as shown in Figure 1, where |u|^{2} is plotted over a {x,y}-plane. Figure 2 shows the contours of the same set of solutions. As a demonstration that stationary solutions have reached, Figure 3 shows that the pulse energy, E, for both u and v remain the same over the distance travelled, where

$\begin{array}{l}{E}_{u}={\displaystyle {\int}_{-20}^{20}{\displaystyle {\int}_{-20}^{20}{\left|u\right|}^{2}\text{d}x\text{d}y}}\\ {E}_{v}={\displaystyle {\int}_{-20}^{20}{\displaystyle {\int}_{-20}^{20}{\left|v\right|}^{2}\text{d}x\text{d}y}}\end{array}$ (20)

It should be noted that Equation (11) does not include a term for transmission loss. In reality, due to transmission loss both u and v would decrease in their intensities along z, needing periodic amplification to restore their pulse energy in long distance transmision.

Figure 1. Waveform of |u|^{2} at z = 5 with σ = −1 and q = 1.

Figure 2. Contour plots of the waveform shown in Figure 1.

Figure 3. The stationary propagation of 2D vortex solitons.

5. Discussions

The stationary solutions for Equation (11) are complex electromagnetic waves. As they are complex quantities they exist only in abstract mathematics. In practice, we can only measure the absolute values of the pulse’s intensities or the pulse energies. But it is much easier to study electromagnetic waves, at least in this particular case, through the use of abstract complex algebra. No one would question whether complex electromagnetic waves do exist because they could not be measured directly. Therefore, a real event may be modelled using one or many abstract mathematical fields with abstract algebraic quantities as long they can be projected back to the real physical system that can be observed and measured.

It should be noted that Equation (11) contains terms associating with 1/r and this equation is used to cover a domain including the point at r = 0. These terms will become infinite unless the fields there have zero value. This is a necessary condition only due to the physical limitation that there is no infinite energy to support such infinite wave intensity. In a mathematical model, there is no such limitation as such a point could be replaced numerically by a nearby point where the field has a very large value. This limitation should be noted when a mathematical model is used to replace the physical event.

Equation (1) is much simpler than Equation (11). But both require integration along one variable, t, the time, in the case of Equation (1) and z, the distance travelled, in the case of Equation (11). We could solve Equation (1) analytically and there no problem to start integration forward from time zero or to start integration backward from a later time. The use of a different starting value could, however, give a completely different set of solutions. Equation (11) needs to be solved numerically. For a nonlinear problem like the one considered here, it is unlikely that using a solution for a later history could be used to recover the initial condition by integrating backward.

It should be noted that both equations are based on energy balance.

As with all other mathematical models, Equation (1) leads to Equation (3) in that two system parameters, c_{1} and c_{2} need to be found from observations or measurements, as we cannot find directly the material properties involved in Equation (1). Therefore, two sets of data are needed. Unless we can find the analytical solution, it is not possible to use one set of data and integrate forward or backward in time to discover the other set as reported in Reference [10].

For a more comprehensive model, more terms, or even equations, could be added to Equation (1). For example, an equation for Q could be developed to account for not just energy production but also energy changes due to expansion or gravity.

6. Conclusions

Our energy balance model has been used to show that the present observed background temperature of the universe need not be the remnant of a very high temperature past. The initial background temperature of the universe could be quite similar to what it is now.

Our model could accommodate an infinite initial temperature. This is done by following standard numerical technique: that is to replace the initial singular point by a very nearby point, for example, using time at 10^{−35} second when temperature has assumed to have cooled to 10^{27} K.

We have explained that an infinite temperature could exist in a mathematical model. In reality, this is not possible as such a temperature needs to be supported by an infinite amount of energy and this is physically not possible. Nevertheless, the universe’s huge mass, if converted, is near to infinite equivalent energy in mathematical terms; we are unable to explain, both mathematically and physically, from where this amount of energy could have come.

Cite this paper

Chen, P. (2019) An Energy Balance Simulation of the Universe.*Applied Mathematics*, **10**, 956-966. doi: 10.4236/am.2019.1011067.

Chen, P. (2019) An Energy Balance Simulation of the Universe.

References

[1] Carroll, A.O. and Ostlie, D.A. (2017) An Introduction to Modern Astrophysics. Cambridge University Press, Cambridge.

https://doi.org/10.1017/9781108380980

[2] NASA/WMAP Science Team. Our Universe? NASA.

https://map.gsfc.nasa.gov/universe/universe.html

[3] NASA/WMAP Science Team. Discovery of the Microwave Background. NASA.

https://map.gsfc.nasa.gov/universe/bb_tests_cmb.html

[4] Malomed, B.A. and Kevrekidis, P.G. (2001) Discrete Vortex Solitons. Physical Review E, 64, 026601.

https://doi.org/10.1103/PhysRevE.64.026601

[5] Neshev, D.N., et al. (2004) Observation of Discrete Vortex Solitons in Optically Induced Photonic Lattices. Physical Review Letters, 92, 123903.

https://doi.org/10.1103/PhysRevLett.92.123903

[6] Swartzlander Jr., G.A. and Law, C.T. (1999) Optical Vortex Solitons Observed in Kerr Nonlinear Media. Physical Review Letters, 69, 2503.

https://doi.org/10.1103/PhysRevLett.69.2503

[7] NASA/WMAP Science Team. Foundation of the Big Band, NASA.

https://map.gsfc.nasa.gov/universe/bb_theory.html

[8] Chen, P.Y.P. (2016) The Lanczos-Chebyshev Pseudospectral Method for Solution of Differential Equations. Applied Mathematics, 7, 927-938.

https://doi.org/10.4236/am.2016.79083

[9] Chen, P.Y.P. and Malomed, B.A. (2012) Lanczos-Chebyshev Pseudospectral Methods for Wave-Propagation Problems. Mathematics and Computers in Simulation, 82, 1056-1068.

https://doi.org/10.1016/j.matcom.2011.05.013

[10] Chen, P. (2019) On Transient Simulation of Field Equations. Applied Mathematics, 10, 719-727.

https://doi.org/10.4236/am.2019.109051