(242i) Algorithm to Deal with Very Large and Very Small Numbers to Enable Utilization of Analytical Solution of Free Radical Polymerization during Gel Effect

Algorithm to Deal
With Very Large and Very Small Numbers to Enable Utilization of Analytical
Solution of Free Radical Polymerization During Gel Effect

Dhiraj Kumar Garg^a,
Christophe A. Serra^b, Yannick Hoarau^c

^aDept. of
Chemical Engineering, Shiv Nadar University, dhiraj.garg@snu.edu.in, ^bICS
(UPR 22 CNRS) & University of Strasbourg, ca.serra@unistra.fr, ^cICUBE
(UMR 7357 CNRS) &  University of Strasbourg,
hoarau@unistra.fr

Introduction

Numerical modelling and simulation is a great tool used
extensively now-a-days. But to get correct results using it, several
precautions/guidelines have been discussed by several researchers. To keep
numerical simulations stable and to get correct results out of them is every
researcher's priority. There are several situations which may be adverse to
this requirement and everyone wants to avoid them and dealing with very large
and very small numbers (inverse of very large numbers) is one of them. From
machine point of view, a very large number would be some which is larger than
the largest number stored and processed by that machine and the machine shows
overflow for this condition. Similarly for the very small number and the
machine shows underflow for the same. The other problems that are faced because
of very large numbers are truncation error, round off error, subtractive
cancellation and problem associated with adding a large and a small number
(Chapra and Canale 2007). Exponential, power and factorial functions are some
functions which are notorious in building up very large numbers rapidly.

We faced one such problem during the application of
analytical solution (AS) of free radical polymerization (FRP) (Garg et al.
2014a) with gel model (Garg et al 2014b; 2014c) to simulate gel effect during
polymerization.

Problem Faced

The AS for FRP worked excellently before gel effect as
expected but there was some problem in getting correct results for one of the
variables while modelling and simulating gel effect at lower temperatures
(<90°C) using AS with different gel models. The variable was identified to
be µ₂, second order moment of dead polymer chain length
distribution whose AS consisted of infinite series solution for case 1 and 2 (Garg
et al. 2014a) as given in eqn. (1).

               (1)

The
important portion of this equation causing problem is given below.

                                                                 (2)

where  ,                                                                                       (3)

n represents n^th
time step and m represents the m^th term.

Considering
only those terms which were causing problem, eqn. 2 reduces to eqn. 4  

                                                              (4)

which,
for the ease of notation, can be rewritten as

                                                                           (5)

where                                                                                                (6)

Three
functions resulting in very large (or small) numbers can be noticed i.e. e^-x,
x^m, and m!. Now, from physical aspects of the problem

x, z
>>1                                                                                                                                  (7)

                                                                                                                             (8)

So x
> z                                                                                                                                   (9)

In our
previous work (Garg et al. 2014) we have shown and proven that the infinite
series present in the eqn. 5 are converging series for any value of x,
with m^th term under the limit . This provides the stopping
criteria to calculate the sum of the series and thus reducing the truncation
error by accommodating large enough number of terms. The physical reality of
the problem states that eqn. 4 should have finite small value irrespective of
how large (or small) values individual functions may or can have.

Modelling
and simulations are usually performed by computational software like Matlab and
they have limitation on maximum and minimum numbers they can work with for a
given computing machine. Just to elaborate the problem further, lets see what
is the practical upper limit of values of x and m that they can have during
computations against physical reality. The maximum value Matlab can
process is 1.7977x10⁺³⁰⁸. Generally
such software have inbuilt library functions for mathematical functions which
restricts maximum value of m = 169 for factorial function, m! =
169! = 4.2691x10⁺³⁰⁴ and lowest value criteria implies m =
745 for exponential function e ^-745 = 4.9406x10^-324.
Power function values corresponding to maximum value of x = 745 and m
= 169 respectively are x^m= 745¹⁰⁷ = 2.0928x10⁺³⁰⁷
and x^m= 66¹⁶⁹ = 3.1837x10⁺³⁰⁷.

