An Operation Reduction Using Fast Computation of an Iteration-Based Simulation Method with Microsimulation-Semi-Symbolic Analysis †

This paper presents a method for shortening the computation time and reducing the number of math operations required in complex calculations for the analysis, simulation, and design of processes and systems. The method is suitable for education and engineering applications. The efficacy of the method is illustrated with a case study of a complex wireless communication system. The computer algebra system (CAS) was applied to formulate hypotheses and define the joint probability density function of a certain modulation technique. This innovative method was used to prepare microsimulation-semi-symbolic analyses to fully specify the wireless system. The development of an iteration-based simulation method that provides closed form solutions is presented. Previously, expressions were solved using time-consuming numerical methods. Students can apply this method for performance analysis and to understand data transfer processes. Engineers and researchers may use the method to gain insight into the impact of the parameters necessary to properly transmit and detect information, unlike traditional numerical methods. This research contributes to this field by improving the ability to obtain closed form solutions of the probability density function, outage probability, and considerably improves time efficiency with shortened computation time and reducing the number of calculation operations.


Introduction
In general, theoretical, experimental, and computational approaches are the basis for the study of observed phenomena. Every scientific and experimental result is expected to be placed into a function for its use, so the commercial use of products and services, and many engineering uses, emanate from a scientific approach that has been translated into an engineering approach.
Emerging developments are posing challenges in information technology [1][2][3] that include searching large databases [1,2], solving complex processes described by mathematical models, analyzing phenomena in communications in the information space, such as the transmission of wireless signals in urban environments [4][5][6], and the continuous high speed delivery of information without A large number of simulations do not guarantee that tolerances will not be exceeded. This is one of the numerous drawbacks of numerical-based tools. Our study had three goals. The first was to solve any analysis in closed form to allow further simplification and manipulation by using an iteration-based simulation method. The second goal was to develop an algorithm to quickly compute the method. Finally, we wanted to reduce the number of operations in the algorithm prior to its implementation. All phases of development and testing were observed by microsimulation-semi-symbolic analysis.
The IBSM was developed using the computer algebra system (CAS) to simplify complex algebraic expressions that offers an acceptable reduced analytic form for further manipulation or simulation as a closed form solution as previously published [12]. As integrals are present in the majority of the analyses, we approached the analyses with elementary calculating when the integrals are presented using Riemann sum. The method converts low-complexity implementation into a high-complexity structure. This approach allows implementation in the hardware environment.
The CAS performs symbolic mathematical operations and is used in the fields of mathematics and computer science. The CAS is based on algebraic calculations and manipulations performed using the same process as manual derivations. The CAS exclusively includes working with symbols, and numerical calculation is a special case for a CAS. Since symbols are used as variables, CAS deals with symbolic processing. Symbolic processing (SP) involves the development, implementation, and application of algorithms that manipulate and analyze mathematical expressions. CAS provides a deeper understanding and helps students to learn and engineers to simulate and design. The Wolfram language (WL) is the programming language suitable for CAS. WL has the ability to manipulate Entropy 2018, 20, 62 3 of 19 symbolic expressions using a method similar to traditional manual derivation [13]. The WL is characterized by high-performance computing and the generation of compact and short program codes.
The goal of IBSM is to introduce a new parameter to obtain a closed form expression. Since the iteration is a new parameter, we used a transformation to change the integral into a sum, i.e., a series. For this purpose, we used Riemann sum transformation for the features of the improper integrals. Using this method, we obtained closed form expressions that can be manipulated and simplified, with a short computation time, while reducing the operations. The results were tested and verified. The general form of the Riemann integral transformation into a series is given as follows [14,15]: By observing the integrals in the previous session, we defined two types of Riemann sums: a single sum and a double sum. We first solved the single sum, then solved the double integrals. So, given Equation (1), the Wolfram language code is shown in Figure 1, where q is the value of the iteration in the defined transformation. The goal of IBSM is to introduce a new parameter to obtain a closed form expression. Since the iteration is a new parameter, we used a transformation to change the integral into a sum, i.e., a series. For this purpose, we used Riemann sum transformation for the features of the improper integrals. Using this method, we obtained closed form expressions that can be manipulated and simplified, with a short computation time, while reducing the operations. The results were tested and verified. The general form of the Riemann integral transformation into a series is given as follows [14,15]: By observing the integrals in the previous session, we defined two types of Riemann sums: a single sum and a double sum. We first solved the single sum, then solved the double integrals. So, given Equation (1), the Wolfram language code is shown in Figure 1, where q is the value of the iteration in the defined transformation. Microsimulation mimics a complex phenomenon by describing its micro-components. Essentially, the system is left free to develop without too many constraints and simplifying assumptions [16]. However, when microsimulation is used with only symbolic content, and the particular numerical values are changed in the final stage, it is called microsimulation semi-symbolic analysis (MSSA). We observed each element of the symbolic calculation through MSSA, which provides faster and better testing and verification as well as a reduction in the operations [17,18]. MSSA also directly calculates in the first run without requiring more simulation attempts. The next step was the development of an algorithm to allow fast computation. To achieve this, we treated the expression as a series. As a reminder, a short explanation of the concept of fast computation follows. A series is said to converge slowly if a large number of members of the series need to be added to determine a sum with the required accuracy. During the addition of series members using the member-by-member technique, the process automatically occurs and is interrupted when a selected criterion for error evaluation is fulfilled. Given the ultimate summation, the absolute value of the relationship between the last member and the calculated sum is most often used. This criterion is not always reliable, especially for the addition of a trigonometric series. An error caused by an interrupted summing is always higher than estimated. Conversely, contemporary computing machines can quickly add a large number of members in the series. However, due to the limitation on the format of the records in the registers, a certain number of decimal places are eliminated, which leads to the accumulation of errors and to completely absurd results in the process of summing. Therefore, procedures exist for accelerating the convergence of a series, such as Kummer, Aitken, Cesar, and Euler. This paper presents an effective method for accelerating the convergence of a series based on Kummer's transformation.
We adhered to two theorems. The first states that if   Microsimulation mimics a complex phenomenon by describing its micro-components. Essentially, the system is left free to develop without too many constraints and simplifying assumptions [16]. However, when microsimulation is used with only symbolic content, and the particular numerical values are changed in the final stage, it is called microsimulation semi-symbolic analysis (MSSA). We observed each element of the symbolic calculation through MSSA, which provides faster and better testing and verification as well as a reduction in the operations [17,18]. MSSA also directly calculates in the first run without requiring more simulation attempts. The next step was the development of an algorithm to allow fast computation. To achieve this, we treated the expression as a series. As a reminder, a short explanation of the concept of fast computation follows. A series is said to converge slowly if a large number of members of the series need to be added to determine a sum with the required accuracy. During the addition of series members using the member-by-member technique, the process automatically occurs and is interrupted when a selected criterion for error evaluation is fulfilled. Given the ultimate summation, the absolute value of the relationship between the last member and the calculated sum is most often used. This criterion is not always reliable, especially for the addition of a trigonometric series. An error caused by an interrupted summing is always higher than estimated. Conversely, contemporary computing machines can quickly add a large number of members in the series. However, due to the limitation on the format of the records in the registers, a certain number of decimal places are eliminated, which leads to the accumulation of errors and to completely absurd results in the process of summing. Therefore, procedures exist for accelerating the convergence of a series, such as Kummer, Aitken, Cesar, and Euler. This paper presents an effective method for accelerating the convergence of a series based on Kummer's transformation.
We adhered to two theorems. The first states that if convergence or divergence occur simultaneously. Kummer's transformation, better known as Kummer's acceleration method, accelerates the convergence of many series. The method subtracts from a given convergent series ∑ a k , and another equivalent series ∑ b k , whose sum C = ∑ k≥0 b k is well known and finite. Kummer's transformation is described as: The convergence of the right hand side of Equation (2) is faster because 1 − ρ·b k /a k tends to 0 as k tends to infinity [19]. The complete procedure is shown in Figure 2. from a given convergent series  k a , and another equivalent series  k b , whose sum is well known and finite. Kummer's transformation is described as: The convergence of the right hand side of Equation (2) is faster because 1 − ρ·bk/ak tends to 0 as k tends to infinity [19]. The complete procedure is shown in Figure 2. The reduction in operations was performed by counting all math operations and functions contained in the final expressions. Wolfram language allows the performance of direct counting. Mathematical operations and functions in WL can be viewed both symbolically and as commands. Operations are recognized using the FullForm command, and the counting is performed using the StringPosition command. Since we had sums where the numbers are repeated q times, the WL code for completely counting the operations is: StringPosition[InnerOperations,{"Times","Power","Plus","Rational", "BesselI","Log","Exp"}]; The orders of Times, Plus, BesselI, Log, and Exp are functions used in close form expressions. Similarly, substituting s for ak[z,q], we obtained the number of operations in the accelerated algorithm.

