Форум программистов, компьютерный форум, киберфорум
Наши страницы
Matlab
Войти
Регистрация
Восстановить пароль
Показать сообщение отдельно
milka-milka
0 / 0 / 0
Регистрация: 19.11.2013
Сообщений: 3
0

Решение системы диф.уравнений

19.11.2013, 11:25. Просмотров 668. Ответов 5
Метки (Все метки)

Всем привет!

Столкнулась с проблемой при решении системы обыкновенных диф.уравнений, прошу помощи.
Сама система представлена в файле
Решение системы диф.уравнений
.

М-файл:
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
function f=tst2(t,x)
 
Rc = 8.31; %универсальная газовая постоянная, Дж/мольК
E = 600*10^3; %энергия активации, Дж/моль
mu_al=26.981*10^(-3); %масса моля алюминия, кг/моль
mu_ox=101.962*10^(-3); %масса моля оксида алюминия, кг/моль
ak=0.5; %окислительный потенциал среды
m2=0.5; %const
c=852.5; %теплоемкость капли, Дж/кгК
delH=459*10^3; %тепловой эффект взаимодействия газа окислителя с металлом, Дж/моль
a=100; %коэффициент теплоотдачи, Вт/м2К
%eps_ox=0.2; %степень черноты оксида алюминия
%eps_sur=0.1; %степень черноты окружающей среды
eps_gen=0.07; %степень черноты приведенная
k1=10^4; %const
ut=10^(-2); %скорость горения топлива, м/с
sigma=5.67*(10^(-8)); %постоянная Стефана-Больцмана, Вт/м2К4
T_sr = 1000; %температура окружающей среды, К
ro_al = 2.6989*10^3 ;%плотность алюминия, кг/м3
ro_ox = 3.99*10^3;%плотность оксида алюминия, кг/м3
 
A=400; %варьируемая переменная 
 
R1 = 10^(-6); %радиус сферической капли, 10 микрон, м
theta_grad = 10; %краевой угол смачивания, градусы
theta = theta_grad*pi/180; %краевой угол смачивания, радианы
R = R1*sin(theta); %радиус сегмента на поверхности с конд.фазой
h = R1*(1-cos(theta)); %высота капли
V = pi*(R*h^2-(1/3)*h^3); %объем капли целиком
Sp = 2*pi*R*h; %площадь поверхности капли 
%V_ox/V_al = 1.45; %отношение объема оксида к объему алюминия
V_ox = 1.45*V/2.45; %объем оксида, м3 (начальное значение оксида)
l = V_ox/Sp; %толщина пленки, м
V_al = V - V_ox; %объем алюминия
 
m_al = V_al*ro_al; %начальная масса алюминия
m_ox = V_ox*ro_ox; %начальная масса оксида
 
koef1 = -mu_al*Sp*A*ak^m2;
koef2 = -E/Rc;
koef3 = mu_ox*A*Sp*ak^m2;
koef4 = A*ak^m2*delH*Sp;
koef5 = R^2*pi*k1*ut;
koef6 = eps_gen*sigma;
koef7 = sin(theta)/(1-cos(theta));
koef8 = pi*(koef7-1/3);
 
%Sp = (2*pi*koef7*((x(1)/ro_al+x(2)/ro_ox)/koef8)^(2/3))
 
f = [koef1*x(3)*exp(koef2/x(3))*(1/x(4)); %dMal
     koef3*x(3)*exp(koef2/x(3))*(1/(2*x(4))); %dMox
    ((koef4*x(3)*exp(koef2/x(3))*(1/x(4))) - (koef5*x(3)/10 + koef6*Sp*(x(3)^4-T_sr^4) + a*Sp*(x(3)-T_sr))-c*x(3)*((koef1*x(3)*exp(koef2/x(3))*(1/x(4)))+(koef3*x(3)*exp(koef2/x(3))*(1/(2*x(4))))))/(c*(x(1)+x(2))); %dT
    (koef3*x(3)*exp(koef2/x(3))*(1/(2*x(4))))/(ro_ox*Sp)]; %l
 
end
Решение:

Matlab M
1
2
3
4
5
6
7
8
9
10
T=[0, 3*10^(-3)];
x0=[1*10^(-19); 3*10^(-19); 2317; 4*10^(-9)];
[t,x]=ode23('tst2',T,x0);
plot(t,x(:,1)),grid,title('dif'),legend('Mal')
figure
plot(t,x(:,2)),grid,title('dif'),legend('Mal2o3')
figure
plot(t,x(:,3)),grid,title('dif'),legend('Tm')
figure
plot(t,x(:,4)),grid,title('толщина'),legend('l')
Результаты получаются неправдоподобными: Температура растет до 10^4 степени (что, конечно, является ложью), + значение Mal переходит в отрицательные значения (что также недопустимо).

Уточните, пожалуйста, где задавать ограничения для значений переменных, если это возможно?
Коэффициент А в данной системе -- варьируемый, но вариации не приводят к желаемому результату.

Заранее спасибо за любой отклик!

Вернуться к обсуждению:
Решение системы диф.уравнений
0
Надоела реклама? Зарегистрируйтесь и она исчезнет полностью.
Similar
Эксперт
41792 / 34177 / 6122
Регистрация: 12.04.2006
Сообщений: 57,940
19.11.2013, 11:25
Готовые ответы и решения:

Решение системы диф.уравнений (Движение снаряда)
Здравствуйте, уважаемые форумчане, помогите пожалуйста разобраться в следующей...

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

Решение диф. Уравнений
Друзья, хелп ми плиз... Найти аналитическое решение уравнения : y′ + y2 =...

Обычные системы диф. уравнений
Доброе время суток. Подскажите пожалуйста как нужно правильно расписать систему...

Решение диф.уравнений Mathlab
Помогите, пожалуйста, написать программу в Mathlab

5
Другие темы раздела
Matlab Написать функцию, описывающую формулу вычисления длины кривой http://www.cyberforum.ru/matlab/thread1011710.html
x= y= sum(abs(diff(complex(x,y)))) не подскажите, правильно ли это?
Matlab Задача на определение коэффициентов полинома Здравствуйте, помогите с задачкой, не знаю как по таким данным определить коэффициенты полинома...Текст задачи прилагаю, задача номер 2 http://www.cyberforum.ru/matlab/thread1011699.html
Matlab Добавление символов в начало строки
Есть массив чисел перевожу его в двоичную систему и записываю в массив строк for i=1:popS popN(i)=randint(1,1,); os{i}=dec2bin(popN(i)); end Нужно чтобы все элементы массива...
Решение СЛАУ по методу Якоби Matlab
У меня проблема. Есть код, по которому должна вычисляться СЛАУ по методу Якоби с заданной точностью, но он вычисляет корень до конца, т.е. делает 1113 подсчётов и не может подсчитать корень с...
Matlab Динамическое формирование анонимной функции http://www.cyberforum.ru/matlab/thread1011523.html
Собственно, узнал я благодаря форуму, что такое динамическая функция, спасибо, пригодилось, но теперь у меня очередная проблема. Я получаю с помощью своей небольшой функции коэффициента ряда Фурье...
Matlab Рассчитать объем тела заданного точками Добрый день. Есть трехмерное тело, примерно такое как на прикрепленном рисунке. Известны координаты его точек, которые заданы отдельными файлами x.txt, y.txt, z.txt для соответствующих координат.... подробнее
КиберФорум - форум программистов, компьютерный форум, программирование
Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2018, vBulletin Solutions, Inc.
Рейтинг@Mail.ru