2014 dxdy logo

Научный форум dxdy

Математика, Физика, Computer Science, Machine Learning, LaTeX, Механика и Техника, Химия,
Биология и Медицина, Экономика и Финансовая Математика, Гуманитарные науки




Начать новую тему Ответить на тему На страницу 1, 2  След.
 
 Начально-краевая задача для уравнения переноса.
Сообщение26.05.2018, 13:40 


25/05/18
8
Нужно решить численно начально-краевую задачу для уравнения переноса $u_t+u_x=0$
Начальное условие: $u(0,x) = 1, x\leqslant5 $ и $u(0,x) = 0, 5\leqslant x \leqslant20$ и граничное условие $u(t,0) = 1$.
Нужно использовать начальную сетку из 41 точки при $h = 0,5$ и провести расчёт до значения $T = 9$. Решить для чисел Куранта $\tau/h = 1.0, 0.6, 0.3$. Рекомендованная разностная схема Лакса-Вендроффа.

Моё решение:
Схема Лакса-Вендроффа: $\frac{u_j^{i+1} - u_j^i}{\tau} + a \frac{u_{j+1}^i - u_{j-1}^i}{2h} - \frac{a^2\tau}{2}\frac{u_{j+1}^i - 2u_j^i + u_{j-1}^i}{h^2} = 0 $
Выражаю $u_j^{i+1}$:
$u_j^{i+1} = u_j^i + \frac{a^2\tau^2}{2h^2}(u_{j+1}^i - 2u_j^i + u_{j-1}^i) - \frac{a\tau}{2h}(u_{j+1}^i - u_{j-1}^i)$
В нашем случае (для начала число куранта будет 1.0) $a = 1$, $\tau = 0.5$:
$u_j^{i+1} = u_j^i + \frac{1}{2}(u_{j+1}^i - 2u_j^i + u_{j-1}^i) - \frac{1}{2}(u_{j+1}^i - u_{j-1}^i)$
Ну и тут получается бред, что $u_j^{i+1} = u_{j-1}^i$ (*)

В программе для начала, чтобы не путаться, вбил статический массив (с учётом шага по вертикали и горизонтали по 0.5):
Код:
double T = 9;
   double L = 20;
   double u[19][41];

Дальше задавал граничные и начальные условия:
Код:
   for (int i = 0; i <= 18; i++){
      u[i][0] = 1;
   }
   for (int j = 0; j <= 10; j++){
      u[0][j] = 1;
   }
   for (int j = 11; j <= 40; j++){
       u[0][j] = 0;
   }

Ну и для следующего момента времени функцию u выражал по вышеприведенной формуле:
Код:
for (int i = 0; i <= 17; i++){
      for (int j = 1; j <= 39; j++) {
         u[i + 1][j] = u[i][j] + 1 / 2 * (u[i][j + 1] - 2 * u[i][j] + u[i][j - 1]) - 1 / 2 * (u[i][j + 1] - u[i][j - 1]);
      }
   }

И вот, собственно, результат. Очевидно, это не то.
При числах куранта 0.6 и 0.3 получается тоже самое, хотя формула (*) уже не выглядит такой тривиальной.
Изображение

Я понимаю, что ошибся, и что ошибка моя в самом начале и связана, как мне кажется, с недопониманием. Прошу навести на верный путь .

 Профиль  
                  
 
 Re: Начально-краевая задача для уравнения переноса.
Сообщение26.05.2018, 19:37 
Заслуженный участник
Аватара пользователя


01/08/06
3134
Уфа
a3ro писал(а):
Ну и тут получается бред, что $u_j^{i+1} = u_{j-1}^i$ (*)
Почему бред? Так и должно было получиться.
a3ro писал(а):
Используется синтаксис C
for (int i = 0; i <= 17; i++){
  for (int j = 1; j <= 39; j++) {
    u[i + 1][j] = u[ i ][j] + 1 / 2 * (u[ i ][j + 1] - 2 * u[ i ][j] + u[ i ][j - 1]) - 1 / 2 * (u[ i ][j + 1] - u[ i ][j - 1]);
  }
}
Не стесняйтесь, пишите, как получилось после упрощений:
Используется синтаксис C
    u[i + 1][j] = u[ i ][j - 1];
Здесь всё верно. Вот результат уже немного смущает. Единички должны на каждом шаге на одну вверх смещаться. Точнее, вдоль прямых $x-t=\rm{const}$ решение должно быть постоянным: где на пересечении этой прямой с границей были нулики, там и будет 0, а где были единички, будет 1. Это верно для точного решения задачи. Но в данном случае (когда число Куранта = 1) это верно и для сеточного решения, ибо здесь сеточное совпадает с точным (в тех точках, в которых вычисляется).
Что-нибудь понятно из моего объяснения? Ключевой вопрос, который нужно здесь понять — вопрос характеристик уравнения переноса. Знаете, что это такое?

 Профиль  
                  
 
 Re: Начально-краевая задача для уравнения переноса.
