From b8636c40de2b8c066bf4ab6cf37c8f0a4b73bad1 Mon Sep 17 00:00:00 2001 From: Peter Stahlecker Date: Mon, 30 Dec 2024 10:50:53 +0100 Subject: [PATCH] car on a racecourse, limited acceleration --- .../car_on_racecourse_solution.npy | Bin 0 -> 48296 bytes examples-gallery/plot_car_on_racecourse.py | 415 ++++++++++++++++++ 2 files changed, 415 insertions(+) create mode 100644 examples-gallery/car_on_racecourse_solution.npy create mode 100644 examples-gallery/plot_car_on_racecourse.py diff --git a/examples-gallery/car_on_racecourse_solution.npy b/examples-gallery/car_on_racecourse_solution.npy new file mode 100644 index 0000000000000000000000000000000000000000..2b2930385ec413624bc9a5ce69ca640978232c8a GIT binary patch literal 48296 zcmbSy`9D=(^tXAq=2^yT2u0G2(q0XQsAw=$l2S?`M2VDy5<)3HQDiDIM5avpkU1pt zoY^&8bLQc>-{<-M1<(0qz4kii?0wGJ`@GlQ>;2wuP8**xx$^_QPFL-(+;zHs?WMSVdh_S~ZAl`3 z#a(c2yCgxJ;+@_*AweWWQUqF4BnZ=8r@w38l^`51U9axnDM83#%RcofaU$P|>j7Pg zIN@VbqD8h5CzQjAcPGe+6Hymn~|&ieVU(`F^Lek7{O0mvgpgD*Wo7?eL51pP4f{8&Ux4DBKU~GN7nI& zb@+&f8*&lKGrUBK^nF{cNM2$g)WiItE-x|edaJu?h=+)cc=m}+i-$;*j4`PRL2W!sfIf&?>(9nb5NW|8zWItz4 z65%-bNGR5Zo$yx-{M22+Mlf98QnO8wjp!{*S$p@KmAHB)>fY{A7Q(Ay&m`+f7UJFh z;#RqD%*01_VaGQ@%mlp*I)D8U6QQV+YZ^PqNDS;k@zOZNf-ZvfmrB$pF4^-aIv&<%HH&K92~{I z@G@=<*H-B`EJv(jZ(kYZ$gmZ>=Va^Aps-~e_B(B5GGYmj$JeyA#x3GuCjMgUj0Ic; zm(h>Xe|Wh=2_`$|aXaUaVTGkR+$xpd*i*y?txj{xQ|ZG#C8V1X5QW{_;DKN zD08~We4WCER@p^?TPN{##=fq&qzT+=zccr+`Z$KR#eJ?-qxg2+#hKeC6kOw*?A+fn zj9-;*aPa94VG&{8Fp19tSkPHT=gndt-YaxDPeG#>FFbX3+27uU&-_RZzqPd;e>>T8 z&)BRPzdl!3a_rt8e9O;~^~q#4e%nl^qgq{pPh8m{|KVRI9=wiXK#YR`qdO zbsv&6^8AhuL3T_`1N70tC6ML6GkE5gBXEi%_Orm@)E#)WX)5!dj z^8Mhr8C1{j`jT;K7AahzdzOE29`SA+uzUaVAF4_H@t;@c0xAq1e{Z6?gz`^UFf4m5 zqdOjk$4<7dptHMr(s%D#LkCiZ17E#dM>5_vbVr9aP=k$+e;An#0@!Z1ivOmAu%MM6 zQjh6jRKA%fWCsJ->-|=>>Sln(`+eTNfs7FLsny&sz9B!wbg)4xN+N;e*)q18D6(KHzq0n>3{I z!IK^1H${OT79MsyIQ^O*YJA6~f<=Iw;eYb+9}yUy zu;X{(6NRu6r+>aiqL4{ z$J|mNy>`MZ&{_&&UoeW~6-mJdHl~|h!qOnVi(aP2UK*G=7+zm3mxev}760aOlVQU# z>CB7<84?fgEGWKChC2T2`(i<4a5|{*eJqy@W>z1=Wc$f5w+3gvlVreV(MBm~w+uv7 z=Gi1(kO9xC<$#4JGO$%#hH@uf29#f1#+y|#aEER~IFTv??iu-^S9oQCTqkEOzFQV9 z@BO;1`K&BVynfu}<0K1gp6WFUA+m4|-}XySm4y#;=l@Fok_ER??x{P&vOtQ`qYtK& z1Ld#x^L)3+!ByUqdYpzL^=#H<`}yX_B)S(Jm5LI=caSLGli?VD5hh8*zzzj+n>@B9k%TB1%%UM0h? zABr~vEXfcPb$ImV6*3sPAKjmPnGCZG4!4XikzvX=cjby189o<${8xFN4DlKfF9MCo zAdocjCdq&dYm_`|2^})1^2L7dL1d^NEci5`P6kD_mA|ly3}X+OlbV&t@cr&eL7X%h zyt}1Od=?_Ze_||4%_K4``?u~tyCDrT9YQR-hor&AG348!YH29`llS?0iZm!HAomZi zq#;L&JF3%88Z?{GK|up)h;ry+(Nd6xx@|Q@2AfjwmAF#7w?zs>Q}02C z>uGvdH~ z;2^IIy*OM=iq1+)5d+(`atS|EF;G@>>nq&Y0{2A2PTR+B0ZRtu=6dZd;OB5kX3wxF zm^o}?9P}22h+kBN6&XEF{f`OhiEDjsEHTlfs~M*rrdwTNqqlItgxw z2t!JX{A6O15Qv|xX#aOe2pIlpX%zn!gr1L9K5bTlpzZTKuwhvMx+?qVM?(a_>`PGi zK@|bm{>-1wl4BB-r|Rwlea2j+4OOgBXhmtG8_E3JlNT=j}6Xs`-gf6uz?}lT2=W1EA*TT)JXly3Q~c5 z%zNrs;rP?Rw(}*dU}#yc$ydM%!R6uExAIxR^J3TCiv%m|Pd9GsuV#feHQ#;dJ6Yks z@U=^s)2y)QZm+P&!3J-lwEde^*dSQUD-NAw16NtrJ-eQ;!M*w?yAr>#fo@WtLEH!% z6vP^5NXxT>+P=zEQEPTcmJrqrNMQ%w;*K2;R@h;m^Q*SzlO%w2*7CSG5@<*!{`*VM z0g5uZ49~7}z}QH*p*po&Th1aKIizF{> z0jg>^$AywD5K|$#d1{*&nA81Cp9mC#e<$|MDNT#PfYH^T%ZB2hy`uNfEJYlexc!C@ zn*@k7u4jFZ|aS6*8yrtk~RF0NZr4)SUifG%-E)Acy&L0cal?Ge6C$G^nY0$T|*t5Gp8p51| zn{8&LfgSJJHOo&1r9?SKBV{tQrEt7Fb(9SLpI!frnv&sxz-KRcTQZ1>Fg;p(Oopc` zcMdoP({Q#lOvo#m3^ECdUaDzi*euAJb15JLFMr9nX(buh@&gwd8_3YC-kmhsNrv*7 zJxSpMWJvaS^{IE14AL8^I>nP@xZd7op*c&2Q*WQ-8qCw$c}}7-t-h+S0i)Mx^|u7SxVoR#UW|_0!RQJypm4X8 zn0zv9d-HiwE0yLyV0zstiuPQ6RBfQe;qJlzmOp!t;gEp;NaAfWjBCPE9TPH~6;nL+ z7-;b&b1mb{HZqX>Y_Asc(&C1#ZgT`x8mQ|!)TwM~5YzbOBqu0 za|UTpK5fL5Oh`dj`B$gIkEEbg|M{%YJ}F3ONcp&TToRl+peFL8By5xyF`qG%gnrI{ zyys~+zEh60_T`HNcdICn+@Dumy=nJkLKFqb6cycdUngj9dGgW_NU^~*e4Pn)4!kG|l4jwu10y2f zDk1MjaT0-~%1cko)`UUJ-a#bTM;KP9ZnqkFgrQp`X8)%MA)p)>z9F+s2-NtR)17h! z;lR5_$+Q!K;H&GsZ?a1Ow(kG=QO8LDoEJ3o>(~Wge>{un`&fQBtfxdRLi}Lj-o|a% z!v`Wu;n&Td@PXfSCflMcAAFnrV!We@7t~%8KRBFtA?A@&NTeh$>`ve0uMgt^Av1ab z$}%^k(*>V?ag`gS7xQ(#*Kt9cpmn$P5iY>~C+|i6;DiyAY0(o0IN?)XVZ-G-4#?jr zF;j7r1CE6jKXR-mff>Jd_7N)*5TjiQlB?`+%5gfA62=Y*H!c|ZXs|Bzigmv zWon!Ek_}?$R~v;j*}&r8?$#PAD~K!R?!ob_AmO&ySAUfiKD5Ry49c*AiFT)3&kzfI zdM=?8ox%eBu4NsZ4_P2w@`eVB9t(6oQ!cg>V}VsA$!hs&X1F-bo5fMZ3}*LFN;AeY zgR|L>U*pf1LI2!?vV!~4KuDh%8h;)-ExL;tZaEWyL1N5cXqr}N#=;D}uS*XZ z&NG48xg`EmgG|snoA!Q33lmhE^@jhcVgg%R{kw9Q34Ro=oZXVm1Y;e62fwE>!Oc#a zhjU+;-~)Hmj>b49cqVh(e_sp}cuecv>xpK9c)_FR9clL6|Jt}-ee*3IlSuZ*H{Cbx zlW2m3uD*4hM0Z#$pX_#?M5(7_tRCA>q8~+T(|$K;eu*iWYgZ;wQlyj>zsV$$c-VjI zv%w_lXqoR)gGnSx`p{IoXA+f(YkzQ1m_&~sNPC!wPNLx2*t#J0Nu<)E_TtbA6`4;S z>GK?;q7KV|M?M`?! z;D1%}2^C!`;k&EiL`8Rv&%8WsO+_26LHwzgsA$Lcqz@ivsi?`vUZGi^isbfDi^a96 zh+}e`)8q&hO;|gBE<#lF1DVdI02RqhUdyvTPDPBa_`wMSDvCc^J^0F)iY{nBcYAL| zMQNFD_BOask@trL{{&k5(IeZ2Uw@&s|2>H$MSI?S@sy{XGgM@oV>e+#oZhZ#bHGQ6V&teui zwGbD5N@tO!=JpI8)j3qA$>z}f-yFI_7iU+*GLMdC`mHiJ&m%Xh58jvi=8;jVD0bI_2lXwB=xsx) zT|6uy**EV=35iSSmNZjKDc>@x9Ig8D&vhC7lbT3T`MZo1QFYYi<145~$b!NB%L>{k zHuQ7hSVd*5YzGe6ts)hD_dA8Pt4PY-QIDdshK8ObzDs$zhD=2dFw%{zAs)S>U5!W9 zQBSq$SE-0~l%rG2+B&0x}3 z-+l5mJ(yhK)r4evXiK`6`lgf~N=vTz)OXTDX3Z1Z=^1)ZCn=qoVqpNc!f2feF$Va) z@1dp>E#JVFt1>=%m;rdlk3Lp3WB`W2EqCQC7$7={>%ZnZ3=mVTYrfBw0kR9b4r@GP zfWs9}7h3`tfZ%X@%J`N6d`9j`eve=P|4an!(G1Y0R}^?Fjse>F544_5V1Nzp5e7e+ zO;B-^TBZ3D=2-&)Jl_Y1>F!nGOwPzj7Gs~zoyScV8)6Yf3FuJGf)Qv@SQ>!%ZKHDN<5pNv*m%M-m zS9u;Zo?SrMxx83->-lJX;EjpGO}j7fsa8 z&!ahab#X(cd9=}dAY>_P4)uz){!z1;Lu#QWwX2eIsDi|8dbw#9`Jd+5eIZ~LMed0t z&uPyhob}&^J`HFn^G_;20I@xf|3oEKQ*N>nih}-zE^rOh3}d zZUS9<^qTdR>I8}oc{yM@KaNTr=a2m^8b^D6k>^ET(ei*xGoQ__jiZ|%%armDj-&1` z4?Jzf#!+$APMHrIW5`tVziCg(82V-$S4ZN-MV;mhrN6T)0qR@|{7a}wCC9>maIPT%-q`^4SkrvVQLQ4)+k>7VsqvL@KN{O$g5gM9TASS0#@naJmVWS!3 z)%L7mBYg(#4z1`oBRGqK$|C>G-=9Tc5nI%o+Gmlo=Ki3bLvv_O`4$Ib^c=d$(0xXY zaURX@dhd1KdLFesiQ9I$Y9770m5&an{6oJ5HPnB<`iISvdq z=U35_OZJYOpH~s-tgh;h#Z^?9$9yJNXAOlG3}`I9Swr?!pNwz*T|+%-u1j-D>*&CF z>YD(Eb;SE~qDMM=9Yu$qSc_aewPkUY%5Xa1zv6fK<0m?3x=80Dn@C5SZy!@PN~43gp3;*FnRM_$S^QjL4ju3@ zQ(T9C&;hl>#9=;<4p#ryiu&uWlCD*wpPo4_#@njV7r#IQfz&D#^)Iv1U$hD-tDjX% zajir>!&fTOdn(ZRKL@fl&sHGsEXk2Gjpb;Qr1HzjtQ?&PlhohYS%#Qqjc(sHDMJD& zU)S_1N|Dx*f~tZ_DOxrRP^^0P6Gd9xgwHi4=o;gRhkuz$5Pq5z*0Q}A9h%EJv*R#< z(iwf~TRn@=_k-yx-ZllOTI3XKiB=v;(=Oca!kCNxmNX^G6=k7%q3Oi%yWi1~^X3jV zqUp%%uYd=KZ3@a;u#PhD_<~e#l(p+8$07yeia{ppNYv=OXksSv8c{ncpM?Z?piM*e zXHg9f_=(Z`OX60+xZ%q7jRJuW_-M+Pl`CrTSWqisi{<`gy!FtJ#+K%E{K4scK=$t} zd^fFaZ#r)tE>A1SSN>Fl*E(~5d@L)*Z@$`DaoCpP^2ARi11BmlQo=bv5CrtJPkA3vmtIHbKFR%mu#pa(bODE>=7`tB*~DE8c@%Ggb_M^9Zd zkz^XhZ#qc5*{c-X%|Q*Sq*8D&lXZf9I|b_*sehy|qG0kR+Z4$N3ckZ^_tVRPg6CRE z{?dmjSZ+^MR39S+TOOSssxBVESg?ti@?ZpCdr-AiM}hVpN|#@()C}XK*;JCW?J!2( z!F(-@!}yQh9dpCbAzUJeIJ+c=u;uH5v=h;TxMt&2^rHA6u4M@(v4;%cD=B|Inb8g4 z<&Zs5r)>JM(_~{vZcQI%^Gi(^R_?=pmmIk6dH3R`$K}z-GXCP>z{ZrBvpx8vmZ|mE zsV=Pj;iaUDS0{}xvphV+)q%w(6lOijo7EZ zL}%`91Gdy9yJUCPW6kM!z6swt?B-VYiF*1M{`z9!qCHO~KF`Bqc&f4#Kb+;eZ2hnp zhaA6nn?bt(&p7##?0vKFgD3IBLbfT`&?9^(M?VUmR8wzVGrohR^tE5*)rOTqH@#dXcBq3cM z+EZAMs)NbPehm%C`2BsJ&iN)(Tfz~|Za{Ol&w)Q?&fboA1N2GPas=8`854j~2Ka(9x*FdCJYI8p6C zjIv7@C7x9bqh|4k8Hd+~QK3fsrl0BviZu*xh`%y|*bLq?jj6z#yF?4y@P!RXIr(L!F6f|Y|pB_%6pxT3Pl(>FU&}oQ2 z+|o)x<#O7yY!nKr91Fq4a};D3dhe>*D$Q>rrS88$LDe21gH`Jk^zQl%gN`K%k_%BO zxH&~ZhRVk{dU`47@Y$LZc{LQ|W4^T~IfH_9O3FKYUQ-aGXt=MUEd?nH{xEY#6tuUq zB>6Tc?ff;lQ}YcY=qmY*o${Lz#Cy1@SL)OV+E-&MH9$9lRGEPf~p;8Bz*BQ(^MN0_{urEk=%-U zb>#LK-D*LD0>Wq1WSY^>ukVxFTN_dM?zeF^0S!pw>)QO(;Xmj^^T_x_&u_F|!Ijc~ zuof-W56MM_RwE-xs(1~h0*O*In^@G#P_5a!t))*(5WOa2Qprng2jU>S05}g^eVFV)+elgc1MkOdJ1-u zU?90mXJK8jn%4 zaoSvcz{S8n*f%G(wKt#Sta z(AtADcFuW32lwJ{_Cb;65BjjR=URgSb3e8_y6640em~YRcW~x$?Z>Wn%!AXT`|&R2 zgxky&{rC#ul(T!HAG=Jv+1Cf zwT@uY+&9_?9l`75A#(*TT71YJ*UpipV8K1w#)*3=n7>N*58W{eeiiE%_1`H9eqg}J zoN%6kr>T)r8>TeYNIaJtNrLx=-n1?FepP{HAg#Z3HjR zH+q@hAHj?IH-0_gAHfkGvR;NRwD=Mqc3kn@5MG$lnmXh%h+iCu6jc`-!1tnpCSE@3 z!?P>4U!L*L71UeRT=cn_k7^#?r*vHSh5Q0P`$pVvLg%-Tbb9x7qOUU^%IW`l5tE^S zrs~N-#5P=z@bu#dx+19hQet8ZJrrQAV%R&09@4T9FY>0*aG>dJ0aMyKOG#>n(9|3{ zkW-kn^7J3_ZxY)t!Lx{#7q2||@pcI%o7-=agjSH{`(~<5z$#jhiJNs_T|?i($NwC? zwt*%x<=ZxhP4xBiNX42z9gO)Fuw?9^2N6!I)Z|v$dU;5cZ{`659Jsj7Xe`JGHWy6W zm693Ztc&j!)}u@iFs1YTY8?|8@Kj34Ut)&4`^=8A4lvX5pX=wPtyrML*5@`uHw%Q> z{H;82hQ=GQtQ(49R#=gBBMy;7o0~r_jO%{}R79$0@Ue?yC1= z(+BL($KUDoJf9uTMIAsNY5C#!t&bZ%Xp#Um#oamckOaFf-0jUk z?#J?U5|=n(!mHMy=^`f_i&YrEWyT3Dybl^1&T&GX8{0bFDNaCT{rt7ZXjZ9M;}tC* zjyzg^2W{ho?~(G(1zR|w#>+9shm8~5itq1rnB@TBHU9XnE)JmEA?D3q&Hq`sJs1Z7K;ToxFAC3m6yo|?~B{qA|2RZdfoN~*LK>t2^6}?)5i+hqle|byk&)} z3VVHqbZPnb@5W)Ji!31P{2Q8CEVl+NGsD&&7mSpnnW4RJgwOFf zGZbt4_FNujf^=M{G~vSp*@0sUf^tl7V))ijMIj?_M7BUJgP$vO?MgD$HwWQ z-*msju9pmevvzGcOJ;=ON<)G3sf?ig^ZA{`y-YBZb2W8dn2Gnf;~uxG)YG}?bMF#6=*FHRU8 zKH{(Jf$A8p@I1u@ z-_{xXSWLJe%Gc+~R&y>e`!B1y`3e_|8=ImUOD>?#q<`i?D`Oog!H7#-p!qpUMb?-L zu3NV6+i{W$cJTjs#zo_oY#amRgZFU3%E#A^KjgUJ{o3hwUVL0I^6;9~<{BqhG~M=4 zpgmu(aBt#F4JXKi-<9cMn|FHR6ENp{jO$5**;MW@l?OLV}_% zR=-Z9ksx-pZRb~K63k|b+Fd@0}0W%LcL6wvf_j`T$f;TzTNZ3aO?MuM0(KxT)`8IsBCc z82Lqy7VKkzHkIveY z%rxBnmHYY?3s^ip$=tJ#6%Owx@#tu0h5o*5# zT!n_O@1HZ@NFhP!hL=m4JO_-*HN7qk;ef5zMcOxL`U_1Z4V6mTxT%@QdLHqX6PizN zdBxJi3Em^q-AT+`U>|gqjJ9*Z220bQAX@)9e>JFlV!#E1x3pZ@^tj-8wc)150WSD# zAIorX3m5oJ1YFmf0npR1}#orN?&Q=`Y z`%|jQl7|27J-TK)ACaIZP2V|?7B4QWZ5=prj~4eH$LRMCutD13!x@&-Z18edsMo#> zRKuX&pJljUHk#UPB!o z%{G$BE6D7sYhB^`68idl_WtCFMP&8wNYRhre<%~IlrA>RA-5*>(W^YOC{W>tLxtWn zI+!G`_HJ_mIiAf8ceyi)vIXfvR9lB=c|KVdtz-SD(5U+S>r35;F>zbQm8~skuHX98 z{w$i7Lgm@j{oRF#j@aT6lNQ`%aIo^vay!s`j5Wi=LqPOakerX?*gg=<(v-PtAX&<1w0+sL-h+r%G7m$tb@ z(GkT0%C%$G^n_CH$y%1(421u2xx<=tj6~c{A}Zz&Bk^Cvi)~H`Oa#Z*v9`j;%*658 z)vr%3vJkf>PkFu9U?o!O|NU+fXCsO}J-pCP$4=b0D5a1n?1Z)%{n1zTBqFuRvV`>q z2l4B!BJ~1IKYb(2*hA5ii!gmR^H|e_o5;zyu;qXl53%XX_3%b14^ba2qQ84PFHwG5 zG||qUmxxm9>8;D=CD!7bZ+I~C5jXC5gzvn-N9^t{WtLCoBlg%@Pu-vABN%u5+ICwpB**5icu_FVP?4BjS6%>{}e*B__<+ z>=ZqDi65TK!&$&d^nWwh=h4bR ze2UljadehM_=|Y?exd2UsXgy4vJ}{fF)yueTEIp;4SB!GaDkP$K%(juIi^mB@Dz<0c8{YW_p7E=P`bb5jx_Hg`ZO8%qHHR z^eiSva09C;xtkgvSi{*(a;H)+tYAy7TPAhxOPIg$vMc_ufP3wiO9Svc_HoUTIodai z!N!(4&MPZkdxt{=t=y3T>OXzPr5w=M>c^d3xpzMyJ! zuobiSdNFVH|Bd&B7dfH?0xudpJ#=-$ACnpQ@1NDnL+V;Fhsrj7A@UK6M0wk0^nGbk zQmCa11+0ENr6JjmZh`{cR@(pP6Te;G({~h2dxWX8@1deuwIFc1JB_4<-0P1YpGBpI z#nqKW=Ml;KJImJbe@NZ^ew2LSBAVeey1?;j87=R>HduOj6>ZMiUp}w0j+jpA+XgLf zpc1bG=aa6{L5Cr^uTqyDI8=7}N~$nG_ZO$r(L9WBjwRMYf0hwS6@&CdnrQ2aeO))# zGnql@o`T70C<|mM8`RM|v%26kvO$E?2^OZEG+kJZEdzx_g7{74jD#T)yiB{$ zxsx^@cAN`$a0ufBjYsr5jO@9fQx=;fALE89spN?BB0L~x)IG*E$^%;4y_NLMdExY+ z)s@;3UJy^HIy0fj2a2R?QPKW( zT;zOGc1r+y!lDPp%xOBWE4e-f(gI+9v{~pw3qLq>8Q)=h&JPOtobN?<^TVQW#mA5~ zKG^+R+*#6{4=!HF*R~Yn1526fB^%ki@UA3l@ZA|+$YkFLvTWvoVv=32IWrG%ePrDk zy@wla&puzjN1IPr-%S6S;Liz<^xEDZ$>M<7hVe6^-6ZfamuROjkl~KSR zx%R6r8=NxG-jiz03dfqhtmQvrf%`HG5yxpd@=Tst`%gcZph|T0EOR3xR2}SCc}1It zwrZ?@ILAoiRqvkoU@b;Z(@BdS=~1Htkvpf8UIlI-yWj6>4RY7en8+>(mi`r#(Dh96 zH)$C?@C;$h-o1$URR3sIUHFGu)wW{O_T6%$X|)Umjuh1!DkP(- zu4o7C&uKVtVDRO1cR8Nm?)$KPv;jYJmnuDP)qz)SlM-G}_Tq7C&bhr`2l2uEdL|EU zQZSqD;D7pZ6Zl|lr{JreDO|De`1PBb87x1SX5^7LhpYXyf|q^%;pp~HE6Z0Gv1xMX z)LpG*oK4-wYbd&k-G};`Ud*mxwt`DNyfqt`qvO|E+p~0prZPkCZzm!$~05ogamnxQJQlysL3TT*P?> zg+E1Q+yuk!(m>@n9wKL=H1jm@5~GA_?fDd5;!f=CU+toNM0J4wuIpZWgz%WP{v( zOy4&U%)7`(NSJ8G8LjXVxvr1cUcKWb7)k`VV|Vfr^BPBWUiE@TZyfe6o3dz{` zj3zF^K14?e&(h+~jm|q4xi|?M8|fQ|w{j4NU!U8ia*#x5*UD~C4A}```R?-vud@*x zZESfXuB?RM+>wHmS1d$1{#kuIftfIwJI^wsF!PL^FO+2w2U#(U{+kYd~`P5Qp4KqLUco{$_!b+OS{=JG7kO2t%p*!8qLgM!(gq@4{5V_^x=KooChHU6uG-Qa3bEDbZ4m*^SacnAJ%Ec#xKHt zmXxzkE;)3M#@brgvO9r?T+J z_?z*^et*Ra5ntMlKZ(TKxGu5!lD#lp^jnKS2C`Q5*`9+Hmb92c$cwA&@x+d*L7jPZ zL0H|=xs(1(7(Nyd%Ncj(BObh*C#m%)1~0vwvDlFkk5hNQ?JvS#a02&drhso>akh_n zPR@f=eB0C2PRREgp5S7B^feYu_YoN93S_@^I{(oBzy0n2Zj@3Mssh_x$EzR#}34xvXlo(g(trvTrVnSmXv|U zI97Dp_vNCd+QRz{CPk?89jR$Bumr`0g*=!%QH~CUwiwXeszzncpA0{oszb8fQ4ALI ze-OK`W~WU`6Vl8yf*0|vNb^G6HiqsFl)oqc-xk9jj(q)HV%+PcS$IosbnL*qU-24F|v*^pq8&ky!vxv(hDf`y;S>zZQPEuEz zMa|_lC+*3!^4F^q-vwsTT5w?dC5BmKa%}E-`~*!$ASSWg(lCQel=4oU%bG#2Tn=xE z3!OpP^+IWv9cK^&N!F!&6RN8?r_kqsiFPN&DfFD*>4s+VBy#kR){&;| z*GeLL>%2^&BHYZv+AL2+!CZVc`y(e%{lcN9QQirpXRene`fMCMFPm8Fn;k<+Lxz^t zmSadH>e)s7YZNj4Y>;H!H;Qy@j7ttiQV=>=apO8Y1$hi5`E9*Ef>sjM{-~7=BQ?>< zqyNbdBe%-!hl0I^5Rq*tJYqJsudohLqgW1t^Zf0z`i+TV}b z->OJ-`}U#b#MhVR8hVk%#gnAfZM{gERYp6buLsE~Zmnd!*^RPyxWccPUumpV5s zI#8Yb-qR~}ZOD%7aVKtnD?0MXaccM5X5=8slJj%45q)>5ll0otfV7%_-L7`5N1I2s zh5ys3Md~~GOgCRup{YcTF8remd6&0o>6R6v;K($XvdBlesw`IwCNfZOX>_u$c^nFV zE}eU6(-A4=fuwPd(g>%e|P2-n6KGmanGkEXUzkdnAS*#MXJ6lU>7JF)1+&;Q<7E9{i@^zruicitEe74VG z<}t1!nR4cgmIa471o}smTjdXRu1{$B5MrZ%H)8-l&+yZs*-^ z*+x!dWxEpv(s!pZ8|AsRJ(|Y)SJj5>xu)^uzDsv2Go~=}T9doH?iBWFQ*}JkHi@5q zGfIDFJ&8jqVwP2Bsdy=Kg>39j#jcZf%|;AV%;*<&&-nQS-n()8j~i(MuZ$Hc4*HMd zOVTge(mBSlkyeJsy5AVyed2L?C(9VFRXls=vezh1=*!t2yhg#-m~OrFaG_w5!ljj%*VZW zxonl8rSLCi{!6;4^RNfsKCsVvNwpiN1tdl|_H^QH)j>02&pWVH$Gd&*O6|C|dgowY zaVvIhHL7@bwgr!QZs?7+H{nj#@*t9FBc9xP*-y3V5ANQVBW1s{9&;Y@?v3TE#f$x` z*FS%$!f#{FMit7G!vszcJ0af0`cn(Z#=9ztW~q+`OFHvD;H9@}5Q^yA(|~O!r`2 zW~cX>PN1j= z<0vmekpDeRC+GioT<`DUF(i=i)aOLfD2jh7HPgj2if*Ob+z7c# z!(o}-A04Vj(4Qk&Ngtzl=+V~=ATd-JTyAdUq`k|EcLFgdP_jWL95Yq1TL66!N*_EVT;5_j`)3Ut>RCIRw7?674 zFn{mkt%s_dy# zmP_DO?Z+W0+d@!%7Sw<7M-Ds<(qbegr$hK3O>O6@M2N~b5p&Kq1Q<|BI^X(<{PUH* zEBq_&>e6`!SA07&C?{)vIAl7NQ-pstNplzdrSP5ltk7jzO(~5RqMbibO>sG`|CcSg zk}{EG(4S!Smm<2n|A4j0UrI{#LAFWTN=mzv`});W)s#T#j)sS)YAI>nox8UO)lsH5 zYbbu}ucz!X?{K+o&_J=J2u<!9r2XZA{*(n0BsQYhm6*Fm}c?ew0vc=f3K7OUCQLFrQ~ zP3O_+q&R(j6`CX7K{;4$6!B}Li}IN5%y;R84$7sWW)t7#4ho&$Rpgy#7e&)DP1R7O zn_?$Py7}yWH$_(Rg>4YMn{rFXE;`%<`=V0nNE?5tl#Cvyl$moql-}~KrG`sAl#NRJ z*rC*3N(#^E))*f=hpiRe)_$j-aDzn`N+E~muj zyvKByV1MLI|4D}`lP5P@Sf`*`Y5z*L#S}>0*10S3a|(ofX&WXxrvYA%**tzY4bN`Z zPIYumgUE}i#VXAiFm%k{rZ%)m^Q`5km3<#w@T<2^3Jg9tLau~^+haXGvF70pU;UIr-bb`%15Xws$=2YfkwTW?SEhB!u)Onw9 z=^Vtpt8P)p^Yx8Wy|<&4=U{f++i-Gh7Mwy1FL-^Jg{)`=Ls^4aknH~ED6ljGoA=M0 zWQdyqKGn@hGACvrJ9=O!XKEUJj$KIl6grJ{Fh(w2^=UYIk1xclX9_eU;(z?}or32R zO0EsEQy~4Jeoy8NIcmMyKm zr%{#?zW?Q98KamvRiBHI9iwnXy4ych!+eMeCBBEcV-$fSdlU>U#wfzPl;HFyV-$HC z(ea$(G0KznxAgKIw zlHy@@-L!aolCsPHXOe|9owDm{`!8M-I^|)TjI+E4o#IH<^0<;jr&MTb2pTofDZNX{ zqEEP{C^B!TBUOe|l)WDxrR@rtqNq!IX3cg@QJT)PX^rE0QT0mnp#S}8$_?}P#q$Nz zlzfii-yZBUl=f1M0Z*eD%Gswnl{Z6YC{?5B-)h=sC~e;a-NkmyQg~n7huL$pl;%q+ z_qzgTDN}{6w!USvl+mx<|xkWHhxDx%ux(NR;BYi72ejPd(wBc{PbU8hs@ z?PMOSe5F%-;&;D!Ye~m?Sr1`;3!SniG;ML&XOdD++$_O8HbI&G`1!(va}yMf>bjJw z@^Q*S5tq_Q)p1H(P2?{M-fz7sWh8>58Cm&64Api%s)4sD4!OrsbkWV*>TjZ(as zAt2silrl{}s&68M*X`7f_0|yz%P4PY-{%p^1ug#4op>L==VjmE8azUgIi4czcLVS1 ztCc!aw!;)Lo2(S?=Y13vo~_!x6;+fY;g{;&(`o^}f1M(h^+Or;P46qV5s0@qFg>(TI0czB99VQf1J1Le82nuDW|yv1kn)659OCNmqg zC*d_Qbc_9&5wJ3_6!e`RhF8`LS1tq$!&W^~yrIxA(3D$QMiPc#!}-|Yu*47;oia-` zjUEIYfqSMg+XkVu?+Sag&j8#UKXoB(t{=LLe;7vH=!d9VE@?Nq`XJ}pkk`CL9|&oL zY7ho{;VNIO(A`J9kS_ht=Z$bLWG5a?{9V@r$n6*Vr&m2N_n_YnPd+F5@8W`J@ZdSHGCqi+4fUt`J&xb|(xm1~-1^>x9rTQIAlr z4*0fd$x84}8z|oYwYy2N1+G4pqxybognSnNU2D(&f$R0-lxO6BkRd_mfBT~mLSOG# z4jStKBfHq-nzVkHa`@tU`t&#)Z_xM1`aT1H8C{!G@)tqGdM;}%c@2UUf~z$r7!W@{ zfti6}8?w*tjMRNM3Z&``vx2*$(CwZQN>CmJ9v(L_W7-HOTUfxGpLA2#i1u z?aIv;p~JA{?hg4H=3%J)%W!V@r6D*dBX;L5=E_J;lOhYo2H-uGH+BmRfLpZo;P0dT z5M$TkJa@Sd_H7CaoptDi4(*?dqMki)!D}{W%$o}Jdel$;?A_3DZcCmmQztAP<%&&N zXamjjT2^6WEs%X{zCX0T3HTROgLd{cz{0J=iXUn9aCs#C!$-zC2yI}R|97eyEcGb4{);y1Z zh{UtSq@qZ;@pD&>gkCfZ%g6C-HHn3)w7s_`|9pf+MR(oBEdt5L((d5vC9vZ6 z`H_~n3?&cLZI3suz^T5GoANks7LYux7qh$uj@05-`o=o25qLXXS2o~+uiaTAoGWQu zynMoSgaO@|FM6|amJt=Sy^Y8zU_|zXmt>uUnUMSYE62p$nUIc<(yJeJOlW8Dt3fUi zW>j}lNj5~A{j8SQZ~HrYAKj7)y+EaeknLCoXt+a-)yP)iEW4#$W1@%!N( zl{gkO;z+%lR>FeH#bf>nbz@&$n2meQfBYL~t~G69BcMrHUk+7X0^-isO-~TUeveQw zjx0$6iZ6T5Q@e+N%s=Rgo68Z<*=(C*WwHb`_~S*`BN+l(*DR@B!oEGNuMx4yy9p@Z z*~^_IK?2%(aO~`h?F95;hk^Ge9s;^*A|=qeg@8gAwzCCoCLm@hy<2+tv9Jb*fE@JZ;K0|8wWKO(m$ z82|r0P9L3q643bgw=wz<0dXA`zCVj|*Y!L;iS=h#5s#6nT3IkFiVRN)t;9aaCg+>Y zPO?Pw!8IwJ-IIvsdwwUT{3D{{@}6Tyh4DNy^2+h9%WP)y%4F#vZU?7dMq5Qe*A?7;F^|0+t4gJD~ z%u_c4NS@etb6N9#ye{^?7zlP{64{Vob^hZIVwc>9r0SP>McjX;^{oaNy!808!D0o#k(A}H`IVW*@RyQ&u?XrK@ zp6D~95I0p9#X=@DBs>}5v7ZU8JmPx08qJ7QW`?zOxET@gj#2lVI|HKh%)dF_xB(oz zCj^?(Fc0fbmTrFF8rW1dFN)t^g*)}ckNeD4;K<4L)j6GI&`bXJ+(&f@#QWS0`xO`A znv~^Jvg!gP8p7?NWAmV(+(SBhX$~BeueFJL&w}c!+Iu`1GeDQ#;Z#eX29Gv_9=hr@ zSWmtWTnm~4GkxQ0)=P8<;V?ISbDIuAaUbu!9GC^^+%8O~cL{tA!6Icv{^Q2~=Dvv7z?I$*Rm3)?9B(y3e^PCBN$?kX{dxHsmYRTDMdKmxwLW5&D z+nG?;q}!FCNk&xD+4TEg2_t$gb+vXQgb}^saEb}GW<*@XySF!~Ga>qJ5V}JOZr;2~yiw*eoY59%#r41m_7PgM7Zotw;{M8q`HsBmXN8~QH4S>(0gLjtK zVVhzVFYD|&T>iiHe|6T~RpssOYDHkVdBj1NR0IP|-mbZYg>YaWG1Bv5A!w$ak^ERg z0iFKP6ST|%&=9lUYwDU0CBsqo<#y)b&mrfmL3}QFO)%BT@5liM%_Nq8-dPZt@*S0r z{D3qamh+p8Ghl6_e^{tI9a2^DIMPn0Vcqt@ZHJL>kY^ph{LkwPTx5O7ZS9^6`>M0_ zI!8Xiy;L2ixb7s-kj?V#x|al(x19YYmXicQynoq`{`v&D)^{U0Gn26nNe{0L{R+E8 zw9m(*G*EO9D;FsH4o^NAY3%FD1lFo71?sz8AokT*>{u&+&-6~d0%9?!XI(M~9xa78 z#k@0i(!Zd6C+&4x*`U?28-R_oG3CQ}6Qm|p?beRM zIV;NI%rBP?*sS_s_?=TXcy=zX^6B-!u@&dbOt_wp(=AL-_3sCoQysgg(jX9@lxqI$ z8-nR_8Hw93M&MHDlDnZS4LEGWyYExRAXrAe=a<$5xTi*3^v|AzOSA5@HJlIDs#o!- zRK+ctV@uS*|VPeOoPUwjLYRRGdLeQI{LK{=fKnSS;wBw!kO32wL{0};MLiq zuRM6?;mxy$?u`@ku#2Fl;a#@?==4n?GOoKH>yEWme_H|+D|tr-a~S00e*b;=aRnku z-9u5_Dm<+8o!lO`23KYK4QFE4Ax9k&Cu26io@y5xYsi4I7h2bOs~Awxft_zAbQw|m zJ);{1#f&H>xr~Qbl?i=5-zenyl?kbEk7fD@F{9pE2fxf0%xJHBvS8XAGs@m;oGXES zj%IGJMBA{><<9G3?f2U-4`tAURqGZ3AtH-nS}p;-u-g=~jhzf$(&Ir~KQQTe`?l;(vN`0~OKxPNdBj{c*h;Jhzd1iGeT;o5nur?r`k7LFiOBDU^_s&?B1-(Iy<_|+5m7dKy!R3! zqSVf`XLfU}NN<?(uvKX1r_gG=h0PRL9$Y1 zQgQRl$iwwc{>*n~WW*t#|KKh&lJnWl@4S~8=^ON5%nB146?+g^`h^Jv@1OJZzs-cU z^zeFr-^YYZ3OnWR&oZJcQ|7RZ97beCjgquEX}dJH=(6ti$`j zcYKF-t;5v8sNC|THApcXKD#2Z2DkEJSU29S!nYaK_Y4fH@HJ=PwG{TXxtA7Pc%8Nk zYVVg`X|H2G16{yQLwgCdR#xZiycgkVDRr)xvH)praW(rE=7AwRTPtPHJUBVh*7sV@ z0jFo1Q$fHigs%{DE*H(f>prnuk=bcj^Jq>h0?RDv3}A z*YLdL!sh)9YgG7Z=+hSZl?vl5Mt`rnQbDWY{^o28D$p-=>J*t%A$OleW6XUjDC}`p zxf)LeuJRqv7wJ^U&|3=CJlg{!TnjlX4Ly*mul%{wwik@HCU!pD)dx?G?hL~{9qfLl zqCQmK4`)7?knR)=K#BS$e)Y;hkk}gG;5s=3y!_J)TA0(iBkxjZ^u1A#@Zt-5UP}YH z->;AKV}6dNEMt`~Wda!N6N`3A(!oCNZn*gz9XcqD$<4o~AWQbZ#lE;{aEY0_UFkIg z8V{nwJFI3Qr2NVaGt9s6=&GB1t~d{D0h4TBg%{vZ^7jZL=OVCbr8?CwEdov9)DGjp zB{+R<87k_QVRkA{;oi>`c$@uv+uPVxcwC~RV(+^K(+`epJ9d2?c4vzbFJTVM^OA$U z9%Kd-k#|1qAg(vwaXQ)DZDv50(!zNTx0w~S~?Ri!5s*BREfS%1Yahb1T_ zgGm*0Pza1J`!3frp{J_nVidPBBY{o)0nrA`=)B%8pG+@iH2Uh7j>!*ZwA&y&e)q36|9|F-Cf3tay*XT*~vsizq;}(?No^nY6tTU&}} zJ;Tn12CTllzs?o`P*?q^)(B3M)R zFd(;fL&n>A_`c;`z9STKA;{t%(k!DeFCy=Ek%}O$iwqY60}|FCdC2qFPvtcrn|tMP zHLk+xsqmskw^i6_eO2F74)ZCjlapIVR-iT2>d4!)6;O4sTTJ&}0Vay#`QuktV6h_L zxt_ra=w^&2dmdSV-QWMa{;$p&+$8vUbFVw-4QX~suX;j+wt4f_I&XM2TI{yC^c*fT zZG=W02?VwIh~T^0Z{c@%S|`J%FxZ);N#XeM0sbBIHhU8r3$F1eatf^zKnoXSHo2dH zccq#;G3py6Z|J`X3Q7mj{9qfun&ozd!5g0Z3CsJ!e!Kq4&Xk2E@+E)H}HY!_#ox#fWBR{|U zjkm8(!3DqPA||h=;oAjY_X`zRhdFT)F0x?%>N)Au4V`(ooFbkP`Emi|ws)lHH!T9i zGb*!2a2ej;Z>o&HwgOanD>nSI6doHgGoYK&wcI<@DglR?b*SAYL3is{7q&+ z;bf82bbS5m(xsMI3mH*K+;O>EM={SRqd`0n_Yvj)B=_p-F{50Gud?ianb8?DYSyqG z3+`9ke!s(c#GU+qY;Wse-C?4~Ex8cuBrczs!j)LjY`x<2%g?MxqPojw6z6H0v|n>d zcoEUJ9|`p}14Kkv%sTN_iw%{gZs~{jY)Dc14SV(&))DA`*T$4^F3~%GQ-eD@%1^)g z=-N+q%#(c86Eee&LRcT_TZxm)vO{RCn4)tt`hwH)l1*AzLQXc z%#LeTQ8>@XQ|NremxQF`tNUDTlF+|lq5Azg*v}l{n-VNULhN!olq+Z1QLG95*J>F% z^4RnBkl`D4v~l+M^^ptgsM)6EB7HZ`^<-tHJsn^}2aFFsZvBXJrHSIq1gxv1Ck?;d zBY^Xmo8P#);r+(qn{m`;Z@k~s`kv;&d7;Hz4VAP$Tt97O{QMKZitYsOHZ4?PMNgL$ zLw@#P-6*5NtK%sFiD{qeJ%)9tAAE8}!Zj?2D(&O@><$Y$7o5F6dJ_vGh78+2`^Joj zt{3-R#{2clJ@z!+Q6|JGJzZzv%Y?*MqFwVvm{8-M9$oJ|MzkfF`Ct{+;e_`&)~C@K z&`^hkj-DUp!fX@R`33W0I7^mJbh~bV+tzm36DjMERIkuef!2 z?N#u&sjAoOjO!D7s}Hl7i}SEn>Edtd5**w5+}e+K2|_L{4_`OHTpcpX4E9-oZJYR5 zO(}TJ>QYy+g1JxugB-b?c)q!J>hK90*IBsfbd5nYe+Kxsc`F98%z%q>N#6cb(@^%I zFI6*o3T|tqJ~qdGLd#sga+)6YlMDs?O*=dZdRJ~1rwoikHLJK*y~h}wrdB#Xz;kZ4 zrC_P6(IZf%bZW+2b{HJv+>7_*3_^FX?MJ8M15hSY`S)9QA8g~JKIg>qiRQJ1QZETy zhddVZdHkar(%QeSdMbCqDt%K@L`xf7xNdet;#~{qFg>-`cWr`0f3M_{9^qWIS)Gw` zQa!xmVKMSwss-MX*_A;Cpf)+;dww2(%Lk z`~{^@Q`s>8v#10*Mq3-4gi7G|@4C~wn~OoS`SBO~;bI759+aCfE&)rvD6Vb1rNAG4 z&F{$mGT`$4$6Ov)4*lrbqh{M*@a90w@kswa5W?5rF0)YyAG3>&dY0Eft=MRq_MJL# z(a^fVCGig^F}26mm>Z$y%Hbnc?9E_5_~+(J#a6h`los^uX*(>rR#nWccfv#K!3Z#* z!a3(Bj$cxGprSnEo;I$d$~awfAAaZuql8$kKYIt^0$+2+tJ)!msEaw`BA2;#U6V;QEloBfE3caDry~@X+}g82WdA zyY4V>q}&TGO<;MoGO%_*<@zF!2{_Mgl<)0QBxc3S1quVolp zb)mOVS74L3wE}5*70=g_n4Gt+gY(ke@EF+*_&Cx|_u0*We0uBUcYMVadFiaW)ySab1X7tL4V=%wU>=7 z=&@C;P7}@n=C#FIRF@M_(|ILVaaC4C99{~2`U&S8c*4{owh_@V&B}}5hVy`opSvw_ z|5b6#Tv=V24PEpnzLR^!h9nNJ8P{~OA%WiFd$Tw<@ui-ivgH;#dNQFRs+olO`6qq4 zfAr#<>X6Yop8yGcn^o@EFd!k<%A33i4@ro$SoBpbiG(Dz+qi<`)3JH@x-DlZxJUkYnRTghnLqT%y%Rx+|--(4^|NkT;uxjcJnNa%5W zkAWQK8^73Pw26M3gpS@3*m${@gv=$acr*st(Y($tF9W>)ZR$oODZGC_rrn=9QIB8G z;9`Kj6&t#BFgCE7PDGrcyBED3iO7ZUk41i#6>;&ftOVGxqS(A?x??W^eO`I?CFL{$ zIkpFV-6&!~+-nQAhvZq1&&rL?&0)-F-(XLo9y2o%+pNVCV8?{W<(|*IS{RWM|3~XQ zZAQc-;IUDW&VW?SnEk}?ewDaRV)Ms;=LZ}=4@qoa2L%tBUn9pF=>LqW1-4aK-^s7A zggK1P-n{AF{L5gMs_^5C{1V974Te?W{pn>=d+<@`1-QLOM~O9g9$Y8(MH-CEfmhg} zUPFaB*z!rYSJ`(KPK%tEksp|Wn#tzGQ-(88oMQ4eAM+)j5ml$dUV|%ZE}Z!Ji+u^v4U?d$ zI_KjXKMD1LzjbiC)Yd+Tu`X;i&16p-><>ZF6bTW{C((kVEnP|@)2$25i~4L##bn*r`QDUo-Xvv6yc zf0(g+4vt!m8MAg{{hgd{W59@coEGQxtYvZjC!{W+;o>sZsXnPsg|EQ=+IwY$&Q)+% zP)}78Ux&QWO`qT2+W_Kn8i^I>!5n7Y%EUSF{jcD%(|jf)B3SudOuCG9)ExfJhGb?m z)2i{bJDnM=3aH;WaE=A(1dp2UTEJX}8czAB*97$E{)&h?<~eMlwVkNRU`48{?o&n@ zSciFVz@;4XCK%{~Tzk~nkp23mKep*Oe@yS{zJ%)pXSV}CGW^-mU1g2l<#BeztEM=V zrAtDOEC%=;-;|8Nk$^mpWK`;k+DfdyZ>Pz87T}=dy=y7I@sO8 z-$X{MV(W`6V`Nlwvg3Bv5*bw{Qa(N1AfsUq8;sN@qX%UN1s7+@D5m6_`I!MSQl}`p zYS)laRL{L)vF~K`Vcg+>G0q)p-Q~*LbBm0UpJ-KV(IlfzwOa$$oMd!LeS6aD77{x7 z?Z&AATrb@_aq~wu<_mW>)9TNOkK2pP z(cz048klN@e3wLN9QW9e#o>K-taq~^5ttX}!}?T|@RLP_N4VcocOPz&#(eiL1>Xf4 zSdp;#hbf+itVra-?2}e8ynmSfD0BFW=R!6f-)HU+(Ef%*9y>foJ8kP8*oV0pGCc*X z7x8>3VMD!sVU-!t$fkK;BAF2@!SdNl4Q9kHwPz%ug9)89NnfEnU_x#&^D?{n@#<^X z+meB~ouYBF1E(0#z>63agE0oIN9{Yx^OON8(1awB2m|V$@8*qh#JY4#1XoZZ=B~0Y ze>vB(29l;t;j_dw=m{uevOTm4-_HvR7vX*>#mg^lK42M?^F6Y-elEc|#q1FBz#_b1 z=p3KKeSQQJrw&DN0cJfKr93Xo!@b8bMoM0DaQvV1j~VP&qnXWfm;IT6Z)^4z#pBb! z5S6P)+B6McPh7Jo-9H7BN7U>@Pt&2UDvDS&JONG}w!E2p#vzR4rIBxsueZhDuZn(* zK;^I_T{Ur4^4J%WruL zDoMYh4Lqaa_=ln+*IkmJ;;OytC_NP-{bx)!oy~+mL2Xz0t9cONU~ypnX%W=6zTGDx zQwEEZTG>ao{f1_@6ULcGs^F~i+Gav)Eo5ao@IT)B58%9nxnW!*SlK>nBpql0Hht-i zt2h_)*4XTbNK7Z3{g87rKok3T>$`4U}U~c}Z+XWj~~I879888iY40 zhlX<5hCxB6@2)cbI}Be@Z1{kEm1W<()1OI=LyXcFt)Phsu(S&E-(yLK;ozR;T-_;H zxpnhq*`aAzVcsjTU1|mds`b*J?8JV*ep`#Wt#cs7e4yzM=RB~FhI6Qq7JzyAN&7I{ zB5>%l*@m$%LEBHk^-UbhusJ)`YTKq2d_GX_pIDGuo{oGX)_bH5jC5E;6OfNPOFsX0 zR@AM1K{yikH~yV3I+K}+=qJ0~T}E6V44G3)P5OvP`X-k3Fh`)(Bek{g2ksNKW8-DF zv!hY@$x!6NjuI{Xx|t4XO`|9cSci?eI{)P&bD|E2cE zHw(N9%wMk4#m^6`e3LoWO-ASDS}N)b$f&BW%k5zh84XB%`Y>~zjOAcrKq7=6x@S6HgU z{fY8_8p6$liqDQu79L?lV!aPd6df5**xBy&LBkCQRMz0_$XbJpB`1Gwu~qm&ood+~ zvJCV%$$|jxC0KdFo)rIj0fbguYoG3#hb)6*i*gy5r%DK7xO#R5Zf(4B_{}m6e`|s) zu{#>d=nF%G+>;>Jl-uz#cML+7X@3Siu;0vi|M``R!*ILcJ<}JnK{&KLz?^5_4_CvZ zxT-(&Lc{Y>C5Is@M2vD;A-k~;`c|mNP%!9^#ZgAXP;s3AytF!uFKWiDs{)2~=#@2yX3FsGfTlo520{W@i z!Pj)IW`E0qY7>GHO7Z6; ze!j%7S`c%cemuOHevu6wRw}R|#o+7Y?YV!py=>@{!ewSVQFbI4lWiA)`BgFX=Un|` z*-=BRwC&$+c9c~gRXC0N2D2SI4WFGRA@ARPH9F5pD9pHbE#)T(-LAK|V>wDfWTAKL z+j+?7@99l*p=xBrJ|e6A-hzzAl&#;>ACS=`)7w&~H)Qnm^|Y!`BIdPlG&D}){CC0Z z)gDZHKno)`&%MFdm!w$Pr7_G`N;7sFeOO9HTb}Q*YR1>k(fgmiCufq8^5Jtw8dAte zid)d$2y=4|@N-9Ly~Fj1N$t*-05W=+6l*2whwGU=mswVx;_LH>0OQPae0^5z&AEf? zpXQsj_95@dh=cpG;Lc>s;kl#yoQQw^);VJbmwGZfF1K&n5Uz8M6eaI0z}Io(=y9ip z6*4OGI~JF|NJh`nD^B;1U~bQ@;uIBJ2OS`O{X2>|Q1?_GYv?=U>wBU`zFw7#R4|a6 zi1{a*|8W>TizA^ZRX*?8<0Ox?48bsn1>Ru<+nZuuIn6%V`JXnb7p!5 zx{bDB{dwx)Yus5>k?27&+>~M_}l_; zo;X@DfP99$a4WE)XnKGEWwv#usZtDu%Im7wGlCTW^LK+g> z3`?-MJR<7LvIx2}Vy7*(E`TYQx;3ZVJP2=fZ(%T=gC}}Th3w+a%1K z2z-4lYXV~RKfT#Qnt9BsxSi1 z_mEg-^7yaOU{U9#qmd&GN;-ePY<`M;=z9%6R)^s{=Ig%`=ReWF>9b|%`wSW+ioEO2 z%AvvR#pbh3nKXDzo(vWJjQ`z1$uISbu8IxA>PgD8j1%6dvPSpMD_8#gry=R-2> zABT;?#0#^>T5_Y{b8A8Gbn*yrxQN{r-!cLxb&-MkrD3R8d2~Sj?GWT7(Juc;8HD7Y zAGXkvFjq6!Ag1bBKO_lT#)_Tn1OM%oime+xAmh?<$pCW*@?NgLuU+eceOqZ6oKqdJ z3su-epmuPVc@UR^`2^HwVqYW{Ti|=^;IQPy7BGEHkvbvL0xKM$ah2*VaNaF^<=eLw z;5Sw-raHGmYUG&9kvDB1{gBBsd7~XN=Q*x7)ZutV$7r*d@jNn zS*HhO>)6-uBz4E14XC1>QeVb>ot!CEm!kvgnO6gr8w?pySzpzb@q9+qmmR7zhy64} zo9WWfr}!L!E$0KC_A()TwsfV(hnP`)Ym>o=C(P)FYyIiWU(86k`D|PSD+?0x@Fn!C zv!EIK%lhs%EciU6PjZ)nSkR7sE%o3uJm0B3^k@&($2$yMWi?uH{o{NkZlZ?;Z5TXM zx!uQtzBRn9oatmiqyt@V3+phy@O4U^AcY0RY)}qgPGCXa-y*AepRyqB_0VC$H5Ozg zrX+Dt4L|QQmD$P3f-0h4x4LvNBa6i0!;cf0kz;e?J!gAnw9d#n+I^52C5M&qX-O;$282q=O576gyj#n|e$WNKpSo_g+fOkdDWhoVkFpGiF-O6k zmzx3Y-uK3Ug}wnHkB@&_s@Z^~={&&$lnt;wUoRY8umQLKd^_q@rAjgji_CebiP^^=?!3|t@)@>Y$+WVi!o73V3Ss2h!$6;RY<_!?`PFB>!Tyv%d z2l6!b{e+$LUoEa*hjeeBD{E%!5a=Mvy{&N#TK(E)BaW=WsTTWP2Ke*wZ|MBu%=!wb zmp#~UG+F_NXf;athh-o?&$-Y$jJfU>zYEIuFTt%_v);CDc%EYMZJAQB0LGrBisAeV z5S5kxGS_|{p6&ArDXgA@!Ltd)5=wLM@sZD8jnA`?tfN0GC5_J!kzi45$9_M;t9BEf zlX!09Ke=&jZW>~&x&MkKPs5GcKZXzOr=f3ozT^54+_x>6*nbh7hHiho8`E2+K~lHR z^gH&|O)IL7>2IF~QO2Q*411?RIb2$E(Qq1~?JF|$-KL?^bAKW+bsE~1^{z?bKJe2k z9a7@{8Av?&BG=h{2KpXdyZf?o1`Z3eskIBwg1NS_t*Hy{AHS+xWT>Bo&cUtu4VaTp zo!TYZ{d5j~Z=7FL=*9D;ny#uz)p>ZqIVhn0dLHI4i5J`V&%;Mo-+N^T7Qo@_&BjlU zFxS)mXG=gS_RW22vzx~Ei*i?%ZUa~Zf=9W-&ijj?e3oF_7PSb|+t|joU>_Zk!7So_ z?;>Puh+gSkScGDej_m=2B^dmpFtS8mg0TNv|5saqn>iuJ1($d=7pfqpZX(+Z~Z)l<=>G)`h=rrt0@w z66+=p;B9`i53YZ@6dRZE&%NJKE!}*Xj0*4ZJmNY|Mxlpqtt;;#qb1MH;(ToQye45U zL2t~Fbe)w>KAT2DRd1?=t+9@)c(UhgtQHB0{ZhVz&h-uifXApBE3$p3GvYXk^$6rMZ555qQGGpuw4;>N%3HPHnDY; z{U;+bQF_bAVa|xM&vYsm&N84ZwW7np&lu3XHM1TKJ_fuGg;HmaV4qCZR`LC()}hd% zC7=ABbCDf&yz|DO*dICf+&dn%)PV# zft&Q1>aFJC;k$~~-M8o9gn-5VCl6+UlH}0v&wmDz^-T&lVyEF%k7Uqn;S{J?^3s9^ z=@3o)R2;btbEOVv@mhsWK-iY@yP^E!pq4u~7W9P%G1nOeXU>kodhg^OHJ%ZOE0Ppz z9v%XH87?Ey`a#fkp0u~A8-RDcAvebc`(dg2z4bpl+WeQX<(MYsVO#J798vq)3u_Vm zw!bmo>NAhV;nuhw2pM^5Wr*j-He;{f%}!HcJ)`rR?@uaN9X?m`J&Xz$D_U;oURPuTj77?!ye6@z9Rl&WHQYlS@?$FJ?h(YAAnjbr!}xUOQ}s=RdXb z`z3hm=3vNvN$>T6c@SyKVY3OGha-g!BML)!o%!g0PiX;e9j)2@)pG%4-#IAWDqeti zA$rE~)dl!E!Ft0)ei0m!sJz$DFM?MYx@ql+_or+@XRXhR@Fmc+eG_#NTof`b&PXi5 zP;Fn^KX;sSu)i_o`BeObT+_fc!{qA3aQ@*`(7Fb`PH%PH|0RDLpM$n3EP(V2`)~{#@Uf}{WU*_4sEqs7bA|~Cm;!;$ z(`9`uNUUg|1(8J)b7K!(F*h=L&5FFm`@GDnSW)7Bmq;-dA{sd+o#BG#(m(!`Zk594 zznS;Z?TJr`NH@QLi=Ke{)TEy?6zsd%mSV|n)_`@}f;|Qam{%*$6}R$wn21Ueb-!|t z5K;a2BtJf^>+GM4-k;b+M31;~I*Uq(==2vxvpb)ONJrgkbCDkrHH_U2{&kIrWYt2R zMJr=Hcp^LOCM!O_GjeL5PBkmC;W1ph8i+Z!^pz#LKK{K|T4F_5@Hum?C9pS>fR<_{ z+s@(s+({@mtdorA<1fu8SyNdMJ%gezc!~vi=*mXj9b?9MhYAH_A7->um`7&Y4rcVx z`uZs@?BCSd>+>rN`=5r2zol=kVnmq{CoAjE<9cZmADbFJ56e#Z*-DZ<-nY)^UtOHu zfQ>xE=W?tYz+LsN{|lZoB_$CGmF3q!ic72NC7y@d9&@7EV85P+dQN$u>oRcr`k;ja={T zdnUnjcyK0w^{pHobb@EjFW^3`@~64rL;nf5Dm^4a^~bz~n!p=WJSW%OBFGzwb;$tj zgLm&6j)EU&@s_^*BOn?oZ|cc242v0KzlPTaA*>=WE^}f42<&-W4zzv{W^B`Zwa^Dt zZ?-C3JQq6kMbuNSVrKEV(i+5=GG z7$4t!%_F`Wc=y%*xrpcNHl`;pO|wuTY@bBoMI$N%O$^<>L&1FVgACwt(i_X+=~ zRh$<*9H~P5-3ycaPMn!(eemdQtVL=@KkQEr=dY|80CK^Z=RfiJ2G3{SY49E$hR@aR zIak9+V2y0IRTQ72w5s-WwBBwEmfk+;Sr{CLk2$96A~$d@YMSD==tu{QH<5XvG6lXx z-p}@}O+h$!ms=ABzwe$(ZElz+SHPj!u}2f1D;JjBy@B)Ov8u6}LSN?KV4lKL3fA%K z{CB#D;eN&_zG6BmX#xI;rLMo$T7*b9mHwO9zk}K@pNzS*1V4@((k$;;f|gV^`#a~C zA(P)DTf1Tzx<5YDnvz+8B`zyXng95&LL)~~?LkTYi_cvr#FFg)Mp~wA>qtXcN_3kRH$9O#+3@#?49h1T7j#%$GS|_k( znm|U3_R&VX!DPf&T6p#GV=^**XJun!Lq-RU<5J>`F*jQTZC*Tt>oS|J4~qP_PP0;b zD>eTgp6p-o2z*ZE_U{S1%o6{@7nk`X^WZ=H?0~zCcn&26eE$rx! zsz6U}7&}ta`gbJ~&y|?wboEYRy+wpm;BIg}8|s}jAG&JIhGw??`0|jO4W-gQX_sRk z^taQN(~5eSNBcQKv80O?DSJK7xryh`!h8WdvhxI#q`re!*cIom^73w4&9k7*z7qHL z;C=MBx?f-#l^N|&(+So+$&AWX)$C2O@Z6{_ugGX86OwE*y;S><5wV@2NqXWwvsQXs ze7_1l7im!R^mpvP)AGD@d|%HRSiO)wVD^3$M7oYQ<>TDQd5Yq5XOU&#T{Ol&FH>84OdN5W9)6F zz&kxIVcePy=^>AVEgB|3&3p3wTz^rYuf50DNJe%<6h$HTRzwL=iIkKiEtO;@84XFKL}p~S&4|y(-q|yd%UeEK_d42fu{oH%bx#!+{&UwEpP7a~#<-W(vA_q|NL~f| zp~`#GQfkq&y91K@W~9J4~)lyWfcKMz$7 zEJj}foDM#LdY0owf1(V36(KfHxvv6>#VCOJMhAagG4kZGA$?0qP@bt_K?+M5x}EgT zr#+$^1&E3k-*SPx&28FW{gKruQ$18Yf$;|-LS<33MLp^;xtAkU-+&C`1CKnkZ9-Q1 znY9_4TF~K`LE3vGttj^3k=sR0?dah3vkR&nok+I0aky=<3pwWJn#U{kAosjM*s zlz$-8``c6>I`Vz+ORB;ElHg+z*1R={f+&hSE#wf=PIZd;x9K+;6p~BgaUDUk%6)Us zdPmW%se-$`M&oGbTZ2`P;tAwsK0egFZxSV?itfFi5B*k6O0rY2(`cS!f6PAH88nzI zA%8>RFB<+Jd`4&XFUnylVE6kzi_DVu-B^S=gZ$K`*S2sEM;pU$87Z~^@k7V)>i7a0 zZmBL$$XP^AT8tUu50{YhgmhEJ@ns~cSSzLl_Lb!AZ~S6%1$poV{0)J6;icE5nWB%^ z(B$(s=h;rIqY&kaShX#{Bl+iUsM9FqF=H46!@q8z{VP4vyZ6yx7MHoTkazI(PaA^!?P7)ZV`Vb!4;i`-c9*i!ZT>miXO9!BIN`lX{gQ&i>kp z8X?X$#N~?kchX@7Mv-SB%5bmo+qUlILwamjpyt2^b;0B1>=({!GGO9C=kDLHp^x#e z0^KO|2}reBb7h(UUqa`gS=yr!+-I=b%ra!@;u780vcc zrLXi(nKR?wzV#~;PI<{CaNTYNUWV|*^knUuEI6yVbbcun_;By_y!z40 zf+s2d+Ru1caUwfYGRr>5A6?^}-eSXw4-=M(nLbcYSbk<9JDL?6dw;qwmcfcu(_i4h zVpeQ={_6|oAAl1-P}KPh>I(z?xXk}+}3FPZSE zAEpBzVITGU`}tt&cIYcH6LjFX!-(ZIs!Y4s7_m@*wBK9|;FqVV*?a9^zdA!dD zePFJ-L_K=U`|{yx)DJq$1JHfNrgYely#XI?reI6|VK=Lv=YXe6AP?9o7S&aHRQlLEI+PY;8&AW%z z&NHkcg)OJnA5Sbp{Yap#RnroxQq&il%Uwi`dr!wz#xEf0m~B5U0I%PR(;|!e+~?4} zjL0M1S7s6I@h8TmXaAyL+aMp06EnzN@>=FoW4LeJdCF-$dJ?(NoVs3>GJ&>VDHu#7 zp>99)&Xcmi5p*(KQuvDCZ{%scDYnmI2sx!i@I3!8fa;E2UKpkJqoL(f4}1fAQBjBV zJ~8cIh#{dEg#`XYy$k%8{dhW%^3;g^nTl2vaHT>0#;0a9{dC7LXF(%!bZMb88LLN& zL4BLOtbQP#%$r&i@oIEqk7IbWOa(f3zn7}(QHtK=8}K&nC`K20wlaJ+E=0#@FRO7@ z63B!>!INh$AN7mKIrysPqgMk>>_e~f(a{5U5~`F4sLSpVJwaQD=0ncDGvX*lOKtrU z3i+jI|JKQsf83QQeAS#&oT(PcmbRV{{85iAD)L7|VEuDvTK+&B0{rrL>pQX6JCI=Q z9>P-?gPrdKO9C{x`8=G>5nvvK1fQ zok!d2XlZZ0S@^FXV?68IB054hzjz1qEVZCh*#UhGz3zNY+h$i#k+!um34Dp~dwyz;|-rZ+)Xt3I!GlL$fz*FdJx6^`+7LR-x z$F4D*kum54`t(ff6APKY2;dVzm1|`iK5x=J>6i(Hti4j1nF?X zH_2<;4e0Qzb&)S#t`O&^B#gNr&+3}E!O~qqhYe)%sgr$lc)H!_!NFxZydtO9+Qv$c zDe2x;0ECO>A`EXl5~0Ttum4y6-#$y&{H@OkS{7_J#LLM6{nDquHxGx;FyqyxqwIeN znena9b$qNX%$Qj_6qLP{q8xRPF=rKLyrZ$ z?zsD320ra?Z_m1J(BUf-ZT~&N(C1MU7*r!lhj}kYzEVwqJTx^x(-rCh?)qJwN{2du zUcXtk4G~&gS^uQfBZvmi&KbV^`fmgE#sR@f(?4|HUdcisb{)NzG%Mz`T|?13pJ}W} zu0sEHv2jV~G7{W0b@Jjf;16Pr9Glf#M4kN1P8q-9p3+x~N_{a0_*3Kflh!OsfA*vM z{OAlCI_f)v-%KN(ms~a~@o*oJdNRmgY62~g6?#73Mp48LUKcLA->CHX%lab4A;izc zu93q!fZX_1whYtsAto>8=SgC}(0z87PvO@((S#c_?eat`@{qb{b~>&JMQRC5$_3S< zWTW>15AIY0?@{WlkNRaOLP31BD6SB3ZCC%S&6$Tv@(TUAXVOv4!Rg6|y=f>&wd_$j ze;g9}%^nmy`x3>_{iJ1D^+9!HNuZ9TJ^Dt&ArQo`gl=ZiywkBWL5A#3p~YNxQTI8m zWAD=Z(b)M|o8q%!Xh)cyO`Tdi@(e4HbWBS_Pug~RKFU_Xp}R#$md$TDm9`9x>$|7iDXB!F-Uk}Ib8FFwT&dR5sd^L@cJE`SU=tEx8#-GA zeGUyClDU80wxJ~ket`v=PSl{aX0j_?Hc|bV7qkzX=TMYwahG_D0%YBeZYYiW)UB? z@&}xXtu`~)&6iQ&0)uJt))myAByY>Mw1UKMYOIG;t|Ftvo*r#DFP*&SbcO}^R|k=o z{)`UXLsse8&uxbD{kkMCH=OfH&*d}XkY7#Hesz_EzMDy;d}>e|a4YzAWFEe!#d8hx zE>W9-pPE*`K@9q)8NNymtt?V-PKRNGs2v@CV0Gc}jaKLjW#3d#qe~C`HW%G^GwHF# zy_nD3fcKF3z%cUp3&3}v6FQ^a2o=<^3 zr$ztAmYmGkAnL4qAl&=LGtpW$)iL9_$6H_TP-MY9D@~WL__1I~{AaehlLg15>rF=} zu;TOEe=})9olMj1#17GXR;=yHSX8?Pd`~UHMx(pgu=YG_h|Fa+e4fi?#1QiPs{yad z4nTcj-u0hBi%n3p|Ln7y1LVtXf9Xr-F0x?{Q-{Rqe{k;mgife#Kpo-Xu+zzqSHJn( zM93Y^leZ}AaX;bwc&^){@?SpG3oOrF5`+8me&1d1U0m33yv?3_#)fQ|t!H<@Fz~%( z$yG&q0pE+w*M*6*z=uLp=-!+Mb;k$wf4*sk_vbeu&-M(m;L*zf&-jJ~ckGJb3pHTD zj(?33qL-PmSF4-vaoG2L%*<%}X#n-dVGS~8p#B&eFefeqGvNk}Q;(~6K%dvhZ?)wb zm{%>u(8x=SIQEvxp@wznOI3O_@CEd6DCzH!|N7K8o2;6;K~J<_`5Y{D=&`Axs>WD7 z9Twls^)k+q4*OPB8$2Zwda768xryJv?@1n<(~Yo5a|H!lP4t^Zfo zem5Gt&fp;}xVVAl(icN6E&M~WI(@WLv+HPD&-l&w)EeL}=yFGlts*Dd&&&G;SJ2R_ zNJYW!W%Q?5kM(NH68gGFn7`QuaQH8a1CzU3dpP8~XAV?T~ecUlHK z_8dj&9FuHO5m2{37U=o3a0vB1sd;j2W&rVb8EC91zt;BPFbnU{Fhjg)CU zRF+LTQFG^?#y4thDA)TTUftJ>QaC=)?Z4W99DeE+J@TnV?IlSgy#5vFpri1($4oIQ zdLkM0(VIXzjW_$Gd%De|EKT~8uYiBrCPkc9=<^bs!;40x z^jFX!#wT~Cf3BjVTR69N+pi=3Gwu6tPX9x6+N%oN)xmBk#e!NhH29!J7u{$AEzVoY zy8q(@+*97m(T!NA;OJLb^*ph3_!uqa5F_NLPaJP}d>06uC*ub9Al#!0nR^^HXTlUKhQE?!Oqh$wsyPGp7k~18FG|C?;!(Lgb2#9- z^dJ8uy34fpw+IDUj|1^hMd+%hP^zjukIyCJ!8b!-`4413Wx% zuk-lQAE9{~;E5|&f5Hw;Kb8U2oDL`0@c2QB+8$>%d^a`C@iE*R(cJXq@cG1snQr`L z)rY#Fy47^bIMg%CTzGcPeh})NpKp;*_yhIH69-<|{e|b@|CNUpj>W&^^ZD-m+r7nh zfQo;2)?C=pPsOszhZ75WsW?GJtxK$jiVvPVUU9RVif_$m9k2LF#n+QxF7S3z@s*Af z@k;Gf+#vel`9=#B$BtF5J2z4BiRm8CuzD))?=Q8uUQ5Mt{1UP3RaAVfc~5nC85L`_ zvj-a&Q}GY_9jV?V6))U8EB>&6iib9Sr#R+QvD<~qQrsXVq&Kt{zf*DG4l}kckW6gv z4Z}eiZf=kj0XbGTtZbP_#o1J;WXD`8&VSGLqUakHYnOHghi1Wej$1vspF_p3ZGTND zr&F;-gn~-nXDasGDP_h5@}T!?o9y>gY_n_@(v(ca6^_&TF7K$g*K@o2yBI1i!Bs|z zZ>ZQ~Ja_*<2#o9H;GggSD()B#$m92=;=9|Posw{);=MACD*JB0xLnFlidj=}&};^~ zf-x22{I?mS8dR*=PS!F?QSpt0A(K2NDmETJwl{N}7hh}{=BaJu#d3FO-b&^0;^2eB zO9Rflc<8s@ak|aC_z(xxRw9N6{~d4Lv9}F)O|-eU7+>ea4olCK7oh)7w&hjO(h z3-8_-cUVTMjmLbFtsu{s8WEWm!;Uw9dv=&6mj_>$Ih$vh&5O+*X=n{lKwnw1HxKTo z;z;w?hmP4%v6VFS>9q$`ye|4H)%+zDYkbStRuu>H?5DS1{v#E;_daT~$fn|4lY)sh z(39f1lCicDD&7+$*>t^EHO>x+~!UyHcuT+QU*Kxn$xl((?!J!s3?B}B=0YM=kH*ruNz|MY=2R)Pxtr~ z`+xR(F~Z0TWSIP=l&>Jo%DhiDfPL%b9yq29vc_m;CvP7WOX+m6cfxZ>?|!Z-cz&a* z(zpkt4M*@DYw(YEfyVRZ@cT%o=XF}(CnXG)mqg(AH)YK=N8tNz8BONP!?S~DY@rW) zEj_q0n8$%1we+npzZV48hUs9wd0h_fT>3%9w3k=b(m-G3K3$x=;FlUte4>TR;qzsO zRWiXZW%OsCgEslZsPcFWkkSz%ByJ<1QtEUKY0WCqz+kOo0od zW-#cpI`GWA7ubhdLf`Ba@I%?^&3y(`{JlMQ?`LsXALwocZ|~>Dp4)!26xQ=#=YgDC zc#sJ{*(d0A^CBbu%6R5PTOtqM*Hh+U_<W0rP*_G}O8o?ANCuS$GBXa{Pow z96#u5q+D;W2`Pa5diGkc~ z9E7Am$}XQ{{_iU@kWY$UU08w7Q%cMpmw_~NvPk*{@{#r0M;nlNoEq8UAiZLw zb=u(jBN~RJTtLb-EnTI8v{YS{jDYu~HJz3g0lA>2)%qC5ahTnhV+_Vs6C^^f1@htf zJ4|jc-l6OJe3D=u=*|X zsQ8{V^V7Qo6@PrS%y&DFimh#h|NVftFnG4;g&o+30%NoRa~iCFET=T=;~_qo2gC)!`e6HSE3?;s>%w;vzUyGeu__UaCtx1M{p8;7fcNpb7#h04i*=4l{rS_x zi?>~N6HN01zA3tx=}K!x?A70Vx9rY8w8H#$XPzhr?o04lF#p1XZ*M#GXwZTeC-Q5E z^T4{C-*DncDa3<(hFZEitg!CbzL2nx06(m!PtQ@M;`Dr)m#+>`aph89WSz2~&Yw3I zQn7D_AJ+%)yA6}N04nI`NXN*SA?ShMDXX)C@IJd9>{|zU@74(G!yzgb z;%ZY^0;$>@+azSt^|<(x6it~lxVE!O2uD%UX^*?r{dd^{dGQ%s91P0|Bc3T zDn2%J$(H676}QX%G5-kqSnvB(d?5nXr(3p{B%;Cpxaf0_f?gI?Buoq5!Mt-fG4sQE z$3x#;GoMPu(OXZBOM-r6^9WOikKjLJpR?bA>}KM%ae|+pO)ejK44(^2@=1#Yy}ppx zp8XB}ui9U}qYvA-YP(U>>gB z4`AH`cKG+^b3R3wulo72^+{OQcdfZZPx0bn|Mz#~t9bG0{)|cve_mYPrf@ggm=~wq zrt6aD=f#(ek9(c013z$gE1PlU!7JEn!+^?zvn-XwEW)|5x7A%sRTgeM*uAAzh2X?i zNg;ytN*uVD_Bfe-jt#GPO6`+7%!02kyqa+W+@vr!uXwX3^f)+k^O+8B3XaTwIBySq z9T)EXe53oG9y1cnfxpkN;&*l`l8-KNW3LkpGJc2`>k$XDwk`7F-#$ODzB~f_Z)@um zdKW5wD|X*%8vG(8iel*uan8~0&MF_Q8xh@)$BYRo?)oEQ#8?7xxpu&-2;}cjhL%AR z{I&2*UW_cOQYKtBq*i!8lBpE@-X1IpeIPc3rK-?s*Tc+mgEdK=jN zl{&ua243tge(1cN0WZ!Dzv`15%Y!}ax()YCbK{=dw?+#OapNN8$z#t?apC7>kLFYl zaALJ)yTW6}9DvL6Cd9&%9a}Iiv^PP&K;NFUvMmBo=f`^H+^`udULITeb=H#wH>Hbf zhzBv_GbN#Xd%~fPp1sz~CXo?;+4`F~wgCE!7LtQ=M(DA$sLMK*p~qXDbq8GT(&4`{ zjq4-O?`J7V?${0Yz?DZIH|!dr#bw`q-aQjQi>qiUf71@q;*iLBD?aE0m5cQKb?g@n z=DxZhKT!?+G>-565&@_0jO^Q_CqHTM<$ovV&lS+%DDiXDH_WtHDOom_1JU9=VvTF8 zP(QnTv8klTkrvAcsXI(~)8dE~I@(ltuAH?GJo$|di$Z_Ia5XK?L>IG`A&;hgjNNsL z1@QPDUy~?@dnA@{uOnrEU#Yujw{`*g(B)c-N8;hVa~i8|yg3c`*uJ_PtftThK&fu6 z)1lyP6YiT1$x`s}3t{8pC0Z5u#>wD|H5%E~*)Q|F{E zCGPcty68&IdmVraBG5g|aO2trYBbu$5*_vr{SY`K{i1dq^_Mse#S5&XiZ5lxQxDe= zPmM#H_3kxf;L>owecu|o9k+!p4C)I`3|*faX)lFy5Py{$Yv20J^u3t}SAL_X z)+Y1fjE@#3-tT#_>H2$*Q)RqRFC^kyAq@Up8aS?|0{hdUJ<$wWV88wKB8_`tUyUbz zv~Gv<&t<)A_015EHAQ_Ymw7ScrZ>~2mAu%3-Z`iBAup~Q?Zw-bd9jzAeWPy)51#UT zw+|iV!4%V?rsf)MtW}Vusb;{9?M~h{Qj6!pSH335KF{LBWwT$4|Fv^q%Cm$s+c`L} zTE4dfBk+)4k~CdbjbX$7G#ja_P`4WPTJxftCM!O$$%0-=M?cl%R-mSvnlT^SpTb6yUZTm@*14qu|?CQ>+5-DfqkYcLD!s3Kkv` zOv{U;U^hGWQ~GdEpW075-IPVaVtlOmr|aRqDx%+I67bRv8x>k;3(#TB+vU=s81fsq z!|ytL&-rt``_>z*9JB@ReHcfs)=}XjvmME)>zyN>({5! zvlCCh(&LjF=c|?=uUqNIRyWANfN%dx4e)|KHi^hqY07>EtghE7@%t15u6g&Y!_1BW zd&->D7IOd|$eX3r*3i$}mwuMR&>nDe*Q!RFPcvXr{!oE0;NV{0Y{fVv$$-nEr{nUM z=<($JA#09OdYlqizApBd9*dc|FkL)Ik97?WZL^z&y8W#zcVz>Uz~w2dIa)>DlYK?=|GL&$Q-ffFIjaPqWJYqrq0!KHT}}Ve&p87C*M;-Xy zkL+aY`nr{g_?TUXiYIS8?G#ss{ka0s5iJMj zuR{L}xlIsf+0}LR7~q_!?J0K=;%9IByRotsUhEnmRuu;OL%yhN;l&VMd|(J%4mbW38Qc;5hYM%L45SD_JPAME$UI)biFXaYyArDl@z!$JbEP^CoN0YV zOW*`KtUIU)&H?#F_Z0YQHYs(>_k5}07S;xv335x8veXQDma}2`7fcy3jcj<+AIG~aaId{clk#ozPIgSMzu-ax{OzqLmZD#VvE#(ZgN|E@R|<|)pO-Q}2T=hQfGs70~# zxgJhDnt6cfoijIfT|FacnIHC@<3_cD__7ok|MgzwcR^`V<|i9xg_<uHJRv}FW##Wlg4oZ`n*~f04r;#BWUw!)G zAS^>hJWx%3x?6^He(LAsXCypHdS~cQddQHh`AR2hgJj6>mVTXk$jDGUS_Z-}%{Q8PZ>Rd;X3=8B+AxmgTm;GUUrK`Km%1S+awEfh0L( z$(U_|lC_&;$rh=mr?13h$#>&=&0Q+8WPDz6)e$~fa=dr{`ge6%axastLdq^#()H=F zXOZyvx%2}=D)h2sU((nJWk!ZPXFR&T6~0e8>vxu4iVUgcVm^DqL55s@aYuZ+lni+! zo#UN#p)^@qS86M#D@_^+jYgwHDe|R{Vo`XHBvrwtBMY5dYMhf;$t;zbtE$2ags4IU!?R~ynUTGb^nPM$0Z7>V9$J~=m;IDc}Oy`{pl*v zYPhN3gTM$OKS_1gE^0;Aja(CRJLZwqn-CpdxpfpB_c|gX1@PVC_Z=_?T+j=aQm?ro zzr&{YFDn7^9u^PhLTyHXU$@8qq(K8G*4s90$DYrH|5eW}l^^E9wv*4~9Ai1KMVk-r z1<2PdwdTJYf_>7~sAOmA6E z;E&%osZ@0lf{br&t7Q%oZ=Li6{a7c7!t7sn^j^#oI>~|I)wD~*?2CNnvWzui_VV<` z@5z5ez54va%XSo!x=AZ~>u2CmO=A~s+rvOQd}yFgW@RP?<~>E9IkNy8`__}UyjaN_ zYHAF_COqWLBKIE^e>q4=HUH^;J}%N<;YsO+KR4+kDoei7<{?dkYS}KV@Q@76Pj4CL z@RDbjlV%L+@n}f7M0qxmA0^5fu8hFv$5(pgQv^w`Q|WWUCW55Ls@1Ee zZUIufkDHb0r~oO)d2u2xl%H&_SBPrJ<|E?*jt~}URI;K`ko))zUh>;b+rzK^agz!s z!r2ZGoTQCW@0UblHu8kdfZL7$W-_kG;!QH})utF04VvVzldhjAyBwE!$zJaPVGU^> z(tU)epWQQ0*k;IQ&{q#4%Og`tVrjF;l{07M^~n{KS`c?=er^N(i%+?9#h-#3=;NC7 z0e6G?@TKuA_X4COjClCUSwF3A-)+8`Px*eA{Qi zbyw2qasTwigo$T#_~Le3uEh%!yi1TaZ8OxFy>@DQhLWKVt*KW)1JhvH9L;-t0~@H$ zEbeuc{RWb3_FB0a`VR#%ZxWsxSwpV94A~xzE2!6@b!hjaMO1fOs{85jIkZ27^!PD7 zjS`P=Ua>D3LsotZYhtRvllk`cFTu7}aB%1gG8k$nYEyd#Ten}uocgv# zk1IEmhd#aw{URbv^7<+Bov{-pkFJg7X1x<8@z=g%td+v#IWKy7iA7;DPO&_5h+Twq zq>nhUz$rpz3Or>H|0hg-aHE`yX%{9hXnyOqih#c#!P|!p2$Mw#@rFCAH zOxg5j-h;$_=O{^PS3B`p+-j%UkrHA;pY6Z|$5)hLwM~k>vk7fEM^2p<>PK;j$3q!@ zjG?Pv^L#u37a{Zh*<;^(7f>Rb(lgHaRkZ35u#^`GxWSLv?B7DYoEg=f`HLVU-ah{A zrTK4W>^Q<5#SMAz^lv^jS6_1Cq;2yN&oj7j?j^wlZ(%r(jBEYnyhp_=6xtiQZG2eY zCcW{XoB+1pd*_Gp4MD7Q=JG1|e+8egqYaZN>xVZj!Z@^U6L67> z;M4@ww?%{~ra$i#ddgP}Z^K)83RSmYsS@sm;>s;}oj&vtb*nhuU|P6;`I z4id*nObnCDx#HOMN8J9uz2f+@^%;c=bQ1UyR!%-4C4oH)9=d(pFM$m{xHg_XCxNeZ zGL^C1lfdsBPca<~m%#FT$ERtsC9r8{tLDBY2{;F3(D_eF;M~uD6T6rtv6vWV;Ri`c z?5_Nh#d)73{#vx&?`SEBv+lT%jh>P?vbr_+Zn7jMr>+f?t&&)M@OCIKgA|tO`l{=o zD20UvZT@KLO5p}6ft&}2LE_rJ0(B|;gI6qpDj)^?6u6$4mBdV)@ixbYC9(7k3rB_x zNxUPQO24Qng^kB97f8L5`u}}aMRui1WbP`VqTQdUw6a0GxxZKb+Lj3-smMp%EchoO zJ22Q<^(PDk$6dHHMAW15#ki7Ukt+1zDtpW6SANLYXsqsudN~m_$b(0z&4gsaufCm; z?F7;|GkPwolTeZ3{A0hRiwHipF%`4WO^DM@a!X405?8NhT72K=CGw6Eozi{EjjHrKuxc{yFrxGPQ|#d+ZWp4QeCy zP3*N~`O!&)9i{pF^!_hmgms{6e_9_Q#yXuJ2E108$iQZk*Dx_Gr&v*Wf0Pj3qDkrC z{6pkUcAs=OFiqG5Tzm9a6!>ddy|isArU>fm^xCv|;NRW8_exkr3GqdcmUpnC0c|;a z`NK24ezg5&!CU%jz)}7+W8cZSju@+-DpubDoKUT6=Up=k0ulQ&wSq^`l+sANszV7n%}@L9+x2gVUG1}x zo~tuD(3!jYXHgt+sI7&8$1;(4wxIf~UgIS3?^-V7m(cHMoQGv~`-4W**6lN%9Mpv_ z55M!@B-ewqycN_qp8Z1XZigAJw*N$iZ}hGe+q5GqDc>z-s2S;7`#Fgk)uFqeA8kL~ zR*o#46ls$9%TS)>>A-}^YNRCja9Pi+30+ayZ#v!IjvAeG)1AdfAN5?w+Crn4$ICH??JE4-MLTM{6devC#nwI>_!Hs zJv775cOj3CO!{2GpNRWoba}>&4s^D!eeX%&iw#o=xIZ)8hRAN6Kjj|+vAl@ z=#C^|m~o;C)$Va}7CF{{(#Kn06(^M;&Xhg(zaIFGj`=a{tgg&O$9Bw!vxj{|uKUyZ zTg;vy)3_rXEvhxd`?ao=5}Pk*i~iM;nfXj&m(fH{jAIcI-SsV~C$@|T|K)n<#`$XE z_{8gLCf@Z#BvCBlnbb%memkC$;Mq)2GhQar?rkM5Cw%EB+|o{r7gt7C%XAQq2^>ER zIXj8_p5Nm-8GaI*A7z-Ip8QG7=9isUkLV)Cp0w%(DRmR)gSAa%vb%{v5gCfqL^lzl z`e{P_%r7G0a)XY!+%IA?)1}H*@n3|Jb_$=za5rHbCfP$qbQ8fF%_Gig-Ncl-`_hw~ zE<$8D`9hy`7jd9j{M37gpF~W(koyalPC_YFvEi9|2k=4n*+1`UBfdO8xXCP}m9Q1K ztRd{tLO4rVhDnRJ5J!ny&TrYQ) zWa-yP5OtO$o7mhQ`KN6qCvK@`+rJVb&0|IV>^nz^PO&kC$A=q;P4f&Lcp;bQHyrqy zTN6(dnmIVH?|DmXTgo|hBQ1?kku{uU_5MZ#w6i`foGu|!>?AK!LTib+zN@j9jy4iL z7k}x51T+&lJ3f|W*0m5WorTweO`3`6TJ6yuqaMk4)nRpWYTJrUcoW9bP; zEm3!|^@iM1HF0`km4DE>n%Ex|MH=}1AQ;RXG=kF_iLz40_Jdzq2%fo~#T{qbhz8ZH zA`R7cBGMz<>}p6mabql*_f1thasCK3z>%|q2=uUGy;j^#xVqcFl$&TH3bv_Ee%RYa zs4QGx&Aix3 z)Svcfukj|#bV^R!&4wXemrtsPQa_+3$z9e*BD0XwymIH79f8t)efwgZi&4*mk0zVg z%h99PohAm=Rp>N3+;eWJLkoxOhkmWsqp@|<4%TOl==yPydW*p(q_7s=dDpN7T}NbG zxm7DNiP)2>VcdqgQ*U_t?`=owieL6U9BD^K;)}y7f;v!{i$#Y4LnnIA|M{w&c_(Td z))SI`-ih2aO1_SM=tQ2Mni~vrJJH(ReM66`JJIQ=G4DC2PDEs{i^p#1Ks8QZ%28S; zI_4J-VjyL0sCr26az51X!??BRc^%{kfwA^vQBq zK(i;JMi}@~=&TK|o{o-eRhRP84>z^WG`Z#EeBYY2S>yVaeAo9C* zC38nL5=s}h(CnsZCPJt2SDS+^MDx5|t^4U#qB!JUTGZ|~;^w1wXWfAx{KMbkoY9ds zVzPnh$@lO!LbT8Lfq8W+(RP`?S;4G@nDqT^nl;fxa608t8g?}hvK*coe#&*kY2&CP>v#sw zlEl!pGVjNjAPh{RpQSMu5i@mWO>b|0SroZdGpV%wpAC#&Gj_aMmQq#h_to)*9% pU5sr>m78(HPz{GH$SZjT%)8Y#V@(m>aV1>=tjqg7_7Gnf`XBpyDmnlF literal 0 HcmV?d00001 diff --git a/examples-gallery/plot_car_on_racecourse.py b/examples-gallery/plot_car_on_racecourse.py new file mode 100644 index 00000000..d18308fa --- /dev/null +++ b/examples-gallery/plot_car_on_racecourse.py @@ -0,0 +1,415 @@ +""" +Car on a Race Course +==================== +A car with rear wheel drive and fron steering must go from start of finish on +a racecourse in minimal time. + +The street is modelled as two function :math:`y = f(x)` and :math:`y = g(x)` , +with :math:`f(x) > g(x), \\forall{x}`. + +The car is modelled as three rods, one for the body, one for the rear axle and +one for the front axle. + +To ensure the car stays on the road, I look at *number* points +:math:`P_i(x_i | y_i)` +distributed evenly along the body of the car and force +:math:`f(x_i) > y_i > g(x_i)` by introducing additional state variables +:math:`py_{upper}` and :math:`py_{lower}` and constraints +:math:`f(x_i) - py_{upper_i} > 0` and :math:`py_{lower_i} - g(x_i) > 0`. + +Also I limit the acceleration at the front and at the rear end of the car to +avoid sliding off the road. + +The is no speed allowed perpendicular to the wheels, realized by two speed +constraints. + +**States** + +- :math:`x, y` : coordinates of fron of the car +- :math:`u_x, u_y` : velocity of the front of the car +- :math:`q_0` : angle of the body of the car +- :math:`u_0` : angular velocity of the body of the car +- :math:`q_f` : angle of the front axis relative to the body +- :math:`u_f` : angular velocity of the front axis relative to the body +- :math:`py_{upper_i}` : distance of point i of the car to :math:`f(x_i)` +- :math:`py_{lower_i}` : distance of point i of the car to :math:`g(x_i)` +- :math:`acc_f` : acceleration of the front of the car +- :math:`acc_b` : acceleration of the back of the car + +**Specifieds** + +- :math:`T_f` : torque at the front axis, that is steering torque +- :math:`F_b` : force at the rear axis, that is driving force + +**Known Parameters** + +- :math:`l` : length of the car +- :math:`m_0` : mass of the body of the car +- :math:`m_b` : mass of the rear axis +- :math:`m_f` : mass of the front axis +- :math:`iZZ_0` : moment of inertia of the body of the car +- :math:`iZZ_b` : moment of inertia of the rear axis +- :math:`iZZ_f` : moment of inertia of the front axis +- :math:`reibung` : friction coefficient between the car and the street +- :math:`a, b, c, d` : parameters of the street + +**Unknown Parameters** + +- :math:`h` : time step + +""" + +import sympy.physics.mechanics as me +import numpy as np +import sympy as sm +from scipy.interpolate import CubicSpline + +from opty.direct_collocation import Problem +from opty.utils import parse_free +import matplotlib.pyplot as plt +from matplotlib.animation import FuncAnimation + +# %% +# Kane's Equations of Motion +# -------------------------- + +# %% + +N, A0, Ab, Af = sm.symbols('N A0 Ab Af', cls= me.ReferenceFrame) +t = me.dynamicsymbols._t +O, Pb, Dmc, Pf = sm.symbols('O Pb Dmc Pf', cls= me.Point) +O.set_vel(N, 0) + +q0, qf = me.dynamicsymbols('q_0 q_f') +u0, uf = me.dynamicsymbols('u_0 u_f') +x, y = me.dynamicsymbols('x y') +ux, uy = me.dynamicsymbols('u_x u_y') +Tf, Fb = me.dynamicsymbols('T_f F_b') + +reibung = sm.symbols('reibung') + +l, m0, mb, mf, iZZ0, iZZb, iZZf = sm.symbols('l m0 mb mf iZZ0, iZZb, iZZf') + +A0.orient_axis(N, q0, N.z) +A0.set_ang_vel(N, u0 * N.z) + +Ab.orient_axis(A0, 0, N.z) + +Af.orient_axis(A0, qf, N.z) +rot = Af.ang_vel_in(N) +Af.set_ang_vel(N, uf * N.z) +rot1 = Af.ang_vel_in(N) + +Pf.set_pos(O, x * N.x + y * N.y) +Pf.set_vel(N, ux * N.x + uy * N.y) + +Pb.set_pos(Pf, -l * A0.y) +Pb.v2pt_theory(Pf, N, A0) + +Dmc.set_pos(Pf, -l/2 * A0.y) +Dmc.v2pt_theory(Pf, N, A0) +prevent_print = 1 + +# %% +# No speed perpendicular to the wheels. +vel1 = me.dot(Pb.vel(N), Ab.x) - 0 +vel2 = me.dot(Pf.vel(N), Af.x) - 0 + +# %% +# The equations of motion. +I0 = me.inertia(A0, 0, 0, iZZ0) +body0 = me.RigidBody('body0', Dmc, A0, m0, (I0, Dmc)) +Ib = me.inertia(Ab, 0, 0, iZZb) +bodyb = me.RigidBody('bodyb', Pb, Ab, mb, (Ib, Pb)) +If = me.inertia(Af, 0, 0, iZZf) +bodyf = me.RigidBody('bodyf', Pf, Af, mf, (If, Pf)) +BODY = [body0, bodyb, bodyf] + +FL = [(Pb, Fb * Ab.y), (Af, Tf * N.z), (Dmc, -reibung * Dmc.vel(N))] + +kd = sm.Matrix([ux - x.diff(t), uy - y.diff(t), u0 - q0.diff(t), + me.dot(rot1- rot, N.z)]) +speed_constr = sm.Matrix([vel1, vel2]) + +q_ind = [x, y, q0, qf] +u_ind = [u0, uf] +u_dep = [ux, uy] + +KM = me.KanesMethod( + N, q_ind=q_ind, u_ind=u_ind, + kd_eqs=kd, + u_dependent=u_dep, + velocity_constraints=speed_constr, + ) +(fr, frstar) = KM.kanes_equations(BODY, FL) +eom = fr + frstar +eom = kd.col_join(eom) +eom = eom.col_join(speed_constr) + +# %% +# Constraints so the car 'respects' the road. +# + +# %% +XX, a, b, c, d = sm.symbols('XX a b c d') +def street(XX, a, b, c): + return a*sm.sin(b*XX) + c + +number = 4 +py_upper = me.dynamicsymbols(f'py_upper:{number}') +py_lower = me.dynamicsymbols(f'py_lower:{number}') +acc_f, acc_b = me.dynamicsymbols('acc_f acc_b') + +park1y = Pf.pos_from(O).dot(N.y) +park2y = Pb.pos_from(O).dot(N.y) +park1x = Pf.pos_from(O).dot(N.x) +park2x = Pb.pos_from(O).dot(N.x) + +delta_x = np.linspace(park1x, park2x, number) +delta_y = np.linspace(park1y, park2y, number) + +delta_p_u = [delta_y[i] - street(delta_x[i], a, b, c) for i in range(number)] +delta_p_l = [-delta_y[i] + street(delta_x[i], a, b, c+d) for i in range(number)] + +eom_add = sm.Matrix([ + *[-py_upper[i] + delta_p_u[i] for i in range(number)], + *[-py_lower[i] + delta_p_l[i] for i in range(number)], +]) + +eom = eom.col_join(eom_add) + +# %% +# Acceleration constraints, so the car does not slide off the road. +accel_front = Pf.acc(N).magnitude() +accel_back = Pb.acc(N).magnitude() +beschleunigung = sm.Matrix([-acc_f + accel_front, -acc_b + accel_back]) + +eom = eom.col_join(beschleunigung) + +print(f'eom too large to print out. Its shape is {eom.shape} and it has ' + + f'{sm.count_ops(eom)} operations') + +# %% +# Set up the Optimization Problem and Solve it +# -------------------------------------------- +h = sm.symbols('h') +state_symbols = ([x, y, q0, qf, ux, uy, u0, uf] + py_upper + py_lower + + [acc_f, acc_b]) +#laenge = len(state_symbols) +constant_symbols = (l, m0, mb, mf, iZZ0, iZZb, iZZf, reibung, a, b, c, d) +specified_symbols = (Fb, Tf) +unknown_symbols = () + +num_nodes = 301 +t0 = 0.0 +interval_value = h +tf = h * (num_nodes - 1) + +# %% +# Specify the known system parameters. +par_map = {} +par_map[m0] = 1.0 +par_map[mb] = 0.5 +par_map[mf] = 0.5 +par_map[iZZ0] = 1. +par_map[iZZb] = 0.5 +par_map[iZZf] = 0.5 +par_map[l] = 3.0 +par_map[reibung] = 0.5 +# below for the shape of the street +par_map[a] = 3.5 +par_map[b] = 0.5 +par_map[c] = 4.0 +par_map[d] = 3.5 + +# %% +# Define the objective function and its gradient. +def obj(free): + return free[-1] + +def obj_grad(free): + grad = np.zeros_like(free) + grad[-1] = 1.0 + return grad + +# %% +# Set up the constraints and the bounds. +instance_constraints = ( + x.func(t0) + 10.0, + ux.func(t0), + uy.func(t0), + u0.func(t0), + uf.func(t0), + + x.func(tf) - 10.0, + ux.func(tf), + uy.func(tf), +) + +grenze = 20.0 +grenze1 = 5.0 +delta = np.pi/4.0 +bounds1 = { + Fb: (-grenze, grenze), + Tf: (-grenze, grenze), + # restrict the steering angle to avoid locking + qf: (-np.pi/2. + delta, np.pi/2. - delta), + x: (-15, 15), + y: (0.0, 25), + h: (0.0, 0.5), + acc_f: (-grenze1, grenze1), + acc_b: (-grenze1, grenze1), +} +bounds2 = {py_upper[i]: (0, 10.0) for i in range(number)} +bounds3 = {py_lower[i]: (0.0, 10.0) for i in range(number)} +bounds = {**bounds1, **bounds2, **bounds3} + +prob = Problem( + obj, + obj_grad, + eom, + state_symbols, + num_nodes, + interval_value, + known_parameter_map=par_map, + instance_constraints=instance_constraints, + bounds=bounds, +) + +# %% +# For the initial guess I use the result of some previous run, to expedite +# execution. +initial_guess = np.ones(prob.num_free) +initial_guess = np.load('car_on_racecourse_solution.npy') + +prob.add_option('max_iter', 1000) +for i in range(1): +# Find the optimal solution. + solution, info = prob.solve(initial_guess) + initial_guess = solution + print(f'{i+1} - th iteration') + print('message from optimizer:', info['status_msg']) + print('Iterations needed',len(prob.obj_value)) + print(f"objective value {info['obj_val']:.3e} \n") +#np.save('car_on_racecourse_solution.npy', solution) +prob.plot_objective_value() + +# %% +# Plot the constraint violations. + +# %% +prob.plot_constraint_violations(solution) + +# %% +# Plot generalized coordinates / speeds and forces / torques +prob.plot_trajectories(solution) + +# %% +# Aminate the Car +# --------------- +fps = 10 +street_lam = sm.lambdify((x, a, b, c), street(x, a, b, c)) + +tf = solution[-1] * (num_nodes - 1) +def add_point_to_data(line, x, y): +# to trace the path of the point. Copied from Timo. + old_x, old_y = line.get_data() + line.set_data(np.append(old_x, x), np.append(old_y, y)) + + +state_vals, input_vals, _ = parse_free(solution, len(state_symbols), + len(specified_symbols), num_nodes) +t_arr = np.linspace(t0, tf, num_nodes) +state_sol = CubicSpline(t_arr, state_vals.T) +input_sol = CubicSpline(t_arr, input_vals.T) + +# create additional points for the axles +Pbl, Pbr, Pfl, Pfr = sm.symbols('Pbl Pbr Pfl Pfr', cls= me.Point) + +# end points of the force, length of the axles +Fbq = me.Point('Fbq') +la = sm.symbols('la') +fb, tq = sm.symbols('f_b, t_q') + +Pbl.set_pos(Pb, -la/2 * Ab.x) +Pbr.set_pos(Pb, la/2 * Ab.x) +Pfl.set_pos(Pf, -la/2 * Af.x) +Pfr.set_pos(Pf, la/2 * Af.x) + +Fbq.set_pos(Pb, fb * Ab.y) + +coordinates = Pb.pos_from(O).to_matrix(N) +for point in (Dmc, Pf, Pbl, Pbr, Pfl, Pfr, Fbq): + coordinates = coordinates.row_join(point.pos_from(O).to_matrix(N)) + +pL, pL_vals = zip(*par_map.items()) +la1 = par_map[l] / 4. # length of an axle +la2 = la1/1.8 +coords_lam = sm.lambdify((*state_symbols, fb, tq, *pL, la), coordinates, + cse=True) + +def init(): + xmin, xmax = -13, 13. + ymin, ymax = -0.75, 12.5 + + fig = plt.figure(figsize=(8, 8)) + ax = fig.add_subplot(111) + ax.set_xlim(xmin, xmax) + ax.set_ylim(ymin, ymax) + ax.set_aspect('equal') + ax.grid() + + XX = np.linspace(xmin, xmax, 200) + ax.fill_between(XX, ymin, street_lam(XX, par_map[a], par_map[b], + par_map[c]-la2), color='grey', alpha=0.25) + ax.fill_between(XX, street_lam(XX, par_map[a], par_map[b], + par_map[c]+par_map[d]+la2), ymax, color='grey', alpha=0.25) + + ax.plot(XX, street_lam(XX, par_map[a], par_map[b], par_map[c]-la2), + color='black') + ax.plot(XX, street_lam(XX, par_map[a], + par_map[b], par_map[c]+par_map[d]+la2), color='black') + ax.vlines(-10.0, street_lam(-10.0, par_map[a], par_map[b], + par_map[c]-la2), street_lam(-10, par_map[a], par_map[b], + par_map[c]+par_map[d]+la2), color='red', linestyle='--') + ax.vlines(10.0, street_lam(10, par_map[a], par_map[b], par_map[c]-la2), + street_lam(10, par_map[a], par_map[b], par_map[c]+par_map[d] + +la2), color='green', linestyle='--') + + line1, = ax.plot([], [], color='orange', lw=2) + line2, = ax.plot([], [], color='red', lw=2) + line3, = ax.plot([], [], color='magenta', lw=2) + line4 = ax.quiver([], [], [], [], color='green', scale=35, width=0.004, + headwidth=8) + + return fig, ax, line1, line2, line3, line4 + +# Function to update the plot for each animation frame +fig, ax, line1, line2, line3, line4 = init() + +def update(t): + message = (f'running time {t:.2f} sec \n The rear axle is red, the ' + + f'front axle is magenta \n The driving/breaking force is green') + ax.set_title(message, fontsize=12) + + coords = coords_lam(*state_sol(t), *input_sol(t), *pL_vals, la1) + + # Pb, Dmc, Pf, Pbl, Pbr, Pfl, Pfr, Fbq + line1.set_data([coords[0, 0], coords[0, 2]], [coords[1, 0], coords[1, 2]]) + line2.set_data([coords[0, 3], coords[0, 4]], [coords[1, 3], coords[1, 4]]) + line3.set_data([coords[0, 5], coords[0, 6]], [coords[1, 5], coords[1, 6]]) + + line4.set_offsets([coords[0, 0], coords[1, 0]]) + line4.set_UVC(coords[0, 7] - coords[0, 0] , coords[1, 7] - coords[1, 0]) + +frames = np.linspace(t0, tf, int(fps * (tf - t0))) +animation = FuncAnimation(fig, update, frames=frames, interval=1000/fps) + +# %% + +# sphinx_gallery_thumbnail_number = 5 + +fig, ax, line1, line2, line3, line4 = init() +update(4.7) + +plt.show()