A Serious Round-Off Problem !

[monirg]

Hello;
The following two sample codes produce two different results !!

SUB TEST1()
' ... my code...
  SUM = 0.0
  FOR I = 1 TO 20000
         SUM = SUM + A(I)*B(I)/C(I)
  NEXT
  RESULT1 = SUM
' ...my code...
END SUB

SUB TEST2()
' ... my code...
  SUM = 0.0
  FOR I = 1 TO 20000
         SUM = SUM + A(I)*B(I)/C(I) /1000.
  NEXT
  RESULT2 = SUM*1000.
' ...my code
END SUB

Variables: Result1, Result2, SUM, A,B,C are DOUBLE data types
Values of A,B,C vary from very small to very large, in no specific order
For a typical set of analytical data, Result1 = 0.9986, while Result2 = 4.2432, a considerable difference.

The different results are most likely associated with accumilated round-off errors known to be associated with floating-point arithmetic precission.

I would very much appreciate your suggestion(s) on how to minimize this serious round-off error.
Can one use a higher precission than DOUBLE ... either in VBA or FORTRAN ??

Thank you kindly.

 

[Damon Ostrander]

Hi monirg,

My first reaction was "no way!" But after examining the code I see that there is a possibility of this being due to a pathological set of numbers in either the A, B, or C arrays (or are they functions? you didn't say). If A or B or C or the combination of A*B/C are such that the values are very close to the same but switch signs a lot such that there are the same number of + and - values, you could legitimately expect such a result. Usually this type of problem is the result of a mathematical formulation that inherently presents computational difficulties, and often can be fixed by reformulating the math while taking this into account.

Of course there is always the possibility that the problem lies somewhere else, such as the assumption that in the two cases the A, B and C arrays are exactly the same.

Regarding the availability of higher precision, the highest floating point precision available built-in is 64-bit (Double). But higher precision for worksheet functions is available via addins (see this link: http://www.informit.com/discussion/index.asp?postid=c36340aa-2058-4892-b6df-9696a2a1c86e&rl=1)
I'm not sure these apply to VBA though.

Keep Excelling.

Damon

 

[monirg]

Hi Damon;
Here's more relevant info.

SUB TEST3() 
' ... my code... 
     mySum1 = 0.0
       mySum2 = 0.0
          mySum3 = 0.0 
  FOR I = 1 TO 20000 
     mySum1 = mySum1 + A(I)*B(I)/C(I) 
       mySum2 = mySum2 + A(I)*B(I)/C(I) /1000. 
          mySum3 = mySum3 + A(I)*B(I)/C(I) *1000. 
  NEXT 
    Result1 = mySum1 
        Result2 = mySum2 *1000.
           Result3 = mySum3 /1000.  
' ...my code... 
END SUB

Variables: Result1, Result2, Result3, mySum1, mySum2, mySum3, A,B,C are DOUBLE data type.
Values in A,B,C vary from very small to very large, and could be +ve or -ve, in no specific order.
For a typical set of analytical data, the above sample code produces:
...Result1 = 0.9986
.......Result2 = 4.2432
...........Result3 = 2.5388

The specified accuracies in computing the arrays' elements A,B,C are 1.D-8, 5.D-8, 1.D-8 respectively, which are quite reasonable realizing the complexity of the computational methods.

The FOR loop in the above sample code involves multiplication, division, and addition of a wide range of a large number of small and large, +ve and -ve numbers. How the processor handles each two operands at each stage, individually and in combinations, and how the intermediate results are represented and stored (in memory), appear to be at the root of the problem here. (These specialized topics are beyond my understanding of computing!!).

Sure enough, some real numbers of "nearly" same value but opposite sign could possibly be encountered in A*B/C, but should not influence the resulting value in mySum = mySum+A*B/C, provided no unnecessary rounding takes place. The situation is considerably different from an ill-posed or ill-conditioned problem

That said, let me re-phrase the question:
Given one-dimensional arrays A, B, C, each of size N (=20000 to 30000), and values of arrays' elements vary between, say, +/- 1.E-5 and +/- 1.E+5 (i.e.; the value of A*B/C could potentially vary between -1.E-15 and +1.E+15), WHAT IS THE BEST WAY TO PERFORM THE FOLLOWING ARITHMETIC OPERATION WITH MINIMUM ROUND-OFF ERROR ???

  mySum = 0.0 
  FOR I = 1 TO 20000 
     mySum = mySum + A(I)*B(I)/C(I) 
  NEXT 
  Result = mySum

Thank you.

 

[Damon Ostrander]

Hi again Monir,

Yes, upon seeing the code it is obvious there is no way the values in A, B or C could be different. But I am wondering since you didn't show the declarations if you are aware that

Dim A(30000), B(30000), C(30000) As Double

actually declares only the C array as Double, and A and B as Variant?

They must be declared

Dim A(30000) As Double
Dim B(30000) As Double
Dim C(30000) As Double

Also, if you think the +/- thing might be the culprit, you can check this by summing the positive and negative values separately. If the positive and negative sums only differ in the 13th or 14th significant digit, this is definitely the problem. And if so, then it must be the result of the nature of the problem you are solving that results in such a balance of +/- values, in which case it might be helpful to let us know what the problem is and how A, B, and C were computed.

I would also like to suggest that instead of multiplying and dividing the mySum2 and mySum3 results by 1000, you use 1024 which represents only a shift of the exponent in binary. And 1000 may cause a bit of a problem for the numerical accuracy of division because I don't believe 1/1000 has an exact representation in binary, whereas 1000 does.

Damon

 

[monirg]

Hello;

I've successfully followed certain guidelines on how additions & subtractions should be done when rounding is a serious issue.
At least for now, I'm getting some consistency in the results, and will post the entire process once it's fully verified.

Here's the latest:

SUB TEST4() 
' ... my code... 
  Factor = 'see below
  mySum = 0.0 
  FOR I = 1 TO 30000 
     mySum = mySum + A(I)*B(I)/C(I) /(Factor) 
  NEXT 
  Result = mySum *(Factor)  
' ...my code... 
END SUB

Factor =.....1.0..........10.0........100.0........1000.0.....5000.0.....10000.0
Result = 0.85933...0.85975....0.31101....0.11709....0.11709....0.11709

Factor =.....1.0........0.100........0.010........0.0010.....0.0001
Result = 0.85933...0.85911....0.85916....0.85933....0.85933

The "make sense" result for the set of data is around the 0.11500 mark.

The remaining inconsistency, I believe, is primarily due to how the term A*B/C is formulated, realizing the wide variation in the individual values. Each element in the arrays A, B & C has a value in the range, say, +/- 1.E-6 to +/- 1.E+6. NO zeros.

The question is basically related to Arithmetic Precision in Multiplication and Division. It should not be confused with which operation is performed first and how the compiler does it, but rather how one would formulate the expression in the first place to minimize the round-off error.

For example:
>would: 1.0/C*A*B be more accurate than A*B/C ?
>how about: A/C*B, or B*A/C ?
>would assigning intermediate variable help: D = A*B, Result = D/C
>or D = 1.0/C, E = A*B, Result = D*E
>etc.

Surely, one may continue trying until a specific formulation is established for each set of data! But this trial & error approach would be cumbersome at best, and would be hardly described as scientific computing!

Your suggestion(s) regarding Arithmetic Precision in Multiplication and Division ... in VBA or FORTRAN ...would be greatly appreciated.

Regards.

 

[PaddyD]

MS's views on the subject:

http://support.microsoft.com/default.aspx?scid=kb;en-us;Q78113

http://support.microsoft.com/default.aspx?scid=kb;en-us;Q48606

& a couple of worked examples:

http://www.cpearson.com/excel/rounding.htm

 

[monirg]

Hi PaddyD;

Thank you for your reply. Will access the links shortly.

Regards.