Доклады Российской академии наук. Математика, информатика, процессы управления, 2020, T. 494, № 1, стр. 17-20

ДВУХЭТАПНЫЙ МЕТОД РЕШЕНИЯ СИСТЕМ НЕЛИНЕЙНЫХ УРАВНЕНИЙ И ЕГО ПРИЛОЖЕНИЕ К ОБРАТНОЙ ЗАДАЧЕ ЗОНДИРОВАНИЯ АТМОСФЕРЫ

Член-корреспондент РАН В. В. Васин^1, 2, *, Г. Г. Скорик^1, 2, **

¹ Институт математики и механики им. Н.Н. Красовского
Екатеринбург, Россия

² Уральский федеральный университет им. Б.Н. Ельцина
Екатеринбург, Россия

^* E-mail: vasin@imm.uran.ru
^** E-mail: skorik@imm.uran.ru

Поступила в редакцию 04.06.2020
После доработки 04.06.2020
Принята к публикации 16.07.2020

DOI: 10.31857/S2686954320050458

Полный текст (PDF)

Аннотация

Для переопределенной системы нелинейных уравнений предлагается двухэтапный метод построения устойчивого к ошибкам приближенного решения. Первый этап состоит в построении регуляризованного семейства приближенных решений для нахождения нормальных квазирешений исходной системы. На втором этапе для аппроксимации регуляризованных квазирешений строится итерационный процесс, основанный на квадратичной аппроксимации тихоновского функционала и привлечении prox-метода. Для сформированного процесса ньютоновского типа доказывается теорема сходимости и устанавливается свойство фейеровости итераций. Обсуждается приложение двухэтапного метода к решению обратной задачи по восстановления относительного содержания тяжелой воды (HDO) в атмосфере по инфракрасным спектрам пропускания солнечного света через атмосферу.

Ключевые слова: нелинейная система, регуляризация, итерационный процесс, инфракрасный спектр, зондирование атмосферы, восстановление HDO

ВВЕДЕНИЕ

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

(1)

$P(z) = f,\quad z \in {{\mathbb{R}}^{n}},\quad f \in {{\mathbb{R}}^{m}}.$

При этом размерность m существенно больше, чем n, что приводит к переопределенности системы (1). В работе [1] был предложен метод решения системы нелинейных уравнений (1), в котором можно выделить следующие этапы построения приближенного решения:

тихоновская регуляризация системы (1)

(2)

$inf\left\{ {{{{\left\| {P(z) - {{f}_{\delta }}} \right\|}}^{2}} + \alpha \left\langle {B(z - \bar {z}),\;z - \bar {z}} \right\rangle {\text{:}}\;z \in {{\mathbb{R}}^{n}}} \right\},$

где B – симметричная положительно определенная матрица, $\bar {z}$ – некоторый фиксированный вектор пространства ${{\mathbb{R}}^{n}}$, $\left\| {f - {{f}_{\delta }}} \right\| \leqslant \delta $;

в предположении, что оператор P дважды непрерывно дифференцируем, квадратичная аппроксимация целевой функции $\Psi (z)$ в задаче (2) и привлечение итерационного prox-метода для решения задачи минимизации с полученной квадратичной функцией, т.е.

(3)

