Application of Minimum Energy Effect to Numerical Reconstruction of Insolation Curves

: Increasing the efﬁciency of the solar energy harvesting system is an urgent need in light of the climate changes we live in nowadays. The most signiﬁcant data to be processed in the photovoltaic harvesters are the curve of solar radiation intensity to achieve the maximum beneﬁts of the solar incident light. This processing contains complicated procedures, and the used algorithms are also high computational power-consuming which makes using special software and high potential hardware essential requirements. An explanation of the Minimum Energy Effect method is presented in this article. Our proposed algorithm uses this method to provide a simple and high-accuracy mathematical tool for generating a simple alternative curve instead of the complicated original nonlinear curve of solar radiation intensity. The produced curve is suitable for further operations, such as derivatives, integrals, or even simple addition/subtraction. Our algorithm provides a gradual procedure to ﬁnd an optimum solution of the equation system, unlike the iterative methods. In addition, the results of analyzing the effect of time-division density showed the relationship between the speed of solving the task and the accuracy of results.


Generally on Solar Insolation
There are currently a number of ways and technologies for using solar energy. In this study, we will focus only on those that come into consideration in terrestrial conditions. One large category is the direct use of the thermal energy of solar radiation falling on the Earth's surface for heating (usually water as an accumulation medium, but also, for example, cooking surfaces) or drying (linen, agricultural products, seawater in salt mining), or lighting areas below the surface or inside buildings using different types of light guides. The second group is the photovoltaic conversion of solar energy into electricity, usually in conjunction with battery storage (either local or central).
Limiting factors for photovoltaics include the following: • Costly construction; • Relatively low efficiency; • Limited service life; • Environmental issues associated with subsequent disposal; • Usability depending on geographical location (especially latitude); • Unstable output power in connection with the pollution of panels, but especially with weather conditions (clouds and precipitation, air pollution of natural origin, and due to human activity).
Photovoltaics is the most important area in the use of solar irradiation/insolation data on the Earth's surface. Operators of individual power plants and the entire energy system are usually primarily interested in the forecast so that they are able to respond well in advance to expected fluctuations in the generated power [1]. An example of a service that provides such predictions is Solcast [2]. Different models use, for example, curve approximations [3][4][5] or neural networks [6,7]. Historical measurement data from meteorological stations and satellites are then used for verification.

Processing of Solar Data
The authors use splice interpolation to set an interval in the solar radiation curve; then, they use the Discrete-Calculus Technique to calculate the space under the curve in the specified interval [8]. The shape produced by this method is a trapezoid, so the corresponding method for calculation is used, and the obtained area is multiplied by solar panel efficiency to get the equivalent output electrical energy [8].
In [9], the authors used a data set of the global solar irradiance in a specific location on the Earth to calculate the solar power generated by solar panels. They calculated the average of the mean hourly solar irradiance.
The authors in [10] used a PV panel connected through a 50-ohm resistor for comparison between the output pattern and the open-circuit voltage. Then they used MATLAB to derive an equation from the insolation data collected from the solar lab and the PV output voltage. The obtained equation gave a regression analysis of 98.9%, which is highly accurate for the calculation of insolation.
The authors used a statistical representation to analyze the insolation graph, which is presented by the means of average and the standard deviation. First, they obtained the irradiation data from meteorological stations; then, they calculated the monthly average of the solar source and the intermediate monthly means for each selected year [11].
In [12], the authors used the D10 wavelets method to analyze the solar radiation data collected at Universiti Teknologi PETRONAS (UTP). The basic equation contains an approximation coefficient, a detail coefficient, a low-pass filter, and a high-pass filter. Then a threshold equation was used to define the global threshold value because it is important to choose the proper threshold value for solar radiation data processing. The authors dropped the value of the reflected solar radiation because it was small compared to the value of the direct solar radiation.
A data fitting method to apply smoothing and filtering operation to the collected global solar radiation data was used in [13]. The authors experiment with several degrees of polynomial curve fitting and calculate the RMSE (root-mean-square error) and R-squared (the coefficient of determination), which present the quality of data fitting. Then, a mathematical model was presented to predict the missing data, being based on the quadratic polynomial fitting.
Machine learning regression algorithms to analyze solar irradiation, which was the Gaussian process regression (GPR), are presented in [14]. Their contribution consists in using the time domain in their model because the previous machine learning models were using just the spatial domain; thanks to that, the computational efficiency of the model has been improved.
Typical applications have also been presented in [15], where the authors used mean and standard deviation on the insolation to find the optimal size of a photovoltaic battery system. In [16], the authors used a new algorithm to analyze the insolation and optimize the most important statistic in the insolation data by using a neural network for applications of Building Integrated Photo Voltaic (BIPV).

