2014 dxdy logo

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

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




Начать новую тему Ответить на тему
 
 Эллиптическая задача, сетка,контурный график, Mathematica
Сообщение26.06.2012, 12:07 


27/10/11
228
Здравствуйте.
Товарищи, помогите пожалуйста разобраться , что показывает выводимый 3х мерный график, после решения эллиптической задачи?
$$$-\Delta u = 0, \,\, 0<x<1, \,\, 0<y<1 $$$
$$$u(x,0)=0,   \,\,  u(x,1)=x     \,\, 0\leq x \leq 1$$$
$$$u(0,y)=0,  \,\,  u(1,y)=y,    \, \,0\leq u\leq 1$$$
Код решения задачи:

Код:
Dirichlet[n_, m_] :=
  Module[{Err = 1.0, i, j, k = 0, \[Omega]},
   \[Omega] = 4/(
    2 + Sqrt[4 - (Cos[\[Pi]/(n - 1)] + Cos[\[Pi]/(m - 1)])^2]); 
   While[ And[Err > 0.001 , k <= 70 ],
    Err = 0.0;
    For[ j = 2, j <= m - 1, j++,
     For[ i = 2, i <= n - 1, i++,
       Relax =
        1/4 (\[Omega] (u[[i, j + 1]] + u[[i, j - 1]] + u[[i + 1, j]] +
              u[[i - 1, j]] - 4 u[[i, j]]));
       u[[i, j]] = u[[i, j]] + Relax;
       Err = Max[Err, Abs[Relax]];  ]; ];
    Print["Max grid change = ", Err];
    k = k + 1;  ]; ];


n = 9;
m = 9;
Print["n = ", n]
Print["m = ", m]

u = Table[70.0, {n}, {m}];
For[i = 2, i < n, i++, u[[i, 1]] = i];
For[i = 2, i < n, i++, u[[i, m]] = 0];
For[j = 2, j < m, j++, u[[1, j]] = 0];
For[j = 2, j < m, j++, u[[n, j]] = j];

u[[1, 1]] = 1/2 (u[[1, 2]] + u[[2, 1]]);
u[[1, m]] = 1/2 (u[[1, m - 1]] + u[[2, m]]);
u[[n, 1]] = 1/2 (u[[n - 1, 1]] + u[[n, 2]]);
u[[n, m]] = 1/2 (u[[n - 1, m]] + u[[n, m - 1]]);


Dirichlet[n, m];


ListPlot3D[Transpose[u], AxesLabel -> {"y(j)", "x(i)", "u "},
ViewPoint -> {4, 2, 3}, Lighting -> False, ColorFunction -> Hue]
Print["Solução numerica para P.D.E:"]
Print["-delta u = 0"]


В Книжке Burden Analisys Numerical
http://ifolder.ru/31217077 на стр 697,701,702 есть пример решения такой задачи, но я не понял если честно также, в таблице на стр 702 (Talbe 12.2) , какой столбец соответсвутют столбцу "max greed change" из программы
http://math.fullerton.edu/mathews/n2003/ellipticpde/EllipticPDEMod/Links/EllipticPDEMod_lnk_3.html

Так же не очень понятно, зачем здесь мы задаём число 70 ? Разве гарница интервала не ограничена 1 ?
Код:
[b]u = Table[70.0, {n}, {m}];[/b]
For[i = 2, i < n, i++, u[[i, 1]] = i];
For[i = 2, i < n, i++, u[[i, m]] = 0];
For[j = 2, j < m, j++, u[[1, j]] = 0];
For[j = 2, j < m, j++, u[[n, j]] = j];

 Профиль  
                  
 
 Re: Эллиптическая задача, сетка,контурный график, Mathematica
Сообщение26.06.2012, 12:41 
Заслуженный участник


25/02/11
1797
Это приближенное решение задачи с нач. условиями $u(x,0)=x$, $u(x,1)=0$, $u(0,y)=0$, $u(1,y)=y$ :-)
Для $u(x,y)=xy$ должно быть
Код:
u = Table[70.0, {n}, {m}];
For[i = 2, i < n, i++, u[[i, 1]] = 0];
For[i = 2, i < n, i++, u[[i, m]] = i];
For[j = 2, j < m, j++, u[[1, j]] = 0];
For[j = 2, j < m, j++, u[[n, j]] = j];

