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.