CODEWITHSUNDEEP
matlabfortranprecisiongfortran
Andrew
Andrew

Reputation: 734

Issues with numerical precision between Fortran and MATLAB

NOTE: Any solutions or work arounds to this issue need to occur in the MATLAB, it is the code where the issues arise. Making the fortran match the MATLAB is counter productive because then neither piece of code will work... I understand that the difference is because of differences in the way the fortran compiler and whatever the hell MATLAB does in place of compiling interpret single versus double precision floats but I am hoping that someone can help me come up with a solution to work around this.

I am working on debugging some code that I have translated from Fortran into MATLAB and have come across something that has me stumped. In both fortran and MATLAB I have the following line

pcnt = 0.9999*(-0.5+z2)

where z2 = 0.51482129528868548. The issue I am having is that there is a difference of 2.4594e-10 in the pcnt calculated in MATLAB and the pcnt calculated in fortran. I have confirmed that z2 is precisely the same (ie z2_matlab-z2_fortran=0) and so I am stumped. both z2 and pcnt are double precision (real*8 for fortran) in both Fortran and MATLAB so as far as I am concerned they should have exactly the same accuracy (since this is being run on the same machine).

Normally I wouldn't care about a difference this small but z2 ends up getting multiplied by a large number, which is then used to calculate an index later on and the small difference ends up being a difference of around 2 for the array index later, leading to huge errors later on in the algorithm (on the order of 1e6 for number that are only at most order 1e7).

Does anyone have any idea why this issue is happening and some way to fix it? I am performing this work on MATLAB R2011a and used the gfortran compiler to compile the fortran all on a 64 bit MacBook pro with an I5 (3rd gen I think) processor.

If someone has any suggestions please let me know because if I can't find a solution to this all of the translation of some 5000 lines of code I performed over the last two weeks would be pretty much worthless otherwise.

Also, any solutions would have to be to the MATLAB code because the Fortran code is the one that currently works.

Thanks in advance,

Andrew

Upvotes: 0

Views: 449

Answers (1)

user2271770
user2271770

Reputation:

The Fortran numeric literals are single precision unless the d modifier used, while MATLAB uses double as default numeric literal type. So maybe you should rewrite your pcnt expression like:

pcnt = 0.9999d+0 * (-0.5d+0 + z2)

Conversely, you should convert to single the numeric literals of MATLAB, in order to emulate the Fortran behaviour:

pcnt = single(0.9999) * (single(-0.5) + z2);

Later edit

In extremis, one should not rely on the different compilers' algorithms of interpreting numeric literals; but instead use the native (binary) representation of the said literals:

  • write unformatted from Fortran all the numeric literals that cause you grief (in this formula, 0.9999 is the most probable suspect, because 0.5 can be represented exactly in both FP precision), in a file (see WRITE).
  • load in MATLAB these values stored in the file (see fread()).
  • use the loaded values instead of numeric literals (which should have a suggestive name, as a good programming practice).

Upvotes: 5

Related Questions