Сообщение26.05.2018, 20:57 


25/05/18
8
Цитата:
Единички должны на каждом шаге на одну вверх смещаться. Точнее, вдоль прямых $x-t=\rm{const}$ решение должно быть постоянным

Если учитывать, что у меня ось $\tau$ направлена "вниз", то наверное в каждой следующей строке должно быть на одну единичку больше. Т.е. от "единиц" идут "единицы" вправо вниз, от нулей должны идти нули. Поправьте, если сказал не так.
Значит ошибка всё-таки была не в начале, а в коде, сейчас попробую пофиксить.

Цитата:
Ключевой вопрос, который нужно здесь понять — вопрос характеристик уравнения переноса. Знаете, что это такое?

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

 Профиль  
                  
 
 Re: Начально-краевая задача для уравнения переноса.
Сообщение26.05.2018, 21:41 
Заслуженный участник
Аватара пользователя


01/08/06
3134
Уфа
Ну да, всё правильно. Кроме ответа, почему-то.
Возможно, компилятор воспринимает 1 / 2 как целочисленное деление и выдаёт в промежуточном результате 0. Это ж Си, он такое любит. Я бы во всех таких местах поставил 0.5. Или 1.0 / 2.0.

 Профиль  
                  
 
 Re: Начально-краевая задача для уравнения переноса.
Сообщение26.05.2018, 22:34 


25/05/18
8
Пофиксил. Я уже забыл о куче подводных камней в языке Си да и в программировании в целом.
Изображение

Получается для числа Куранта $q=1$ сгущение сеток для данных граничных и начальных условий не будет влиять на погрешность численного решения (т.к. погрешности нет). Но для чисел Куранта $q_2=0.6$ и $q_3=0.3$ такая погрешность появится.
И дальше немного глупый вопрос, который вытекает из условия задачи: "просят оценить порядок сходимости схемы сравнением изменения погрешности численного решения по сравнению с точным на последовательности сгущающихся сеток при фиксированном числе Куранта (применить сгущение сеток вдвое)".
Т.е. порядок такой:
1) Мы берем число Куранта. (А дальше у меня в голове путаница...)
2) Для каждого параметра h (0.5, 0.25 или 0.125) сетки мы берем разность $([u]_h - u)$, где u точное решение, и сравниваем его c $h^k$ $(ch^k ?)$ где k и есть наш порядок сходимости.
Т.е. $([u]_h - u) \leqslant  h^kc$.
У нас такое соотношение должно выполняться для каждой точки? Как находить число c?

 Профиль  
                  
 
 Re: Начально-краевая задача для уравнения переноса.
Сообщение26.05.2018, 22:49 
Заслуженный участник
Аватара пользователя


01/08/06
3134
Уфа
a3ro писал(а):
Как находить число c?
Тут главное оценить порядок (как вас и просят в условии), т.е. $k$. Например, находите максимум $([u]_h - u)$ по всем точкам, а потом то же самое для уменьшенного в 2 (3, 4, 10, ...) раз шага. Максимум погрешности должен при этом уменьшиться примерно в $2^k$ ($3^k$, $4^k$, $10^k$, ...) раз. После этого и $c$ можно оценить (хотя этого и не просят).

 Профиль  
                  
 
 Re: Начально-краевая задача для уравнения переноса.
Сообщение27.05.2018, 17:45 


25/05/18
8
Собственно решил для примера вывести решение для сетки $h = 0.5, q = 0.6$
ИзображениеИзображение
(Пришлось в 2 картинки из-за ограничения по ширине).
Где приведенная численная схема выглядит: $u_j^{i+1} = 0.64u_j^i + 0.48u_{j-1}^i - 0.12u_{j+1}^i$

Ну и вот, что я подумал. Тут у нас полный диапазон значений от 0 до 1.2. А точные решения в данной задаче либо 1, либо 0. Максимум $([u]_h-u)$ по всем точкам будет 0.5 (как видно из скрина). Я сделал тоже самое для вдвое сгущенной сетки ($h = 0.25$), там тоже встречаются значения $u_h = 0.5$, а значит и максимум $([u]_h-u)$ будет тоже 0.5.

Поправьте, пожалуйста, если я ошибаюсь.

 Профиль  
                  
 
 Re: Начально-краевая задача для уравнения переноса.
Сообщение28.05.2018, 00:21 
Заслуженный участник
Аватара пользователя


