Основной идеей является использование итераций Гаусса-Ньютона для решения нелинейных уравнений. Скажите, что вы пытаетесь решить уравнение

В настройке КЭМ вы решаете слабую форму r (u) = 0. Установите как обычно

где x представляет 2-D или 3-D точку. Затем умножьте уравнение на произвольную тестовую функцию _ϕi, интегрируйте в области, и используйте формулу Грина и граничные условия, чтобы получить

который должен храниться для всех индексов i.

Вектор невязок ρ (U) можно легко вычислить как

где матрицы K, M, Q и векторы F и G получаются путем сборки задачи

Предположим, что у вас есть догадка U⁽ⁿ⁾ решения. Если U⁽ⁿ⁾ достаточно близко к точному решению, улучшенному приближению U⁽ⁿ⁺¹⁾ получается путем решения линеаризированной задачи

где $$\alpha$$ - положительное число. (Не обязательно, чтобы ρ (U ) = 0 имели решение, даже если ρ ( u) = 0 имеет.) В этом случае итерация Гаусса-Ньютона имеет тенденцию быть минимизатором остаточного, то есть решения _min U$$\Vert \rho (U)\Vert$$.

Хорошо известно, что для достаточно малых $$\alpha$$

и

называется направлением спуска для $$\Vert \rho (U)\Vert$$, где $$\Vert \cdot \Vert$$ является L 2-нормой. Итерация

где $$\alpha$$ ≤ 1 выбирают как можно большим, так что шаг имеет разумное спуск.

Метод Гаусса-Ньютона является локальным, и сходимость гарантируется только при U⁽⁰⁾ достаточно близко к решению. В целом первое предположение может находиться вне области сходимости. Чтобы улучшить сходимость из плохих начальных догадок, реализуется стратегия демпфирования для выбора α, поиска линии Армиджо-Гольдштейна. Он выбирает самый большой коэффициент демпфирования α из последовательности 1, 1/2, 1/4,. таким образом обеспечивается следующее неравенство:

который гарантирует уменьшение нормы невязки не менее чем на 1 - $$\alpha$$/2. Каждый шаг алгоритма линейного поиска требует оценки остаточного $$\rho \left(U^{(n)}+\alpha p_n\right)$$.

Важной точкой этой стратегии является то, что когда U⁽ⁿ⁾ приближается к решению, затем $$\alpha$$→1 и, таким образом, скорость сходимости увеличивается. Если существует решение ρ (U) = 0, схема в конечном счете восстанавливает квадратичную скорость сходимости стандартной итерации Ньютона.

С предыдущей задачей тесно связан выбор исходного предположения U⁽⁰⁾. По умолчанию решатель устанавливает U⁽⁰⁾ а затем собирает матрицы КЭМ K а также F и вычисляет

Демпфированная итерация Гаусса-Ньютона начинается с U⁽¹⁾, что должно быть лучшим предположением, чем U⁽⁰⁾. Если граничные условия не зависят от u решения, то U⁽¹⁾ удовлетворяет их, даже если U⁽⁰⁾ Кроме того, если уравнение линейное, то U⁽¹⁾ является точным решением КЭМ, и решатель не входит в цикл Гаусса-Ньютона.

Бывают ситуации, когда U⁽⁰⁾ = 0 не имеет смысла, или сходимость невозможна.

В некоторых ситуациях вы можете уже иметь хорошее приближение и нелинейный решатель может быть запущен с ним, избегая режима медленной сходимости. Эта идея используется в генераторе адаптивной сетки. Он вычисляет решение $$\tilde{U}$$ на mesh оценивает ошибку и может уточнить некоторые треугольники. Интерполяция $$\tilde{U}$$ очень хорошее начальное приближение для решения на рафинированном mesh.

В целом точный якобиан

недоступен. Приближение _Jn конечными различиями следующим образом дорого, но допустимо. i-й столбец _Jn может быть аппроксимирован

что подразумевает сборку матриц КЭМ для треугольников, содержащих i точек сетки. Очень простое приближение к _Jn, которое дает фиксированную итерацию точек, также возможно следующим образом. По существу, для заданного U⁽ⁿ⁾, вычислите матрицы КЭМ K и F и установите

Это эквивалентно аппроксимации якобиана с матрицей жесткости. Действительно, с ρ года (U⁽ⁿ⁾) = KU⁽ⁿ⁾ - F, помещая _Jn = K выражения

Во многих случаях темп сходимости медленный, но затраты на каждую итерацию малы.

Нелинейный решатель Partial Differential Equation Toolbox™ производными также обеспечивает компромисс между двумя крайностями. Чтобы вычислить производную от U отображения → KU, продолжите следующим образом. Термин a был опущен для ясности, но снова появляется в конечном результате.

Первый интегральный термин не более чем _Ki,j.

Второй член является «комом», т.е. заменен диагональной матрицей, которая содержит суммы строк. Так как и j ϕj = 1, второй член аппроксимируется

который является i-м компонентом K^(c')U, где K^(c') - матрица жесткости, связанная с коэффициентом ∂ c/ ∂ u а не с c. Те же аргументы могут быть применены к производной от отображения U → MU. Производная от отображения U − F в точности

который является большой матрицей, сопоставленной с коэффициентом ∂ f/ ∂ u. Таким образом, якобиан невязки ρ (U) аппроксимируется

где дифференцирование относительно u, K и M обозначают жесткость и большие матрицы, а их индексы обозначают коэффициенты, относительно которых они собраны. На каждой итерации Гаусса-Ньютона нелинейный решатель собирает матрицы, соответствующие уравнениям

а затем производит приблизительный якобиан. Дифференциации коэффициентов выполняются численно.

В общей настройке эллиптических систем граничные условия добавляются к матрице жесткости, чтобы сформировать полную линейную систему:

где коэффициенты $$\tilde{K}$$ и $$\tilde{F}$$ может зависеть от решения $$\tilde{U}$$. «Ограниченный» подход аппроксимирует производное отображения невязки

Нелинейности граничных условий и зависимости коэффициентов от производных $$\tilde{U}$$ не были правильно линеаризированы этой схемой. Когда такие нелинейности сильны, схема уменьшается до итерации с фиксированной точкой и может сходиться медленно или вообще не сходиться. Когда граничные условия линейны, они не влияют на свойства сходимости схем итерации. В случае Неймана они невидимы (H является пустой матрицей), и в случае Дирихле они просто заявляют, что невязка равна нулю на соответствующих граничных точках.