matlab - Matlab 和 fortran 精度
问题描述
我在 Matlab 和 fortran 90 中运行代码,但尽管代码相同,但我得到不同的结果。这是由于语言的精度不同吗?我的代码发布在下面
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.1682185947500601;
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);
对于 ORIGINAL_FE,我在 matlab 中得到 -82.670010385176070,但在 fortran 90 中得到 -82.670007683885615,我不知道为什么。我的 fortran 代码发布在下面(您仍然可以得到我使用隐式双精度 (AH,OZ) 的结果
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
解决方案
非常感谢你们的帮助。我想出了我的问题的答案。下面是我的 matlab 和 fortran 代码以获得相同的结果。第一个是fortran代码,第二个是matlab代码。如果您运行这两个代码,您应该得到相同的结果(-82.670010385176070)。如果您在 matlab 中取消注释双(单)表达式,那么您必须在 fortran 中取出所有D0以获得与 matlab 相同的结果(-82.670007683885615)。Fortran 代码
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,F1
XSTART = 2.D0
EPSA = 1.D0
EPSW = 80.D0
BULK_STRENGTH = 9.42629D0
KAPPA = 8.486902807D0*BULK_STRENGTH/EPSW
VK = sqrt(KAPPA)*1.D0
EC2_KBT = (332.06364D0/0.5921830D0)
F1 = 1.1682217268107287D0
UNC1 = F1 - ((EC2_KBT*1.D0)/(EPSA*XSTART))
FREE_ENERGY = (0.5D0)*(1.D0)*(UNC1)*(0.5921830)
ORIGINAL_FE = (0.5D0)*(1.D0**2)*(332.06364D0)*(0.5D0)* &
(1.D0/((VK*XSTART+1.D0)*EPSW) - 1.D0/EPSA)
print *, original_fe
end program
Matlab代码
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)
推荐阅读
- php - 即使为空,数组也充满了空白值
- mysql - eZ平台本地安装无法连接MySQL数据库
- amazon-s3 - Cloud Storage 客户访问最佳做法
- amazon-web-services - AWS 节点的 Docker 上是否有时间同步过程?
- python - pdoc 生成的 HTML 中的方法参数格式不正确
- go - 如何在 EC2 实例的 2 个端口上运行简单的 Go 服务器?
- ios - 多个命令在存档上产生 gRPCCertificates.bundle 错误(iOS Xcode 11.1,Firebase SDK)
- c# - 无法访问 C# 中 Activator 返回的对象方法
- swift - 如何在 Swift 中转换 @Binding
- excel - 对列求和,但仅包括在同一行的多列中具有匹配条件的单元格