How to solve the same numerical Problem in 7 different Programming Languages

Most of the Fortran discussion is spent explaining how weird implicit none is.

4 Likes

I find it weirder to use subroutines instead of a functions.

1 Like

Well, he has a point. It is weird to have such an old feature sticking around for such a long time. IMHO, implicit none should be the default, at least for free form source code.

9 Likes

Totally agree.

It is not a matter of simplicity. It is a matter of marketing, for situations such like that in the link.

The genuine disappointment for me is that people dare to advise others about topics they barely know (here, Fortran) because, in their perceived (potentially biased) view of the world, they are the epicenter of justice, fairness, and the peak of the mountain of knowledge with no room for any doubts. This problem goes beyond a simple language comparison and is endemic to all inexperienced. Also, this article seems to be a Julia marketing post.

13 Likes

Ah, yes. The magazine article author who is “often wrong but never uncertain”. After running the Python,…,Perl, C++, Fortran versions, goes on to say, " I wasn’t satisfied in seeing such old languages win. There has to be some progress made in language design. So I took a second look at the Julia version, as it was the one that showed the biggest potential. Using all the tricks the C++ guys showed me, I settled on this optimized Julia version."

Let us also note that the author is a graduate student and the article is one way of “having some fun” as @macneacail sometimes does.

3 Likes

I must be misreading the problem definition, but where is N or dt defined ?
A function is much better and the initial condition can be provided as u(1).
Also no discussion of accuracy with first derivitive extrapolation ? (or is that part of the solution for other languages ?)
An alternative is to also use the second derivitive, with a larger dt
I needed N=500 for stability and N=125 with second derivitive, assuming I interpreted the rest of the problem description.

The analytical solution is in terms of the tanh function.

When reading such articles, I wonder if more people are influenced by the author’s opinions

nobody should ever use FORTRAN outside of existing legacy projects

or the actual code. Looking at the C++ code,

#include <vector>

static double f(double u){
    return u*(1.0-u);
}

void Cplusplus_logistic(int N, vector<double>& u, double dt){
    // Parameters
    int T = 25;
    double u0 = 1e-5;
    // Discretization
    u[0] = u0;
    for(auto i = u.begin(); i != u.end(); ++i){
        *(i+1) = *i + dt*f(*i);
    }
    return (u);    
}

it took me a while to understand that

*(i+1) = *i + dt*f(*i);

sets a value of array element u[i+1].

I think the Fortran code is clearer, and it could be improved, as discussed above.

subroutine f(u, r)
implicit none
  real(kind(0.d0)), intent(out) :: r
  real(kind(0.d0)), intent(in) :: u
  r = u*(1.d0-u)
end subroutine f

subroutine  fortran_logistic(N, dt, u)
  integer, intent(in) :: N
  real(kind(0.d0)), intent(in) :: dt
  real(kind(0.d0)), intent(out) :: u(N)
  
  integer :: i
  real(kind(0.d0)) :: temp
  
  u(1) = 1.d-5
  do i=1, N-1
    call f(u(i), temp)
    u(i+1) = u(i) + dt*temp
  end do
end subroutine  fortran_logistic

I just purchased the book High Performance Python: Practical Performant Programming for Humans, 2nd Ed. by Gorelick and Ozsvald. Along with using NumPy, Cython, and calling C, they do consider using Fortran with f2py for speed.

3 Likes

I tried to improve his code
Compiler: ifort
Flags: -O3

  1. I didn’t understand why he needed two subroutines
  2. He also adds some overhead with the subroutine call inside a loop
  3. I added some comments and spaces to make the code readable
  4. I also graph the data at the end
    If possible, how can I improve his code further ?