An Overview of Numerical Methods for Solving Solar "Field" Distributions
Numerical methods for solving technical problems can be categorized according to many (often inconsistent) criteria. For the research we are presenting in this article, the concentration (or vice versa: distribution) of the parameters of the investigated events is an important factor. When concentrating parameters in the "center of gravity" of the whole task-or its part-we get into the domain of equivalent circuits. Here, solutions are sought by computer simulations, and the results are very often presented as a function of time. Typical examples of such simulations are [17] where a photovoltaic harvester equivalent circuit, quartz piezoelectric harvester equivalent circuit, equivalent circuit of an electromagnetic force harvester, etc., are presented.
The importance of simulations based on the principle of equivalent circuits certainly cannot be questioned. It is a powerful tool available to technicians and scientists. Nevertheless, our approach is focused on simulations of events distributed in "space" and "time". These can be characterized as field resolution methods. That is why the next part of this section will focus on a more detailed overview of how methods analyzing the distribution of fields in the energy harvesting (EH) area are applied.
The basic and historically oldest numerical method based on differential calculus (and subsequent development over the Taylor series) is the Finite Difference Method (FDM). Countless examples can be found of using this method in the field of EH. A typical example may be [18] where, like many others (e.g., [19,20]), the authors deal with the analysis of piezoelectric systems.
The Finite Element Method (FEM) is also an extremely widespread method. It is undoubtedly possible to declare this method to be the most frequently used today. Its principles are sufficiently clearly described in [21]. The penetration of FEM into the EH area is well illustrated by more than 530,000 records offered by one of the major web search engines. The most cited articles include [22][23][24]. However, it should be noted here that the research of these authors is directed to the field of chemistry rather than to electrical engineering or electronics. A typical example of the FEM application in the field of electrical engineering is [25]. The authors of the article analyzed the piezoelectric harvester. The responses of piezoelectric actuators were analyzed by FEM in [26]. Special application of an energy harvester device was analyzed in [27]. The issues of solar radiation were discussed in [28].
In addition to the above two methods of solving arrays, it is advisable to add a third one: Boundary Element Method. This, unlike the former two, is based on integration principles. The required mathematical complexity of the method ranks it among the least used. The method is in principle, particularly suitable for open-type tasks [29,30]. Quite often, this method is used in a combined mode with FEM [31].
All three of the above methods are robust and computationally intensive. The further use of the results (postprocessing) provided by these methods is also complicated as the used software solutions are closed.
The numerical reconstruction of the solar radiation intensity curve was dealt with in [32] (including its theoretical analysis) or [33] (design of the solar simulator). The authors rather analyzed the effect of altitude on solar radiation, including altitudes corresponding to Low Earth Orbit (LEO) for satellites. In [34], the numerical simulation was focused on the analysis of solar radiation for the purpose of identifying possible heat or sound sources. The vast majority of articles published in this field focus on the technologies for or ecological properties of solar cells. An overview article summarizing the situation in the area is [35].

