QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
qrm_apply.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 
47 subroutine _qrm_apply2d(qrm_mat, transp, b)
48 
49  use _qrm_spmat_mod
50  use qrm_common_mod
51  use qrm_string_mod
52  use _qrm_solve_mod, protect => _qrm_apply2d
53 #if defined(_OPENMP)
54  use omp_lib
55 #endif
56  implicit none
57 
58  type(_qrm_spmat_type), intent(in) :: qrm_mat
59  _qrm_data, intent(inout) :: b(:,:)
60  character(len=*), intent(in) :: transp
61 
62  integer :: nb, rhs_nthreads, nrhs, i, keeph
63  ! error management
64  integer :: err_act
65  character(len=*), parameter :: name='qrm_apply'
66 
67  call qrm_err_act_save(err_act)
68 
69  ! immediately check if the facto was done. Otherwise push an error and return
70  if(.not. qrm_mat%fdata%ok) then
71  call qrm_err_push(14, 'qrm_apply')
72  goto 9999
73  end if
74 
75  call qrm_get(qrm_mat, 'qrm_keeph', keeph)
76  if(keeph .eq. qrm_no_) then
77  call qrm_err_push(30, 'qrm_apply')
78  goto 9999
79  end if
80 
81  ! blocking to deal with multiple rhs
82  call qrm_get(qrm_mat, 'qrm_rhsnb', nb)
83  call qrm_get(qrm_mat, 'qrm_rhsnthreads', rhs_nthreads)
84  nrhs = size(b,2)
85  if(nb.le.0) nb = nrhs
86 
87 #if defined(_OPENMP)
88  call omp_set_nested(.true.)
89 #endif
90  if(qrm_str_tolower(transp(1:1)) .eq. 't') then
91 
92  !$omp parallel do num_threads(rhs_nthreads) private(i)
93  do i=1, nrhs, nb
94  call qrm_apply_qt(qrm_mat, b(:,i:min(nrhs,i+nb-1)))
95  end do
96  !$omp end parallel do
97  __qrm_check_ret(name,'qrm_apply_qt',9999)
98 
99  else
100  !$omp parallel do num_threads(rhs_nthreads) private(i)
101  do i=1, nrhs, nb
102  call qrm_apply_q(qrm_mat, b(:,i:min(nrhs,i+nb-1)))
103  end do
104  !$omp end parallel do
105  __qrm_check_ret(name,'qrm_apply_q',9999)
106  end if
107 
108 #if defined(_OPENMP)
109  call omp_set_nested(.false.)
110 #endif
111 
112  call qrm_err_act_restore(err_act)
113  return
114 
115 9999 continue ! error management
116  call qrm_err_act_restore(err_act)
117  if(err_act .eq. qrm_abort_) then
118  call qrm_err_check()
119  end if
120  return
121 
122 end subroutine _qrm_apply2d
123 
124 
135 subroutine _qrm_apply1d(qrm_mat, transp, b)
136 
137  use _qrm_spmat_mod
138  use qrm_common_mod
139  use qrm_string_mod
140  use _qrm_solve_mod, protect => _qrm_apply1d
141  use _qrm_utils_mod
142  implicit none
143 
144  type(_qrm_spmat_type), intent(in) :: qrm_mat
145  _qrm_data, intent(inout) :: b(:)
146  character(len=*), intent(in) :: transp
147 
148  _qrm_data, pointer :: pnt(:,:)
149  integer :: n, keeph
150  ! error management
151  integer :: err_act
152  character(len=*), parameter :: name='qrm_apply1d'
153 
154  call qrm_err_act_save(err_act)
155 
156  ! immediately check if the facto was done. Otherwise push an error and return
157  if( qrm_mat%fdata%done .eq. 0) then
158  call qrm_err_push(14, 'qrm_apply1d')
159  goto 9999
160  end if
161 
162  call qrm_get(qrm_mat, 'qrm_keeph', keeph)
163  if(keeph .eq. qrm_no_) then
164  call qrm_err_push(30, 'qrm_apply')
165  goto 9999
166  end if
167 
168 
169  n = size(b,1)
170 
171  call _qrm_remap_pnt(b, pnt, n)
172 
173  if(qrm_str_tolower(transp(1:1)) .eq. 't') then
174  call qrm_apply_qt(qrm_mat, pnt)
175  __qrm_check_ret(name,'qrm_apply_qt',9999)
176  else
177  call qrm_apply_q(qrm_mat, pnt)
178  __qrm_check_ret(name,'qrm_apply_q',9999)
179  end if
180 
181  call qrm_err_act_restore(err_act)
182  return
183 
184 9999 continue ! error management
185  call qrm_err_act_restore(err_act)
186  if(err_act .eq. qrm_abort_) then
187  call qrm_err_check()
188  end if
189  return
190 
191 end subroutine _qrm_apply1d
subroutine _qrm_remap_pnt(arr1d, pnt2d, n)
This function makes a 2D pointer point to a 1D array.
This module contains the interfaces of all non-typed routines.
subroutine qrm_err_act_restore(err_act)
Restores the value of the qrm_err_act variable.
Generic interface for the ::_qrm_apply_qt routine.
subroutine _qrm_apply2d(qrm_mat, transp, b)
This function applies Q or Q^T to a set of vectors.
Definition: qrm_apply.F90:47
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_apply_q routine.
subroutine _qrm_apply1d(qrm_mat, transp, b)
This function applies Q or Q^T to a single vector.
Definition: qrm_apply.F90:135
This module contains all the interfaces for the typed routines in the solve phase.
This type defines the data structure used to store a matrix.
This module contains the definition of the basic sparse matrix type and of the associated methods...
Generif interface for the ::_qrm_pgeti, ::_qrm_pgetr and.
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.