[MITgcm-devel] vectorization of layers

[Ryan Abernathey]

I wrote this code. What is going on here is the most important and
expensive part of the LAYERS package: the sorting of each point in the
water column into the appropriate density class. This is essentially a
search problem (it has to hunt for the correct bin), and I based the code
on one of the search algorithms from Numerical Recipes. I did some basic
testing and optimization of the code in MATLAB before putting it into the
model, but I am sure there is lots of room for improvement.

The "hunting" takes place inside the while loops. kgu is an index for
layers_bounds that tracks which target density bin is being tested. If the
bin is too cold for the current point (whose density value is TatU), kgu
is augmented. The reason for using a while-loop is that it only runs as
long as it has to until the search condition is met. This seemed like the
fastest thing at the time.

This search has to run on each point in the refined vertical grid (Nr *
FineGridFact), so any improvement to its execution time would really
multiply.

I have only a vague understanding of what vectorization is, but I would be
happy to work with someone who understands such things better in order to
make the code more efficient.

-Ryan


On Fri, Jan 4, 2013 at 10:58 AM, Martin Losch <Martin.Losch at awi.de> wrote:

> Hi there,
> I am not sure how much effort one still should invest into vectorizing the
> code, since our vector machine will only run for another year or so, but
> can you think of a good way of rewrite layers_fluxcalc.F so that it
> vectorizes? Currently these lines (and the corresponding for TatV):
>
>          ELSE IF (TatU .GE. layers_bounds(kgu(i,j),iLa)) THEN
> C        have to hunt for the right bin by getting hotter
>           DO WHILE (TatU .GE. layers_bounds(kgu(i,j)+1,iLa))
>            kgu(i,j) = kgu(i,j) + 1
>           ENDDO
> C         now layers_bounds(kgu(i,j)+1,iLa) < TatU <=
> layers_bounds(kgu(i,j)+1,iLa)
>          ELSE IF (TatU .LT. layers_bounds(kgu(i,j)+1,iLa)) THEN
> C        have to hunt for the right bin by getting colder
>           DO WHILE (TatU .LT. layers_bounds(kgu(i,j),iLa))
>            kgu(i,j) = kgu(i,j) - 1
>           ENDDO
> C         now layers_bounds(kgu(i,j)+1,iLa) <= TatU <
> layers_bounds(kgu(i,j)+1,iLa)
>          ELSE
>
> inhibit vectorization (because the while-loops don't vectorize), and thus
> making it impossible to use the layers-pkg on a vector computer (This
> routines has a vector operation ratio of 5% and takes 86% of the computer
> time). Currently I have no clue what to do here, but I don't quite
> understand, what's going on in thie package/routine.
>
> Martin
> _______________________________________________
> MITgcm-devel mailing list
> MITgcm-devel at mitgcm.org
> http://mitgcm.org/mailman/listinfo/mitgcm-devel
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mitgcm.org/pipermail/mitgcm-devel/attachments/20130104/677afc84/attachment.htm>