<div dir="ltr"><div>Another potential interesting project regarding SIMD is xsimd: <a href="https://github.com/xtensor-stack/xsimd">https://github.com/xtensor-stack/xsimd</a><br></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Fri, 17 Apr 2020 at 21:36, Even Rouault <<a href="mailto:even.rouault@spatialys.com">even.rouault@spatialys.com</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><u></u>
<div style="font-family:"monospace";font-size:9pt;font-weight:400;font-style:normal">
<p style="margin:0px;text-indent:0px">On vendredi 17 avril 2020 13:29:17 CEST Charles Karney wrote:</p>
<p style="margin:0px;text-indent:0px">> I wonder whether it's possible to get the benefits of vectorization</p>
<p style="margin:0px;text-indent:0px">> without massive of changes to the code. </p>
<p style="margin:0px;text-indent:0px"> </p>
<p style="margin:0px;text-indent:0px">That's exactly what I had in mind.</p>
<p style="margin:0px;text-indent:0px"> </p>
<p style="margin:0px;text-indent:0px">> Perhaps the basic projection</p>
<p style="margin:0px;text-indent:0px">> functions could be templated to allow Eigen arrays of PJ_LP's to be</p>
<p style="margin:0px;text-indent:0px">> passed.  Eigen already has overloads for componentwise arithmetic</p>
<p style="margin:0px;text-indent:0px">> operations and handles SIMD vectorization automatically.</p>
<p style="margin:0px;text-indent:0px"> </p>
<p style="margin:0px;text-indent:0px">I didn't know Eigen, but from a quick look & try, it seems I've done something similar.</p>
<p style="margin:0px;text-indent:0px">At least for types like Vector2d, Vector4d. But Eigen is quite limited in the operations it supports. Just +, -, *, /. No comparisons, sqrt, sin, cos, etc.</p>
<p style="margin:0px;text-indent:0px">Which my prototype supports (through sleef for sin, cos, etc).</p>
<p style="margin:0px;text-indent:0px"> </p>
<p style="margin:0px;text-indent:0px">Adding support for functions that take Eigen types, like sqrt(Vector4d) should be possible, but comparison operators require to be able to modify the class definition.</p>
<p style="margin:0px;text-indent:0px"> </p>
<p style="margin:0px;text-indent:0px">Another pain point I found with Eigen, is that it refuses to do things like</p>
<p style="margin:0px;text-indent:0px">a * b on vectors because of potential ambiguity of the operator in the context of linear algebra (matrix multiplication, cross product, member-to-member multiplication ?). But in a PROJ world where we'd port from double to multiple-double, we just want member-to-member multipication, and with Eigen you have to do it explictly with a.cwiseProduct(b). Would require more changes and would clutter the readability of the code.</p>
<p style="margin:0px;text-indent:0px"> </p>
<p style="margin:0px;text-indent:0px">> A major advantage of this approach is that PROJ doesn't need to get into</p>
<p style="margin:0px;text-indent:0px">> the weeds with SIMD instructions.  So when a new instruction set comes</p>
<p style="margin:0px;text-indent:0px">> along we can (I hope) rely on the maintains of Eigen to do the hard</p>
<p style="margin:0px;text-indent:0px">> work.  (I see that there's already some interoperability with Eigen and</p>
<p style="margin:0px;text-indent:0px">> CUDA.)</p>
<p style="margin:0px;text-indent:0px"> </p>
<p style="margin:0px;text-indent:0px">Admitedly I've only prototyped SSE and AVX. But supporting other archs would only require to add a specialization in the generic template stuff (not a big deal. It is of the order of 10 intrinsincs per arch), not projection code. And for non-supported archs, the templated operations would just expand to a 1-coordinate at a time implementation, like the current serial code. But I guess if Intel is covered, we have already > 90% use cases that matter (nobody ever paid me to optimize something for non-Intel archs).</p>
<p style="margin:0px;text-indent:0px"> </p>
<p style="margin:0px;text-indent:0px">One potential downside (in the context of PROJ licencing) of Eigen is its license. It is MPL 2, weak copyleft, which would require binary distributions of PROJ to provide a way of accessing Eigen source code. Not a huge deal, but this would put additional requirements for proprietary shops using PROJ.</p>
<p style="margin:0px;text-indent:0px"> </p>
<p style="margin:0px;text-indent:0px">> If this approach is followed, I would also recommend that the basic</p>
<p style="margin:0px;text-indent:0px">> floating point type, double, be either templated or typedef'ed. </p>
<p style="margin:0px;text-indent:0px"> </p>
<p style="margin:0px;text-indent:0px">Yep, exactly my idea of template<typename T, int N> struct VF</p>
<p style="margin:0px;text-indent:0px"> </p>
<p style="margin:0px;text-indent:0px">-- </p>
<p style="margin:0px;text-indent:0px">Spatialys - Geospatial professional services</p>
<p style="margin:0px;text-indent:0px"><a href="http://www.spatialys.com" target="_blank">http://www.spatialys.com</a></p></div>_______________________________________________<br>
PROJ mailing list<br>
<a href="mailto:PROJ@lists.osgeo.org" target="_blank">PROJ@lists.osgeo.org</a><br>
<a href="https://lists.osgeo.org/mailman/listinfo/proj" rel="noreferrer" target="_blank">https://lists.osgeo.org/mailman/listinfo/proj</a><br>
</blockquote></div></div>