А.А. Булычев, И.В. Лыгин, В.Р. Мелихов - Численные методы решения прямых задач грави- и магниторазведки, страница 13

С учетом того, что∑ν =1σ ν +1 − σ ν, и окончание N–ой стороны совпадает с началом 1–ойσ ν +1 − σ νстороны, т.е. σ N +1 = σ 1 , получим:αν =Nσ−σν +1ν(σ ν +1 − σ ν ) = (σ 2 − σ 1 ) + (σ 3 − σ 2 ) + L + (σ 1 − σ∑−σν =1 σν +1N) = 0.νОкончательное выражение функции G(s) для многоугольника спостоянной плотностью, будет иметь вид:NG ( s ) = Gδ ∑ (αν s + βν − s ) lnν =1σ ν +1 − s.σν − sПолучим теперь выражение для комплексной напряженностимагнитного поля, создаваемого этой же моделью с постояннойнамагниченностью I. Напомним, что под величиной I понимаетсякомплекснаявеличина,образованнаяизкомпонентвекторанамагниченности: I = Ix + iIz.

При выводе выражения для поля H(s)воспользуемся 2–ой формулой О.–Г. и явным заданием уравнений сторонмногоугольника:85NNdσ11= I∑H ( s ) = 2iI ∫dS = I ∫dσ = I ∑ ∫2ν =1 Гν σ − sν =1∂ Dσ − sD (σ − s )Nσ ν +1ν =1σν= I ∑ ανОкончательноNH ( s ) = I ∑ α ν lnν =1∫σ ν +1∫σνd (αν σ + βν )=σ −sNdσσ= I ∑ αν ln(σ − s ) σν +1 .νσ −sν =1σ ν +1 − s.σν − s17. Теперь перейдем к рассмотрению более сложных моделей. Одна изтаких моделейбудет представлять собой многоугольник скриволинейными сторонами. Для получения аналитических выраженийполей необходимо задать уравнение дуг, описывающих сторонымногоугольника. Одной из форм представления таких кривых может бытьполином n–ой степени:σ ν ( t ) = aνn t n + aνn−1t n−1 + L + aν0 ,где aνj – заданные коэффициенты, 0 ≤ t ≤ 1 .18.

Получим выражение для поля G(s) от такого многоугольника спостоянной плотностью δ. Для простоты будем решать только внешнююзадачу. В этом случае мы можем записать:1σdS = Gδ ∫dσ = Gδ ∫ σ d ln(σ − s ) =Dσ − s∂ Dσ − s∂DG ( s ) = 2iGδ ∫N 1N= Gδ ∑ ∫ σ d ln(σ − s ) = Gδ ∑ ∫ σ ν ( t ) d ln(σ ν ( t ) − s ) .ν =1 Гνν =1 0Выражение σ ν ( t ) − s для случая задания σ ν ( t ) в виде полинома n–ой степени от t можно представить в виде:σ ν ( t ) − s = aνn (t − t1,ν ( s ))(t − t 2,ν ( s ))L (t − t n ,ν ( s )) = aνn ∏ (t − t j ,ν ( s )),nj =1где t j ,ν ( s ) – корни уравнения σ ν ( t ) − s = 0 .

В результате такогопредставления можно записать:86⎞ n⎛ ν ndt.d ln(σ ( t ) − s ) = d ln⎜⎜ a n ∏ (t − t j ,ν ( s ))⎟⎟ = ∑t−t(s)j=1j=1j ,ν⎠⎝νВ результатеN 1N ⎛ n 1⎞σ ν (t )G ( S ) = Gδ ∑ ∫ σ ν ( t ) d ln(σ ν ( t ) − s ) = Gδ ∑ ⎜ ∑ ∫dt ⎟ .⎜⎟ν =1 0ν =1 ⎝ j =1 0 t − t j ,ν ( s )⎠19. Для поля H(s), создаваемого такой же моделью с постояннойкомплексной намагниченностью I, можно, сделав следующиепреобразования, получить:Ndσd ln(σ − s )1H ( s ) = 2iI ∫dS = I ∫= I∑ ∫dσ =2dσν =1 ГνD (σ − s )∂ Dσ − sN σ ν +1= I∑∫ν =1 σ νN 1dσd ln(σ − s ) = I ∑ ∫dσν =1 0dσ ν (t )dt d ln(σ ν ( t ) − s ) =dσ ν ( t )dt⎛⎞dσ ν (t )⎜⎟Nn 1dtdt⎜⎟.= I∑ ∑∫⋅ν⎜⎟ν =1 j =1 0 dσ ( t ) t − t j ,ν ( s )⎜⎟dt⎝⎠20. В качестве примера получим выражение для H(s) от многоугольникас постоянной намагниченностью. В этом случае уравнение ν–ой стороныбудет представляться в параметрическом виде следующим образом:σ ν ( t ) = aν t + bν .Уравнение σ ν ( t ) − s = 0 имеет единственный корень t1,ν ( s ) =NdσdσH ( s) = I ∫=I∫d ln(σ − s ) = I ∑ν =1∂ Dσ − s∂D dσσ ν +1∫σνdσd ln(σ − s ) =dσNad (aν t + bν )= I∑∫d ln( aν t + bν − s ) = I ∑ νν =1 0 d ( aν t + bν )ν =1 aνN 187s − bν.

