DISPATCH
list_mod.f90
1 !*******************************************************************************
2 !> Module with list handling for generic class task_t objects
3 !*******************************************************************************
4 MODULE list_mod
5  USE io_mod
6  USE omp_mod
7  USE omp_lock_mod
8  USE omp_timer_mod
9  USE mpi_mod
10  USE mpi_mesg_mod
11  USE trace_mod
12  USE bits_mod
13  USE task_mod
14  USE patch_mod
15  USE link_mod
16  !USE load_balance_mod
17  USE timer_mod
18  implicit none
19  private
20  type, public:: list_t
21  character(len=64):: name=''
22  class(link_t), pointer:: head => null()
23  class(link_t), pointer:: tail => null()
24  class(link_t), pointer:: queue => null()
25  class(link_t), pointer:: active => null()
26  integer:: n=0 ! total
27  integer:: nq=0 ! ready
28  integer:: nv=0 ! virtual
29  integer:: nb=0 ! boundary
30  integer:: ni=0 ! internal
31  integer:: ne=0 ! external
32  integer:: na=0 ! active = internal+boundary
33  integer:: nac=0 ! active tasks
34  real(8):: lc(3) = +huge(0d0) ! lower corner
35  real(8):: uc(3) = -huge(0d0) ! upper corner
36  real(8):: llc(3) = +huge(0d0) ! lower left corner
37  real(8):: urc(3) = -huge(0d0) ! upper right corner
38  real(8):: position(3) = 0.0_8 ! center offset
39  real(8):: size(3) = 1.0_8 ! full extent = uc-lc
40  integer:: n_behind=0, n_check=0, n_ready=0, n_nbor=0
41  integer:: debug=0
42  integer:: n_tasks=1
43  integer:: dims(3)=0
44  integer:: verbose=0
45  type(lock_t):: lock
46  character(len=4):: kind='list'
47  logical:: face_nbors=.false.
48  logical:: detailed_timer=.false.
49  contains
50  procedure:: init
51  procedure:: init_bdries
52  procedure:: reset
53  procedure:: append_task
54  procedure:: append_link
55  generic:: append => append_task, append_link
56  procedure:: prepend_link
57  procedure:: remove_task
58  procedure:: remove_link
59  generic:: remove => remove_task, remove_link
60  procedure:: remove_and_reset
61  procedure:: update_counts
62  procedure:: count_status
63  procedure:: reset_status
64  procedure:: check_nbors
65  procedure:: check_ready
66  procedure:: check_all
67  procedure:: resend_bdry
68  procedure:: check_oldest
69  procedure:: check_queue
70  procedure:: add_new_link
71  procedure:: init_nbors
72  procedure:: set_init_nbors
73  procedure:: init_nbor_nbors
74  procedure:: check_nbor_nbors
75  procedure:: refresh_nbors
76  procedure:: init_all_nbors
77  procedure:: remove_parents
78  procedure:: intpos
79  procedure:: consistency
80  procedure:: qconsistency
81  procedure:: aconsistency
82  procedure:: statistics
83  procedure:: append_list
84  procedure:: add_by_quality
85  procedure:: queue_by_time
86  procedure:: queue_active
87  procedure:: remove_active
88  procedure:: print_queue
89  procedure:: print_queue_until
90  procedure:: print_queue_times
91  procedure:: print => print_list
92  procedure:: print_tasks
93  procedure, nopass:: check_nbor_list
94  procedure:: give_to
95  procedure:: update_nbor_status
96  procedure:: update_status
97  procedure:: send_to_vnbors
98  procedure:: test_nbor_status
99  procedure:: info
100  end type
101  integer, dimension(:), allocatable:: sequence
102  integer, parameter:: mnbor=100
103  type(list_t), public:: list
104  public link_t, task2patch
105 CONTAINS
106 
107 !===============================================================================
108 !> Initialize list
109 !===============================================================================
110 SUBROUTINE init (self, name)
111  class(list_t):: self
112  character(len=*), optional:: name
113  !.............................................................................
114  call trace%begin ('list_t%init')
115  if (present(name)) then
116  self%name = name
117  write (io_unit%log,*)'init_list: ', name
118  call self%lock%init (name(1:4))
119  else
120  write (io_unit%log,*)'init_list: list_t'
121  call self%lock%init ('list_t')
122  end if
123  !$omp critical (sequence_alloc_cr)
124  if (.not.allocated(sequence)) then
125  allocate (sequence(0:mpi%size-1))
126  sequence(:) = 0
127  end if
128  !$omp end critical (sequence_alloc_cr)
129  call trace%end()
130 END SUBROUTINE init
131 
132 !===============================================================================
133 !> Reset list, by deallocating memory used by tasks, and deallocating tasks and
134 !> links
135 !===============================================================================
136 SUBROUTINE reset (self)
137  class(list_t):: self
138  class(link_t), pointer:: link, old
139  class(task_t), pointer:: task
140  !.............................................................................
141  call trace%begin ('list_t%reset '//self%name)
142  link => self%head
143  do while (associated(link))
144  task => link%task
145  if (associated(link)) then
146  call task%dealloc
147  deallocate (task)
148  end if
149  old => link
150  link => link%next
151  deallocate (old)
152  end do
153  nullify (self%head)
154  nullify (self%tail)
155  self%n = 0
156  call trace%end()
157 END SUBROUTINE reset
158 
159 !===============================================================================
160 !> Allocate a new link and append to the list, optionally setting attributes
161 !===============================================================================
162 SUBROUTINE append_task (self, task, nbor, needed, needs_me)
163  class(list_t):: self
164  class(link_t), pointer, optional:: nbor
165  logical, optional:: needed, needs_me
166  class(task_t), pointer:: task
167  class(link_t), pointer:: link
168  !.............................................................................
169  call trace%begin ('list_t%append_task '//self%name)
170  allocate (link)
171  call link%init
172  link%task => task
173  if (present(nbor)) link%nbor => nbor
174  if (present(needed)) link%needed = needed
175  if (present(needs_me)) link%needs_me = needs_me
176  call self%append_link (link)
177  call trace%end()
178 END SUBROUTINE append_task
179 
180 !===============================================================================
181 !> Append to the task list, counting active tasks only
182 !===============================================================================
183 SUBROUTINE append_link (self, link)
184  class(list_t):: self
185  class(link_t), pointer:: link
186  class(task_t), pointer:: task
187  class(*), pointer:: anon
188  logical, save:: first_time=.true.
189  !.............................................................................
190  call trace%begin ('list_t%append_link ')
191  call link%init
192  if (self%verbose > 2) &
193  write (io%output,*) wallclock(),' thread',omp%thread,' wait for tasklist(2)'
194  call self%lock%set ('list_t%append_link')
195  if (self%verbose > 2) &
196  write (io%output,*) wallclock(),' thread',omp%thread,' locked tasklist(2)'
197  nullify(link%next)
198  if (associated(self%tail)) then
199  link%prev => self%tail
200  link%prev%next => link
201  self%tail => link
202  else
203  self%head => link
204  self%tail => link
205  end if
206  if (associated(link%task)) then
207  task => link%task
208  select type (task)
209  class is (patch_t)
210  task%link => link
211  end select
212  else
213  call io%abort('append_link: no task associated')
214  end if
215  !-----------------------------------------------------------------------------
216  ! IMPORTANT: Some procedures depend on knowing the size of the list!
217  !-----------------------------------------------------------------------------
218  self%n = self%n + 1
219  if (self%verbose > 2) &
220  print *,'list_t%append: (1) n =', trim(self%name), self%n, link%task%id
221  call self%update_counts (link, +1)
222  if (self%verbose > 1) &
223  print *,'list _t%append:(2) n =', trim(self%name), self%n, link%task%id
224  call self%lock%unset ('list_t%append_link')
225  if (self%verbose > 2) &
226  write (io%output,*) wallclock(),' thread',omp%thread,' unlocked tasklist(2)'
227  call trace%end()
228 END SUBROUTINE append_link
229 
230 !===============================================================================
231 !> Prepend to the task list, counting active tasks only
232 !===============================================================================
233 SUBROUTINE prepend_link (self, link)
234  class(list_t):: self
235  class(link_t), pointer:: link
236  class(task_t), pointer:: task
237  class(*), pointer:: anon
238  logical, save:: first_time=.true.
239  !.............................................................................
240  call trace%begin ('list_t%prepend_link '//self%name)
241  call link%init
242  call self%lock%set ('list_t%prepend_link')
243  nullify(link%prev)
244  if (associated(self%head)) then
245  link%next => self%head
246  link%next%prev => link
247  self%head => link
248  else
249  self%head => link
250  self%tail => link
251  end if
252  if (associated(link%task)) then
253  task => link%task
254  select type (task)
255  class is (patch_t)
256  task%link => link
257  end select
258  call self%update_counts (link, +1)
259  else
260  call io%abort('prepend_link: no task associated')
261  end if
262  !-----------------------------------------------------------------------------
263  ! IMPORTANT: Some procedures depend on knowing the size of the list!
264  !-----------------------------------------------------------------------------
265  call self%lock%unset ('list_t%prepend_link')
266  call trace%end()
267 END SUBROUTINE prepend_link
268 
269 !===============================================================================
270 !> Remove a task from consideration, making sure that virtual copies do the same.
271 !> Set bits%remove, which excludes the task from nbor lists in check_ready(),
272 !> and signals rank nbors to repeat this remove_and_reset().
273 !===============================================================================
274 SUBROUTINE remove_and_reset (self, link)
275  class(list_t):: self
276  class(link_t), pointer:: link
277  integer:: count
278  !-----------------------------------------------------------------------------
279  call trace%begin ('list_t%remove_and_reset')
280  !-----------------------------------------------------------------------------
281  ! By setting bits%remove immediately, we make sure that the task is no longer
282  ! considered in check_ready() calls, and also that the task is no longer
283  ! counted in check_ready() calls to nbors of it.
284  !-----------------------------------------------------------------------------
285  call link%task%set (bits%remove) ! cf. check_ready() and unpack()
286  !-----------------------------------------------------------------------------
287  ! Send a suicide note to virtual nbors, if this is a boundary patch
288  !-----------------------------------------------------------------------------
289  if (link%task%is_set (bits%boundary)) then
290  call self%send_to_vnbors (link)
291  end if
292  !-----------------------------------------------------------------------------
293  if (self%verbose > 0) then
294  write (io_unit%log,'(f12.6,2x,a,i6,2l3)') &
295  wallclock(), 'list_t%remove_and_reset: id =', link%task%id, &
296  link%task%is_set(bits%boundary), link%task%is_set(bits%virtual)
297  end if
298  !-----------------------------------------------------------------------------
299  call self%lock%set ('remove_and_reset') ! lock task list while changing
300  call self%remove_link (link) ! remove link from task list
301  call self%lock%unset ('remove_and_reset')! release the task list
302  call self%check_nbors (link) ! check nbors, now link is removed
303  call self%set_init_nbors (link) ! trigger new nbor lists on nbors
304  !$omp atomic update
305  timer%levelpop (link%task%level) = & ! decrement level population count
306  timer%levelpop (link%task%level) - 1
307  call link%garbage_collect (link) ! add to garbage bin
308  call trace%end()
309 END SUBROUTINE remove_and_reset
310 
311 !===============================================================================
312 !> Remove a link from the list, taking special care with the head and tail
313 !===============================================================================
314 SUBROUTINE remove_link (self, link)
315  class(list_t):: self
316  class(link_t), pointer:: link, next
317  class(task_t), pointer:: task
318  !.............................................................................
319  call trace%begin ('list_t%remove_link '//self%name)
320  next => link%next
321  !-----------------------------------------------------------------------------
322  ! Check for tail task
323  !-----------------------------------------------------------------------------
324  if (associated(next)) then
325  !call next%qlock%set ('remove_link')
326  next%prev => link%prev ! point backwards, possibly null
327  !call next%qlock%unset ('remove_link')
328  else
329  !call self%tail%qlock%set ('remove_link')
330  self%tail => link%prev ! new tail link
331  !call self%tail%qlock%unset ('remove_link')
332  end if
333  !-----------------------------------------------------------------------------
334  ! Check for head task
335  !-----------------------------------------------------------------------------
336  if (associated(link%prev)) then ! not head task
337  !call link%prev%qlock%set ('remove_link')
338  link%prev%next => next ! make prev point to next
339  !call link%prev%qlock%unset ('remove_link')
340  else ! link was head
341  !call self%head%qlock%set ('remove_link')
342  self%head => next ! make head point to next
343  !call self%head%qlock%unset ('remove_link')
344  end if
345  !-----------------------------------------------------------------------------
346  ! IMPORTANT: Some procedures depend on knowing the size of the list!
347  !-----------------------------------------------------------------------------
348  call self%update_counts (link, -1)
349  call trace_end
350 END SUBROUTINE remove_link
351 
352 !===============================================================================
353 !> Updating the status counts using the delta of a single task, which is faster
354 !> but less robust than counting all tasks. Should be thoroughly validated
355 !> before put into use, which is also not motivated unless count_status() takes
356 !> a noticeable time.
357 !===============================================================================
358 SUBROUTINE update_counts (self, link, delta)
359  class(list_t):: self
360  class(link_t), pointer:: link
361  integer:: delta, n, na
362  !-----------------------------------------------------------------------------
363  ! The counter updates below do not need to be atomic if the list is locked,
364  ! but the penalty is minimal
365  !-----------------------------------------------------------------------------
366  if (link%task%is_set (bits%virtual)) then
367  !$omp atomic update
368  self%nv = self%nv + delta
369  else if (link%task%is_set (bits%boundary)) then
370  !$omp atomic update
371  self%nb = self%nb + delta
372  !$omp atomic update
373  self%na = self%na + delta
374  else if (link%task%is_set (bits%internal)) then
375  !$omp atomic update
376  self%ni = self%ni + delta
377  !$omp atomic update
378  self%na = self%na + delta
379  end if
380  !-----------------------------------------------------------------------------
381  ! Validate using the incremental method, when task_list_t%verbose > 0
382  !-----------------------------------------------------------------------------
383  if (self%verbose > 0) then
384  n = self%n
385  na = self%na
386  call self%count_status ()
387  if (n /= self%n .or. na /= self%na) then
388  write (stdout,'(a,2(2i6,2x))') &
389  'list_t%update_counts WARNING: inconsistent counts', &
390  n, self%n, na, self%na
391  self%n = n
392  self%na = na
393  end if
394  end if
395 END SUBROUTINE update_counts
396 
397 !===============================================================================
398 !> Count the number of internal, boundary, virtual, frozen, and external tasks
399 !===============================================================================
400 SUBROUTINE count_status (self, label)
401  class(list_t):: self
402  character(len=*), optional:: label
403  class(link_t) , pointer:: link, nbor
404  integer:: n, ni, nb, nv, ne, nf
405  integer, save:: itimer=0
406  !.............................................................................
407  call trace%begin ('list_t%count_status', itimer=itimer)
408  call self%lock%set ('count_status')
409  nb=0; ni=0; nv=0; ne=0; nf=0; n = 0
410  link => self%head
411  do while (associated(link))
412  n = n+1
413  if (link%task%is_set(bits%frozen )) nf = nf+1
414  if (link%task%is_set(bits%boundary)) nb = nb+1
415  if (link%task%is_set(bits%virtual )) nv = nv+1
416  if (link%task%is_set(bits%internal)) ni = ni+1
417  if (link%task%is_set(bits%external)) ne = ne+1
418  link => link%next
419  end do
420  self%ni = ni
421  self%nb = nb
422  self%nv = nv
423  self%n = n
424  if (self%verbose > 0) then
425  if (present(label)) then
426  write (io_unit%log,'(1x,a,6i6,2x,a)') &
427  'list_t%count_status: n[ibveaf] = ', &
428  self%ni, self%nb, self%nv, self%ne, self%na, nf, label
429  else
430  write (io_unit%log,'(1x,a,6i6,2x,a)') &
431  'list_t%count_status: n[ibveaf] = ', &
432  self%ni, self%nb, self%nv, self%ne, self%na, nf
433  end if
434  end if
435  self%na = ni + nb - nf
436  !-----------------------------------------------------------------------------
437  call self%lock%unset ('count_status')
438  call trace%end (itimer)
439 END SUBROUTINE count_status
440 
441 !===============================================================================
442 !> Clear and then set the virtual and boundary bits
443 !===============================================================================
444 SUBROUTINE reset_status (self, check)
445  class(list_t):: self
446  logical, optional:: check
447  logical:: flags(3,2)
448  class(link_t), pointer:: link, nb
449  integer:: n_internal, n_boundary, n_virtual, n_external, n_frozen
450  !-----------------------------------------------------------------------------
451  ! Clear virtual and boundary bits
452  !-----------------------------------------------------------------------------
453  call trace%begin ('list_t%reset_status')
454  if (self%verbose > 1) &
455  write (io_unit%output,*) 'reset_status'
456  n_internal = 0
457  n_boundary = 0
458  n_virtual = 0
459  n_external = 0
460  n_frozen = 0
461  call self%lock%set ('reset_status')
462  !-----------------------------------------------------------------------------
463  ! Set bit
464  !-----------------------------------------------------------------------------
465  link => self%head
466  do while (associated(link))
467  if (link%task%rank == mpi%rank) then
468  !-------------------------------------------------------------------------
469  ! Assume internal -- this may change in set_status_and_flags()
470  !-------------------------------------------------------------------------
471  call link%task%set(bits%internal)
472  call link%task%clear(bits%boundary+bits%virtual+bits%external)
473  else
474  !-------------------------------------------------------------------------
475  ! Assume external -- this may change in set_status_and_flags()
476  !-------------------------------------------------------------------------
477  call link%task%set(bits%external)
478  call link%task%clear(bits%internal+bits%boundary+bits%virtual)
479  end if
480  !---------------------------------------------------------------------------
481  ! Run through nbors and set their status, optionally checking
482  !---------------------------------------------------------------------------
483  nb => link%nbor
484  do while (associated(nb))
485  if (present(check)) then
486  flags(:,1) = [nb%needed, nb%needs_me, nb%download]
487  end if
488  call link%task%nbor_relations (nb%task, nb%needed, nb%needs_me, nb%download)
489  if (present(check)) then
490  flags(:,2) = [nb%needed, nb%needs_me, nb%download]
491  if (.not.all(flags(:,2) .eqv. flags(:,1))) then
492  write (stderr,*) link%task%id, flags(:,1)
493  write (stderr,*) nb%task%id, flags(:,2)
494  flush (stderr)
495  call io%abort ('nbor_relations() returned value inconsistent with initial values')
496  end if
497  end if
498  nb => nb%next
499  end do
500  if (self%verbose > 1) &
501  write(io_unit%log,'(a,i6,i4,3x,3l1)') &
502  'reset_status: task, rank, IBV =', link%task%id, link%task%rank, &
503  link%task%is_set(bits%internal), &
504  link%task%is_set(bits%boundary), &
505  link%task%is_set(bits%virtual)
506  n_internal = n_internal + merge(1,0,link%task%is_set(bits%internal))
507  n_boundary = n_boundary + merge(1,0,link%task%is_set(bits%boundary))
508  n_virtual = n_virtual + merge(1,0,link%task%is_set(bits%virtual ))
509  n_external = n_external + merge(1,0,link%task%is_set(bits%external))
510  n_frozen = n_frozen + merge(1,0,link%task%is_set(bits%frozen ))
511  link => link%next
512  end do
513  self%nb = n_boundary
514  self%nv = n_virtual
515  self%ni = n_internal
516  !-----------------------------------------------------------------------------
517  ! Check that this does not interfere with end_time logics
518  !-----------------------------------------------------------------------------
519  self%na = n_internal+n_boundary-n_frozen
520  self%n = n_internal+n_boundary+n_virtual
521  if (self%verbose > 1) then
522  write (io_unit%log,'(5(a,i7,5x))') &
523  'rank:',mpi%rank, &
524  'n_internal:', n_internal, &
525  'n_boundary:', n_boundary, &
526  'n_virtual:' , n_virtual, &
527  'n_frozen:' , n_frozen
528  flush (io_unit%log)
529  end if
530  !-----------------------------------------------------------------------------
531  call self%lock%unset ('reset_status')
532  call trace%end()
533 END SUBROUTINE reset_status
534 
535 !===============================================================================
536 !> Check the nbor nbors for tasks ready to update, and the task itself
537 !===============================================================================
538 SUBROUTINE check_nbor_nbors (self, link)
539  class(list_t):: self
540  class(link_t), pointer:: link
541  !.............................................................................
542  class(link_t), pointer:: nbor
543  !-----------------------------------------------------------------------------
544  call trace%begin ('list_t%check_nbor_nbors')
545  nbor => link%nbor
546  do while (associated(nbor))
547  call self%check_nbors (nbor%link)
548  nbor => nbor%next
549  end do
550  call self%check_ready(link)
551  call trace%end()
552 END SUBROUTINE check_nbor_nbors
553 
554 !===============================================================================
555 !> Among a task and its neighbor tasks, move local tasks to ready_queue if they
556 !> are ready. If the ready bit is already set it means the patch has already
557 !> been put in the ready queue, and should not be checked again.
558 !>
559 !> The link pointer and everything it points to are private to this task, and
560 !> are not at this point in time accessible from the ready queue.
561 !>
562 !> To prevent AMR from changing the nbor lists while they are being accessed,
563 !> it is necessary (and sufficient) to lock the current link
564 !===============================================================================
565 SUBROUTINE check_nbors (self, link)
566  class(list_t):: self
567  class(link_t), pointer:: link
568  class(link_t), pointer:: nbor, nbors
569  class(task_t), pointer:: task
570  integer:: itr
571  integer, save:: itimer=0
572  !.............................................................................
573  task => link%task ! main task
574  itr = 1
575  call trace%begin('list_t%check_nbors', itimer=itimer)
576  !-----------------------------------------------------------------------------
577  ! Lock the other task link while accessing its nbor list, to prevent changes
578  !-----------------------------------------------------------------------------
579  nbor => link%nbor ! first nbor
580  do while (associated (nbor)) ! keep going until end
581  if (nbor%task%id==io%id_debug) print *, &
582  'task', task%id,' needs task ', nbor%task%id, nbor%needs_me
583  if (nbor%needs_me) then ! needs checking?
584  if (self%verbose > 1) &
585  write (io_unit%log,*) wallclock(), task%id, 'checks nbor', nbor%task%id
586  call self%check_ready (nbor%link, lock=.true.) ! pointer back
587  else
588  if (self%verbose > 1) &
589  write (io_unit%log,*) wallclock(), task%id, 'no need to check', nbor%task%id
590  end if
591  nbor => nbor%next ! next nbor
592  end do
593  !-----------------------------------------------------------------------------
594  ! Since we may be called upon to check nbors of a patch that should not itself
595  ! be checked, we check if any bits are set indicating this to be the case.
596  ! Note that this call to check_ready(), which concerns the task being updated,
597  ! does NOT need locking of the task link.
598  !-----------------------------------------------------------------------------
599  if (task%is_clear (bits%virtual+bits%external)) then
600  call self%check_ready (link) ! finally check link task
601  end if
602  call trace%end (itimer)
603  return
604 END SUBROUTINE check_nbors
605 
606 !===============================================================================
607 !> Among a task and its neighbor tasks, move local tasks to ready_queue if they
608 !> are ready, and send task data to non-locals.
609 !===============================================================================
610 SUBROUTINE check_all (self, repair)
611  class(list_t):: self
612  logical, optional:: repair
613  !.............................................................................
614  class(link_t), pointer:: link, nbor
615  integer, save:: itimer=0
616  integer:: nq
617  !-----------------------------------------------------------------------------
618  call trace%begin('list_t%check_all', itimer=itimer)
619  write (io_unit%log,*) wallclock(), 'check_all: phase 1', &
620  self%nq, associated(self%queue)
621  if (present(repair)) then
622  nullify (self%queue)
623  end if
624  nq = self%nq
625  link => self%head
626  do while (associated (link)) ! keep going until end
627  if (link%task%is_clear (bits%virtual+bits%external)) then
628  call self%check_ready (link) ! link ready?
629  if (self%nq > nq) then
630  write (io_unit%log,*) 'check_all found', link%task%id, link%task%time
631  nbor => link%nbor
632  do while (associated(nbor))
633  write (io_unit%log,*) 'nbor, time =', nbor%task%id, nbor%task%time, &
634  nbor%needed, nbor%needs_me, nbor%download
635  nbor => nbor%next
636  end do
637  nq = self%nq
638  end if
639  end if
640  link => link%next
641  end do
642  write (io_unit%log,*) wallclock(), 'check_all: done ', &
643  self%nq, associated(self%queue)
644  if (.not.associated(self%queue)) then
645  write (io_unit%log,*) 'check_all: phase 2, clearing ready bits'
646  link => self%head
647  do while (associated (link)) ! keep going until end
648  if (link%task%is_clear (bits%virtual+bits%external)) then
649  if (link%task%is_set (bits%ready)) then
650  write (io_unit%log,*) 'clearing ready bit on', link%task%id
651  call link%task%clear (bits%ready)
652  end if
653  call self%check_ready (link) ! link ready?
654  end if
655  link => link%next
656  end do
657  end if
658  call trace%end (itimer)
659 END SUBROUTINE check_all
660 
661 !===============================================================================
662 !> Resend all boundary patches
663 !===============================================================================
664 SUBROUTINE resend_bdry (self)
665  class(list_t):: self
666  class(link_t), pointer:: link
667  !.............................................................................
668  call trace%begin ('list_t%resend_bdry')
669  call self%lock%set ('resend_bdry')
670  link => self%head
671  do while (associated(link))
672  if (link%task%is_set (bits%boundary)) then
673  call self%send_to_vnbors (link)
674  end if
675  link => link%next
676  end do
677  call self%lock%unset ('resend_bdry')
678  call trace%end()
679 END SUBROUTINE resend_bdry
680 
681 !===============================================================================
682 !> Find the oldest local and virtual patches
683 !===============================================================================
684 SUBROUTINE check_oldest (self)
685  class(list_t):: self
686  class(link_t), pointer:: link, oldest, oldestv, oldestnb, nbor
687  real(8):: told, toldv
688  integer, save:: phase=1
689  !.............................................................................
690  call trace%begin('list_t%check_oldest ')
691  write (io_unit%log,*) 'check_oldest: phase 1', wallclock()
692  flush (io_unit%log)
693  told = 1d30
694  toldv = told
695  nullify (oldest)
696  nullify (oldestv)
697  nullify (oldestnb)
698  link => self%head
699  do while (associated (link)) ! keep going until end
700  if (link%task%is_clear (bits%frozen)) then
701  if (link%task%is_set(bits%virtual)) then
702  if (link%task%time < toldv) then
703  toldv = link%task%time
704  oldestv => link
705  end if
706  else
707  if (link%task%time < told) then
708  told = link%task%time
709  oldest => link
710  end if
711  end if
712  end if
713  link => link%next
714  end do
715  if (associated(oldest)) then
716  write (io_unit%log,*) 'oldest active is ', oldest%task%id , told, &
717  oldest%task%is_set(bits%ready), oldest%task%is_set(bits%busy)
718  else
719  write (stderr,*) 'STALLED: oldest active not found '
720  end if
721  flush (io_unit%log)
722  if (associated(oldestv)) then
723  write (io_unit%log,*) 'the oldest virtual task is ', oldestv%task%id, toldv
724  write (io_unit%log,*) 'it was last unpacked at', oldestv%task%unpack_time
725  flush (io_unit%log)
726  end if
727  if (associated(oldest)) then
728  nbor => oldest%nbor
729  write (io_unit%log,2) oldest%task%id, oldest%task%time, oldest%task%rank, oldest%task%level, &
730  oldest%needed, oldest%needs_me, oldest%task%is_set(bits%boundary), oldest%task%is_set(bits%virtual)
731  2 format(i6,f12.6,2i4,5l3,f12.6)
732  write (io_unit%log,*) 'nbor list:'
733  do while (associated(nbor))
734  write (io_unit%log,2) nbor%task%id, nbor%task%time, nbor%task%rank, nbor%task%level, &
735  nbor%needed, nbor%needs_me, nbor%task%is_set(bits%boundary), nbor%task%is_set(bits%virtual), &
736  nbor%task%is_ahead_of (oldest%task), nbor%task%unpack_time
737  nbor => nbor%next
738  end do
739  call oldest%task%clear (bits%ready+bits%busy)
740  call self%queue_by_time (oldest) ! add task to queue
741  write (io_unit%log,*) 'task queued'
742  flush (io_unit%log)
743  write (io_unit%log,*) 'check_oldest: phase 1 done'
744  flush (io_unit%log)
745  link => oldest
746  else
747  call trace%end ()
748  return
749  end if
750  !
751  if (phase > 1) then
752  if (associated(oldest)) then
753  write (io_unit%log,*) 'check_oldest: phase 2', associated(oldest)
754  write (io_unit%log,1) 'the oldest local task is', oldest%task%id, &
755  'rank', oldest%task%rank, oldest%task%time, &
756  oldest%task%is_set (bits%virtual), oldest%task%is_set (bits%external)
757  told = 1d30
758  link => oldest%nbor
759  do while (associated(link))
760  if (link%task%time < told) then
761  told =link%task%time
762  oldestnb => link
763  end if
764  write (io_unit%log,1) &
765  'nbor', link%task%id, 'rank', link%task%rank, link%task%time, &
766  link%task%is_set(bits%boundary), link%task%is_set (bits%virtual)
767  1 format (1x,a,i9,3x,a,i6,1p,e16.6,2l3)
768  link => link%next
769  end do
770  end if
771  write (io_unit%log,'(f12.6,2(3x,a,i8,i6,f12.6,3x,4l1),3i5)') wallclock(), &
772  'oldest:', oldest%task%id, oldest%task%rank, oldest%task%time, &
773  oldest%task%is_set(bits%internal), oldest%task%is_set(bits%boundary), &
774  oldest%task%is_set(bits%virtual), oldest%task%is_set(bits%external), &
775  'nbor:', oldestnb%task%id, oldestnb%task%rank, oldestnb%task%time, &
776  oldestnb%task%is_set(bits%internal), oldestnb%task%is_set(bits%boundary), &
777  oldestnb%task%is_set(bits%virtual), oldestnb%task%is_set(bits%external), &
778  self%nq, sent_list%n, recv_list%n, unpk_list%n
779  flush (io_unit%log)
780  if (associated(oldestv)) then
781  write (io_unit%log,*) 'the oldest virtual task is ', oldestv%task%id, oldestv%task%time
782  write (io_unit%log,*) 'check_oldest: phase 3', associated(oldestv)
783  write (io_unit%log,1) 'the oldest virtual task is', oldestv%task%id, &
784  'rank', oldestv%task%rank, oldestv%task%time, &
785  oldest%task%is_set (bits%virtual), oldest%task%is_set (bits%external)
786  link => oldestv%nbor
787  do while (associated(link))
788  write (io_unit%log,1) &
789  'nbor', link%task%id, 'rank', link%task%rank, link%task%time, &
790  link%task%is_set(bits%boundary), link%task%is_set (bits%virtual)
791  link => link%next
792  end do
793  end if
794  end if
795  call trace%end()
796 END SUBROUTINE check_oldest
797 
798 !===============================================================================
799 !> Check if the link task is ready, and if so, add it to self (the ready_queue)
800 !===============================================================================
801 SUBROUTINE check_ready (self, link, lock)
802  class(list_t):: self
803  class(link_t), pointer:: link
804  logical, optional:: lock
805  !.............................................................................
806  class(task_t), pointer:: task
807  class(link_t), pointer:: nbor
808  integer:: ignore, unit
809  integer, save:: itimer=0
810  logical:: ok
811  !-----------------------------------------------------------------------------
812  ! Since we may be called upon to check nbors of a patch that should not itself
813  ! be checked, we first check if any bits are set indicating this to be the case
814  !-----------------------------------------------------------------------------
815  task => link%task
816  ignore = bits%ready + bits%busy + bits%remove + bits%virtual + bits%frozen
817  if (task%is_set (ignore)) &
818  return
819  if (self%detailed_timer) &
820  call trace%begin ('list_t%check_ready', itimer=itimer)
821  !-----------------------------------------------------------------------------
822  ! Since only this thread can change the nbor list we are free to use it w/o
823  ! lock in cases where task is the active task, but in other cases we need to
824  ! lock the nbor list, so it isn't changed during the loop
825  !-----------------------------------------------------------------------------
826  if (present(lock) .and. omp_lock%links) &
827  call link%lock%set ('check_ready')
828  !-----------------------------------------------------------------------------
829  ! Tasks that are frozen or removed should not be used in the loop
830  !-----------------------------------------------------------------------------
831  ignore = bits%frozen + bits%remove
832  ok = .true.
833  nbor => link%nbor ! use original nbor list
834  do while (associated (nbor)) ! keep going until end
835  !---------------------------------------------------------------------------
836  ! Check nbors that are needed, not frozen, and have support
837  !---------------------------------------------------------------------------
838  if (nbor%needed .and. nbor%task%is_clear(ignore)) then
839  if (.not. nbor%task%is_ahead_of(task)) then ! and is ahead in time
840  !$omp atomic update
841  task%n_failed = task%n_failed+1
842  if (self%verbose > 1 .or. task%n_failed > 10000) then
843  if (self%verbose > 1) then
844  unit = io_unit%log
845  else
846  unit = stderr
847  end if
848 ! write (unit ,'(a,i4,1p,2(2x,a,i6,2x,a,g14.5),2(2x,a,g11.3))') &
849 ! 'list_t%check_ready: on rank', task%rank, &
850 ! 'task' , task%id , 'at t=', task%time, &
851 ! 'failed on nbor task', nbor%task%id, 'at t=', nbor%task%time, &
852 ! 'grace =', task%grace, ' dt=', nbor%task%dtime
853  end if
854  ok = .false.
855  exit
856  end if
857  end if
858  nbor => nbor%next ! next nbor
859  end do
860  if (present(lock) .and. omp_lock%links) &
861  call link%lock%unset ('check_ready')
862  if (ok) then
863  !$omp atomic write
864  task%n_failed = 0
865  call self%queue_by_time (link) ! add task to queue
866  end if
867  if (self%detailed_timer) &
868  call trace%end (itimer)
869 END SUBROUTINE check_ready
870 
871 !===============================================================================
872 !> Check if the number of nbor tasks is consistent
873 !===============================================================================
874 SUBROUTINE consistency (self, link, i)
875  class(list_t):: self
876  class(link_t), pointer:: link
877  class(task_t), pointer:: task
878  class(link_t), pointer:: nbor
879  integer:: n, id1=0, id2=0, i
880  !.............................................................................
881  if (self%debug < 2) return
882  call trace%begin ('list_t%consistency',1)
883  call self%lock%set ('list_t%consistency')
884  n = 0
885  task => link%task
886  nbor => link%nbor ! start on nbor list
887  if (associated(nbor)) id1 = nbor%task%id
888  do while (associated (nbor)) ! keep going until end
889  id2 = nbor%task%id
890  !write (io_unit%log,*)id2
891  nbor => nbor%next ! next nbor
892  n = n+1
893  end do
894  if (n /= link%task%n_nbors) then
895  write (io_unit%log,*) task%id, 'ERROR: n, nbors', n, id1, id2, i
896  else
897 ! write (io_unit%log,*) task%id, 'consistent', id1, id2, i
898  end if
899  call self%lock%unset ('list_t%consistency')
900  call trace%end()
901 END SUBROUTINE consistency
902 
903 !===============================================================================
904 !> Check if the number of queued tasks is consistent; i.e., of the number of
905 !> tasks linked at self%queue is the same as nq
906 !===============================================================================
907 SUBROUTINE qconsistency (self, ident)
908  class(list_t):: self
909  class(link_t), pointer:: link
910  class(task_t), pointer:: task, prev
911  class(link_t), pointer:: nbor
912  integer:: n, id1=0, id2=0, ident
913  real(8):: time
914  !.............................................................................
915  !if (self%debug < 1) return
916  call trace%begin ('list_t%qconsistency',4)
917  n = 0
918  if (io%verbose>2) then
919  link => self%head
920  do while (associated(link))
921  print *,'task list: id, time =', link%task%id, link%task%time
922  link => link%next
923  end do
924  print *,'nq:', self%nq
925  end if
926  link => self%queue
927  prev => link%task
928  do while (associated(link))
929  n = n+1
930  if (n > self%nq) exit
931  task => link%task
932  time = task%time
933  if (io%verbose>2) then
934  print *,'queue: n, id, time =', n, link%task%id, link%task%time, &
935  link%task%is_set (bits%ready), link%task%is_set (bits%busy)
936  end if
937  if (task%time < time) then
938  write (io_unit%log,'(a,i4,2(i9,1p,g15.6,3x))') &
939  'ERROR: queue out of order ident, prev%id, time, task%id', ident, prev%id, time, task%id, task%time
940  end if
941  prev => task
942  link => link%next_time
943  end do
944  if (n /= self%nq) then
945  write (io_unit%log,*) &
946  'ERROR: qconsistency ident, n, nq, ident =', ident, n, self%nq
947  n = 0
948  link => self%queue
949  do while (associated(link))
950  n = n+1
951  write (io_unit%log,*) n, link%task%id, link%task%time
952  link => link%next_time
953  if (n > self%nq+1) exit
954  end do
955  !call mpi%abort ('ERROR: qconsistency size')
956  end if
957  call trace%end()
958 END SUBROUTINE qconsistency
959 
960 !===============================================================================
961 !> Check if the number of queued tasks is consistent; i.e., of the number of
962 !> tasks linked at self%queue is the same as nq
963 !===============================================================================
964 SUBROUTINE aconsistency (self, ident)
965  class(list_t):: self
966  class(link_t), pointer:: link
967  class(task_t), pointer:: task, prev
968  class(link_t), pointer:: nbor
969  integer:: n, id1=0, id2=0, ident
970  real(8):: time
971  !.............................................................................
972  !if (self%debug < 1) return
973  call trace%begin ('list_t%aconsistency',4)
974  n = 0
975  link => self%active
976  if (associated(link)) then
977  time = link%task%atime
978  prev => link%task
979  link => link%next_active
980  n = n+1
981  end if
982  do while (associated(link))
983  task => link%task
984  n = n+1
985  if (n > self%nac) exit
986  if (task%atime < time) then
987  write (io_unit%log,'(a,i4,2(i9,1p,g15.6,3x),i5)') &
988  'ERROR: active queue out of order ident, prev%id, time, task%id, task%atime =', &
989  ident, prev%id, time, task%id, task%atime, self%nac
990  n = 0
991  link => self%active
992  do while (associated(link))
993  n = n+1
994  write (io_unit%log,*) n, link%task%id, link%task%atime
995  link => link%next_active
996  end do
997  exit
998  end if
999  time = task%atime
1000  prev => task
1001  link => link%next_active
1002  end do
1003  if (n /= self%nac) then
1004  write (io_unit%log,*) &
1005  'ERROR: aconsistency ident, n, nac =', ident, n, self%nac
1006  n = 0
1007  link => self%active
1008  do while (associated(link))
1009  n = n+1
1010  write (io_unit%log,*) n, link%task%id, link%task%atime
1011  link => link%next_active
1012  if (n > self%nac+2) exit
1013  end do
1014  end if
1015  call trace%end()
1016 END SUBROUTINE aconsistency
1017 
1018 !===============================================================================
1019 !> Set non_leaf flag on non-leaf patches
1020 !===============================================================================
1021 SUBROUTINE init_nonleaf (self)
1022  class(list_t):: self
1023  class(link_t), pointer:: link, scan
1024  class(link_t), pointer:: nbor
1025  integer:: n_add
1026  !.............................................................................
1027  return ! FIXME: remove call?
1028  call trace%begin('list_t%init_nonleaf ')
1029  call self%lock%set ('list_t%init_nonleaf')
1030  !-----------------------------------------------------------------------------
1031  ! First find and mark all non-leaf patches; these are patches that the link
1032  ! patch is fully contained inside -- we *assume* that finer level patches
1033  ! cover the whole parent patch!
1034  !-----------------------------------------------------------------------------
1035  link => self%head ! start checking links
1036  do while (associated(link))
1037  scan => self%head ! scan all links (FIXME)
1038  do while (associated(scan)) ! until end
1039  if (.not.associated(scan,link).and. & ! skip same patch
1040  .not.scan%task%is_set(bits%not_leaf)) then ! skip non-leaf patch
1041  if (scan%task%overlaps(link%task)) then ! if patches overlap
1042  if (scan%task%level==link%task%level-1) then ! and level is one below
1043  if (all(abs(scan%task%position-link%task%position) < scan%task%size*0.55d0)) then
1044  call scan%task%set (bits%not_leaf) ! mark non-leaf
1045  if (io%verbose>0) &
1046  write (io_unit%log,*) 'patch', link%task%id, ' parent:', scan%task%id
1047  end if
1048  end if
1049  end if
1050  end if
1051  scan => scan%next ! continue search
1052  end do
1053  link => link%next ! continue list
1054  end do
1055  call self%lock%unset ('list_t%init_nonleaf')
1056  call trace%end()
1057 END SUBROUTINE init_nonleaf
1058 
1059 !===============================================================================
1060 !> Add a new task (e.g. originating from AMR or load balancing)
1061 !===============================================================================
1062 SUBROUTINE add_new_link (self, link)
1063  class(list_t):: self
1064  class(link_t), pointer:: link
1065  !.............................................................................
1066  class(link_t), pointer:: nbor
1067  class(task_t), pointer:: task
1068  !-----------------------------------------------------------------------------
1069  call trace%begin ('list_t%add_new_task')
1070  !-----------------------------------------------------------------------------
1071  ! Make sure class(patch_t) tasks have a pointer back to the task link
1072  !-----------------------------------------------------------------------------
1073  task => link%task
1074  select type (task)
1075  class is (patch_t)
1076  task%link => link
1077  end select
1078  !-----------------------------------------------------------------------------
1079  ! Add an nbor list, and set status bits and flags. Add the task to the
1080  ! task_list, and if it turns out to be a boundary task, send it to rank nbors.
1081  !-----------------------------------------------------------------------------
1082  call self%init_nbors (link)
1083  call self%append (link)
1084  if (link%task%is_set (bits%boundary)) &
1085  call self%send_to_vnbors (link)
1086  !-----------------------------------------------------------------------------
1087  ! Increment task counters on the rank
1088  !-----------------------------------------------------------------------------
1089  if (task%is_clear (bits%virtual)) then
1090  !$omp atomic update
1091  timer%levelpop (task%level) = & ! increment level population count
1092  timer%levelpop (task%level) + 1
1093  end if
1094  !-----------------------------------------------------------------------------
1095  ! Set the flag init_nbors in the nbor tasks, so the next time they are being
1096  ! updated, they obtain updated nbor lists
1097  !-----------------------------------------------------------------------------
1098  call self%set_init_nbors (link)
1099  !-----------------------------------------------------------------------------
1100  ! Check the ready state of the new task, which already knows its nbors.
1101  ! The nbors do not have the new task in their nbor lists yet, but will check
1102  ! correspondingly when the new task has been added.
1103  !-----------------------------------------------------------------------------
1104  call self%check_ready (link)
1105  call trace%end ()
1106 END SUBROUTINE add_new_link
1107 
1108 !===============================================================================
1109 !> The very first initialization of a neighbor list. After this, it can be
1110 !> maintained without searching the entire task list on the rank.
1111 !===============================================================================
1112 SUBROUTINE init_nbors (self, link)
1113  class(list_t):: self
1114  class(link_t), pointer:: link
1115  !.............................................................................
1116  class(link_t), pointer:: nbor, scan, new_head, new_sort, old_head, old_sort, nbors
1117  integer:: n_add
1118  logical:: overlaps
1119  integer, save:: itimer=0
1120  !-----------------------------------------------------------------------------
1121  call trace%begin('list_t%init_nbors ', itimer=itimer)
1122  nullify (new_head, new_sort)
1123  n_add = 0
1124  call self%lock%set ('init_nbors')
1125  scan => self%head ! scan all links (FIXME)
1126  do while (associated(scan)) ! until end
1127  !---------------------------------------------------------------------------
1128  ! Add ranks that differ from the process rank to the list of load balance
1129  ! neighbors (each rank is only added once)
1130  !---------------------------------------------------------------------------
1131  !if (scan%task%rank /= mpi%rank .and. load_balance%on) then
1132  !call load_balance%add (scan%task%rank)
1133  !end if
1134  if (.not.associated(scan,link) .and. associated(scan%task)) then
1135  overlaps = scan%task%overlaps(link%task) &
1136  .and. scan%task%level <= link%task%level+1
1137  !-------------------------------------------------------------------------
1138  ! Potential nbors are other patches that overlap
1139  !-------------------------------------------------------------------------
1140  if (io%verbose>2 .or. link%task%id==io%id_debug) &
1141  write (io_unit%log, '(a,3(i6,2(3f7.4,1x)),2x,2l2)') &
1142  'init_nbors:', &
1143  link%task%id, link%task%position, link%task%size, &
1144  scan%task%id, scan%task%position, scan%task%size, &
1145  0, scan%task%distance (link%task), self%size, &
1146  overlaps, scan%task%is_set(bits%virtual)
1147  if (overlaps) then ! if its task overlaps
1148  !-----------------------------------------------------------------------
1149  ! The ifs select to add nbors to a patch if either the patch is not a
1150  ! virtual patch, or the nbor is not a virtual patch, thus excluding
1151  ! adding virtual patches as nbors to virtual patches
1152  !-----------------------------------------------------------------------
1153  allocate (nbor)
1154  nbor%task => scan%task ! pointer to task
1155  nbor%link => scan ! pointer to task link
1156  call link%add_nbor_by_rank (new_head, nbor) ! add in rank order
1157  call link%task%nbor_relations (nbor%task, nbor%needed, &
1158  nbor%needs_me, nbor%download) ! set flags and status bits
1159  n_add = n_add+1
1160  end if
1161  end if
1162  scan => scan%next ! continue search
1163  end do
1164  call self%lock%unset ('init_nbors')
1165  !-----------------------------------------------------------------------------
1166  if (self%verbose > 1 .or. link%task%id == io%id_debug) &
1167  write(io_unit%log,'(a,i6,i4,3x,3l1)') &
1168  'init_nbors: link, rank, IBV =', link%task%id, link%task%rank, &
1169  link%task%is_set(bits%internal), &
1170  link%task%is_set(bits%boundary), &
1171  link%task%is_set(bits%virtual)
1172  if (self%verbose > 2 .or. link%task%id == io%id_debug) &
1173  write (io_unit%log,*) &
1174  'added ', n_add, ' neighbors to', link%task%id, &
1175  link%task%is_set(bits%boundary), link%task%is_set(bits%virtual)
1176  !-----------------------------------------------------------------------------
1177  ! The links only needs to be locked the brief moment while the nbor pointers
1178  ! are changed.
1179  !-----------------------------------------------------------------------------
1180  nullify (new_sort)
1181  call link%sort_nbors_by_level (new_head, new_sort)
1182  call link%lock%set ('init_nbors')
1183  old_head => link%nbor
1184  old_sort => link%nbors_by_level
1185  link%nbor => new_head ! switch nbor list
1186  link%nbors_by_level => new_sort
1187  call link%lock%unset ('init_nbors')
1188  !-----------------------------------------------------------------------------
1189  ! Remove the two previour nbor lists
1190  !-----------------------------------------------------------------------------
1191  call link%remove_nbor_list (old_head) ! recover nbor list mem
1192  call link%remove_nbor_list (old_sort)
1193  !-----------------------------------------------------------------------------
1194  link%task%n_nbors = n_add ! remember n_nbors
1195  !-----------------------------------------------------------------------------
1196  ! Finally, clear the bits%init_nbors, in case it was set, and perform a
1197  ! check_nbors() on the link, in order to detect if an nbor that was waiting
1198  ! for an update on this task, which did not yet have it in its nbor list,
1199  ! now should be added to the ready queue.
1200  !-----------------------------------------------------------------------------
1201  call link%task%clear (bits%init_nbors)
1202  call trace%end (itimer)
1203 END SUBROUTINE init_nbors
1204 
1205 !===============================================================================
1206 !> Set bits%init_nbors on nbor tasks, if the task itself is new (and hence has
1207 !> a new nbor list, which should be matched by the nbors nbor list). For safe
1208 !> measure, we do this twice, even though one time should in principle be eno
1209 !===============================================================================
1210 SUBROUTINE set_init_nbors (self, link)
1211  class(list_t):: self
1212  class(link_t), pointer:: link, nbor
1213  !-----------------------------------------------------------------------------
1214  call trace%begin ('list_t%set_init_nbors')
1215  if (self%verbose > 0) &
1216  write (io_unit%log,*) &
1217  'set_init_nbors: id =', link%task%id
1218  nbor => link%nbor
1219  do while (associated(nbor))
1220  call nbor%task%set (bits%init_nbors)
1221  if (self%verbose > 1) then
1222  write (io_unit%log,*) 'bits%init_nbors set ', nbor%task%id, nbor%task%istep
1223  end if
1224  nbor => nbor%next
1225  end do
1226  call trace%end ()
1227 END SUBROUTINE set_init_nbors
1228 
1229 !===============================================================================
1230 !> Initialize the nbor lists of all nbors of link. Note that the nbor list of
1231 !> an nbor starts at nbor%link%nbor. To prevent other threads from modifying
1232 !> nbor list we need to loop over, we lock the link in the meantime, and to
1233 !> avoid a possible deadlock caused by having more than one lock active, we use
1234 !> a copy of the nbor list (cf. check_nbors())
1235 !===============================================================================
1236 SUBROUTINE init_nbor_nbors (self, link)
1237  class(list_t):: self
1238  class(link_t), pointer:: link
1239  !.............................................................................
1240  class(link_t), pointer:: nbor, nbors
1241  !-----------------------------------------------------------------------------
1242  ! Re-generate the link nbor list, to make sure a correct nbor loop below
1243  !-----------------------------------------------------------------------------
1244  call trace%begin('list_t%init_nbor_nbors')
1245  call link%lock%set ('init_nbnbs') ! lock the link while modifying
1246  call self%init_nbors (link) ! initialize for loop below
1247  call self%check_nbors (link)
1248  call link%copy_nbor_list (link%nbor, nbors) ! copy to temporary nbor list
1249  call link%lock%unset ('init_nbnbs') ! unlock the link
1250  !-----------------------------------------------------------------------------
1251  ! When resetting the nbors of nbors, notify rank nbors to do the same, by
1252  ! setting bits%init_nbors
1253  !-----------------------------------------------------------------------------
1254  nbor => nbors ! loop over temporary list
1255  do while (associated(nbor))
1256  if (nbor%link%task%is_set (bits%boundary)) &
1257  call nbor%link%task%set (bits%init_nbors)
1258  if (self%verbose > 0) &
1259  write (io_unit%log,*) 'init_nbor_nbors: id =', nbor%link%task%id
1260  call self%init_nbors (nbor%link)
1261  call self%check_nbors (nbor%link)
1262  nbor => nbor%next
1263  end do
1264  call link%remove_nbor_list (nbors) ! remove the nbor list copy
1265  call trace%end()
1266 END SUBROUTINE init_nbor_nbors
1267 
1268 !===============================================================================
1269 !> Check if the nbor list is still OK. A first test is to check if any of the
1270 !> nbors no longer overlaps. One could / should also check if any nbor nbor
1271 !> now overlaps.
1272 !===============================================================================
1273 SUBROUTINE refresh_nbors (self, link)
1274  class(list_t):: self
1275  class(link_t), pointer:: link
1276  !.............................................................................
1277  class(link_t), pointer:: nbor1, nbor2
1278  logical:: ok
1279  integer, save:: itimer=0
1280  !-----------------------------------------------------------------------------
1281  call trace%begin('list_t%refresh_nbors ', itimer=itimer)
1282  ok = .true.
1283  nbor1 => link%nbor
1284  do while (ok .and. associated(nbor1))
1285  ok = ok .and. nbor1%task%overlaps(link%task)
1286  !nbor2 => nbor1%link%nbor
1287  !do while (ok .and. associated(nbor2))
1288  ! ok = ok .and. .not. nbor2%task%overlaps(link%task)
1289  ! nbor2 => nbor2%next
1290  !end do
1291  nbor1 => nbor1%next
1292  end do
1293  if (.not.ok) then
1294  print *, link%task%id, 'list_t%refresh_nbors: init_nbprs'
1295  call self%init_nbors (link)
1296  call self%check_nbors (link)
1297  end if
1298  call trace%end (itimer)
1299 END SUBROUTINE refresh_nbors
1300 
1301 !===============================================================================
1302 !> The very first initialization of the neighbor lists. After this, it can be
1303 !> maintained without a global search. Distribute the work over available
1304 !> threads.
1305 !===============================================================================
1306 SUBROUTINE init_all_nbors (self)
1307  class(list_t):: self
1308  class(link_t), pointer:: link
1309  integer:: i, i0, i1, i2, n=0, nv=0, nb=0
1310  real(8):: time
1311  integer, save:: itimer=0
1312  !.............................................................................
1313  call trace%begin('list_t%init_all_nbors', itimer=itimer)
1314  call timer%tic (time)
1315  call init_nonleaf (self) ! set non_leaf bit
1316  call self%reset_status ! update task bits
1317  call self%count_status ! update task counts
1318  i2 = self%n/omp%nthreads + 1 ! tasks per thread
1319  !$omp parallel private(i,i0,i1,link) shared(self,i2) default(none)
1320  i0 = omp%thread*i2
1321  i1 = i0 + i2
1322  link => self%head ! start checking links
1323  i = 0 ! task counter
1324  do while (associated(link))
1325  if (i >= i0 .and. i<i1) then ! thread interval
1326  call self%init_nbors (link)
1327  !write (io_unit%log,*) 'init_all_nbors:', wallclock(), i0, i, i1
1328  end if
1329  i = i+1 ! increment counter
1330  link => link%next ! continue list
1331  end do
1332  !$omp end parallel
1333  !call self%reset_status ! set bits and count
1334  call timer%toc ('init_all_nbors', 1, time)
1335  call trace%end (itimer)
1336 END SUBROUTINE init_all_nbors
1337 
1338 !===============================================================================
1339 !> Remove parent patches, which have higher level patches inside
1340 !===============================================================================
1341 SUBROUTINE remove_parents (self)
1342  class(list_t):: self
1343  class(link_t), pointer:: link, scan, parent
1344  !.............................................................................
1345  call trace%begin('list_t%remove_parents')
1346  call self%lock%set ('remove_parent')
1347  link => self%head ! start checking links
1348  do while (associated(link))
1349  scan => self%head ! scan all links (FIXME)
1350  do while (associated(scan)) ! until end
1351  !---------------------------------------------------------------------
1352  ! Find parents, as patches on the next lower level that have the
1353  ! current patch position inside
1354  !---------------------------------------------------------------------
1355  nullify (parent)
1356  if (scan%task%level==link%task%level-1) then
1357  if (all(abs(scan%task%position-link%task%position) < scan%task%size*0.55d0)) then
1358  parent => scan
1359  end if
1360  end if
1361  scan => scan%next ! continue search
1362  if (associated(parent)) then
1363  write (*,'(2(1x,a,i9,5x))') &
1364  'removing parent task', parent%task%id, 'level', parent%task%level
1365  io%nwrite = io%nwrite - 1
1366  io%ntask = io%ntask - 1
1367  call parent%remove_from_nbors (parent)
1368  call self%remove_link (parent)
1369  end if
1370  end do
1371  link => link%next ! continue list
1372  end do
1373  call self%lock%unset ('remove_parent')
1374  call trace%end()
1375 END SUBROUTINE remove_parents
1376 
1377 !===============================================================================
1378 !===============================================================================
1379 FUNCTION intpos (self, task1, task2) RESULT (out)
1380  class(list_t):: self
1381  class(task_t), pointer:: task1
1382  class(task_t), pointer, optional:: task2
1383  real:: dist(3)
1384  integer:: out(3)
1385  !.............................................................................
1386  if (present(task2)) then
1387  dist = task1%position - task2%position
1388  dist = modulo(dist + 0.5*self%size, self%size) - 0.5*self%size
1389  out = floor(dist/task1%size + 0.5)
1390  else
1391  out = floor(task1%position/task1%size + 0.5)
1392  end if
1393 END FUNCTION intpos
1394 
1395 !===============================================================================
1396 !===============================================================================
1397 SUBROUTINE statistics (self)
1398  class(list_t):: self
1399  write (io_unit%log,'(a,3(5x,a,i6))') 'output_experiment:', &
1400  'n_ready:', self%n_ready, &
1401  'n_check:', self%n_check, &
1402  'n_nbor:' , self%n_nbor
1403 END SUBROUTINE statistics
1404 
1405 !===============================================================================
1406 !===============================================================================
1407 SUBROUTINE append_list (self, task_list)
1408  class(list_t):: self
1409  class(list_t), pointer:: task_list
1410  class(link_t), pointer:: link
1411  integer:: n
1412  !-----------------------------------------------------------------------------
1413  call trace%begin ('list_t%append_list')
1414  if (associated(self%tail)) then
1415  self%tail%next => task_list%head
1416  else
1417  self%head => task_list%head
1418  end if
1419  self%tail => task_list%tail
1420  self%n = self%n + task_list%n
1421  link => self%head
1422  n = 0
1423  do while (associated(link))
1424  link => link%next
1425  n = n+1
1426  end do
1427  if (n==self%n) then
1428  write (io_unit%log,*) trim(self%name)//': OK append of '//task_list%name
1429  else
1430  write (io_unit%log,*) trim(self%name)//': inconsistent append of '//task_list%name
1431  end if
1432  call trace%end()
1433 END SUBROUTINE append_list
1434 
1435 !===============================================================================
1436 SUBROUTINE remove_task (self, task)
1437  class(list_t):: self
1438  class(task_t), pointer:: task
1439  class(link_t), pointer:: link, prev
1440  !.............................................................................
1441  call trace%begin('list_t%remove_task')
1442  if (self%verbose > 2) &
1443  write (io%output,*) wallclock(),' thread',omp%thread,' wait for tasklist(3)'
1444  call self%lock%set ('remove_task')
1445  if (self%verbose > 2) &
1446  write (io%output,*) wallclock(),' thread',omp%thread,' locked tasklist(3)'
1447  link => self%head
1448  nullify (prev)
1449  do while (associated(link))
1450  if (link%task%id == task%id) then
1451  if (associated(prev)) then
1452  prev%next => link%next
1453  else
1454  self%head => link%next
1455  end if
1456  goto 9
1457  end if
1458  prev => link
1459  link => link%next
1460  end do
1461  print*,'ERROR: failed to remove from task list, task', task%id
1462 9 continue
1463  !-----------------------------------------------------------------------------
1464  ! IMPORTANT: Some procedures depend on knowing the size of the list!
1465  !-----------------------------------------------------------------------------
1466  call self%count_status()
1467  call self%lock%unset ('remove_task')
1468  if (self%verbose > 1) &
1469  write (io%output,*) wallclock(),' thread',omp%thread,' locked tasklist(3)'
1470  call trace%end()
1471 END SUBROUTINE remove_task
1472 
1473 !===============================================================================
1474 !> Add a link node into the queue chain, in ascending time order. When a link
1475 !> to a task with time larger than this%task%time is found, or the end of the
1476 !> list is fount, this link is inserted before.
1477 !>
1478 !> -> link0 -next_time-> link1 -next_time-> link2
1479 !> | | |
1480 !> t0 t1 t2
1481 !> tnew0 tnew1 tnew2
1482 !>
1483 !> In the sketch above the worst case is when two tasks are being queued at
1484 !> nearly the same time, into the same interval. This may be handled by having
1485 !> a lock on each link, so when the two tasks have identified the interval
1486 !> they should be in, which means they want to modify the next_time of the
1487 !> previous link (link0), they need to compete about who gets the lock. The one
1488 !> that succeeds modifies the next_time of link0 and its own next_time and then
1489 !> unlocks the link.
1490 !===============================================================================
1491 SUBROUTINE queue_by_time (self, this)
1492  class(list_t):: self
1493  class(link_t), pointer:: this
1494  logical:: error
1495  integer:: verbose
1496  !-----------------------------------------------------------------------------
1497  if (self%verbose > 2) &
1498  write (io%output,*) wallclock(),' thread',omp%thread,' waitfor tasklist(4)'
1499  call self%lock%set ('queue_by_time')
1500  if (self%verbose > 2) &
1501  write (io%output,*) wallclock(),' thread',omp%thread,' locked tasklist(4)'
1502  call real_queue_by_time (self, this, error)
1503  if (error) then
1504  io%do_trace = .true.
1505  verbose = io%verbose
1506  io%verbose = 4
1507  call self%qconsistency (999)
1508  call self%check_all (repair=.true.)
1509  call real_queue_by_time (self, this, error)
1510  if (error) then
1511  print *,'ERROR: queue repair failed'
1512  call io%abort()
1513  end if
1514  io%verbose = verbose
1515  io%do_trace = .false.
1516  end if
1517  call self%lock%unset ('queue_by_time')
1518  if (self%verbose > 2) &
1519  write (io%output,*) wallclock(),' thread',omp%thread,' unlocked tasklist(4)'
1520 END SUBROUTINE queue_by_time
1521 
1522 SUBROUTINE real_queue_by_time (self, this, error)
1523  class(list_t):: self
1524  class(link_t), pointer:: this
1525  logical:: error
1526  class(link_t), pointer:: next, prev, link
1527  class(task_t), pointer:: task
1528  integer, save:: itimer=0
1529  integer:: nit, i, j
1530  !-----------------------------------------------------------------------------
1531  ! Refuse to queue virtual tasks and finished tasks
1532  !-----------------------------------------------------------------------------
1533  error=.false.
1534  task => this%task
1535  if (task%is_set (bits%virtual + bits%frozen)) then
1536  write(stderr,*) mpi%rank, omp%thread, &
1537  'ERROR: tried to queue virtual or frozen task, VF =', &
1538  task%is_set(bits%virtual), task%is_set(bits%frozen)
1539  write(io_unit%log,*) mpi%rank, omp%thread, &
1540  'ERROR: tried to queue virtual or frozen task, VF =', &
1541  task%is_set(bits%virtual), task%is_set(bits%frozen)
1542  return
1543  end if
1544  if (task%time > io%end_time) then
1545  write (io_unit%mpi,*) mpi%rank, omp%thread, &
1546  'ERROR: tried to queue a task beyond end_time, VF =', &
1547  task%id, task%time, io%end_time, &
1548  task%is_set(bits%virtual), task%is_set(bits%frozen)
1549  link => self%head
1550  i = 0
1551  j = 0
1552  write (io_unit%mpi,*) 'task list:'
1553  do while (associated(link))
1554  i = i+1
1555  task => link%task
1556  if (task%time > io%end_time) then
1557  j = j+1
1558  write (io_unit%mpi,*) i, j, task%id, task%time, &
1559  task%is_set(bits%ready), task%is_set(bits%busy), &
1560  task%is_set(bits%boundary), task%is_set(bits%virtual), &
1561  task%is_set(bits%frozen), associated(link%next_time)
1562  end if
1563  link => link%next
1564  end do
1565  j = 0
1566  link => self%queue
1567  write (io_unit%mpi,*) 'queue:'
1568  do while (associated(link))
1569  task => link%task
1570  j = j+1
1571  write (io_unit%mpi,*) j, task%id, task%time, &
1572  task%is_set(bits%ready), task%is_set(bits%busy), &
1573  task%is_set(bits%boundary), task%is_set(bits%virtual), &
1574  task%is_set(bits%frozen), associated(link%next_time)
1575  link => link%next_time
1576  end do
1577  write (io_unit%mpi,*) 'na =', self%na
1578  call self%reset_status()
1579  call self%count_status('overtime')
1580  flush (io_unit%mpi)
1581  !call mpi%delay (5e3)
1582  !call io%abort ('task overtime')
1583  return
1584  end if
1585  !-----------------------------------------------------------------------------
1586  if (self%detailed_timer) &
1587  call trace%begin ('list_t%queue_by_time ', itimer=itimer)
1588  !$omp atomic
1589  mpi_mesg%n_ready = mpi_mesg%n_ready+1
1590  !-----------------------------------------------------------------------------
1591  ! Lock the queue while inserting the task, and mark the task as being in the
1592  ! ready queue (or executing -- the bit is cleared only after updating).
1593  ! Double check that no other thread has set the ready bit on this task.
1594  !-----------------------------------------------------------------------------
1595  if (task%is_clear (bits%ready+bits%busy)) then
1596  call task%set(bits%ready)
1597  if (io%verbose>1) &
1598  write (io_unit%log,'(f12.6,i5,2x,a,i9,a,f12.6,2x,a,i5)') wallclock(), omp%thread, &
1599  'list_t%queue_by_time: adding task', task%id, ' at time', task%time, &
1600  'nq', self%nq
1601  nullify (prev)
1602  next => self%queue
1603  nit = 0
1604  do while (associated(next))
1605  nit = nit+1
1606  if (nit > self%nq) then
1607  print *,'ERROR: hang in queue_by_time, nit =', nit
1608  error = .true.
1609  if (self%detailed_timer) &
1610  call trace%end (itimer)
1611  return
1612  end if
1613  if (associated(next%task, task)) then
1614  write (io_unit%log,*) omp_mythread, ' WARNING: task', task%id, ' is already in ready queue'
1615  go to 9
1616  else if (next%task%time > task%time) then
1617  this%next_time => next
1618  if (associated(prev)) then
1619  if (io%verbose > 1) then
1620  write (io_unit%output,'(i4,2x,a,i6,a,i6,i5,1p,3g16.6)') &
1621  omp%thread, 'task',task%id,' in ready queue between',prev%task%id,next%task%id, &
1622  prev%task%time, this%task%time, next%task%time
1623  write (io_unit%log,'(i4,2x,a,i6,a,i6,i5,1p,3g16.6)') &
1624  omp%thread, 'task',task%id,' in ready queue between',prev%task%id,next%task%id, &
1625  prev%task%time, this%task%time, next%task%time
1626  end if
1627  !call prev%qlock%set ('queue_by_time')
1628  prev%next_time => this
1629  !call prev%qlock%unset ('queue_by_time')
1630  else
1631  !call self%queue%qlock%set ('queue_by_time')
1632  self%queue => this
1633  !call self%queue%qlock%unset ('queue_by_time')
1634  if (io%verbose > 1) then
1635  write (io_unit%output,*) 'task',task%id,' at ready queue head',self%nq
1636  write (io_unit%log,*) 'task',task%id,' at ready queue head',self%nq
1637  end if
1638  end if
1639  !$omp atomic
1640  self%nq = self%nq+1
1641  go to 9
1642  end if
1643  prev => next
1644  next => next%next_time
1645  end do
1646  self%nq = self%nq+1
1647  if (associated(prev)) then
1648  if (io%verbose > 1) then
1649  write (io_unit%output,*) 'task',task%id,' in ready queue after',prev%task%id,self%nq
1650  write (io_unit%log,*) 'task',task%id,' in ready queue after',prev%task%id,self%nq
1651  end if
1652  prev%next_time => this
1653  else
1654  if (io%verbose > 1) then
1655  write (io_unit%output,*) 'task',task%id,' at ready queue head',self%nq
1656  write (io_unit%log,*) 'task',task%id,' at ready queue head',self%nq
1657  end if
1658  self%queue => this
1659  end if
1660  nullify (this%next_time)
1661  9 continue
1662  call self%remove_active (this) ! remove from active_queue
1663  !-----------------------------------------------------------------------------
1664  ! Since we clear the ready bit in task_list%update, outside of any critical
1665  ! region, it is perfectly normal that several threads compete about putting
1666  ! a given task onto the queue, so there is no reason to print a warning here
1667  !-----------------------------------------------------------------------------
1668  else if (io%verbose > 1) then
1669  if (io%verbose > 2) &
1670  print '(i6,a,i5,a,i6,1p,e15.6)', mpi%rank,' INFO: thread',omp%thread, &
1671  ' found ready bit in queue_by_time for task', task%id, task%time
1672  write (io_unit%log,*) wallclock(),' INFO: thread',omp%thread, &
1673  ' found ready bit in queue_by_time for task', task%id, task%time
1674  flush (io_unit%log)
1675  end if
1676  if (self%detailed_timer) &
1677  call trace%end (itimer)
1678 END SUBROUTINE real_queue_by_time
1679 
1680 !===============================================================================
1681 !===============================================================================
1682 SUBROUTINE queue_active (self, this)
1683  class(list_t):: self
1684  class(link_t), pointer:: this
1685  class(link_t), pointer:: next, prev
1686  class(task_t), pointer:: task
1687  integer, save:: itimer=0
1688  !-----------------------------------------------------------------------------
1689  if (self%detailed_timer) &
1690  call trace%begin ('list_t%queue_active ', itimer=itimer)
1691  !-----------------------------------------------------------------------------
1692  task => this%task
1693  task%atime = task%time
1694  if (io%verbose>1) &
1695  write (io_unit%log,'(1x,i4,2x,a,i9,a,f10.5,i5,l4)') omp%thread, &
1696  'list_t%queue_active: adding task', task%id, ' at time', task%atime, &
1697  self%nac, associated(self%active)
1698  nullify (prev)
1699  !$omp critical (active_cr)
1700  self%nac = self%nac+1
1701  next => self%active
1702  do while (associated(next))
1703  if (associated(next%task, task)) then
1704  write (io_unit%log,*) omp_mythread, ' WARNING: task', task%id, ' is already in active queue'
1705  exit
1706  else if (next%task%atime > task%atime) then
1707  this%next_active => next
1708  if (associated(prev)) then
1709  if (io%verbose > 1) then
1710  write (io_unit%log,*) 'task',task%id,' in active queue after',prev%task%id,self%nac
1711  end if
1712  prev%next_active => this
1713  else
1714  self%active => this
1715  if (io%verbose > 1) then
1716  write (io_unit%log,*) 'task',task%id,' at ready queue head',self%nac
1717  end if
1718  end if
1719  exit
1720  end if
1721  prev => next
1722  next => next%next_active
1723  end do
1724  if (.not.associated(next)) then
1725  if (associated(prev)) then
1726  if (io%verbose > 1) then
1727  write (io_unit%log,*) 'task',task%id,' at active queue tail',prev%task%id,self%nac
1728  end if
1729  prev%next_active => this
1730  else
1731  if (io%verbose > 1) then
1732  write (io_unit%log,*) 'task',task%id,' at active queue head',self%nac
1733  end if
1734  self%active => this
1735  end if
1736  nullify (this%next_active)
1737  end if
1738  call self%aconsistency (1)
1739  !$omp end critical (active_cr)
1740  if (self%detailed_timer) &
1741  call trace%end (itimer)
1742 END SUBROUTINE queue_active
1743 
1744 !===============================================================================
1745 !> Remove a task link from the active queue
1746 !===============================================================================
1747 SUBROUTINE remove_active (self, this)
1748  class(list_t):: self
1749  class(link_t), pointer:: this
1750  class(link_t), pointer:: next, prev
1751  !-----------------------------------------------------------------------------
1752  if (.not.associated(self%active)) return
1753  call trace%begin ('list_t%remove_active')
1754  !$omp critical (active_cr)
1755  nullify (prev)
1756  next => self%active
1757  do while (associated(next))
1758  if (associated(next%task, this%task)) then
1759  if (associated(prev)) then
1760  prev%next_active => next%next_active
1761  else
1762  self%active => next%next_active
1763  end if
1764  self%nac = self%nac-1
1765  if (io%verbose > 1) &
1766  write (io_unit%mpi,*) 'remove_active: id, nac =', this%task%id, self%nac
1767  exit
1768  end if
1769  prev => next
1770  next => next%next_active
1771  end do
1772  !$omp end critical (active_cr)
1773  call trace%end()
1774 END SUBROUTINE remove_active
1775 
1776 !===============================================================================
1777 !> Add a link node into the task list, in increasing quality order
1778 !===============================================================================
1779 SUBROUTINE add_by_quality (self, this)
1780  class(list_t):: self
1781  class(link_t), pointer:: this
1782  class(link_t), pointer:: next, prev
1783  !.............................................................................
1784  call trace%begin ('list_t%add_by_quality')
1785  call self%lock%set ('add_by_quality')
1786  nullify (prev)
1787  next => self%head
1788  do while (associated(next))
1789  if (associated(next%task,this%task)) then
1790  write (io_unit%log,*)'ERROR: Trying to add task', this%task%id, ' which is already there'
1791  goto 9
1792  end if
1793  if (next%task%level >= this%task%level) then
1794  !-------------------------------------------------------------------------
1795  ! Found a task we should be ahead of, so insert here
1796  !-------------------------------------------------------------------------
1797  this%next => next
1798  if (associated(prev)) then
1799  prev%next => this
1800  else
1801  self%head => this
1802  end if
1803  goto 9
1804  end if
1805  prev => next
1806  next => next%next
1807  end do
1808  !-----------------------------------------------------------------------------
1809  ! This task has the highest quality, so append it
1810  !-----------------------------------------------------------------------------
1811  if (associated(prev)) then
1812  prev%next => this
1813  else
1814  self%head => this
1815  end if
1816  nullify (this%next)
1817  self%tail => this
1818  self%n = self%n+1
1819 9 continue
1820  !-----------------------------------------------------------------------------
1821  ! IMPORTANT: Some procedures depend on knowing the size of the list!
1822  !-----------------------------------------------------------------------------
1823  call self%count_status()
1824  call self%lock%unset ('add_by_quality')
1825  call trace%end()
1826 END SUBROUTINE add_by_quality
1827 
1828 !===============================================================================
1829 !> Print a table with patch id, positions, and all neighbors
1830 !===============================================================================
1831 SUBROUTINE check_queue (self)
1832  class(list_t):: self
1833  class(link_t), pointer:: link
1834  real:: previous
1835  integer:: nq
1836  !.............................................................................
1837  call self%lock%set ('check_queue')
1838  link => self%queue
1839  previous = -1.0
1840  nq = 0
1841  do while (associated(link))
1842  if (link%task%time < previous) &
1843  write (io_unit%log,*)'ERROR: queue out of order at', link%task%id, link%task%time, &
1844  previous
1845  previous = link%task%time
1846  link => link%next_time
1847  nq = nq+1
1848  end do
1849  if (nq /= self%nq) write (io_unit%log,*)'ERROR: inconsistent number of tasks in queue', nq, self%nq
1850  call self%lock%unset ('check_queue')
1851 END SUBROUTINE check_queue
1852 
1853 !===============================================================================
1854 !> Insert 'new' task before 'old' task
1855 !===============================================================================
1856 SUBROUTINE insert (self, old, new)
1857  class(list_t):: self
1858  class(link_t), pointer:: old, new
1859  !.............................................................................
1860  call trace%begin ('list_t%insert '//self%name)
1861  call self%lock%set ('insert')
1862  if (associated(old%prev)) then
1863  new%prev => old%prev
1864  new%prev%next => new
1865  else
1866  self%head => new
1867  end if
1868  new%next => old
1869  new%next%prev => new
1870  !-----------------------------------------------------------------------------
1871  ! IMPORTANT: Some procedures depend on knowing the size of the list!
1872  !-----------------------------------------------------------------------------
1873  call self%count_status()
1874  call self%lock%unset ('insert')
1875  call trace%end()
1876 END SUBROUTINE insert
1877 
1878 !===============================================================================
1879 !> Print a table with task id, positions, and all neighbors
1880 !===============================================================================
1881 SUBROUTINE print_queue (self)
1882  class(list_t):: self
1883  class(link_t) , pointer:: link, nbor
1884  !.............................................................................
1885  call self%print_queue_times ('print_queue')
1886 END SUBROUTINE print_queue
1887 
1888 !===============================================================================
1889 !> Print the queue until and including a given task
1890 !===============================================================================
1891 SUBROUTINE print_queue_until (self, task)
1892  class(list_t):: self
1893  class(task_t):: task
1894  class(link_t) , pointer:: link, nbor
1895  !.............................................................................
1896  call self%lock%set ('print_queue_until')
1897  link => self%queue
1898  if (io%master) write (io_unit%log,*)'print_queue_until:'
1899  do while (associated(link))
1900  write (io_unit%log,'("queue: id, t0, t, t0-t =",i7,1p,3e12.3)') &
1901  link%task%id, task%time, link%task%time, task%time-link%task%time
1902  if (link%task%id==task%id) exit
1903  link => link%next_time
1904  end do
1905  call self%lock%unset ('print_queue_until')
1906 END SUBROUTINE print_queue_until
1907 
1908 !===============================================================================
1909 !> Print the queue
1910 !===============================================================================
1911 SUBROUTINE print_queue_times (self, label)
1912  class(list_t):: self
1913  character(len=*), optional:: label
1914  class(link_t) , pointer:: link, nbor
1915  !.............................................................................
1916  call self%lock%set ('print_queue_times')
1917  link => self%queue
1918  if (present(label)) then
1919  write (io_unit%log,*)'print_queue_times: '//trim(label)
1920  else
1921  write (io_unit%log,*)'print_queue_times:'
1922  end if
1923  do while (associated(link))
1924  write (io_unit%log,'("queue: id, time =",i7,1p,g14.5,l5)') link%task%id, link%task%time, &
1925  link%task%is_set(bits%ready)
1926  link => link%next_time
1927  end do
1928  call self%lock%unset ('print_queue_times')
1929  write (io_unit%log,*)'done'
1930 END SUBROUTINE print_queue_times
1931 
1932 !===============================================================================
1933 !===============================================================================
1934 SUBROUTINE print_list (self, label)
1935  class(list_t):: self
1936  character(len=*), optional:: label
1937  class(link_t), pointer:: link, prev
1938  class(task_t), pointer:: task
1939  character(len=32):: type
1940  !.............................................................................
1941  call trace%begin ('list_t%print')
1942  if (io%verbose>0) &
1943  write (io_unit%log,*) 'list_t%print: ',trim(self%name), self%n
1944  call self%lock%set ('print_list')
1945  link => self%head
1946  do while (associated(link))
1947  task => link%task
1948  type = ''
1949  !select type (task)
1950  !type is (task_t)
1951  ! type = 'task_t'
1952  !type is (particle_t)
1953  ! type = 'particle_t'
1954  !type is (planet_t)
1955  ! type = 'planet_t'
1956  !type is (star_t)
1957  ! type = 'star_t'
1958  !type is (pebble_t)
1959  ! type = 'pebble_t'
1960  !class is (patch_t)
1961  type = 'patch_t'
1962  !end select
1963  if (io%verbose>1) &
1964  write (io_unit%log,'(a,i9,f10.6,3x,3f10.6,i4,i6,l5,2x,a)') &
1965  'list_t%print:', task%id, task%time, task%position, task%status, &
1966  task%n_nbors, associated(link%nbor), type
1967  link => link%next
1968  end do
1969  call self%lock%unset ('print_list')
1970  call trace%end()
1971 END SUBROUTINE print_list
1972 
1973 !===============================================================================
1974 !> Print a table with task id, positions, and all neighbors
1975 !===============================================================================
1976 SUBROUTINE print_tasks (self, label)
1977  class(list_t):: self
1978  character(len=*), optional:: label
1979  class(link_t) , pointer:: link, nbor
1980  !.............................................................................
1981  !if (io%verbose < 2) return
1982  call trace%begin ('list_t%print_tasks')
1983  if (present(label)) then
1984  write (io_unit%log,'(a,i8,":")') 'task_list: '//trim(label), self%n
1985  else
1986  write (io_unit%log,'(a,i8,":")') 'task_list', self%n
1987  end if
1988  link => self%head
1989  do while (associated(link))
1990  if (io%verbose < 3) then
1991  write (io_unit%log,'(i8,g15.6,3x,3l1)') link%task%id, link%task%time, &
1992  link%task%is_set(bits%ready), link%task%is_set(bits%boundary), &
1993  link%task%is_set(bits%virtual)
1994  else
1995  call link%print_nbors
1996  end if
1997  link => link%next
1998  end do
1999  write (io_unit%log,*) ''
2000  call trace%end()
2001 END SUBROUTINE print_tasks
2002 
2003 !===============================================================================
2004 !> Initialize boundary bits, which may (optionally) be used to apply boundary
2005 !> conditions in experiment_t%update() procedures.
2006 !===============================================================================
2007 SUBROUTINE init_bdries (self)
2008  implicit none
2009  class(list_t):: self
2010  class(link_t), pointer:: link
2011  class(task_t), pointer:: patch
2012  real(8):: position(3), limit(3)
2013  !-----------------------------------------------------------------------------
2014  call trace%begin('task_list_t%init_bdries')
2015  link => self%head
2016  do while (associated(link))
2017  patch => link%task
2018  select type (patch)
2019  class is (patch_t)
2020  call patch%init_bdries
2021  end select
2022  self%lc = min(self%lc, link%task%position)
2023  self%uc = max(self%uc, link%task%position)
2024  self%llc = min(self%llc, link%task%position-0.5000000001_8*link%task%size)
2025  self%urc = max(self%urc, link%task%position+0.5000000001_8*link%task%size)
2026  link => link%next
2027  end do
2028  call trace%end()
2029 contains
2030  !-----------------------------------------------------------------------------
2031  function formatted (task) result (out)
2032  class(patch_t):: task
2033  character(len=6):: out
2034  out = '......'
2035  if (task%boundaries%is_set(bits%xl)) out(1:1) = 'T'
2036  if (task%boundaries%is_set(bits%xu)) out(2:2) = 'T'
2037  if (task%boundaries%is_set(bits%yl)) out(3:3) = 'T'
2038  if (task%boundaries%is_set(bits%yu)) out(4:4) = 'T'
2039  if (task%boundaries%is_set(bits%zl)) out(5:5) = 'T'
2040  if (task%boundaries%is_set(bits%zu)) out(6:6) = 'T'
2041  end function formatted
2042 END SUBROUTINE init_bdries
2043 
2044 !===============================================================================
2045 !> Check nbor list of self, verifying that %task and %link are consistent
2046 !===============================================================================
2047 SUBROUTINE check_nbor_list (link, label, verbose)
2048  class(link_t):: link
2049  character(len=*):: label
2050  integer, optional:: verbose
2051  class(link_t), pointer:: nbor
2052  class(task_t), pointer:: task
2053  !.............................................................................
2054  nbor => link%nbor
2055  do while (associated(nbor))
2056  if (associated(nbor%link)) then
2057  if (associated(nbor%task)) then
2058  task => nbor%task
2059  select type (task)
2060  class is (patch_t)
2061  if (.not.associated(task%link,nbor%link)) then
2062  print *, label, link%task%id, task%id, &
2063  'ERROR: nbor%link, task%link not associated'
2064  end if
2065  end select
2066  else
2067  print *, label, link%task%id, 'ERROR: nbor%task not associated'
2068  end if
2069  else
2070  print *, label, link%task%id, 'ERROR: nbor%link not associated'
2071  end if
2072  if (present(verbose)) then
2073  task => nbor%link%task
2074  select type (task)
2075  class is (patch_t)
2076  print *, link%task%id, nbor%task%id, task%id, associated(task%link,nbor%link)
2077  end select
2078  end if
2079  nbor => nbor%next
2080  end do
2081 END SUBROUTINE check_nbor_list
2082 
2083 !===============================================================================
2084 !> Give a task to another rank, making sure to update the nbor relations on
2085 !> this rank -- doing the same on the receiver rank is handled when unpacking,
2086 !===============================================================================
2087 SUBROUTINE give_to (self, link, rank)
2088  class(list_t):: self
2089  class(link_t), pointer:: link
2090  integer:: rank
2091  !-----------------------------------------------------------------------
2092  ! Set the bits and rank as they should be on the local rank, and mark
2093  ! the task with a swap_request bit
2094  !-----------------------------------------------------------------------
2095  call link%task%clear (bits%boundary)
2096  call link%task%set (bits%virtual+bits%swap_request)
2097  link%task%rank = rank
2098  !-----------------------------------------------------------------------
2099  ! With the new rank and status settings, update the status of the nbors
2100  ! and their nbor lists, which may all be affected by the change of rank.
2101  ! The ones that are new boundary tasks will be sent to nbor ranks by
2102  ! the update_nbor_status procedure. The rank itself will be sent by
2103  ! task_list_mod::update, and should not be sent from here.
2104  !-----------------------------------------------------------------------
2105  call self%update_nbor_status (link)
2106  if (io%verbose>1) &
2107  call self%test_nbor_status (link)
2108 END SUBROUTINE give_to
2109 
2110 !===============================================================================
2111 !> Update the nbor/vnbor status of a links nbor list, including the self link.
2112 !> If the boundary bit gets set, send the task to its vnbors. Before doing so,
2113 !> the nbor list of the nbors need to be updated, to reflect the new task rank
2114 !===============================================================================
2115 SUBROUTINE update_nbor_status (self, link)
2116  class(list_t):: self
2117  class(link_t), pointer:: link, nbor, this
2118  integer:: status
2119  !-----------------------------------------------------------------------------
2120  call trace%begin ('link_mod::update_nbor_status')
2121  call link%lock%set ('update_nbor_status')
2122  !
2123  nbor => link%nbor ! then help others
2124  do while (associated(nbor))
2125  if (io%verbose>1) &
2126  write (io_unit%log,*) 'update_nbor_status:', nbor%task%id, &
2127  nbor%task%is_set(bits%boundary), nbor%task%is_set(bits%virtual)
2128  status = nbor%task%status ! current status
2129  call nbor%link%remove_nbor (link) ! remove out-of-order
2130  allocate (this) ! new nbor link
2131  this%link => link ! attach task link
2132  this%task => link%task ! attach task
2133  call nbor%link%add_nbor_by_rank (link%nbor, this) ! add in rank order
2134  call self%update_status (nbor%link) ! new status
2135  if (nbor%task%is_set (bits%boundary).and.iand(status,bits%boundary)==0) then
2136  if (io%verbose>0) &
2137  write (io_unit%log,*) 'sending new bdry task',nbor%task%id,' to vnbors'
2138  if (.not.nbor%task%is_clear (bits%virtual)) then
2139  call nbor%task%clear (bits%virtual)
2140  write (io_unit%log,*) 'WARNING: needed to clear virtual bit(1)'
2141  end if
2142  call nbor%task%set (bits%swap_request)
2143  nbor%task%rank = mpi%rank ! since it's a boundary
2144  call self%send_to_vnbors (nbor%link)
2145  call nbor%link%task%clear (bits%swap_request + & ! sometimes set -- make
2146  bits%init_nbors) ! sure clear after send
2147  end if
2148  nbor => nbor%next ! check next nbor
2149  end do
2150  !
2151  if (io%verbose>0) &
2152  write (io_unit%log,*) 'update_nbor_status: id =', link%task%id, associated(link)
2153  status = link%task%status ! current status
2154  call self%update_status (link) ! new status
2155  if (link%task%is_set (bits%boundary).and.iand(status,bits%boundary)==0) then
2156  if (io%verbose>0) &
2157  write (io_unit%log,*) 'sending new bdry task',link%task%id,' to vnbors'
2158  if (.not.nbor%task%is_clear (bits%virtual)) then
2159  write (io_unit%log,*) 'WARNING: needed to clear virtual bit(2)'
2160  call nbor%task%clear (bits%virtual)
2161  end if
2162  call nbor%task%set (bits%swap_request)
2163  nbor%task%rank = mpi%rank ! since it's a boundary
2164  call self%send_to_vnbors (link) ! need to send ...
2165  call nbor%task%clear (bits%swap_request)
2166  end if
2167  call link%lock%unset ('update_nbor_status')
2168  call trace%end()
2169 END SUBROUTINE update_nbor_status
2170 
2171 !===============================================================================
2172 !> Update the nbor/vnbor status of a link. If self belongs to the local rank
2173 !> it is either an internal or boundary rank, if not it is virtual or external,
2174 !> depending on what nbors it has.
2175 !===============================================================================
2176 SUBROUTINE update_status (self, link)
2177  class(list_t):: self
2178  class(link_t):: link
2179  class(link_t), pointer:: nbor
2180  !-----------------------------------------------------------------------------
2181  call link%lock%set ('update_status')
2182  nbor => link%nbor
2183  if (link%task%rank == mpi%rank) then
2184  do while (associated(nbor))
2185  if (nbor%task%rank /= mpi%rank) then
2186  call link%task%set (bits%boundary)
2187  call link%task%clear (bits%internal+bits%external+bits%virtual)
2188  if (io%verbose>0) &
2189  write (io_unit%log,*) ' set to boundary id =', link%task%id
2190  return
2191  end if
2192  nbor => nbor%next
2193  end do
2194  call link%task%set (bits%internal)
2195  call link%task%clear (bits%boundary+bits%external+bits%virtual)
2196  if (io%verbose>0) &
2197  write (io_unit%log,*) ' set to internal id =', link%task%id, associated(link%nbor)
2198  else
2199  do while (associated(nbor))
2200  if (nbor%task%rank == mpi%rank) then
2201  call link%task%set (bits%virtual)
2202  call link%task%clear (bits%internal+bits%external+bits%boundary)
2203  if (io%verbose>0) &
2204  write (io_unit%log,*) ' set to virtual id =', link%task%id
2205  return
2206  end if
2207  nbor => nbor%next
2208  end do
2209  call link%task%set (bits%external)
2210  call link%task%clear (bits%internal+bits%boundary+bits%virtual)
2211  !call link%remove_from_nbors
2212  if (io%verbose>0) &
2213  write (io_unit%log,*) ' set to external id =', link%task%id
2214  end if
2215  call link%lock%unset ('update_status')
2216 END SUBROUTINE update_status
2217 
2218 !===============================================================================
2219 !> Send a package with updated time slice info to neighbor ranks.
2220 !>
2221 !> 1) assemble task data into a message
2222 !> 3) send the message to the nbors rank (only once per rank)
2223 !> 4) add the message to a sent_list
2224 !> 5) check the sent list for completed sends
2225 !> 6) remove the message from the sent_list and disassemble the message
2226 !===============================================================================
2227 SUBROUTINE send_to_vnbors (self, link)
2228  class(list_t):: self
2229  class(link_t), pointer:: link
2230  class(link_t), pointer:: nbor
2231  class(task_t), pointer:: task, nbtask
2232  class(mesg_t), pointer:: mesg
2233  character(len=24):: label
2234  integer:: ierr, rank, tag, seq
2235  integer, save:: itimer=0
2236  !.............................................................................
2237  task => link%task
2238  if (task%is_clear (bits%boundary)) then
2239  write (stdout,*) mpi%rank, omp%thread, &
2240  'ERROR: trying to send non-boundary task', task%id
2241  return
2242  end if
2243  if (task%is_set (bits%virtual)) then
2244  write (stdout,*) mpi%rank, omp%thread, &
2245  'ERROR: trying to send virtual task', task%id
2246  return
2247  end if
2248  !-----------------------------------------------------------------------------
2249  call trace%begin('list_t%send_to_vnbors', itimer=itimer)
2250  if (task%id == io%id_debug) &
2251  write(io%output,*) 'DBG link_t%send_to_vnbors: id, rank =', &
2252  task%id, mpi%rank
2253  !-----------------------------------------------------------------------------
2254  ! Pack task into mesg
2255  !-----------------------------------------------------------------------------
2256  select type (task)
2257  class is (patch_t)
2258  call task%pack (mesg)
2259  end select
2260  if (mpi_mesg%uniq_mesg .and. mpi_mesg%tag_type == 1) then
2261  task%seq = task%seq + 1
2262  tag = mod(task%seq,100) + 100*task%id
2263  else
2264  tag = task%id
2265  end if
2266  !-----------------------------------------------------------------------------
2267  ! Same id for all steps, and increment mesg%seq only once, for all ranks
2268  ! If the task is new, force the sequence number to be 1
2269  !-----------------------------------------------------------------------------
2270  mesg%id = task%id
2271  !-----------------------------------------------------------------------------
2272  ! Send to all unique ranks in the nbor list (which is sorted by rank).
2273  !-----------------------------------------------------------------------------
2274  rank = -1 ! previous rank
2275  nbor => link%nbor ! first nbor
2276  do while (associated(nbor)) ! until end
2277  nbtask => nbor%task
2278  if (nbtask%rank/=mpi%rank .and. nbtask%rank/=rank) then ! new rank?
2279  rank = nbtask%rank ! current rank
2280  if (mpi_mesg%uniq_mesg .and. mpi_mesg%tag_type == 2) then
2281  !$omp atomic capture
2282  sequence(rank) = sequence(rank)+1
2283  seq = sequence(rank)
2284  !$omp end atomic
2285  tag = mod(seq,100) + 100*task%id
2286  end if
2287  if (task%logging > 1) then
2288  write (label,'(a,i4,i8)') 'vnbor ', rank, tag
2289  call task%log (label)
2290  end if
2291  call mesg%send (rank, tag=tag) ! send it
2292  !write (stdout,*) 'send: ', task%id, mesg%id, tag, mesg%tag, rank
2293  if (self%verbose > 0 .or. task%id == io%id_debug) then
2294  write (io_unit%mpi,'(f12.6,2x,a,i9,1p,e18.6,i9,2x,a,i5,2x,5l1)') &
2295  wallclock(), 'send_to_vnbors: sent', &
2296  mesg%id, task%time, tag, 'to', rank, &
2297  task%is_set (bits%internal), &
2298  task%is_set (bits%boundary), &
2299  task%is_set (bits%virtual), &
2300  task%is_set (bits%external), &
2301  task%is_set (bits%swap_request)
2302  flush (io_unit%mpi)
2303  end if
2304  end if
2305  nbor => nbor%next ! next nbor
2306  end do
2307  !-----------------------------------------------------------------------------
2308  ! Add the message to mesg%sent_list, if it was sent, and clear the swap bit
2309  !-----------------------------------------------------------------------------
2310  call mpi_mesg%sent (mesg) ! add to sent_list
2311  call task%clear (bits%swap_request)
2312  !-----------------------------------------------------------------------------
2313  call trace%end (itimer)
2314 END SUBROUTINE send_to_vnbors
2315 
2316 !===============================================================================
2317 !> Check the status of a link and all nbor links
2318 !===============================================================================
2319 SUBROUTINE test_nbor_status (self, link)
2320  class(list_t):: self
2321  class(link_t):: link
2322  class(link_t), pointer:: nbor
2323  !.............................................................................
2324  write (io_unit%log,*) 'test_nbor_status: link', link%task%id
2325  call test_status (link)
2326  nbor => link%nbor
2327  do while (associated(nbor))
2328  call test_status (nbor%link)
2329  nbor => nbor%next
2330  end do
2331 END SUBROUTINE test_nbor_status
2332 
2333 !===============================================================================
2334 !> Check the status of a link task with respect to its nbors
2335 !===============================================================================
2336 SUBROUTINE test_status (self)
2337  class(link_t):: self
2338  class(link_t), pointer:: nbor
2339  integer:: status, test, n
2340  !.............................................................................
2341  write (io_unit%log,*) 'test_status: link', self%task%id
2342  if (self%task%rank == mpi%rank) then
2343  status = bits%internal
2344  else
2345  status = bits%external
2346  end if
2347  nbor => self%nbor
2348  do while (associated(nbor))
2349  if (nbor%task%rank /= mpi%rank .and. status==bits%internal) then
2350  status = bits%boundary
2351  else if (nbor%task%rank == mpi%rank .and. status==bits%external) then
2352  status = bits%virtual
2353  end if
2354  nbor => nbor%next
2355  end do
2356  test = iand(self%task%status,bits%internal+bits%boundary+bits%virtual+bits%external)
2357  if (status /= test) then
2358  write (io_unit%log,'(i9,2x,a,i6,2(2x,4l1))') self%task%id, 'inconsistent:', &
2359  self%task%rank, &
2360  self%task%is_set(bits%internal), &
2361  self%task%is_set(bits%boundary), &
2362  self%task%is_set(bits%virtual) , &
2363  self%task%is_set(bits%external), &
2364  iand(status,bits%internal)/=0, &
2365  iand(status,bits%boundary)/=0, &
2366  iand(status,bits%virtual )/=0, &
2367  iand(status,bits%external)/=0
2368  n = 0
2369  nbor => self%nbor
2370  do while (associated(nbor))
2371  n = n+1
2372  write (io_unit%log,'(2i9,2x,a,2i6)') self%task%id, nbor%task%id, &
2373  'inconsistent ranks', self%task%rank, nbor%task%rank
2374  nbor => nbor%next
2375  end do
2376  if (n==0) write (io_unit%log,*) self%task%id, ' inconsistent bits: has no nbors'
2377  end if
2378 END SUBROUTINE test_status
2379 
2380 !===============================================================================
2381 !> Dump nbor relations, flags and status bits for the list
2382 !===============================================================================
2383 SUBROUTINE info (self)
2384  class(list_t):: self
2385  class(link_t), pointer:: link, nbor
2386  class(task_t), pointer:: task, nbtask
2387  !-----------------------------------------------------------------------------
2388  call trace%begin ('task_list_t%info')
2389  link => self%head
2390  do while (associated(link))
2391  call link%info
2392  link => link%next
2393  end do
2394  call trace%end()
2395 END SUBROUTINE info
2396 
2397 END MODULE list_mod
Each thread uses a private timer data type, with arrays for start time and total time for each regist...
Definition: timer_mod.f90:11
Support tic/toc timing, as in MATLAB, and accurate wallclock() function. The timing is generally much...
Module with list handling for generic class task_t objects.
Definition: list_mod.f90:4
Template module for patches, which adds pointers to memory and mesh, and number of dimensions and var...
Definition: patch_mod.f90:6
Definition: io_mod.f90:4
The lock module uses nested locks, to allow versatile use of locks, where a procedure may want to mak...
Template module for tasks.
Definition: task_mod.f90:4