Applications of the Accelerating Procedure and Operation Reduction with Microsimulation Semi-Symbolic Analysis
In this section, the operation reduction using fast computation of an iteration-based simulation method with microsimulation-semi-symbolic analysis was applied to two processing problems to illustrate the shorter computation time of the algorithm, and to demonstrate the variety of applications for which the operation may be used. A case with complex calculation is illustrated in the example with non-coherent Amplitude-Shift Keying (ASK) with shadowing, interference, and correlated noise. The second example treats second-order statistics in the SC macrodiversity system operating over Gamma shadowed Nakagami-m fading channels [20].

Non-Coherent Amplitude Shift Keying (ASK) with Shadowing, Interference, and Correlated Noise
Non-coherent ASK is a modulation scheme used to send digital information between digital equipment and it is shown on Figure 3. Similar part of the system, where real-time estimation is Operations are recognized using the FullForm command, and the counting is performed using the StringPosition command. Since we had sums where the numbers are repeated q times, the WL code for completely counting the operations is: StringPosition[InnerOperations, {"Times","Power","Plus","Rational", "BesselI", "Log","Exp"}]; The orders of Times, Plus, BesselI, Log, and Exp are functions used in close form expressions. Similarly, substituting s for ak[z,q], we obtained the number of operations in the accelerated algorithm.

