using maxloc 
Author Message
 using maxloc

Could someone tell me how to correct the following program so that it
will behave as desired?

program using_maxloc
   implicit none

   real :: A(3,2), B(3,2)

   A(1,1) = 1.0
   A(2,1) = 7.0
   A(3,1) = 8.0

   A(1,2) = 4.0
   A(2,2) = 9.0
   A(3,2) = 5.0

   call random_number(B)

   write(*,*) maxloc(A(:,:))
   ! Is there any efficient way to correct the following statement
   ! so that it will actually print the the value of A(2,2)/B(2,2)
   !   where (2,2) is the location of maximum of A
   write(*,*) A(maxloc(A(:,:)))/B(maxloc(A(:,:)))
end program using_maxloc

$gfortran using_maxloc.f90
  In file using_maxloc.f90:20

   write(*,*) A(maxloc(A(:,:)))/B(maxloc(A(:,:)))
                                1
Error: Rank mismatch in array reference at (1) (1/2)

As far as I understand, maxloc returns an array. But this cannot be
passed as array subscripts. Any other workarounds?

thanks
raju

--
Kamaraju S Kusumanchi
http://www.*-*-*.com/
http://www.*-*-*.com/



Tue, 22 Jul 2008 10:30:37 GMT  
 using maxloc

Quote:

> Could someone tell me how to correct the following program so that it
> will behave as desired?

> program using_maxloc
>   implicit none

>   real :: A(3,2), B(3,2)

>   A(1,1) = 1.0
>   A(2,1) = 7.0
>   A(3,1) = 8.0

>   A(1,2) = 4.0
>   A(2,2) = 9.0
>   A(3,2) = 5.0

>   call random_number(B)

>   write(*,*) maxloc(A(:,:))
>   ! Is there any efficient way to correct the following statement
>   ! so that it will actually print the the value of A(2,2)/B(2,2)
>   !   where (2,2) is the location of maximum of A
>   write(*,*) A(maxloc(A(:,:)))/B(maxloc(A(:,:)))
> end program using_maxloc

> $gfortran using_maxloc.f90
>  In file using_maxloc.f90:20

>   write(*,*) A(maxloc(A(:,:)))/B(maxloc(A(:,:)))
>                                1
> Error: Rank mismatch in array reference at (1) (1/2)

> As far as I understand, maxloc returns an array. But this cannot be
> passed as array subscripts. Any other workarounds?

> thanks
> raju

First, a suggestion: rather than writing A(:,:), you should write just A. On
some compilers, it may make a difference to performance, especially if A is large.

Now on to your question. I'm sure you'll get lots of suggestions on how to do
this, each with its own advantages. Some that spring to my mind are:

i = MAXLOC(A)
write(*,*) A(i(1),i(2))/B(i(1),i(2))

(where i is a default integer with DIMENSION(2)). If you don't want to use a
temporary, you might try this:

write(*,*) A(MAXLOC(A,DIM=1),MAXLOC(A,DIM=2))/B(MAXLOC(A,DIM=1),MAXLOC(A,DIM=2))

...but that's 4 calls to MAXLOC, which is {*filter*}.

cheers,

Rich



Tue, 22 Jul 2008 10:38:18 GMT  
 using maxloc

Quote:
> write(*,*)

A(MAXLOC(A,DIM=1),MAXLOC(A,DIM=2))/B(MAXLOC(A,DIM=1),MAXLOC(A,DIM=2))

Quote:
> ...but that's 4 calls to MAXLOC, which is {*filter*}.

And it gives the wrong answer, possibly even an array subscript
out of bounds for non-square A, which is even nastier.  How
about:

write(*,*) A(MAXLOC(MAXVAL(A,DIM=2)),MAXLOC(MAXVAL(A,DIM=1)))/ &
           B(MAXLOC(MAXVAL(A,DIM=2)),MAXLOC(MAXVAL(A,DIM=1)))