Motivation
Relatively extensive research by many authors is linked to the processing of solar radiation intensity (a general overview, including the methods used, is provided above). With respect to the stated goals, these authors used various more or less sophisticated mathematical methods and operations to process the curves of solar radiation-from the analytical calculation (see Section 2), which is relatively complex due to multiple nonlinearities, to the various statistically-based interpretations of the individual parts of the curve. Our objective was to provide these authors with another tool that will be mathematically simple, and which will replace-with high accuracy-the original nonlinear curve of solar radiation intensity with a curve that will be more suitable for further processing. Due to the partial linearization of the curve, the subsequent use of the data provided by us is relatively simple. This applies to practically all required operations, such as derivatives, integrals, or even simple addition/subtraction.

Structure of the Article
The structure of the article is as follows. Section 1: Introduction: an overview of the current state in the field of solar energy with special regard to the further processing of the curve of solar radiation intensity. An overview of numerical methods commonly used in the field. Section 2: Introduction to the Minimum Energy Effect method. Section 3: Introduction to the basic mathematical apparatus used in analytical calculations of solar radiation intensity.
Section 4: Application of the Minimum Energy Effect method to the field of solar radiation intensity analysis. Quadratization of the function (Section 4.1). Division of the quadratization region into three time periods (Section 4.2). Design of an algorithm for solving the problem divided into several elements (Section 4.3).
Section 5: A case study of the MEE application for a specific day of the year and for a specific location on the surface of the Earth.
Section 6: Discussion on results that we have achieved and challenges for future work. Section 7: Conclusions.

Minimum Energy Effect
Many different simulation tools and computer-aided design software products are currently available to technicians and designers. Some of them are focused on real-time simulations, where computer software evaluates data based on a mathematical model. A characteristic feature of this type of simulation is that it is possible to enter into ongoing calculations and change some of the input variables. Sufficient computing capacity then provides calculations of the system response in real-time. A typical example of this logic is the MATLAB-Simulink software. These are mostly 1D problems-mainly from the point of view of the dependence of the investigated quantities on time.
On the other hand, there are software simulation products that are based on different physical interpretations of the problem. Here we usually talk about solving fields (for example, electromagnetic ones). In general, the currently available computing capacity for simulation does not allow real-time calculations of this type to be performed. Typical representatives are, for example, the Boundary Element Method (a method based on transforming a problem through integral calculus), the Finite Element Method (based on differential calculus), or the Finite Difference Method (based on the differential calculus of Taylor expansion). All these methods have one thing in common: they look for the optimal solution for the distribution of fields in the area of interest and offer the solution that best meets the initial conditions of the task on average. The popular least squares method, which minimizes the sums of the squares of the deviations of the proposed solution from the assignment, behaves similarly in this respect.
A completely different principle applies to the often-neglected method of Minimum Energy Effect (MEE). The principles of this method were described, for example, by Richard Feynman in [36]. Feynman illustrated the reformulation of the physical interpretation of the problem using the example of a circle: the mathematical definition of a circle is clearly given-it is a set of points in a plane, the distance of which from the center point is constant. The reformulation of the problem into the optimization plane consists of the definition of a circle so that it is a "curve" of a given (understand: constant) length, which has the task of encircling the maximum area.
According to [36], the distribution of a physical field (electromagnetic, gravitational, etc.) will always be such that it meets the condition of minimum energy. However, the concept of minimum energy is often used incorrectly, as it does not concern energies from the physical point of view. Therefore, we will call this condition the energy effect, and the method based on it is the Minimum Energy Effect method. A specific example is demonstrated in [36], focusing on the motion of a body in the gravitational field of the Earth.

