Planet
navi homePPSaboutscreenshotsdownloaddevelopmentforum

source: downloads/boost_1_34_1/libs/numeric/ublas/doc/operations_overview.htm @ 29

Last change on this file since 29 was 29, checked in by landauf, 17 years ago

updated boost from 1_33_1 to 1_34_1
File size: 9.3 KB
Line 
1<!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Transitional//EN"
2    "http://www.w3.org/TR/xhtml1/DTD/xhtml1-transitional.dtd">
3<html xmlns="http://www.w3.org/1999/xhtml">
4<head>
5<meta name="generator" content=
6"HTML Tidy for Linux/x86 (vers 1st March 2004), see www.w3.org" />
7<meta name="GENERATOR" content="Quanta Plus" />
8<meta http-equiv="Content-Type" content=
9"text/html; charset=us-ascii" />
10<link href="ublas.css" type="text/css" />
11<title>uBLAS operations overview</title>
12</head>
13<body>
14<h1><img src="../../../../boost.png" align="middle" />
15Overview of Matrix and Vector Operations</h1>
16
17<dl>
18<dt>Contents:</dt>
19<dd><a href="#blas">Basic Linear Algebra</a></dd>
20<dd><a href="#advanced">Advanced Functions</a></dd>
21<dd><a href="#sub">Submatrices, Subvectors</a></dd>
22<dd><a href="#speed">Speed Improvements</a></dd>
23</dl>
24
25<h3>Definitions:</h3>
26
27<table style="" summary="notation">
28<tr><td><code>A, B, C</code></td>
29<td> are matrices</td></tr>
30<tr><td><code>u, v, w</code></td> 
31<td>are vectors</td></tr>
32<tr><td><code>i, j, k</code></td> 
33<td>are integer values</td></tr>
34<tr><td><code>t, t1, t2</code></td> 
35<td>are scalar values</td></tr>
36<tr><td><code>r, r1, r2</code></td> 
37<td>are <a href="range.htm">ranges</a>, e.g. <code>range(0, 3)</code></td></tr>
38<tr><td><code>s, s1, s2</code></td> 
39<td>are <a href="range.htm#slice">slices</a>, e.g. <code>slice(0, 1, 3)</code></td></tr>
40</table>
41
42<h2><a name="blas">Basic Linear Algebra</a></h2>
43
44<h3>standard operations: addition, subtraction, multiplication by a
45scalar</h3>
46
47<pre><code>
48C = A + B; C = A - B; C = -A;
49w = u + v; w = u - v; w = -u;
50C = t * A; C = A * t; C = A / t;
51w = t * u; w = u * t; w = u / t;
52</code></pre>
53
54<h3>computed assignements</h3>
55
56<pre><code>
57C += A; C -= A;
58w += u; w -= u;
59C *= t; C /= t;
60w *= t; w /= t;
61</code></pre>
62
63<h3>inner, outer and other products</h3>
64
65<pre><code>
66t = inner_prod(u, v);
67C = outer_prod(u, v);
68w = prod(A, u); w = prod(u, A); w = prec_prod(A, u); w = prec_prod(u, A);
69C = prod(A, B); C = prec_prod(A, B);
70w = element_prod(u, v); w = element_div(u, v);
71C = element_prod(A, B); C = element_div(A, B);
72</code></pre>
73
74<h3>transformations</h3>
75
76<pre><code>
77w = conj(u); w = real(u); w = imag(u);
78C = trans(A); C = conj(A); C = herm(A); C = real(A); C = imag(A);
79</code></pre>
80
81<h2><a name="advanced">Advanced functions</a></h2>
82
83<h3>norms</h3>
84
85<pre><code>
86t = norm_inf(v); i = index_norm_inf(v);
87t = norm_1(v);   t = norm_2(v);
88t = norm_inf(A); i = index_norm_inf(A);
89t = norm_1(A);   t = norm_frobenius(A);
90</code></pre>
91
92<h3>products</h3>
93
94<pre><code>
95axpy_prod(A, u, w, true);  // w = A * u
96axpy_prod(A, u, w, false); // w += A * u
97axpy_prod(u, A, w, true);  // w = trans(A) * u
98axpy_prod(u, A, w, false); // w += trans(A) * u
99axpy_prod(A, B, C, true);  // C = A * B
100axpy_prod(A, B, C, false); // C += A * B
101</code></pre>
102<p><em>Note:</em> The last argument (<code>bool init</code>) of
103<code>axpy_prod</code> is optional. Currently it defaults to
104<code>true</code>, but this may change in the future. Set the
105<code>init</code> to <code>true</code> is equivalent to call
106<code>w.clear()</code> before <code>axpy_prod</code>. Up to now
107there are some specialisation for compressed matrices that give a
108large speed up compared to <code>prod</code>.</p>
109<pre><code>
110w = block_prod&lt;matrix_type, 64&gt; (A, u); // w = A * u
111w = block_prod&lt;matrix_type, 64&gt; (u, A); // w = trans(A) * u
112C = block_prod&lt;matrix_type, 64&gt; (A, B); // w = A * B
113</code></pre>
114<p><em>Note:</em> The blocksize can be any integer. However, the
115total speed depends very strong on the combination of blocksize,
116CPU and compiler. The function <code>block_prod</code> is designed
117for large dense matrices.</p>
118<h3>rank-k updates</h3>
119<pre><code>
120opb_prod(A, B, C, true);  // C = A * B
121opb_prod(A, B, C, false); // C += A * B
122</code></pre>
123<p><em>Note:</em> The last argument (<code>bool init</code>) of
124<code>opb_prod</code> is optional. Currently it defaults to
125<code>true</code>, but this may change in the future. This function
126may give a speedup if <code>A</code> has less columns than rows,
127because the product is computed as a sum of outer products.</p>
128
129<h2><a name="sub">Submatrices, Subvectors</a></h2>
130<p>Accessing submatrices and subvectors via <b>proxies</b> using <code>project</code> functions:</p>
131<pre><code>
132w = project(u, r);         // the subvector of u specifed by the index range r
133w = project(u, s);         // the subvector of u specifed by the index slice s
134C = project(A, r1, r2);    // the submatrix of A specified by the two index ranges r1 and r2
135C = project(A, s1, s2);    // the submatrix of A specified by the two index slices s1 and s2
136w = row(A, i); w = column(A, j);    // a row or column of matrix as a vector
137</code></pre>
138<p>Assigning to submatrices and subvectors via <b>proxies</b> using <code>project</code> functions:</p>
139<pre><code>
140project(u, r) = w;         // assign the subvector of u specifed by the index range r
141project(u, s) = w;         // assign the subvector of u specifed by the index slice s
142project(A, r1, r2) = C;    // assign the submatrix of A specified by the two index ranges r1 and r2
143project(A, s1, s2) = C;    // assign the submatrix of A specified by the two index slices s1 and s2
144row(A, i) = w; column(A, j) = w;    // a row or column of matrix as a vector
145</code></pre>
146<p><em>Note:</em> A range <code>r = range(start, stop)</code>
147contains all indices <code>i</code> with <code>start &lt;= i &lt;
148stop</code>. A slice is something more general. The slice
149<code>s = slice(start, stride, size)</code> contains the indices
150<code>start, start+stride, ..., start+(size-1)*stride</code>. The
151stride can be 0 or negative! If <code>start >= stop</code> for a range
152or <code>size == 0</code> for a slice then it contains no elements.</p>
153<p>Sub-ranges and sub-slices of vectors and matrices can be created directly with the <code>subrange</code> and <code>sublice</code> functions:</p>
154<pre><code>
155w = subrange(u, 0, 2);         // the 2 element subvector of u
156w = subslice(u, 0, 1, 2);      // the 2 element subvector of u
157C = subrange(A, 0,2, 0,3);     // the 2x3 element submatrix of A
158C = subslice(A, 0,1,2, 0,1,3); // the 2x3 element submatrix of A
159subrange(u, 0, 2) = w;         // assign the 2 element subvector of u
160subslice(u, 0, 1, 2) = w;      // assign the 2 element subvector of u
161subrange(A, 0,2, 0,3) = C;     // assign the 2x3 element submatrix of A
162subrange(A, 0,1,2, 0,1,3) = C; // assigne the 2x3 element submatrix of A
163</code></pre>
164<p>There are to more ways to access some matrix elements as a
165vector:</p>
166<pre><code>matrix_vector_range&lt;matrix_type&gt; (A, r1, r2);
167matrix_vector_slice&lt;matrix_type&gt; (A, s1, s2);
168</code></pre>
169<p><em>Note:</em> These matrix proxies take a sequence of elements
170of a matrix and allow you to access these as a vector. In
171particular <code>matrix_vector_slice</code> can do this in a very
172general way. <code>matrix_vector_range</code> is less useful as the
173elements must lie along a diagonal.</p>
174<p><em>Example:</em> To access the first two elements of a sub
175column of a matrix we access the row with a slice with stride 1 and
176the column with a slice with stride 0 thus:<br />
177<code>matrix_vector_slice&lt;matrix_type&gt; (A, slice(0,1,2),
178slice(0,0,2));
179</code></p>
180
181<h2><a name="speed">Speed improvements</a></h2>
182<h3><a name='noalias'>Matrix / Vector assignment</a></h3>
183<p>If you know for sure that the left hand expression and the right
184hand expression have no common storage, then assignment has
185no <em>aliasing</em>. A more efficient assignment can be specified
186in this case:</p>
187<pre><code>noalias(C) = prod(A, B);
188</code></pre>
189<p>This avoids the creation of a temporary matrix that is required in a normal assignment.
190'noalias' assignment requires that the left and right hand side be size conformant.</p>
191
192<h3>Sparse element access</h3>
193<p>The matrix element access function <code>A(i1,i2)</code> or the equivalent vector
194element access functions (<code>v(i) or v[i]</code>) usually create 'sparse element proxies'
195when applied to a sparse matrix or vector. These <em>proxies</em> allow access to elements
196without having to worry about nasty C++ issues where references are invalidated.</p>
197<p>These 'sparse element proxies' can be implemented more efficiently when applied to <code>const</code>
198objects.
199Sadly in C++ there is no way to distinguish between an element access on the left and right hand side of
200an assignment. Most often elements on the right hand side will not be changed and therefore it would
201be better to use the <code>const</code> proxies. We can do this by making the matrix or vector
202<code>const</code> before accessing it's elements. For example:</p>
203<pre><code>value = const_cast&lt;const VEC&gt;(v)[i];   // VEC is the type of V
204</code></pre>
205<p>If more then one element needs to be accessed <code>const_iterator</code>'s should be used
206in preference to <code>iterator</code>'s for the same reason. For the more daring 'sparse element proxies'
207can be completely turned off in uBLAS by defining the configuration macro <code>BOOST_UBLAS_NO_ELEMENT_PROXIES</code>.
208</p>
209
210<hr />
211<p>Copyright (&copy;) 2000-2004 Joerg Walter, Mathias Koch, Gunter
212Winkler, Michael Stevens<br />
213Permission to copy, use, modify, sell and distribute this document
214is granted provided this copyright notice appears in all copies.
215This document is provided ``as is'' without express or implied
216warranty, and with no claim as to its suitability for any
217purpose.</p>
218</body>
219</html>
Note: See TracBrowser for help on using the repository browser.