CC BY
401
УДК 528.34: 629.783 Ю.В. Сурнин СГГ А, Новосибирск
МАТЕМАТИЧЕСКАЯ МОДЕЛЬ ВНЕШНЕГО ГРАВИТАЦИОННОГО ПОЛЯ НА ОГРАНИЧЕННОМ УЧАСТКЕ ЗЕМНОЙ ПОВЕРХНОСТИ ПО СПУТНИКОВЫМ И ТРАДИЦИОННЫМ ГЕОДЕЗИЧЕСКИМ ДАННЫМ
Введение. В соответствии с теорией М.С. Молоденского [1], традиционные астрономо-геодезические методы определения плановых и высотных координат (триангуляция, трилатерация, полигонометрия, геометрическое нивелирование, астрономические определения), а также гравиметрические измерения дают возможность, определять поверхность Земли лишь совместно с определением ее внешнего гравитационного поля. В этой теории поверхность Земли, как граничная поверхность краевой задачи, на которой заданы результаты измерений, является неизвестной. Появление современных методов спутниковых координатных определений позволяет с высокой точностью определять координаты точек земной поверхности в единой геоцентрической системе отсчета. Таким образом, независимо от решения краевой задачи в постановке М.С. Молоденского, стало возможным получение физической поверхности Земли в виде множества дискретных точек в общей земной системе координат. В связи с этим, включение спутниковых высокоточных координатных данных в состав астрономогеодезической и гравиметрической информации с целью определения внешнего гравитационного поля Земли (ГПЗ) позволяет по иному подойти к решению задачи. Такой подход с одной стороны существенно облегчает ее решение, с другой повышает точность определения гравитационного поля Земли.
Кроме того, в современной практике астрономо-гравиметрического нивелирования квазигеоида (по методу М.С. Молоденского [2]) карты аномалий высот £ и составляющих отклонений отвеса £, ц строятся, главным образом, по аномалиям ускорения силы тяжести Дg посредством уточненных интегральных формул Стокса и Венинг-Мейноса. Астрономические долготы X, широты ф
и азимуты а на пунктах Лапласа, как дополнительная и независимая высокоточная информация о ГПЗ, традиционно непосредственно не вкладываются
в определение параметров единственной геопотенциальной функции -потенциала ускорения силы тяжести, которая порождает все разнообразие возможных трансформант гравитационного поля Земли. В интегральных уравнениях М.С. Молоденского [1] астрономические координаты X, ф выполняют лишь роль аргументов геопотенциальной функции (наряду или взамен геодезических долгот L и широт В). Но X и ф непосредственно характеризуют направление силовых линий ГПЗ на поверхности Земли. Астрономический азимут косвенно определяет составляющую отклонения
отвеса в плоскости первого вертикала ц через уравнение Лапласа. Поэтому не включение в состав измерительной информации астрономических величин X, ф и а, относительная погрешность которых одного порядка с первоклассным нивелированием и гравиметрическими измерениями, приводит к потере полезной информации и к снижению точности определения гравитационного поля Земли.
Метод М.С. Молоденского [2] определения фигуры геоида при совместном использовании астрономо-геодезических отклонений отвеса ц, £ и карт аномалий силы тяжести Дg, не решает задачу одновременного вложения трех видов измерений ц, £ и Дg в определение ГПЗ. Уравнения наблюдений для разностей между астрономо-геодезическими (индекс - АГ) и гравиметрическими (индекс - Гр) отклонениями отвеса формируются и решаются независимо друг от друга [3], [4]:
£аг - £гр = аф + ЬХ + с, (1)
Цаг - Цгр = аф + рх + у. (2)
Это приводит к тому, что два неизвестных параметра Ь и а (которые для единой потенциальной функции должны быть равными друг другу) определяются порознь из двух систем уравнений (1) и (2), ухудшая обусловленность и точность решения задачи. В числовом примере Н. П. Макарова [3, с. 281] Ь = -0,15",
а = +0,01". Различие в оценках Ь и а в этом примере лежит в доверительном интервале, определяемом точностью исходных данных.
Корректировка систем уравнений (1) и (2) - с учетом того, что составляющие отклонения отвеса есть производные одной потенциальной функции - возмущающего потенциала Т, приводит к решению объединенной системы из уравнений (1), (2) с пятью неизвестными параметрами {а, Ь, с, в, у}. Теперь в уравнениях вида (2) на месте параметра а должен стоять параметр Ь.
В современном комбинированном методе [5] развития сети нормальных высот на территории России с помощью геометрического нивелирования (в виде нормальных высот Н7) в комплексе со спутниковыми координатными GPS/ГЛОНАСС-определениями (в виде геодезических высот Н) и детальной гравиметрической съемкой (в виде аномалий силы тяжести Дg), информация об астрономических трансформантах ГПЗ X, ф, а, по-прежнему, не используется.
В усовершенствованной методике создания цифровых карт уклонений отвесных линий (ЦК УОЛ) [6] основными видами исходной информации являются:
- Гармонические коэффициенты до 36 степени модели «Параметры Земли 1990 года» для учета планетарной части ГПЗ;
- Средние значения аномалий силы тяжести (АСТ) по трапециям 1 х 1° и 5 х 7,5' для учета региональных особенностей ГПЗ;
- Гравиметрические карты масштаба 1 : 200 000 и топографические карты масштаба 1 : 50 000 для учета тонкой структуры ГПЗ.
В этой методике информация об астрономических долготах X и широтах Ф на 16 астропунктах для района Карпат площадью 40 х 60' использована только для контроля экспериментального варианта ЦК УОЛ.
Таким образом, в теории М.С. Молоденского [1] и ее практическом применении [2 - 6], имеющаяся высокоточная информация о ГПЗ в виде [Ь,
В, Н, Н7, X, ф, а, А^} вкладывается в определение ГПЗ не полностью и не всегда комплексно.
Первое предложение о комплексном использовании такого рода астрономо-геодезической и гравиметрической информации [Н7, X, ф, а, А^} для определения геопотенциала было высказано, по-видимому, Г. Морицем [7] в методе коллокации. Однако идея применения коллокации к задаче определения ГПЗ,
с нашей точки зрения, полезна только в том случае, когда «запросы потребителя» (количество, состав определяемых параметров, конструкция модели ГПЗ
и требуемая точность) не соизмеримы с имеющейся информацией о ГПЗ (по количеству, по распределению в пространстве, по составу, по точности). Другими словами, мы имеем дело с так называемыми, некорректно поставленными задачами. Такой, например, является задача определения глобальной модели ГПЗ в виде бесконечного ряда шаровых функций по информации [Ну, X, ф, а, Ag}, которая распределена крайне неравномерно из-за неполной гравиметрической изученности. В такой ситуации в методе коллокации Г. Морица предлагается недостаток информации восполнить введением дополнительной - эргодической - модели возмущающего потенциала Т. Кроме двух естественных моделей: детерминированной модели ГПЗ - (АХ) и стохастической модели (в виде ковариационной матрицы Кп для вектора п погрешностей измерений), в уравнение наблюдений (для вектора измерений I) включается стохастическая модель возмущающего потенциала Т, (в виде ВТ) [7, с. 187]:
I = АХ + ВТ + п. (3)
Но получение достоверной ковариационной функции Т требует знания аномальной части ГПЗ (которая неизвестна!). Для выхода из этого «заколдованного круга» в методе коллокации предлагается принять гипотезу об однородности и изотропности аномального, но не случайного поля Т. Таким образом, вынужденная необходимость использования в методе коллокации, в общем случае, неадекватной объекту (ГПЗ) математической модели сигнала Т, неизбежно будет приводить и к смещенным оценкам параметров и к искаженной оценке точности получаемого решения.
Регулярная - «правильная» - стратегия решения таких обратных задач геодезии должна строиться на трех условия регулярности [8]:
- Адекватности математической модели объекту (процессу); (4)
- Наблюдаемости математической модели по имеющейся измерительной информации;
- Состоятельности решения задачи.
Первое условие регулярности требует конструирования адекватной математической модели, состоящей из двух частей: детерминированной части (в виде модели АХ, в которой X включает все параметры ГПЗ, в том числе и аномальной части поля, влияние погрешности априорного знания которых превосходит погрешности используемых измерений) и стохастической части модели (в виде ковариационной матрицы Кп) для погрешностей п измерений и остаточных случайных погрешностей, не учитываемых детерминированной частью модели.
Проблемы, возникающие при выполнении второго условия регулярности
- наблюдаемости - должны решаться не на основе каких-либо гипотез - типа эргодичности поля аномалий силы тяжести, а путем, во-первых, привлечения достоверной информации о гравитационном поле, как в виде результатов измерений, так и иного характера, во-вторых, путем приведения в соответствие требований к количественному и качественному составу оцениваемых параметров детерминированной модели ГПЗ согласно количеству и качеству имеющейся
в распоряжении измерительной информации (по распределению измерений в пространстве-времени, по видам измерений и т. п.).
Третье условие регулярности - состоятельность решения - будет выполняться автоматически в случае выполнения первых двух условий регулярности и применения подходящего метода математико-статистической обработки результатов измерений - например, метода максимального правдоподобия или как частный случай метода наименьших квадратов.
Поэтому замечания Г. Морица [7, с. 174] о том, что «Операционный подход (метод коллокации, прим. авт.) к физической геодезии проявился относительно недавно - после того как появились новые возможности для геодезических измерений и, в особенности, после того как классический подход (теория М.С. Молоденского, прим. авт.) обнаружил неспособность дать полный ответ (определить глобальное ГПЗ, прим. авт.) при неполной гравиметрической изученности.», не совсем справедливо по отношению к глубокой и строгой теории М.С. Молоденского, имеющей большое практическое значение.
Теорию М.С. Молоденского определения глобального ГПЗ [1] и метод астрономо-гравиметрического нивелирования квазигеоида [2-6] на локальных участках земной поверхности необходимо дополнить включением современных видов измерений, и лучше не подменять ее коллокацией, используя неадекватные - эргодичные - модели реальных гравитационных аномалий. Кроме того, на основе теории М.С. Молоденского целесообразно построить такую модель ГПЗ, которую можно было бы сравнительно просто адаптировать с одной стороны - к «рельефу» аномалий геопотенциала, с другой - к имеющейся измерительной информации (распределению
измерений в пространстве, их составу
и точности).
В качестве такой модели наиболее целесообразно использовать конечноэлементную модель ГПЗ. Первые результаты применения метода конечных элементов для математического описания гравитационного поля получены в работах Я. Райнера [9] и А.В. Елагина [10].
В первой работе [9] создана законченная модель в виде опорной поверхности геоида, легко адаптируемая к измерениям и рельефу поверхности и основанная на использования почти всего спектра измерительной информации
{С = Н - И7, %, щ, Л^} (без астрономических азимутов а). Однако модель является двумерной (тогда как исследуемый объект - трехмерный) и содержит на границах конечных элементов разрывы первого рода, как самой функции, так и ее производных. Но величина разрывов по высотам геоида практически не ощутима (по-видимому, одного порядка с погрешностями исходных данных), по наклонам же поверхности величины разрывов должны превышать погрешности «измеренных» составляющих отклонений отвеса.
Во второй оригинальной работе [10] по определению параметров конечно-элементной модели возмущающего потенциала на локальном участке 5 земной поверхности (внутри нивелирного полигона) используется несколько иная измерительная информация (но без астрономических данных X, ф, а). По периметру полигона 5 требуются результаты нивелирных и гравиметрических измерений, внутри полигона 5 необходимо располагать гравиметрическими, градиентометрическими измерениями и данными для вычисления наклонов физической поверхности Земли по отношению к поверхности эллипсоида. В статье [10], к сожалению, не рассматривается вопрос о необходимой точности градиентометрических измерений и знания наклонов земной поверхности. Из рассмотрения двух работ [9], [10] следует, что создание эффективной математической модели ГПЗ с помощью метода конечных элементов требует дальнейших исследований. Поэтому совершенствование теории, методики и практики определения конечноэлементной модели ГПЗ является актуальной задачей.
В данной работе рассматривается один из этапов построения математической модели ГПЗ в соответствии с условиями регулярности (4) на локальном участке 5 земной поверхности. Область 5 рассматривается как одна ячейка
в конечно-элементной модели локального, регионального или глобального гравитационного поля Земли. Такое поле в области 5, в дальнейшем для краткости, будем называть локальным ГПЗ. Отличительными особенностями рассматриваемой модели локального ГПЗ является ее трехмерность и использование расширенного состава высокоточной измерительной информации, по сравнению
с работами [2 - 6, 9, 10], путем привлечения нивелирных, астрономических, гравиметрических и спутниковых данных к решению задачи определения
параметров одной геопотенциальной функции - возмущающего потенциала T.
Постановка задачи. Исходной информацией для отыскания параметров модели локального ГПЗ служат результаты спутниковых и традиционных (классических) средств измерений в узлах хаотичной пространственной решетки, расположенных на локальном участке 5 поверхности Земли.
В качестве результатов спутниковых измерений в рассматриваемом способе построения модели используются геодезические долготы, широты и высоты
{и, Bl, И}, г = 1, 2,..., k (5)
относительно общего земного эллипсоида с параметрами {аЕ, еЕ}. Не представляет никаких принципиальных трудностей вместо общеземных координат (5) использовать геодезические координаты {и, В'г, И'} в системе референц-эллипсоида с 8-ю параметрами {а'Е, е'Е, 5Ха, дУ0, дZo, юх, юу, ю2,} или квазигеоцентрические геодезические координаты В*г, И*г}
относительно квазигеоцентрического WGS-эллипсоида с параметрами {а*Е, е*Е, дХ*, ЗУ*, дZ*}. Чтобы не усложнять изложение методики определения модели ГПЗ, будем далее применять геодезические координаты (5), которые можно получить на основе радиотехнических средств наблюдений космических аппаратов, входящих в системы глобального позиционирования ГЛОНАСС и GPS, методом относительных спутниковых координатных определений (название метода заимствованы из работы [5]).
В качестве результатов традиционных геодезических измерений для определения параметров модели ГПЗ применяются нормальные высоты, астрономические долготы, широты и азимуты, а также ускорения силы тяжести
{И7г, Хг, фг, а, gi}, 1=1, 2,..., k, (6)
где I, k - порядковый номер узла и количество узлов решетки.
Каждый узел может содержать от максимального количества измерительной информации ^, В, Иг, И7, Хг, ф, аг, gi} до минимального { Ll, В, И, И1}.
Помимо этой информации необходимо знать ковариационную матрицу погрешностей измеренных величин (5) и (6), хотя бы диагонального вида.
Требуется решить - по терминологии Ю.М. Неймана [11] - задачу глобальной аппроксимации ГПЗ на локальном участке земной поверхности 5 в регулярной постановке [8]. Это означает, что требуется найти такую аналитическую функцию Т, которая имела бы простейшую структуру, наилучшим образом (в соответствии с принятым критерием) приближала возмущающий потенциал Т на 5 и отвечала трем условиям регулярности (4).
Построение модели состоит из трех основных частей: математического описания локального ГПЗ, методики решения обратной задачи (определение параметров модели по измерениям) и практической обработки результатов измерений с целью получения числовых значений параметров модели ГПЗ. Здесь рассмотрим только первые два этапа создания модели локального гравитационного поля: коротко - модель и методику.
Общее описание математической модели локального ГПЗ дается
формулами:
Ж = Ж(а, Ь, В, Н), (7)
Ж = и + Т, (8)
и = и(Ь, Ь, В, Н), (9)
Т = Т(с, Ь, В, Н), (10)
где Ж, и, Т - потенциалы, соответственно, ускорения силы тяжести, нормальный и возмущающий, как функции геодезических координат L, В, И текущей точки на физической поверхности Земли и векторов а, Ь, с, содержащих в качестве своих компонент параметры локальной модели ГПЗ. В частности,
а = {Ь, с} (11)
- Вектор параметров модели потенциала ускорения силы тяжести (7),
(8);
Ь = {ц, аЕ, еЕ,со], (12)
- Вектор параметров модели нормальной Земли (9), в котором ц -гравитационный параметр Земли (произведение универсальной постоянной тяготения И. Ньютона на массу Земли), со - угловая скорость вращения Земли, аЕ, еЕ - параметры уровенного общего земного эллипсоида, принятые при вычислении координат (5);
С = {С/}, 1 = 1,2,..., Ч (13)
- Вектор параметров модели локального возмущающего потенциала (10) в области 5.
Математическая модель нормальной Земли. В качестве нормального гравитационного поля Земли (9) принимается потенциал и нормального ускорения силы тяжести у уровенного двухосного эллипсоида с параметрами (12), вращающегося вместе с Землей. Ниже даются замкнутые формулы М.С. Молоденского [1], адаптированные автором для текущей точки ^, В, И), расположенной на или вне поверхности уровенного эллипсоида:
и = (р,/е)жсХ£(с/Ь) + (ш2а2/2)[(аЕ/ЬЕ)2(Ь/а)2(д/дЕ)(2/3 - б1п2м) + б1п2м],
(14)
у = &/а2 - (2/3)ю2Ь [(®(1 - (3/2)б1п2м + (3/2)б1п2м)]}(1 - (с/а)^1п2и)-112,
(15)
где
q = [3 + (с/Ь)2]агС£(с/Ь) - 3(с/Ь), дЕ = [3 + (с/ЬЕ)2]аг^(с/ЬЕ) - 3(с/Ье),
/ = 3[агС£(с/Ь) - (с/Ь)] + (с/а)2(е/Ь),
fE = {3[arctg(c/b) - (c/b)] + (c/bE)2arctg(c/bE)}(bE/aE)2,
_ и 2 _ 2 2
c = aEeE, bE = aE - c ,
a = (R1 + R2)/2, R1 = R2 + c2 + 2cD, R2 = R2 + c2 - 2cD, b2 = a2 - c2, sin2u = (D/a)2, D2 = X2 + Y2, R2 = D2 + Z2,
X = (N + H)cosBcosL, Y = (N + H)cosBsinL, Z = (N + H - NeE2)sinB, (lб) N = aE(1 - Єе sin B) .
Математическая модель возмущающего потенциала, ускорения силы тяжести и его основных трансформант. Аналитический вид модели T = T(c, L, B, H), которая аппроксимирует спутниковые (З) и наземные данные (б) в ограниченной области S, выбирается из класса гармонических функций. Для описания модели выбирается прямоугольная горизонтальная система координат (oxyz), начало которой помещается в точку д с координатами (Lo, Bo, Ho), расположенную, примерно, в средине локального участка земной поверхности S. Ось абсцисс ox направляется в геодезический зенит начальной точки d(Lo, Bo, Ho), ось аппликат oz - по линии пересечения плоскостей геодезического горизонта и вертикала с геодезическим азимутом Ao, ось ординат oy дополняет систему до правой. Азимут Ao задает ориентировку короткой оси симметрии сферической трапеции, охватывающей локальную область S, где лежит информация. В прямоугольной системе (oxyz) для текущей точки на поверхности Земли вводятся сферические координаты (u, v, r), связанные с прямоугольными (x, y, z) прямыми и обратными соотношениями:
.X = rcosvcosu, y =rcosvsinu, z =rsinv,
d2 = x2 + y2, r2 = d2 + z2, cosu = x/d, sinu = y/d, cosv = d/r, sinv = z/r. (17) Переход от геодезических координат (L, B, H) текущей точки, заданных в общеземной системе, к прямоугольным координатам (x, y, z) в локальной системе (oxyz) осуществляется по формуле
[x-y- z]T = {[X YZ]T- [Xo YoZo]T}R3(-Lo)R2(Bo)R1(Ao), (18)
где {X, Y, Z} и {Xo, Yo, Zo} - прямоугольные координаты в общеземной системе отсчета текущей (L, B, H) и начальной (Lo, Bo, Ho) точек, вычисляемые по формулам вида (1б). Ri(Z) - обозначение элементарной матрицы вращения вокруг координатной оси с номером і на угол Z. Положительный поворот Z считается происходящим против часовой стрелки для наблюдателя, смотрящего с положительного конца оси поворота на начало координат. Номер координатной оси устанавливается по правилу: для оси абсцисс і = 1, для оси ординат і = 2, для оси аппликат і = 3.
Возмущающий потенциал Т и его частные производные Ти, Ту, Тг по касательным к координатным линиям и, V, г, как функции гармонических коэффициентов {спт, snm} и сферических координат текущей точки {и, V, г} аппроксимируются следующими разложениями в ряды шаровых функций:
Т= | Д.
и+1 Z(cnm cos mu + snm slnmu)Pnm(sln v)-n=О rn+1 m=О
1
rcosv 1
■T =
Tu
dT
N 1 n / • ^Pnmfanv)
-----—= І-Г77 I m' (snm cos mu~ cnm sin mu)-------------,
rcosvou n=1 rn +2 m=О cosv
^ 8T N 1 n / , • fiPnm(sinv)
-■Tv= — = I—- X (cnmcosmu + snmsinmu)----------------Г-------
r rdv n=1 r 2 m=О ov
ST
N n +1 n
dr n=О rn+2 m=О
I(cnm cos mu + snm sln mu)Pnm (sln v)-
(19)
Для вычисления сферических гармоник применяются следующие рекуррентные соотношения:
m
для n=m, Рпт (sin v) = cos" vn (li -1);
7=1
для n = m + 1, Pm+1,m(sinv) = (2m + 1)sinvPnm(sinv);
для n Ф m, (n - m)Pnm(sinv) = (2n - 1)sinvPn-1m(sinv) - (n + m - 1)Pn-2,m(sinv);
dPnm(sinv)/dv = (1 - Snm)Pn,m+i(sinv) - mtgvPnm(sinv); cosmu = cos[(m - 1)u]cosu - sin[(m -1)u]sinu; sinmu = sin[(m-1)u]cosu + cos[(m - 1)u]sinu;
dnm = {1, если n = m, иначе 0}. (20)
Связь измеряемых величин с определяемыми параметрами ГПЗ.
Основные «измеряемые» трансформанты возмущающего потенциала T, такие
как:
С = н-н7,
щ = (X - L)cosB,
£ = ф - B, (21)
Ла = а - A,
= g - Y
- Аномалии высоты
- Отклонения отвеса в плоскости первого вертикала
- Отклонения отвеса в плоскости меридиана
- Разность между геодезическим и астрономическим азимутами
- Чистые аномалии ускорения силы тяжести
связываются с возмущающим потенциалом Т и его частными производными Ти,Ту, Тг равенствами:
ЛТТ = [С ц(Ы + Н) Ш + Н) А&(-унЯ\
(22)
Т = [т -Ти -Ту -Тг ]т, лтт = ТТк1сТдоТиТ,
(23)
^ = ^{у, гсоъу/(Ы+Н), г/(М + И), -ун}, Уи = -2у/Ят,
СТ = diag{1, Я\(п/2 - у)Яз(п/2+и)}, доТ = diag{1, Яі(-Ао)Я2(-Во)Яз(Ьо)},
(24)
И = ё1ав{1, Яі(ж/2 - Б)Яз(ж/2 + £)}, £ = 4 + ^,
5% = [>(Н 11т)81п2П, (> (ур уе)/уе,
да = ЯуАа, А а = ц^Ві + (&іпАу - яь- = Д/ш0г},
(25)
где лтт, тТ - векторы, составленные из трансформант возмущающего потенциала Т; в - гравиметрическое сжатие уровенного (нормального) эллипсоида; ур, уе, у - ускорения силы тяжести на полюсе (В = п/2), на экваторе (В = 0) и на широте В нормального эллипсоида, Ы, М, Ят - радиусы кривизн эллипсоида на широте В, соответственно: меридиана, первого вертикала и среднее интегральной значение; ун - приближенное значение вертикального градиента нормального ускорения силы тяжести на широте В и высоте И = 0.
Методика составления и решения системы уравнений наблюдений.
Вектор опорных параметров модели локального ГПЗ (не уточняемых по результатам измерений) формируется из четырех констант нормального поля и и четырех параметров ориентирования локальной системы координат
(оху£)
Ь = {и, ав, ев, т, Lo, Во, Но, Л}. (26)
Вектор оцениваемых параметров модели локального ГПЗ образуется из коэффициентов разложения в ряд шаровых функций возмущающего потенциала Т, упорядоченных по правилу
с = {C00, Cl0, С1Ь C20, С21 C22,■■■ , СЫЫ s11, Я2Ь s22,■■■, SNN}, (27)
где индекс N - максимальный порядок разложений (19).
Неизвестный вектор с размера q х1 (д = (Ы + 1) ) определяется под условием метода наименьших квадратов (МНК) из решения линейной системы уравнений
Ас = / + V, К„ (28)
/={-.ЛлЛЛ А-.-Ь
/= Q = Hi - fi= ПгЩ + Щ= (К - Ц(Ыг + Hi)cOSBi, /,■ = UMr + Щ, /ёг = Ag/(-Ун),
= лaijDijSin&ij = (Оу - AiJDijSinGij'
где f - вектор свободных членов размераp х 1, A - матрица размераpxq коэффициентов {Ay}, вычисляемых на основе формул (19) - (25), v - вектор случайных погрешностей измерений и модели размера pxl; Kv -ковариационная матрица вектора v.
Оценка векторов с и v под условием МНК производится по формулам:
С = Af v = (aA - I)f, A = K;mA, f = Kf2f (29)
Aj
Псевдообратную матрицу A рекомендуется вычислять посредством сингулярного разложения [12], [13]
A# = WZ#UT, Z = diag{ou ..., aq}, (ЗО)
где а1, ..., aq - сингулярные числа матрицы A; U и W- ортогональные матрицы левых и правых сингулярных векторов матрицы A.
Ковариационные матрицы K£ и Ky найденных оценок с и v параметров модели ГПЗ находятся по формулам:
kc_ = я2 wz#2 WT, kv=я2(АА# -1) = M2(UrUrT -1)) = я22u0uJ, (з 1)
о о
Я = iivl \/(p-r), r = rang (A),
где Uo - ортонормированный базис для ортогонального дополнения области значений оператора A, Ur - ортонормированный базис для области
значений A. Прямоугольные матрицы Uo и Ur имеют размеры p х (p - r) и
2
p х r и являются блоками матрицы U = [Ur|Uo], д_ - масштабный множитель, приводящий в соответствие априорную Kv и апостериорную Кг ковариационные матрицы.
Декомпозиция и регуляризация решения. Анализ сингулярного спектра {а1, ..., aq} матрицы A, расположенного в порядке убывания сингулярных чисел, и анализ вектора свободных членов G = {g1, ..., gp} диагональной системы уравнений
ZY = G, Y = WTc, G = UTf (32)
где Y = {y1, ..., yq} = преобразованный вектор неизвестных параметров c, позволяют произвести регуляризацию [14, 15] решения обратной задачи (28)
- (31) и обоснованно определить количественный состав вектора
оцениваемых параметров с на основе информации о погрешностях вектора f и матрицы A. Для установления качественного состава вектора с, т. е. для указания конкретных параметров C/, оценки которых Cj получаются достоверными, необходимо выполнить физическую декомпозицию [14, 16] системы уравнений наблюдений (28), таким образом, чтобы правые
сингулярные векторы были близки к единичным векторам канонического вида. Иначе говоря, чтобы физическое пространство оцениваемых параметров (в котором ищутся параметры с) практически совпадало с алгебраическим пространством области определения оператора A. Более детальное описание процесса регуляризации и декомпозиции обратной задачи (28) - (31) является предметом отдельного исследования, которое выходит за рамки данной работы.
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Молоденский М.С., Еремеев В.Ф., Юркина М.И. Методы изучения внешнего гравитационного поля Земли // Тр. ЦНИИГАиК. - 1960. - Вып. 131. - С. 251.
2. Молоденский М.С. Определение фигуры геоида при совместном использовании астрономо-геодезических уклонений отвеса и карты аномалий силы тяжести // Тр. ЦНИИГАиК, вып. 17. - Редбюро ГУГСК, 1937.
3. Макаров Н.П. Геодезическая гравиметрия // М.: Недра, 1968. С. - 408.
4. Бровар В.В., Магницкий В.А., Шимберев Б.П. Теория фигуры Земли. - М.: Геодезиздат, 1961. - С. 256.
5. Непоклонов В.Б., Чугунов И.П., Яковенко П.Э., Орлов В.В. Новые возможности развития сети нормальных высот на территории России // Геодезия и картография. - 1996.
- № 7. - С. 20 - 22.
6. Непоклонов В.Б., Яковенко П.Э., Кузьмин Ю.А., Переверткин С.В. О создании цифровых карт уклонений отвесных линий / Геодезия и картография. - 1996. - № 9. - С. 1
- 5.
7. Мориц Г. Современная физическая геодезия // М.: Недра, 1983. - С. 393.
8. Брандин В.Н., Разоренов Г.Н. Определение траекторий космических аппаратов // «Машиностроение». - М., 1978. - С. 216.
9. Reiner J. Precise Transformation of Classical Networks to ITRF by COPAG and Precise Vertical Reference Surface Representation by DFHRS - General Concepts and Realisation of Databases for GIS, GNSS and Navigation Applications. // Hochschule Karlsruhe Technik und Wirtschaft - University of Applied Sciences. - Siberian State Academy of Geodesy (SSGA). - Novosibirsk, 17-26 Februar 2006. - www.dfhbf.de
10. Елагин А.В. Применение метода конечных элементов для решения краевой задачи теории потенциала внутри нивелирного полигона // Вестник СГГА, вып. 10. -Новосибирск. - 2005. - С. 25 - 30.
11. Журкин И. Г., Нейман Ю. М. Методы вычислений в геодезии // М.: Недра. -1988. - С. 304.
12. Форсайт Дж., Малькольм М., Моулер К. Машинные методы математических вычислений // М.: Мир, 1980. - С. 280.
13. Годунов С.К., Антонов А.Г., Кирилюк О.П., Костин В.И. Гарантированная точность решения систем линейных уравнений в евклидовых пространствах // Новосибирск: ВО «Наука», 1992. - С. 360.
14. Sournin Yu. The regular approach to the estimation of parameters of the mathematical model of the Earth’s crust motion and displacements using satellite data // PROCEEGINGS OF THE INTERNATIOAL SEMINAR. “On the Use of Space Techniques for Asia-Pacific Regional Crustal Movements Studies” (Project: Asia-Pacific Space Geodinamcs) APSG-IRKUTSK-2002. Irkutsk, 5-10 August, 2002. Moscow GEOS. - 2002. - С. 206 - 212.
15. Сурнин Ю.В. Сравнительный анализ непрерывной и дискретной регуляризации решений некорректных задач космической геодезии // ТРЕТИИ СИБИРСКИМ КОНГРЕСС ПО ПРИКЛАДНОЙ И ИНДУСТРИАЛЬНОЙ МАТЕМАТИКЕ (ИНПРИМ-98). Тез. докл., часть III, Изд. Института математики СО РАН. - 1998. - С. 122.
16. Сурнин Ю.В., Гиенко Е.Г. Алгебраическая и физическая декомпозиция математических моделей при решении плохо обусловленных обратных задач геодезии // ЧЕТВЕРТЫЙ СИБИРСКИЙ КОНГРЕСС ПО ПРИКЛАДНОЙ И ИНДУСТРИАЛЬНОЙ МАТЕМАТИКЕ (ИНПРИМ-2000). Тез. докл., часть IV, изд. Института математики СО РАН.
- 2000. - С. 73 - 74.
© Ю.В. Сурнин, 2006