А 70 здесь начальное приближение. Это типа мы решаем уравнение теплопроводности. Нач. т-ра 70, а в процессе стремится к стационарному распределению, определяемому граничными условиями.

 Профиль  
                  
 
 Re: Эллиптическая задача, сетка,контурный график, Mathematica
Сообщение26.06.2012, 12:47 


27/10/11
228
Спасибо :-)
Только Вы не могли бы пояснить пожалуйста, а что тут означает u=xy, точнее откуда мы это получили :-)

И как отразить в коде (и надо ли это делать) знак минус?

$$$-\Delta u = 0, \,\, 0<x<1, \,\, 0<y<1 $$$

 Профиль  
                  
 
 Re: Эллиптическая задача, сетка,контурный график, Mathematica
Сообщение26.06.2012, 21:06 
Заслуженный участник


25/02/11
1797
Поскольку уравнение однородное, то знак без разницы. А решение я просто угадал. Раз на границе полиномы, то...

 Профиль  
                  
 
 Re: Эллиптическая задача, сетка,контурный график, Mathematica
Сообщение27.06.2012, 00:04 


27/10/11
228
Ааа., понятно :-) Спасибо большое :-)

Вы не могли бы пожалуйста ещё пояснить, а что значат здесь значения n,m ? И почему при различных их значениях мы получаем разные графики трёхмерные ?

n=100,m=100
n=21,m=21
n=9,m=9
дают нам совершенно разные графики...

Код:
n = 100;
m = 100;
Print["n = ", n]
Print["m = ", m]

u = Table[70.0, {n}, {m}];
For[i = 2, i < n, i++, u[[i, 1]] = 0];
For[i = 2, i < n, i++, u[[i, m]] = 1];
For[j = 2, j < m, j++, u[[1, j]] = 1];
For[j = 2, j < m, j++, u[[n, j]] = 0];

u[[1, 1]] = 1/2 (u[[1, 2]] + u[[2, 1]]);
u[[1, m]] = 1/2 (u[[1, m - 1]] + u[[2, m]]);
u[[n, 1]] = 1/2 (u[[n - 1, 1]] + u[[n, 2]]);
u[[n, m]] = 1/2 (u[[n - 1, m]] + u[[n, m - 1]]);


Dirichlet[n, m];



п.с. т.е. можно сказать, что например, на графике изображено искривление металлического листа стали,задаваемое эллиптическим уравнением, верно ?

-- 27.06.2012, 01:12 --

Vince Diesel в сообщении #589238 писал(а):
Это приближенное решение задачи с нач. условиями $u(x,0)=x$, $u(x,1)=0$, $u(0,y)=0$, $u(1,y)=y$ :-)
Для $u(x,y)=xy$ должно быть
Код:
u = Table[70.0, {n}, {m}];
For[i = 2, i < n, i++, u[[i, 1]] = 0];
For[i = 2, i < n, i++, u[[i, m]] = i];
For[j = 2, j < m, j++, u[[1, j]] = 0];
For[j = 2, j < m, j++, u[[n, j]] = j];

А 70 здесь начальное приближение. Это типа мы решаем уравнение теплопроводности. Нач. т-ра 70, а в процессе стремится к стационарному распределению, определяемому граничными условиями.



Но в этом случае, получается что просто нет ограничений для x,y ?
$0\leq x \leq 1$
$0\leq y \leq 1$

-- 27.06.2012, 01:17 --

И можно у Вас ещё пожалуйста уточнить, а что изображает (охначает контурный график в этой задаче и как мы для него записываем эти параметры входные ? 0,180,15 ?
http://math.fullerton.edu/mathews/n2003/ellipticpde/EllipticPDEMod/Links/EllipticPDEMod_lnk_3.html

И как написать его же для моего номера ?

Код:
c = Table[x, {x, 0, 180, 15}];
ListContourPlot[Reverse[u], Contours -> c, ColorFunction -> Hue];
Print["Solução numerical para P.D.E."];
Print["- delta u = 0 "];

 Профиль  
                  
 
 Re: Эллиптическая задача, сетка,контурный график, Mathematica
Сообщение29.06.2012, 10:03 


27/10/11
228
Товарищи, ну так что Вы думаете?

 Профиль  
                  
Показать сообщения за:  Поле сортировки  
Начать новую тему Ответить на тему  [ Сообщений: 6 ] 

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



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

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


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

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