? Answer: No, because it could fail if two elements of A
have the same maximum value.  I guess I have to bite the
bullet with the more pedestrian:

write(*,*) A(DOT_PRODUCT(MAXVAL(A),(/1,0/)), &
             DOT_PRODUCT(MAXVAL(A),(/0,1/)))/ &
           B(DOT_PRODUCT(MAXVAL(A),(/1,0/)), &

             DOT_PRODUCT(MAXVAL(A),(/0,1/)))

--
write(*,*) transfer((/17.392111325966148d0,6.5794487871554595D-85, &
6.0134700243160014d-154/),(/'x'/)); end



Tue, 22 Jul 2008 12:14:18 GMT  
 using maxloc

Quote:


>> Could someone tell me how to correct the following program so that it
>> will behave as desired?

>> program using_maxloc
>>   implicit none

>>   real :: A(3,2), B(3,2)

>>   A(1,1) = 1.0
>>   A(2,1) = 7.0
>>   A(3,1) = 8.0

>>   A(1,2) = 4.0
>>   A(2,2) = 9.0
>>   A(3,2) = 5.0

>>   call random_number(B)

>>   write(*,*) maxloc(A(:,:))
>>   ! Is there any efficient way to correct the following statement
>>   ! so that it will actually print the the value of A(2,2)/B(2,2)
>>   !   where (2,2) is the location of maximum of A
>>   write(*,*) A(maxloc(A(:,:)))/B(maxloc(A(:,:)))
>> end program using_maxloc

>> $gfortran using_maxloc.f90
>>  In file using_maxloc.f90:20

>>   write(*,*) A(maxloc(A(:,:)))/B(maxloc(A(:,:)))
>>                                1
>> Error: Rank mismatch in array reference at (1) (1/2)

>> As far as I understand, maxloc returns an array. But this cannot be
>> passed as array subscripts. Any other workarounds?

>> thanks
>> raju

> First, a suggestion: rather than writing A(:,:), you should write just
> A. On some compilers, it may make a difference to performance,
> especially if A is large.

That is cool. I did not know that before. One book I read before
suggested to always write A(:,:) instead of A so that, the reader will
know for sure that it is a 2D array etc., From then on, I think I have
fallen into the habit of writing A(:,:).

Also, Thanks for all the suggestions on the actual problem.

raju

--
Kamaraju S Kusumanchi
http://www.people.cornell.edu/pages/kk288/
http://malayamaarutham.blogspot.com/



Wed, 23 Jul 2008 04:33:09 GMT  
 using maxloc

Quote:

> That is cool. I did not know that before. One book I read before
> suggested to always write A(:,:) instead of A so that, the reader will
> know for sure that it is a 2D array etc., From then on, I think I have
> fallen into the habit of writing A(:,:).

Hmm. I suspect I know what book that was. I don't like many of its other
recommendations either. :-(

I strongly disagree with the book's recommendation in this case. If you
start getting into th ehabit of writing A(:) as a synonym for A, you
will start introdusing bugs. They are *NOT* the same thing. One is a
named array. The other is a slice of that array. Although they have the
sameelements, their propertties are often different. For example:

1. A(:) is never allocatable. Nor is it ever a pointer. Even if A is.
 Thus, for example, with f2003 automatic allocation of allocatables on
assignment, A=something vs A(:)=something can do completely different
things. The A=something form might reallocate; the A(:) form never will.

2. The bounds may be different. The lower bound of A(:) is *ALWAYS* 1,
regardless of the bounds of A. This can cause errors in applications
like the one here, where you obtained index values for A(:) and then
applied them to A.

3. Various other differences.

--
Richard Maine                    | Good judgement comes from experience;
email: last name at domain . net | experience comes from bad judgement.
domain: summertriangle           |  -- Mark Twain



Wed, 23 Jul 2008 05:31:17 GMT  
 using maxloc

Quote:


>>That is cool. I did not know that before. One book I read before
>>suggested to always write A(:,:) instead of A so that, the reader will
>>know for sure that it is a 2D array etc., From then on, I think I have
>>fallen into the habit of writing A(:,:).

When I was first exposed to Fortran90 way back when, the first piece of advice I can
remember from an f90-experienced colleague was to use A(:) instead of just A so readers of
the code can tell that "A" is an array, not a scalar. I have since disabused myself of
that advice (and the subsequent practice I also adopted). Maybe the novelty of array
syntax in Fortran discombobulated people? :o)  To reiterate what Richard Maine said,
"A(:)" and "A" are not the same thing. Things will start working normally once this
particular perception anchor has been aweighed.

