r2 = х2 + у2. (3.7)

Тогда дифференциальное уравнение теплопроводности для бесконечного цилиндра можно преобразовать так:

(3.8)

для бесконечного цилиндра можно преобразовать так:

(3.9)

(3.10)

Дифференцируя (3.8) по х, а (3.10) по у, получаем

(3.11)

(3.12)

Складывая уравнения (3.11) и (3.12) и принимая во внимание (3.7), получим для уравнения теплопроводности следующее выражение:

В общем случае, когда температура зависит от всех трех координат (х, у, г), дифференциальное уравнение теплопроводности конечного ци­линдра имеет вид

; (3.13)

4 СОСТАВЛЕНИЕ АЛГОРИТМА

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

; (4.1)

Заменим частный дифференциал разностным отношением:

; (4.2)

Осуществим следующее преобразование функции:

; (4.3)

; (4.4)

; (4.5) (4.6)

; (4.7)

; (4.8)

Подготовим уравнение (4.8) для рекуррентного вычисления в MatLab V6.0

Произведём переобозначения:

; (4.9)

; (4.10)

; (4.11)

; (4.12)

; (4.13)

Имеем формулу:

T(i+1,j+1)=T(i,j+1)+(a*dt/dr)*(((T(i,j+2)-2*T(i,j+1)+T(i,j))/dr)+((1/r)*(T(i,j+2)-T(i,j+1)))); (4.14)

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

1.Программа начинается c задание переменных: начального и конечного момента времени, числа дискретных отсчётов по времени, радиус цилиндра и число его разбиений, констант характеризующих тепло-физические свойства полимера.

2.Следующим этапом является вычисление шага аргументов, по которым будет вычисляться исходная функция.

3.Краевые условия: значения искомой функции в начальный момент времени t0 = 0 в зависимости от радиуса, и температуры стенки литникового канала в любой момент времени задаются циклом For.

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

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

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

7.При переходе к внешнему циклу отсчёт по времени увеличивается на единицу. Значение производной температуры по радиусу в любой момент времени равно нулю и поэтому, чтобы учесть ещё одно краевое условие при переходе от внешнего цикла к внутреннему значение последней температуры копируется два раза.

8.После получения матрицы температур надо построить график. Чтобы координатные оси были проградуированные удобно для использования в матрице температур переставляют столбцы. Осуществляется это с использованием двух вложенных циклов.

9.Далее следует вывод графика и градуировка его осей.

5 СОСТАВЛЕНИЕ ПРОГРАММЫ

Программа для MatLab v6.0 R12 начинается очищения переменных графических окон функций и окна вывода результата. Осуществляют это с помощью: clear, clc, clf, clg

Чтобы программа была легка в использовании и проста в конфигурировании под любые задачи разработаем её используя понятные обозначения:

Задаём переменные:

начальный момент времени выбираем как t0=0;

конечный момент времени tk=120;

число дискретных отсчётов времени nt=120;

температура стенки Tc=30;

максимальная температура материала в середине цилиндра Tpol=170;

число дискретных отсчетов длинны цилиндра nR=10;

радиус цилиндра R=0.01 м;

температуропроводность полистирола a = 0.00000056 град/м с

Рассчитаем интервалы изменения температуры и радиуса

dr=R/(nR-1);

dt=(tk-t0)/(nt-1);

Присвоим начальные значения температуры стенки в цикле For:

for i=1:nt+1

T(i,1)=Tc;

end

Присвоим начальные значения температурного поля полимера в цикле:

for j=1:nR+2

T(1,j)=Tpol*exp(-2000*(R-dr*(j-1))^2);

end

Рассчитаем матрицу температурного поля T во вложенном цикле For:

for i=1:nt

for j=1:nR

r=R-dr*(j-1)+0.0001*dr;

T(i+1,j+1)=T(i,j+1)+(a*dt/dr)*(((T(i,j+2)-2*T(i,j+1)+T(i,j))/dr)+((1/r)*(T(i,j+2)-T(i,j+1))));