. Компьютерное моделирование (стр. 6 )

Пренебрегая при этом изменением на участке скольжения ψ, ε, имеем:

Это выражение и можно использовать для определения z’. Однако поскольку в программе расчёта имеется процедура интегрирования расширенной системы уравнений, то целесообразнее её использовать для нахождения z’, т. е. для z’ на основании (1.37) получим уравнение:

Так как = sin (ψτ) (1 + z/R0) (см. рис. 1.5), то

Из геометрии оси антиволновода (рис. 1.5) находим

,

здесь определяется из следующего уравнения, полученного из (1.28).

Значение ψ′ в точке Bτ получим из третьего уравнения расширенной системы с учётом (1.39):

Используя найденные таким образом значения ψ, z′, ψ′ в качестве начальных, проводим далее геометрооптический расчёт поля вне окрестности оси антиволновода (обозначенной пунктиром на рис. 1.4). Подбирая на оси точку τ, определяем значение поля в заданной точке земной поверхности (рис. 1.4).

Точность рассматриваемого метода зависит от того, насколько совпадают геометрооптическое решение и первая направляемая волна (1.24) в месте их сшивания, что, в свою очередь, зависит от близости к оси скользящих траекторий. Эту близость характеризует расстояние d (см. выше) между высотами точки отражения траектории и оси. Тестовые расчёты показали, что при величине d порядка 1 км результаты вычислений данным способом совпадают с результатами, получающимися с использованием только геометрооптического метода. Это значение для d и закреплено при сшивке двух решений.

Алгоритм, реализующий данный комбинированный способ, показал его высокую эффективность. На рис. 1.6 представлен график, показывающий рост относительных затрат времени расчёта квазикритического распространения при увеличении дальности трассы для обычного геометрооптического метода (∆t0) по сравнению с комбинированным (∆t2).

Рис. 1.6. Зависимость относительных затрат времени расчёта квазикритического распространения на ЭВМ от дальности трассы для двух способов

Как видно, на больших расстояниях (около 3,5 тыс. км) использование одного геометрооптического метода для расчёта квазикритического распространения затруднительно (время, затрачиваемое на вычисления, резко растёт). Отметим, что при моделировании, в ряде случаев (см. разд. 2) требуется расчёт скользящих лучей и на трассах с большей дальностью, именно в таких случаях и целесообразно применение комбинированного способа.

1.5. Аналитическое представление дискретных данных
ионосферной модели

Особенности построения полуэмпирической модели ионосферы (как и других глобальных моделей, учитывающих физические процессы, протекающие в ионосфере) приводят к дискретной выдаче значений электронной концентрации в определённых точках области, охватывающей радиотрассу. В то же время в программе расчёта распространения радиоволн требуется аналитическое представление зависимости N(x, y, z). Поэтому использование современных сложных моделей среды вызвало необходимость решения задачи их сопряжения с методом характеристик. Впервые автором была применена [87] для этих целей кубичная сплайн-интерполяция [96]. Были показаны её удобство и хорошая точность. Вопрос об интерполяции сплайнами данных ионосферных моделей рассматривался в дальнейшем в работах [34, 49, 102]. Детальное исследование этого вопроса было проведено нами в [9, 46, 57, 81]. Изложим кратко основные моменты разработанной методики интерполяции на примере двухмерно-неоднородной ионосферы.

Пусть некоторая прямоугольная область [x0, xn; z0, zm] разделена линиями, параллельными осям координат и пересекающими их в точках:

x0< x1<…< xj<…<xn – 1<xn,

z0 < z1<…< zj <…< z m – 1< zm.

Кубическим сплайном называется функция y (x, z), если она обладает следующими свойствами:

1)  в каждом прямоугольнике [x0, xn; z0, zm] является бикубическим полиномом вида

2)  непрерывна и имеет непрерывные частные производные до включительно во всей заданной области.

Задача интерполирования функции ε (x, z) формулируется так. Пусть в области в узлах прямоугольной сетки заданы значения. Требуется построить бикубический сплайн на той же сетке, интерполирующий в заданных узлах, т. е. , где , и удовлетворяющий дополнительным граничным условиям.

Такая постановка определяет единственным образом и приводит к системе линейных алгебраических уравнений для коэффициентов.

Матрица системы при этом получается трёхдиагональной, для которой имеются разработанные варианты решения системы методом прогонки со слабой чувствительностью к ошибкам округления в промежуточных вычислениях. В работе [64] приводится процедура кусочно-полиномиальной кубичной интерполяции, которая, как подтвердили проведённые вычисления, с успехом может быть применена для аналитического представления ионосферных данных. При этом достаточно эффективен вариант выбора так называемых естественных граничных условий, когда вторая производная интерполируемой функции на концах отрезка равна нулю (это соответствует предположению о том, что функция вне отрезка интерполирования является линейной).

Выбор степени интерполяции связан, в первую очередь, с необходимостью обеспечить при интерполяции непрерывность вторых производных от ε, входящих в систему лучевых уравнений. Кроме того, сплайны третьей степени обладают наилучшими гладкими свойствами [96]. Вместе с тем для интерполяции зависимости основных ионосферных параметров вдоль трассы, используемой в некоторых случаях, или для интерполяции формы ДУХ (см. подраздел 1.3), вполне достаточно применения квадратичных сплайнов, при вычислении коэффициентов которых легко могут быть получены рекуррентные выражения, не требующие хранения полного массива коэффициентов на всем интервале интерполяции. Так, записывая сплайн на каком-то k-м интервале и вводя для удобства коэффициенты a, b, с имеем:

На основании свойств сплайна можем записать:

Задание граничных условий допускает несколько вариантов, стандартным из которых считается тот, когда задано значение производных на концах отрезка

Учитывая, что gk(xk) = ak, получаем выражение для определения коэффициентов ak:

Далее, в узле k выполняется

Предположим для простоты, что шаг интерполяции является равномерным. Обозначая его через h, имеем:

Затем:

Отсюда:

и

Подставляя теперь значение ck – 1h, получим

тогда

и

Таким образом, выражения (1.45)–(1.48) определяют коэффициенты сплайна на любом интервале. В крайних интервалах для вычисления (1.46) и (1.48) используются условия: bn+1 = f ’(xn), b0 = f ’(x). Случай, когда известны значения производной, является довольно редким, поэтому обычно для них берут первые разделённые разности. Учитывая, что от точности задания граничных условий существенно зависит ошибка интерполяции в ближних интервалах (далее к середине она уменьшается), можно при необходимости задать производные и точнее, применив, например, экстраполяцию полиномом Лагранжа [42].

При организации алгоритма интерполяции массивы значений аргумента и функции расширяются в зависимости от выбранного типа граничных условий. В нашем случае массивы будут