Feynmann's Approach to the Problem
In his lectures on physics, Feynman focused on motion in the gravitational field of the Earth. Feynman processed the kinetic (W k ) and potential (W p ) energy of motion in such a way of evaluating its difference W diff [36]: In this step, it is important to note that the mass m occurs in both members of the difference. Therefore, its value does not play a role in finding the extreme (minimum or maximum) of the investigated function.
According to [36], an important element is the effect S of the differential energy W diff in the time interval from t 1 (start time of the event) to t 2 (end time of the event) [37]: The analytical solution of this integral is [37]: After substituting specific values, the value of the energy effect is obtained. What is important is the fact that this value is minimal for the actual curve of the physical process. The effect of energy intensity (the principle of the least action) will logically manifest itself in such a way that in any other curve of the investigated event (i.e., in this case, the path of movement, would be other, different from the real one), the value of the energy effect will be higher. At the same time, the greater the deviation from the original function, the greater the value of the energy effect. We will also use this principle in the numerical reconstruction of the investigated curve. Instead of the original quadratic function, we will use linear functions, which are simpler and less demanding for possible further processing. This means that if we use a linear approximation, the subsequent mathematical operations (such as integrals or derivatives) can be replaced by the operations of multiplication (division) and addition (subtraction) that are time-saving and more suitable for digital computational techniques. The main point of this method is that the optimization of linear (or any other) replacement of the original event is possible without knowledge of the original course of the sought quantity, i.e., without knowledge of the global minimum value of the energy effect. We get a typical example of an optimization task: finding a new global or local minimum for several variables. Given the theory presented in Section 4.3 of this article, it is also possible to circumvent this step and significantly simplify the calculation procedures.
A detailed analysis of the characteristic properties of the Method is presented in [37]. The situation is illustrated in Figure 1, where three points are defined: The coordinates of these two points (t 0 , h 0 ; t 2 , h 2 ) are clearly given-these are the boundary conditions of the task. Let us define a single internal point C, which we place on the time axis to the value t 1 : Our task is to find the optimal value of h 1 for point C. It is obvious that the time axis is divided into two intervals (elements): <t 0 , t 1 > and <t 1 , t 2 >. At each of these intervals, the required linear approximation is expressed as follows: At the same time, for reasons of physical continuity of the sought quantity, it is The energy effect for linear waveforms in combination with Equation (2) is expressed as: Through a detailed analysis of the described event, the authors concluded in [37] that the method of Minimum Energy Effect tends to reconstruct the value of the original quantity at the selected point. This means that when moving point C up or down, the value of the energy effect is always higher than it would be if the height of point C was set to the original (blue) curve. From the mathematical and system point of view, this is a very interesting fact, which will be used to present a possible scenario for mapping the curve with divisions on a larger number of time intervals (elements), which bypasses computationally intensive optimization of a problem with a large number of degrees of freedom (see Section 4.3).

Analytical Calculation of Solar Radiation
The basic tool for assessing the usability of solar energy are well-known relationships for calculating the intensity of solar radiation falling on the Earth's surface. If we neglect disturbing influences (typically current atmospheric conditions or reflections and shading caused by surrounding objects), we must take into account only latitude (or longitude if we want to work with local time) and the day of the year. Furthermore, we can assume the ideal impact of solar radiation is always perpendicular to the measured area, regardless of the angular height of the Sun above the horizon (zenith angle) and the slope of the terrain or other ground (for example, when using automatic positioning of the solar panel).
The fundamental analytical formula expressing the whole daily extraterrestrial radiation from the sunrise to the sunset is as follows [38]: π I sc 1 + 0.033 cos 360n 365 cos L cos δ cos h s + 2πh s 360 sin L sin δ , where: I SC the solar constant (1353 ± 21) [W/m 2 ]; n the day of the year that counts from 1 January; L local latitude; δ declination; h s the time of sunrise (sunset) beyond the horizon.
On the other hand, there exist already implemented solutions that transform Equation (5) as well as other necessary equations from [38] to the required insolation waveforms. One of the solutions can be found on the website [39]. Here, after entering the input conditions (latitude, day of the year, etc.), it is possible to obtain not only a graphical representation of the ideal curve of insolation but also tabular data for individual parts of the day in sufficient resolution. These materials served as input data for subsequent numerical solutions.

