QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
qrm_solve_mod.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 
38 
40  interface qrm_apply_q
41  subroutine _qrm_apply_q(qrm_mat, b)
42  use _qrm_spmat_mod
43  type(_qrm_spmat_type), target :: qrm_mat
44  _qrm_data, intent(inout) :: b(:,:)
45  end subroutine _qrm_apply_q
46  end interface
47 
49  interface qrm_apply_qt
50  subroutine _qrm_apply_qt(qrm_mat, b)
51  use _qrm_spmat_mod
52  type(_qrm_spmat_type) :: qrm_mat
53  _qrm_data :: b(:,:)
54  end subroutine _qrm_apply_qt
55  end interface
56 
57 
62  interface _qrm_apply
63  module procedure _qrm_apply2dw, _qrm_apply1dw
64  end interface _qrm_apply
65 
66 
67  ! !> @brief Generic interface for the
68  ! !! @link ::_qrm_apply @endlink and
69  ! !! @link ::_qrm_apply1d @endlink routines
70  ! !!
71  interface qrm_apply
72  subroutine _qrm_apply2d(qrm_mat, transp, b)
73  use _qrm_spmat_mod
74  type(_qrm_spmat_type), intent(in) :: qrm_mat
75  _qrm_data, intent(inout) :: b(:,:)
76  character(len=*), intent(in) :: transp
77  end subroutine _qrm_apply2d
78  subroutine _qrm_apply1d(qrm_mat, transp, b)
79  use _qrm_spmat_mod
80  type(_qrm_spmat_type), intent(in) :: qrm_mat
81  _qrm_data, intent(inout) :: b(:)
82  character(len=*), intent(in) :: transp
83  end subroutine _qrm_apply1d
84  end interface qrm_apply
85 
87  interface qrm_solve_r
88  subroutine _qrm_solve_r(qrm_mat, b, x)
89  use _qrm_spmat_mod
90  type(_qrm_spmat_type), target :: qrm_mat
91  _qrm_data, intent(inout) :: b(:,:)
92  _qrm_data, intent(out) :: x(:,:)
93  end subroutine _qrm_solve_r
94  end interface
95 
97  interface qrm_solve_rt
98  subroutine _qrm_solve_rt(qrm_mat, b, x)
99  use _qrm_spmat_mod
100  type(_qrm_spmat_type), target :: qrm_mat
101  _qrm_data, intent(inout) :: b(:,:)
102  _qrm_data, intent(out) :: x(:,:)
103  end subroutine _qrm_solve_rt
104  end interface
105 
110  interface _qrm_solve
111  module procedure _qrm_solve1dw, _qrm_solve2dw
112  end interface
113 
118  interface qrm_solve
119  subroutine _qrm_solve2d(qrm_mat, transp, b, x)
120  use _qrm_spmat_mod
121  type(_qrm_spmat_type), target :: qrm_mat
122  _qrm_data, intent(inout) :: b(:,:)
123  _qrm_data, intent(out) :: x(:,:)
124  character(len=*) :: transp
125  end subroutine _qrm_solve2d
126  subroutine _qrm_solve1d(qrm_mat, transp, b, x)
127  use _qrm_spmat_mod
128  type(_qrm_spmat_type), target :: qrm_mat
129  _qrm_data, intent(inout) :: b(:)
130  _qrm_data, intent(out) :: x(:)
131  character(len=*) :: transp
132  end subroutine _qrm_solve1d
133  end interface
134 
137  subroutine _qrm_solve_sing_front(qrm_mat, b, x, trans)
138  use _qrm_spmat_mod
139  use _qrm_fdata_mod
140  type(_qrm_spmat_type), target :: qrm_mat
141  _qrm_data, intent(inout) :: b(:,:)
142  _qrm_data, intent(inout) :: x(:,:)
143  character :: trans
144  end subroutine _qrm_solve_sing_front
145  end interface
146 
147 contains
148 
149  subroutine _qrm_apply2dw(qrm_mat, transp, b)
150  use _qrm_spmat_mod
151  type(_qrm_spmat_type), intent(in) :: qrm_mat
152  _qrm_data, intent(inout) :: b(:,:)
153  character(len=*), intent(in) :: transp
154  call _qrm_apply2d(qrm_mat, transp, b)
155  return
156  end subroutine _qrm_apply2dw
157 
158  subroutine _qrm_apply1dw(qrm_mat, transp, b)
159  use _qrm_spmat_mod
160  type(_qrm_spmat_type), intent(in) :: qrm_mat
161  _qrm_data, intent(inout) :: b(:)
162  character(len=*), intent(in) :: transp
163  call _qrm_apply1d(qrm_mat, transp, b)
164  return
165  end subroutine _qrm_apply1dw
166 
167  subroutine _qrm_solve2dw(qrm_mat, transp, b, x)
168  use _qrm_spmat_mod
169  type(_qrm_spmat_type), target :: qrm_mat
170  _qrm_data, intent(inout) :: b(:,:)
171  _qrm_data, intent(out) :: x(:,:)
172  character(len=*) :: transp
173  call _qrm_solve2d(qrm_mat, transp, b, x)
174  return
175  end subroutine _qrm_solve2dw
176 
177  subroutine _qrm_solve1dw(qrm_mat, transp, b, x)
178  use _qrm_spmat_mod
179  type(_qrm_spmat_type), target :: qrm_mat
180  _qrm_data, intent(inout) :: b(:)
181  _qrm_data, intent(out) :: x(:)
182  character(len=*) :: transp
183  call _qrm_solve1d(qrm_mat, transp, b, x)
184  return
185  end subroutine _qrm_solve1dw
186 
187 
188 end module _qrm_solve_mod
Generic interface for the ::_qrm_solve_sing_front routine.
This module contains the definition of all the data related to the factorization phase.
Generic interface for the ::_qrm_solve and ::_qrm_solve1d routines.
subroutine _qrm_apply_q(qrm_mat, b)
This function applies Q to a vector/matrix.
Definition: qrm_apply_q.F90:44
subroutine _qrm_solve_r(qrm_mat, b, x)
This function solves for R against multiple vectors.
Definition: qrm_solve_r.F90:45
subroutine _qrm_solve2dw(qrm_mat, transp, b, x)
Generic interface for the ::_qrm_solve and ::_qrm_solve1d routines.
subroutine _qrm_apply2dw(qrm_mat, transp, b)
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
subroutine _qrm_solve_sing_front(qrm_mat, b, x, trans)
This function handles the front containing the singletons during the solve for R or R'...
Generic interface for the ::_qrm_apply_q routine.
subroutine _qrm_solve_rt(qrm_mat, b, x)
This function solves for R' against multiple vectors.
subroutine _qrm_apply_qt(qrm_mat, b)
This function applies Q' to a vector/matrix.
subroutine _qrm_solve1d(qrm_mat, transp, b, x)
This function solves for R or R' against a single vector.
Definition: qrm_solve.F90:128
subroutine _qrm_apply1d(qrm_mat, transp, b)
This function applies Q or Q^T to a single vector.
Definition: qrm_apply.F90:135
subroutine _qrm_apply1dw(qrm_mat, transp, b)
Generic interface for the ::_qrm_solve_r routine.
Generic interface for the ::_qrm_solve_rt routine.
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.
subroutine _qrm_solve1dw(qrm_mat, transp, b, x)
This module contains the definition of the basic sparse matrix type and of the associated methods...
Generic interface for the ::_qrm_apply and ::_qrm_apply1d routines.
subroutine _qrm_solve2d(qrm_mat, transp, b, x)
This function solves for R or R' against multiple vectors.
Definition: qrm_solve.F90:50