Решение смешанной задачи для уравнения теплопроводности методом конечных разностей


Решение смешанной задачи для уравнения теплопроводности методом конечных разностей

1. Цель работы

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

2. Описание метода

Рассмотрим стержень из теплопроводящего материала с коэффициентом теплопроводности k. Предположим, что температура на концах стержня задана, а боковая поверхность стержня теплоизолирована. Пусть ось x направлена вдоль оси стержня, а его концы расположены в точках x=0 и x=L. Тогда задача сводится к определению зависимости от времени температуры u в точках стержня, то есть функции двух переменных u(x, t). Функция u(x, t) должна удовлетворять уравнению теплопроводности

(0<x<L), (1)

Начальному условию

U(x,0)=f(x), (0<x<L), (2)

И условиям на концах стержня

U(0,t)=j1(t), u(L, t)=j2(t), (tV0). (3)

Значения u(0,0) и u(L,0), полученные из (2) и (3), должны совпадать. Это будет если j1(0)=f(0), j2(0)=f(L).

Следует отметить, что путем замены переменных t^=a2T уравнение (1) можно преобразовать к виду

. (4)

Это означает, что решение задачи (1)-(3) путем замены переменных сводится к решению задачи (4),(2),(3). Далее будем полагать а=1.

Построим на плоскости (x, t) сетку с шагом h по переменной x (xI = (i-1)h, i=1,..,n+1, h=L/n) и с шагом t по переменной t (tJ = (j-1)t). Обозначим uIj = u(xI,tJ). Производные в уравнении (1) аппроксимируем следующим образом:

, (5)

. (6)

Подставляя (5) и (6) в (1) при a=1, получим разностное уравнение:

(7)

В соответствии с (2) и (3) значения

UI0 = f(xI), u0j = j1(tJ), uNj = j2(tJ) (8)

Являются известными. Тогда, подставляя в (7) j=0, получим систему n-1 линейных уравнений, решив которую можно определить uI1, i=1,..,n-1. При этом, поскольку u01=j1(t1), uN1=j2(t1), известными оказываются все значения временного слоя j=1, (t=t1). Затем, подставляя в (7) j=2, решаем систему уравнений относительно uI2 и т. д. для всех j=2,..,m.

Из (7) следует, что в каждое i-тое уравнение (i=1,..,n-1) с ненулевыми коэффициентами входят только три неизвестных uI-1,j; uIj; uI+1,j. Величина uI, j-1 к этому моменту является известной и потому отнесена в правую часть уравнения. Таким образом, матрица системы уравнений является трехдиагональной и эту систему можно решить методом прогонки. Для этого представим ее в стандартном виде:

.(9)

Для данной задачи xI=uIj, aI = l, gI = l, bI = 1-2l, b0 = 1, g0 = 0, j0 = u0j = j1(tJ), jN = uNj = j2(tJ), jI = - uI, j-1 (i=1,..,n-1).

Пусть на j-том шаге заданными являются параметры uI, j-1 (i=1,..,n-1), u0j, uNj, l. Все неизвестные значения uIj можно разместить в массиве xI (xI=uIj, i=0,..,n). Ищем связь xI-1 с xI в виде рекуррентного соотношения

XI-1=cI-1XI+nI-1, i=1,..,n. (10)

Подставляя (10) в (7), получаем

Lci-1xi-(1+2l)xi+lxi+1 = - ui, j-1-lni-1.

Отсюда

(11)

Сравнивая (11) с (10), находим рекуррентные соотношения

,

, (12)

C0= 0, n0 = u0j .

Таким образом, алгоритм определения значений uIj по известным uI, j-1 состоит из двух этапов: прямого хода прогонки по формулам (12) при i=1,..,n-1 и обратного хода прогонки по формуле (10) при i=n,..,2.

А)

Б)

Рис. 5.1

Теплопроводность смешанный конечный разность

Необходимо отметить, что разностное уравнение (7) связывает одно известное значение UI, j-1 (из предыдущего j-1 временного слоя) и три неизвестных UI, j, UI-1,j, UI+1,j. Поэтому найти значения UI, j (i=1,...,n-1) можно только все сразу путем решения системы уравнений. Такая схема связи переменных в разностном уравнении называется неявной. Шаблон неявной разностной схемы представлен на рис. 5.1а.

Наряду с неявной возможна организация явной разностной схемы. Для этого вместо выражения (5) для первой разностной производной по времени используют формулу

, (13)

Тогда разностное уравнение запишется в виде

(14)

В этом случае связываются три неизвестные значения, относящиеся к предыдущему временному слою (здесь j-тому) и только одно неизвестное UI, j+1. Шаблон явной разностной схемы представлен на рис. 5.1а.

При использовании этой схемы неизвестные параметры определяются путем последовательного применения формулы (2.14) при i=1,...n-1. Поскольку при этом не надо решать системы уравнений, то процесс определения параметров одного временного слоя требует меньших затрат времени, чем в случае неявной схемы.

Однако, неявная схема устойчива (ошибка не возрастает от шага к шагу) при любых значениях l=t/h2. Явная схема является устойчивой только при l<1/2. В противном случае развивается экспоненциальный рост погрешности так, что обычно происходит аварийная остановка ЭВМ по переполнению порядка. Поэтому при использовании явной схемы вычисления приходится вести с очень малым шагом по времени.

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

При выполнении данной работы будем предполагать, что температура на концах стержня поддерживается постоянной, то есть

