From bcfd50f1e326717836485f2228c268ecc06cd849 Mon Sep 17 00:00:00 2001 From: wyq Date: Thu, 7 Aug 2025 15:46:16 +0800 Subject: [PATCH] add lambertw function --- .../java/org/meteoinfo/jython/JythonUtil.java | 18 ++ meteoinfo-lab/milconfig.xml | 38 ++--- .../pylib/mipylib/numeric/core/_ndarray.py | 5 +- .../mipylib/numeric/core/fromnumeric$py.class | Bin 12880 -> 13638 bytes .../pylib/mipylib/numeric/core/fromnumeric.py | 29 +++- .../mipylib/numeric/core/numeric$py.class | Bin 155059 -> 155153 bytes .../pylib/mipylib/numeric/core/numeric.py | 5 +- .../pylib/mipylib/numeric/lib/_util$py.class | Bin 18242 -> 18191 bytes .../pylib/mipylib/numeric/lib/_util.py | 2 +- .../pylib/mipylib/numeric/special/__init__.py | 2 + .../mipylib/numeric/special/_lambertw.py | 155 ++++++++++++++++++ .../org/meteoinfo/math/special/LambertW.java | 103 ++++++++++++ .../org/meteoinfo/ndarray/math/ArrayUtil.java | 149 ++++++++++------- 13 files changed, 417 insertions(+), 89 deletions(-) create mode 100644 meteoinfo-lab/pylib/mipylib/numeric/special/_lambertw.py create mode 100644 meteoinfo-math/src/main/java/org/meteoinfo/math/special/LambertW.java diff --git a/meteoinfo-jython/src/main/java/org/meteoinfo/jython/JythonUtil.java b/meteoinfo-jython/src/main/java/org/meteoinfo/jython/JythonUtil.java index 5eb5e8a4..15f06dd1 100644 --- a/meteoinfo-jython/src/main/java/org/meteoinfo/jython/JythonUtil.java +++ b/meteoinfo-jython/src/main/java/org/meteoinfo/jython/JythonUtil.java @@ -14,6 +14,24 @@ import java.util.List; public class JythonUtil { + /** + * Convert jython complex to java complex + * @param v Jython complex + * @return Java complex + */ + public static Complex toComplex(PyComplex v) { + return new Complex(v.real, v.imag); + } + + /** + * Convert java complex to jython complex + * @param v Java complex + * @return Jython complex + */ + public static PyComplex toComplex(Complex v) { + return new PyComplex(v.real(), v.imag()); + } + /** * Convert PyComplex value to ArrayComplex * @param data PyComplex value diff --git a/meteoinfo-lab/milconfig.xml b/meteoinfo-lab/milconfig.xml index 41067650..b30e3faf 100644 --- a/meteoinfo-lab/milconfig.xml +++ b/meteoinfo-lab/milconfig.xml @@ -1,32 +1,32 @@ - - + + + - - - - - - - + + - + + + - - + + + + - - - + + + - - - + + + @@ -34,5 +34,5 @@
- + diff --git a/meteoinfo-lab/pylib/mipylib/numeric/core/_ndarray.py b/meteoinfo-lab/pylib/mipylib/numeric/core/_ndarray.py index db84cc85..98b23b28 100644 --- a/meteoinfo-lab/pylib/mipylib/numeric/core/_ndarray.py +++ b/meteoinfo-lab/pylib/mipylib/numeric/core/_ndarray.py @@ -100,10 +100,7 @@ class NDArray(object): #deal with Ellipsis if Ellipsis in indices: - n = 0 - for ii in indices: - if ii is not None: - n += 1; + n = self.ndim - len(indices) + 1 indices1 = [] for ii in indices: diff --git a/meteoinfo-lab/pylib/mipylib/numeric/core/fromnumeric$py.class b/meteoinfo-lab/pylib/mipylib/numeric/core/fromnumeric$py.class index 1af77dbf324d09254767eb5aac30c78b972e64cc..4c8a574a39cbdf3b6c9335d84b7381ebf3d4bfe7 100644 GIT binary patch delta 4327 zcmb_f32;-_6}{JztZ!Ka`8zg?u@Sae>@S-oE&)RcTS7upLuLWQt@ut2tgWQyHR zt1dV&EEC{mCbTo@umv~MmI7{s%T}Ol*;?9^76>JfF6ojsHA`q~H}t&sWmzzUnNH)G zt9Spq_uPBlz3;vMeWv#064z@V?LSOJqj=#S*Q{P?+#weBS}7j27sP794Yk3>V7S3s z>0KCzL;}tJP;jZ~^DYX9LknHr#NV7cGh%vYMV77CXx27%}j+*5%XtOuaQ0wiX)nJ;n(q)n9 zT@(oiYHI?~#=xSGDd(c$dea-N3p`_rUVR6~EdD*YCN+4~io7d=jdfl#WY(Jvjfvb% z$S#-ElPDNXtWV~@-MJahQBwiR)e1G*$7Bnc4I2Fzir-x3rUPT?yDloFgG`xz|6;Q- z*l5=K{cbvJ&HseSE{m$sQKodiVy4kgnN|(kSe7p~xGKcs`q=UlGR8?a#p&4$YN69W zh0LcLHF}N7iKq>MdecpB&|DYU>CLvU-#5WF?hKq{Z#4QPljD&9$ZGT}CU?{fL~81) z!;wa_R-<2oYOn#J&bjI5R`?5T-_gf$-fx*)R`XXx0?%mlj+-u84p-XT_fO@_-@EB| z*2Ev%=I4BaGym+SKUouhZM&Ykk~9DAroUMe{}A8oJ3i|gk_{`F?WQSo!$sHWLvg0> z$lQO@zvQa^&7_U557#aanG-bX5`FsZihV@4Tr?e|iB%joeei42*n^>2V;#d6HNGF_Xt37QI2Xl_%7P$PUkPdo0@Rd+BZc7fJMl1Uide z0NgI&-5}soK=108k|3TaCJpKnL*Uze-qZh(^r=h+K7H@N@9Td_f+~`rEB6L<=pVF; zZG-apDRE>_NnR6*^a7*}avE|4p#zHC~;}H~_MRnn=qG^n%d4u~S#^}N0rp!mg8LtGW7#or>Y;aMWq5AGi`JR#{W z{iIKvj#SzWM?XBj$k06X-Bx?jn%-k=P^mVkR5pY#^kdND^lg=7xZw#)2kTeb)pVxe zfZ7T}dOSw{DkJR`3YYCRJj+zH`;yT>929Il=^`l(GKS@l!Zw=Xke+l%A>61OijhK} z9H)pcW=LB4ap;PbrAesL;&ciOoIMPdE4)`xbXp0z`J-cLESuV!>|Hd|(M4Z&Om$MJ zQ?++A^@6Dm5`-*)gdp|P@Y4VZCuAAM2qX$=ge=D%z~vzK0K%^!dGAwNo{>c>!9<;_ zIe2IhBZ?{<<`1qvzXg7+ zlws&+p;Q@8iwa;^1a`FQ+HWLNMOy7y!)n5*Dnm-+xT*V4HhU|Y{lIXR_1u=S zRN2VH5x%_-`{B4egrM;j1)UA`4-zgp} zDJpI4#wFtGmsJ<;H~?^eS*$E6=I@D@OA7h>;*FBwIX{q4>DXOMc>Pd#N(X1|?hdM7 z786SgaJ%_S2N?1yw92>z+E&Om$acsM$V-WFC&pG$;L9tu$6#oOMepK)wk$d#CK?0zsMu=ci?0~Hb8u;C>iDOSJPvC`=7|)NmksX%dA>6! zt&{4yXrV(jy%4WgqJxk_5c#Ej1acH|400TD0)53h+~SG%x%Vhz>{qd^|% z9EI+q(ER|`U$6FQi3k%sO)9;X65;Kh^hk*Sc`5*cfpmU#ESY{HKGzi5lgj>GM>2mZ z9BVf8$}DQL5UG^k>p@tIt%M3ukwvF1L^s8lEdKlVfcUcLiUr9?Syjg#>!a4~1MOt$XP|MtG#`jxsoQJe?g|D=oM_AD^cqETX`0L7> zC|z%9=P_11=CNEU?PAr4>EgCm60HeD#;MCrvei&XV3R;9 zGHV%WJHtcZb5LPqrjIGK2nM!pQeL%{N42!2+SXPfJnX|_rPNaU2sQnE-`y<1L1ygq zkKcah-tYX*`Q3BQckdosbh6F+>WBLe6VViYV6V41-9aCSh3QWiwM#>h(2~N|&W^VB z?)J_OF)rGtwPy)96*PZOn$86O!>$ zL}>3omY?hNlAZY%v6Gq8c3*Vv1~b+xyJHP#`vOV$aTPhqoHV*{I^ z-u5nyO$6!1oln~>nHy^F5UqI)#zil!;XIwQd91X?i~hUFQ+1k04jJNRoz~gbtvY?vwr-Qdp?0Uq zW-04-S$$?tl3Pt;y`|GiihuFDij-iS)$Qb>t6r*rOUg;df&JB zb7D$COZp{Z{tTs&Mku-1QV^Joslf1#lTB7DdaG8nQLX46)-fnIiUW32>~cXrNB5fp z12m=>Xlg!1^9F2lmlx z=8ux>`UE?JR)g-6^br#Hqe!RBUnNPRqR1QFeq;Z0=9#2_Z8Ayh#-w-5-z7=aNz#ey zk4AaC*fv5(!4X6d`5EO%>I`V$WU5u!;Wx0qkuWVXp z6C7rPB1z?z5rllzuNKQ#2-`Ai;`E^3a>fZ;t&NuBMd}Rn!I@PyuhJC6^+II9D9S8E zj+G=|PUO%nqf>^IO+?70V(-LiSHLosi3=0UN*gRAqLOS(Cb7&y*m2q#r-z@VV$0|f zx0d@HMvr))Trro6P32`SFjt6!<>OuO=oN$IzNqo2jAZ}1p>|}l;&hZ7RWRK$k0bUD zO0~>aU^H4DwCy=Ljw_Xl65KR@5~qm$yp{`LERTBfpEX+Ai1TBkBx};eJQ=6Lq~>hP zSS3B=)=0f9BPzaLQ4p2A$|27XbqvvrAqwI53i;i$JybD7?W&^nI6*|6z*1lt5COV? z9$-1J0_X)W5(+D^J_f7;qUZu%1BosyMi<;@)v%R8YoW=j&P5BZ)}awot5}1%e$?a0 zHaiFCi*fqRwT;*h#|>ham6P|RbhpQ0PvZ`%v{G_Vr{}Q0QJk@gGvz)i6{~gSxQ-7W z!~IPn%U?Wc^9U~y-+WUQ+AeQ$u;Bib(qo9|gmL0=e+q9AU-cL9)8f1Ss;p<^o}?rT zn?YN}MSp49wvq7Wo5C0HVa9I{OodBIXkWsxiWIcDv!RAu&kUv|@C{1@}&@SPtEa2xwv~s+7uyRb+ zZYfX=X5tKi^#XRQ)A~k9KU`V(c{wDhgTPtfJoL??i-|s8pmcKL;0#vu)88s$^be== z5b1OMmiaHI)431V%+3sh2IxOA8g}7RgiU8A-hFVn#Yj^=F6k3BSUG*|BTjY36y@Yr zfoz{|xICi*xiZ<=kNir$`5}YU-S5bzhus0cwC+pbn@9W&ksR2H?&< zmQ;MY@l!}?!}zGXJotC$!NT{G6Si#*L*FObiGu%$N&iPq^ncWSLkZ;I`KXtJ@8cG% ascjU!l z0s6H)57#5~>nwRag_H9sk_XD`?|Ezw`Ed~V+)bw zb2FK=sW#2kv&XeU_Q(45^YS{(*U_Q;wfrHUm$!nD9(+nBL)p<9 z%EoY-Op3U(Og8c-t?}#^UCh4LXo>UswWSTNZ*puKN)%3o^V=Ysli!p{8qbhP2i_%< zSNM)hhHzY4B;RmBTXfZ89wF~aRxEByBWPF!MNtz#6VL_>0AGUnpadKS*A>NV1+gF< zA-7YqiY!DO%qK(^&4MHyx}kBd8sVvPc+pa8rI z-UZ{qHAS(RK@3O%^*}y=0k%0{1Na$0j_sbJ*!@8yNCdsW8vqru&jpVa#SsQj-T~ts zLzTLUV+<~ez((*3fFQ?BMR6(sf1GKc19$}t0pEbdpkzhg_B6Yho%|HLq9XSWc7I2e zs)zx_G{2aBET-Zwl?t9fHK16rdi;uMaiK>n7nbB7@LLd1oJFD3RTOGG#|@=er&`%-QqV+hYLR;DfWNJ!Nt797+tymiAenZn8$)cMfIKGrVjL#g}UU>#{&L zXvD|nkh^?K&7?twNv+saM7`9Oc;1R%E5h?m+VdgC^Ugf;dy1znJYz2A=k%Y*&L@h< zLS1=xln4`Qx6kKdqIT!5(IS{r=2C#G8_Fja)42?j+qCWR7BoT+ZZ(%GQBQtrE@te@ zQk3)zyXI14WN%po+3Ur0D#KKyBzY`uwV2kSSRW}gl-={l5iGO0O;W1rR<}LghR1!m z&OA!M?Cmm-qUk?nheto2HIF>hpMRc5P4haTZk<6F&=qt8-HjB(TkQ#62ABqFZ_vj` zF&|Zo7xh1&AL!3r=A(Vy)7od0(LO^sMu;e(4&zKABGf7{?0w!gA3gQ~Uz(5Qa3nhy zP%?eU84K{^DDJp`BDiD$mV^-t(7mI127aKiyk-H_qj7vs-c96?h4dVK#;-0!(aHSP zLdwB{b9x~h6LZ)iD(Leb48k;1hk^IO2Vf-l5R3w!05}qLEEoqSg3kaPiwY;CGKTC| z-XW_T(}%1w59ZI(>NnM>-)wEjViAxtLa0UTSxkxW;F;f3wZK^rG#ktT@D1u*{&X=~ zXt^d~rXitNm#{*UP)=QmhgyfbXgJ@zOD=~Qqh%cxLiH+5Z2EGDDh4aSO0WvH&hpwS zj$`_oF}>QMz-CR;97EF7cHSGbddLLCpTYRLMg*GqpDlLRepzZ;CjTCNO z-3RuA1K=QtuR+GtOzr#aO+lt`|59w z;NHckf25%_nmex4>+j&86%>*d4iymq>#I2$ow3^286O*+(ExQx(X_8Iv{yH@ryAO8 z7}}rFwb#_O%QkqHAFrSc#4uXxn`=w$S#c@QRUM>)8sHgF6FdvBG$Y%(m8fzfP3Zel?#xC;c3*P#x6j00D2m%^|CLjkD+U%{6Qm9u*c?O2H)?{xr zWVexa#q-0}5d!d!wsL)TFGkfC_& z!4gy^4l%#aUzAW99M&Er^$~NK;eXBWzviI;{?|MlHvQ(a=^r$UijlV&)f@MeI)vg< z-AwsMQ#avX)P2TX)*uL2zXqZHBu!@{_;B-h@F|!8;LFXEp!1|xr-)Zq~OZ&5F?iYttmY4lZ%V} zo-Op4^vp$yAD++Qduu6)udO9J;&EKp!pe}kj*@5|ca+I|9=(oYXaUcc$wJ;Ke_X^} zdco%flu}UO3{+z#_!fKzW^sd33aP1>7vN?gScE#B_tkNhVo|K{{h?643dXtoE*f=Q z&F>GO=myb%5JeOg1o3}kKMTU)1g&ZVW_*E zO-LIMXPP%NZ9sgom6JD6gljFb`c|Ds)A4*8zqkPr&4CRR*1QyrupVqM(v3(rfz1Fb zg?TO5hDN&TYouHHKoBD#=^$jg{H_~~bjWZrhq;w>tiN*fpHxY$iF`+R>PD*0ft#@F zJF$@h>6q?!PH@&Bic#wz%SnD^6IJBrH&GCsme(0{hUaXe%5;`@BXOPgwuu(f4?>X( zY}$;rzR1sRrp9!M-`Pwpvk$`HLjYSQ^RM6tI0}w|6W}B`4bFhG085D(ONkjviTM&L z`X}GmLcsy%+gjEBGOBjRsM_C}+weEvmBs|M#!BI$j&iApXnD)ZmqF1W{6w!lFk_&}rG(cHD ziw$Nz^7Y3{@O5}X_dJ1?d@0(b7^0i;{jHc61-$QFLiD_i9MqgMwjt6eYi|W^Z9fh&$vN zWXS(a$`9s~yQy*VM3nms! zOsM%D@)a@dMW{KKtL}xsdB~vqY<*XiYIrwfoX=hNQsrQoql}16c(Oo`MVIcS>i_u~ zm7WTwfp5Tc0AFqS7JLV0f!Sb=abJXVE|>@A8|eZxT8g*Pl>J7eFaw&FX^oa@G}>}w zf)yJLwt}bpjBO7#1Hb4qH86b{L@o!|66l%n+J@*& zoVt(VqzZtD@ zjE9e;Xl*hc=dAshj3@Q2W@8jR#UJdas+CSdf~=Wuz0;Ow%Nbc8-3OW)^A8oVLqD^h zg5j@Y&Y@)^4q$nK2e!ZiTYdxZz?S0x9@ugUzzbW>fU_{Ysn2wO#MEeQWheAs*G$hd zOuxbJA0T(cO{vxx&Ddt8`?gt@KY7^!N^!{yPaIeY~-8 ze;Vm6q<_JV0ur(&N<3FN zMDA4S;q?ZShG?5ls~NwQC?h)jNRL*~_(D&BZ(nM)$W<`fOlzUa%99SkR2#23gjmMT ze;lGHZ00S8vBPt6mBZN7mE-1zscNVb0-sc~RO=I~KaV|(9e#P9dl>uu06z9B7RW#j z`IRb$`ae~mLy_DQTV34zS8Otaxi_-l$Y{*3C?3j1zv5W6f?hm~oku9l<$9`k63V){ z_7O^@a6MZDx9X1TACJJDSNV;S^Bg!gbb@k5>W{QM2mpaV-n&3B2mzs>0tf^0cQ?RZ z#~Q&Ce!~)w0R1kLw-XYLPN=GPLN&b;JX$;vQxC=^A-fhK>5j_ON(&ogZzglcqiD?( ze(xwXib{ZnL{Jq}10IkBl0ga|J4)`V^)=0%jPB3GUD=K_LRJWFz|qHWT9w6}u2D4n zuQi)JC#a$-+o~ZCKZbs7%qLG!h3ZYDaWd%jhHo;I)1xM!4FR2?atjEJ&6G_NUK?IqnqC0kW*Bt`8(2( zC*Go(M>P&*tlAUV2T7iw<=PLTv{gn4Q$OGC z30gnl1*h>UWHg^S4VONa)!*s8yiqXvqkrKO{5Belg#&!cTM2PVcYtCz#2wQrVI1mO zkMqvosWY9By%x-wXDAknM&~nBNT+0|St(A=GP!7YTz*qz9bBfts6T43(r;oZp8c+u zw4B9iafTz$!sVRfI%hE{&vW~;REI9`gtLfTFY*ua?h^lT7Te=1a$U*c0p}<<4{L|@ zcYyW7dJdci7r;eu30y%{KJry%jJ|@4Sm?f^Rpk?-KmOJuu|wx*6y4>H=V8fxUNRYb z61)_!J}9#uTOZ0J6aDu`T>6;eg-v0{1!`)ukycl{;+culCro^s9}l{K4W3C4;*VTF zRCb457b&i~V)Mf*5);fD?=_Ep)MFRI@rx};GjD=!p3SAlF1BEOwYG&y8{*};)J7ba z+T7gvA)NKYi&PGt*YO9sI)an`K#oXmg}ZW5(lG5n%LVJaCtChoYY~=+rGHRF^B`CQ zXKRBKw!s10;DBv#yf!#q8xEOlksu1zPVreQ_Zub*?r2zBMYEO-YZLU?!iJY~L#26D zx%%rk|FKmw%=Pf$OX$ZWy;bmbCD>gB;u8QounivA1`lk52eu_~$IF<~X_|cFu-I0c zKfH|fCY_5fQ<6DD*3HCwAE2_PE96c~L%G@@9b~}dMP8GYPWpUQ8p7nJn#k`BlXJQC z6>5;&6h(4Tq}W#^7gwkyUKj-0w38>(N~1^{UU!A6s;yC^4PSxc^4f8TsnmEr(`bcn zjaHb&uV1Ck@|4XWt?#hUAtU*a9sfTj|`L@fg;&qp874N%7ajvDv{ltaTz=IN2 zuT!P8HUDLS1-Z1~#U{@Xw8d0E+LrM?JSR2~wskVzuN181cp`K&o|N+E*Rj#Oah<{{ zBdE5m1gpSOPy*HfoG#eb8EGlU-9Tf?LvojCuhA$wIPV6Hq@BFy2Hq9#V(U$u3+(1b zH?bMn!>`<=%#iIUk0{u-6YK)J!5&_B6SDNs_kbbm2%o+99;H9RUb!Z8%`CA-uIA zytN^`wVeWRQ8u_J8(fr4hPY?FV+$J!?f8yrp_Q&_)je;RcwIj-!D}Ad4L#=Z(w=Wt-l{RJzSWZX?qjo^czV>u-HxoxDx%ylc>kK-Y#q*M=b1b_@IkZi74E zZ&>6z)1h68zc!~o^aR>53D5z?nT}nR9w;6!^|dK&x9I1_cB?+6?Kb}D4$i9WZ2KFD zgPZ(Koe-yLZyn_8)92ehP?~73Ppe=#QS6vfcC3tcGq3|uJh+1-X+;XLRrrI&;EvLk%7!-LtwdEH%%fy!FB2S&LneEu%g#2Zz` zUWKFY;dv|hTs4kN>}@&c9;I5I!;kVv^X5JZgFm;o=P&P}ojdRvJf#=-s!U$=t$#MI zd!H&q$x zcoDn=I{zOa;<9%GGVK9+Dx?Pcx}`bBiu8?e9ro8?flD)(AYz?+!01+4AsZhFryEL%hfC%FE(+69igT}_y~@9goEsjSmf}AQ+r`# zA7zZZkL0tuAsX zyZpCWJH{IA_?_O4vy66}&4WqQ#$k@O4VOoe8Fs8X_V2(fFdNK4hsJw5RQ*M9$=Rto z(E?{R(UlF+=j8im8J6Moa}2_=PD11$9-A#h7rdY|`-xfsS5Tlt^>kRFCUrphD_r0w z+R{~?=_i`wo1p7{A}#pp)6F^o`EGC>lXwkZ(9AQ5rf_pNO(HY+7V>!lZ0H?rx3&`b zZgV45EDM&K`X?RpNcPoTc9=y)9J!{LMHa%g{$}w4zKGgs7S)39J(c}2UMt!kvehEq z8g>@7KL^f(3*ZvC46YblMGEI>mnod59j63#I8Zwrs2$EzyT6O{9=H!47}pQcqRHMC zrRta<@Nvq4B1Kk0wvoOI3De#UkM`i}1#JP+|VH^M>Le;pX|5Ne31e&eBcn(Tr$HGKZkzAgHi`3O(Q<+5G`Ddr2wD(DR|z5 z)19JPlnmUSc-L;Q$>&`ix%|FUEbf>A%j~kRO!zUWX$2;gZ z{o|c0(E9=uy3}4qod7(sqYLQ9GyKu;12x_K4Bc-Uy5BN%zpd-mHoRqY zYdhNirF$?xUmm*O({&GojyJ(u;B7DnyaNV<_xN~u;SSe!tg5kNb!gj5{hjz0{Y6Jo zfUx5f9Pfu?`1n6Ks_Gmd%riYoLqX9tQzn#w=)q;TPw#(Q1J^9^5?$#h;92wyalj|bxJJH8hv z@@W>23l^C?K1jsYp9b~cfazc+fWvdlLiY~ycCV<5jewjzC(*r2Gz*6tP5fgS3zzB^ zF4HYsu3K2FTew2EaHV14D&4}>x`idWg==*SOSxe%EL^W!csN+ZhOB}`t3e4^3rfNI z&qL5&J2h)Y8P@EQhx%^L4iO26cKd{gy7*R8dm-=G%iBXl20oPZ3l*v92zng|dL6sL zPhcS>L;7SyT;sT_nLfp6pL_gTsEGED!-M;}?bz#Z z*YV^Mjt83gn$7?Gi_MSZ917;C^3V5SqFe}^xdYDJaUZ~$JN^NWcwZRIwP;Ge(v>=` znkM}{w{Y4x%q^-T^2K{Pd_wyVnO)kKdrk|!N<>)e#DeOy^H4WB*R9E!X~+rJhuUg4 zCPM`8bz?9`@?|#~GK&4fk>Km)aQqm<&BL+i#_2=NSy2wP@Mveew^AyO!HA1-YEn`fCsD!U(8_T3B_mN389v^`@=i%iM!h>D;=?D=^ z$t)tpOL$k?IZ|}MSHG(wF&}I2?~$;pCaY1Rq1%m$g@XtX38FwWh~ZbGgeTjB%OsEt zs)JNe1Jp#{{^;#nzl%6Uzyj%)2CL~gK?9Zby z9~yGKXjGz+{E9EmC38PnReC{a6P_Zo@7Rth3PuhcwF zjS&??WsWC;wiYSBLsm1rbivc5hn6YboTtU0p$d6xj7XvuQdUT~eq`sA#h!@7XQ9rP zvIF8QMGDSPoUJ)OR)pg-^q#RI4KHE8jur9vLb@~-zOFr=h!rpP%tZ4w09ezU*#PUB z6YH3Ek9EwM1F)7kb3vYwHbYtfnu9{n!nkjV6lWLCHUL`#XFJdyU9*b2#EX#fa!27( zON_qgCbwCUJR=_Cusc7BLreDHOYw*6JvcC~PF0-9UHH1H24+f!?4GzZ8#2@`k2)qoH|V8O?9<`gmx5OV>O|*ZhvA z8JjcbVEH(N->rmpdskQY2E<~eb7G})Vx@Cpm2+a1a}EaYa&%>=(^om;vmob3`kcf& za(k@vV_6kBO3J>^(*n{rq3#nNSs7hDn%}P?8u5+FLX8>YDJh!ca`Sz(SpNl<(w$wR4(`L7DPGU+xrpa760Z%{YehK0+*2;c~A}RPwsYcto zXtL$h`AVKeq{vp7%BK=VZH})hB9q6Vn&ZKz06PfhX8^kh=Vb6X_yT+hz5-LZOI0+- zT&?ww7_C2#Kdy={pU-ovijnxfCb=30<3cW|Cb9zOqR>1rA1nY1d0I6H(9WRLKMesZ zc}q1>xA_V@Ux~qW&o}S;Dt>YW7W-lHdQGRXHFs{%?jAxwnS*TSMp z+|`5ktDAYSM?{5h(F#2-Tc}8+LiSc(g6Fi2_sZPc`J&9dL*_R7J(Zg_BKJ;?NkZ;j zTt5l9ck^>eB0BshEw|!#qwFl*d|nn*=!Zf=4usuBQv)6DsFX z`7Wgv$EJvv${j{${0feMqbr7_h;%znl2_CbO|U*ZtRw2tu*b{O8Eq8R8jOGCxDqz`wv8lXq4wUl(_=oLv{*FP;0<#Q@LcujE}j z-X`yQ@GY5T09USuWCXXb2MtqrYCX{g|KD<>o@j>uZ>d)w!me=t`pAN1ygr^e_;!6P zQqOZ_CX&9~NG5~1uS`bs#7sz;%!e|eVG)Z4&{x9M8^GM3xOW58=seG9Ao}5eEF=p* zMsZ;ld{-)ul1Txt%YwdF`KG*km#bzY8P8p_q5T?9%7!pAug^x77``R%QaH9Dl6u^{ zA++c7z=mj$cD%G9+6>>SH53K-q${TpO!}F}G=h}#ytWaHxW^9~p`1S_HiocBezCFW zinVxcWBk~QFEvKFH@J2at*u{cB8K2Ii)&5r>+Y_wm@6*c6)s*KOYKH15Ce1AN-Gkfnn=iGMMolIM6TTyHqy^`FC zRFA_FDVAz-gG3rdwRvqK^`~bzF^S$0_0v{=n?xt9R(^Qc`U{z&xK=F+Uj0=qTJNGD zuGNB~cvv=tQaqo^rn+33O_iyMd~QSSd0BIEaaMDxK?C*YiMiz9$<3(>&5%Ddq5{4l ziy}#>9v$FUT2M8*!k@OFWKwvYe2e6Z@~xiEtY%9v$kw0xwFKji`tw)v_n!Lm0r~qt z{n?&_=h6CemVEw%Q(I94|4%+I;*WB`VH>Z`p)9(j*95f!rF-)E8M1QER*>sC{=601 z_2bR5c$=@vVgkEcqnO6cT0^`AJWv)z{JH#bhh!MZr&^Pn&T{oOU>MQ{3~$NjXDCen zX=EF4N@Z??jyb$h7O(QL^nI3^cG z3vMlojy$G4rX16pMss+x{AoR3&83F4n^W6SCWQHKJCM4|(`Dh%#XQ~)V!HWmJ2X$= z*!C#;aGUnvJV?@Og~i&{p6b&)y?Ja0H2+6`enCD5dlx!{zmPw);8h(!q$i(}#c+0X z1hWsgt}Nzqo-9guVn@{N(aG%Vh?zL6KU+HC`6|bBq9oyX$gMh|nT_9)MQxrfi|2W# zEMDQevKYp(&!YH(+dPY#Kl{eZWD3BWX90RY;T9VXK#%UL|!GMTJl zKq`;}yb25jCIOdCCYuis0eFBkpalQ{Y;%D1z-|C^Yd47xRE?4op`)Y9letJQdu(m8n^`SM4F*to6dQ6>Ii3I!_h=vtqD zo(B;Y$m);BQD9a1W4e^ZEXTO8(DZsrArTxvr+W?>jZ%7$SYZ?b4p-5-70YQq?6 z3?r9&6;MB=ExP3L&jsZ6dr7N**QoEx%RZwhl`f-aH-7XTRpy7g$zG*!rfD61>W;?6 zrm1>PgmSZ>Elw= z^y3-dQzB0A`tK=<`t#^J6epDbvRFVK8o)jJK^?VTPStZ>!hpL1-GJ^u9?-)mvBi~M zKyRQA&=+{wDEpzrj#l~u*m%kSUcUg7_6|gG`HVLvZ5V$wm?DKTf;a4^aE@6>X6L(Q z!*M|sVoL{^A?_${xeyj|4EJA1sp@;OiJ4E`p$MUj<>G}DPUHFTLaGXz@JPN*T9D%MVel-ZjYU4wNOCCIAzG4}ro(VTD4+*D!zZZ6-yR%8EeVgEr?IYQYg zr(WtCY&M*}8Xi?{EV>3D^v50k#6$jItEvcHk!f=0e#C>@rGdedT9hFR%~z1yi}OY$|D@u9YbZ0-uI- zE}tKasXT)>kBbWAzQbJcl$8_=+bZ=no#R@s!#@~QNPa>&&tnQHMkv2?{uGK*dV$gf zp1FcD=py^gq*$T+!IM`}9Zp$E&YaT_;S6vVI0u{ueg`f97lA(@-cql49~j~_0@3@L zc*_m(9_WhxA6=k_6$N@^2=rJNsHZ1DHF3>Vluo1zr20sKvSROp-va;!Pk|j%9s!SG zm0@920(|(|Dr{S6om@Vvb;hboGOi-jAdVf0VXMJj*8UVCImc3{cP7-zysZ#3;WikD zv;U72$4f$m%DJm4G$#;jTmT#l6*gB50YZVw0PM3Gj_LUEY15&4G^y4aQl%JDr5aLI zuOL+oU8i&u=Z*lj#z`) z$9S}1*~Cfhc4~W%A1&<$5jg*H~9jn&}*w6Qu8qxqEw6=O6XKEZHO zIm3@;tTXrpfBbKLkY@_2>(bP6M!`_y>LlPJ;9~&lT%F3J*O5E;D~-u%V^m*Da)AqZ z<7^7f1gUT2yHtxV`U^fvu16`e@9Jz>CB2tipCbA)HLsBFQ*I$aug0Nc{*3h$jYy&A zAq=+ZdUC-%*Pe&T=XpGNJ*CimE|x_BUzCl%V`~YD?>VUiKQ7R1HRegFHeRs~zzM8p z@`c=A0bc{(0JDL709HkXRZ-{hv=aDU%S))TU!gXzi^jlK^XNAysxF;?Y(-^*syLiA z__17ZS8!z1A30?M0?lIXx&aZ^TAsRrLTDWqY(Uhrp3iTfq@Z=^EPa#n>8e{TeGEga zmT>Gwiiz4Fd8u*&g^zUdJpLBEz+S^ORX1|pMz|%LdD2D-cm0T_UJLW+8`N*%H5=(Q zcsV9@H2~LNT?70G6a#C4b-;R~EJ3*e*a*NTXzjLOJN;R9p2RJaP4i~Qf%q?A=kmE` z%+G#9TO5#wB~s~%1_ycHCfMIYe18*$f24x?II5`+PTmZSKjj$tD*r+A#Up#c{H?5G8A3Tajz{PeV#wrLd|R!P+v@2yC}Tx1_0?^GhavdLBJcpo4{b;EeIR#6*k-u_B8sA(l|sK z97by##^@a0E9WqlhwlQ1_j$=K$_yR_ilc!sz}M_b2l|h zg%53+0x1&9q!2eCg#{k#07&t<)-BoSHWS@y^Qqlbi)QJfe4*RvPH6ficm5d_Uvca{ z3TuELq=qinTTeQi7EMPRUabtF-iBL?69~7~lFui9#@YOql|67`=jbkJ!#!|GWv2=* z=?FANYGn^4rp@~|FH5;|YnhKq%^n!xT1}6p5+$z!F4;o~L9#^!!kba^o$k&C?4=s7 zeU3Sr3Cse%0KNpE>@86Amf1i)@GUUMc!$Eb%md~F1xEQD=D2#<9GmtVN?jmIEsQ1R)j#Am+MdtY{3D z6DO6@(V;|RRo`H>L0X>Z1Q! zKz%{R0}HKhE!Q8QezGZ+vP$C09OyuF?_k+Egut zVL8RI2axVN&Hh5f@jC~oa$D$e%OT(}a1{6zI1c;X1Q8ET{3)2^M-GhYaDoxs=H)^a%EQxLf5&|L72B2JmMf_HogL;e*#y5 zYesn;;*`hj6S+-0%=`(IkDI1~sK2 z(C0?Vfb@2;dEsG74zZze`N&JOWwAQ=kHbhOI9WYH4Um&cKMFfsi3cB{ zDj|+1I&>+}#P#QeM-aLO@HRA|K)!SY-Ce94#o`8O-IIbfsT`0hT`|r+;VHoC{{+Dl zrlauD<{zbKPCJV1`}4i=>i;XfVz)fRj6VV%^N?SuMkca|RtMk&;QCoB8Snlm1Asul z1t3(iB2==<-w`fZ5iVIP^VMIWk>kPJWj0pHnjlvxlB(*9l&CLKlAP{n1WhSq^N>Q-l)`mRs1XWQ|@Ysk?9>?N0=FF2w@i*a^-{2-Tcw#o2 zp;=Sza1M*WCn<=U@uHIyrZz_-WRmii;P-{U(KD&(5H|IA z!RC|nALRYCk@eO`CKHd|o%-k$n)RtLvcCGr`c)X&sgqQ%ZWl=WBJdK>73co<(Z&rMg|AcXf-n8N$~`((vl#$ zH@i`H=e2sT{Ewlm%6gaw%G1yJp(2RQrzt1%EyytB{{flb@NT84a1NmI$ObB7H8BSp zVvgevPg7*CasL*xEMaBUYNTWkbA=|Rj6WcwA*QTA0nhX`-jo_nnicG_;<+hm&`FnOS5dfr<-N#eQ&trJcI81(;FjIA5Z2r14V++s|9X6Mjw;vz7Xk18GY$1{rm|wz_h4QPHu&R}L z(k0xDxOs_u3*&v4=%oZdTnIsv+Mr==(5^OUR~xjd4O-QPB$Uk!gmH&Iu+VXu95Mrq zaMl(Ne<5Ai5)3h_YF8q{mZ&G+Y)M@AGFBv6AKQ|hWKHO%#r09r5w^tv@c>k<4XV}# zRcnK)wIy@$Ww5HPu`*K4wmN+FGS!6f2)Y6lQ=c1Oq4c`7(Yg*$52%lY{Gn_iP2Kg= zWNHE_n`wlW7&2$`(ks*`wHdl(qf4Q;OG`YN;d4L_+qBe^xyb0!o+JLG2(*HI_F=ZUt9I zwnDCPjbdF(QCqGJYN4WtdtamYy8ljySkOvyU2HN@q2)C6T)A}gp%vli+lr;5A0MQ( zC?~oFHEa3yHOjRu#!#04OS#c?3U@EXBeEK{HNcNRu~Dw&QP;8NGUx6x?=fa$8y8%s zaa77FH*g8PonO9z6vI!Pe*?+#o&3`cY8bo~{kH+7z;*!HZ`)3exCvSZHChJ@T8FsV zO(-|`wws`KX5S<;5*+j7Pe=8o`qfyfH2z@r^}j{ zzZ+s+(d1LAqT)|Zy93Ybnpf4>u1i%D4+nPM9UQBhd`#Bg;zxI&k#6gUEb}h8b1s7+ z+*=#mTN_+k+jZaua1(%XvE7C+-V_n-G*y!Op5FHa*hv${NDCvUho>M80x~y?qgcJ^6vZa zAG@;iF}{|tJU~9TD~CRyDoteOS-YG?hR-m!oi$h6E^Erq!Xj`p+DdclJyiDYJp2Lf zdh>Wc$ho>d(P=I^^oqy&|1L{nR=ywrKn0`}J+flD!RFZO|QczPbVHPV<|?IFb|uYtfoZubzgJWOjH zZmhy^W2cNTR&OM)dk8f(N@lR5@y!f9!iIgX;`SM9Y@czmE}E}jLWr^W5!r*^myBh~ zv^)s5PvG5;CTz!Pam?#$`g+Ao8$CR8gNw(KkAYDZeh^1Ih>1u%qMhjV4`$SY% zNI`s~FVV-clXmawRvM$#WGO`)pLmRjXR23Br!?^PX?&`Rh(f{v7xeZ|yyE&!A1FfU zv+3i-lu^SlFT;Tmz(`;eFdBFdfCFV82fPnV06qXF0?5$YCmGKlqnr#(1*QR?0JwBj zqPz?B8^R2^PzxY~OaQr*7-Q;@yy3AcaSLV7HzxmEK0=~CGCYXD>!iWf*k=Rzz_-90 z(2Xmji<>|fB3zT=1YPa!Pf0N7p5b?ekgk@z&{kTZ{v59rqB;C5A0N>J7khntL_PmY z5T;o1bXl#WbU~X-yue31OMmcVAJGJt3u|{OBW9WWgcZt%?`%tyPjmR13TN zMDxeEj#B04@TT4E%u-n$#{!Xj1!C0NT_JZEA-m z)!y%*ybJsV+%ulSzOmC z+PNA_0$%y6qrM5VQ&fwT{+u^lPphJtRHV!Fqb}20o*e|4*6A`Wh5$=|rCc>w#7SS)u@+dz-GeblrJ6kB4SBZ9 zY`&Yn2^NWPUrU2U2EOXlK8tYd;;JDc9iLS83=uUm;l(=O#X5chb^yD8-5BLW@2Rf` zzY(Y2bREaxu*PJP!Q=@49D?aS$~QtpEMjl%LyQ;E?6`i79EwK--Iaw4-faI+G(F89 zhl;G$hrtBytK%qe3^)$JeRZ4!;I`sxpt4Wy+8WH^EjTo{RWSlmhs@DR_o`7bC02Jp zYdWr==cnGD$BmxwYaDkpv8Nki-<4uVDY2;dONt)E116JIIPU2KkW&AxOZ`BeiXi?{ z{_%AS-(Vmo!Kd<+~KDUzr< zXGMxCm1=5J&_hcO`CEC~C?o_XM2aYD2!Of5_ z=<>1!^|^@17qjG37sAL?r%PFBOt#EcRt#aiS#GDa6BR7N$N0KWH+{~A+&v11q7l!E z!Z5Ssw|*g>1h}1g0!!;zZYuhLXA?=Q?8_vDJFKB(m2h$8Y|1sFMRGl9BPvQEEsYk- zICoCythT!4rmV`uyA$+ zp2ZrkDO+P&sb@FnEyyjy6PH~nHWoRLW2#`mdT`e&A{m#xE#u*E_U6&?P(6L5zZ~oA zD~G7r7nj1PKgU`4GhX=7%lt52^vdoIZh1ft0AY!|Hat~|K~f53lOe<# zx)5)2R6K+jtPAm$F2oQ`2;a9W2(g7Xp2Kn1q!_A65jhClU^tyHoK6@{Cyb^OM$`E= zFqB6m;QZ@G)A+#1IZi+SxNx_}INz7VmBdoG;xjx5S};@o;W80 z9|8zZoF4&`fhhpO6X!Glw%9oxm;o>_lh-F>s^)>OOF3$+&V2qo87i=V)oNlKKJ}SY z4g2_eUQkVB1Irs3zRuS{_TeX%Nt}SS_KlT*__se+6IGWqX9X=}L51 zg)Q=T*%tBX4>tv>TkOa9)@ec**#YNz?d>5bR!DR^OE}*nB54Dcc(4K+`3T;^Hfc>C zS8Q4!YJ%cs_DezQEu5MnD*tCInrTWb*YdA5j6mCMl8EZ_1QFVVic)?TL^^KQ>hbRa zb~m`+1k7dqATb@ML`k z>SfYSpI(3GQ|KPj`njI$hq&B%gd0>B^=#{~H`W6s{CahfM@Lp4t1i;)a2SHqL^iDH zOKBn#*J?B6TYEm4CYs@&d*ai@f8ecllSKi~NEZ$9zc&ZdMQ!{Oku^iq!aotEXW&iY z*W_Cyf0F?fUypZX;CyEDJ^9v=6B^*HCwFLoCIfk_EJpLH2H@}s-)|r~;Xj0OGev9s zhtRA{5WB>?GtuN8Kh8v*l{+vPk2FO~Abc|I!2#(vbt2V$SAp=cb|!K7yOu6e8{BOHD!L zEQdFP5O=vtGxT%t_+}vH<~7YkPk1ij+4!*!x6DSrK|DQMo9mym#k=?>uuXH+UE+6| zqr*SEs5v+U^G*3<3yyAqq9?a$0S?1?WDD^7fH%uxE?ZloLkVZJ#DwhOezG{rpSBbU zNGPmnDQd}o_PNp$`2ZWY$U!=yHh-BTol3!m|S)M diff --git a/meteoinfo-lab/pylib/mipylib/numeric/core/numeric.py b/meteoinfo-lab/pylib/mipylib/numeric/core/numeric.py index 8e3147fd..b0b2e533 100644 --- a/meteoinfo-lab/pylib/mipylib/numeric/core/numeric.py +++ b/meteoinfo-lab/pylib/mipylib/numeric/core/numeric.py @@ -1292,7 +1292,10 @@ def any(x, axis=None): :returns: (*array_like*) Any result """ - if isinstance(x, list): + if isinstance(x, bool): + return x + + if isinstance(x, (list, tuple)): x = array(x) return x.any(axis) diff --git a/meteoinfo-lab/pylib/mipylib/numeric/lib/_util$py.class b/meteoinfo-lab/pylib/mipylib/numeric/lib/_util$py.class index 27b25ca90600b9126c7b6f551204b682049881b9..33295bc3ea52254e64613128153b2f8a7a256400 100644 GIT binary patch delta 3111 zcmcgtYiwLc6+Y+e!xcPZ;~dCN!|VQD z+q>I+v`O>Wh(Gj48xZ#&Dp7@0AyrkI#P!DDJ}3bs5G4did6ZIGC~bhYyrFrSGjrDu zq#6mS!diLG%zWp3-#KSy?l+I)=Z~ZRd)L1HJOH%Qu~YtQs}1*9x=#@Gf^hdWoU&TZ@sth9N zMvlrfH`$GA`vb4|VUXfFwNr<*Hrm0kQF9B|b5}OO&KlipKeqNRv4y(*w3)Wr*Vgum zkJC0j19XeMW1WfH?I+eHMIY_(($=YAQ|qHI+%_%punwN)iQBfUJB(r${DS~ zaXS{?=`J=jykQ>;e;Mz$Z)mw557_aRPuIK(NErx#fFizTpJ*us`8x3pT*~K(Pap3d znrLc0f$o>5@E7)XA|X6%za8noL-va3@%kf-Wa(ESra>?NEItUb#*KG*&d9kszF{fH)1|{%L3(iaY~1ww=avcm_VxA2INJ6;;8XoaFe(AvSK3 z)lF=Aghjl!DHMMIG&lz^!_8uV2IuK88-6uC!PaAd^4Q2ua#fD=j+@6F6Weco2cuog z{{r(f)r%Fdq7hdGwi(Dv_eLxDLA>#wI z*r&ReG3LA77a23Lwdn)9(7>_Y-O8uKuIowHOZ|Lm`2Vx4olmmVZhxWYI(*UoYR{VT z{QkEXHH5D?xw%SE11ohMQx`F|i!nY>yqT8!> zT?YPPAnG>ooq?#wz@r1vNZmf6;BEy+35QR)Jym!a-&6QMOXOz*0R{CefF2Wpn%U5Q zRnWhy(El8WdJTMkAgUR7c_1q4as$x-uUzGUmu7en1|bhVgIqlL4f61?#2}3a-5?)R z8fEka!iv=EUO*!=Q&LfzlS)e@J=I7z`bc1#^3rgj8+(ORf)^I^+OE9BoL3n+saZGl z1iS%XRfbf{bU~Rd88Ta{m)X)X&eovioCGg@lm}j0$pbI0lG-p9*riwk-uXJ>diShQ z3Mo`ypvISE1*GazQ{D-ug`X#_V z*SpVq$H@xadD4#D+F3^Do#fXTzXE|7xI?|h{4UCGi~OimKQl8dLcY@<@=p2xMYIS{ u`9sfW5u!zS_WjRA1)`(eP@)o*d0YuD%JQ?nzmIZspEPz*o{G{aQ0d=NY=6xF delta 3158 zcmb_ddvKgp6+h>E-|pVsehJ@Z`+k_Dh19%kFnBqat*TL_h%fVC~(YKzT&y}OSWYOoWbc# zviID3f9IUv@0@$@_tk^={6P%-@bWXy0ze0S{)xcln$_1?$y6e5jij>^@qFIS3j8jD zHe{;y^zl?KFYw0(z6^KzVHkgc;F|{s{3*hs zQ9D0kji(dw@tlD_gHHuO$DbQG!(sg1j6rE=V9~=|WymK)FTR zt&VnfGdP~rw1=%q1~3OfCTHF-YERmGGTA;sRR~7n9y>8+jie@7prFMJDCRQpY|ggQ zV}h1)KhnGJwG;WAnlERO6Y<=bnylo_c+N_3xH#8>Rx{6BB9+-I=wqC5A~!0ih8^;0 zRRBS2nNiNRlY;8lnv73y#s=2FLld`>8BAEY^j3K=1hp`g^h7FejoDl~iXc>~Y1U2_ zCX?|=PWW-jE}OKo^Ozersth8igQIeqE(hy}%+moFrW+0Frl_=f9jm??`<-I_o#G}M z2+&6Qr1N_Hptza%^wiNNr@6txtItg=b-6-@8hAX_B}p&(8HCTa^7iAGmnvum+&E{tFva=d}tr-2lT=B-UIDCtYF+& zXEoy<={&`_U0pTN55Y0zD$x0N)BFBWHc!>*sjl6`A(lmJxJc`wJ6P=3q8qVJid_Tx z9BgKR3o<>D!(ikB+;IU0n8ds5+L%Ps`WZawoLRqm)hl=$b`FI@_-bnj7UNL3gugh( zfsnDwdt`w(_I!o2o$RTN?gJmZ3^2z$7~q3bG{r_hO%Jnm1W=wFNs}HyJb4U-yBnS;?_x`}xQ?r~6Y8 zsfdpf|9+d*@bPxlz)jo|*_39M@uDl^lvCZi)cMlJisBWg7uj_?o_3wI0#twopru-l zk5t%4xiDvgd9u@74csR|8I#rV44c*xe7WZ%!-B3|54EQROkz_h&P0nVi8{m=VB32 zm5D`6Z-|fwUx+jwM2I{*_(SC5ArK-z4~s&yhzBD?6`Ydu`^_em)GM1ODIuqk$`n=f zm8sHoGCc*Q8->gP<>kYx-Po&E34Ha#yf!H>QT8gzYzuDaF}6S>fv`ng{ z6;ds&VznV!!yhulp5a)@cXM>pT^5>^eHRT5g2gf=CiT?IKA zgrjbdfKN&<$`&w~wKJHtGnj=k7<~&(`Oy_Iohry|5FU4fEGm~)x|G}aSh*qiEfrh& zQCkq6cKwV6zoN=pZpl+uQhh?EBECeku1f_joMOq54?VLkCp@v>C+jBPs35%Hx>a6P zt5=@jdnsQm{-n;|)_YC9hw>8ZPvSZGoB{Z09j%|=uH{=PcbCk<@b(b(9H$K>>YJg# zmKnOmjpnD#v}GP%rNI0yR1MG2Ha8s7r|Ay4?;IMu+o78~iY1ynR-()?be~Utjq#gc z&cQbI8uN20zdG`BGyegaV-oV6zD7I6Ya8w(M=R34Ypm#_Q5aG6{eLI`(E&aTM0@FB b9@l_IIezBsIq_KtWr#;fh(cb?8df$Tm diff --git a/meteoinfo-lab/pylib/mipylib/numeric/lib/_util.py b/meteoinfo-lab/pylib/mipylib/numeric/lib/_util.py index 452cb7e2..477da26f 100644 --- a/meteoinfo-lab/pylib/mipylib/numeric/lib/_util.py +++ b/meteoinfo-lab/pylib/mipylib/numeric/lib/_util.py @@ -188,7 +188,7 @@ class _RichResult(dict): return sorted(omit_redundant(d.items()), key=key) if self.keys(): - return _dict_formatter(self, sorter=item_sorter) + return list(self.keys()).__repr__() else: return self.__class__.__name__ + "()" diff --git a/meteoinfo-lab/pylib/mipylib/numeric/special/__init__.py b/meteoinfo-lab/pylib/mipylib/numeric/special/__init__.py index 2d48dd54..8c64a629 100644 --- a/meteoinfo-lab/pylib/mipylib/numeric/special/__init__.py +++ b/meteoinfo-lab/pylib/mipylib/numeric/special/__init__.py @@ -2,9 +2,11 @@ from ._gamma import * from ._basic import * from ._erf import * from ._airy import * +from ._lambertw import * __all__ = [] __all__.extend(_basic.__all__) __all__.extend(_gamma.__all__) __all__.extend(_erf.__all__) __all__.extend(_airy.__all__) +__all__ += ['lambertw'] diff --git a/meteoinfo-lab/pylib/mipylib/numeric/special/_lambertw.py b/meteoinfo-lab/pylib/mipylib/numeric/special/_lambertw.py new file mode 100644 index 00000000..d79b01de --- /dev/null +++ b/meteoinfo-lab/pylib/mipylib/numeric/special/_lambertw.py @@ -0,0 +1,155 @@ +from ..core import numeric as np + +from org.meteoinfo.math.special import LambertW +from org.meteoinfo.jython import JythonUtil + + +def lambertw(z, k=0, tol=1e-8): + r""" + lambertw(z, k=0, tol=1e-8) + + Lambert W function. + + The Lambert W function `W(z)` is defined as the inverse function + of ``w * exp(w)``. In other words, the value of ``W(z)`` is + such that ``z = W(z) * exp(W(z))`` for any complex number + ``z``. + + The Lambert W function is a multivalued function with infinitely + many branches. Each branch gives a separate solution of the + equation ``z = w exp(w)``. Here, the branches are indexed by the + integer `k`. + + Parameters + ---------- + z : array_like + Input argument. + k : int, optional + Branch index. + tol : float, optional + Evaluation tolerance. + + Returns + ------- + w : array + `w` will have the same shape as `z`. + + See Also + -------- + wrightomega : the Wright Omega function + + Notes + ----- + All branches are supported by `lambertw`: + + * ``lambertw(z)`` gives the principal solution (branch 0) + * ``lambertw(z, k)`` gives the solution on branch `k` + + The Lambert W function has two partially real branches: the + principal branch (`k = 0`) is real for real ``z > -1/e``, and the + ``k = -1`` branch is real for ``-1/e < z < 0``. All branches except + ``k = 0`` have a logarithmic singularity at ``z = 0``. + + **Possible issues** + + The evaluation can become inaccurate very close to the branch point + at ``-1/e``. In some corner cases, `lambertw` might currently + fail to converge, or can end up on the wrong branch. + + **Algorithm** + + Halley's iteration is used to invert ``w * exp(w)``, using a first-order + asymptotic approximation (O(log(w)) or `O(w)`) as the initial estimate. + + The definition, implementation and choice of branches is based on [2]_. + + References + ---------- + .. [1] https://en.wikipedia.org/wiki/Lambert_W_function + .. [2] Corless et al, "On the Lambert W function", Adv. Comp. Math. 5 + (1996) 329-359. + https://cs.uwaterloo.ca/research/tr/1993/03/W.pdf + + Examples + -------- + The Lambert W function is the inverse of ``w exp(w)``: + + >>> from mipylib.numeric.special import lambertw + >>> w = lambertw(1) + >>> w + (0.56714329040978384+0j) + >>> w * np.exp(w) + (1.0+0j) + + Any branch gives a valid inverse: + + >>> w = lambertw(1, k=3) + >>> w + (-2.8535817554090377+17.113535539412148j) + >>> w*np.exp(w) + (1.0000000000000002+1.609823385706477e-15j) + + **Applications to equation-solving** + + The Lambert W function may be used to solve various kinds of + equations. We give two examples here. + + First, the function can be used to solve implicit equations of the + form + + :math:`x = a + b e^{c x}` + + for :math:`x`. We assume :math:`c` is not zero. After a little + algebra, the equation may be written + + :math:`z e^z = -b c e^{a c}` + + where :math:`z = c (a - x)`. :math:`z` may then be expressed using + the Lambert W function + + :math:`z = W(-b c e^{a c})` + + giving + + :math:`x = a - W(-b c e^{a c})/c` + + For example, + + >>> a = 3 + >>> b = 2 + >>> c = -0.5 + + The solution to :math:`x = a + b e^{c x}` is: + + >>> x = a - lambertw(-b*c*np.exp(a*c))/c + >>> x + (3.3707498368978794+0j) + + Verify that it solves the equation: + + >>> a + b*np.exp(c*x) + (3.37074983689788+0j) + + The Lambert W function may also be used find the value of the infinite + power tower :math:`z^{z^{z^{\ldots}}}`: + + >>> def tower(z, n): + ... if n == 0: + ... return z + ... return z ** tower(z, n-1) + ... + >>> tower(0.5, 100) + 0.641185744504986 + >>> -lambertw(-np.log(0.5)) / np.log(0.5) + (0.64118574450498589+0j) + """ + LambertW.setEPS(tol) + if isinstance(z, (int, float)): + z = complex(z) + + if isinstance(z, complex): + r = LambertW.Wk(JythonUtil.toComplex(z), k) + return JythonUtil.toComplex(r) + else: + r = LambertW.Wk(z, k) + return np.array(r) diff --git a/meteoinfo-math/src/main/java/org/meteoinfo/math/special/LambertW.java b/meteoinfo-math/src/main/java/org/meteoinfo/math/special/LambertW.java new file mode 100644 index 00000000..252c8632 --- /dev/null +++ b/meteoinfo-math/src/main/java/org/meteoinfo/math/special/LambertW.java @@ -0,0 +1,103 @@ +package org.meteoinfo.math.special; + +import org.meteoinfo.ndarray.Array; +import org.meteoinfo.ndarray.Complex; +import org.meteoinfo.ndarray.DataType; +import org.meteoinfo.ndarray.IndexIterator; + +/** + * Implementation of an algorithm for the Lambert W + */ +public class LambertW { + + public static int MAXIT = 15; + public static double EPS = 1e-15; + + public static void setEPS(double eps) { + EPS = eps; + } + + /** main branch W₀(z) */ + public static Complex W0(Complex z) { return eval(z, 0); } + + /** main branch W₀(z) */ + public static Array W0(Array z) { + Array r = Array.factory(DataType.COMPLEX, z.getShape()); + IndexIterator iterR = r.getIndexIterator(); + IndexIterator iterZ = z.getIndexIterator(); + while (iterR.hasNext()) { + iterR.setComplexNext(eval(iterZ.getComplexNext(), 0)); + } + + return r; + } + + /** other branch Wₖ(z) */ + public static Complex Wk(Complex z, int k) { return eval(z, k); } + + /** other branch Wₖ(z) */ + public static Array Wk(Array z, int k) { + Array r = Array.factory(DataType.COMPLEX, z.getShape()); + IndexIterator iterR = r.getIndexIterator(); + IndexIterator iterZ = z.getIndexIterator(); + while (iterR.hasNext()) { + iterR.setComplexNext(eval(iterZ.getComplexNext(), k)); + } + + return r; + } + + /* ---------- core iteration ---------- */ + private static Complex eval(Complex z, int k) { + if (z.real() == Double.NEGATIVE_INFINITY || z.imag() == Double.NEGATIVE_INFINITY + || Double.isNaN(z.real()) || Double.isNaN(z.imag())) + return new Complex(Double.NaN, Double.NaN); + + /* deal special points */ + if (k == 0 && z.real() == -1.0 / Math.E && z.imag() == 0.0) + return new Complex(-1, 0); + if (k == 0 && z.real() == 0 && z.imag() == 0) + return new Complex(0, 0); + + Complex w = initialGuess(z, k); + + for (int iter = 0; iter < MAXIT; iter++) { + Complex e = w.exp(); + Complex f = w.multiply(e).subtract(z); + Complex df = e.multiply(w.add(new Complex(1, 0))); + Complex ddf = e.multiply(w.add(new Complex(2, 0))); + Complex delta = f.multiply(df).divide(df.multiply(df).subtract(f.multiply(ddf).multiply(0.5))); + w = w.subtract(delta); + if (delta.abs() < EPS * w.abs()) + break; + } + return w; + } + + /* ---------- initial guess ---------- */ + private static Complex initialGuess(Complex z, int k) { + if (k == 0) return initialGuess0(z); + else return initialGuessK(z, k); + } + + /* main branch W₀(z) initial value */ + private static Complex initialGuess0(Complex z) { + double x = z.real(), y = z.imag(); + if (Math.abs(y) < 1e-15 && Math.abs(x) <= 0.5) { + /* from SciPy cephes/lambertw.c Pade approximation */ + if (x == 0) return new Complex(0, 0); + double L1 = Math.log(1 + x); + // first order correction: w ≈ L1 / (1 + L1) + return new Complex(L1 / (1 + L1), 0); + } + // approximate: w ≈ log(z) (k=0) + return z.log(); + } + + /* other branch Wₖ(z), k≠0 initial value */ + private static Complex initialGuessK(Complex z, int k) { + Complex logZ = z.log(); + Complex twoPiK = new Complex(0, 2 * Math.PI * k); + return logZ.add(twoPiK); + } +} diff --git a/meteoinfo-ndarray/src/main/java/org/meteoinfo/ndarray/math/ArrayUtil.java b/meteoinfo-ndarray/src/main/java/org/meteoinfo/ndarray/math/ArrayUtil.java index 4adc83f3..265e8b57 100644 --- a/meteoinfo-ndarray/src/main/java/org/meteoinfo/ndarray/math/ArrayUtil.java +++ b/meteoinfo-ndarray/src/main/java/org/meteoinfo/ndarray/math/ArrayUtil.java @@ -1761,66 +1761,6 @@ public class ArrayUtil { } } - /** - * Return the indices of the elements that are non-zero. - * - * @param a Input array - * @return Indices - */ - public static List nonzero(Array a) { - List> r = new ArrayList<>(); - int ndim = a.getRank(); - for (int i = 0; i < ndim; i++) { - r.add(new ArrayList()); - } - - IndexIterator iterA = a.getIndexIterator(); - int[] counter; - double v; - while (iterA.hasNext()) { - v = iterA.getDoubleNext(); - if (!Double.isNaN(v) && v != 0) { - counter = iterA.getCurrentCounter(); - for (int j = 0; j < ndim; j++) { - r.get(j).add(counter[j]); - } - } - } - - if (r.get(0).isEmpty()) { - return null; - } - - List ra = new ArrayList<>(); - for (int i = 0; i < ndim; i++) { - ra.add(ArrayUtil.array(r.get(i), null)); - } - return ra; - } - - /** - * Return the flat indices of the elements that are non-zero. - * - * @param a Input array - * @return Flat indices - */ - public static Array flatNonZero(Array a) { - List r = new ArrayList<>(); - IndexIterator iterA = a.getIndexIterator(); - int[] counter; - double v; - int i = 0; - while (iterA.hasNext()) { - v = iterA.getDoubleNext(); - if (!Double.isNaN(v) && v != 0) { - r.add(i); - } - i += 1; - } - - return ArrayUtil.array(r, DataType.INT); - } - // // @@ -3021,6 +2961,95 @@ public class ArrayUtil { return r; } + /** + * Return the indices of the elements that are non-zero. + * + * @param a Input array + * @return Indices + */ + public static List nonzero(Array a) { + List> r = new ArrayList<>(); + int ndim = a.getRank(); + for (int i = 0; i < ndim; i++) { + r.add(new ArrayList()); + } + + IndexIterator iterA = a.getIndexIterator(); + int[] counter; + double v; + while (iterA.hasNext()) { + v = iterA.getDoubleNext(); + if (!Double.isNaN(v) && v != 0) { + counter = iterA.getCurrentCounter(); + for (int j = 0; j < ndim; j++) { + r.get(j).add(counter[j]); + } + } + } + + if (r.get(0).isEmpty()) { + return null; + } + + List ra = new ArrayList<>(); + for (int i = 0; i < ndim; i++) { + ra.add(ArrayUtil.array(r.get(i), null)); + } + return ra; + } + + /** + * Return the flat indices of the elements that are non-zero. + * + * @param a Input array + * @return Flat indices + */ + public static Array flatNonZero(Array a) { + List r = new ArrayList<>(); + IndexIterator iterA = a.getIndexIterator(); + int[] counter; + double v; + int i = 0; + while (iterA.hasNext()) { + v = iterA.getDoubleNext(); + if (!Double.isNaN(v) && v != 0) { + r.add(i); + } + i += 1; + } + + return ArrayUtil.array(r, DataType.INT); + } + + /** + * Return elements chosen from x or y depending on condition. + * + * @param a Input condition array + * @param x X array for true condition + * @param y Y array for false condition + * @return An array with elements from x where condition is True, and elements from y elsewhere. + */ + public static Array where(Array a, Array x, Array y) { + DataType dataType = ArrayMath.commonType(x.getDataType(), y.getDataType()); + Array r = Array.factory(dataType, a.getShape()); + IndexIterator iterA = a.getIndexIterator(); + IndexIterator iterX = x.getIndexIterator(); + IndexIterator iterY = y.getIndexIterator(); + IndexIterator iterR = r.getIndexIterator(); + while (iterR.hasNext()) { + double v = iterA.getDoubleNext(); + if (!Double.isNaN(v) && v != 0) { + iterR.setObjectNext(iterX.getObjectNext()); + iterY.next(); + } else { + iterR.setObjectNext(iterY.getObjectNext()); + iterX.next(); + } + } + + return r; + } + // // /**