QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
qrm_vecnrm.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 
50 subroutine _qrm_vecnrm2d(vec, n, ntype, nrm)
51 
52  use qrm_string_mod
53  use qrm_error_mod
54  implicit none
55 
56  _qrm_data, intent(in) :: vec(:,:)
57  _qrm_real :: nrm(:)
58  integer, intent(in) :: n
59  character :: ntype
60 
61  integer :: i, j
62  _qrm_real :: _rxnrm2
63 
64  ! error management
65  integer :: err_act
66  character(len=*), parameter :: name='qrm_vecnrm'
67 
68  call qrm_err_act_save(err_act)
69 
70  nrm = _qrm_rzero
71 
72  if(qrm_str_tolower(ntype) .eq. 'i') then
73  do j=1, size(vec,2)
74  nrm(j) = maxval(abs(vec(:,j)))
75  end do
76  else if(qrm_str_tolower(ntype) .eq. '1') then
77  do j=1, size(vec,2)
78  nrm(j) = _qrm_zero
79  do i=1, n
80  nrm(j) = nrm(j) + abs(vec(i,j))
81  end do
82  end do
83  else if(qrm_str_tolower(ntype) .eq. '2') then
84  do j=1, size(vec,2)
85  nrm(j) = _rxnrm2(n, vec(1,j), 1)
86  end do
87  else
88  call qrm_err_push(15, '_qrm_vecnorm')
89  end if
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  return
100 
101 end subroutine _qrm_vecnrm2d
102 
103 
104 
118 subroutine _qrm_vecnrm1d(vec, n, ntype, nrm)
119 
120  use qrm_string_mod
121  use qrm_error_mod
122  implicit none
123 
124  _qrm_data, intent(in) :: vec(:)
125  _qrm_real :: nrm
126  integer, intent(in) :: n
127  character :: ntype
128 
129  integer :: i
130  _qrm_real :: _rxnrm2
131 
132  ! error management
133  integer :: err_act
134  character(len=*), parameter :: name='qrm_vecnrm'
135 
136  call qrm_err_act_save(err_act)
137 
138  nrm = _qrm_rzero
139 
140  if(qrm_str_tolower(ntype) .eq. 'i') then
141  nrm = maxval(abs(vec))
142  else if(qrm_str_tolower(ntype) .eq. '1') then
143  nrm = _qrm_zero
144  do i=1, n
145  nrm = nrm + abs(vec(i))
146  end do
147  else if(qrm_str_tolower(ntype) .eq. '2') then
148  nrm = _rxnrm2(n, vec, 1)
149  else
150  call qrm_err_push(15, '_qrm_vecnorm')
151  end if
152 
153  call qrm_err_act_restore(err_act)
154  return
155 
156 9999 continue ! error management
157  call qrm_err_act_restore(err_act)
158  if(err_act .eq. qrm_abort_) then
159  call qrm_err_check()
160  end if
161  return
162 
163 end subroutine _qrm_vecnrm1d
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_vecnrm1d(vec, n, ntype, nrm)
This subroutine computes the norm of a vector.
Definition: qrm_vecnrm.F90:118
subroutine qrm_err_check()
This subroutine checks the errors stack. If something is found all the entries in the stack are poppe...
subroutine _qrm_vecnrm2d(vec, n, ntype, nrm)
This subroutine computes the norm of multiple vectors.
Definition: qrm_vecnrm.F90:50
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.