Ускорение сходимости ньютоновских методов

Полный текст

Введение

Известным средством ускорения сходимости ньютоновских методов к особым решениям нелинейных уравнений является экстраполяция [1, 2]. Однако, как демонстрируется в этой статье, ускоряющий эффект этой процедуры для разных ньютоновских методов может быть разным, и целью работы является теоретическое объяснение этой разницы для линейно-квадратичных уравнений, а также численный анализ эффекта от экстраполяции для глобализованных одномерным поиском методов Ньютона, Левенберга–Марквардта, а также ньютоновского метода с подзадачами линейного программирования (соответствующие ссылки на литературу приводятся ниже).

Используемые в статье обозначения вполне традиционны. Через $L^{\perp }$ обозначается ортогональное дополнение линейного подпространства L. Если не указано иное, будем считать, что $\parallel \cdot \parallel$ — евклидова норма; в некоторых случаях будет также использоваться $\infty$-норма, которая для $u\in {ℝ}^p$ определяется как $\parallel u{\parallel }_{\infty }=\max \left\{\right|u_j\left|\right|j\in \{\mathrm{1,}\dots ,p\}\}.$ Для заданных $\overline{u},\overline{v}\in {ℝ}^p,$$\epsilon >0$ и $\delta >0$ определим множество

$K_{\epsilon ,\delta }(\overline{u};\overline{v})=\{u\in {ℝ}^p|\parallel u-\overline{u}\parallel \le \epsilon ,\parallel \parallel \overline{v}\parallel (u-\overline{u})-\parallel u-\overline{u}\parallel \overline{v}\parallel \le \delta \parallel u-\overline{u}\parallel \parallel \overline{v}\parallel \}.$

Через I обозначается единичная матрица, a через $i{m}^A$ и $\ker A$ — образ (множество значений) и ядро (множество нулей) матрицы (линейного оператора) A. Для симметричного билинейного отображения $B:{ℝ}^p\times {ℝ}^p\to {ℝ}^p$ и для заданного $v\in {ℝ}^p$ через $B\left[v\right]$ обозначается матрица линейного оператора $B[v,\cdot ].$ Наконец, запись $O\left(t\right)u$ используется для произведения вида $Au$ с матрицей $A=O\left(t\right)$ в случаях, когда другая специфика этой матрицы не имеет значения.

1. Постановка задачи и ньютоновские методы

Рассматриваются ньютоновские методы для решения нелинейного уравнения

$\Phi \left(u\right)=0$ (1.1)

с достаточно гладким отображением $\Phi :{ℝ}^p\to {ℝ}^p.$ Такие методы по имеющемуся приближению $u\in {ℝ}^p$ генерируют следующее приближение как $u+v,$ где $v\in {ℝ}^p$ есть решение подзадачи, характеризующей конкретный метод. Например:

• подзадача метода Ньютона (NM; см., например, [3, разд. 2.1.1]) — это линеаризованное уравнение (1.1), т. е.

