AU2021101629A4

AU2021101629A4 - Boundary truncation layer method and apparatus for three-dimensional forward modelling of low-frequency

Info

Publication number: AU2021101629A4
Application number: AU2021101629A
Authority: AU
Inventors: Da LEI; Juntao LU; Hao Ren; Liangyong YANG
Original assignee: Institute of Geology and Geophysics of CAS
Current assignee: Institute of Geology and Geophysics of CAS
Priority date: 2021-02-25
Filing date: 2021-03-30
Publication date: 2021-05-20
Anticipated expiration: 2029-03-30
Also published as: CN113505516A; AU2021204490A1; CN113505516B

Abstract

The patent discloses a method and an apparatus for low-frequency three-dimensional magnetotelluric forward modelling, the method including causing a low-frequency electromagnetic wave to be incident from a target region into an anisotropic boundary layer without reflection and determining a tuning parameter expression of the anisotropic BTL; performing gradual processing on a tuning parameter of the BTL to obtain an absorption parameter of a transmitted wave, and determining an optimal tuning parameter value of the BTL of a finite thickness; and introducing the tuning parameter of the BTL of a finite thickness into an electromagnetic field Maxwell equation according to a method adopted in low-frequency magnetotelluric forward modelling to obtain parameter information of low-frequency magnetotelluric forward modelling. With the method and an apparatus, a target region of a low-frequency electromagnetic method such as the magnetotelluric method can be effectively truncated. Drawings 110 Causing a low-frequency electromagnetic wave to be incident from a target region into an anisotropic boundary layer without reflection, and determining a tuning parameter expression of the anisotropic boundary truncation layer 120 Performing gradual processing on a tuning parameter of the boundary truncation layer to obtain an absorption parameter of a transmitted wave, and determining an optimal tuning parameter 130 value of the boundary truncation layer of a finite thickness Introducing the tuning parameter of the boundary truncation layer of a finite thickness into an electromagnetic field Maxwell equation according to a method adopted in low-frequency magnetotelluric forward modelling to obtain parameter information of the low-frequency magnetotelluric forward modelling FIG. 1

Description

Drawings

110

Causing a low-frequency electromagnetic wave to be incident from a target region into an anisotropic boundary layer without reflection, and determining a tuning parameter expression of the anisotropic boundary truncation layer 120

Performing gradual processing on a tuning parameter of the boundary truncation layer to obtain an absorption parameter of a transmitted wave, and determining an optimal tuning parameter 130 value of the boundary truncation layer of a finite thickness

Introducing the tuning parameter of the boundary truncation layer of a finite thickness into an electromagnetic field Maxwell equation according to a method adopted in low-frequency magnetotelluric forward modelling to obtain parameter information of the low-frequency magnetotelluric forward modelling

FIG. 1

EDITORIAL NOTE 2021101629

There are 27 pages of description only.

Description

BOUNDARY TRUNCATION LAYER METHOD AND APPARATUS FOR THREE-DIMENSIONAL FORWARD MODELLING OF LOW-FREQUENCY MAGNETOTELLURIC METHOD TECHNICAL FIELD

[0001] This study relates to the field of geophysical electromagnetic detection and, more

particularly, to a boundary truncation layer (BTL) method and apparatus for three-dimensional

forward modelling of a low-frequency magnetotelluric method.

BACKGROUND

[0002] A target region of an electromagnetic wave currently needs to be truncated artificially.

Limited computing capability means that it is impossible to perform calculations on an

electromagnetic field of an infinite region. In previous simulations concerning a boundary

condition in the electromagnetic wave field a grid of a gradually increasing size was used to

extend a skin depth several times and to attenuate a quadratic field to 0, setting the Dirichlet

boundary condition or Neumann boundary condition at an outer boundary of a simulation region.

These procedures may allow for an external travelling wave without reflection and make it decay

quickly. However, this method needs to have an adequate extension. Regarding a wider frequency

range, the size of the extended grid needs to increase slowly. The number of grids that should be

extended is, therefore, large, which requires significant calculation time.

SUMMARY

[0003] The objective of this paper is to address the problems described above and to offer

partial solutions or mitigations.

[0004] This study proposes a BTL formula with respect to a diffusion field. This research

indicates that a BTL method for three-dimensional forward modelling of a low-frequency magnetotelluric method may effectively truncate a target region with respect to a low-frequency electromagnetic method such as the magnetotelluric method. The method includes

[0005] causing a low-frequency electromagnetic wave to occur from a target region into an

anisotropic BTL without reflection, and determining a tuning parameter expression of the

anisotropic BTL;

[0006] performing gradual processing on a tuning parameter of the BTL to obtain an absorption

parameter of a transmitted wave and determining an optimal tuning parameter value of the BTL of

a finite thickness; and

[0007] introducing the tuning parameter of the BTL of a finite thickness into an electromagnetic

field Maxwell equation according to a method adopted in low-frequency magnetotelluric forward

modelling. This process enables the collection of parameter information on the low-frequency

magnetotelluric forward modelling.

[0008] Determining a tuning parameter expression of the anisotropic BTL would preferably

include

[0009] setting electrical conductivity of an anisotropic medium to C' A, a dielectric constant

to &.A, and magnetic conductivity to W~A,where A is a key tuning parameter and satisfies

the formula below:

sz 0 0 A= 0 s. 0

0 0

[0010] sz]

[0011] where s.is an anisotropic principal axis parameter;

[0012] setting a time factor to e and a Max equation expression of a frequency domain

being

