File: modomp.f90

package info (click to toggle)
elkcode 5.4.24-2
  • links: PTS, VCS
  • area: main
  • in suites: buster, sid
  • size: 12,840 kB
  • sloc: f90: 48,415; fortran: 22,457; perl: 965; makefile: 384; sh: 369; python: 105; ansic: 67
file content (148 lines) | stat: -rw-r--r-- 3,326 bytes parent folder | download
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148

! Copyright (C) 2015 J. K. Dewhurst, S. Sharma and E. K. U. Gross.
! This file is distributed under the terms of the GNU General Public License.
! See the file COPYING for license details.

module modomp

! number of processors available for OpenMP
integer ncpu
! maximum number of OpenMP threads available
integer maxthd
! maximum OpenMP nesting level
integer maxlvl
! number of active OpenMP threads for each nesting level
integer, allocatable :: nathd(:)

interface

integer function omp_get_num_procs()
end function

integer function omp_get_max_threads()
end function

integer function omp_get_num_threads()
end function

integer function omp_get_thread_num()
end function

logical function omp_get_nested()
end function

integer function omp_get_max_active_levels()
end function

logical function omp_get_dynamic()
end function

integer function omp_get_level()
end function

subroutine omp_set_num_threads(num_threads)
integer, intent(in) :: num_threads
end subroutine

subroutine omp_set_nested(nested)
logical, intent(in) :: nested
end subroutine

subroutine omp_set_max_active_levels(max_levels)
integer, intent(in) :: max_levels
end subroutine

subroutine omp_set_dynamic(dynamic_threads)
logical, intent(in) :: dynamic_threads
end subroutine

end interface

contains

subroutine omp_init
implicit none
! get number of processors
ncpu=omp_get_num_procs()
if (maxthd.lt.0) then
! set the number of threads equal to the number of processors
  maxthd=ncpu
  call omp_set_num_threads(maxthd)
else if (maxthd.eq.0) then
! use the system default number of threads
  maxthd=omp_get_max_threads()
else
! use the number of threads specified in the input file
  call omp_set_num_threads(maxthd)
end if
! switch off dynamic allocation of threads
call omp_set_dynamic(.false.)
! allow nested parallelism
call omp_set_nested(.true.)
! set the maximum nesting level
call omp_set_max_active_levels(maxlvl)
! allocate the number of active threads array
if (allocated(nathd)) deallocate(nathd)
allocate(nathd(0:maxlvl))
! initialise the number of active threads
call omp_reset
return
end subroutine

subroutine omp_reset
implicit none
! number of active threads at each nesting level
nathd(0)=1
nathd(1:)=0
return
end subroutine

subroutine omp_hold(nloop,nthd)
implicit none
! arguments
integer, intent(in) :: nloop
integer, intent(out) :: nthd
! local variables
integer lvl,na,n
! current nesting level
lvl=omp_get_level()
if ((lvl.lt.0).or.(lvl.ge.maxlvl)) then
  nthd=1
  return
end if
!$OMP CRITICAL(omp_hold_)
! determine number of active threads at the current nesting level
na=nathd(lvl)
na=max(min(na,maxthd),1)
! number of threads allowed for this loop
nthd=maxthd/na
if (mod(maxthd,na).gt.0) nthd=nthd+1
nthd=max(min(nthd,maxthd,nloop),1)
! add to number of active threads in next nesting level
n=nathd(lvl+1)+nthd
n=max(min(n,maxthd),0)
nathd(lvl+1)=n
!$OMP END CRITICAL(omp_hold_)
return
end subroutine

subroutine omp_free(nthd)
implicit none
! arguments
integer, intent(in) :: nthd
! local variables
integer lvl,n
! current nesting level
lvl=omp_get_level()
if ((lvl.lt.0).or.(lvl.ge.maxlvl)) return
!$OMP CRITICAL(omp_free_)
! subtract from the number of active threads in next nesting level
n=nathd(lvl+1)-nthd
n=max(min(n,maxthd),0)
nathd(lvl+1)=n
!$OMP END CRITICAL(omp_free_)
return
end subroutine

end module