Как указать начальные условия в частях домена в Mathematica? - PullRequest
Новая
       42

Как указать начальные условия в частях домена в Mathematica?

1 голос
/ 31 марта 2011

Я пытаюсь сравнить решение моего метода конечных объемов с решением Mathematica для простого уравнения теплопроводности U_t = U_xx. Для этого я должен указать начальное и граничное условие для функции NDSolve в Mathematica. Я хотел бы иметь на границах U = 90. В качестве начального условия я хотел бы иметь 100 во всех квадратных областях, кроме границ. Как я могу это сделать? Очевидно, мой код не работает. 

FUNC = NDSolve[{D[T[x, y, t], t] == (D[T[x, y, t], x, x] + D[T[x, y, t], y, y]), 
    T[x, y, 0] == 100, T[0, y, t] == 90, T[9, y, t] == 90, 
    T[x, 0, t] == 90, T[x, 9, t] == 90}, 
   T, {x, 0, 9}, {y, 0, 9}, {t, 0, 6}];

Попытка установить начальное условие для всех доменов, равное 100.

1 Ответ

6 голосов
/ 31 марта 2011

На самом деле, с небольшой модификацией ваш код работает. 

FUNC = T /. 
  NDSolve[{
     D[T[x, y, t], t] == (D[T[x, y, t], x, x] + D[T[x, y, t], y, y]), 
     T[x, y, 0] == 100, T[0, y, t] == 90, T[9, y, t] == 90, 
     T[x, 0, t] == 90, T[x, 9, t] == 90}, 
    T, {x, 0, 9}, {y, 0, 9}, {t, 0, 10}][[1]]

Модификации - это часть T/. (ReplaceAll) впереди и [[1]] (Партия) в конце;Вы можете посмотреть эти операции в документации.Они необходимы, чтобы преобразовать вывод в правильную форму.

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

В качестве альтернативы вы можете изменить начальное условие так, чтобы оно читалось как 

FUNC = T /. 
  NDSolve[{
     D[T[x, y, t], t] == (D[T[x, y, t], x, x] + D[T[x, y, t], y, y]), 
     T[x, y, 0] == If[x == 0 || x == 9 || y == 0 || y == 9, 90, 100], 
     T[0, y, t] == 90, T[9, y, t] == 90, T[x, 0, t] == 90, 
     T[x, 9, t] == 90}, T, {x, 0, 9}, {y, 0, 9}, {t, 0, 10}][[1]]

Теперь вы можете делать множество вещей с помощью этой функции, например, например, Manipulate:

Manipulate[ContourPlot[FUNC[x, y, t], {x, 0, 9}, {y, 0, 9}], {t, 0, 10}]

Manipulate ContourPlot

или анимированный GIF: 

anim = Table[DensityPlot[FUNC[x, y, t], {x, 0, 9}, {y, 0, 9}, 
          ColorFunctionScaling -> False, PlotPoints -> 50,
          ColorFunction -> (Hue[Rescale[#, {50, 100}], 1, 1] &)], 
        {t, 0, 10, .2}];
Export[ToFileName[$UserDocumentsDirectory, "anim.gif"], anim, "GIF"]

animated gif

...