QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
dqrm_min_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 
49 subroutine dqrm_min_norm2d(qrm_mat, b, x)
50  use dqrm_spmat_mod
51  use qrm_mem_mod
52  use qrm_error_mod
55  use dqrm_solve_mod
56  implicit none
57 
58  type(dqrm_spmat_type) :: qrm_mat
59  real(kind(1.d0)) :: b(:,:), x(:,:)
60 
61  ! error management
62  integer :: err_act
63  character(len=*), parameter :: name='qrm_min_norm'
64 
65  call qrm_err_act_save(err_act)
66 
67  __qrm_prnt_dbg('("Entering the min-norm driver")')
68 
69  call dqrm_check_spmat(qrm_mat)
70  __qrm_check_ret(name,'qrm_check_spmat',9999)
71 
72  if(qrm_mat%n .lt. qrm_mat%m) then
73  call qrm_err_push(31, name,ied=(/qrm_mat%m,qrm_mat%n,0,0,0/))
74  goto 9999
75  end if
76 
77  ! analysis
78  call dqrm_analyse(qrm_mat, 't')
79  __qrm_check_ret(name,'qrm_analyse',9999)
80  ! factorization
81  call dqrm_factorize(qrm_mat, 't')
82  __qrm_check_ret(name,'qrm_factorize',9999)
83 
84  call dqrm_solve2d(qrm_mat, 't', b, x)
85  __qrm_check_ret(name,'qrm_solve',9999)
86  call dqrm_apply2d(qrm_mat, 'n', x)
87  __qrm_check_ret(name,'qrm_apply',9999)
88 
89  call qrm_err_act_restore(err_act)
90  return
91 
92 9999 continue ! error management
93  call qrm_err_act_restore(err_act)
94  if(err_act .eq. qrm_abort_) then
95  call qrm_err_check()
96  end if
97  return
98 
99 end subroutine dqrm_min_norm2d
100 
101 
102 
103 
105 
117 subroutine dqrm_min_norm1d(qrm_mat, b, x)
118  use dqrm_spmat_mod
119  use qrm_mem_mod
120  use qrm_error_mod
123  use dqrm_solve_mod
124  implicit none
125 
126  type(dqrm_spmat_type) :: qrm_mat
127  real(kind(1.d0)) :: b(:), x(:)
128 
129  integer :: info
130  ! error management
131  integer :: err_act
132  character(len=*), parameter :: name='qrm_min_norm'
133 
134  call qrm_err_act_save(err_act)
135 
136  __qrm_prnt_dbg('("Entering the least-squares driver")')
137 
138  call dqrm_check_spmat(qrm_mat)
139  __qrm_check_ret(name,'qrm_check_spmat',9999)
140 
141  if(qrm_mat%n .lt. qrm_mat%m) then
142  call qrm_err_push(31, name,ied=(/qrm_mat%m,qrm_mat%n,0,0,0/))
143  goto 9999
144  end if
145 
146  ! analysis
147  call dqrm_analyse(qrm_mat, 't')
148  __qrm_check_ret(name,'qrm_analyse',9999)
149 
150  ! factorization
151  call dqrm_factorize(qrm_mat, 't')
152  __qrm_check_ret(name,'qrm_factorize',9999)
153 
154  call dqrm_solve1d(qrm_mat, 't', b, x)
155  __qrm_check_ret(name,'qrm_solve',9999)
156  call dqrm_apply1d(qrm_mat, 'n', x)
157  __qrm_check_ret(name,'qrm_apply',9999)
158 
159  call qrm_err_act_restore(err_act)
160  return
161 
162 9999 continue ! error management
163  call qrm_err_act_restore(err_act)
164  if(err_act .eq. qrm_abort_) then
165  call qrm_err_check()
166  end if
167  return
168 
169 end subroutine dqrm_min_norm1d
This module contains all the error management routines and data.
subroutine dqrm_solve2d(qrm_mat, transp, b, x)
This function solves for R or R' against multiple vectors.
Definition: dqrm_solve.F90:50
This module contains all the interfaces for the typed routines in the solve phase.
subroutine dqrm_apply1d(qrm_mat, transp, b)
This function applies Q or Q^T to a single vector.
Definition: dqrm_apply.F90:135
subroutine qrm_err_act_restore(err_act)
Restores the value of the qrm_err_act variable.
subroutine dqrm_apply2d(qrm_mat, transp, b)
This function applies Q or Q^T to a set of vectors.
Definition: dqrm_apply.F90:47
This module contains the definition of the basic sparse matrix type and of the associated methods...
subroutine dqrm_min_norm2d(qrm_mat, b, x)
This routine computes the min-norm solution of a problem.
subroutine qrm_err_check()
This subroutine checks the errors stack. If something is found all the entries in the stack are poppe...
subroutine dqrm_analyse(qrm_mat, transp)
This is the driver routine for the analysis phase.
This module contains all the generic interfaces for the typed routines in the factorization phase...
This type defines the data structure used to store a matrix.
This module contains the generic interfaces for all the analysis routines.
subroutine dqrm_check_spmat(qrm_spmat, op)
Check the compatibility and correctness of icntl and rcntl parameters.
subroutine dqrm_min_norm1d(qrm_mat, b, x)
This routine computes the min-norm solution of a problem.
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38
subroutine dqrm_solve1d(qrm_mat, transp, b, x)
This function solves for R or R' against a single vector.
Definition: dqrm_solve.F90:128
subroutine qrm_err_push(code, sub, ied, aed)
This subroutine pushes an error on top of the stack.
subroutine dqrm_factorize(qrm_mat, transp)
This routine is the main factorization driver.
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.