35 #include "qrm_common.h"
56 real(kind(1.d0)) :: r(:,:)
57 real(kind(1.d0)) :: nrm(:)
59 real(kind(1.d0)),
allocatable :: tmp(:,:)
60 real(kind(1.d0)),
allocatable :: rnrm(:)
61 real(kind(1.d0)) :: anrm
65 character(len=*),
parameter :: name=
'qrm_residual_orth'
70 __qrm_check_ret(name,
'qrm_check_spmat',9999)
76 call
qrm_matmul(qrm_mat,
't', 1.d0, r, 0.d0, tmp)
91 if(err_act .eq. qrm_abort_)
then
117 real(kind(1.d0)) :: r(:)
118 real(kind(1.d0)) :: nrm
120 real(kind(1.d0)),
allocatable :: tmp(:)
121 real(kind(1.d0)) :: rnrm
122 real(kind(1.d0)) :: anrm
126 character(len=*),
parameter :: name=
'qrm_residual_orth'
131 __qrm_check_ret(name,
'qrm_check_spmat',9999)
136 call
qrm_matmul(qrm_mat,
't', 1.d0, r, 0.d0, tmp)
142 nrm = nrm/(rnrm*anrm)
150 if(err_act .eq. qrm_abort_)
then
181 real(kind(1.d0)) :: b(:,:), x(:,:)
182 real(kind(1.d0)) :: nrm(:)
184 real(kind(1.d0)),
allocatable :: tmp(:,:)
185 real(kind(1.d0)),
allocatable :: rnrm(:)
189 character(len=*),
parameter :: name=
'qrm_residual_orth'
194 __qrm_check_ret(name,
'qrm_check_spmat',9999)
197 call
qrm_matmul(qrm_mat,
'n', -1.d0, x, 1.d0, b)
203 call
qrm_matmul(qrm_mat,
't', 1.d0, b, 0.d0, tmp)
216 if(err_act .eq. qrm_abort_)
then
244 real(kind(1.d0)) :: b(:), x(:)
245 real(kind(1.d0)) :: nrm
247 real(kind(1.d0)),
allocatable :: tmp(:)
248 real(kind(1.d0)) :: rnrm
252 character(len=*),
parameter :: name=
'qrm_residual_orth'
257 __qrm_check_ret(name,
'qrm_check_spmat',9999)
260 call
qrm_matmul(qrm_mat,
'n', -1.d0, x, 1.d0, b)
265 call
qrm_matmul(qrm_mat,
't', 1.d0, b, 0.d0, tmp)
278 if(err_act .eq. qrm_abort_)
then
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.
Generic interface for the ::dqrm_matmul2d and ::dqrm_matmul1d routines.
subroutine dqrm_residual_and_orth1d(qrm_mat, b, x, nrm)
This routine computes the scaled norm of the product A'*r.
subroutine qrm_err_act_restore(err_act)
Restores the value of the qrm_err_act variable.
This module contains the definition of the basic sparse matrix type and of the associated methods...
Generic interface for the ::dqrm_vecnrm2d and ::dqrm_vecnrm1d routines.
This module contains generic interfaces for a number of auxiliary tools.
subroutine dqrm_residual_orth1d(qrm_mat, r, nrm)
This routine computes the scaled norm of the product A'*r.
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_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.
subroutine dqrm_residual_and_orth2d(qrm_mat, b, x, nrm)
This routine computes the scaled norm of the product A'*r.
This type defines the data structure used to store a matrix.
subroutine dqrm_check_spmat(qrm_spmat, op)
Check the compatibility and correctness of icntl and rcntl parameters.
subroutine dqrm_residual_orth2d(qrm_mat, r, nrm)
This routine computes the scaled norm of the product A'*r for multiple residuals. ...
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.