QR_MUMPS
dqrm_compute_graph.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 
35 #include "qrm_common.h"
36 
38 
46 subroutine dqrm_compute_graph(qrm_mat, graph)
47 
48  use dqrm_spmat_mod
49  use qrm_error_mod
50  use qrm_mem_mod
52 
53  implicit none
54 
55  type(dqrm_spmat_type) :: qrm_mat
56  type(dqrm_spmat_type), intent(out) :: graph
57 
58  integer :: i, j, pnt, savepnt, dups, ii
59  integer, allocatable :: dupsmap(:)
60  ! error management
61  integer :: err_act
62  character(len=*), parameter :: name='qrm_compute_graph'
63 
64  call qrm_err_act_save(err_act)
65 
66  call qrm_aalloc(dupsmap, qrm_mat%m)
67 
68  ! At this moment, if the matrix is centralized, we just need to
69  ! convert the structure from COO to CSC. Otherwise, we also need
70  ! to gather it from the other procs.
71 
72  call dqrm_spmat_convert(qrm_mat, graph, 'csc', values=.false.)
73  __qrm_check_ret(name,'alloc/convert',9999)
74 
75  ! make a pass in order to remove duplicates
76  dupsmap = 0
77  pnt = 0
78  dups = 0
79  savepnt = 1
80  do j=1, graph%n
81  do ii=graph%jptr(j), graph%jptr(j+1)-1
82  i = graph%irn(ii)
83  if(dupsmap(i) .eq. j) then
84  ! duplicate entry, skip
85  dups = dups+1
86  else
87  ! flag the entry as visited
88  dupsmap(i) = j
89  pnt = pnt+1
90  graph%irn(pnt) = i
91  end if
92  end do
93  graph%jptr(j) = savepnt
94  savepnt = pnt+1
95  end do
96  graph%jptr(graph%n+1)=savepnt
97  graph%nz = pnt
98  graph%icntl = qrm_mat%icntl
99  graph%rcntl = qrm_mat%rcntl
100 
101  __qrm_prnt_dbg('("Number of duplicates in the matrix: ",i10)')dups
102 
103 
104  call qrm_adealloc(dupsmap)
105  __qrm_check_ret(name,'qrm_pdealloc',9999)
106 
107  call qrm_err_act_restore(err_act)
108  return
109 
110 9999 continue ! error management
111  call qrm_err_act_restore(err_act)
112  if(err_act .eq. qrm_abort_) then
113  call qrm_err_check()
114  end if
115 
116  return
117 
118 end subroutine dqrm_compute_graph
Generic interface for the qrm_adealloc_i, qrm_adealloc_2i, qrm_adealloc_s, qrm_adealloc_2s, qrm_adealloc_3s, qrm_adealloc_d, qrm_adealloc_2d, qrm_adealloc_3d, qrm_adealloc_c, qrm_adealloc_2c, qrm_adealloc_3c, qrm_adealloc_z, qrm_adealloc_2z, qrm_adealloc_3z, routines.
This module contains the generic interfaces for all the analysis routines.
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.
This module contains all the error management routines and data.
This module contains the definition of the basic sparse matrix type and of the associated methods...
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
integer, parameter qrm_abort_
Possible actions to be performed upon detection of an error.
subroutine qrm_err_check()
This subroutine checks the errors stack. If something is found all the entries in the stack are poppe...
This type defines the data structure used to store a matrix.
subroutine dqrm_compute_graph(qrm_mat, graph)
Computes the adjacency graph of a matrix.
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38
subroutine qrm_err_act_restore(err_act)
Restores the value of the qrm_err_act variable.
subroutine dqrm_spmat_convert(in_mat, out_mat, fmt, values)
This subroutine converts an input matrix into a different storage format. Optionally the values may b...