From e1a62a0b5abdbfc91df7a6d84b069d734115eee3 Mon Sep 17 00:00:00 2001 From: Nicolas Bertin Date: Fri, 11 Sep 2026 10:13:34 -0700 Subject: [PATCH] Fix bug, extend format support, and optimize performance --- README.md | 4 +- examples/ReaderExample_01.ovito | Bin 43960 -> 17261 bytes examples/ReaderExample_02.py | 63 +-- pyproject.toml | 6 +- src/OpenDiSFileReader/__init__.py | 659 ++++++++++++++++++------------ 5 files changed, 415 insertions(+), 317 deletions(-) diff --git a/README.md b/README.md index c3c02d3..b94245a 100644 --- a/README.md +++ b/README.md @@ -6,7 +6,7 @@ An example for this type of data file can be found [here](https://github.com/Ope ![Example of an OpenDIS file imported into OVITO Pro](examples/OpenDISFileReader.png) ## Description -This file reader imports the "nodes" in the OpenDIS file as particles and the "dislocation segments" or "arms" as [lines](https://docs.ovito.org/reference/pipelines/data_objects/lines.html) into OVITO. +This file reader imports the "nodes" in the OpenDIS file as particles and the "dislocation segments" as [lines](https://docs.ovito.org/reference/pipelines/data_objects/lines.html) into OVITO. For more information and an example see this [discussion](https://github.com/OpenDiS/OpenDiS/issues/3). @@ -25,7 +25,7 @@ For more information and an example see this [discussion](https://github.com/Ope ``` ## Technical information / dependencies -- Tested on OVITO version 3.10.6 +- Tested on OVITO version 3.16.0 ## Contact For questions or support regarding this file reader, please contact: diff --git a/examples/ReaderExample_01.ovito b/examples/ReaderExample_01.ovito index 01c9b8787d6cb8192680f35dd8cd3a42ff0c1dc3..bc13e41047506767aa6d63d154bbed0967e65818 100644 GIT binary patch literal 17261 zcmd5@dypLEaUY%TPP_MhNjgb~bP@u@OD7=^LI&=nJ4q<$ap+D0fo(CjJ9j&1cV|5_ zE9r-lR!vTGQd9TsCB|RXieij zY)x9e<)L&rum?!!uy$Bmt?kxMT({wQnf$jlJrvPaGlc+6K+?RkKIvAV0-p!ebN4V}k*z|go^D--l8&@!~uN>C$GvL+I1x&!yM33^>bi%N$K zO;ue4)XK50ery0dh4mn)F?OO;Aa4V3P_^242X7e2dX2CV1v21GAlM++TH#_4EaQ@G zHuTN7Uy`6N1+DtwQr!T&5lA_O%Rq2)gr=0iZ9zHv;Xp12awTbZ9Pn=76+rTnu+%o$ zPO?p!*Dm1s+Ds}*$HjixhK@n&a=r>zf>8cd3HoY7(|Wfi%11!6ek**{vnrsr0c9z+ zD{&17)f>n~umdG`t3DTTO~FsfM>E#~Gh=St!*!s37X9u-KWZ1+`SlGy;Fc}|)u z$q5sZHcTs^XA!6rj+uArfj}R@Uh$~cwDzOlC&_&-S%J??DbgC}mP7X^Fl=6dlw}t} zWX>vBw&lVl1M8GYW45E*fm9xx5Eu=Uj!qbf76e=D(KBx7?cO~(9~Qjw(V2UkTzJfx zb$q9sqkmg13(dDmXD9H$7=I}ixg!2Qfi9=4VR(RR6;U#b1$NL?9#lXJU3TWy38=BL`l1=Sif{bM$<%@0UF3?;Ss-MpIZXh+dm5fkvl5W0P9qdiz-@;H z&y;XAa1h8Ga4nEIABW8xZ#{4Wa2U7=$h4vaRBe?>n`i@4z{52Qq$3J6R$wQqV4t2< zg*9@=+K;=9hQ+Q;EjIH#U>l#`fFFhUSAo)Y!$(G>!daOA98?>ETg56?!KH!cM7<&C zGI(3%xr93P?y>mUM~zvuF@p9YcOhNE(CEPbo4QpD;z>@3Q|5D0K>BN1Qr6)$mCaZV~QL3uJS`)PqIYS~JsusFyJhu7?ZQk(wn(~S2TY#E9c?t+^BA?5t{ zep3CxV)V&D_FOA}9l$c!C~(>;_`+pkG?3?3iw24-f zH6gXBfw8>hZ=i^lw*|O4#cC?BVrCL_GR<(`+JzxYq=L_=c+Tmuel z4EXc#t!vTmwhNAPSa!oX7#5>RRrW$uGQhFZ4AaN{)6{6RKXQ}zWkf~m|5~5Fk9l!^ zgg*e1?gsv$p}BwF1{8Zqj|%Ehi?XrLYQ4{apkK4`O~CaML(-?XRo>q`3n}UH^Gf`Y z5`V13b0F9v?fp-H$$0h!KEL!`TqAfMH_BcE{32HU^KJq%t>Qfk zY&pa?KOipjnq?60Xqi$m>)VYbtUX!Qmr=K_#}`n*q0vg2y_3EU_-BURiMs;rgZv~c zwd~JPNmTk5BuM=wkgbvaE8a=)$ft!bgU*a?#x;UhaMSKzP0(LS(2-=69jj}wY)eAs z@&v7RaRY@6F#c;l+y|8Aen7e%NQcz&5ztY}aR@q}g$uvq4@Rgdkgy#ss=MolqqmZ) zs}wzOp%;BHQn|Nk?&QmAx6!YzJ##;W zx4TX^Yen4OiHEL5y;Vr#qeOw94x$W(3Bx<~*U>izcOCMBmq-8pi%nf?%TJ6x`K|VW zZ~oKrmtLOBteyVyQ={Mb*p9UydG6`aH&3p-BKQ2<=xdkWeEp->4!!i9Xa4&~gKy4^ zzIO8^Klu90D`Rhbd-P1p7YeT|d;Y~I7oNCy?At2_ZtVTfy%?FS+2}t%a^0SX@A@I0 z8JvHXA9xeK#0Wj`59ptfx)b@`j9rYr)o@F?hozSoZWq4XP>U7RH{``@!u%FcpL$u$X6V)#&?~2Ng+2uKBUwuD43D#7I<9@O|GkzF| zRGpr$)OFVvoOUbw&ba|i=wQj7lguX`XIo&Oai+cfm`JkGfjOt_KR^hXo7fl&g`wRr)6iZx#4Hxf@i9quAyt?;kdyIX z>}iaaAO|~w@_N`%R>j%F5icJzmuQvzp^X=o6qz+>4K*W{hp7s2aZIl!Xw9mRuF@JRS1e!*Rae3Pu26NaB zefM0lIKv?r2NYS(xSK!TIO@;1unPQwgGOHCdA0K7&{;1VI7MBy=3A)0HkJZmhOmE* zHL)|uQf&s8%IwIySb2M1BPNO>MyD(Yb$=ugmPG$Gdx#u}g~F(~6^Hu>+u$%d45O#$ zd8h4g+Vv0U$H`APWskao0$80up;dK4*sB=T>qXik zr5yK*bZoIr*37si<`!bG6P)VRz^%X`8uTr8z+GVu9pl$Su(Q}!WmQD5i$5$%vQ}i+ zPCt+uTJp;AXxJGtB%C@G6zsf5Ulao%rx%M7d5zXGRRX8ysq*-wh5|i=oSv2|}l0c2TKGnSOD7`5n%9(GG%S%jEau zVE5AweyRoI#`+gf7FX_I)>Dcz!;F~XC6L(0nrA)#G;BNyYp{8g36R+uTg(c@0_1l^<~-?p6~_-_ z#YC+6Qpur5N~&Zki5@c8`A;M>1B>@dhcJRk&tZT({0FHCs57rdrI7>#SrQTA)}$ zCOkRLm zXr>m-v#73j(d?7lbr+a5Q^Oho0 z+853DV1Oi;p?^k>(IUZ0tr1&m&@j54^c|DDD7rOjot%y`a>g`9r*{^o#m0!sDl1G! zAO%>~BH4goN`+p0L9 zG3s+0ikQ98{HNdp&VC%RYA-|B=*qk^Yhx3iG;JRf+d=&W)?Wwg%vI(6%Vg1xJ{yoB z-!3~#US5Ls0qkqogj$BVe7g;|#VzB6j0ZXVFyrC(&WzRpDGa`b>n zC*>>SOdh?@A^8FQDeCSfAQe&{#wczAx>45+yot46DF4<)DD)1)x;z)T&cFA0#;`@3i7tStP1q-D!oL zw9UGHN^wUn%Kucp8KQhYv*I%*`C*p3g3E)GDoCWSvb}NF?mE`^a#Okd8n%i%zU`lo z6J4E_x0p8Qn*9-f#~wWp-|A4gmm85m&^ zrXP=$L84EE4#^uyW-2Jv@iNq8R3cH%Q&h`T^ag~7Vww5?3Oiz_t5*wi&a703j??u$ zOxDif+H)T&wDHt2WzMxccS$!g>gW#v$Uj&+%i-W65;FNnVxJyxVkg>Wgk)-5iU8+g zbcTdgDl?J|Xuy}mpmw&({ zbn5(EMWlh6o0+vyE;qp5dT#wpMCTg(IYne-(zQ%or-*WVI+1c&K&@b2_F&Ii%xh4Y zSc9lC(TP|nj$PPAa$)IIS*PuJc;u$iiT_@##$6P8bu3 z!4w{ezetY4n)RI+1xbYz63mACs5AjtjFdBH;*-+c=SMFHi*K$y)Y3+_}{g0Zn`<-C~f> zdvIwr$mjRsqMsM|jO;+(bFmE2^k2BJ4fDL5{qwzW1_)UGFfQ0XMFzL#f5)W*NEtxO z3%Jk_C#EZC0P31x>Tj=O2Tn zx0&azPk@#k7HtfT3^TavhL50}Oo&mz_ZB3@}9dMr;zD z;dT-2ZWqwl`$e>$8-T{%FM1L5?izms+d;5U`3+}5b4*5l!}FkTs_{2++*}jG-*^{j Wz@~5F@)aSi>se|nZNYgY#vgJpd_#rIW3E^R{?p;ZnT-|%S zcduDtB3y)r&AAji|^8hn@O>r-NU0U2w)X@ZTOHPd>+ z@A`SINH0ZkADcJDWiyhXknPk1+RH!_-ew~c=%co)ALu)aPwchStJXc1|;;&nVE zK5h$8Sl!DYAwUg-6*7RJ?}H&MT~WI5%|36alb<+uAZ$S!bkGXRRtfj4=p2<|T4py!y?9&%SHRx#HDlzI*hmkN)??yS&ehePw3f?~Z?Q&#yhy|Cuje zed_y9f9YfY_pR}7|KlINVc+0mSHJS~oxj)frKRz&ef~2qyma!z}kw*QnQLKZeELtUj#jFo455(P6eRNvhM)?J$y%YcT<$))keDd?+-A}tdWV<4d zG{5TA&vNaB?XuT8iFPt)`Jzl6)(ENMbv0B+!Bg#D&G|;@;yDGSi8Z3rH-ei2e%65Q zYUp(AuX_Pp54w*1zuQHX^Qj&A41h&DIX=VZa^<>0A!BE>d(CWlPc=n_@6Gi}Y2~ z;_2z%`1`+^I;w$+G=c0+qtA7$Le%5NDde*@hOUocjJFi|tmCAd`PC(`CzVL&(KPgm z^zFjvO94iI3NVt4OU+O6!Z2*q3xRamkG&InmgO#iF=+0_pK0dai@-atyD{>)10$cT zwP4w>#tN{_F5unmYE(1^X&JS0_~IjEtbI!R+14W7XOQbf^~jDoDP6l;!j$@1)0gy2 z^Pgk1tq5%A+hcim2&|FpBb~a&pZUU@KAAuJghc;a`+|xG<`RS{5r4VVI4e# z5iucN3Hn)^3rWw9h2_ z;O(Oa$DGle?b``V+bIfbNU7$C=GOCMLj94YG+1Ih$W*ezR4&aZNJUo$V0e)+e5O|a znK6ijgqM;A>jl5{fu63v`vTu}*tEiSY;BoZ@)wGQGqVplS^umv=Xg#b>kz@nqdo7s zjuBk<;oqAdA0^)f?{zO3q!E8b8lxinx+d(W9{Jv#&-_~DpE_l zkr+s%Y0#_94S}|I{C{p$mw0DjHy#Y=LY@st)H1tEfMf!h+2E0Zq9&!x7AdFkXTHbl z`{N(qp)P@{cjA$Mc;@w!FTGW~GcXtr5_;WMqD!O`>*u%*zOgBk?c>jc@HPJD$DVhG z&uay@Wt_a-kFXcsJNgR>))Dkj(^dw(kZa{tA9NOfZS+v^D#W1Q6h5oU)(V8A+=f7< zu?OJ=2zwD;h_DagMF>RsVuYI!UV<=+a3BUe&RYE(7YA8K#4XdYn`8PYjU-qmOtOaJ z?=+9c1P~HV1iKM-S*J0w_t10O;9Q1Y7v+E9yHQICpU7#{c@qyB`#}LCvH4#8%F+jq zzGLEem`>-rMdXXv5mFzuJ)KyyrqMHtPE+StOnkKP2+^_bY@DVJiKK#+LG`Zx-_74!Z}DQ%M(Dt{wWFh+9BP0Tr_t>GNq2GOtD#tekQ) zq-f3F2`L#AX(!OKjwZnaD#e_7;=sqLOUP->lKM8PADGKL{d68T}`2(1dP-p z-y;wLDvm`cnHg(|)<`q5shv-q$R@3lJyDC5gd7J^B*{@4W~mNAzuZ_xw?Ibqd?;k} zzApYVD$y_jx}zSM@R6((42iK0u6^>LbxMVK-pF30&c5Qz@X8U zpgueZl?$1tm0eE{8i_ND3dn6WJNGeQ@__b zgc2<_SCZEt0~vglW*e?u2m>%oe6+$dup}PBb3Z0#*jrSY6QRJ@ZPl+RMyC42VhGk7 z3^bV4)`RQ{P0;eoPrHFh`HR<${>JxG=E;8fjYF}0m^(J*GJ`7M$d9&y3**d19A?Zi zvLf*(HA0D#!@#`>I8W5#6c{xGn(_w``hGQEEX*_DftPM4`p(`fx&lDz6$fukqDM+; zkY(_^Gl6hF=(r39X@#IiB$_GXqDYcy@i^PYjs#rwJ4rHzqw6(p!N-l40xItHSh5>w zmO+NE(tWfK&wCO`E1}b*sY{vRP#`UN*$jL`z-K@gY?P8Qfwb###ma2nL5!Cc2ko&V zn;5#S2h`tJsJ~d-6kI}3n4<C>g!4;Rv8JcdbG*ufa-LRpduAEE+CfDJMT((d z!9#YEq`10}WVTu);@hn3Wln!o$7@N}W(URM@}@*Ic)wCb`$dt;5wLIrI%&PCD0v{u z%Iiog(?E)$K~wvW6_U#gubuGu^hBGS*Ywq$%A@L=|8OIW}MsM$Yl3 z9pA?SnwC&^y+fm>>XxW?D%1t#U+@ARs!?hIWT;VUr~z60w~&ZSu5-Cm^!&+UVa}az zrTkSQ=w_2d?8k-~I#Mpf$<28XKO7g)y)v&XQC4v@QxuW}AVHrWRyaw#BP6KJVs=Al zuKRPqA*^!Y+SqG$>u+QLe1GU2G z%Hpi!9k$C(`&=;P zjDVi1QhXao-#%++AD;J$l|oKwbaZCESIp;~j&!^3WA^}q?B2BhjZVHIX16%MrB;p! zhCaePff`MR2u>wPB6gs24;SJaU4AVoN_V*wbBbgkMcs8TXzEjtPZZ{g?aLy-gsej> zhzHNgR6G^}PQ}WgMbPa52`z&6B5Vly)3e2*Kk4LkgX$0$;96RG=zx+Hx2MK9Itepz z36fM%hLS^xV)2DEy&W)y7ws$VqWfAWr#t0z_yW)an$}E=)}c79F}i*e={ML}-@W8K zcZhlnMD=O2R%`aoD|ECSCV`Y1dA=vCH%>}R!MHh^#4Z8mG(1Nu5f}fQQc6Y zt@j{MqTh(n^RTm|i(M;x_W|w_#wOvayn3+FN)GD+?*mXt=r+vQ&{;~ish$B-V1;xb z%d!REpXhR&wwn%X67?svtEx$f?d%GR&H_rXHO|@TU~|n%iK`8qL)fcR zbjz@M%|T=89DOyA4V2s~PJRVT1spXhW+G%?t=@}!B$n-Tsam>z0Pt++G`((Z50G=U#_v3Lxxs-SPBT9jxnRQx>IEr`TMXxFYY|<2fck8~PLq;nGzn;YD z4!d*WdLz>-L%_l>x12eHq&3oXb4W}EpPGcDIwfv(qU2mb=);)A@nB?XvB}_NLtsb~ zu^28>Bg3w^QHB*|+=}7BaD*XQ?I0Bmmzlc=&&AJI15+G(VjP&etrMSjwtHdWj z#9zOs|Y@TS@<8p7y}rip4<%jq>bu7U%UIKGTcxG7$|0eFY) z0(?Bg8=I7WRn()eYy=beG>W%>MnFVVPL568kLC1uZgt z%)o0JZFQ*MWEO$A0>nf4V)0?ypRv7p#}^%MlhB}BCv>@+%pwS&r*#`nBM{w2jr+Gr zsDSjSCZtuTFY?x?}oCHz)^q#h|DtSvxJwLrB55P}hjRg;Ve*$cEaWx6|! z7R}PF83AyXj?0oqZ1@!ps1JX8Oi>#YcIRkPWO2U|KN+9O*5WC!I4%A)L-`eVNm9CR z(WG#_pX70Mg%ihu@ijRfQf~?g;P9_a?02?@{UKmqkMjW9!)o2$9#+~eJv?V9)Rh@L zOODl(%E-EGx3*VpN=Tqsrp1OJ8pEnh`5d$<$HGQZYoF?q7{u`x)e2?bu#VOr`-XcE zHhK=u@N%<0R;q~gOihCY2&3<1zrK}L-0^Uk3Ktt-c_xfG?Qm`ZwrQ1coxyNq^>DS( zOx0`90#uV?S_|7EsI;YyV0whTCAza#xeQnaJGlI6ivb*q(ztx5bI{aU{XzG5&dck5 zsTDRvh6jb=jR?K&BFC$dpQ9CEk7=b}I~9;1kvK@52OtEiLTlHma4K#meJ}UeWiD}= zN2TB#ISlws*GDa`!|;8Vq0Hj_e5^G9zG2L8T(F5p)m2f*=8zhm-vOQtCZ15lIzX{u zNb|mSU8Y_}&9wu;5~aw~Q`HJe8*s9&S$nRfo;J>FHse9eS?8S#7s?BEuBd&owNMX* zeksmt1pDQcrz-56j2fh&I}DR0SO%LWWHz*QnCSK^x{-B6=wX0ctBN1ie5=q51w(IT zUsUW)Ap%$cJ>ak!fk=mP{Op3=@Un3MNJY&#=|v!us4dD`%&=9XeKM2kW(J`zi4r<*y0;Bop?TItIk92h3E|QWa}e7id^ekSrWfpJ%{84fZJu%ehSK>gnwM$Fi$xYvUrrR6k#ArOk%b~ zm=taW@Un+eIyv^A+P-Q5pxIiui)K@Kp?LXRX{IJ#HKkY2|WV1>uIaDU2v05 z1yAi$aB`D&spE3frc}pxN^7cZE=JJg2936%t!l7qN`u`MYMU{Ll=Qx4tW!`BJj#;W zH?vR^+T7>i4dL{9pe&T3|5wOQZT41#!Np?E#Tp#%IAp%m2?9Asygx?a4ue9fw&P08 zRf1n35{it=8P>W5>~iXi7AvC%-nC=6xC5ewONu+(0Nks&(-|>C5iyw0^+p;C9##XEE}hLG!X0Sq@R1klF zFkFV6j=r^K#YT##fhMww*qO-^HM5oPy~s1-72yZM5~$}qp`~KkeJ&Hq+E5l&*kw=> zTlgJ%+Vq=VWnyKJswW zJnje;sGiWB0Ef`QeoRfXvZ5{Rljdoj+^wmM>3m{hQ zw?C%V>xxBmDR39HBRcL60ac%K1wU28O-P}7N&(AIREM;ElLTVFl1D9jN3o_qIQKm~ z={)Y87Bs&K5I4EyQ*acx9E+5FrzAFo*ZTdNoGwv>hIkOKX-)Od%xG#(K(!v8uF2NJ z9SqkM)z8N0u}aHR4FD1CWvJD^nx>eSU+|F{V9?3AKF_(B@`_iq*QjY=ARw}lb_1I1 z0UDpxQb3M`MG-Av@7t|4MQe*e*a;41JqHVSa)+11!patK1BBqEnmcr{NYz(ndDxKJ z0tdh~q9+x3e7xO;^FfK?gFr)3RN5w(_Jk0W953=7QV(>54fzc&1*n!57ah)!t1~Es z%1IR!7eDsm9H8@JiGAwjc_onaV&OiklFxIlGLV*=rzP!U9c>)`yhwTm+tBCAnt``K16^ zw6r#5T<>Fb*Hg&^$WQP$kt`Uq{sB zHGexM@NjVxNRegT8?&6oihxWEV4?#DX^eJ%Hy&$N z0%Y(~lGYFs(g#oj`qT=5dqZXsg9J7jO=xg}ks}-aAXqTwX*mABHqNRz85VsaR`fBW z=x5R2@T@}Z5&gL@{4_@7;T~9wi;QFHRXg}OIHlxtSo5h42x$?l*#cRjUVR8*liFaY z23C{U4N(?@=adiwi7|xck4Fq;v59EL!LPVN6N`X}!6SQoG;eFT7A)8QNr(bU@bhoED50FtYIhL0|oWgQ1<`IaNG zg@0pw4|wvoo*KeO9%BYd0WZyFA6!Y`uUJXo53B@mo)pnTo`;H}U0xe)SOF!kR|4oO zRs!e)X#(X|g?_FsH4e8(g|n;L@tV<&a@HVuwP(UxQ!}o810;!>@z)V>E@(-)Y2~~@ z0th5ccpDlamUd`x(jg(4dK=6GR@b$~w^S->a_)kY^XzE5%yIy0*3Gs;%@On?{B_?P zj1`WGet$(7S@ z9mXJ$r_smRpq8{J?FFIf7w^QGlBab65@&M3NqLA9Qjorrw#~hNyyn~MOI1zehNs|1xV7Ro(A4VC-Fd%7R?p_GWr+R zWw3ZGrb;@5OBFJ}pN2{otRlB~8NXW=Te90axrXR*Vf1w&gTt!p5F*DkYt^c}cHx<- zh6}BG$k2AdA-`yA_u#Rwz!Rk0>}_z-Fblsa3KRjbWy?3n&tI21X#Gd zy1t4(VZ!cP5%1w!d_ImpYgN?h{xR`3AcOZTp0jQdC!|4ZCA1@X<#$9;1cC!Fhp z`1Pb^OBnw$?pwq7f8xF^jQ<4p?LmAC=_0>yZDsxF_@mYMwy#y=+e;$;qiTG|r>gNA zcyd?>XD7?5 bool: + def detect_format(filename: str) -> DDDFileFormat: # OpenDiS data files always start with "dataFileVersion = ..." + # ExaDiS restart files always start with "# ExaDiS restart file" try: with open(filename, "r") as f: - line = f.readline() - return line.strip().startswith("dataFileVersion =") + line = f.readline().strip() + if line.startswith("dataFileVersion ="): + return DDDFileFormat.DATAFILE + elif line.startswith("# ExaDiS restart file"): + return DDDFileFormat.EXADIS_RESTART + else: + return DDDFileFormat.UNKNOWN + except OSError: + return DDDFileFormat.UNKNOWN + + @staticmethod + def detect(filename: str) -> bool: + try: + return __class__.detect_format(filename) != DDDFileFormat.UNKNOWN except OSError: return False def scan(self, filename: str, register_frame: Callable[..., None]) -> None: # OpenDiS files contain a single static snapshot - register_frame(frame_info=(0, 0)) - + file_format = __class__.detect_format(filename) + if file_format == DDDFileFormat.EXADIS_RESTART: + try: + header, *_ = __class__.parse_exadis_header(filename) + step_number = header["step"] + register_frame(frame_info=(step_number), label=f"Timestep {step_number}") + except OSError: + register_frame(frame_info=(0)) + else: + register_frame(frame_info=(0)) + @staticmethod def skip_line(line: str) -> bool: return line.startswith("#") or not line.strip() - + @staticmethod def parse_number(number: str) -> int | float: try: return int(number) except ValueError: return float(number) - + @staticmethod - def read_array(f: TextIOWrapper) -> list[int | float]: - # Reads one value per line until a line containing "]" is encountered. + def read_array(f: TextIOWrapper) -> list[int | float | list[int | float]]: + # Reads scalar values or vector rows until a line containing "]" is encountered. # Called after the opening "[" has already been consumed by parse_header. array = [] - line = f.readline().strip() - while "]" not in line: - if line.startswith("#") or not line: + while line := f.readline(): + array_end = "]" in line + line = line.split("#", 1)[0].replace("]", "").strip() + if not line: + if array_end: + break continue - array.append(__class__.parse_number(line)) - line = f.readline().strip() + tokens = line.split() + values = [__class__.parse_number(token) for token in tokens] + array.append(values[0] if len(values) == 1 else values) + if array_end: + break return array @staticmethod def parse_header(f: TextIOWrapper) -> dict[str, Any]: key = None - header_end = "END OF DATA FILE PARAMETERS" + header_end = "nodalData" header = {} while line := f.readline(): if header_end in line: @@ -94,6 +106,9 @@ def parse_header(f: TextIOWrapper) -> dict[str, Any]: if "=" in line: key = line.split("=")[0].strip() + if key == "domainDecomposition": + continue + if "[" in line: # Value is a multi-line bracketed array; read until closing "]" array = __class__.read_array(f) @@ -118,71 +133,67 @@ def parse_domain_decomposition( return data @staticmethod - def parse_primary_line(line: str) -> Node: - # Lines may be prefixed with a domain tag: "domain,tag x y z ..." - # Split on "," and take the last part to strip any domain prefix. - line = line.strip().split(",")[-1] - tokens = line.split() - return Node( - int(tokens[0]), - [float(t) for t in tokens[1:4]], - int(tokens[4]), - [], - int(tokens[5]), - ) + def parse_nodal_data(f: TextIOWrapper, node_count: int | None = None) -> tuple[np.ndarray, np.ndarray]: + # Node records alternate between a primary line (tag, position, num_arms, constraint) + # and secondary lines (one entry per arm with bvec and nvec). + nodes = np.empty((node_count, 6), dtype=float) if node_count else [] + segs = [] + nodes_map = {} + pending_segs: dict[tuple[int, int], list[int]] = {} + node_index = 0 + + content = re.sub(r"(?m)#.*$", "", f.read()).replace(",", " ") + values = np.fromstring(content, sep=" ") + cursor = 0 + + while cursor < values.size: + domain = int(values[cursor]) + tag = int(values[cursor + 1]) + num_arms = int(values[cursor + 5]) + node = [ + domain, tag, + values[cursor + 2], values[cursor + 3], values[cursor + 4], + int(values[cursor + 6]) + ] + cursor += 7 + + if node_count: + nodes[node_index] = node + else: + nodes.append(node) + node_tag = (node[0], node[1]) + nodes_map[node_tag] = node_index - @staticmethod - def parse_secondary_line(f: TextIOWrapper, num_arms: int) -> list[Arm]: - # Arms are stored in "secodary" lines, one per line arm and *num_arms* entries. - # Each arm's data consists of: arm_tag, bvec (3 floats), nvec (3 floats). - # This function reads all arms belonging to a node and returns them as a list - arms = [] - current_arm = [] - while len(arms) < num_arms: - line = f.readline() - if __class__.skip_line(line): - continue - if "," in line: - line = line.strip().split(",")[-1] + for seg_index in pending_segs.pop(node_tag, []): + segs[seg_index][1] = node_index - tokens = line.split() + for _ in range(num_arms): + arm_tag = (int(values[cursor]), int(values[cursor + 1])) + if arm_tag in nodes_map: + cursor += 8 + continue - if len(current_arm) == 0: - # First token on the first line of an arm entry is the neighbor tag - current_arm.append(int(tokens[0])) - tokens = tokens[1:] - assert len(tokens) == 3 or len(tokens) == 6 - while tokens: - current_arm.append([float(t) for t in tokens[:3]]) - tokens = tokens[3:] - - assert len(current_arm) <= 3 - if len(current_arm) == 3: - arms.append(Arm(*current_arm)) - current_arm = [] - return arms + pending_segs.setdefault(arm_tag, []).append(len(segs)) + segs.append([ + node_index, -1, + *values[cursor + 2:cursor + 5], + *values[cursor + 5:cursor + 8] + ]) + cursor += 8 - @staticmethod - def parse_nodal_data(f: TextIOWrapper) -> list[Node]: - # Node records alternate between a primary line (tag, position, num_arms, constraint) - # and secondary lines (one entry per arm with bvec and nvec). - data = [] - primary_line = True - while line := f.readline(): - if __class__.skip_line(line): - continue + node_index += 1 - if primary_line: - data.append(__class__.parse_primary_line(line)) - primary_line = False - if not primary_line: - node = data[-1] - node.arms = __class__.parse_secondary_line(f, node.num_arms) - primary_line = True - return data + if pending_segs: + missing_tag = next(iter(pending_segs)) + raise AssertionError(f"Nodal data references missing node tag {missing_tag}") + + if node_count: + assert node_index == node_count + + return np.asarray(nodes), np.asarray(segs, dtype=float).reshape((-1, 8)) @staticmethod - def parse_body(f: TextIOWrapper, num_domains: int) -> dict[str, Any]: + def parse_body(f: TextIOWrapper, num_domains: int, node_count: int | None = None) -> dict[str, Any]: key = None body = {} while line := f.readline(): @@ -196,210 +207,314 @@ def parse_body(f: TextIOWrapper, num_domains: int) -> dict[str, Any]: body[key] = __class__.parse_domain_decomposition(f, num_domains) if key == "nodalData": - body[key] = __class__.parse_nodal_data(f) + nodes, segs = __class__.parse_nodal_data(f, node_count) + body["nodes"] = nodes + body["segs"] = segs return body @staticmethod - def point_in_cell(cell: SimulationCell, point: np.ndarray) -> bool: - if np.any(point < cell[:, 3]): - return False + def parse_data_file(filename: str): + with open(filename, "r") as f: + header = __class__.parse_header(f) + num_domains = int(np.prod(header["dataDecompGeometry"])) + with open(filename, "r") as f: + body = __class__.parse_body(f, num_domains, header.get("nodeCount")) + + header["segmentCount"] = len(body["segs"]) + pbc = (True, True, True) + cell = np.zeros((3, 4)) + cell[:, 3] = header["minCoordinates"] for i in range(3): - if point[i] >= cell[i, 3] + cell[i, i]: - return False - return True + cell[i, i] = header["maxCoordinates"][i] - header["minCoordinates"][i] + return header, num_domains, pbc, cell, body + @staticmethod - def get_next_start_node( - nodes: list[Node], start: int - ) -> tuple[int, Node] | tuple[None, None]: - while start < len(nodes): - # Any node that is not in a chain, either junction or end point - if nodes[start].num_arms != 2 and not nodes[start].processed: - return start, nodes[start] - start += 1 - return None, None + def parse_exadis_header(filename: str): + header: dict[str, Any] = {"dataDecompGeometry": [1, 1, 1]} + pbc = (True, True, True) + cell = np.zeros((3, 4)) - @staticmethod - def trace_path(start: Node, arm: Arm, node_dict: dict[int, Node]) -> list[Node]: - # Follow a chain of degree-2 nodes from start through arm until reaching - # a junction (num_arms != 2) or an already-processed node. - segment = [start] - prev_tag = start.node_tag - curr = node_dict[arm.arm_tag] - while curr.num_arms == 2 and not curr.processed: - segment.append(curr) - curr.processed = True - prev = curr - for next_arm in curr.arms: - if next_arm.arm_tag != prev_tag: - prev_tag = curr.node_tag - curr = node_dict[next_arm.arm_tag] + def parse_values(tokens: list[str]) -> int | float | str | list[int | float]: + values = [] + for token in tokens: + try: + values.append(__class__.parse_number(token)) + except ValueError: + return " ".join(tokens) + return values[0] if len(values) == 1 else values + + line_num = 0 + with open(filename, "r") as f: + while line := f.readline(): + line_num += 1 + + if __class__.skip_line(line): + continue + tokens = line.split() + if not tokens: + continue + + key = tokens[0] + if key == "Nnodes": + header["nodeCount"] = int(tokens[1]) break - if curr is prev: - # Degenerate two-node loop: both arms point back where we came - # from, so close the segment by stepping back there. - curr = node_dict[prev_tag] - break - segment.append(curr) - return segment + elif key == "pbc": + pbc = tuple(bool(int(t)) for t in tokens[1:4]) + header["pbc"] = [int(t) for t in tokens[1:4]] + elif key == "H": + h = [float(t) for t in tokens[1:10]] + cell[:, :3] = np.asarray(h, dtype=float).reshape((3, 3)) + header["cellVectors"] = [ + h[0:3], + h[3:6], + h[6:9], + ] + elif key == "origin": + origin = [float(t) for t in tokens[1:4]] + cell[:, 3] = origin + header["origin"] = origin + header["minCoordinates"] = origin + elif key == "crystal" and len(tokens) > 2: + header[f"crystal_{tokens[1]}"] = parse_values(tokens[2:]) + else: + header[key] = parse_values(tokens[1:]) + + return header, pbc, cell, line_num @staticmethod - def walk_lines(nodes: list[Node]) -> Generator[float, None, list[Line]]: - # Returns the collected Line objects; callers must use "yield from" to - # receive the return value while forwarding progress yields upstream. - node_dict = {node.node_tag: node for node in nodes} - - lines = [] - start, start_node = __class__.get_next_start_node(nodes, 0) - while start_node is not None and start is not None: - yield 0.0 - start_node.processed = True - segments = [] - for arm in start_node.arms: - if not node_dict[arm.arm_tag].processed: - segments.append(__class__.trace_path(start_node, arm, node_dict)) - if segments: - lines.append(Line(segments)) - start, start_node = __class__.get_next_start_node(nodes, start) - - # Any node left unprocessed belongs to a closed loop, which has no junction - # or endpoint to start from. Pick an arbitrary node on each remaining loop. - for node in nodes: - if node.processed or not node.arms: - continue - yield 0.0 - # Marking the start node first makes trace_path stop when it comes back - # around to it and append it, so the returned segment is closed. - node.processed = True - lines.append(Line([__class__.trace_path(node, node.arms[0], node_dict)])) - return lines + def parse_exadis_file(filename: str): - @staticmethod - def walk_line( - segment: list[Node], - ref_point: np.ndarray, - cell: SimulationCell, - positions: list[np.ndarray], - sections: list[int], - bvecs: list[np.ndarray], - nvecs: list[np.ndarray], - counter: int, - ) -> Generator[float, None, np.ndarray]: - # ref_point carries the last unwrapped position across calls so that - # delta_vector can resolve PBC images consistently along the full path. - for node_id in range(1, len(segment)): - yield 0.0 - n0 = segment[node_id - 1] - n1 = segment[node_id] - - # delta_vector returns the shortest-image displacement, respecting PBC - p0 = ref_point + cell.delta_vector(ref_point, n0.pos) - p1 = p0 + cell.delta_vector(p0, n1.pos) - ref_point = p1 - - # Avoid duplicating the shared start point when appending the next - # segment from the same junction (sections[-1] already equals counter). - new_segment = len(sections) == 0 or sections[-1] != counter - - if new_segment: - positions.append(p0) - sections.append(counter) - positions.append(p1) - sections.append(counter) - - # Burgers vector and normal come from the arm in n0 that points to n1 - matching_arm = None - for arm in n0.arms: - if arm.arm_tag == n1.node_tag: - matching_arm = arm - break - if matching_arm is None: - raise ValueError(f"Could not find matching arm for node {n1.node_tag}") - if new_segment: - bvecs.append(matching_arm.bvec) - nvecs.append(matching_arm.nvec) - bvecs.append(matching_arm.bvec) - nvecs.append(matching_arm.nvec) + header, pbc, cell, line_num = __class__.parse_exadis_header(filename) - return ref_point + node_count = header["nodeCount"] + nodes = np.loadtxt(filename, skiprows=line_num, max_rows=node_count) + if nodes.shape[1] == 5: + nodes = np.hstack((np.zeros((node_count, 1)), nodes)) + elif nodes.shape[1] == 7: + nodes = nodes[:,1:] - def parse(self, data: DataCollection, filename: str, **kwargs: Any) -> Generator[str | float, None, None] | None: # type: ignore[override] + line_num = line_num + node_count + 2 + segs = np.loadtxt(filename, skiprows=line_num) + header["segmentCount"] = segs.shape[0] - with open(filename, "r") as f: - header = __class__.parse_header(f) - num_domains = np.prod(header["dataDecompGeometry"]) - body = __class__.parse_body(f, num_domains) - assert len(body["nodalData"]) == header["nodeCount"] + if "minCoordinates" not in header: + header["minCoordinates"] = [0.0, 0.0, 0.0] + if "cellVectors" not in header: + raise ValueError("ExaDiS restart file is missing H cell matrix") + + cell_vectors = np.asarray(header["cellVectors"], dtype=float) + lower = np.asarray(header["minCoordinates"], dtype=float) + upper = lower + np.sum(cell_vectors, axis=0) + header["maxCoordinates"] = upper.tolist() + + num_domains = 1 + body = {"nodes": nodes, "segs": segs} + + return header, num_domains, pbc, cell, body + + @staticmethod + def generate_connectivity( + nodes_constraint: np.ndarray, segs_indices: np.ndarray + ) -> list[list[tuple[int, int, int]]]: + conn: list[list[tuple[int, int, int]]] = [[] for _ in nodes_constraint] + for i, (n1, n2) in enumerate(segs_indices): + conn[n1].append((n2, i, 1)) + conn[n2].append((n1, i, -1)) + return conn + + @staticmethod + def build_links( + cell, + nodes_pos: np.ndarray, nodes_constraint: np.ndarray, + segs_indices: np.ndarray, segs_bvec: np.ndarray, segs_nvec: np.ndarray, + ): + conn = __class__.generate_connectivity(nodes_constraint, segs_indices) + is_discretization = [ + len(node_conn) == 2 and constr == 0 + for constr, node_conn in zip(nodes_constraint, conn) + ] + visited = [-1] * len(nodes_constraint) + try: + seg_deltas = np.asarray( + cell.delta_vector( + nodes_pos[segs_indices[:,0]], nodes_pos[segs_indices[:,1]] + ) + ) + if seg_deltas.shape != segs_bvec.shape: + raise TypeError + except (TypeError, ValueError): + seg_deltas = None + + def next_connection(node_index: int, prev_index: int) -> tuple[int, int, int]: + node_conn = conn[node_index] + neighbor, seg_index, order = node_conn[0] + if neighbor != prev_index: + return neighbor, seg_index, order + return node_conn[1] + + num_physical_nodes = 0 + nl = 0 + + link_positions = [] + link_sections = [] + link_bvecs = [] + link_nvecs = [] + + # Links connected to physical nodes: junctions, endpoints, and constrained nodes. + for n, node_conn in enumerate(conn): + if is_discretization[n]: + continue + if visited[n] == -1: + visited[n] = num_physical_nodes + num_physical_nodes += 1 + + for nn, il, order in node_conn: + if visited[nn] == -1: + position = nodes_pos[n] + link_positions.append(position) + if seg_deltas is None: + position = position + cell.delta_vector(position, nodes_pos[nn]) + else: + position = position + order * seg_deltas[il] + link_positions.append(position) + link_sections.append(nl) # first node + link_sections.append(nl) # second node + link_bvecs.append(order * segs_bvec[il]) + link_nvecs.append(order * segs_nvec[il]) + prev = n + + if not is_discretization[nn]: + visited[nn] = num_physical_nodes + num_physical_nodes += 1 + nl += 1 + continue + + visited[nn] = 1 + while is_discretization[nn]: + prev, (nn, ilp, order) = nn, next_connection(nn, prev) + if seg_deltas is None: + position = position + cell.delta_vector(position, nodes_pos[nn]) + else: + position = position + order * seg_deltas[ilp] + link_positions.append(position) + link_sections.append(nl) + + if not is_discretization[nn]: + if visited[nn] == -1: + visited[nn] = num_physical_nodes + num_physical_nodes += 1 + nl += 1 + else: + visited[nn] = 1 + elif not is_discretization[nn] and nn > n: + position = nodes_pos[n] + link_positions.append(position) + if seg_deltas is None: + position = position + cell.delta_vector(position, nodes_pos[nn]) + else: + position = position + order * seg_deltas[il] + link_positions.append(position) + link_sections.append(nl) # first node + link_sections.append(nl) # second node + link_bvecs.append(order * segs_bvec[il]) + link_nvecs.append(order * segs_nvec[il]) + nl += 1 + + # Closed loops made only of discretization nodes. + for n, node_conn in enumerate(conn): + if visited[n] != -1 or not is_discretization[n]: + continue + visited[n] = num_physical_nodes + num_physical_nodes += 1 + prev = n + nn, il, order = node_conn[0] + position = nodes_pos[n] + link_positions.append(position) + if seg_deltas is None: + position = position + cell.delta_vector(position, nodes_pos[nn]) + else: + position = position + order * seg_deltas[il] + link_positions.append(position) + link_sections.append(nl) # first node + link_sections.append(nl) # second node + link_bvecs.append(order * segs_bvec[il]) + link_nvecs.append(order * segs_nvec[il]) + while nn != n: + visited[nn] = 1 + prev, (nn, ilp, order) = nn, next_connection(nn, prev) + if seg_deltas is None: + position = position + cell.delta_vector(position, nodes_pos[nn]) + else: + position = position + order * seg_deltas[ilp] + link_positions.append(position) + link_sections.append(nl) + nl += 1 + + link_positions = np.asarray(link_positions) + link_sections = np.asarray(link_sections) + link_bvecs = np.asarray(link_bvecs) + link_bvecs = link_bvecs[link_sections] + link_nvecs = np.asarray(link_nvecs) + link_nvecs = link_nvecs[link_sections] + + return link_positions, link_sections, link_bvecs, link_nvecs + + def parse(self, data: DataCollection, filename: str, frame_info: Any, **kwargs: Any): + + # Read data + file_format = __class__.detect_format(filename) + if file_format == DDDFileFormat.DATAFILE: + header, num_domains, pbc, cell, body = __class__.parse_data_file(filename) + elif file_format == DDDFileFormat.EXADIS_RESTART: + header, num_domains, pbc, cell, body = __class__.parse_exadis_file(filename) + else: + raise ValueError(f"OpenDiSFileReader does not support file type") for k, v in header.items(): data.attributes[k] = v - cell = np.zeros((3, 4)) - cell[:, 3] = header["minCoordinates"] - for i in range(3): - cell[i, i] = header["maxCoordinates"][i] - header["minCoordinates"][i] - cell = data.create_cell(cell, pbc=(True, True, True)) + cell = data.create_cell(cell, pbc=pbc) # Scale line/node width to ~0.1 % of the cell diagonal for visual clarity self.lines_vis.width = 1 / 1000 * np.linalg.norm(cell[:3, :3].diagonal()) + # Create nodes + nodes = body["nodes"] + assert len(nodes) == header["nodeCount"] + nodes_pos, nodes_constraint = nodes[:,2:5], nodes[:,5] + nodes_pos = cell.wrap_point(nodes_pos) + particles = data.create_particles(count=header["nodeCount"], vis_params={'title': 'Nodes'}) - identifier = particles.create_property("Particle Identifier") + tags = particles.create_property("Node Tag", dtype=int, components=('domain', 'index'), data=nodes[:,0:2].astype(int)) particle_type = particles.create_property("Particle Type") - positions = particles.create_property("Position") - num_arms = particles.create_property("Num Arms", dtype=int) - constraint = particles.create_property("Constraint", dtype=int) + positions = particles.create_property("Position", data=nodes_pos) + constraints = particles.create_property("Constraint", dtype=int, data=nodes_constraint) node_type = particle_type.add_type_name("Node", data.particles) node_type.radius = self.lines_vis.width / 2 particle_type[:] = node_type.id - for i, node in enumerate(body["nodalData"]): - identifier[i] = node.node_tag - positions[i] = node.pos - num_arms[i] = node.num_arms - constraint[i] = node.constrain - yield i / len(body["nodalData"]) - - positions = [] - sections = [] - bvecs = [] - nvecs = [] - - # Line objects returned by the generator via its StopIteration value. - lines = yield from self.walk_lines(body["nodalData"]) - - ref_point = np.asarray(data.cell[:, 3]) - for counter, line in enumerate(lines): - ref_point = data.cell.wrap_point(ref_point) - for segment in line.segments: - # Reverse so the segment walks from the far end back to the junction, - # keeping ref_point continuous across consecutive segments of the line. - segment = list(reversed(segment)) - ref_point = yield from self.walk_line( - segment, - ref_point, - cell, - positions, - sections, - bvecs, - nvecs, - counter, - ) - self.lines_vis.color = node_type.color - # Empty lists would be interpreted as shape (0,) arrays, which cannot be - # assigned to the (0, 3) vector properties, so give them the right shape. - positions = np.asarray(positions, dtype=float).reshape((-1, 3)) - bvecs = np.asarray(bvecs, dtype=float).reshape((-1, 3)) - nvecs = np.asarray(nvecs, dtype=float).reshape((-1, 3)) + # Create lines + segs = body["segs"] + assert len(segs) == header["segmentCount"] + segs_indices = segs[:,0:2].astype(int) + segs_bvec, segs_nvec = segs[:,2:5], segs[:,5:8] - lines = data.lines.create("Arms", count=len(positions), vis=self.lines_vis) - lines.create_property("Position", data=positions) - lines.create_property("Section", data=sections) - lines.create_property("Burgers vector", data=bvecs, components=["X", "Y", "Z"]) + link_positions, link_sections, link_bvecs, link_nvecs = \ + __class__.build_links(cell, nodes_pos, nodes_constraint, segs_indices, segs_bvec, segs_nvec) + + lines = data.lines.create("Dislocations", count=len(link_positions), vis=self.lines_vis) + lines.create_property("Position", data=link_positions) + lines.create_property("Section", data=link_sections) + lines.create_property("Burgers vector", data=link_bvecs, components=["X", "Y", "Z"]) lines.create_property( - "Burgers vector magnitude", data=np.linalg.norm(bvecs, axis=1) + "Burgers vector magnitude", data=np.linalg.norm(link_bvecs, axis=1) ) - lines.create_property("Normal vector", data=nvecs, components=["X", "Y", "Z"]) + lines.create_property("Normal vector", data=link_nvecs, components=["X", "Y", "Z"]) + + print(f"Number of nodes: {header['nodeCount']}") + print(f"Number of segments: {header['segmentCount']}") + print(f"Number of links: {link_sections[-1]+1 if len(link_sections) > 0 else 0}")