QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
qrm_matnrm.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 
47 subroutine _qrm_matnrm(qrm_mat, ntype, nrm)
48 
49  use _qrm_spmat_mod
50  use qrm_string_mod
51  use qrm_error_mod
52  implicit none
53 
54  type(_qrm_spmat_type), intent(in) :: qrm_mat
55  _qrm_real :: nrm
56  character :: ntype
57 
58  _qrm_real, allocatable :: tmp(:)
59  integer :: r, c, i
60  _qrm_real :: _rxnrm2
61 
62  ! error management
63  integer :: err_act
64  character(len=*), parameter :: name='qrm_matnrm'
65 
66  call qrm_err_act_save(err_act)
67 
68  if(qrm_str_tolower(ntype) .eq. 'i') then
69  call qrm_aalloc(tmp, qrm_mat%m)
70  __qrm_check_ret(name,'qrm_aalloc',9999)
71  tmp = _qrm_zero
72  do i=1, qrm_mat%nz
73  r = qrm_mat%irn(i)
74  tmp(r) = tmp(r)+abs(qrm_mat%val(i))
75  end do
76  nrm = maxval(tmp)
77  else if(qrm_str_tolower(ntype) .eq. '1') then
78  call qrm_aalloc(tmp, qrm_mat%n)
79  __qrm_check_ret(name,'qrm_aalloc',9999)
80  tmp = _qrm_zero
81  do i=1, qrm_mat%nz
82  c = qrm_mat%jcn(i)
83  tmp(c) = tmp(c)+abs(qrm_mat%val(i))
84  end do
85  nrm = maxval(tmp)
86  else if(qrm_str_tolower(ntype) .eq. 'f') then
87  nrm = _rxnrm2(qrm_mat%nz, qrm_mat%val, 1)
88  else
89  call qrm_err_push(15, 'qrm_matnrm')
90  goto 9999
91  end if
92 
93  call qrm_adealloc(tmp)
94  __qrm_check_ret(name,'qrm_adealloc',9999)
95 
96  call qrm_err_act_restore(err_act)
97  return
98 
99 9999 continue ! error management
100  call qrm_err_act_restore(err_act)
101  if(err_act .eq. qrm_abort_) then
102  call qrm_err_check()
103  end if
104  return
105 
106 end subroutine _qrm_matnrm
This module contains all the error management routines and data.
subroutine qrm_err_act_restore(err_act)
Restores the value of the qrm_err_act variable.
subroutine qrm_err_check()
This subroutine checks the errors stack. If something is found all the entries in the stack are poppe...
subroutine _qrm_matnrm(qrm_mat, ntype, nrm)
This subroutine computes the matrix norm. The return value is a real scalar.
Definition: qrm_matnrm.F90:47
This type defines the data structure used to store a matrix.
This module contains the definition of the basic sparse matrix type and of the associated methods...
This module contains various string handling routines.
int i
Definition: secs.c:40
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.