Fortran 77 common blocks in multithreading C++ application

1.4k Views Asked by At

I develop one C++ program which calls a Fortran 77 routine. The main C++ program may run multithreaded. However, it happens that the Fortran 77 routine hides several common blocks which are modified on each call depending on its arguments.

I am afraid that all common blocks may be shared between multiple threads and that concurrent accesses to these blocks will probably mess everything.

  • First question : Am I right? Would common blocks be shared between multiple threads?

  • Second question : Is there a simple way to avoid it? Rewriting the Fortran routines seems unaffordable, I am rather looking for a way so that each thread has its own copy of all common blocks (which are not large, should be fast to copy). I do not know if a compiling option would help or if OpenMP could help me.

4

There are 4 best solutions below

0
On BEST ANSWER

Yes, common blocks are shared.

In OpenMP it is possible to specify a common block as THREADPRIVATE. Each thread than makes a new instance of the common block dynamically. To copy the data from the original one use the COPYIN specifier. See also Difference between OpenMP threadprivate and private

The basic syntax is

!$OMP THREADPRIVATE (/cb/, ...)  

where cb is the name of a common block. See https://computing.llnl.gov/tutorials/openMP/#THREADPRIVATE

2
On

Yes, you can't use common areas with multithreading. And no, there is no way to avoid this. All common areas are actually conflated by the linker into single block, and there is no way to copy it between threads. It is a known problem everywhere where legacy Fortran code exists. Most common solution is to use multiprocessing instead of multithreading.

0
On

You are correct that common blocks are not threadsafe. They are global data that lets you declare variables in any scoping unit that all share the same storage association. The effect is essentially the same if you were writing to global variables in C++ with all of the thread synchronization issues that would cause.

Unfortunately, I don't think there is a simple way to avoid it. If you need to maintain a multi-threaded approach, one idea I've seen thrown around in the past is to move all of the variables from a common block into a user defined type and to pass instances of that type to any procedure needing access to them (one instance per thread). This would involve potentially expensive changes to the code to implement though.

You would also need to look at other thread safety issues with the Fortran code (this is not an exhaustive list):

  • IO units should be unique per thread, otherwise file input/output would not be reliable
  • Any variables with the SAVE attribute (implicit in module variables and in variables initialized when declared) are problematic (these variables are persistent between procedure calls). The implitness of this attribute is also compiler/standard dependent, making this an even bigger potential issue.
  • declare procedures with the RECURSIVE attribute -- this implies the function is re-entrant. This can also be satisfied by compiling with your compilers openmp option rather than changing code.

Another route you could explore is to use multi-processing or message-passing to parallelize your code rather than mutli-threading. This avoids the thread-safety problems with your Fortran code, but presents another potentially expensive code architecture change.

Also see:

0
On

Thanks for your answeers, especially the hint about OpenMP, it is indeed doable. I made a small program in order to be completely sure. It consists of one fortran 77 part called in one main C++ program (which is my concern) :

the fortran 77 routines func.f :

  subroutine set(ii, jj)
  implicit none

  include "func.inc"
  integer ii, jj
  integer OMP_GET_NUM_THREADS, OMP_GET_THREAD_NUM

  i = ii + 1
  j = jj

  !$OMP CRITICAL
  print *, OMP_GET_THREAD_NUM(), OMP_GET_NUM_THREADS(), i, j
  !$OMP END CRITICAL
  return
  end


  subroutine func(n, v)
  implicit none

  include "func.inc"

  integer n, k
  integer v(n)

  do k = i, j
     a = k + 1
     b = a * a
     c = k - 1
     v(k) = b - c * c
  enddo

  return
  end

with the include file func.inc

  integer i, j
  integer a, b, c

  common /mycom1/ i, j
  !$OMP THREADPRIVATE(/mycom1/)
  common /mycom2/ a, b, c
  !$OMP THREADPRIVATE(/mycom2/)

and finally the C++ program main.cpp :

#include<iostream>
#include<sstream>
#include<vector>
using namespace std;

#include<omp.h>

extern "C"
{
  void set_(int*, int*);
  void func_(int*, int*);
};


int main(int argc, char *argv[])
{
  int nthread;
  {
    istringstream iss(argv[1]);
    iss >> nthread;
  }

  int n;
  {
    istringstream iss(argv[2]);
    iss >> n;
  }

  vector<int> a(n, -1);

#pragma omp parallel num_threads(nthread) shared(a)
  {
    const int this_thread = omp_get_thread_num();
    const int num_threads = omp_get_num_threads();

    const int m = n / num_threads;
    int start = m * this_thread;
    int end = start + m;

    const int p = n % num_threads;
    for (int i = 0; i < this_thread; ++i)
      if (p > i) start++;
    for (int i = 0; i <= this_thread; ++i)
      if (p > i) end++;

#pragma omp critical
    {
      cout << "#t " << this_thread << " : [" << start
           << ", " << end << "[" << endl;
    }

    set_(&start, &end);
    func_(&n, a.data());
  }

  cout << "[ " << a[0];
  for (int i = 1; i < n; ++i)
    cout << ", " << a[i];
  cout << "]" << endl;

  ostringstream oss;
  for (int i = 1; i < n; ++i)
    if ((a[i] - a[i - 1]) != int(4))
      oss << i << " ";

  if (! oss.str().empty())
    cout << "<<!!  Error occured at index " << oss.str()
         << " !!>>" << endl;

  return 0;
}
  • Compilation steps (gcc version 4.8.1) :

    gfortran -c func.f -fopenmp
    g++ -c main.cpp  -std=gnu++11 -fopenmp
    g++ -o test main.o func.o -lgfortran -fopenmp
    
  • You can launch it as follows :

    ./test 10 1000
    

    where

    • the first integer (10) is the number of threads you want,
    • the second one (1000) is the length of one vector.

    The purpose of this program is to split this vector between threads and to let each thread of fill one portion of it.

    The filling of the vector is made within fortran 77 :

    • set routine first sets the lower and upper bound managed by the thread,
    • func routine then fills the vector between previous bounds.

Normally, if there are no errors and if common fortran 77 blocks are not shared, the final vector should be filled with 4 * k values, k going from 1 to 1000.

I could not catch the program out. Conversely, if I remove fortran 77 OMP directives in func.inc, then common blocks are no more private and lots of errors arise.

So to conclude, the only thing I need to do to solve my initial problem is to add OMP directives just behind any common blocks, which hopefully is not to complicated as they are all gathered in one include file (like my test).

Hopes this helps.

Best regards.