Partial Quadratization
In order to be able to use the MEE method in the field of insolation analysis, it is necessary to look for the shape similarity of the waveforms with those presented in Section 2. The similarity with the quadratic function is obvious here, but the theoretical nature of the curves is completely different (see Section 3). The question is, what is the degree of similarity between these curves? For this reason, we sought suitable approximations of insolation curves through a trend line. The analysis was not performed for each day of the year-we selected days number 35 (4 February), 80 (21 March), 126 (6 May), 172 (21 June), 219 (7 August), 266 (23 September), 310 (6 November), 355 (21 December). Due to their distribution in the calendar year, they sufficiently represent all seasons and relevant possible states. In Figure 2, it is possible to compare the original distribution of insolation intensity (bold blue curve) with the approximation (thin curve) for a specific day No. 126 (6 May). The points of the curve correspond to the values of insolation according [39]. Figure 2 is completed with real measured data (red curve) provided by the Czech Hydrometeorological Institute; however, due to the actual weather conditions, day 131 (measured in 2015) is represented in the illustration. It should be noted here that although this approximation led for some days to a relatively high value of reliability R 2 (in this case, R 2 = 0.967 for day No. 35), this is not the case for all days of the year. For some days, this value decreased to 0.93, and similar values (specific values of R 2 for such an approximation are provided in Table 1). Therefore, this form of this step is not applicable universally to the subsequent procedure. We decided to improve the accuracy of the approximation (the criterion is the value of the reliability coefficient R 2 ) by dividing each day into three intervals: 1. from sunrise to 9:00 a.m.; 2.
from 15:00 p.m. to sunset. Note: The reliability coefficient R 2 states how many percent of the variance is explained by the approximation and how many is not. For example, at a value of R 2 equal to 0.99, it means that 99% of the variance of the investigated quantity falls within the applied approximation. For each of the three intervals, it is possible to approximate a concave quadratic function in the following form: Another possible problem is the unidirectional component in a given course, both in the direction of the x-axis and in the direction of the y-axis. This problem can be evaded by defining a new relative coordinate system with the center at a point corresponding to sunrise, or more precisely, the first non-zero value of solar radiation intensity. An example of such offset is shown in Figure 3. By setting this point to the beginning of the coordinate system, Equation (6) is simplified to the form: which is equivalent to the equations provided in [36]. We must always keep in mind that-from the point of view of accuracy-it is more advantageous to divide the day into three independent intervals and to find the numerical solution of the equation virtually independently of each of these intervals.

Application to Insolation Curve
As described in Section 4.1, there are two procedures for solving the problem: (A) A universal solution covering the period from sunrise to sunset: (B) The solution by dividing the period of sunshine into three intervals.
The procedure (A) is seemingly simpler; on the other hand, the inaccuracy of the quadratic approximation (R 2 ranging from 0.93 to 0.98) of the solved curve can lead to local inaccuracies, which can fundamentally affect the results of the analysis. Therefore, we preferred procedure (B), i.e., dividing the approximated curves into three intervals. Here, the value of R 2 was always higher than 0.99, so the approximation can be considered sufficiently accurate. The complexity of applying the mathematical apparatus to the type of problem formulated in this way consists in the fact that it is necessary to determine the parameters a 1 and b 1 in Equation (7) for essentially each interval. In principle, it is then possible to search for a numerical solution for each of the three intervals independently of each other and thus reduce the complexity of the algorithm used. Each interval is approached so that the curve is shifted (see Figure 3 in Section 4.1) to the origin of the coordinate system (the "DC" component of the waveform is always zero).

