/usr/share/doc/libblitz-doc/html/User-et.html

<!DOCTYPE html PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN" "http://www.w3.org/TR/html4/loose.dtd">
<html>
<!-- Created by GNU Texinfo 6.0, http://www.gnu.org/software/texinfo/ -->
<head>
<title>Blitz++: User et</title>

<meta name="description" content="Blitz++: User et">
<meta name="keywords" content="Blitz++: User et">
<meta name="resource-type" content="document">
<meta name="distribution" content="global">
<meta name="Generator" content="makeinfo">
<meta http-equiv="Content-Type" content="text/html; charset=utf-8">
<link href="index.html#Top" rel="start" title="Top">
<link href="Keyword-Index.html#Keyword-Index" rel="index" title="Keyword Index">
<link href="index.html#SEC_Contents" rel="contents" title="Table of Contents">
<link href="Array-Expressions.html#Array-Expressions" rel="up" title="Array Expressions">
<link href="Where-expr.html#Where-expr" rel="next" title="Where expr">
<link href="Math-functions-2.html#Math-functions-2" rel="prev" title="Math functions 2">
<style type="text/css">
<!--
a.summary-letter {text-decoration: none}
blockquote.indentedblock {margin-right: 0em}
blockquote.smallindentedblock {margin-right: 0em; font-size: smaller}
blockquote.smallquotation {font-size: smaller}
div.display {margin-left: 3.2em}
div.example {margin-left: 3.2em}
div.lisp {margin-left: 3.2em}
div.smalldisplay {margin-left: 3.2em}
div.smallexample {margin-left: 3.2em}
div.smalllisp {margin-left: 3.2em}
kbd {font-style: oblique}
pre.display {font-family: inherit}
pre.format {font-family: inherit}
pre.menu-comment {font-family: serif}
pre.menu-preformatted {font-family: serif}
pre.smalldisplay {font-family: inherit; font-size: smaller}
pre.smallexample {font-size: smaller}
pre.smallformat {font-family: inherit; font-size: smaller}
pre.smalllisp {font-size: smaller}
span.nocodebreak {white-space: nowrap}
span.nolinebreak {white-space: nowrap}
span.roman {font-family: serif; font-weight: normal}
span.sansserif {font-family: sans-serif; font-weight: normal}
ul.no-bullet {list-style: none}
-->
</style>


</head>

<body lang="en">
<a name="User-et"></a>
<div class="header">
<p>
Next: <a href="Where-expr.html#Where-expr" accesskey="n" rel="next">Where expr</a>, Previous: <a href="Math-functions-2.html#Math-functions-2" accesskey="p" rel="prev">Math functions 2</a>, Up: <a href="Array-Expressions.html#Array-Expressions" accesskey="u" rel="up">Array Expressions</a> &nbsp; [<a href="index.html#SEC_Contents" title="Table of contents" rel="contents">Contents</a>][<a href="Keyword-Index.html#Keyword-Index" title="Index" rel="index">Index</a>]</p>
</div>
<hr>
<a name="Declaring-your-own-math-functions-on-arrays"></a>
<h3 class="section">3.10 Declaring your own math functions on arrays</h3>

<a name="index-math-functions-declaring-your-own"></a>
<a name="index-Array-declaring-your-own-math-functions-on"></a>

<p>There are four macros which make it easy to turn your own scalar functions
into functions defined on arrays.  They are:
</p>
<a name="index-BZ_005fDECLARE_005fFUNCTION"></a>

<div class="example">
<pre class="example">BZ_DECLARE_FUNCTION(f)                   // 1
BZ_DECLARE_FUNCTION_RET(f,return_type)   // 2
BZ_DECLARE_FUNCTION2(f)                  // 3
BZ_DECLARE_FUNCTION2_RET(f,return_type)  // 4
</pre></div>

<p>Use version 1 when you have a function which takes one argument and returns
a result of the same type.  For example:
</p>
<div class="example">
<pre class="example">#include &lt;blitz/array.h&gt;

using namespace blitz;

double myFunction(double x)
{ 
    return 1.0 / (1 + x); 
}

BZ_DECLARE_FUNCTION(myFunction)

int main()
{
    Array&lt;double,2&gt; A(4,4), B(4,4);  // ...
    B = myFunction(A);
}
</pre></div>

<p>Use version 2 when you have a one argument function whose return type is
different than the argument type, such as
</p>
<div class="example">
<pre class="example">int g(double x);
</pre></div>

