e14dafd3001bdf4ad6f945216242e6f6a53ed50f
[scilab.git] / scilab / modules / sparse / help / en_US / iterativesolvers / qmr.xml
1 <?xml version="1.0" encoding="UTF-8"?>
2 <!--
3  * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
4  * Copyright (C) XXXX-2008 - INRIA
5  * 
6  * This file must be used under the terms of the CeCILL.
7  * This source file is licensed as described in the file COPYING, which
8  * you should have received as part of this distribution.  The terms
9  * are also available at    
10  * http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
11  *
12  -->
13 <refentry xmlns="http://docbook.org/ns/docbook" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:svg="http://www.w3.org/2000/svg" xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:db="http://docbook.org/ns/docbook" xmlns:scilab="http://www.scilab.org" xml:lang="en" xml:id="qmr">
14     <refnamediv>
15         <refname>qmr</refname>
16         <refpurpose>quasi minimal resiqual method with preconditioning  </refpurpose>
17     </refnamediv>
18     <refsynopsisdiv>
19         <title>Calling Sequence</title>
20         <synopsis>[x,flag,err,iter,res] = qmr(A,Ap,b,x0,M1,M1p,M2,M2p,maxi,tol)
21             [x,flag,err,iter,res] = qmr(A,b,x0,M1,M2,maxi,tol)
22         </synopsis>
23     </refsynopsisdiv>
24     <refsection>
25         <title>Arguments</title>
26         <variablelist>
27             <varlistentry>
28                 <term>A</term>
29                 <listitem>
30                     <para>
31                         matrix of size n-by-n or function.                 
32                     </para>
33                     <itemizedlist>
34                         <listitem>
35                             <para>
36                                 <literal>matrix.</literal>If A is a matrix, it can be
37                                 dense or sparse
38                             </para>
39                         </listitem>
40                         <listitem>
41                             <para>
42                                 <literal>function.</literal>If A is a function which returns <literal>A*x</literal>, it must
43                                 have the following header :
44                             </para>
45                             <programlisting role=""><![CDATA[ 
46 function y = A ( x )
47  ]]></programlisting>
48                             <para>
49                                 If A is a function which returns <literal>A*x</literal> or <literal>A'*x</literal> depending t. 
50                                 If <literal>t = "notransp"</literal>, the function returns <literal>A*x</literal>. 
51                                 If <literal>t = "transp"</literal>, the function returns <literal>A'*x</literal>. It must
52                                 have the following header :
53                             </para>
54                             <programlisting role=""><![CDATA[ 
55 function y = A ( x, t )
56  ]]></programlisting>
57                         </listitem>
58                     </itemizedlist>
59                 </listitem>
60             </varlistentry>
61             <varlistentry>
62                 <term>Ap</term>
63                 <listitem>
64                     <para>
65                         function returning <literal>A'*x</literal>. It must have the followinf header :
66                     </para>
67                     <programlisting role=""><![CDATA[ 
68 function y = Ap ( x )
69  ]]></programlisting>
70                 </listitem>
71             </varlistentry>
72             <varlistentry>
73                 <term>b</term>
74                 <listitem>
75                     <para>right hand side vector</para>
76                 </listitem>
77             </varlistentry>
78             <varlistentry>
79                 <term>x0</term>
80                 <listitem>
81                     <para>initial guess vector (default: zeros(n,1))</para>
82                 </listitem>
83             </varlistentry>
84             <varlistentry>
85                 <term>M1</term>
86                 <listitem>
87                     <para>
88                         left preconditioner : matrix or function (In the first case, default: eye(n,n)). If <literal>M1</literal> is a function, she returns either,
89                     </para>
90                     <itemizedlist>
91                         <listitem>
92                             <para>
93                                 only <literal>M1*x</literal>
94                             </para>
95                         </listitem>
96                         <para> or
97                         </para>
98                         <listitem>
99                             <para>
100                                 <literal>M1*x</literal> or <literal>M1'*x</literal> depending <literal>t</literal>.
101                             </para>
102                         </listitem>
103                     </itemizedlist>
104                 </listitem>
105             </varlistentry>
106             <varlistentry>
107                 <term>M1p</term>
108                 <listitem>
109                     <para>
110                         must only be provided when <literal>M1</literal> is a function returning  <literal>M1*x</literal>. 
111                         In this case <literal>M1p</literal> is the function which returns <literal>M1'*x</literal>. 
112                     </para>
113                 </listitem>
114             </varlistentry>
115             <varlistentry>
116                 <term>M2</term>
117                 <listitem>
118                     <para>
119                         right preconditioner : matrix or function (In the first case, default: eye(n,n)). If <literal>M2</literal> is a function, she returns either
120                     </para>
121                     <itemizedlist>
122                         <listitem>
123                             <para>
124                                 only <literal>M2*x</literal>
125                             </para>
126                         </listitem>
127                         <para> or
128                         </para>
129                         <listitem>
130                             <para>
131                                 <literal>M2*x</literal> or <literal>M2'*x</literal> depending <literal>t</literal>.
132                             </para>
133                         </listitem>
134                     </itemizedlist>
135                 </listitem>
136             </varlistentry>
137             <varlistentry>
138                 <term>M2p</term>
139                 <listitem>
140                     <para>
141                         must only be provided when <literal>M2</literal> is a function returning  <literal>M2*x</literal>. 
142                         In this case <literal>M2p</literal> is the function which returns <literal>M2'*x</literal>
143                     </para>
144                 </listitem>
145             </varlistentry>
146             <varlistentry>
147                 <term>maxi</term>
148                 <listitem>
149                     <para>maximum number of iterations (default: n)
150                     </para>
151                 </listitem>
152             </varlistentry>
153             <varlistentry>
154                 <term>tol</term>
155                 <listitem>
156                     <para>error tolerance (default: 1000*%eps)</para>
157                 </listitem>
158             </varlistentry>
159             <varlistentry>
160                 <term>x</term>
161                 <listitem>
162                     <para>solution vector</para>
163                 </listitem>
164             </varlistentry>
165             <varlistentry>
166                 <term>flag</term>
167                 <listitem>
168                     <variablelist>
169                         <varlistentry>
170                             <term>0 =</term>
171                             <listitem>
172                                 <para>
173                                     <literal>gmres</literal> converged to the desired tolerance within <literal>maxi</literal> iterations
174                                 </para>
175                             </listitem>
176                         </varlistentry>
177                         <varlistentry>
178                             <term>1 =</term>
179                             <listitem>
180                                 <para>
181                                     no convergence given <literal>maxi</literal>
182                                 </para>
183                             </listitem>
184                         </varlistentry>
185                     </variablelist>
186                 </listitem>
187             </varlistentry>
188             <varlistentry>
189                 <term>res</term>
190                 <listitem>
191                     <para>residual vector</para>
192                 </listitem>
193             </varlistentry>
194             <varlistentry>
195                 <term>err</term>
196                 <listitem>
197                     <para>final residual norm</para>
198                 </listitem>
199             </varlistentry>
200             <varlistentry>
201                 <term>iter</term>
202                 <listitem>
203                     <para>number of iterations performed</para>
204                 </listitem>
205             </varlistentry>
206         </variablelist>
207     </refsection>
208     <refsection>
209         <title>Description</title>
210         <para>
211             Solves the linear system <literal>Ax=b</literal> using the Quasi Minimal Residual Method with preconditioning.
212         </para>
213     </refsection>
214     <refsection>
215         <title>Examples</title>
216         <programlisting role="example"><![CDATA[ 
217         // If A is a matrix
218 A=[ 94  0   0   0    0   28  0   0   32  0  
219      0   59  13  5    0   0   0   10  0   0  
220      0   13  72  34   2   0   0   0   0   65 
221      0   5   34  114  0   0   0   0   0   55 
222      0   0   2   0    70  0   28  32  12  0  
223      28  0   0   0    0   87  20  0   33  0  
224      0   0   0   0    28  20  71  39  0   0  
225      0   10  0   0    32  0   39  46  8   0  
226      32  0   0   0    12  33  0   8   82  11 
227      0   0   65  55   0   0   0   0   11  100];
228 b=ones(10,1);
229 [x,flag,err,iter,res] = qmr(A, b)
230
231 [x,flag,err,iter,res] = qmr(A, b, zeros(10,1), eye(10,10), eye(10,10), 10, 1d-12)
232
233         // If A is a function
234         function y = Atimesx(x,t)
235         A=[ 94  0   0   0    0   28  0   0   32  0  
236      0   59  13  5    0   0   0   10  0   0  
237      0   13  72  34   2   0   0   0   0   65 
238      0   5   34  114  0   0   0   0   0   55 
239      0   0   2   0    70  0   28  32  12  0  
240      28  0   0   0    0   87  20  0   33  0  
241      0   0   0   0    28  20  71  39  0   0  
242      0   10  0   0    32  0   39  46  8   0  
243      32  0   0   0    12  33  0   8   82  11 
244      0   0   65  55   0   0   0   0   11  100];
245          if (t == 'notransp') then
246         y = A*x;
247     elseif (t ==  'transp') then
248         y = A'*x;
249     end
250         endfunction
251          
252          [x,flag,err,iter,res] = qmr(Atimesx, b)
253          
254          [x,flag,err,iter,res] = qmr(Atimesx, b, zeros(10,1), eye(10,10), eye(10,10), 10, 1d-12)
255          
256          // OR
257          
258          function y = funA(x)
259         A = [ 94  0   0   0    0   28  0   0   32  0  
260      0   59  13  5    0   0   0   10  0   0  
261      0   13  72  34   2   0   0   0   0   65 
262      0   5   34  114  0   0   0   0   0   55 
263      0   0   2   0    70  0   28  32  12  0  
264      28  0   0   0    0   87  20  0   33  0  
265      0   0   0   0    28  20  71  39  0   0  
266      0   10  0   0    32  0   39  46  8   0  
267      32  0   0   0    12  33  0   8   82  11 
268      0   0   65  55   0   0   0   0   11  100];
269          y = A*x
270         endfunction
271         
272          function y = funAp(x)
273         A = [ 94  0   0   0    0   28  0   0   32  0  
274      0   59  13  5    0   0   0   10  0   0  
275      0   13  72  34   2   0   0   0   0   65 
276      0   5   34  114  0   0   0   0   0   55 
277      0   0   2   0    70  0   28  32  12  0  
278      28  0   0   0    0   87  20  0   33  0  
279      0   0   0   0    28  20  71  39  0   0  
280      0   10  0   0    32  0   39  46  8   0  
281      32  0   0   0    12  33  0   8   82  11 
282      0   0   65  55   0   0   0   0   11  100];
283          y = A'*x
284         endfunction
285          
286          [x,flag,err,iter,res] = qmr(funA, funAp, b)
287          
288          [x,flag,err,iter,res] = qmr(funA, funAp, b, zeros(10,1), eye(10,10), eye(10,10), 10, 1d-12)
289          
290          // If A is a matrix, M1 and M2 are functions
291          function y = M1timesx(x,t)
292          M1 = eye(10,10);
293     if(t=="notransp") then
294         y = M1*x;
295     elseif (t=="transp") then
296         y = M1'*x;
297     end
298         endfunction
299         
300         function y = M2timesx(x,t)
301          M2 = eye(10,10);
302     if(t=="notransp") then
303         y = M2*x;
304     elseif (t=="transp") then
305         y = M2'*x;
306     end
307         endfunction
308         
309         [x,flag,err,iter,res] = qmr(A, b, zeros(10,1), M1timesx, M2timesx, 10, 1d-12)
310         
311         // OR
312         
313         function y = funM1(x)
314         M1 = eye(10,10);
315         y = M1*x;
316         endfunction
317         
318         function y = funM1p(x)
319         M1 = eye(10,10);
320         y = M1'*x;
321         endfunction
322         
323         function y = funM2(x)
324         M2 = eye(10,10);
325         y = M2*x;
326         endfunction
327         
328         function y = funM2p(x)
329         M2 = eye(10,10);
330         y = M2'*x;
331         endfunction
332         
333         [x,flag,err,iter,res] = qmr(A, b, zeros(10,1), funM1, funM1p, funM2, funM2p, 10, 1d-12)
334         
335         // If A, M1, M2 are functions
336         [x,flag,err,iter,res] = qmr(funA, funAp, b, zeros(10,1), funM1, funM1p, funM2, funM2p, 10, 1d-12)
337         [x,flag,err,iter,res] = qmr(Atimesx, b, zeros(10,1), M1timesx, M2timesx, 10, 1d-12)
338
339  ]]></programlisting>
340     </refsection>
341     <refsection role="see also">
342         <title>See Also</title>
343         <simplelist type="inline">
344             <member>
345                 <link linkend="gmres">gmres</link>
346             </member>
347             <member>
348                 <link linkend="pcg">pcg</link>
349             </member>
350         </simplelist>
351     </refsection>
352     <refsection>
353         <title>History</title>
354         <revhistory>
355             <revision>
356                 <revnumber>5.4.0</revnumber>
357                 <revdescription>
358                     Calling qmr(A, Ap) is deprecated. qmr(A) should be used instead. The following function is an example :
359                     <programlisting role=""><![CDATA[ 
360 function y = A ( x, t )
361 Amat = eye(2,2);
362 if ( t== "notransp") then
363 y = Amat*x;
364 elseif (t == "transp") then
365 y = Amat'*x;
366 end
367 endfunction
368  ]]></programlisting>
369                 </revdescription>
370             </revision>
371         </revhistory>
372     </refsection>
373 </refentry>