QR_MUMPS
dqrm_do_metis.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_metis(graph, cperm)
50 
51  use dqrm_spmat_mod
52  use qrm_error_mod
53  use dqrm_analysis_mod, savesym => dqrm_do_metis
54  use qrm_mem_mod
55  use iso_c_binding
56 
57  implicit none
58 
59  type(dqrm_spmat_type) :: graph
60  integer :: cperm(:)
61 
62  interface
63  subroutine qrm_metis(n, iptr, jcn, cperm, iperm) bind(c, name='qrm_metis')
64  use iso_c_binding
65  integer(c_int) :: n
66  integer(c_int) :: iptr(*), jcn(*), cperm(*), iperm(*)
67  end subroutine qrm_metis
68  end interface
69 
70  integer :: i, idx, cnt, tmp, alen
71  type(dqrm_spmat_type) :: ata_graph
72  integer, allocatable :: iperm(:)
73  ! error management
74  integer :: err_act
75  character(len=*), parameter :: name='qrm_do_metis'
76 
77  call qrm_err_act_save(err_act)
78 
79  call dqrm_ata_graph(graph, ata_graph)
80  __qrm_check_ret(name,'qrm_ata_graph',9999)
81 
82  call qrm_aalloc(iperm, graph%n)
83  __qrm_check_ret(name,'qrm_ata_graph/alloc',9999)
84 
85  call qrm_metis(graph%n, ata_graph%iptr, ata_graph%jcn, cperm, iperm)
86 
87  call qrm_adealloc(iperm)
88  call dqrm_spmat_destroy(ata_graph, all=.true.)
89  __qrm_check_ret(name,'dealloc/destroy',9999)
90 
91  call qrm_err_act_restore(err_act)
92  return
93 
94 9999 continue ! error management
95  call qrm_err_act_restore(err_act)
96  if(err_act .eq. qrm_abort_) then
97  call qrm_err_check()
98  end if
99 
100  return
101 end subroutine dqrm_do_metis
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 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...
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 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_metis(graph, cperm)
Please refer to:
This type defines the data structure used to store a matrix.
void qrm_metis(int *n, int *iptr, int *jcn, int *cperm, int *iperm)
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.