Applications of the Accelerating Procedure and Operation Reduction with Microsimulation Semi-Symbolic Analysis
In this section, the operation reduction using fast computation of an iteration-based simulation method with microsimulation-semi-symbolic analysis was applied to two processing problems to illustrate the shorter computation time of the algorithm, and to demonstrate the variety of applications for which the operation may be used. A case with complex calculation is illustrated in the example with non-coherent Amplitude-Shift Keying (ASK) with shadowing, interference, and correlated noise. The second example treats second-order statistics in the SC macrodiversity system operating over Gamma shadowed Nakagami-m fading channels [20].

Non-Coherent Amplitude Shift Keying (ASK) with Shadowing, Interference, and Correlated Noise
Non-coherent ASK is a modulation scheme used to send digital information between digital equipment and it is shown on Figure 3. Similar part of the system, where real-time estimation is needed, can be found in [21]. The data is transmitted by the non-coherent system without a carrier in a binary manner. needed, can be found in [21]. The data is transmitted by the non-coherent system without a carrier in a binary manner. Shadowing with interference is one of the most common models used in wireless communications to describe the phenomenon of multiple scattering [21][22][23][24]. The basic components of the system are shown in Figure 1. Both shadowing and interference cause strong fluctuations in the amplitude of the useful signal. This occurs in urban areas and is described as a log-normal distribution. In our analysis, we performed an outage probability. Transmitting signals using two symbols were observed in the non-coherent ASK system [25,26]. The noise, as a narrow-band stochastic process, is correlated and the coefficient of correlation is denoted by R (R ≠ 1). Mathematically, the noise can be described as The receiver is sheltered, and no optical visibility exists toward the transmitter, but interference i1(t) = A1·cos(ωt) is present. If the system sends logical zero, then the signal s0(t) = a0·cos(ωt) has been sent, but if the system sends a logical unit, then the signal s1(t) = a1·cos(ωt) has been sent. The parameters a0 and a1 are the signal elements from which the code words are formed. The receiver detects information signal b0·cos(ωt) and b1·cos(ωt) with envelops z0 and z1 after passing through a transmitting channel. The bm (m = 0, 1) are the elements of the detected signals. The receiver system includes a filter and detector envelope. In the receiver input, the signal is: with envelopes z0 and z1, and phases φ0 and φ1, respectively. The general form of the condition joint probability density function is: where R is the coefficient of correlation and σ is variance. To ensure the set of expressions is solved continuously, using the polar coordinates is necessary, as follows: The next step was determining the condition joint probability density function (JPDF). Substituting Equation (5) into Equation (4), we obtained: where |J| is Jacobian. A joint probability density function has a log-normal distribution described as [22]: Shadowing with interference is one of the most common models used in wireless communications to describe the phenomenon of multiple scattering [21][22][23][24]. The basic components of the system are shown in Figure 1. Both shadowing and interference cause strong fluctuations in the amplitude of the useful signal. This occurs in urban areas and is described as a log-normal distribution. In our analysis, we performed an outage probability. Transmitting signals using two symbols were observed in the non-coherent ASK system [25,26]. The noise, as a narrow-band stochastic process, is correlated and the coefficient of correlation is denoted by R (R = 1). Mathematically, the noise can be described as The receiver is sheltered, and no optical visibility exists toward the transmitter, but interference i 1 (t) = A 1 ·cos(ωt) is present. If the system sends logical zero, then the signal s 0 (t) = a 0 ·cos(ωt) has been sent, but if the system sends a logical unit, then the signal s 1 (t) = a 1 ·cos(ωt) has been sent. The parameters a 0 and a 1 are the signal elements from which the code words are formed. The receiver detects information signal b 0 ·cos(ωt) and b 1 ·cos(ωt) with envelops z 0 and z 1 after passing through a transmitting channel. The b m (m = 0, 1) are the elements of the detected signals. The receiver system includes a filter and detector envelope. In the receiver input, the signal is: with envelopes z 0 and z 1 , and phases φ 0 and φ 1 , respectively. The general form of the condition joint probability density function is: where R is the coefficient of correlation and σ is variance. To ensure the set of expressions is solved continuously, using the polar coordinates is necessary, as follows: The next step was determining the condition joint probability density function (JPDF). Substituting Equation (5) into Equation (4), we obtained: where |J| is Jacobian. A joint probability density function has a log-normal distribution described as [22]: For i = j = 0, the code word 00 was sent; for i = 1 and j = 0, the code word 01 was sent; for i = 0 and j = 1, the code word 10 was sent; and for i = 1 and j = 1, the code word 11 was sent. So: The last expression can be transformed using a modified Bessel function [27] before derivation of the closed form expression: and applying trigonometric transformation: Equation (8) becomes: where: Using the Bessel identity I n (x) = I -n (x), it follows that: The present interference is described with the Rayleigh distribution over the probability density function (PDF) [23,24] as: Entropy 2018, 20, 62 To eliminate the interference, performing averaging is necessary for all values of interference A 1 .
The integral in Equation (15) is solved using integral: The distribution is obtained by averaging ϕ 0 and ϕ 1 for all values between −π and π.
For all code word combinations, distributions of envelopes are obtained by integrating all values between b 0 and b 1 . So, when the code word |ij| (i = 0, 1; j = 0, 1) has been sent, marked with H i H j in Equation (18), and when the same is detected in the input of the receiver marked with D i D j , the detection of the signals is described as: The outage probability is: where P(H i H j ) = P(H i ) · P(H j ) = 1 2 · 1 2 = 1 4 , i = 0, 1, and j = 0, 1. For the outage probability, Equation (19) represents closed form expression and is often not present in the closed form solution. Closed form expression represents an implicit solution that is contained in a mathematical expression [12]. A closed form solution provides a solved problem in terms of functions and mathematical operations from a given and generally-accepted set [28]. In other words, a closed form solution provides an explicit solution to an observed problem, whereas closed form expression shows an implicit or insufficient solution.
From Equation (7), the joint probability density function is shown in Figure 4.
The distribution is obtained by averaging φ0 and φ1 for all values between −π and π.
For all code word combinations, distributions of envelopes are obtained by integrating all values between b0 and b1. So, when the code word |ij| (i = 0, 1; j = 0, 1) has been sent, marked with HiHj in Equation (18), and when the same is detected in the input of the receiver marked with DiDj, the detection of the signals is described as: The outage probability is: where 4 , i = 0, 1, and j = 0, 1. For the outage probability, Equation (19) represents closed form expression and is often not present in the closed form solution. Closed form expression represents an implicit solution that is contained in a mathematical expression [12]. A closed form solution provides a solved problem in terms of functions and mathematical operations from a given and generally-accepted set [28]. In other words, a closed form solution provides an explicit solution to an observed problem, whereas closed form expression shows an implicit or insufficient solution. From Equation (7), the joint probability density function is shown in Figure 4.  Figure 5 shows the manipulating of Equation (8) and substituting into Equation (12) for the changing of coefficients for simplification.      Interference A 1 is present per Equation (14) and described by Wolfram language code in Figure 6, Interference A1 is present per Equation (14) and described by Wolfram language code in Figure 6, Averaging all A1 values is necessary, according to Equation (15). The general form of the condition joint probability density function is defined in Equation (14), and is described in Figure 7.  The closed form of PDFoutage in Figure 8 provides the next parameters: iteration q, h0, and h1 are the resolution of the iteration, z0 and z1 are envelopes, R is the coefficient of correlation, and σ is variance. This expression cannot be manually obtained by using numerical tools. The resultant closed form solution of Poutage is an expression that is ready for further processing. Accordingly, the viewpoint is an insight into the parameters and variables that participate in obtaining all the features of this case study. Drawing the characteristics is now possible, but this calculation would take too long, regardless of the chosen accuracy. On the other hand, for greater accuracy, a number of iterations is required, which is not beneficial for this form of expression. Averaging all A 1 values is necessary, according to Equation (15). The general form of the condition joint probability density function is defined in Equation (14), and is described in Figure 7. Interference A1 is present per Equation (14) and described by Wolfram language code in Figure 6, Averaging all A1 values is necessary, according to Equation (15). The general form of the condition joint probability density function is defined in Equation (14), and is described in Figure 7.  The closed form of PDFoutage in Figure 8 provides the next parameters: iteration q, h0, and h1 are the resolution of the iteration, z0 and z1 are envelopes, R is the coefficient of correlation, and σ is variance. This expression cannot be manually obtained by using numerical tools. The resultant closed form solution of Poutage is an expression that is ready for further processing. Accordingly, the viewpoint is an insight into the parameters and variables that participate in obtaining all the features of this case study. Drawing the characteristics is now possible, but this calculation would take too long, regardless of the chosen accuracy. On the other hand, for greater accuracy, a number of iterations is required, Interference A1 is present per Equation (14) and described by Wolfram language code in Figure 6, Averaging all A1 values is necessary, according to Equation (15). The general form of the condition joint probability density function is defined in Equation (14), and is described in Figure 7.  The closed form of PDFoutage in Figure 8 provides the next parameters: iteration q, h0, and h1 are the resolution of the iteration, z0 and z1 are envelopes, R is the coefficient of correlation, and σ is variance. This expression cannot be manually obtained by using numerical tools. The resultant closed form solution of Poutage is an expression that is ready for further processing. Accordingly, the viewpoint is an insight into the parameters and variables that participate in obtaining all the features of this case  The closed form of PDF outage in Figure 8 provides the next parameters: iteration q, h 0 , and h 1 are the resolution of the iteration, z 0 and z 1 are envelopes, R is the coefficient of correlation, and σ is variance. This expression cannot be manually obtained by using numerical tools. The resultant closed form solution of P outage is an expression that is ready for further processing. Accordingly, the viewpoint is an insight into the parameters and variables that participate in obtaining all the features of this case study. Drawing the characteristics is now possible, but this calculation would take too long, regardless of the chosen accuracy. On the other hand, for greater accuracy, a number of iterations is required, which is not beneficial for this form of expression.
Finally, the closed form solution of P outage is shown in Figure 9. In our case, a member ak represents a general member of the series in Poutage, from the closed form solution in Figure 9.
Convergence testing of the ak verified that: Convergence testing was performed with assumptions that 0 ≤ R < 1, σ > 0, z ≥ 0, and q ≥ 1. The selection of the auxiliary function is one of the most important aspects of the MSSA [29]. In testing many series, the authors of this paper highlighted the series that shows the best performance to accelerate convergence, meaning a shorter computation time with the optimum number of iteration. Comparative analysis of different auxiliary series can be the subject of particular surveys, and the reader(s) are encouraged to do so. Therefore, in our case, the auxiliary series is: The series converges to 2log2. To fully use Equation (2), we made a minor modification to the member bk, with respect to the convergence theorems that have been mentioned above. The new In our case, a member a k represents a general member of the series in P outage , from the closed form solution in Figure 9.
Convergence testing of the a k verified that: Convergence testing was performed with assumptions that 0 ≤ R < 1, σ > 0, z ≥ 0, and q ≥ 1. The selection of the auxiliary function is one of the most important aspects of the MSSA [29]. In testing many series, the authors of this paper highlighted the series that shows the best performance to accelerate convergence, meaning a shorter computation time with the optimum number of iteration. Comparative analysis of different auxiliary series can be the subject of particular surveys, and the reader(s) are encouraged to do so. Therefore, in our case, the auxiliary series is: The series converges to 2log2. To fully use Equation (2), we made a minor modification to the member b k , with respect to the convergence theorems that have been mentioned above. The new member becomes b k → a k + c k , so: where c k is general term in Equation (21). We obtain the general member of P outage marked as a k in Figure 10, separating it from Figure 9. Following the next step in MSSA, we derived the term ρ ( Figure 11).
where ck is general term in Equation (21). We obtain the general member of Poutage marked as ak in Figure  10, separating it from Figure 9. Following the next step in MSSA, we derived the term ρ ( Figure 11).  We checked that the value ρ tends to 1 after convergence testing. The quicker computation was performed by assuming how much iteration is required to calculate the outage probability Poutage obtained by the IBSM. Otherwise, a large number of iterations are required to calculate the closest exact values of Poutage, but the computation is time consuming. Then, the resulting Poutage equalizes with a new series obtained by the Kummer's transformation, and performs point matching for the various values of the envelopes, followed by a new reduced number of iterations. After that, the verification where ck is general term in Equation (21). We obtain the general member of Poutage marked as ak in Figure  10, separating it from Figure 9. Following the next step in MSSA, we derived the term ρ ( Figure 11).  We checked that the value ρ tends to 1 after convergence testing. The quicker computation was performed by assuming how much iteration is required to calculate the outage probability Poutage obtained by the IBSM. Otherwise, a large number of iterations are required to calculate the closest exact values of Poutage, but the computation is time consuming. Then, the resulting Poutage equalizes with a new series obtained by the Kummer's transformation, and performs point matching for the various values of the envelopes, followed by a new reduced number of iterations. After that, the verification We checked that the value ρ tends to 1 after convergence testing. The quicker computation was performed by assuming how much iteration is required to calculate the outage probability P outage obtained by the IBSM. Otherwise, a large number of iterations are required to calculate the closest exact values of P outage , but the computation is time consuming. Then, the resulting P outage equalizes with a new series obtained by the Kummer's transformation, and performs point matching for the various values of the envelopes, followed by a new reduced number of iterations. After that, the verification of the obtained results was performed by checking the relative error, which determined the degree of adjustability of the algorithm [29]. Finally, we checked the number of operations of calculations in the expression in Figure 9, and then obtained a reduced number of operations with a new decreased number of iterations.
After all symbolic derivations, we used closed form solutions to directly obtain results in the first attempt. To obtain concrete numerical results, we needed to set the initial parameters. We supposed that the closest exact value was obtained after 500 iterations by using the outage probability P outage in Figure 9, and the resolution of the iteration was h 0 = 0 and h 1 = 1. We also used z 0 = z 1 = z to simplify the analysis. The next step was calculating the new numbers of iterations that are reduced for various values of the envelope z. This was performed using the command FindRoot[s==Poutage,{q,1}]. s is a new expression obtained by Kummer's transformation in Equation (22), and P outage is a closed form solution in Figure 8. We took the range of values z = {1, 15} for a concrete case [29]. Experiments were performed for various values of the coefficient of correlation R (R = 7/10,8/10) and the variance σ (σ = 2, 3). All calculations were performed with a precision of 10 −6 . All tests were performed on a PC with: Intel ® Core™ i5-6500 CPU@ 3.2 GHz, 8 GB RAM, 64-bit Operating System, Windows 10, and Mathematica Wolfram 11.1. The reduced number of iterations are shown in Table 1.   1  25  32  21  30  2  20  27  16  25  3  19  25  14  22  4  18  24  12  21  5  17  23  11  20  6  17  23  10  20  7  18  23  9  19  8  18  23  9  19  9  20  24  9  20  10  21  25  9  19  11  23  26  10  20  12  27  27  11  21  13  29  28  12  21  14  32  30  14  22  15 36 31 16 23 In Figure 12, the changing iteration values in terms of the envelope z are shown for the accelerated algorithm. Notably, the reduced number of iterations is not the same for each envelope value. The minimum number of iterations is z = 10, where the value provides a true detection. However, the number of iterations is in range of 9 to 35 if we observe the total range of the envelope, which is a significant reduction compared to the original 500 iterations.
In Figure 12, the changing iteration values in terms of the envelope z are shown for the accelerated algorithm. Notably, the reduced number of iterations is not the same for each envelope value. The minimum number of iterations is z = 10, where the value provides a true detection. However, the number of iterations is in range of 9 to 35 if we observe the total range of the envelope, which is a significant reduction compared to the original 500 iterations. Since the absolute error is not precisely characterized by accuracy, the relative error is used as: Since the absolute error is not precisely characterized by accuracy, the relative error is used as: Relative errors do not exceed more than 10% of the value as it is shown in Figure 13. This indicates that the algorithm is quite accurate. In Figure 14, the comparative characteristics of P outage and s are shown. The accelerated algorithm s is marked as P e,approx .
Relative errors do not exceed more than 10% of the value as it is shown in Figure 13. This indicates that the algorithm is quite accurate. In Figure 14, the comparative characteristics of Poutage and s are shown. The accelerated algorithm s is marked as Pe,approx.
Relative errors do not exceed more than 10% of the value as it is shown in Figure 13. This indicates that the algorithm is quite accurate. In Figure 14, the comparative characteristics of Poutage and s are shown. The accelerated algorithm s is marked as Pe,approx.   The total calculation of formula P outage required 1193.97 s, or 19 min and 54 s, so the average time per iteration was 70.2335 s. The sped up algorithm's total calculation time for the accelerated formula was 1.25 s, so the average time per iteration was 0.0735294 s. Wolfram language code for time consumed is:  Table provides a calculation for any value of envelope z, and command Timing provides the exact time of calculation. Command Total summarizes total time per envelope. Similarly, changing the parameter Poutage with s for the accelerated algorithm in the previous WL command line provides the time consumed for fast computation. Our algorithm is accelerated as: Figure 15 shows the number of operations in terms of the number of iterations q for fast computation. The number of iterations is fixed at q = 500 for P outage because we initially assumed that this number of iterations was satisfied for the closest exact value of P outage . The number of operations for fast computation of IBSM is less than P outage . For 500 iterations, we counted 120,000 math operations for P outage . The number of math operations changes in the range of 9000 to 34,000, which is the result of variety in the number of iterations for fast computation. Poutage. The number of math operations changes in the range of 9000 to 34,000, which is the result of variety in the number of iterations for fast computation.

