From c3f88ed36039337de1d0e79458f5eb6ac1b0e579 Mon Sep 17 00:00:00 2001 From: Daniel Schwen Date: Mon, 29 Aug 2022 17:32:05 -0600 Subject: [PATCH 1/2] Add weight functions (#21960) --- .../source/userobject/RadialAverage.md | 10 ++++- framework/include/userobject/RadialAverage.h | 12 ++++-- framework/src/userobject/RadialAverage.C | 13 +++--- .../userobject/ThreadedRadialAverageLoop.C | 20 ++++++++- .../ad_nonlocal_scalar_damage.i | 4 +- .../nonlocal_scalar_damage.i | 4 +- .../radial_average/gold/constant.e | Bin 0 -> 50216 bytes .../userobjects/radial_average/gold/cosine.e | Bin 0 -> 50216 bytes .../userobjects/radial_average/gold/linear.e | Bin 0 -> 50216 bytes .../gold/time_changing_test_out.e | Bin 49988 -> 0 bytes .../{time_changing_test.i => test.i} | 9 ++-- test/tests/userobjects/radial_average/tests | 39 ++++++++++++++---- 12 files changed, 83 insertions(+), 28 deletions(-) create mode 100644 test/tests/userobjects/radial_average/gold/constant.e create mode 100644 test/tests/userobjects/radial_average/gold/cosine.e create mode 100644 test/tests/userobjects/radial_average/gold/linear.e delete mode 100644 test/tests/userobjects/radial_average/gold/time_changing_test_out.e rename test/tests/userobjects/radial_average/{time_changing_test.i => test.i} (89%) diff --git a/framework/doc/content/source/userobject/RadialAverage.md b/framework/doc/content/source/userobject/RadialAverage.md index c8de8cc8a4a7..e4e031a6b4d6 100644 --- a/framework/doc/content/source/userobject/RadialAverage.md +++ b/framework/doc/content/source/userobject/RadialAverage.md @@ -6,9 +6,10 @@ Given a material property and a radius for averaging, the RadialAverage object computes the spatial average value of the property over the radius. ## Applications + This can be used for nonlocal damage models in the `TensorMechanics` module where the damage_index that is used for computing the damage stress is average -over a certain radius. This can help alleviate mesh sensitivity in certain +over a certain radius $r_0$. This can help alleviate mesh sensitivity in certain cases. This can be accomplished by running the RadialAverage object on a local damage material property. Then using the `NonlocalDamage` model in conjunction with the `ComputeDamageStress` the damage index used for updating the stress is @@ -20,7 +21,7 @@ The RadialAverage user object is derived from `ElementUserObject` and works in two stages. 1. In the element loop (in the `execute()` method) a list of all quadrature - points, their locations, indices, and selected variable value is compiled. + points, their locations, indices, and selected material property value is compiled. 2. In the `finalize()` method @@ -36,6 +37,11 @@ works in two stages. 4. the results from the search are used to spatially averaged +The [!param](/UserObjects/RadialAverage/weights) parameter determines the distance +based weight function to be used in the averaging process. `constant` assigns an equal weight +to all material points, `linear` weights each point by $r_0-r$ (a linear fall-off with distance), +and `cosine` weights each point by $1+\cos (\frac r{r_0}\pi}$. + !syntax parameters /UserObjects/RadialAverage !syntax inputs /UserObjects/RadialAverage diff --git a/framework/include/userobject/RadialAverage.h b/framework/include/userobject/RadialAverage.h index 1dd11f482fce..b0cd32df8c5e 100644 --- a/framework/include/userobject/RadialAverage.h +++ b/framework/include/userobject/RadialAverage.h @@ -72,10 +72,16 @@ class RadialAverage : public ElementUserObject protected: void updateCommunicationLists(); - /// material name to get averaged - std::string _averaged_material_name; + /// distance based weight function + enum class WeightsType + { + CONSTANT, + LINEAR, + COSINE + } _weights_type; + /// material to be averaged - const MaterialProperty & _averaged_material; + const MaterialProperty & _prop; /// cut-off radius const Real _radius; diff --git a/framework/src/userobject/RadialAverage.C b/framework/src/userobject/RadialAverage.C index 9fff6c736aad..a676612916e4 100644 --- a/framework/src/userobject/RadialAverage.C +++ b/framework/src/userobject/RadialAverage.C @@ -34,8 +34,11 @@ InputParameters RadialAverage::validParams() { InputParameters params = ElementUserObject::validParams(); - params.addClassDescription("Perform a radial equal weight average of a material property"); - params.addRequiredParam("material_name", "Name of the material to average."); + params.addClassDescription("Perform a radial average of a material property"); + params.addRequiredParam("prop_name", + "Name of the material property to average"); + MooseEnum weights_type("constant linear cosine", "linear"); + params.addRequiredParam("weights", weights_type, "Distance based weight function"); params.addRequiredParam("radius", "Cut-off radius for the averaging"); params.addRangeCheckedParam( "padding", @@ -57,8 +60,8 @@ RadialAverage::validParams() RadialAverage::RadialAverage(const InputParameters & parameters) : ElementUserObject(parameters), - _averaged_material_name(getParam("material_name")), - _averaged_material(getMaterialProperty(_averaged_material_name)), + _weights_type(getParam("weights").getEnum()), + _prop(getMaterialProperty("prop_name")), _radius(getParam("radius")), _padding(getParam("padding")), _update_communication_lists(false), @@ -89,7 +92,7 @@ RadialAverage::execute() // collect all QP data for (const auto qp : make_range(_qrule->n_points())) - _qp_data.emplace_back(_q_point[qp], id, qp, _JxW[qp] * _coord[qp], _averaged_material[qp]); + _qp_data.emplace_back(_q_point[qp], id, qp, _JxW[qp] * _coord[qp], _prop[qp]); // make sure the result map entry for the current element is sized correctly auto i = _average.find(id); diff --git a/framework/src/userobject/ThreadedRadialAverageLoop.C b/framework/src/userobject/ThreadedRadialAverageLoop.C index 6df37150f4f5..d6b76a619191 100644 --- a/framework/src/userobject/ThreadedRadialAverageLoop.C +++ b/framework/src/userobject/ThreadedRadialAverageLoop.C @@ -26,6 +26,7 @@ ThreadedRadialAverageLoop::operator()(const QPDataRange & qpdata_range) const auto radius = _radavg._radius; const auto & qp_data = _radavg._qp_data; const auto & kd_tree = _radavg._kd_tree; + const auto & weights_type = _radavg._weights_type; // tree search data structures std::vector> ret_matches; @@ -52,11 +53,26 @@ ThreadedRadialAverageLoop::operator()(const QPDataRange & qpdata_range) std::size_t n_result = kd_tree->radiusSearch(&(local_qp._q_point(0)), radius * radius, ret_matches, search_params); Real total_vol = 0.0; + Real weight = 1.0; for (std::size_t j = 0; j < n_result; ++j) { const auto & other_qp = qp_data[ret_matches[j].first]; - sum += other_qp._value * other_qp._volume; - total_vol += other_qp._volume; + switch (weights_type) + { + case RadialAverage::WeightsType::CONSTANT: + break; + + case RadialAverage::WeightsType::LINEAR: + weight = radius - std::sqrt(ret_matches[j].second); + break; + + case RadialAverage::WeightsType::COSINE: + weight = std::cos(std::sqrt(ret_matches[j].second) / radius * libMesh::pi) + 1.0; + break; + } + + sum += other_qp._value * other_qp._volume * weight; + total_vol += other_qp._volume * weight; } sum /= total_vol; } diff --git a/modules/tensor_mechanics/test/tests/scalar_material_damage/ad_nonlocal_scalar_damage.i b/modules/tensor_mechanics/test/tests/scalar_material_damage/ad_nonlocal_scalar_damage.i index ed2968c11989..f90b38195853 100644 --- a/modules/tensor_mechanics/test/tests/scalar_material_damage/ad_nonlocal_scalar_damage.i +++ b/modules/tensor_mechanics/test/tests/scalar_material_damage/ad_nonlocal_scalar_damage.i @@ -60,9 +60,9 @@ [UserObjects] [ele_avg] type = RadialAverage - material_name = local_damage_reg + prop_name = local_damage_reg + weights = constant execute_on = "INITIAL timestep_end" - block = 0 radius = 0.55 [] [] diff --git a/modules/tensor_mechanics/test/tests/scalar_material_damage/nonlocal_scalar_damage.i b/modules/tensor_mechanics/test/tests/scalar_material_damage/nonlocal_scalar_damage.i index c015692d12b5..4affe0574718 100644 --- a/modules/tensor_mechanics/test/tests/scalar_material_damage/nonlocal_scalar_damage.i +++ b/modules/tensor_mechanics/test/tests/scalar_material_damage/nonlocal_scalar_damage.i @@ -59,9 +59,9 @@ [UserObjects] [ele_avg] type = RadialAverage - material_name = local_damage + prop_name = local_damage + weights = constant execute_on = "INITIAL timestep_end" - block = 0 radius = 0.55 [] [] diff --git a/test/tests/userobjects/radial_average/gold/constant.e b/test/tests/userobjects/radial_average/gold/constant.e new file mode 100644 index 0000000000000000000000000000000000000000..3beb60434cba8598aa17c4894128633479ce9818 GIT binary patch literal 50216 zcmeHQYmi*UbsjB3F9hacV;f`pf_Mmt)f<6i4YP~1!UDZO5(Z;4x;t}cchKzIVII;- zRY^G}<`4K$N-E_zq+;VgHiex^uq&^v+Ko*mkeD*K^207yewm_-U1j_u3!?b+@=x9>NnyHB6)+qWM(ckLcZrBaJQD?nF+mX*A6w(L%Msnj9_ zQ_JhcDKA^AdliQCwY1!r%I1qx4DtEp443`9SJU@T!dXV1SMnGs4_cwa+40gb6^Y=N z%5&M8S2xd67R9d>^Z4C1{$&QESSHLyg+uK38JujztKv7+4dSdMU*pusWExv}o4W8e zcj0X@yi4SF$_g%Z3HU8M!i9N+3-buq&O7c_$zz@<`Iv86%_L9dWwBiF5#>t^vy4i` zGRi1*AD;m&1zqA+irM2{wN~`Y2&OIuEgD_)gD=Bf8o~ciq2#;uVtFE4@r&hp_Jm)} zBhO9|n|Y{yAMWx(v4rx_H(KuT{#8@%ba9G00KLgVxQY5C*unD-{3`*zbqf^9tv6WI zf(>u#;(D=O;t%4YK3$UY%ZOO6uhUHC^Yx$=N*7N|`X(}en!v}YGM`dNAY48qIHwT|DTOmZ>q z3@*@q8GJy#Ryy-mzJPs+{;*A%y(ul^@8OuWhAn#$nGwdN)NcDhSmqhk#dJ$yjQ$J zEaYh`mgIk*Sf)5rS1K&zS1Xp}ttHmg#9G0;3(#~rnGRlk;z{&th_{mVHcNvv z&oQpSrSWB4+6uW^asLZU!HtW{x-7Uf;KXG<$<=~;^9^7MZd}~GjBDg&wi5SJkX)^} zbua}tF7DHetMwts-^&CPE^%+w;hSH8D{&tFLW|EAn0u z5%-heuVlqyr@^#y)X6JNeCD+(_x|erYUb6Y$K4WYtHz^}<02+;TI)rNo%&V<_byd( zs6sKVuzv0;+~Z{_%Po}U3MGqRZuG`O^04m7MR^<@cbg%h^uk8WW+UcWW+Uc zWJDgr59AKZWqwa&==k06?_9uT0eo&b&<;=tMNy4oLZvN5ANS=!EEl=!WQp=!oct=nBVLsVCGK%0xLtXXvLK)E)BZ zrySHH^5~}=)G6}lrySHT^5~}=)HU+xrySHf^5~}=)Isv-ryLH5Jo+iybsz^sJ*FJj zgH{~Q= z?N$(Zl!>xQ-JqYcv5t^OKV{nhB9DH`!MZ~p{gh)9h&=i!2kR7h^iz&4AoA#^9IR{P z(N8(Hg2-|<89!wkGd%h!+it_7pR(O;c=S`YPZ}Qmlx>gU z(N8(kqn~o@1(8QTW!ncLk1|oV{h$NJPuUI{9{rT>*9LEihe#&ve@aU%;(}qVsWjkqji{K7{E&?qE zT?|?RS_)bQS`JzPx&*WmbSdaE(EC7_gRTHw3AzgOe$WR%SA#wX`Vi>Dpld)M0bL9F zDClFLkAoc0b)f4(t3Wq^ZUo%~x*4<@bPH$=Xf5bg&@gBnXgz2HXd`G7XftRF=r+(+ z(Cwf*K-)l4$9KS=2HgqDfJQ*0pu0djLAyX>pxvOmL7xQe0o?<-7ql0&540b|dLgZc;}@4@KWZuUy<42T^CJC(LW+BlC#qXWA$?W#aDflEl>RYZzeGIZP=~3%Oh5IOI!k?}t};JF zSFeMgx<#G23&b)O8zXha2d@Hp3{(ZxKy^?9#JYC^#5%@$MLC%k>i;Sb^_69Kj)S(ai#86i*5*6b-7Z%yz}`RPo192{BqT`uRM}_ zJ@aPeyDuFb+L3wlw724g##~AIf`G5J;mJ|Oy`^JhppZe&tPt1ICPTk1p@E0;~@_rWqSEbYs1~`Hz?we2k>*2q+ z;m5OYEM9Z(J>UP*%=eBOToneOv(ulNa?5#VFQ#rfqt%JVlvl3T9Na=oBj;J;z8lD` zabGsLPtJVva|+j4Kj}|->ub5m6Z!QAtNujQovN*$^8K1OTy^uPSHmYfeh;XvuY0w+ zzLTYQ*ss)`;i5zTx}tNSQLn69w{E6^m}-*exLfnmG6mB)QuC?@#vjGhvD$hUlSthO zZ@mF-Mj9UWs#U-076k|EspQ=YLk+Zh z9kiObjDG*x==IwnInQL;^5Njv;X}DjILFz%ZsXQ6H*dcaUdCBu( zSKq~dN*D&jyu7B6X=l&=2S)dr$>70rksil!OC|qAw(L#7avk@wj~P2EkZFfK`~EM- z$u(;D;EmM~e(qRJV47|?596Q%Kujcrd2RxU&mNGoAYmS|3Q=I1ZrT|=a+v+YJ}<|K z01`-hKsDJo8r2?uk2a?7^{QpBRC|xg>NpRt zTOYi63_wVhYfl4aUsC{B(;cl**hMs+?U%iAaGai+n0D@0-NaGdsY4m|Rx*m?>bTh> z2ik>Gp|srYjT7AEL&|b>!7_lh1{Q#RqnvlEC%aV83S=Xz63G=%rIxQ@$O;DGDjOyW`5MEhzmgw!gl1S*O9+;Qy&_TkdotyXU-*b4+ z=w9dWo_%Of#|~!i8oPVX{+qjt$!ECeC_CZR@xGs}H425|G`Z$VQd4T4&VRvZbJNa- z?xI9+`EaA=VMSxDQo`FwK8s$%WD8^xK=c%fk;RgB0$bGAC~ZUlw)s~fP-(OOs5!TERXJSyNyU< zUFvCdBsr#C-p4^4QU1Fy8q-r!iq!e0 z((>vPdhzKA#+%&%TV8#VB-yQ2-6*|Aq1apOf%hCR-c#|*YDA7X(0(OpM&lcLDPs<- z$6OQ#!+Bn~TNsF$Bo7DPb0{1SXvIo7^AKdg_+31aleW|mU6b{moTzaaEomkZm-yGF6UJgA~tuxXxcQo#d?Q`B6{&^ z)u?+}Y3}hdTFk?kWyigOU-ejtuwzPgJXh!t05*O)PTm$EC5TT8i4h1D%-;L zH6>q%VQs@`3$ReuZnjYy02U*AhpLhChBt1OLX34|#+8z^dR{r7<$Fn6LYslK)4y~S z9}@DTGh4mxX}2$w)^@H-g;q7isxG$ymL~jau|7G~nq;k@v87=mV1S^dC5FwN%8q}2TmOnt$d3c)*Am*@`20lP(@uVH|<+e{qNZ16- zRxnWn^GsUsVEHT93oh3u9E}5WaJqnv)EbpMrUtcw(Vup zzal^GrJWO%zQsdC$0@5)#!xWRjvZerc=5pygbXyHP%O$=17fyBwcea{BKbR?H?p=v$5(#K)zI-GP(*h^_ zGIz>C!>q<#q{&O5hdsyUs5Z7=Z3jaXkJP?-XRl9*rl(jveMps84060I5yQz*OcHLuX`uQ-Z;S2xYKHg z3j=AnI_9;uUG@?xK#YI{JZ8NdQ8e~71yivx*R5Tw zmQdO$PE|_9T(LgvHA)^ARP?;?-SJqllJ%xJ5X%Ni4}mlh?Mt~<_G(iF%wgeIkZyTs z+eqor=Ey#Af;>bH)wvyUGnuWv1b{W2Z{2#EK(;py*b}zbKudx`QEDAIE|LV$bMq(F zPQ;kh;?*%!2lWZEae;uk7PgAJKtOce@j;<=;Mu-VC=Pf3E@J;rj538g5B6%)3aC=t z4!Iz}b|_sg$pEE;9m+pKsywtB?!(Fq?q6-{OQ0^|u`gu>htJQon zvDm)VYX$7E*G#f)Gpa`h^eBN_S7VbZOUJzl4yD>BnRayw^VN$k#>P9FXL~@*w;ac0 z5GjPPYYTV|%=_An*IM-gFp#gyJzLmhH`Tl^zI+c#Ko)AxPuBK|9rMZBz6TpSm-wJt zt$BGVyttfYl4GXF2vi?)WAJMP9QbUOn0c^!_v{}Vy+1Cga-R}%Uue2rFZ(7y0d{_V z6LiefEf!W!J9nlx+`4+>+WIZG4lC!D9y4`eRtG8J`8PdgpIOxW-kh_lY6_M@v@gJ3 zm?5LOo%81FZV&UBia9YlQi@V}LoSm%@cDJj9FF-?G!LjO%b-piF-eA??+WWZlb(u@wgO_4DG3e<5WDHJ^l ztIl=ts`aARX1LoQLf)mh>#Foog?5HN{9R7t)(JXkBUfZ(sF#Hlj()SIGTvl223fxtu-x^~81C)*t_kMxID#i>DFb}!61$Rf0M|!(hlUXiu^iYx()H&m^ zmr!rNfZa+OQ{(9Ehd(md+e5chnRK&NABP&YjDM@O|c%xX# zb4s@{4*m%b^Z0`Wv-X9O!y<>OU&QJH(P(J1o?y&8wa>$*`3mfnHeZ1?18Dy(U(8!d zk~~z$C`QjDOFIuf#O>1hz1S6LnNup-*K}=?%O8yyN=sr7VVtl`0n698b|m~kJxLx; zxY$c;kX28*GSI_17A^3k(`F7Pkx-w5nuQhXI(6J$u{mR+ z)wN=mR~`Yh!h_nN3d1sRnL*cCl4)nx*zVCIdu0orI0aXLx|PF}fi9Q*EYt3V8$e@H zs8yONA+^PPNW#_-82+1RS?7wC{vN*J#HBR2?>N;{S#2M1rf!#_NjuzKv;WI9^tyYH zJu!)sgZWy3oeqpz+8YPVmT21NAnO8;tz!3|mOf*u4>T|}1AE|YDd^afKk=X00t`2kB86HCb2ah53*$-d?m@lT4Owq^}Jkx)G;9_!E1%% zOrrx97Ay!S(ZlIp3-;2^aQ`nUrY@fxw1g)IE#av_OL(fcgjN+Lbw2L1#O}brdP@CN zcH20AfbQId&mXph@CT27p_CV%#Sb6QyX6e5I%AWKN`@qr_a1$iC z)@$FU1f`t!3sY{y{uPkXUc#PJ@d((tDdpJyT@e%A)Ci}V#?V?0ZGQXiT(EV18%_{N z^UY1op6O3}@YQT;z=bkVZd6KsdwmP!OHwG4+}%I=vVW86Imki@Z6|p1-L8s>9T7s> zm-3^vZgZakczeYj_!zADq+Tqj1wI%{#Qr(O=5k8=Lcy|%0^VueO6!)WL5dU9u(>}@5 zd3Vr&T|XVkk0g+ADylhAn6`lDz>=JURUeykO%3q1V(<-pudVZ0Q>@0)QW!7|?G+&U zVqA1M_2EG#Zz6k|K@$))A|3$l9y3^dwJqE*XS|Ay?$wutkaP^%8wUZzC(Uk6y*;3> zj*U_RHk+oVsP+c7igkHdpocBL^np}rUdKWJ2R}VVJ|J15yj1n%AN}TOa$0OGk$VTc2qP$v~L2;HgPOhwO@`Dwn4O)@R&Ly2mOC@ib`lh2YWy$ zi#kCw0@6Oo0bt(axSv*Zwg*)Enw1-zE{ILtc*X{%%h!U*lJ-fC#xfGxe(?$B zERAD}`XbsEyaUx-tpC=p81@q4BmeQD->7A=PnT3a*oDU$*fbUM((M6F9f{ug#MU_c zVWt?ld3_!P=J|=l7&c~a-ag4{=M(1NNi_Ge`KzgS$pE=6b`v1%vzJiqMNq-g{A^c^ zv_^bn0N5AGIq;rxqWudrNqyupXJ`#XLh4D+ zcPaZ$`~Usc0d2>S#bzawtS9p4F`SLv&afZ4{yMy!wK(Lf_xBX+OssYzhDoe<+40`q9tMEuQg_)we?fr79J93j zN&-0W#k>_ZqTzKZB4-i_)l%6(p1h_2qXtKkti9z~)s0ZwBq?ScUT3o@I3Ikkuf27g z(6O&UYL&QBVWr|f8yunY~dl}E;6`H97(6Rz0Z*U zqB-OI)ehX?d3CmV-)G2Y^jY;^Hn^`C+?0`1@Z8v-zT??B1)MSAjK2R{xx@P`KyBU3dE~D?ip3%?m z$mr+a8qv{$`j^piyldRE#(l@QXXZ@&*?E06!{+_&ZPlNDV&-{$C!@n($js~h{)m2l z$heCP?);qQci4;BdG)^Z!|9p7pHsjgd4CjxD|!{=x5 zHQ^wA?cpHb&XxR=bj;~=TxQ%Q#$6h?rW{Rp@!ZlsfBb{X7JYke>39C~#j76tr zOFr|F-&{U7bj^{ccmKz)cKlG+I~mU$(f6+yG51IG^Usdx=RZHvo{yK9{F8dO*o1#? z+%rZFkyG>>ay_&1PbMCEeQxOAzx&43nb&8&ujLrg;nzm={a=~iTN$|R=?mjtJh${) zugq8e>aXX9e)#c!oL+Xz%s+fi0gveLnh|~fnE+3xXStDA?9kFV{X6d(`4LvQUk`A< z(}ug+;C|iUt_WO>Be+Wp?(YrmGK2dKgS$q?2aUXf`!@lemiJ;K@9zxGMF#h*3FCL# z<=ydblD{fn)*5*Q_hN%D^+M!bVep2GyWHS@)5v@44%4z3d2cniPaF9}FP0j7DX&Wm zp6J60qYu9^`tZwUK4|?g`tS>b`>x;_cTosWJ*aqftT7w5*|$o)7R zxi9hKT!|m&#^K2QI2^ez@#9>HALqv5$o)7Rxi9hKT!|m&#^K2QI2^ez@#9>HALqv5 K$o)7Rx&OalmkI3v literal 0 HcmV?d00001 diff --git a/test/tests/userobjects/radial_average/gold/cosine.e b/test/tests/userobjects/radial_average/gold/cosine.e new file mode 100644 index 0000000000000000000000000000000000000000..9055172a2fe86272438cbb4964f4898cbf218795 GIT binary patch literal 50216 zcmeHQX^b4lbskc(Na~g)+ma>Qt(!WQw@6Xm(e8>|O0;-MB;`Z4M!hpVyDiR4kLFk| z1t=DfLk48T3F0J{92m4^2S#8!Kw=vRu*wp$fmps`8wLXZ2x1@sjQof$$ks)2^Lgeh^wr$_BJd??s4LTR}UeL;tSIJe}vX{v$!)0bw z!!LWedc&*oO1l-6W;s{z%e>lF^MFETi}s#n8rsvE>PpL~r|CzEN6@-}qgZS2C^ zWO(Pv?^G0A<~;BtJi^6!gp2bC*Up=DYvfTUNry z_wgCf3eb6O)z8g(wYncvaG7}*XxYTF=RXI3MGyYv#Zurl{K|B$8u*n)?r=~mKxdB; zn>tj#4}VqBFCjnljmSOTzqIVm`DNMw>?Q}{rW-S0htE3*pAYcOn_x(8qsgomYy>m! zYWR&3e-KXDQ$DC8b!!Yu>FOHXJ6GA_;hDfh=5Le4orAzieerb)u*#ci{UELxmpTzg z8Ax!e-oU?=h+nM*RoaTgrEVEMTbBtX-&s4igQK67=a1uhRi)nXdz48W<6f*UXo|?! zdfV?YZqo0`^YmBMJAaRMATIiQ4wx$uKbk(q9rX9sA?`W(0M;>wkk-jrDy08G_OpJV zgXfOJ*Kt`B^o@HV^NV_?dTAHTf5{_}?cCfBzf^j_Ej6`8FplIq!V)RYFq_SJZHI)#}4U|YRetVUwa~emgrqD6YIrkl$*sf|TW$WjRg9rD)Z{y$j zDt`O#i}DzWzyD_pWHQ&ihJPzn{JN?*DXS^~ybH(YR;y6vhChP8GE!&vi~9}6HMp#N z8JBgkxG3%~!4%x2xGc+p%L<&h)RVXf?v2-hDY!{-cQLM!mue;M1t4)z+y$EAHKvDJzS*kssGWqFGu1^+Ai)MclA|~xFYW*J>uR2{`o9e z>@=8mj`sF$6Q8hdP?=&{Vg1~t zxJPFx%T1K!5+#dZuJ@*5^04fQqdbnj_{Qzwqp~!=v%1`DK+|B$w zw5`P;SM&0yo$I*fGnawh0}=<`pBt515dT*up7<}Q>n|)ApO{fMX8576pZmg?_!@W9nD|DHG4YKYVbo_4och2B41D`pgU#j60)`Qrw_)<>9kMsY(uch1$+y9l$$Aw}0re2YHZdkt5 zQ^WSdde6u`v}>)`R{cG!YpIuq?Hl_K+mm5koAM*|=CJ%>+rg=ZHx)D*ff|-N)g&K| zyfuG%h8@Dw)ir|Z>6&zv?vhtLowGl4y4(Nhp8uzN`9Iy;ztg(%Ka+TBfMZlLcEvF= zj-Bxv0LSiS%$sA+94F@39>-rfR?0EZi$NT-4E8K^*h^ z0O*6D4}l!e)u3xY*MhDCT@ShebR&pk+&6)42HgVUnCu8>4QMTB9f(7P9GBe)+61}{ zv>CJo#PPwcpgTYx291KUpgTchpmERy=q}JU(00%yXb0$S5XW72f<6Me7qkns8?*n8yrDK9k3=4%#4j3{wsVL>|MG?P`z%q8(F?Ye3hUFy*+; z@EE2X*Bc(gl;Z}&W0-Q#rpaTNa;yfC$1vrft&_(v<+vF{9>Wgk77%$1Q?^?{bJv9f&-JDF@3Ac??sI4IuIurW`C&2a(4x z<(L4Cn=oa&3q;;l5M|p2+HS&>ZPM@bI^XLt-#wueFFF-+N9(8DH7*``3`F-+O=Ao3^^Wh;O@6Q*oM!(*7TO&cD=lx@cF z7^ZB#;W11(9x*(I9|j!)k;gFQD1pdhm~xasv=zP!xpbJ44fi4DJ0(v*-QqX0f_ki9DdLQU=&=sKfgRTT!1^NK!gP;$A z9MILEYe3h6t^-{Ux&d?}Xf^02(9NJ*K(~TMKx;s2LF+*4K^s6DL7PCgfi{D-fNlqE z1xXnnMK}w(6Ep@I2Tg$P0&N3r2Tg)@fbIs}1KJ7t2Glo1vXh*VH3*#Pl<5l$$c~_lU)K%Lk-oh{wERzEJ1XCv`D($CHdb=cb*|l9<=EYa5V^K%X0n7 z{i}|D@ka|U{`R3yzy0#Fqfb7eZj700{2uOi5%AiK`oR!K@Wg-VmB;USXvg<~`<_iR z*FQ7*%(&q48i3BWpj>t<1!or~H#!rw>1NrhH0lm+A*PY@r14LLa%n7HYHr@6A|B9lRGpoex1y52xOjX_=hIuQbSrg7 z^KA!`=wX8DA@yJl^#>-p-jPxI{vz^?`hUUL$+FO1I`Y4Tt1ZPybBEJ%(|*IrBe&pf zS--)(eRF=}ZR?$&;%r#AVZ+F}+ebER%?Pgck1M%%`BS^S`ivy71(|tj{vTAfWsH44 zW$gOwn4B>^?3+BWKi>)GI2+fj-@ImX=ID1uX9QQ3pYNJ-Q;Xs1a9mHM=}nO^t{lXG6a ziOEkTPhn!laJj00+X2^axMSny2w+bZIu&Z(M=^ToMzK>=r-T@wc;_p`%FC(p-yeHZ^JVH^)V>ot=9gnAl||gNM)c3^l0f*Q}R&)YPN`nRVE+AKW=kzFEfyZ>)w0 z@`vgI)AHjwOo9#pF_92ebPGsA_JEuPiRvMz5Cx{?XPt?I2iQLx&^b;IAc3?8RFjQ+ zquK*<619la#m_o-PaNF8f4~IeO@*&d?`ZenM04(5uU7F&^>?VOj`QG}wc(q`5QJpD z_N~C|O9}w%x~nyTdJ&ChyR)wx9H*}(W}W+0H*rFD>yU^2jf~>BI&SXZ-ge>CRaS2I z#|iK9AyspA!7_lh6)XUQW~JcPj&#Meudz`h&UpDlIf2D2>w2T+=lQAumnKuN5Tj8G zN|@h0-C4t;wh7BpsHi>x3gJ~$XNkd1mqtQQ^}xJbhYb=w>)g0=&&~roCw4gpcJ4-d zI=L@**W}$h_i*&8i{mp~bd;O+8hGE&)tkklKS!={X-Z1n)A|>THb3jE>mo;Q79Vcb zJ*;S~S4((9Ddf;=n2A890Yp#1@0l%GC#;JG8Un}GR3~e9(vCK<2}`&qDJJXfp(Xs8BK2JC6766^NU zvE!wp*$A*M9m^L@i_%`#+84dPZ63|0ePk&jOxIw-ix~lJ$8GvH00cZxb1QWpl~GrT z;!NmpDyU(QrIM=$C5$W;UBAu3?G4+_=e?>z#O4kdO`AivSnm+wqZhA6-XZ3|)ixw- z1FZUa7aL8a!cwU>X5$uS-AcaOaEpEAKs)3&g?z4SgKrP0Ume|M^lc|n?^vkMB$-WA zv*G2WxyQ?>UqEG+oAQc5&0`_Ljw!jRe0`w2!ON`yGzo?QwdvB{u-WEeSXb>%g56u( zT7D*{Mm=FZCAU0PaB~>kr+ECxpXexICt#@&h;^sPlzy+@X(|UZ2J24Mb1k`A*$f?e zBexA)P%_~dY*i0iw~$OIMUwR>am~Y^L#VkHu99%w?TMKo5ge!Bqqc73s#PrB$EF-D zC>^plY`l8&Wx1`M35~lF)XF(kmIktFh)D@Gps8n;fTdCmz(zckZP)fCrO-gd+J@0~ z!A!Yjqm9`B(C^thRE?B3vVNmvVxpfgu9T+L^D2cL-%HvO+PcU(gG)#8At66Jv(@XK zcKb|;wsTzyv??i9Ww{NoG#%9Z#!NYyWNn~{xnUw;fFRNmqt4HB4~<$m>>%40y(l1? zdcVE5y)+%dB8Lmur4Cb>VpVU1vIu6HMDS4a zSFsmden2=@4%FaW5gVyBs|8FAih?nmbvA9-x@qg?EgQBDlx_ePi>SYPe%#AChpPjN zhlq|-QMrtvV5S{AzEttzgCPhRXhNn~8T_zb~rc~%?oCuYlmX9 zf_{FdQ)XjfI3ZwH@qnLoEV5W6WJNvHT<%iyhiSKO0@!F`Q1>y^VL_?MVQl1jv8Wo; z+|&_tD3I*I-|C@5?8G)*^9$W_MCovFFq6yNyAs8uC!lTh6(D(-7C0SLxKkDyW;O02 zO^n9`v$6ept7`ZhW+x&$x+HsvfJic?4~90QGG$2;cE0?ikWi0<^y)q+!?q}W z_|&Qa%}Z&=J7j7JyYw}cD($N;BVn%&j!m)CB_({v>)y$+uN+X-xN~ZV3j=BS2IjW5 zIeQ5eAVxsK9%maKN{a1(q#oz#`A%w|De46|H(Sp=;)cUt=!Rj!ZOU8@-?%a1nV`Xx z?^e*{V-~=SS4&ab0FZqfJ+Y5dl$WYQ8S|R11?Vt!WQ~1E!BlL_b!!(Za+P)ba<$~= z{lg6Klu<$ELHyzqGQhKzxXP-EQ z4xyn2w%*4qTKedT~XVfziVq{tLy);-5X(g1pH;fUIa7?WDO28QZj zKA{?y2&hY8RNN&3qU(+i3ekaQ`%IxY-2L0f{+}3SigzCD%_a({Qrr%?B*1nkoiFJC zrGp*QA0d?vt%m!sI)nRHoAMHxiv;X*SplNoP*$Ig>J7{S%X{^DU?vvZw|Y^)j{2HO zwrxiB$bcRta2slDQmJ&xo90leeUe$PmN8%5cQH2J**x0=V!q`RCWA;O#C==9Q()fL z4!qW?7l5IBUGCe$rn{-;eevacSOPLr`+l;vH|&H@*7iNv#JR-#+*;i$NaiJJmPt;S z9wSf#%#FdX5pd|USrX>K?%KI$a^n7^q{@9t$bF&dcXjqnfCB9N{3hs_sp}V4&pCHy z*WJ2${Vk1~ZXHqHO?_tSqFNoMgs0#1nB8Vk^E-3Sswyd13emm*dubIJ&Fz%0zV7xg zpUKZlr6aj0g*WCi$wQxiC`$|o-r>?z4rX#k2bN-_X#!AVqc65IyVc*}G>l zrl>xbd@!L2kX#tN2`G5^l8aA!9g7V3YiXL%p|>eA20(#2T$)VLqp<2+7q8avy|(H* zjoI`a*3CNJTu@La`1Hx!)-5J+=&L^zpgNePLxap^yB!un3Ds1VHtlKjti#S;D@jQq zKJ`vtpK=7M3{^u$NG8zkuCJn&m0< z_TwKJ?CV3fRGo2iwE%}2Mn~91QhPg`@)~{cQW$pMx0Ey`Wvc0y3Y^kyDu>{(hk5+r zf?4}a$zzd2E%33rKr9;8tS1;#r}jE*RbOF!rK_*7HUsG3EnmW0N*W!iV-$POB+EJ< z`#85t8`QBY5}8vf*4OfFlFJ`W7)nc1AEI)?G6gJOw+v;7r{vS(A0RyXN4|H1xWAkUKnslY{wMfSnGE zS=v_)m@Uz&pF`D!0b9c!JS_vpR3B(i)eP-*k4fwMa+YRfsc)ClK2s9OLzIL;h8U~$31>oj&a z*KfgI))^W6CB>BGBg2;P$gm}RY}gV$)?Y%af|5BM_gP|hU}!z1K?=KVoZmxt?#$;8 z+f4Yq$3K(HOIPuOd-ZNPLo3eMWTRT;s_kwIDB1Mcx62VgE{nr{8(|Zw0d0rw)vS26 zGf5keT-ZyP4@#W_R05|uOIA6hy6l)I-!_`mf=da=9VEMq4k+9dNv`+Xw<$#~7lLBh z?XiCaWVDyC?^HYjc5X^Jxo3NiiEe6yQ%z%t=0jW9zB?Cgo!^ENUS#p+R>v2bL2U2;+0L0fuy(tow`6Ej(Bs`4|;yANJ*;9%n zM*(H6U8_Mz3&`(ySTE&PC|&er+}^Ntl3SY_bfnr&M$w_lno`iYKFnT1jySaJxFcj= zQo^O=GngY0bnHQ$Lc$}Jyk3>kAsrd^^&!@R4YvuHZ5-+4CkoRR@Dx~*v#%ClbFT6bUn_>+(D&OqpC!d=EG>fp)7V}CVlT!=hf^OOWOUQB zr)e6BgwHdSaL1&<8mMjIhB;F;Y;>=_EX1VarG4cffcT`@t*N&M4Aj9a#@(z`h@l(wNQgXS3)te|pUdwXI+iJ|Ly39KwC)!i^2k&w0$N1?qWQ zLqhdM3MWEv8+n*b`z(HxW!q7a7NB%W+DhP30B5nDT-NZ0WzQZ^Zzy7RbqI07fPF3( zXSwldM?n)r`zF9_6UTg2`{hVx8$^2v52-V9(Em55sKhpOum^;(s1r1MK-woc1T1(Q z_wy=wsjpO`?E%%kW|bzV3u03@p0R=H^7U}Cqpy^pp9??81I z>%a9YhP{ON$Uo}`&3X>|bV=buz3^xgo2Fu3x;>yNBfWP%v8^0_KU0j{f<6xd^ZZ0& z3>&jIZ=Ynf^9l8L1kJr{{%XoyIzVoV-2#Y*>?Krt5md1>KiAbF?EyI}D|}my%a^8+ zGFu0P`w_G?0QQ-33cRNrY#mlQw0D((I^%;zQXjd@8CpY;kap7dUCO@G{(rx9K-)27 ziK=9h^+f&xhO@ET8TLchUkA7t0)J+6^nq&~0``3muPSc85YGQ&pAeqpIwQrZ% z<>WD|7KeQG|DJ-KiPdhzs1h4pcD%Q*4?{q7sk`i;za&6EPnEX2G=M{2%v)t68ef+p za;A|`EtMVS$!iL**WgH+wYNO0x*2PmB*iSlYiu?J=Y#L`wYQEFJN7k9tqRzGC_Hs$ zvzW8J<&1BAW;TxFHvHE0%eq!v&}YC03+fE``Nai&M*O!H^cnG=8U3b?H~N>yRA6i= zXUH!uo*Xl0kQ@KFoGonRi=IV3z!~yi{INM(SkA6rW^h|Ll1}mapCK=2od2K$H+){5 zZ9ecB@;LLJXVV`wXQzw286&6YL(Za?{E_@azO|ouKmUhMn0!#UUop5p3C~zIXQdl_ z!AZh>di2lRaKAbruGm#e-&)=uF6eYWJ^HK@l~?p@b%s0vpalX{cwlGMj%+B@BDCAFhsS7+Kg(~mRd z4|-Wzdv{D%&zU5D+{8=rlW^pI5{}%L_({IRPx6y+EX7f_nWx|5f|PjlUSb7W{=pRW5^tMfP*? zes}WK>ev-WzqN4kmE%`_-ucYva~b7EzowrU)AzqUrtg1tOh5nVn0o%X<>M-LCL{is zmg9u+Pa6L%;~zg};x8;}*Sc2y^&c&qRQ>s9NB{Z>6?4p7UKk-1y6kpD}WtH2w+WL%$05dj|LE z(SQ7PhR)y~g}nS-kvB88sGR8I#Dcz-{F3|=+)&>6zjz{w`xS%xlX1nx?-?B93r-U5 z)1zN*!~N=jxMEkKzKtG!*wTZeB#Y_yg>l8rgnBjdMsRdOj_G*NyLrFU`Tw=~zkNd8 z0j_?()AxTmuH$#$E@-%JdeHB8+Mb`9e_jI^`+f>~;Qg}=u7!S0d#B6GcSh#~SCy~t znsQT%;YvAq>)1*i|JVwnhq!O+q&%Xel&o3B#{m%Huja+Ame1`wFxjx&-A@Yh`$B(Ue?ZlIh zH-2$!#qH<($MM%@M?a^_|G2*XgK>TT-$PuT&gF7#a0Tyo=K4*f^8>%z@tyOIt=RIf zzu9xuy7?ysS6y#LI(YxmP)@C{69!kx^BV^5_15pG=gv2{Yv;fFumT>}zq8TcUTkpR zF#cP{m+~$?(y!_{=tbdPVsJN%{zAW)8{BmUca^~ry}up975xe>e4S2%yLt47UFAK` z;Jj}9H)FWar;_tE;|uO(26yAUhRqmx*BE&PSLA)e;EGSw@ZdHUn|Z$2Toy#5p5s&Z%WB>hremLI!VMb9Mp<0f8`pM)d#lW^p|#82`iev+Sr zBlnYV9#E1g;ha3cO0xYtO>?mPi;#bYHxQ%V% zHnoY{Y~ohP?-UhWW(DFK;s_UyBV0U=aLsYkUYT)}iQ>gjWK5%VsE&SFN$A z1slQ4mDQkH;19wnePZ6c@>&z6WOWPfEmyjDVlp&|`P(FMixF6=FJ2%6%lvHp2XW1` zl!-XjN1ntn?J9ji zRYbW~ntqRIlYUR0r@ySy`g_y^agD#{AaW_vH|CFNd;PujNV_-}LOW&<(lS}3Li!(g zU-VBpxO@b@PD@MBHSN_bFUpP4y-eAh>dBI+bWF4O;X|W4l(ka2zIX(2aI60s z|5mE>&A%_gVtrn=Gyo_veTd+@!Ru%YsV-PF%`KTm$YMw;@t+lhW>HTEj2JO5CeJ;u>+Qh!otU zwBKS{tq+lWuNF}F#Jxvf-|++Zl2+t!jlQS;H>Q27A+6+XRr|DSKGcv__`SA6+K(fC zB`X#?4W^x=zFlY1Q`V~7d#m?rD66Zcy#i{h#-ozM5tBHL^`gN}y->ovs}vuqP)sYV zpSvFS=q!1;le}E3coEDEem2Gr>z+9BB ze8>4Rzb8Dj{BHbrF5|KQpT8_$is2%bgXppNQcuK>$N#afrQY`2KdZ}euHU|~D`MyR z<%^x_w;$O(!*}1VwOs4=yIh0epZTv51o)X{~m5g0+ zjErMv{06|Wdl~cQ*fYn8Ikv~~SB{l(40IKUW0o9;BfNloe0=gA+8|WjT+d+4LIL3V^=q}LRAdblnfYyT6f!2dKRLF7JO`y%7 zdqG=3_klP*xE1sO=%b)P&=BZjpkdGmXcV*!v>mhqGzQuUdJx2M*Il4bfF1(v2JHdu z1?>av2R#gW1ats&5OfIiDCm=*$3TyRo&X&NJqdaWbOhvrj)FW;7L)_!K|W|4Gy$3f z1)!%v$3O)Td+S9|2x8~$IH(NbI41itHPAHZ1c>97{Gzk~v=Foi#4=`B%9(KtOL|caSW3W>KEe}CLh!_#xYDjsCSHGn0!zN8OJdBa6pV>n7rK#azNB$@^K63Rud*4 zx0yJG$;U@b9K+<}b`!@i`JhfSj$!h#2E;gq$p>|vaSW4>yFiR%*a6)QVjRQd?H&;0 z$P;;!y1_7cV;x}}!{lu}h;a;)57r&VF-$%-f*8jz`Cy%59K+;eGl+2vlMmK4#xYDj zwtyJNFnPNV#JIH}@*(w+Ve;{SiDQ_2eAL7-Og;up9K+;e$iy*BK0apR7$zUXCXQk9 zF#=*7!{lQWG-AT!Z5xPjTS4S)J7|XqleaMw$1r)@Y2p|rZx5O{hRNHEaSW4>M@<~Vwr&x2lO_n0(Ys9K+;e+QczTK2Df8hRMf_iDQ_&oiuR^;4cI% z0$l;R60{h!1hf>i4741y0<;oz73gZvHK0|XYeCn6t_Qsb^j^^WK<@|L0J;(M0ni6Q z9|CJ31$N}99x&?GA=r+(tK(~YL0IdPt3Azh(H|QSF0B9{}9cVph185^?6KFH& zUeFfMeW3e6TR~FC2N50ueGD`V8Uc-hwt=>Tc7Vn}J3$YEJ`UOi`UL19&~DHk&|VPh zh3LNMxaczVSM-#+NgWkE6nzuD@^CK;qRvoHC~L}*GGhLjH}XxM_c$-GeiZv|0bsfRm2 z)Ekze4=NQy%v}p}UI7zVs(v}}3W&XgoXYzpl;1X4U2tm~mwBV~z1NRC_qE|Sj{ok) zZ$AF1;Wy9jc3$4LbokBJj;#3dt8U9Tk~in+&JWYS)!AP5p-y=xd(zpnZ{LA2he}8{ zr_;^Py(QOYe){7xE3cV*zd#!+9d9)6Se+c0=*M*U!fBY5J!dG@!{ z-+un7xi_ADaqlnOGlS0_HMr^;fX?=CYRW6-o!yw;=!}*pYEyo(T5)g-DGi_R82^04 zx5hmzxU;7|t8krlli`%Vu9BNPkzcpJ98Q$Ismi*kFs%3kWiO9z9msmQW7V>k^C^f&^&C}a7v`;~?9b1i z8`kCYc0J!p+zsV+b|#4EjR(KWX%2L8LD&8=SD{Y+S!_4v2YqbtZ?LU3(uL-EAg=N6&Q(IF45+geTmhKLN`%?YqZ~B^Aga zhduk=o#W(c6@2i!cH51$a!$|(18OzW*BcOe0_RXyThZknTPyx(Jxfq zrMx=M;kE0cH;+Cr$#U&#!0c-Z04utyH43|k#Pb0v_X*hlM#^-GVaD?4 z9yi;FGzPZstzPH89ONdvcUv1*4T%AeH*StboM0*Y| zy>?j-&%y5d<@8NmgR$`&z)agnHu&~{`qj~^ zHNNdM*c}V;nJlx3s#SegntQyA26-4WH|vjwWuKJ@JEpkVT&1V7!ON`yG#Q2hHObOm zvDxNfSXb>%g56uZa&FR9qn=Qof;W}TdoBj|$sRxQCpvQ23RtLyqTMM9rQ7Rw8slKb zV7+X`t@BkcX2{SVxVPtul8MG(OM2M4j$}e9vaCmm%RUAjBFWWp6^HBZPb>^c;5hjJ z*1GDJN?5#)O*!gNI%Kcdc=hJXa+95jl)D&~r(9K+deUmJNe(rju`?HdrE(3xYCM-s z*Y-6fUxi_9!f3i+q1?U6Mr{BXbnG3fM#>-9ut^Fr(N7pxO4I83#k|Y+lBSH-FNU1n zxg-1FkRP7e>b1|ieW5hAb8RZLswq}=xe2f^5tf7MrY>N^oWz8>!Vwc}xvz1Y>x}*}QS<=B-=q+qku-d;@TN z9QLc@$GsuvM5$-?kkE08s+2JltS8&JCA|1x2to#$P$(8Zr=p3(ZrzcDcG>E zP~$K*N_~7>HKuO%BsvtxcJFWXkRf(rn9Z2-u=jh@)YDascrLmBg$D2M1UwG@qgO~F)b%ynyK ztKn+M38qSgAQw~z{93{1f{LyV-yM$yCD)(fKr9<5Jp@umv@hjK(XUL6V-5?yg0z#N zZ6l>en>+T2Q^*h!s&YHxdNvzFDF7>4ZtZ%TK(;py*b}zfKud~3QEJ_BTqF&k@8wUb zorp21#jj$h4(bz$@d5#LEo@}>0s+x=#|MSRfoJRb!J%q*;H0L#g&zrd^%FeD%P?*m!I6Y!8U}mRU>& zkwS?3wtyGGysw>jtyM1oefhfFwS`T0Q_cJ0%XhyRWTAHbWNok537@R(d$5UfiT8Wu zil3LlOOhS?3(}u*!uZR z&@odt7+*8vd~9g_J!>}HUA^<30p;DM6r%x9MSKL1db7!bV0r7;d>az_i6ib&H0z+$5>%*9;w zw}}XS$l1GZ?^sMwT`u`xLK7gR(0db*_j3ggpY~c78SvNAG^2fQQ)CQ)47Ip4g`!7c z)wwQyxf=LQ<~xnrbRE_ma{QSvuTJpkQnvM5%wpeHe<(n8FiVRBS;}@hEJPftsVq&> z)96DEJA3skC58CZyM2AiCD<-wNza9~bIs+hrWmQ^LV6Pitt$=x?Wbj_j zfVRVZK`~dT;n*SFeZVKAiOH~1AFZ=}8^vN?HN86qQx~E=plX%fKrPOz<3u6*LeYTn zFGm*wXcj1TU@`@V0@$8OdqCC9mur~drtWWfcEL-b*j?XJ($th}EhywU zrP~;X@PvU6F=4rJ{Xx-zK~K(S)J2H1-h23Ck3)e2r^I;vdx0$Z*2L zURwREdeW7F9#*kvfhV2TC76aoeGaM%E7o@ExV>U?#zevQ##?yPe4>=CW|Cpf593M-T0mEqIa?TmkA<4pRo&IQv{Ll%BNgyA8@8_o1z(VxVvWW&NTG8d*Gg! z#L2;YEx=9(MlJ1)17=Iq?QOQh; z*+`+O?^Nvwn+mW~wKIT`O7s%oGnGkfjmLv*c@cMMWLT+X^H|Tz6-X@;f>OLzIL-_@ zU~$Hx>oj^e(`~`tkTcNxONyz>C;QFe$$oQqs^1)*>dv86Maf)>`z*0L(6^paFO}Uk z%^#pUclq;&Z6W-@(_b#-3m5U>eR{W?zEx*zvQa8=)polDlx+I!+vNs8E{nr{8&MOg z0d0rw)vS26Gf5MWT-bA%3k$6SR05|uOI9VNy6l)I-!z(3hf4{_9VFX~4k+9dS*~>3 zw<$#_=fm+SufzTo;L)DLu2b;{*t#j@*xnr-Cc3E+PBo38u^gJh_T9N?>-;91=mOIc z@O{SX_`)66DSS1X>T{t?6l3CL)!`7e77rOVn+ui z?MwOTO1rsF5xl)(_k9f3d{Pg_)dC-kC1U@apuU{azEH5NVjST(?r z+XS`l|88F>F+hBcngTKpZ%__^7R#Y&9O8g#HAo#13!>xIwd*wi@wHKJisY49Bp4E& zzz1=hS)%kQMV1=@WvyMkK}ZY8?|4`*wH=QlLsc||uyuWy zJ%=1|sOxb{$iAjTOUWlOM@0)G8iz8?G+&U zVtjNs_2EHAHXVDKreR3a#@(z^)&V(x^@UXS3)te>%+ywKcE_ACS@*hiKoqXk!EPb58l? zkUct_@src@!FpA2aHzgW;Y0{-BM-G{p2d&4Y&$Aa2b4}pQw_Wjz-g=}mo>am-LnVO z8;Y1+9U`14U|-7P)7<#9rJ@OMBjRAL)C*aJdY z)Crm$Anmgp0p@*<`}xJ3*ejK3dq56(xWyW$3u03@p0R=H@|9?^qUEihb zJMI7fTL&~9LzXCKW?4_<&to_nyPaV_bp3UJiy`o5Hb?K->JYH&dw_l4CUKfpgQg~C zm;;gteOUW;nO#mEvubh3SNHEJ*qK=EMhuf!ZL{OOy*-Qo(WP#)gZ>2p204nf*`)#O z`(oY_8`1c>6yY-shia*8KTlp$fSm?M(yYDZS=F^z+ax(=9bRj*DL5Z|x39f*oY=9i zQEpYh{zKu#Gn>Vn?JZ|~>oc=)9Jh(DU%#cT;+#GMe#M+V13o*a&$xbOPM@Lu>A^p7 zRKmf()alI`*y+xYpPzqc*qlLb{Bv@)u;H)aEb?B?kY6cha4$4x3@hjk79x?o$nSDVE&fumxL!OMyy{YtJ)!ZA;jqdr( zD>n^#id~$^uC7-PoBNR-nDUDBq45m)`S}LiU~ z&yat6PV_CokMv3M+rht#r}ZT+hs*La^(dw{aXEFm{$t5ce%G`cQlG`niCwB+OFg>W z-l+mN_V{vp*S&p?x1*w0m)kqjkBjOZp@~2}kZH;mCbSpX5vWBtHp9?kC~MeMz6>OZp@~2}kZH;l%IH zsn-wmU$uYE`1A2=!JnHy%YNvJxp~#k&CZ=Keed-n&pk7D{`l{1{O04I9{dX@LpS)h z%F7Jv`(GQ@_dhkPpIwe`@ahvoG%bg?nc3&yT8E!G5@5ucY zoDTRS=68~OgZEB^tML{_*Yf-A{rMPv948ayKgMV5xp~Q-@fXVd=y{W$7+>>pKlXek zx)#3W_r-s=;>G!j^J}i7{54-ka=TdMC;2et`uo$8U+rIEuEk$$uCF}397Tq@c%ipMqH~X@-TO;QFh`xWzh<^U7 zBl`KD2_ME^VElzKKZ-xC`CDXs;V%QezV&=|*~-(4*T4AoO*j4P;E!Kae2nPpGx%NJ zzslfV5yO}Ive4iHPuUIBtJiC1&TV?>^x}1&|LE=i_{+h+UvF^L^}mnk`~MKb6MW%M z{3UgKg?rlI{$lXk+Fp(5-x1uG4Q`Zwb1it{FEjj}H~jvc!F_Eo--P>bhTk6;J$=Xc zVh03I?+ literal 0 HcmV?d00001 diff --git a/test/tests/userobjects/radial_average/gold/time_changing_test_out.e b/test/tests/userobjects/radial_average/gold/time_changing_test_out.e deleted file mode 100644 index b07f78d1ae6c161e3a965ef6868d562919feb646..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 49988 zcmeHQX^b4lbsmbew4`qNlJAye-6BO&6eUVr?=HzDO=)>ZBxTF8-Rha?nQd`qdNjv! zsRW3C$Up!Gv0=-wVkLor7;uo-PMkRZ7TC3sA5M(qR{+C^4H$?KAc76hh64$LP3HTm z>h<*WOixesj$sjGW8YTw>#A>FRlRyuU03fvaPW#kp|A*r_TOp)i;d{sC9^sNa!u9i(ycT)P6D1$>&DBiuR9*&+N{AqzxN2t+>?pni|1|*Ly$g!ub-FBS!G^c+ zwocHg@dt5HpWdT%uC(YisuZl_w>xgw?QAnaRNmc+a9dQkGmD{#Oj$^BT?x;WzWMq~ za5wwKo&dOKSmu#jl!=6E3*${8e6tlcsXG#uc}@Q&T}PC>8xJ10~;1jfp40-}_U>q5XAf zpz2-|_soQNc_1Mheg+<{_x}|Crb=O>;g>qpZ30o>H-qH2H!B;|IomXqj$y8R;?&#$ zWxbTHUpWgNt`GhT{#~QO_y0bAGrHj4@OgcAJNZ*8d|Me(%1YZT7oXpyZ|-;*?o`Z% zj|#iOum+dLm|L$$#D^rIPrQBvr}LHQ#`v? z$)O6xw8}c(+wqK-r7U+-mK&8Ug1N&lCgfq=lZ)~=`sVLW!qxQ(qPQKehs?+<^YZ|?aMxXc$`bmJsEVldz0&m@lU zuUI|D;Dy5K|Nc3zuSBl!|N8-Tt&0026HeTJRrjx~8dszB{?8bGI-we-y4WP2jJ&O&y4aZv)J2KX37Z$C+>;Wz=2tG|;xAn4_5V^Y|I7C8(zbt- zvZI%AF4`R67?q4&ag2;(XPg7z*gePHIOff(rw94qA*=nWu_S#lhb<83mA z#(C0rgE&|59?*M1?*qLb^a0QZK_3D+pj$w?nm0njecgP`4@hd>X5_JAG%Jqp?jng)FYGy|Fi&4KoT_Ja<9=0OKR zhd>_%9R@uHdK|O>Is!ThItDrpdIIz$=mh8_=oIKF(8oX=w|yG)3D7ejj%R-g^l6X_ zdKTn?il7px4Dvx0P!+TY3P8_+PJ?P7_SPGq5JcnkJg5cYI41itUCI`L~9HKMy zQx575dGu2b>JfSLQx57BdGu2b>KA$RQx57HdGu2b>K%FXQx57NdGu2b2SgtIl~Q=gSt*0{gmS#5P9@F zAl3o$=%;M=fykpwluhae{gjP$ggp8w+jbCn^ivMj9rEa>96Lee(N8&8r^ut9a_j<; zM?d9YT_cZv%CQ?n9{rTxE4f=%fQ?_RekABMb zNyDR`vVF?%=%;L-Haz+%n`?OVQ?_S8f49{rTZc;}@4@KWZuRJ^}f~YgplMSFNKw0I%aN`ooj^)BUXMQrT zm=8=N>nmkr*^fii1(q|*kmbX?XTCDem>*0l)4?>1L*wQP^J5&6{G%RIhpE3zKlPS6 zOMRuTGCxFDZ-Jk>MV;9PVi}8#k-8Fs*91KeYJu9I4yX%a-8%zf9b>(s49rXF|7{@Z zEA{XIh^gJZs4=pb~65M3GwENX_ z)cnJZX1C*D8L4loTZCheH#~iUBW{{JNj}3rgC3H^OLLR)((YA!#|aulu{iC#cZsa0h&_QGER$5`*Rg1%EWw=p50Z%&%b0_CdoJ_?Ld6j(V z2M+lG zm*YIM?f$4=KLH_Gt|JYYeN6#i`&jWgtbm_`T|{4_|1!p5?965Ck(xN|Jn5BTFXyli z$E*4%!_n9T^mJJF)G^*Au4xCA-Sw8Lqj93AjIJzK0w=mn8rUdiibny4-A37KosB@^ zo;-l1MZa{~EgD#?VQzO?L5UrCTzLuw+Z;Nru!hys)jHKkb-~6G&NAK(tEV^MFZs~EWeg|*U zZo6Bl1k2=_D^E>n`#S#xqtBgownx7{j)#%Ect0rz?Pkp@`DGWqhQ%0U9zgUIg4AM} zb|PCe)))iL*S7m7V74SX*qwW8Uu;vJGI{vcY*qT(l&Bo@=PWg@;?;xNaO1_QcjVa7 zV@FW+#w$AlPQH;cy=<7VGWy&7HX@IKop`I4-j_puOgje;KZ%!saWCd!8)Jdk_DS#H zl_$v=?eY;0mx%K3^V%GhOus17$XB1}#U~x3@5yk$*sD*TBzvuvm!{V!6nl$3@ty<5 zdzxWGjdu+zW%X*-mrtW^h>+dTOBr+6I^!ZB-WL62w~%~w<;lZ|_Z%uqBF+#ukA&D{ zjbV|&ev(OY{|k6V78V#uHS}GtrT~)v>?OSP+GRaF2fOdL@;7y{lEO`#%PEF>GOwrQx>2 z8b+2XUeKrE_KNM5N`6xz;wS}-rY)metj}Wz(2G|i?>1%|67~T$gOZ119dco*^c%JD z%1d6OG;FxVzLBCHa!#S-4sG!50rjh+*NwmJ#7vHTp|FYSc6?Wwd%TPWWf(KJ=vTs) z&q{>zN!((oJyzMUiPC^3!BC*SytG$rws{!VRcC?VtQ4&`*Y7;W!HmIr#kSj%t5?j(Lx1aokdbmo$?GN|o{D?cpDF?9D9k`bKG zjj&aBRb3iOt06fh)PTm$tOHAB1OPk9RQBE5*OYPxhP4l)?}mkP?*le!1Hd4)cc>aE zf9sA1MwDy+ql|H-JguJJD7&1J=}TzuX4)BFI?4(OTP++}f>v*M+U*M^-p&oF(5j|b z)#X0GS~YA1oyB@Q$y!1Cvp)R}0|c>_7&gDe6DnHusDo@@_2PhR>Vy8?_WJS=8aY|O zAvR1Eik0~as>m~+$9kwCy0(jqblLh#yKiS=9jE3FQI>cq+cp7QdNfcJeSq6>rRz!J z2_+n6`6EP~hgj1v<^@2k5HSsWgc9>iJ7(rK51`xhN~bX;5yOM!Z{p;&(wJ~G4$Q&j z3XU!7Hp^Ie69=Pz+S#>p&#pbYAKJNRtaJmgQi1(So!pyt&NRms4*?ygp-LG;!OTpY zr_#iW4~8IQpb3RyQN|h&t0h{UZc`nP7a0Xfy1j8|Xw5DMOGEXCYICoZJv4xvlplte zTf#K$R*G#aVOnyv+91h0Y^olTtPC~Zb_R-sdN@89)eCB)twW{TKtF%*QfgzNI1%7b z@qj<=SY)6`NJBl*QXXRYqqM71wrT>{Xkt+J1=z5#*5xoZrhh6G)tI`)v*=JD+2g;} zBM))TS+y0Ehn0xR!|}mPDf6UARFhOdTk{nnd6@65h7F!Wg@##;yGWClLk~xe%~5S^ zN9UAjXrSv3xMn8C59jl5p_lCx?W393-v#U~j(63y*mz9e9V#E#Z*9 zrpiit^JOGV>)_Z_J6%(vcf8@99DCycQ{yeGAubH0l{#T;%@@;bFQEd&2uRf9R9P4O zvmcPu<2nH4S%56rHwI)ego4d;zs!MHHc)y9q=#u=DD8&du2-;zg|i^T^3b-C(xc6(ec~K>h#czh@ZVlC zfVvjO z#a$;Ly6%`zhz~s57YfDUS=#|l*Tg7Ovh!fCHgQ0e;(o|=0ro@ba>)lM9qdH@5mM!$ z)o@>jQ%)|_J85S`LfiKg1&Dq_T}>UeJ6Hu)^4smuEG)Kf_2Pg7_L@nyZ36VjfF33A zI%;fEWogmJG2aBVO)~9j9qZKt4`btl&9glq)>{^_7(^Gr;D1}di(uW?LA=(g7l4V( zE{|+s^F37SzWDMzDFIomBd4tG6+2_f+P(*yxt932*J}G^DZH$lWs)&2`|0nF-Odn;kV|RRaH~46{3H8^ZE=K)$O7;Ux$3GX9`MUbfgre z@+Mp+dE)C2Ws3p92d+HE!7T0=z*3AnO#m!5`oi4I)q9(W(5IcF$By!xhA}Tj%xD6n z6vl4?%6_TlVQOz+lL3D%Pcu64Hbuq&DA2%_r%?1Ltomw--|7T@pZU(CHY101vBzgQ zEUWM3j8)D|;>34(>Ca#7J!y(nlm{U5Wj-uX-aJ(LCmBhj@APxI|+vinL|PW94@3iplaq@ zT`cfZTcc|MK6a+|)T#~9 z^3;^JTFsue>)k5DG(r?I$N0?CqgfYc6_j zE5s)R<0I@Msl6T+{muxy9Ev^iEhSG)DRzTenM=BjaR|@&SjQi2n6)pI5;i%s!T`Gq zM5Ce2dVw+X)IJY;<}0#S+I&UU44~t;d>L;kdGb(wq%eIY*|hWVr+FCJxEH%3u{EWl zeLdGEx%|9+vB&(iuWuS*0Y+B%VMSD4zM?(G8 zQ5RNh=+bd}#pa5Mn(tKxWYcn@l&z*Kux1I<)$L%&D~|wL;X!Rsg<%=E&0y#%$!X`n z{K2_X3*+v-Re-vc!;*m^m;EZ!;e{JPV^OG8nkgZ*#Y`mOXb24dRb%^HvD4qj9L~5^ z8F|=x!PWKwU$7lgG}8{xt{nd|4ZZFjxMvpeonF2c;G_eimiER0t0j8&Im)`oW2@NX zr)A8LYJvu)W@4{z?+3yU7U?xmIiGb7BEM{OECB@vjOizHh%QfTVasvTi-0j5o!Vm+mCD!Xl%-$8fo^7MyoA^gt6UoPeKXYn(~^l3R0 ztIjxNquJ!H?O_`z+4R}B%P~N1i^F*vQ4^{GZHMmFYE^YINgt5h*h^RnYl8z+0;f7l zRwbpn>{ut?H=5Lg%L&L6B!`R+DBK)LZjU;*DMu-n!%E#tIllrj+DkZcDINg_52c(x zdLU(?n;PL%%NXM2(C4@B&Zz^Z2lrG6s{Lz4wLalOsW!UJTG(ITL_E0mg|ZlwQ?*(F z+ScWfgDjNLc7ZqF?OIsa5ld#=YhTLGwTG>JiZJc-aN=XIW=cJ%s0}_COT>XZL2o;y zy&hs)MFsD)UNe3qYLwy}YSRP6F+qdpzuOl|0uZxNbs+QO4eA1*#dc^Kha{ld4bnrz zhG@$spo!T=eJB#P`J+fMBwWRWIKBv?^eIP@Z2>QWEjhT(j3;*+XIH`YX(ZdUenaeNP7d@!oECg(8G~m3K<7C;ZrZ`0i68u7@6U7yMSp9v!~%Q_}0f##L}gsAi<_MF<`v54GuE#Si0P`zF-`luk)s z4O|c4684kJ9^R<#*#qhWMa-!V5l-Z>FXhS-4?Z2JXo5Iih3nZ~O!@0$;#jWgyc{WP zgJ>_|Y4wd9^#9FQR1zDl?E#@I>id)_koHNA0Lwnd{rpBr?3GHiJ)k<*tkLCiK^*GF zZ){+>d^=h!X`keDEF+=qj8C*?sfZ)$1GFu82WorR|E*s!>?On`|56Zk+b+)OlFA3W z@O&4Crea;XJ)o&0=_jAq8i(J_5+kpyzXySJej+i3joF*GPqI4sg!y+C&AlA{YU*7+ zKpuX<+j3IAJdKn&Iv_fapsxY2FO-YmJ>^u}Z|N^+ zy+wE5s!1;Mk*tZlO;JxqzDwD6+5tG~UCO=*$awWhe1@#o0H^qvWWA8TjNxn?pN8|$ z_3QvQL*UOGjy|^4A>hdO0QT5QL&#-!O6txG-8;<&X5!D z^W>ocM3;KV3Hs{-3`)$>{wojQ#5eOs{jXs!^|^MovQhGn$ zyz(_2f90z$D$mS%zNZcZ^S#~;#vLM${=HT8z4WZ_sW*QvI7xT)oZuu~`HuL|@2Thd z-Zy_?)qE#B`cAv~p8AaWZoB1o#AD@`koRSs9)o)!12=}F;>GS0xY76Ob-q~Q#>z2) z`xkxkt|NczKRo!|NFN?s`J&){d6V_`)C+}JKAbvm>W$C6ZuI?&D_?t20o#2~J^Ide z4A;~vU0)Vf{$AvL_32k{^6~xmk>68~;plj2`kd7JwEML9Nq$=VBtIE1?LHYV?LIAj zlAjho$xp^hyHCbTyHAUs?vwG-?$hEY`DyW!{A9ee`((Vd`?UB;ep>uA zzKK6}eOWd2?s9uKmYt8;^`sq2*!j!t-T7Y~di4Ct<@QdsZ?^5-npUqumOE?0Ww}{6 z@;nPio=f;FSHfqxSvc}M3rC(y_$*h#XSrE8@;nPio=f;FSHfqxSvc}M3rC(y_$*h# zXSrE8@;nPC`P{_YAn~(Z-OpSrsDADnt8c3Q=h4+SuV32t2fuuI(AdleRyo;uLW1X-x+=Q#jIuz=--;*8#k5@=KW6V^A}gX z^P&oA^!p1LdL6@2(PH=J{Z7}Hg_TwHykP1N^jANh)$F9a8`h@ugDiK}gv)ZXaO8Ow zjy#v}S+0c7aW}{C#b5vMS8E%8{l~kG{NT;C zjh}tx+duxlZF|3}>)qae((bH&{+n6zeAdLDHSuTl@BP-e7n1TN`9_Wl#uYg)tWACG zd#~MkcX4gwFZX|Y*MGdR@@@Hjb^qq9d7i|R-y2AOWItZLdE@r)ytKCI%9#iMxpsc< ziw0NS*9)^#>iK6PT%Dc^=3e9*bagzzy|D619Jh8p!WH{u$|b^$zAoFKS<&k`8TcIpu=amvnE`Yn}s9KvvB0OgwJv%e3qMqBhRyN Date: Mon, 29 Aug 2022 18:37:25 -0600 Subject: [PATCH 2/2] Fix MPI bug (#21960) --- framework/include/userobject/RadialAverage.h | 3 +++ framework/src/userobject/RadialAverage.C | 26 +++++++++----------- 2 files changed, 15 insertions(+), 14 deletions(-) diff --git a/framework/include/userobject/RadialAverage.h b/framework/include/userobject/RadialAverage.h index b0cd32df8c5e..4770ac0ee563 100644 --- a/framework/include/userobject/RadialAverage.h +++ b/framework/include/userobject/RadialAverage.h @@ -117,6 +117,9 @@ class RadialAverage : public ElementUserObject std::vector> _communication_lists; bool _update_communication_lists; + /// processors to send (potentially empty) data to + std::vector _candidate_procs; + processor_id_type _my_pid; //@{ PerfGraph identifiers diff --git a/framework/src/userobject/RadialAverage.C b/framework/src/userobject/RadialAverage.C index a676612916e4..6d29aad389e0 100644 --- a/framework/src/userobject/RadialAverage.C +++ b/framework/src/userobject/RadialAverage.C @@ -123,24 +123,18 @@ RadialAverage::finalize() libMesh::n_threads() > 1) updateCommunicationLists(); - // sparse send data (processor ID,) - std::vector non_zero_comm; - for (auto i = beginIndex(_communication_lists); i < _communication_lists.size(); ++i) - if (!_communication_lists[i].empty()) - non_zero_comm.push_back(i); - // data structures for sparse point to point communication - std::vector> send(non_zero_comm.size()); - std::vector send_requests(non_zero_comm.size()); + std::vector> send(_candidate_procs.size()); + std::vector send_requests(_candidate_procs.size()); Parallel::MessageTag send_tag = _communicator.get_unique_tag(4711); std::vector receive; const auto item_type = TIMPI::StandardType(&(_qp_data[0])); // fill buffer and send structures - for (auto i = beginIndex(non_zero_comm); i < non_zero_comm.size(); ++i) + for (const auto i : index_range(_candidate_procs)) { - const auto pid = non_zero_comm[i]; + const auto pid = _candidate_procs[i]; const auto & list = _communication_lists[pid]; // fill send buffer for transfer to pid @@ -153,8 +147,12 @@ RadialAverage::finalize() } // receive messages - we assume that we receive as many messages as we send! - for (auto i = beginIndex(non_zero_comm); i < non_zero_comm.size(); ++i) + // bounding box overlapp is transitive, but data exhange between overlapping procs could still + // be unidirectional! + for (const auto i : index_range(_candidate_procs)) { + libmesh_ignore(i); + // inspect incoming message Parallel::Status status(_communicator.probe(Parallel::any_source, send_tag)); const auto source_pid = TIMPI::cast_int(status.source()); @@ -279,14 +277,14 @@ RadialAverage::updateCommunicationLists() bbs.emplace_back(pp.first - rpoint, pp.second + rpoint); // get candidate processors (overlapping bounding boxes) - std::vector candidate_procs; + _candidate_procs.clear(); for (const auto pid : index_range(bbs)) if (pid != _my_pid && bbs[pid].intersects(mypp)) - candidate_procs.push_back(pid); + _candidate_procs.push_back(pid); // go over all boundary data items and send them to the proc they overlap with for (const auto i : _boundary_data_indices) - for (const auto pid : candidate_procs) + for (const auto pid : _candidate_procs) if (bbs[pid].contains_point(_qp_data[i]._q_point)) _communication_lists[pid].insert(i);