-
Notifications
You must be signed in to change notification settings - Fork 0
/
RootFindingTests.cpp
63 lines (44 loc) · 2.24 KB
/
RootFindingTests.cpp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
#define BOOST_TEST_MODULE gk_dispersion_tests
#include <boost/test/included/unit_test.hpp>
#include <boost/mpl/list.hpp>
#include <boost/math/constants/constants.hpp>
#include "RootFinder.h"
BOOST_AUTO_TEST_SUITE( root_finding_test_suite, * boost::unit_test::tolerance( 1e-14 ) )
using namespace RootFinder;
BOOST_AUTO_TEST_CASE( winding_tests )
{
Real M_2PI = boost::math::double_constants::two_pi;
Path Circle = [=]( Real x ) { return std::complex<Real>( std::sin( M_2PI * x ), std::cos( M_2PI * x ) );};
Path DoubleCircle = [=]( Real x ) { return std::complex<Real>( std::sin( M_2PI * x * 2 ), std::cos( M_2PI * x * 2) );};
Path Squircle = [=]( Real x ) { return ( ( x - .5 )*( x-.5 ) + 0.5 )*std::complex<Real>( std::cos( M_2PI * x * 2 ), std::sin( M_2PI * x * 2) );};
BOOST_TEST( WindingNumber( Circle ) == -1.0 );
BOOST_TEST( WindingNumber( DoubleCircle ) == -2.0 );
BOOST_TEST( WindingNumber( Squircle ) == 2.0 );
Complex a( -1,-1 ),b( 1,1 );
BOOST_TEST( WindingNumber( Rectangle( a, b ) ) == 1.0 );
BOOST_TEST( WindingNumber( Rectangle( b, a ) ) == 1.0 );
Complex c( 2,2 );
BOOST_TEST( WindingNumber( Rectangle( b, c ) ) == 0.0 );
}
BOOST_AUTO_TEST_CASE( polynomial_root_tests, * boost::unit_test::tolerance( 1e-14 ) )
{
// z^2 - 4i == 0, with roots +/- sqrt(2)*(1+i)
Func P = []( Complex z ){ return z*z - std::complex<Real>( 0.0, 4.0 );};
Complex a( 0.0, 0.0 ), b( 1.0, 1.0 ), c( 2.0, 2.0 ),d( -2.0, -2.0 );
BOOST_TEST( WindingNumber( RectangleImage( a, b, P ) ) == 0 );
BOOST_TEST( WindingNumber( RectangleImage( b, c, P ) ) == 1 );
BOOST_TEST( WindingNumber( RectangleImage( d, c, P ) ) == 2 );
std::list<Complex> roots = FindWithin< std::list<Complex> >( RootBoundingBox( b, c, 1 ), P, 1e-13 );
Complex exact = ::sqrt( 2.0 )*std::complex<double>( 1.0, 1.0 );
BOOST_TEST( roots.size() == 1 );
BOOST_TEST( std::abs( roots.front() - exact ) == 0.0 );
// And with higher-order roots
auto Q = [=]( Complex z, unsigned int n ){ return pow( z-exact, n );};
for ( unsigned int n=3; n < 10; n++ )
{
roots = FindWithin< std::list<Complex> >( RootBoundingBox( b, c, 1 ), std::bind( Q, std::placeholders::_1, n ), 1e-15 );
BOOST_TEST( roots.size() == 1 );
BOOST_TEST( std::abs( roots.front() - exact ) == 0.0 );
}
}
BOOST_AUTO_TEST_SUITE_END()