<p>Use version 3 for a function which takes two arguments and returns a result
of the same type, such as:
</p>
<div class="example">
<pre class="example">double g(double x, double y);
</pre></div>

<p>Use version 4 for a function of two arguments which returns a different
type, such as:
</p>
<div class="example">
<pre class="example">int g(double x, double y);
</pre></div>

<a name="Tensor-notation"></a>
<h3 class="section">3.11 Tensor notation</h3>

<a name="index-tensor-notation"></a>
<a name="index-Array-tensor-notation"></a>

<p>Blitz++ arrays support a tensor-like notation.  Here&rsquo;s an example of
real-world tensor notation:
<pre>
 ijk    ij k
A    = B  C
</pre></p>
<p><em>A</em> is a rank 3 tensor (a three dimensional array), <em>B</em> is a rank
2 tensor (a two dimensional array), and <em>C</em> is a rank 1 tensor (a one
dimensional array).  The above expression sets 
<code>A(i,j,k) = B(i,j) * C(k)</code>.
</p>
<p>To implement this product using Blitz++, we&rsquo;ll need the arrays and some
index placeholders:
</p>
<a name="index-index-placeholders-used-for-tensor-notation"></a>

<div class="example">
<pre class="example">Array&lt;float,3&gt; A(4,4,4);
Array&lt;float,2&gt; B(4,4);
Array&lt;float,1&gt; C(4);

firstIndex i;    // Alternately, could just say
secondIndex j;   // using namespace blitz::tensor;
thirdIndex k;
</pre></div>

<p>Here&rsquo;s the Blitz++ code which is equivalent to the tensor expression:
</p>
<div class="example">
<pre class="example">A = B(i,j) * C(k);
</pre></div>

<p>The index placeholder arguments tell an array how to map its dimensions onto
the dimensions of the destination array.  For example, here&rsquo;s some
real-world tensor notation:
<pre>
 ijk    ij k    jk i
C    = A  x  - A  y
</pre></p>
<p>In Blitz++, this would be coded as:
</p>
<div class="example">
<pre class="example">using namespace blitz::tensor;

C = A(i,j) * x(k) - A(j,k) * y(i);
</pre></div>

<p>This tensor expression can be visualized in the following way:
</p>
<div align="center"><img src="tensor1.gif" alt="tensor1">
</div><div align="center">Examples of array indexing, subarrays, and slicing.
</div>
<p>Here&rsquo;s an example which computes an outer product of two one-dimensional
arrays:
<a name="index-outer-product"></a>
<a name="index-kronecker-product"></a>
<a name="index-tensor-product"></a>
</p>
<div class="smallexample">
<pre class="smallexample">#include &lt;blitz/array.h&gt;

using namespace blitz;

int main()
{
    Array&lt;float,1&gt; x(4), y(4);
    Array&lt;float,2&gt; A(4,4);

    x = 1, 2, 3, 4;
    y = 1, 0, 0, 1;

    firstIndex i;
    secondIndex j;

    A = x(i) * y(j);

    cout &lt;&lt; A &lt;&lt; endl;

    return 0;
}

</pre></div>

<p>And the output:
</p>
<div class="smallexample">
<pre class="smallexample">(0,3) x (0,3)
[ 1 0 0 1 
  2 0 0 2 
  3 0 0 3 
  4 0 0 4 ]

</pre></div>

<p>Index placeholders can <em>not</em> be used on the left-hand side of an
expression.  If you need to reorder the indices, you must do this on the
right-hand side.
</p>
<p>In real-world tensor notation, repeated indices imply a contraction (or
summation).  For example, this tensor expression computes a matrix-matrix
product:
<pre>
 ij    ik  kj
C   = A   B
</pre></p>
<p>The repeated k index is interpreted as meaning
<pre>
c    = sum of (a   * b  ) over k
 ij             ik    kj
</pre></p>
<a name="index-contraction"></a>
<a name="index-tensor-contraction"></a>

<p>In Blitz++, repeated indices do <em>not</em> imply contraction.  If you want
to contract (sum along) an index, you must use the <code>sum()</code> function:
</p>
<div class="example">
<pre class="example">Array&lt;float,2&gt; A, B, C;   // ...
firstIndex i;
secondIndex j;
thirdIndex k;

C = sum(A(i,k) * B(k,j), k);
</pre></div>

