1!-------------------------------------------------------------------------------
2! Copyright (c) 2020 FrontISTR Commons
3! This software is released under the MIT License, see LICENSE.txt
4!-------------------------------------------------------------------------------
5
6!----------------------------------------------------------------------
7!> @brief HECMW_ORDERING_METIS is a program for the Nested Dissection
8!         ordering provided by Metis
9!----------------------------------------------------------------------
10module hecmw_ordering_metis
11  use hecmw_util
12  implicit none
13
14  private
15  public :: hecmw_ordering_metis_nodend
16
17contains
18
19  !======================================================================!
20  !> @brief hecmw_ordering_metis_NodeND
21  !======================================================================!
22  subroutine hecmw_ordering_metis_NodeND(Neqns,Xadj,Adj0,Perm,Invp)
23    implicit none
24    !------
25    integer(kind=kint), intent(in):: NEQns
26    integer(kind=kint), intent(in):: Adj0(:)
27    integer(kind=kint), intent(in):: Xadj(:)
28    integer(kind=kint), intent(out):: Perm(:)
29    integer(kind=kint), intent(out):: Invp(:)
30    !------
31#ifdef HECMW_WITH_METIS
32
33#  if HECMW_METIS_VER == 5
34
35    integer(kind=kint), allocatable:: vwght(:)
36    integer(kind=kint):: options(40)
37    integer(kind=kint):: ierror
38
39    allocate(vwght(Neqns),stat=ierror)
40    if ( ierror/=0 ) stop "ALLOCATION ERROR, vwght: SUB. gennd"
41
42    vwght(:)=1
43
44    call METIS_SetDefaultOptions(options)
45    ! set fortran numbering
46    options(18)=1
47
48    call METIS_NodeND(Neqns,Xadj,Adj0,vwght,options,Perm,Invp)
49
50    deallocate(vwght)
51
52#  elif HECMW_METIS_VER == 4
53
54    integer(kind=kint):: numflag
55    integer(kind=kint):: options(8)
56
57    numflag=1
58    options(:)=0
59
60    call METIS_NodeND(Neqns,Xadj,Adj0,numflag,options,Perm,Invp)
61
62#  else
63#    error unknown HECMW_METIS_VER
64#  endif
65
66#else
67    stop "METIS not available"
68#endif
69  end subroutine hecmw_ordering_metis_NodeND
70
71end module hecmw_ordering_metis
72