Matlab and Fortran Precision Issue
4 次查看(过去 30 天)
显示 更早的评论
I am running a code in Matlab and fortran 90 but I get different results althought the codes are the same. Is this due to different precisions in the languages?
My code is posted below
XSTART = 2.0;
EPSA = 1.0;
EPSW = 80.0;
BULK_STRENGTH = 9.42629*1.0;
KAPPA = 8.486902807*BULK_STRENGTH/EPSW;
VK = sqrt(KAPPA);
EC2_KBT = (332.06364/0.5921830)*1.0;
F1 = 1.1682217268107287;
UNC1 = F1 - ((EC2_KBT*1.0)/(EPSA*XSTART));
FREE_ENERGY = (0.50)*(1.0)*(UNC1)*(0.5921830*1.0);
ORIGINAL_FE = (0.50)*(1.0^2)*(332.06364)*(0.50)* ...
(1.0/((VK*XSTART+1.0)*EPSW) - 1.0/EPSA)
abs(FREE_ENERGY-ORIGINAL_FE);
for ORIGINAL_FE, I get -82.670010385176070 in matlab but -82.670007683885615 in fortran 90 and I am not sure why. My fortran code is posted below (you still get the results I had using implicit double precision (A-H,O-Z)
PROGRAM MIB_HDM
IMPLICIT real*8 (A-H,O-Z)
REAL*8 :: EPSW,VK,XSTART
REAL*8 :: EC2_KBT,KAPPA
REAL*8 :: UNC1,BULK_STRENGTH
REAL*8 :: ORIGINAL_FE
REAL*8 :: EPSA
XSTART = 2.0
EPSA = 1.0
EPSW = 80.0
BULK_STRENGTH = 9.42629*1.0
KAPPA = 8.486902807*BULK_STRENGTH/EPSW
VK = sqrt(KAPPA)
EC2_KBT = (332.06364/0.5921830)*1.0
F1 = 1.1682217268107287
UNC1 = F1 - ((EC2_KBT*1.0)/(EPSA*XSTART))
FREE_ENERGY = (0.50)*(1.0)*(UNC1)*(0.5921830)
ORIGINAL_FE = (0.50)*(1.0**2)*(332.06364)*(0.50)* &
(1.0/((VK*XSTART+1.0)*EPSW) - 1.0/EPSA)
print *, original_fe
end program
Any explantion will be greatly appreciated. Thanks
0 个评论
采纳的回答
David Goodmanson
2021-10-15
Hello SA
It's the constants. First of all for simplicity I commented out everythng not concerned with calculation of ORIGINAL_FE. I don't have a fortran complier but I tried to simulate what is going on by replacing three constants with double(single(constant)). If you uncomment the three double(single)) lines in the code below, you get exactly the fortran result. I was surprised to get exactly that result to all decimal places, but that's what happens.
XSTART = 2.0;
EPSA = 1.0;
EPSW = 80.0;
BULK_STRENGTH = 9.42629*1.0;
% BULK_STRENGTH = double(single(BULK_STRENGTH));
c1 = 8.486902807;
% c1 = double(single(c1));
KAPPA = c1*BULK_STRENGTH/EPSW;
VK = sqrt(KAPPA);
% EC2_KBT = (332.06364/0.5921830)*1.0;
% F1 = 1.1682217268107287;
% F1 = double(single(F1))
% UNC1 = F1 - ((EC2_KBT*1.0)/(EPSA*XSTART));
% FREE_ENERGY = (0.50)*(1.0)*(UNC1)*(0.5921830*1.0);
c2 = 332.06364;
% c2 = double(single(c2));
ORIGINAL_FE = (0.50)*(1.0^2)*c2*(0.50)* ...
(1.0/((VK*XSTART+1.0)*EPSW) - 1.0/EPSA)
Incidentally it's not needed for this calculation but your fortran version there is no declaration of F1 as REAL*8.
更多回答(2 个)
John D'Errico
2021-10-14
No.
In the most common case, such a disagreement is because you copied over some numbers from Fortran to MATLAB (or the other way around) but only used a limited number of significant digits.
For example, we see computations like this: (332.06364/0.5921830)*1.0. Are those the EXACT values used in both cases? Do you know those to be EXACTLY the numbers that were used in both cases?
The reason i suggest this, is because the results that you show as different are the same, down to about the same number of significant digits. It is a very common mistake that we see.
The precision of MATLAB and FORTRAN is the same. They will both be using the same IEEE standard form for double precision variables. So assuming that those variables were declared in Fortran to be doubles, the precision will be the same. In MATLAB, the default is a double. And this brings up a second possibility. Could you have you declared something to be a single in the Fortran code? We don't know.
3 个评论
John D'Errico
2021-10-14
I don't have your code, as written in BOTH languages to see if you are telling the truth, or just what you THINK to be true. However, when I see things written like:
BULK_STRENGTH = 9.42629*1.0;
or
ORIGINAL_FE = (0.50)*(1.0^2)*(332.06364)*(0.50)* ...
it very much implies you are a novice, and therefore easily subject to the type of mistakes I have alluded to. The presence of many mutiplications by 1.0 or even 1 squared, suggests a novice to programming at work.
My expectation is that one (or more) of those numbers is not known to be exactly what you represent it as.
The fact is, this is NOT a question of precision carried, IF both languages are using double precision. Nothing you have done is even remotely close to creating a precision problem.
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Fortran with MATLAB 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!