01/08/06
3134
Уфа
a3ro писал(а):
$u_j^{i+1} = 0.64u_j^i + 0.48u_{j-1}^i - 0.12u_{j+1}^i$
Ниже неправильные рассуждения:
Глядя только на это, без предыдущих выкладок, можно сразу сказать, что условие Куранта тут нарушено. Число Куранта > 1. В результате возник отрицательный коэффициент (-0.12) там, где (при выполнении условия Куранта) должны были быть только положительные.
Это приводит к нарушению устойчивости схемы и делает её бесполезной чуть менее, чем полностью. По результату это тоже видно, как только значения начинают превышать максимум краевых значений 1.0: при устойчивой схеме такое немыслимо Думаю, дальше там и отрицательные значения повылезают. Ошибки округления растут экспоненциально, и довольно быстро забивают всё решение. Нет смысла говорить о порядке точности в таких условиях.

 Профиль  
                  
 
 Re: Начально-краевая задача для уравнения переноса.
Сообщение28.05.2018, 01:25 


25/05/18
8
worm2 в сообщении #1315432 писал(а):
a3ro писал(а):
$u_j^{i+1} = 0.64u_j^i + 0.48u_{j-1}^i - 0.12u_{j+1}^i$
Глядя только на это, без предыдущих выкладок, можно сразу сказать, что условие Куранта тут нарушено. Число Куранта > 1. В результате возник отрицательный коэффициент (-0.12) там, где (при выполнении условия Куранта) должны были быть только положительные.
Это приводит к нарушению устойчивости схемы и делает её бесполезной чуть менее, чем полностью. По результату это тоже видно, как только значения начинают превышать максимум краевых значений 1.0: при устойчивой схеме такое немыслимо Думаю, дальше там и отрицательные значения повылезают. Ошибки округления растут экспоненциально, и довольно быстро забивают всё решение. Нет смысла говорить о порядке точности в таких условиях.


Так, давайте распишу тогда подробнее:
Цитата:
Вот наша формула:
$u_j^{i+1} = u_j^i + \frac{a^2\tau^2}{2h^2}(u_{j+1}^i - 2u_j^i + u_{j-1}^i) - \frac{a\tau}{2h}(u_{j+1}^i - u_{j-1}^i)$

$a = 1, q=\tau/h = 0.6$. Тогда $\tau^2/h^2 = 0.36. Подставляю в уравнение выше:
$u_j^{i+1} = u_j^i + 1/2 \cdot 0.36 \cdot (u_{j+1}^i - 2u_j^i + u_{j-1}^i) - 1/2 \cdot 0.6 \cdot(u_{j+1}^i - u_{j-1}^i)$
$u_j^{i+1} = u_j^i + 0.18 \cdot (u_{j+1}^i - 2u_j^i + u_{j-1}^i) - 0.3 \cdot (u_{j+1}^i - u_{j-1}^i)$
И вот тут у нас вылазит минус перед членом $u_{j+1}^i$. Для числа Куранта 0.3 по таким же рассуждениям тоже появится отрицательный член.

 Профиль  
                  
 
 Re: Начально-краевая задача для уравнения переноса.
Сообщение28.05.2018, 13:27 
Заслуженный участник
Аватара пользователя


01/08/06
3134
Уфа
Да, я наврал в предыдущем сообщении, поспешил :oops: Хотел спинным мозгом обойтись, а тут надобно головным думать :-)

Условие Куранта выполняется. И формулы у вас правильные. И отрицательные коэффициенты для схемы Лакса–Вендроффа неизбежны. Но устойчивость всё-таки есть, просто она более хитро проверяется, не на основе принципа максимума, а по фон Нейману.

Соответственно, и норму невязки придётся считать не как максимум $([u]_h - u)$ по всем точкам, а как-то усредняя погрешность по всем точкам.
Например, в сеточном аналоге нормы $L_2$ как-то так:
$$||u_h - u||=\sqrt{\frac{1}{(N_x+1) (N_t+1)}\sum\limits_{i=0}^{N_t}\sum\limits_{j=0}^{N_x}\left(u_j^i-u_{\text{точн}}(x_j, t_i)\right)^2}.$$
И вот эта погрешность должна убывать как $k$-я степень шага при уменьшении шагов с сохранением условия Куранта.

 Профиль  
                  
 
 Re: Начально-краевая задача для уравнения переноса.
Сообщение29.05.2018, 16:57 


25/05/18
8
Благодарю, вчера еще увидел, но уснул на 15 часов :oops:

Не будет ли в данном случае халтурой, если я буду просто каждый элемент численных решений сравнивать с 0 и 1 и брать наименьшую разность, чтобы подставить в фон Неймана?

