34 #include "qrm_common.h"
99 integer,
parameter :: qrm_task_exit_ = 0
100 integer,
parameter :: qrm_task_pnl_ = 1
101 integer,
parameter :: qrm_task_upd_ = 2
102 integer,
parameter :: qrm_task_act_ = 3
103 integer,
parameter :: qrm_task_asm_ = 4
104 integer,
parameter :: qrm_task_free_ = 5
105 integer,
parameter :: qrm_task_cln_ = 6
106 integer,
parameter :: qrm_task_app_ = 7
107 integer,
parameter :: qrm_task_sol_ = 8
125 integer,
parameter :: max_tasks = 300
127 integer,
private :: qrm_task_thn, qrm_task_nth
135 integer(kind=omp_lock_kind) :: lock
143 integer,
allocatable :: proxy_list(:,:)
144 integer :: stolen=0, ntsk=0, ncores, nnodes, cnode
161 #if defined (_OPENMP)
162 qrm_task_nth = omp_get_num_threads()
163 qrm_task_thn = omp_get_thread_num()
178 h%ncores = qrm_task_nth
180 h%cnode = qrm_task_nth
187 allocate(h%queues(0:qrm_task_nth-1))
188 do i=0, qrm_task_nth-1
209 integer ::
i, j, ii, jj, cnt, id, c, n
210 integer,
allocatable :: topo(:,:), core_to_node(:)
215 allocate(topo(0:h%cnode-1,0:h%nnodes-1))
216 allocate(h%proxy_list(0:qrm_task_nth-1,0:qrm_task_nth-1))
219 h%cnode = min(qrm_task_nth,1)
220 h%ncores = qrm_task_nth
221 h%nnodes = (qrm_task_nth-1)/h%cnode+1
222 allocate(topo(0:h%cnode-1,0:h%nnodes-1))
223 allocate(h%proxy_list(0:qrm_task_nth-1,0:qrm_task_nth-1))
233 allocate(core_to_node(0:h%nnodes*h%cnode-1))
236 core_to_node(topo(
i,j)) = j
240 do id=0, qrm_task_nth-1
245 if(c*h%ncores +id .ge. qrm_task_nth)
exit
246 h%proxy_list(id,cnt) = c*h%ncores+id
253 n = core_to_node(mod(id,h%ncores))
254 n = mod(n+j,h%nnodes)
258 if(c*h%ncores +topo(
i,n) .ge. qrm_task_nth)
exit
259 if(c*h%ncores +topo(
i,n) .eq. id)
exit
260 h%proxy_list(id,cnt) = c*h%ncores+topo(
i,n)
269 deallocate(topo, core_to_node)
293 integer,
optional :: q
305 if(h%queues(iq)%n .lt. max_tasks)
then
310 if(h%queues(iq)%n .eq. 0)
then
315 if(pol .eq.
't')
then
316 i = mod(h%queues(iq)%h+h%queues(iq)%n-1,max_tasks)+1
317 else if(pol .eq.
'h')
then
318 h%queues(iq)%h = h%queues(iq)%h-1
319 if (h%queues(iq)%h .le. 0) h%queues(iq)%h=max_tasks
322 h%queues(iq)%n = h%queues(iq)%n+1
324 h%queues(iq)%q(
i) = tsk
353 do i=0, qrm_task_nth-1
354 queue = h%proxy_list(qrm_task_thn,
i)
357 if(omp_test_lock(h%queues(queue)%lock))
then
359 if(h%queues(queue)%n .gt. 0)
then
360 tsk = h%queues(queue)%q(h%queues(queue)%h)
361 if(tsk%bind .and.
i .ne. 0) goto 10
363 h%queues(queue)%q(h%queues(queue)%h) =
qrm_task_type(0, 0, 0, 0, .false.)
364 h%queues(queue)%h = mod(h%queues(queue)%h, max_tasks)+1
365 h%queues(queue)%n = h%queues(queue)%n-1
366 if (h%queues(queue)%n .eq. 0) h%queues(queue)%h = 0
368 h%stolen = h%stolen+1
378 call omp_unset_lock(h%queues(queue)%lock)
412 integer,
optional :: who
417 if(present(who))
then
446 do i=0, qrm_task_nth-1
450 __qrm_prnt_dbg(
'("Queues cleaned -- ntsk:",i6," stolen:",i6)')h%ntsk,h%stolen
453 deallocate(h%queues, h%proxy_list)
484 allocate(h%queues(0:0))
516 integer,
optional :: q
522 if(h%queues(0)%n .lt. max_tasks)
then
525 if(h%queues(0)%n .eq. 0)
then
530 if(pol .eq.
't')
then
531 i = mod(h%queues(0)%h+h%queues(0)%n-1,max_tasks)+1
532 else if(pol .eq.
'h')
then
533 h%queues(0)%h = h%queues(0)%h-1
534 if (h%queues(0)%h .le. 0) h%queues(0)%h=max_tasks
537 h%queues(0)%n = h%queues(0)%n+1
539 h%queues(0)%q(
i) = tsk
568 call omp_set_lock(h%queues(0)%lock)
570 if(h%queues(0)%n .gt. 0)
then
571 tsk = h%queues(0)%q(h%queues(0)%h)
573 h%queues(0)%q(h%queues(0)%h) =
qrm_task_type(0, 0, 0, 0, .false.)
574 h%queues(0)%h = mod(h%queues(0)%h, max_tasks)+1
575 h%queues(0)%n = h%queues(0)%n-1
576 if (h%queues(0)%n .eq. 0) h%queues(0)%h = 0
584 call omp_unset_lock(h%queues(0)%lock)
615 integer,
optional :: who
634 __qrm_prnt_dbg(
'("Queues cleaned -- ntsk:",i6," stolen:",i6)')h%ntsk,h%stolen
This module contains all the error management routines and data.
subroutine qrm_init_task_queue(h)
Inititalizes a set of queues attached to a family of threads referenced through the handle h...
logical function qrm_task_queue_empty(h, who)
Tells whether one, or all, queues are empty.
Generic interface for the ::qrm_hwloc_topo routine.
subroutine qrm_clean_task_queue(h)
Destroyes a set of queues.
This type defines the handle for the queues attached to a family of threads.
logical function qrm_get_task(h, tsk)
Pops a task from a queue. Tasks are always popped from the head of the queue. The return value is ...
This module contains the interfaces of all non-typed routines.
This type defines a computational task.
Generic interface for the ::qrm_hwloc_info routine.
integer function qrm_task_queue_card(h)
Returns the number of tasks present on a set of queues referenced by a handle.
subroutine qrm_task_proximity(h)
Defines the order in which queues have to be visited by each thread.
logical function qrm_sched_task(h, tsk, pol, q)
Pushes a task on a queue.
This module contains the definition of a task type that is used for scheduling tasks during the facto...
This type defines the task queue attached to a thread.