$\Phi \left(u\right)+{\Phi }^{\text{'}}\left(u\right)v=0;$ (1.2)

• подзадача метода Левенберга–Марквардта (LM; происходит из работ [4, 5]; см. также, например, [6, разд. 10.2], а также обзор современных теорий сходимости в [7]) имеет вид

$\frac{1}{2}\parallel \Phi \left(u\right)+{\Phi }^{\text{'}}\left(u\right)v{\parallel }^2+\frac{1}{2}\sigma \left(u\right)\parallel v{\parallel }^2\to \min ,$ (1.3)

с функцией $\sigma :{ℝ}^p\to {ℝ}_{+},$ определяющей значения параметра регуляризации, что равносильно линейному уравнению

${\left({\Phi }^{\text{'}}\right(u\left)\right)}^{\top }\Phi \left(u\right)+\left({\left({\Phi }^{\text{'}}\right(u\left)\right)}^{\top }{\Phi }^{\text{'}}\right(u)+\sigma (u\left)I\right)v=0;$ (1.4)

• подзадача ньютоновского метода с подзадачами линейного программирования (LPN; [8]) имеет вид

$\gamma \to \min ,\parallel \Phi \left(u\right)+{\Phi }^{\text{'}}\left(u\right)v\parallel \le \gamma \parallel \Phi \left(u\right){\parallel }^2,\parallel v\parallel \le \gamma \parallel \Phi \left(u\right)\parallel ,$ (5)

относительно $(v,\gamma )$ с $\gamma \in ℝ\mathrm{,}$ и если используется $\infty$-норма, то это задача линейного программирования, что и дает название методу.

Как обсуждается ниже, все перечисленные методы могут интерпретироваться как возмущенный метод Ньютона (pNM) с подзадачей

$\Phi \left(u\right)+\left({\Phi }^{\text{'}}\right(u)+\Omega (u\left)\right)v=\omega \left(u\right),$ (6)

где отображения $\Omega :{ℝ}^p\to {ℝ}^{p\times p}$ и $\omega :{ℝ}^p\to {ℝ}^p$ характеризуют различные возмущения базового NM. Свойствах базового NM и pNM вблизи неособых решений излагаются в [3, разд. 2.1.1]. Здесь решение $\overline{u}$ уравнения (1.1) называется неособым, если ${\Phi }^{\text{'}}\left(\overline{u}\right)$ — невырожденная матрица, и особым в противном случае.

Данная работа посвящена именно случаю особого решения $\overline{u}\mathrm{.}$ Ключевую роль при этом играет условие 2-регулярности отображения $\Phi$ в точке $\overline{u}$ по направлению $v\in {ℝ}^p$ [9, разд. 1.3.3], сводящееся в рассматриваемом здесь случае отображения из ${ℝ}^p$ в ${ℝ}^p$ к невырожденности линейного оператора $u\mapsto \Pi {{\Phi }^{\text{'}}}^{\text{'}}\left(\overline{u}\right)[v,u]:\ker {\Phi }^{\text{'}}\left(\overline{u}\right)\to {\left(im{\Phi }^{\text{'}}\right(\overline{u}\left)\right)}^{\perp },$ где $\Pi$ — оператор ортогонального проектирования на ${\left(im{\Phi }^{\text{'}}\right(\overline{u}\left)\right)}^{\perp }.$ Важное свойство, которое легко проверить, состоит в том, что из 2-регулярности $\Phi$ в точке $\overline{u}$ по направлению v вытекает существование $\epsilon >0$ и $\delta >0$ таких, что для $u\in K_{\epsilon ,\delta }(\overline{u},v)$ матрица ${\Phi }^{\text{'}}\left(u\right)$ невырождена, причем

${\left({\Phi }^{\text{'}}\right(u\left)\right)}^{-1}=O(\parallel u-\overline{u}{\parallel }^{-1})$ (1.7)

при $u\to \overline{u}\mathrm{.}$ С использованием этого факта, следующий результат о локальной сходимости pNM с линейной скоростью был получен в [10, теорема 1] (см. также [11, теорема 3.1], где используются более слабые предположения гладкости), и является развитием результатов, полученных [1, 2] для базового NM. Для всякого  введем однозначное разложение $u=u_1+u_2,$ $u_1\in {\left(\ker {\Phi }^{\text{'}}\right(\overline{u}\left)\right)}^{\perp },$ $u_2\in \ker {\Phi }^{\text{'}}\left(\overline{u}\right).$

Теорема 1.1. Пусть $\Phi {ℝ}^p\to {ℝ}^p$ дважды дифференцируемо вблизи $\overline{u}\in {ℝ}^p\mathrm{,}$ и его вторая производная липшицева относительно $\overline{u}\mathrm{,}$ т. е.

${{\Phi }^{\text{'}}}^{\text{'}}\left(u\right)-{{\Phi }^{\text{'}}}^{\text{'}}\left(\overline{u}\right)=O(\parallel u-\overline{u}\parallel )$

при $u\to \overline{u}\mathrm{.}$ Пусть $\overline{u}$ является решением уравнения (1.1) причем $\Phi$ 2-регулярно в точке $\overline{u}$ по направлению $\overline{v}\in \ker {\Phi }^{\text{'}}\left(\overline{u}\right)\backslash \left\{0\right\}.$ Пусть для $\Omega {ℝ}^p\to {ℝ}^{p\times p}$ и $\omega {ℝ}^p\to {ℝ}^p$ существует $\delta >0$ такое, что оценки

$\Omega \left(u\right)=O(\parallel u-\overline{u}\parallel ),\omega \left(u\right)=O(\parallel u-\overline{u}{\parallel }^2),$

$\Pi \Omega \left(u\right)=O(\parallel u_1-{\overline{u}}_1\parallel )+O(\parallel u-\overline{u}{\parallel }^2)$

$\Pi \omega \left(u\right)=O(\parallel u-\overline{u}\parallel \parallel u_1-{\overline{u}}_1\parallel )+O(\parallel u-\overline{u}{\parallel }^3)$

выполняются для $u\in K_{\epsilon ,\delta }(\overline{u};\overline{v})$ при $\epsilon \to 0+\mathrm{.}$

Тогда для любых $\widehat{\epsilon }>0$ и $\widehat{\delta }>0$ найдутся $\epsilon =\epsilon \left(\overline{v}\right)>0$ и $\delta =\delta \left(\overline{v}\right)>0$ такие, что для всякой начальной точки $u^0\in K_{\epsilon ,\delta }(\overline{u};\overline{v})$ существует единственная последовательность $\left\{u^k\right\}\subset {ℝ}^p$ такая, что для каждого k выполняется $u^{k+1}=u^k+v^k,$ где $v^k$ удовлетворяет (1.6), и для этой последовательности и для каждого k имеет место $u_2^k\cancel{=}{\overline{u}}_2\mathrm{,}$ $u^k\in K_{\widehat{\epsilon },\widehat{\delta }}(\overline{u};\overline{v}),$  сходится к $\overline{u}\mathrm{,}$ $\{\parallel u^k-\overline{u}\parallel \}$ сходится к нулю монотонно,

$\frac{\parallel u_1^{k+1}-{\overline{u}}_1\parallel }{\parallel u_2^{k+1}-{\overline{u}}_2\parallel }=O(\parallel u^k-\overline{u}\parallel )$

при $k\to \infty \mathrm{,}$ и

$\underset{k\to \infty }{lim}\frac{\parallel u_2^{k+1}-{\overline{u}}_2\parallel }{\parallel u_2^k-{\overline{u}}_2\parallel }=\frac{1}{2}.$

Характер сходимости в теореме 1.1 служит основой для техник ускорения сходимости ньютоновских методов к особым решениям [1, 2], одной из которых является экстраполяция, которая в своем простейшем варианте состоит в генерировании, параллельно с основной последовательностью $\left\{u^k\right\},$ вспомогательной последовательности $\left\{{\widehat{u}}^k\right\}$ посредством удвоения шага (p)NM:

${\widehat{u}}^{k+1}=u^k+2v^k.$ (8)

Как показано в [2, Theorem 4.1], для базового NM можно ожидать $\left\{{\widehat{u}}^k\right\}$ линейной сходимости с асимптотическим общим частным 1/4 в то время как для основной итерационной последовательности $\left\{u^k\right\}$ оценки в теореме 1.1 гарантируют такое свойство только с 1/2. Заметим, что (1.8) не связано ни с какими дополнительными вычислительными затратами, кроме одного лишнего вычисления $\Phi$ на каждой итерации, что требуется для оценки качества полученного приближения ${\widehat{u}}^{k+1}\mathrm{.}$ Кроме того, экстраполяция никак не влияет на основную последовательность $\left\{u^k\right\},$ и поэтому легко комбинируется с различными реализациями pNM, и в частности, с глобализациями сходимости методов этого класса, в том числе с рассматриваемыми в разд. 4. алгоритмами, использующими одномерный поиск для адаптивного выбора параметра длины шага. Однако, эффект от экстраполяции можно ожидать только при асимптотическом принятии полного ньютоновского шага, поскольку данная процедура основана на характере сходимости именно полношагового pNM; см. [12, 13].

Возвращаясь к базовому pNM, для метода LM, использующего $\sigma \left(u\right)=O(\parallel \Phi (u\left){\parallel }^{\theta }\right)$ с фиксированным $\theta >\mathrm{0,}$ на основе (7) в [10, разд. 3.1] показано, что при 2-регулярноcти $\Phi$ в точке $\overline{u}$ по направлению $\overline{v}\in \ker {\Phi }^{\text{'}}\left(\overline{u}\right)\backslash \left\{0\right\}$ найдутся $\epsilon >0$ и $\delta >0$ такие, что для $u\in K_{\epsilon ,\delta }(\overline{u},\overline{v})$ (1.4) записывается в виде (1.6) с $\omega (\cdot )\equiv 0$ и 1.

$\Omega \left(u\right)=\sigma \left(u\right){\left({\left({\Phi }^{\text{'}}\right(u\left)\right)}^{-1}\right)}^{\top }=O(\parallel u-\overline{u}{\parallel }^{-1}\parallel \Phi (u\left){\parallel }^{\theta }\right)$ (9)

при $u\to \overline{u}\mathrm{.}$ В [14, разд. 3.2] предложена альтернативная интерпретация подзадачи (1.3) метода LM как подзадачи (1.6) pNM, с $\Omega (\cdot )\equiv 0$ и соответствующей оценкой на $\omega \mathrm{.}$ В отличие от приведенной выше интерпретации, явно основанной на форме (1.4) итерационной подзадачи метода LM, альтернативная интерпретация применима и в случае уравнения с дополнительным ограничением, но получаемые при этом оценки несколько слабее и не подходят для использования в разд. 3 ниже.

Для метода LPN, снова с использованием (1.7), в [10, разд. 3.2] показано, что в тех же предположениях найдутся $\epsilon >0$ и $\delta >0$ такие, что для $u\in K_{\epsilon ,\delta }(\overline{u},\overline{v})$ любое решение подзадачи (1.5) удовлетворяет (1.6) с $\Omega (\cdot )\equiv 0$ и некоторым $\omega \mathrm{,}$ для которого

$\parallel \omega \left(u\right)\parallel \le \gamma \left(u\right)\parallel \Phi \left(u\right){\parallel }^2=O(\parallel u-\overline{u}\parallel \parallel \Phi (u)\parallel )$ (10)

при $u\to \overline{u}\mathrm{,}$ где $\gamma \left(u\right)$ — оптимальное значение задачи (1.5).

2. Сравнительный эффект от экстраполяции

Начнем с простого примера, хорошо демонстрирующего эффекты, о которых идет речь в этом разделе.

Пример 2.1. Пусть $p=\mathrm{1,}$ $\Phi \left(u\right)=u^2.$ Единственным решением уравнения (1.1) является $\overline{u}=\mathrm{0,}$ оно является особым, и $\Phi$ 2-регулярно в точке 0 по любому ненулевому направлению.

Итерационное уравнение (1.2) NM имеет вид $u^2+2uv=\mathrm{0,}$ и при $u\cancel{=}0$ его единственным решением является $v=-u/2.$ Поэтому $u+2v=\mathrm{0,}$ т.е. NM с экстраполяцией за один шаг попадает в точное решение из любого начального приближения.

Далее, итерационное уравнение (1.4) метода LM имеет вид $2u^3+(4u^2+\sigma (u\left)\right)v=\mathrm{0,}$ и при $u\cancel{=}0$ его единственным решением является $v=-2u^3/(4u^2+\sigma (u\left)\right).$ Если $\sigma \left(u\right)=|u{|}^{2\theta },$ то

$u+2v=\frac{\sigma \left(u\right)u}{(4u^2+\sigma (u\left)\right)}=\frac{|u{|}^{2\theta +1}}{(4u^2+|u{|}^{2\theta })},$

что при $\theta >1$ дает оценки

$u+2v=\frac{1}{4}|u{|}^{2\theta -1}+o(\left|u{|}^{2\theta -1}\right)=O\left(\right|u{|}^{2\theta -1}),$

что при $\theta \ge 2$ превращается в оценку $u+2v=O\left(\right|u{|}^3).$

Наконец, итерационная подзадача (1.5) метода LPN имеет вид

$\gamma \to \min ,-\gamma u^4\le u^2+2uv\le \gamma u^4,-\gamma u^2\le v\le \gamma u^2,$

и, скажем, при $u>0$ единственным решением этой задачи линейного программирования является $(v,\gamma )=(-u/(2+u),1/(2u+u^2\left)\right).$ Тогда $u+2v=u^2/(2+u),$ что дает оценки

$u+2v=\frac{1}{2}{u}^2+o\left(u^2\right)=O\left(u^2\right).$

Функция в примере 2.1 является квадратичной скалярной функцией одной переменной. Обратимся теперь к случаю общего линейно-квадратичного отображения $\Phi$ пусть

$\Phi \left(u\right)=Au+\frac{1}{2}B[u,u],$ (2.1)

где A — $(n\times n)$-матрица, а $B{ℝ}^p\times {ℝ}^p\to {ℝ}^p$ — симметричное билинейное отображение. Тогда ${\Phi }^{\text{'}}\left(u\right)=A+B\left[u\right],$ и уравнение (1.1) имеет решение $\overline{u}=\mathrm{0,}$ в котором ${\Phi }^{\text{'}}\left(0\right)=A,$ и, соответственно, это решение является особым, если матрица A вырождена. Условие 2-регулярности $\Phi$ в точке 0 по направлению $v\in {ℝ}^p$ при этом состоит в невырожденности линейного оператора $\widehat{\mathcal{B}}\left(v\right):\ker A\to {\left(imA\right)}^{\perp },$ $\widehat{\mathcal{B}}\left(v\right)u=\Pi B[v,u],$ где $\Pi$ — оператор ортогонального проектирования на ${\left(imA\right)}^{\perp }.$

Итерационная подзадача (1.2) NM принимает вид

$(A+B[u\left]\right)v=-Au-\frac{1}{2}B[u,u].$ (2.2)

Если $u\in \ker A\mathrm{,}$ то правую часть (2.2) можно заменить на $-(A+B[u\left]\right)u/\mathrm{2,}$ и тогда уравнение (2.2) имеет решение $v=-u/\mathrm{2,}$ причем это решение единственно в случае невырожденности $A+B\left[u\right].$ При этом $u+2v=\mathrm{0,}$ т. е. экстраполяция немедленно дает точное решение. В чисто квадратичном случае, когда $A=\mathrm{0,}$ условие $u\in \ker A$ выполняется автоматически, и значит, экстраполяция дает точное решение за один шаг из любого начального приближения $u\mathrm{,}$ в котором матрица $B\left[u\right]$ невырождена.

Далее, итерационная подзадача (1.6) pNM принимает вид

$(A+B[u]+\Omega (u\left)\right)v=-Au-\frac{1}{2}B[u,u]+\omega \left(u\right),$ (2.3)

и при $u\in \ker A$ правую часть можно заменить на $-(A+B[u]+\Omega (u\left)\right)u-\Omega \left(u\right)u)/2+\omega (u).$ Если матрица $A+B\left[u\right]+\Omega \left(u\right)$ обратима, то единственным решением (2.3) является

$v=-\frac{1}{2}u+{(A+B[u]+\Omega (u\left)\right)}^{-1}\left(\frac{1}{2}\Omega \left(u\right)u+\omega \left(u\right)\right).$ (2.4)

При $u\in \ker A$ имеет место оценка $\Phi \left(u\right)=O(\parallel u{\parallel }^2),$ и, предполагая 2-регулярность $\Phi$ в точке 0 по направлению $\overline{v}\in \ker A\backslash \left\{0\right\},$ для метода LM с $\theta >\mathrm{1,}$ в силу (1.9), для $u\in K_{\epsilon ,\delta }\left(\mathrm{0,}\overline{v}\right)$ при достаточно малых $\epsilon >0$ и $\delta >0$ имеем

$\Omega \left(u\right)=O(\parallel u{\parallel }^{2\theta -1}),$ (2.5)

и матрица $A+B\left[u\right]+\Omega \left(u\right)$ обратима, причем, аналогично (1.7),

${(A+B[u]+\Omega (u\left)\right)}^{-1}=O(\parallel u{\parallel }^{-1}).$ (2.6)

Тогда из (2.4) (и того, что $\omega (\cdot )\equiv 0$) получаем

$v=-\frac{1}{2}u+O(\parallel u{\parallel }^{2\theta -1})$ (2.7)

при $u\to 0.$ Поэтому экстраполяция дает приближение

$u+2v=O(\parallel u{\parallel }^{2\theta -1}),$ (2.8)

и, в частности, при $\theta \ge 2$ имеет место оценка

$u+2v=O(\parallel u{\parallel }^3)$

при $u\to 0.$

Для метода LPN аналогичным образом при $u\in \ker A$ из (10) выводится оценка

$\omega \left(u\right)=O(\parallel u{\parallel }^3),$ (2.9)

и в силу (2.4), (2.6) (и того, что $\Omega (\cdot )\equiv 0$)

$v=-\frac{1}{2}u+O(\parallel u{\parallel }^2),$

и следовательно,

$u+2v=O(\parallel u{\parallel }^2)$ (2.10)

при $u\to 0.$

Тем самым доказано

Предложение 2.1. Пусть A — $(n\times n)$-матрица, а $B:{ℝ}^p\times {ℝ}^p\to {ℝ}^p$ — симметричное билинейное отображение, и пусть для некоторого $\overline{v}\in \ker A\backslash \left\{0\right\}$ линейный оператор $\widehat{\mathcal{B}}\left(\overline{v}\right)$ невырожден.

Тогда найдутся $\epsilon >0$ и $\delta >0$ такие, что для $u\in K_{\epsilon ,\delta }\left(\mathrm{0,}\overline{v}\right)\cap \ker A,$ для уравнения (1.1) с линейно-квадратичным отображением $\Phi$ из (11) справедливо следующее: 

• NM с подзадачей (1.2) однозначно определяет направление $v\mathrm{,}$ и для него экстраполяция дает $u+2v=0;$
• метод LM с подзадачей (1.3), использующей $\sigma \left(u\right)=O(\parallel \Phi (u\left){\parallel }^{\theta }\right)$ с фиксированным $\theta >\mathrm{1,}$ однозначно определяет направление v, и для него экстраполяция дает приближение $u+2v\mathrm{,}$ для которого выполняется оценка (2.8);
• метод LPN с подзадачей (1.5) определяет некоторое направление v, и для всякого такого направления v, экстраполяция дает приближение $u+2v\mathrm{,}$ для которого выполняется оценка (2.8).

В частности, при $\theta >1$ оценка (2.8) для метода LM лучше, чем оценка (2.10) для метода LPN, и это может служить объяснением большего эффекта экстраполяции для метода LM, чем для LPN.

В общем линейно-квадратичном случае (когда $A\cancel{=}\mathrm{0,}$ и, соответственно, выполнение $u\in \ker A$ не является автоматическим), и даже просто в случае общего достаточно гладкого $\Phi \mathrm{,}$ согласно [15, теорема 6.1] (см. также [13, лемма 7]), за один шаг NM «приближенно» попадает в $\ker A$ из асимптотически плотного звездного множества начальных точек $u^0\in {ℝ}^p\mathrm{,}$ а именно, попадает в $K_{\epsilon \mathrm{,}\delta }\overline{v}$ с соответствующим $\overline{v}=\overline{v}(u^0/\parallel u^0\parallel ),$при любых наперед заданных $\epsilon >0$ и $\delta >0.$ (Точного попадания в $\ker A$ может не быть даже для NM, даже в линейно-квадратичном случае, что демонстрируется примером 2.2 ниже, в котором $\mathrm{dimker}A=1;$ аналогичный эффект при $\mathrm{dimker}A=2$ наблюдается в тесте Misc 10 в разд. 4.) Однако, как указано в комментариях в конце [10, разд. 2], гарантировать указанное свойство для общего pNM не удается, если только не накладывать жесткие предположения на  глобальное (не только вблизи ядра) поведение $\Omega$ и $\omega \mathrm{,}$ а именно, $\Pi \Omega \left(u\right)=O(\parallel u{\parallel }^2)$ и $\Pi \omega \left(u\right)=O(\parallel u{\parallel }^3)$ для u из рассматриваемого звездного множества начальных приближений. Для метода LM эти предположения сводятся к оценке $\sigma \left(u\right)=O(\parallel u{\parallel }^3),$ что можно гарантировать только при $\theta \ge 3.$ Для метода LPN гарантировать требуемое можно при выполнении оценки $\gamma \left(u\right)=O(\parallel u\parallel ),$ а эта оценка ниоткуда не следует. Тем не менее, на практике попадание LM и LPN в $K_{\epsilon ,\delta }\left(\mathrm{0,}\overline{v}\right)$ с достаточно малыми $\epsilon >0$ и $\delta >0$ обычно наблюдается, если не за одну, то за небольшое количество начальных итераций.

 

Рис. 1. Пример 1: области попадания одного шага NM в ${\mathit{K}}_{\mathbf{\epsilon }\mathbf{,}\mathbf{\delta }}\mathbf{(}\mathbf{0}\mathbf{,}\overline{\mathbf{v}}\mathbf{)}$ с $\mathit{\epsilon }\mathbf{=}\mathbf{0}\mathbf{.}\mathbf{1}$

 

Пример 2.2. Пусть $p=\mathrm{2,}$ $\Phi \left(u\right)=(u_1+u_2^2/\mathrm{2,}{u}_1{u}_2+u_2^2/2).$ Единственным решением уравнения (1.1) является $\overline{u}=\mathrm{0,}$ оно является особым, $\ker {\Phi }^{\text{'}}\left(0\right)$ натянуто на $\overline{v}=\left(\mathrm{0,}1\right),$ и $\Phi$ 2-регулярно в точке 0 по направлению v.

На рис. 1 показана линия, задаваемая равенством $\det {\Phi }^{\text{'}}\left(u\right)=\mathrm{0,}$ а также области начальных точек, один шаг рассматриваемых методов из которых приводит в $K_{\epsilon ,\delta }\left(\mathrm{0,}\overline{v}\right)$ с $\epsilon =\mathrm{0.1,}$ при двух разных значениях $\delta \mathrm{.}$ Эти области уменьшаются при уменьшении $\delta ,$ но, по-видимому, остаются асимптотически плотными в точке 0. 

В соответствии со сказанным выше, далее будем рассматривать случай, когда $u^0\in K_{\epsilon ,\delta }\left(\mathrm{0,}\overline{v}\right).$ По-прежнему будем использовать разложение всякого $u\in {ℝ}^p$ в сумму $u=u_1+u_2,$ где в линейно-квадратичном случае $u_1\in {\left(\ker A\right)}^{\perp },$ $u_2\in \ker A\mathrm{.}$ 

Согласно [11, лемма 3.1] для общего (не обязательно линейно-квадратично) случая, и даже всего лишь при сильной полугладкости ${\Phi }^{\text{'}}$ в решении $\overline{u}=\mathrm{0,}$ для $u\in K_{\overline{\epsilon },\overline{\delta }}\left(\mathrm{0,}\overline{v}\right)$ при малых $\overline{\epsilon }>0$ и $\overline{\delta }>0$ направление v NM однозначно определяется и удовлетворяет оценкам

$u_1+v_1=O(\parallel u_1\parallel \parallel u\parallel )+O(\parallel u{\parallel }^3),$ (2.11)

$u_2+v_2=\frac{1}{2}{u}_2+O(\parallel u_1\parallel )+O(\parallel u{\parallel }^2)$ (2.12)

при $u\to 0.$ При этом $\parallel u+v\parallel \ge \gamma \parallel u\parallel$ с некоторым фиксированным $\gamma >\mathrm{0,}$ и, поскольку  $\parallel u_1\parallel =O(\parallel u\parallel ),$

$\frac{\parallel u_1+v_1\parallel }{\parallel u+v{\parallel }^2}=O\left(\frac{\parallel u_1+v_1\parallel }{\parallel u{\parallel }^2}\right)=O\left(\frac{\parallel u_1\parallel }{\parallel u\parallel }\right)+O(\parallel u\parallel )=O\left(1\right),$ (2.13)

т. е. $u_1+v_1=O(\parallel u+v{\parallel }^2).$ Далее, согласно [11, теорема 3.1], если $\epsilon >0$ и $\delta >0$ достаточно малы, то для любого $u^0\in K_{\epsilon ,\delta }\left(\mathrm{0,}\overline{v}\right)$ последовательность $\left\{u^k\right\}$ с шагами $v^k$ NM однозначно определяется, сходится к 0, причем $u^k\in K_{\overline{\epsilon },\overline{\delta }}\left(\mathrm{0,}\overline{v}\right)$ для любого k, откуда и из оценки (2.13) следует, что $u_1^k=O(\parallel u^k{\parallel }^2)$ при $k\to \infty \mathrm{.}$ Тогда, снова привлекая (2.11),

$\frac{\parallel u_1^k+v_1^k\parallel }{\parallel u^k+v^k{\parallel }^3}=O\left(\frac{\parallel u_1^k+v_1^k\parallel }{\parallel u^k{\parallel }^3}\right)=O\left(\frac{\parallel u_1^k\parallel }{\parallel u^k{\parallel }^2}\right)+O\left(1\right)=O\left(1\right),$

т. е. имеет место оценка

$u_1^k=O(\parallel u^k{\parallel }^3)$ (2.14)

при $k\to \infty \mathrm{,}$ которая согласуется с [2, теорема 2.1, (iii)], где эта оценка приведена без доказательства.

Для pNM, согласно [11, лемма 3.1], оценки (2.11), (2.12) заменяются на

$u_1+v_1=O(\parallel u_1\parallel \parallel u\parallel )+O(\parallel u\parallel \parallel \Omega (u)\parallel )+O(\parallel \omega (u)\parallel )+O(\parallel u{\parallel }^3),$ (2.15)

$u_2+v_2=\frac{1}{2}{u}_2+O(\parallel u_1\parallel )+O(\parallel \Pi \Omega (u)\parallel )+O(\parallel u{\parallel }^{-1}\parallel \Pi \omega (u)\parallel )+O(\parallel u{\parallel }^2)$ (2.16)

при $u\to 0.$ Для метода LM, согласно (2.5), поскольку в этом случае $\omega (\cdot )\equiv \mathrm{0,}$ эти оценки принимают вид

$u_1+v_1=O(\parallel u_1\parallel \parallel u\parallel )+O(\parallel u{\parallel }^{2\theta })+O(\parallel u{\parallel }^3),$

$u_2+v_2=\frac{1}{2}{u}_2+O(\parallel u_1\parallel )+O(\parallel u{\parallel }^{2\theta -1})+O(\parallel u{\parallel }^2),$

что при $\theta \ge 3/2$ дает (2.11), (2.12). Для метода LPN, согласно (2.9), поскольку в этом случае $\Omega (\cdot )\equiv \mathrm{0,}$ оценки (2.15), (2.16) сразу принимают вид (2.11), (2.12). Значит, и для этих методов в указанных случаях проходят рассуждения, проведенные выше для NM, и приводящие к оценке (2.14).

Следующая лемма уточняет в линейно-квадратичном случае оценку (1.7).

Лемма 2.1. Пусть A — $(n\times n)$-матрица, а $B:{ℝ}^p\times {ℝ}^p\to {ℝ}^p$ — симметричное билинейное отображение, и пусть для некоторого $\overline{v}\in {ℝ}^p$ линейный оператор $\widehat{\mathcal{B}}\left(\overline{v}\right)$ невырожден.

Тогда найдутся $\overline{\epsilon }>0$ и $\overline{\delta }>0$ такие, что для $u\in K_{\overline{\epsilon },\overline{\delta }}\left(\mathrm{0,}\overline{v}\right)$ матрица $A+B\left[u\right]$ невырождена, причем для $w\in imA$ справедлива оценка

${(A+B[u\left]\right)}^{-1}w=O(\parallel w\parallel )$ (2.17)

при $w\to 0.$

Доказательство. Для заданных $u\in {ℝ}^p$ и $w\in imA\mathrm{,}$ умножая левую и правую части уравнения

$(A+B[u\left]\right)v=w$ (2.18)

на $I-\Pi$ и $\Pi \mathrm{,}$ получаем уравнения

$(A+(I-\Pi \left)B\right[u\left]\right)v_1=w-(I-\Pi )B\left[u\right]v_2,\Pi B\left[u\right]v_2=-\Pi B\left[u\right]v_1.$ (2.19)

Обозначим через $\widehat{A}:{\left(\ker A\right)}^{\perp }\to imA$ сужение линейного оператора A на ${\left(\ker A\right)}^{\perp },$ а через $\mathcal{A}\left(u\right):{\left(\ker A\right)}^{\perp }\to imA$ — сужение линейного оператора $A+(I-\Pi )B\left[u\right]$ на ${\left(\ker A\right)}^{\perp }.$ Тогда

$\mathcal{A}\left(u\right)=\widehat{A}+O(\parallel u\parallel ),$

причем оператор $\widehat{A}$ обратим, откуда и, например, из [3, лемма A.6] следует, что для $u$ достаточно близких к 0 оператор $\mathcal{A}\left(u\right)$ тоже обратим, причем

${\left(\mathcal{A}\right(u\left)\right)}^{-1}={\widehat{A}}^{-1}+O(\parallel u\parallel )$

при $u\to 0.$ Тогда первое уравнение в (2.19) при любом $v_2\in {ℝ}^p$ имеет единственное решение

$v_1={\left(\mathcal{A}\right(u\left)\right)}^{-1}(w-(I-\Pi \left)B\right[u\left]v_2\right)=O(\parallel w\parallel )+O(\parallel u\parallel )v_2$ (2.20)

при $u\to 0$ и $w\to 0.$

Подставляя выражение из (2.20) во второе уравнение в (2.19), получаем

$\mathcal{B}\left(u\right)v_2=\Pi O(\parallel u\parallel \parallel w\parallel ),$ (2.21)

где линейный оператор $\mathcal{B}\left(u\right):\ker A\to {\left(imA\right)}^{\perp }$ в левой части имеет вид

$\mathcal{B}\left(u\right)=\widehat{\mathcal{B}}\left(u\right)+O(\parallel u{\parallel }^2).$

Из невырожденности $\widehat{\mathcal{B}}\left(\overline{v}\right),$ снова применяя [3, лемма A.6], получаем существование $\overline{\epsilon }>0$и $\overline{\delta }>0$ таких, что для $u\in K_{\overline{\epsilon },\overline{\delta }}\left(\mathrm{0,}\overline{v}\right)$ оператор $\mathcal{B}\left(u\right)$ невырожден, причем

${\left(\mathcal{B}\right(u\left)\right)}^{-1}={\left(\widehat{\mathcal{B}}\right(u\left)\right)}^{-1}+O\left(1\right)=O(\parallel u{\parallel }^{-1})$

при $u\to 0.$ Отсюда следует, что (2.21), а значит, и второе уравнение в (2.19) имеет единственное решение

$v_2=O(\parallel w\parallel )$ (2.22)

при $u\to 0$ и $w\to 0.$

Тем самым показано, что для $u\in K_{\overline{\epsilon },\overline{\delta }}\left(\mathrm{0,}\overline{v}\right)$ уравнение (2.18) имеет единственное решение $v=v_1+v_2,$ т. е. матрица $A+B\left[u\right]$ невырождена, причем из (2.20) и (2.22) вытекает оценка (2.17).

Для NM, согласно (2.2),

$(A+B[u\left]\right)v=-\frac{1}{2}(A+B[u\left]\right)u-\frac{1}{2}Au,$

и в предположениях леммы 2.1 для $u\in K_{\overline{\epsilon },\overline{\delta }}\left(\mathrm{0,}\overline{v}\right)$

$v=-\frac{1}{2}u-\frac{1}{2}{(A+B[u\left]\right)}^{-1}Au=-\frac{1}{2}u+O(\parallel Au\parallel )=-\frac{1}{2}u+O(\parallel u_1\parallel )$  

при $u\to 0.$ Поэтому, с учетом (2.14), если $u^0\in K_{\epsilon ,\delta }\left(\mathrm{0,}\overline{v}\right)$ с достаточно малыми $\epsilon >0$ и $\delta >\mathrm{0,}$то

$u^k+2v^k=O(\parallel u_1^k\parallel )=O(\parallel u^k{\parallel }^3)$

при $k\to \infty \mathrm{.}$

Для метода LM, аналогично построениям выше, для $u\in K_{\overline{\epsilon },\overline{\delta }}\left(\mathrm{0,}\overline{v}\right)$ выводится следующее обобщение оценки (2.7):

$v=-\frac{1}{2}u+O(\parallel u{\parallel }^{2\theta -1})+O(\parallel Au\parallel )=-\frac{1}{2}u+O(\parallel u{\parallel }^{2\theta -1})+O(\parallel u_1\parallel )$

при $u\to \mathrm{0,}$ откуда и из (2.14) вытекает оценка

$u^k+2v^k=O(\parallel u^k{\parallel }^{2\theta -1})+O(\parallel u_1^k\parallel )$

при $k\to \infty \mathrm{.}$ Тогда, согласно (2.14),

$u^k+2v^k=O(\parallel u^k{\parallel }^{\min \{2\theta -\mathrm{1,}3\}})$ (2.23)

и, в частности, при $\theta \ge 2$ 

$u^k+2v^k=O(\parallel u^k{\parallel }^3)$ (2.24)

при $k\to \infty \mathrm{,}$ как и для NM.

Наконец, для метода LPN аналогичным образом для $u\in K_{\overline{\epsilon },\overline{\delta }}\left(\mathrm{0,}\overline{v}\right)$ выводится оценка

$v=-\frac{1}{2}u+O(\parallel u{\parallel }^2)+O(\parallel Au\parallel )=-\frac{1}{2}u+O(\parallel u{\parallel }^2)+O(\parallel u_1\parallel )$

при $u\to \mathrm{0,}$ а значит, согласно (2.14), и оценка

$u^k+2v^k=O(\parallel u^k{\parallel }^2)$ (2.25)

при $k\to \infty \mathrm{.}$

Предложение 2.2. В предположениях предложения 2.1 найдутся $\epsilon >0$ и $\delta >0$ такие, что для $u^0\in K_{\epsilon ,\delta }\left(\mathrm{0,}\overline{v}\right),$ для уравнения (2.1) с линейно-квадратичным отображением $\Phi$ из (2.1) справедливо следующее: 

• NM с подзадачей (1.2) однозначно определяет последовательности $\left\{v^k\right\}$ шагов и $\left\{u^k\right\}$ приближений, и для соответствующих экстраполированных приближений выполняется оценка (2.24) при $k\to \infty$;
• метод LM с подзадачей (1.3), использующей $\sigma \left(u\right)=O(\parallel \Phi (u\left){\parallel }^{\theta }\right)$ с фиксированным $\theta \ge 3/\mathrm{2,}$ однозначно определяет последовательности $\left\{v^k\right\}$ шагов и $\left\{u^k\right\}$ приближений, и для соответствующих экстраполированных приближений выполняется оценка (2.23) при $k\to \infty$;
• метод LPN с подзадачей (1.5) определяет последовательности $\left\{v^k\right\}$ шагов и $\left\{u^k\right\}$ приближений, и для любых таких последовательностей, для соответствующих экстраполированных приближений выполняется оценка (2.25) при $k\to \infty \mathrm{.}$ 

Полученные в этом разделе результаты существенно используют линейно-квадратичный характер отображения $\Phi \mathrm{.}$ За пределами линейно-квадратичного случая оценка (2.24) для метода LM с $\theta =\mathrm{2,}$ и даже для NM, может нарушаться, в том числе и при $u\in \ker A\mathrm{,}$ и эффект от экстраполяции для MN и методов LM и LPN при этом может выравниваться, что будет продемонстрировано в следующем разделе.

3. Глобализованные алгоритмы и численные результаты

Алгоритм 1 ниже, глобализующий сходимость NM, использует страховочные шаги градиентного метода для функции $\phi :{ℝ}^p\to ℝ,$ 

$\phi \left(u\right)=\frac{1}{2}\parallel \Phi \left(u\right){\parallel }^2,$ (3.1)

градиент которой имеет вид

${\phi }^{\text{'}}\left(u\right)={\left({\Phi }^{\text{'}}\right(u\left)\right)}^{\top }\Phi \left(u\right).$ (3.2)

Такие шаги используются в тех случаях, когда направление NM либо не удается определить, либо оно оказывается «слишком длинным» в смысле нарушения теста (38). Этот алгоритм был предложен в [12, алгоритм 3.1], где были установлены свойства его глобальной сходимости.

Алгоритм 1. Фиксируем параметры $C>\mathrm{0,}$ $\tau >\mathrm{0,}$ $\rho \in \left(\mathrm{0,}1\right)$ и $\varkappa \in \left(\mathrm{0,}1\right).$ Выбираем $u^0\in {ℝ}^p$ и полагаем $k=0.$  

1. Если $\Phi \left(u^k\right)=\mathrm{0,}$ стоп.
2. Вычисляем $v^k$ как решение линейного уравнения (.12). Если $v^k$ не удается вычислить, или $v^k$ нарушает неравенство

$\parallel v^k\parallel \le \max \{C,1/\parallel \Phi (u^k\left){\parallel }^{\tau }\right\},$ (3.3)

переходим к шагу 4.

3. Полагаем $\alpha =1.$ Если выполняется неравенство

$\parallel \Phi (u^k+\alpha v^k)\parallel \le (1-\rho \alpha )\parallel \Phi \left(u^k\right)\parallel ,$ (3.4)

полагаем ${\alpha }_k=\alpha$ и переходим к шагу 6. В противном случае заменяем a на $\varkappa \alpha$ и проверяем снова неравенство (3.4) до тех пор, пока оно не выполнится, после чего полагаем ${\alpha }_k=\alpha$ и переходим к шагу 6.

4. Полагаем $v^k=-{\phi }^{\text{'}}\left(u^k\right),$ где функция $\phi$ определена в (3.1) (см. (3.2)). Если $v^k=\mathrm{0,}$ стоп.

5. Полагаем $\alpha =1.$ Если выполняется неравенство Армихо

$\phi (u^k+\alpha v^k)\le \phi \left(u^k\right)-\rho \alpha \parallel v^k{\parallel }^2,$ (3.5)

полагаем ${\alpha }_k=\alpha .$ В противном случае заменяем a на $\varkappa \alpha$ и проверяем снова неравенство (3.5) до тех пор, пока оно не выполнится, после чего полагаем ${\alpha }_k=\alpha .$

6. Полагаем $u^{k+1}=u^k+{\alpha }_k{v}^k,$ увеличиваем k на l и переходим к шагу 1.

В отличие от NM, метод LM может быть глобализован «чистым» образом, без каких-либо гибридных элементов и, в частности, без необходимости использования каких-либо страховочных шагов.

Алгоритм 2. Фиксируем параметры $\theta >\mathrm{0,}$ $\rho \in \left(\mathrm{0,}1\right)$ и $\varkappa \in \text{\hspace{0.05em}(0,1)}$ Выбираем $u^0\in {ℝ}^p$ и полагаем $k=0.$  

1. Если $\Phi \left(u^k\right)=\mathrm{0,}$ стоп.
2. Полагаем $\sigma \left(u^k\right)=\parallel \Phi \left(u^k\right){\parallel }^{\theta },$ и вычисляем $v^k$ как решение линейного уравнения (4). Если $v^k=\mathrm{0,}$ стоп.
3. Полагаем $\alpha =1.$ Если выполняется неравенство

$\phi (u^k+\alpha v^k)\le \phi \left(u^k\right)-\frac{1}{2}\rho \sigma \left(u^k\right)\alpha \parallel v^k{\parallel }^2,$ (3.6)

где функция $\phi$ определена в (3.1), полагаем ${\alpha }_k=\alpha .$ В противном случае заменяем a на $\varkappa \alpha$ и проверяем снова неравенство (3.6) до тех пор, пока оно не выполнится, после чего полагаем ${\alpha }_k=\alpha .$

4. Полагаем $u^{k+1}=u^k+{\alpha }_k{v}^k,$ увеличиваем k на l и переходим к шагу 1.

Алгоритм 2 (с $\theta =2$) и анализ его глобальной сходимости были предложены в [16]; см. также [7, алгоритм 6.2, теорема 6.2].

Метод LPN также допускает «чистую» глобализацию сходимости, разработанную в [17, алгоритм 1, теорема 3.1]. Алгоритм 3 ниже основан на одномерном поиске для функции  

$f:{ℝ}^p\to ℝ,$

$f\left(u\right)=\parallel \Phi \left(u\right){\parallel }_{\infty }.$ (3.7)

Алгоритм 3. Фиксируем параметры $\rho \in \left(\mathrm{0,}1\right)$ и $\varkappa \in \text{\hspace{0.05em}(0,1)}$ Выбираем $u^0\in {ℝ}^p$ и полагаем $k=0.$  

1. Если $\Phi \left(u^k\right)=\mathrm{0,}$ стоп.
2. Вычисляем $(v^k,{\gamma }_k)$ как решение задачи линейного программирования (5), в которой используется $\infty$-норма. Если для величины

${\Delta }_k=-f\left(u^k\right)(1-{\gamma }_{k}f(u^k\left)\right),$

где функция f определена в (3.7), выполняется ${\Delta }_k=\mathrm{0,}$ стоп.

3. Полагаем $\alpha =1.$ Если выполняется неравенство

$f(u^k+\alpha v^k)\le f\left(u^k\right)+\rho \alpha {\Delta }_k,$ (3.8)

полагаем ${\alpha }_k=\alpha .$ В противном случае заменяем a на $\varkappa \alpha$ и проверяем снова неравенство (3.8) до тех пор, пока оно не выполнится, после чего полагаем ${\alpha }_k=\alpha .$

4. Полагаем $u^{k+1}=u^k+{\alpha }_k{v}^k,$ увеличиваем k на l и переходим к шагу 1.

В экспериментах, результаты которых представлены ниже, алгоритм 3 был также снабжен всеми предложенными в [17, разд. 5] модификациями, улучшающими его поведение, кроме немонотонного одномерного поиска. В процедурах одномерного поиска использовались параметры $\rho =0.01$и $\kappa =\mathrm{0.5.}$ В алгоритме 1 использовались $C={10}^7$ и $\tau =2.$

Запуск завершался с успехом в случае выполнения для некоторого $k\le 100$ условия

$\parallel \Phi \left(u^k\right)\parallel \le {10}^{-8},$

либо условия

$\min \{\parallel \Phi (u^k)\parallel ,\parallel \Phi (u^{k-1}+2v^{k-1})\parallel \}\le {10}^{-8}$ (3.9)

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

$\parallel {\left({\Phi }^{\text{'}}\right(u^k\left)\right)}^{\top }\Phi \left(u^k\right)\parallel \le {10}^{-20}$

для алгоритмов 1 и 2, и

$\left|{\Delta }_k\right|\le {10}^{-16}$

для алгоритма 3, либо условия $\alpha \parallel v^k\parallel \le {10}^{-16}$ в процессе одномерного поиска на соответствующих шагах алгоритмов.

В алгоритме 2 правило выбора параметра регуляризации было модифицировано следующим образом, исключающим слишком большие значения этого параметра вдали от решений: $\sigma \left(u\right)=\min \{\mathrm{1,}\parallel \Phi (u\left){\parallel }^{\theta }\right\},$ причем использовалось значение $\theta =2.$

Эксперимент проводился на всех задачах из тестового набора MGH [18], для которых известны точные решения, за исключением трех задач с линейными отображениями $\Phi \mathrm{.}$ Некоторые из выбранных таким образом задач имеют больше уравнений, чем переменных, и в таких случаях лишние уравнения удалялись (всегда линейные, за исключением теста  Beale, в котором выбор уравнения для удаления не влияет существенным образом на поведение алгоритмов). Информация о полученном в результате тестовом наборе содержится в [12, таблица 1]. Далее, выбранные задачи трансформировались посредством техники из [19] таким образом, чтобы известное решение $\overline{u}$ оставалось решением модифицированного уравнения, но при этом выполнялось $rank{\Phi }^{\text{'}}\left(\overline{u}\right)=p-1$ (за исключением трех задач Powell singular, Extended Powell singular, и Variably dimensioned, которые имеют известные особые решения в своей оригинальной форме). Задачи  Variably dimensioned и Brown almost-linear использовались в двух версиях: с $p=10$ и с $p=500.$

Запуски сравниваемых алгоритмов осуществлялись из 100 одинаковых случайных начальных точек, равномерно распределенных в $l_{\infty }$-шаре (гиперкубе) радиуса l с центром в решении.

 

Рис. 2. Сравнение на наборе MGH

 

Рис. 3. Сравнение на наборе Misc

 

Результаты представлены в форме так называемых performance profiles, изначально предложенных в [20], и адаптированных в [21] на случай многократных запусков для каждой тестовой задачи, с использованием средних характеристик эффективности и процентов успешных запусков. На рис. 2 для каждого алгоритма показан график функции, значение которой в  по существу соответствует доле задач в наборе, для которых средний показатель эффективности для данного алгоритма был хуже наилучшего не более, чем в t раз. При этом считается, что результат любого неудачного запуска в бесконечное число раз хуже любого другого. Значения функции на правом конце (при больших t) соответствует доле успешных запусков алгоритма, а на левом — доле задач, для которых его средний показатель эффективности был наилучшим. В данном сравнении в роли показателей эффективности использовались количество итераций и время запуска.

Представленные результаты демонстрируют, что эффект от экстраполяции для метода LM слабее, чем для NM, а для метода LPN слабее, чем для LM.

Эксперимент на тестовом наборе MGH был дополнен аналогичным экспериментом на наборе Misc(ellaneous), состоящим из небольших систем нелинейных уравнений с особыми решениями, собранных из различных источников (первые 11 тестов взяты из [22, разд. 5]). А именно, набор состоит из следующих задач.

Misc 1. Это уравнение (1.1) из примера 2.1, с $p=\mathrm{1,}$ $\Phi \left(u\right)=u^2;$ изолированное особое решение $\overline{u}=0.$

Misc 2. Уравнение (1.1) с $p=\mathrm{2,}$ $\Phi \left(u\right)=(u_1^2-u_2^2,u_1{u}_2);$ изолированное особое решение $\overline{u}=\left(\mathrm{0,}0\right).$

Misc 3. Уравнение (1.1) с $p=\mathrm{2,}$ $\Phi \left(u\right)=(u_1^2-u_2^2,u_1{u}_2);$ изолированное особое решение $\overline{u}=\left(\mathrm{0,}0\right).$

Misc 4. Уравнение (1.1) с $p=\mathrm{2,}$$\Phi \left(u\right)=(u_1+u_2,-u_1-u_2+u_1{u}_2);$ изолированное особое решение $\overline{u}=\left(\mathrm{0,}0\right).$

Misc 5. Уравнение (1.1) с $p=\mathrm{2,}$ $\Phi \left(u\right)=(u_1^2,u_2^2);$ изолированное особое решение $\overline{u}=\left(\mathrm{0,}0\right).$

Misc 6. Уравнение (1.1) с $p=\mathrm{2,}$ $\Phi \left(u\right)=\left(2\right(u_1-u_2^2),u_2^2);$ изолированное особое решение $\overline{u}=\left(\mathrm{0,}0\right).$

Misc 7. Уравнение (1.1) с $p=\mathrm{2,}$ $\Phi \left(u\right)=\left(u_1\right(u_1^2+u_2),u_2(1+u_2\left)\right);$ изолированное особое решение $\overline{u}=\left(\mathrm{0,}0\right)$ (имеются также неособые решения $(\mathrm{0,}-1),$ $(\mathrm{1,}-1),$ $(-\mathrm{1,}-1)$).

Misc 8. Уравнение (1.1) с $p=\mathrm{2,}$ $\Phi \left(u\right)=(u_1+u_2^2,3u_1{u}_2/2+u_2^2(1+u_2\left)\right);$ изолированное особое решение $\overline{u}=\left(\mathrm{0,}0\right)$ (имеется также неособое решение $(-\mathrm{4,}2)$).

Misc 9. Уравнение (1.1) с $p=\mathrm{3,}$ $\Phi \left(u\right)=(u_1+u_2+u_3-\mathrm{1,}{u}_1^3/5+u_2^2/2-u_3+u_3^2/2+1/\mathrm{2,}{u}_1+u_2+u_3^2/2-1/2);$ изолированное особое решение $\overline{u}=\left(\mathrm{0,}\mathrm{0,}1\right)$ (имеется также другое особое решение $(-5/\mathrm{2,}5/\mathrm{2,}1)$).

Misc 10. Уравнение (1.1) с $p=\mathrm{3,}$ $\Phi \left(u\right)=(u_1+u_1{u}_2+u_2^2,-2u_1+u_1^2+u_2^2,u_1+u_3^2);$ изолированное особое решение $\overline{u}=\left(\mathrm{0,}\mathrm{0,}0\right).$

Misc 11. (дискретизованное H-уравнение Чандрасехара). Уравнение (1.1) с выбираемым значением p (в экспериментах использовалось $p=5$),

${\Phi }_i\left(u\right)=u_i-\frac{i}{2n}{u}_i\sum _{j=1}^n\frac{u_j}{i+j}-\mathrm{1,}i\in \{\mathrm{1,}\dots ,p\};$

для $p=5$ приближенное изолированное особое решение

$\overline{u}=\left(\mathrm{1.359753,}\mathrm{1.688205,}\mathrm{2.005894,}\mathrm{2.318350,}2.627810\right).$

Misc 12. Уравнение (1.1) с $p=\mathrm{2,}$ $\Phi \left(u\right)=(u_1+a{u}_2^2/\mathrm{2,}{u}_2^2/2)$ с вещественным параметром a (в экспериментах использовалось $a=\sqrt{15}$); изолированное особое решение $\overline{u}=\left(\mathrm{0,}0\right).$

Misc 13. Это уравнение (1.1) из примера 2.2, с $p=\mathrm{2,}$ $\Phi \left(u\right)=(u_1+a{u}_2^2/\mathrm{2,}{u}_1{u}_2+u_2^2/2)$ с вещественным параметром a (в экспериментах использовалось $a=\mathrm{1,}$ как и в примере 1); изолированное особое решение $\overline{u}=\left(\mathrm{0,}0\right)$ (имеется также неособое решение $(-1/(2a),1/a)$).

Misc 14. Уравнение (1.1) с $p=\mathrm{2,}$ $\Phi \left(u\right)=(u_1^2+u_2^3,u_1{u}_2);$ изолированное особое решение $\overline{u}=\left(\mathrm{0,}0\right).$

Misc 15. Уравнение (1.1) с $p=\mathrm{2,}$ $\Phi \left(u\right)=(u_1+u_1{u}_2+u_2^2,u_1^2-2u_1+u_2^2);$ изолированное особое решение $\overline{u}=\left(\mathrm{0,}0\right).$

Misc 16. Уравнение (1.1) с $p=\mathrm{2,}$ $\Phi \left(u\right)=(u_1^2-u_2,u_1^2+u_2^2);$ изолированное особое решение $\overline{u}=\left(\mathrm{0,}0\right).$

Misc 17. Уравнение (1.1) с $p=\mathrm{2,}$ $\Phi \left(u\right)=(u_1^2-u_2^2,3u_1^2-3u_2^2);$ неизолированное решение $\overline{u}=\left(\mathrm{0,}0\right).$

Misc 18. Уравнение (1.1) с числом переменных $p=\mathrm{5,}$ и с числом уравнений, равным 4, неизолированное решение $\overline{u}=\left(\mathrm{1,}\mathrm{1,}\mathrm{0,}\mathrm{0,}0\right).$ $\Phi \left(u\right)=(u_1+u_2+u_3^2+u_4^2+u_5^2-\mathrm{2,}{u}_1-u_2+u_3^2+u_4^2+u_5^2,-u_3^2+u_4^2+u_5^2,u_3^2+u_4^2-u_5^2);$ Для этой задачи NM работает как метод Гаусса–Ньютона, использующий решение минимальной нормы итерационного уравнения (2).

Misc 19. Уравнение (1.1) с $p=3$ из [1, приложение] с параметрами k=1, m=1; изолированное особое решение $\overline{u}=\left(\mathrm{0,}\mathrm{0,}0\right).$

Misc 20. Уравнение (1.1) с p=2, $\Phi \left(u\right)=(u_1^2+u_2^2-\mathrm{2,}{e}^{u_1-1}+u_2^2-2);$ неособое решение $\overline{u}=(\mathrm{1,}-1),$ которое трансформировалось в особое посредством упомянутой выше в связи с набором MGH процедуры из [19] (у исходного уравнения имеется также другое неособое решение (1,1)).

Misc 21. Уравнение (1.1) с $p=\mathrm{2,}$  неособое решение  которое трансформировалось в особое процедурой из [19].

Misc 22. Уравнение (1.1) с $p=\mathrm{2,}$ $\Phi \left(u\right)=(e^{u_1^2+u_2^2}-\mathrm{1,}{u}_1+u_2-\sin (3(u_1+u_2)\left)\right);$ изолированное особое решение $\overline{u}=\left(\mathrm{0,}0\right).$

Misc 23. Уравнение из Misc 22, к решению которого дополнительно применялась процедура из [19].

Misc 24. Уравнение (1.1) с $p=\mathrm{2,}$ $\Phi \left(u\right)=(u_1,1-u_1/(u_1+0.1)+2u_2^2);$ неособое решение $\overline{u}=\left(\mathrm{5000,}10\right),$ которое трансформировалось в особое процедурой из [19].

Misc 25. Уравнение (1.1) с $p=\mathrm{2,}$ $\Phi \left(u\right)=(u_1+u_2^2,2(u_1-1\left)u_2\right);$ неособое решение $\overline{u}=\left(\mathrm{0,}0\right),$ которое трансформировалось в особое процедурой из [?].

Результаты эксперимента на наборе Misc представлены на рис. 2, и выводы из них аналогичны выводам для MGH: эффект от экстраполяции для NM больше, чем для метода LM, а для метода LM больше, чем для метода LPN.

 

Рис. 4. Сравнение на задачах, не являющихся линейно-квадратичными

 

Дальнейший эксперимент был направлен на демонстрацию того, что за пределами линейно-квадратичного случая эффект от экстраполяции для NM и методов LM и LPN выравнивается, что согласуется с рассмотрениями в разд. 3. Для этого из наборов MGH и Misc были отобраны задачи, в которых отображение  не является линейно-квадратичным, а именно, Freudenstein and Roth, Beale, Helical valley, Gulf research and development, Box three-dimensional, Biggs EXP6, Trigonometric, Brown almost-linear; Misc 7–9, = Misc 14, Misc 20, Misc 22–24. Результаты на рис. 4 показывают, что на данной выборке тестов сравниваемые алгоритмы имеют схожий эффект от экстраполяции, хотя для метода LPN эффект все же меньше, чем для других алгоритмов.

 

Рис. 5. Сравнение на полном тестовом наборе

 

На рис. 5 приведенные выше результаты дополнены сравнением трех алгоритмов 1, 2 и 3 между собой. Сравнение по времени здесь не приводится, поскольку метод LPN для рассматриваемых здесь уравнений без дополнительных ограничений неизбежно проигрывает по этому показателю, так как решение его подзадач (задач линейного программирования) является гораздо более затратным, чем решение подзадач NM и метода LM (систем линейных уравнений). Вместо этого приводится сравнение по другой важной характеристике, а именно, по количеству вычислений $\Phi \mathrm{.}$ Важное наблюдение состоит в том, что по этой характеристике методы с экстраполяцией могут быть существенно менее эффективны, чем по количеству итераций, поскольку на каждой итерации требуют по крайней мере одного лишнего вычисления  для проверки условия остановки (3.9).

 

Таблица 1. Процент среднего количества последних полных шагов

Тест	NM	LM	LPN
Rosenbrock	98.67	49.70	70.87
Freudenstein and Roth	99.44	98.95	99.75
Brown badly scaled	32.14	NaN	100.00
Beale	99.08	100.00	99.93
Helical valley	95.77	33.95	49.63
Gulf research and development	24.05	13.64	17.39
Box three-dimensional	90.71	28.46	44.85
Powell singular	100.00	100.00	100.00
Wood	99.23	94.97	98.16
Biggs EXP6	83.20	4.20	6.91
Extended Rosenbrock	93.48	26.68	87.02
Extended Powell singular	100.00	100.00	100.00
Variably dimensioned,  $p=10$	100.00	100.00	99.82
Variably dimensioned,  $p=500$	100.00	100.00	62.91
Trigonometric	NaN	53.00	NaN
Brown almost-linear,  $p=10$	98.24	100.00	95.34
Brown almost-linear,  $p=500$	100.00	98.27	65.41
Misc 1	100.00	100.00	100.00
Misc 2	100.00	100.00	100.00
Misc 3	100.00	100.00	100.00
Misc 4	98.65	100.00	100.00
Misc 5	100.00	100.00	100.00
Misc 6	100.00	100.00	100.00
Misc 7	97.89	99.74	99.18
Misc 8	98.57	100.00	100.00
Misc 9	97.75	100.00	98.77
Misc 10	98.50	100.00	100.00
Misc 11	100.00	100.00	100.00
Misc 12	100.00	100.00	100.00
Misc 13	96.75	99.33	99.50
Misc 14	99.86	100.00	100.00
Misc 15	98.43	100.00	100.00
Misc 16	97.94	100.00	100.00
Misc 17	87.50	100.00	100.00
Misc 18	99.76	100.00	100.00
Misc 19	89.61	83.42	86.06
Misc 20	89.69	55.38	79.58
Misc 21	100.00	NaN	13.21
Misc 22	98.07	100.00	99.07
Misc 23	13.29	18.83	13.12
Misc 24	96.56	52.57	57.29
Misc 25	76.92	73.04	84.40

 

Для алгоритмов с одномерным поиском эффект от экстраполяции во многом определяется асимптотическим принятием полного шага, поскольку эта техника основана на характере сходимости к особым решениям именно полношаговых ньютоновских методов. В этих экспериментах единственный тест, на котором иногда наблюдались случаи неполного последнего шага (для всех трех алгоритмов), — это Misc 23, в котором ${\Phi }^{\text{'}}\left(\overline{u}\right)=\mathrm{0,}$ но 2-регулярности не имеет места ни по какому направлению. Таблица 1 содержит информацию о проценте среднего количества последних полных шагов на один успешный запуск, от общего количества итераций.

Интересное наблюдение состоит в том, что по времени алгоритм 2 показывает лучший результат на примерно 70% задач (по указанным выше причинам эти результаты здесь не приводятся), а по количеству вычислений $\Phi$ и особенно по количеству итераций такой процент гораздо ниже; см. рис. 5а, 5б. Возможное объяснение состоит в том, что по средним значениям принимаемых алгоритмами параметров длины шага алгоритм 2 также оказывается лучшим на сравнимой доле задач, как и по среднему количеству последних полных шагов; см. таблицу 1.

В заключение кратко опишем результаты эксперимента с методом LM, в котором в правиле для параметра регуляризации $\sigma \left(u\right)=\min \{\overline{\sigma },\parallel \Phi (u\left){\parallel }^2\right\}$ использовались разные значения $\overline{\sigma }\mathrm{,}$ a именно, наряду с l, рассматривались значения ${10}^{-1}\mathrm{,}$ ${10}^{-5}$ и ${10}^{-10}\mathrm{.}$ На данных тестовых наборах версии метода LM с меньшими значениями $\overline{\sigma }$ (что в принципе позволяет алгоритму 2 вдали от решений генерировать и принимать более длинные шаги, оказываются несколько более эффективными, но при этом несколько менее робастными, однако разница невелика. Вместе с тем, выбор $\overline{\sigma }=1$значительно проигрывает по времени рассматриваемым альтернативам.

Об авторах

Алексей Феридович Измаилов

ФГБОУ ВО «Московский государственный университет им. М. В. Ломоносова»

Автор, ответственный за переписку.
Email: izmaf@cs.msu.ru
ORCID iD: 0000-0001-9851-0524

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

Россия, 119991, ГСП-2, Москва, Ленинские горы, 1

Евгений Иванович Усков

ФГБОУ ВО «Тамбовский государственный университет им. Г. Р. Державина»

Email: euskov@cs.msu.ru
ORCID iD: 0000-0002-3639-0317

кандидат физико-математических наук, научный сотрудник

Россия, 392000, Тамбов, ул. Интернациональная, 33

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

Дополнительные файлы