QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
qrm_task_mod.F90
Go to the documentation of this file.
1 !! ##############################################################################################
2 !!
3 !! Copyright 2012 CNRS, INPT
4 !!
5 !! This file is part of qr_mumps.
6 !!
7 !! qr_mumps is free software: you can redistribute it and/or modify
8 !! it under the terms of the GNU Lesser General Public License as
9 !! published by the Free Software Foundation, either version 3 of
10 !! the License, or (at your option) any later version.
11 !!
12 !! qr_mumps is distributed in the hope that it will be useful,
13 !! but WITHOUT ANY WARRANTY; without even the implied warranty of
14 !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 !! GNU Lesser General Public License for more details.
16 !!
17 !! You can find a copy of the GNU Lesser General Public License
18 !! in the qr_mumps/doc directory.
19 !!
20 !! ##############################################################################################
21 
22 
23 !! ##############################################################################################
33 
34 #include "qrm_common.h"
35 
39 
90 
91 #if defined (_OPENMP)
92  use omp_lib
93 #endif
94  use qrm_common_mod
95  use qrm_error_mod
96 
97  ! parameter variables follow the convention of having an underscore at the end
98  ! 0=exit, 1=panel, 2=update, 3=activate, 4=assemble, 5=free, 6=clean
99  integer, parameter :: qrm_task_exit_ = 0 ! exit
100  integer, parameter :: qrm_task_pnl_ = 1 ! panel
101  integer, parameter :: qrm_task_upd_ = 2 ! update
102  integer, parameter :: qrm_task_act_ = 3 ! activate
103  integer, parameter :: qrm_task_asm_ = 4 ! assemble
104  integer, parameter :: qrm_task_free_ = 5 ! free
105  integer, parameter :: qrm_task_cln_ = 6 ! clean
106  integer, parameter :: qrm_task_app_ = 7 ! apply (for Q or Q')
107  integer, parameter :: qrm_task_sol_ = 8 ! solve (for R or R')
108 
112  integer :: id
113 
114  integer :: front
115 
116  integer :: pnl
117 
118  integer :: col
119 
120  logical :: bind
121  end type qrm_task_type
122 
123 
125  integer, parameter :: max_tasks = 300
126 
127  integer, private :: qrm_task_thn, qrm_task_nth
128  !$omp threadprivate(qrm_task_thn, qrm_task_nth)
129 
132  integer :: h, t, n
133  type(qrm_task_type) :: q(max_tasks)
134 #if defined(_OPENMP)
135  integer(kind=omp_lock_kind) :: lock
136 #endif
137  end type qrm_task_queue
138 
142  type(qrm_task_queue), allocatable :: queues(:)
143  integer, allocatable :: proxy_list(:,:)
144  integer :: stolen=0, ntsk=0, ncores, nnodes, cnode
145  end type qrm_task_queue_handle
146 
147 contains
148 
149 #ifndef NOLOC
150 
155  subroutine qrm_init_task_queue(h)
156  implicit none
157  type(qrm_task_queue_handle) :: h
158 
159  integer :: i
160 
161 #if defined (_OPENMP)
162  qrm_task_nth = omp_get_num_threads()
163  qrm_task_thn = omp_get_thread_num()
164 #else
165  qrm_task_nth = 1
166  qrm_task_thn = 0
167 #endif
168 
169 #if defined (hwloc)
170  !$omp master
171  call qrm_hwloc_info(h%ncores, h%nnodes, h%cnode)
172  !$omp end master
173  !$omp barrier
174  ! bind the thread but only if we are at level 1
175  !$ if (omp_get_level() .eq. 1) call qrm_hwloc_bind(mod(qrm_task_thn,h%ncores))
176 #else
177  !$omp master
178  h%ncores = qrm_task_nth
179  h%nnodes = 1
180  h%cnode = qrm_task_nth
181  !$omp end master
182 #endif
183 
184  !$omp master
185  call qrm_task_proximity(h)
186 
187  allocate(h%queues(0:qrm_task_nth-1))
188  do i=0, qrm_task_nth-1
189  h%queues(i)%h = 0
190  h%queues(i)%t = 0
191  h%queues(i)%n = 0
192  h%queues(i)%q = qrm_task_type(0, 0, 0, 0, .false.)
193  !$ call omp_init_lock(h%queues(i)%lock)
194  end do
195  !$omp end master
196  !$omp barrier !TODO maybe remove
197  return
198  end subroutine qrm_init_task_queue
199 
200 
204  subroutine qrm_task_proximity(h)
205  use iso_c_binding
206  implicit none
207  type(qrm_task_queue_handle) :: h
208 
209  integer :: i, j, ii, jj, cnt, id, c, n
210  integer, allocatable :: topo(:,:), core_to_node(:)
211 
212 #if defined (hwloc)
213  call qrm_hwloc_info(h%ncores, h%nnodes, h%cnode)
214  ! write(*,'("ncores: ",i2," nnodes: ",i2," cnode: ",i2)')h%ncores,h%nnodes,h%cnode
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))
217  call qrm_hwloc_topo(h%nnodes, topo)
218 #else
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))
224  cnt = 0
225  do j=0, h%nnodes-1
226  do i=0, h%cnode-1
227  topo(i,j)=cnt
228  cnt = cnt+1
229  end do
230  end do
231 #endif
232 
233  allocate(core_to_node(0:h%nnodes*h%cnode-1))
234  do j=0, h%nnodes-1
235  do i=0, h%cnode-1
236  core_to_node(topo(i,j)) = j
237  end do
238  end do
239 
240  do id=0, qrm_task_nth-1
241  cnt = 0
242  ! add thread id in first position
243  c = 0
244  do
245  if(c*h%ncores +id .ge. qrm_task_nth) exit
246  h%proxy_list(id,cnt) = c*h%ncores+id
247  cnt = cnt+1
248  c = c+1
249  end do
250 
251  ! add all the other threads
252  do j=0, h%nnodes-1
253  n = core_to_node(mod(id,h%ncores))
254  n = mod(n+j,h%nnodes)
255  do i = 0, h%cnode-1
256  c = 0
257  do
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)
261  cnt = cnt+1
262  c = c+1
263  end do
264  end do
265  end do
266 
267  end do
268 
269  deallocate(topo, core_to_node)
270 
271  return
272  end subroutine qrm_task_proximity
273 
288  function qrm_sched_task(h, tsk, pol, q)
289  implicit none
290  type(qrm_task_queue_handle) :: h
291  type(qrm_task_type) :: tsk
292  character :: pol
293  integer, optional :: q
294  logical :: qrm_sched_task
295 
296  integer :: i, iq
297 
298  if(present(q)) then
299  iq = q
300  else
301  iq = qrm_task_thn
302  end if
303 
304  !$ call omp_set_lock(h%queues(iq)%lock)
305  if(h%queues(iq)%n .lt. max_tasks) then
306  qrm_sched_task = .true. ! scheduling succeeded
307 ! !$omp critical(numtsk)
308  h%ntsk = h%ntsk+1
309 ! !$omp end critical(numtsk)
310  if(h%queues(iq)%n .eq. 0) then
311  i = 1
312  h%queues(iq)%n = 1
313  h%queues(iq)%h = 1
314  else
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
320  i = h%queues(iq)%h
321  end if
322  h%queues(iq)%n = h%queues(iq)%n+1
323  end if
324  h%queues(iq)%q(i) = tsk
325  else
326  qrm_sched_task = .false. ! sched failed
327  end if
328  !$ call omp_unset_lock(h%queues(iq)%lock)
329 
330  return
331  end function qrm_sched_task
332 
333 
342  function qrm_get_task(h, tsk)
343  implicit none
344  type(qrm_task_queue_handle) :: h
345  type(qrm_task_type) :: tsk
346  logical :: qrm_get_task
347 
348  integer :: i, queue
349 
350  tsk = qrm_task_type(-10, 0, 0, 0, .false.)
351  qrm_get_task = .false.
352 
353  do i=0, qrm_task_nth-1
354  queue = h%proxy_list(qrm_task_thn,i)
355  ! !$ call omp_set_lock(h%queues(queue)%lock)
356 #if defined(_OPENMP)
357  if(omp_test_lock(h%queues(queue)%lock)) then
358 #endif
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
362  qrm_get_task = .true.
363  h%queues(queue)%q(h%queues(queue)%h) = qrm_task_type(0, 0, 0, 0, .false.) !TODO: REMOVE!!!
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
367  if(i .gt. 0) then
368  h%stolen = h%stolen+1
369  end if
370  !$ call omp_unset_lock(h%queues(queue)%lock)
371  return
372  else
373  tsk = qrm_task_type(-10, 0, 0, 0, .false.)
374  qrm_get_task = .false.
375  end if
376 10 continue
377 #if defined(_OPENMP)
378  call omp_unset_lock(h%queues(queue)%lock)
379  end if
380 #endif
381  end do
382 
383 
384  return
385  end function qrm_get_task
386 
387 
393  type(qrm_task_queue_handle) :: h
394  integer :: qrm_task_queue_card
395 
396  ! qrm_task_queue_card = h%queues(qrm_task_thn)%n
397  qrm_task_queue_card = sum(h%queues(:)%n)
398 
399  return
400  end function qrm_task_queue_card
401 
410  function qrm_task_queue_empty(h, who)
411  type(qrm_task_queue_handle) :: h
412  integer, optional :: who
413  logical :: qrm_task_queue_empty
414 
415  integer :: i
416 
417  if(present(who)) then
418  !$ call omp_set_lock(h%queues(who)%lock)
419  qrm_task_queue_empty = h%queues(who)%n .eq. 0
420  !$ call omp_unset_lock(h%queues(who)%lock)
421  else
422  qrm_task_queue_empty = .true.
423  ! !$ do i=0, qrm_task_thn-1
424  ! !$ call omp_set_lock(h%queues(i)%lock)
425  ! !$ end do
426  qrm_task_queue_empty = qrm_task_queue_empty .and. (maxval(h%queues(0:qrm_task_nth-1)%n) .eq. 0)
427  ! !$ do i=qrm_task_thn-1, 0, -1
428  ! !$ call omp_unset_lock(h%queues(i)%lock)
429  ! !$ end do
430  end if
431 
432  return
433  end function qrm_task_queue_empty
434 
439  subroutine qrm_clean_task_queue(h)
440 
441  type(qrm_task_queue_handle) :: h
442 
443  integer :: i
444 
445  !$omp master
446  do i=0, qrm_task_nth-1
447  !$ call omp_destroy_lock(h%queues(i)%lock)
448  end do
449 #if defined(debug)
450  __qrm_prnt_dbg('("Queues cleaned -- ntsk:",i6," stolen:",i6)')h%ntsk,h%stolen
451 #endif
452 
453  deallocate(h%queues, h%proxy_list)
454  !$omp end master
455  return
456  end subroutine qrm_clean_task_queue
457 
458 
459 
460 
461 #else
462  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
463  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
464  !!
465  !! the following code is for a basic scheduling with no locality at
466  !! all. mostly used for testing
467  !!
468  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
469  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
470 
471 
472 
473 
479  subroutine qrm_init_task_queue(h)
480  implicit none
481  type(qrm_task_queue_handle) :: h
482 
483  !$omp master
484  allocate(h%queues(0:0))
485  h%queues(0)%h = 0
486  h%queues(0)%t = 0
487  h%queues(0)%n = 0
488  h%queues(0)%q = qrm_task_type(0, 0, 0, 0, .false.)
489  !$ call omp_init_lock(h%queues(0)%lock)
490  !$omp end master
491  !$omp barrier !TODO maybe remove
492  return
493  end subroutine qrm_init_task_queue
494 
495 
496 
511  function qrm_sched_task(h, tsk, pol, q)
512  implicit none
513  type(qrm_task_queue_handle) :: h
514  type(qrm_task_type) :: tsk
515  character :: pol
516  integer, optional :: q
517  logical :: qrm_sched_task
518 
519  integer :: i
520 
521  !$ call omp_set_lock(h%queues(0)%lock)
522  if(h%queues(0)%n .lt. max_tasks) then
523  qrm_sched_task = .true. ! scheduling succeeded
524  h%ntsk = h%ntsk+1
525  if(h%queues(0)%n .eq. 0) then
526  i = 1
527  h%queues(0)%n = 1
528  h%queues(0)%h = 1
529  else
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
535  i = h%queues(0)%h
536  end if
537  h%queues(0)%n = h%queues(0)%n+1
538  end if
539  h%queues(0)%q(i) = tsk
540  else
541  qrm_sched_task = .false. ! sched failed
542  end if
543  !$ call omp_unset_lock(h%queues(0)%lock)
544 
545  return
546  end function qrm_sched_task
547 
548 
557  function qrm_get_task(h, tsk)
558  implicit none
559  type(qrm_task_queue_handle) :: h
560  type(qrm_task_type) :: tsk
561  logical :: qrm_get_task
562 
563  tsk = qrm_task_type(-10, 0, 0, 0, .false.)
564  qrm_get_task = .false.
565 
566 #if defined(_OPENMP)
567  ! if(omp_test_lock(h%queues(0)%lock)) then
568  call omp_set_lock(h%queues(0)%lock)
569 #endif
570  if(h%queues(0)%n .gt. 0) then
571  tsk = h%queues(0)%q(h%queues(0)%h)
572  qrm_get_task = .true.
573  h%queues(0)%q(h%queues(0)%h) = qrm_task_type(0, 0, 0, 0, .false.) !TODO: REMOVE!!!
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
577  goto 10
578  else
579  tsk = qrm_task_type(-10, 0, 0, 0, .false.)
580  qrm_get_task = .false.
581  end if
582 10 continue
583 #if defined(_OPENMP)
584  call omp_unset_lock(h%queues(0)%lock)
585  ! end if
586 #endif
587 
588  return
589  end function qrm_get_task
590 
591 
597  type(qrm_task_queue_handle) :: h
598  integer :: qrm_task_queue_card
599 
600  qrm_task_queue_card = h%queues(0)%n
601 
602  return
603  end function qrm_task_queue_card
604 
613  function qrm_task_queue_empty(h, who)
614  type(qrm_task_queue_handle) :: h
615  integer, optional :: who
616  logical :: qrm_task_queue_empty
617 
618  qrm_task_queue_empty = h%queues(0)%n .eq. 0
619 
620  return
621  end function qrm_task_queue_empty
622 
627  subroutine qrm_clean_task_queue(h)
628 
629  type(qrm_task_queue_handle) :: h
630 
631  !$omp master
632  !$ call omp_destroy_lock(h%queues(0)%lock)
633 #if defined(debug)
634  __qrm_prnt_dbg('("Queues cleaned -- ntsk:",i6," stolen:",i6)')h%ntsk,h%stolen
635 #endif
636 
637  deallocate(h%queues)
638  !$omp end master
639  return
640  end subroutine qrm_clean_task_queue
641 
642 
643 #endif
644 
645 
646 end module qrm_task_mod
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.
int i
Definition: secs.c:40
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.