Need help converting berkeley madonna code to matlab code
5 Ansichten (letzte 30 Tage)
Ältere Kommentare anzeigen
Not sure where to start and how to convert the berkeley code to matlab
Here is the Berkeley Madonna code:
{Top model}
METHOD RK4
STARTTIME = 0
STOPTIME = 24
DT = 0.02
{Reservoirs Blood, Fat, NonFat and Liver}
d/dt (NF) = + Jnf
INIT NF = 0
d/dt (F) = + Jf
INIT F = 0
d/dt (B) = - Jnf - Jf - Jl + Jresp
INIT B = 0
d/dt (L) = + Jl - Jmetab
INIT L = 0
{Flows}
Jnf = Qnf*(Cb-Cnf)
Jf = Qf*(Cb-Cfv)
Jl = Ql*(Cb - Clv)
Jmetab = vmax*Cl/(Km + Cl)
{Replace squarepulse with " IF t>6h Jresp=0 "}
Jresp = Qp*(Ci - (Cb/Pb))*squarepulse(0,6)
{Functions}
Vnf = 1
Cnf = NF/Vnf
Qnf = 1
Vf = 1
Cf = F/Vf
Qf = 1
Cb = B/Vb
Vb = 1
Cfv = Cf/Pf
Pf = 20
Vl = 1
Cl = L/Vl
Ql = 1
vmax = 1
Km = 1
Pl = 2
Clv = Cl/Pl
Pb = 18
{Benzene conc in inhaled air}
Ci = 0.32
{alveolar ventilation rate}
Qp = 5.74
3 Kommentare
Cris LaPierre
am 18 Mai 2021
Berkley Madonna is a differential equation solver. You can likely implement your code using MATLAB's main ode solver, ode45.
Walter Roberson
am 18 Mai 2021
Bearbeitet: Walter Roberson
am 18 Mai 2021
SimBiology produces and solves differential equations -- it is one of the several options, one that might look closer to the Madonna code.
ode45() and related functions are other useful options for numeric solving. The ode functions whose name include 's' are useful for stiff systems.
It can also be useful to set up differential equations using the Symbolic Toolbox, and see if dsolve() can happen to solve them exactly, and if not then use odeFunction() to convert into functions for use with the numeric solvers.
Antworten (2)
Cris LaPierre
am 18 Mai 2021
Here's equivalent MATLAB code.
STARTTIME = 0;
STOPTIME = 24;
% initial values
NF = 0;
F = 0;
B = 0;
L = 0;
y0 = [NF; F; B; L];
% solve ode
[t,y] = ode45(@odefun,[STARTTIME STOPTIME],y0);
% Plot results
yyaxis left
plot(t,y(:,[1,3]))
ylabel("NF,B")
yyaxis right
plot(t,y(:,[2,4]))
ylabel("F,L")
legend(["NF","B","F","L"])
function ddt = odefun(t,y)
NF = y(1);
F = y(2);
B = y(3);
L = y(4);
% Constants
Vnf = 1;
Cnf = NF/Vnf;
Qnf = 1;
Vf = 1;
Cf = F/Vf;
Qf = 1;
Vb = 1;
Cb = B/Vb;
Pf = 20;
Cfv = Cf/Pf;
Vl = 1;
Cl = L/Vl;
Ql = 1;
vmax = 1;
Km = 1;
Pl = 2;
Clv = Cl/Pl;
Pb = 18;
% Benzene conc in inhaled air
Ci = 0.32;
% alveolar ventilation rate
Qp = 5.74;
% Flows
Jnf = Qnf*(Cb-Cnf);
Jf = Qf*(Cb-Cfv);
Jl = Ql*(Cb - Clv);
Jmetab = vmax*Cl/(Km + Cl);
% Replace squarepulse with " IF t>6h Jresp=0 "
if t<=6
Jresp = Qp*(Ci - (Cb/Pb));
else
Jresp = 0;
end
% Reservoirs Blood, Fat, NonFat and Liver
ddt(1,1) = Jnf;
ddt(2,1) = Jf;
ddt(3,1) = -Jnf - Jf - Jl + Jresp;
ddt(4,1) = Jl - Jmetab;
end
Compare those figures to the plots generated by BM
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/621893/image.png)
6 Kommentare
Jeremy Huard
am 19 Mai 2021
I'm attaching the SimBiology implementation of your model.
You can either open justing.sbproj in SimBiology and use the Model Analyzer App to run simulations (the file contains a pre-configured program to replicate Cris' plot), or you can use SimBiology commands to write scripts.
For instance:
% load model
sbioloadproject justin.sbproj
getequations(m1)
% run simulation
sd = sbiosimulate(m1);
[t,y] = selectbyname(sd, ["NF","B","F","L"]);
% Plot results
figure;
ax = axes(XLimitMethod='padded', YLimitMethod='padded');
yyaxis(ax, 'left');
plot(t,y(:,[1,2]))
ylabel("NF,B")
yyaxis(ax, 'right');
plot(t,y(:,[3,4]))
ylabel("F,L")
grid(ax,'on');
legend(["NF","B","F","L"],Location='eastoutside',Box='off')
I hope this helps.
Walter Roberson
am 20 Mai 2021
Those sound like useful features of SimBiology.
In Answers, we get a fair number of questions about ode45() and ode23s() in which people use if/else in discontinuous ways, and we have to keep talking to them about how that is (usually) mathematically incompatible with the ode*() functions. Is SimBiology an alternative we should be talking up for general ODE systems that have events and impulses?
... and would there just happen to be a conversion function that could take symbolic ODE with dirac or heaviside or piecewise and build and run an appropriate SimBiology system?
Arthur Goldsipe
am 18 Mai 2021
If you'd like to convert a Berkeley Madonna model to the equivalent SimBiology model, you can try using this converter from the File Exchange.
0 Kommentare
Siehe auch
Kategorien
Mehr zu Extend Modeling Environment finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!