QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
qrm_residual_orth.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 !! ##############################################################################################
34 
35 #include "qrm_common.h"
36 
38 
48 subroutine _qrm_residual_orth2d(qrm_mat, r, nrm)
49  use _qrm_spmat_mod
50  use qrm_error_mod
51  use qrm_mem_mod
52  use _qrm_utils_mod
53  implicit none
54 
55  type(_qrm_spmat_type) :: qrm_mat
56  _qrm_data :: r(:,:)
57  _qrm_real :: nrm(:)
58 
59  _qrm_data, allocatable :: tmp(:,:)
60  _qrm_real, allocatable :: rnrm(:)
61  _qrm_real :: anrm
62 
63  ! error management
64  integer :: err_act
65  character(len=*), parameter :: name='qrm_residual_orth'
66 
67  call qrm_err_act_save(err_act)
68 
69  call _qrm_check_spmat(qrm_mat)
70  __qrm_check_ret(name,'qrm_check_spmat',9999)
71 
72  call qrm_aalloc(tmp, qrm_mat%n, size(r,2))
73  call qrm_aalloc(rnrm, size(r,2))
74 
75  ! compute A'*r
76  call qrm_matmul(qrm_mat, 't', _qrm_one, r, _qrm_zero, tmp)
77 
78  call qrm_vecnrm(r , qrm_mat%m, '2', rnrm)
79  call qrm_vecnrm(tmp, qrm_mat%n, '2', nrm)
80  call qrm_matnrm(qrm_mat, 'f', anrm)
81  nrm = nrm/(rnrm*anrm)
82 
83  call qrm_adealloc(tmp)
84  call qrm_adealloc(rnrm)
85 
86  call qrm_err_act_restore(err_act)
87  return
88 
89 9999 continue ! error management
90  call qrm_err_act_restore(err_act)
91  if(err_act .eq. qrm_abort_) then
92  call qrm_err_check()
93  end if
94  return
95 
96 end subroutine _qrm_residual_orth2d
97 
99 
109 subroutine _qrm_residual_orth1d(qrm_mat, r, nrm)
110  use _qrm_spmat_mod
111  use qrm_error_mod
112  use qrm_mem_mod
113  use _qrm_utils_mod
114  implicit none
115 
116  type(_qrm_spmat_type) :: qrm_mat
117  _qrm_data :: r(:)
118  _qrm_real :: nrm
119 
120  _qrm_data, allocatable :: tmp(:)
121  _qrm_real :: rnrm
122  _qrm_real :: anrm
123 
124  ! error management
125  integer :: err_act
126  character(len=*), parameter :: name='qrm_residual_orth'
127 
128  call qrm_err_act_save(err_act)
129 
130  call _qrm_check_spmat(qrm_mat)
131  __qrm_check_ret(name,'qrm_check_spmat',9999)
132 
133  call qrm_aalloc(tmp, qrm_mat%n)
134 
135  ! compute A'*r
136  call qrm_matmul(qrm_mat, 't', _qrm_one, r, _qrm_zero, tmp)
137 
138  call qrm_vecnrm(r, qrm_mat%m, '2', rnrm)
139  call qrm_vecnrm(tmp, qrm_mat%n, '2', nrm)
140  call qrm_matnrm(qrm_mat, 'f', anrm)
141 
142  nrm = nrm/(rnrm*anrm)
143  call qrm_adealloc(tmp)
144 
145  call qrm_err_act_restore(err_act)
146  return
147 
148 9999 continue ! error management
149  call qrm_err_act_restore(err_act)
150  if(err_act .eq. qrm_abort_) then
151  call qrm_err_check()
152  end if
153  return
154 
155 end subroutine _qrm_residual_orth1d
156 
157 
158 
159 
160 
162 
173 subroutine _qrm_residual_and_orth2d(qrm_mat, b, x, nrm)
174  use _qrm_spmat_mod
175  use qrm_error_mod
176  use qrm_mem_mod
177  use _qrm_utils_mod
178  implicit none
179 
180  type(_qrm_spmat_type) :: qrm_mat
181  _qrm_data :: b(:,:), x(:,:)
182  _qrm_real :: nrm(:)
183 
184  _qrm_data, allocatable :: tmp(:,:)
185  _qrm_real, allocatable :: rnrm(:)
186 
187  ! error management
188  integer :: err_act
189  character(len=*), parameter :: name='qrm_residual_orth'
190 
191  call qrm_err_act_save(err_act)
192 
193  call _qrm_check_spmat(qrm_mat)
194  __qrm_check_ret(name,'qrm_check_spmat',9999)
195 
196  ! compute the residual
197  call qrm_matmul(qrm_mat, 'n', -_qrm_one, x, _qrm_one, b)
198 
199  call qrm_aalloc(tmp , qrm_mat%n, size(x,2))
200  call qrm_aalloc(rnrm, size(x,2))
201 
202  ! compute A'*r
203  call qrm_matmul(qrm_mat, 't', _qrm_one, b, _qrm_zero, tmp)
204 
205  call qrm_vecnrm(b, qrm_mat%m, '2', rnrm)
206  call qrm_vecnrm(tmp, qrm_mat%n, '2', nrm)
207  nrm = nrm/rnrm
208 
209  call qrm_adealloc(tmp)
210 
211  call qrm_err_act_restore(err_act)
212  return
213 
214 9999 continue ! error management
215  call qrm_err_act_restore(err_act)
216  if(err_act .eq. qrm_abort_) then
217  call qrm_err_check()
218  end if
219  return
220 
221 end subroutine _qrm_residual_and_orth2d
222 
223 
225 
236 subroutine _qrm_residual_and_orth1d(qrm_mat, b, x, nrm)
237  use _qrm_spmat_mod
238  use qrm_error_mod
239  use qrm_mem_mod
240  use _qrm_utils_mod
241  implicit none
242 
243  type(_qrm_spmat_type) :: qrm_mat
244  _qrm_data :: b(:), x(:)
245  _qrm_real :: nrm
246 
247  _qrm_data, allocatable :: tmp(:)
248  _qrm_real :: rnrm
249 
250  ! error management
251  integer :: err_act
252  character(len=*), parameter :: name='qrm_residual_orth'
253 
254  call qrm_err_act_save(err_act)
255 
256  call _qrm_check_spmat(qrm_mat)
257  __qrm_check_ret(name,'qrm_check_spmat',9999)
258 
259  ! compute the residual
260  call qrm_matmul(qrm_mat, 'n', -_qrm_one, x, _qrm_one, b)
261 
262  call qrm_aalloc(tmp, qrm_mat%n)
263 
264  ! compute A'*r
265  call qrm_matmul(qrm_mat, 't', _qrm_one, b, _qrm_zero, tmp)
266 
267  call qrm_vecnrm(b, qrm_mat%m, '2', rnrm)
268  call qrm_vecnrm(tmp, qrm_mat%n, '2', nrm)
269  nrm = nrm/rnrm
270 
271  call qrm_adealloc(tmp)
272 
273  call qrm_err_act_restore(err_act)
274  return
275 
276 9999 continue ! error management
277  call qrm_err_act_restore(err_act)
278  if(err_act .eq. qrm_abort_) then
279  call qrm_err_check()
280  end if
281  return
282 
283 end subroutine _qrm_residual_and_orth1d
This module contains all the error management routines and data.
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_residual_orth2d(qrm_mat, r, nrm)
This routine computes the scaled norm of the product A'*r for multiple residuals. ...
subroutine qrm_err_act_restore(err_act)
Restores the value of the qrm_err_act variable.
subroutine _qrm_residual_and_orth1d(qrm_mat, b, x, nrm)
This routine computes the scaled norm of the product A'*r.
subroutine _qrm_residual_and_orth2d(qrm_mat, b, x, nrm)
This routine computes the scaled norm of the product A'*r.
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.
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
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...
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38
subroutine _qrm_residual_orth1d(qrm_mat, r, nrm)
This routine computes the scaled norm of the product A'*r.
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.