QR_MUMPS
dqrm_do_colamd.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 
37 
39 
52 subroutine dqrm_do_colamd(graph, cperm)
53 
54  use dqrm_spmat_mod
55  use qrm_error_mod
56  use dqrm_analysis_mod, savesym => dqrm_do_ordering, savesym2 => dqrm_do_colamd
57  use qrm_mem_mod
58  use iso_c_binding
59 
60  implicit none
61 
62  type(dqrm_spmat_type) :: graph
63  integer, target :: cperm(:)
64 
65  interface
66  subroutine qrm_colamd(n_row, n_col, Alen, A, p, err) bind(c, name='qrm_colamd')
67  use iso_c_binding
68  integer(c_int), value :: n_row, n_col, Alen
69  integer(c_int) :: A(*), p(*), err
70  end subroutine qrm_colamd
71  end interface
72 
73  interface
74  subroutine qrm_colamd_recommended(Alen, nnz, n_row, n_col) bind(c, name='qrm_colamd_recommended')
75  use iso_c_binding
76  integer(c_int), value :: nnz, n_row, n_col
77  integer(c_int) :: Alen
78  end subroutine qrm_colamd_recommended
79  end interface
80 
81  type(dqrm_spmat_type) :: gcopy
82  integer :: i, idx, cnt, tmp, alen, err
83  ! error management
84  integer :: err_act
85  character(len=*), parameter :: name='qrm_do_colamd'
86 
87  call qrm_err_act_save(err_act)
88 
89 
90  ! at this point we have to make a copy of the graph.
91  ! this is a huge waste of mem but we don't have a choice
92  ! since ccolamd destroys the graph which, instead, we want to
93  ! save for successive computations
94 
95  ! compute the memory required by ccolamd (a lot) and allocate
96  call qrm_colamd_recommended(alen, graph%nz, graph%m, graph%n)
97  call qrm_palloc(gcopy%irn, alen)
98  __qrm_check_ret(name,'qrm_palloc',9999)
99  gcopy%jptr => cperm
100 
101  ! make the copy
102  call dqrm_spmat_copy(graph, gcopy, values=.false.)
103  __qrm_check_ret(name,'qrm_spmat_copy',9999)
104 
105  ! ccolamd wants 0 based indices (argh!)
106  gcopy%irn(1:gcopy%nz) = gcopy%irn(1:gcopy%nz)-1
107  gcopy%jptr(1:gcopy%n+1) = gcopy%jptr(1:gcopy%n+1)-1
108 
109  ! call ccolamd
110  call qrm_colamd(gcopy%m, gcopy%n, alen, gcopy%irn, gcopy%jptr, err)
111  if(err .eq. 0) then
112  call qrm_err_push(18,name)
113  goto 9999
114  end if
115 
116  nullify(gcopy%jptr)
117  ! ccolamd return a 0-based permutation (re-argh!)
118  cperm = cperm+1
119 
120  call dqrm_spmat_destroy(gcopy, all=.true.)
121  __qrm_check_ret(name,'qrm_spmat_destroy',9999)
122 
123  call qrm_err_act_restore(err_act)
124  return
125 
126 9999 continue ! error management
127  call qrm_err_act_restore(err_act)
128  if(err_act .eq. qrm_abort_) then
129  call qrm_err_check()
130  end if
131 
132  return
133 
134 end subroutine dqrm_do_colamd
void qrm_colamd(int n_row, int n_col, int Alen, int *A, int *p, int *err)
subroutine qrm_err_push(code, sub, ied, aed)
This subroutine pushes an error on top of the stack.
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.
subroutine dqrm_spmat_copy(in_mat, out_mat, values)
This subroutine makes a copy of a matrix. Optionally the values may be ignored (this comes handy duri...
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...
integer, parameter qrm_abort_
Possible actions to be performed upon detection of an error.
subroutine dqrm_spmat_destroy(qrm_spmat, all)
This subroutine destroyes a qrm_spmat instance.
subroutine qrm_err_check()
This subroutine checks the errors stack. If something is found all the entries in the stack are poppe...
subroutine dqrm_do_colamd(graph, cperm)
This subroutine computes the fill reducing ordering using COLAMD.
This type defines the data structure used to store a matrix.
subroutine dqrm_do_ordering(graph, cperm, cperm_in)
This routine computes (through different methods) a column permutation of the input matrix in order t...
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
void qrm_colamd_recommended(int *alen, int nnz, int n_row, int n_col)
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.