On 3/5/2009 10:11 AM, Dag Sverre Seljebotn wrote:
> Cython can relatively easily transform things like
>
> cdef int[:,:] a = ..., b = ...
> c = a + b * b
Now you are wandering far into Fortran territory...
> If a and b are declared as contiguous arrays and "restrict", I suppose
> the C compiler could do the most efficient thing in a lot of cases?
> (I.e. "cdef restrict int[:,:,"c"]" or similar)
A Fortran compiler can compile a vectorized expression like
a = b*c(i,:) + sin(k)
into
do j=1,n
a(j) = b(j)*c(i,j) + sin(k(j))
end do
The compiler do this because Fortran has strict rules prohibiting
aliasing, and because the instrinsic function sin is declared
'elemental'. On the other hand, if the expression contains functions not
declared 'elemental' or 'pure', there may be side effects, and temporary
copies must be made. The same could happen if the expression contained
variables declared 'pointer', in which case it could contain aliases.
I think in the case of numexpr, it assumes that NumPy ufuncs are
elemental like Fortran intrinsics.
Matlab's JIT compiler works with the assumption that arrays are
inherently immutable (everything has copy-on-write semantics). That
makes life easier.
> However if one has a strided array, numexpr could still give an
> advantage over such a loop. Or?
Fortran compilers often makes copies of strided arrays. The trick is to
make sure the working arrays fit in cache. Again, this is safe when the
expression only contains 'elemental' or 'pure' functions. Fortran also
often does "copy-in copy-out" if a function is called with a strided
array as argument.
S.M.
_______________________________________________
Numpy-discussion mailing list
Numpy-discussion@scipy.org
http://projects.scipy.org/mailman/listinfo/numpy-discussion
