* Bug #7927 fixed - Output "flag" in qmr function was not well documented.
[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 residual 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                     <itemizedlist>
169                         <listitem>
170                             <para><varname>flag</varname>=0: <code>qmr</code> converged to the desired tolerance within <varname>maxi</varname>
171                             iterations,</para>
172                         </listitem>
173                         <listitem>
174                             <para><varname>flag</varname>=1: no convergence given <varname>maxi</varname>,</para>
175                         </listitem>
176                         <listitem>
177                             <para><literal>-7 &lt; flag &lt; 0</literal>: A breakdown occurred because one of the scalar quantities calculated during
178                             <code>qmr</code> was equal to zero.</para>
179                         </listitem>
180                     </itemizedlist>
181                 </listitem>
182             </varlistentry>
183             <varlistentry>
184                 <term>res</term>
185                 <listitem>
186                     <para>residual vector</para>
187                 </listitem>
188             </varlistentry>
189             <varlistentry>
190                 <term>err</term>
191                 <listitem>
192                     <para>final residual norm</para>
193                 </listitem>
194             </varlistentry>
195             <varlistentry>
196                 <term>iter</term>
197                 <listitem>
198                     <para>number of iterations performed</para>
199                 </listitem>
200             </varlistentry>
201         </variablelist>
202     </refsection>
203     <refsection>
204         <title>Description</title>
205         <para>
206             Solves the linear system <literal>Ax=b</literal> using the Quasi Minimal Residual Method with preconditioning.
207         </para>
208     </refsection>
209     <refsection>
210         <title>Examples</title>
211         <programlisting role="example"><![CDATA[ 
212         // If A is a matrix
213 A=[ 94  0   0   0    0   28  0   0   32  0  
214      0   59  13  5    0   0   0   10  0   0  
215      0   13  72  34   2   0   0   0   0   65 
216      0   5   34  114  0   0   0   0   0   55 
217      0   0   2   0    70  0   28  32  12  0  
218      28  0   0   0    0   87  20  0   33  0  
219      0   0   0   0    28  20  71  39  0   0  
220      0   10  0   0    32  0   39  46  8   0  
221      32  0   0   0    12  33  0   8   82  11 
222      0   0   65  55   0   0   0   0   11  100];
223 b=ones(10,1);
224 [x,flag,err,iter,res] = qmr(A, b)
225
226 [x,flag,err,iter,res] = qmr(A, b, zeros(10,1), eye(10,10), eye(10,10), 10, 1d-12)
227
228         // If A is a function
229         function y = Atimesx(x,t)
230         A=[ 94  0   0   0    0   28  0   0   32  0  
231      0   59  13  5    0   0   0   10  0   0  
232      0   13  72  34   2   0   0   0   0   65 
233      0   5   34  114  0   0   0   0   0   55 
234      0   0   2   0    70  0   28  32  12  0  
235      28  0   0   0    0   87  20  0   33  0  
236      0   0   0   0    28  20  71  39  0   0  
237      0   10  0   0    32  0   39  46  8   0  
238      32  0   0   0    12  33  0   8   82  11 
239      0   0   65  55   0   0   0   0   11  100];
240          if (t == 'notransp') then
241         y = A*x;
242     elseif (t ==  'transp') then
243         y = A'*x;
244     end
245         endfunction
246          
247          [x,flag,err,iter,res] = qmr(Atimesx, b)
248          
249          [x,flag,err,iter,res] = qmr(Atimesx, b, zeros(10,1), eye(10,10), eye(10,10), 10, 1d-12)
250          
251          // OR
252          
253          function y = funA(x)
254         A = [ 94  0   0   0    0   28  0   0   32  0  
255      0   59  13  5    0   0   0   10  0   0  
256      0   13  72  34   2   0   0   0   0   65 
257      0   5   34  114  0   0   0   0   0   55 
258      0   0   2   0    70  0   28  32  12  0  
259      28  0   0   0    0   87  20  0   33  0  
260      0   0   0   0    28  20  71  39  0   0  
261      0   10  0   0    32  0   39  46  8   0  
262      32  0   0   0    12  33  0   8   82  11 
263      0   0   65  55   0   0   0   0   11  100];
264          y = A*x
265         endfunction
266         
267          function y = funAp(x)
268         A = [ 94  0   0   0    0   28  0   0   32  0  
269      0   59  13  5    0   0   0   10  0   0  
270      0   13  72  34   2   0   0   0   0   65 
271      0   5   34  114  0   0   0   0   0   55 
272      0   0   2   0    70  0   28  32  12  0  
273      28  0   0   0    0   87  20  0   33  0  
274      0   0   0   0    28  20  71  39  0   0  
275      0   10  0   0    32  0   39  46  8   0  
276      32  0   0   0    12  33  0   8   82  11 
277      0   0   65  55   0   0   0   0   11  100];
278          y = A'*x
279         endfunction
280          
281          [x,flag,err,iter,res] = qmr(funA, funAp, b)
282          
283          [x,flag,err,iter,res] = qmr(funA, funAp, b, zeros(10,1), eye(10,10), eye(10,10), 10, 1d-12)
284          
285          // If A is a matrix, M1 and M2 are functions
286          function y = M1timesx(x,t)
287          M1 = eye(10,10);
288     if(t=="notransp") then
289         y = M1*x;
290     elseif (t=="transp") then
291         y = M1'*x;
292     end
293         endfunction
294         
295         function y = M2timesx(x,t)
296          M2 = eye(10,10);
297     if(t=="notransp") then
298         y = M2*x;
299     elseif (t=="transp") then
300         y = M2'*x;
301     end
302         endfunction
303         
304         [x,flag,err,iter,res] = qmr(A, b, zeros(10,1), M1timesx, M2timesx, 10, 1d-12)
305         
306         // OR
307         
308         function y = funM1(x)
309         M1 = eye(10,10);
310         y = M1*x;
311         endfunction
312         
313         function y = funM1p(x)
314         M1 = eye(10,10);
315         y = M1'*x;
316         endfunction
317         
318         function y = funM2(x)
319         M2 = eye(10,10);
320         y = M2*x;
321         endfunction
322         
323         function y = funM2p(x)
324         M2 = eye(10,10);
325         y = M2'*x;
326         endfunction
327         
328         [x,flag,err,iter,res] = qmr(A, b, zeros(10,1), funM1, funM1p, funM2, funM2p, 10, 1d-12)
329         
330         // If A, M1, M2 are functions
331         [x,flag,err,iter,res] = qmr(funA, funAp, b, zeros(10,1), funM1, funM1p, funM2, funM2p, 10, 1d-12)
332         [x,flag,err,iter,res] = qmr(Atimesx, b, zeros(10,1), M1timesx, M2timesx, 10, 1d-12)
333
334  ]]></programlisting>
335     </refsection>
336     <refsection role="see also">
337         <title>See Also</title>
338         <simplelist type="inline">
339             <member>
340                 <link linkend="gmres">gmres</link>
341             </member>
342             <member>
343                 <link linkend="pcg">pcg</link>
344             </member>
345         </simplelist>
346     </refsection>
347     <refsection>
348         <title>History</title>
349         <revhistory>
350             <revision>
351                 <revnumber>5.4.0</revnumber>
352                 <revdescription>
353                     Calling qmr(A, Ap) is deprecated. qmr(A) should be used instead. The following function is an example :
354                     <programlisting role=""><![CDATA[ 
355 function y = A ( x, t )
356 Amat = eye(2,2);
357 if ( t== "notransp") then
358 y = Amat*x;
359 elseif (t == "transp") then
360 y = Amat'*x;
361 end
362 endfunction
363  ]]></programlisting>
364                 </revdescription>
365             </revision>
366         </revhistory>
367     </refsection>
368 </refentry>