Principle of the Method Application in a Single Time Interval
The simplest example of MEE application in a single time interval is illustrated in Figure 1. For a universal solution, according to Section 4.1, point A corresponds to sunrise and point B to sunset. For a solution by dividing the period into three intervals, points A and B may represent, for example, the values for times 9:00 (point A) and 15:00 (point B). The value of the sought quantity is known for these points, and therefore these two points are considered as boundary conditions of the problem. A solution (value) is sought for the central point C. Its position in time is at the choice of the operator-this basically creates the division of the task into time intervals-elements. However, the value at the selected time is already calculated by the MEE system so that the height of point C tends to be placed on the original (blue) curve. Subsequently, the values at any other time are approximated by the linear approximation according to Figure 1.
The course of the effect of energies depending on the value at point C is a simple task that can be performed by virtually any mathematical apparatus. In [37], the authors of the article proved that a MEE system really tends to reconstruct the original curve and thus set the value of the approximation in accordance with reality. Linear interpolation, of course, introduces considerable inaccuracy into the system-this manifests itself (see Section 5) in a significantly different value of the energy effect compared to the original curve. The logical step to achieve more accurate results with a given access logic is to divide the timeline into a larger number of intervals-elements. We present the logic of solving the problem for such a case in Section 4.3.

Possible Approaches to Solution
Taking into account the stricter requirements for accuracy, the logical step is to divide the timeline into a larger number of elements. Consequently, the demands on the mathematical apparatus capable of solving the optimization problem formulated in this way increase. Possible approaches are presented through the example in Figure 4, where the problem has boundary conditions (fixed points A and B) and three internal points (C 1 , C 2 , and C 3 ). The solution algorithm for a large number of elements can be based on either of the following principles: (A) Calculation of a system of equations "at once". This is a "classical" approach to solving an optimization problem. A relatively sophisticated mathematical apparatus is needed. The "optimal" solution of the system of equations is sought either by direct calculation (with a smaller number of variables) or various iterative methods are used (Newton-Raphson and the like), which do not find an "exact" solution to the problem, but at least a solution close enough to the optimal one. This procedure is illustrated in Figure 4, where the system tries to set the values of all internal points C 1 -C 3 at once, based on the defined criteria.
(B) Recursive solution. This solution is based on the already mentioned feature of the MEE system, the trend of which is to adhere to the values of the original investigated curve. The principle of this solution is that the value of the sought quantity at point C 1 should remain the same regardless of the presence of the other points C 2 and C 3 . Thus, an identical solution at point C 1 can be achieved by solving an equivalent problem according to Figure 5a, which is based (except for the boundary conditions in points A and B) on a single internal point C 1 (so it can be solved directly according to Section 4.2). After solving this problem, the assignment itself will change: points C 1 and B will be considered as boundary conditions. For this reformulated problem, a solution is sought again only one internal point C 2 (Figure 5b). After finding the solution, the boundary conditions are redefined to points C 2 and B, and a solution is sought for C 3 (Figure 5c). This procedure is continued as necessary until a solution is found for all internal points of the system. The advantage of the recursive solution is the simplicity of the apparatus used (and consequently of finding an exact-i.e., not approximated by iterations-value of the solution). On the other hand, a possible error is "chained" for the following points and can more fundamentally affect the overall solution results.

