QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
qrm_residual_norm.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 #include "qrm_common.h"
35 
37 
50 subroutine _qrm_residual_norm2d(qrm_mat, b, x, nrm)
51  use _qrm_spmat_mod
52  use qrm_error_mod
53  use _qrm_utils_mod
54  implicit none
55 
56  type(_qrm_spmat_type) :: qrm_mat
57  _qrm_data :: b(:,:), x(:,:)
58  _qrm_real :: nrm(:)
59 
60  _qrm_real, allocatable :: nrmb(:), nrmx(:)
61  _qrm_real :: nrma
62  integer :: nrhs
63 
64  ! error management
65  integer :: err_act
66  character(len=*), parameter :: name='qrm_residual_norm'
67 
68  call qrm_err_act_save(err_act)
69 
70  call _qrm_check_spmat(qrm_mat)
71  __qrm_check_ret(name,'qrm_check_spmat',9999)
72 
73  nrhs = min(size(b,2),size(x,2))
74 
75  call qrm_aalloc(nrmb, nrhs)
76  call qrm_aalloc(nrmx, nrhs)
77 
78  call qrm_vecnrm(b(:,1:nrhs), qrm_mat%m, 'i', nrmb)
79  call qrm_vecnrm(x(:,1:nrhs), qrm_mat%n, 'i', nrmx)
80 
81  ! compute the residual
82  call qrm_matmul(qrm_mat, 'n', -_qrm_one, x, _qrm_one, b)
83  call qrm_matnrm(qrm_mat, 'i', nrma)
84  call qrm_vecnrm(b, qrm_mat%m, 'i', nrm)
85 
86  nrmb = nrmb + nrma*nrmx
87 
88  nrm = nrm/nrmb
89 
90  call qrm_adealloc(nrmx)
91  call qrm_adealloc(nrmb)
92 
93  call qrm_err_act_restore(err_act)
94  return
95 
96 9999 continue ! error management
97  call qrm_err_act_restore(err_act)
98  if(err_act .eq. qrm_abort_) then
99  call qrm_err_check()
100  end if
101  return
102 
103 end subroutine _qrm_residual_norm2d
104 
105 
106 
108 
121 subroutine _qrm_residual_norm1d(qrm_mat, b, x, nrm)
122  use _qrm_spmat_mod
123  use qrm_error_mod
124  use _qrm_utils_mod
125  implicit none
126 
127  type(_qrm_spmat_type) :: qrm_mat
128  _qrm_data :: b(:), x(:)
129  _qrm_real :: nrm
130 
131  _qrm_real :: nrmb, nrmx
132  _qrm_real :: nrma
133 
134 
135  ! error management
136  integer :: err_act
137  character(len=*), parameter :: name='qrm_residual_norm'
138 
139  call qrm_err_act_save(err_act)
140 
141  call _qrm_check_spmat(qrm_mat)
142  __qrm_check_ret(name,'qrm_check_spmat',9999)
143 
144  call qrm_vecnrm(b, qrm_mat%m, 'i', nrmb)
145  call qrm_vecnrm(x, qrm_mat%n, 'i', nrmx)
146 
147  ! compute the residual
148  call qrm_matmul(qrm_mat, 'n', -_qrm_one, x, _qrm_one, b)
149  call qrm_matnrm(qrm_mat, 'i', nrma)
150  call qrm_vecnrm(b, qrm_mat%m, 'i', nrm)
151 
152  nrmb = nrmb + nrma*nrmx
153  nrm = nrm/nrmb
154 
155  call qrm_err_act_restore(err_act)
156  return
157 
158 9999 continue ! error management
159  call qrm_err_act_restore(err_act)
160  if(err_act .eq. qrm_abort_) then
161  call qrm_err_check()
162  end if
163  return
164 
165 end subroutine _qrm_residual_norm1d
This module contains all the error management routines and data.
subroutine _qrm_residual_norm1d(qrm_mat, b, x, nrm)
This routine computes the scaled norm of the residual.
subroutine qrm_err_act_restore(err_act)
Restores the value of the qrm_err_act variable.
subroutine _qrm_residual_norm2d(qrm_mat, b, x, nrm)
This routine computes the scaled norm of multiple residuals.
This module contains generic interfaces for a number of auxiliary tools.
subroutine qrm_err_check()
This subroutine checks the errors stack. If something is found all the entries in the stack are poppe...
Generic interface for the ::_qrm_vecnrm2d and ::_qrm_vecnrm1d routines.
This type defines the data structure used to store a matrix.
subroutine _qrm_check_spmat(qrm_spmat, op)
Check the compatibility and correctness of icntl and rcntl parameters.
This module contains the definition of the basic sparse matrix type and of the associated methods...
Generic interface for the ::_qrm_matmul2d and ::_qrm_matmul1d routines.
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.