• 大小: 1KB
    文件类型: .zip
    金币: 1
    下载: 0 次
    发布日期: 2021-06-03
  • 语言: Matlab
  • 标签:

资源简介

均匀长方体重力异常正演模拟Matlab代码

资源截图

代码片段和文件信息

function gg=cftmx()
%     均匀长方体重力异常     %
%质心坐标(x0y0z0)
%埋深H(m)
%长为a(m),宽为b(m),高为c(m),剩余密度p(kg/m^3)
%xyz为采样点
%重力异常gg (m/s^2)
%G=6.67e-11(m^3/kg.s^2) 万有引力常数
%数据保存在‘均匀长方体重力异常.txt’

clear
clc

%长方体模型参数%
a=2000;       %长
b=200;       %宽
c=100;       %高
x0=0;        %质心坐标9x0y0z0)
y0=0;
H=1000;       %质心埋深z0
ph=2*10^3;   %剩余密度

%采样区间%
x=(-2000:50:2000);
y=(-2000:50:2000);
z=0;

%常数%
G=6.67e-11;

%计算异常%
[x1y1]=meshgrid(xy);   %生成网线节点矩阵
r1=sqrt((x0+a/2-x1).^2+(y0+b/2-y1).^2+(H+c/2-z).^2);
r2=sqrt((x0+a/2-x1).^2+(y0+b/2-y1).^2+(H-c/2-z).^2);
r3=sqrt((x0+a/2-x1).^2+(y0-b/2-y1).^2+(H+c/2-z).^2);
r4=sqrt((x0+a/2-x1).^2+(y0-b/2-y1).^2+(H-c/2-z).^2);
r5=sqrt((x0-a/2-x1).^2+(y0+b/2-y1).^2+(H+c/2-z).^2);
r6=sqrt((x0-a/2-x1).^2+(y0+b/2-y1).^2+(H-c/2-z).^2);
r7=sqrt((x0-a/2-x1).^2+(y0-b/2-y1).^2+(H+c/2-z).^2);
r8=sqrt((x0-a/2-x1).^2+(y0-b/2-y1).^2+(H-c/2-z).^2);
g1=G*ph*((x0+a/2-x1).*log((y0+b/2-y1)+r1)+(y0+b/2-y1).*log((x0+a/2-x1)+r1)+(H+c/2-z).*atan((x0+a/2-x1).*(y0+b/2-y1))./((H+c/2-z).*r1));
g2=G*ph*((x0+a/2-x1).*log((y0+b/2-y1)+r2)+(y0+b/2-y1).*log((x0+a/2-x1)+r2)+(H-c/2-z).*atan((x0+a/2-x1).*(y0+b/2-y1))./((H-c/2-z).*r2));
g3=G*ph*((x0+a/2-x1).*log((y0-b/2-y1)+r3)+(y0-b/2-y1).*log((x0+a/2-x1)+r3)+(H+c/2-z).*atan((x0+a/2-x1).*(y0-b/2-y1))./((H+c/2-z).*r3));
g4=G*ph*((x0+a/2-x1).*log((y0-b/2-y1)+r4)+(y0-b/2-y1).*log((x0+a/2-x1)+r4)+(H-c/2-z).*atan((x0+a/2-x1).*(y0-b/2-y1))./((H-c/2-z).*r4));
g5=G*ph*((x0-a/2-x1).*log((y0+b/2-y1)+r5)+(y0+b/2-y1).*log((x0-a/2-x1)+r5)+(H+c/2-z).*atan((x0-a/2-x1).*(y0+b/2-y1))./((H+c/2-z).*r5));
g6=G*ph*((x0-a/2-x1).*log((y0+b/2-y1)+r6)+(y0+b/2-y1).*log((x0-a/2-x1)+r6)+(H-c/2-z).*atan((x0-a/2-x1).*(y0+b/2-y1))./((H-c/2-z).*r6));
g7=G*ph*((x0-a/2-x1).*log((y0-b/2-y1)+r7)+(y0-b/2-y1).*log((x0-a/2-x1)+r7)+(H+c/2-z).*atan((x0-a/2-x1).*(y0-b/2-y1))./((H+c/2-z).*r7));
g8=G*ph*((x0-a/2-x1).*log((y0-b/2-y1)+r8)+(y0-b/2-y1).*log((x0-a/2-x1)+r8)+(H-c/2-z).*atan((x0-a/2-x1).*(y0-b/2-y1))./((H-c/2-z).*r8));
gg=-(g1-g2-g3+g4-g5+g6+g7-g8)*10^5;    %单位mGal

%成图%
figure(1)
mesh(x1y1gg)%三维
xlabel(‘‘)
ylabel(‘‘)
title(‘均匀长方体重力异常‘)
figure(2)
contourf(x1y1gg)%二维
title(‘均匀长方体重力异常‘)


%数据生成文本%
%t=[x1(:)‘
%   y1(:)‘
%   gg(:)‘];
%fid=fopen(‘均匀长方体重力异常.txt‘‘wt‘); %wt以文本格式写入
%fprintf(fid‘%4.2f %4.2f %.2e\n‘t);
%fclose(fid);


 属性            大小     日期    时间   名称
----------- ---------  ---------- -----  ----
     文件        2487  2015-12-01 21:13  cftmx.m

评论

共有 条评论

相关资源