Матричная прогонка по разностной схеме для УЧП

[kitnone]

Добрый день, есть УЧП, но есть вопрос - как работать с $f(x,y,t)A$?
$iA_t=A_{xx}+A_{yy}+f(x,y,t)A$

К нему составил разностную схему метода покоординатного расщепления, неявная шеститочечная, первый слой вдоль x, второй слой вдоль y:
$i\frac{A^{n+0.5}_{k,l}-A^{n}_{k,l}}{\tau}=\sigma\frac{A^{n+0.5}_{k+1,l}-2A^{n+0.5}_{k,l}+A^{n+0.5}_{k-1,l}}{h^2}+(1-\sigma)\frac{A^{n}_{k+1,l}-2A^{n}_{k,l}+A^{n}_{k-1,l}}{h^2}+f^{n}_{k,l}A^{n}_{k,l}$

$i\frac{A^{n+1}_{k,l}-A^{n+0.5}_{k,l}}{\tau}=\frac{A^{n+1}_{k,l+1}-2A^{n+1}_{k,l}+A^{n+1}_{k,l-1}}{h^2}+(1-\sigma)\frac{A^{n+0.5}_{k,l+1}-2A^{n+0.5}_{k,l}+A^{n+0.5}_{k,l-1}}{h^2}$

Корректно ли она составлена?
Позволительно ли $f^{n}_{k,l}A^{n}_{k,l}$ держать в первом слое? Или лучше раскинуть на оба или, только на второй? Исходил из того что раз $A^{n}$ используется на первой половинке то лучше туда отправить и все слагаемые, содержащие это.

Перегруппировал слои по разным сторонам:
$ \frac{1-\sigma}{h^2}A^{n}_{k+1,l}+(\frac{i}{\tau}-\frac{2(1-\sigma)}{h^2}+f^{n}_{k,l})A^{n}_{k,l}+\frac{1-\sigma}{h^2}A^{n}_{k-1,l}=-\frac{\sigma}{h^2}A^{n+0.5}_{k+1,l}+(\frac{i}{\tau}+\frac{2\sigma}{h^2})A^{n+0.5}_{k,l}-\frac{\sigma}{h^2}A^{n+0.5}_{k-1,l}$
$ \frac{1-\sigma}{h^2}A^{n+0.5}_{k,l+1}+(\frac{i}{\tau}-\frac{2(1-\sigma)}{h^2})A^{n+0.5}_{k,l}+\frac{1-\sigma}{h^2}A^{n+0.5}_{k,l-1}=-\frac{\sigma}{h^2}A^{n+1}_{k,l+1}+(\frac{i}{\tau}+\frac{2\sigma}{h^2})A^{n+1}_{k,l}-\frac{\sigma}{h^2}A^{n+1}_{k,l-1} $

Подхожу к вопросу - при переходе к самим матрицам не понимаю что делать с $f^{n}_{k,l}A^{n}_{k,l}$ - без нее были бы обычные трехдиагональные матрицы вида
$M=\begin{pmatrix} \frac{i}{\tau}-\frac{2(1-\sigma)}{h^2} & \frac{1-\sigma}{h^2} & 0 & 0 & \ldots\\ \frac{1-\sigma}{h^2}& \frac{i}{\tau}-\frac{2(1-\sigma)}{h^2} & \frac{1-\sigma}{h^2} & 0 & \ldots\\ 0 & \frac{1-\sigma}{h^2}& \frac{i}{\tau}-\frac{2(1-\sigma)}{h^2} &\frac{1-\sigma}{h^2} & \ldots\\ 0 & 0 & \frac{1-\sigma}{h^2}& \frac{i}{\tau}-\frac{2(1-\sigma)}{h^2} & \ddots\\ \vdots & \vdots & \vdots & \ddots & \ddots \end{pmatrix}$
$N=\begin{pmatrix} \frac{i}{\tau}+\frac{2\sigma}{h^2} & -\frac{\sigma}{h^2} & 0 & 0 & \ldots\\ -\frac{\sigma}{h^2}& \frac{i}{\tau}+\frac{2\sigma}{h^2} & -\frac{\sigma}{h^2} & 0 & \ldots\\ 0 & -\frac{\sigma}{h^2}& \frac{i}{\tau}+\frac{2\sigma}{h^2} &-\frac{\sigma}{h^2}& \ldots\\ 0 & 0 &-\frac{\sigma}{h^2}& \frac{i}{\tau}+\frac{2\sigma}{h^2} & \ddots\\ \vdots & \vdots & \vdots & \ddots & \ddots \end{pmatrix}$

И решаемое в два этапа матричное уравнение $A^{n+0.5} = N^{-1}MA^{n} \qquad A^{n+1} = N^{-1}MA^{n+0.5}$ (которое по первому слою для всех l, по второму слою для всех k).

А как $f^{n}_{k,l}$ включить в это матричное уравнение? Закинуть в матрицу М на главную диагональ? Тогда получается на каждой(!) итерации необходимо переопределять матрицу? Или отдельно складывать как столбец какой-нибудь?

[waxtep]

