Тема 6. Геометрия земного эллипсоида.
1. Земной эллипсоид.
1.2 Методы определения его параметров.
2. Исходные геодезические даты.
3. Системы координат на поверхности земного эллипсоида.
4. Линейный элемент поверхности.
5. Дифференциалы и длины дуг меридианов и параллелей.
6. Площадь сферической трапеции.
7. Уравнение геодезической линии на поверхности эллипсоида.
Вопрос 1.
Поскольку земной эллипсоид является поверхностью относительной, на
которой решаются задачи по определению координат точек из которого точки и
линии проецируются на плоскость, то необходимо задать (округлить) параметры
земного эллипсоида.
Земной эллипсоид образуется вращением эллипса вокруг его малой оси.
Другими словами эллипсоид называется сфероидом.
Эллипсоид определяется параметрами:
1. Большой полуосью a и малой b
𝑎−𝑏
2. 𝛼 =
– сжатие
𝑎
3. Параметрами a и эксцентриситет e
𝑒=
√𝑎2 −𝑏2
𝑎
или 𝑒 ′ =
√𝑎2 −𝑏2
𝑏
Вопрос 1.2.
Для решения всех геодезических задач необходимо подобрать такие
параметры, при которых эллипсоид наилучшим образом описывает фигуру
реальной земли.
Существуют следующие методы определения параметра:
Геологический
Физический
Геологический – заключается в измерении двух меридианов. Из сравнении этих
измерений можно определить параметры. Такой метод называется методом дуг.
Вторым геологическим методом является метод площадей. В этом методе
определяются астрономические координаты пунктов геодезической сети
расположенных на большой территории. Далее задается приближение
параметрами земного эллипсоида а0, α0. По этим параметрам определяем
геодезические координаты (В0 и L0) пунктов этой сети.
В результате обработки этих данных надо получить поправку в эллипсоиде δа,
δα, δB, δL. Поправки в координаты пункты будет столько каким является
двойным числом пунктов. Чтобы упростить и решить коррекцию задачи
поправки δВi и δLi выражают через поправки в исходных пунктов δВi и δLi, δAi,
где Ai – азимут на исходном пункте. В результате будет только 5 неизвестных:
δВi, δLi, δAi, δa, δα.
Физический метод округления параметров – он заключен в измерении силы
тяжести на пунктах геодезической сети. По расхождениям значений силы тяжести
можно так же определить параметры земного эллипсоида. Например, сила
тяжести на полосе и на экваторе различаются сильно между собой. Можно по
разности таких измерений определить схожесть.
Геометрический метод развивается дополнением GPS измерений. На больших
расстояниях 600 и больше измерений приращения координат между пунктами. Их
тоже можно использовать для определения параметров земного эллипсоида.
Вопрос 2.
Под исходными геодезическими датами понимаем значение:
1. 2х параметров земного эллипсоида
2. Координат исходными пункта на эллипсоиде, а именно геодезической
широты и геодезической долготы и азимута, которым равны астрономические
координаты. Таким образом названные астрономических координат определяют
ориентировку эллипсоида в теле земли.
Вопрос 3.
Вообще раздел геодезия в котором решаются задачи на поверхности земного
эллипсоида называется сфероидоической геодезией. Все задачи необходимо
решать в определенной системе координат.
Различают следующие системы координат:
1. Прямоугольная геоцентричная система координат. В ней ось X направлена
в пересечение экватора и меридиана Гринвича. Ось Z совподает с осью
вращения земли. Начало в центре масс земли. Ось Y дополняет систему
правой.
2. Прямоугольная топоцентрическая. Ось Z направлена в зенит по нормали, X
– на север, Y – дополняет систему до правой на восток.
3. Геодезическая. Положение точки в геодезической системе координат
определяется геодезической широтой В (угол между нормалью и
плоскостью экватора) и долготой L и геодезической высотой H.
4. Геоцентрическая угловая. Положение точки А определяется: Ф –
геоцентрической широтой, геодезической долготой и высотой Н по радиус
вектору направленному с центра эллипсоида.
5. Система координат с приведенной широтой. Если через точку А провести
радиус вектор длинной «а» (большая полуось), таким образом чтобы он
косался оси вращения в точки С, то угол между этим радиусом вектором и
плоскостью экватора называется приведенной и обозначается U.
6. Система прямоугольных сферических координат. Одной из осей системы
является меридиан, а второй кривая линия перпендикулярная к меридиану в
точке О. Положение точки D будет определяться криволинейными
координатами p и q.
7. Плоские прямоугольные координаты. Примером является система координат
Гаусса –Крюгера.
Связь между системами координат.
a) Связь между геоцентричными координатами и геодезической широтой.
𝑑𝑧
𝑑𝑦
b) Связь между геодезичеой широтой и геоцентричной широтой.
𝑧
𝑡𝑔𝜑 =
𝑦
𝑑𝑧
с) Исходя из формулы 𝑡𝑔(90° − 𝐵) = и уравнения эллипса меридианного
𝑡𝑔(90° − 𝐵) =
𝑑𝑦
𝑦2
𝑧2
𝑎
𝑏
2𝑦𝑑𝑦
сечения. 2 + 2 = 1 найдем завсимость между геоцентричной и геодезической
широтами.
𝑎2
+
2𝑧𝑑𝑦
𝑏2
𝑎2
= 0;
𝑦𝑑𝑦
2𝑑𝑧
𝑦𝑑𝑦
𝑎
𝑏
2𝑑𝑧𝑎
=
2
,
2
=
2
−1 𝑑𝑦
𝑏
,
2
𝑑𝑧
𝑎2
𝑧
𝑏
𝑦
=− 2∗ ,
𝑐𝑡𝑔 (90° − 𝐵) = − 2 ∗ 𝑐𝑡𝑔Ф.
𝑏
d) связь между приведеной и геодизическами широтами.
Чтобы установить эту связь найдем чуму равен отрезок AD. Z = AD*sin U.
Y = AC*sin(90-U) = a*cos U. Подставим это значение в уравнение эллипса.
(𝐴𝐶∗sin(90°−𝑈))
𝑎2
+
(𝐴𝐷∗𝑠𝑖𝑛𝑈)2
𝑏2
=
𝑎2 𝑠𝑖𝑛2 𝑈
𝑎2
+
𝐴𝐷2 ∗𝑠𝑖𝑛2 𝑈
𝑏2
, следовательно AD=d.
С помощью b можно установить зависимость между геодезической и
приведенной широтой. Для вычисления dz и dy выразим z и y через
приведенную широту. Z = b*sin U, y = a*cos U, dz =b*cosUdU, dy = -a*sinUdU.
𝑏 𝑐𝑜𝑠𝑈𝑑𝑈
𝑏
tg(90°-B) =
= − 𝑐𝑡𝑔𝑈.
−𝑎𝑠𝑖𝑛𝑈𝑑𝑈
𝑎
Данная зависимость является базовой сферической геодезией при вычислении
длинны дуг и координатных точек. Так зная приведенную широту и
геодезическую долготу легко вычислить прямоугольные геоцентрические
координаты.
Очевидно отрезок r перейдет на n следовательно у = r*sin L. Для вычисления z
воспроизводят отрезок AD*b x = r*cos L, r = b*sin U. Для упрощения точку А
перенесем в плоскость чертежа
r = a * cos n, следовательно
y = a*cosU*sinL
x=a*cosU*cosL
z=b*sinU
Эти формулы являются базовыми вычисления геоцентричных координат по
геодезическим.
Вопрос 4.
При решении многих задач сферической геодезии используют
дифференциальные и интегральные исчисления. Поэтому необходимо знать
дифференциалы дуг вдоль меридиана и параллели.
Пусть имеется меридиан и параллель и элементарный отрезок ds. Обозначим
азимут, обозначим элементарные вдоль меридиана и параллели. Поскольку это
элементарный отрезок можно применить формулы тригонометрии.
dx = ds *cos A
dy = ds*sin A
С другой стороны можно выразить dx и dy через дифференциал долготы и
широты.
dx = M*dB , где М – радиус меридиана на элементе dx.
Выразим радиус параллели через радиус меридиана в первой вертикали.
r= N cosB, dy = rdL, dx = MdB, следовательно можно записать соотношение
между линейными элементами.
dx = ds*cosA=MdB. dy = ds*sinA=NcosBdL.
ds = √(𝑑𝑥)2 + (𝑑𝑦)2 = √(𝑀𝑑𝐵)2 + (𝑁𝑐𝑜𝑠𝐵𝑑𝐿)2
Эта формула называется первой квадратичной формулой поверхности
эллипсоида.
Вопрос 5.
dx и dy располагаются соответственно вдоль меридиана и параллели. Это
есть дифференциалы меридиана и параллели, а именно
dx=M·dB
dy=N·cosB·dL
На основе этих формул выведем длину дуги меридиана и параллели.
Длина дуги меридиана определяется интегрированием
B2
dx MdB
B1
L2
dy N cos BdL .
L1
Рассмотрим длину дуги меридиана.
В подынтегральном выражении значение М является радиусом меридиана в
плоскости меридиана. Это значение зависит от широты.
Запишем без вывода
M
a(1 l 2 )
.
2
2
3
( 1 l sin B )
Тогда длина дуги
B2
dB
.
x a(1 l 2 )
2
2
3
(
1
l
sin
B
)
B1
Обычно выполняют численное интегрирование, в результате которого
получают формулу для конкретных значений параметров земного эллипсоида.
например,
для
эллипсоида
Красовского
х=6367558,5·ВsinB·cosB·(32005,6+134,6·sin²B)/
Длина дуги параллели легко вычисляется
L2
L2
y N cos BdL N cos B dL N cos B( L2 L1 ) .
L1
L1
Вопрос 6.
dp=dx·dy
(1l 2 )
dxMdba
( 1l 2 sin 2 B )3
dy=NcosBdL
L2 B2 N cos B(1 l 2 )a
P
dBdL ,
2
2
3
(
1
l
sin
B
)
L1 B1
Без вывода
a cos B
N
,
1 l 2 sin 2 B
L2 B2 a 2 cos B(1 l 2 )
P
dBdL .
2
2
2
(
1
l
sin
B
)
L1 B1
Сделаем интегрирование по В с заменой переменной или расположением в
ряд подынтегрального выражения.
Поскольку имеется компьютерная техника, то можно выбрать любой
вариант. Мы возьмем интегрирование с заменой переменной.
esinB=sinψ
Тогда
ecosBdB=cosψdψ,
1
L2 B2 cosd
P a 2 (1 l 2 ) dL e
.
4
cos
L1 B1
Окончательно
B2 d
2 L2
2 1 e
.
Pa (
) dL
e L B cos3
1
1
1
cos 3 d является табличным интегралом.
1 e2
sin B
1 1 e sin B B2
P ( L2 L1 )a 2 (
)
ln
.
2
2
e 1 e sin B 2 1 e sin B B1
7. Геодезическая линия – это кратчайшее расстояние между точками на
эллипсоиде.
Геодезическая линия на плоскости – это отрезок.
Если на плоскости записывают уравнение прямой, то такое же уравнение
можно записать и для геодезической линии на эллипсоиде.
В основу вывода положим линейные элементы поверхности
dx=ds·cosA=M·dB,
dy=ds·sinA=N·cosB·dL.
Отсюда найдем cosA, sinA
MdB
cos A
,
dS
N cos BdL
sin A
,
dS
cos(90-B)=tgdL·ctg(90-dA),
sinB·tgdL=tgdA,
MdB
rdA cos A r
sin BdL .
dS
cosA умножим на величину rdA, а sinA – dr.
dB
cos ArdA M
rdA
dS
dL
sin Adr N cos B dr
dS
rdA sin BdL
dB
cos ArdA N
sin BdL
dS
Для простоты выводов вернемся к радиусу параллели
NcosB=r
dL
dL
sin Adr r dr (MdB sin B)
(1)
dS
dS
rdr MdB sin B
Правые части двух записанных равенств 1 равны по модулю и отличаются
лишь по знаку.
Тогда сумма этих выражений равняется нулю.
cos ArdA sin Adr 0
r cos AdA sin Adr 0
Это есть производная следующего выражения
(rsinA)'=0
rsinA=const.
Это и есть уравнение геодезической линии.
Тема 7. Решение сферических треугольников
1. Общие сведения о решении сферических треугольников. Порядок
решения по способу Лежандра
2. Способ аддитаментов
Вопрос 1.
Суть способа Лежандра заключается в том, чтобы зная длину исходной
стороны по формулам плоской тригонометрии вычислить остальные стороны.
Рассмотрим этот способ
Из сферической тригонометрии известно, что сумма безошибочно известных
углов сферического треугольника равна 180+ε, где ε – сферический избыток.
Без вывода приведем формулы вычисления сферического избытка
a
b
sin 2 sin 2 sin c
sin
,
c
2
cos
2
где соответствующие величины представлены на чертеже
Углы А, В, С измерены, а а, b, с – длины сторон, выраженные в угловой мере.
Решение треугольников на сфере осуществляется по формуле
sin a sin b sin c
.
sin A sin B sin C
Если sina известно, то
sin a sin B
.
sin b
sin A
По синусу найти значение дуги в градусной мере, а потом, зная радиус шара
в данной точке, определяем значение дуги в линейной мере.
R·b=Sb.
Раньше когда вычисления осуществляются по таблицам использовались
формулы геометрии. Для этого из каждого угла измеряемого на эллипсоиде
вычислялась величина, равна 1/3 сферического избытка.
Если задана исходная сторона Sa и исправленные за сферический избыток
углы
А А
3
B B
C C
3
3
Записывают теорему синусов плоской геометрии
Sa
S
b ,
sin A sin B
S sin B
.
Sb a
sin A
В заключении покажем простейший вывод формулы сферического избытка
sin
а
2
sin sin
b
2 sin c
c
2
Поскольку сферический избыток малая величина то его можно заменить
2
величиной
cos
2
2
а
2
sin sin
cos
или
c
2
b
2 sin c
аb
sin c
2
Он получается в радианах.
И чтобы перевести в секунды
а
Sa
S
;b b
R
R
S a Sb
sin c
2R 2
И тогда
или
2R 2
S a Sb sin c
Иногда дуги заменяются хордами. Но в этом случае длина каждой хорды
вычисляется так
b
2 sin b
R
2
b 2 R sin
b
2
Но тогда углы исправляются на 1/4 избытка
А А
В В
С С
4
4
4
а sin B
b
sin A
Если длина дуги b измерена, то
b
Sb
R
а
Sa
R
А если а
тогда
а 2R sin
Sa
2
Вопрос 2.
Аддитамент (от англ. слова – «добавить»)в этом случае вводятся поправки не
в углы как в способе Лежандра, а в стороны, и решение выполняют по формулам
плоского треугольника.
Запишем в качестве исходной формулы, формулу вычисления сторон на
сфере
sin b
sin b sin a
sin a
Sa
R
S
S sin B
sin b sin a
R
R sin A
а
С тем чтобы решалось по правилам плоской тригонометрии sin сторон
разлаживают в ряд Тейлора ограничиваясь первыми членами разложения
3
1 Sa
R
2
S
S
S
S3
sin a a a a3
R
R
R 6R
3
3
Sb S 3b Sa S a sin B
R 6 R 3 R 6 R3 sin A
Можно считать, что приведенные члены в 3-й степени являются какими-то
поправками.
Тема 8. Решение главных геодезических задач на эллипсоиде.
1.Общие сведенья о главных геодезических задачах на плоскости
эллипсоида.
2. Пути и методы решения главных геодезических задач.
3. О точности решения главных геодезических задач.
4. Вывод формул путем разложения в ряд широт, долгот и азимутов.
5. Решение обратной геодезической задачи со средними аргументами.
6. Решение прямой геодезической задачи методом вспомогательной
точки.
7. Решение обратной геодезической задачи.
8. О способе Бесселя решения главных геодезических задач.
9. Переход на эллипсоид способа Бесселя.
10. Дифференциальные уравнения способа Бесселя.
11. Численный метод решения прямых геодезических задач.
Вопрос 1.
В геодезии на плоскости существует 2 главные геодезические задачи: прямая
и обратная.
В прямой геодезической задаче даны координаты Х У первой точки,
дирекционный угол направления с точки 1 на точку 2. Необходимо определить
координаты Х У второй точки.
В обратной геодезической задаче по координатам конечных точек отрезка
или прямой необходимо вычислить ее длину S и дирекционный угол α.
На эллипсоиде в сфероедической геодезии также рассматривается 2
геодезические задачи: прямая и обратная.
Прямая геодезическая задача
Даны B1 L1 – широта и долгота точки А, S – длина геодезической линии, А12
– азимут линии.
Необходимо определить B2 L2 – широту и долготу точки В и обратный
азимут А21.
Обратная геодезическая задача
Даны B1 L1, B2 L2. Необходимо найти S – длину геодезической линии и
азимуты А12, А21.
Вопрос 2.
Решение главных геодезических задач на эллипсоиде зависит в основном от
длины геодезической линии.
Длина геодезической линии бывает:
1) малой – до 45 км;
2) средней – до 600 км;
3) большой – до 5000 км;
4) сверхбольшой – до 20000 км.
Решение главных геодезических задач на малые и средние расстояния
осуществляется способами Гаусса и Шрейбера.
Решение на большие расстояния осуществляется способом Бесселя.
Сверхбольшие – численный метод.
Все эти задачи для всех расстояний могут решаться следующими путями:
- прямым
- косвенным.
Вопрос 3.
Все задачи решаются вне зависимости от длины геодезической линии с
высокой точностью. Так широты и долготы вычисляют с точностью 0,0001'',
азимуты – 0,001''.
Длина S должна быть задана с точностью до 1 мм.
Вопрос 4.
Прямая геодезическая задача
Можно записать, что B2=B1+ΔB=B1+(B2-B1)=B1+b
L2=L1+ΔL=L1+(L2-L1)=L1+l
А21= А12±180°+α.
Необходимо найти b, l, α.
Можно отметить, что увеличение широты, долготы и азимута
осуществляется с увеличением длины геодезической линии.
Если знать значение приращения одной из координат в зависимости от
B L A
B L A
длины геодезической линии
или производные
,
,
, , , то в
S S S
S S S
соответствии с разложением в ряд Тейлора можно записать
i B S i
B2 B1 b
S i i!
i B
- производная широта по S i-того порядка.
S i
Аналогичные выражения запишем для долгот и азимута
i L S i
L2 L1
S i i!
A21 A12 180 0
i A S i
.
S i i!
Исходя из этого чертежа, установили зависимость между дифференциалами
dB=dScosA,
dL=dSsinA.
Установим зависимость между dS и dА.
Очевидно, что дугу dB можно вычислить дважды
dB= dL=dSsinA,
а второй раз dА, задаваясь радиусом первого вертикала N (радиусом
эллипсоида в точке DB).
DB=N·dA.
Мы знаем, что N можно выразить через радиус параллели.
Очевидно
r
.
N
cos B
Но такая запись в дальнейших преобразованиях неудобна.
DB TD dA
r
cos( 90 B)
TD
r
r
TD
cos( 90 B) sin B
Таким образом
r
DB TD dA
dA .
sin B
Очевидно
r
dS sin A
dA ,
sin D
dA sin A sin A sin B
,
r
dS
r
sin B
dB
Для определения
выразим dB через дугу, приращение широты и радиус
dS
меридиана.
ED=dS·cosA
ED=M·dB,
где
М
–
радиус
кривизны.
MdB=dScosA
dB cos A
dS
M
DB=dL=dSsinA
Поскольку dL есть
приращение долготы, то
его удобно вычислить через
радиус параллели
DB rdL
dS sin A rdL
dL sin A
dS
r
Сведем эти производные вместе
dB cos A
dS
M
dL sin A
dS
r
dA sin A sin B
dS
r
Для того чтобы можно дифференцировать значение H и r выразим через
широту.
c
Поскольку M ,
v
2
a
c
где
b
v 1 c 2 cos 2 B
r N cos B
c
N
v
dB v
cos A
dS c
dL
v
sin A
dS c cos B
dA v sin A sin B
dS
c cos B
Все производные зависят от 2-х аргументов B и А.
Очевидно, что последующие производные необходимо брать по двум
dB
аргументам В и А. посмотрим на примере
dS
dB
dB
d 2 B d ( dS ) dB d ( dS ) dA
.
dB dS
dA dS
dS 2
На основе этого примера можно составить рекуррентные формулы взятия
производных, а именно
d i 1 B
d i 1 B
d
(
)
d
(
)
i 1
i 1
diB
dB
dA
dS
dS
.
dB
dS
dA
dS
dS i
d i 1 B
d i 1 L
Остальные значения определяются также, только вместо
берется
dS i 1
dS i 1
d i 1 A
и
.
dS i 1
Полученные производные подставляют в суммы для вычисления b, l и α. И
этим решается прямая геодезическая задача.
Вопрос 5.
Решение обратной геодезической задачи со средними аргументами неразрывно
связано с решением прямой геодезической задачей со средними аргументами.
Такое решение называется решением способом Гаусса..
Пусть известна точка D с координатами B;L;A12 , дана длинна геодезической
линии DQ – S. Необходимо определить координаты точки Q B2; L2; A21. Решение.
По Гауссу выбирается строго по середине геодезической линии DQ и условно
вводятся координаты этой точки Bm; Lm; Am.
Записываем выражение разности широт B1 – Bm, B2 – Bm. Аналогично
предыдущему случаю будем записывать B1 – Bm*
𝑑𝐵
𝑑𝑆
𝑆
∗ . Поскольку значение B1
2
меньше чем значение Bm то перед первой производной ставится знак минус.
𝑑𝐵 𝑆 𝑑 2 𝐵 𝑆 2 1 𝑑 3 𝐵 𝑆 3 1
𝐵1 − 𝐵𝑚 = −
∗ +
∗( ) ∗ − 3 ∗( ) −
𝑑𝑆 2 𝑑𝑆 2
2
2! 𝑑𝑆
2
3!
𝑑𝐵 𝑆 𝑑 2 𝐵 𝑆 2 𝑑 3 𝐵 𝑆 3
𝐵1 − 𝐵𝑚 = −
∗ +
∗ −
∗
𝑑𝑆 2 𝑑𝑆 2 8 𝑑𝑆 3 48
𝑑𝐵 𝑆 𝑑 2 𝐵 𝑆 2 𝑑 3 𝐵 𝑆 3
𝐵2 − 𝐵𝑚 =
∗ +
∗ +
∗
𝑑𝑆 2 𝑑𝑆 2 8 𝑑𝑆 3 48
Из этих выражений найдем величину B2 = B1+(B2-B1) очевидно что
𝑑𝐵
𝑑2𝐵
𝑆3
𝐵2 − 𝐵1 = ( ) ∗ 𝑆 + ( 2 ) ∗
𝑑𝑆 𝑚
𝑑𝑆 𝑚 48
Особенность данной формулы является то, что производные берутся в точке М.
Для вычисления производных берутся широта, долгота и азимут геодезической
линии в точке М. Но в реальности мы не можем это сделать, т.к координаты точки
М неизвестны.
Гауссом предложено вместо координат точки М брать среднее координат точки
D и Q.𝐵 =
𝐵1 −𝐵2
2
;𝐿=
𝐿1 −𝐿2
2
;𝐴 =
𝐴12 −𝐴21
2
.
Но эти формулы тоже не известны, следовательно задачу решают методом
приближений. В нулевом приближении задаются грубыми координатами точки Q
(B20, L20, A21) с точностью до 1°. Выполняют первое приближение и вычисляют
первое приближение B2, L2, A21. Во втором приближении подставляют эти
значения, в третьем значения со второго приближения и т.д. Приближения
повторяются до тех пор, пока результат предыдущего и последнего приближения
не будут одинаковыми.
𝑑3𝐵
𝑑𝐵
𝑆3
Замечание: в выражении 𝐵2 − 𝐵1 = ( ) ∗ 𝑆 + ( 3 ) ∗ производные даны в
𝑑𝑆
𝑑𝑆
24
𝑚
𝑚
точке M, а их выражения необходимы в точке В с координатами B, L, A. Эта
точка находится на каком-то расстоянии от M до M`. Расстояние между точками
M` и M в координатах будет B – Bm, L – Lm, A – Am, поэтому функцию
𝑑𝐵
𝑑𝑆
разлаживаем в ряд Тэйлора.
Поскольку рассмотренное произведение широты по расстоянию, то она будет
зависеть от самой широты. Азимута и геодезической линии. Тогда разложение
вряд Тэйлора будет иметь вид:
𝑑𝐵
𝑑𝐵
( 𝑑𝑆 ) = 𝑑𝑆 +
𝑚
𝑑(𝑑𝐵⁄𝑑𝑆)
𝑑(𝑑𝐵⁄𝑑𝑆)
∗ (𝐵𝑚 − 𝐵) +
∗ (𝐴𝑚 − 𝐴)
𝑑𝐵
𝑑𝐵
𝑑3𝐵
𝑑𝐵
𝑆3
Тогда можно записать 𝐵2 − 𝐵1 = ( ) ∗ 𝑆 + ( 3 ) ∗
𝑑𝑆
𝑑𝑆
24
𝑚
Но для вычисления самого
𝑑𝐵
𝑑𝑆
𝑚
необходимо знать Bm-B, Am-A.
Найдем Bm-B исходя из записанных разностей широт относительно точки
М.Чтобы найти Bm-B нужно
знак минус. 𝐵𝑚 − 𝐵 = −
𝑑𝐵
записать 𝐵2 − 𝐵1 = ((
𝑑𝑆
(𝐵1 −𝐵𝑚 )+(𝐵2 −𝐵𝑚 )
2
𝑑2𝐵
𝑆2
𝑑𝑆
8
+
∗
2
и вычеслить его присвоив ему
. Тогда окончательно для B2-B1 можно
𝑑𝐵
𝑑2𝐵
𝑆2
) /𝑑𝐵) ∗ ( 𝑑𝑆 2 ∗ 8 ) следовательно 𝐵2` = 𝐵20 +
𝑑𝑆
(𝐵2 − 𝐵1 )0 ; 𝐵2`` = 𝐵2` + (𝐵2 − 𝐵1 )` и т.д.
Аналогичное приближение по L и A:
𝐿𝑛2 = 𝐿𝑛−1
+ (𝐿2 − 𝐿1 )𝑛−1
2
𝐴𝑛2 = 𝐴𝑛−1
+ (𝐴2 − 𝐴1 )𝑛−1
2
Вопрос 6.
Способ Шрейдора.
Дано B1, L1,A1,S.
Найти B2,L2,A.
Решение.
Берется вспомогательная точка на меридиане, делается прямое сечение из точки С
в точку Q, таким образом чтобы с точки С углы были прямые.
Тогда Шрэйдором было предложено вычислить координаты точки С
относительно точки D, а потом относительно точки С вычеслить координату
точки Q.
В основу решения положены известные формулы:
𝑑𝐵
𝑑2𝐵 𝑆 2 𝑑3𝐵 𝑆 3
𝑏 = 𝐵2 − 𝐵1 =
∗𝑆+ 2 ∗ + 3 ∗ +⋯
𝑑𝑆
𝑑𝑆
2! 𝑑𝑆
3!
𝑑𝐿
𝑑2𝐿 𝑆 2 𝑑3𝐿 𝑆 3
𝑙 = 𝐿2 − 𝐿1 =
∗𝑆+ 2∗ + 3∗ +⋯
𝑑𝑆
𝑑𝑆
2! 𝑑𝑆
3!
𝑑𝐴
𝑑2𝐴 𝑆 2 𝑑3𝐴 𝑆 3
𝛼 = 𝐴2 − 𝐴1 =
∗𝑆+ 2 ∗ + 3 ∗ +⋯
𝑑𝑆
𝑑𝑆
2! 𝑑𝑆
3!
Прежде чем приступить к алгоритму решения сделаем некоторые общие
земечания:
Обозначим дугу DC – х, CQ – y. Эти дуги нужны для того чтобы вычислить
координаты точек C и Q.
Решение будет на основании треугольника DCQ. Запишем теорему синусов.
𝑥
𝑦
𝑆
=
=
𝜀
𝜀
𝜀
sin (𝛽 − 3) sin (𝐴1 − 3) sin (900 − 3)
x,y,s – геодезические линии.
𝑥=
𝜀
3
𝜀
sin(900 −3)
𝑆∙sin(𝛽− )
;
𝑦=
𝜀
3
𝜀
sin(900 −3)
𝑆∙sin(𝐴1 − )
Угол β находится из треугольника.
1800 + 𝜀 = 900 + 𝛽 + 𝐴1 => 𝛽 = 1800 − 900 + 𝜀 − 𝐴1 = 900 + 𝜀 − 𝐴1
Можем записать что 𝑥 =
𝜀
3
𝑆∙cos(𝜀−𝐴1 − )
𝜀
sin(900 −3)
=
𝜀
3
𝑆∙sin(900 +𝜀−𝐴1 − )
𝜀
sin(900 −3)
𝜀
𝜀
3
3
В дальнейших выводах будем полагать что sin (90° − ) = cos .
𝑥=
2
3
𝑆∙cos( 𝑆−𝐴1 )
𝜀
3
cos
; 𝑦=
2
3
𝑆∙sin(𝐴1 − )
cos
𝜀
3
Для упрощения выражений разложим cos и sin разности углов через разность
произведений. cos(𝛼 ± 𝛽) = 𝑐𝑜𝑠𝛼 ∗ 𝑐𝑜𝑠𝛽 ∓ 𝑠𝑖𝑛𝛼 ∗ 𝑠𝑖𝑛𝛽
sin(𝛼 ± 𝛽) = 𝑠𝑖𝑛𝛼 ∗ 𝑐𝑜𝑠𝛽 ± 𝑐𝑜𝑠𝛼 ∗ 𝑠𝑖𝑛𝛽
cos (
𝑥=
𝑆
𝜀
cos
3
𝜀
2𝜀
2𝜀
− 𝐴1 ) = cos(𝐴1 − )
3
3
2
𝑆
3
𝜀
cos
3
∙ (cos 𝐴1 cos 2 + sin 𝐴1 sin 𝜀) ; 𝑦 =
3
𝜀
𝜀
3
3
∙ (sin 𝐴1 cos − sin cos 𝐴1 )
Дальнейшее решение можно проводить при этих значениях x и y. Переход к
координатам точки С.
1.Найдем разность широт точки С и точки D. В данном случае S=x, A2=A1
следовательно 𝑥 =
𝑆
𝜀
cos
3
2
2
3
3
∙ (cos ( 𝜀) cos 𝐴1 + sin 𝐴1 sin ( 𝜀))
Далее алгоритм вычислений соответствует выше описанному, а именно значение
b вычисленного по формуле. 𝑙 = 𝐵𝑐 − 𝐵𝑝 =
𝑑𝐵𝑥
𝑑𝑆
𝑥+
𝑑2𝑆
𝑥2
𝑑𝑆
2!
∙
2
+⋯
Lc=L1, A = 0.
2. Вычисление координат точки Q относительно точки С.
S=y, ACQ=90°.
Снова используем для вычисления известные величины l,b,α.
Особенностью вычислений l,b,α заключается в значении производных
Поскольку А=0, то
𝑑𝐵
𝑑𝑆
𝑉
𝑑𝑙
𝐶
𝑑𝑆
= ;
= 0;
𝑑𝐴
𝑑𝑆
𝑑𝐵
𝑑𝑆
, 𝑑𝑆.
=0.
Дуга в точке Q параллельна меридиану DP – называется параллелью этого
меридиана.
Вопрос 7.
Из всех способов решения обратной геодезической задачи самой
распространенной является способ Щрейбера.
Даны координаты B1, L1; B2, L2. Необходимо вычислить длину геодезической
линии и азимута A1 и А2.
В основу решения положена формула b, l, α. В каждой этой формуле
присутствуют значения S только в различной степени. Левая часть b известна.
Покажем это на примере разности широт.
𝑑𝐵
𝑑2𝐵 𝑆 2 𝑑3𝐵 𝑆 3
𝑏 = 𝐵2 − 𝐵1 =
∗𝑆+ 2 ∗ + 3 ∗ +⋯
𝑑𝑆
𝑑𝑆
2! 𝑑𝑆
3!
Левая часть известна. Можно найти производные
𝑑𝐵
𝑑𝑆
,
𝑑2𝐵
𝑑𝑆 2
и т.д и решить
настоящее уравнение. Значение производной берется в точке 𝐵 =
𝐵1 +𝐵2
2
. Если
ограничится разложением до второй степени. Можно решить уравнение второго
порядка и найти производные S. Если ограничить уравнение n-й степени, то
необходимо решить уравнение n-й степени и найти S.
Приближенность решения определяется тем, что производная в данном ряде в
точке В, поэтому данный способ пригоден лишь для малых расстояний, не больше
40 км. Точность определения длинны составляет порядка
0,008 мм *40 000 000мм/206 000 = 1,5 – точность определения длинны. Где 0,008
– ошибка определения азимута.
Для вычисления азимутов можно использовать дифференциальные формулы:
𝑑𝐿 𝑑𝐴 𝑑𝐵
,
.
𝑑𝑆 𝑑𝑆 𝑑𝑆
𝑣
=
𝑐
∗ cos 𝐴,
𝑑𝐿
𝑑𝑆
Очевидно, что разделив
=
𝑑𝐵
𝑑𝑆
𝑣∗𝑠𝑖𝑛𝐴
𝑐∗𝑐𝑜𝑠𝐵
/
𝑑𝐿
𝑑𝑆
𝑑𝐵
𝑑𝑆
,
.
можно найти
𝑑𝐿
tan 𝐴1 = (𝑑𝑆⁄𝑑𝐵) ∙ cos 𝛽
𝑑𝑆
Чтобы получилась замкнутая формула снова необходимо воспользоваться рядами
𝑑𝐿
1
𝑑2𝐿
𝑆2
𝑑3𝐿
𝑆3
𝑑𝐵
1
𝑑2𝐵
𝑆2
𝑑3𝐵 𝑆 3
= (𝑙 − 2 ∙ − 3 ∙ )
𝑑𝑆
𝑆
𝑑𝑆
2!
𝑑𝑆
3!
= (𝑏 − 2 ∙ − 3 ∙ )
𝑑𝑆
𝑆
𝑑𝑆
2!
𝑑𝑆
3!
𝛼=
𝑑𝐴
𝑑𝑆
∙𝑆+
𝑑2𝐴 𝑆 2
𝑑𝑆
∙
2!
+⋯
Производные берутся в начальной точки.
Вопрос 8.
Этот способ предназначен для решения на любые расстояния. В начале решается
прямая и обратная задачи на шаре, чтобы при введение поправок получилось
решение на сфероиде.
Решение главной геодезической задачи на шаре.
Дано: φ1, σ, λ1, α1.
Найти: φ2, λ2, α2.
Решение:
Для вычисления широты φ2 записывается формула косинуса стороны
сферического треугольника Q1PQ2. Для этого используется теория геометрии на
сфере.
cos(90° − 𝜑2 ) = cos(90° − 𝜑1 ) ∗ 𝑐𝑜𝑠𝜎 + sin(90° − 𝜑1 ) ∗ 𝑠𝑖𝑛𝜎 ∗ 𝑐𝑜𝑠𝜎𝑚
В правой части все параметры известны. Чтобы найти долготу λ2= λ1+λ.
Для вычисления разности долгот записывают 2 теории.
1. Теорема синусов по элементам λ, σ, α, (90°-φ2)
2. Теорема пяти элементов (90°-φ2), λ, α1, (90°-φ1), σ.
Теорема синусов.
sin 𝜆
sin 𝛼1
=
sin 𝜎 sin(90° − 𝜑2 )
sin 𝜆 ∗ 𝑐𝑜𝑠𝜑2 = sin 𝜎 ∗ sin 𝛼
Теорема пяти элементов.
sin(90° − 𝜑2 ) 𝑐𝑜𝑠𝜆 = sin(90° − 𝜑1 ) 𝑐𝑜𝑠𝜎 − 𝑠𝑖𝑛𝛿 cos(90° − 𝜑1 ) 𝑐𝑜𝑠𝛼1
𝑐𝑜𝑠𝜑2 𝑐𝑜𝑠𝜆 = 𝑐𝑜𝑠𝜑2 𝑐𝑜𝑠𝜎 − 𝑠𝑖𝑛𝜎𝑠𝑖𝑛𝜑1 𝑐𝑜𝑠𝛼1
Если выражение из теоремы синусов разделить на формулу пяти элементов.
𝑡𝑔𝜆 =
𝑠𝑖𝑛𝜎 ∗ 𝑠𝑖𝑛𝛼1
sin(90° − 𝜑1 ) 𝑐𝑜𝑠𝜎 − 𝑠𝑖𝑛𝜎 cos(90° − 𝜑1 ) 𝑐𝑜𝑠𝛼1
Модуль второй части находится тоже из двух формул.
1. Формула синусов относительно элементов: α2, (90°-φ1), σ1, (90°-φ2).
2. Формула пяти элементов относительно: (90°-φ2), α2, σ, α1, (90°-φ1).
Формула синусов.
− sin 𝛼2
sin 𝛼1
=
sin(90° − 𝜑1 ) sin(90° − 𝜑2 )
sin 𝛼2 ∗ sin(90° − 𝜑2 ) = 𝑠𝑖𝑛𝛼1 ∗ sin(90° − 𝜑1 )
Формула пяти элементов.
sin(90° − 𝜑2 ) 𝑐𝑜𝑠𝛼2 = 𝑠𝑖𝑛𝜎 ∗ cos(90° − 𝜑1 ) − 𝑐𝑜𝑠𝜎 ∗ sin(90° − 𝜑1 ) 𝑐𝑜𝑠𝛼1
Разделить выражение формулы синусов на формулу пяти элементов.
𝑡𝑔𝛼2 =
𝑠𝑖𝑛𝛼1 ∗ 𝑐𝑜𝑠𝜑1
𝑠𝑖𝑛𝜎 ∗ cos(90° − 𝜑1 ) − 𝑐𝑜𝑠𝜎 ∗ sin(90° − 𝜑1 ) 𝑐𝑜𝑠𝛼1
Обратная геодезическая задача на шаре.
Дано φ1 λ1, φ2 λ2.
Найти α1 α2 σ
Решение.
Для решения задачи тоже используют формулы сферической тригонометрии. Для
вычисления угла α записываем снова 2 формулы.
1. Теорема синусов относительных элементов λ2 σ α1 (90°-φ2)
2. Формула пяти элементов: σ α1 (90°-φ2) (90°-φ1) λ
Теорема синусов.
𝑠𝑖𝑛𝜎 sin(90° − 𝜑2 )
=
𝑠𝑖𝑛𝜆
𝑠𝑖𝑛𝛼1
𝑠𝑖𝑛𝜎 ∗ 𝑠𝑖𝑛𝛼1 = 𝑠𝑖𝑛𝜆 ∗ sin(90° − 𝜑2 )
Формула пяти элементов.
𝑠𝑖𝑛𝜎 ∗ 𝑐𝑜𝑠𝛼1 = sin(90° − 𝜑1 ) ∗ cos(90° − 𝜑2 ) − sin(90° − 𝜑2 ) ∗ cos(90° − 𝜑1 ) ∗ 𝑐𝑜𝑠𝛼
Разделив эти формулы получим:
𝑡𝑔𝛼 =
𝑠𝑖𝑛𝐴 ∗ sin(90° − 𝜑2 )
sin(90° − 𝜑1 ) cos(90° − 𝜑2 ) − sin(90° − 𝜑2 ) cos(90° − 𝜑1 ) 𝑐𝑜𝑠𝜆1
Для вычисления α2 запишем:
1. Теорему синусов: α2 (90°-φ2) λ σ
2. Формула пяти элементов σ α2 (90°-φ1) (90°-φ2) λ
Теорема синусов.
𝑠𝑖𝑛𝜎 sin(90° − 𝜑2 )
=
𝑠𝑖𝑛𝜆
𝑠𝑖𝑛𝛼2
𝑠𝑖𝑛𝜎 ∗ 𝑠𝑖𝑛𝛼1 = sin(90° − 𝜑2 ) ∗ 𝑠𝑖𝑛𝜆
Формула пяти элементов.
𝑠𝑖𝑛𝜎 ∗ 𝑠𝑖𝑛𝛼2 = sin(90° − 𝜑2 ) cos(90° − 𝜑1 ) − sin(90° − 𝜑1 ) cos(90° − 𝜑1 ) 𝑐𝑜𝑠𝜎
Разделив второе уравнение выравнивания на первое получим.
𝑡𝑔𝛼0 =
𝑠𝑖𝑛𝛼 ∗ sin(90° − 𝜑1 )
sin(90° − 𝜑2 ) cos(90° − 𝜑1 ) − sin(90° − 𝜑1 ) cos(90° − 𝜑2 ) ∗ 𝑐𝑜𝑠𝜆
σ можно вычеслить из любой формулы.
Вопрос 9.
В способе Бесселя переход на эллипсоид осуществляется при следующих
предположениях.
1. σ =S где σ – дуга большого круга на шаре
S – геодезическая линия
Т.е когда решают прямую геодезическую задачу на шаре, то берут σ измеренное
на эллипсоиде.
2. Широты на шаре (φ) принимаются равными приведенным широтам на
эллипсоиде.
Исходные данные B1 L1 , то по формулам связи приведенных широт и
геодезических вычисляют элемент.
𝑠𝑖𝑛𝑈1 =
𝑠𝑖𝑛𝐵1
𝑣
или 𝑐𝑜𝑠𝑈1 =
𝑐𝑜𝑠𝐵1
𝑤
. Полагаем что φ1 = U1.
2. Азимут линии на шаре равен азимуту на эллипсоиде, т.е B1 L1 A1 –
исходные данные. α1 = А1 Таким образом если даны B1 L1 A1 S, то от В
переходим к 𝑈1 , следовательно φ1 = U1, α1 = А1, σ =S и приступает к
решению задач на сфере. При этих предположениях из некоторых формул
пяти элементов уже можно найти некоторые решения на эллипсоиде.
Так записывая теорему косинусов для элемента 90° − 𝜑2
cos(90° −
𝜑2 ) = cos(90° − 𝜑1 ) 𝑐𝑜𝑠𝜎 + sin(90° − 𝜑1 ) 𝑠𝑖𝑛𝜎𝑐𝑜𝑠𝛼1
И учитывается предположение Бесселя найдем решение для широты
cos(90° − 𝑈2 ) = cos(90° − 𝑈1 ) 𝑐𝑜𝑠𝑆 + sin(90° − 𝑈1 ) 𝑠𝑖𝑛𝑆𝑐𝑜𝑠𝐴1
Вторым элементом который можно легко получить прямым путем является
азимут А2. Для этого записываем формулу tgα2 при решении прямой
геодезической задачи на шаре.
𝑠𝑖𝑛𝛼sin(90° − φ1 )
𝑠𝑖𝑛𝛼𝑐𝑜𝑠(90° − 𝜑1 ) − 𝑐𝑜𝑠𝜎 sin(90° − 𝜑1 ) 𝑐𝑜𝑠𝛼1
𝑠𝑖𝑛𝐴1 ∗ sin(90° − 𝑈1 )
𝑡𝑔𝐴2 =
𝑠𝑖𝑛𝑆 ∗ cos(90° − 𝑈1 ) − 𝑐𝑜𝑠𝑆 ∗ sin(90° − 𝑈1 ) 𝑐𝑜𝑠𝐴1
𝑡𝑔𝛼2 =
Исходя из постановления геодезической задачи надо найти L2 = L1 + l.
На шаре находим λ, следовательно нужно l выразить через λ. Эта задача решается
с помощью дифференциальных уравнений способа Бесселя.
Вопрос 10.
l=f(λ). Необходимо вырезать dL=f(dx). Необходимо выразить дифференциалы
долготы на шаре. Попутно найдем все дифференциалы и выразим
дифференциалы всех элементов на эллипсоиде через диференциалы элементов на
шаре.
Другими словами задача заключается в установлении зависимости
𝑑𝐵 𝑑𝐿 𝑑𝐴
;
;
𝑑𝜑 𝑑𝜆 𝑑𝛼
.
Запишем дифференциальные уравнения для эллипсоида:
dx=dS*cosA; dy = dS*sinA. В этих уравнениях нет еще дифференциалов широт и
долгот, они появятся в следующих формулах: dx=MdB; dy=rdL. Приравняв
правые и левые равенства в результате получится MdB=dScosA; rdL=dSsinA,
следовательно 𝑑𝐵 =
𝑑𝑆𝑐𝑜𝑠𝐴
𝑀
; 𝑑𝐿 =
𝑑𝑆𝑠𝑖𝑛𝐴
𝑟
. здесь М - радиус меридиана в точке на
отрезке, r – радиус параллели, которую выражают через радиус первого
вертикала. r = NcosB. Следовательно 𝑑𝐿 =
Без вывода известно что dA=dLsinB, dB =
𝑠𝑖𝑛𝐴
𝑀𝑐𝑜𝑠𝐵
𝑐𝑜𝑠𝐴
𝑀
𝑑𝑆.
𝑑𝑆.
Запишем дифференциалы на шаре по аналогии с эллипсоидом.
dφ; dλ; dα
dφ = cosα*dσ; dλ =
𝑠𝑖𝑛𝐴
𝑐𝑜𝑠𝜑
𝑑𝜎; dα=sinσ*tgφ*dσ.
Решим главную геодезическую задачу (выразим dL через dλ)
𝑑𝐿
𝑠𝑖𝑛𝐴 ∗ 𝑐𝑜𝑠𝜑 𝑑𝑆
=
∗
𝑑𝜆 𝑁𝑐𝑜𝑠𝐵 ∗ 𝑠𝑖𝑛𝛼 𝑑𝜎
В учебниках по высшей геодезии находят
𝑑𝐴 =
𝑑𝐿
𝑑𝜆
=
𝑠𝑖𝑛𝐴 𝑡𝑔𝐵
𝑛∗𝑠𝑖𝑛𝛼∗𝑡𝑔𝜑
∗
𝑑𝐿
𝑑𝜆
𝑠𝑖𝑛𝐴 𝑡𝑔𝐵
∗ 𝑑𝑆
𝑛
𝑑𝑆
𝑑𝜎
По условию Бесселя можно прнять dA=dα, тогда
1=
Тогда
𝑑𝑆
𝑑𝜎
=
𝑁∗𝑠𝑖𝑛𝛼∗𝑡𝑔𝜑
𝑠𝑖𝑛𝐴 𝑡𝑔𝐵
Без вывода запишем
1
dL= 𝑑𝜆
𝑣
-
𝑠𝑖𝑛𝐴 𝑡𝑔𝐵
𝑑𝑆
∗
𝑁 ∗ 𝑠𝑖𝑛𝛼 ∗ 𝑡𝑔𝜑 𝑑𝜎
.
𝑑𝐿
𝑑𝜆
=
𝑄2 1
dL = ∫𝑄1
𝑣
𝑠𝑖𝑛𝜑 𝑑𝐿
;
𝑠𝑖𝑛𝐵 𝑑𝜆
𝑑𝜆
=
𝑠𝑖𝑛𝜑
𝑠𝑖𝑛𝑈∗𝑣
=
1
𝑣
U=φ по условию Бесселя.
Вопрос 11.
Численный метод решения прямых геодезических задач
Из предыдущего материала нами были выведены формулы зависимости
приращений широт от приращений длинны геодезической линии, другими
словами, были установлены значения производных
𝑑𝐵 𝑑𝐵 𝑑𝐵
; ; .
𝑑𝑆 𝑑𝑆 𝑑𝑆
𝑑𝐵 𝑉 3
=
𝑐𝑜𝑠𝐴
𝑑𝑆
𝐶
𝑑𝐿 𝑉 𝑠𝑖𝑛𝐴
=
𝑑𝑆 𝐶 𝑐𝑜𝑠𝐵
𝑑𝐴 𝑉 𝑠𝑖𝑛𝐴𝑠𝑖𝑛𝐵
=
𝑑𝑆 𝐶 𝑐𝑜𝑠𝐵
Эти уравнения являются дифференциальными уравнениями 1 порядка. Для них
можно записать общий вид:
𝑑𝑦
= 𝑓(𝑦)
𝑑𝑥
Настоящее дифференциальное уравнение решается любым известным в
математике методом. В данном случае мы
рассмотрим решение методом Рунге-Кутта.
В данном методе ось x разбивается на интервалы с
шагом h. Приближённый вариант решения будет
следующий:
𝑦𝑖+1 − 𝑦𝑖
= 𝐻(𝑦𝑖 )
𝑥𝑖+1 − 𝑥𝑖
𝑦𝑖+1 − 𝑦𝑖 = 𝑓(𝑦𝑖 )(𝑥𝑖+1 − 𝑥𝑖 )
𝑦𝑖+1 − 𝑦𝑖 = 𝑦𝑖 − 𝑓(𝑦𝑖 )(𝑥𝑖+1 − 𝑥𝑖 )
В нашем случае
𝑑𝐵 𝑉 3
=
𝑐𝑜𝑠𝐴
𝑑𝑆
𝐶
Очевидно тогда, если h принять равным S, то
𝑉𝑖3
𝐵𝑖+1 = 𝐵𝑖 + ( 𝑐𝑜𝑠𝐴) 𝑆
𝐶
По этой формуле вычислим приближённое значение широты
𝑉𝑖3
𝐵2 = 𝐵1 + ( 𝑐𝑜𝑠𝐴) 𝑆
𝐶
Данное решение является приближённым, т.к. мы ограничились производной
первой степени
𝑓(ℎ) = 𝑦𝑥′ = 𝑓𝑥′
Однако решение можно уточнить как точно но для этого необходимо
использовать разложение функции в ряд Тейлора до как угодно большой степени
Если записать 𝑦 = 𝑓(𝑥 + ℎ), то разложение в ряд Тейлора даст
𝑦 = 𝑓 ′ (𝑥) + 𝑓(𝑥)ℎ = 𝑓(𝑥) +
𝑦 = 𝑓(𝑥) +
𝑑𝑓
ℎ
𝑑𝑥
𝑑𝑓
ℎ
𝑑𝑥
𝑓(𝑥) = 𝑦𝑖
𝑦 = 𝑦𝑖+1
𝑓(ℎ𝑖 ) = 𝑓(𝑦)
𝑦𝑖+1 = 𝑦𝑖 + 𝑓(𝑦𝑖 )ℎ
Рассмотрим разложение в ряд Тейлора до второй степени, для обобщения будем
считать, что y0 известно, необходимо найти
𝑦1 = 𝜑(𝑥0 + ℎ)
𝑑𝜑
𝑑2𝜑 2 1
𝑦1 = 𝜑(𝑥0 ) +
ℎ+ 2ℎ
𝑑𝑥
𝑑𝑥
2!
Поскольку 𝜑(𝑥0 ) = 𝑦0
𝑑𝜑
= 𝑓(𝑦0 )
𝑑𝑥
Поскольку 𝑓(𝑥0 ) = 𝑦0 ; dφ/dx=f(y). Тогда согласно Рунге-Кутта можно записать
𝑦1 = 𝑦0 + 𝑓(𝑦0 ) ∗ ℎ +
𝑑2𝜑
ℎ2
𝑑𝑥
2
∗
2
.
Очевидно, что
𝑑2𝜑
𝑑𝑥 2
=
𝑑𝜑
)
𝑑𝑥
𝑑(
=
𝑑𝑥
𝑑(𝑓(𝑦))
𝑑𝑥
=
𝑑(𝑓(𝑦))
𝑑𝑦
=
𝑑𝑦
𝑑𝑥
.
Тогда запишем
𝑦1 = 𝑦0 + 𝑓(𝑦0 ) ∗ ℎ +
𝑑(𝑓(𝑦))
𝑑𝑦
∗
ℎ2
= 𝑦0 + 𝑓(𝑦0 ) ∗ ℎ +
2
𝑑𝑓(𝑦)
𝑑𝑥
∗
𝑑𝑦
𝑑𝑥
∗
ℎ2
2
.
Теперь запишем это выражение в приращениях, а не в дифференциалах.
𝑦1 = 𝑦0 +∗ ℎ +
∆𝑓
∆𝑦
ℎ2
∗ 𝑓(𝑦) .
2
Очевидно что
∆𝑓
𝑑𝑦
=
𝑓(𝑦1 )−𝑓(𝑦0 )
𝑦1 −𝑦0
.
В качестве 𝑦1 − 𝑦0 возьмем 𝑦1 − 𝑦0 =f(y0)*h, следовательно
𝑦1 = 𝑦0 + 𝑓(𝑦0 ) ∗ ℎ +
𝑓(𝑦1 )−𝑓(𝑦0 )
∗
ℎ2
𝑓(𝑦0 )∗ℎ
2
𝑦0 +𝑓(𝑦0 )∗ℎ+𝑓(𝑦0 )+𝑓(𝑦0 )∗ℎ−𝑓(𝑦0 )
𝑓(𝑦0 )∗ℎ
или 𝑦1 =
+ 𝑓(ℎ1 )
∗ 𝑓(𝑦1 ) −
ℎ2
2
.
Введем обозначения 𝑘0 = 𝑓(𝑦0 ) ∗ ℎ, 𝑘1 = 𝑓(𝑦0 + 𝑓(𝑦0 )) ∗ 𝛼; 𝑘 = 𝑓(𝑦0 + 𝑘0 ) ∗ ℎ.
После сокращения на h и принятых обозначений запишем
𝑘 −𝑓(𝑦0 )
𝑦 = 𝑦0 + 𝑘0 + 1
𝑓(𝑦0 )
ℎ
∗ 𝑓(𝑦1 ) ∗ .
2
Если дoпустим в данном выражении y1 = y0, следовательно
𝑘1 − 𝑓(𝑦0 )
ℎ
1
ℎ
𝑦 = 𝑦0 + 𝑘0 +
∗ 𝑓(𝑦0 ) = 𝑦0 + 𝑘1 + ∗ 𝑘 − 𝑓(𝑦0 ) ∗
𝑓(𝑦0 )
2
𝜆
2
1
1
𝑦 = 𝑦0 + 𝑘0 + ∗ (𝑘1 − 𝑘0 )или 𝑦 = 𝑦0 + ∗ (𝑘1 + 𝑘0 )
2
2
Выражение Рунте-Кутта для вычисления функций при решении данного
дифференциального уравнения с разложением в ряд Тэйлора до второй степени.
Аналогичный метод можно применить для разложений по третей, четвертой и
выше степеням. Без вывода приводим запись Рунге-Кутта для вычисления
функций при решении дифференциальных уравнений с разложением до пятой
1
степени. 𝑦1 = 𝑦0 + ∗ (𝑘1 + 4𝑘4 + 𝑘5 )
2
где 𝑘1 = ℎ ∗ 𝑓(𝑥0 )
В этом выражении в вычислении функции принимают участие элементы K1,
K4, K5, которые равны следующим значениям:
𝐾1 =
ℎ
𝑓(𝑦0 )
3
ℎ
𝑓(𝑦0 + 𝐾1 )
3
ℎ
𝐾1 𝐾2
𝐾3 = 𝑓 (𝑦0 + + )
3
2
2
ℎ
3
9
𝐾4 = 𝑓(𝑦0 + 𝐾1 + 𝐾3 )
3
8
8
ℎ
3
9
𝐾5 = 𝑓(𝑦0 + 𝐾1 − 𝐾3 + 6𝐾4 )
3
2
2
Вернёмся к нашей задаче. В ней роль y0 выполняет широта точки 1, B1
В учебниках 𝐾𝑖 = ∆𝐵𝑖
𝐾2 =
Для широты 𝑓(𝑦) =
𝑉3
𝐶
𝑐𝑜𝑠𝐴
𝑦=𝐵
𝑦0 = 𝐵1
ℎ
Тогда 𝐾1 = 𝑓(𝐵1 )
3
Обратим внимание на одно обстоятельство, вдоль геодезической линии
измеряется не только широта, но и азимут, поэтому Мерсоном было предложено
вычисление коэффициентов вести зависимости от B и A, тогда
ℎ
𝐾1 = 𝑓(𝐵1 𝐴1 )
3
ℎ
𝑓(𝐵1 + 𝐾1 , 𝐴1 + 𝐾1 )
3
ℎ
𝐾1 𝐾2
𝐾1 𝐾2
𝐾3 = 𝑓 (𝐵1 + + , 𝐴1 + + )
3
2
2
2
2
ℎ
3
9
3
9
𝐾4 = 𝑓(𝐵1 + 𝐾1 + 𝐾3 , 𝐴1 + 𝐾1 + 𝐾3 )
3
8
8
8
8
ℎ
3
9
3
9
𝐾5 = 𝑓(𝐵1 + 𝐾1 − 𝐾3 + 6𝐾4 , 𝐴1 + 𝐾1 − 𝐾3 + 6𝐾4 )
3
2
2
2
2
𝐾2 =
1
𝑦 = 𝑦0 + ∗ (𝑘1 + 4𝑘4 + 𝑘5 )
2
h=s
𝑆 𝑉13
𝐾1 = ∗
𝑐𝑜𝑠𝐴1
3 𝐶
𝑆 𝑉43
𝐾4 = ∗
𝑐𝑜𝑠𝐴4
3 𝐶
где
3
9
𝐵4 = 𝐵1 + 𝐾1 + 𝐾3
8
8
3
9
𝐴4 = 𝐴1 + 𝐾1 + 𝐾3
8
8
Аналогично вычисляются долготы
𝑆 𝑉53
𝐾5 = ∗
𝑐𝑜𝑠𝐴5
3 𝐶
3
9
𝐵5 = 𝐵1 + 𝐾1 − 𝐾3 + 6𝐾4
2
2
3
9
𝐴5 = 𝐴1 + 𝐾1 − 𝐾3 + 6𝐾4
2
2
На основе общей формулы
1
𝐵2 = 𝐵1 + (𝑘1 + 4𝑘4 + 𝑘5 )
2
Аналогично делается вычисления по долготе и азимуту, только выражение для
азимута будет
1
𝐴2 = 𝐴1 + (𝑘1 + 4𝑘4 + 𝑘5 ) ± 180°
2
Метод Рунге- Кутта был предложен в 1900 году немецкими геодезистами Рунге и
Куттом, а теоретическое применение его к вычислению геодезических координат
нашло только в 60-х годах ХХ века, к этому времени уже во всех странах
геодезические работы были выполнены и не было необходимости вычисления
данным методом, тем более он требовал применение электронных
вычислительных машин. Поэтому данных о точности этого метода нет.
Известно, что если ограничиваться рядом Тейлора 5 степени, то данный метод
даёт такую же точность на расстоянии 50 км, как и способ Гаусса
