ode15s struggle: Implement ODE solver within GUI (without GUIDE)

2 Ansichten (letzte 30 Tage)
Kristyna
Kristyna am 2 Jan. 2022
Beantwortet: Kristyna am 3 Jan. 2022
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 Kommentare
Kristyna
Kristyna am 3 Jan. 2022
Bearbeitet: Kristyna am 3 Jan. 2022
Ok, update: the parameters are passed inside "model" function just fine, last thing that is calculated correctly is F(1) inside this function. The rest is not calculated & the ode15s is not working properly, which is weird as it works just fine standing alone - I uploaded one more file (reaktor_BR15K.m) which shows what the code should be doing (without major user interaction).
My problem seems similar to this question, except the nesting does not work for me - may be there's a nuance I cannot see due to lack of experience..?
Kristyna
Kristyna am 3 Jan. 2022
@Cris LaPierre I uploaded the "reaktor_BR15K.m" to show how the whole thing should work.
With th GUI, I want to allow the user to change values of some of the parameters and recalculate the integration accordingly.

Melden Sie sich an, um zu kommentieren.

Antworten (1)

Kristyna
Kristyna am 3 Jan. 2022
Found the error, thanks everyone for their time: I lost my head in the code!
Anyway, your inputs were incredibly useful and I will apply them to make my code more efficient.
Thanks again!!

Produkte


Version

R2020a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by