瞬变电磁一维解析解MATLAB代码

% 中心回线中心点感应电动势解析解 % H=0是地面瞬变电磁法 均匀大地解析解 clear;clc; a=70/(pi^0.5); I0=1; u0=4*pi*10^-7; p1=50; ConE=1./p1; HzA=zeros(5000,1); t=logspace(-8,-1,5000); for m=1:5000 % u1=a*(u0*ConE/2/t(m))^0.5; u=a*(u0*ConE/t(m))^0.5/2.; Hzt(m)=(I0*p1/a.^3)*(3*erf(u)-2*u*(3+2*u^2)*exp(-u^2)/sqrt(pi)); %均匀大地磁场对时间的导数 Hz(m)=I0*(3/sqrt(pi)/u*exp(-u^2)+(1-3/2/u^2)*erf(u))/2/a; %均匀大地磁场 end figure(1); loglog(t,Hz,'-b','LineWidth',2) % 解析解 hold on; % tt=load('ttt.dat'); % 有限差分数值解 % hz=load('hhh.dat'); % loglog(tt,hz,'-r','LineWidth',2)