<p>The <code>sum()</code> function is an example of an <em>array reduction</em>,
described in the next section.
</p>
<p>Index placeholders can be used in any order in an expression.  This example
computes a kronecker product of a pair of two-dimensional arrays, and
permutes the indices along the way:
</p>
<div class="example">
<pre class="example">Array&lt;float,2&gt; A, B;   // ...
Array&lt;float,4&gt; C;      // ...
fourthIndex l;

C = A(l,j) * B(k,i);
</pre></div>

<p>This is equivalent to the tensor notation
<pre>
 ijkl    lj ki
C     = A  B
 </pre></p>
<p>Tensor-like notation can be mixed with other array notations:
</p>
<div class="example">
<pre class="example">Array&lt;float,2&gt; A, B;  // ...
Array&lt;double,4&gt; C;    // ...

C = cos(A(l,j)) * sin(B(k,i)) + 1./(i+j+k+l);
</pre></div>

<a name="index-tensor-notation-efficiency-issues"></a>
<p>An important efficiency note about tensor-like notation: the right-hand side
of an expression is <em>completely evaluated</em> for <em>every</em> element in
the destination array.  For example, in this code:
</p>
<div class="example">
<pre class="example">Array&lt;float,1&gt; x(4), y(4);
Array&lt;float,2&gt; A(4,4):

A = cos(x(i)) * sin(y(j));
</pre></div>

<p>The resulting implementation will look something like this:
</p>
<div class="example">
<pre class="example">for (int n=0; n &lt; 4; ++n)
  for (int m=0; m &lt; 4; ++m)
    A(n,m) = cos(x(n)) * sin(y(m));
</pre></div>

<p>The functions <code>cos</code> and <code>sin</code> will be invoked sixteen times each.
It&rsquo;s possible that a good optimizing compiler could hoist the <code>cos</code>
evaluation out of the inner loop, but don&rsquo;t hold your breath &ndash; there&rsquo;s a
lot of complicated machinery behind the scenes to handle tensor notation,
and most optimizing compilers are easily confused.  In a situation like the
above, you are probably best off manually creating temporaries for
<code>cos(x)</code> and <code>sin(y)</code> first.
</p>
<a name="Array-reductions"></a>
<h3 class="section">3.12 Array reductions</h3>
<a name="index-Array-reductions"></a>
<a name="index-reductions"></a>

<p>Currently, Blitz++ arrays support two forms of reduction:
</p>
<ul>
<li> Reductions which transform an array into a scalar (for example,
summing the elements).  These are referred to as <strong>complete
reductions</strong>.

</li><li> Reducing an N dimensional array (or array expression) to an N-1
dimensional array expression.  These are called <strong>partial reductions</strong>.

</li></ul>

<a name="index-Array-reductions-complete"></a>
<a name="index-complete-reductions"></a>
<a name="index-reductions-complete"></a>

<a name="Complete-reductions"></a>
<h3 class="section">3.13 Complete reductions</h3>

<p>Complete reductions transform an array (or array expression) into 
a scalar.  Here are some examples:
</p>
<div class="example">
<pre class="example">Array&lt;float,2&gt; A(3,3);
A = 0, 1, 2,
    3, 4, 5,
    6, 7, 8;
cout &lt;&lt; sum(A) &lt;&lt; endl          // 36
     &lt;&lt; min(A) &lt;&lt; endl          // 0
     &lt;&lt; count(A &gt;= 4) &lt;&lt; endl;  // 5
</pre></div>