Method Verification
To verify the system, we used day No. 126 (6 May) as an example, because in the conditions of Central Europe, it has a relatively long period of daylight: sunrise is about 4:51 and sunset about 19:04. The value of the reliability coefficient for the approximation of the all-day curve by a single quadratic function was R 2 = 0.933. Therefore, it was appropriate to use the division of the day into several time periods-in accordance with the description provided in Section 4.1. We used three time periods (from 4:51 to 9:00, from 9:00 to 15:00, and from 15:00 to 19:04). We will present the response of the MEE system to the middle of the given time periods, i.e., from 9:00 to 15:00. The boundary conditions for this period were as follows: Following Section 4.1, we adjusted these conditions so that the "DC" component (with respect to the global coordinate system) was removed, i.e., we shifted the time by −9 h and the radiation intensity by −0.92595772 kW/m 2 (so that point A was located at the beginning of the coordinate system): For this case, the approximation using the quadratic equation results as follows: From the above equation, it is clear that y stands for the insolation intensity and x for time. The accuracy of this approximation, expressed by the value of R 2 , was approximately 0.9991 (see Table 1 in Section 4.1). The following example will show what the response of the system will be at time x = 3 (absolute time corresponds to the high noon, i.e., 12:00). At the beginning, it should be noted that the original curve according to Equation (8) showed the energy effect S = −0.002093936, and for the set x = 3, the value y = 0.06863949. It was thus a division of the time axis into two elements. The first element E 1 was defined on the X-axis from 0 to 3, the second element E 2 from 3 to 6. On the edges of the task defined in this way, the boundary conditions are firmly defined-points A and B, and inside the task, there is a single internal point C 1 -its position in time is given by setting x = 3. We are therefore looking for optimal values of radiation intensity at point C 1 :C 1 (y).
For the assumed linear approximation within each of the two elements, the following equations were needed: The coefficients e 1 , f 1 , e 2 , f 2 was calculated as follows: The situation is shown in Figure 6. Only the value of C 1 (y) was unknown in the previous set-we knew the values of all other quantities (boundary conditions were defined, and we had chosen the time for C 1 (x) equal to 3. We can experiment with the value of C 1 (y). We evaluated the energy effect for the waveforms according to Equation (9)-see Section 2. For example, for the value C 1 (y) = 0.05, we obtained the value of the energy effect S = −0.0014546. The whole dependency is provided in Table 2.   Table 2 shows that for C 1 (y) values approaching gradually from 0.03 to 0.06, the value of S continued to decrease, and after a further increase in C 1 (y), the value of S began to increase. The solution sought is, therefore, the global minimum of the function, which was close to C 1 (y) = 0.07 (see Figure 7). We implemented a simple solution using the MS Excel Solver add-on. The Solver found the value C 1 (y) = 0.0686395 (marked with an asterisk in Table 2). For this result, the minimum value of the energy effect was S min = −0.0015704. Two remarks should be made here. Firstly: of course, it is not necessary to calculate all the values in Table 2-these are listed only to illustrate the principle of the method. In our calculations, we only used the Solver that immediately delivered the required result.
Secondly: we chose to use the Solver as an add-on to MS Excel with respect to its universal availability and ease of use. It is possible to use any other commonly available software and obtain the same results. Alternatively, it is possible to program your own application to search for the minimum of the function (or more correctly: to search for the minimum of the functional, since the function defined by us is composed of nested functions). Graphical representation of the data from Table 2 is provided in Figure 8.

The Effect of the Number of Elements on the Calculation Accuracy
To demonstrate the accuracy of the calculation, we used the knowledge of the actual course, or more precisely, the value of the energy effect S for the actual course (S = −0.002093936), which acquires its global minimum for this course. The influence of the number of elements on the accuracy of the calculation presented by us took into account the deviation of the value of the energy effect for a specific solution from the global minimum valid for the actual course. It follows from the above procedure that the accuracy of the energy effect calculation was approximately 75% for the two elements (−0.0015704/−0.002093936 × 100%). The logical solution was to increase the number of elements (i.e., divide the timeline into multiple intervals). The results are presented in Table 3 and Figure 9.  The solution procedure for a larger number of elements is described in Section 4.3.

Discussion on Results and Challenges
We added basic information about the speed of calculation when using analytical formulas and when using our method. We performed a series of calculations in MATLAB 2020a installed on a machine with Intel ® Core™ i7-8550U CPU @ 1.80 GHz (8 CPUs), 2.0 GHz, and 16,384 MB of RAM. We made a calculation for the equivalent of analytical Equation (5) and compared it with the calculation of our proposed linear system of Equation (4). For a better assessment of the calculation properties (calculation time), we always performed 20 calculation cycles, from which the average calculation time for one cycle was subsequently derived. When calculating the analytical equation, the time of one cycle was equal to 0.449 ms. The calculation of one cycle of a linear system took 0.095 ms, which is an almost 80% savings in computational time. Here comes another question concerning the number of points needed for the reconstruction of the solar curve. Our proposed numerical solution requires five computational cycles (for the required accuracy of 97.22% achieved by dividing the investigated area into six elements, i.e., five nodes-see Table 3 in the article; the dependence of the accuracy of the solution on the number of used nodes/elements can be seen in Figure 9; the decision on the satisfactory precision is up to the designer himself). These five calculation cycles took 5 × 0.095 = 0.475 ms. For five calculations of the analytical equation, the system needed 5 × 0.449 = 2.245 ms. An equation based on the sine and cosine functions was used, for example, also by [40].
Nevertheless, there also exists an approach based on a simplified exponential equation that uses only one trigonometric function [41]. For such a case, the calculation time for one cycle was up to 5.769 ms. So, compared to this approach, the calculation we propose saves up to 98% of the computational time.
Our proposed procedure was based on the principles of the MEE numerical method. In addition to the previous results presented in this article, it is possible to extend the method, for example, by adding one of the following functions.

Influence of Adaptive Meshing
The research presented by us was based on an even distribution of the X-axis (time axis). The MEE is, in principle, independent of the uniformity of this division. It is, therefore, possible to proceed to adaptive meshing without any further modifications. This can be done on the basis of two facts: the decision of the designer or the evaluation of the weight function (energy effect function S). The first option assumes the knowledge of the designer in the given field, on the basis of which he decides to adjust the density of the network. The challenge consists in the second option, where it is necessary to build a SW solution for evaluating the S effect and the subsequent change in network density.

Higher-Order Elements
We used linear interpolations within an element, i.e., first-order elements. In principle, it is possible to insert higher-order elements (e.g., polynomial) into the system, but also elements based on other functions (sine or cosine, etc.). We believe that, unlike solutions for electromagnetic and other fields, such an extension lacks its significance. By increasing the order of elements, we get into a situation where the computational complexity of the numerical method will approach the complexity of the analytical solution. In the Introduction section, we gave an overview of authors who performed further processing of insolation curves in various ways. This was our motivation: to provide authors of similar focus with a new tool that will make the required calculations easier. The aim of the article was not to offer a ready-made "out-of-the-box" solution but rather to define a calculation procedure that will be more efficient than previous solutions.

Nonlinear Weighting Functions
Further extensions in the form of nonlinear weighting functions-it aims at a possible reduction in the number of elements while maintaining accuracy, but at the expense of increased calculation demands. The same applies here, as we stated, for higher-order elements.

Conclusions
The novelty of the presented approach consists mainly in offering new tools for numerical reconstruction of insolation curves. We consider this a useful aid for other authors working in the field of solar radiation. An overview of the authors' work, where such a tool could contribute to the efficiency of their solution, is given in the Introduction section.
Another significant contribution/novelty consists in the design of solution procedures according to Section 4.3 (part B). Based on the MEE behavioral analysis, we demonstrated that, unlike the FEM or BEM methods, which offer a solution that meets weighting criteria on average across the observed area, our method tends to find the original value in each defined (and also in a single) node.
Among the novelties, we also include the analysis of the course of the energy effect in a given specific insolation environment (Figure 8), from which it follows that there is a global minimum of the solved task. This proves the previous statement that MEE tends to reach the original value in each defined (even in a single) node.
We have determined the basic possible and applicable criteria for the "accuracy" of numerical reconstruction of insolation curves-this is the energy/weight function S.
As another basic novelty, we presented an analysis of the course of the energy effect S depending on the number of elements (Figure 9). A perfect match occurred for a theoretically infinite number of nodes in the system but matches over 90% of the system were reached using five nodes (six elements in Table 3).
The original numerical solutions presented in this article are intended to help in the further processing of complicated functions of solar radiation. They can use relatively simple procedures to replace the original functions and allow other authors to apply the required mathematical operations to the daily course of solar radiation in a computationally simpler form.
In addition, the practical use of the presented solution is not limited only to solar radiation, as there are many natural and technological processes of similar nature that require precise determination of values in time.