![]() |
Blitz Support : |
From: Julian Cummings (cummings_at_[hidden])
Date: 2004-06-16 12:46:45
Derrick Bass reminded me that he had sent in a patch for a similar
problem with computing the proper ordering values that was occurring
when constructing an Array from an Array expression object. So I
borrowed generously from this patch and made a similar change for the
Array reduction object defined in <blitz/array/reduce.h>. This patch is
now checked into the blitz repository and it fully resolves the problem
with Patrik's tensor expressions.
-- Julian C.
Julian C. Cummings wrote:
>After thinking about this a bit more tonight, I believe that the ordering
>computation for expressions involving a reduction needs to be corrected. We
>need to add some checks that (a) the ordering value is less than the reduced
>Array rank, and (b) the ordering value has not already appeared earlier in
>the ordering array. If either condition is violated, we take the ordering
>value for the next highest dimension. These fixes would eliminate the
>mismatches in ordering values evident in Patrik's examples.
>
>Regards, 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
>>Julian C. Cummings
>>Sent: Monday, June 14, 2004 9:04 PM
>>To: 'Support list for Blitz++'
>>Subject: RE: [Blitz-support] tensor expression weirdness
>>
>>
>>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
>>>
>>>
>>>
>>_______________________________________________
>>Blitz-support mailing list
>>Blitz-support_at_[hidden]
>>http://www.oonumerics.org/mailman/listinfo.cgi/blitz-support
>>
>>
>>
>
>
>_______________________________________________
>Blitz-support mailing list
>Blitz-support_at_[hidden]
>http://www.oonumerics.org/mailman/listinfo.cgi/blitz-support
>
>
-- Dr. Julian C. Cummings E-mail: cummings_at_[hidden] California Institute of Technology Phone: 626-395-2543 1200 E. California Blvd., Mail Code 158-79 Fax: 626-584-5917 Pasadena, CA 91125