|
|
Blitz Support : |
From: Julian C. Cummings (cummings_at_[hidden])
Date: 2004-06-14 23:03:42
Hi Patrik,
The ordering array specifies the order in which a blitz array is stored in
memory. Thus, ordering[0] gives the rank of the array with unit stride in
memory and ordering[N_rank-1] gives the rank with the largest stride. What
is happening is that you are mixing two arrays with different storage
orders, which is not allowed when using tensor notation.
I think what is wrong here is that the ordering is not correctly computed
for an array subexpression containing a reduction. According to the code in
<blitz/array/reduce.h>, the ordering for a reduction is the same as the
ordering for the initial expression being reduced, with the rank reduced by
one. It is not obvious to me why this should be the case. The behavior of
the array ordering is also somewhat convoluted for arrays that are indexed
with tensor notation. This is handled by the code in <blitz/array/map.h>,
which indicates that the ordering for a given dimension d is equal to the
ordering of the underlying array for the dimension r that has been mapped to
d, as long as r is less than the rank of the underlying array.
A couple of simple examples may help. By default, blitz arrays have C-style
ordering, so the ordering array will contain [N_rank-1, N_rank-2,...,1,0].
The ordering is unchanged if we index an array with tensor indices in the
normal ascending ordering. Thus, the ordering of A and
A(tensor::i,tensor::j) are the same for a 2d array A. If we reverse the
tensor indices, the ordering will be reversed. If we use a tensor index
beyond the rank of the array, the ordering will be set to a flag value
INT_MIN for each index that is skipped over. So if B is a 3d array, the
ordering of B(i,j,l) is [2,1,INT_MIN,INT_MIN].
When combining arrays in an expression, the ordering of each array must
match for each rank, but the INT_MIN flag value is ignored.
In your first expression, the ordering for the summed expression is
[2,1,1,0], which then reduces to [2,1,1] after the summation. The ordering
for the denominator will be [INT_MIN,INT_MIN,0], so it does not match in the
final dimension. In your second expression, the ordering for the summed
expression is [1,0], which reduces to [1] after summation. This obviously
doesn't match up with the ordering of the denominator, which must be [0].
In neither case does the ordering of the denominator match the ordering of
the numerator.
So why does your first expression work and your second expression fail? The
first one uses the assignment operator while the second one uses an array
constructor. The array constructor sets the ordering of the new array to
that of the expression, so it has to compute this ordering, and that is
where the mismatch is detected. The assignment operator simply assumes that
the array on the left-hand side is filled in its own default ordering
scheme, so the ordering of the right-hand side is not checked.
I'm not sure if the behavior of the ordering is incorrect for summations or
not. If you change your second expression to use assignment, it should work
OK but I would check the results carefully. When in doubt, store the
summation in an array, and then perform the division separately and make
sure your results do not change.
Hope this helps,
Julian C.
Dr. Julian C. Cummings
Staff Scientist, CACR/Caltech
(626) 395-2543
cummings_at_[hidden]
> -----Original Message-----
> From: blitz-support-bounces_at_[hidden]
> [mailto:blitz-support-bounces_at_[hidden]] On Behalf Of Patrik
> Sent: Sunday, June 13, 2004 1:48 AM
> To: Support list for Blitz++
> Subject: [Blitz-support] tensor expression weirdness
>
>
> hi all,
> i'm using tensor expressions to do some reductions:
>
> broadband_flux =
> (sum (image (tensor::i, tensor::j, tensor::l)*
> passbands (tensor::k, tensor::l)*delta_lambda (tensor::l),
> tensor::l)/ewidth (tensor::k));
>
> Array<double,1> broadband_L_integrated
> (sum (L_lambda_integrated (tensor::j)*
> passbands_integrated (tensor::i, tensor::j)*
> delta_lambda_integrated (tensor::j), tensor::j)/
> ewidth_integrated(tensor::i));
>
> The first one works, while the second one gives the following
> error (it's
> the division that's the offending part, the sum by itself works):
>
> /u7/patrik/include/blitz/array/expr.h:276 Two array operands
> have different
> orders: for rank 0, the orders are 1 and 0
> Assertion failed: 0, file
> /u7/patrik/include/blitz/array/expr.h, line 279
>
> I don't even know what the error message means... The two
> expressions are
> very similar (the only difference is that the first one has a
> 3d array and
> the second one has a 1d one, so I don't understand why one
> works but not
> the other. Any hints?
>
> Regards,
>
> /Patrik
>
> _______________________________________________
> Blitz-support mailing list
> Blitz-support_at_[hidden]
> http://www.oonumerics.org/mailman/listinfo.cgi/blitz-support
>