Second-Order Statistics in Wireless Channels
The level crossing rate (LCR) and the average duration of fade (ADF) are important second-order statistical characteristics describing the fading channel in mobile communications. These values are suitable for designing mobile radio communication systems and for analyzing their performance. In digital telecommunications, a sudden drop in the value of the received signal directly leads to a drastic increase in the probability of error. For optimizing the coding system required to correct errors, the number of times the received signal passes through the given level in time and how long, on average, the signal is below the specified level must be known. The LCR and ADF are the appropriate measures closely related to the quality of the received signal [24].
The LCR of signal Z(t), marked as NZ(z), is defined as the signal speed crossing through level z with a positive derivative at the intersection point z. The ADF, marked as TZ(z), represents the mean time for which the signal overlay is below the specified z level.
The LCR at envelope z is mathematically defined by [22]: where z is the envelope of the received signal,  z is its derivative in time, and ) , ( is the joined probability density function. The average fade duration (AFD) is determined as [22]:

Second-Order Statistics in Wireless Channels
The level crossing rate (LCR) and the average duration of fade (ADF) are important second-order statistical characteristics describing the fading channel in mobile communications. These values are suitable for designing mobile radio communication systems and for analyzing their performance. In digital telecommunications, a sudden drop in the value of the received signal directly leads to a drastic increase in the probability of error. For optimizing the coding system required to correct errors, the number of times the received signal passes through the given level in time and how long, on average, the signal is below the specified level must be known. The LCR and ADF are the appropriate measures closely related to the quality of the received signal [24].
The LCR of signal Z(t), marked as N Z (z), is defined as the signal speed crossing through level z with a positive derivative at the intersection point z. The ADF, marked as T Z (z), represents the mean time for which the signal overlay is below the specified z level. The LCR at envelope z is mathematically defined by [22]: where z is the envelope of the received signal, • z is its derivative in time, and p z • Z (z, • z) is the joined probability density function. The average fade duration (AFD) is determined as [22]: where F Z (z ≤ Z) represents the probability that the signal level Z(t) is less than the level z. Evaluation and calculation of LCR and ADF are trivial in an environment where no large reflections exists with a large number of transmission channels and shadowing, which simplifies the mathematical description of the distribution of the signal. However, in complex environments, obtaining LCR and ADF characteristics is time-consuming. An example of a complex environment is described in Stefanovic et al. [20]. In this example, the LCR and ADF expressions were obtained. Their analytical shapes are closed forms, but the complexity shows a long computation time. Thus, the LCR value is normalized by the Doppler shift frequency f d [20] through Equation (15): Ώ 0i is related to the average powers of the Gamma long-term fading distributions, and K v (x) is the modified Bessel function of the second order. Similarly, the AFD is obtained as [20] per Equation (16): As in the previous example, we defined a general term a k from Equation (27), shown in Figure 16.   Using the expression in Figure 17, we derived the term ρ that tends to 1 when q → ∞. Using the expression in Figure 17, we derived the term ρ that tends to 1 when q → ∞. In this case, Equations (27) and (28) have already been provided in advance in a closed form where the iteration parameter q is present, so applying the IBSM would be excessive. To compute the closest exact values of LCR and AFD, 100 iterations were required in Stefanovic et al. [20]. Using Kummer's transformation, both LCR and AFD were calculated in the first iteration. All computations were performed using the values of m = 1, L = 2, Ώ = 1, c = 2, and R = 1/5. An auxiliary series was used: In this case, Equations (27) and (28) have already been provided in advance in a closed form where the iteration parameter q is present, so applying the IBSM would be excessive. To compute the closest exact values of LCR and AFD, 100 iterations were required in Stefanovic et al. [20]. Using Kummer's transformation, both LCR and AFD were calculated in the first iteration. All computations were performed using the values of m = 1, L = 2, Ώ = 1, c = 2, and R = 1/5. An auxiliary series was used: The series C converges to (1/2)·(ϑ 3 (0, e −1 )-1), where ϑ a (u, x), (a = 1, . . . ,4) is the theta function, defined as [30]: Figure 18 shows the comparative characteristics of LCR and accelerated LCR. The deviation of the accelerated series is small in relation to the original series, and the relative error is shown in Figure 19, in the specified range of envelope -35 ≤ z ≤ 30.    Figure 18 shows the number of operations of LCR (NZ) in terms of the number of iterations q for fast computation. The number of iterations was fixed at q = 100 for LCRorig because we initially assumed that this number of iterations satisfied the closest exact value of LCRorig. For 100 iterations, we counted 20,200 math operations for LCRorig. The number of math operations was 1184 for LCRaccelerated calculated in the first iteration using fast computation. Using the same method, the AFD was obtained by applying Equation (22). Figure 20 shows the comparative characteristics of AFD and accelerated AFD. A small deviation in the range of −35 ≤ z ≤ −28 was observed, perceived through the relative error in Figure 21.     Figure 18 shows the number of operations of LCR (NZ) in terms of the number of iterations q for fast computation. The number of iterations was fixed at q = 100 for LCRorig because we initially assumed that this number of iterations satisfied the closest exact value of LCRorig. For 100 iterations, we counted 20,200 math operations for LCRorig. The number of math operations was 1184 for LCRaccelerated calculated in the first iteration using fast computation. Using the same method, the AFD was obtained by applying Equation (22). Figure 20 shows the comparative characteristics of AFD and accelerated AFD. A small deviation in the range of −35 ≤ z ≤ −28 was observed, perceived through the relative error in Figure 21.  Figure 18 shows the number of operations of LCR (N Z ) in terms of the number of iterations q for fast computation. The number of iterations was fixed at q = 100 for LCR orig because we initially assumed that this number of iterations satisfied the closest exact value of LCR orig . For 100 iterations, we counted 20,200 math operations for LCR orig . The number of math operations was 1184 for LCR accelerated calculated in the first iteration using fast computation. Using the same method, the AFD was obtained by applying Equation (22). Figure 20 shows the comparative characteristics of AFD and accelerated AFD. A small deviation in the range of −35 ≤ z ≤ −28 was observed, perceived through the relative error in Figure 21.
The total calculation of formula AFD orig required 19,553.1 s, or 5 h and 25 min, so, the average time per iteration was 279.33 s, or 4 min and 19.33 s. The sped up algorithm total calculation time with the accelerated formula was 1.29688 s, so, the average time per iteration was 0.0185268 s. An obvious difference in the time calculation exists because the number of sums for AFD increased in Equation (27), where we have sums for k, l, and q. In this case, our algorithm is accelerated as: For 100 iterations, we counted the 344 × 10 6 math operations for AFD orig . The number of math operations was 5619 for LCR accelerated calculated in the first iteration for fast computation. For 100 iterations, we counted the 344 × 10 6 math operations for AFDorig. The number of math operations was 5619 for LCRaccelerated calculated in the first iteration for fast computation.

Conclusions
This paper presents a new method to accelerate the computation and reduce the number of calculation operations in the iteration-based simulation method. The method was developed to simulate the systems and processes when obtaining mathematical formulas in the final closed form is not possible. Often, many phenomena show that closed form expressions and simulations are executed with numerical based tools. In these cases, the users do not have insight into the phenomena  For 100 iterations, we counted the 344 × 10 6 math operations for AFDorig. The number of math operations was 5619 for LCRaccelerated calculated in the first iteration for fast computation.

Conclusions
This paper presents a new method to accelerate the computation and reduce the number of calculation operations in the iteration-based simulation method. The method was developed to simulate the systems and processes when obtaining mathematical formulas in the final closed form is not possible. Often, many phenomena show that closed form expressions and simulations are executed with numerical based tools. In these cases, the users do not have insight into the phenomena

Conclusions
This paper presents a new method to accelerate the computation and reduce the number of calculation operations in the iteration-based simulation method. The method was developed to simulate the systems and processes when obtaining mathematical formulas in the final closed form is not possible. Often, many phenomena show that closed form expressions and simulations are executed with numerical based tools. In these cases, the users do not have insight into the phenomena that affect the flow of processes, which can lead to incorrect assumptions and results. The method provides insight into processes and systems using symbolic processing, with significant acceleration and reduction in the number of computation operations required. For symbolic derivation, the computer algebra system was used, and Kummer's transformation was used to shorten the computation time. The complete method to reduce the number of operations and shorten the computation time was illustrated in two examples. Both cases require complex and time-consuming calculations. Due to the large number of operations, the memory resources can also play a significant role in the speed of the calculation. The acceleration of the algorithm and the reduction of the number of operations significantly affected efficiency in terms of time savings and the rapid production of results. The method can be used in many fields where fast computation in one-step simulation runs is required.