Тогдаaν1dt∫0 t − t1 ν ( s ) =,NN(a + bν ) − saaν 1 − t1,ν ( s )aν1= I ∑ ln(t − t1,ν ( s ))0 = I ∑ ln.= I ∑ ν ln νbν − s− t1,ν ( s )ν =1 aνν =1 aνν =1 aνNЕсли заданы координаты вершин многоугольника, то коэффициенты aν иbν определяются из следующих условий:σ ν (0) = bν = σ ν , σ ν (1) = aν + bν = σ ν +1 , aν = σ ν +1 − σ ν .Окончательно для поля H(s) можно записатьaν σ ν +1 − sln,σν − sν =1 aνNH ( s) = I ∑что совпадает с выражением, полученным нами ранее.21.

Поучим выражение для поля G(s) от параболическогомногоугольника, т.е. от многоугольника ограниченного дугами общимчислом N, и каждая из которых представлена полиномом второго порядкаσ ν ( t ) = aν t 2 + bν t + cν , 0 ≤ t ≤ 1 .Уравнение σ ν ( t ) − s = 0 имеет два корня:− bν + bν2 − 4aν (cν − s ),t1,ν ( s ) =2aν− bν − bν2 − 4aν (cν − s ).t 2 ,ν ( s ) =2aνТогдаN2 1G ( s ) = Gδ ∑∑ ∫ν =1 j =1 0σ ν (t )t − t j ,ν ( s )dt .Для того чтобы аналитически взять этот интеграл, преобразуемвыражение для σ ν ( t ) следующим образом:σ ν ( t ) = aν t 2 + bν t + cν = aν ( t − t j ,ν ( s )) 2 + bν ( t − t j ,ν ( s )) + cν ++ 2aν t ⋅ t j ,ν ( s ) + bν t j ,ν ( s ) − aν t 2j ,ν ( s ) =()= aν ( t − t j ,ν ( s )) 2 + bν + 2aν t j ,ν ( s ) ( t − t j ,ν ( s )) + aν t 2j ,ν ( s ) + bν t j ,ν ( s ) + cν .88Теперь для G(s) можно записать:⎧ 1G ( s ) = Gδ ∑∑ ∫dt = Gδ ∑∑ ⎨aν ∫ ( t − t j ,ν ( s ))dt +ν =1 j =1 0 t − t j ,ν ( s )ν =1 j =1 ⎩02 1Nσ ν (t )+ (bν+ 2aν t j ,νN2( s ))∫ dt + (a t1ν2j ,ν( s ) + bν t j ,ν ( s ) + cν0)∫ t − tdt ( s ) ⎫⎪⎬ =10j ,ν⎪⎭N 2 ⎧⎛1 ⎞⎟ ⎫⎪⎪a= Gδ ∑∑ ⎨ ν + bν + aν t j ,ν ( s ) + (aν t 2j ,ν ( s ) + bν t j ,ν ( s ) + cν )ln⎜⎜ 1 −⎟ ⎬⎪ .2t(s)ν =1 j =1 ⎪j ,ν⎝⎠⎭⎩22.

Следующая модель – многоугольник с полиномиальнымраспределением плотности δ (ξ ,ζ ) = Pn (ξ ,ζ ) , n – степень полинома:δ (ξ , ζ ) = Pn (ξ , ζ ) =Используя подстановку ξ =σ +σnm∑0 ∑0 Am ,km= k=, ζ =ξ kζ m −k .σ −σполучим выражение для22iплотности как функции переменных σ и σ :nmδˆ (σ ,σ ) = ∑∑m =0 k =0Am ,k2m i m−knm(σ + σ ) k (σ − σ ) m − k = ∑∑ a m ,k σ k σ m − k = Pˆn (σ ,σ ) .m =0 k =0При переходе от площадного интеграла к контурному по 1–ойформуле О.–Г. в комплексной форме вычисляется формальный интегралΦ(σ ,σ ) = ∫ Pn (σ ,σ )dσ =nm∑∑ am =0 k =0m ,kσkσ m − k +1m − k +1.Получим выражение поля для всех точек пространства, т.е. как вне, так ивнутри области D. В этом случае выражение для G(s) будет иметь вид:NΦ(σ ,σ ) − Φ( s , s )G( s) = G ∫dσ = G ∑σ −sν =1∂Dσ ν +1∫σνΦ(σ ,σ ) − Φ( s , s )dσ .σ −sПридальнейшихпреобразованияхбудемиспользоватьпредставление уравнения прямой в явном виде σ = αν σ + βν .