program simple_num_prob
  implicit none

  integer :: i
  integer,parameter :: tf=25,n=10000
  real,allocatable :: t(:),u(:)
  real :: dt
  real(kind=8) :: t1,t2
  character(len=*),parameter :: plt_file='simple_num_prob.gp'

  call cpu_time(t1)
  allocate(t(n))
  allocate(u(n))

  ! Assigning step value
  dt=real(tf)/real(n)


  ! Initial u value assigned
  u(1)=1.0E-5


  ! Population of t vector
  do i=1,n
     t(i)=real(i)
  end do

  ! Calculating and populating the function vector
  do i=1,n-1
     u(i+1) = u(i) + (dt*(u(i)*(1.0-u(i))))
  end do

  call cpu_time(t2)

  print *, "Elapsed time",t2-t1

  !writing to file
  open (unit=1,file='simple_num_prob.dat')
  do i=1,n
     write(1,*) t(i),' ',u(i)
  end do
  close(1)
  print *, "File written !"

  call cpu_time(t2)
  print *, "Elapsed time with writing to file",t2-t1

  !deallocation
  deallocate(t)
  deallocate(u)

  call execute_command_line('gnuplot ' // plt_file)

end program simple_num_prob

Gnuplot file:

plot 'simple_num_prob.dat' using 1:2
 pause -1

2 Likes

Me too! I don’t see any value being added in articles by making such broad generalisations. Especially when they are followed by performance metrics that prove that very sentence false.

Also, I chuckled when I read that quote. Someone should probably let the rest of the world know…

1 Like

What’s a void function doing with a return statement?

void Cplusplus_logistic(int N, vector<double>& u, double dt){
...
    return (u);    
}

I find it weird to use functions all the time, as C does. Both subroutines and functions have their uses. The concept of a void “function” is not only misleading, but plain wrong from a Mathematical point of view, if you ask me. Furthermore, a function that returns two or more separate entities (by its name and by altering its arguments) is also weird - but still allowed in Fortran. If a “function” has to return multiple things, either wrap them in a derived type and return just that, or make it a subroutine. It’s way more clear and easily understood. Compilers supporting the F-standard actually issue a syntax error if a function has any intent(out) or intent(inout) arguments, recommending using a subroutine instead.

Totally agree. Again, the F-standard makes implicit none obsolete since it is implied by default (it is a pity popular compilers lack -std=F.) However I have to tell people still bashing Fortran because of implicit none that backwards compatibility is the reason it exists, and I could find way weirder “features” in the programming languages they praise 24/7, be it C, C++ or, even worse, Python, Java, C#, whatever…

The quoted statement is correct if and only if the author means FORTRAN and not Fortran - but I doubt he does. Most people bashing Fortran aren’t even aware of the difference, but that doesn’t stop them from having a strong opinion and hate anyway. I really hope that article would influence as less people as possible - preferably nobody.

2 Likes

Totally agree with you.

Just to add my 1 and 1/2 cents, I understand void functions in C/C++ as the analog of fortran subroutines with no intent(out) arguments. Otherwise it also seems weird to me.

1 Like

Someone give me a double “like” please. Oh wait, make it triple like. Perhaps quad like. Ok, 100 like. The quoted post definitely deserves it. :laughing:
The author doesn’t hide under the carpet the fact “FORTRAN” is fast though… you have to grant him that. Although his sorrow for admitting Fortran is fast is more than obvious.

1 Like

Totally agree, both concepts are useful. My comment was directly related to the original code from Fortran logistic · GitHub

subroutine f(u, r)
implicit none
  real(kind(0.d0)), intent(out) :: r
  real(kind(0.d0)), intent(in) :: u
  r = u*(1.d0-u)
end subroutine f

subroutine  fortran_logistic(N, dt, u)
  integer, intent(in) :: N
  real(kind(0.d0)), intent(in) :: dt
  real(kind(0.d0)), intent(out) :: u(N)
  
  integer :: i
  real(kind(0.d0)) :: temp
  
  u(1) = 1.d-5
  do i=1, N-1
    call f(u(i), temp)
    u(i+1) = u(i) + dt*temp
  end do
end subroutine  fortran_logistic

where for me the subroutine f should be a function:

function f(u)
  implicit none
  real(kind(0.d0)) :: f
  real(kind(0.d0)), intent(in) :: u
  f = u*(1.d0-u)
end function f

subroutine  fortran_logistic(N, dt, u)
  integer, intent(in) :: N
  real(kind(0.d0)), intent(in) :: dt
  real(kind(0.d0)), intent(out) :: u(N)
  
  integer :: i
  
  u(1) = 1.d-5
  do i=1, N-1
    u(i+1) = u(i) + dt*f(u(i))
  end do
end subroutine  fortran_logistic

This reads IMHO much better and is not too different from the Python variant.
I would even make fortran_logistic a function here because it has a single return value.

3 Likes

I agree. But it does not have to appear in functions/subroutines because a single implicit none at the module level is sufficient.

@Carltoffel , @Pap , @conradoat , @jacobwilliams , and all other readers who would like to see implicit none be the default in a future revision of the Fortran standard:

Please see this thread at GitHub J3-Fortran site for language proposals:

If you agree with suggestion, will you please consider upvoting it?

You will know things change slowly, perhaps even more so than the wheels of justice (!!), in the “world of Fortran”. Should one begin chartering a new direction now, one will see begin sighting the light at the end of the tunnel only possible a decade plus from now.

Hence I personally sense a pressing need and urgency to seek support from the wider community and this is where your upvote can help gain traction and influence on a simple change like this that can commence with Fortran 202Y.

7 Likes

@MarDie , but the current situation is pernicious with INTERFACE blocks that are effectively stand-alone constructs uninfluenced by a host and its IMPLICIT mapping - please see the link in my prior post:

interface
   function foo(..) result(..)
      [import .. ]
      implicit none !<-- Per current standard, forget this and face peril
      ..
   end function
   function bar(..) result(..)
      [import .. ]
      implicit none !<-- Per current standard, forget this and face peril
      ..
   end function
   subroutine foobar(..)
      [import .. ]
      implicit none !<-- Per current standard, forget this and face peril
      ..
   end subroutine
end interface
2 Likes

Yes this is just terrible. Why oh why doesn’t the implicit none in the module apply to the interface blocks?