175 / 141 / 50
Регистрация: 07.02.2014
Сообщений: 480
1

На выходе программы получаю искривление теплового фронта

03.05.2015, 10:19. Показов 397. Ответов 5
Метки нет (Все метки)

Здравствуйте!

Помогите найти ошибку в решении... На выходе получаю искривление теплового фронта (причем с одной стороны все в порядке T(:,Ny) = T(:,Ny-1), а вот с другой проблемы):
Matlab M
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
clc, clear all
% Дано
Lx = 10; Ly = 10; hx = 0.1; hy = 0.1; Nx = Lx/hx; Ny = Ly/hy;
lamda_fld = 1.35; lamda_sld = 2.25; rho_fld = 1250; rho_sld = 1550;
C_fld = 3325e03; C_sld = 2500e03; a_sld = lamda_sld/(rho_sld*C_sld);
a_fld = lamda_fld/(rho_fld*C_fld); T0 = -2.5; tau = 3600; t_end = 3600*24*31; Lq = 155e06; Tbf = 0;
delta = 0.1;
% Начальная температура
T = ones(Nx,Ny)*T0;
% Граничные условия
T(1,:) = 10;
T(Nx,:) = -2;
 
time = 0;
 
while time<t_end
    
    time = time + tau;
    % Учитываем нелинейность теплоемкости и теплопроводности
    for i = 1:Nx
        for j = 1:Ny
            if T(i,j) < (Tbf-delta)
                C(i,j) = C_sld;
                lamda(i,j) = lamda_sld;
            elseif T(i,j) >= (Tbf-delta) && T(i,j) <= (Tbf+delta)
                C(i,j) = 0.5*(C_fld+C_sld) + 0.5*Lq/delta;
                lamda(i,j) = lamda_sld + 0.5*(lamda_fld-lamda_sld)/delta*(T(i,j) - (Tbf - delta));
            else
                C(i,j) = C_fld;
                lamda(i,j) = lamda_fld;
            end
        end
    end
    % Здаем нулевой теплопоток на границах
    T(:,1) = T(:,2);
    T(:,Ny) = T(:,Ny-1);
    % Считаем уравнение теплопроводности
    for i = 2:Nx-1
        for j = 2:Ny-1
            T(i,j) = T(i,j) + tau/hx^2*(lamda(i,j)/C(i,j))*(T(i+1,j) - 2*T(i,j) + T(i-1,j))+...
                tau/hy^2*(lamda(i,j)/C(i,j))*(T(i,j+1) - 2*T(i,j) + T(i,j-1));
        end
    end
end
X = 0:Lx/(Nx-1):Lx; Y = 0:Ly/(Ny-1):Ly;
contourf(X, Y, T) % визуализируем
% colorbar, 
grid on % закрепляем сетку и цвет.шкалу
title('Температурное поле массива грунта на 30 год эксплуатации (на 15 сентября)') %название
xlabel('Ширина, м'); ylabel('Глубина, м'); zlabel('Температура, deg[C]') %название осей
axis auto, caxis([-5 5]) %X и Y оси автоматически, T - в заданных пределах
На выходе программы получаю искривление теплового фронта
0

Помощь в написании контрольных, курсовых и дипломных работ здесь.

Programming
Эксперт
94731 / 64177 / 26122
Регистрация: 12.04.2006
Сообщений: 116,782
03.05.2015, 10:19
Ответы с готовыми решениями:

На выходе функции получаю None
Добрый день! Подскажите, пожалуйста, по такому вопросу. Есть функция, давление насыщенных паров...

Enum - почему на выходе получаю 0, а не 1 ?
начал писать на шарпе, вот появился вопрос есть перечисление enum days { ...

На выходе получаю лист одинаковых элементов
Доброго времени суток. Сразу к сути: в цикле создаются некие объекты и заносятся в лист. При дебаге...

Мусор на выходе программы
Доброго времени суток. Дана программа: определённое поле, которое я задаю как двоичный массив,...

5
175 / 141 / 50
Регистрация: 07.02.2014
Сообщений: 480
08.08.2015, 19:19  [ТС] 2
Уважаемые форумчане, вопрос до сих пор актуален, потому что ошибку я так и не нашел. Есть ли какие-то соображения по этому поводу?
0
kRosis
09.08.2015, 02:12
  #3

Не по теме:

А что такое hx hy?

0
175 / 141 / 50
Регистрация: 07.02.2014
Сообщений: 480
09.08.2015, 12:14  [ТС] 4
kRosis, Lx, Ly - длины осей, hx, hy - шаги сетки, Nx, Ny - кол-во узлов сетки

Добавлено через 7 минут
Вот код с небольшими пояснениями что к чему:
Кликните здесь для просмотра всего текста
Matlab M
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
clc, clear all, close all
 
%% Дано
Lx = 10; Ly = 10; %Длины осей x и y
hx = 0.1; hy = 0.1; %Шаги сетки по x и y
Nx = Lx/hx; Ny = Ly/hy; %Кол-во узлов сетки
lamda_fld = 1.35; lamda_sld = 2.25; %Теплопроводность в талом и мерзлом состоянии
rho_fld = 1250; rho_sld = 1550; %Плотность в талом и мерзлом состоянии
C_fld = 3325e03; C_sld = 2500e03; %Теплоемкость в талом и мерзлом состоянии
a_sld = lamda_sld/(rho_sld*C_sld); %Температуропроводность в твердом состоянии
a_fld = lamda_fld/(rho_fld*C_fld); %Температуропроводность в талом состоянии
T0 = -2.5; %Начальная температура массива
tau = 3600; %Шаг по времени
t_end = 3600*24*31; %Конечное время расчета
Lq = 155e06; %Теплота фазового перехода
Tbf = 0; %Температура фазового перехода
delta = 0.1; %Зона фазового перехода в методе эффективной теплоемкости
 
%% Начальные условия
T = ones(Nx,Ny)*T0; 
lamda = zeros(Nx,Ny);
C = zeros(Nx,Ny);
 
%% Граничные условия
T(1,:) = 10;
T(Nx,:) = -2;
 
time = 0;
 
%% Решаем
while time<t_end
    
    time = time + tau;
    % Учитываем нелинейность теплоемкости и теплопроводности
    for i = 1:Nx
        for j = 1:Ny
            if T(i,j) < (Tbf-delta)
                C(i,j) = C_sld;
                lamda(i,j) = lamda_sld;
            elseif T(i,j) >= (Tbf-delta) && T(i,j) <= (Tbf+delta)
                C(i,j) = 0.5*(C_fld+C_sld) + 0.5*Lq/delta;
                lamda(i,j) = lamda_sld + 0.5*(lamda_fld-lamda_sld)/delta*(T(i,j) - (Tbf - delta));
            else
                C(i,j) = C_fld;
                lamda(i,j) = lamda_fld;
            end
        end
    end
    % Здаем нулевой теплопоток на границах
    T(:,1) = T(:,2);
    T(:,Ny) = T(:,Ny-1);
    % Считаем уравнение теплопроводности
    for i = 2:Nx-1
        for j = 2:Ny-1
            T(i,j) = T(i,j) + tau/hx^2*(lamda(i,j)/C(i,j))*(T(i+1,j) - 2*T(i,j) + T(i-1,j))+...
                tau/hy^2*(lamda(i,j)/C(i,j))*(T(i,j+1) - 2*T(i,j) + T(i,j-1));
        end
    end
end
 
%% Визуализируем
X = 0:Lx/(Nx-1):Lx; Y = 0:Ly/(Ny-1):Ly;
contourf(X, Y, T) % визуализируем
colorbar, grid on % закрепляем сетку и цвет.шкалу
title('Температурное поле массива грунта') %название
xlabel('Ширина, м'); ylabel('Глубина, м'); zlabel('Температура, deg[C]') %название осей
axis equal, caxis([-5 5]) %X и Y оси автоматически, T - в заданных пределах


Выяснил, что не зависимо от длины осей, искривление происходит на крайних левых 3 метрах (см.рис.)

На выходе программы получаю искривление теплового фронта
0
125 / 124 / 10
Регистрация: 09.11.2010
Сообщений: 200
09.08.2015, 14:17 5
С уменьшением шага сетки hx и hy растёт нелинейность полученного графика. Попробуйте задать эти величины в районе 0.08.

При шаге 0.5 получается практически пригодный график. Возможно выражение для T(i,j) оперирует числами на грани точности допустимой для Matlab.

При шаге 0.05 решение вообще не определено.
0
175 / 141 / 50
Регистрация: 07.02.2014
Сообщений: 480
10.08.2015, 19:36  [ТС] 6
Понял, в чем была ошибка! Решаю явным методом, а шаг по времени взял такой, что решение не устойчивое. Сделав шаг по времени меньше, все сошлось. Спасибо всем за внимание
0
IT_Exp
Эксперт
87844 / 49110 / 22898
Регистрация: 17.06.2006
Сообщений: 92,604
10.08.2015, 19:36

Подтверждение при выходе из программы
Ребят как здесь сделать подтверждение выхода из программы? private void...

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

Ошибка при выходе из программы
Помогите пожалуйста найти ошибку, которая возникает при выходе из программы. Error -...

Ошибка при выходе из программы
Работаю с экселем.Открываю книгу,обрабатываю.Если открываю документ из моей программы,то после...


Искать еще темы с ответами

Или воспользуйтесь поиском по форуму:
6
Ответ Создать тему
Опции темы

КиберФорум - форум программистов, компьютерный форум, программирование
Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2021, vBulletin Solutions, Inc.