改实验一：二度直立柱体正演程序设计实验报告(202,24-李龙坤,2016-11-18)

《重磁资料处理与解释》实验一

二度直立柱体正演程序设计

专业名称：勘查技术与工程

学生姓名：李龙坤

学生学号：201326020224

指导老师：王万银、纪新林、纪晓琳、邱世灿提交日期：2016-11-18

目录

1 基本原理 (1)

2 输入/输出数据格式设计 (1)

2.1 场源参数数据格式设计 (1)

2.2 计算点坐标数据格式设计 (2)

2.3 计算结果输出数据格式设计 (2)

2.4 参数文件数据格式设计 (2)

3 总体设计 (2)

4 测试结果 (3)

4.1 测试参数 (3)

4.2 测试结果 (4)

五结论及建议 (5)

附录：源程序代码 (5)

1基本原理

在直角坐标系o-xz 中，形体（二度直立柱体）模型如图1所示。设该二度直立柱体x 方向的坐标范围为21~ξξ，z 方向（铅垂向下为正）坐标为21~ζζ；又设该二度直立柱体剩余密度为i σ。

根据正演理论得知，其在空间任意一点),(z x 处产生的重力异常为

()()

()()()()[]

()2

1

21

12

22

22ln 2),(i i i i z x tg z z x x G d d z x z

G z x g i S

i i ζζξξ

ζξζζξξσζξζξζσδ⎭⎬⎫⎩

⎨⎧⎪⎪⎭⎫ ⎝⎛---+-+--=-+--=-⎰⎰

（1-1）

式中，G 为万有引力常数，在国际单位制中其值为()

2311-m 10676s kg ⋅⨯/.。

若有n 个二度直立柱体，其引起的重力异常符合叠加原理，由此得到

()()∑==n

i i z x g z x g 1,,δ∆（1-2）

2输入/输出数据格式设计

2.1 场源参数数据格式设计

场源参数按照一个二度直立柱体为一个记录进行设计，在数据文件中占一行。第一列为剩余密度density_source(g/cm 3)；第二列～第三列为x 坐标的起点1ξ和终点2ξ(m)；第四列～第五列为z 坐标的起点1ζ和终点2ζ(m ，向下为正)。以上各量均为实型变量，各量的意义见图1所示。例如：

0.2 -100 -50 50 200 0.3 -50 50 50 200

2.2 计算点坐标数据格式设计

计算点坐标数据格式设计为规则测线，采用一个计算点为一个记录的方式设计。第1列保存计算点x坐标x_coord（m），第z列保存计算点z坐标z_coord（m）。以上各量均为实型变量。保存在Input_file_source文件名变量中。例如：

-400 -20

-398 -20

2.3计算结果输出数据格式设计

计算结果输出数据格式与输入格式对应，设计为规则测线，采用一个计算点为一个记录的方式设计。第1列保存计算点x坐标x_coord（m），第2列保存计算点z坐标z_coord （m），第3列保存计算点计算结果field（mGal）。以上各量均为实型变量。保存在Output_file_source文件名变量中。例如：

-400 -20 0.24

-398 -200.35

2.4 参数文件数据格式设计

将以上部分量保存在一个文件中，该文件名变量为cmdfile，字符串变量，长度不超过80，全路径名。在该文件中保存的参数如下：

场源参数个数：N_SOURCE，整型数据。

计算点参数个数：N_COORD，整型数据。

场源参数文件名input_file_source，字符串变量，长度不超过80；

计算点坐标文件名input_file_coordinate，字符串变量，长度不超过80；

计算结果输出文件名output_file_field，字符串变量，长度不超过80

3总体设计

此次程序采用IPO结构设计，首先通过读取cmd文件，得到相关输入参数：输入场源参数个数、计算点参数个数、场源文件名、计算点坐标文件名、输出结果文件名；再通过场源参数文件，得到场源密度及其X、Y坐标；然后用过计算点参数文件，得到计算点X、Y坐标。然后利用重力异常计算公式得到每个计算点所对应的重力异常，然后将其输出到指定文件，利用GRAFER画图。

4测试结果

4.1 测试参数

（1）有关参数保存在cmdfile.txt文件中，如下:

3

401

source.dat

coord.txt

field.txt

（2）场源参数保存在“source.dat”中。第一列为剩余密度(g/cm3)；第二列～第三列为x坐标的起点和终点(m)；第四列～第五列为z坐标的起点和终点(m，向下为正)。

0.2 -100 -50 50 200

0.3 -50 50 50 200

0.2 50 100 50 200

（3）计算点为规则网数据保存在“coord.txt”文件中（txt格式）。形式如下：-400 -20

-398 -20

-397 -20

∙∙

∙∙

398 20

400 20

4.2 测试结果

利用上面的测试参数，通过重力异常程序计算并利用GRAFER 画图，如下：

通过上图3可知，地下密度异常体的在关于0轴呈对称分布，最大异常在0.6mgal,最小异常在0.08mgal 。

-400

400

X/M 0

0.2

0.4

0.6

0.8

1

∆g /m g a

l

∆g