kitnone, делюсь очень старым и позабытым опытом, так что не взыщите, если советы покажутся тривиальными. Возможно кто-то из активно действующих специалистов найдет время дать больше. Из учебных материалов, мне в свое время показалось очень внятным описание прогонки в книге Р.П. Федоренко "Введение в вычислительную физику".

kitnone в сообщении #1551998 писал(а):

Корректно ли она составлена?
Позволительно ли $f^{n}_{k,l}A^{n}_{k,l}$ держать в первом слое? Или лучше раскинуть на оба или, только на второй? Исходил из того что раз $A^{n}$ используется на первой половинке то лучше туда отправить и все слагаемые, содержащие это.

Вы решаете уравнение Шредингера, или параболическое волновое уравнение (PWE; в этом случае $t$ имеет смысл третьей пространственной координаты, а $fA$ - наличие оптически плотной среды). Для ответа на вопросы хорошо ли делать так или этак, в целом есть два пути: 1. доказывать устойчивость схемы и оценивать зависимость величины невязки от различных параметров или 2. экспериментально подбирать параметры схемы на близких к Вашей задаче вариантах, имеющих точное решение (+ полезно сверяться с физическим смыслом задачи, что не полезли какие-то численные артефакты - тут обычно крайне полезна визуализация). Первый вариант теоретически нагружен (это непросто), во втором сможете получить правдоподобные результаты, но без теоретического обоснования. Из общих рукомахательных соображений "размазывание" $f$ между двумя слоями выглядит симпатичнее, а то Вы сейчас по сути на одном слое решаете неоднородную задачу, а на другом однородную - вряд ли это положительно влияет на устойчивость/невязку схемы (и Вам же это еще с граничными условиями сопрягать, всякие отражения от границ чувствительны к ммм чересчур простым вариантам схемы). Пробуйте, смотрите, что получается именно в Вашей задаче (т.е. я бы рекомендовал второй подход). Тоже из старинных воспоминаний: хорошую устойчивость и медленный рост невязки обеспечивают энергетические (интегральные) подходы: как определяется энергия в Вашей схеме, сколько-куда ее должно втекать-вытекать и т.п. - это может дать неочевидные с ходу малые поправки к коэффициентам схемы, повышающие устойчивость и уменьшающие невязку разностной схемы.

kitnone в сообщении #1551998 писал(а):

А как $f^{n}_{k,l}$ включить в это матричное уравнение? Закинуть в матрицу М на главную диагональ? Тогда получается на каждой(!) итерации необходимо переопределять матрицу? Или отдельно складывать как столбец какой-нибудь?

Да, на главную диагональ, а куда же деваться Если вы работаете с каким-нибудь матпакетом, там как правило масса эффективных алгоритмов работы с разреженными матрицами (sparse matrices). Если же реализуете прогонку руками, тоже наличие неоднородности $fA$ не должно вызывать проблем (кроме, может быть, случаев, когда $f$ имеет разрывы/резкие границы - здесь тоже требуется аккуратность при построении схемы, чтоб не нахватать численных артефактов)

[kitnone]

waxtep в сообщении #1552085 писал(а):

kitnone, делюсь очень старым и позабытым опытом

Огромное спасибо, хоть в целом успокоили. Раньше с Шрёдингером работал, правда с одномерным и однородным. С точностью до размерности численное решение двумерной однородной совпадает с численной одномерной (модуль лазерного пучка не изменяется, а вещественная и комплексная составляющие колеблются), а уже одномерная задача решена аналитически - на нестрогое доказательство сойдет)

Насчёт того что одна половинка однородная, а другая - нет -- резонно. Тогда попробую к обеим половинкам прибавить половинки от неоднородности.

Вы упомянули матпакет для матриц, а случаем не знаете есть ли что нибудь для подобных краевых задач? Вроде даже в Maple не встречал что то красивое чтоб аналитику/численное решение получить . А матлабовские приложения крайне не гибкие. Вот и пишу ручками.

А то отдельная проблема будет когда в качестве f возьму эффект Керра или покельса - типа $f(x,y,t)=|A(x,y,t)|$ или квадрат модуля. Там для решения сопряжённой краевой задачи придётся запоминать решение прямой задачи и с этим всем идти в обратный ход. Ухх.

[waxtep]

kitnone в сообщении #1552273 писал(а):

Тогда попробую к обеим половинкам прибавить половинки от неоднородности.

Давайте, удачи Раз у Вас нелинейная задача, могут быть вообще самые неожиданные сюрпризы, тут только численный эксперимент рулит. У меня в основном были более банальные линейные прямые и немного обратные (оптика, радиоволны), с которыми боролся как Вы пишете ниже.

kitnone в сообщении #1552273 писал(а):

А матлабовские приложения крайне не гибкие. Вот и пишу ручками.

Очень активно пользовался матлабом, шикарная среда для возни с УЧП на мой вкус. Его негибкость можно в значительной степени компенсировать дописыванием/переписыванием его встроенных функций, тот же PDE Toolbox в основном на матлабе и написан, удобно. Для быстрой проверки и экспериментов самое оно; но основные расчетные программы тоже ручками, на, господи прости, дельфях