search

Open-Systems-Pharmacology / OSPSuite-R

R package for the OSPSuite
https://www.open-systems-pharmacology.org/OSPSuite-R/
Other
29 stars 12 forks source link

Export a model as (R-)ODE system: define the export format #296

Open Yuri05 opened 4 years ago

Yuri05 commented 4 years ago

Define export format for exporting of a Model to the R-ODE (similar to the current Matlab ODE export)

Here an example of Matlab ODE-Export: ODEMatlabExportExample.zip

All important ingredients of Matlab ODE export are below. For R the export might be structured differently.

Initial Values ``` function y = ODEInitialValues y(1) = 0; % S4|Organism|VenousBlood|Plasma|C1 y(2) = 0; % S4|Organism|VenousBlood|BloodCells|C1 ... y(154) = 0; % S4|Organism|Saliva|Saliva|C1 y(155) = 0; % S4|Organism|Saliva|SalivaGland|C1 y=y'; ```
RHS function ``` function dy = ODERHSFunction(Time,y) P_0 = 9500; % S4|Organism|SA proportionality factor P_1 = 1; % S4|Organism|Surface area scaling exponent ... global P_29; P_30 = 1; % S4|Organism|Plasma protein scale factor ... global P_46; ... global P_272; if isempty(P_272) P_272 = 1; % S4|Organism|Lumen|OralApplicationsEnabled end P_273 = 1; % S4|Organism|Lumen|Paracellular absorption sink condition ... P_1219 = 0; % S4|Neighborhoods|Lumen_coltrans_ColonTransversum_cell|C1|Sum of active process rates P_1220 = 0; % S4|Neighborhoods|Lumen_coldesc_ColonDescendens_pls|C1|Sum of active process rates P_1221 = 0; % S4|Neighborhoods|Lumen_coldesc_ColonDescendens_pls|C1|Sum of passive process rates P_1222 = 0; % S4|Neighborhoods|Lumen_coldesc_ColonDescendens_cell|C1|Sum of active process rates P_1223 = 0; % S4|Neighborhoods|Lumen_colsigm_ColonSigmoid_pls|C1|Sum of active process rates P_1224 = 0; % S4|Neighborhoods|Lumen_colsigm_ColonSigmoid_pls|C1|Sum of passive process rates P_1225 = 0; % S4|Neighborhoods|Lumen_colsigm_ColonSigmoid_cell|C1|Sum of active process rates ... P_1284 = (IIf(((P_39-P_40) >= 0),(P_39-P_40),0)); % S4|Organism|Specific clearance (endosome) P_1285 = (P_1114+(Time*P_1115)); % S4|Organism|Age P_1286 = (EvalParameter(P_46, Time, y)*P_43); % S4|Organism|VenousBlood|Weight (tissue) P_1287 = ((1-P_1113)*EvalParameter(P_46, Time, y)); % S4|Organism|VenousBlood|Plasma|Volume P_1288 = (P_1113*EvalParameter(P_46, Time, y)); % S4|Organism|VenousBlood|BloodCells|Volume P_1289 = (EvalParameter(P_51, Time, y)*P_47); % S4|Organism|ArterialBlood|Weight (tissue) P_1290 = ((1-P_1113)*EvalParameter(P_51, Time, y)); % S4|Organism|ArterialBlood|Plasma|Volume P_1291 = (P_1113*EvalParameter(P_51, Time, y)); % S4|Organism|ArterialBlood|BloodCells|Volume P_1292 = (EvalParameter(P_86, Time, y)*P_52); % S4|Organism|Bone|Weight (tissue) P_1293 = (1-(P_68+P_63)); % S4|Organism|Bone|Fraction intracellular P_1294 = (P_63*(1-P_1113)*EvalParameter(P_86, Time, y)); % S4|Organism|Bone|Plasma|Volume P_1295 = (P_63*P_1113*EvalParameter(P_86, Time, y)); % S4|Organism|Bone|BloodCells|Volume P_1296 = (P_68*EvalParameter(P_86, Time, y)); % S4|Organism|Bone|Interstitial|Volume P_1297 = P_3; % S4|Organism|Bone|Interstitial|pH P_1298 = EvalParameter(P_86, Time, y); % S4|Organism|Bone|Intracellular|Volume of protein container P_1299 = (6*P_89); % S4|Organism|Gallbladder|Time to complete gallbladder refilling P_1300 = (IIf((P_88 > 0),(y(9)/P_88),0)); % S4|Organism|Gallbladder|C1|Concentration P_1301 = (EvalParameter(P_123, Time, y)*P_91); % S4|Organism|Brain|Weight (tissue) P_1302 = (1-(P_106+P_101)); % S4|Organism|Brain|Fraction intracellular P_1303 = (P_101*(1-P_1113)*EvalParameter(P_123, Time, y)); % S4|Organism|Brain|Plasma|Volume P_1304 = (P_101*P_1113*EvalParameter(P_123, Time, y)); % S4|Organism|Brain|BloodCells|Volume P_1305 = (P_106*EvalParameter(P_123, Time, y)); % S4|Organism|Brain|Interstitial|Volume P_1306 = P_3; % S4|Organism|Brain|Interstitial|pH P_1307 = EvalParameter(P_123, Time, y); % S4|Organism|Brain|Intracellular|Volume of protein container ... P_1345 = (P_676*0.08699999999999999); % S4|Organism|Lumen|Duodenum|Volume (gut wall) P_1346 = (IIf((P_1238 ~= P_1243),(1/(1+(10^(P_1238*(P_1249-P_302))))),1)); % S4|Organism|Lumen|Duodenum|C1|pKa_pH_WS_lum_F1 P_1347 = (IIf((P_1239 ~= P_1243),(1/(1+(10^(P_1239*(P_1250-P_302))))),1)); % S4|Organism|Lumen|Duodenum|C1|pKa_pH_WS_lum_F2 P_1348 = (IIf((P_1240 ~= P_1243),(1/(1+(10^(P_1240*(P_1251-P_302))))),1)); % S4|Organism|Lumen|Duodenum|C1|pKa_pH_WS_lum_F3 ... P_1619 = (IIf((P_1238 ~= P_1243),(1/(1+(10^(P_1238*(P_1249-P_332))))),1)); % S4|Neighborhoods|Lumen_lje_LowerJejunum_cell|C1|pKa_pH_WS_lum_F1 P_1620 = (IIf((P_1239 ~= P_1243),(1/(1+(10^(P_1239*(P_1250-P_332))))),1)); % S4|Neighborhoods|Lumen_lje_LowerJejunum_cell|C1|pKa_pH_WS_lum_F2 P_1621 = (IIf((P_1240 ~= P_1243),(1/(1+(10^(P_1240*(P_1251-P_332))))),1)); % S4|Neighborhoods|Lumen_lje_LowerJejunum_cell|C1|pKa_pH_WS_lum_F3 P_1622 = (y(35)+y(38)+y(41)+y(44)+y(47)+y(50)+y(53)+y(56)+y(59)+y(62)+y(65)); % S4|Neighborhoods|Lumen_sto_Lumen_duo|C1|Oral mass absorbed P_1623 = (((P_488*909.090909090909)^0.75)*100000); % S4|Neighborhoods|Stomach_int_Stomach_cell|Surface area (interstitial/intracellular) P_1624 = (P_1113*P_22*0.6*P_488*P_466); % S4|Neighborhoods|Stomach_pls_Stomach_bc|Surface area (blood cells/plasma) ... P_2588 = (y(36)/P_1807); % S4|Organism|Lumen|UpperJejunum|Fraction of geometric volume filled with liquid P_2589 = ((P_1256*(P_2563/P_2292))+P_1358); % S4|Organism|Lumen|UpperJejunum|C1|Solubility P_2590 = (((P_1806*P_305*P_1807)+(P_2294-(P_1817*P_320*P_1818)))/(P_320*P_1818)); % S4|Organism|Lumen|LowerJejunum|Absorption of liquid P_2591 = (y(39)/P_1818); % S4|Organism|Lumen|LowerJejunum|Fraction of geometric volume filled with liquid P_2592 = ((P_1256*(P_2563/P_2297))+P_1367); % S4|Organism|Lumen|LowerJejunum|C1|Solubility P_2593 = (((P_1817*P_320*P_1818)+(P_2299-(P_1828*P_335*P_1829)))/(P_335*P_1829)); % S4|Organism|Lumen|UpperIleum|Absorption of liquid P_2594 = (y(42)/P_1829); % S4|Organism|Lumen|UpperIleum|Fraction of geometric volume filled with liquid P_2595 = ((P_1256*(P_2563/P_2302))+P_1376); % S4|Organism|Lumen|UpperIleum|C1|Solubility P_2596 = (((P_1828*P_335*P_1829)+(P_2304-(P_1839*P_350*P_1840)))/(P_350*P_1840)); % S4|Organism|Lumen|LowerIleum|Absorption of liquid ... P_3203 = P_3192; % S4|Neighborhoods|Lumen_cae_Caecum_cell|C1|Sum of passive process rates P_3204 = P_3193; % S4|Neighborhoods|Lumen_colasc_ColonAscendens_cell|C1|Sum of passive process rates P_3205 = P_3194; % S4|Neighborhoods|Lumen_coltrans_ColonTransversum_cell|C1|Sum of passive process rates P_3206 = P_3195; % S4|Neighborhoods|Lumen_coldesc_ColonDescendens_cell|C1|Sum of passive process rates P_3207 = P_3196; % S4|Neighborhoods|Lumen_colsigm_ColonSigmoid_cell|C1|Sum of passive process rates dy(1) = ... (0-(P_2440*(1-P_1113)*P_1716)) ... +(0-(P_1590*P_1699*((P_2558*P_1716)-(P_2557*(P_1719/P_2806))))) ... +(P_1725*(1-P_1113)*P_1728) ... +(P_1736*(1-P_1113)*P_1738) ... +(P_1746*(1-P_1113)*P_1748) ... +(P_1756*(1-P_1113)*P_1758) ... +(P_1766*(1-P_1113)*P_1768) ... +(P_1777*(1-P_1113)*P_1779) ... +(P_1978*(1-P_1113)*P_1980) ... +(P_2003*(1-P_1113)*P_2005) ... +(IIf(P_896,((P_2432+P_2446)*(1-P_1113)*P_2434),0)) ... +(IIf(P_896,0,((P_2423+P_2446)*(1-P_1113)*P_2425))); dy(2) = ... (0-(P_2440*P_1113*P_1719)) ... +(P_1590*P_1699*((P_2558*P_1716)-(P_2557*(P_1719/P_2806)))) ... +(P_1725*P_1113*P_1731) ... +(P_1736*P_1113*P_1741) ... +(P_1746*P_1113*P_1751) ... +(P_1756*P_1113*P_1761) ... +(P_1766*P_1113*P_1771) ... +(P_1777*P_1113*P_1782) ... +(P_1978*P_1113*P_1983) ... +(P_2003*P_1113*P_2008) ... +(IIf(P_896,((P_2432+P_2446)*P_1113*P_2436),0)) ... +(IIf(P_896,0,((P_2423+P_2446)*P_1113*P_2427))); dy(3) = ... (0-(P_1591*P_1699*((P_2558*P_1721)-(P_2557*(P_1724/P_2806))))) ... +(P_2440*(1-P_1113)*P_1970) ... +(0-(P_1725*(1-P_1113)*P_1721)) ... +(0-(P_1736*(1-P_1113)*P_1721)) ... +(0-(P_1746*(1-P_1113)*P_1721)) ... +(0-(P_1756*(1-P_1113)*P_1721)) ... +(0-(P_1766*(1-P_1113)*P_1721)) ... +(0-(P_1777*(1-P_1113)*P_1721)) ... +(0-(P_1937*(1-P_1113)*P_1721)) ... +(0-(P_1978*(1-P_1113)*P_1721)) ... +(0-(P_1988*(1-P_1113)*P_1721)) ... +(0-(P_2003*(1-P_1113)*P_1721)) ... +(0-(P_1925*(1-P_1113)*P_1721)) ... +(0-(P_2013*(1-P_1113)*P_1721)) ... +(0-(P_1916*(1-P_1113)*P_1721)) ... +(0-(P_2423*(1-P_1113)*P_1721)); dy(4) = ... (P_1591*P_1699*((P_2558*P_1721)-(P_2557*(P_1724/P_2806)))) ... +(P_2440*P_1113*P_1973) ... +(0-(P_1725*P_1113*P_1724)) ... +(0-(P_1736*P_1113*P_1724)) ... +(0-(P_1746*P_1113*P_1724)) ... +(0-(P_1756*P_1113*P_1724)) ... +(0-(P_1766*P_1113*P_1724)) ... +(0-(P_1777*P_1113*P_1724)) ... +(0-(P_1937*P_1113*P_1724)) ... +(0-(P_1978*P_1113*P_1724)) ... +(0-(P_1988*P_1113*P_1724)) ... +(0-(P_2003*P_1113*P_1724)) ... +(0-(P_1925*P_1113*P_1724)) ... +(0-(P_2013*P_1113*P_1724)) ... +(0-(P_1916*P_1113*P_1724)) ... +(0-(P_2423*P_1113*P_1724)); dy(5) = ... (0-(P_1593*P_1699*((P_2558*P_1728)-(P_2557*(P_1731/P_2806))))) ... +(0-(P_1699*P_1118*P_1594*(P_1728-(P_1733/P_2022)))) ... +(0-(P_1725*(1-P_1113)*P_1728)) ... +(P_1725*(1-P_1113)*P_1721); dy(6) = ... (0-(P_1725*P_1113*P_1731)) ... +(P_1593*P_1699*((P_2558*P_1728)-(P_2557*(P_1731/P_2806)))) ... +(P_1725*P_1113*P_1724); dy(7) = ... (0-((P_2453*P_2760*P_1733)-(P_3098*P_2761*P_2262))) ... +(P_1699*P_1118*P_1594*(P_1728-(P_1733/P_2022))); dy(8) = ... ((P_2453*P_2760*P_1733)-(P_3098*P_2761*P_2262)); dy(9) = ... (0-(IIf(P_87,(0.6931471805599453*y(9)*(P_90/P_89)),0))); dy(10) = ... (0-(P_1596*P_1699*((P_2558*P_1738)-(P_2557*(P_1741/P_2806))))) ... +(0-(P_1699*P_2457*P_1597*(P_1738-(P_1743/P_2023)))) ... +(0-(P_1736*(1-P_1113)*P_1738)) ... +(P_1736*(1-P_1113)*P_1721); dy(11) = ... (0-(P_1736*P_1113*P_1741)) ... +(P_1596*P_1699*((P_2558*P_1738)-(P_2557*(P_1741/P_2806)))) ... +(P_1736*P_1113*P_1724); dy(12) = ... (0-((P_2456*P_2762*P_1743)-(P_3099*P_2763*P_2264))) ... +(P_1699*P_2457*P_1597*(P_1738-(P_1743/P_2023))); dy(13) = ... ((P_2456*P_2762*P_1743)-(P_3099*P_2763*P_2264)); dy(14) = ... (0-(P_1599*P_1699*((P_2558*P_1748)-(P_2557*(P_1751/P_2806))))) ... +(0-(P_1699*P_1123*P_1600*(P_1748-(P_1753/P_2024)))) ... +(0-(P_1746*(1-P_1113)*P_1748)) ... +(P_1746*(1-P_1113)*P_1721); dy(15) = ... (0-(P_1746*P_1113*P_1751)) ... +(P_1599*P_1699*((P_2558*P_1748)-(P_2557*(P_1751/P_2806)))) ... +(P_1746*P_1113*P_1724); dy(16) = ... (0-((P_2460*P_2764*P_1753)-(P_3100*P_2765*P_2266))) ... +(P_1699*P_1123*P_1600*(P_1748-(P_1753/P_2024))); dy(17) = ... ((P_2460*P_2764*P_1753)-(P_3100*P_2765*P_2266)); dy(18) = ... (0-(P_1602*P_1699*((P_2558*P_1758)-(P_2557*(P_1761/P_2806))))) ... +(0-(P_1699*P_1126*P_1603*(P_1758-(P_1763/P_2025)))) ... +(0-(P_1756*(1-P_1113)*P_1758)) ... +(P_1756*(1-P_1113)*P_1721); dy(19) = ... (0-(P_1756*P_1113*P_1761)) ... +(P_1602*P_1699*((P_2558*P_1758)-(P_2557*(P_1761/P_2806)))) ... +(P_1756*P_1113*P_1724); dy(20) = ... (0-((P_2463*P_2766*P_1763)-(P_3101*P_2767*P_2268))) ... +(P_1699*P_1126*P_1603*(P_1758-(P_1763/P_2025))); dy(21) = ... ((P_2463*P_2766*P_1763)-(P_3101*P_2767*P_2268)); dy(22) = ... (0-(P_1605*P_1699*((P_2558*P_1768)-(P_2557*(P_1771/P_2806))))) ... +(0-(P_1699*P_1129*P_1606*(P_1768-(P_1773/P_2026)))) ... +(0-(P_1766*(1-P_1113)*P_1768)) ... +(P_1766*(1-P_1113)*P_1721); dy(23) = ... (0-(P_1766*P_1113*P_1771)) ... +(P_1605*P_1699*((P_2558*P_1768)-(P_2557*(P_1771/P_2806)))) ... +(P_1766*P_1113*P_1724); dy(24) = ... (0-((P_2466*P_2768*P_1773)-(P_3102*P_2769*P_2270))) ... +(P_1699*P_1129*P_1606*(P_1768-(P_1773/P_2026))); dy(25) = ... ((P_2466*P_2768*P_1773)-(P_3102*P_2769*P_2270)); dy(26) = ... (0-(P_1608*P_1699*((P_2558*P_1779)-(P_2557*(P_1782/P_2806))))) ... +(0-(P_1699*P_1132*P_1609*(P_1779-(P_1784/P_2027)))) ... +(0-(P_1699*P_2272*P_1133*P_1329*P_1779)) ... +(0-(P_1777*(1-P_1113)*P_1779)) ... +(P_1777*(1-P_1113)*P_1721); dy(27) = ... (0-(P_1777*P_1113*P_1782)) ... +(P_1608*P_1699*((P_2558*P_1779)-(P_2557*(P_1782/P_2806)))) ... +(P_1777*P_1113*P_1724); dy(28) = ... (0-((P_2469*P_2770*P_1784)-(P_3103*P_2771*P_2273))) ... +(P_1699*P_1132*P_1609*(P_1779-(P_1784/P_2027))); dy(29) = ... ((P_2469*P_2770*P_1784)-(P_3103*P_2771*P_2273)); dy(30) = ... (P_1699*P_2272*P_1133*P_1329*P_1779); dy(31) = ... (IIf(EvalParameter(P_272, Time, y),(0+((P_1786*P_2582)-P_2947)),0)); dy(32) = ... (0-(IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1338*y(32)*P_2821),0))); dy(33) = ... (IIf(EvalParameter(P_272, Time, y),(P_2947+((P_2284*P_2821)-(P_2949+(P_2584*y(33))))),0)); dy(34) = ... (IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1338*y(32)*P_2821),0)) ... +(0-(IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1795*y(34)*P_2823),0))) ... +(IIf(P_87,(0.6931471805599453*y(9)*(P_90/P_89)),0)) ... +(0-P_3013) ... +(0-P_3174); dy(35) = ... (P_3064+(P_3202-(P_1209+P_1210))); dy(36) = ... (IIf(EvalParameter(P_272, Time, y),(P_2949+((P_2289*P_2823)-(P_2950+(P_2587*y(36))))),0)); dy(37) = ... (0-P_2964) ... +(0-P_3164) ... +(IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1795*y(34)*P_2823),0)) ... +(0-(IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1806*y(37)*P_2825),0))); dy(38) = ... (P_3060+(P_3197-(P_1134+P_1135))); dy(39) = ... (IIf(EvalParameter(P_272, Time, y),(P_2950+((P_2294*P_2825)-(P_2951+(P_2590*y(39))))),0)); dy(40) = ... (0-P_2965) ... +(0-P_3166) ... +(IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1806*y(37)*P_2825),0)) ... +(0-(IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1817*y(40)*P_2827),0))); dy(41) = ... (P_3061+(P_3198-(P_1136+P_1137))); dy(42) = ... (IIf(EvalParameter(P_272, Time, y),(P_2951+((P_2299*P_2827)-(P_2952+(P_2593*y(42))))),0)); dy(43) = ... (IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1817*y(40)*P_2827),0)) ... +(0-(IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1828*y(43)*P_2829),0))) ... +(0-P_2966) ... +(0-P_3168); dy(44) = ... (P_3062+(P_3199-(P_1141+P_1142))); dy(45) = ... (IIf(EvalParameter(P_272, Time, y),(P_2952+((P_2304*P_2829)-(P_2953+(P_2596*y(45))))),0)); dy(46) = ... (IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1828*y(43)*P_2829),0)) ... +(0-(IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1839*y(46)*P_2831),0))) ... +(0-P_2975) ... +(0-P_3170); dy(47) = ... (P_3063+(P_3200-(P_1146+P_1147))); dy(48) = ... (IIf(EvalParameter(P_272, Time, y),(P_2953+((P_2309*P_2831)-(P_2954+(P_2599*y(48))))),0)); dy(49) = ... (IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1839*y(46)*P_2831),0)) ... +(0-(IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1850*y(49)*P_2833),0))) ... +(0-P_3176); dy(50) = ... (P_1212+(P_3203-(P_1211+P_1213))); dy(51) = ... (IIf(EvalParameter(P_272, Time, y),(P_2954+((P_2314*P_2833)-(P_2955+(P_2602*y(51))))),0)); dy(52) = ... (IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1850*y(49)*P_2833),0)) ... +(0-(IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1861*y(52)*P_2835),0))) ... +(0-P_3178); dy(53) = ... (P_1215+(P_3204-(P_1214+P_1216))); dy(54) = ... (IIf(EvalParameter(P_272, Time, y),(P_2955+((P_2319*P_2835)-(P_2956+(P_2605*y(54))))),0)); dy(55) = ... (IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1861*y(52)*P_2835),0)) ... +(0-(IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1872*y(55)*P_2837),0))) ... +(0-P_3180); dy(56) = ... (P_1218+(P_3205-(P_1217+P_1219))); dy(57) = ... (IIf(EvalParameter(P_272, Time, y),(P_2956+((P_2324*P_2837)-(P_2957+(P_2608*y(57))))),0)); dy(58) = ... (IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1872*y(55)*P_2837),0)) ... +(0-(IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1883*y(58)*P_2839),0))) ... +(0-P_3182); dy(59) = ... (P_1221+(P_3206-(P_1220+P_1222))); dy(60) = ... (IIf(EvalParameter(P_272, Time, y),(P_2957+((P_2329*P_2839)-(P_2958+(P_2611*y(60))))),0)); dy(61) = ... (IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1883*y(58)*P_2839),0)) ... +(0-(IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1894*y(61)*P_2841),0))) ... +(0-P_3184); dy(62) = ... (P_1224+(P_3207-(P_1223+P_1225))); dy(63) = ... (IIf(EvalParameter(P_272, Time, y),(P_2958+((P_2334*P_2841)-P_2616)),0)); dy(64) = ... (IIf((EvalParameter(P_272, Time, y) & P_1246),(P_1894*y(61)*P_2841),0)) ... +(0-(IIf(EvalParameter(P_272, Time, y),(P_1905*y(64)),0))) ... +(0-P_3172); dy(65) = ... (P_1167+(P_3201-(P_1166+P_1168))); dy(66) = ... (IIf(EvalParameter(P_272, Time, y),P_2616,0)); dy(67) = ... (IIf(EvalParameter(P_272, Time, y),(P_1905*y(64)),0)); dy(68) = ... (0-(P_1916*(1-P_1113)*P_1918)) ... +(0-(P_1624*P_1699*((P_2558*P_1918)-(P_2557*(P_1921/P_2806))))) ... +(0-(P_1699*P_1140*P_1625*(P_1918-(P_1923/P_2060)))) ... +(P_1916*(1-P_1113)*P_1721); dy(69) = ... (0-(P_1916*P_1113*P_1921)) ... +(P_1624*P_1699*((P_2558*P_1918)-(P_2557*(P_1921/P_2806)))) ... +(P_1916*P_1113*P_1724); dy(70) = ... (0-((P_2476*P_2772*P_1923)-(P_3106*P_2773*P_2339))) ... +(P_1699*P_1140*P_1625*(P_1918-(P_1923/P_2060))); dy(71) = ... ((P_2476*P_2772*P_1923)-(P_3106*P_2773*P_2339)); dy(72) = ... (0-((1-P_495)*P_1925*(1-P_1113)*P_2847)) ... +(0-(P_2775*P_1699*((P_2558*P_2847)-(P_2557*(P_2849/P_2806))))) ... +(0-(P_1699*P_1145*P_2776*(P_2847-(P_2852/P_2077)))) ... +(0-P_2967) ... +(0-P_2969) ... +(0-P_2971) ... +(0-P_2973) ... +(0-P_2976) ... +(P_1925*(1-P_1113)*P_1721); dy(73) = ... (0-((1-P_495)*P_1925*P_1113*P_2849)) ... +(P_2775*P_1699*((P_2558*P_2847)-(P_2557*(P_2849/P_2806)))) ... +(0-P_2968) ... +(0-P_2970) ... +(0-P_2972) ... +(0-P_2974) ... +(0-P_2977) ... +(P_1925*P_1113*P_1724); dy(74) = ... (0-((P_2481*P_2893*P_2852)-(P_3108*P_2894*P_2853))) ... +(P_1699*P_1145*P_2776*(P_2847-(P_2852/P_2077))); dy(75) = ... ((P_2481*P_2893*P_2852)-(P_3108*P_2894*P_2853)); dy(76) = ... P_2967 ... +(0-(P_2511*P_1699*((P_2558*P_2629)-(P_2557*(P_2631/P_2806))))) ... +(0-(P_1699*P_1178*P_2512*(P_2629-(P_2633/P_2125)))) ... +P_3013 ... +(0-P_2923); dy(77) = ... P_2968 ... +(P_2511*P_1699*((P_2558*P_2629)-(P_2557*(P_2631/P_2806)))) ... +(0-P_2924); dy(78) = ... (0-((P_2510*P_2991*P_2633)-(P_3119*P_2992*P_2635))) ... +(P_1699*P_1178*P_2512*(P_2629-(P_2633/P_2125))); dy(79) = ... ((P_2510*P_2991*P_2633)-(P_3119*P_2992*P_2635)) ... +P_3174; dy(80) = ... P_2964 ... +(0-P_2887) ... +P_2969 ... +(0-(P_2514*P_1699*((P_2558*P_2639)-(P_2557*(P_2641/P_2806))))) ... +(0-(P_1699*P_1181*P_2515*(P_2639-(P_2643/P_2126)))); dy(81) = ... (0-P_2888) ... +P_2970 ... +(P_2514*P_1699*((P_2558*P_2639)-(P_2557*(P_2641/P_2806)))); dy(82) = ... (0-((P_2513*P_2993*P_2643)-(P_3120*P_2994*P_2645))) ... +(P_1699*P_1181*P_2515*(P_2639-(P_2643/P_2126))); dy(83) = ... P_3164 ... +((P_2513*P_2993*P_2643)-(P_3120*P_2994*P_2645)); dy(84) = ... P_2965 ... +(0-P_2889) ... +P_2971 ... +(0-(P_2517*P_1699*((P_2558*P_2649)-(P_2557*(P_2651/P_2806))))) ... +(0-(P_1699*P_1184*P_2518*(P_2649-(P_2653/P_2127)))); dy(85) = ... (0-P_2890) ... +P_2972 ... +(P_2517*P_1699*((P_2558*P_2649)-(P_2557*(P_2651/P_2806)))); dy(86) = ... (0-((P_2516*P_2995*P_2653)-(P_3121*P_2996*P_2655))) ... +(P_1699*P_1184*P_2518*(P_2649-(P_2653/P_2127))); dy(87) = ... P_3166 ... +((P_2516*P_2995*P_2653)-(P_3121*P_2996*P_2655)); dy(88) = ... P_2966 ... +(0-P_2891) ... +P_2973 ... +(0-(P_2520*P_1699*((P_2558*P_2659)-(P_2557*(P_2661/P_2806))))) ... +(0-(P_1699*P_1187*P_2521*(P_2659-(P_2663/P_2128)))); dy(89) = ... (0-P_2892) ... +P_2974 ... +(P_2520*P_1699*((P_2558*P_2659)-(P_2557*(P_2661/P_2806)))); dy(90) = ... (0-((P_2519*P_2997*P_2663)-(P_3122*P_2998*P_2665))) ... +(P_1699*P_1187*P_2521*(P_2659-(P_2663/P_2128))); dy(91) = ... P_3168 ... +((P_2519*P_2997*P_2663)-(P_3122*P_2998*P_2665)); dy(92) = ... P_2975 ... +P_2976 ... +(0-P_2895) ... +(0-(P_2523*P_1699*((P_2558*P_2669)-(P_2557*(P_2671/P_2806))))) ... +(0-(P_1699*P_1190*P_2524*(P_2669-(P_2673/P_2129)))); dy(93) = ... P_2977 ... +(0-P_2896) ... +(P_2523*P_1699*((P_2558*P_2669)-(P_2557*(P_2671/P_2806)))); dy(94) = ... (0-((P_2522*P_2999*P_2673)-(P_3123*P_3000*P_2675))) ... +(P_1699*P_1190*P_2524*(P_2669-(P_2673/P_2129))); dy(95) = ... P_3170 ... +((P_2522*P_2999*P_2673)-(P_3123*P_3000*P_2675)); dy(96) = ... (0-(P_2778*P_1699*((P_2558*P_2862)-(P_2557*(P_2864/P_2806))))) ... +(0-(P_1699*P_1150*P_2779*(P_2862-(P_2867/P_2094)))) ... +(0-((1-P_685)*P_1937*(1-P_1113)*P_2862)) ... +(0-P_2978) ... +(0-P_2980) ... +(0-P_2982) ... +(0-P_2984) ... +(0-P_2986) ... +(0-P_2988) ... +(P_1937*(1-P_1113)*P_1721); dy(97) = ... (0-((1-P_685)*P_1937*P_1113*P_2864)) ... +(P_2778*P_1699*((P_2558*P_2862)-(P_2557*(P_2864/P_2806)))) ... +(0-P_2979) ... +(0-P_2981) ... +(0-P_2983) ... +(0-P_2985) ... +(0-P_2987) ... +(0-P_2989) ... +(P_1937*P_1113*P_1724); dy(98) = ... (0-((P_2486*P_2897*P_2867)-(P_3110*P_2898*P_2868))) ... +(P_1699*P_1150*P_2779*(P_2862-(P_2867/P_2094))); dy(99) = ... ((P_2486*P_2897*P_2867)-(P_3110*P_2898*P_2868)); dy(100) = ... P_2978 ... +(0-(P_2526*P_1699*((P_2558*P_2685)-(P_2557*(P_2687/P_2806))))) ... +(0-(P_1699*P_1193*P_2527*(P_2685-(P_2689/P_2130)))) ... +(0-P_2925); dy(101) = ... P_2979 ... +(P_2526*P_1699*((P_2558*P_2685)-(P_2557*(P_2687/P_2806)))) ... +(0-P_2926); dy(102) = ... (0-((P_2525*P_3001*P_2689)-(P_3124*P_3002*P_2691))) ... +(P_1699*P_1193*P_2527*(P_2685-(P_2689/P_2130))); dy(103) = ... ((P_2525*P_3001*P_2689)-(P_3124*P_3002*P_2691)) ... +P_3176; dy(104) = ... P_2980 ... +(0-(P_2529*P_1699*((P_2558*P_2695)-(P_2557*(P_2697/P_2806))))) ... +(0-(P_1699*P_1196*P_2530*(P_2695-(P_2699/P_2131)))) ... +(0-P_2927); dy(105) = ... P_2981 ... +(P_2529*P_1699*((P_2558*P_2695)-(P_2557*(P_2697/P_2806)))) ... +(0-P_2928); dy(106) = ... (0-((P_2528*P_3003*P_2699)-(P_3125*P_3004*P_2701))) ... +(P_1699*P_1196*P_2530*(P_2695-(P_2699/P_2131))); dy(107) = ... ((P_2528*P_3003*P_2699)-(P_3125*P_3004*P_2701)) ... +P_3178; dy(108) = ... P_2982 ... +(0-(P_2532*P_1699*((P_2558*P_2705)-(P_2557*(P_2707/P_2806))))) ... +(0-(P_1699*P_1199*P_2533*(P_2705-(P_2709/P_2132)))) ... +(0-P_2929); dy(109) = ... P_2983 ... +(P_2532*P_1699*((P_2558*P_2705)-(P_2557*(P_2707/P_2806)))) ... +(0-P_2930); dy(110) = ... (0-((P_2531*P_3005*P_2709)-(P_3126*P_3006*P_2711))) ... +(P_1699*P_1199*P_2533*(P_2705-(P_2709/P_2132))); dy(111) = ... ((P_2531*P_3005*P_2709)-(P_3126*P_3006*P_2711)) ... +P_3180; dy(112) = ... P_2984 ... +(0-(P_2535*P_1699*((P_2558*P_2715)-(P_2557*(P_2717/P_2806))))) ... +(0-(P_1699*P_1202*P_2536*(P_2715-(P_2719/P_2133)))) ... +(0-P_2931); dy(113) = ... P_2985 ... +(P_2535*P_1699*((P_2558*P_2715)-(P_2557*(P_2717/P_2806)))) ... +(0-P_2932); dy(114) = ... (0-((P_2534*P_3007*P_2719)-(P_3127*P_3008*P_2721))) ... +(P_1699*P_1202*P_2536*(P_2715-(P_2719/P_2133))); dy(115) = ... ((P_2534*P_3007*P_2719)-(P_3127*P_3008*P_2721)) ... +P_3182; dy(116) = ... P_2986 ... +(0-(P_2538*P_1699*((P_2558*P_2725)-(P_2557*(P_2727/P_2806))))) ... +(0-(P_1699*P_1205*P_2539*(P_2725-(P_2729/P_2134)))) ... +(0-P_2933); dy(117) = ... P_2987 ... +(P_2538*P_1699*((P_2558*P_2725)-(P_2557*(P_2727/P_2806)))) ... +(0-P_2934); dy(118) = ... (0-((P_2537*P_3009*P_2729)-(P_3128*P_3010*P_2731))) ... +(P_1699*P_1205*P_2539*(P_2725-(P_2729/P_2134))); dy(119) = ... ((P_2537*P_3009*P_2729)-(P_3128*P_3010*P_2731)) ... +P_3184; dy(120) = ... P_2988 ... +(0-P_2899) ... +(0-(P_2541*P_1699*((P_2558*P_2735)-(P_2557*(P_2737/P_2806))))) ... +(0-(P_1699*P_1208*P_2542*(P_2735-(P_2739/P_2135)))); dy(121) = ... P_2989 ... +(0-P_2900) ... +(P_2541*P_1699*((P_2558*P_2735)-(P_2557*(P_2737/P_2806)))); dy(122) = ... (0-((P_2540*P_3011*P_2739)-(P_3129*P_3012*P_2741))) ... +(P_1699*P_1208*P_2542*(P_2735-(P_2739/P_2135))); dy(123) = ... P_3172 ... +((P_2540*P_3011*P_2739)-(P_3129*P_3012*P_2741)); dy(124) = ... (0-(P_2096*P_1699*((P_2558*P_2425)-(P_2557*(P_2427/P_2806))))) ... +(0-(P_1699*P_1153*P_2098*(P_2425-(P_2429/P_2097)))) ... +(0-(IIf(P_896,((P_2423+P_2446)*(1-P_1113)*P_2425),0))) ... +(P_2423*(1-P_1113)*P_1721) ... +(P_2446*(1-P_1113)*P_1998) ... +(0-(IIf(P_896,0,((P_2423+P_2446)*(1-P_1113)*P_2425)))); dy(125) = ... (P_2096*P_1699*((P_2558*P_2425)-(P_2557*(P_2427/P_2806)))) ... +(0-(IIf(P_896,((P_2423+P_2446)*P_1113*P_2427),0))) ... +(P_2423*P_1113*P_1724) ... +(P_2446*P_1113*P_2001) ... +(0-(IIf(P_896,0,((P_2423+P_2446)*P_1113*P_2427)))); dy(126) = ... (0-((P_2489*P_2780*P_2429)-(P_3111*P_2781*P_2745))) ... +(P_1699*P_1153*P_2098*(P_2425-(P_2429/P_2097))); dy(127) = ... ((P_2489*P_2780*P_2429)-(P_3111*P_2781*P_2745)); dy(128) = ... (0-(P_2100*P_1699*((P_2558*P_2434)-(P_2557*(P_2436/P_2806))))) ... +(0-(P_1699*P_1156*P_2102*(P_2434-(P_2438/P_2101)))) ... +(IIf(P_896,((P_2423+P_2446)*(1-P_1113)*P_2425),0)) ... +(0-(IIf(P_896,((P_2432+P_2446)*(1-P_1113)*P_2434),0))); dy(129) = ... (P_2100*P_1699*((P_2558*P_2434)-(P_2557*(P_2436/P_2806)))) ... +(IIf(P_896,((P_2423+P_2446)*P_1113*P_2427),0)) ... +(0-(IIf(P_896,((P_2432+P_2446)*P_1113*P_2436),0))); dy(130) = ... (0-((P_2492*P_2782*P_2438)-(P_3112*P_2783*P_2748))) ... +(P_1699*P_1156*P_2102*(P_2434-(P_2438/P_2101))); dy(131) = ... ((P_2492*P_2782*P_2438)-(P_3112*P_2783*P_2748)); dy(132) = ... (P_2440*(1-P_1113)*P_1716) ... +(0-(P_2440*(1-P_1113)*P_1970)) ... +(0-(P_1642*P_1699*((P_2558*P_1970)-(P_2557*(P_1973/P_2806))))) ... +(0-(P_1699*P_1159*P_1643*(P_1970-(P_1975/P_2103)))); dy(133) = ... (P_2440*P_1113*P_1719) ... +(0-(P_2440*P_1113*P_1973)) ... +(P_1642*P_1699*((P_2558*P_1970)-(P_2557*(P_1973/P_2806)))); dy(134) = ... (0-((P_2495*P_2784*P_1975)-(P_3113*P_2785*P_2441))) ... +(P_1699*P_1159*P_1643*(P_1970-(P_1975/P_2103))); dy(135) = ... ((P_2495*P_2784*P_1975)-(P_3113*P_2785*P_2441)); dy(136) = ... (0-(P_1645*P_1699*((P_2558*P_1980)-(P_2557*(P_1983/P_2806))))) ... +(0-(P_1699*P_1162*P_1646*(P_1980-(P_1985/P_2104)))) ... +(0-(P_1978*(1-P_1113)*P_1980)) ... +(P_1978*(1-P_1113)*P_1721); dy(137) = ... (0-(P_1978*P_1113*P_1983)) ... +(P_1645*P_1699*((P_2558*P_1980)-(P_2557*(P_1983/P_2806)))) ... +(P_1978*P_1113*P_1724); dy(138) = ... (0-((P_2498*P_2786*P_1985)-(P_3114*P_2787*P_2443))) ... +(P_1699*P_1162*P_1646*(P_1980-(P_1985/P_2104))); dy(139) = ... ((P_2498*P_2786*P_1985)-(P_3114*P_2787*P_2443)); dy(140) = ... (0-(P_1648*P_1699*((P_2558*P_1990)-(P_2557*(P_1993/P_2806))))) ... +(0-(P_1699*P_1165*P_1649*(P_1990-(P_1995/P_2105)))) ... +(0-(P_1988*(1-P_1113)*P_1990)) ... +(P_1988*(1-P_1113)*P_1721); dy(141) = ... (0-(P_1988*P_1113*P_1993)) ... +(P_1648*P_1699*((P_2558*P_1990)-(P_2557*(P_1993/P_2806)))) ... +(P_1988*P_1113*P_1724); dy(142) = ... (0-((P_2501*P_2788*P_1995)-(P_3115*P_2789*P_2445))) ... +(P_1699*P_1165*P_1649*(P_1990-(P_1995/P_2105))); dy(143) = ... ((P_2501*P_2788*P_1995)-(P_3115*P_2789*P_2445)); dy(144) = ... P_2887 ... +P_2889 ... +(P_1916*(1-P_1113)*P_1918) ... +P_2891 ... +((1-P_495)*P_1925*(1-P_1113)*P_2847) ... +P_2895 ... +((1-P_685)*P_1937*(1-P_1113)*P_2862) ... +(P_1988*(1-P_1113)*P_1990) ... +P_2899 ... +(0-(P_1656*P_1699*((P_2558*P_1998)-(P_2557*(P_2001/P_2806))))) ... +(P_2013*(1-P_1113)*P_2015) ... +P_2923 ... +P_2925 ... +P_2927 ... +P_2929 ... +P_2931 ... +P_2933 ... +(0-(P_2446*(1-P_1113)*P_1998)); dy(145) = ... P_2888 ... +P_2890 ... +(P_1916*P_1113*P_1921) ... +P_2892 ... +((1-P_495)*P_1925*P_1113*P_2849) ... +P_2896 ... +((1-P_685)*P_1937*P_1113*P_2864) ... +(P_1988*P_1113*P_1993) ... +P_2900 ... +(P_1656*P_1699*((P_2558*P_1998)-(P_2557*(P_2001/P_2806)))) ... +(P_2013*P_1113*P_2018) ... +P_2924 ... +P_2926 ... +P_2928 ... +P_2930 ... +P_2932 ... +P_2934 ... +(0-(P_2446*P_1113*P_2001)); dy(146) = ... (0-(P_1658*P_1699*((P_2558*P_2005)-(P_2557*(P_2008/P_2806))))) ... +(0-(P_1699*P_1171*P_1659*(P_2005-(P_2010/P_2122)))) ... +(0-(P_2003*(1-P_1113)*P_2005)) ... +(P_2003*(1-P_1113)*P_1721); dy(147) = ... (0-(P_2003*P_1113*P_2008)) ... +(P_1658*P_1699*((P_2558*P_2005)-(P_2557*(P_2008/P_2806)))) ... +(P_2003*P_1113*P_1724); dy(148) = ... (0-((P_2506*P_2790*P_2010)-(P_3117*P_2791*P_2448))) ... +(P_1699*P_1171*P_1659*(P_2005-(P_2010/P_2122))); dy(149) = ... ((P_2506*P_2790*P_2010)-(P_3117*P_2791*P_2448)); dy(150) = ... (0-(P_2013*(1-P_1113)*P_2015)) ... +(0-(P_1661*P_1699*((P_2558*P_2015)-(P_2557*(P_2018/P_2806))))) ... +(0-(P_1699*P_1174*P_1662*(P_2015-(P_2020/P_2123)))) ... +(P_2013*(1-P_1113)*P_1721); dy(151) = ... (0-(P_2013*P_1113*P_2018)) ... +(P_1661*P_1699*((P_2558*P_2015)-(P_2557*(P_2018/P_2806)))) ... +(P_2013*P_1113*P_1724); dy(152) = ... (0-((P_2509*P_2792*P_2020)-(P_3118*P_2793*P_2450))) ... +(P_1699*P_1174*P_1662*(P_2015-(P_2020/P_2123))); dy(153) = ... ((P_2509*P_2792*P_2020)-(P_3118*P_2793*P_2450)); dy(154) = ... (P_1699*(P_1721-(P_3033/P_2124))*P_1111); dy(155) = ... (0-(P_1699*(P_1721-(P_3033/P_2124))*P_1111)); dy = dy'; ```
Events ``` function [yOut switchUpdate] = PerformSwitches(Time, y) switchUpdate = false; yOut = y; global SwitchUpdateTimePoints; global P_29; global P_46; global P_51; global P_86; global P_123; global P_124; global P_157; global P_192; global P_226; global P_227; global P_270; global P_271; global P_272; global P_489; global P_677; global P_894; global P_895; global P_933; global P_934; global P_966; global P_1000; global P_1035; global P_1036; global P_1040; global P_1074; global P_1108; global P_1109; global P_1253; global P_1272; global P_3021; global P_2946; switchConditionApplies = (Time == P_1272); if switchConditionApplies && isempty(intersect(SwitchUpdateTimePoints, Time)) newFormula = @(Time,y) (EvalParameter(P_1253, Time, y)+P_3021); if isnumeric(P_1253) P_1253 = newFormula; switchUpdate = true; else if ~strcmp(func2str(P_1253), func2str(newFormula)) P_1253 = newFormula; switchUpdate = true; end end newFormula = @(Time,y) 1; if isnumeric(P_272) P_272 = newFormula; switchUpdate = true; else if ~strcmp(func2str(P_272), func2str(newFormula)) P_272 = newFormula; switchUpdate = true; end end newValue = (y(31)+P_2946); if newValue ~= y(31) yOut(31) = newValue; switchUpdate = true; end newValue = (y(32)+P_3021); if newValue ~= y(32) yOut(32) = newValue; switchUpdate = true; end end if switchUpdate SwitchUpdateTimePoints = [SwitchUpdateTimePoints Time]; end ```
Table Parameters ``` function SetupTableParameters global P_29; P_29 = @(Time,y) P_29_Table(Time,y); global P_46; P_46 = @(Time,y) P_46_Table(Time,y); global P_51; P_51 = @(Time,y) P_51_Table(Time,y); global P_86; P_86 = @(Time,y) P_86_Table(Time,y); global P_123; P_123 = @(Time,y) P_123_Table(Time,y); global P_124; P_124 = @(Time,y) P_124_Table(Time,y); global P_157; P_157 = @(Time,y) P_157_Table(Time,y); global P_192; P_192 = @(Time,y) P_192_Table(Time,y); global P_226; P_226 = @(Time,y) P_226_Table(Time,y); global P_227; P_227 = @(Time,y) P_227_Table(Time,y); global P_270; P_270 = @(Time,y) P_270_Table(Time,y); global P_271; P_271 = @(Time,y) P_271_Table(Time,y); global P_489; P_489 = @(Time,y) P_489_Table(Time,y); global P_677; P_677 = @(Time,y) P_677_Table(Time,y); global P_894; P_894 = @(Time,y) P_894_Table(Time,y); global P_895; P_895 = @(Time,y) P_895_Table(Time,y); global P_933; P_933 = @(Time,y) P_933_Table(Time,y); global P_934; P_934 = @(Time,y) P_934_Table(Time,y); global P_966; P_966 = @(Time,y) P_966_Table(Time,y); global P_1000; P_1000 = @(Time,y) P_1000_Table(Time,y); global P_1035; P_1035 = @(Time,y) P_1035_Table(Time,y); global P_1036; P_1036 = @(Time,y) P_1036_Table(Time,y); global P_1040; P_1040 = @(Time,y) P_1040_Table(Time,y); global P_1074; P_1074 = @(Time,y) P_1074_Table(Time,y); global P_1108; P_1108 = @(Time,y) P_1108_Table(Time,y); global P_1109; P_1109 = @(Time,y) P_1109_Table(Time,y); function yout = P_29_Table(Time, y) if Time <= 0 yout = 1; return; end if Time >= 31395600 yout = 0.986596986265264; return; end X_Values = [ 0 5.60575585950573e-09 31395600 ] Y_Values = [ 1 1 0.986596986265264 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_46_Table(Time, y) if Time <= 0 yout = 0.964017053589008; return; end if Time >= 36817200 yout = 0.73571357788497; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 0.964017053589008 0.959787566895977 0.945182619104018 0.922028360566593 0.892912374931973 0.8638237273463319 0.835927847643979 0.73571357788497 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_51_Table(Time, y) if Time <= 0 yout = 0.419106180714442; return; end if Time >= 36817200 yout = 0.320330995758924; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 0.419106180714442 0.418026557875077 0.411658846073675 0.401563661384098 0.388869182227036 0.376186622374279 0.364024105607944 0.320330995758924 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_86_Table(Time, y) if Time <= 0 yout = 11.8188930262479; return; end if Time >= 36817200 yout = 8.301639682783989; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 11.8188930262479 11.8490318812716 11.4826371518858 11.0925844718037 10.3559302892674 9.637984748548901 8.96981221566644 8.301639682783989 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_123_Table(Time, y) if Time <= 0 yout = 1.50913809221829; return; end if Time >= 36817200 yout = 1.3445199164005; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 1.50913809221829 1.50731333183525 1.50312848636575 1.46358674340526 1.43462527585336 1.40035715341436 1.37243746537672 1.3445199164005 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_124_Table(Time, y) if Time <= 0 yout = 0.51675; return; end if Time >= 36817200 yout = 0.406727314773112; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 0.51675 0.4960804488543 0.484331152796907 0.473498747384986 0.453596455943061 0.445553419601616 0.4331690347028 0.406727314773112 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_157_Table(Time, y) if Time <= 0 yout = 14.6567054773171; return; end if Time >= 36817200 yout = 16.8394985174644; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 14.6567054773171 16.0779516653242 16.4394854328631 17.6660429431114 19.1710006654368 19.6642652815604 20.3969475767517 16.8394985174644 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_192_Table(Time, y) if Time <= 0 yout = 0.0402694788859839; return; end if Time >= 36817200 yout = 0.0294584103768549; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 0.0402694788859839 0.0407581168035073 0.0350254245116991 0.0335367986988893 0.0318958021534651 0.0305862728081995 0.0303467427702173 0.0294584103768549 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_226_Table(Time, y) if Time <= 0 yout = 0.417642184690095; return; end if Time >= 36817200 yout = 0.425536369520118; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 0.417642184690095 0.434893559202714 0.439696885224335 0.455154645600953 0.465348529409829 0.441845878673725 0.438368059642297 0.425536369520118 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_227_Table(Time, y) if Time <= 0 yout = 0.62335; return; end if Time >= 36817200 yout = 0.758275609546323; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 0.62335 0.6242228869999999 0.633113186 0.658058594688688 0.664250454152935 0.758275609546323 0.758275609546323 0.758275609546323 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_270_Table(Time, y) if Time <= 0 yout = 0.438501357013879; return; end if Time >= 36817200 yout = 0.356578867622123; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 0.438501357013879 0.476784485647609 0.468965997211384 0.456371150969536 0.44247447331563 0.396414759459036 0.367331673182402 0.356578867622123 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_271_Table(Time, y) if Time <= 0 yout = 3.02705; return; end if Time >= 36817200 yout = 0.846482685136696; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 3.02705 2.81404023345742 2.38459219750622 2.00479432050997 1.65343691084705 1.43363143935299 1.15714878467975 0.846482685136696 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_489_Table(Time, y) if Time <= 0 yout = 0.3861; return; end if Time >= 36817200 yout = 0.12926396886118; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 0.3861 0.388271282920694 0.386967740037282 0.330358733226773 0.254208230630086 0.190893176829743 0.160353559331436 0.12926396886118 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_677_Table(Time, y) if Time <= 0 yout = 0.8976499999999999; return; end if Time >= 36817200 yout = 0.300527846796784; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 0.8976499999999999 0.90269805002269 0.899667422544589 0.7680562467780701 0.591012738215738 0.443810567679924 0.372808527671234 0.300527846796784 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_894_Table(Time, y) if Time <= 0 yout = 0.412581159732335; return; end if Time >= 36817200 yout = 0.412464901777145; return; end X_Values = [ 0 36817200 ] Y_Values = [ 0.412581159732335 0.412464901777145 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_895_Table(Time, y) if Time <= 0 yout = 0.6304999999999999; return; end if Time >= 36817200 yout = 0.211087625918089; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 0.6304999999999999 0.6340456976987759 0.631917016559197 0.5394746990403529 0.415121184698961 0.311727915025001 0.26185682247726 0.211087625918089 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_933_Table(Time, y) if Time <= 0 yout = 2.37726035978654; return; end if Time >= 36817200 yout = 1.21556700014631; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 2.37726035978654 2.34563732370137 2.2263141786317 2.05898412852037 1.73079486771463 1.43093626903022 1.33792181280061 1.21556700014631 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_934_Table(Time, y) if Time <= 0 yout = 0.1794; return; end if Time >= 36817200 yout = 0.117814239620869; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 0.1794 0.182996245052804 0.192180754561188 0.177307451573124 0.16243414858506 0.147560845596997 0.132687542608933 0.117814239620869 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_966_Table(Time, y) if Time <= 0 yout = 1.21442383944443; return; end if Time >= 36817200 yout = 0.917569006530751; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 1.21442383944443 1.25225965825849 1.25484693437177 1.26567574930572 1.19329111133064 1.06772488177619 0.945238736880048 0.917569006530751 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_1000_Table(Time, y) if Time <= 0 yout = 32.6404056170836; return; end if Time >= 36817200 yout = 19.716229778381; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 32.6404056170836 32.6198123190286 31.4715808157702 29.8762536053942 28.82867098645 26.7723074443089 23.8777207725935 19.716229778381 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_1035_Table(Time, y) if Time <= 0 yout = 0.190430290456271; return; end if Time >= 36817200 yout = 0.139612113607712; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 0.190430290456271 0.190996559839795 0.183769059499443 0.178318424484128 0.165554752063014 0.15004683392541 0.143822183378515 0.139612113607712 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_1036_Table(Time, y) if Time <= 0 yout = 0.3419; return; end if Time >= 36817200 yout = 0.156131205584736; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 0.3419 0.342803347940327 0.355089391282635 0.312409982697301 0.258930545105473 0.214535367462637 0.188013076428367 0.156131205584736 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_1040_Table(Time, y) if Time <= 0 yout = 1.03743257442277; return; end if Time >= 36817200 yout = 0.683742388704714; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 1.03743257442277 0.8919876747757201 0.878414428079786 0.856895798388897 0.829836581098742 0.802802770683589 0.776877470409073 0.683742388704714 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_1074_Table(Time, y) if Time <= 0 yout = 3.76319586580084; return; end if Time >= 36817200 yout = 3.12452723880335; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 3.76319586580084 3.79247892178574 3.74737632131139 3.69422174914786 3.63737127536962 3.53914497306394 3.44175171110928 3.12452723880335 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_1108_Table(Time, y) if Time <= 0 yout = 0.206933891921373; return; end if Time >= 36817200 yout = 0.0921722919939123; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 0.206933891921373 0.188664456442222 0.177006478847315 0.167475634828187 0.158658237757135 0.136359411816881 0.113442012564629 0.0921722919939123 ] yout = interp1(X_Values, Y_Values, Time); function yout = P_1109_Table(Time, y) if Time <= 0 yout = 0.801125; return; end if Time >= 36817200 yout = 0.602156263312432; return; end X_Values = [ 0 5259600 10519200 15778800 21038400 26298000 31557600 36817200 ] Y_Values = [ 0.801125 0.883643914390751 0.938680307865232 0.846966568639193 0.687953509419963 0.6010870826449241 0.606927464410815 0.602156263312432 ] yout = interp1(X_Values, Y_Values, Time); ```
ODE Options ``` function opt = ODEoptions global eventsFunctionHandle; opt = odeset('AbsTol',1e-10, ... 'BDF','on', ... 'InitialStep',1e-10, ... 'MaxOrder',5, ... 'MaxStep',100000, ... 'RelTol',1e-05, ... 'Events',eventsFunctionHandle ... ); ```
Main function/script ``` function [tout, yout] = ODEMain tic clear global; global eventsFunctionHandle; eventsFunctionHandle = @events; SetupTableParameters; y0 = ODEInitialValues; tstart = 0; tout = tstart; %call ODE RHS function in order to init all parameters ODERHSFunction(tstart, y0); %perform initial switches update [y0 switchUpdate] = PerformSwitches(tstart, y0); yout = y0.'; outtimes = [ 0 3 6 9 12 15 18 21 24 27 30 33 36 39 42 45 48 51 54 57 60 63 66 69 72 75 78 81 84 87 90 93 96 99 102 105 108 111 114 117 120 135 150 165 180 195 210 225 240 255 270 285 300 315 330 345 360 375 390 405 420 435 450 465 480 495 510 525 540 555 570 585 600 615 630 645 660 675 690 705 720 735 750 765 780 795 810 825 840 855 870 885 900 915 930 945 960 975 990 1005 1020 1035 1050 1065 1080 1095 1110 1125 1140 1155 1170 1185 1200 1215 1230 1245 1260 1275 1290 1305 1320 1335 1350 1365 1380 1395 1410 1425 1440 ]; global switchtimes; switchtimes = [ 0 ]; NextRestartTime = tstart; while NextRestartTime < outtimes(end) switchtimes = switchtimes(switchtimes > NextRestartTime); if isempty(switchtimes) NextRestartTime = outtimes(end); else NextRestartTime = switchtimes(1); end while (1) outtimes_step_i = [tstart outtimes(outtimes > tstart & outtimes <= NextRestartTime)]; if (outtimes_step_i(end) < NextRestartTime) outtimes_step_i = [outtimes_step_i NextRestartTime]; %#ok end [t, y] = ode15s(@ODERHSFunction, outtimes_step_i, y0, ODEoptions); nt = length(t); tout = [tout; t(2:nt)]; %#ok yout = [yout; y(2:nt,:)]; %#ok y0 = y(nt,:); tstart = t(nt); %update start vector by switches (if applies) [y0 switchUpdate] = PerformSwitches(tstart, y0); if (tstart==outtimes_step_i(end)) break end end end [tout idx_out] = intersect(tout, outtimes); yout = yout(idx_out, :); toc function [value,isterminal,direction] = events(Time,y) value = 1; %no event per default isterminal = 1; % stop the integration direction = 0; % negative direction global switchtimes; if (min(abs(switchtimes - Time)) == 0.0) return end [yOut switchUpdate] = PerformSwitches(Time,y); if ~switchUpdate %no switch fired return end value = 0; % force system restart ```
Observers Are not part of the current Matlab ODE Export, but should be exported as well
Jacobian of the RHS Is not part of the current Matlab ODE Export, but should be exported as well
msevestre commented 4 years ago

Isn't it a simmodel job?

Yuri05 commented 4 years ago

Isn't it a simmodel job?

Sure. But it should know how to export :) This issue is about the definition of an export format. Once the format is clear, export can be implemented relatively easily.

PavelBal commented 4 years ago

How is the ODE sysmtem to be solved? With the sundialr package?

Yuri05 commented 4 years ago

Whatever the best choice for stiff ODE solver in R is 😄

msevestre commented 4 years ago

But it should know how to export :)

Ok sorry. Got it. Well I can't help here

PavelBal commented 4 years ago

I am not an expert here, maybe we should ask someone who directly solved ODE in R?

msevestre commented 4 years ago

@abdullahhamadeh @pchelle Any idea here? Stiff ODE Solver in R? Or do you know anyone who knows?

PavelBal commented 4 years ago

@StephanSchaller RxODE? However, I don't like RxODE so would prefer not to...

StephanSchaller commented 4 years ago

Well, we already have discussions with the RxODE crew (via nlmixr)… maybe we can start with one format (RxODE) and then expand. I'll also ask an R-Guru I know.