(j3.2006) Would there be a technical problem if...

[Damian Rouson]

Bill,

Yes, of course the type and operator would be implemented in a module or at least the procedure interfaces would be in a module.  


I would share your concerns if the example were exactly as I wrote it. I was attempting to be brief and apparently erred on the side of not providing enough context and detail so I?ll have to be much more long-winded to explain. ?I generally recommend against such tiny data structures as the one that might have been assumed from my brief email. ?In fact, my very first publication on scientific programming focused largely on issues related to data structure granularity. :) ?I probably should have explained that the derived type name ?vector_field? was chosen to imply that there?s a lot of data inside each element. ?For example, u(1) might have an array component that stores all of the nodal values of a field quantity sampled on a grid at many millions or billions of points. ?In that context, any overhead incurred in the type-bound procedure invocations is swamped by the computation that takes place inside the procedure.?

Moreover, Jim Xia also worked out a nice way to write such operators using coarray components so that the operators might be used to encapsulate parallel algorithms operating on distributed data. ?It works as long as the operators produce non-coarray results. ?In this context, the operators actually make it much easier to write high-performing, asynchronous code. The reason relates to Fortran?s requirement that operator arguments have INTENT(IN), which lends itself naturally to making the operators PURE. ?I make the point in lectures all the time that it is this requirement that makes this much more than just syntactic sugar.

With PURE operators, one can compose reasonably long and complicated expressions that execute asynchronously. ?For example, the right-hand side of the Navier-Stokes equations can be evaluated in a form

du_dt = - rho * (.grad. p) + nu*(.laplacian. u) - u .dot. (.grad. u)

wherein the RHS might represent parallel operators with distributed-memory operands and the PURE attribute guarantees that even those operators that require communication (e.g., .grad. and .laplacian.) can be executed asynchronously. ?In this regard, the notation actually aids performance. ?And the gain in clarity is priceless. ?If you took most students who had just taken their first fluid dynamics course and showed them most CFD codes, they would be hard pressed to find anything on the screen that is recognizable. By contrast, students who don?t even know Fortran look at the above code and know precisely what physics is being represented. ?

Jim wrote the first example like the one above for our book at a time when no compiler on the planet could handle it because of the coarray components. ?(Cray came the closest but had a bug.) ?That code, a Burgers PDE solver, was the only code in our book that we weren?t able to test before publication. ?When I finally got my hands on a compiler that could compile the code (I think it was Cray but might have been Intel), it compiled and ran correctly the first time without modification and scaled to 16,384 cores with nearly 90% parallel efficiency in weak scaling. ?Try writing a parallel PDE solver with any more traditional program structure and tell me if you can get it to compile, run correctly, and scale to 16K cores on your first attempt. The code in Chapter 12 of our book did exactly that.

Q.E.D.

Damian

On January 17, 2018 at 7:11:56 AM, Bill Long (longb at cray.com) wrote:

I don?t find this example persuasive. Firstly, the obvious way to do something like this is to define the operator generic in the module where the type is defined and supply the implementation function in the same module. I don?t see how this example has anything to do with OOP or type-bound procedures. For simple computations, tbp are almost always the wrong choice for performance. I hope people are not being taught to do things like this. Finally, I would expect an operator named .dot. to have a scalar result.  

Cheers,  
Bill  

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.j3-fortran.org/pipermail/j3/attachments/20180117/e4f037f9/attachment-0001.html>