So we have practical limitations on the upper limits of
values of x and m. It limits the number of terms m ≤
169 for x ≥ 66 before going overflow for factorial and
power function and limits the value of x ≤ 745 for m ≥
107 to avoid underflow for exponential function. This limitation will have
serious truncation errors during the computation of eqn. 4. The problem is
further compounded by the physical aspect of FRP where, for lower temperatures
< 90°C, the value of x
= C_n can
go upto the order of 10³ during gel effect. It is quite clear that
computations will simply fail for such large values for x = C_n cannot be processed. Hence
the approach to evaluate various functions like exponential, power and
factorial, separately and then using them for the final computation in the
equation is neither an efficient nor effective strategy.  Thus we need to
devise some new algorithm to calculate each term in such a way that the no term
or number overflows or underflows at any stage while keeping the truncation
error minimum and that too without increasing the computational load by
including unnecessarily large number of terms. 

Solution of the problem

A new efficient and effective algorithm is proposed to avoid
dealing with very large (and small) numbers during computation for any value of x = C_n. Following
measures were taken in that:

1)      Exponential,
power and factorial function were not evaluated independently in eqn. 4 for
each term.

2)      Exponential
term was clubbed with each term of infinite series instead of calculating it
separately.

3)      Each
term of infinite series was calculated using natural log to avoid very large
numbers.

4)      Calculation
of term of second series (RHS) is calculated from term of first series to avoid
the duplication of calculations.

                                                                                                       (10)

5)      Successive
terms were obtained from previous terms using a recursive relationship, thus again
avoiding duplication of calculations as well as avoiding the need to use power
and factorial functions.

         (11)

6)      Instead
of calculating the sum of series and then taking the difference, sum of the
series of differences of terms were calculated with the help of a logic loop.

7)      To
ensure that sufficient numbers of terms are taken for the summation to keep the
truncation error low as well as to avoid unnecessarily large number of terms
without much gain towards summation, following criteria MUST BE USED to stop
the above mentioned summation loop (step 6).

       
i.           
The
relative change in the sum of the series of differences to be less than
pre-specified tolerance value, AND

      ii.           
n^th (current) term
becomes lesser than (n-1)^th (previous) term of the first
series in eqn. 4., AND

    iii.           
n^th (current) term
becomes lesser than some pre-specified tolerance value.

The
criteria mentioned in step 7 is the most important step to make this algorithm both
efficient with running the loop only for required number of terms instead of
fixed number to keep the truncation error low, and effective to accommodate
large values of x = C_n as well as m without ever
reaching overflow or underflow.

This
algorithm has been applied successfully in the AS of FRP. Just to give an idea
as to how much large numbers this algorithm could handle, a successful
modelling and simulation was run for step change in FRP batch reactor from 70°C
to 50°C, at t = 45 min after starting the reaction. During computations at
different times, the maximum value of m and C_n reached
was 1408 and ≈1327 are way higher than 169 and 745 respectively as shown earlier.
This gives the numbers as e¹³²⁷, 1327¹⁴⁰⁸and
1408! which are clearly very large numbers and whose calculations were
easily bypassed by using above algorithm.

References

1.     
S.C.
Chapra, R.P. Canale, Numerical Methods for Engineers, (2007) Tata McGraw Hill,
New Delhi.

2.     
D.K.
Garg, C.A. Serra, Y. Hoarau, D. Parida, M. Bouquey, R. Muller, Analytical
Solution of Free Radical Polymerization: Derivation and Validation.
Macromolecules,2014, 47 (14), 4567?4586

3.     
D.K.
Garg, C.A. Serra, Y. Hoarau, D. Parida, M. Bouquey, R. Muller, Analytical
Solution of Free Radical Polymerization: Applications- Implementing Gel Effect
Using CCS Model. Macromolecules, 2014, 47 (23), 8178-8189.

4.     
D.K.
Garg, C.A. Serra, Y. Hoarau, D. Parida, M. Bouquey, R. Muller, Analytical
Solution of Free Radical Polymerization: Applications- Implementing Gel Effect
Using AK Model. Macromolecules, 2014, 47 (21), 7370?7377.