QR_MUMPS
dqrm_do_scotch.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 
49 subroutine dqrm_do_scotch(graph, cperm)
50 
51  use dqrm_spmat_mod
52  use qrm_error_mod
53  use dqrm_analysis_mod, savesym => dqrm_do_scotch
54  use qrm_mem_mod
55  use iso_c_binding
56  implicit none
57 
58 #if defined(have_scotch)
59  include "scotchf.h"
60 #endif
61 
62  type(dqrm_spmat_type) :: graph
63  integer :: cperm(:)
64 
65 #if defined(have_scotch)
66  integer :: i, info, cblknbr
67  type(dqrm_spmat_type) :: ata_graph
68  integer, allocatable :: iperm(:)
69  real(kind(1.d0)) :: grafdat(scotch_graphdim), stradat(scotch_stratdim), orderdat(scotch_orderdim)
70  ! error management
71  integer :: err_act
72  character(len=*), parameter :: name='qrm_do_scotch'
73 
74  call qrm_err_act_save(err_act)
75 
76  call dqrm_ata_graph(graph, ata_graph)
77  __qrm_check_ret(name,'qrm_ata_graph',9999)
78 
79  info = 0
80  call scotchfgraphinit(grafdat, info)
81  call scotchfstratinit(stradat, info)
82  if(info .ne. 0) then
83  call qrm_err_push(19,name)
84  goto 9999
85  end if
86 
87  call scotchfgraphbuild(grafdat, 1, ata_graph%n, ata_graph%iptr(1), &
88  & ata_graph%iptr(2), ata_graph%iptr, ata_graph%iptr, ata_graph%nz, &
89  & ata_graph%jcn, ata_graph%jcn, info)
90  if(info .ne. 0) then
91  call qrm_err_push(19,name)
92  goto 9999
93  end if
94 
95  call scotchfgraphorder(grafdat, stradat, grafdat, cperm, cblknbr, &
96  & grafdat, grafdat, info)
97  if(info .ne. 0) then
98  call qrm_err_push(19,name)
99  goto 9999
100  end if
101 
102  call scotchfgraphexit(grafdat)
103  call scotchfstratexit(stradat)
104 
105  call dqrm_spmat_destroy(ata_graph, all=.true.)
106  __qrm_check_ret(name,'qrm_spmat_destroy',9999)
107 
108  call qrm_err_act_restore(err_act)
109  return
110 
111 9999 continue ! error management
112  call qrm_err_act_restore(err_act)
113  if(err_act .eq. qrm_abort_) then
114  call qrm_err_check()
115  end if
116 #endif
117 
118  return
119 end subroutine dqrm_do_scotch
subroutine dqrm_do_scotch(graph, cperm)
Please refer to:
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.
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...
This type defines the data structure used to store a matrix.
subroutine dqrm_ata_graph(g_csc, ata_graph)
This subroutine computes the fill reducing ordering using METIS.
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.