QR_MUMPS
qrm_queue_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 
36 
48 
49 #if defined(_OPENMP)
50  use omp_lib
51 #endif
52 
54  integer, parameter :: qrm_fifo_=0
55 
57  integer, parameter :: qrm_lifo_=1
58 
60  type qrm_queue
61  integer, allocatable :: elems(:)
62  integer :: nelems
63  integer :: h, t
64  integer :: maxelems
65  integer :: pol
66 #if defined (_OPENMP)
67  integer(omp_lock_kind) :: lock
68 #endif
69  end type qrm_queue
70 
71  interface qrm_queue_init
72  module procedure qrm_queue_init
73  end interface
74 
75 
76 contains
77 
89  subroutine qrm_queue_init(q, nelems, pol)
91  implicit none
92 
93  type(qrm_queue) :: q
94  integer :: nelems, pol
95 
96  call qrm_aalloc(q%elems, nelems)
97  q%elems = 0
98  q%h = 0
99  q%t = 0
100  q%nelems = 0
101  q%maxelems = nelems
102  q%pol = pol
103  !$ call omp_init_lock(q%lock)
104  return
105 
106  end subroutine qrm_queue_init
107 
108 
113  subroutine qrm_queue_free(q)
115  implicit none
116 
117  type(qrm_queue) :: q
118 
119  call qrm_adealloc(q%elems)
120 
121  q%nelems = 0
122  q%h = 0
123  q%t = 0
124  !$ call omp_destroy_lock(q%lock)
125 
126  return
127  end subroutine qrm_queue_free
128 
129 
130 
137  subroutine qrm_queue_push(q, elem)
138  implicit none
139 
140  type(qrm_queue) :: q
141  integer :: elem
142 
143  !$ call omp_set_lock(q%lock)
144  if (q%nelems .eq. 0) then
145  q%h = elem
146  q%t = elem
147  q%elems(elem) = -1
148  else if (q%nelems .eq. q%maxelems) then
149  write(*,'("Cannot push anymore", i4,x,i4)')elem, q%nelems
150  goto 10
151  else
152  if(q%pol .eq. qrm_fifo_) then
153  q%elems(q%t) = elem
154  q%elems(elem) = -1
155  q%t = elem
156  else if(q%pol .eq. qrm_lifo_) then
157  q%elems(elem) = q%h
158  q%h = elem
159  end if
160  end if
161 
162  q%nelems = q%nelems+1
163  ! write(*,'("Pushed ",i4,x,i4,x,i4)')elem,q%nelems,q%h,q%t
164 10 continue
165  !$ call omp_unset_lock(q%lock)
166 
167  return
168 
169  end subroutine qrm_queue_push
170 
171 
172 
173 
178  subroutine qrm_queue_prnt(q)
179  implicit none
180 
181  type(qrm_queue) :: q
182 
183  integer i, elem
184 
185  !$ call omp_set_lock(q%lock)
186  elem = q%h
187  do i=1, q%nelems
188  write(*,'(i3,"->")',advance='no')elem
189  elem = q%elems(elem)
190  end do
191  write(*,'(" ")')
192  !$ call omp_unset_lock(q%lock)
193  return
194  end subroutine qrm_queue_prnt
195 
196 
201  function qrm_queue_pop(q)
202  implicit none
203 
204  type(qrm_queue) :: q
205  integer :: qrm_queue_pop
206 
207  !$ call omp_set_lock(q%lock)
208  if (q%nelems .eq. 0) then
209  qrm_queue_pop = 0
210  else
211  qrm_queue_pop = q%h
212  q%h = q%elems(q%h)
213  q%elems(qrm_queue_pop) = 0
214  q%nelems = q%nelems-1
215  end if
216  !$ call omp_unset_lock(q%lock)
217 
218  return
219 
220  end function qrm_queue_pop
221 
222 
229  subroutine qrm_queue_rm(q, n)
230  implicit none
231 
232  type(qrm_queue) :: q
233  integer :: n
234 
235  integer :: tmp, i, s
236 
237 
238  !$ call omp_set_lock(q%lock)
239  if (q%nelems .eq. 0) then
240  goto 10
241  end if
242 
243  if(n .eq. q%h) then
244  q%h = q%elems(n)
245  q%nelems = q%nelems-1
246  else
247  tmp = q%h
248  do i =1, q%nelems-1
249  if (q%elems(tmp) .eq. n) then
250  q%elems(tmp) = q%elems(n)
251  q%nelems = q%nelems-1
252  if(n .eq. q%t) q%t = tmp
253  exit
254  end if
255  tmp = q%elems(tmp)
256  end do
257 
258  end if
259  q%elems(n) = 0
260 
261 10 continue
262  !$ call omp_unset_lock(q%lock)
263 
264  return
265 
266  end subroutine qrm_queue_rm
267 
268 
286  function qrm_queue_next(q, n)
287  implicit none
288 
289  type(qrm_queue) :: q
290  integer :: qrm_queue_next, n
291 
292  integer :: i
293 
294 #if defined(_OPENMP)
295  call omp_set_lock(q%lock)
296 #endif
297  if (q%nelems .eq. 0) then
298  qrm_queue_next = 0
299  else if (n .eq. q%t) then
300  qrm_queue_next = 0
301  else
302 
303  if(n.eq.0) then
304  qrm_queue_next = q%h
305  else
306  qrm_queue_next = q%elems(n)
307  end if
308  end if
309 #if defined(_OPENMP)
310  call omp_unset_lock(q%lock)
311 #endif
312  return
313  end function qrm_queue_next
314 
315 end module qrm_queue_mod
Generic interface for the qrm_adealloc_i, qrm_adealloc_2i, qrm_adealloc_s, qrm_adealloc_2s, qrm_adealloc_3s, qrm_adealloc_d, qrm_adealloc_2d, qrm_adealloc_3d, qrm_adealloc_c, qrm_adealloc_2c, qrm_adealloc_3c, qrm_adealloc_z, qrm_adealloc_2z, qrm_adealloc_3z, routines.
This module contains all the facilities for front queues.
subroutine qrm_queue_prnt(q)
Prints the content of a queue.
A data type meant to to define a queue.
subroutine qrm_queue_rm(q, n)
Removes (without returning it) an element from a queue.
integer function qrm_queue_pop(q)
Pops an element from a queue.
Generic interface for the qrm_aalloc_i, qrm_aalloc_2i, qrm_aalloc_s, qrm_aalloc_2s, qrm_aalloc_3s, qrm_aalloc_d, qrm_aalloc_2d, qrm_aalloc_3d, qrm_aalloc_c, qrm_aalloc_2c, qrm_aalloc_3c, qrm_aalloc_z, qrm_aalloc_2z, qrm_aalloc_3z, routines.
Definition: qrm_mem_mod.F90:78
subroutine qrm_queue_push(q, elem)
Pushes an element on a queue.
integer, parameter qrm_fifo_
parameter to define the policy of the queue: FIFO
subroutine qrm_queue_free(q)
Frees a queue.
integer, parameter qrm_lifo_
parameter to define the policy of the queue: LIFO
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38
integer function qrm_queue_next(q, n)
Returns the element that follows n in the queue q. Very useful for sweeping through a queue...