QR_MUMPS
qrm_do_ordering.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 
38 
41 
51 subroutine _qrm_do_ordering(graph, cperm, cperm_in)
52 
53  use _qrm_spmat_mod
54  use qrm_error_mod
55  use _qrm_analysis_mod, savesym => _qrm_do_ordering
56  use qrm_mem_mod
57  use qrm_const_mod
58  use qrm_common_mod
59  implicit none
60 
61  type(_qrm_spmat_type) :: graph
62  integer :: cperm(:)
63  integer, pointer :: cperm_in(:)
64 
65  integer :: iord=0
66  integer :: i
67  ! error management
68  integer :: err_act
69  character(len=*), parameter :: name='qrm_do_ordering'
70 
71  call qrm_err_act_save(err_act)
72 
73  call qrm_get(graph, 'qrm_ordering', iord)
74 
75  if(iord .eq. qrm_auto_) then
76  iord = qrm_choose_ordering()
77  end if
78 
79 
80  select case(iord)
81  case(qrm_natural_)
82  ! perm has to be allocated and built explicitly because
83  ! we are going to use an equivalent permutation
84  do i=1, graph%n
85  cperm(i) = i
86  end do
87  case(qrm_given_)
88  ! in this case we just need to check that
89  ! the given permutation makes sense
90  if(.not. associated(cperm_in)) then
91  call qrm_err_push(8, name)
92  goto 9999
93  else
94  call qrm_check_cperm(cperm_in, graph%n)
95  __qrm_check_ret(name,'qrm_check_perm',9999)
96  cperm(1:graph%n) = cperm_in(1:graph%n)
97  end if
98  case(qrm_colamd_)
99  ! COLAMD
100 #if defined (have_colamd)
101  call _qrm_do_colamd(graph, cperm)
102 #else
103  call qrm_err_push(16, 'qrm_do_ordering',aed='colamd')
104  goto 9999
105 #endif
106  case(qrm_metis_)
107  ! METIS
108 #if defined (have_metis)
109  call _qrm_do_metis(graph, cperm)
110 #else
111  call qrm_err_push(16, 'qrm_do_ordering',aed='metis')
112  goto 9999
113 #endif
114  case(qrm_scotch_)
115  ! SCOTCH
116 #if defined (have_scotch)
117  call _qrm_do_scotch(graph, cperm)
118 #else
119  call qrm_err_push(16, 'qrm_do_ordering',aed='scotch')
120  goto 9999
121 #endif
122  case default
123  call qrm_err_push(9, 'qrm_do_ordering',ied=(/iord, 0, 0, 0, 0/))
124  goto 9999
125  end select
126 
127  call qrm_err_act_restore(err_act)
128  return
129 
130 9999 continue ! error management
131  call qrm_err_act_restore(err_act)
132  if(err_act .eq. qrm_abort_) then
133  call qrm_err_check()
134  end if
135 
136  return
137 
138 contains
139 
140  function qrm_choose_ordering()
141  ! Function: qrm_choose_ordering
142  ! This function sets the ordering for the case where an automatic
143  ! choice is requested
144  !
145 
146  integer :: qrm_choose_ordering
147  integer :: iord
148 
149  iord = 1
150 
151 #if defined(ccolamd)
152  iord = qrm_colamd_
153 #endif
154 
155 #if defined(metis)
156  iord = qrm_metis_
157 #endif
158 
159 #if defined(scotch)
160  iord = qrm_scotch_
161 #endif
162 
163  qrm_choose_ordering = iord
164  return
165 
166  end function qrm_choose_ordering
167 
168 
169 end subroutine _qrm_do_ordering
integer, parameter qrm_given_
This module contains the interfaces of all non-typed routines.
subroutine qrm_err_push(code, sub, ied, aed)
This subroutine pushes an error on top of the stack.
subroutine _qrm_do_colamd(graph, cperm)
This subroutine computes the fill reducing ordering using COLAMD.
subroutine _qrm_do_ordering(graph, cperm, cperm_in)
This routine computes (through different methods) a column permutation of the input matrix in order t...
integer function qrm_choose_ordering()
integer, parameter qrm_scotch_
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.
integer, parameter qrm_natural_
integer, parameter qrm_abort_
Possible actions to be performed upon detection of an error.
subroutine _qrm_do_metis(graph, cperm)
Please refer to:
This module contains the generic interfaces for all the analysis routines.
subroutine qrm_err_check()
This subroutine checks the errors stack. If something is found all the entries in the stack are poppe...
subroutine _qrm_do_scotch(graph, cperm)
Please refer to:
This type defines the data structure used to store a matrix.
integer, parameter qrm_metis_
This module contains the definition of the basic sparse matrix type and of the associated methods...
Generic interface for the ::qrm_check_cperm routine.
integer, parameter qrm_auto_
constants to set control parameters
Generif interface for the ::_qrm_pgeti, ::_qrm_pgetr and.
integer, parameter qrm_colamd_
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.