Вопрос о89том, как найти коэффициенты этого уравнения нами уже обсуждалсяранее. Используя это представление, получим:Φ(σ ,σ ) = Φ(σ ,αν σ + β ν ) = Ρn+1,ν (σ ) .Таким образом для каждой ν–ой стороны функция Φ(σ ,σ ) будетпредставляться в виде нового полинома от переменной σ степени (n+1).Далее можно осуществить следующие преобразования:N σ ν +1G( s) = G ∑ν =1∫σΡn+1,ν (σ ) − Ρn+1,ν ( s ) + Ρn+1,ν ( s ) − Φ( s , s )σ −sνdσ =σ ν +1⎧⎪σν +1 Ρn+1,ν (σ ) − Ρn+1,ν ( s )dσ ⎫⎪dσ + (Ρn+1,ν ( s ) − Φ( s , s )) ∫= G∑ ⎨ ∫⎬.σsσs−−⎪⎭ν =1 ⎪σσ⎩ ννNЧислитель первого интеграла может быть представлен как некоторыйполином степени (n+1) от разности (σ – s) – Ρˆ n+1,ν (σ − s ) . В результатеподинтегральное выражение будет представлять собой полином n–ойстепени от (σ – s).

Первообразная от такого полинома вновь будетполиномом степени (n+1), и в результате:σ ν +1∫σΡn+1,ν (σ ) − Ρn+1,ν ( s )σ −sνdσ =σ ν +1∫σνΡˆ n+1,ν (σ − s )σ −sdσ = Qn+1,ν ( s ) .Окончательно:N⎧σ − s⎫G ( s ) = G ∑ ⎨Qn+1,ν ( s ) + (Ρn+1,ν ( s ) − Φ( s , s ))ln ν +1⎬.σs−ν =1 ⎩ν⎭23. Получим выражение для гравитационного поля для многоугольникас линейно изменяющейся плотностью δ ( x , z ) = a x x + a z z + a 0 . Дляупрощения выкладок будем решать только внешнюю задачу.На первом шаге выразим закон изменения плотности черезкомплексные координаты σ и σ :δˆ (σ ,σ ) = a xσ +σ2=+ azσ −σ2i+ a0 =aσaσa xσ a x σ+− i z + i z + a0 =2222a x − ia za + ia zσ+ xσ + a 0 = cσ + c σ + a 0 .2290На втором шаге получим выражение для функции Φ(σ ,σ ) :Φ(σ ,σ ) = ∫ (cσ + c σ + a 0 )dσ = cσσ + cσ22+ a 0σ .На следующем шаге для каждой ν–ой стороны получим:(αν σ + βν ) 2Φ(σ ,σ ) = Φ(σ ,αν σ + βν ) = Ρ2 ,ν (σ ) = (cσ + a 0 )(αν σ + βν ) + c=2= cαν σ + cβν σ + a 0αν σ + a 0 βν + c2= (cαν + cαν22αν22σ + c αν β ν σ + c2)σ + (cβν + a 0αν + c αν βν )σ + a 0 β ν + c2βν22βν22=== K 2 ,ν σ 2 + K 1,ν σ + K 0 ,ν .Далее преобразуем это выражение:Ρˆ 2 ,ν (σ − s ) = K 2 ,ν (σ − s ) 2 + 2 K 2 ,ν σ s − K 2 ,ν s 2 + K 1,ν σ + K 0 ,ν == K 2 ,ν (σ − s ) 2 + ( 2 K 2 ,ν s + K 1,ν )(σ − s ) + ( 2 K 2 ,ν s + K 1,ν ) s − K 2 ,ν s 2 + K 0 ,ν == K 2 ,ν (σ − s ) 2 + ( 2 K 2 ,ν s + K 1,ν )(σ − s ) + ( K 2 ,ν s 2 + K 1,ν s + K 0 ,ν ) .Теперь можем записать окончательное выражение для G(s):⎧KG ( s ) = G ∑ ⎨ 2 ,ν (σ ν2+1 − σ ν2 ) + (K 2 ,ν s + K 1,ν )(σ ν +1 − σ ν ) +ν =1 ⎩ 2N⎛ σ − s ⎞⎫⎟⎟ ⎬ .+ (K 2 ,ν s 2 + K 1,ν s + K 0 ,ν )ln⎜⎜ ν +1σs−⎝ ν⎠⎭24.

