1 2! Copyright (C) 2015 J. K. Dewhurst, S. Sharma and E. K. U. Gross. 3! This file is distributed under the terms of the GNU General Public License. 4! See the file COPYING for license details. 5 6module modomp 7 8! maximum number of OpenMP threads available 9integer maxthd 10! maximum number of OpenMP threads for the first nesting level 11integer maxthd1 12! maximum number of threads available to MKL 13integer maxthdmkl 14! maximum OpenMP nesting level 15integer maxlvl 16! number of active OpenMP threads for each nesting level 17integer, allocatable, protected :: nathd(:) 18 19interface 20 21integer function omp_get_num_procs() 22end function 23 24integer function omp_get_max_threads() 25end function 26 27integer function omp_get_num_threads() 28end function 29 30integer function omp_get_thread_num() 31end function 32 33integer function omp_get_max_active_levels() 34end function 35 36logical function omp_get_dynamic() 37end function 38 39integer function omp_get_level() 40end function 41 42subroutine omp_set_num_threads(num_threads) 43integer, intent(in) :: num_threads 44end subroutine 45 46integer function mkl_set_num_threads_local(num_threads) 47integer, intent(in) :: num_threads 48end function 49 50subroutine omp_set_nested(nested) 51logical, intent(in) :: nested 52end subroutine 53 54subroutine omp_set_max_active_levels(max_levels) 55integer, intent(in) :: max_levels 56end subroutine 57 58subroutine omp_set_dynamic(dynamic_threads) 59logical, intent(in) :: dynamic_threads 60end subroutine 61 62end interface 63 64contains 65 66subroutine omp_init 67implicit none 68! determine the maximum number of available threads 69select case(maxthd) 70case(:-1) 71! set the number of threads equal to a fraction of the number of processors 72 maxthd=omp_get_num_procs()/abs(maxthd) 73 maxthd=max(maxthd,1) 74 call omp_set_num_threads(maxthd) 75case(0) 76! use the system default number of threads 77 maxthd=omp_get_max_threads() 78case default 79! use the number of threads specified in the input file 80 call omp_set_num_threads(maxthd) 81end select 82! determine the maximum number of threads available at first nesting level 83select case(maxthd1) 84case(:-1) 85 maxthd1=maxthd/abs(maxthd1) 86 maxthd1=max(maxthd1,1) 87case(0) 88 maxthd1=maxthd 89case default 90 maxthd1=min(maxthd1,maxthd) 91end select 92! switch off dynamic allocation of threads 93call omp_set_dynamic(.false.) 94! allow nested parallelism (deprecated in OpenMP version 5) 95call omp_set_nested(.true.) 96! set the maximum nesting level 97call omp_set_max_active_levels(maxlvl) 98! allocate the number of active threads array 99if (allocated(nathd)) deallocate(nathd) 100allocate(nathd(0:maxlvl)) 101! initialise the number of active threads 102call omp_reset 103end subroutine 104 105subroutine omp_reset 106implicit none 107! number of active threads at each nesting level 108nathd(0)=1 109nathd(1:)=0 110end subroutine 111 112subroutine holdthd(nloop,nthd) 113implicit none 114! arguments 115integer, intent(in) :: nloop 116integer, intent(out) :: nthd 117! local variables 118integer lvl,na,n 119! current nesting level 120lvl=omp_get_level() 121if ((lvl.lt.0).or.(lvl.ge.maxlvl)) then 122 nthd=1 123 return 124end if 125!$OMP CRITICAL(modomp_) 126! determine number of active threads at the current nesting level 127na=nathd(lvl) 128na=max(min(na,maxthd),1) 129! number of threads allowed for this loop 130nthd=maxthd/na 131if (mod(maxthd,na).gt.0) nthd=nthd+1 132if (lvl.eq.0) nthd=min(nthd,maxthd1) 133nthd=max(min(nthd,maxthd,nloop),1) 134! add to number of active threads in next nesting level 135n=nathd(lvl+1)+nthd 136n=max(min(n,maxthd),0) 137nathd(lvl+1)=n 138!$OMP END CRITICAL(modomp_) 139end subroutine 140 141subroutine freethd(nthd) 142implicit none 143! arguments 144integer, intent(in) :: nthd 145! local variables 146integer lvl,n 147! current nesting level 148lvl=omp_get_level() 149if ((lvl.lt.0).or.(lvl.ge.maxlvl)) return 150!$OMP CRITICAL(modomp_) 151! subtract from the number of active threads in next nesting level 152n=nathd(lvl+1)-nthd 153n=max(min(n,maxthd),0) 154nathd(lvl+1)=n 155!$OMP END CRITICAL(modomp_) 156end subroutine 157 158end module 159 160