Колиденкова (1211099), страница 4
Текст из файла (страница 4)
; 
; 
;
;
;
;
;
;
В результате получим
.
Вычислив интегралы, получим локальную матрицу теплопроводности (Рис.2). Для наглядности пусть шаг по всем осям будет одинаковым.
Рисунок2 - Локальная матрица теплопроводности
Правую часть вычисляем по формуле:
Далее строим матрицу связности, для ансамблирования правой части и локальной матрицы теплопроводности в глобальные матрицы.
Где 
– число узлов в сетке. Здесь знак равенства означает совпадение глобального и локального узлов.
Таким образом, мы получаем систему линейных алгебраических уравнений
,
которую решаем итерационным методом Зейделя, код программы представлен в приложении А.
В результате работы программы, было получено наглядное представление распределения температуры. На рисунке 3 взят срез области по центру параллельно граням, 
, 
. На рисунке 4 - 
, .
Рисунок 3 - Распределение температуры при ,
Рисунок 4 - Распределение температуры при ,
Видно, что температура у источника самая высокая и постепенно уменьшается к границе области. При отсутствии источника тепла, т.е. 
, движения не наблюдается, результат представлен на рисунке 5. Это свидетельствует о том, что результат расчетов был верный.
Рисунок 4 - Распределение температуры при
Заключение
В результате работы были выполнены в полном объёме все поставленные задачи, была рассмотрена задача о распространения температуры в однородном теле с внутренними источниками тепла и терморегулятором. Был проведен анализ модели данной задачи с помощью вариационного метода.
Были изучены необходимые теоретические сведения из теории вариационных неравенств. Но несмотря на эффективность вариационного метода, зачастую им пренебрегают при рассмотрении определенных прикладных задач. Такие задачи просто решаются широко известными и глубоко изученными численными методами. При этом существование самого решения и его единственность остаются не рассмотренными. Такой подход является в корне не верным, поскольку ставит под сомнение сам результат численного решения задачи. Именно поэтому в данной работе приведены доказательства существования и единственности решения рассматриваемой задачи с помощью теории вариационных неравенств.
После чего была с помощью метода конечных элементов решена исследуемая задача и построен соответствующий вычислительный алгоритм. С помощью разработанного комплекса прикладных программ на объектно-ориентированном языке C++, осуществлена численная реализация этого алгоритма, в результате котрого было получено наглядное представление распределения температуры в однородном теле с внутренними источниками тепла в трехмерной области.
Значимость исследования сложно переоценить, так как в практической деятельности моделирование играет немаловажную роль и является одним из универсальных методов познания.
Список использованных источников
1 Араманович И.Г., Левин В.И. Уравнения математической физики.-М.: Наука, 1964. -288 с.
2 Арсенин В.Я. Методы математической физики и специальные функции. - М.: Наука, 1974. -421 с.
3 Бицадзе А.В. Уравнения математической физики. - М.: Наука, 1976. -296 с.
4 Джеффрис Г., Свирлс Б. Методы математической физики. - М.: Мир, 1972. -344 с.
5 Тихонов А.Н., Самарский А.А. Уравнения математической физики. -5-е изд. - М.: Наука, 1977. -736 с.
6 Галлагер Р. Метод конечных элементов. Основы. М.: Мир, 1984 г. – 428 с.
7 Бате Н., Вилсон Е. Численные методы анализа и метод конечных элементов. – М.: Стройиздат, 1982. – 448 с.
8 Бабенко К. И. Основы численного анализа. М.: Наука, 1986. 744 с.
Приложение А
Листинг программы
-
#include <iostream>
-
#include <cmath>
-
#include <clocale>
-
-
using namespace std;
-
-
double hx, hy, hz, h;
-
double lx = 1, ly = 1, lpx = 0.5, lpy = 0.5;
-
double k1, k2, sigma, g0, g1;
-
int Ni, Nj, Nk, ni, nj, nk;
-
double K[8][8], Fm[8];
-
double **A, *F, **Omega, **tmp1, *u, *v;
-
#define eps 1e-5
-
-
//создание матрицы связности
-
void makeOmega1(int i, int j, int k){
-
for (int _i = 0; _i < 8; _i++){
-
for (int _j = 0; _j < nj*nj*nk; _j++){
-
Omega[_i][_j] = 0;
-
}
-
}
-
Omega[0][i + ni*j + ni*nj*k] = 1;
-
Omega[1][i+1 + ni*j + ni*nj*k] = 1;
-
Omega[2][i+1 + ni*(j+1) + ni*nj*k] = 1;
-
Omega[3][i + ni*(j+1) + ni*nj*k] = 1;
-
Omega[4][i + ni*j + ni*nj*(k+1)] = 1;
-
Omega[5][i+1 + ni*j + ni*nj*(k+1)] = 1;
-
Omega[6][i+1 + ni*(j+1) + ni*nj*(k+1)] = 1;
-
Omega[7][i + ni*(j+1) + ni*nj*(k+1)] = 1;
-
}
-
-
//умножение матрицы связности на локальную матрицу теплопроводности
-
void mul1(){
-
for (int i = 0; i < nj*ni*nk; i++){
-
for (int j = 0; j < 8; j++){
-
double sum = 0;
-
for (int l = 0; l < 8; l++){
-
sum += Omega[l][i]*K[l][j];
-
}
-
tmp1[i][j] = sum;
-
}
-
}
-
}
-
-
void mul2(){
-
for (int i = 0; i < nj*ni*nk; i++){
-
for (int j = 0; j < nj*ni*nk; j++){
-
double sum = 0;
-
for (int l = 0; l < 8; l++){
-
sum += tmp1[i][l]*Omega[l][j];
-
}
-
A[i][j] += sum;
-
//printf("%.3f\t", sum);
-
} //cout << endl;
-
} //cout << endl;
-
}
-
-
//умножение матрицы связности на правую часть
-
void mulF(){
-
for (int i = 0; i < nj*ni*nk; i++){
-
for (int j = 0; j < 1; j++){
-
double sum = 0;
-
for (int l = 0; l < 8; l++){
-
sum += Omega[l][i]*Fm[l];
-
}
-
F[i] += sum;
-
}
-
}
-
}
-
-
// Условие окончания метода Зейделя
-
bool converge()
-
{
-
double norm = 0;
-
for (int i = 0; i < ni*nj*nk; i++)
-
{
-
norm += (v[i] - u[i])*(v[i] - u[i]);
-
}
-
cout << norm << endl;
-
if(sqrt(norm) >= eps)
-
return false;
-
return true;
-
}
-
-
int main(){
-
setlocale(LC_CTYPE, "rus");
-
cout << "Введите количество узлов по X: ";
-
cin >> Ni;
-
cout << " Введите количество узлов по Y: ";
-
cin >> Nj;
-
cout << " Введите количество узлов по Z: ";
-
cin >> Nk;
-
ni = Ni;
-
nj = Nj;
-
nk = Nk;
-
hy = 1.0/(Ni-1);
-
hx = 1.0/(Nj-1);
-
hz = 1.0/(Nk-1);
-
double h = min(hx,min(hy,hz));
-
g0 = 0;
-
A = new double*[ni*nj*nk];
-
F = new double[ni*nj*nk];
-
Omega = new double*[8];
-
tmp1 = new double*[ni*nj*nk];
-
for (int i = 0; i < ni*nj*nk; i++){
-
A[i] = new double[ni*nj*nk];
-
for (int j = 0; j < ni*nj*nk; j++){
-
A[i][j] = 0;
-
}
-
F[i] = 0;
-
tmp1[i] = new double[8];
-
}
-
for (int i = 0; i < 8; i++){
-
Omega[i] = new double[ni*nj*nk];
-
for (int j = 0; j < ni*nj*nk; j++){
-
Omega[i][j] = 0;
-
}
-
Fm[i] = h*h*h/8/0.05*0;
-
}
-
-
K[0][0] = h/3;
-
K[0][1] = 0;
-
K[0][2] = -h/12;
-
K[0][3] = 0;
-
K[0][4] = 0;
-
K[0][5] = -h/12;
-
K[0][6] = -h/12;
-
K[0][7] = -h/12;
-
K[1][0] = 0;
-
K[1][1] = h/3;
-
K[1][2] = 0;
-
K[1][3] = -h/12;
-
K[1][4] = -h/12;
-
K[1][5] = 0;
-
K[1][6] = -h/12;
-
K[1][7] = -h/12;
-
K[2][0] = -h/12;
-
K[2][1] = 0;
-
K[2][2] = h/3;
-
K[2][3] = 0;
-
K[2][4] = -h/12;
-
K[2][5] = -h/12;
-
K[2][6] = 0;
-
K[2][7] = -h/12;
-
K[3][0] = 0;
-
K[3][1] = -h/12;
-
K[3][2] = 0;
-
K[3][3] = h/3;
-
K[3][4] = -h/12;
-
K[3][5] = -h/12;
-
K[3][6] = -h/12;
-
K[3][7] = 0;
-
K[4][0] = 0;
-
K[4][1] = -h/12;
-
K[4][2] = -h/12;
-
K[4][3] = -h/12;
-
K[4][4] = h/3;
-
K[4][5] = 0;
-
K[4][6] = -h/12;
-
K[4][7] = 0;
-
K[5][0] = -h/12;
-
K[5][1] = 0;
-
K[5][2] = -h/12;
-
K[5][3] = -h/12;
-
K[5][4] = 0;
-
K[5][5] = h/3;
-
K[5][6] = 0;
-
K[5][7] = -h/12;
-
K[6][0] = -h/12;
-
K[6][1] = -h/12;
-
K[6][2] = 0;
-
K[6][3] = -h/12;
-
K[6][4] = -h/12;
-
K[6][5] = 0;
-
K[6][6] = h/3;
-
K[6][7] = 0;
-
K[7][0] = -h/12;
-
K[7][1] = -h/12;
-
K[7][2] = -h/12;
-
K[7][3] = 0;
-
K[7][4] = 0;
-
K[7][5] = -h/12;
-
K[7][6] = 0;
-
K[7][7] = h/3;
-
-
//return 0;
-
-
freopen("1.txt", "w", stdout);
-
for (int i = 0; i < 8; i++){
-
for (int j = 0; j < 8; j++){
-
printf("%.3f\t", K[i][j]);
-
} cout << endl;
-
}
-
cout << endl << endl;
-
for (int i = 0; i < 8; i++){
-
printf("%.3f\t", Fm[i]);
-
cout << endl;
-
}
-
cout << endl << endl;
-
-
//return 0;
-
-
//заполняем глобальную матрицу
-
for (int i = 0; i < Ni-1; i++){
-
for (int j = 0; j < Nj-1; j++){
-
for (int k = 0; k < Nk-1; k++){
-
makeOmega1(i, j, k);
-
mul1();
-
mul2();
-
mulF();
-
}
-
}
-
}
-
-
//return 0;
-
-
for (int ii = 0; ii < nj*ni*nk; ii++){
-
for (int jj = 0; jj < nj*ni*nk; jj++){
-
printf("%.3f\t", A[ii][jj]);
-
} cout << endl;
-
} cout << endl;
-
-
//return 0;
-
-
-
u = new double[ni*nj*nk];
-
v = new double[ni*nj*nk];
-
for (int i = 0; i < ni*nj*nk; i++)
-
u[i] = v[i] = 0;
-
for (int i = 0; i < ni; i++){
-
for (int j = 0; j < nj; j++){
-
u[i + ni*j + ni*nj*0] = v[i + ni*j + ni*nj*0] = g0;
-
u[i + ni*j + ni*nj*(nk-1)] = v[i + ni*j + ni*nj*(nk-1)] = g0;
-
}
-
}
-
for (int i = 0; i < ni; i++){
-
for (int k = 0; k < nk; k++){
-
u[i + ni*0 + ni*nj*k] = v[i + ni*0 + ni*nj*k] = g0;
-
u[i + ni*(nj-1) + ni*nj*k] = v[i + ni*(nj-1) + ni*nj*k] = g0;
-
}
-
}
-
for (int j = 0; j < nj; j++){
-
for (int k = 0; k < nk; k++){
-
u[0 + ni*j + ni*nj*k] = v[0 + ni*j + ni*nj*k] = g0;
-
u[ni-1 + ni*j + ni*nj*k] = v[ni-1 + ni*j + ni*nj*k] = g0;
-
}
-
}
-
do
-
{
-
for (int i = 0; i < ni*nj*nk; i++)
-
u[i] = v[i];
-
for (int i = 0; i < ni; i++){
-
for (int j = 0; j < nj; j++){
-
for (int k = 0; k < nk; k++){
-
if (i == 0 || i == ni-1 || j == 0 || j == nj-1 || k == 0 || k == nk-1){
-
continue;
-
}
-
int ii = i + ni*j + ni*nj*k;
-
double var = 0;
-
for (int jj = 0; jj < ii; jj++)
-
var += (A[ii][jj] * v[jj]);
-
for (int jj = ii + 1; jj < ni*nj*nk; jj++)
-
var += (A[ii][jj] * u[jj]);
-
v[ii] = (F[ii] - var) / A[ii][ii];
-
cout << v[ii] << endl;
-
}
-
}
-
}
-
}
-
while (!converge());
-
-
//Вывод в файл
-
freopen("MKE.txt", "w", stdout);
-
for (int i = 0; i < ni; i++){
-
for (int j = 0; j < nj; j++){
-
for (int k = 0; k < nk; k++){
-
k = nk/2;
-
cout << h*i << "\t" << h*j << "\t" << h*k << "\t" << u[i + ni*j + ni*nj*k] << endl;
-
break;
-
}
-
} cout << endl;
-
}
-
}
22
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