<p>Here are the available complete reductions:
</p>
<dl compact="compact">
<dt><code>sum()</code></dt>
<dd><a name="index-sum_0028_0029-reduction"></a>
<p>Summation (may be promoted to a higher-precision type)
</p>
</dd>
<dt><code>product()</code></dt>
<dd><a name="index-product_0028_0029-reduction"></a>
<p>Product 
</p>
</dd>
<dt><code>mean()</code></dt>
<dd><a name="index-mean_0028_0029-reduction"></a>
<p>Arithmetic mean (promoted to floating-point type if necessary) 
</p>
</dd>
<dt><code>min()</code></dt>
<dd><a name="index-min_0028_0029-reduction"></a>
<p>Minimum value 
</p>
</dd>
<dt><code>max()</code></dt>
<dd><a name="index-max_0028_0029-reduction"></a>
<p>Maximum value 
</p>
</dd>
<dt><code>minmax()</code></dt>
<dd><a name="index-minmax_0028_0029-reduction"></a>
<p>Simultaneous minimum and maximum value (returns a value of type MinMaxValue&lt;T&gt;)
</p>
</dd>
<dt><code>minIndex()</code></dt>
<dd><a name="index-minIndex_0028_0029-reduction"></a>
<p>Index of the minimum value (<code>TinyVector&lt;int,N_rank&gt;</code>)
</p>
</dd>
<dt><code>maxIndex()</code></dt>
<dd><a name="index-maxIndex_0028_0029-reduction"></a>
<p>Index of the maximum value (<code>TinyVector&lt;int,N_rank&gt;</code>)
</p>
</dd>
<dt><code>count()</code></dt>
<dd><a name="index-count_0028_0029-reduction"></a>
<p>Counts the number of times the expression is logical true (<code>int</code>)
</p>
</dd>
<dt><code>any()</code></dt>
<dd><a name="index-any_0028_0029-reduction"></a>
<p>True if the expression is true anywhere (<code>bool</code>)
</p>
</dd>
<dt><code>all()</code></dt>
<dd><a name="index-all_0028_0029-reduction"></a>
<p>True if the expression is true everywhere (<code>bool</code>)
</p></dd>
</dl>

<p><strong>Caution:</strong> <code>minIndex()</code> and <code>maxIndex()</code> return TinyVectors, 
even when the rank of the array (or array expression) is 1.
</p>
<p>Reductions can be combined with <code>where</code> expressions (<a href="Where-expr.html#Where-expr">Where expr</a>)
to reduce over some part of an array.  For example, <code>sum(where(A &gt; 0,
A, 0))</code> sums only the positive elements in an array.
</p>
<a name="Partial-Reductions"></a>
<h3 class="section">3.14 Partial Reductions</h3>

<a name="index-Array-reductions-partial"></a>
<a name="index-partial-reductions"></a>
<a name="index-reductions-partial"></a>

<p>Here&rsquo;s an example which computes the sum of each row of a two-dimensional
array:
</p>
<div class="example">
<pre class="example">Array&lt;float,2&gt; A;    // ...
Array&lt;float,1&gt; rs;   // ...
firstIndex i;
secondIndex j;

rs = sum(A, j);
</pre></div>

<p>The reduction <code>sum()</code> takes two arguments:
</p>
<ul>
<li> The first argument is an array or array expression.

</li><li> The second argument is an index placeholder indicating the
dimension over which the reduction is to occur.  

</li></ul>

<p>Reductions have an <strong>important restriction</strong>: It is currently only
possible to reduce over the <em>last</em> dimension of an array or array
expression.  Reducing a dimension other than the last would require Blitz++
to reorder the dimensions to fill the hole left behind.  For example, in
order for this reduction to work:
</p>
<div class="example">
<pre class="example">Array&lt;float,3&gt; A;   // ...
Array&lt;float,2&gt; B;   // ...
secondIndex j;

// Reduce over dimension 2 of a 3-D array?
B = sum(A, j);
</pre></div>

<p>Blitz++ would have to remap the dimensions so that the third dimension
became the second.  It&rsquo;s not currently smart enough to do this.
</p>
<p>However, there is a simple workaround which solves some of the problems
created by this limitation: you can do the reordering manually, prior to the
reduction:
</p>
<div class="example">
<pre class="example">B = sum(A(i,k,j), k);
</pre></div>

<p>Writing <code>A(i,k,j)</code> interchanges the second and third dimensions,
permitting you to reduce over the second dimension.  Here&rsquo;s a list of the
reduction operations currently supported:
</p>
<dl compact="compact">
<dt><code>sum()</code></dt>
<dd><p>Summation
</p>
</dd>
<dt><code>product()</code></dt>
<dd><p>Product 
</p>
</dd>
<dt><code>mean()</code></dt>
<dd><p>Arithmetic mean (promoted to floating-point type if necessary)
</p>
</dd>
<dt><code>min()</code></dt>
<dd><p>Minimum value
</p>
</dd>
<dt><code>max()</code></dt>
<dd><p>Maximum value
</p>
</dd>
<dt><code>minIndex()</code></dt>
<dd><p>Index of the minimum value (int)
</p>
</dd>
<dt><code>maxIndex()</code></dt>
<dd><p>Index of the maximum value (int)
</p>
</dd>
<dt><code>count()</code></dt>
<dd><p>Counts the number of times the expression is logical true (int)
</p>
</dd>
<dt><code>any()</code></dt>
<dd><p>True if the expression is true anywhere (bool)
</p>
</dd>
<dt><code>all()</code></dt>
<dd><p>True if the expression is true everywhere (bool)
</p>
</dd>
<dt><code>first()</code></dt>
<dd><p>First index at which the expression is logical true (int); if the expression
is logical true nowhere, then <code>tiny(int())</code> (INT_MIN) is returned.
</p>
</dd>
<dt><code>last()</code></dt>
<dd><p>Last index at which the expression is logical true (int); if the expression
is logical true nowhere, then <code>huge(int())</code> (INT_MAX) is returned.  
</p></dd>
</dl>

