[Date Prev][Date Next] [Thread Prev][Thread Next] [Date Index] [Thread Index]

Re: Is Linux Unix?



Erik Steffl <steffl@bigfoot.com> wrote in message news:<2lnqj-1oX-1@gated-at.bofh.it>...
> Ryo Furue wrote:
[...]
> > That's out of question! :)  For numerical computation, Fortran 95 is "better"
> > than C++. (I have plenty of experience in C++ programming.)  Sure I greatly
> > miss C++'s objects and templates and other goodies.  But codes for numerical
> > computation can be much more "easily" written in Fortran 95.  That's the most
> > important point for a non-proffessional programmer like me.
> 
>    I find that hard to believe. is it because of fortran libs that are 
> not available for c/c++ or is performance much better or does fortran 
> have some features that make it easy to write numerical computation code 
> (I forgot most of what I had to learn about fortran but I don't remember 
> ANYTHING being easy/readable in fortran).

First of all, thanks for asking this question.  Second, let us be clear
about which "Fortran" we are talking about.  *Many* people (I don't know
if you are included) have very bad impression on Fortran because all they
know about it is Fortran 77 or 66.  A major revision, called Fortran 90,
and its minor revision called Fortran 95, are a very different language.
Fortran 90/95 is a superset of Fortran 77, but the addition over Fortran
77 is a huge improvement.  Fortran 90/95 still lacks "modern" features
such as object-oriented programming support, parametric polymorphism,
functional programming support, etc.  But it's "usable" to me, whereas
I would call Fortran 77 unusable.

Third, and most important for answering your question, I use Fortran 95
mostly because I can write a certain type of code easily. I think this
example is convincing:

   forall (i=1:nx)
      where (a(i,:,:) /= b(i))
         a(i,:,:) = cos(c(:,:) + e(:,:));
      else where
         a(i,:,:) = sum(d(:,:,:), dim=3) / count(mask(:,:))
      end where
   end forall

Fortran 95 has lots of convenient intrinsic capabilities to handle arrays.
And I use arrays a *lot* in numerical codes.   The ":" notation means
"all indices".  In a simple-minded C++, the above code is roughly equivalent
to

    double s[ny][nz];

    for (int j=0; j < ny; ++j) {
       for (int k=0; k < nz; ++k) {
          s[j][k] = 0;
          for (int m=0; m < nt; ++m) {
             s[j][k] = s[j][k] + d[j][k][m];
          }
       }
    }

    int c = 0;
    for (int i=0; i < nx; ++i) {
       for (int j=0; j < ny; ++j) {
          if (mask[i][j]) ++c;
       }
    }

    for (int i=0; i < nx; ++i) {
      for (int j=0; j < ny; ++j) {
         for (int k=0; k < nk; ++k) {
            if (a[i][j][k] != b[i]) {
               a[i][j][k] = cos(c[j][k] + e[j][k]);
            } else {
               a[i][j][k] = s / c;
            }
         }
      }
   }

You can see that the Fortran version is easier to write and read.  You may
point out that in C++ version, the computation of s and c is outside of the
main loop, so C++ version is more efficient.  But, any decent Fortran compiler
will move loop-constants out of the loop.  You don't need to worry about it
too much.  There's another problem with the above C++ code:  If ny and nz
aren't constant, you can't write

    double s[ny][nz];

Instead, either you allocate the array in two stages:

    double** s = new (double*)[ny];
    for (j=0; j < ny; ++j) s[j] = new double[nz];

or you write a multidimensional array class, something like:

    template<class T>
    class ArrayTwoDim {
       T* a_;  int nx_;
    public:
       ArrayTwoDim(int nx, int ny) { a_ = new T[nx*ny]; nx_ = nx; }
       ~ArrayTwoDim() { delete[] a_; }
      // operator[]
      // operator+, operator-, operator*, operator/,
    };
    // Also, define sin(), cos(), log(), etc. for a whole array.

The first is cumbersome (Don't forget to deallocate those arrays!) and the
second involves signigicant work.  The dynamic allocation of multidimensional
arrays is a piece of cake in Fortran,

   real(8):: a(nx,ny)  !! nx or ny need not be constant.
                       !! a is allocated on the stack

or

   real(8), allocatable:: a
   . . . .
   allocate(a(nx,ny))  !! allcated on the heap.

and the deallocation is automatic at the end of the subroutine (or you
can deallocate for yourself if you want an earlier deallocation).

But, you might say, once we've written multidimensional array classes, we
can happily use them.  No.  Unless the array classes are part of the C++
standard, you'll have problems when you share programs with other people.
And, such classes will be still significantly less powerful than the Fortran
counterpart.  For example, sections of array are directly supported by the
core language:

    a(:, 3:10, :) = c(1:7, :, :) + sin(d(1:7, :, :))

whereas in C++, either you write

    for (int i=0; i < nx; ++i) {   // code (1)
       for (int j=2; j < 10; ++j) {
          for (int k=0; k < nz; ++k) {
             a[i][j][k] = c[j-2][i][k] + sin(c[j-2][i][k]);
          }
       }
    }

or you define a Section method for your array class and say something like

      a.Section(0,nx,  2,10,  0,nz)  // code (2)
         = c.Section(0,7,  0,nx,  0,nz) + sin(d.Section(0,7, 0,nx, 0,nz));

(or you write a proxy class to represent sections).  Either of these is far less
readable than the Fortran counterpart.  Also, will code (2) be optimzed to an
eauivalent of code (1)?

This convenince in handling multidimensional arrays alone is reason enough to
use Fortran, for me.  For other applications where arrays play a less significant
role, I would use other languages.

Also, *It is said* that Fortran code tends to be faster.  I emphasize
"it is said" because I haven't myself compared speeds of Fortran and C++
programs. All I can point out is that the recent C standard introduced a
keyword which says arrays pointed to by pointers don't overlap.  When arrays
overlap, a significant type of optimization can't be applied.  And in C/C++,
the compiler can't always tell whether arrays overlap or not.  The purpose of
the new C keyword is specifically to allow this type of optimization, which
Fortran compilers have been doing since (or before?) the days of Fortran 77.
Why does C introduced such a keyword?  Because it wanted to match the
performance of Fortran, so I heard. (See
http://www.lysator.liu.se/c/restrict.html for a discussion of this new keyword.
This page does mention Fortran.) I admit I heard these things from people
biased toward Fortran. You may want to verify them for yourselves.
By the way, C++ doesn't have the keyword.

Finally, overall I prefer C++ to Fortran 95.  C++ is generally more powerful
than Fortran 95.  But, for numerical computations, no.

Cheers,
Ryo



Reply to: