QR_MUMPS
qrm_common_mod.F90
Go to the documentation of this file.
1 !! ##############################################################################################
2 !!
3 !! Copyright 2012 CNRS, INPT
4 !!
5 !! This file is part of qr_mumps.
6 !!
7 !! qr_mumps is free software: you can redistribute it and/or modify
8 !! it under the terms of the GNU Lesser General Public License as
9 !! published by the Free Software Foundation, either version 3 of
10 !! the License, or (at your option) any later version.
11 !!
12 !! qr_mumps is distributed in the hope that it will be useful,
13 !! but WITHOUT ANY WARRANTY; without even the implied warranty of
14 !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 !! GNU Lesser General Public License for more details.
16 !!
17 !! You can find a copy of the GNU Lesser General Public License
18 !! in the qr_mumps/doc directory.
19 !!
20 !! ##############################################################################################
21 
22 
23 !! ##############################################################################################
33 
34 #include "qrm_common.h"
35 
38  use qrm_const_mod
39 
40 
44  interface qrm_print_tree
45  subroutine qrm_print_nsteps_tree(file, adata, weight)
47  type(qrm_adata_type) :: adata
48  character :: file*(*)
49  real(kind(1.d0)), optional :: weight(:)
50  end subroutine qrm_print_nsteps_tree
51  subroutine qrm_print_elim_tree(file, parent, n)
52  character :: file*(*)
53  integer :: parent(:)
54  integer :: n
55  end subroutine qrm_print_elim_tree
56  subroutine qrm_print_asm_tree(file, parent, rc, n)
57  character :: file*(*)
58  integer :: parent(:), rc(:)
59  integer :: n
60  end subroutine qrm_print_asm_tree
61  end interface
62 
64  interface qrm_postorder
65  subroutine qrm_postorder(parent, n, porder, weight)
66  integer :: n
67  integer :: parent(:), porder(:)
68  integer, optional :: weight(:)
69  end subroutine qrm_postorder
70  end interface
71 
73  interface
74  subroutine qrm_amalg_tree(n, parent, rc, porder, nvar, min_var, fill_thresh)
75  integer :: n, min_var
76  integer :: parent(:), rc(:), porder(:), nvar(:)
77  real(kind(1.d0)) :: fill_thresh
78  end subroutine qrm_amalg_tree
79  end interface
80 
83  subroutine qrm_compress_data(adata, porder, parent, rc, stair, n)
85  use qrm_adata_mod
86  type(qrm_adata_type) :: adata
87  integer :: porder(:), parent(:), rc(:), stair(:)
88  integer :: n
89  end subroutine qrm_compress_data
90  end interface
91 
93  interface qrm_reorder_tree
94  subroutine qrm_reorder_tree(adata)
96  type(qrm_adata_type) :: adata
97  end subroutine qrm_reorder_tree
98  end interface
99 
104  interface qrm_prnt_array
105  subroutine qrm_prnt_iarray(a, lab, unit)
106  integer :: a(:)
107  character :: lab*(*)
108  integer, optional :: unit
109  end subroutine qrm_prnt_iarray
110  subroutine qrm_prnt_sarray(a, lab, digits, unit)
111  real(kind(1.e0)) :: a(:)
112  character :: lab*(*)
113  integer :: digits
114  integer, optional :: unit
115  end subroutine qrm_prnt_sarray
116  subroutine qrm_prnt_darray(a, lab, digits, unit)
117  real(kind(1.d0)) :: a(:)
118  character :: lab*(*)
119  integer :: digits
120  integer, optional :: unit
121  end subroutine qrm_prnt_darray
122  subroutine qrm_prnt_carray(a, lab, digits, unit)
123  complex(kind(1.e0)) :: a(:)
124  character :: lab*(*)
125  integer :: digits
126  integer, optional :: unit
127  end subroutine qrm_prnt_carray
128  subroutine qrm_prnt_zarray(a, lab, digits, unit)
129  complex(kind(1.d0)) :: a(:)
130  character :: lab*(*)
131  integer :: digits
132  integer, optional :: unit
133  end subroutine qrm_prnt_zarray
134  end interface
135 
137  interface qrm_swtime
138  function qrm_swtime() bind(c, name='qrm_swtime')
139  real(kind(1.d0)) :: qrm_swtime
140  end function qrm_swtime
141  end interface
142 
144  interface qrm_uwtime
145  function qrm_uwtime() bind(c, name='qrm_uwtime')
146  real(kind(1.d0)) :: qrm_uwtime
147  end function qrm_uwtime
148  end interface
149 
151  interface qrm_msleep
152  subroutine qrm_msleep(n) bind(c, name='qrm_msleep')
153  use iso_c_binding
154  integer(c_int), value :: n
155  end subroutine qrm_msleep
156  end interface
157 
158 
160  interface qrm_check_cperm
161  subroutine qrm_check_cperm(cperm, n)
162  integer :: cperm(:)
163  integer :: n
164  end subroutine qrm_check_cperm
165  end interface
166 
167 #if defined(hwloc)
168 
169  interface
170  subroutine qrm_hwloc_bind(id) bind(c, name='qrm_hwloc_bind')
171  use iso_c_binding
172  integer(c_int), value :: id
173  end subroutine qrm_hwloc_bind
174  end interface
175 
177  interface
178  subroutine qrm_hwloc_topo(nnodes, topo) bind(c, name='qrm_hwloc_topo')
179  use iso_c_binding
180  integer(c_int) :: topo(nnodes,*)
181  end subroutine qrm_hwloc_topo
182  end interface
183 
185  interface
186  subroutine qrm_hwloc_info(ncores, nnodes, cnode) bind(c, name='qrm_hwloc_info')
187  use iso_c_binding
188  integer(c_int) :: ncores, nnodes, cnode
189  end subroutine qrm_hwloc_info
190  end interface
191 
192 #endif
193 
194 
198  interface qrm_count_flops
199  module procedure qrm_count_realflops, qrm_count_pureflops
200  ! function qrm_count_realflops(m, n, k, op)
201  ! real(kind(1.d0)) :: qrm_count_realflops
202  ! integer :: m, k, n
203  ! character :: op*(*)
204  ! end function qrm_count_realflops
205  ! function qrm_count_pureflops(stair, n, j, nb)
206  ! implicit none
207  ! real(kind(1.d0)) :: qrm_count_pureflops
208  ! integer :: stair(:)
209  ! integer :: n, nb, j
210  ! end function qrm_count_pureflops
211  end interface
212 
213  interface qrm_set
214  module procedure qrm_gseti
215  end interface
216 
217  interface qrm_get
218  module procedure qrm_ggeti, qrm_ggetii
219  end interface
220 
221 contains
222 
223  subroutine qrm_gseti(string, ival)
225  use qrm_error_mod
226  use qrm_mem_mod
227  implicit none
228 
229  character(len=*) :: string
230  integer :: ival
231 
232  character(len=len(string)) :: istring
233  ! error management
234  integer :: err_act
235  character(len=*), parameter :: name='qrm_gseti'
236 
237  call qrm_err_act_save(err_act)
238 
239  istring = qrm_str_tolower(string)
240 
241  if(index(istring,'qrm_eunit') .eq. 1) then
242  qrm_eunit = ival
243  else if(index(istring,'qrm_ounit') .eq. 1) then
244  qrm_ounit = ival
245  else if(index(istring,'qrm_dunit') .eq. 1) then
246  qrm_dunit = ival
247  else if(index(istring,'qrm_max_mem') .eq. 1) then
248  qrm_max_mem(0) = ival
249  else if(index(istring,'qrm_tot_mem') .eq. 1) then
250  qrm_tot_mem(0) = ival
251  else if(index(istring,'qrm_exact_mem') .eq. 1) then
252  qrm_exact_mem = ival
253  else if(index(istring,'qrm_error_action') .eq. 1) then
254  call qrm_err_act_set(ival)
255  else
256  call qrm_err_push(23, name, aed=string)
257  goto 9999
258  end if
259 
260  call qrm_err_act_restore(err_act)
261  return
262 
263 9999 continue ! error management
264  call qrm_err_act_restore(err_act)
265  if(err_act .eq. qrm_abort_) then
266  call qrm_err_check()
267  end if
268 
269 
270  return
271  end subroutine qrm_gseti
272 
273  subroutine qrm_ggeti(string, ival)
275  use qrm_error_mod
276  implicit none
277 
278  character(len=*) :: string
279  integer :: ival
280 
281  integer(kind=8) :: iival
282 
283  ! error management
284  integer :: err_act
285  character(len=*), parameter :: name='qrm_ggeti'
286 
287  call qrm_err_act_save(err_act)
288 
289  call qrm_ggetii(string, iival)
290  __qrm_check_ret(name,'qrm_ggetii',9999)
291  ival = iival
292 
293  call qrm_err_act_restore(err_act)
294  return
295 
296 9999 continue ! error management
297  call qrm_err_act_restore(err_act)
298  if(err_act .eq. qrm_abort_) then
299  call qrm_err_check()
300  end if
301 
302  return
303  end subroutine qrm_ggeti
304 
305  subroutine qrm_ggetii(string, iival)
307  use qrm_error_mod
308  use qrm_mem_mod
309  implicit none
310 
311  character(len=*) :: string
312  integer(kind=8) :: iival
313 
314  character(len=len(string)) :: istring
315  integer :: ival
316  ! error management
317  integer :: err_act
318  character(len=*), parameter :: name='qrm_ggeti'
319 
320  call qrm_err_act_save(err_act)
321 
322  istring = qrm_str_tolower(string)
323 
324  if(index(istring,'qrm_error') .eq. 1) then
325  call qrm_err_get(ival)
326  iival = ival
327  else if(index(istring,'qrm_error_action') .eq. 1) then
328  iival = qrm_err_act
329  else if(index(istring,'qrm_max_mem') .eq. 1) then
330  iival = qrm_max_mem(0)
331  else if(index(istring,'qrm_tot_mem') .eq. 1) then
332  iival = qrm_tot_mem(0)
333  else if(index(istring,'qrm_ounit') .eq. 1) then
334  iival = qrm_ounit
335  else if(index(istring,'qrm_eunit') .eq. 1) then
336  iival = qrm_eunit
337  else if(index(istring,'qrm_dunit') .eq. 1) then
338  iival = qrm_dunit
339  else
340  call qrm_err_push(23, name, aed=string)
341  goto 9999
342  end if
343 
344  call qrm_err_act_restore(err_act)
345  return
346 
347 9999 continue ! error management
348  call qrm_err_act_restore(err_act)
349  if(err_act .eq. qrm_abort_) then
350  call qrm_err_check()
351  end if
352 
353  return
354  end subroutine qrm_ggetii
355 
356 
357 
358 
360 function qrm_count_realflops(m, n, k, op)
361  real(kind(1.d0)) :: qrm_count_realflops
362  integer :: m, k, n
363  character :: op*(*)
364 
365  real(kind(1.d0)) :: rk, rm, rn
366 
367  rk = real(k, kind(1.d0))
368  rm = real(m, kind(1.d0))
369  rn = real(n, kind(1.d0))
370 
371  qrm_count_realflops = 0.d0
372 
373  select case(op)
374  case('panel')
375  if(m .gt. k) then
376  qrm_count_realflops = 2*rk*rk*(rm-rk/3.d0)
377  else
378  qrm_count_realflops = 2*rm*rm*(rk-rm/3.d0)
379  end if
380  ! qrm_count_realflops = qrm_count_realflops + rk*rk*(rm-rk/3.d0)
381  case('update')
382  qrm_count_realflops = rk*rn*(4*rm-rk)
383  end select
384 
385  return
386 end function qrm_count_realflops
387 
388 
391 function qrm_count_pureflops(stair, n, j, nb)
392  implicit none
393  real(kind(1.d0)) :: qrm_count_pureflops
394  integer :: stair(:)
395  integer :: n, nb, j
396 
397  integer :: i
398 
399  qrm_count_pureflops=0
400  do i=j, min(j+nb-1, n)
401  qrm_count_pureflops = qrm_count_pureflops+(stair(i)-i+1)*(3 + 4*(n-i))
402  end do
403 
404  return
405 end function qrm_count_pureflops
406 
407 end module qrm_common_mod
subroutine qrm_print_asm_tree(file, parent, rc, n)
This subroutine prints on a file the assembly tree described by a parent and a postorder arrays...
Generic interface for the ::qrm_hwloc_topo routine.
subroutine qrm_err_get(info)
This subroutine return the code of the first error on the stack or zero if the stack is empty...
subroutine qrm_print_nsteps_tree(file, adata, weight)
prints an assembly tree in compressed format
subroutine qrm_prnt_carray(a, lab, digits, unit)
Generic interface for the ::qrm_prnt_iarray, ::qrm_prnt_sarray, ::qrm_prnt_darray, ::qrm_prnt_carray and ::qrm_prnt_zarray routines.
This module contains the interfaces of all non-typed routines.
Generic interface for the ::qrm_msleep routine.
subroutine qrm_err_push(code, sub, ied, aed)
This subroutine pushes an error on top of the stack.
integer, save qrm_dunit
The unit for debug messages.
real(kind(1.d0)) function qrm_count_pureflops(stair, n, j, nb)
Used for counting the real flops (i.e., it ignores the zeros included by the blocking) ...
Generic interface for the ::qrm_uwtime routine.
integer, save qrm_eunit
The unit for error messages.
integer qrm_exact_mem
This module contains the definition of the analysis data type.
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.
subroutine qrm_ggetii(string, iival)
subroutine qrm_print_elim_tree(file, parent, n)
This subroutine prints on a file the elimination tree described by a parent array. The tree is written in "dot" format (see graphviz)
This module contains all the error management routines and data.
Generic interface for the ::qrm_reorder_tree routine.
subroutine qrm_prnt_darray(a, lab, digits, unit)
integer(kind=8), dimension(0:qrm_maxthreads-1) qrm_max_mem
a counter to keep track of the peak memory, per thread
subroutine qrm_prnt_iarray(a, lab, unit)
Generic interface for the ::qrm_swtime routine.
The main data type for the analysis phase.
integer, parameter qrm_abort_
Possible actions to be performed upon detection of an error.
Generic interface for the ::qrm_compress_data routine.
integer qrm_err_act
Default action.
Generic interface for the ::qrm_amalg_tree routine.
subroutine qrm_err_check()
This subroutine checks the errors stack. If something is found all the entries in the stack are poppe...
integer(kind=8), dimension(0:qrm_maxthreads-1) qrm_tot_mem
a counter to keep track of the currently allocated memory, per thread
Generic interface for the ::qrm_hwloc_info routine.
subroutine qrm_prnt_sarray(a, lab, digits, unit)
Generic interface for the ::qrm_count_realflops ::qrm_count_pureflops.
Generic interface for the ::qrm_print_nsteps_tree, ::qrm_print_elim_tree and ::qrm_print_asm_tree rou...
Generic interface for the ::qrm_postorder routine.
Generic interface for the ::qrm_hwloc_bind routine.
integer, save qrm_ounit
The unit for output messages.
subroutine qrm_prnt_zarray(a, lab, digits, unit)
Generic interface for the ::qrm_check_cperm routine.
real(kind(1.d0)) function qrm_count_realflops(m, n, k, op)
Used for counting the actual flops.
subroutine qrm_gseti(string, ival)
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38
subroutine qrm_err_act_set(err_act)
Sets the default error action.
This module contains various string handling routines.
subroutine qrm_err_act_restore(err_act)
Restores the value of the qrm_err_act variable.
subroutine qrm_ggeti(string, ival)