Quote:
> Hmm. I suspect I know what book that was. I don't like many of its other
> recommendations either. :-(

I went through the obvious parts of my copy of "that book" and I don't think it's the one
the OP was referring to. There was no recommendation (that I could find) to use things
like A(:,:) instead of just A -- and none of the code examples used it either. BUT, I did
come across the recommendation to "[n]ever use /assumed-size dummy arrays/. They are
undesirable and are likely to be eliminated from a future version of the Fortran
language." (pg315). Crikey.   :o)

cheers,

paulv

Quote:

> I strongly disagree with the book's recommendation in this case. If you
> start getting into th ehabit of writing A(:) as a synonym for A, you
> will start introdusing bugs. They are *NOT* the same thing. One is a
> named array. The other is a slice of that array. Although they have the
> sameelements, their propertties are often different. For example:

> 1. A(:) is never allocatable. Nor is it ever a pointer. Even if A is.
>  Thus, for example, with f2003 automatic allocation of allocatables on
> assignment, A=something vs A(:)=something can do completely different
> things. The A=something form might reallocate; the A(:) form never will.

> 2. The bounds may be different. The lower bound of A(:) is *ALWAYS* 1,
> regardless of the bounds of A. This can cause errors in applications
> like the one here, where you obtained index values for A(:) and then
> applied them to A.

> 3. Various other differences.

--
Paul van Delst



Wed, 23 Jul 2008 07:12:35 GMT  
 using maxloc
Thank you for all the information on using A against A(:).
But could you please talk more on why this is slow?
Also, does it matter if I specify the actual range of A, such as
A(1:N)?
Thank you very much.

Shi



Thu, 24 Jul 2008 05:44:02 GMT  
 using maxloc

Quote:

> Thank you for all the information on using A against A(:).
> But could you please talk more on why this is slow?
> Also, does it matter if I specify the actual range of A, such as
> A(1:N)?

There isn't any real reason the A(:) must be slower
than A, or that A(1:N) would be faster/slower, or
that any of them must even differ.  It's just an observation
that some implementations do better with one than another.

If you passed A with stride, so the actual argument is say
A(1:N:5), some implementations really have to make a copy
on the way into the procedure so that the data accessed in
the procedure was guaranteed contiguous.  Hence, some
implementation adopt a habit of copying the argument
whenever a slice notation is present at all.  Well, making the
copy costs time, and in the cases you mention (rank 1 array
with no stride) the copy *shouldn't* be necessary even if
the implementation requires contiguous arguments.  (It would
be desirable if compilers were better at making copies only
when they're required.)

Now, why would an implementation require contiguous
arguments?  Well, if the dummy argument is declared
as explicit shape or assumed size, the old rules of F77
apply, and it's really tough to have those work correctly
without contiguous arguments.  But, even for assumed
shape dummy arguments (which need explicit interfaces:
just to remind you), having contiguous arguments is often
faster since it maintains cache coherency better if you're
really working those arrays hard - this may save so much
time that copy-in/copy-out argument passing is still the
most efficient way to go.

--
J. Giles

"I conclude that there are two ways of constructing a software
design: One way is to make it so simple that there are obviously
no deficiencies and the other way is to make it so complicated
that there are no obvious deficiencies."   --  C. A. R. Hoare



Thu, 24 Jul 2008 06:58:48 GMT  
 using maxloc

Quote:

> Thank you for all the information on using A against A(:).
> But could you please talk more on why this is slow?

James has discussed that in his followup.

Quote:
> Also, does it matter if I specify the actual range of A, such as
> A(1:N)?

No. In fact, that means exactly the same thing as A(:) (assuming that A
is dimensioned 1 to N). It is still an array slice rather than the whole
array.

This is somewhat like, for example, writing x%y when y is the only
component of x - say a real component. It has all the same values, but
it is a different kind of thing. X is of derived type; x%y is of real
type. You can't use them interchangeably.

--
Richard Maine                    | Good judgement comes from experience;
email: last name at domain . net | experience comes from bad judgement.
domain: summertriangle           |  -- Mark Twain



Thu, 24 Jul 2008 08:23:42 GMT  
 using maxloc

Quote:

> > First, a suggestion: rather than writing A(:,:), you should write just
> > A. On some compilers, it may make a difference to performance,
> > especially if A is large.

> That is cool. I did not know that before. One book I read before
> suggested to always write A(:,:) instead of A so that, the reader will
> know for sure that it is a 2D array etc., From then on, I think I have
> fallen into the habit of writing A(:,:).

In most cases, I would also prefer to have it obvious that it is a 2D
array etc.  On what compilers does this actually make a difference, and
why?


Sat, 26 Jul 2008 03:35:41 GMT  
 using maxloc


Quote:


>>>First, a suggestion: rather than writing A(:,:), you should write just
>>>A. On some compilers, it may make a difference to performance,
>>>especially if A is large.

>>That is cool. I did not know that before. One book I read before
>>suggested to always write A(:,:) instead of A so that, the reader will
>>know for sure that it is a 2D array etc., From then on, I think I have
>>fallen into the habit of writing A(:,:).

> In most cases, I would also prefer to have it obvious that it is a 2D
> array etc.  On what compilers does this actually make a difference, and
> why?

It for sure makes a difference with allocatable arrays with
F2003 compilers.  Array assignment to a whole allocatable
array requires the compiler to detect size mismatches across
the = sign and deallocate the left and reallocate it with
the new correct size.  Assignment to a section is still
an error if the sizes don't match; there's no magic
de/reallocation.  Because new compilers will (sometimes)
have to treat
        array = different_sized_array
differently from
        array(:,:) = different_sized_array
it's probably a good idea to avoid the (:.:) notation
all of the time.

Dick Hendrickson



Sat, 26 Jul 2008 03:41:05 GMT  
 using maxloc
Phillip Helbig---remove CLOTHES to reply

on writing things like A(:,:) vs just A.

Quote:
> In most cases, I would also prefer to have it obvious that it is a 2D
> array etc.  On what compilers does this actually make a difference, and
> why?

See my previous post. This is not a matter of compiler-dependence. The
fact that one notation is a slice and the other is not is specified in
the standard. Any compiler that doesn't distinguish them in the
situations I mentioned is broken and will probably be fixed.

Specifically, the things I mentioned were differences in the lower
bounds and differences in attributes - notably allocatable and pointer.
If you start ignoring the allocatable attribute difference, you'll find
yourself in major trouble with allocatable string length, which I'm
guessing will be a widely used feature of f2003.

The performance differences are compiler-dependent matters. The
difference in meaning is not.

--
Richard Maine                     | Good judgment comes from experience;
email: my first.last at org.domain| experience comes from bad judgment.
org: nasa, domain: gov            |       -- Mark Twain



Sat, 26 Jul 2008 04:05:33 GMT  
 
 [ 12 post ] 

 Relevant Pages 

1. maxloc with mask

2. maxloc/minloc functions

3. MAXLOC 'problem'

4. Strange behavior of MAXLOC ???

5. maxloc returns??

6. Use of MAXLOC in CVF v6.0

7. maxloc

8. MAXLOC

9. maxloc()

10. Maxloc - multiple maxima

11. return value maxloc of an empty array!?

12. maxloc problem

 

 
Powered by phpBB® Forum Software