QR_MUMPS
qrm_rowperm.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 !! ##############################################################################################
34 
35 
36 #include "qrm_common.h"
37 
41 
68 subroutine _qrm_rowperm(graph, cperm, rperm, nvar, stair)
69 
70  use _qrm_spmat_mod
71  use qrm_mem_mod
72  use qrm_common_mod
73  implicit none
74 
75  type(_qrm_spmat_type) :: graph
76  integer :: cperm(:), rperm(:), stair(:), nvar(:)
77 
78 
79  integer, allocatable :: mark(:)
80  integer :: ii, jj, i, j, node, nv, pnt, rpnt
81  ! error management
82  integer :: err_act
83  character(len=*), parameter :: name='qrm_rowperm'
84 
85  call qrm_err_act_save(err_act)
86 
87  call qrm_aalloc(mark, graph%m)
88  __qrm_check_ret(name,'qrm_aalloc',9999)
89  mark = 0
90  stair = 0
91  ! scan the matrix by columns
92 
93  rpnt = 0
94  pnt=1
95  do
96  if (pnt .gt. graph%n) exit
97  node = cperm(pnt)
98  stair(node) = rpnt
99  nv = nvar(node)
100  do jj = pnt, pnt+nv-1
101  j = cperm(jj)
102  do ii=graph%jptr(j), graph%jptr(j+1)-1
103  i = graph%irn(ii)
104  if(mark(i) .eq. 0) then
105  stair(node) = stair(node)+1
106  rpnt = rpnt+1
107  rperm(rpnt) = i
108  mark(i) = j
109  end if
110  end do
111  end do
112  pnt = pnt+nv
113  if (pnt .gt. graph%n) exit
114  end do
115 
116  ! look for empty rows and put them at the end
117  do i=1, graph%m
118  if(mark(i) .eq. 0) then
119  rpnt = rpnt+1
120  rperm(rpnt) = i
121  mark(i) = i
122  end if
123  end do
124 
125  call qrm_adealloc(mark)
126  __qrm_check_ret(name,'qrm_adealloc',9999)
127 
128  call qrm_err_act_restore(err_act)
129  return
130 
131 9999 continue ! error management
132  call qrm_err_act_restore(err_act)
133  if(err_act .eq. qrm_abort_) then
134  call qrm_err_check()
135  end if
136 
137  return
138 
139 end subroutine _qrm_rowperm
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 interfaces of all non-typed routines.
subroutine _qrm_rowperm(graph, cperm, rperm, nvar, stair)
This routine computes a row permutation that puts the input matrix into a staircase format...
Definition: qrm_rowperm.F90:69
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 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...
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38