Диссертация (Алгоритмы оценки квантильного критерия с заданной точностью в задачах стохастического программирования с кусочно-линейными и квадратичными функциями потерь), страница 12

Вычислениеэтой матрицы описано ниже.63A = ∂Q/∂q — матрица частных производных размера 2 × 5 .Элементы матрицы A . Первая строка:a)A11 = ∂S/∂xa = 1, A12 = ∂S/∂za = 0 , A13 = ∂S/∂vxa = Lv cos θa +Lθ (1+tg2 (θtg(θ,a ))|va | cos θa1A14 = ∂S/∂vya = Lv sin θa − Lθ (1+tg2 (θa ))|v, A15 = ∂S/∂vza = 0 .a | cos θaВторая строка:A21 = ∂B/∂xa = 0, A22 = ∂B/∂za = 1 , A23 = ∂B/∂vxa = 0 , A24 = ∂B/∂vya = 0,A25 = ∂B/∂vza =L.|va | cos θaλmin, λmax — наименьшее и наибольшее собственные значения (2 × 2)-матрицыACK0 C T AT .Разработанный алгоритм включает в себя следующие три этапа.Этап 1. Подготовка исходных данных для вероятностного анализа в результатемоделирования дискретных пучков траекторий ЛА на ПУТ, соответствующих различным задаваемым с шагом 500 км значениям сферической дальности полета иразличным значениям угла бросания в начале ПУТ.

Подробно этот этап описан ниже. Моделирование проводится с помощью разработанной для этой цели программыPrTrass2. На этом же этапе используется программа PassAtmProject для определения значений функции L (|va | , −θa ) для смоделированного набора значений скоростии угла входа в атмосферу.Этап 2. Определение параметров λmin , λmax для смоделированных на этапе 1 траекторий. Этот этап подробно описан ниже и реализован в виде программы Ballist.Этап 3. Для каждой смоделированной траектории искомое КВО κ определяетсяпо формуле:(115)κ=σpλmax f (γ) ,где f (γ) — КВО для двумерного нормального закона с плотностью11y22(116)p (x, y) = √ exp −x +,2π γ2γаγ=λmin.λmaxДля нахождения f (γ) в [43] разработана специальная процедура ее вы-числения с заданной точностью, реализованная как функция программы Probably.Параметр σ определяется по известному КВО κmax для траектории максимальной дальности.

Для этого сначала выполняются первые два этапа алгоритма. Затемвычисляем(117)σ=√κmax.λmax f (γ)643.4Процедура МНК оценки производной функции .Для численной оценки ряда частных производных используется следующий численный метод, основанный на МНК. Пусть требуется оценить производную f ′ (x0 )функции f (x) в точке x0 , имея возможность лишь вычислять значения f (x) в различных точках xk с некоторыми ошибками ξk (не обязательно случайными). Вычисляем значения функции в трех точках:y1 = f (x0 − h) + ξ1 , y0 = f (x0 ) + ξ0 , y2 = f (x0 + h) + ξ2 ,и с помощью МНК проводим на плоскости (x, y) прямую y = ax + b через точки(xk , yk ) , k = 0, 2 .

Получившееся при этом значение параметра a и является оценкойискомой производной:(118)a=2Pk=0xk yk −2Pk=0x2k13−2Pk=013xk ×2Pxkk=02Pykk=02.Отметим, что использование МНК обусловлено желанием уменьшить погрешности численного интегрирования дифференциальных уравнений движения ЛА ватмосфере.3.5Оценивание матрицы баллистических производных .В этом разделе описывается алгоритм определения матрицы C баллистическихпроизводных. Для реализации алгоритма раздела 3.3 необходимо рассчитать эту матрицу в СК x , y , z .

Обозначим −r→ и −v→ — вектор положения и скорости фрагментаaaavhvhпри входе в атмосферу. Сначала определим матрицы∂−r→vhTr = −=→∂ v0∂xvh∂vx0 ∂y vh ∂vx0∂zvh∂vx0∂xvh∂vy0∂yvh∂vy0∂zvh∂vy0∂xvh∂vz0∂yvh ,∂vz0 ∂zvh∂vz0∂−v→vhTv = −=→∂ v0∂vxvh∂vx0 ∂v yvh ∂vx0∂vzvh∂vx0∂vxvh∂vy0∂vyvh∂vy0∂vzvh∂vy0∂vxvh∂vz0∂vyvh ∂vz0 ∂zzvh∂vz0в АГСК . Это описано в подразделах 3.3.1 – 3.3.3. Затем вычисляем эти матрицы вСК xa , ya , za в соответствии с правилом преобразования тензора 2-го ранга:(119)Tra = T Tr T T ,Tva = T Tv T T ,→→→где строками матрицы T являются вектора −e1 , −e2 , −e3 . Искомая матрица C являетсясоставной:TraC = Tva1  ,Tva365где Tva1 и Tva3 — 1-я и 3-я строки матрицы Tva .3.5.1Дифференцирование в АГСК по vx0 .Процедура определения производных основана на дифференцировании следующихинвариантов задачи: векторного интеграла площадей−→→−→r0 × −v0 = −r→vh × vvh(120)(три уравнения), условия входа в атмосферу22|−r→vh | = const = (R3 + 90km)(121)(одно уравнение), где R3 — радиус сферической Земли, и интеграла энергии22−|−v→|→v0 |vh |−→→− Π (|rvh |) =− Π (|−r0 |)22(122)(одно уравнение), где Π — силовая функция центрального поля тяготения Земли.Дифференцирование по vx0 векторного соотношения (120) приводит к следующимтрем соотношениям:(123)∂yvh∂vzvh ∂zvh∂vyvhvzvh + yvh−vyvh − zvh= 0,∂vx0∂vx0∂vx0∂vx0(124)∂zvh∂vxvh ∂xvh∂vzvhvxvh + zvh−vzvh − xvh= z0 ,∂vx0∂vx0∂vx0∂vx0(125)∂xvh∂vyvh ∂yvh∂vxvhvyvh + xvh−vxvh − yvh= −y0 .∂vx0∂vx0∂vx0∂vx0Дифференцируя по vx0 соотношение (121), получаем(126)∂yvh∂zvh∂xvhxvh +yvh +zvh = 0.∂vx0∂vx0∂vx0Дифференцируя по соотношение (122), с учетом (121) имеем:(127)∂vxvh∂vyvh∂vzvhvxvh +vyvh +vzvh = vx0 .∂vx0∂vx0∂vx0Таким образом, получены 5 линейных уравнений (123) – (127) относительно шестинеизвестных производных∂xvh ∂yvh,∂vx0 ∂vx0,∂zvh∂vx0xvh, ∂v,∂vx066∂vyvh∂vx0,∂vzvh∂vx0.

Первую из них, т.е.∂xvh∂vx0, определяем численно по методике раздела 3.4. Обозначим ее значение через a.Тогда для остальных пяти производных получаем 5 линейных уравнений:(128)∂vzvh ∂zvh∂vyvh∂yvhvzvh + yvh−vyvh − zvh= 0,∂vx0∂vx0∂vx0∂vx0(129)∂zvh∂vxvh∂vzvhvxvh + zvh− xvh= z0 + avzvh ,∂vx0∂vx0∂vx0(130)xvh∂vyvh ∂yvh∂vxvh−vxvh − yvh= −y0 − avyvh ,∂vx0∂vx0∂vx0(131)∂yvh∂zvhyvh +zvh = −axvh ,∂vx0∂vx0(132)∂vyvh∂vzvh∂vxvhvxvh +vyvh +vzvh = vx0 .∂vx0∂vx0∂vx0Правые части уравнений (128) – (132) обозначим посредством b1÷5 . Полученнаясистема пяти линейных уравнений решается аналитически следующим образом. Из(129) получаем∂zvh1=∂vx0vxvh(133)∂vxvh∂vzvhb2 − zvh+ xvh.∂vx0∂vx0А из (131) имеем∂yvh1=−∂vx0vxvh(134)∂vyvh∂vxvhb3 − xvh+ yvh∂vx0∂vx0.Подставляя (133) и (134) в (128), (131) и (132), получаем систему из трех линейныхуравненийvyvhzvh −vxvh −v2yvhvxvh+∂vxvhxvh ∂vx0vzvhyvxvh vh2zvhvxvh+ vyvh∂vxvh∂vx0∂vxvh∂vx0∂vyvh∂vx0++vzvhxvxvh vhyvh xvh ∂vyvhvxvh ∂vx0+− zvh∂vyvh∂vx0xvh zvh ∂vzvhvxvh ∂vx0+ yvh −= b4 +vyvhxvxvh vhyvhbvxvh 3−∂vzvh∂vx0=vzvhbvxvh 3+vyvhb,vxvh 2zvhb,vxvh 2zvh+ vzvh ∂v= b5 ,∂vx0относительно трех неизвестных производных∂vxvh ∂vyvh, ∂vx0∂vx0,∂vzvh∂vx0.

Система трехлинейных уравнений решается аналитически, поскольку обратная матрица к матрицеразмера 3 × 3 аналитически выписывается в соответствием с правилом Крамера.Производные∂yvh∂vx0,∂zvh∂vx0определяются после этого из (133) и (134).673.5.2Дифференцирование в АГСК по vy0 .∂v∂yvhvhxvhzvh, ∂z, ∂v, ∂vyvh, ∂v∂vy0∂vy0∂vy0∂vy0y0vhопределяем численно повыводятся аналогично подразделу 3.5.1.

Производную ∂y∂vy0vhметодике раздела 3.4. Обозначим ее значение через a. Тогда для производных ∂x,∂vy0∂zvhполучаются следующие соотношения:∂vy0Соотношения для определения производных(135)∂zvh1=−∂vy0vyvh(136)∂xvh1=∂vy0vyvh∂xvh∂vy0,∂vyvh∂vzvhb1 + zvh− yvh∂vy0∂vy0∂vxvh∂vyvhb3 + yvh− xvh∂vy0∂vy0,,где b1 = −z0 − avzvh , b3 = x0 + avxvh .∂vxvh∂vy0Система же трех линейных уравнений относительно производных∂vzvh∂vy0∂vyvh∂vy0,получается следующей:zvh −∂vxvh∂vy0+ vyvh∂vyvh∂vy0vzvhyvyvh vhxvh yvh ∂vxvh vyvh ∂vy0v,∂vxvhxvh ∂vy0−x2vhvyvh++2zvhvyvhvzvhxvyvh vh∂vyvh∂vy0∂vyvh∂vy0−vxvhzvyvh vh+zvh yvh ∂vzvhvyvh ∂vy0+= b4 +vzvhyvyvh vhzvhbvyvh 1− xvh−∂vzvh∂vy0vxvhbvyvh 1=+vzvhb,vyvh 3xvhb,vyvh 3zvh+ vzvh ∂v= b5∂vy0где b4 = −ayvh , b5 = vy0 .3.5.3Дифференцирование в АГСК по vz0 .∂xvh ∂yvhvhxvh,, ∂z, ∂v∂vz0 ∂vz0∂vz0∂vz0vhПроизводную ∂zопределяем∂vz0∂vyvh∂vz0∂vzvh∂vz0Соотношения для определения производных,дятся аналогично подразделу 3.5.1.численно по мето-,дике раздела 3.4.

Обозначим ее значение через a . Тогда для производныхвыво-∂yvh∂vz0vh, ∂x∂vz0получаются следующие соотношения:∂vzvh∂vyvhb1 − yvh+ zvh,∂vz0∂vz0(137)∂yvh1=∂vz0vzvh(138)∂xvh1=−∂vz0vzvh∂vxvh∂vzvhb2 − zvh+ xvh∂vz0∂vz0где b1 = y0 + avyvh , b2 = −x0 − avxvh .,Система же трех линейных уравнений относительно производных∂vzvh∂vz0получается следующей:68∂vxvh∂vz0,∂vyvh∂vz0,vyvh∂vxvhzvh − yvh ∂vz0 + xvh − vzvhxvh zvh ∂vxvh vzvh ∂vz0v∂vxvhxvh ∂vz0+yvh zvh ∂vyvhvzvh ∂vz0+ vyvh∂vyvh∂vz0−x2vhvzvhvxvhzvzvh vh+2yvhvzvhzvh+ vzvh ∂v= b5∂vz0∂vyvh∂vz0+∂vzvh∂vz0= b4 +vxvhyvzvh vh−xvhbvzvh 2vyvhxvzvh vh−∂vzvh∂vz0=vyvhbvzvh 2+vxvhb,vzvh 1yvhb,vzvh 1где b4 = −azvh , b5 = vz0 .3.63.6.1Моделирование траекторий.Расчет начальных условий баллистического полета на заданнуюдальность.Отметим, что поставлена задача исследования некоторых свойств баллистическогополета (БП) с использованием сферической модели Земли без учета влияния атмосферы и вращения Земли.

Исходными данными являются требуемая сферическаядальность полета lсф и определенные ниже в подразделе 3.6.3 параметры конца АУТ- дальность и высота полета на АУТ: lк и hк . Требуется определить в векторнойформе радиус вектор и скорость полета на момент начала БП, а также координатыточек падения в целевой системе координат, соответствующие выбранным начальным условиям (НУ). Предполагается, что величина и наклон вектора скорости намомент окончания АУТ будут варьироваться с помощью датчика случайных чиселв определенных пределах.В сферической системе координат дальность баллистического полета lб определяется параметрами конца АУТ, выступающими в качестве НУ БП: lб = lб (Vк , θк , hк ),где θк ) — наклон вектора скорости Vк в конце АУТ (в точке К); hк — высота конечной точки АУТ.

Если требуется определить значения перечисленных параметров дляслучая полета на заданную дальность (что требуется в настоящем исследовании), товозникает краевая задача, имеющая аналитическое решение (исчерпывающе подробно эта задача исследована Д.А. Погореловым [72]). С учетом того обстоятельства, чтов данной задаче lк , известная дальность полета на АУТ, требуемая дальность БП lби соответствующая ей сферическая дальность (F ) БП вычисляются по формулам:(139)lб = lсф − lк ,F = lб /Rз .Здесь радиус Земли Rз = 6371 км.Уравнение, связывающее НУ БП при заданной угловой дальности F (уравнениегодографа скорости) имеет вид [72]:(140)V2 =2b0hi ,F2F2r sin 2θк ctg 2 + cos θк 2 + ĥ + ĥ ctg 269где относительная высота полета (параметр ĥ ) рассчитывается по формуле:(141)ĥ =hк,Rз + hцгравитационный параметр b0 = 398600.44 км3 /с2 , а V –– требуемая для полета назаданную угловую дальность скорость в начале БП при заданном же угле наклонаэтой скорости.Основным вариантом расчетов при использовании рассматриваемого алгоритмаявляется случай, когда параметр θк задается в исходных данных (варьируется отминимального до максимального).