# How to speed up the simulation time of a big system of ODEs?

4 views (last 30 days)
Gabriele Galli on 10 Sep 2020
Edited: Gabriele Galli on 11 Sep 2020
Hi Everyone,
I have a big system of equations (378 diff. equations) and I have to simulate it using an ODE solver (I am using ODE23s) but the simulation time is too long.
Is there a way to improve my code in order to have a small simulation time?
Attached you can find the code.
Thank you!
Gabri

J. Alex Lee on 11 Sep 2020
don't use syms

Gabriele Galli on 11 Sep 2020
Hi Alex,
Unfortunately using a paper-and-pencil approach for this problem is very difficult due to its dimension.
These 378 equations are obtained starting from a system of just 18 equations (and 18 variables) and 21 parameters and performing partial differentiation of each equation w.r.t. each of the 21 parameters that is why the final system is composed by 21x18=378 equations.
So far I solved this problem manually but that system was much smaller and easy to handle than this one. That is why, for this one I need to use symbolic calculations (or, maybe, there are other approaches that I don't know).
By the way, I though to convert my symbolic matrix into a double one (see code attached here) but, unfortunately, it contains the variables x1 through x378 that are not numerical values hence I cannot convert it.
Due you have any idea either on how I can perform this conversion or other approaches in order to obtain a faster simulation?
Thank you so much for your help!
J. Alex Lee on 11 Sep 2020
well, if you insist on the symbolics approach, i'm afraid i can't help because i don't know how to use the symbolic toolbox (and don't have it).
if you want to step back, the only anyone can help you is if you show the original 18 equations, and why/how you are taking partial derivatives w.r.t. parameters.
Gabriele Galli on 11 Sep 2020
Yes, sure. I mean maybe we can find a better solution. I am not familiar either with the symbolic toolbox. The equations are contained in the f matrix below:
r1 = k1*TGFb_TGFbR;
r2 = k2*TGFbR*TGFb;
r3 = k3*TGFb_TGFbR*(1-exp(-((t-k20)/(k21))^10));
r4 = k4*TGFb_TGFbR_P;
f=[r1 - r2;
r1 - r2;
- r1 + r2 - r3 + r4 + r6;
r3 - r4 - r5;
r5 - r6;
- r7 - r8 + r9;
r7 - r10 + r11 - r12 + r13 - r17 + r18;
-r12 + r13 - r14 + r15;
r10 - r11 - r16;
r12 - r13 - r19;
r8 - r9 + r20;
r16 + r21 - r22;
r17 - r18 - r20 - r21 + r22 - r23 + r24;
r19 + r23 - r24;
r14 - r15 - r23 + r24;
r25 - r26;
r26 - r27;
r28 - r29 - r5 + r6];
In particular, k1 through k21 are the parameters with nominal values known and TGFb_TGFbR through I_Smad are the variables, for semplicity I called them x1 thorough x18.
What I want to find is the sensitivity matrix Sxd:
Sxd=dfdx*Sx+dfdk;
%dfdx is the Jacobian of f wrt x(=vector contained the 18 parameters), dim(dfdx=18x18)
%dfdk= is the Jacobian of f wrt k(=vector contained the 21 parameters), dim(dfdk=18x21)
%Sx=sensitivity matrix=dx(t,k)/dk. This is what I want to calculate using the ODE...
... I for simplicity I created I called each of these x1, x2, ..., x378 and reshape...
... them in a matrix 18x21. So basically, each line of the matrix contains the sensitivity...
... of a variable x wrt to each of the 21 parameters.
%Sxd is the first derivative wrt t of Sx. dim(Sxd)=18x21. For this reason, I reshape...
... this matrix in a column vector (dim=378x1) in order to be a suitable input for the ODE solver.
The time span I want to simulate is [0, 180] and we can assume all the 378 initial conditions = 0 except for:
IC(1)=1; IC(2)=1; IC(6)=40.98; IC(8)=34.15; IC(11)=19.02; IC(15)=15.85;