QR_MUMPS
qrm_check_cperm.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 
44 subroutine qrm_check_cperm(cperm, n)
46  use qrm_mem_mod
47  implicit none
48 
49  integer :: cperm(:)
50  integer :: n
51 
52  integer, allocatable :: tmp(:)
53  integer :: i
54  ! error management
55  integer :: err_act
56  character(len=*), parameter :: name='qrm_check_cperm'
57 
58  call qrm_err_act_save(err_act)
59 
60  call qrm_aalloc(tmp, n)
61  __qrm_check_ret(name,'qrm_aalloc',9999)
62  tmp = 0
63 
64  do i=1, n
65  if (cperm(i) .gt. n .or. cperm(i) .lt. 1) then
66  call qrm_err_push(8, name)
67  call qrm_adealloc(tmp)
68  goto 9999
69  end if
70  if(tmp(cperm(i)) .gt. 0) then
71  call qrm_err_push(8, name)
72  call qrm_adealloc(tmp)
73  goto 9999
74  else
75  tmp(cperm(i)) = 1
76  end if
77  end do
78 
79  call qrm_adealloc(tmp)
80 
81 9999 continue ! error management
82  call qrm_err_act_restore(err_act)
83  if(err_act .eq. qrm_abort_) then
84  call qrm_err_check()
85  end if
86 
87  return
88 
89 end subroutine qrm_check_cperm
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.
subroutine qrm_err_push(code, sub, ied, aed)
This subroutine pushes an error on top of the stack.
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.
subroutine qrm_check_cperm(cperm, n)
This routine simply checks whether a column permutation provided by the user makes sens...
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 qrm_err_check()
This subroutine checks the errors stack. If something is found all the entries in the stack are poppe...
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.