← Back to branch summary

~ubuntu-branches/ubuntu/utopic/blitz++/utopic

« back to all changes in this revision

Viewing changes to examples/erf.cpp

  • Committer: Package Import Robot
  • Author(s): Christophe Trophime
  • Date: 2012-07-06 09:15:30 UTC
  • mfrom: (11.1.5 sid)
  • Revision ID: package-import@ubuntu.com-20120706091530-vzrb8zf0vpbf8tp9
Tags: 1:0.10-1
http://bugs.debian.org/679407
* New upstream release
  Closes: #679407
* debian/rules:
  - update for new release
  - add override_dh_auto_test target
  - regenerate configure and Makefile.am
* debian/control:
  - add libtool, automake to BuildDepends
* debian/libblitz-doc.install
  - modify path for html files
* remove uneeded patches
* add examples.patch

Show diffs side-by-side

added added

removed removed

Lines of Context:
examples/erf.cpp
1
 
/*****************************************************************************
2
 
 * erf.cpp        Blitz++ Vector<T> example
3
 
 *****************************************************************************/
4
 
 
5
 
#include <blitz/blitz.h>
6
 
 
7
 
#ifdef BZ_HAVE_IEEE_MATH
8
 
 
9
 
#include <blitz/vector-et.h>
10
 
 
11
 
// This program uses erf(), which is not part of ANSI C/C++.  In order
12
 
// to compile this example, you must have a Posix or X/Open-compliant 
13
 
// compiler, and #define an appropriate source level (-D_XOPEN_SOURCE or
14
 
// -D_POSIX_SOURCE).  You will need to run compiler/bzconfig with
15
 
// this flag defined as well, so that the Blitz++ flag BZ_HAVE_IEEE_MATH
16
 
// will be set properly in <blitz/config.h>.
17
 
 
18
 
BZ_USING_NAMESPACE(blitz)
19
 
 
20
 
int main()
21
 
{
22
 
    // Integrate the expression
23
 
    //            x
24
 
    //   2        /
25
 
    // -------    | exp(-t^2) dt
26
 
    // sqrt(Pi)   /
27
 
    //           0
28
 
    // to estimate the error function erf(x) on the interval [0,1].
29
 
    //
30
 
    // Three methods are compared:
31
 
    //  1. Naive summation
32
 
    //  2. Two-point trapezoid
33
 
    //  3. Three-point Simpson's
34
 
 
35
 
    const double m_2_sqrtpi = 1.12837916709551257389615890312154517;
36
 
    const int numSamples = 1024;       // Number of samples
37
 
    double dt = 1.0 / numSamples;      // Distance between samples
38
 
    Range R(0, numSamples - 1);        // Index set for sample vectors
39
 
    Range displayRange(0, 900, 100);   // Indices to output to screen
40
 
 
41
 
    // For comparison purposes, use the built-in erf(x) function.
42
 
    Vector<double> z = erf(R * dt);
43
 
    cout << "  Built-in: " << z(displayRange+1) << endl;
44
 
 
45
 
 
46
 
    // Naive summation
47
 
    Vector<double> z1 = accumulate(m_2_sqrtpi * exp(-sqr(R * dt)) * dt);
48
 
    cout << "     Naive: " << z1(displayRange+1) << endl;
49
 
 
50
 
    // Calculate the rms error
51
 
    double error1 = sqrt(mean(sqr(z1 - z)));
52
 
    cout << "RMS Error: " << error1 << endl;
53
 
 
54
 
 
55
 
    // For the following rules, it's easier to work with a
56
 
    // sampled integrand.
57
 
    Vector<double> a = m_2_sqrtpi * exp(-sqr(R * dt));
58
 
 
59
 
 
60
 
    // Two-point trapezoid
61
 
    Range J(1, numSamples-1);
62
 
    Vector<double> z2 = accumulate(dt / 2.0 * (a(J) + a(J-1)));
63
 
 
64
 
    cout << " Trapezoid: " << z2(displayRange) << endl;
65
 
    cout << "RMS Error: " << sqrt(mean(sqr(z2-z(J)))) << endl;
66
 
 
67
 
 
68
 
    // Three-point Simpson's (parabolic)
69
 
    Range I(1, numSamples-2);
70
 
    Vector<double> z3 = accumulate(dt / 6.0 * (a(I-1) + 4 * a(I) + a(I+1)));
71
 
 
72
 
    cout << " Parabolic: " << z3(displayRange) << endl;
73
 
    cout << "RMS Error: " << sqrt(mean(sqr(z3 - z(I)))) << endl;
74
 
 
75
 
 
76
 
    return 0;
77
 
}
78
 
 
79
 
#else
80
 
 
81
 
#include <iostream.h>
82
 
 
83
 
int main()
84
 
{
85
 
    cout << "This example requires IEEE math functions." << endl;
86
 
    return 0;
87
 
}
88
 
 
89
 
#endif  // BZ_HAVE_IEEE_MATH