From a9ad3ea70ee0d1b8185f6d7358611e73c445cff0 Mon Sep 17 00:00:00 2001 From: JuliusMartensen Date: Thu, 1 Feb 2024 14:56:04 +0100 Subject: [PATCH] Initial draft for JOSS publication Co-authored-by: chplate christoph.plate@ovgu.de Co-authored-by: SebastianSager sager@ovgu.de --- .gitignore | 4 ++ paper/bibliography.bib | 101 +++++++++++++++++++++++++++ paper/figures/Lotka.pdf | Bin 0 -> 35520 bytes paper/paper.md | 148 ++++++++++++++++++++++++++++++++++++++++ 4 files changed, 253 insertions(+) create mode 100644 paper/bibliography.bib create mode 100644 paper/figures/Lotka.pdf create mode 100644 paper/paper.md diff --git a/.gitignore b/.gitignore index 6232a26..9b9a713 100644 --- a/.gitignore +++ b/.gitignore @@ -304,3 +304,7 @@ src/__old/** test/__old/** .DS_Store Manifest.toml + +paper/media/* +paper/*.jats +paper/*.pdf \ No newline at end of file diff --git a/paper/bibliography.bib b/paper/bibliography.bib new file mode 100644 index 0000000..2cc126f --- /dev/null +++ b/paper/bibliography.bib @@ -0,0 +1,101 @@ +@misc{ma2021modelingtoolkit, + title={ModelingToolkit: A Composable Graph Transformation System For Equation-Based Modeling}, + author={Yingbo Ma and Shashi Gowda and Ranjan Anantharaman and Chris Laughman and Viral Shah and Chris Rackauckas}, + year={2021}, + eprint={2103.05244}, + archivePrefix={arXiv}, + primaryClass={cs.MS} +} + +@article{bezanson2017julia, + title={Julia: A fresh approach to numerical computing}, + author={Bezanson, Jeff and Edelman, Alan and Karpinski, Stefan and Shah, Viral B}, + journal={SIAM review}, + volume={59}, + number={1}, + pages={65--98}, + year={2017}, + publisher={SIAM}, + url={https://doi.org/10.1137/141000671} +} + +@software{vaibhav_kumar_dixit_2023_7738525, + author = {Vaibhav Kumar Dixit and Christopher Rackauckas}, + month = mar, + publisher = {Zenodo}, + title = {Optimization.jl: A Unified Optimization Package}, + version = {v3.12.1}, + doi = {10.5281/zenodo.7738525}, + url = {https://doi.org/10.5281/zenodo.7738525}, + year = 2023} + +@Article{Waechter2006, + Title = {{O}n the {I}mplementation of an {I}nterior-{P}oint {F}ilter {L}ine-{S}earch {A}lgorithm for {L}arge-{S}cale {N}onlinear {P}rogramming}, + Author = {A. W\"{a}chter and L.T. Biegler}, + Journal = {{M}athematical {P}rogramming}, + Year = {2006}, + Number = {1}, + Pages = {25--57}, + Volume = {106}, + doi = {10.1007/s10107-004-0559-y}, + File = {:2006/Waechter2006.pdf:PDF}, + Journal-iso = {Math. Program.}, + Keywords = {interior point}, +} + +@Article{Sager2013, +author = {Sager, Sebastian}, +title = {Sampling Decisions in Optimum Experimental Design in the Light of Pontryagin's Maximum Principle}, +journal = {SIAM Journal on Control and Optimization}, +volume = {51}, +number = {4}, +pages = {3181-3207}, +year = {2013}, +doi = {10.1137/110835098}, +URL = {https://doi.org/10.1137/110835098}, +eprint = {https://doi.org/10.1137/110835098}, +abstract = { Optimum experimental design (OED) problems are optimization problems in which an experimental setting and decisions on when to measure---the so- called sampling design---are to be determined such that a follow-up parameter estimation yields accurate results for model parameters. In this paper we use the interpretation of OED as optimal control problems with a very particular structure for the analysis of optimal sampling decisions. We introduce the information gain function, motivated by an analysis of necessary conditions of optimality. We highlight differences between problem formulations and propose to use a linear penalization of sampling decisions to overcome the intrinsic ill-conditioning of OED. The results of this paper are independent from the actual numerical method to compute the solution to the OED problem and of the question of local and global optima. } +} + + +@article{DifferentialEquations.jl-2017, + author = {Rackauckas, Christopher and Nie, Qing}, + doi = {10.5334/jors.151}, + journal = {The Journal of Open Research Software}, + keywords = {Applied Mathematics}, + note = {Exported from https://app.dimensions.ai on 2019/05/05}, + number = {1}, + pages = {}, + title = {DifferentialEquations.jl – A Performant and Feature-Rich Ecosystem for Solving Differential Equations in Julia}, + url = {https://app.dimensions.ai/details/publication/pub.1085583166 and http://openresearchsoftware.metajnl.com/articles/10.5334/jors.151/galley/245/download/}, + volume = {5}, + year = {2017} +} + +@PhdThesis{Koerkel2002, + Title = {{N}umerische {M}ethoden f\"ur {O}ptimale {V}ersuchsplanungsprobleme bei nichtlinearen {DAE}-{M}odellen}, + Author = {K{\"o}rkel, S.}, + School = {Universit\"at {H}eidelberg}, + Year = {2002}, + + Address = {Heidelberg}, + Keywords = {agbock experimental design DAE}, +} + +@article{Li2000SensitivityAnalysisDifferential, + title = {Sensitivity Analysis of Differential\textendash Algebraic Equations: {{A}} Comparison of Methods on a Special Problem}, + shorttitle = {Sensitivity Analysis of Differential\textendash Algebraic Equations}, + author = {Li, Shengtai and Petzold, Linda and Zhu, Wenjie}, + year = {2000}, + month = feb, + journal = {Applied Numerical Mathematics}, + volume = {32}, + number = {2}, + pages = {161--174}, + issn = {0168-9274}, + doi = {10.1016/S0168-9274(99)00020-3}, + urldate = {2023-06-14}, + abstract = {We compare several methods for sensitivity analysis of differential\textendash algebraic equations (DAEs). Computational complexity, efficiency and numerical conditioning issues are discussed. Numerical results for a chemical kinetics problem arising in model reduction are presented.}, + langid = {english}, + keywords = {Automatic differentiation,Differential\textendash algebraic equations,Sensitivity analysis}, +} \ No newline at end of file diff --git a/paper/figures/Lotka.pdf b/paper/figures/Lotka.pdf new file mode 100644 index 0000000000000000000000000000000000000000..83511c45c40d70f69fe96de91f2d22d476868eb9 GIT binary patch literal 35520 zcmc$_Wm_CU*Dag`8GLXJ6Wn2N9ozy0cXxLmECdD%PH=(~AV~1w5G=R_cXxLP`cCfq zdEV=sUvR$AGTptqs+y{`*IJuaMM{R1gN+M~wrBfd`>^;jcf6+;jf;Yv!r9CY4FaJ6 z%h)=2Sh~Tt4yGQKQkLe<7M2v?50*~W9ySy_926oVXzm_vmZpwq(AUfwC>jhJjEcjm41-7TP;kP3R-*u`S-LxWx|v(LQ*i(Hgru{R2mE(;ihpkozd+H_ z!q!y68A_qg4nM)m&Q8J2&TfbXf7$nNMqS(}|3niJ`PU2o^WcA{{`bNEZNUHD1XlMn^Y~Ab z{u6)#tZi!n@56sD2Yzp9Yi;B4?zLDzO*mEHgy!~r0~{=t?l+B zW5exA{Y9Jq`8&x4t3)~87Qc%n03E&c=y2=6(bx3hO!s|U@{m`QVkp>e{@-q1i&1cK z)heX2y=F*@q(MY^E%6DP0@2}9M<_bcN!zVa41sIQa%gBMX|i05Y;JCD$Kmxwang9O zn4H|_$^EH(IYG!9P+Yh#LtNMm;d|Za6?a7iP?xnoQiWd4tmwnz1FF+p6x9 z&Ts>wi#MV}2ZJ95ma(&H{p~uluVbT%A-|Cv5@}JYI zh={B?>iu1t{dBFWZ~I0nQV9i4<2+iMd;V_gX_p;W{E;=h9rkK%8Wy;W{tE7@3E^U; z^ghnq68}YWW0Tfq;Wk>DPBQsG4}9inmuw3W{WY2)I&*l=UV?F+i$&a>4(Rky-uiIT0w2Yn7@Df$g1#g( z({+pcOrj$k=OK{!AIfL`wJV-aG_teHcH{;cc(7uU+Hjg{r_45sFVrV>Q%OR6T!Fi# zjWeg)jz$b^2bCevapRucciE6e>fbOT>8h0}^^+}$^H){<+she8B>w$RgJ#yTlK#mE zy0Y+|7COm{j!s8Nynj@Ivdd`Io^_;5B{(MD^9?Q5OK1g!Xi??H0rcx^IGmRFnGbEN z%7LPf>=96xHSt+ZytfFzL9SfvNZbaZ%1&pNpJs-ZYO8e8^fzb<+Va-Yajy*tv*rrE zI(K2x`8ehmwqEP+>}fq!r3A$epR`k zts0ba{e}v07&1Kkweal{b=fPgmh~_)-`kTm|G#d-rdkySB;-V41Rq`<-85}K)% zAyE3R4mmimI4RowLS0~m6DQUn)b`zV>wz{p&{^SNQp!+vpNTGTqW7i5L(>r6VX<*X z(`W0Qu>I`59f8s9!7=y2;BNHd9It?y`$$IoW8ow|s;~V<{370xNk_Dvd9&L+Z*4nI zejeC)XhwD%OqZ2u{a{L#tg#sQVmDptntJBvHm#}sySI1HZUU-8NLVky+J<4wHmKkzl~4PwaQRFve3yEubD)=U~g|UDeu=|l^A>?sAFI6{3E*@!)GXN@oe@} zU?67O+MkV8I7wSe-U1xeZz!z4UrH#THDywmPN+Dsj-#^P+0)bGwmlg0PE>@4=ci;i z*6nKGuP^pk8{#14OI#7Z`9I}D#%vcl(K}V;+&|xorlW~>$e$)9IELZ|?0>oJJ8v3a zZS8?QQz1XOtaWnqEs6F;(DQnl@H{uR%%=n%-yP+N9GkCm-(dgMxGafpgPN`M{uP;i zA0#Ul3drMfT&T0JXhNnI$KFU`)xR^yN9_GPYv{ea8Ek%cw(c@Q*#@p(csaQj?0bi0_dDkUd9><9L9O2xMv@AD7- z#bZeogM+m-HEwQh7T^RG)S$;l|B7%XM#L8XyP*6cvMAl3`M--cIv=7;6x+h@S=h|C z%j@zJA^wgL?g(~^kY)njD@hXOWtz~l!X4$TnEG-6X72z4$QmTG6`#Z0}GF=wurcQimEl?$lVO?&kq%mR-4Ee=^~>xs)Pdk zt`BGpL?I{=8nxcL3OV&`*?fi?+f6po*st`r zrqT8uhqNLw0%J*e4myGmZ1`seI%d)9wJBnPjP@rA$FhZ)(NWHGSy&uGseMveH+?=q0+MTWCmp2v*H;Y^=8ny_`yQl4RJ;~5!$>#m;xJA#Nc|hhA@Z!t}xzk z!^hUNx;=s&nwIP}_@D{UwYv%$Q!s3xy)exNjg4WP;|s(X#IX@L88zyq$kX z5>y*#d88iq!qV0|7o1gT({pZwT;><6jIFBwE+;@oD_I1S)-7?c!l$`#fu)<&empcj zRT~(xkv;ku4ws(j9OjqMq|Nt5fv1HoYDpuro=!7kPSq7td;otYXr529{1eB6!f3l|| zUyO}TCDS6*IXB5WRIBgK7j0Og=PADnM$fvxJ{(#Rc{p1e)liNJ{iv)x$SWU(?P|T| zW*kN5arUw7l=-vF1J%zd4su~yX)4M<`Ts^WgP^>-XuGu!h?V}7R>&lD+sRTh1^x>Rs#g4#JFq; zZATbrn0pmREPcbTlhBEfvzH2$j+mG;uH3K-gQ_5=rcA1AbRz8vUsqgxEJAu6Aovz2 zxU_OeXX|U*L#kSiKdqh8Za`9tOmL+40}S52_$+2~wzlxoo`X6v+-GMvU57=(z<|ZB zQ4xB!IXx=9`}BZJX%`Wz4t0Rll3qtwluz}Xtd<_8x$;;~M`vzNTCm~)2fOmwFi1er zW4Hbm3vW~z!`=cVL@pv88XQI1buCzlzy1<_)%gY2FW^+3sM0^aK9oif8|5`91R(~E zfXQL8KF0o;KF5Ee#!>+X-)pT?Rd!i$=G$k|{8RiKYQJN3zW~LYse-d{t$<$49f@(0 zyTw#eN>Gy7hB_`M=n?i>GSK|&P`p=@#=nEUB7V_8NCpi;l*WjL{v`K3>k3 z@wLv;7Q+D6zbyTmj9Z&Jg*!H=cWP~0nNyDDYuzx7gz!Lf1TxQ8T;|3^KQcwd~F z6Z`Sv@!BaaJr^GLlZeM33ZY0iK(s&OX;`m8ETSU`x(?T51rY$85p~|GpCv2i+$(0^ zJ13DI4Q*!W6?0Iqf%hT(xe`|bWAx0PzK)zabqqo&dnS3Y?}=CU)sq=r=vleFFwn$Fot-8CZps7E{2x@gd z(AuEqW~=Q1hciDGzkSxrbqNXSn7QCvffiUROUK3bz^qnvp%bHHv-`hOec-s=MPVXv z2y%+Bph4o$)*c8Aju_!eq- zIzOEkeRWcc6kr&JPn*z#DaKaZvse`~Gc#;DojI9r2sK7YfxnWEE%QnBydr~ zGKqhq0;+*BT^F8?!^GHV)DQAdZbIUfrTe^!pC+Y|Zo~`+jRXa}6)ke*8|0dpAam^j zej)nz8&Tk$(X5Ag>l{j6;g6aYW#nXJQD+}i*5ZN?02?Lo!fx*FlM@pzs0eJTko)0u z*qVKRhN!hf`FB>yFpN$qjt~I>0n5Rp&QaflH&AvaCi!|*I@H4-EZMHD@myHr*tARE zf<4j#ph^P?o@Huu*Ha0z!7R3yFLK5rtXy1pc%+2*MDo%dLu|8>uG-xIz-FmyLgTUZjWzLsfT`5wg?VCfXhUI!X*Gb?&btmy$`G0qK<=b20r^gYowvX_^8&QY>L;BqK9i$rku;zZ?`fqoW>Mt zmve6#fkLVPyoMeQ+7=C?R(&2^6Efwh%B166`>k(esFAyAoX1C^E>_18{9qcsdnFvg zdP;Rv?O&V_YY3kc8z}X>Z%oD*wmaeJN*2HCo8zt?lldCUrmJf6KKWEO%=Sb_=l-)b zI6!8vYCCOzID4jcEJ=pkT?-QR!ny>%UvcJTSlaqu<+SgplpO~>Qf<~&rmSetF(U##9DF|%uZ4z0J1d+$fm*Quf?arPoeIN*pZw>i@ zt&NJs@_9W&w?wTl6#eCC|I_L|P~`5_p#^r#ds{X9CW$b>pOlBwU+)xP7N?DtkO&@U zoxgOYFf}kN$_g`668U&nRZjOhTh426cv@XDKDIN4OkK4A-2nTP0(! za{umozO8GL5Q&@st)>px-@e_oF-(wxoA=;r(?BBqucwAxAzoUV0C zQCzo+cs~|bB5C*%XtEo4jq?nzKX{q{DIcD2>#OJy{KSK*C4sO3{vpEtNVLD{& zc6L2|ohZISb+8OByvkkcG%U#`db2B`c%von$FM`pG0C8D(!4mmxbQuqYl5e!mN1$b z(yU&n>D?A3W*%Ot;77K`QsX8vq4#5&(5R4pXU z1iw1ZB}|&~9NmG6;P6+bFe&(@YxKR5gXeFMP@nwZ8CvLy9-U|w^F-a9Mn5b&Klm>8 znx#UFozOZEbX2>#Y}s&Z)g&!kx3{Mx^*T}H6((O}QI3}Xa`FK{DwsbzVuvI?AJ(7( z2N}0m>$Nz2;=g0RrhV(o^B4ZM%xslsQ0fQbUZBin0>C?d zWurBl@$yrWT-4!&ImzdxZKBJU!~R0;+!w%CsU3QaU)=koBm_0Z#rS}+BwWZ!nZeua zrw?n5<|jKGA}h|)Z6CJ1(!3dc&MM)=n#J*p`_j`S;XB3n5t~w#9vB0&<@=#8=&8$6 z<9Q*UMtsh$X-5@<&r;=oR2i`y<&AEk%RTkHw!dR)tM6%=Pf=X6k3giVr*Mi=k+Ahd zygE1eq4!fhnF$H7-023nTZ37TW~zou^l8TE14G(A$TIfUu5HJ!vQgTPg!vs5LLyFi z?(C@Ht#d<}O9xnIdo`@uw#>3(@tFdDf1@fQOwRrux*On+s5)p<4YB?^m8r?a%Ek?i zn(TvhE*jPIykM{ezkt<>BP~tA>Z)9)xHI}@(}j^w^+E>s)8OH2(D^$?xFHv4*EgYN_ReYfMSis&_j(aF>6p>yi1_Bo1xdVksk1`$7g z&k$fR_5r|9kdqTPA^GF;Uu%9<)sn~a&Ayz#XCHnl;eh+P{{kJ6?HaZz3-|N34>zUr zT}OpDP4-p~EJzAQ9J`z9Bb-RIAQq9IqMjarp;89iAV?Ggnl@E?7CacjfK7yywOziLx`>?@^=j; zR(FAHTyU}7Y{k>GwlQV3{Vdz@T#BFd_WzeX#uZ zJwQC5dij;erRD)$>lHNV@_mtESO^m&La^^=DQmevp#Y*oKg`jq5qa@wvi8b$-Iw7s z*r=4^OYYa)+K7-F2O4kY=+Gma=q5@`Jv7`$JxP<2bXzO{*Bw^G)@K z&*csT*nYPW6+`yY@~6Okcl2GrR=C1j0+9dhviCVA&UNPzLoe)L70{WsC-m6aMoL+5 z^~$)+06r&$N^Kf$PCtw6oviI%IAd|I1Y3UoGFQOZt^c5DFcLb7+rSzYaqg?JHtT_Eu8D`=gRXDmWzh zHa+EO|0aD194ZoU4|zQ6hMM6pg3PPuYhMhzKUIlilJQPC$%g_Cr^~`E^Os(O%2f+` zBMF!-SSMS*v#`i_^}KQ3PX7pg!NLbxeK^#v@`kDDTWsH?|GpT2p>)KsFV*}&(a-w< zb=ks6b>zc{%*MVNhx6Nh!p{L0W$eodcBiM={#@>SQkLMiV>t_y2}wyB6$VY#**--c z`xD+PKJ>_W20yX`Zsvczq}=+`>N_$5&A6lumI&HS14g58frDmcG$G7N{978$r8Gbk zU{^B3+eju&Ik}H)vwp_kR;4N}Mk1EzfnCs(d!P9zI$}oJo{NQWz+_Zg!{Q1vPv9-r zjO^jrz`T79S!cNqHLo(w%kokJv3;5R37nUPBdw7Pp6*?N3Ns?GQP0^ioUhV5AR!TR@=U0* z5)K!ic9_uC{Xh4@NKZxRhAT62FVo>IH9LiPjoUk9sbJFV&D!ah8im=$K=$7Ae}z;b zb-Te9k`$pXy|%XvWLK>m$cT>{lYqpq264cwYV*hwLWEvYt(&jWFm3{SRk2Pgr}xRy z^V5U5Lqvv&?NV)7#y;9Cu;H&ue(bpqNkO3b%S>>^GyZ z21R!ZVyHRyMaSL{*=lumF$#H|*lH_-w$1)PELH0BtwjJ2-xuZzN5e|44si1sJp$Oy zfH?#*moJU}al1Hy_dYolszT0&m?7Kf<4K;Qjd{bFa}$cTH~PUp*7SjkR~(&HKl3IP zpenrBtX7B7+;L%Yck_d$LXr{a) zkzVN_=YlKG4FsX%SZ#~Y?id{u3gF3y>nR=Lej)Cj)e8Y z^rUOzS425>JB%GeBq=@qa0V9@74l5GGY}MVkS$oNugoWHg#}ChYetsN# zp9m9l&V^PjE6sM~W)Sxnj;0Aznaq-&c#)mvxo{y_t9M?X;UDQ`aoDSv9A>7a%~hMr zMx*`Xd}=I*o*r%i-$Z&~@>QaoW3G^>kp}0@whtv=sC*XgM5AI)@KIwjz5XHc)XoSU zS@NrM3N(qeKhd^No=L=Sbl-Hr{3an|tSFplVp2*y+gqpniJ-wVB3NbGJZj_GP4k;p zliv4aI7%CR?*e=3N?C&Pg)m6DQ9AcHCdxI-Qu^mJsth*7DM|f(JH)kDpN5MN)2}$& zo#auqW!%ihWaPE42e_izwh^{y)%bwLS;%ntztV! zOKBKOwVf_DfakDB+Y0UCv4QfY0tv?fP_`j&8oIU*dddkj_CbL`=2G*?N-_Kssx>9j zP9Xp|Y0M}YIIAWfv9B`*%kpZvuf7X&I*)k|x{6{JiUB8S0^5A^yN^w|KKIPBYT(I3 z?q>9zDUmCTfAls_ULVdr46FULo3JIuFn+uG5wL>c8#bVGG z!SpGJ?9HUClvKplK{0&ANJy(_IT2Ka(5mbWJdxrY02awn_I= z28|YD60XI1>$C@TQ2_J*?k6@AfmV940S{km(g3TvCy*@iG)oPc6Y0SSURrA)I&0NI zUfR1}mg;fN!a|fK|H~2oDjM52L$cLSnG6N#jYNREo7-73Ho99gwyFDNVY0ogphR?xR`w^29~dv4kJ1OC$2c zS&KM1IoJG+qKYE14A8zCV74hVH8rJNkI!&Ow|NmDlZVY0W6?e@8)!I&3D@0BZ0xBc zQ9b0jmP%1Q58@cUF@AM8RpS2gbW2Ss$c_9CpRjZk1|2BTt)=t-S<;Bn1;r-gi$n?^ zD!eNp`h$Uvjvh*YQe7wKm#fp}$dre9no(koZ|s@m_Lr7#&vd%7hTn2?Q+zIA^wa)# zIc&>n2S|U=0?JNDXEm6__#J}$>9qB7w5+oI!2M1J(!YF|W--Pb zHL`&4yFT?ZNJ*RS;TSic$|Ym-Yz6Xr0{E^Rik;l2#BoHF$jfv^Gjp~U_tg*&;Q(fM z8#AOU0+$Ou*m~6qdpx>&*Rn>|+D%M!>$<^NQ5UUBqqdORf|{=s7OMZ;7QDaqaqruzPD$l}A>8>gmZx$g!+-t}RoA=lY?v(fjqK-;V#fY4km*Jt1qS8=+u_mLrzO%Q4ZMfE(SZv&X-eYvQpFXRl_ztD0}UudV;Y_!@TyVl%Qi*qVCo z%?0<5oV2q1i@1TF9jB$sO=B?%8=boQCL)@KscLfkbr6wPgK8Jk3*7jerKe^`X zA4vqzfL=*8Zvvm=0>0576%vQIH$XFhuy!?EG7R{;ol2<9IqY+`{d@GsK+&(&j6XwX znTS1>Xpqs)(@Z2~Mt6`@)+H}^I1*$ozKx(6;sZ%jS!=T|Dnu`wGfU%Lugs?wN^;2f zBe`DqF1sX2o<@MT)Li*pIKl?odF&sNhof5OUr#e;{$1&Sbw1=3bCy?2(MDl*^w-Wi z9|WDnpF*XjGClZKERbxh&L~K~+$g8!o416rq7~NpK<@S+eH=3D+U0LOZ}h25V%M?UzY*3!~MQ zijL9OH6Nk5U|hYkBsXaj<#!Qr?sur(`Z1(DZ&+8nR|AOYOZh)?l=zws!B8^N)olKh zWstP7cxtr3L9?PS)6@u-$d*%W_Z4NH4@AeqOnG5$^uQ5gqRi@S?IJ2080=MmctS04 zMC%|X*`>;mopXwg0bOC{z&`m})tVI%caVo~>#P#SNu0rd;q{M?;Ivznc7BJZ-@W?J zNw|P_^S)W$YA4>ZrIBCT{@m3Bpq~ z6T4#oaTdOyxEh;n(ZIO58%}L2d9%2%*7Hv^EH1THyg=Wq*v9(yKuUx+l_wJO7i<00 zPQyWm7189JJLPF9wQgZ!8CA!5`5R}p1n+dDB?OaWL{W6A0$t{+O!dy@e+7!pOjAUx z^R4<9d;?|`XZxZh?5*2y3km7Vu`WT+W=dkqaoy>!zNyq6Tk5J7T-39qq?lDBg7DX6 z|FlYs?1)c6t{#7IY;}5~g8U`QRs9VIvsN-zuLEXdJ9=R&Q05;;#BGM_!lMPz*}GWb z7ibhk-!QVdt?#J{9~42DVPQ>VT?p&acWOeipck%azqgBCe)OcP_1IUV1sT*k%%g)_ zJ19Sm2g5;yqOx)-5{*BU{ZkxOw^2%)jsM9#!@Rz2+gU2{?3YU`Ow=Fi;-8TDs9ku? zgng}&`Xd235pJMSj1Dqt{SVk?s6i=pa~5-w$&uZDNGohf6=8V>`wzLlRKA~?@FiRm3@h{JA8r7a z%GN^SM!3*5hoIQd-k&9#p6S!r%%?H;P`R z2h%ay9W4Q<$AIAt*Evs!(Byor^<{CFl!F6HZ}BWjY=vP9RaoUwLSka#KZ<5xw5xpS z(Vjj9u_~;jq-2i4uJuicyUe%BSZ=5AywEXi0!0D6{&5OtH1&HINPkX7M!$$BzmGqk zdQ?>}f94I@QQjsN!N->R{)~kTJ$Z=+XSauN?t$@e#J^_2E zv=mN+Bc(lL%jfLtbw#@M@O>C9dWlB=w^U2%n9;?W_jZsqTtxE>Gl?siC?9Ng`1<^{ z!B)KEoLlWR45H&(Q@}veZc6L*GLxKKN-jqytfAma$e%l|7W3+$!+Q^tZ|RdOSR3FU z)&jRI!q@ypyJ@D7%)R9nzB#nv?}cG$7TvKAAD`}*J!cKyn?xa|mcb8}yu}{FI-nmu zU~ItvK)8u7E_pbN7l*bZ zF4!2W%ZVcTR^8`&^P(rakHZtgKGb!Gl9N#2w9eb$ODd^KneBKD`^EZW=E~MCzL=<9 zA)9duql>LX<>&LhcaJ6V8(Ra3&M2{~A0J}WcYLsm<>cfdj8yX^kV4$Mlq0|V#NH8q zaQMWzGQ_zX+rQeZu}_C=@H|O4D_VYwiyQ(7qVZSlqPFayoYsE544Jbju+Kn+`3C%8 zPh)PE6jpW(GM0(HRXP4`$U}*5mHo3p`x|O~SBjSrzePVb#;dru0$mDcc`B9#?<#wF zaY`X5vNF?-mxbLXVE0`Mnn%r?;{a=W*0wpn2kzcBAr**CV4dbkmIw%=urRvGo~kxg zNh|36&6iFPlfAke|5|R@So;|{t{L7_+wIEGsXP=8a-G}%S{+HPtxLEVt#1$-IoVNZ z_?*?-#zvd99L;dd4f)H&a9pGQvB~tY$rd}OSkU)(3*2c~68JGb^VMN^;<^+D?ujHmqi*AmqhC!5Z28u;e~wdXi-N zS?@zA4RGODAw0ooXHco(GXlwVALHsiD_N)8dtxz=l)}?OKv8UTvi8$u{$2w#yjg@1Na$ zT&->WC+0^&jI@mJ_ja^!AA_!+_Oo3AH*kki^m-n4M??flxqqjOXWn?QN5naWR=8Qr zccjL@njt4fTK+o}8^`yWS|-VrIa>^_wdt`~{J=K@a5;zQD8ai2UnetJ{x%9wHtg23 zUE8E9e_gGXEnuzlT)%xUKea4*O%QHyhzD^-{zc}#Jkp@w3zvru5zy^gzrtGO-E*XOK)O+)EDc)<=YIXcQ>ZP$Ds68#(*}K_oAv#&Oz?LC`TLK zb#VyBU67|G9akXJF#Goez~K)sHvV^NzdCegHJ}#ZPVcX;RhECdk|^&&-(_}RBsO>6 zhq&Ls1{a`wuyZT2n^dn3`!x;RDjyzR%rhUHg8GqUKw=eL_( z!=Kr;tBRvykCh&ZqY||C4)bK>w=Q_&Tw~CG_?f;<{G9P^Rp_#39Dc((aa&j5F6@AgY4b|OVnq*8=$v7cVNP=`bwMo1?i;@dxOJ{G zl{`@Xvlu1z5rjDBw9PMN`5jH5UwUU8d5~|)om-UQ3K_<&lm0nH1~NUD3r}~$pKMf5 z&|9CLr^R{|YH_-Rq^ILNjML}6f{Gx3GJ=3A^La*bI6U%lMJWA zcTVTDH1a$C+2l@eJB7#t+~IHKQ?3?5-4MB79D>_)bv-doRWS zUh}u(lIljAcCx~!W&U>9gAn?Gi6u@&c^%xAMGXG43bJdJ$_s@Xp=5(a52~xAXTIUs z=BpLT&-uJ2EYq#Ef};(oZ4kf9+Wlc=c#dPpNEU}V@R zS0-O9x$80dRnyP(FfD1xN{j;;Zg`gje1FIVqC_0M^$nUhLT>P@{r-wF#0uC7O`4~P6e=VrinxmXv>Uh81aATE0F4&MmWFGe2RMe2hG&$iH`+7GsgvhhksmkBL06@V%FcYlr8so!&LmSQV!*9r8#q+#1a4T9wU3 zEB}W$2<=6V>+zD{b+e=5O)K2R)c_QV12oYIsc+)8aOix?!^$<*4&4+Vu6=40;Fb-B zb3nz>Hk(A82ZN-3-Y?92sX60s0FlCjbP21ykxa6{5{Xg+M7(&QAR% zY7A703@3gKjy^6r5P%;txIa6C$6)a$Ej^RUq5+#Rd=idN_t&kl z3^@3Zh^F8M*gRj$m=V01Aa<6ZoCWNgZ=npdgksys+uszU3!}b3fiaSz2~c zmTqTj;IB~|)oiDF_O4tkt0 zhuJGuBu5xjmy4FEK9|k%B33M90@h1vc11)#LXLf8HGFmdtA@?hE)3gidKcn>{_tdA z3WB6s-~!KfMcl-kSiK$Gf1zJ2==E*R%o=a8=zQi2+$@{`jTN;_Q`yhMQ265%`3QH^ zM1C&H@~N7)noUv{R}OAely57)^N0_f=Xvf={_JA%P)sP)*$=bPzR5pN|6^Xq9iKNE ziw~SpwN8fq_2k5eS?tlc`69Q|X3lUu{;EDeGqd;x?L~wZv8uk)sLmL&=O*}LUN(?& zO#fqsxqJRhHEiTG*a3WrPh>(}HDW*7>TB)R3bcu8BK!tH>A9$9QzV~H5%7>lw$5># zHMlwjT7AEupv&Z#(URzBG^h_J<+C4`*j=t@e~>h#Fu?|(*F5yWeut)$eE=edpya`e zi!2E3NqKZ$C7e@Pz?b*a-qAY5&0g53ghu)9<uuMGGG zs_7LmO!m$xAUA2}>eRu8Qio11VR>)D zU}#_Rm?;-=ZgQ}9gCa@i&?0AR3gjt4p|%>^_C1ZFVd8$17uR|Zjpxb``RRN!7}#xe z#~^dfX*Ix(8C%#@sweKBaB4Y6pC(e`husk3InDgR`y1iF>TUeGxXI4@2)N7WW9O(Q z^7)WEbk*=gEKs!uPUffaR{K5A6w0$&GVl5^=+|QWIR2+OF}z%KS0z>PScU*z>pfTS z?5)8|*Y*(DYv*9jZ$(52_PLyL3=F`M9V#K1Kb3ztsNo&?b_WOAT1tEdG|gBJlRD5C ztqgOWxT(J-!<*ju0QbMu5~SNTaD11?I2+r$j+o<8GDuQ*oPmv68azHGZ&s zox}XZdOQi|rlb*kcEPcIqniG2xqRU?PUq~0v4^4NBvG3n%^v@3%3elG+tc?dthSb&jMpnFS`PrZ(JPo5#M038Qz}~S*WZyHZ0SX9>h~h zVq0y%t2ZV8j@rHJ zr`urGm|tm_G`?54Zpu_}rv2?95?hitU4|DRAyD>4!^xDj9=KRoR9uWsaL2tdW$Zq} zhH3KFJ9zVa7$zrOT&3@Xa4ZIuWhVa^i#YOR9$wW@IROfsrJ z6py?BVi0ZiY$#XsWJdhFZCzJ|rP;&A)92^sqstd&Kg_Hs$`uHG1E&vb;{NLE%RS9*g*z$ zHaObpSwh|`zK1{ZBd<@EAKD*IHbdbaAKiKfl0c>p?3qLa8Pd&ra}`NAj_iu}x8i&^ zK{Ar?i(%HQ^+X`zdj;v71T3z70~Y6LMSKxm>|D8tQXt)s=E1inEv>Tp!= zwA}1kNWBEH5|QRf=dvFCm%O9HEA;Gy;l2dv-Q{C|bnnkn$Ls=Vd(4fT2@xt?3!F2( z&_tUs>YZJ^M_kD;(%1N6$WV*Da6hP{Znit3yXki$G_y;z+!C)c1}_MMyXhu-z3P!o zlhuNo!#~z6%uB49YZazbdf6^Fdu|VmMd@VB7`;g4~@qYZsJRgCkid)O`HkQXJ)X% z-1poQ%dRk;^w`hfPHqVs#1*wZ$1c}fcwHzAXqe6J3RlB9%*7TO+3x4yC8OhquHEXc zqZcIyYdY>cbAEes<|hx-_2$b}qJk7{)Bj$mM1O_d{Ls zW`KUl@zDu#n)Ru}i(lsisJ`+c`ldUyd?2`31!(v9NI01299>YaY7sYszNU+E!KLX6@sHPf^Cq6kEy)Ln;IIHNp1jdi1(X+SptK-s7xl;nTwbHvUn+( zbRC3-Dm{uR&ir&Cr3~QO#po-e77;6B`BBFTG<&9BuSC`%v>)@s_Z7xx9&yTyjvL)z z^!D3$cc!<$CXR?g>@xmzW&n81F=v~C;)w<4N+7xoE)ySrh~aY@cR)`o4(9UbACdXy zz-k$41fHC*g9)mEwlJvd*6~!8j=&Kkn8DSQ$W;+U#;g&lTf@($&rpf~^tTOC`{yHAJ*NEz5pXB9 zv1){*+HK*~BdgyEgC*b544O8z8@?4JN^RQxSFPM-Xt5QWrw5{WZY$cAA8H5Yy7tl9 z5U&{5wO%47rovdwti@RAN#jGLo{4QSaEXJxK7y+=RK=>u_|sKeOB@)$&BXu|<*?B3 z2EXN(V3T*Co?HU^|F62Qj*9Z_)~1J&E{S14TDn1MXlW@Ckdp3H8f5530TGduR>2~q z(*P6@1r?Q&5~M>!`1btbJ?}ZM?>XyR>wJHFu33w5=9wq%+SlIu+SlDbbF_~N+fX+= z`=ILR@|hFUNK(pO#a)|pQ@p5GycdpLKdve*-_NjB$TMHG46e|=a)DVwe&)5^p5nZ& zgfqp<3sP!d@ch%9+)Pim`yOuIHT|~tsPx$)(}Q?lkvGJ}6jDLhHWU1L{+o59UQpn~ zlDrRQdK!rJYagu*HRVJdR--zF0&mpd;|g7%>ck^+n>h+rUM;%M#3BFa&2*%#uIBZ? z&gT|X6udM{kLSHocb{S|*$A_S@-covoc^)lP%$f!dzLa)fAF9Wvu@++ByzghSKQ@? zq#!Dxg{0S4M7Cbz>p`t z56i;O-(Jf~qHy-9r|0yT4i*)EnmU|BK>O2~g`^1ilz${GDOAv-*ZO7CC2zM@iAN<3 zl@VvQP#JlZqJX@K?oxH?Izv&?hM3ruos}7hgnHhrargooZl7stMTvV|%l3)sO>G*gd~{S^GL;(J1Fg@{*tZE50OtMunjL5<&roD9^} zp*HlDEnBz@k*Y@byeHUU!Lzy&$BI^R649a$zrgQw!x z9z4Q-fcLlZ=Bj(e|sp)@sX_TR`Q zt9fa|kfqe796n}Ol_-V3h-qe7IF)?jIrY!?UAL4@FZt!#oj~f@f3w}~-ng*!jCF)RV*J^bM@2w&N10_p&_wUm16H`wB{sM|~9L=g5EU&(we?7a*C`bTJio_-pRbDz2|9!vQov9+tvvm$bZI^7hXz#S}(Nt{>P1^wRUa9=c=%^D|LqJ zt`ExEKJ-gX8$|UbA0L{^UVpW)W{L6RJkSnr$$)CbBB?0y+`G|yT^G0+4WIfyY87T5 zd0#yBQm1brsk`T+=wgqrh|}Vl#X$@k6?l=u5!nqx%s%N(3{#WaLtAt_d9{Z<<0EON6zsoZe+dBenO3V(Iq#yjoG4kw^oifiYmScuBj;}?n| zl)T3^>|aWd{h`&DP@B3{hgJ8VafSML(`Im+PsCJh7 zs+-^ad2>S{Rnk4W7_roFO6gorqt6!#851a9UDh#dJ)CLuQ{piI%?Z7ehgL7wm4{8% z6S__X3jrmQmQg1hpQgUEUjlQC->_PZo_Hy13M!v@JSPTo#JkbeoDO{4N; zdc%{(c1GGU(Qu4PUEq&mEmG5g2fowYP|@Qt!o|jcr^4iy8Aqe@Q{4_Ka)J1*^j6PU zCDHJ@;4OCs#lX=!XDb+}uk}C6VDA9!g-lQ89?d124k1o9dG5bOTK>{%`M>$o!a+q# zE8#n(68ssnW%K?i2sxy)%F%p#GxgTUoP#{>V}loqr0eiZU&?KRlZ|N_8X7EwTH)bT z1Z>hjd7NEb$WZU|AYW=VlP~2#M9x$XIpp`X(>$>)mK?6yyj)y~DxJA}CFuk7?Z})5 z0?0F`DcR(?VqAZ0e%(*9dg7j={i+D56dnvrSP0Ydu9Ye9*yjM*t{=P|D1z=#gz}Jy z5kapr+Ek>_2cZkMuOOoiH(C|G+`m_@sj9d6Ic%*BC{e(Cj-1QVLWVAKEPd+&wp3l) zuH}m9OnOyS6_uKbRsa-^8$o9K7&Diqe@psoS}K#+3?$ck11uoJ3bh;B^xgW`3JL~K z?Lvs$Ht4Ncur2-7VMy!6h|d9iDVh>l^D^i3Kv^dne8|}vdaDnJd?&RkuBOx^OfJm>!P6Y-AqLVhiB!)v!4oq^%xQn z8Yb16#e|RL*rZcST;LI|!afipV_BF-cX@btW~iB^J@YdE;LJUWpGI!)t$swbU)PNN zF@>S*%rR2gg&Zp-MyuNITC6A@l!^9RP!c?Qg5~W@kO;erBwmXOVpj~6EIgT|npgA= z>K?DGedeSRMYjH!za?fBM?hJp%{yrzh76y&iijq$K77#-v9tVkpMFLJBufGS?9UqV zZw$$NIXwEgX3W}&v5oxk6&MZ=A3nTuzDJ*Prt8Png7VIoi30K4+VhP4wK1%Z6VFK> z>oxd%If4?#`5k_03SyMa=Svxi{f%@=74m_Lokg5ieNL3??7x&jNA0+p&}1A+BuAgdra1+e(in9_%X-g-86mErLZf!g_vOf z0LP(qeKte6KSA|Ka*&p3ypwmL3xg<^jYf#^Z*(00T3fEIBjT80xoOH9&c%sj6?aTr zfBWlR<7b*7*V~F=3<}RX?m9Tz$QF>ZPj)w#?|bhomSt=bA3_G&>tn-L!>mp&awJV5 z;(;JR5p^8tbU}ZHF%0)<20cOooPBsiEPY8pe)=*=O5%F|&I{h6MCSD7#PhZhIgF|$ zpxtxmmd6Xc_k}Hf)fL^;hCo#4f}-R4xeR3uJZ10wq}dz9GXe~% z#Pob@o5XhIeC^NIhkEfb8YYMQH#)Bf>oxnj?c4S(1ZSxR7Z%F;m8!;(PPz~V5j3SU zysF*5vv21KUXVtX^>*$Kkj4yeN$cOCZK!10NQ&%qD(Tc6hnwOK5SSFz~SN zH(au3an|`^!nL`z)rC3CHUQBa7`I_1c;#wHTm^h$)h_*pDHbuFjdlJDWA3fFTc`K1 z&kKHP5?%K?-B^6l1`IVfw=?l_#MQ+b65)aC0lzi}L{dJayw$Q`>G{QJ&rIb;%JvE! zxHJ??Thav?U)zs{_d;YPo|VEE9p z3}i0HbPJ^y2YHc;=7g~zqL0T!f%mZgq{6xGI`aMO!PPQDeSLao@{f|8XsuaF7<%{d z&+y=oy>put=pTX5UXFJ-h2f-V3x%FaNL~%i_om)?r%1GTW|zNGXJVVe;b`1?tnP&S z-8V|UQ_)5-aTbT(Y&gVr8$IcS$DJ_peh&eaJI|jgYAO~SIPBM28L5NE3BGn|$%tHyI^VS2wFC=L4=IDY;9o~EadyU; zef^{z&w)4Q_=naR@l{_P@c0#b}6aa9l+= zJ5qNi<{2fb9w$oG;WF~#_`}3OQBMI**CeKTJ6iH%-B!#*|-J zoLIcSzklXBC+C_hE5$1cPFa<14dq^eH&#Av^dqGEoK$acd-~tqv)ZE3`G=I8G+N~p zr4p>q;c<#cO4kjrFtgks6@Wv@DM`r%^nEqv6tDWN7-iAeyC#(sb%rE=HUzF}Cl9gf zG+N!MZx3xBC$^eo)X?g<_PK3JDO@s^fIRW>w#v4t**34ETgb>)86N>ocqahQ}4RF zRkJb;A~=nVi~!$6)MUj_UA=aVrt%Gcw$TUPKV$)`472Z$j;WRBO&D!HKb2wU%WL%j zhJ3*G;>ct<)a=^<6x4?(M(FIUP%KSD zzl`Dt=X}>%Q-=^`_?3-s(%ltLKQuBH$|K16g!!iYXt_ZlY;$u{U3I{qt2;NRuO6kB z9Ir{K!OJW*lN>p6{}X+skQFIA_m$cMl8-e~R7w*pMf(?mpH9c^o-b&hI2-)*d)(z< z>yQ`mE|yzvT+(L8=Yo?T+xBZHnmUHh;w#fd2<<;`9ee(^w2$7)BS3@K|H>k&b=m#n zySGgnQBiF2bS^Fu;ia(^vAxupkL4*hzYug`AUm(%PK@-NY!b6ITSrZ#L8HHd@%MAJ z(U68UY{YT5rP_7c#@fBc+ApH@0;~pGQ4w z)&)kzgeu6mJ=1)=r6WvciNA}~c&jh`dXC68FbEH4nF2oo@Zu1*d8k{(y|_DaC>uKmj21%OM|(aa&IIe6^^eg z{aE;v!f-QiFAm#aG?e~ZrJ(*b=l-_Q$JSe(ZmSd~+kEVB^VeLf2Ct7iRsj)Pr9-!5 z*2e2>!!y}^+IHT-Kz6bh;~ad^95^6!o@w3Zozv(YT~H6xIGR#;-@Ii!lwV_3XZ+PT zy!9H

KLY87JkNmKvoi!HhsqO*4_36FM+?wLz1fk;esRFr;xezwx# z`siqRz9h!zY%eA$DT(Qc-)a?rho8HvS^W+MsXqaeCP=>PpK9T7C0#1lz2Si7x1p1Z z80k0&uGdImutsxwbOjv%?5Cy7zcP@B80s;9A&ynuS@xen2ZiHBJ~68`2X~riS6>16 zn%b<2p>$?YYP2HpTA}CdDzGV`s>=2{kG>7aFZyzRE|m)bhA=XGy$dLot5AB>Em>i# zhAk)`5vtrkA=oyGziVk<(VRg?41;iJ&(K3_1O@8^pRDMQj8Y(&6gY4tR-X%Brg=+h zJNV$LYJ@WNMlrGEY-KNo6@35xBJLYNOaN(UCq8;bc$3Zcn_9noR8xV!bjzRVyYg zeofb~KsIYnML0~&jlE4!@GhMCZswSGLa@;pb5(t8Vcg)Gb?5%oeOk+Sf-`t=*KjJ- zq7lVI7(6RBB9j`n%+0nmJ=zNUedvki1?+E`i2)dmT;@z z+SipQXoZ9%M&0$XT?6xJXf_q4IDttSKj`*VA?dkF*u4*e~GAB_U2kS4y1TC{~(EMC?pg z6q5N(i)~f=r?Sc!e~NKBuZBTvyrAEZ5;%W9^>vIveXE^!N=(E17_AW}|I zPNoE6?E@t~yx3J3FIEaqEWifFTr-6}P9{-%0R`c+Gc!!gj@|-}5XKV68p$_>ao7QH zR}wDBarU3a`KVIpZ$@NX&hWP#u;}YhmwfMHi*)HCDW$mv(_=LSS1J{^%P@**XV6QJ z>>FZ(jo7<;F$cbH@@ip`un;nVm3*Nv%s@Lv+XFi+Hc$bvznv}&JgRpuwtj#5O_XV{ zOQh@~C5-Yb@xLAA_MCbK2K{#Hzyi>_d($yCI8(0vx3i3lw&de#J5*!tf11?u6qmtD z7T%_Yo}0~dJ3giq{sj&{QX$xgP|Q-K1}}X>7pA}1_Bc%KjzqtMKOe&b@DioscIG{4 z@?>`;oWNy2p2fTKe#l*&n(bpsF5l#P+i|{o72#Br0+)KGNH4jJFoE!13+_ zyb}LW(xEqT4SeWtbsh10=SD7`Lek-O@4}Q9Dzv;yP3ew4Tcf7H*B4ttpuktEe3Wg^ zuc;@-LK&)%7~}|oQ(e9fKY{2Hhgy86U)PtYW5vKlqrr^b#SzSTb-rYgj>BupUuk+V zoGKHAG?WkfOj1D9R*3Ix5OxB_=ER<)B%+fdxL9+pOcT<%5G-nhph||)0kQYHseoIx zbk15pUhRR{T*|Vx!fyv&(kHIO+H`GF$fa+H?m~<#7S-xX7fFg0l9r|W%DFZF_3Kpzr}L|zG^R+eHwSapZhHU}H*I}4UiM-X z=yA;L5GG$FF*q}!yN&>R z)gQI3jxT5guT|T$$#c$J?Q2@0<G4pd9AQDxV`j7&Z&1SB{9+HWukDj7c*!^ z&(d{b4ttA@8J??MQm7hL$tjQ2ENFICh8x+%^{W^#=GoOf5^eZU+%q012&LAEnC$GifPn zqRZB3R?8NAVe0Q4GJ$}06{X@ICx1cc+2vr-L+Qjp>F}Ot4xrYC+Zo zufN~38gFx@j(RE0F@Um+%`52Dly)V3hOrnq2StWL?>9_sBX-_`-?oq)1!a@qPz*Tt zrpqbewV2|z>@U2!eL6=?b*uDYJrnh4D_;U~D1|vBdHZX5$hnB7k=81uNy+)-B0`QU32{O3+=o>vO7~|vRc96W)O%95K`L8CE${Oj zMJfS0{R3IIAot}Sdn8SYM#ugPnP4;8#*-^nN~wTWNh(+PB6H(488VUvH{W&-0rR`E zq1&!U*Z9z|e)R5x;Y#4xDBf^+VM~S;VeVhTuAM`5))px-S*Y;B0>!Xg_6oPgXhYro z=S+-@yo3HL6JkenMa<20nXLr69`$;3ahy{zp-E3ee6f<~UyGGozYrXE=4GLflK5sp zw|X(l`+OIPR1`qS9KP2tN zJggU`gLDRX%Irmia^u2>9povZ4m#|V%=}_=TL20h^MZ#Dv~GjhoBz{ZoR*hoF4;I* z`#42^{6APcbJe@aAs__wpuuFBG24uq2%46ez$pKAaApq}(gV*tAr636qJSzUAJ6H4 z%_H`OSdwKnemPV>vxfW@JjNwawM`+FM8$@Z(cL+%>Dz})Aqr!(b?5Q>Ns41DzhSDI&i6Gl3UYWJ8$%Kao;p{(bJ`}C( zt0&&O{919ZAd?{LookmT7?ycr+TX%wz1T!rT%tqBoP{7<^#&hGsVY2YFa-26hvG&$ z9=ZOI)K*xe<-9;RgT7#qN8G$qs)L)6M`u3?mcWl(1K_8e3cy3kctQa02P!;R9GYvr z-j>8f3o!9*DJvfvT>Sa(FMc-`e2UQZE(Eq_>Ln`R3}9{!?%!V8f&}D+u6*=rw67-T zoF0`wX9>>~(2LWHb?iNU^?Zljob+Yvh~>%mn}LROH% zt##@Cd20+JHmMro;D;Yp7DbY-xm(`{@63K9Rg=uQmYPyrEARLH&wbbY9y8}i&j;}H zLXdj;`h;?Y%29=8{F*B;TL0ZhxJLGUtn&&M?mx^#O#7VJjQF|hCg#upjrQ;H-~RRH z-47mJ(e*rYkb0qvL4ewMUbe$9^fnr&{5`-%J%A^#6GzVy0cpRGw zi#wspkW@%OhH94)K6;yoxOqzNw{yCl#$oXaOuU|(gCMAlYg-1HH%J)ni85I?`*YWF zl7HYpN)i`Ko`FS)WTSZdy{Hv7^=to!J$?NoV8Nh?iEz&kVV>+jUA?;^e52N5{62g5 zIvR}fhW)vGKY${wmHF)vKQ~E%c8JY{+?Efj*T-+e+(`KYb)^pidwI4`eGM-K zc*bVh^1v^;F>H!N!$8CDc}}js(<6v$_O>n(9((#^#&WIJiKe5so4*6 zke!h+p%!#6z=uyyOQUy-K_r}D;i38dkDy!p;om_whU)2Ta4wS2B(MyC|Jj8&s|8*F zTQLW$3E0I-fS(^X1*_)=r>AgW#O|hWSPp)lhPwKCrS%M?OM8N;PIZCtxrC!0{dIt}AU3lGfjFS>0havCTNqqgRYmNdokC(}vjwWNK-d%pZf3!mE8qKB ze0uk^j`HM#1X7weLb?#NBGw(h%JxErq zxBk%DJQv5PAn;`8bHtU6n-V&t`<-d&Wr1%!<>n0sg}cdr?kiHJA@d7uOdD3{Z)>v& zpa1smRphU+Ak+I$hh(<0*fe(kcI|nQ%JkHS7be-mjcC|}!d^5PX|qO_PwNqpUGg1+ z>?$4#FmSz1T@`tCP){Z}6H zq1oeID<0z95Ejw2`=Yt#3^(R+#AgioO9*i*-<)Jtl$04ifgBrkG>>69KqRQKZRRt5 z*MR!U%29%L`pYX5)SKX&v>z- zK=y{nDAS@@wK1WGxt!EU2iDd3cau(VK+5BtU_Dbi4Pb(xx#wBd%}n1Fw+1DwnLW^c z^jdJ%jhxInCi7QQDI9YaT#vvNuo2dw0!0K=8>z=2>h`VsLZpo6U>ij6)vb1_ySZt(#DHkHX9)VMO~ueE{|r68R_A7#}G zqG26$$Iir)F9;KgrqG8Ty~MR({qJVZ0gCW90)f?xSSJ!V!sw3$+E!R_X|moPzwrjH z;aw%*)lh@CN&5~g_rFHW7|d1Q@05a}NIt){iqlA@bw~E+mIP0wn~bPC3Nst{anse*M2SHZS59wV8f#^YW^lY?;q0 zl<`51g14p4jO+V^Nz5wYzq@d^V)R-%QjydHW{75f3No=OMtcjF;AVFDM=+?~i`n8} zf}TUvLXl(s1s1~mg%j$dsP3za;e>;8V1o@ajRVQ93&SMjF6mH%4c`R+K73=|4I5NM zZ~9KDO~ZtjeUz9<7U^gfiw0Olgq!Wmi<}b$oqi=2XI%)8ar_bET5YCam5r8dhT^^X zN71YJM&S|`c*mjeSH6ZIT>Nk-ds{}cSdqSJwFeynsl?2Le+I>pJU1bHqfRWQgCz4||a z_HKy_PXpav#eIOyXudr0I4?YY?6GD&Y#c}%!w-0L7{AvE4hDm`O?cv zoXcBTT&m#}*EqN`5})tr-Pz=!kN+k{G<)+ITi>E)`G&Ut$)3Tcj_*%}zOA$^tbp+3 zg5;~C+pohdJXVh_!`9X={5;$ZJ$lgRy87YCy=*E1mgVBv%xVlCf6(wojguuRXqYam z1UG1yD&pI-yU^75NZ7dxe~}VDp}O>o!~Jn zJbYJF)EuzzGEhHSz6y6HDqEh9$r$mgPPaW*g>5@0ufCGv!>qFyBb_^}bKSCHV(0Wu z#t~cEYIWFk#lWUS)3%No;1m;>4pr0A&`f=B1NsdFP<4QLWnyd`i-x!$$8A-TY7Db0 zZkjQI&}P|ArS_)%^cY%Aup(&;ZK6)WTp$-fbfRZv#q*mJR+b0C!L+-19za0S6A9r6 z&npoiFUAuZR0)|;S#apnBp&LQkY}?&$K65M=Pie9}J|A`XLWAwGR))N5#;qOe z^(V1s#j$&Wm>4Q@`LwBTwECvSOg7_B+^iv&Zr2 zsDEZtF)}i;Id;$rj$}{g*PfRM4;GQHTU~?fl5Sga&*zS5WkWWNZ;FqbvhEuT+q70y z1e3TJ-p+s1*20Kk`3k|Z8!VE~cl@7(WMy5@N}qFy9#?Kf%wNe~iYHhPLoypKF1@+Y z-J^))?D?6KbpGc~i@1As>H+im63NaJ-?d$Qgbkb0Ve%E0Z(0GCPoB00KL`pd65QFe z7b6tTE!LNEA5*#E?C5ysIbjJ_sWFEp_Lw#Vw zJoJy|GqGo;faDO=wka37&G#WAo_%}KA$oM3ik2I|K|nD^JQ*kQH4yV^;87NiX4yj& zLxaVd&QW-M%6|hA1&E;>kg}Cq_;UP?>m^@};{3H2u9^6*7KA@i(+eot_0?uc6r+?7 z@~*efevZa%zSM8gIBs#U9qYb%aFR_y)@M%pnSj;j>xT39Ev6T^zuvegBFtctFXd-iurFQX*4+2(osa@r^c$^rH(4>v@3aWp~=_&K`#0&`i>f6%cIy$Cp{wr68f?< zBDU?}kQGR&7Qu?ZgKdcg$a4C$X;t)FT6fJryHZ?EV`=Zn#E?6#mNi7;dTiNz1*Ap@ z!n4c8CQ>g$<;Kw8?waZ@7wh!A7oAXDn{vbJ^PZ?L0i}S_@D(o1Ho$BvIktF4a7QNY zUXYL!s?@v5`OuDoLnoS<)3PHnBK9kS_EIKsR^L)48G99%Icolj@UmT`qi3ZFf>3_B zcvZmh%qf@w@K({AA)rLfAJ{UY<{vQ-d*mTHA(V?66um`wPg`oB+N5%66KZK`iN=T%5B9HXv(ar=VUET)!$PVtyG`d1tA3C2JYf0)@a-NOx^w3a;6t@M z%{1n?SED4fJQji+!_H&qXC!bF#;`7#4Lnv^F6QS{!(T4uU_j6UcM~%$_ueN7M>?`= z@Fc=;+VI>MI%@Icd^rWF;264RnR-|q?znXf9X1*MG7cAg_atN|B`;$@Bq+l$1m8=K zEwDf|$Ou(VA#d|UkL(UbKBu^3F$J@OpCKo#wrlzFe3`}wwmQoSLSuwJ^+{$Uc#tOE zEcwO#EBb77+SQo!^FnX%sxcLb-w{lZ=QA=w zIOPyqM$5ZWN%G5Owh1j!x;K%w)|r7TXCUGdDrCU!*GqnzASZZ@`819K+j%kE00b@QB4qz zXI*0xjFn$R>7AnsGV(3{mnRmbwra8a5f11^H&-yp%aO;RV9{>+h*5)WyQrjsV7VAa zm2J#46{QY+hcbb_+md3uOv@(jNW7c-;e*56B)rr~XC?h}HhOJ78eE0jMkaFrhbdsGZSqTLx^BDSxIRUR8 z+v8_9v8(W&wQn8^1HpHL3#(qU2U%0Rjk>FC&HSLNS^Mb z){@=oq{zbV>Xk8crS?nzDhwWBkdb*(F<}*E*9hAvWH+thmmNC4eWpSfjk8~*{x=8Q zB+RP7!5HlCvVnjVU@D)tHE$t! z3?{}k*dCiu1wuhf&^O@wD!KW7pKONUNq1Wx6xwbCKL2kmKEq>c2x^7Hwr9D57e-6e z5VBiui}q4P+N-RaH(NutimyzJYy*9|d;Tf{xd<}`5(L1+O!hL^IjKP!B&sal@CF`( zEok^M{lnW?FvZ|4vTmtRYB?DN&$tE;LBX0TUl8payzh-3f=?SBnWX%?r+Re z8|dn$Bq!4(tq1tM(U);@a&i+-=EI^ZWxx#R4EaSvp@%naz<2>9fIY2cv%eFJ-mlG; z4W4q~cb?S1nIZntJ}jQK4>FLCy*@HbgZYP*_MGKszh`Xx{QTlIMB$@dgM)+E4Jn6S z)fYzFAT~;`eF;lc`uXIhMBUDcTV`Y#_m64}gf#@G+AT?^7-M7>Aj$U-$0o1t#&*%M zH5ug|_nXeE`){|nF$t?656lsVJlmB_e$wA-+jYT!f$@shTY-h-(^TC0MM`W|0Ry@r z^EX|tUAyK;cNMtWgcOC5JtAwzN@R?}sYyv_aM+I)jp&p1q1f~(h)Ui}t-+NPO`$*i z&gh5UZFze9toR})6NbuPt4H|d%r82I%!+dC2^xYmO->aVXFvgRGqXmp%XG+5jX#p|oWa5?cc+(bZ+K*Q z<_~OlmTxc$KmHZnsRjw^1W%Ve-r2`)??*<)*UG#~6Ag}N5*Djlsg(S*f*%h^6 z`}ET_DXeZhL7ehMJ#p0&LjG?oXzuNci#1Qqh@LV# zVk>A-iwCntIr#p&UtojQHf7jrI{nLQ=e(vL9lh*+v#`29zGC|0%KSIBkm~BqI;Zy! zD?w>)WODVR(|$%YPkw&>tsJN4$ML!ffqdt?G8-=?Xa(-&5B_@eX6{Q@NXyuVYZI+& z`{`%hUv|AO%<)Qh8S^kxs1lQ|5-KcrTZ!O5?7Zu>>~vgL_|v#)zIvmT2GgArc0p-w zcE853#${YdJ(BL|Sz54-o(h@rPWdynI}xA`yE)eJspVIQ{P3mA+l%L-bn0)*@baE- zql}g6K;&0y(tc|HG<)+y^G3|UDtlgaZL*HnRj2;*OYiEe>C^4vfhT1P`fDNm5aVSu zQFAgS0N*<&KOZ)yRDOJtTJzR-n0Rzakchu3(F^>M6na5(dN_lO!C2yjFO?{d%l;hL7VseVaS%V_1tt% z0yvAX8APmTyaDXs$8*0G94!mutG)&XuRdM*g5?dzig`R%PojZ9xeAjI`pzY;a}yp% ziato zVD7$qZacB2t>VRt7ux9lXRy)P^Q==pE{cARDb6#r_-RK+ZRt6bYd%h5M<6XNy&GD< zK}I0zX!6KepDZ3%Y-yxSOXcMW-x5fK3!7BvMH|w@^g*zJRXK+v_)_|3nm}GYM%dntsuVYG2LU$y)pXw^VYYktpSRFfO0+hV@DLex{o z?&$DG-ky^?0}8)Ylc;JgnjwXk?m7WpXmX|?s|1ahoHCcsP}i-QL$9<@aUDD^mH)dJ zhngK^+$|EX-rf^k8(Xia=luwbpo|}G6K{uyJL8fpKE@1PeBey`k#W(7eL(rRI=mOp zu865SWQ2vkOy2kPr4i&)2@b< z=Fw;5Xo8^LNH^sAJKps(9~sKSOm`H-G%+K{I|J70JHy+=ZGRh&{q?5|u|gQxF|LZG4B= zTp;OnsYZ%*-?MdR8P5qSB=U&&nlp=x_t-5URgknCmmnRik{|@UbYe2se{=5a*_G)8 z_Wt?gVB5DU(O-;%9j`mzk5HYCz`5lO!Pl7$Q|hjBZ2QH=ST0KZ8OBx6AubHH37(FF zl|#zsGb>9@)%%kaA$YDh<~61yWbVF|u8-_aHno3^iSX<-|7yQM(?Rl`NdR=Xrq_ae{dp_O8U}V}v$`X%pF_T76u}t(oE|pquG`y_h(G zY#}V|xQLG%77?AN7n`6@PWQ5E0N-`o&afVGWlL~hBmCHMvhFxCH7YU9ibg-fWKG0$ zopnUev(GvCNTH^pXo{CN%(e$>htBT#25jASm|O-n;&j3WjQEipZ4*taR6NwvX|<5I&1ccFjSVb}Q4?j7a?kGX ze{yc6(Z@1n>KlPqQ4%fh(_9*F>oN=zRqwl7G{o`8qF)U)sXH@RB4_64(5!vx#BZag zGkK$LmNI2!H71c_&YFk|;4v-Gq>EDZ9n~_Wk25t%hJ%+k5npIb z-xo*5Sd4k7>#{!81^gtF6HG8{tMGkQIpnQGHB$`~$hkefGuM6fpSMY%V2h_$q+o%h3|4u3!|WI6Wo zE=SV|4Vh}9yS1gvY`vZA0;Ayh;o|BVP?_^KlZfZFOS!yarPkcrTzUUUBk%ZWLMKEb z;*J=|=>W&tTqX|Re8CIx7Ccp8nP7{(A{`*yyk5s^TQG^yj;h>0fgmlWWhcKy%JTSj zU$;!Ej$NU2dTG(gfaV>WoE@k5aMI6@#Yk}Vupt8(;r&pC#z|Na?%d&`?&ly!&7`rJ zVNE~B>)x3Y*lY^ZzDe5jI|7MGbscMK*3SmG40U(&Tm$C^Pre;hv~RkPPmucEULA?s zG)J{}wDRb+g9TL%Y1!+%U&1saUn$rxbEYm?Fkvf;Dmd&X6Lg6`nYfId=?gRx6xUG9 zB6Be+8#_lt&Ycpp4GbW6&GC3#1ZDoAHd}8i#oas{hQT-e=N^{73z968AgX%u<=w$0r#W!S{f)nw*V77j_~|JF-|k2+dc zX}`JtA-gW$C#khK0_|u5?M{&RBDM@Id5f5`rB641kiaD^6+}YyndgrSN7oYDQJ!s1 zA-+5wJ7=HaeKY3Q#Ruh-rLcCi)z+^&pM^6-6Ig4jG@0dm)P9|F+hcq7kE17A{lPnQ+7akEIg<*3}I@q~~ZlyC-P(X7(Hij3^^CpGzsfqCez>>@7 z;{QCWc5yUSeHA8CLbj(na%p0_1)^i??gUY%eNR37vvhdY2$u9Y%`d*av+LP<_y9t+;!cnPA4WILqjmiK z#k=!7F*%Hz7#>Fml$I;v|7F&VRERA`&`|fbUmV)CxEp)g>gGe(tjm?1=tLGCp}NKK z)jfB&*Z*{t^lzww95sTymG>qp-kCShe}MxVV%htcrKP1|*!o>vg!e1QU(dvcB48ay zPmb$~w8fH?pVCb29Annw>0X5tv(n6hF%SvgBybYf!ezOcqg43($4?uClYA0Yq7&_V&4%pVG+taIJU+c3+RD3WtdP24lkkPz zq0Qg2Bkw5-{~TFHNu5woPM}uVbz7Yamm0{3&pJ3gb1K+N`so!)w?XYv_F|zRMZcX? zVW62jMBOZ@{ve5WSU+dOeyu?DR;WX%0psNAf%h#z;Wb2$s_t1hr~&T+DRkIUaB zw5Gs|SSyhZ#n~>z8E@2k^0t`mn8ar6eL*IsZTKR-d`84^03p?^%Iy2SDEIa+viR~g z(JX6}<9<%E|JrxbfAyUeO607XnyP((lOqb-pT?+*e|N9=KXtEZh%a=ik-dvk07^vP z$7MAUqOyaQ2UC=p#Z#{iTSy!+#P*be?{#K`{3Db?%) z?LB>5{ysPMn#TY11G#_tfswzDW6&iSJ$zpFmkdo%=1%?rZa&^9F-c)@VKJ07pK6es zrz1*SLQO`@njdv8$lgCN6lJ8YspjJy=oB2t4|fQ|HxPcK4x^8hPhDPIR9sS2TvSX< zLKOa5ii+~W3y#zGal{^i`uB(Y$1?>(UWrr`B_>KL`j3AoX(=fQDU>tn@6W^}B~gET zi}L>KGYLsqX=o1YzkhrN@5P}7{>x|LlK*;~gamZfe>x7%Cn51KpGn9<3;vhm`_C9gjJLB73fnc<_v+Bo3J#L? z;+Gu7q(#LYo$Z~)B^~9Y9qh%OrN!+XB^)nF*(;I$&xfEc{<#!DxBbzk5>jH)5~SSR J8ityr{|EZal>Pt! literal 0 HcmV?d00001 diff --git a/paper/paper.md b/paper/paper.md new file mode 100644 index 0000000..e0cbaa2 --- /dev/null +++ b/paper/paper.md @@ -0,0 +1,148 @@ +--- +title: 'NeuralOED.jl: A Julia package for solving dynamic optimum experimental design problems' +tags: + - Julia + - optimization + - experimental design + - parameter estimation +authors: + - name: Julius Martensen + #equal-contrib: true # (This is how you can denote equal contributions between multiple authors) + orcid: 0000-0003-4143-3040 + corresponding : true + affiliation: 1 + - name: Christoph Plate + orcid: 0000-0003-0354-8904 + #equal-contrib: true + affiliation: 1 # (Multiple affiliations must be quoted) + - name: Sebastian Sager + orcid : 0000-0002-0283-9075 + affiliation: 1 +affiliations: + - name: Otto von Guericke University Magdeburg, Germany + index: 1 +date: 14 February 2023 +bibliography: bibliography.bib + +# Optional fields if submitting to a AAS journal too, see this blog post: +# https://blog.joss.theoj.org/2018/12/a-new-collaboration-with-aas-publishing +#aas-doi: 10.3847/xxxxx <- update this with the DOI from AAS once you know it. +#aas-journal: Astrophysical Journal <- The name of the AAS journal. +--- + +# Summary + +Optimum experimental design (OED) problems are typically encountered when unknown or uncertain +parameters of mathematical models are to be estimated from experimental data. +In this scenario, OED helps to decide on an experimental setup, i.e., deciding on how to control the process and +when to measure in order to maximize the amount of information gathered such that the parameters can be accurately estimated. + +Solving OED problems is of interest for several reasons. First, all model-based optimization strategies rely on the knowledge of the accurate values of the model's parameters. Second, +it allows to reduce the number of experiments or measurements in practical applications, which is helpful when measuring quantities of interest is only possible to a limited extent, e.g., due to high costs of measurements. +Following ideas presented in [@Sager2013], we cast the OED problem into an optimal control +problem. This can be done by augmenting the user-provided system of ordinary differential equations (ODE) or differential algebraic equations (DAE) with their variational differential (algebraic) equations and the differential equation governing the evolution of the Fisher information matrix (FIM). A suitable criterion based on the FIM is then optimized in the resulting control problem employing a direct *first discretize, then optimize* single shooting approach with standard NLP solvers. +For more information on optimal experimental design for DAEs and their sensitivity analysis, we refer to [@Li2000SensitivityAnalysisDifferential; @Koerkel2002]. + +# Statement of need + +`DynamicOED.jl` is a Julia [@bezanson2017julia] package for solving optimum experimental +design problems. It was designed to be used by researchers and modelers to easily investigate and analyze models and be able to collect insightful data for parameter estimation of dynamical systems. To our knowledge, it is the first package for solving dynamic optimal experimental design problems with differential equation models in the Julia programming language. + +\pagebreak + +# Problem statement + +The problem we are interested in solving reads + +$$ +\begin{array}{clcll} +\displaystyle \min_{x, G, F, z, w, u, x_0} +& \multicolumn{3}{l}{ \phi(F(t_f))} & \\[1.5ex] +\mbox{s.t.} & 0 & = & f(\dot x(t), x(t), u(t), p) \\ + & 0 & = & f_{\dot x}(\dot x(t), x(t), u(t), p) \dot{G}(t) + f_x(\dot x(t), x(t), u(t), p)G(t) \\ + & & & + f_p(\dot x(t), x(t), u(t), p), \\ + & \dot{F}(t) & = & \sum_{i=1}^{n_h} w_i(t) (h^i_x(x(t)) G(t))^\top (h^i_x(x(t)) G(t)), \\ + & \dot{z}(t) & = & w,\\ + & x(0) & = & x_0, \quad G(0) = 0, \quad F(0) = 0, \quad z(0) = 0,\\ + & x_0 & \geq & \underbar{$x$}_0,\\ + & x_0 & \leq & \bar{x}_0, \\ + & u(t) & \in & \mathcal{U},\\ + & w(t) & \in & \mathcal{W},\\ + & z(t_f) - M & \leq & 0, +\end{array} +$$ +where $\mathcal{T} = [t_0, t_f]$ is the fixed time horizon and $x : \mathcal{T} \mapsto \mathbb{R}^{n_x}$ are the differential states. Matrix- or vector-valued inequalities and equalities are to be understood componentwise. The first and second constraint denote the dynamical system and the sensitivities of the solution of the dynamical system with respect to the uncertain parameters, respectively, and are given in an implicit form. Here, $f_{\dot x}$ ($f_x$) denotes the partial derivative of $f$ with respect to $\dot x$ ($x$). The objective $\phi(F(t_f))$ of Bolza type is a suited objective function, e.g., the A-criterion $\phi(F(t_f)) = \frac{1}{n_p} \textrm{trace}(F^{-1}(t_f))$. The evolution of the symmetric FIM $F : \mathcal{T} \mapsto \mathbb{R}^{n_p \times n_p}$ is influenced by the measurement function $h: \mathbb{R}^{n_x} \mapsto \mathbb{R}^{n_h}$, the sensitivities +$G : \mathcal{T} \mapsto \mathbb{R}^{n_x \times n_p}$ and the sampling decisions $w(t) \in \{0,1\}^{n_h}$. The latter are the main optimization variables and represent the decision whether to measure at a given time point or not. In our direct approach, these variables are discretized, hence we write $w(t) \in \mathcal{W} := \{0,1\}^{N \times n_h}$, where $N$ is the number of discretization intervals on an equidistant time grid on $\mathcal{T}$. The sampling decisions are then accumulated in the variables $z$ and constrained by $M \in \mathbb{R}^{n_h}_{+}$. For now the controls $u \in \mathcal{U}$ are assumed to be given and fixed. The initial conditions may also be considered variables $x_0 \in \mathbb{R}^{n_x}$ with their lower and upper bounds $\underbar{x}_0$ and $\bar{x}_0$, respectively. + +# Example : Lotka-Volterra Equations + +The functionality in this package integrates into Julia's [`SciML`](https://sciml.ai/) ecosystem. The model is provided in symbolic form as an `ODESystem` using `ModelingToolkit.jl`[@ma2021modelingtoolkit] with additional frequency information for the observed and control variables. Both ODE and DAE system can be provided. `DynamicOED.jl` augments the given system symbolically with the its sensitivty equations and the dynamics of the FIM. The resulting system together with a sufficient information criterion defines an `OEDProblem`. Here, all sampling and control decisions are discretized in time and can be used to model additional constraints. At last, the `OEDProblem` can be transformed into an `OptimizationProblem` as a sufficient input to `Optimization.jl` [@vaibhav_kumar_dixit_2023_7738525]. Here, a variety of optimization solvers for nonlinear programming and mixed-integer nonlinear programming available as additional backends, e.g. `Ipopt` [@Waechter2006]. A simple example demonstrates the usage of `DynamicOED.jl` for the Lotka-Volterra system [@Sager2013]. + +\pagebreak + +```julia +using DynamicOED +using ModelingToolkit +using Optimization, OptimizationMOI, Ipopt + +@variables t +@variables x(t)=0.5 [description = "Biomass Prey"] y(t)=0.7 [description ="Biomass Predator"] +@variables u(t) [description = "Control"] +@parameters p[1:2]=[1.0; 1.0] [description = "Fixed Parameters", tunable = false] +@parameters p_est[1:2]=[1.0; 1.0] [description = "Tunable Parameters", tunable = true] +D = Differential(t) +@variables obs(t)[1:2] [description = "Observed", measurement_rate = 96] +obs = collect(obs) + +# Define the system +@named lotka_volterra = ODESystem( + [ + D(x) ~ p[1]*x - p_est[1]*x*y; + D(y) ~ -p[2]*y + p_est[2]*x*y + ], tspan = (0.0, 12.0), + observed = obs .~ [x; y] +) + +# Extend the system +@named oed_system = OEDSystem(lotka_volterra) + +# Define a OED problem +oed_problem = OEDProblem(structural_simplify(oed_system), DCriterion()) + +# Define constraints +optimization_variables = states(oed_problem) + +constraint_equations = [ + sum(optimization_variables.measurements.w₁) ≲ 32, + sum(optimization_variables.measurements.w₂) ≲ 32, +] + +@named constraint_system = ConstraintsSystem( + constraint_equations, optimization_variables, [] +) + +# Solve the OED problem using Optimization.jl +optimization_problem = OptimizationProblem( + oed_problem, AutoForwardDiff(), constraints = constraint_system, + integer_constraints = false +) + +optimal_design = solve(optimization_problem, Ipopt.Optimizer(); hessian_approximation="limited-memory") +``` + +\pagebreak + +![State trajectory and optimal sampling design for Lotka-Volterra system. \label{fig:lotka}](figures/Lotka.pdf) + + +The package comes with several plotting functionalities with which the solutions can be quickly visualized and analyzed. \autoref{fig:lotka} shows the solution of the example above including the state trajectory $x(t), y(t)$ and the sampling decisions $w$. + +More examples can be found at the documentation. + +# Acknowledgements + +The work was funded by the German Research Foundation DFG within the priority +program 2331 'Machine Learning in Chemical Engineering' under grants KI 417/9-1, SA +2016/3-1, SE 586/25-1 + +# References \ No newline at end of file