Use of OpenMP in filling the array elements

Question

Program test
    !$ use omp_lib
    
    implicit none
    !$ call omp_set_num_threads(4)
    
    integer :: i, j
    integer, parameter :: dimen = 5
    real (kind = 8), dimension(dimen, dimen) :: tstarr = 0.0

    !.....................................................................
    ! Later I will have to fill the elements of the following matrix in the same way
    ! real (kind = 8), dimension(:,:), allocatable :: tstarr2
    ! allocate(tstarr2(2 * dimen, 2 * dimen))
    !.....................................................................

    !$ omp parallel do default(shared) private(i,j)
    do i = 1, dimen
        tstarr(i,i) = 3
    end do
    
    do i = 1, dimen - 1
        tstarr(i + 1, i) = 2
        tstarr(i, i + 1) = 2
    end do

    do i = 1, dimen - 2
        tstarr(i + 2, i) = 4
        tstarr(i, i + 2) = 4
    end do

    do i = 1, dimen - 3
        tstarr(i + 3, i) = 5
        tstarr(i, i + 3) = 5
    end do
    !$ omp end parallel do
    do i = 1, dimen
        do j = 1, dimen
            print *, tstarr(i,j)
        end do
    end do


end program test

This code is working when compiling serially that means with gfortran test.f90 but giving errors when compiling with gfortran -fopenmp test.f90

It would very much help, as @francescalus says, if you edit the question to include the text of the errors you are encountering. But I'm guessing that it's misuse of `!$omp parallel do`. If you can confirm this by editing the question and showing the errors I (or somebody else) can write an answer, but IMO the first step to solving this is forgetting you ever knew about `parallel do` at all. — Ian Bush, Jul 17 '23 at 10:47
Welcome, please take the [tour] and read [ask]. Your question is really missing the errors you are getting. The error is important, because it helps to diagnose the problem and also because other people may encounter the same problem and may be able to find the answer to your question as a solution to their problem. As Ian mentions, I also suggest avoiding `parallel do` altogether until you are confident in using `omp parallel` and `omp do` independently. — Vladimir F Героям слава, Jul 17 '23 at 11:51

score 1 · Answer 1 · answered Jul 17 '23 at 14:09

1

The directive parallel do is used for a single do loop. Your compile error indicates that you can not have multiple loops in this directive.

So I would use an outer parallel and then do for each loop.

answered Jul 17 '23 at 14:09

Victor Eijkhout

5,088
2
22
23

score 1 · Answer 2 · answered Jul 18 '23 at 10:32

You have a number of problems with your code:

In Fortran you can't mix executable statements such as call omp_set_num_threads(4) with the variable declarations. Actually this is probably the cause of most of your (unreported) errors
Don't use kind = 8. It is not portable, not required to be supported by your compiler, and may not do what you expect. I suggest learning about Fortran kinds, and see the code below for a better way - and also how you should deal with constants.
You are misusing omp parallel do. In fact I suggest you never use this construct, I think it is a major source of confusion for beginners in the language, and in fact a mistake it is included in OMP at all. Why do I think this? Well it is important to understand that there are two stages when parallelizing via OpenMP. First you create the (non-master) threads with !$omp parallel. The second stage is to divide the work up between the threads with a worksharing construct such as !$omp do. The method you are using conflates the two, and blurs this important separation which in my experience leads to confusion and hence incorrect and inefficient code. I strongly recommend you always separate the two stages in your code - as I have done below.

I would also recommend you use default( none ) for scoping variables in the parallel region, it forces you to think about whether a variable should be shared or private, and also to avoid call omp_set_num_threads and use the environment variable OMP_NUM_THREADS instead, as you can then change the number of threads without having to recompile.

Fixing these leads me to the following code, which I think does what you want:

ijb@ijb-Latitude-5410:~/work/stack$ cat omp.f90
Program test

  Use iso_fortran_env, Only : wp => real64, stdout => output_unit

  !$ use omp_lib

  Implicit None

  Integer :: i
  Integer, Parameter :: dimen = 5
  ! Note
  ! 1) Kind = 8 is not portable and should not be used
  ! 2) Real constants are single precision unless otherwise specified
  Real (kind = wp), Dimension(dimen, dimen) :: tstarr = 0.0_wp

  !.....................................................................
  ! Later I will have to fill the elements of the following matrix in the same way
  ! real (kind = wp), dimension(:,:), allocatable :: tstarr2
  ! allocate(tstarr2(2 * dimen, 2 * dimen))
  !.....................................................................

  ! Note in Fortran executable statements can not be mixed with declarations
  ! However I wouldn't do it this way - use the environment variable OMP_NUM_THREADS, you don't have
  ! to recompile just to change the number of threads. Hence I have commented it out.
  !!!$ call omp_set_num_threads(4)

  ! First create the threads. I also strongly recommend default( none )
  !$omp parallel default( none ) shared( tstarr ) private( i )

  ! Report how many threads we are using
  !$omp single
  !$ Write( stdout, '( "Using ", i0, " threads" )' ) omp_get_num_threads()
  !$omp end single
  

  ! Now share out the work for each loop

  !$omp do
  Do i = 1, dimen
     tstarr(i,i) = 3.0_wp
  End Do
  !$omp end do

  !$omp do
  Do i = 1, dimen - 1
     tstarr(i + 1, i) = 2.0_wp
     tstarr(i, i + 1) = 2.0_wp
  End Do
  !$omp end do
  
  !$omp do
  Do i = 1, dimen - 2
     tstarr(i + 2, i) = 4.0_wp
     tstarr(i, i + 2) = 4.0_wp
  End Do
  !$omp end do
  
  !$omp do
  Do i = 1, dimen - 3
     tstarr(i + 3, i) = 5.0_wp
     tstarr(i, i + 3) = 5.0_wp
  End Do
  !$omp end do

  ! And end the parallel region
  !$omp end parallel

  ! Report the matrix - prettily
  Do i = 1, dimen
     Write( stdout, '( *( f3.1, 2x ) )' ) tstarr(i,:)
  End Do


End Program test
ijb@ijb-Latitude-5410:~/work/stack$ gfortran-12 --version
GNU Fortran (GCC) 12.2.0
Copyright © 2022 Free Software Foundation, Inc.
This is free software; see the source for copying conditions.  There is NO
warranty; not even for MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.

ijb@ijb-Latitude-5410:~/work/stack$ gfortran-12 -fopenmp -Wall -Wextra -fcheck=all -std=f2018 -Werror -g omp.f90
ijb@ijb-Latitude-5410:~/work/stack$ export OMP_NUM_THREADS=4
ijb@ijb-Latitude-5410:~/work/stack$ ./a.out
Using 4 threads
3.0  2.0  4.0  5.0  0.0
2.0  3.0  2.0  4.0  5.0
4.0  2.0  3.0  2.0  4.0
5.0  4.0  2.0  3.0  2.0
0.0  5.0  4.0  2.0  3.0
ijb@ijb-Latitude-5410:~/work/stack$ ./a.out
Using 4 threads
3.0  2.0  4.0  5.0  0.0
2.0  3.0  2.0  4.0  5.0
4.0  2.0  3.0  2.0  4.0
5.0  4.0  2.0  3.0  2.0
0.0  5.0  4.0  2.0  3.0
ijb@ijb-Latitude-5410:~/work/stack$ export OMP_NUM_THREADS=12
ijb@ijb-Latitude-5410:~/work/stack$ ./a.out
Using 12 threads
3.0  2.0  4.0  5.0  0.0
2.0  3.0  2.0  4.0  5.0
4.0  2.0  3.0  2.0  4.0
5.0  4.0  2.0  3.0  2.0
0.0  5.0  4.0  2.0  3.0
ijb@ijb-Latitude-5410:~/work/stack$ ./a.out
Using 12 threads
3.0  2.0  4.0  5.0  0.0
2.0  3.0  2.0  4.0  5.0
4.0  2.0  3.0  2.0  4.0
5.0  4.0  2.0  3.0  2.0
0.0  5.0  4.0  2.0  3.0
ijb@ijb-Latitude-5410:~/work/stack$ ./a.out
Using 12 threads
3.0  2.0  4.0  5.0  0.0
2.0  3.0  2.0  4.0  5.0
4.0  2.0  3.0  2.0  4.0
5.0  4.0  2.0  3.0  2.0
0.0  5.0  4.0  2.0  3.0
ijb@ijb-Latitude-5410:~/work/stack$

Use of OpenMP in filling the array elements

2 Answers2