QR_MUMPS
qrm_print_tree.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 
45 subroutine qrm_print_elim_tree(file, parent, n)
46 
47  character :: file*(*)
48  integer :: parent(:)
49  integer :: n
50 
51  integer :: i
52 
53  open(4, file=file, action='write')
54 
55  write(4,'("graph G {")')
56  write(4,'("node [color=black,")')
57  write(4,'("fillcolor=white,")')
58  write(4,'("shape=circle,")')
59  write(4,'("style=filled")')
60  write(4,'("];")')
61  do i=1, n
62  if(parent(i) .ne. 0) then
63  write(4,'(i6," -- ",i6)')parent(i),i
64  else
65  write(4,'(i6)')i
66  end if
67  end do
68 
69  write(4,'("}")')
70 
71  close(4)
72  return
73 
74 end subroutine qrm_print_elim_tree
75 
76 
77 
78 
97 subroutine qrm_print_asm_tree(file, parent, rc, n)
98 
99  use qrm_mem_mod
100 
101  character :: file*(*)
102  integer :: parent(:), rc(:)
103  integer :: n
104 
105  integer :: i, maxnv
106  integer, allocatable :: nvar(:)
107  real(kind(1.d0)) :: s
108 
109  call qrm_aalloc(nvar, n)
110  nvar=1
111  maxnv = 0
112  do i=1, n
113  if(parent(i) .lt. 0) then
114  nvar(-parent(i)) = nvar(-parent(i))+1
115  if(nvar(-parent(i)) .gt. maxnv) maxnv=nvar(-parent(i))
116  end if
117  end do
118 
119 
120  open(4, file=file, action='write')
121 
122  write(4,'("graph G {")')
123  write(4,'("node [color=black,")')
124  write(4,'("fillcolor=white,")')
125  write(4,'("shape=circle,")')
126  write(4,'("style=filled")')
127  write(4,'("];")')
128 
129  do i=1, n
130  if(parent(i) .ge. 0) then
131  s = 5.d0/(real(maxnv,kind(1.d0))/real(nvar(i),kind(1.d0)) )
132  write(4,'("node",i6.6,"[label="" node:",i6,"\n rc:",i6,'//&
133  &'"\n nv:",i6,""", shape=circle, height=",i6.6,",'// &
134  &'fontsize=",i6.6,"];")')i,i,rc(i),nvar(i),nvar(i),20*nvar(i)
135  end if
136  end do
137 
138  do i=1, n
139  if(parent(i) .gt. 0) then
140  write(4,'("node",i6.6," -- node",i6.6)')parent(i),i
141  else if (parent(i) .eq. 0) then
142  write(4,'("node",i6.6)')i
143  end if
144  end do
145 
146  write(4,'("}")')
147 
148  close(4)
149 
150  call qrm_adealloc(nvar)
151 
152  return
153 
154 end subroutine qrm_print_asm_tree
155 
156 
163 subroutine qrm_print_nsteps_tree(file, adata, weight)
165  use qrm_adata_mod
166  use qrm_mem_mod
167  implicit none
168 
169  type(qrm_adata_type) :: adata
170  character :: file*(*)
171  real(kind(1.d0)), optional :: weight(:)
172 
173  integer :: i, p, n
174  integer, allocatable :: tmp(:)
175 
176 
177  open(4, file=file, action='write')
178 
179  write(4,'("graph G {")')
180  write(4,'("node [color=black,")')
181  write(4,'("fillcolor=white,")')
182  write(4,'("shape=circle,")')
183  write(4,'("style=filled")')
184  write(4,'("];")')
185 
186  if(present(weight)) then
187  ! if(adata%nnodes .gt. 1000) then
188  if(.false.) then
189  call qrm_aalloc(tmp, adata%nnodes+1, lbnd=0)
190  tmp = 0
191 
192  do i=1, adata%nnodes
193  if(adata%small(i) .eq. 1) then
194  write(4,'("node",i6.6,"[fillcolor=red, label="" node:",i6,"\n m:",i6,"\n n:",'//&
195  &'i6,"\n np:",i6,"\n pv:",i6,"\n w:",f4.1,"% ""];")')i,i,adata%nfrows(i),adata%rc(i),&
196  &adata%cp_ptr(i+1)-adata%cp_ptr(i),adata%cperm(adata%cp_ptr(i)), weight(i)
197  tmp(i) = 1
198  p = i
199  do
200  p = adata%parent(p)
201  if((p .gt. 0) .and. (tmp(p) .eq.0)) then
202  write(4,'("node",i6.6,"[label="" node:",i6,"\n m:",i6,"\n n:",'//&
203  &'i6,"\n np:",i6,"\n pv:",i6,"\n w:",f4.1,"% ""];")')p,p,adata%nfrows(p),adata%rc(p),&
204  &adata%cp_ptr(p+1)-adata%cp_ptr(p),adata%cperm(adata%cp_ptr(p)), weight(p)
205  tmp(p) = 1
206  else
207  exit
208  end if
209  end do
210  end if
211  end do
212 
213  tmp = 0
214 
215  do i=1, adata%nnodes
216  if(adata%small(i) .eq. 1) then
217  n = i
218  do
219  p = adata%parent(n)
220  if (tmp(n) .eq. 0) then
221  if(p .gt. 0) then
222  write(4,'("node",i6.6," -- node",i6.6)')p,n
223  tmp(n) = 1
224  n = p
225  else if (p .eq. 0) then
226  write(4,'("node",i6.6)')n
227  tmp(n) = 1
228  exit
229  end if
230  else
231  exit
232  end if
233  end do
234  end if
235  end do
236 
237  call qrm_adealloc(tmp)
238 
239  else
240  do i=1, adata%nnodes
241  if(adata%small(i) .eq. 1) then
242  write(4,'("node",i6.6,"[fillcolor=red, label="" node:",i6,"\n m:",i6,"\n n:",'//&
243  &'i6,"\n np:",i6,"\n pv:",i6,"\n w:",f4.1,"% ""];")')i,i,adata%nfrows(i),adata%rc(i),&
244  &adata%cp_ptr(i+1)-adata%cp_ptr(i),adata%cperm(adata%cp_ptr(i)),weight(i)
245  else
246  write(4,'("node",i6.6,"[label="" node:",i6,"\n m:",i6,"\n n:",'//&
247  &'i6,"\n np:",i6,"\n pv:",i6,"\n w:",f4.1,"% ""];")')i,i,adata%nfrows(i),adata%rc(i),&
248  &adata%cp_ptr(i+1)-adata%cp_ptr(i),adata%cperm(adata%cp_ptr(i)),weight(i)
249  end if
250  end do
251 
252  do i=1, adata%nnodes
253  if(adata%parent(i) .gt. 0) then
254  write(4,'("node",i6.6," -- node",i6.6)')adata%parent(i),i
255  else if (adata%parent(i) .eq. 0) then
256  write(4,'("node",i6.6)')i
257  end if
258  end do
259  end if
260  else
261  if(adata%nnodes .gt. 10) then
262  call qrm_aalloc(tmp, adata%nnodes+1, lbnd=0)
263  tmp = 0
264 
265  do i=1, adata%nnodes
266  if(adata%small(i) .eq. 1) then
267  write(4,'("node",i6.6,"[fillcolor=red, label="" node:",i6,"\n m:",i6,"\n n:",'//&
268  &'i6,"\n np:",i6,"\n pv:",i6,"""];")')i,i,adata%nfrows(i),adata%rc(i),&
269  &adata%cp_ptr(i+1)-adata%cp_ptr(i),adata%cperm(adata%cp_ptr(i))
270  tmp(i) = 1
271  p = i
272  do
273  p = adata%parent(p)
274  if((p .gt. 0) .and. (tmp(p) .eq.0)) then
275  write(4,'("node",i6.6,"[label="" node:",i6,"\n m:",i6,"\n n:",'//&
276  &'i6,"\n np:",i6,"\n pv:",i6,"""];")')p,p,adata%nfrows(p),adata%rc(p),&
277  &adata%cp_ptr(p+1)-adata%cp_ptr(p),adata%cperm(adata%cp_ptr(p))
278  tmp(p) = 1
279  else
280  exit
281  end if
282  end do
283  end if
284  end do
285 
286  tmp = 0
287 
288  do i=1, adata%nnodes
289  if(adata%small(i) .eq. 1) then
290  n = i
291  do
292  p = adata%parent(n)
293  if (tmp(n) .eq. 0) then
294  if(p .gt. 0) then
295  write(4,'("node",i6.6," -- node",i6.6)')p,n
296  tmp(n) = 1
297  n = p
298  else if (p .eq. 0) then
299  write(4,'("node",i6.6)')n
300  tmp(n) = 1
301  exit
302  end if
303  else
304  exit
305  end if
306  end do
307  end if
308  end do
309 
310  call qrm_adealloc(tmp)
311 
312  else
313  do i=1, adata%nnodes
314  if(adata%small(i) .eq. 1) then
315  write(4,'("node",i6.6,"[fillcolor=red, label="" node:",i6,"\n m:",i6,"\n n:",'//&
316  &'i6,"\n np:",i6,"\n pv:",i6,"""];")')i,i,adata%nfrows(i),adata%rc(i),&
317  &adata%cp_ptr(i+1)-adata%cp_ptr(i),adata%cperm(adata%cp_ptr(i))
318  else
319  write(4,'("node",i6.6,"[label="" node:",i6,"\n m:",i6,"\n n:",'//&
320  &'i6,"\n np:",i6,"\n pv:",i6,"""];")')i,i,adata%nfrows(i),adata%rc(i),&
321  &adata%cp_ptr(i+1)-adata%cp_ptr(i),adata%cperm(adata%cp_ptr(i))
322  end if
323  end do
324 
325  do i=1, adata%nnodes
326  if(adata%parent(i) .gt. 0) then
327  write(4,'("node",i6.6," -- node",i6.6)')adata%parent(i),i
328  else if (adata%parent(i) .eq. 0) then
329  write(4,'("node",i6.6)')i
330  end if
331  end do
332  end if
333  end if
334 
335  write(4,'("}")')
336 
337  close(4)
338 
339 
340  return
341 
342 end subroutine qrm_print_nsteps_tree
343 
344 
345 
346 
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_adealloc_i, qrm_adealloc_2i, qrm_adealloc_s, qrm_adealloc_2s, qrm_adealloc_3s, qrm_adealloc_d, qrm_adealloc_2d, qrm_adealloc_3d, qrm_adealloc_c, qrm_adealloc_2c, qrm_adealloc_3c, qrm_adealloc_z, qrm_adealloc_2z, qrm_adealloc_3z, routines.
subroutine qrm_print_nsteps_tree(file, adata, weight)
prints an assembly tree in compressed format
This module contains the definition of the analysis data type.
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)
The main data type for the analysis phase.
Generic interface for the qrm_aalloc_i, qrm_aalloc_2i, qrm_aalloc_s, qrm_aalloc_2s, qrm_aalloc_3s, qrm_aalloc_d, qrm_aalloc_2d, qrm_aalloc_3d, qrm_aalloc_c, qrm_aalloc_2c, qrm_aalloc_3c, qrm_aalloc_z, qrm_aalloc_2z, qrm_aalloc_3z, routines.
Definition: qrm_mem_mod.F90:78
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38