[00131 VxE=iopH

[0014] V x H = YE - imeE

[0015] where E is an electric field vector, H is a magnetic field vector, o is an angular

frequency, c is electrical conductivity of a target region medium, P is magnetic conductivity

of the target region medium, E is a dielectric constant of the target region medium, i= is an imaginary number, 5E denotes a diffusion current, and i(oE denotes a displacement current.

[0016] It is derived that

V x V x E+ -io -yo(2 psE =0

[0017] (_jF0

[0018] which allows 6=-io 02

[0019]

[0020] It is calculated that

6= -, 1+II a 2 2 oec)o)

[0021]

[0022] setting the tuning parameter of the BTL to si,= xyz an expression being

s, =I1+ a,1

[00231 6C,

[0024] where Ci is a key parameter of the BTL, which controls attenuation of an

electromagnetic wave in the BTL. It is deduced that

s, =1+-,271

o) 1i+ _-1_u - i - 1+ +1 C

[0025]

[0026] After determining a tuning parameter expression of the anisotropic BTL, the method

preferably also includes simplifying the tuning parameter of the BTL.

[0027] Simplifying the tuning parameter of the BTL preferably includes

[0028] at a high-frequency limit o0 ? G , taking an expression of si, =X, y, z as

s,=1

[00301 simplifying,atalow-frequencylimit 0). s ,

0)8 ',)8

as

S =1+

[0031]

[0032] After determining a tuning parameter expression of the anisotropic BTL, the method

preferably also includes

[0033] determining a transmitted wave expression of the BTL according to an electromagnetic

wave expression of the target region and the tuning parameter.

[0034] Determining a transmitted wave expression of the BTL according to an electromagnetic

wave expression of the target region and the tuning parameter preferably includes

[0035] in a diffusion equation, making a wave number equal to

[0036]1G kk

[0037] according to a phase tuning principle, making the wave number k x of an incident wave

in an x direction equal to the wave number kt- of the transmitted wave in the x direction; and in a

z direction, the wave number ki of the incident wave and the wave number ktz of the

transmitted wave to satisfy the relationship below

[00381 k, = sk,

[0039] in the target region, expressing the electromagnetic wave as

[00401 1 = eEoe*'ek'= where represents the polarisation direction of the electric field and

EO represents an incident electric field amplitude;

[0041] expressing the transmitted wave of the BTL as

[0042] E,=eE0 e

[0043] Assuming an included angle between an incident direction of the incident wave and a

decomposition surface is 0 , then a wave number expression in the z direction would be

k= -iopac cosO and derived as

[0044] E = eEeEekz=eEoek' ek'=ze C

[0045] Performing gradual processing on a tuning parameter of the BTL to obtain an absorption

parameter of a transmitted wave and determining an optimal tuning parameter value of the BTL of a finite thickness preferably includes s =1+

[0046] modifying the tuning parameter - of the BTL

S.=K+ Vi+

[0047] into a form of i ) C) to achieve a preset effect on absorption of the

incident wave;

[0048] where Ki and are two control parameters, K controls the propagation direction

of the transmitted wave and absorption of an evanescent wave, a controls an attenuation speed

of the transmitted wave at different depths in the BTL; and Ci, Ki and a work together to

jointly control attenuation of the transmitted wave;

[0049] setting the BTL of a finite thickness, making the electromagnetic wave attenuate to 0 at

the outer boundary of the BTL, and setting a Direchlet boundary condition on an outermost side of

the BTL;

[0050] performing gradual processing on the BTL parameter si, starting a gentle rise from 1 at

an interface between the target region and the BTL, accelerating a rising rate at a distance away

from the interface; and concluding settings of parameters C7 , Ki and ai in the BTL parameter

si as follows,

K =1+max exp m- -1

-dz -exp m- -1 ai = Umax/

m 1- -1 a ex

[0051]

[0052] where d is a total thickness of the BTL; x is a distance from different units of the

BTL to the interface between the target region and the BTL; m is a constant for adjusting a

change rate of the parameter; dz is a thickness of a unit of the BTL; Kmax, Ymaxand (max are,

respectively, key parameters for controlling sizes of parameters , Ki and ai; optimal values

of Kmax,rmax and 1max are °Pt O°t and O°t, which are selected by numerical experiments.

[0053] Introducing the tuning parameter of the BTL of a finite thickness into an

electromagnetic field Maxwell equation according to a method adopted in low-frequency

magnetotelluric forward modelling to obtain parameter information of low-frequency

magnetotelluric forward modelling includes

[0054] using the BTL for grid truncation and absorption with respect to a quadratic field;

[0055] taking the tuning parameter A as a constitutive parameter matrix of the medium and

introducing A into Maxwell equations,

[00561 VxE=iopA-H

[00571 V x H = YA- E - imEA E,

[0058] a double curl expression of the electric field being

Vx VxE (iop02ps)A-E=0

[0059]

[0060] in the target region, making the tuning parameter A equal to an identity matrix;

expressing, in the BTL, the tuning parameter A as

0 0 S 0 0 s 0 , s, 0 0 A= 0 s, 0 0 0 0 s, 0 s 0 0 s ' 1 0 0 S 0 0

s 0 0 SY = 0 0

0 0 5

[00611 sz j

[0062] where an electromagnetic wave primary field of an earth background model and a total

field when the background model contains an anomalous body both conform to the double curl

expression of the electric field. By subtracting a primary field formula from a total field formula, a

magnetotelluric quadratic field forward formula may be obtained as follows:

Ix E ~(VC+(2 A -~~!cAE = io~(c7-c7)A.EP

[00631 Vp Vx E,

[0064] where P is electrical conductivity distribution of an earth background medium, ( is electrical conductivity distribution after adding an anomalous body into the earth background medium, EP is a primary field of the earth background model, and a one-dimensional analytical solution is used for calculation;

[0065] solving a partial differential equation of a quadratic field forward formula by using a

finite element method; obtaining the quadratic field and, thereafter, adding the primary field

E P, to obtain a total field value of the electric field; and

[0066] calculating through the electric field and the magnetic field to obtain impedance

information, to further calculate apparent resistivity and phase.

[0067] According to a second aspect, this research provides a BTL device for three-dimensional

forward modelling of a low-frequency magnetotelluric method, including

[0068] a parameter determining module, configured to cause a low-frequency electromagnetic

wave to be incident from a target region into an anisotropic BTL without reflection, and determine

a tuning parameter expression of the anisotropic BTL;

[0069] a gradual processing module, configured to perform gradual processing on a tuning

parameter of the BTL, to obtain an absorption parameter of a transmitted wave, and determine an

optimal tuning parameter value of the BTL of a finite thickness; and

[0070] a forward modelling processing module, configured to introduce the tuning parameter of

the BTL of a finite thickness into an electromagnetic field Maxwell equation, according to a

method adopted in low-frequency magnetotelluric forward modelling, to obtain parameter

information of the low-frequency magnetotelluric forward modelling.

[0071] According to a third aspect, this research provides a computing device, including a

memory, a processor, and a computer program stored in the memory and can be run by the

processor, wherein the processor implements the above-described method when executing the

computer program.

[0072] This research adopts a BTL based on a diffusion field, which can effectively truncate the

target region of a low-frequency electromagnetic method such as the magnetotelluric method. The

BTL is applied to finite element method magnetotelluric forward modelling to solve field value

distribution in the target region and to further calculate apparent resistivity and phase information

required for forward modelling.

[0073] Specific outcomes of this research are described below in conjunction with the

accompanying drawings. The objectives, advantages and features of the study are discussed.

BRIEF DESCRIPTION OF THE DRAWINGS

[0074] Specific findings of this research are considered in detail in an exemplary rather than a

restrictive manner with reference to the accompanying drawings. The same reference signs in the

drawings indicate the same or similar components or parts. These drawings are not necessarily

drawn to scale.

[0075] FIG. 1 is a schematic flow chart of a BTL method for three-dimensional forward

modelling of a low-frequency magnetotelluric method according to the findings of this research.

[0076] FIG. 2 is a structural schematic diagram of a BTL apparatus for three-dimensional

forward modelling of a low-frequency magnetotelluric method according to the findings of this

research.

[0077] FIG. 3 is a schematic diagram of a target region and a BTL according to the findings of

this research.

[0078] FIG. 4 is a three-dimensional COMMEMI3D-1A model diagram truncated by the BTL

according to the findings of this research.

[0079] FIG. 5 is a schematic diagram of a skin depth at different frequencies according to the

findings of this research.

[0080] FIG. 6 is a comparison diagram of electric field values when a frequency is 0.1 Hz

according to the findings of this research, in which (a) is a comparison diagram of field values of a

total field in a z direction, at x=125, y--0 m; (b) is a comparison diagram of field values of a

quadratic field Exin an x direction, at x=125, z=1250 m; (c) is a comparison diagram of field

values of the quadratic field in ay direction, atx=125,z=1250 m; and (d) is a comparison diagram

of field values of the quadratic field in a z direction, at x=125, y=O m.

[0081] FIG. 7 is a comparison diagram of electric field values when a frequency is 10 Hz

values of the quadratic field in ay direction, atx=125,z=1250 m; and (d) is a comparison diagram of field values of the quadratic field in a z direction, at x=125, y=O m.

[0082] FIG. 8 is a forward modelling result diagram of apparent resistivity and phase when the

frequency is 0.1 Hz according to the findings of this research.

[0083] FIG. 9 is a forward modelling result diagram of apparent resistivity and phase when the

frequency is 10 Hz according to the findings of this research.

[0084] FIG. 10 is a structural schematic diagram of a computing device according to the

findings of this research.

[0085] FIG. 11 is a structural schematic diagram of a computer readable storage medium

according to the findings of this research.

DETAILED DESCRIPTION

[0086] FIG. 1 is a schematic flow chart of a BTL method for three-dimensional forward

The BTL method for three-dimensional forward modelling of the low-frequency magnetotelluric

method according to the findings of this research includes

[0087] S110: causing a low-frequency electromagnetic wave to be incident from a target region

into an anisotropic boundary layer without reflection, and determining a tuning parameter

expression of the anisotropic BTL;

[0088] S120: performing gradual processing on a tuning parameter of the BTL to obtain an

absorption parameter of a transmitted wave, and determining an optimal tuning parameter value of

the BTL of a finite thickness; and

[0089] S130: introducing the tuning parameter of the BTL of a finite thickness into an

magnetotelluric forward modelling.

[0090] In the findings of this research, step S110, determining a tuning parameter expression of

the anisotropic BTL includes

[0091] setting electrical conductivity of an anisotropic medium to CF' A, a dielectric constant

to &-A, and magnetic conductivity to p A, where A is a key tuning parameter and satisfies the formula below: sz 0 0 A= 0 s. 0

0 0

!

[0092] sj

[0093] where s,is an anisotropic principal axis parameter;

[0094] setting a time factor to e , a Max equation expression of a frequency domain being:

[00951 VxE=io)pH

[00961 V x H = YE - imeE:,

[0097] where E is an electric field vector, H is a magnetic field vector, o is an angular

frequency, C is electrical conductivity of a target region medium, t is magnetic conductivity

of the target region medium, 8 is a dielectric constant of the target region medium, i= is

an imaginary number, 7E denotes a diffusion current, and i(oE denotes a displacement

current.

[0098] It is derived that

VxVx E+ -io jj,= Fxx -o_ 02 psE=0

[00991

[00100] Let: 6=- -o2

[00101]

[00102] It is calculated that

6= oI 1+IIa12 i2 2 W ) o7

[001031

[00104] setting the tuning parameter of the BTL to si, ixy,zan expression being

si =I 11

[00105] 6 .

[00106] Where Gi is a key parameter of the BTL, which controls attenuation of an

electromagnetic wave in the BTL, it is deduced that s, =1+-,271 o) 1i+ _-1__7 - i - 1+ +1 I

[001071

[00108] In the findings of this research, after determining a tuning parameter expression of the

anisotropic BTL, the method further includes simplifying the tuning parameter of the BTL.

[00109] In the findings of this research, simplifying the tuning parameter of the BTL includes

. ,Y si = _

[00110] taking an expression of si, = x, as iO at a high-frequency limit, s = 1+ '1

[00111] obtaining ios , when the time factor is e";

[00112] simplifying, at a low-frequency limit(a,; , s, =1+ ori+rj-2-i-i 1+Ka+1}8 as

si =I+--r

[001131 (1- This formula may be applied to truncation of a low-frequency

electromagnetic field.

[00114] In the findings of this research, s.is an anisotropic principal axis parameter and selection

of the formula is very important. After determining a tuning parameter expression of the

anisotropic BTL, the method further includes

[00115] determining a transmitted wave expression of the BTL according to an electromagnetic

wave expression of the target region and the tuning parameter, specifically including

[00116] in a diffusion equation, making a wave number equal to

[00117] k -lj;

[00118] according to a phase tuning principle, making the wave number k x of an incident wave

z direction, the wave number kz of the incident wave and the wave number kz of the

transmitted wave to satisfy the relationship below:

[00119] k,, = szki,

[00120] where represents the polarisation direction of the electric field, and EO representsan

incident electric field amplitude. The transmitted wave of the BTL obtained is expressed as

[00121] E,=eEe' e'

[00122] Supposing an included angle between an incident direction of the incident wave and a

decomposition surface is 0 , then a wave number expression in the z direction is k~~ ~ ~ .~k 4cosjo kz= acosO toderive E,=eEoe ke =eiE ee

.

[00123] In the findings of this research, step S120, performing gradual processing on a tuning

parameter of the BTL to obtain an absorption parameter of a transmitted wave, and determining an

optimal tuning parameter value of the BTL of a finite thickness includes

[00124] modifying the tuning parameter of the BTL to improve an absorption effect of the BTL

on the low-frequency electromagnetic wave

Si =K+

[00125] into a form of to achieve a preset effect on absorption of the

incident wave, which may achieve a better absorption effect on the incident wave.

[00126] Where Ki and ai are two control parameters, K controls the propagation direction

of the transmitted wave and absorption of an evanescent wave, 0a controls an attenuation speed

of the transmitted wave at different depths in the BTL, and yi, Ki and ai work together to

jointly control attenuation of the transmitted wave.

[00127] Due to limited computing capacity, it is necessary to set up a BTL of a finite thickness.

Making the electromagnetic wave attenuate to 0 at an outer boundary of the BTL, and setting a

Direchlet boundary condition on an outermost side of the BTL;

[00128] performing gradual processing on the BTL parameter si, starting a gentle rise from 1 at

from the interface, and concluding settings of parameters C77, Ki and a in the BTL parameter

si as follows:

K =1+ x exp m 1d

a, =(max/( C- dz -exp m- -1

LI rn)Caxx mI -XkJI

[00129]

[00130] where d is a total thickness of the BTL; x is a distance from different units of the

change rate of the parameter; dz isathicknessofaunitoftheBTL; Kmax, 3maxand Omax are

respectively key parameters for controlling sizes of parameters ', Ki and ai; and optimal

values of Kmax, are Opt G° and apt, which are selected by numerical

experiments;

[00131] and where the parameter a1 in si in the three-dimensional magnetotelluric forward

modelling, conforms to requirements that a low-frequency quadratic field is stronger than a

high-frequency quadratic field.

[00132] In the findings of this research, step S130, introducing the tuning parameter of the

boundary truncation layer of a finite thickness into an electromagnetic field Maxwell equation

according to a method adopted in low-frequency magnetotelluric forward modelling, to obtain

parameter information of low-frequency magnetotelluric forward modelling includes

[00133] using the BTL for grid truncation and absorption with respect to a quadratic field;

[00134] taking the tuning parameter A as a constitutive parameter matrix of the medium, and

introducing A into Maxwell equations,

[001351 VxE=iopA-H

[001361 V x H = YA- E - ioEA-E ;

[00137] with a double curl expression of the electric field being:

A ;XvxE

[001381

[00139] making, in the target region, the tuning parameter A equal to an identity matrix;

expressing, in the BTL, the tuning parameter A as

1 0 0 S 0 0 sY s, 0 0 A= 0 s, 0 0 0 0 S, 0 s 0 0 Sr ' 1 0 0 s 0 0

S5 0 0 Sx

= SSz S

0 0 S

[00140] Sz

[00141] where an electromagnetic wave primary field of an earth background model and a total

expression of the electric field; and by subtracting a primary field formula from a total field

formula, a magnetotelluric quadratic field forward formula may be obtained as follows:

Vx E VL-(iww A 0mwktc)EE=iwoc-c~)A-E, V A Vx E,

[00142]

[00143] where P is electrical conductivity distribution of an earth background medium, 7 is

electrical conductivity distribution after adding an anomalous body into the earth background

medium, EP is a primary field of the earth background model, and a one-dimensional analytical

solution is used for calculation;

[00144] solving a partial differential equation of a quadratic field forward formula by using a

finite element method; obtaining the quadratic field and thereafter adding the primary field EE EP, to obtain a total field value of the electric field; and

[00145] calculating through the electric field and the magnetic field to obtain impedance

information, to further calculate apparent resistivity and phase.

[00146] As shown in FIG. 2, an embodiment of the present disclosure provides an apparatus for

three-dimensional forward modelling, including

[00147] a parameter determining module 100, configured to cause a low-frequency

electromagnetic wave to be incident from a target region into an anisotropic BTL without

reflection, and determine a tuning parameter expression of the anisotropic BTL;

[00148] a gradual processing module 200, configured to perform gradual processing on a tuning

optimal tuning parameter value of the BTL of a finite thickness; and

[00149] a forward modelling processing module 300, configured to introduce the tuning

parameter of the BTL of a finite thickness into an electromagnetic field Maxwell equation

parameter information of the low-frequency magnetotelluric forward modelling.

[00150] As shown in FIG. 3, on an x-z plane, the target region is on an upper side, an anisotropic

medium is on a lower side, and their interface is perpendicular to a z axis. The electrical

conductivity of the target region is C7, the dielectric constant is 8, and the magnetic conductivity is t . The electromagnetic wave is incident from the target region into the anisotropic medium on

the lower side, and the incident angle is 0. If the electrical conductivity of the anisotropic medium on the lower side is F , A, the dielectric constant is 8, A, the magnetic conductivity is

,and A satisfies a form below:

s 0 0 A= 0 s. 0

0 0 1

[00151] s ).

[00152] The electromagnetic wave may, therefore, be incident from the target region on the

upper side into the anisotropic medium on the lower side without reflection. The target region

medium on the upper side and the anisotropic medium on the lower side are impedance-tuned, so

that the anisotropic medium on the lower side may be referred to as a boundary truncation tuning

layer. Theoretically, Sz may be set arbitrarily, and is not likely to affect the electromagnetic wave

entering the BTL without reflection. However, to make the electromagnetic wave entering the

BTL attenuate quickly, in the embodiment of the present disclosure, the formula of sz needs to

be carefully defined.

[00153] First, in the embodiment of the present disclosure, it is assumed that the time factor is ,e e -i( , and the Max equation in the frequency domain may be written as

[00154] V x E = iopH (2),

[00155] V x H = aE - icEE (3),

[00156] where E is an electric field vector, H is a magnetic field vector, 0 is an angular

frequency, 7 is electrical conductivity of the medium, E is magnetic conductivity of the

medium, 8 is a dielectric constant of the medium, and = is an imaginary number, and

where 7E may be regarded as a diffusion current, and ioE may be regarded as a

displacement current. Combined with formula (3), formula (2) may be written as

V x V x E+ -io -7o(2 psE =0

[00157] (_K F ic(4).

[00158] It is defined that 6=-io 2

[00159] (5).

[00160] It may be calculated that

1 2 = (a2

2 o()o

[001611 (6).

[00162] The parameter s = x, Y, z )of the BTL according to the embodiment of the present

disclosure is defined as

si =1+a'

[001631 6C (7).

[00164] i is a key parameter of the BTL, which may control attenuation of the

electromagnetic waves in the boundary truncation tuning layer. By substituting formula (6) for

formula (7), it is obtained that

N2( s, =1+

o) 1,+ -1 - i - 1+ +1 I

[00165] (8).

[00166] C in a denominator of formula (8) refers to electrical conductivity of the medium, and

has a different meaning from Gi in a numerator. This BTL form is suitable for truncation of the

electromagnetic wave target region in a full frequency band.

[00167] At the high-frequency limit (OE ? a ), formula (8) may be simplified as

[001681 "oF (9).

[00169] When the time factor is e , formula (9) should be

s = 1+--~

[00170] koc (10).

[00171] At the low-frequency limit (o« ),formula (8) maybe simplified as

si =1+

[00172] - (11).

[00173] In magnetotellurics, since a researched frequency is relatively low, the diffusion current

is usually much larger than the displacement current, which meets the low-frequency limit

) `« 3 . Therefore, the BTL of formula (11) may be used for grid truncation in magnetotelluric

forward modelling simulation. It can be seen that the BTL formula (11) of the diffusion field

proposed by the embodiment of the present disclosure differs from the traditional BTL formula (9)

for the high-frequency electromagnetic wave. It is related to the electrical conductivity of the

target region and has different frequency dependence as compared with the high-frequency BTL.

[00174] Ina diffusion equation, a wave number is equal to

[00175] -1 07 (12).

[00176] As shown in FIG. 3, according to a phase tuning principle, the wave number k. of the

incident wave in the x direction is equal to the wave number k- of the transmitted wave in the x

direction. In the z direction, the wave number k. of the incident wave and the wave number k,

of the transmitted wave satisfy the relationship below:

[00177] ktz = szkiz (13).

[00178] In the target region, the electromagnetic wave may be expressed as

[00179] E,=eEe'xekz (14).

[00180] The transmitted wave in the BTL may, therefore, be expressed as

[001811 E, = eEe'e' (15).

[00182] Assuming an included angle between the incident direction of the incident wave and the

decomposition surface is 0 , then the wave number in the z direction may be written as k= -iOpa cos, and by substituting formula (11) into formula (15), it may be obtained that

[001831 = eEoekxekt=z = eEge**ekie (16).

[00184] Compared with the electromagnetic wave formula (14) in the target region, it can be

seen that the first two exponential terms of the electromagnetic wave formula (16) in the BTL are

the natural attenuation terms of electromagnetic waves, and the last exponential term is an

additional attenuation term of the BTL. It can be seen that the attenuation ability of the BTL is

related to the parameter yz and the incident angle 0. When the electromagnetic wave is incident

vertically, the attenuation ability is the highest. When the electromagnetic wave is incident

obliquely, the attenuation ability is relatively weak.

[00185] In consideration of boundary truncation layers in three directions of x, y and z, a more

general formula is obtained as follows:

0 0 S 0 00 s1 s, 0 0 A= 0 s, 0 0 0 0 s, 0 s 0 0 s ' 1 0 0 S 0 0

___ 0 0 S

= 0 0 S, S S 0 0 x

[001861 S (17),

[00187] where s s and s are the BTL parameters in the x, y and z directions, all of which

S may be written as formula (11). When the interface of the BTL is perpendicular to the x-axis, s

takes a form of formula (11), and s and sz in formula (17) are regarded as 1. When the

interface of the BTL is perpendicular to the y or z axis, a similar approach is adopted.

[00188] To improve the absorption ability of the parameter si of the BTL with respect to the

transmitted wave, in the embodiment of the present disclosure, formula (11) is modified as

S =K .+

[001891 (18),

[00190] where Ki and ai are two newly added parameters with no specific physical meaning.

Ki mainly controls a propagation direction of the transmitted wave and absorption of the

evanescent wave, and al mainly controls an attenuation speed of the transmitted wave at

different depths in the BTL. , i and at 1 work together to jointly control attenuation of the

transmitted wave.

[00191] Due to limited computing capacity, the embodiment of the present disclosure can only

set a BTL of a finite thickness, make the electromagnetic wave quickly attenuate to 0, and set a

Direchlet boundary condition on the outermost side of the BTL. To reduce emission caused by a

sudden change of physical properties at the interface between the target region and the BTL, the

embodiment of the present disclosure also needs to perform gradual processing on the BTL

parameter si, start a gentle rise from 1 at the interface between the target region and the BTL, and

accelerate a rising rate at a distance away from the interface. Settings of "i, Ki and at in the

BTL parameter s are concluded as follows:

K =1+x exp m- -1

S= max/( -dz -exp m.- 1

(19),

[00192] a =

[00193] where d is a total thickness of the BTL; x isadistance from different units of the

change rate of the parameter; dz is a thickness of a unit of the BTL; Kmax, Ymaxand Omax are

respectively key parameters for controlling sizes of parameters a7 , Ki andci ; and optimal

values of OP O' and a°ptof Kmax, Ymax and Omax are selected by numerical experiments.

[00194] Magnetotelluric forward modelling usually includes two algorithms: a total field method

and a quadratic field method. In the total field method, an earth background medium and an

anomalous body medium are considered at the same time to calculate the total field. In the

quadratic field method, a primary field value of the earth background medium is first calculated, then a quadratic field value generated by the anomalous body is calculated and, finally, the primary field value and the quadratic field value are added to obtain the total field value. With respect to the quadratic field, it may be regarded as an external travelling wave, so the BTL is used for grid truncation and absorption. The matrix A may be regarded as a constitutive parameter matrix of the medium, and A is introduced into Maxwell equations:

[00195] VxE=iopA-H (20),

[001961 VxH=c7A-E-ioEA-E (21).

[00197] Combined with formula (21) and formula (22), the equations are written into a double

curl form of the electric field

VxVx (opa+o02p)A- E= 0

[001981 (22).

[00199] In the target region, the matrix A is equal to an identity matrix. In the BTL, the matrix A takes the form of formula (17). Both the electromagnetic wave primary field of the earth

background model and the total field when the background model contains an anomalous body

may be written into the form of formula (23). By subtracting a primary field formula from a total

field formula, a magnetotelluric quadratic field forward formula may be obtained as follows:

VxI ~~(OV +(2L - &ICAE = io~(c7-c)A. EP

[00200] V Vx E - A (23).

[00201] P is electrical conductivity distribution of an earth background medium, and C is

electrical conductivity distribution after adding an anomalous body to the earth background

medium. EP is a primary field of the earth background model, and a one-dimensional analytical

solution is used for calculation. A partial differential equation of formula (24) is solved by using a

finite element method, and after the quadratic field Es is obtained, the primary field EP is

further added, to obtain a total field value of the electric field. By calculating through the electric

field and the magnetic field, impedance information is obtained, to further calculate apparent

resistivity and phase.

[00202] In magnetotelluric forward modelling, an air layer above the ground needs to be

considered. The resistivity of the air layer is usually greater than 106Qm, so the resistivity of the air layer is usually much greater than common earth dielectric resistivity. When the BTL is applied to magnetotelluric forward modelling, existence of the air layer is a significant obstacle.

Since the denominator in the BTL formula has electrical conductivity of a medium surrounded

usually only when the electrical conductivity of a medium surrounded is consistent or has small

disparity, a better electromagnetic wave absorption effect can be achieved. In application of a

traditional BTL, the key parameters are the dielectric constant and the magnetic conductivity;

usually, different media are not greatly different in dielectric constant and magnetic conductivity,

so influence of the difference in physical properties of the medium surrounded on the use effect of

BTL is limited. However, in magnetotelluric forward modelling, difference in electrical

conductivity between the air layer and the earth is usually multiplied thousands of times or even

greater. Thus, this brings great difficulty in application of the BTL. Fortunately, after a large

number of numerical experiments, it is found that if the electrical conductivity of the medium in

the denominator of the BTL parameters in all directions is set to the electrical conductivity of the

earth background medium, after carefully searching for the parameters OPt, Opt and Opt,

even if there is an air layer, a better grid truncation effect may also be achieved. The parameters

OPt Opt and UOpt obtained by selection may be applied to the wide frequency band of the

magnetotelluric method and the commonly used earth background medium resistivity range; and

after changing the frequency or model background medium resistivity, there is no need to find best

parameters again.

[00203] To test the grid truncation effect of the BTL in the forward modelling of the

three-dimensional magnetotelluric method, the embodiment of the present disclosure establishes a

COMMEMI3D-1A model shown in FIG. 4. The COMMEMI3D-1A model was proposed in 1997,

and forward modelling results by different researchers were collected, giving an approximate

range of the forward modelling results. Air resistivity of the model is set to 108 Q- m, and

resistivity of the earth background medium is set to 100 0- m. There is a low-resistance prism with a maximum buried depth of 250 m in a background medium, and resistivity is 0.5 Q- m. Lengths of the low-resistance prism in the x, y and z directions are respectively 1,000 m, 2,000 m,

and 2,000 m. In the embodiment of the present disclosure, a simulation region is discretised into

cube units which has lengths of 250 m in both the x and y directions. The air layer is dissected into four layers in the z direction, with each layer having a thickness of 250 m, so a thickness of the air layer is 1,000 m. In an underground medium, to calculate reliable apparent resistivity and phase information, a z-direction thickness of three layers of underground grids is set to 5 m, and a thickness of the fourth and fifth layers of underground grids is 117.5 m. Starting from the sixth layer of underground grids, grid thicknesses in the z direction are all 250 m. The number of grids of the underground medium in the z direction is 15, and the number of grids of the effective simulation region in the x and y directions is 16, so an effective observation range on the earth's surface is -2,000 m to 2,000 m in the x direction, and -2,000 m to 2,000 m in the y direction. In the embodiment of the present disclosure, all outer sides of the simulation region are each loaded with a BTL layer of six layers of grids, with each layer having a thickness of 50 m, so a total thickness of the BTL on each side is 300 m, which is much smaller than the skin depth of the low-frequency electromagnetic wave in the medium. Therefore, the total dissection number of grids in the model including the simulation region and the BTL layer region is 282,831, and the number of unknown numbers in the finite element equations is 78,039. In three-dimensional simulation, the BTL layer parameters used in the embodiment of the present disclosure are m= opt 8.57

=0.2357 and . In the embodiment of the present disclosure, seven frequency points

are calculated, respectively 0.001 Hz, 0.01 Hz, 0.1 Hz, 1 Hz, 10 Hz, 100 Hz, and 1,000 Hz.

[00204] To compare with the traditional method, the embodiment of the present disclosure

adopts a traditional grid extension method to calculate a reference solution. When the earth

background resistivity is 100 *m, the skin depths of the seven frequency points within the

frequency range of 0.001 Hz to 1,000 Hz are shown in FIG. 5. It can be seen that the skin depths

of electromagnetic waves of different frequencies in the medium vary greatly, the skin depth may

reach 159 km at 0.001 Hz, while it is only 0.16 km at 1,000 Hz. In magnetotelluric simulation, a

traditional grid truncation solution is to load a certain number of grids with increasing sizes on the

periphery of the simulation region. In consideration of the inconsistency of different skin depths at

different frequencies, at least two to three sampling points are required within a skin depth of each

frequency. Therefore, in simulation of electromagnetic methods with a wide band, such as a

magnetotelluric method, a size of the extended grid needs to be gradually increased to obtain

better accuracy. With respect to the quadratic field method, a Dirichlet boundary condition is usually set on an outermost boundary of the extended grid, that is, a field value is set to 0. Thus, the total length of the extended grid needs to be two to three times the skin depth to meet the demand. In this model, the embodiment of the present disclosure allows the number of extended grids on each side to be 16, and the size of the extended grid to gradually increase at a rate of 1.55 times from 250 m. A total length of grid extension on each side is 780 km, which is much larger than the skin depth of 159 km in the medium when the frequency is 0.001 Hz. Therefore, the traditional method, as a reference solution, is reliable and accurate. The total number of grids in the grid extension method is 484,851, and the number of unknown numbers in the finite element equations is 367,059, which is much larger than the unknown numbers in the boundary truncation tuning layer method. In terms of time calculation, the boundary truncation tuning layer method takes 13 minutes, while the traditional grid extension method takes 4 hours and 36 minutes. The traditional grid extension method requires more grids for extension, occupies more memory, and takes longer calculation time.

[00205] In the following, the embodiment of the present disclosure only shows the forward

modelling results at 0.1 Hz and 10 Hz. FIG. 6 shows comparison of field values obtained by using

the finite element method when the frequency is 0.1 Hz. FIG. 6(a) is a comparison diagram of

total field values between the BTL method and the traditional grid extension method. It can be

seen that the two are highly consistent, wherein a purple dashed line is a primary field value. FIG.

6(b), FIG. 6(c) and FIG. 6(d) are comparison diagrams of quadratic field values of the electric

fields in the x, y and z directions, in which a blue dashed line with an 'o' mark is data of the BTL

layer method, and a red solid line with a solid dot mark is data of the traditional grid extension

method. It can be seen from the diagrams that, the quadratic electric field obtained by using the

BTL method fits well with the quadratic electric field of the traditional large-distance grid

extension method, indicating that the boundary truncation tuning layer proposed by the

embodiment of the present disclosure is effective and has a truncated boundary with high

calculation accuracy, and that the boundary truncation tuning layer method may effectively save

calculation amount and calculation time. FIG. 7(a), FIG. 7(b), FIG. 7(c) and FIG. 7(d) are

comparison diagrams of field values at a frequency of 10 Hz, which have the same meaning as

FIG. 6(a), FIG. 6(b), FIG. 6(c) and FIG. 6(d), with only the frequencies being inconsistent.

[00206] FIG. 8 and FIG. 9 show, respectively, the forward modelling results of apparent resistivity and phase at 0.1 Hz and 10 Hz. FIG. 8(a) and FIG. 8(c) are, respectively, plan views of apparent resistivity and phase. It can be seen from the diagrams that there is an obvious low-resistance anomaly at a corresponding earth's surface observation position above the anomalous body. FIG. 8(b) and FIG. 8(d) are apparent resistivity and phase curves in the x direction at y=0. A green solid line is a result in an article by Ren Zhengyong (2013), and a magenta error bar is an apparent resistivity and phase range given by different researchers in literature in 1997. A blue dashed line with an 'o' mark is data of the boundary truncation tuning layer method, and a red solid line with a solid dot mark is data from the traditional grid extension method. It can be seen that the forward modelling results of the boundary truncation tuning layer method are highly consistent with the forward modelling results of the traditional grid extension method, indicating that the boundary truncation tuning layer method can reach a very accurate forward modelling result when applied to grid truncation of magnetotelluric forward modelling.

The results calculated by using the boundary truncation tuning layer method and the grid

extension method differ slightly from Zhengyong's data. This difference may be caused by use of

various impedance calculation modes, and the cause of the difference is understandable and

acceptable.

[00207] The boundary truncation tuning layer proposed by the present disclosure may be applied

to magnetotelluric forward modelling. The traditional boundary truncation tuning layer mainly

applies to the high-frequency electromagnetic wave, which is based on the displacement current.

This boundary truncation tuning layer is no longer applicable to forward modelling of the

low-frequency electromagnetic wave that is mainly based on a diffusion current. The boundary

truncation tuning layer proposed by the present disclosure probably truncates the simulation

region of the electromagnetic wave of the low-frequency diffusion field with stable absorption

performance of the external travelling wave and with more accurate forward modelling results.

The traditional boundary truncation tuning layer parameters are usually discrete in a polynomial

manner in the boundary truncation tuning layer, while the boundary truncation tuning layer

proposed in the present disclosure is exponentially discrete to achieve a better result. The air layer

is a significant obstacle for the boundary truncation tuning layer to be applied to magnetotelluric

forward modelling. Therefore, when setting the boundary truncation tuning layer parameters at the

outer boundary of the air layer, approximate setting is performed by using the earth background electrical conductivity. Numerical values have demonstrated that this procedure may be one of the best forms of compromise. In three-dimensional magnetotelluric forward modelling, a better grid truncation effect may be achieved by using the six-layer boundary truncation tuning layer. The more layers that the boundary truncation tuning layer has, the better the truncation effect is, however, the calculation required also increases. Different numbers of layers of the boundary truncation tuning layer require different optimal boundary truncation tuning layer parameters.

Once the best boundary truncation tuning layer parameters have been found, they may be applied

to the frequency range and the background resistivity range commonly used in the magnetotelluric

method. When the frequency parameters and background resistivity change during forward

modelling, there is no need to find the optimal boundary truncation tuning layer parameters again.

Therefore, in the embodiment of the present disclosure, a fixed number of layers of the boundary

truncation tuning layer may be used in future magnetotelluric forward modelling, and by always

using the same optimal boundary truncation tuning layer parameters workloads may be reduced

and a satisfactory grid truncation effect may be achieved. Usually, the traditional grid extension

method is aimed at the frequency range of the magnetotelluric method with a wide frequency band

and requires a long grid extension distance, which may be hundreds of kilometres long. The

extended grid size must be gradual, which causes the grid extension method to have a larger

number of extended grids. Thus, to obtain a reliable result, a substantial amount of memory and

calculation time may be required. When the boundary truncation tuning layer is applied to

magnetotelluric forward modelling, the six-layer boundary truncation tuning layer may obtain

very high-precision forward modelling results, with less memory occupied and less time for

calculation. Therefore, the boundary truncation tuning layer method is a highly recommended

method in magnetotelluric forward modelling.

[00208] An embodiment of the present disclosure further provides a computing device. Referring

to FIG. 10, the computing device includes a memory 1120, a processor 1110, and a computer

program that is stored in the memory 1120 and can be run by the processor 1110. The computer

program is stored in space 1130 for program codes in the memory 1120, and the computer

program is configured to execute any method of step 1131, according to the present disclosure

when executed by the processor 1110.

[00209] An embodiment of the present disclosure further provides a computer readable storage medium. Referring to FIG. 11, the computer readable storage medium includes a storage unit for program codes. The storage unit is provided with a program 1131 for executing the method steps according to the present disclosure, and the program is executed by a processor.

[00210] An embodiment of the present disclosure further provides a computer program product

containing instructions. When running on the computer, the computer program product causes the

computer to execute the method steps according to the present disclosure.

[00211] The above-described embodiments may be implemented completely or partly by

software, hardware, firmware, or any combination of these. When implemented by software, the

above-described embodiments may be implemented in a form of a computer program product

completely or partly. The computer program product includes one or more computer instructions.

When the computer loads and executes the computer program instructions, the flows or functions

described in the embodiments of the present disclosure are generated completely or partly. The

computer may be a general-purpose computer, a special-purpose computer, a computer network,

or any other programmable apparatus. The computer instructions may be stored in a computer

readable storage medium, or transmitted from one computer readable storage medium to another

computer readable storage medium. For example, the computer instructions may be transmitted

from one website, computer, server, or data centre to another website site, computer, server or data

centre in a wired manner (e.g., by coaxial cable, optical fibre, digital subscriber line) or by

wireless (e.g., infrared, radio, microwave, etc.). The computer readable storage medium may be

any medium that can be accessed by a computer or a data storage device including a server or a

data centre, integrated by one or more available media. The usable medium may be a magnetic

medium (e.g., a floppy disk, a hard disk, a magnetic tape), an optical medium (e.g., a DVD), or a

semiconductor medium (e.g., a Solid State Disk).

[00212] Professionals should be aware that the units and algorithm steps of the respective

examples described with the embodiments disclosed may be implemented by electronic hardware,

computer software, or a combination of the two, to clearly illustrate interchangeability of

hardware and software. In the above description, the composition and steps of the respective

examples have been generally described according to functions. Whether these functions are

executed by hardware or software depends on the specific application and design constraint

conditions of the technical solution. Professionals and technicians may use different methods to implement the described functions with respect to each specific application, but such implementation should not be considered beyond the scope of the present disclosure.

[00213] People who are usually skilled in their professionmay understand that all or part of the

steps in the method according to the foregoing embodiments may be completed as a program

instructs a processor. The program may be stored in a computer-readable storage medium, and the

storage medium is a non-transitory medium, for example, a random access memory, a read only

memory, a flash memory, a hard disk, a solid state drive, a magnetic tape, a floppy disk, an optical

disc and any combination of these.

[00214] The above are the preferred embodiments of the present disclosure, but the scope of the

disclosure is not limited to these embodiments. Any changes or substitutions that may readily

occur to one skilled in the profession without departing from the scope of the present disclosure

are intended to be encompassed by the present disclosure. The scope of the present disclosure is

subject to the scope of the claims below.

EDITORIAL NOTE 2021101629

There are 6 pages of claims only.

Claims

1. A BTL method for three-dimensional forward modelling of a low-frequency magnetotelluric

method, comprises

causing a low-frequency electromagnetic wave to be incident from a target region into an

anisotropic boundary layer without reflection, and determining a tuning parameter expression of

the anisotropic BTL;

performing gradual processing on a tuning parameter of the BTL to obtain an absorption

parameter of a transmitted wave, and determining an optimal tuning parameter value of the BTL

of a finite thickness; and

introducing the tuning parameter of the BTL of a finite thickness into an electromagnetic field

Maxwell equation according to a method adopted in low-frequency magnetotelluric forward

modelling, to obtain parameter information of the low-frequency magnetotelluric forward

modelling.

2. The method according to claim 1, wherein determining a tuning parameter expression of the

anisotropic BTL includes

setting electrical conductivity of an anisotropic medium to C7 A, a dielectric constant to 8, A, and magnetic conductivity to p *A, where A is a key tuning parameter and satisfies the formula

below:

s, 0 0 A= 0 s. 0

0 0 1

where s.is an anisotropic principal axis parameter;

setting a time factor to e , a Max equation expression of a frequency domain being

VxE=iopH

Vx H= aE - iwcE:,

where E is an electric field vector, H is a magnetic field vector, 0 is an angular frequency,

C is electrical conductivity of a target region medium, t is magnetic conductivity of the target

region medium, 8 is a dielectric constant of the target region medium, isan imaginary i=

number, 5E denotes a diffusion current, and ioE denotes a displacement current.

It is derived that

VxVxE+ -io - o2 sE= 0

let:

6= - io 02

It is calculated that

6= of 10W+ 2 2

setting the tuning parameter of the BTL to siI= xz an expression being

s, =I1+ a,1 6C

where Ci is a key parameter of the BTL, which controls attenuation of an electromagnetic wave

in the BTL, and it is obtained that

N2( si =1+

o) 1i+ _-1_u - i - 1+ +1 C

3. The method according to claim 2, wherein after determining a tuning parameter expression of

the anisotropic BTL, the method further comprises simplifying the tuning parameter of the BTL.

4. The method according to claim 3, wherein simplifying the tuning parameter of the BTL

includes

at a high-frequency limit o ? taking an expression of s, as , s:

si = 1+ -- ~ obtaining iOs , when the time factor is e ;

simplifying, at a low-frequency limit( ` (3,

-F s, =1+

o + - -1 - Fi + 8+1 as

5. The method according to claim 2, wherein after determining a tuning parameter expression of

the anisotropic BTL, the method further comprises determining a transmitted wave expression of

the BTL according to an electromagnetic wave expression of the target region and the tuning

parameter.

6. The method according to claim 5, wherein determining a transmitted wave expression of the

BTL according to an electromagnetic wave expression of the target region and the tuning

parameter includes

making, in a diffusion equation, a wave number equal to

k = 4-iopacy kk making, according to a phase tuning principle, the wave number k. of an incident wave in an x

direction equal to the wave number kt- of the transmitted wave in the x direction; and in a z

k k direction, the wave number kz of the incident wave and the wave number kz of the transmitted

wave to satisfy the relationship below:

k, = ski

expressing, in the target region, the electromagnetic wave as

E, = eEge*ek , where ei represents a polarisation direction of the electric field, and E0

represents an incident electric field amplitude;

expressing the obtained transmitted wave of the BTL as

E, = eE0 e'e'

supposing an included angle between an incident direction of the incident wave and a

decomposition surface is 0 , then a wave number expression in the z direction being

k, = -iopa~cosO, to derive

7. The method according to claim 4, wherein performing a gradual processing on a tuning

optimal tuning parameter value of the BTL of a finite thickness includes

s1+ modifying the tuning parameter of the BTL - into a form of

s.=K I +r ( aS- i)N to achieve a preset effect on absorption of the incident wave;

where K and at are two harmonic control parameters, K controls a propagation direction of

the transmitted wave and absorption of an evanescent wave; ai controls an attenuation speed of

the transmitted wave at different depths in the BTL; and yi, Ki and at work together to

jointly control attenuation of the transmitted wave;

setting the BTL of a finite thickness, making the electromagnetic wave attenuate to 0 at the outer

boundary of the BTL, and setting a Direchlet boundary condition on the outermost side of the

BTL;

performing gradual processing on the BTL parameter si, starting a gentle rise from 1 at an

interface between the target region and the BTL, accelerating a rising rate at a distance away from

the interface; and concluding settings of parameters C77, Ki and a' in the BTL parameter s,

as follows:

K =l+Kmax exP m

S= (Umax/( 7- / dz -)) exp m. -1

a = amax exP k -m -1

where d is a total thickness of the BTL; x is a distance from different units of the BTL to the

interface between the target region and the BTL; m is a constant for adjusting a change rate of

the parameter; dz is a thickness of a unit of the BTL; Kmax, Ymaxand Omax are respectively key parameters for controlling sizes of parameters C7 , and ai; and optimal values of Kmax , (max and max are OPt O°t and O°, which are selected by numerical experiments.

8. The method according to claim 7, wherein introducing the tuning parameter of the BTL of a

finite thickness into an electromagnetic field Maxwell equation, according to a method adopted in

low-frequency magnetotelluric forward modelling, to obtain parameter information of

low-frequency magnetotelluric forward modelling includes

using the BTL for grid truncation and absorption with respect to a quadratic field;

taking the tuning parameter A as a constitutive parameter matrix of the medium, and introducing

A into Maxwell equations VxE=iwoPA.H

V x H = c7A -E - i)E:A -E

a double curl expression of the electric field being

Vx iVxE (iopaco2ps)A-E=0

making, in the target region, the tuning parameter A equal to an identity matrix; expressing, in

the BTL, the tuning parameter A as

0 0 S 0 07 sY s, 0 0 A= 0 s 0 0 0 0 s 0 xS 0 0 s 'S1 0 0 S 10 0

0 0 Sx

= 0 0 0 0 xS Sz j

where an electromagnetic wave primary field of an earth background model and a total field when

the background model contains an anomalous body both conform to the double curl expression of

the electric field; and by subtracting a primary field formula from a total field formula, a

magnetotelluric quadratic field forward formula is obtained as follows:

Vx i V x ES -(iopaLo2pc)A -E =iops(c-c )A-E,

where y is electrical conductivity distribution of an earth background medium, C iselectrical

conductivity distribution after adding an anomalous body into the earth background medium, EP

is a primary field of the earth background model, and a one-dimensional analytical solution is used

for calculation;

solving a partial differential equation of a quadratic field forward formula by using a finite

element method; obtaining the quadratic field and thereafter adding the primary field to

obtain a total field value of the electric field; and

calculating through the electric field and the magnetic field to obtain impedance information to

further calculate apparent resistivity and phase.

9. A BTL apparatus for three-dimensional forward modelling of a low-frequency magnetotelluric

method, comprising

a parameter determining module, configured to cause a low-frequency electromagnetic wave to be

incident from a target region into an anisotropic BTL without reflection, and determine a tuning

parameter expression of the anisotropic BTL;

a gradual processing module configured to perform gradual processing on a tuning parameter of

the BTL to obtain an absorption parameter of a transmitted wave and determine an optimal tuning

parameter value of the BTL of a finite thickness; and

a forward modelling processing module, configured to introduce the tuning parameter of the BTL

of a finite thickness into an electromagnetic field Maxwell equation according to a method

adopted in low-frequency magnetotelluric forward modelling, to obtain parameter information of

low-frequency magnetotelluric forward modelling.

10. A computing device comprising a memory, a processor and a computer program stored in the

memory and run by the processor, wherein the processor implements the method according to any

of claims 1 to 8 when executing the computer program.