From b6a7eb9cedf84575cdc71b091ab84782ec1bb038 Mon Sep 17 00:00:00 2001 From: Jorge Gonzalez Date: Tue, 15 Jun 2021 10:40:13 +0200 Subject: [PATCH] Adding a time step for collisions A new option has been added in which MCC are computed with its own time step. If no time is provided, then the minimum time step of the simulation is employed. --- doc/user-manual/fpakc_UserManual.pdf | Bin 169499 -> 169553 bytes doc/user-manual/fpakc_UserManual.tex | 10 ++- src/fpakc.f90 | 2 +- .../mesh/inout/0D/moduleMeshOutput0D.f90 | 4 +- .../inout/gmsh2/moduleMeshOutputGmsh2.f90 | 2 +- src/modules/mesh/moduleMesh.f90 | 81 ++++++++++-------- src/modules/moduleCaseParam.f90 | 1 + src/modules/moduleCollisions.f90 | 3 + src/modules/moduleInput.f90 | 38 ++++++-- 9 files changed, 88 insertions(+), 53 deletions(-) diff --git a/doc/user-manual/fpakc_UserManual.pdf b/doc/user-manual/fpakc_UserManual.pdf index 99e12e789f42bdd8eb3e357819b7d27cf3d5c2eb..6b3b20a422160a8ab4b068d4c7dcc3440bc83499 100644 GIT binary patch delta 21663 zcma)icRbbo`@d0y>{*CpkDPOyX+-u&*<@w!y+<4wB|ByB9kMqOviBZYr6{`y3HiO< z_vap+b3Y!x@A>1-eYsw**R`+bbsZ_Mh(})$$9zNgIxrz9nh&+V&M=!<(e*C7kWw&! z-Uo>hyIgVs(L`RBFY$V0&;B4oo4u(Z;mXC{PO9J}e}jFS>-RM$Vo#D;ku_@z5<>39gwNUtAQr&VwjQ!=8q6246PhM! zQIT%Z!P^nbWC^Td#c?dL(8Ve|h{kUDZ6?Q$EZO>WyI8FWFm%MKlKsJ#r%p zW3Ga0MW*E?R&#W2rryoHCbGe8?e!prflJ!#Ax&{0&X^&^@VhG3V8Jx0DXF4D<;g0eW0TBBAigqhRm*1FnRiD+B82eli%8t4F2j@H505g0xLT>V7Ss|dvLRUZYl* z(JyJ>H|-`qXOF&4UTx%Ucbz+*hgF+hIcTUXO<1@2aezf)-C_ce=`>axS!llRnD1 zCQb9vaT#A>?d9c~$vm~DlFnRJ-^o=?IUTvONivlGz4i0`w*%`LmJ54K3Af)}n$+O2;^pjBq>kLpj#qF4^OC-n25Vg? z2+fPAs~}WvdF3qKt-LQ?8&6tN{&a+XPU&5F&I+^r;;34yA43JggX;!EVC6u;~I6luKjnW)#Ezn9+hVfvJ~)JA-@Z?=1=Jhf8I zkkywKID|{2?T4DJZXNY-Kvs8N0~aUmyHD%pI*=A`S~GZSze>T9pPf?6*Nh+eJWl=P zWW+z&9l=}ObT8;1Bs#HgF~tdrNIw=#_eCQUoV5wwuFA*y1o4M}NY{GQX%!NJ_DsG% zuzsug4%u0M{auI7)P$9GRDFcN?Stp+9m+Wo&tu76CGx1%BvVY)ll@=?U6gl$_YEv4 z-i+xlK0LEN4fD8HO+4Z~zG!kbb)f1)H?;2k*llwoF@jMplNYvd`+bdU~Zf$X5Oiqymj%a*9xsup+4` z!1V<$TRi{yhS08?#&y4>#s76m>(6*~3a-3Q%p=!?jI|bjd|w%;XQh`LVF-CXB=r5^ zub(-WdIZ=H4kD8GrXLg3k`Oo()G#g@MGh}MWOg!D>irfgQ8w9k#4ZPb-H5=)cjx8^bDd0o+oitsfnkD-A^a}S$1-^JJ{QeuG;>-`is5!R$M{ZJ%7<> zD=gWUg-_q2y1bo>IyO^X3O)k$4P=`x`EDS_{^xH=x#Nt776eSF+Z%h81)uX#IW= z^GVa*&~R=RsiF!D5xte;*I@Ow`i(vQf%tf@om*v;``BpC1uOn*3p`P#_}4$MfR4yKv|tYK)Ru^@7lr;GXnV#@x79U){U{a}t{= z+In*zGKm`Ls-gWSsAf{vjY7X2TPN*s#FH%I5PGgfw24T5OG=ZvbgeyPCbwOrcKfiY z(h;$>V?Bl3;2zgj!>cuJoX5;Z>8~MleV+(}1dZf;5-+0pgBBsYNi__}@DQ4hIC((Q}K)LDT=xM%inxl^h#JCx!s z_aSxI*t-edyizjPSe+ft?c4KY;HflXHRGebFZ1#=LLXi2A{l*(d;IdK?Q~`4c;R_- zzR=Va69Ul_`0LgDT)qoS#<{^7S<@czuar)H_Ae)wTFYMSrbZYO3Jfc9-Ahj``|M-k zv#zPTGVqY4ljCXYb^H(z#MNuo_E(dO7xV(`q^Nn@-$H3&%0NsY6v3eQTwvSoI6zO` z`odNw6I>asqe`h&kMr#Mde+fx3ubL8PoK+^Zl(E4%DM3?^to{E{bxvPSp3Sk(&Gur z*WQ`wkIYIexwz}wr|+B7c7MoGG4z!Rgoq||9&IP=84W|4!N0 zZWr?E(g8X?E7>^ec@zDewI_Za{!@G|3w^Dh_^#6bvcRVzV3Sh3exp*x!~1ccUdUC! zJ9-l!$wFk=!kGi%v_((!E#h4QK4WTn4Z0-KDBr4AowBv0OdqQ1k3uf6KP32cndH{= zHfZA#nSva7Xk+R_8s_HZ384)Mf^8A;JB%?|1eW0gXns?eU}L@N5GW}JUrN@2=4T~>ibaXjqgYGRO`7S z7K6HkjObg9_JU_Zx@+b>m%CIP6j*plekd3gBpCAGDT}E4d9!Cg)zxnoeBrVYcF48D z?{c_24x5+%?^q60t@Mm()^0y*1%LWuL`OcN90|&C9-JNTXubbJG z&4%-2B|3`5km|xVn07^97wb8H!lliMPR$ylW&y2w*G6$<(8_S#bWIX<{)9Gup0jIe7tT4!yBKx9BY$YpK0#hCYO3VojSnRvNmwH{Z-Lk zL1X%!E1Fez-otbH^}__O>q)AuzLPJubyU6`ss^_!TP3xugCg!8oiJB$$S`y;)!lw% z8#dzl46WWwIz`tmX)E;}pCK?h)95*AuE|pZjr5>S8tKNw2b;ZX?rqmM92X;tVs*u} zVl+xvk7xy?MO5?BUmAZrHaDwe&+NXi8$kTvrA8N8H!Q~crJaRwAc`-pgwdt7q2XoQppIy;btgT8*SNmqueYny8gZtEzaqyt+=d?pAk0%q4^Mg>NqpSrq0+clS{JwwNWH&can}Z zGLNxj2*O*+k*k&ER}*PM>GJ5g$nfCAj8&tjeIdF+ofP%Rl0aq*USw!>KgA;Dqn7<1 z#0Sp~S~0#tDD}o=lH$n3IU|@rMTF5eMsxS8Lzb`9a+d{{vig1bhBd+dmCw#_ie&=l zhebXuMIY&2L86W*4N9fiiRL`TbHf)e>1jUBGFze!)<)Y9S97@%t-?4&om+ZaWZ7LGVtHvRESpf5f#_v;enI@&rA0E+{ltscCdnEaUY0`=mPIpMM5NxR zNwAw}uZn_ROk7rW4;d(?I<(mZX}F76hzL*ipg(YFQZ72P&VObrb%SfpGuK31?mKGe zdfUnXf>^d_(HFOF%B8BrePir&4NTzc?U|9kDHq9Q@TXw>Ms z;z<4MhBMOgiEZVq$ zo5g5ynJB-}_pO)wI&;6M-QgU&y6We}Z>_9{Zq^97sps~oV;ymW;}Uu3y)ki@x>zxD zMaRgvLXOX_R(h^w+VbWNq6JPN$d;G23xlh0s~y`56WS zgPDQO&el^X>dh&RulgU6nho~aU(34D*8fWj%&TI}N^L05@~vC>fw3QB zFgf!_29fNkxm}G%kxar6a88{g4&ehZZn^1xorkEiz%r*eH8rouNNs*-xEiWvb8*u^ z3avG1Z+HEztCkLNc5M^RhP8Vix+X2%P3~3L%Qa)ObId?4tXllJ@5~M6m{tHr^H~%bp>Z>NvH{Y+D7&jGN&CTetPT?{wED zJoy;#`s%<>`k!7rcX;o-xP>RDQtG8AGc;C!cG4(1be{EcaYsm{n62yBWL?Sg9=zqU zYS{dV=?eFg>9WDd(v+&HX&G3dYpjUIW`A4NEEaXoNj4+nc1YZv|Wg` zk7|?Ug*)|))qM*I2AvLbo{Vj;cPllS(Be8k(sTiDH(p1?&V@}pX2Kb;@U|d)#!Y0t<$1BN~H~p4R!OUQZ2PHT^@*% za<|S@^d`U%D}8v(@b8%6gEv&`@2(O&ePaT#K;LRU5C2Pc-`a zG3ni~pw*DX-iSd&`HRR0)1EoSbb_tUtOw{y%WGw7#B6RXmG}lyIJnW1FI{|hMxtt! zG;-NWx}P=fSl?|fQ?xILD|DXRpM3Sop;)`uLXt>d|KsqLxF{TEy6tKWT^I9S~oI=8fjb%nQ*9-#N)a-6D{Ckpqh!-q)))orcWzub;7LFpSy zHL%{uxX-jUw@IgG)lm_LA5~sUAd%EYuWfi`g2A??X`JtV+P2TwNxfy+6?{JDs=MDb zrtiAgJz7Rj7dV0043(Y?djz@kMVy=wzetZ#ZTyImD!7#t^OGJKWPI6AW&zqjYg&MI zw*Q!a$SXsT!ii&P5zVMGRY>k=bl;vYI?rD9;qG?tNpf!}$7NIc+EgQ+im>_Hfe8~^ z&Z>D8`j5+IekvDweRx6YwV{`~yODX%FUNy8C#;YtfzrDXXY$r5#~iq(a*UbY!<#Y^ z_99$Bw$Co%X0A*w&t~NTge*WdpL>l%g9&^r=~Di42Jg#ovQR##l7zU=yXis~A&@zN01)tI@!fC;rX?P3|enWAq&*ekaL9fyEg4 zwnYS|X}?0_eU$I}xg{4j+U={LmoB*Y*+s3s?139NWh?yf*BS{%pR}3=1C;j3_s~a- z#H5rB+fipxH};)}(HSDQ)~+6W*-4C7Q{)lY53ZR%7BkN}ma^v6Y!^K(upb_d;jmFJ zDmlHYVObN@5)&(GT1N)OTYA~)ZEJ{9@Y`HUJ)w6Y9vOK~$9WAGRR!Zj%%FQ?=esB; zyPpd0M^V8U>6zd6lg**JuF~|DTj6!uCi{D{{X*cA9F}Xkem+s3gPn>; zF7lNarir7d-;zkK50#~k8GJ<-?|rr;Pab~fzBEhm3Nn~bvYz`T&*9cnZOOS+|EAis zihRT%%J;iN%c*^nmYRE8k8GXoxHn6nO{2$_el+^W>b2&AJ7t=Cht`jc^(~Wk3VahK zc%EoJSq;~EuU*EEJ}&kU5YBrjICFA#>x}VWgMX=1t;m&1< z?ptNs?q&?F@SZOc-4Ntx+?Ut?GHv9mw__bfT<>5l9`{$nlo^bD?5AN5^ zU6U*(?|L6NDSMC|x6bE$y@h_ns>@qDd)PhVrhg>gW%M*?$=)-FCTyg`BOD8ND?M7feV4%JfcwGfpnYk6 z&TDq{X`LwfH&RbKy~fPGi+9Blj(rwjj9-nt7E5tQ;IXm}<2^+M#tOExmIV_-g|Ncp z8_~&TZ=t!iuqA&hjCiePO%S6@r<}r)vlv=hOI7y=OtYb z?p+>VWu}8hSyX(7rIkT-FLRd2J&#|LPx>&l>G1f{eulGrm(iWiBwT0wDMCryFW%I( z9D!IaTxr-PtDW4uMoyH_8+C2hNqUR15km0qr-Mw0-Tv;0m;>G zKD{I(7Cv?e=0dH;-uHFQ=$8&)+YKY1ivjYuOV)Dzrc=IyuMU-!s%438zY|_yPIqd5 z@_`6ncGR3-9x9b56ZLT$5{Yb0wSH;jPRjcVcW9s5J(=b$DcO}|g|SDOt;5x2z9q_I z{0Spg0%DsI7UrjUFfzV?XW;3I=9_V(iMT;wH;Cri6!ce@%64rXzo0$2?aJyu-}F=J z^qF1gk57EG*>Lk$oGi!FK67PCqcMpQQzuEBvlmwDKi>`%k~!1Az0`=)mRw0x?^D)Y z<(nL5g7o8tYR5Ga_7}zTtSR&P>}EbLUb;odut0Yh@WAK;#4Gz?AceF~mrPZ?ABUFV z5VjWVGU%?=sPLI!7PXAl&e;s8P~skp@D#v(r8}Ic7VJoFry*0`Q6?fuSD$Wl_<;-8 zT>NH7?iZtk%jO|{zv8Rse+(U-q>U2a7A{t7GzR$wB78k1?jF-WsYnpuQ_7r?a3j*K z%2;D23DR(|{P9e#`$N$C0`&b`5;>Byf?-Wc7Q(kAe57>1w=Mj|&_x4kn|O);Im=sSJgD>#~IUoT-g-(uG_a6bMYy`X-dW9u_(8hCjNF z(yaR|{OuhC?l;);feN9n*dd-MFeWSz8d4Dvr|~KCt-MOtr7Hzj zt}k6oAlsKA*v3t|-3-QkhT1$$cqN|_j`vZ{W^g-7KRYzZGX3@_`hI_!mE}nI=MI&? zkdq|0jB#56US6-4hV8j)rZ=||)eByxFW>Xgk~(P(*`WIDnKW*lRr_dj>)YcrmVz%` zt*)`b@i~@wHO-?!3(OCeNYi8jH<~@jkT|ihGu9e1v$ao4cd9QhoIqtR?}LPtT5?j! ztLg9VJtrb$z%@!B&p~&2TYXn13DM>ZyZXqS^^*X92$wZmmb>Vaw|AZu(;eI_c``>I zZg`(yjNrc6;;%sYMcX9S8mgMACTlRyxA-rbSJJC;rQ%2<9lbUCI2=8E&a4-2l+3A3 zuv#@$^PnFE89~n@RNh=o?JMc*`JRyL)^c zQ%EgSUjckn_>=Y@2zsY$~VV}6jsJ27MBrE z736$R1jz&fesO8Tg&#j^SA3xI_KTggW4Iarar2#;MKST?)y)+}RrGy?Oa{20|pl(EgUNV)@3R=4ZY1^*pG z-O$Xqm9d}JcQPkK^2E^ke0yhCZ_beH`Nyt0SK16c9CeX3Y_=NdAWSnwp4sy87?)BB zL6C|&ca_qIog(jCdmV~ec%YkF6w^mP#J!ZZVw2OE;+<{FEI;at>L1wrJRT^;fGA4@ za=Xw#ox2*gmzwh3_?X*tjTVAyTtf%l{N$M!GN=ega){*0+lh&E4Eqv5&XM~S60kx{x)M6E?DWwTh< z>2kU81+7oWN0Z_@-Z!|&i!O7yCC|3|ej7Silw+psh3>s*qh!Z&fMFdhv=a zjQ$Yfn7jgq@0*Ou91ut~_t9$kG@Fgy;1+Zjj-Y<{8ezVC16?ImJ|B_@IkIx>3!XRi zPh`q*NvmzwbUHL@r=yg-_v05K4Cnss##Opl#w&9@;qjD58tusIF?mEB&SZm3IC@SC zqeA{FgD**qi$boZz7&12pElOz`fSy*Sck7r9(j#X3P;P__Y#6u!!ZoW&Cyy|9Q9~B zxI^>pTN3r8>}zQGR~#H)S8J-}q%Y-Gn)_WZCH@ua*_Jo<1z${B4!~^qjm)_hr4=T(^LfZZbhSKh8t*+@c1-2MN=^ z@CI5r%N3lVDGV&DD+w2fyHE~=D7TMyL$q}~FWriCXUCzi(_oD2NX_fmjUe86lGo8I zLF07*F30TB-G`PlkUl5@^5dj!ICEJK~rcf~d?- z^AtK1kQI2VJ&lZmZ%n?{(Az=d5e+*+car-12S1!6w@FMDAc5dOOn(&VEuY&su4r&@)#yW+S z3G+ibQmRoI-QsgKgn5}!K`gpHc#Vx6-T7Aso)TX-aan~E-!3I!_c0bxXS@vot#$OZ zPvH8Ds#3KbKnAJEtP)hDqZ)C zOT8vNfmV0zNsnG$ks6n+IB4#(1R0&^ z(q*(U$rYn*LhF<{j#V4YBN``8T-A0nmc6dWs;VY6*kze9fH_$;REdoG7Z-Cs$d62# zwXd=bS4gxz883P4*UELu6=t%vm8~;NpV7!f+etP!@2aZl;SJMR}bvs71m@5f5t=!lj1-cQZz>Mg@%c;;33 zN`w}7!$-8c7jKmn`B2&Gd04-+Ny2TAe{%iIm1@EZC5mwa;^de6^o-tT2tBVe>%B`a zu+x-x_yH}XJ?WY@N9cdVSLYMef!lO}Y`%v_HB#vWL7Vgu0u|8xF1C;Bu0!>bs4#xM zVavRn=JKdL?ufO|9p$^UCJ#r|s^YG{4Gcuhbq7L>6xCd~mgZN9#yBoW_-4{PnT$O( zv@b8HjrX@<)BAxF@`JKEQ=OV7ZBdl-wuTu` z#?0TD8gJS|Hkapk4&U3wqDY6#nqo4GYv@o&X4H>6jnupzFZ*>1|Idoy5{3yuE*YVe zdq;a{>q?NpacWU| z{l!a#Tsx@?V(iGBDTtqaV8%zzu0l5~XkRZ%{(X-yguT{9SPgWfy z)-tc9w#CVaBYB!o|G@xMhN|B>aebEx8g7?LXz4|23=2i!kRwt=uCr+Mc z-7*Z5o|40JZyqteK$kM4jOU#tJ=Cm(bI_huNii!^m?$5-<{&8$pHwLlvn7{!t$EDH zJ6ouqS1DLCJ@Mu}Gs~URlV*~N43_Ktl7#4A6Ov`pJku4i4>RMF;!1zm*qywvQ_xh~ z(jDq6P0JY(v+nkv+NG$W$a zJU+bJEbfj#%Wp3Et>>?$|2$mij*4e-AX;bHq)5=;2~tL>zVV$KrMj2c@Atje<%?&* zc0+?WMRp*$MKmQz6#`9)&Z!D9;aL?$GIBq=xBB*dI4zmoFJ@Use9j=$T8tvm>hq|2 zqtnjq3eHx)3XY!)P~0|wOBIey z=yq+pBep#0j%4WeWQW64xoyHs((a9(hvdrSI`Rqq~|m-5D~ zf}T~-3WEVz7nx~H@Oo7isyj6zLpH=|vl*oqjlb^G_{tIYR%vMD;{9%p2LXBk=pqWp z?&sGRk+h6EI4AnXCuq?)b;Y68nL8OTp2P+vYWJjyLg@9qo-SzjUA@bDO~^N?pSa%S zO$@c^52{sK&Jb|u?s%(P86>Tib5FS8=`!)#&@UI;>%KhK36yr#yDrU)dj(~@=Wi5p zE$1|VPoK!2#E5f4+l^&e{;ABa+HAyFRu}`Q-01p>)^@kW!*)t^Lx)@5Xaft9C!2eQ zFLuN?a6gYPt1RRG zb~CAx_WqT;jV>pOYRBo7CPqtNF?3S}=duR5IRv`0Ln-lbGZukZzev;TyT_+km#!qA z@|ZxSnveb>(^b=YbP1K!&*HET`LzQ_vEM&8zLE9cVBabY-%;s&ICm%BuB_kxGHxD~ ziS=M^;oM2Km4v^iTR@A9gP6_4FQW7jg%<)V6a(V;v3p*Wx(PG2&$=4O>jY@3Hl-%2 zrT4?H-%$wLzv6O6&GCr_S$(&v?56?wi0+gFO8a!e!Q0z7yB+Oo#m3mJYptoutDlf1^KRFOvkS%X%4$ITt8+f2|OA7UarQG_+vpMU6#Kq>p z?3kg@A4;#9!cMJ;rca}W&-N(xyxEArFyi70QGbn*P; zRQ5g7)z^FL6CwrM+J{Yt{!2HA9`4 zhvVNOXE9cnpZT+w+m-F75!H281#d&|M79E*q-XJ^u>%U}d_?-ScGdyxR-@VGv6Hh9 zbJ{8cXwf~fc3SL~JsO>{z~A!nY7Vi`-Q?%z*Jm9Z(#urIFG*x%bQ_}r5lHBrap|j!}C-vx95dRAZvaPZQs|Ez)RLo z`2}p-vRgAdaW7Onn~}ETT%};4m=x2d`2q9nTH zI3(C}C55F|y7)XeX?g?BI&)d9{SRjk$uC-=)P&wPqtQkX3FYScY#&G3Pa`gt_frNv zS|^h}F4jo5l)qHe54ym1GIeIO1=mLh-*sv;aEb!+ z{YvzO9R4sBgzXN-Y}bXl=4$S{=m z=U&feKcbmUq^^9?4R8n%W(E5 zdB0DuzkzLM_f+$Vtf-;Y6%T`8a0l6)d!(`%E?EQjg4+u@p4w^&jHe3=%05u@iceUl z_!6ii8Xd?jzcFApyu~L*Vu8FYCwk3t4vZrlBWd@1L@Gpklf{Tq?#Pu8>f!MCNwpI% zE`30Y$yMf6YJ`n7ImdMIbrNDD7BqW=YJSv4JFN}HDz#XikaLeH`S zQab#0O7;iAFF-MS+;1j+@;R0~;8MN)Ev1&u4q`Aw=V>D%Frr+yBlaSzPZM2GdL_bV zccrShM{YAsj*++G+Rx>pNKtS4t!qp!m+BGEB0tF0Aqu%polI=0M*7+6WP3_%ESqNM z7ZV7Mbv}O^k2lN|Z(t)ad(M-Zsn%z!Ee+k-HaK$m- zWRG_Of<&t@fm;Vi2j^+)nvT)Xz>^ll5167_3wtlgK*c^Ru#D^ zKc-|unE2g7ueQ|qN62SV3bc!rEz!|2wqBPdm$BzoOBCeo&W=s<5*_gh3?-aEszlL- zM6SJBZVX5XI{=g6>pXDqi{dg;{o3S3m42`#$}aKon9|sBWJ+7}YvdTZszGVGFjqp! zvT4+2N^w&63N`UV8|Sx&BLS9^7GAUwc^=jK%S)R*MqFhJvc`431YS(l>8Q9#l0gZy z-`*6yc4vjyliBn|b)#UPohjR*nycJ7lgfGaoZ+L!s!ii{d7HBZbe1CHpV9DIw?rwRO}|9 zae8fJs|2pavyc3Vt28^yo7kTS#9G=C2Y48F43IPbm@xZ- znrai)W&A)k68(%~ajqKTe|;nRU7HZ^CwFOoPKjdK8w>-#HZPZLq2D0i;_O^bb0U~O zTY0kh-gvjS97l!>bTiI?M^J;Oikuc@K-M%T;&`!kw9d|Ki30NNSKURwp;ptaoa(HF z=(kMH6!^BE3U~O}V5$K`iFM-DzRJUWbKald6Q!a)H_ZDjeZ{*u?nU3yz4_>p5+z%I zX*)&6FM6l4t+1oj*-|w1YXw8C`}8cKArzOMxR+P*H(olHrpxm!Q61DXE&I%lCggS? zQCZ@Ci_EH3NlmD_v2T#M^5VIrYekU1eb$&;UfVK%LytMb2|o7#xf^bylYs7l%MSa&fN>zE6> z6Xz}}Izp>blhA4T5^eI*20S%Y3w^raSMYAQ)%#1RP+X$A%}0r@G94Kcf88+1Q}U=U zg1_cpI=%$s*?Yj4o#|_M!d?5-ZtbSbu;;N&&w0O|jM5@~{`%_$ zD?ONNb0(LJ;fstE)4?6lnb%7*tfj@pA-cc_Ipxr^kocliYXffc!_B%x)z*%w@dbkN|xWE)@0z<$p zOi%~`5X=lL00vtCKLH;BAut$b4l;#+Q5L8KUT3_r{~PD2Gu~}dxBwWY0R{;Rp@hMZ z1SJ=|n;4|WE_i|exBK{A@#Ow*oC&UY9RDZIy7#Vlj(B`wEiDWILBOEDITCz4@Y2yRIEwY$G$2+7R^s!2!Emq;48VlvpTB4zB(MdC zW6%&N5Crx+4Gcp-P}pY?NQ^gt`(R+9^RwuGKLP~Si(RdFyz1ZfcHTlBvc57p$R~PB7hP4{PWL+P#6OHLO247#fO4p2?GV}g~f+N zBd|1vfx*}`2(~6LAt9_s;6f;@vv9!Wm=DALQUv`k?O^c#U_`*NtN=s8u^xe;u&n?G zqyF&mzbyzC0{wR1x&4s{#BV-1w>c6HsD~*v0*+0C!_Rjb=bwMv1V@0d0uGM+GvI#X zLxIs)njk>f7b3tAY(5C&-+^>40TLmEEg%B=A2irs!HWJzYNXKlaCZLrhX!O0mKBg7 zU=G9J1C+xO5DCF{Cy-3P3(W{R&sh-mVK@kuz~_s{iT)ddO%_S^0j|do8U(Zy`0w`t_D2Dh#%6>GVbkCUY#I`)7N9^VD5gtLU?4;=6al&@7>vdu z3PG_R0i29YgJWw4M!>Q80Doe>0RjCjAo_o88H57DA4@{WAA8>O0QtLSpde5vrkhYe z)5DB(6tD}zNLk%~10mA>z2h=Jd%$Nm~!7iXEA>?04@Q-nY zP=6L=AT0kO0=NZBcqjy`Ap%i{QGJ1e2ti_c0gwc{rU420cSxK^#2;k|1%+c<9E!w( z4=5-KhUI>s6JWR+1aLvH`5?a+M8K#(VaLWTD3}m-z{7;F(jElxVYv>4fMF#84Ei5J z{|`0-`u}Bqz=Qx)!|<_?5F9)0VQ?%Qih?09dk5221?xU2rZ4~>{BPGK z5Wm19MuV`BItnm6R(K)c|8YMYtBoK*Fig;M9_d(s={)=X@W#K(E`SuU#y|jn!7#0i z0#F`i!2xL4l@kfVif4f9w{-ra>;g3G28IM(QE_QhKVQv8Z1FI4A2pd*@J*z!LepW07qa( zGr$KVFxFGZzl+dcJ1QUwFk=8X3)IQ;v;Qg*fFD7BHITmn13<&<3jiO`%`hJZ@IQvX z|A-4LX94i<-vKED$N(@xF+Bs&U|6LXpkXCCfIndH1lw?6`o`+gfP;QVBEVz7h5_3d zz_fr-1AzCjnjb(z{fCl&&9VRFJAf+y8~tt}Aael3_&e>K4{J_=3IXf6-}!*}5yISU z0B6B}Ek^$Imd?2)Q84Jg&3kTLU@*gu2?zwU5djxMFf9%U4M!pqjwA5c0HB2gmJb+$ z!5~fD-u z5Ljp$zc0F1$=VcQIb z{bdIr%ySRC|m|GA7-vI-wEn*H1nR~Lt F{|`$TZixT@ delta 21444 zcma)kbyQXB*0+QbA}I~h2Dv+E>6R1$=|;MHgQSRbBPmFObPGr+-4ZsEf=CMpe#?80 z-gEXoV|*|E0n9zu^UR*VITyq6a?;Fl(ztDmvJ*3s^7&TP2bad;JLKr@{_68l=kQz= zKIcc|M|9R_ao)}#%E&$4KMnt^)T;kV*6+b3>3JvX=RR?d8IMC;Qlcl{pZWMPJon=n z7xiZnCwun*dH1>=A}It)g<;zo8eM+-VTgdccs=9H#SbL$rZ^BnUEvlpG<-T1`N~@m zfzbqiiV}D{JZ$>ivJLa|cvFV_g?M!1WEFBWg3T#Bl|Y5=N16~pCOznYY(%Dphlzx9 z#>b68w9sCRF`m-?>MrWLuLhvhSO>%@)lT^rssBSk+npq($5lSmzAsB&ygZ%dEfk)5p( zZgW!9DwfI2PO1rO#M3aWCue1acCs^*D;Lzv8kA|?tUoVq#gWd&uLM)vnF z@;eu9eyTT0-VvNz)zNY%xPtIiwQEnTQgGOn3jNk>upz~_nxSjpcR$*soOJ%QQHQ%c#_aytqTt?#$%G>;aKgL^0 z%K0*fLb`Pn>lzP}sM05t?&nTHy={UX$Mr+G)034-pXzopTGskoEuGzCd`D4HTURL5 zGT5AJ)ZFkThIpYkRPbKS_J#v>K%ZXT;%K5gDi3-6#}tP17WZ&; z-)h!yl+VP*ZR>;?Joa6dLf`yT*0(AH(MjSdNe{u|s)M)2!g)4RY*;xbxb-V8G6^@z z(YaqkRn|ZC64@tX@VzU3&{Hm%b5}2dr~hqPaV?%vAL7)}++e5oXAkmL@Cw;01L%G&sHC*0u2>ZXnR7WVYdU=30Vk|0y64~HT%6f@9MLE zCLvH|(L+&Z4c-s0SHl3qgWo$0id(j-2 zy%@e4i@7x$gHB#H@H_QJhHO@RE>_Er(Y)bIDp2#mVE2pP_oXLZ(NnzlWM$S{3|_#_ zVkGNdFOrH*ypj5ERcizpOF@A8Zj+^~gi0n#h-n{|g5t2Xadigwv zVm|ry!bUkaZegC@ z@}zgAyPXMYe`nZ9#-5x}>1jsy?oQd~y*M_aw7_?-^rWAewsaZaB>D0*$CG}bnXxU0 zd%rLVx}KAY8F3q{-Bp`laqr9%>JV-HvFJ5P|9zalU2F0!`3w8bvRU-^++uk!k*##v zR~LH;gE5fdVp3hIX2jZ!(>i@Mp=~aVrR3Vdk<@Oz5phV4x%65&8HzlBtGwF7QFBzu zAh($j=|{W1$9-V=IK?-^3&ecm{=60kgA80kw8$I?A;7t~_UVQlI z(!3)^dh*$ff=GL>pnvXuO+=HM6&P`$LKT!ZpM;b!~oNSxB%od6e z=>g1~b2J(-fy4cetoyvqR>2iwRfEAf`h#NiGJGqZq*CwHMDNGxQgf?uGc#lEk1m*P zUpAFc;^r%)aP$miE#y~Yymht&zp@{oA9FGLo#E%yNV~_uaV4ro)t*#OyyY$tK3@}B z7+rdrUZWoix}ExC>_Z>fbhyIqYa1d+Gqe^S_ zp2zee+!r&u9UAhNuWfo{`k<}%xv2+sdSFac*TSDGk-D2o^xU{Z|NVv0Rdtex>b#)R zuCwQ+N&D#FAG`-yon)aDpxmm z>J+y#ayD&ueSfG>1h1^#qon>YesA>q+mynKgzrt2AKdN@aJZ-F_rM8L60;WWOBkJE zO5@oy2rXXifAYn{KW9#S-#~w3NR&;M=Yc6f!0JQyq4Bx5-wsW-i>?-bG+5F;%(&#C zYvnOd=~R8g(D)ARV*FGyv=i;I?t=GL07v`k7D_j%NxzLl}N z@|?5&W=U+fi>~w)ad8U;)~?KmHOBbzDSI$Grd`_Dy7OHBM~&-M?Ho55#ZtTW#bB*u zZ%q%*=wf-ln6M1!a*dg~EPrbn*G0D4Ac|OGE9_gIc>9tr+W)wx&ymf4uOAb^-KXdd z@~*5JxFy#9I=bsyjr-zL!z6>?r<>+psChDKlX&iF=QsIYTbx>vR}qzjniz{NWwfP& zw`Is8(bWC6i1b|X;jw6f-m$Ou7l{gTFVuCjNXBaT+1%X;t@FuHF+-)LUUI#|T*q)2 z?M+T6aP>AvOl!zRP?}HZvo$VE=Z9+A8unSnaJz7{+NRMn)%QZLju?qPx7=Ke7qThe z81`17k0ir05x1z>)T9Wj=gho9uA-&wP0ZtbyY+`@M5SsXn;7ZN1wKpp8fZ7VudYh|#WuN#m@Zhl;O68m+&`)9PAFD=Hm;zn&jtR$ZCtIY3X zc$7Bt83Ww}Cq3EobneUnQsdQ69?iU*t+(=f6k>c9B&yGK-$Uf~f^abAMbgb#@;x5u z-kjF~24$e!S80NL*49p!etxxOkm_53`+zAs{l8?7D?okVqwgy__5C~t8CaaNG#wQ9~q792E8!}lw$LVM1()^ z^CNY5UQTKIM(&{{+PWe}NZ%v`1iWQ>SHV0$sdzt)R>2+Iq1eQ)cGYx_VK9ulVOw#E zCRp3#*ye*F<&BP>w};-8MV94spvvRukm9`!@0yusCuu+8B6MB^uS9BRWkb1day};z z6K^;P*tUP>@{F;Nf!FXJ=0#QpWK^MiAi>J%Hknt!t-cg11-GbHbw?4eBE#eiXbpFX z@SF6aQhptZB-(Dt5+|ZyL9LS7LI; zfJ=6^2m@DHrIttLdc<7M2T|f*Q>BmUv~K8+dK!#H(?C17B}OuYNdkDL<>SUpVm{xh zaS!EM#qYCIAAB~uBI@%vRCugA)GAZgbYM=T5p1iJ9!&U|mg58EEF*DpU;L)|8x1b9 ztO@f&>4I~`)~Z6tN2wTYVHuyLDwuvJYy=YOInJD!U)E#Vv6--_R`BIJNky^)*Q2r_ zGEW1_M_lK;6x(O&sQJ0@{N!_v=01&oI`R@ZyTp2z*@CG{|9z~32vzCx{zuM7Kb#&{ z*lE%9aZb>BdB0_qJ+r$QHy@gg!h_m&)ytI8+88@Mx$~Xd*>$BElg@Ncwy{zG-Y^!o z=wlTp;DEm;r`%f7h<2QF5PTxvd4aX1-QQ_dx00FDSo=`XM$*n&ilRp9$~&DFB)4O^ z*Su{)Sf?7ML^mg5E{KA;XFJLF`GgDfHe>&fAVFRHt0W6v5_R;a9QLJT(yu0HLJK6i z)85o_?WCEb3?+3gVb%ge-$e^f@)Y0Iyvr&V7;IHH6alWc>|A%OWAu96H^R-!-+Qg= z;F7d3%)aq`m5GQM$0Ow~yuOi@P|9Zao-gVUYvN6raj{F$B*m9nic(xiM##+f)6QS9 ziTKeJWXzhnA$qGhE6?~9J9S<69W!s2L6*on9g(GIrUisQ1mQ!wvMb-;1 zcpfxOWE3|U@G>yIJW^Z~qCMBG3~jKERkUy(an{Auk0e6Dx3lrrcx^!pM+KrPO!Gzp z-FNo{S4zhQ2hG*$TBxTyJ?&UV=00UOq+AkfFT-5F+)kj8WD+`0sG&qRFT24@Nq0K+ zB6s;F(_?uTNxp3iN5{r=gnhpeY!)hI47$;BIWfI1)Mq7@r^033sg!rGMpj;}JP@gQ z|0CXFqPzUp+X+P~w4iVZK0F|%RRB~|U3}oR9yJh}+rCMOiZWglS^FgF^47R2tTSc${rq+hJTvo$N@Ved%W^z64CoeS1U0$WvZb|6U@rIvg8#%PkK3*h195?IXEu$3A%m^bBT^G3%rq?{v~e~e zlRJxL{oy!S;bFKsQCRpPe?}r+fWdw5 z{TA;-Wh91ibs%Ai?L4n$DW|^4b@fx;WMh#BKUk)%eUmmIAoF?yX!!!f+(DyUQbqYuftw?l$qB<<*s} zo^;jykNEPx!%f?-k20{stuE%L=p$%6B;>7%!*0kl3YTX$l}v-FJ@LuPRh+0 zV`v;Paj??Ln#2(NlPFBJi~V52>tddtRIFXSU&gnKhgB!fDHOMDEZWtiylt9i;v?ZB z9bV#;8f8jcJNaV&LV76W@a~h+YQYEO%&BbKuFs}A`cAFi-P0LVR1+c=IW=D~7g&)? z|2&cXe6?Mi!Z1ymtNooINhGbDY3o$2LXp_+JVL@vW~Ku zw&fUyHwQ*M3nsh7xlm}tdHyd`5%kZ>pC->w)Lj=x`$0Oa5N0IZvtiS;+7EOJNU8l5 zrmyJsGVMeBduN<#KZTXqcN)C^^tmxac>mJoWeS~0*{0~1ljn5%8MegH1v>^wEcobqf&vDJ28Qe~`$V+N)63e=?2mt#aK&!yt4a-`kgPM> z>T5L*f|pD;o<)3orsBN&HCMGBsbBieD3wPm`&;zDsy;_oV$Vi5OGo=!u~!Xu^qG{< z;{)ec?als6-EJ?PS}(HP_Q`t4bFk$5rfcfwCD|LbSLoTVl)SuEDB7za??$)3;Fhw= zD2I9Rx{5QCXEm{7%uY6k5YjwdL3m;KPS^n%aUV$~Uyy6g(0)t^A4%ZLGWgYN#y!VE z1%-w-#zpnEVLfbNg)Q&RkDnGuff(EhGV3idPrg!{%X;zLyd9@#w=vynvyj$S#C$83 zyy0+_V`4+x?nBB9*~Tn8EQIZ+J}NT!V=swS$f0i7N@(5^86^X&T;-SK34n2(Gam|J*@d zCRMd)6*(h{nOe;L3D)tSfNylY2kYG6QC1c+9X<>!Qd_MVl`Ul7Ta$3`3k=te(o`*B z_u|Z&CR>+s-HEQPSK#`_Ts&gnvcx@8jd43LmQWeu0+b_$gz#7v@akN2U?yKN(c@X7 zW^`HXi8rKGR-jWJmg5`JU;GnFSJ;U>wp&+T=t$qISwDz;J2BMukpH>2V0 zwj*;-ymmL87U~m9FUlcW3^PYgD77c+@2&X$=tUl`idp(rWOUquk}V2P`b)g+CR9R~ z3&V`-@H>1hj~0yU?~UJB^6Ogud3f~xnX&beXR>Bi;_UX(`G)JCO7CMRpSnz1`l@=% z?j-BU-6K{!Pp!GtVw4b<@w0wKVw|OvSCdtED4Pv4CffDMVXP0KJF9WczP4~l>7^Ee z-bLb_Wc2W&r$Flm{LsCqg53IXcg}p8cGAGJmhEt3$DO4wo)sdmiP~ykXHO4uP7!@cZ(En1=U)A2ilanXAxg4a31?5h@DzV%ArwDByviWc9dim9^@tZjb~o`}R-R42A6UeyJh>nK^GYLo zyYfnI@>)Ud@Wb`K?*S$> z!h&>|5N_xNlwb^f%@*vdIG+2Whwh!lvbv0*#2-vF=cGR-k@{dTa$I-#hq=cT>2eTT z&vIR04-2c-{0qm;K(=BTbt>`$WL{QYP?_SN(bA9$=EBCPBkQ@xuVnImwC@i`9I?NA zne6qf!s^glp-IJa2xMy}e;}w#rX%P+VbG}{3blS}i#OCkPYycQN#{fbiG8hX6oi;G zZ(UQUpHRKfu_=r3c&KMUG|6{PhV;9di2jw@clWH?o%HkJO2iX0YhS;G(ukAYVgN-6 zn%~T4Db0PXpfGb4aji4wbvN zW=aI%`yy1~atALux;bXkqr|?Gj6s`4ah6=|6Xa)Lvaz3Z2S(b|uUE`#&471i@T*Zn z51GZN$95HrS5aH=Uevk73-_;t(BOYBlF@q2W~Js4*~%a{cAj`9Nwss@Pb~ZGWclUp zHR-(O(kZb!x9e3Do%F{OiQson573>f4NXKp@Gquz(bQYNvsa&Sjm*%$P~6#6rFbpr zUS;issxb292vrA$$&x7AP|DV-?H-BAv~OP$r@x=c6W;qgfp4j5I`;`%DI`d)pSW;b zl_hGmxX!12&{(k4=kQj}F_WbE|GMM~8JX=}Wm#bI0`Z26VdHx{3u+-=F0F5Eb4zi$jHEuWI^%jY5Nw>oV%BQK&t6O_tIT4yY42qnV)^r5s@M)s8$vV8fSOk zos{8Uw};(x!pP@wOyM=`ze?gTMRf_XHS|D4DVLs}=e$l6!uD-cl1MvF?aH8R$BK7x z|LkQg?fKTaE352uk904nm|9Tf{yZQ0{@#wihwJFqr@aN=Mf03rebTenJU_fGP^HxT zd01^rV%V?1!+?|ieuFZfdIr&_eM=~Z*g*Exj{0=2p~z$m6-6%Zcze@_guWkRD)D<% z*~<~0zP`ey-NtPr{WX?gQ7(Czc^BSX?gXf9!}_LYqIn49$Ln??dzq$dU18^vGe4L| z*1L;ugJlSa4<`62y*@1k*ZAbQGq}WSzQ2#`o%Yo#+`XI|EIU?nMK=|0VJjJs^XOvI z#kIU)hRyUBm=4jr8+WZ;&go^FOX{!FQz*!_Io?UaKG*4{aB(*CYd^cv%Kc}OxXSGKj(uzqma<*m>i&=GuZ!+yjwM}nb zv3~!!{q7B=F8LgJ`kF>(zd+0Kr}k*eBZ+n|&y{zt%aB!-Y&ma_UMoF*dz5f?xOO`F zt-pO#m+kS}_Qqo?@V@b3>NmCz1Zi>}p5Ic_BX!+YFLTF?6LRK92{rtPC-Up1jpcZt zEGZUVVz}jvv3HGbil#Douw`^NcVU3;)3GgC!m6YI`UZu3*QY7wi(KhUANh((;ogQ?M&bFZ^_?vwHY?vVNKPnC z8FDzuipN4Z?wFp}{K&s9!vZ#mUz6h}yy)o^7So>>h-AK1=p6{`z-jV4VzJt%UYQixcP2Or7Qh*W_d^empVccU|Nmx0#?XC$= zFm&hf=4z(4JYV-QtC`y2(%lE0q8IW}efyuDnG)U}qSoLy7)Z#YpR`O~3{``Mcnb+K zX(ft2nqq&lH6}}7E}%-FcQk1TOD}jzGd|qt{U|OEbFL-TXYw@}@0XILha1%48x4Fa z`@XNra#AI8x5w>74$(KLQBT+v$p{G`2vW+!&Tb|-egW1dtu_@cLuc-n9OaUU(BsWrFgx92Da`eo+oh5+db*# z4r6pXH8X}wR7QwxFNF5{#5dpdeg@8MxT`h)g`$#D@0uLD;!_r> zZ1ah_{;kn%3Lf2pgRi;+VT&AF*@p|3YN)+jMG=p(*@={>8lUp-0z04bq}Y!0mPzhX zc&4oTX%041>b(B^)tcw>tO|wsNu}L!EQ(V#8`X7ERx}thb({4@=le>Sj;+TIk?X{R zkr?3ZKmSQ41^rWS)DL)$HHQ$rG(B4ury*%|62riX)>h2Cc5oy0fDe?f{HAo$dx5LN zJBChN3#_&FmRj9-QCCEy?oejDq4!oF4mfIQJ2;f%5yz973Nv5p*#lK|b6bE$T-=7K zSiE@OgFnChcm$h?zb-ke;b2ecmO}jn#@oi0f}y@fpz|`Ig82oDqPGGo&B?_RZ}hVo zli;){VMLPYPWqv*jLnF1M*I<>{EOIifri}ZEEngqi(DHso2xss%T1{5v*o7xF-JNZ zok>~#V_}>6H=Pu(XTmR}%9iKNymT#vyza64KrkARRrAHu?qS-swAnkhwMN4ib?Y0F z?HVz(c`yb%GT+;Au?wU>yuhbt=*^2)OXvn-fwkb!1;(Ua+G#DYIi5wM$iqjy&77$W zhAkRAh+Dm<6NM`i%4==?B{^27O^PmAp{9lFMy?{}tg`69%RlX3Uuk9(D7J*(G_hUCK6uL*p+PW05`T8yX3Vz)W}<;h2>9*x?U#9iHbe3T14oK&b z**z!v?t$m44rP$&PCs{lS?Huwa<)kfI{rMYv&cQX& zpEtE%E9Ym7yb5>5-+z7j_*mo_^6GcvJamZiRY3lBXjO&lr5cOJV&W0xaBzO^qOFUx}r3^&N5tm z!uH`)qqKz`!P8?Rh3x?Nz4!>4FD!5FwC6p0=h*Ug=K-=~y7&#?j(8>>ZH3Fc6*y1Yx+V%hyTR&3u| zuhU@vv@)7`+)%Eh#;WFh5kuRZ?yr^PukFnn{WgO4w^N_q*nL znvpV3v}xhS79MP)rTg(F^9R)`2CaPkKn$@{ZTe}$CD6m`AN*5aj&EN`HvQn@FukAa zwGl24Y7c_!b!F;Bh0bZFYa35x++!fEA>TT&#g8_K&0(`b^pb74IO=pL^ec6p-95K` z=1&?~&$5kK`s=$JPc!L-r-DpkG#F;`ZvS{l)T-*KA`MkFb|hvk7a(HY=nh(e`G2#@ zD|!e+rF}QmTN4oVHX$zJ_TG7wZ2YJ0H58?(Es zrKtlc@S~HN9VrwHVF$7Q`UwU`gBuGe@RI2;hE+63%7FLA)v@oQ_V>L-&%EvV;v@?t zGsW}7$JEu8qy(3wOTX_vo{ODro_|;>+Ed+sFod>_+0FAx@h!?(^24k<1}4sSoX)CF zitTT;P_PO}yrG`@StothwVeX7By;*zRbE-tRPb>8+F;)Mtd`W}{$NU4e~@4FowW8N zj2AH0Hj!lBpYuoE*QV)CJuU3}(1ET~G)Yjp zprpeW!}!{nL8ixVP_f0$848wQczDGudsUmbI|L$uWnKHNd zQ`!~nr;{@lM%R!UcTE|n+9khV2%uJuvvK4LE#_k9x@Q0GQu&y99h-9)UzZc!ki}`# zS1y(t)e1(unHs{{R0J33_zUc>+-|M7x|B2UkRucB7tkrJVjIO`eT*F@1MaJhC}da5q2m7eUQqGtK1`iaDK zj-t*yct$|-mT%BZ4Pvl)y$CT_c98}{`%ao>j$vk=dV;qdTp&>B#HWGx`1*E?HA~Jq zIf}gCPF5geaHlHMx1=hVsB8*8hqa?>pbHhMFPWoaHkRh&`1UlH2A*D1yyPTv*E4OR zjwW+1F<$!ONz*pe%FqmfW4x?Q%RTA`p<^Fcdl!T#83SfP}(QT`TbPcpykHKnJn^_d6U0 z{l!E2>pO@Yiu&~l@UQV}1zwrJ@6&P?@K~;bVTM9-_NLZu>`(})&`mdUOGkHhBp602 zByH-fU}o|cY6?p{u*rf=|!Fc1(}T-@!yb;|p(hmNiGWnH?;{JgUKAOpLe9_{+X z=tU_RnI~$ZOkGqVL6Klb{)Vv{1zn_Str8JaL0)qmqUK^dr+bK;S%^qq=*Ed^tZM`-~p3}K8gS1i3r50h~ z*N8t2FSRLEkc>Wvqra{>0a+fpl-&r5V-Hj<=P-+`?-4e;P#v768cJ;tncTz6c|Ve~ zkBTPpK@T@)&^IoIYRrpFP(olqhK$j+(I`(zg7{YQyG!rhpO7hF2nuZ{Tn8ee@!tM|{*=ad|y%J{-|WAiuMaD6+^w-5|{?W7${I(6e$wl7IY>gvS&1%x}*_ z{O4jbh`3)y(nG7p!ZqnsZ|%_LQtjJQf^8;FZh;u7>Gh1d3=2WXiGxS$h7aC29}Csh zs1%!blIrXav-bN{U?PH_2e=3x^qF_B2n~N=u4sNdqh~{qN^JgvufaI*wA+y6m~t(x z>nr9>3HavHm4-ThYsO;RcY1H^YvB!E#fAnFA3}BpW|g?(KQ*3F90_W%xl4+?`LP!6 zR-$emK))iPow+EN(;UX=xU<3avD+_!gWSPN)Ne@l4l>#HE(u0y7H{{a19ixEwe^o* z^DjHX!wS7!ZI`LrmM`z*a@lNFGgBDJeR|(m{nV`3=3EPj=K`vrh~H!;RflQca$}jt zU)SB3U8Rd{k&-2P2}Kh2q51l!WwDiHA=F?@$QoSowWJ4hMF_UO`C`Fwg42daiv2jI z#G{U-T~oI8#=Bxn2-(YEsYtf-GOXZ?WOXtgxsbe7WRv(gzbDe+PnFK?CN*@MQ2C89 zAcL4DdTnC^X634)s9%$1huO>%r9_KlcO7yg+e)zM@w3;lS&x|p&y_Z zJ=5d-rkxMo!w+FPR8Wd;8tDHM# z^5HY#q~#|a%{E-IZ)FyaiCCyV>b~rYlpfDoRTPUGOCsL;^m3o7O=rd4!Pa(DF<5HE z=W|ja^C}Jf<*vKbVT*hz0-=R_2^c%U({1m=8+@(TDxpi)#}0N?F&{@PtDczEDtVB8 zw_{@-9~>f6#k_CvF?fCBa!?AwJza+E;|I#8voAfnlsQ!>0|%#OLP!qA_FHp*?b)jX*_Q~WQr`?1w*DdBAq{fj{>wQ+=x#A1m3>P<1X^h&%na0e9J$l$ zkLZ=~@@a+KXVbs?jX-zhL((fvw|NgLoHGSiA;g{vkQR_@YbW0(%# zO*Ox@kE}$fTZf9E?(F#fVj<*eUfRqwE&TPJA@YYK$t@f{@_8{j!dlF-E<{zN#!P-T zeugBv+G(R0X-K(PK?%4ALzjckyqX+qK8}JNtm!r$ZnZpnn)r6MfrL&HwsGIT?0IDm zsi#eeW68bzob_}{@8vqlTgoEV9$r&(!GXk+^+(StB!`!_l>*4eHW|RR(tbbH&39NF z7g<#6oe5S?UJ@<$^FO@!>4!#Hu~Q4PDib%=Q9WT8hO}G7AVw=WJ?qKshUI!s{dpFC zD;HR|mag5P40sD8u&kAzH`VR){l~RdL8Kq5E0al8`19l`@1ndk)exUD6FkK9*U!EJRySK?Cu3M(@b=0jvWTP^kW2E|0yCs+31@%X(Z8^CT zQ0Mo~LKsVvg;P+!)t!rc-!~s>#~)wsBpjM&p09sVyWu;-=zXqev}x4qEkZ1K7xpg7 z2`s|?A|c}Cr-vYR*PTmm?4FOg?hlSU$rsCctxpf@dqGVcq4TU9B~0NbyIdHPR-{kUg@oR@XQr!A4ER^|3l_bWu-sn*C`CfVSq+oMIS$AGv zgFPWjgy5EJcaER@U8wu_==036U5O%XeJe(W2*Ky^e~onjYh@wh4aMWBg2f8)g51R2XUUZ=*w16T6^)qMvW z>ZMfPtZ;(qSCiYEonWmfj*=B(Z>G7)%emWndbWRtm8lS}WtdZ0zj)l7}yT)8g z9~~UgHmDUjnBIcx*0y@--lQG7I>vb0x>7NEWdOm)`$%1ldFX=MW~YEm^%-Zwvj7Dh zUGr2kAF-t?lM*Sci|u2}C*c({uNyn&UdOVZU5xphX5F~-8BaHw8)^x+fWWL!0v1Rl zTmTHVLI{9>F98?=28Dx>NHEykBGtqNze?!;!g-C|+ri$^P1xJPIn~Jx|4QmzSA2#` zFaa=98w?VLfrTOP)L>V9K`c_6EB^ien%5e7(tzwS-7@xy1s;2IAFLt=$G7^oH4kqL$((GaZs1cN~k5SFfB z2n?sFfI(0m2DM6yl#X z;6IfCjD%rn3-i`YX`^`=f!N!DWQPaA`0UE)A%K|2Xj1c}O_6JJEoMID9~~La>t% z4F1Z$-S3(FFf*rwV2n0JYL1+jR zjmrnePEUXjfxrq@3_wI-2Nys?!*R?7fP;Vf4$v41fnw!5K!YN2gaNP^c3`5RKqq4V zPyYQ51|T}nFbIyh0AzuE36P#hEKNWlAXc$50>~U7cm8R9Fcb)eKj&fo3x@!y5x-O7 zPYXf=s0r&~fbZ8z5XtxXv8nPV)nN*1zZgW;|4s?jG|!SE5Bv>)hTg+6*@H@ zz;BFb1khEm&IM>lSn9`LD&T$u!45wV8i6|n0Ak>n8X*14RZuuih#`Q!^k0!ut3vQE zVQ_qbfMe$m7|7CJEAYRE8;~kU+@^p){2L8{Q!o&J`tKj#0-6HNJJlr=FbvjRz(>H0 z^UwQ1C=dj@m;iie81}P(VQ_K}0)k_gS-{<(-$m?yxEGBCP!!h101b^57+??rh20GS zG=W57@d20MCJrzZVkHv>1Xvdr!K0Bd+!~DpN(+vFaGZVt)JdF)1ju<1wr2p#ph4E4A>Lf-@j<$HVP$+ia|HX$$9g4=MzyLqOez)vDjsrF{xRD6K zZHRza2&|$2v>p)lDi4TPU^2qCFEHQ%|5L+13I!kxP;zk?5r1t_{tQq69}=ff01GH= z;0N%baGEI~;QwYI*k9lo0aPfQ9EG5M_le(@0S*LeXevhxJ||!~C9hSd3z|7T{T&9W+4uHF{!O4hp;p_~(HDBMb+jp#Z>#dkGx776KIx zdxrrC3=F-0gf{qH-Jl5E#sEcPuUG+LP`Dus{WG5YwlFZX0-^aV7XPvajygb~V%heW zNI2*m2`mY5$^$UGqEqqW@M&mq;tQBzA*mX`M>r$euL|}@;otfK8xbr60y6*SZ(RP~ z^8i2!!$!RT4UU~_zi7Y1@qgVL=*c+y6eJKHI5hy6lW_$^U>7bxKqNMr1un*3AOGs! zuzxARz)JzF>;-mDpx^cWkLLg+{4WCH%yqyWXl$Gh2#7+$v9?rT05D?@-oI$LlRPjm zVY>u?#kge|4$QvTIRXd&=^K~+LYF|v1fXOpIN{%#DIAC-oHKy+C^ij=15kl72*(a! zKp!}Eg#u{E-+B8dBmg4t-U|B?;57_RqywNHCv$+!08W?yl@>=41nl>^;4iQK3JJuY zdHuJm0UA){afCtQ1Uwjy^S%UFfP=B84uFQj*&9QE_Kn022?zog`9ToCXp1Eb1O;?a zEE*IFtgNwV2sAdt2apl?&vNDW;u2V=0gF7Whk+#?PI-U>g&oHwC?sy$124v~&jYp~ z*s~-Cn1BIThIR2T`(tYipl+N40ZU_7N5CIA?^=MhJ_Jxaj|v!;n;5oYc+T)YaYF)zXR-crgnsp-6dn JWN*ol{y*(BDE0sV diff --git a/doc/user-manual/fpakc_UserManual.tex b/doc/user-manual/fpakc_UserManual.tex index 6477cff..479efb2 100644 --- a/doc/user-manual/fpakc_UserManual.tex +++ b/doc/user-manual/fpakc_UserManual.tex @@ -29,8 +29,7 @@ \makeglossaries %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - -\newacronym{fpakc}{fpakc}{Finite Element PArticle Code} +\newacronym{fpakc}{fpakc}{Finite element PArticle Code} \newacronym{mpi}{MPI}{Message Passing Interface} \newacronym{gpu}{GPU}{Graphics Processing Unit} \newacronym{cpu}{CPU}{Central Processing Unit} @@ -308,14 +307,13 @@ make \end{lstlisting} to compile the code. If everything is correct, an executable named \textit{fpakc} will be generated. - %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \section{Running the code} To run a case, simply execute: \begin{lstlisting} ./fpakc path/to/input-file.json \end{lstlisting} - in a command line from the root \acrshort{fpakc} folder. + in a command line from the root \acrshort{fpakc} folder. Substitute \lstinline|path/to/input-file.json| with the path to the input file of the case you want to run. The examples in the run directory are presented in Chapter \ref{ch:exampleRuns}. @@ -659,6 +657,10 @@ make The file needs to be located in the folder \textbf{output.folder}. If this value is not present, the mesh defined in \textbf{geometry.meshFile} is used for \acrshort{mcc}. The format of this mesh needs to be the same as the one defined in \textbf{geometry.meshType}. + \item \textbf{timeStep}: Real. + Units of $\unit{s}$. + Time step to calculate MCC. + If no time is provided, the minimum time step is used. \item \textbf{collisions}: Object. Array. Contains the different short range interactions (\acrshort{mcc}). diff --git a/src/fpakc.f90 b/src/fpakc.f90 index 308a0b5..8bd909b 100644 --- a/src/fpakc.f90 +++ b/src/fpakc.f90 @@ -70,7 +70,7 @@ PROGRAM fpakc tColl = omp_get_wtime() !$OMP END SINGLE - IF (ASSOCIATED(meshForMCC)) CALL meshForMCC%doCollisions() + IF (ASSOCIATED(meshForMCC)) CALL meshForMCC%doCollisions(t) !$OMP SINGLE tColl = omp_get_wtime() - tColl diff --git a/src/modules/mesh/inout/0D/moduleMeshOutput0D.f90 b/src/modules/mesh/inout/0D/moduleMeshOutput0D.f90 index db52de4..13e25ff 100644 --- a/src/modules/mesh/inout/0D/moduleMeshOutput0D.f90 +++ b/src/modules/mesh/inout/0D/moduleMeshOutput0D.f90 @@ -42,7 +42,7 @@ MODULE moduleMeshOutput0D USE moduleOutput IMPLICIT NONE - CLASS(meshGeneric), INTENT(in):: self + CLASS(meshGeneric), INTENT(inout):: self INTEGER, INTENT(in):: t CHARACTER(:), ALLOCATABLE:: fileName @@ -56,7 +56,7 @@ MODULE moduleMeshOutput0D END IF OPEN(20, file = path // folder // '/' // fileName, position = 'append', action = 'write') - WRITE(20, "(ES20.6E3, I20)") REAL(t)*tauMin*ti_ref, mesh%vols(1)%obj%nColl + WRITE(20, "(ES20.6E3, I20)") REAL(t)*tauMin*ti_ref, self%vols(1)%obj%nColl CLOSE(20) END SUBROUTINE printColl0D diff --git a/src/modules/mesh/inout/gmsh2/moduleMeshOutputGmsh2.f90 b/src/modules/mesh/inout/gmsh2/moduleMeshOutputGmsh2.f90 index 4159873..66fe345 100644 --- a/src/modules/mesh/inout/gmsh2/moduleMeshOutputGmsh2.f90 +++ b/src/modules/mesh/inout/gmsh2/moduleMeshOutputGmsh2.f90 @@ -95,7 +95,7 @@ MODULE moduleMeshOutputGmsh2 USE moduleOutput IMPLICIT NONE - CLASS(meshGeneric), INTENT(in):: self + CLASS(meshGeneric), INTENT(inout):: self INTEGER, INTENT(in):: t INTEGER:: numEdges INTEGER:: n diff --git a/src/modules/mesh/moduleMesh.f90 b/src/modules/mesh/moduleMesh.f90 index 86343f0..77170da 100644 --- a/src/modules/mesh/moduleMesh.f90 +++ b/src/modules/mesh/moduleMesh.f90 @@ -310,7 +310,7 @@ MODULE moduleMesh SUBROUTINE printColl_interface(self, t) IMPORT meshGeneric - CLASS(meshGeneric), INTENT(in):: self + CLASS(meshGeneric), INTENT(inout):: self INTEGER, INTENT(in):: t END SUBROUTINE printColl_interface @@ -637,7 +637,7 @@ MODULE moduleMesh END FUNCTION findCellBrute !Computes collisions in element - SUBROUTINE doCollisions(self) + SUBROUTINE doCollisions(self, t) USE moduleCollisions USE moduleSpecies USE moduleList @@ -646,6 +646,7 @@ MODULE moduleMesh IMPLICIT NONE CLASS(meshGeneric), INTENT(inout), TARGET:: self + INTEGER, INTENT(in):: t INTEGER:: e CLASS(meshVol), POINTER:: vol INTEGER:: nPart !Number of particles inside the cell @@ -657,49 +658,57 @@ MODULE moduleMesh REAL(8):: sigmaVrelMaxNew TYPE(pointerArray), ALLOCATABLE:: partTemp(:) - !$OMP DO SCHEDULE(DYNAMIC) - DO e=1, self%numVols - vol => self%vols(e)%obj - nPart = vol%listPart_in%amount - !Calculates number of collisions if there is more than one particle in the cell - IF (nPart > 1) THEN - !Probability of collision - pMax = vol%totalWeight*vol%sigmaVrelMax*tauMin/vol%volume + IF (MOD(t, everyColl) == 0) THEN + !Collisions need to be performed in this iteration + !$OMP DO SCHEDULE(DYNAMIC) + DO e=1, self%numVols + vol => self%vols(e)%obj + nPart = vol%listPart_in%amount - !Number of collisions in the cell - vol%nColl = NINT(REAL(nPart)*pMax*0.5D0) + !Resets the number of collisions + vol%nColl = 0 - IF (vol%nColl > 0) THEN - !Converts the list of particles to an array for easy access - partTemp = vol%listPart_in%convert2Array() + !Calculates number of collisions if there is more than one particle in the cell + IF (nPart > 1) THEN + !Probability of collision + pMax = vol%totalWeight*vol%sigmaVrelMax*tauColl/vol%volume - END IF + !Number of collisions in the cell + vol%nColl = NINT(REAL(nPart)*pMax*0.5D0) - DO n = 1, vol%nColl - !Select random numbers - rnd = random(1, nPart) - part_i => partTemp(rnd)%part - rnd = random(1, nPart) - part_j => partTemp(rnd)%part - ij = interactionIndex(part_i%species%n, part_j%species%n) - sigmaVrelMaxNew = 0.D0 - DO k = 1, interactionMatrix(ij)%amount - CALL interactionMatrix(ij)%collisions(k)%obj%collide(vol%sigmaVrelMax, sigmaVrelMaxNew, part_i, part_j) - - END DO - - !Update maximum cross section*v_rel per each collision - IF (sigmaVrelMaxNew > vol%sigmaVrelMax) THEN - vol%sigmaVrelMax = sigmaVrelMaxNew + IF (vol%nColl > 0) THEN + !Converts the list of particles to an array for easy access + partTemp = vol%listPart_in%convert2Array() END IF - END DO + DO n = 1, vol%nColl + !Select random numbers + rnd = random(1, nPart) + part_i => partTemp(rnd)%part + rnd = random(1, nPart) + part_j => partTemp(rnd)%part + ij = interactionIndex(part_i%species%n, part_j%species%n) + sigmaVrelMaxNew = 0.D0 + DO k = 1, interactionMatrix(ij)%amount + CALL interactionMatrix(ij)%collisions(k)%obj%collide(vol%sigmaVrelMax, sigmaVrelMaxNew, part_i, part_j) - END IF + END DO - END DO - !$OMP END DO + !Update maximum cross section*v_rel per each collision + IF (sigmaVrelMaxNew > vol%sigmaVrelMax) THEN + vol%sigmaVrelMax = sigmaVrelMaxNew + + END IF + + END DO + + END IF + + END DO + !$OMP END DO + + END IF END SUBROUTINE doCollisions diff --git a/src/modules/moduleCaseParam.f90 b/src/modules/moduleCaseParam.f90 index c8df0e4..8df3210 100644 --- a/src/modules/moduleCaseParam.f90 +++ b/src/modules/moduleCaseParam.f90 @@ -4,6 +4,7 @@ MODULE moduleCaseParam INTEGER:: tFinal, tInitial = 0 REAL(8), ALLOCATABLE:: tau(:) REAL(8):: tauMin + REAL(8):: tauColl END MODULE moduleCaseParam diff --git a/src/modules/moduleCollisions.f90 b/src/modules/moduleCollisions.f90 index 6557fa7..2f22204 100644 --- a/src/modules/moduleCollisions.f90 +++ b/src/modules/moduleCollisions.f90 @@ -2,6 +2,9 @@ MODULE moduleCollisions USE moduleSpecies USE moduleTable + !Integer for when collisions are computed + INTEGER:: everyColl + !Abstract type for collision between two particles TYPE, ABSTRACT:: collisionBinary REAL(8):: rMass !Reduced mass diff --git a/src/modules/moduleInput.f90 b/src/modules/moduleInput.f90 index 34e7507..19255dd 100644 --- a/src/modules/moduleInput.f90 +++ b/src/modules/moduleInput.f90 @@ -40,6 +40,11 @@ MODULE moduleInput CALL readSpecies(config) CALL checkStatus(config, "readSpecies") + !Reads case parameters + CALL verboseError('Reading Case parameters...') + CALL readCase(config) + CALL checkStatus(config, "readCase") + !Read interactions between species CALL verboseError('Reading interaction between species...') CALL readInteractions(config) @@ -55,10 +60,10 @@ MODULE moduleInput CALL readGeometry(config) CALL checkStatus(config, "readGeometry") - !Reads case parameters - CALL verboseError('Reading Case parameters...') - CALL readCase(config) - CALL checkStatus(config, "readCase") + !Read initial state for species + CALL verboseError('Reading Initial state...') + CALL readInitial(config) + CALL checkStatus(config, "readInitial") !Read injection of particles CALL verboseError('Reading injection of particles from boundaries...') @@ -233,11 +238,6 @@ MODULE moduleInput WRITE(tString, '(I1)') iterationDigits iterationFormat = "(I" // tString // "." // tString // ")" - !Read initial state for species - CALL verboseError('Reading Initial state...') - CALL readInitial(config) - CALL checkStatus(config, "readInitial") - END SUBROUTINE readCase !Reads the initial information for the species @@ -558,6 +558,8 @@ MODULE moduleInput USE moduleCollisions USE moduleErrors USE moduleMesh + USE moduleCaseParam + USE moduleRefParam USE OMP_LIB USE json_module IMPLICIT NONE @@ -595,6 +597,24 @@ MODULE moduleInput END IF + !Reads collision time step + CALL config%info('interactions.timeStep', found) + IF (found) THEN + CALL config%get('interactions.timeStep', tauColl, found) + tauColl = tauColl / ti_ref + + ELSE + tauColl = tauMin + + END IF + + IF (tauColl < tauMin) THEN + CALL criticalError('Collisional time step below minimum time step.', 'readInteractions') + + END IF + + everyColl = NINT(tauColl / tauMin) + !Inits the MCC matrix CALL initInteractionMatrix(interactionMatrix)