$\begin{gathered} {{z}^{{k + 1}}} = arg\,min\left\{ {\Psi ({{z}^{k}}) + \Psi {\text{'}}({{z}^{k}})(z - {{z}^{k}}){{ + }_{{_{{_{{_{{_{{_{{_{{}}}}}}}}}}}}}}}} \right. \\ + \;\frac{1}{2}\langle \Psi ''({{z}^{k}})(z - {{z}^{k}}),z - {{z}^{k}}\rangle + \left. {\frac{1}{2}\beta {\text{||}}z - {{z}^{k}}{\text{|}}{{{\text{|}}}^{2}}{\text{:}}\;z \in {{\mathbb{R}}^{n}}} \right\}. \\ \end{gathered} $

Для численного нахождения очередного приближения в (3) привлекался метод наискорейшего спуска.

Как показано в [2], при условии разрешимости задачи в смысле наименьших квадратов

(4)

$inf\{ {{\left\| {P(z) - f} \right\|}^{2}}{\text{:}}\;z \in {{\mathbb{R}}^{n}}\} $

регуляризованная задача (2) имеет решение z_α и при $\alpha (\delta ) \to 0$, ${{\delta }^{2}}{\text{/}}\alpha (\delta ) \to 0$, δ → 0 все предельные точки последовательности ${\text{\{ }}{{z}_{{\alpha (\delta )}}}{\text{\} }}$ являются решениями задачи (4), ближайшими к вектору $\bar {z}$, т.е. $\bar {z}$-нормальными квазирешениями. Кроме того, была установлена сходимость процесса (3) к некоторому регуляризованному решению u_α и его {z_α}-фейеровость.

Предположим, что множества решений зада-чи (2) и уравнения

(5)

$\Phi (z) + \alpha B(z - \bar {z}) = 0,$

где $\Phi (z) = P\,{\text{'}}{{(z)}^{{\text{т}}}}(P(z) - {{f}_{\delta }})$, совпадают при всех данных f_δ, $\left\| {f - {{f}_{\delta }}} \right\| \leqslant \delta $, δ, α.

Запишем процесс (3) в эквивалентном виде

(6)

$\begin{gathered} {{z}^{{k + 1}}} = {{z}^{k}} - {{[\Phi {\text{'}}({{z}^{k}}) + \alpha B + \beta I]}^{{ - 1}}} \times \\ \times \;[\Phi ({{z}^{k}}) + \alpha B(z - \bar {z})]. \\ \end{gathered} $

При построении итераций согласно (6) наиболее трудоемкой операцией является вычисление и обращение симметричной матрицы

$\begin{gathered} A({{z}^{k}}) = \Phi {\text{'}}({{z}^{k}}) + \alpha B + \beta I = \frac{1}{2}({\text{||}}P({{z}^{k}}) - {{f}_{\delta }}{\text{||}}_{{{{\mathbb{R}}^{n}}}}^{2} + \\ + \;\alpha \langle B({{z}^{k}} - \bar {z}),{{z}^{k}} - \bar {z}\rangle + \beta {\text{||}}{{z}^{k}} - \bar {z}{\text{|}}{{{\text{|}}}^{2}})'' \\ \end{gathered} $

в каждой итерационной точке ${{z}^{k}}$. Поэтому естественно рассмотреть модифицированный вариант метода (6)

(7)

$\begin{gathered} {{z}^{{k + 1}}} = {{z}^{k}} - \gamma {{[\Phi {\text{'}}({{z}^{0}}) + \alpha B + \beta I]}^{{ - 1}}} \times \\ \times \;[\Phi ({{z}^{k}}) + \alpha B(z - \bar {z})], \\ \end{gathered} $

в котором матрица вторых производных A(z) вычисляется один раз в точке z⁰ в течение всего процесса итераций, а также дополнительно введен параметр шага γ. Считаем, что точка z⁰, в которой вычисляется матрица A(z), совпадает со стартовой точкой процесса (7).

Заметим, что при формировании как основного процесса (3), так и его модифицированного варианта (7) наряду с квадратичной аппроксимацией (см. например [3, гл. 5 §10]) применяется prox-метод, что позволяет не только получить сходимость итераций, но и оценить скорость сходимости.

1. ИССЛЕДОВАНИЕ СХОДИМОСТИ МОДИФИЦИРОВАННОГО МЕТОДА

Предположим, что при заданном уровне погрешности δ правой части уравнения (1) выбран параметр $\alpha (\delta )$, при котором ${{z}_{{\alpha (\delta )}}}$ аппроксимирует с приемлемой точностью некоторое $\bar {z}$-нормальное квазирешение $\hat {z}$ задачи (4). Введем обозначения:

$\begin{gathered} \mathcal{D} = (\Phi {\text{'}}({{z}^{0}}) + \alpha B + \beta I), \\ F(z) = {{\mathcal{D}}^{{ - 1}}}\left[ {\Phi (z) + \alpha B(z - \bar {z})} \right], \\ {{\mathcal{D}}_{{\text{1}}}} = (\Phi {\text{'}}({{z}^{0}}) + \alpha B),\quad m = \mathop {inf}\limits_{\left\| z \right\| = 1} \left\langle {{{\mathcal{D}}_{1}}z,z} \right\rangle , \\ M = \mathop {sup}\limits_{\left\| z \right\| = 1} \left\langle {{{\mathcal{D}}_{1}}z,z} \right\rangle . \\ \end{gathered} $

Установим сильную фейеровость оператора шага T процесса (7) и погрешность аппроксимации регуляризованного решения ${{z}_{{\alpha (\delta )}}}$ итерационным процессом (7).

Теорема 1. Пусть выполнены условия

$\begin{gathered} \left\| {\Phi ({{z}_{1}}) - \Phi ({{z}_{2}})} \right\| \leqslant {{N}_{1}}\left\| {{{z}_{1}} - {{z}_{2}}} \right\|, \\ \left\| {\Phi {\text{'}}({{z}_{1}}) - \Phi {\text{'}}({{z}_{2}})} \right\| \leqslant {{N}_{2}}\left\| {{{z}_{1}} - {{z}_{2}}} \right\| \\ \end{gathered} $

в шаре $S({{z}^{0}},d)$, содержащем ${{z}_{{\alpha (\delta )}}}$, $\hat {z}$, а также

${\text{||}}{{z}_{{\alpha (\delta )}}} - {{z}^{0}}{\text{||}} \leqslant {{r}^{0}},\quad {{r}^{0}} = m{\text{/}}5{{N}_{2}},\quad m > 0.$

Тогда при $\gamma < m(m + \beta ){\text{/}}{{({{N}_{1}} + \alpha \left\| B \right\|)}^{2}}$ оператор T есть сильно ${\text{\{ }}{{z}_{{\alpha (\delta )}}}{\text{\} }}$-фейеровский в шаре $S({{z}^{0}},{{r}^{0}})$ и z^k сходится к ${{z}_{{\alpha (\delta )}}}$, а при $\gamma = m(m + \beta ){\text{/}}(2{{({{N}_{1}} + \alpha \left\| B \right\|)}^{2}})$ справедлива оценка

(8)

$\begin{gathered} {\text{||}}z_{{\alpha (\delta )}}^{k} - {{z}_{{\alpha (\delta )}}}{\text{||}} \leqslant {{q}^{k}}{\text{||}}{{z}^{0}} - {{z}_{{\alpha (\delta )}}}{\text{||}}, \\ q = \sqrt {1 - {{m}^{2}}{\text{/}}{{{({{N}_{1}} + \alpha \left\| B \right\|)}}^{2}}} , \\ \end{gathered} $

где z⁰– начальная точка в методе (7).

Доказательство. Используя введенные обозначения, условия теоремы и самосопряженность оператора ${{\mathcal{D}}_{1}}$, имеем для $z \in S({{z}_{{\alpha (\delta )}}},{{r}_{0}})$ соотношения

$\begin{gathered} \left\langle {F(z),z - {{z}_{\alpha }}} \right\rangle = \left\langle {F(z),z - {{z}_{\alpha }}} \right\rangle - \left\langle {F({{z}_{\alpha }}),z - {{z}_{\alpha }}} \right\rangle = \\ = \;\left\langle {{{\mathcal{D}}^{{ - 1}}}\left[ {\Phi (z) - \Phi ({{z}_{\alpha }}) + \alpha B(z - {{z}_{\alpha }})} \right],z - {{z}_{\alpha }}} \right\rangle = \\ = \;\left\langle {{{\mathcal{D}}^{{ - 1}}}} \right.\left[ {\Phi (z)} \right. - (\Phi (z) + \Phi {\text{'}}(z)({{z}_{\alpha }} - z) + \xi ) + \\ + \;\left. {\left. {\alpha B(z - {{z}_{\alpha }})} \right],z - {{z}_{\alpha }}} \right\rangle = \\ = \left\langle {{{\mathcal{D}}^{{ - 1}}}} \right.\left[ {\Phi {\text{'}}(z)} \right.(z - {{z}_{\alpha }}) + \alpha B(z - {{z}_{\alpha }}) + \\ \end{gathered} $

$\begin{gathered} + \;(\Phi {\text{'}}({{z}^{0}})(z - {{z}_{\alpha }}) - \left. {\left. {\Phi {\text{'}}(z)(z - {{z}_{\alpha }})) - \xi } \right],z - {{z}_{\alpha }}} \right\rangle \geqslant \\ \geqslant \;\frac{m}{{m + \beta }}{{\left\| {z - {{z}_{\alpha }}} \right\|}^{2}} - \frac{1}{{m + \beta }} \times \\ \times \;\left[ {{{N}_{2}}{{{\left\| {z - {{z}_{\alpha }}} \right\|}}^{2}}{\text{||}}z - {{z}^{0}}{\text{||}} + \frac{1}{2}{{N}_{2}}{{{\left\| {z - {{z}_{\alpha }}} \right\|}}^{3}}} \right] \geqslant \\ \geqslant \;\frac{m}{{m + \beta }}{{\left\| {z - {{z}_{\alpha }}} \right\|}^{2}} - \frac{1}{{m + \beta }}{{N}_{2}} \times \\ \times \;\left( {\left\| {z - {{z}_{\alpha }}} \right\| + \,{\text{||}}{{z}_{\alpha }} - {{z}^{0}}{\text{||}} + \frac{1}{2}\left\| {z - {{z}_{\alpha }}} \right\|} \right){{\left\| {z - {{z}_{\alpha }}} \right\|}^{2}}. \\ \end{gathered} $

Подставляя в соотношения ${{z}_{\alpha }} = {{z}_{{\alpha (\delta )}}}$ и учитывая в правой части последнего неравенства $z\, \in \,S({{z}_{{\alpha (\delta )}}}$, r⁰), ${\text{||}}z - {{z}^{0}}{\text{||}} \leqslant {\text{||}}z - {{z}_{{\alpha (\delta )}}}{\text{||}} + \,{\text{||}}{{z}_{{\alpha (\delta )}}} - {{z}^{0}}{\text{||}} \leqslant 2r$, $r = m{\text{/}}5{{N}_{2}}$, окончательно получаем оценку снизу:

(9)

$\left\langle {F(z),z - {{z}_{{\alpha (\delta )}}}} \right\rangle \geqslant \frac{m}{{2(m + \beta )}}{{\left\| {z - {{z}_{{\alpha (\delta )}}}} \right\|}^{2}}.$

В силу условия Липшица для $\Phi (z)$, справедлива оценка сверху для ${{\left\| {F(z)} \right\|}^{2}}$:

$\begin{gathered} {{\left\| {F(z)} \right\|}^{2}} = {{\left\| {F(z) - F({{z}_{{\alpha (\delta )}}})} \right\|}^{2}} = \\ = \;{{\left\| {{{\mathcal{D}}^{{ - 1}}}\left[ {\Phi (z) - \Phi ({{z}_{{\alpha (\delta )}}}) + \alpha B(z - {{z}_{\alpha }})} \right]} \right\|}^{2}} \leqslant \\ \end{gathered} $

(10)

$ \leqslant \frac{1}{{{{{(m + \beta )}}^{2}}}}{{({{N}_{1}} + \alpha {\text{||}}B{\text{||}})}^{2}}{\text{||}}z - {{z}_{{\alpha (\delta )}}}{\text{|}}{{{\text{|}}}^{2}}.$

Объединяя (9) и (10), приходим к оценке

(11)

${{\left\| {F(z)} \right\|}^{2}} \leqslant \frac{{2{{{({{N}_{1}} + \alpha \left\| B \right\|)}}^{2}}}}{{m(m + \beta )}}\left\langle {F(z),z - {{z}_{\alpha }}(\delta )} \right\rangle .$

Сильная {u_α}-фейеровость оператора шага T в процессе (7) означает выполнение неравенства (см. определение 1.5 в [4])

(12)

${{\left\| {T(z) - {{z}_{{\alpha (\delta )}}}} \right\|}^{2}} - {{\left\| {z - {{z}_{{\alpha (\delta )}}}} \right\|}^{2}} + \nu {{\left\| {T(z) - z} \right\|}^{2}} \leqslant 0$

для некоторого ν > 0, что эквивалентно неравенству

(13)

${{\left\| {F(z)} \right\|}^{2}} \leqslant \frac{2}{{\gamma (1 + \nu )}}\left\langle {F(z),z - {{z}_{{\alpha (\delta )}}}} \right\rangle .$

Сравнивая неравенства (11) и (13), убеждаемся, что при $\gamma < m(m + \beta ){\text{/}}{{({{N}_{1}} + \alpha \left\| B \right\|)}^{2}}$ оператор T является сильно ${\text{\{ }}{{z}_{{\alpha (\delta )}}}{\text{\} }}$-фейеровским, что влечет сходимость z^k к ${{z}_{{\alpha (\delta )}}}$ [4].

Если теперь подставить оценки (8) и (9) в соотношение

(14)

$\begin{gathered} {{\left\| {T(z) - {{z}_{\alpha }}} \right\|}^{2}} = \\ = {{\left\| {z - {{z}_{\alpha }}} \right\|}^{2}} - 2\gamma \left\langle {F(z),z - {{z}_{\alpha }}} \right\rangle + {{\gamma }^{2}}{{\left\| {F(z)} \right\|}^{2}} \\ \end{gathered} $

и проминимизировать правую часть полученного неравенства по γ, то получим значение γ = m(m + + β)/ $2{{({{N}_{1}} + \alpha \left\| B \right\|)}^{2}}$, при котором при z = z^k (14) перейдет в оценку (8).

Замечание 1. В работе [2] для основного процесса (6) доказана сходимость итераций к регуляризованному решению и свойство фейеровости оператора шага. Переход к модифицированному варианту (7), введение дополнительного параметра γ и использование иной схемы доказательства позволили установить сильную фейеровость оператора T, что влечет не только сходимость, но и оценку скорости сходимости итераций.

2. ОБРАТНАЯ ЗАДАЧА ЗОНДИРОВАНИЯ АТМОСФЕРЫ

В качестве примера реализации представленного метода выбрана задача восстановления высотных профилей $\delta {{D}_{{{\text{HDO}}}}}(h)$ относительной концентрации молекул HDO в водяном паре по инфракрасным спектрам пропускания солнечного света. Количественные данные об относительной концентрации тяжелой воды являются важным маркером скрытых трендов в гидрологических процессах, происходящих в атмосфере.

Относительная концентрация $\delta {{D}_{{{\text{HDO}}}}}$ определяется формулой

$\delta {{D}_{{{\text{HDO}}}}} = \left( {\frac{R}{{{{R}_{0}}}} - 1} \right) \times 1000\% ,\quad R = \frac{{{{N}_{{{\text{HDO}}}}}}}{{{{N}_{{{{{\text{H}}}_{{\text{2}}}}{\text{O}}}}}}},$

где ${{N}_{{{\text{HDO}}}}}$ и ${{N}_{{{{{\text{H}}}_{{\text{2}}}}{\text{O}}}}}$ – соответственно число молекул HDO и H₂O в единице объема; ${{R}_{0}} = 3.1069 \times {{10}^{{ - 4}}}$ – среднее значение R в океанической воде. Величина $\delta {{D}_{{{\text{HDO}}}}}$ в атмосфере лежит в пределах – 1000–0‰.

Исследуемая задача сводится к решению переопределенной системы нелинейных уравнений

(15)

$P(z) = {{\bar {P}}_{\delta }},$

где ${{\bar {P}}_{\delta }}$ – измеренный с ошибкой $\delta $ вектор значений ИК-спектра солнечного света, на сетке частот из микроокон, содержащих выбранные спектральные линии HDO/H₂O. Спектр измерялся фурье-интерферометром FTIR (Fourier Transform InfraRed “Bruker” IFS 125M), установленным на Уральской атмосферной станции, расположенной на территории Коуровской астрономической обсерватории УрФУ, находящейся в окрестности Екатеринбурга. Вектор $z$ составлен из значений $ln{{N}_{{{{{\text{H}}}_{{\text{2}}}}{\text{O}}}}}({{h}_{i}})$ и $\delta {{D}_{{{\text{HDO}}}}}({{h}_{i}})$ на заданной сетке высот {h_i}, P – нелинейный оператор, вычисляющий модельный спектр по заданному z. Оператор P строится на основе данных, полученных в результате работы программы FIRE-ARMS, описание которой содержится в работе [5], решающей прямую задачу вычисления модельного спектра по известным профилям концентрации газов, давления и температуры, полученных из реанализа метеоданных на момент измерения спектра. Построение нелинейного оператора $P(z)$ подробно описано в [2].

Решение задачи (15) сильно чувствительно к ошибкам измерения спектра, поэтому вначале применяется тихоновская регуляризация со стабилизатором, полученным на основе метода максимального правдоподобия,

(16)

$min\left\{ {{{{\left\| {P(z) - {{{\bar {P}}}_{\delta }}} \right\|}}^{2}} + \alpha \langle {{C}^{{ - 1}}}(z - \bar {z}),z - \bar {z}\rangle {\text{:}}\;z \in {{R}^{n}}} \right\},$

где $C = (K + \varepsilon I)$, $\bar {z}$ и K – среднее значение и ковариационная матрица априорного распределения вектора z, $\varepsilon $ – малый положительный параметр, выбираемый из условия положительной определенности матрицы $C$.

Для аппроксимации решения z_α задачи (16) использовался итерационный процесс (7) с γ = 1 и несколько видоизмененный, что позволило исключить обращение матрицы C,

(17)

$\begin{gathered} {{z}^{{k + 1}}} = {{z}^{k}} - {{[C\Phi {\text{'}}({{z}^{0}}) + \alpha I + \beta C]}^{{ - 1}}} \times \\ \times \;(C\Phi ({{z}^{k}}) + \alpha ({{z}^{k}} - \bar {z})), \\ \end{gathered} $

где $\Phi (z) = P\,{\text{'}}{{(z)}^{{\text{т}}}}(P(z) - {{\bar {P}}_{\delta }})$, $\beta $ – параметр prox-метода.

Для сравнения, задача минимизации (16) решалась также с помощью немодифицированного итерационного процесса с вычислением матрицы вторых производных на каждом шаге итерации

(18)

$\begin{gathered} {{z}^{{k + 1}}} = {{z}^{k}} - {{[C\Phi {\text{'}}({{z}^{k}}) + \alpha I + \beta C]}^{{ - 1}}} \times \\ \times \;(C\Phi ({{z}^{k}}) + \alpha ({{z}^{k}} - \bar {z})). \\ \end{gathered} $

Метод (17) был протестирован на множестве наборов (около 260) измеренных спектральных данных, описание которых содержится в [6]. Также было произведено сравнение времени обработки спектров методами (17) и (18). Для данной задачи новый метод дает прирост производительности на несколько процентов.

Замечание 2. Незначительная разница по затратам машинного времени при использовании модифицированного (17) и немодифицированного (18) процессов в описанном эксперименте обусловлена невысокой размерностью задачи (dim(z) = = 100). Однако в задачах высокой размерности преимущество модифицированных методов по времени является существенным (см. [7]).

Список литературы

Skorik G.G. // Evrasian Journal Mathematical and Computer Applications. 2018. V. 6. Iss. 1. P. 56–64.
Skorik G.G., Vasin V.V. // Evrasian J. Math. and Comp. Appl. 2019. V. 7. Iss. 2. P. 79–88.
Васильев Ф.П. Методы оптимизации. Кн. 1 М.: МЦНМО, 2011. 619 с.
Vasin V.V., Eremin I.I. Operators and Iterative Processes of Fejér Type. Theory and Application. B., N.Y.: Walter de Gruyter, 2009. 155 p.
Gribanov K.G., Zakharov V.I., Tashkun S.A. et al. // Quantitative Spectroscopy and Radiative Transfer. 2001. V. 68. Iss. 4. P. 435–451.
Skorik G.G., Vasin V.V., Gribanov K.G. et al. // Doklady Earth Sciences. 2014. V. 454. Iss. 2. P. 208–212.
Васин В.В., Скурыдина А.Ф. // Тр. Ин-та математики и механики УрО РАН. 2019. Т. 23. № 1. С. 57–74.

Дополнительные материалы отсутствуют.

Инструменты

следующая статья выпуска предыдущая статья выпуска содержание выпуска

Доклады Российской академии наук. Математика, информатика, процессы управления

Архивы выпусков Информация о журнале Отправить рукопись в журнал