Ну и вопрос (возможно не совсем по теме) по Лаксу-Вендроффу: как посроить эту схему для стационарного уравнения Эйлера:
$v\frac{dv}{dx}$ = -$\frac{1}{\rho}$$\frac{dp}{dx}$.
Моё непонимание в том, что тут 2 функции, а дифференцируются по 1 переменной, а не наоборот.
Или в данных случаях переворачивают равенство (хотя можно ли это делать, когда у нас скорость и давление функции нескольких переменных для 2D случая), накладывая условия для области определения?

 Профиль  
                  
 
 Re: Начально-краевая задача для уравнения переноса.
Сообщение29.05.2018, 18:11 
Заслуженный участник
Аватара пользователя


01/08/06
3134
Уфа
a3ro в сообщении #1315917 писал(а):
Не будет ли в данном случае халтурой, если я буду просто каждый элемент численных решений сравнивать с 0 и 1 и брать наименьшую разность, чтобы подставить в фон Неймана?
Не знаю. Боюсь, не получится.
a3ro в сообщении #1315917 писал(а):
Моё непонимание в том, что тут 2 функции, а дифференцируются по 1 переменной, а не наоборот.
Когда одна переменная, а функций две, это гораздо проще. Пара обычных дифуров, не в частных производных. Можно обойтись без разностных схем, а методами вроде Рунге-Кутта. Только проблема, что неизвестных больше, чем уравнений :-) Ну то есть одного уравнения Эйлера мало. Или вам давление дано, и вы по нему можете вычислить скорость, или наоборот. Или ещё что-то.
Когда переменных две или больше, тут уже будут уравнения в частных производных, и схема понадобится. Только та же проблема остаётся: не хватает одного этого уравнения, нужно ещё что-то, либо давление, либо скорость, либо ещё какое-то уравнение. Нужно поработать над постановкой, а потом уже задумываться над методами решения.
И потом, тут вы перескакиваете от модельной математической задачи к более сложной задаче гидродинамики, не слишком ли быстро?

 Профиль  
                  
 
 Re: Начально-краевая задача для уравнения переноса.
Сообщение29.05.2018, 19:18 


25/05/18
8
Цитата:
Не знаю. Боюсь, не получится.

Да просто не хочу "вручную" забивать точные данные для каждой точки в сетках, где число куранта не 0, чтобы сравнивать (хочу обойтись спинным мозгом :D)
Цитата:
Только проблема, что неизвестных больше, чем уравнений :-) Ну то есть одного уравнения Эйлера мало. Или вам давление дано, и вы по нему можете вычислить скорость, или наоборот. Или ещё что-то.Нужно поработать над постановкой, а потом уже задумываться над методами решения.

Я опечатался в предыдущем сообщении, там p - это импульс.
Задача заключается в том, что нужно в рамках краевой задачи для уравнения Эйлера решить задачу об обтекании кругового цилиндра поперек потока. Ничего больше не дали :D
Тут нужно еще уравнение $divV = 0$ - уравнение "неразрывности" для жидкости (плотность $\rho$ const).
Цитата:
И потом, тут вы перескакиваете от модельной математической задачи к более сложной задаче гидродинамики, не слишком ли быстро?

Ну гидродинамику я немного знаю, и подобную задачу решал через комплексный потенциал ,находил распределение скорости на поверхности (закон от угла $V = 2V_\infty \sin\theta$, где $\theta$ - угол "от центра"). Но с уравнением Эйлера дел на конкретных примерах не имел.

 Профиль  
                  
 
 Re: Начально-краевая задача для уравнения переноса.
Сообщение29.05.2018, 20:22 
Заслуженный участник
Аватара пользователя


01/08/06
3134
Уфа
Ну это на порядок более сложная задача. Боюсь, здесь у меня кончился бензин. Нужно читать учебники, тут ещё до разностных схем пилить и пилить. Постановка, построение сетки (видимо, не обойтись без нерегулярной)...

Кстати, вроде бы $p$ — всё-таки давление.

 Профиль  
                  
 
 Re: Начально-краевая задача для уравнения переноса.
Сообщение29.05.2018, 22:15 
Заслуженный участник
Аватара пользователя


01/08/06
3134
Уфа
По разностной схеме для задачи гидродинамики я бы порекомендовал создать новую тему.
От заголовка этой темы она отошла.

 Профиль  
                  
Показать сообщения за:  Поле сортировки  
Начать новую тему Ответить на тему  [ Сообщений: 16 ]  На страницу 1, 2  След.

Модераторы: Karan, Toucan, PAV, maxal, Супермодераторы



Кто сейчас на конференции

Сейчас этот форум просматривают: нет зарегистрированных пользователей


Вы не можете начинать темы
Вы не можете отвечать на сообщения
Вы не можете редактировать свои сообщения
Вы не можете удалять свои сообщения
Вы не можете добавлять вложения

Найти:
Powered by phpBB © 2000, 2002, 2005, 2007 phpBB Group