<p>The reductions <code>any()</code>, <code>all()</code>, and <code>first()</code> have
short-circuit semantics: the reduction will halt as soon as the answer is
known.  For example, if you use <code>any()</code>, scanning of the expression
will stop as soon as the first true value is encountered.
</p>
<p>To illustrate, here&rsquo;s an example:
</p>
<div class="example">
<pre class="example">Array&lt;int, 2&gt; A(4,4);

A =  3,   8,   0,   1,
     1,  -1,   9,   3,
     2,  -5,  -1,   1,
     4,   3,   4,   2;

Array&lt;float, 1&gt; z(4);
firstIndex i;
secondIndex j;

z = sum(A(j,i), j);
</pre></div>

<p>The array <code>z</code> now contains the sum of <code>A</code> along each column:
</p>
<div class="example">
<pre class="example">[ 10    5     12    7 ]
</pre></div>

<p>This table shows what the result stored in <code>z</code> would be if
<code>sum()</code> were replaced with other reductions:
</p>
<div class="example">
<pre class="example">sum                     [         10         5        12         7 ]
mean                    [        2.5      1.25         3      1.75 ]
min                     [          1        -5        -1         1 ]
minIndex                [          1         2         2         0 ]
max                     [          4         8         9         3 ]
maxIndex                [          3         0         1         1 ]
first((A &lt; 0), j)       [ -2147483648        1         2 -2147483648 ]
product                 [         24       120         0         6 ]
count((A(j,i) &gt; 0), j)  [          4         2         2         4 ]
any(abs(A(j,i)) &gt; 4, j) [          0         1         1         0 ]
all(A(j,i) &gt; 0, j)      [          1         0         0         1 ]
</pre></div>

<p>Note: the odd numbers for first() are <code>tiny(int())</code> i.e. the smallest
number representable by an int.  The exact value is machine-dependent.
</p>
<a name="index-Array-reductions-chaining"></a>
<a name="index-partial-reductions-chaining"></a>
<a name="index-reductions-chaining"></a>

<p>The result of a reduction is an array expression, so reductions
can be used as operands in an array expression:
</p>
<div class="example">
<pre class="example">Array&lt;int,3&gt; A;
Array&lt;int,2&gt; B;
Array&lt;int,1&gt; C;   // ...

secondIndex j;
thirdIndex k;

B = sqrt(sum(sqr(A), k));

// Do two reductions in a row
C = sum(sum(A, k), j);
</pre></div>

<p>Note that this is not allowed:
</p>
<div class="example">
<pre class="example">Array&lt;int,2&gt; A;
firstIndex i;
secondIndex j;

// Completely sum the array?
int result = sum(sum(A, j), i);
</pre></div>

<p>You cannot reduce an array to zero dimensions!  Instead, use one of the
global functions described in the previous section.
</p>

<hr>
<div class="header">
<p>
Next: <a href="Where-expr.html#Where-expr" accesskey="n" rel="next">Where expr</a>, Previous: <a href="Math-functions-2.html#Math-functions-2" accesskey="p" rel="prev">Math functions 2</a>, Up: <a href="Array-Expressions.html#Array-Expressions" accesskey="u" rel="up">Array Expressions</a> &nbsp; [<a href="index.html#SEC_Contents" title="Table of contents" rel="contents">Contents</a>][<a href="Keyword-Index.html#Keyword-Index" title="Index" rel="index">Index</a>]</p>
</div>



</body>
</html>
libblitz-doc 1:0.10-3.3 / usr / share / doc / libblitz-doc / html / User-et.html