Следующая модель – криволинейный многоугольник со сторонами,представляемыми в виде полинома степени pν от параметра t, и сплотностью δ, описываемой полиномом степени n. Осуществимследующие преобразования.–Представим закон изменения плотности как функцию двухпеременных σ и σ :91nnmδˆ (σ ,σ ) = ∑∑m=0 k =0–Am ,k2m i m−km∑∑ Aδ (ξ , ζ ) = Pn (ξ , ζ ) =m ,km =0 k =0nξ kζ m −k ,m(σ + σ ) k (σ − σ ) m − k = ∑∑ a m ,k σ k σ m − k = Pˆn (σ ,σ ) ;m =0 k =0Определим следующий интегралΦ(σ ,σ ) = ∫ Pn (σ ,σ )dσ =–nm∑∑ am =0 k =0m ,kσσ m − k +1km − k +1;Выпишем выражение для функции G(s) во внешнем пространстве:NΦ(σ ,σ )Φ(σ ,σ )G( s) = G ∫dσ = G ∑ ∫dσ =ν =1 Гν σ − s∂D σ − sN 1= G ∑ ∫ Φ(σ ν ( t ),σ ν ( t ))d ln(σ ν ( t ) − s ) .ν =1 0Поскольку уравнение ν–ой стороны представляется полиномом степени pν(pν)( σ ( t ) = ∑ bνl t l ), то функция Φ σ ν ( t ),σ ν ( t ) будет представлять собойνl =0полином степени qν = n+pν+1, т.е.()qνΦ σ ( t ),σ ( t ) = ∑ cνj t j = Ρqν ( t ) .ννj=0pνПредставив σ ( t ) − s = α pν ∏ (t − t k ,ν ( s )) , можем записать выражениеννk =1для G(s):Npν 1G ( s ) = G ∑∑ ∫ν =1 k =1 0Ρqν ( t )t − t k ,ν ( s )dt .Дальнейшие действия по вычислению таких интегралов нами ужерассматривались.25.

Подведем некоторый итог проделанной работе. Нами были полученывыражения для гравитационного и магнитного поля от тел достаточнопростой формы и от более сложных конфигураций источников. В92частности, нами были рассмотрены тела, стороны которых описывалисьполиномом некоторой степени. В общем случае мы не выводили явныхвыражений полей, создаваемых такими моделями, а определили ходполучения решений. По сути дела нами описан алгоритм вычисленияполей от таких объектов.С этой точки зрения рассмотрим еще раз задачу определения поляG(s) от криволинейного многоугольника с N вершинами и с постояннойплотностью δ.

Если дуга задана уравнениемσ ( t ) = a n t + a n−1tνν nνn −1n+ L + a 0 = ∑ aνk t k ,νk =0то корни уравнения σ ν ( t ) − s = 0 могут быть найдены числено с помощьюитеративных процедур (метод Ньютона и т.п.). После того как корниnнайдены, степенной ряд σ ( t ) = ∑ aνk t k преобразуется к виду:νk =0∑ aνk t k = ∑ Akν ( t j ,ν ( s )) ⋅ (t − t j ,ν ( s )) .nnk =0k =0kПроцедура преобразований коэффициентов полинома aνk в коэффициентыAνk ( t j ,ν ( s )) может быть запрограммирована.

В результате:nNn 1G ( s ) = Gδ ∑∑ ∫ν =1 j =1 0σ (t )νt − t j ,ν ( s )∑ aν tn 1Ndt = Gδ ∑∑ ∫ν =1 j =1 0k =0kkt − t j ,ν ( s )dt =⎛ n 1 n ν⎞k −1= Gδ ∑ ⎜⎜ ∑ ∫ ∑ Ak ( t j ,ν ( s )) ⋅ (t − t j ,ν ( s )) dt ⎟⎟ .ν =1 ⎝ j =1 0 k = 0⎠NВычисление определенных интегралов от степенных функций такжелегко программируемая операция, поскольку1∫ (t − tj ,ν( s )) dt =p( t − t j ,ν ( s )) p+1p +101=0(1 − t j ,ν ( s )) p+1p +1и при k = 0⎛1 ⎞⎟dt⎜1=−ln∫0 t − t j ,ν ( s ) ⎜ t j ,ν ( s ) ⎟ .⎝⎠193−( − t j ,ν ( s )) p+1p +1,Таким образом, для вычисления аномальных полей от таких тел нетнеобходимости получать формулы в окончательном виде. По сути дела всявычислительная процедура сводится к преобразованию полиномов.Литература.1. Дополнительные главы курса гравиразведки и магниторазведки.

–Новосибирск. Изд-во НГУ. 1966. 560 с.2. Гордин В.М., Розе Е.Н., Углов Б.Д. Морская магнитометрия. – М.:Недра. 1986. 232 с.3. Жданов М.С. Аналоги интеграла типа Коши в теории геофизическихполей. – М.: Наука. 1984. 327 с.4. Страхов В.Н. Методы интерпретации гравитационных и магнитныханомалий. – Пермь.: ПГУ.