ode15s struggle: Implement ODE solver within GUI (without GUIDE)
2 ビュー (過去 30 日間)
古いコメントを表示
Hi, I am building a GUI wrapper for my chemical engineering project.
I wrote the calculation functions separately to be sure they work correctly and everything worked fine. Now, when I wrapped everything in a GUI, I started having trouble with ODE solver which apparently is not working. I am quite a beginner using Matlab, only use it occasionally for study purposes, never for a GUI before but I suspect that the trouble is with passing the parameters between app data and the model function (code attached) as they change with user inputs.
I tried multiple ways to solve this already - external function file for "model", write it as internal function as well as nested function inside ''chlazeni" function and nothing worked neither have I found any useful answer among the Community.
Thanks in advance for any help!
*****************************************
I am running the "chlazeni" function from "reaktordatabase" function, which is a ValueChangedFcn of a list box
function reaktordatabase(hObject,eventData,handles,fig)
%'BR15K','BR4K','BR2K','BR200','BR300','SHT15','BR1440'
mdot = get(handles.edt1,'Value'); %[kg/s], hm.prutok chl.vody
setappdata(fig,'mdot',mdot)
tIN = get(handles.edt2,'Value'); %[°C], teplota vstupujici chl.vody
setappdata(fig,'tIN',tIN)
otacky = get(handles.edt3,'Value');
setappdata(fig,'otacky',otacky);
lbxValue = get(handles.lbx,'Value');
if strcmp(lbxValue,'BR15K')
disp(lbxValue)
V = 15; %[m3], objem reaktoru
D = 2.0; %[m], prumer reaktoru
H = 5.1; %[m], vyska reaktoru celk.
h = 4.5; %[m], vyska duplikace
A = 4.115; %[m2], teplosmenna plocha
VB = 6; %[m3], objem vsadky
Q0 = 198595.470/1000; %[kW], teplo k odebrani
setappdata(fig,'V',V)
setappdata(fig,'D',D)
setappdata(fig,'H',H)
setappdata(fig,'h',h)
setappdata(fig,'A',A)
setappdata(fig,'VB',VB)
setappdata(fig,'Q0',Q0)
pomocne_vypocty(handles,fig)
chlazeni(handles,fig)
the "chlazeni" function itself looks like:
function chlazeni(handles,fig)
data = getappdata(fig);
N = data.N;
tstep = data.tstep;
% pocatecni podminky
tOUT0 = 25*ones(N,1);
setappdata(fig,'tOUT0',tOUT0)
tB0 = get(handles.edt2,'Value');
setappdata(fig,'tB0',tB0)
y0 = [tB0; tOUT0];
setappdata(fig,'y0',y0)
pars = [data.mdot data.cpW data.alfaW data.d data.lambda316L data.A,...
data.V data.alfaB data.rhoB data.mB data.cpB data.tIN,...
data.N data.tstep data.Vdot]';
assignin('base','pars',pars);
% integrace
opt = odeset('Events',@EventsFcn);
tEnd = 16*3600; %integracni cas, v sekundach
tspan = 0:tstep:tEnd;
[tt, yy] = ode45(@(t,y) model(t,y,pars),tspan,y0,opt);
tt = tt/3600; %integracni cas, v hodinach
[~,Q] = cellfun(@(t,y) model(t,y.',pars), num2cell(tt), num2cell(yy,2),'uni',0);
Q = cell2mat(Q);
Qcumul = trapz(tt,Q/1000) %[kW]
setappdata(fig,'Qcumul',Qcumul)
set(handles.edt5,'Value',Qcumul)
M = size(yy);
M = M(2);
t_stred = mean(yy(:,M))
setappdata(fig,'t_stred',t_stred)
t_max = max(yy(:,M))
setappdata(fig,'t_max',t_max)
% plotting results
for k = 2:4:M
plot(handles.ax1,tt,yy(:,k))
hold(handles.ax1,'on')
end
hold(handles.ax1,'off');
axis(handles.ax1,'tight');
xlabel(handles.ax1,'Čas [h]');
ylabel(handles.ax1,'Teplota výstupní chl.vody [°C]');
txt = {['t_max: ' num2str(t_max)],['t_střední: ' num2str(t_stred)]};
handles.str = uilabel('Parent',fig,'Position',[700 390 70 45],'Text',txt,...
'HorizontalAlignment','right','VerticalAlignment','center');
plot(handles.ax2,tt,yy(:,1));
axis(handles.ax2,'tight');
xlabel(handles.ax2,['Čas [h] - doba chlazení: ',num2str(tt(end)),' h'])
ylabel(handles.ax2,'Teplota vsádky [°C]');
axis(handles.ax2,'tight');
plot(handles.ax3,tt,Q/1000);
axis(handles.ax3,'tight');
xlabel(handles.ax3,'Čas [h]');
ylabel(handles.ax3,'\Sigma Q [kW]');
axis(handles.ax3,'tight');
end
and the original "model" function that is the input in ode15s, as was before using GUI (=working code) is:
function [dTdt,Q] = model(~,Tvec,pars)
% rozbaleni parametru
% pars = [mdot cpW alfaW d lambda A V alfaB rhoB mB cpB tIN N tstep Vdot]';
mdot = pars(1);
cpW = pars(2);
alfaW = pars(3);
d = pars(4);
lambda316L = pars(5);
A = pars(6);
V = pars(7);
alfaB = pars(8);
rhoB = pars(9);
mB = pars(10);
cpB = pars(11);
tIN = pars(12);
N = pars(13);
tstep = pars(14);
Vdot = pars(15);
tRef = 0;
mW = ((0.0024*4.5)/(Vdot/3600))*mdot; % Amax(m2)*výška duplikace(m)*rho(kg/m3)*čas(s) "zádrž vody v plášti"
K = 1/(d/lambda316L+1/alfaW+1/alfaB);
F = zeros(N+1,1);
Q = zeros(N+1,1);
F(1) = mdot*cpW*(tIN-tRef);
Q(1) = K*A/N*(Tvec(1)-tIN);
for i=2:N+1
Q(i) = K*A/N*(Tvec(1)-Tvec(i));
F(i) = mdot*cpW*(Tvec(i)-tRef);
end
for j = 2:N+1
dTdt(j) = (F(j-1)-F(j)+Q(j))/((mW/N)*cpW);
end
Q = sum(Q);
dTdt(1) = -Q/(cpB*mB);
dTdt = dTdt';
6 件のコメント
回答 (1 件)
参考
カテゴリ
Help Center および File Exchange で Ordinary Differential Equations についてさらに検索
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!