[zt] 浮点数的存储格式

本文详细介绍了IEEE 754浮点数标准,包括单精度(float)与双精度(double)的数据格式。通过实例展示了如何将浮点数转换为二进制形式存储,并探讨了浮点数在计算机中的具体实现。

http://wjbwjb143.blog.163.com/blog/static/2628831620078217539735/

 

IEEE floating point standard

上面我们说了,浮点数的小数点是不固定的,如果每个人都按照自己的爱好存储在电脑里,那不就乱套了吗?那么怎么在计算机中存储这种类型的数字呢?象这类古老的问题前人早都为我们做好了相应的规范,无规矩不成方圆吗。我们平时所说的浮点数的存储规范,就是由IEEE指定的,具体的规范文件是:IEEE Standard 754 for Binary Floating-Point Arithmetic。大家可以很容易的从网络上下载到这篇文档。

下面,偶就大致的描述一下,感兴趣的“同志”们可以阅读原文。

--------------------------------------------------------

声明:

  此文为原创,欢迎转载,转载请保留如下信息

  作者:afreez 北京-中关村

  联系方式:afreez.gan@gmail.com (欢迎与作者交流)

  初次发布时间:2006-12-19

初次发布在: http://blog.youkuaiyun.com/ganxingming/ 

不经本人同意,不得用语商业或赢利性质目的,否则,作者有权追究相关责任!

---------------------------------------------------------

 

       在c语言中,单精度(float)数据类型为32bits,具体的如下图所示:

 

整个32bits分三部分,即

       Sign:符号位,1 bit,0为正,1为负;

       Exponent(bias):指数部分,8 bits,存储格式为移码存储(后面还会说明),偏移量为127;

       Mantissa(fraction):尾数部分。

 

       对应的双精度(double)类型的格式为:

同样,64位也被分为了三部分,对照单精度,不用我说就可以理解各个部分的含义了吧?

       是不是有点迷糊了,不要怕,理论这个东西最能忽悠人了,看起来很高深,其实也就是个屁大的事,举个例子就很容易明白了。

举例说明,如3.24x103,则对应的部分为,Sign为0,3为指数部分(注意计算机里面存储的不是3,这里仅仅为了说明),3.24为尾数。我们知道,计算机“笨”的要死,只认识0和1,那么到底一个浮点数值在计算机存储介质中是如何存储的呢?

例如,我们要想偷窥浮点类型的值4.25在计算机硬盘中存储的庐山真面目,请跟我来:首先把4.25转换成二进制的表达方式,即100.01,在详细点,变成1.0001x22,好了,对号入座把。

Sign=0;

Exponent(bias)=2+127=129 (偏移量为127,就是直接加上个127了);

Mantissa=1.0001-1.0=0001(规格化后,小数点前总是整数1,全世界人都知道前面是1不是0,所以省略不写了,即尾数部分不包括整数部分;当别人问你,为什么23 bit的尾数部分可以表示24位的精度,知道怎么回答了吧。 靠,什么,没有看懂,再仔细读两便就知道了)。

 

对照上面的图示,相信你已经看明白了吧?相信你的智商。为了加深认识,再来一个。如果给定你一个二进制数字串,01000000100010000000000000000000,并告诉你这是一个float类型的值,让你说出它是老几,知道怎么算了吧?如果不知道,看下面的图,我就不废话解释了。

2.2深入理解浮点存储格式

为了更深入的理解浮点数的格式。我们使用C语言来做一件事。在C语言的世界里,强制类型转换,大家应该都很熟悉了。例如:

float f=4.6;

int i;

i = (int)(f+0.5); // i=5

..

下面我们不使用强制类型转化,我们自己来计算f转换成整形应该等于几?

把主要代码帖出来,如下:

 

//取23+1位的尾数部分

int ival= ((*(int *)(&fval)) & 0x07fffff) | 0x800000;

// 提取指数部分

