QR_MUMPS
qrm_ata_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 
43 subroutine _qrm_ata_graph(g_csc, ata_graph)
45  use qrm_error_mod
46  use qrm_mem_mod
47 
48  implicit none
49 
50  type(_qrm_spmat_type), intent(in) :: g_csc
51  type(_qrm_spmat_type), intent(out) :: ata_graph
52 
53  type(_qrm_spmat_type) :: g_csr
54  integer :: i, j, row1, col1, row2, col2
55  integer, allocatable :: mark(:)
56  ! error management
57  integer :: err_act
58  character(len=*), parameter :: name='qrm_ata_graph'
59 
60  call qrm_err_act_save(err_act)
61 
62  call _qrm_spmat_convert(g_csc, g_csr, 'csr', .false.)
63 
64  call qrm_palloc(ata_graph%iptr,g_csc%n+2)
65  ata_graph%iptr = 0
66  ata_graph%iptr(1:2) = 1
67 
68  call qrm_aalloc(mark, g_csc%n)
69  __qrm_check_ret(name,'convert/alloc',9999)
70 
71  mark = 0
72 
73  do col1=1, g_csc%n
74  ! loop over all the columns of A
75  do i=g_csc%jptr(col1), g_csc%jptr(col1+1)-1
76  ! for each nnz in the column, we go through the corresponding
77  ! row and count all the i,j pairs
78  row1 = g_csc%irn(i)
79 
80  do j=g_csr%iptr(row1), g_csr%iptr(row1+1)-1
81  col2 = g_csr%jcn(j)
82  ! now, the element (col1,col2) is present in A'*A. Check if
83  ! it was already found, otherwise add it.
84  ! skip the diagonal
85  if(col1 .eq. col2) cycle
86  if(mark(col2) .lt. col1) then
87  mark(col2) = col1
88  ata_graph%iptr(col1+2) = ata_graph%iptr(col1+2)+1
89  end if
90  end do
91  end do
92  end do
93 
94  do i=3, g_csc%n+2
95  ata_graph%iptr(i) = ata_graph%iptr(i)+ata_graph%iptr(i-1)
96  end do
97 
98  ata_graph%nz = ata_graph%iptr(g_csc%n+2)
99  call qrm_palloc(ata_graph%jcn,ata_graph%nz)
100  __qrm_check_ret(name,'qrm_palloc',9999)
101 
102  ! second pass to fill up
103  mark = 0
104  do col1=1, g_csc%n
105  ! loop over all the columns of A
106  do i=g_csc%jptr(col1), g_csc%jptr(col1+1)-1
107  ! for each nnz in the column, we go through the corresponding
108  ! row and count all the i,j pairs
109  row1 = g_csc%irn(i)
110 
111  do j=g_csr%iptr(row1), g_csr%iptr(row1+1)-1
112  col2 = g_csr%jcn(j)
113  ! now, the element (col1,col2) is present in A'*A. Check if
114  ! it was already found, otherwise add it.
115  ! skip the diagonal
116  if(col1 .eq. col2) cycle
117  if(mark(col2) .lt. col1) then
118  mark(col2) = col1
119  ata_graph%jcn(ata_graph%iptr(col1+1)) = col2
120  ata_graph%iptr(col1+1) = ata_graph%iptr(col1+1)+1
121  end if
122  end do
123  end do
124  end do
125 
126  ata_graph%n = g_csc%n
127  ata_graph%m = g_csc%n
128 
129  call _qrm_spmat_destroy(g_csr, all=.true.)
130  call qrm_adealloc(mark)
131  __qrm_check_ret(name,'destroy/dealloc',9999)
132 
133  call qrm_err_act_restore(err_act)
134  return
135 
136 9999 continue ! error management
137  call qrm_err_act_restore(err_act)
138  if(err_act .eq. qrm_abort_) then
139  call qrm_err_check()
140  end if
141 
142  return
143 
144 end subroutine _qrm_ata_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.
subroutine _qrm_spmat_destroy(qrm_spmat, all)
This subroutine destroyes a qrm_spmat instance.
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.
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.
This module contains the definition of the basic sparse matrix type and of the associated methods...
Generic interface for the qrm_palloc_i, qrm_palloc_2i, qrm_palloc_s, qrm_palloc_2s, qrm_palloc_d, qrm_palloc_2d, qrm_palloc_c, qrm_palloc_2c, qrm_palloc_z, qrm_palloc_2z, routines.
Definition: qrm_mem_mod.F90:57
subroutine _qrm_spmat_convert(in_mat, out_mat, fmt, values)
This subroutine converts an input matrix into a different storage format. Optionally the values may b...
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38
subroutine _qrm_ata_graph(g_csc, ata_graph)
This subroutine computes the fill reducing ordering using METIS.
subroutine qrm_err_act_restore(err_act)
Restores the value of the qrm_err_act variable.