J1(t)Tf(0), j2(t)Tf(L).

4. Варианты заданий

Решить смешанную задачу для уравнения теплопроводности с начальным u(x,0)=f(x) и граничными условиями u(0,t)=f(0), u(1,t)=b (L=1).

А

B

C

F(x)

1

1.1

3.0

0.05

f

b

A

0 X

c 1

Листинг программы:

#include <iostream. h>

#include <stdio. h>

#include <math. h>

#include <stdlib. h>

#include <windows. h>

#include <iomanip. h>

Float funkf(float x);

Float L, tau, lamda, h,T, a,b, c,betta[9999],gamma[9999],alfa[9999],X[9999],V[9999];

Int i, j,n, m;

Float x[9999],t[9999],u[9999][9999];

Void main()

{

L=1;

Cout<<"Vvedite koli4estvo zna4enii po x n"<<endl;

Cin>>n;

Cout<<"Vvedite maximal'nu'u temperaturu t"<<endl;

Cin>>T;

Cout<<"Vvedite shag po temperature tau"<<endl;

Cin>>tau;

Cout<<"Vvedite parametr a, b,c"<<endl;

Cin>>a>>b>>c;

M=T/tau+1;

Cout<<" "<<n<<endl;

H=L/n;

Lamda=tau/(h*h);

For(i=1;i<n+1;i++)

{

X[i]=(i-1)*h;

}

For(j=1;j<m;j++)

{

T[j]=(j-1)*tau;

}

For(j=0;j<m+1;j++)

{

U[0][j]=funkf(x[0]);

U[n][j]=b;

//cout<<"u[0][j]="<<u[0][j]<<" u[n][j]"<<u[n][j]<<endl;

}

/*for(i=1;i<n-1;i++)

{

For(j=1;j<m-1;j++)

{

U[i][j+1]=lamda*u[i-1][j]-(1-2*lamda)*u[i][j]+lamda*u[i+1][j];

}

}*/

Betta[0]=1;

Gamma[0]=0;

For(j=1;j<m-1;j++)

{

For(i=1;i<n;i++)

{

Alfa[i]=lamda;

Gamma[i]=lamda;

Betta[i]=1-2*lamda;

X[0]=0;

V[0]=u[0][j];

X[i]=lamda/(1+2*lamda-lamda*X[i-1]);

V[i]=(u[i][j-1]+lamda*V[i-1])/(1+2*lamda+lamda*X[i-1]);

U[i][j]=x[i];

}

For(i=n;i>1;i--)

{

U[i-1][j]=X[i-1]*u[i][j]+V[i-1];

}

}

For(i=0;i<n+1;i++)

{

For(j=0;j<m;j++)

{

Cout<<"u["<<i<<"]["<<j<<"]="<<u[i][j]<<endl;

}

}

/*for(j=0;j<m;j++)

{

Cout<<"u["<<n<<"]["<<j<<"]="<<u[n][j]<<endl;

}*/

System("PAUSE");

}

Float funkf(float x)

{

Float m;

If(x>=0)

{

If(x<c)

{

M=a;

}

Else

If(x<=1)

{

M=-2*x+3.1;

}

}

Else

Cout<<"x ne vxodit v zadannii promegutok"<<endl;

Return(m);

}

Результаты работы программы:

Vvedite koli4estvo zna4enii po x n

10

Vvedite maximal'nu'u temperaturu t

20

Vvedite shag po temperature tau

4

Vvedite parametr a, b,c

    1.1 3.0 0.05 10

U[0][0]=1.1

U[0][1]=1.1

U[0][2]=1.1

U[0][3]=1.1

U[0][4]=1.1

U[0][5]=1.1

U[1][0]=0

U[1][1]=0.984706

U[1][2]=0.989311

U[1][3]=0.989335

U[1][4]=0.989335

U[1][5]=0

U[2][0]=0

U[2][1]=0.871874

U[2][2]=0.878632

U[2][3]=0.87867

U[2][4]=0.878671

U[2][5]=0

U[3][0]=0

U[3][1]=0.980508

U[3][2]=0.988619

U[3][3]=0.988662

U[3][4]=0.988662

U[3][5]=0

U[4][0]=0

U[4][1]=1.2011

U[4][2]=1.21003

U[4][3]=1.21008

U[4][4]=1.21008

U[4][5]=0

U[5][0]=0

U[5][1]=1.46945

U[5][2]=1.47864

U[5][3]=1.47868

U[5][4]=1.47868

U[5][5]=0

U[6][0]=0

U[6][1]=1.75851

U[6][2]=1.76731

U[6][3]=1.76734

U[6][4]=1.76734

U[6][5]=0

U[7][0]=0

U[7][1]=2.05823

U[7][2]=2.06595

U[7][3]=2.06597

U[7][4]=2.06597

U[7][5]=0

U[8][0]=0

U[8][1]=2.36535

U[8][2]=2.37126

U[8][3]=2.37128

U[8][4]=2.37128

U[8][5]=0

U[9][0]=0

U[9][1]=2.67919

U[9][2]=2.68253

U[9][3]=2.68254

U[9][4]=2.68254

U[9][5]=0

U[10][0]=3

U[10][1]=3

U[10][2]=3

U[10][3]=3

U[10][4]=3

U[10][5]=3

Похожие статьи




Решение смешанной задачи для уравнения теплопроводности методом конечных разностей

Предыдущая | Следующая