int exponent = 150 - (((*(int *)(&fval)) >> 23) & 0xff);

if (exponent < 0)

ival = (ival<< -exponent);

else

ival = (ival >> exponent);

// 如果小于0,则将结果取反

if ((*(int *)&fval) & 0x80000000)

ival = -ival;

class RequirementPool(BaseModel): """ Product需求明细主模型:存储需求池的基本信息 """ id: Annotated[int | None, Field(description="主键(自动生成)")] = None project_no: Annotated[str | None, Field(alias="projectNo", max_length=20, description="项目编号")] = None req_name: Annotated[str | None, Field( alias="reqName", description="需求名称")] = None equip: Annotated[str | None, Field( max_length=100, description="产品名称")] = None req_sub_date: Annotated[datetime | None, Field( alias="reqSubDate", description="需求提交日期")] = None zt_id: Annotated[int, Field(alias="ztId", description="需求编号")] = None estimated_finish_date: Annotated[datetime | None, Field(alias="estimatedFinishDate", description="预估实现时间")] = None req_dev_leader: Annotated[str | None, Field(alias="reqDevLeader", max_length=100, description="开发负责人")] = None rd_req_datas: Annotated[List[RDReqData] | None, Field(alias="rdReqDatas", description="关联的 RDReqData 列表")] = None associated_version: Annotated[str | None, Field(alias="associatedVersion", max_length=255, description="关联版本号")] = None associated_version_id: Annotated[int | None, Field(alias="associatedVersionId", max_length=255, description="关联版本id")] = None release_date: Annotated[datetime | None, Field(alias="releaseDate",description="发布日期(未发布 或 具体日期)")] = None requirement_status: Annotated[str | None, Field(alias="requirementStatus", max_length=20, description="需求状态")] = None software_project: Annotated[str | None, Field(alias="softwareProject", max_length=1050, description="所属软件项目")] = None rd_req_id: Annotated[int | None, Field( alias="rdReqid", description="研发需求ID")] = None class Config: from_attributes = True extra = 'allow' 处理这个模型 时间字段要求是 北京时间
08-26
import numpy as np import os def generate_final_fortran_subroutine(input_txt, output_for): """ 生成最终优化的Fortran子程序(包含DLOAD和UTRACLOAD) 读取多圈数据(255°-285°范围内的所有数据) """ try: print(f"开始处理文件: {input_txt}") # 读取数据:第1列(时间)和第4列(法向力Fz) print("正在读取数据文件...") data = np.loadtxt(input_txt, usecols=(0, 3), dtype=np.float32) print(f"读取到 {len(data)} 行数据") # 基本参数 R = 0.457497 # 接触斑中点所在半径(m) velocity = 400 / 3.6 # 速度(m/s) omega = velocity / R # 角速度(rad/s) # 计算每个数据点对应的角度(度) - 使用连续角度 print("计算每个数据点对应的角度(度)...") angles_deg = np.degrees(data[:, 0] * omega) # 连续角度值 # 定义关键角度范围 (255°-285°) start_angle = 255.0 end_angle = 285.0 angle_range = end_angle - start_angle print(f"提取关键角度范围({start_angle}°-{end_angle}°)内的数据...") # 计算每个点相对于255°的角度偏移 angles_offset = angles_deg - start_angle # 找到所有在255°-285°范围内的点(包括多圈) mask = (angles_offset % 360) <= angle_range # 提取关键角度和载荷值 key_angles = angles_deg[mask].astype(np.float32) key_times = data[mask, 0].astype(np.float32) key_fz = data[mask, 1].astype(np.float32) if len(key_angles) == 0: print("错误: 未找到关键角度范围内的数据") return False nPoints = len(key_angles) print(f"找到 {nPoints} 个关键角度范围内的数据点") print("\n详细数据点对比 (时间, 角度(度), 载荷(N)):") print("=" * 65) print(f"{'时间(s)':<12}{'角度(度)':<12}{'载荷(N)':<12}") print("-" * 65) # 打印每个数据点的详细对比信息 for i in range(nPoints): print(f"{key_times[i]:<12.6f}{key_angles[i]:<12.6f}{key_fz[i]:<12.6f}") print("=" * 65) # 打印关键数据预览 print("\n关键角度数据点预览:") print(f" 起始角度: {key_angles[0]:.2f}°, Fz: {key_fz[0]:.6f}N") print(f" 结束角度: {key_angles[-1]:.2f}°, Fz: {key_fz[-1]:.6}N") # 生成Fortran代码(包含两个子程序) fortran_code = f""" SUBROUTINE DLOAD(F,KSTEP,KINC,TIME,NOEL,NPT,LAYER,KSPT, 1 COORDS,JLTYP,SNAME) C INCLUDE 'ABA_PARAM.INC' C DIMENSION TIME(2),COORDS (3) CHARACTER*80 SNAME **************** 直接使用KINC作为索引 **************** * 1. 直接使用增量步号KINC作为数组索引 * 2. 无需索引计算或循环 * 3. 数组从1开始与KINC匹配 *************************************************** ! ==== 声明所有变量 ==== integer nPoints parameter (nPoints = {nPoints}) real*8 angleData(nPoints) ! 存储角度(度) real*8 FzData(nPoints) ! 法向力 real*8 Fz, Tforce, omega real*8 x0,y0,z0,check,x,y,z,xt,yt,zt real*8 R, b, a, velocity real*8 pi parameter (pi = 3.141592653589793d0) integer idx ! 直接使用KINC作为索引 logical firstCall save firstCall real*8 t real*8 current_angle_deg ! 当前数据点角度(度) real*8 mapped_angle_rad ! 当前角度的弧度值 integer lastKINC ! 用于跟踪上一个增量步 save lastKINC ! 保存lastKINC值 ! ==== 初始化数组 ==== data angleData /""" # 添加角度数据 (度) # 第一行特殊处理(不加&) for i in range(min(10, nPoints)): fortran_code += f"{key_angles[i]:.6f}" if i < min(9, nPoints - 1): # 不是行内最后一个且不是数组最后一个 fortran_code += "," elif i < nPoints - 1: # 行内最后一个但不是数组最后一个 fortran_code += ",\n &," # 行末加逗号,新行以&,开头 # 后续行 for i in range(10, nPoints): if i % 10 == 0: # 每10个元素换行 fortran_code += f"{key_angles[i]:.6f}" # 新行第一个元素 else: fortran_code += f",{key_angles[i]:.6f}" # 行内元素加逗号前缀 # 行末处理(不是数组最后一个) if (i + 1) % 10 == 0 and i < nPoints - 1: fortran_code += ",\n &," # 行末加逗号,新行以&,开头 fortran_code += " /" fortran_code += "\n data FzData /" # 添加法向力数据 # 第一行特殊处理(不加&) for i in range(min(10, nPoints)): fortran_code += f"{key_fz[i]:.6f}" if i < min(9, nPoints - 1): # 不是行内最后一个且不是数组最后一个 fortran_code += "," elif i < nPoints - 1: # 行内最后一个但不是数组最后一个 fortran_code += ",\n &," # 行末加逗号,新行以&,开头 # 后续行 for i in range(10, nPoints): if i % 10 == 0: # 每10个元素换行 fortran_code += f"{key_fz[i]:.6f}" # 新行第一个元素 else: fortran_code += f",{key_fz[i]:.6f}" # 行内元素加逗号前缀 # 行末处理(不是数组最后一个) if (i + 1) % 10 == 0 and i < nPoints - 1: fortran_code += ",\n &," # 行末加逗号,新行以&,开头 fortran_code += " /" fortran_code += "\n ! ==== 结束嵌入数据 ====\n\n" fortran_code += f""" ! ==== 初始化状态标志 ==== data firstCall /.true./ data lastKINC /0/ ! 初始化lastKINC为0 ! ==== 参数设置 ==== x0 = -457.497d-3 y0 = 87.503d-3 z0 = 0.0d0 b = 8.7036d-3 a = 5.4657d-3 R = 457.497d-3 ! ==== 蠕滑区尺寸 ==== a1 = 3.1851d-3 ! 蠕滑区y方向半轴长 b1 = 4.8659d-3 ! 蠕滑区z方向半轴长 ! ==== 蠕滑区中心 ==== xn0 = x0 yn0 = y0 zn0 = z0 - (a - a1) ! 蠕滑区中心位置 ! ==== 载荷速度 ==== velocity = 400.0d0 / 3.6d0 omega = velocity / R ! ==== 获取当前时间 ==== t = TIME(1) ! ==== 直接使用KINC作为索引 ==== idx = KINC ! ==== 检查索引是否在有效范围内 ==== if (idx < 1 .or. idx > nPoints) then write(*,*) '错误: 索引超出范围! KINC=', idx, '有效范围:1-', nPoints F = 0.0d0 return endif ! ==== 获取当前载荷点角度和载荷值 ==== current_angle_deg = angleData(idx) Fz = FzData(idx) ! ==== 输出调试信息 ==== if (firstCall) then firstCall = .false. write(*,*) '直接使用KINC索引载荷模式' write(*,*) '数据点数: ', nPoints write(*,*) '起始角度: ', angleData(1), ' Fz: ', FzData(1) write(*,*) '结束角度: ', angleData(nPoints), ' Fz: ', FzData(nPoints) endif ! ==== 输出当前信息(每个增量步只输出一次) ==== if (KINC .ne. lastKINC) then write(*,*) '增量步号: ', KINC, 1 ' 角度: ', current_angle_deg, '°' write(*,*) 'Fz = ', Fz, 'N' lastKINC = KINC endif ! ==== 计算接触应力 ==== Tforce = 3.0d0 * Fz / (2.0d0 * pi * a * b) ! ==== 检查点是否在接触斑内 ==== x = COORDS(1) y = COORDS(2) z = COORDS(3) ! 坐标转换(使用当前数据点角度 - 转换为弧度) mapped_angle_rad = current_angle_deg * pi / 180.0d0 xt = x * dcos(mapped_angle_rad) + z * dsin(mapped_angle_rad) yt = y zt = -x * dsin(mapped_angle_rad) + z * dcos(mapped_angle_rad) ! 检查点是否在椭圆内 check = ((yt - y0) / a)**2 + ((zt - z0) / b)**2 ! 施加载荷或零 if (check .le. 1.0d0 .and. dabs(xt - x0) .le. 0.1d0 * R) then F = Tforce * dsqrt(1.0d0 - check) else F = 0.0d0 endif RETURN END C C ====================================================================== C SUBROUTINE UTRACLOAD(ALPHA,T_USER,KSTEP,KINC,TIME,NOEL,NPT, 1 COORDS,DIRCOS,JLTYP,SNAME) C INCLUDE 'ABA_PARAM.INC' C DIMENSION T_USER(3), TIME(2), COORDS(3), DIRCOS(3,3) CHARACTER*80 SNAME *******************头部文件************************ ! ==== 声明所有变量 ==== integer nPoints parameter (nPoints = {nPoints}) real*8 angleData(nPoints) ! 存储角度(度) real*8 FzData(nPoints) ! 法向力 real*8 Fz, Tforce, omega real*8 x0,y0,z0,check,x,y,z,xt,yt,zt real*8 R, b, a, velocity real*8 pi parameter (pi = 3.141592653589793d0) integer idx ! 直接使用KINC作为索引 logical firstCall save firstCall real*8 t real*8 current_angle_deg ! 当前数据点角度(度) real*8 mapped_angle_rad ! 当前角度的弧度值 integer lastKINC ! 用于跟踪上一个增量步 save lastKINC ! 保存lastKINC值 real*8 alpha1, E, v, fv, m, n, xp, yp, zp, check1 real*8 xn0, yn0, zn0, xck, F_hua, temp1, temp2 real*8 a1, b1 ! 蠕滑区尺寸 ! ==== 初始化数组 ==== data angleData /""" # 添加角度数据 (度) # 第一行特殊处理(不加&) for i in range(min(10, nPoints)): fortran_code += f"{key_angles[i]:.6f}" if i < min(9, nPoints - 1): # 不是行内最后一个且不是数组最后一个 fortran_code += "," elif i < nPoints - 1: # 行内最后一个但不是数组最后一个 fortran_code += ",\n &," # 行末加逗号,新行以&,开头 # 后续行 for i in range(10, nPoints): if i % 10 == 0: # 每10个元素换行 fortran_code += f"{key_angles[i]:.6f}" # 新行第一个元素 else: fortran_code += f",{key_angles[i]:.6f}" # 行内元素加逗号前缀 # 行末处理(不是数组最后一个) if (i + 1) % 10 == 0 and i < nPoints - 1: fortran_code += ",\n &," # 行末加逗号,新行以&,开头 fortran_code += " /" fortran_code += "\n data FzData /" # 添加法向力数据 # 第一行特殊处理(不加&) for i in range(min(10, nPoints)): fortran_code += f"{key_fz[i]:.6f}" if i < min(9, nPoints - 1): # 不是行内最后一个且不是数组最后一个 fortran_code += "," elif i < nPoints - 1: # 行内最后一个但不是数组最后一个 fortran_code += ",\n &," # 行末加逗号,新行以&,开头 # 后续行 for i in range(10, nPoints): if i % 10 == 0: # 每10个元素换行 fortran_code += f"{key_fz[i]:.6f}" # 新行第一个元素 else: fortran_code += f",{key_fz[i]:.6f}" # 行内元素加逗号前缀 # 行末处理(不是数组最后一个) if (i + 1) % 10 == 0 and i < nPoints - 1: fortran_code += ",\n &," # 行末加逗号,新行以&,开头 fortran_code += " /" fortran_code += "\n ! ==== 结束嵌入数据 ====\n\n" fortran_code += f""" ! ==== 初始化状态标志 ==== data firstCall /.true./ data lastKINC /0/ ! 初始化lastKINC为0 ! ==== 参数设置 ==== x0 = -457.497d-3 y0 = 87.503d-3 z0 = 0.0d0 b = 8.7036d-3 a = 5.4657d-3 R = 457.497d-3 ! ==== 蠕滑区尺寸 ==== a1 = 3.1851d-3 ! 蠕滑区y方向半轴长 b1 = 4.8659d-3 ! 蠕滑区z方向半轴长 ! ==== 蠕滑区中心 ==== xn0 = x0 yn0 = y0 zn0 = z0 - (a - a1) ! 蠕滑区中心位置 ! ==== 载荷速度 ==== velocity = 400.0d0 / 3.6d0 omega = velocity / R ! ==== 获取当前时间 ==== t = TIME(1) ! ==== 直接使用KINC作为索引 ==== idx = KINC ! ==== 检查索引是否在有效范围内 ==== if (idx < 1 .or. idx > nPoints) then write(*,*) '错误: 索引超出范围! KINC=', idx, '有效范围:1-', nPoints ALPHA = 0.0d0 return endif ! ==== 获取当前载荷点角度和载荷值 ==== current_angle_deg = angleData(idx) Fz = FzData(idx) ! ==== 输出调试信息 ==== if (firstCall) then firstCall = .false. write(*,*) 'UTRACLOAD: 直接使用KINC索引载荷模式' write(*,*) '数据点数: ', nPoints write(*,*) '起始角度: ', angleData(1), ' Fz: ', FzData(1) write(*,*) '结束角度: ', angleData(nPoints), ' Fz: ', FzData(nPoints) endif ! ==== 输出当前信息(每个增量步只输出一次) ==== if (KINC .ne. lastKINC) then write(*,*) 'UTRACLOAD: 增量步号: ', KINC, 1 ' 角度: ', current_angle_deg, '°' write(*,*) 'Fz = ', Fz, 'N' lastKINC = KINC endif ! ==== 计算接触应力 ==== Tforce = 3.0d0 * Fz / (2.0d0 * pi * a * b) ! ==== 蠕滑区参数 ==== E = 210e9 ! 弹性模量 v = 0.29 ! 泊松比 fv = 0.4 ! 摩擦系数 ! ==== 获取当前点坐标 ==== x = COORDS(1) y = COORDS(2) z = COORDS(3) ! ==== 坐标转换 ==== mapped_angle_rad = current_angle_deg * pi / 180.0d0 alpha1 = mapped_angle_rad ! 使用当前角度作为alpha1 xt = x * dcos(alpha1) + z * dsin(alpha1) yt = y zt = -x * dsin(alpha1) + z * dcos(alpha1) ! ==== 检查点是否在接触斑内 ==== check = ((yt - y0) / a)**2 + ((zt - z0) / b)**2 check1 = ((yt - yn0) / a1)**2 + ((zt - zn0) / b1)**2 xck = dabs(xt - x0) ! ==== 计算切向力大小 ==== F_hua = fv * Tforce * dsqrt(1.0d0 - check) ! 滑动区切向力 ! ==== 初始化切向力 ==== ALPHA = 0.0d0 ! ==== 检查点位置并计算切向力 ==== if (xck .le. R / 9.0d0) then ! 滑动区切向力 if (check .le. 1.0d0 .and. check1 .gt. 1.0d0) then ALPHA = F_hua endif ! 蠕滑区切向力(衰减) if (check1 .le. 1.0d0) then temp1 = ((xt - a + a1) / a1)**2 + (yt / b1)**2 temp2 = a1 / a * dsqrt(1.0d0 - temp1) ALPHA = F_hua - temp2 endif endif ! ==== 设定切向力方向 ==== T_USER(1) = dsin(alpha1) T_USER(2) = 0.0d0 T_USER(3) = -dcos(alpha1) RETURN END """ # 写入Fortran文件 with open(output_for, 'w') as f: f.write(fortran_code) print(f"\n成功生成包含DLOAD和UTRACLOAD的Fortran子程序: {output_for}") print(f"关键角度数据点数: {nPoints}") print(f"角度范围: {min(key_angles):.2f}° 到 {max(key_angles):.2f}°") print(f"法向力范围: {min(key_fz):.6f}N 到 {max(key_fz):.6f}N") return True except Exception as e: print(f"\n错误: {str(e)}") import traceback traceback.print_exc() return False if __name__ == "__main__": # 输入文件路径 input_txt = r"D:\Desktop\FEIWENTAI\波磨结果\output_V400_u-0.05已修正_W0.08\Normalforce4-Longitudinalforce1-creepage右轮.txt" # 输出文件路径 - 修改为 DONG_all1.for output_for = r"D:\Desktop\DONG_all1.for" print("=" * 50) print("数据对比工具 + FORTRAN子程序生成器") print("=" * 50) print(f"输入文件: {input_txt}") print(f"输出文件: {output_for}") print("=" * 50) # 执行转换 success = generate_final_fortran_subroutine(input_txt, output_for) if success: print("\n" + "=" * 50) print("数据对比完成!FORTRAN代码已成功生成") print("=" * 50) else: print("\n" + "=" * 50) print("处理失败,请检查错误信息") print("=" * 50) 这是生成fortran代码的python代码,针对你之前检查出来的错误对代码进行修改
07-30
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值