Category Archives: C++ Code

C++17: Lambdas

C++ has supported lambda functions since C++11, and post-C++14 even supported generic lambdas. C++17 adds a quirky feature to enable ‘copy capture’ of ‘this’. Here’s an example of it in action:

struct MyStruct
{
    auto lambda_with_this() // Before C++17, copies this pointer
    {
        auto f = [this]{ return value; };
        return f;
    }
    
    auto lambda_with_star_this() // C++17 - takes a local copy of *this
    {
        auto f = [*this]{ return value; };
        return f;
    }
    
    int value;
};

TEST( Cpp17, lambda_capture_this )
{
    MyStruct s{ 42 };
    auto f = s.lambda_with_this();
    s.value = 101;
    EXPECT_EQ( 101, f() );
}

TEST( Cpp17, lambda_capture_star_this )
{
    MyStruct s{ 42 };
    auto f = s.lambda_with_star_this();
    s.value = 101;
    EXPECT_EQ( 42, f() );
}

Notice that in the second case, we capture a copy of our object – so the lambda returns the value held at the point of capture (42) rather than the value when it is called (101). This can be very important if ‘this’ has been destroyed between the creation of the lambda and the time at which it’s called.

Now, C++14 also supported generalised lambda capture, which meant you could (re-)name variables when capturing (and provided a neat way to capture-by-move):

    auto f = [tmp = *this]{ return tmp.value; };

But the C++17 code is more concise. See this useful StackOverflow post for more discussion.

Another advance in C++17 is that lambdas are implicitly constexpr – so you can now use them in compile-time contexts, like declaration of std::array:

// lambda explicitly constexpr since C++17
auto square = []( auto v ){ return v*v; }; 

TEST( Cpp17, lambda_implicitly_constexpr )
{
    // std::array calls 'square()' at compile time
    std::array<int, square(4)> values; 
    EXPECT_EQ( 16, values.size() );
}

See also the previous C++17 post.

Leave a comment

Filed under C++, C++ Code, Programming

C++17: Nested Namespaces

Another crowd-pleaser in C++17 is the ability to declared nested namespaces without literally nesting them. In the past, you had to do this, which involves a lot of wasted whitespace:

namespace NestedNamespaces
{
    namespace Really
    {
        namespace Work
        {
            auto return_int(){ return 42; };
        }
    }
}

Happily, you can now do this instead:

namespace NestedNamespaces::Really::Work
{
    auto return_int(){ return 42; };
}

TEST( Cpp17, nested_namespaces )
{
    EXPECT_EQ( 42, NestedNamespaces::Really::Work::return_int() );
}

See also previous C++17 post.

2 Comments

Filed under C++, C++ Code, Programming

C++17: Structured Bindings

This is one of a series of posts on C++17 features – see also previous post on if initialisers.

Structured bindings are a convenient way of handling multiple return values from functions. Whilst F# has been able to do this:

let f() = 42, "Hello, World" // return a pair of values
let a, b = f() // assign a and b to the values returned by f

in C++, we’ve had to declare the variables up front and use std::tie to assign values (so not only does this take more lines, we also have to default initialise the variables then throw away the defaults).

auto t = std::make_tuple( 42, "Hello, World" );
int a, b;
std::tie( a, b ) = t;

The new structured bindings are much more concise, even if the use of square brackets came as a surprise. Even better, you can use structured bindings with structs and std::array.

int my_int{ 42 };
std::string my_string{ "Hello, World" };
bool my_bool{ true };

auto return_pair()
{
    return std::make_pair( my_int, my_string );
}

auto return_tuple()
{
    return std::make_tuple( my_int, my_string, my_bool );
}

struct MyStruct
{
    int a;
    double b;
    int c;
    
    static MyStruct Expected;
};

MyStruct MyStruct::Expected = { 1, 2.2, 3 };

auto return_struct()
{
    return MyStruct::Expected;
}

auto return_array()
{
    return std::array<int,3>{ 1, 2, 3 };
}

auto return_map()
{
    return std::map<int, std::string>{ {1, "a"}, {2, "b"}, {3, "c"} };
}

TEST( Cpp17, structured_bindings_for_pair )
{
    auto [i, s] = return_pair();
    
    EXPECT_EQ( my_int, i );
    EXPECT_EQ( my_string, s );
}

TEST( Cpp17, structured_bindings_for_tuple )
{
    auto [i, s, b] = return_tuple();
    
    EXPECT_EQ( my_int, i );
    EXPECT_EQ( my_string, s );
    EXPECT_EQ( my_bool, b );
}

TEST( Cpp17, structured_bindings_for_struct )
{
    auto [i1, d, i2] = return_struct();
    
    EXPECT_EQ( MyStruct::Expected.a, i1 );
    EXPECT_EQ( MyStruct::Expected.b, d );
    EXPECT_EQ( MyStruct::Expected.c, i2 );
}

TEST( Cpp17, structured_bindings_for_array )
{
    auto [i1, i2, i3] = return_array();
    
    EXPECT_EQ( 1, i1 );
    EXPECT_EQ( 2, i2 );
    EXPECT_EQ( 3, i3 );
}

TEST( Cpp17, structured_bindings_for_iterating_over_map )
{
    for ( const auto& [key,value] : return_map() )
    {
        switch (key)
        {
            case 1: EXPECT_EQ( "a", value ); break;
            case 2: EXPECT_EQ( "b", value ); break;
            case 3: EXPECT_EQ( "c", value ); break;
            default: break;            
        };
    }
}

For me, the best examples come when combining features – the range-based for loop with structured bindings is a thing of beauty.

See also next C++17 post.

2 Comments

Filed under C++, C++ Code, Programming

C++17: if initialiser

I attended a C++17 presentation by Nicolai Josuttis last year, but at the time, my laptop’s compiler didn’t support any of the features to try them out. After a recent update, it turns out that many are now supported, so I’ve written a few unit tests using GTest.

The first feature I tried was the if initialiser. This feature looks a bit odd at first, because C++ programmers are so conditioned to seeing if statements containing a single condition. Allowing an initialiser statement as well

if ( initialiser; condition )

means that you can initialise a variable and test it on the same line. It also prevents the variable being used outside the scope of the if statement – this prevents accidental re-use if you subsequently mis-type a variable name.

auto return_int()
{
   return 101;
}

TEST( Cpp17, if_initialiser )
{
    // NB we can use i in the body of the if or the else
    // Also, must have a variable name for the object to live in the whole statement
    // (so must name locks taken, even if not used in the body, otherwise it's a temporary).
    if ( auto i = return_int(); i < 100 )
    {
        EXPECT_TRUE( i < 100 );
    }
    else
    {
        EXPECT_TRUE( i >= 100 );
    }
}

TEST( Cpp17, if_initialiser_with_map_insert)
{
    std::map<int, std::string> my_map{ {42, "Hi" } };
    
    if ( auto[it, inserted] = my_map.insert( std::make_pair(42, "Bye" ) ); !inserted )
    {
        // See also StructuredBindings for iterating over a map
        auto& [key,value] = *it; // iterator is pair of key and value
        EXPECT_EQ( 42, key );
        EXPECT_EQ( "Hi", value );
        
        value = "Bye"; // update the value
        EXPECT_EQ( "Bye", my_map[42] );
        
        // key = 43; // compile error! - cannot assign to const key-type
    }
}

See also next post on C++17 Structured Bindings.

1 Comment

Filed under C++, C++ Code, Programming

Higher-Order functions in C++

The other day, I was writing some C++ and found that I was thinking about how to manipulate the data I had as if I was writing F#. It would have been convenient to turn a std::map into an array of tuples, which I could do in F# like this:

let f (xs : map<int,string) =
  let xs = xs |> Map.toArray
  // now treat xs as array of tuples...

There’s no function in STL to do this off the bat – instead, you have to roll your own (not that it’s much code, but it does break your through processes if you have to stop to code this sort of thing every time).

Of course, this is just one of many helpful F# higher-order functions that are provided in the F# Map module – and there are counterparts for each of the collection classes i.e. Array, Set, List etc. In C++, the nearest equivalent is the STL which provides both collection classes and a number of algorithms that operate on them. Better still, from C++11 onwards we have lambdas, which make using STL algorithms much easier. Even so, in most cases, the F# operations seem much more tailored to the sort of data transformation I see at work – our codebase is littered with map/filter/fold operations as people transform/select and accumulate data. Conversely, our C++ codebase is full of … for loops, evidence in my eyes that STL algorithms aren’t as immediately applicable. In fact, the ease of use of higher-order functions was one of the reasons that F# was quickly adopted in my workplace (along with immutability, strong-typing, conciseness, type inference and syntax checking).

I’ve written one-to-one C++ equivalents of the F# module functions that I use the most for Map and Vector – see below. Interestingly, I found that I really did have to ‘engage brain’ to write some of these, particularly Map.filter. For that one, you can’t use the erase-remove idiom because map keys are both const and strictly ordered (whereas for Vector, erase-remove_if implements filter neatly). A library of functions as per my code below would definitely be a productivity boost.

First, I’ve factored some common utilities into namespace Collection:

namespace MusingStudio
{
    namespace Collection
    {
        template <typename C, typename F>
        C& filter( C& collection, F keep_predicate )
        {
            auto erase_predicate = [&pred=keep_predicate]( auto&& x ){ return !pred( std::forward<decltype(x)>(x) ); };
            collection.erase( std::remove_if( collection.begin(), collection.end(), erase_predicate ), collection.end() );
            return collection;
        }
        
        // This form of filter always takes a copy and applies the filter to it
        // - sometimes you want to preserve the original collection
        template <typename C, typename F>
        C filter_copy( const C& collection, F keep_predicate )
        {
            C target;
            std::copy_if( collection.begin(), collection.end(), std::inserter( target, target.end() ), keep_predicate );
            return target;
        }
        
        template <typename C, typename F, typename A>
        A fold( const C& items, F f, A&& init )
        {
            A acc{ std::forward<A>(init) };
            for ( const auto& item : items )
            {
                f( acc, item );
            }
            return acc;
        }
                
        // F( T ) -> T and collection C is mutated
        template <typename C, typename F>
        C& transform( C& items, F f )
        {
            for ( auto& t : items )
            {
                t = f( t );
            }
            
            return items;
        }

    }
}

Next, here are the higher-order functions that I use for Map:

namespace MusingStudio
{
    namespace Map
    {
        // filter_copy takes a copy of the original collection then applies the filter
        template< typename K, typename V, typename F>
        std::map<K,V> filter_copy( const std::map<K,V>& items, F predicate )
        {
            return Collection::filter_copy( items, predicate );
        }
        
        template< typename K, typename V, typename F>
        std::map<K,V>& filter( std::map<K,V>& items, F predicate )
        {
            // NB the erase-remove_if idiom does not work for std::map
            // because the nodes must remain ordered by key.  This is enforced
            // by std::map<K,V> holding keys as const K.  So any assignment
            // to the key (to effecfively re-order the binary tree) fails to compile.
            // http://stackoverflow.com/questions/9515357/map-lambda-remove-if
            
            // instead, manually iterate over the collection, erasing items
            // for which predicate() returns false
            for ( auto it = items.begin(), itEnd = items.end(); it != itEnd; )
            {
                if ( predicate( *it ) )
                {
                    ++it; // ok - keep this item
                }
                else
                {
                    it = items.erase( it );
                }
            }
            
            return items;
        }
        
        template <typename K, typename V>
        auto to_vector( const std::map<K,V>& collection )
        {
            std::vector< std::pair<K,V> > items;
            
            for ( const auto& item : collection )
            {
                items.push_back( std::make_pair( item.first, item.second ) );
            }
            
            return items;
        }
        
        template <typename K, typename V>
        auto keys( const std::map<K,V>& collection )
        {
            std::set< K > items;
            
            for ( const auto& item : collection )
            {
                items.insert( item.first );
            }
            
            return items;
        }
        
        template <typename K, typename V>
        auto values( const std::map<K,V>& collection )
        {
            std::vector< V > items;
            
            for ( const auto& item : collection )
            {
                items.push_back( item.second );
            }
            
            return items;
        }
        
        template<typename K, typename V, typename F, typename A>
        A fold( const std::map<K,V>& items, F f, A&& init )
        {
            return Collection::fold( items, f, std::forward<A>(init) );
        }
        
        // F( std::pair<K,V> ) -> std::pair< L, U >
        // Construct a new std::map<L,U> mapping from (K,V) to (L,U)
        template <typename K, typename V, typename F>
        auto map( const std::map<K,V>& items, F f )
        {
            using KVP = typename std::map<K,V>::value_type;
            using RVP = decltype( f( KVP() ) );
            
            std::map< decltype( RVP().first ), decltype( RVP().second ) > result;
            
            for ( const KVP& kvp : items )
            {
                result.insert( f( kvp ) );
            }
            
            return result;
        }
        
        // F( K, V ) -> V and std::map<K,V> is mutated
        template <typename K, typename V, typename F>
        std::map<K,V>& transform( std::map<K,V>& items, F f )
        {
            using KVP = typename std::map<K,V>::value_type;
            
            for ( const KVP& kvp : items )
            {
                items[kvp.first] = f( kvp.first, kvp.second );
            }
            
            return items;
        }
    }
}

And here are the higher-order functions that I use for Vector:

namespace MusingStudio
{
    namespace Vector
    {
        template< typename T, typename F>
        std::vector<T>& filter( std::vector<T>& items, F predicate )
        {
            Collection::filter( items, predicate );
            return items;
        }
        
        template< typename T, typename F>
        std::vector<T> filter_copy( const std::vector<T>& items, F predicate )
        {
            return Collection::filter_copy( items, predicate );
        }
        
        // Requires F to have signature void( A&, T )
        template< typename T, typename F, typename A>
        A fold( const std::vector<T>& items, F f, A&& init )
        {
            return Collection::fold( items, f, std::forward<A>(init) );
        }
        
        template< typename T, typename P = std::less<T> >
        std::vector<T>& sort( std::vector<T>& items, P compare = P() )
        {
            std::sort( items.begin(), items.end(), compare );
            return items;
        }
        
        template< typename T, typename P = std::less<T> >
        std::vector<T> sort_copy( const std::vector<T>& items, P compare = P() )
        {
            std::vector<T> result( items );
            std::sort( result.begin(), result.end(), compare );
            return result;
        }
        
        // F( T ) -> U, construct a new vector<U>, mapping from T to U
        template <typename T, typename F>
        auto map( const std::vector<T>& items, F f )
        {
            using U = decltype( f(T()) );
            
            std::vector< U > result;
            std::transform( items.begin(), items.end(), std::inserter( result, result.end() ), f );
            return result;
        }
        
        // F( T ) -> T and std::vector<T> is mutated
        template <typename T, typename F>
        std::vector<T>& transform( std::vector<T>& items, F f )
        {
            return Collection::transform( items, f );
        }
    }
}

Here are some unit tests that show how much easier it is to use the Map/Vector functions instead of going directly to STL – I’d argue that this code is comparable to F# for conciseness (although F# code would still benefit from pipelining subsequent operations).

#include <iostream>

#include <gmock/gmock.h>
#include <Vector.hpp>
#include <Map.hpp>

using namespace testing;
using namespace MusingStudio;

TEST( Map, to_vector )
{
    using Mapped = std::map<int, std::string>;
    using Tuples = std::vector<std::pair<int,std::string> >;
    
    Mapped items{ { 1, "Hi" }, { 2, "Bye" } };
    
    EXPECT_EQ( (Tuples{ { 1, "Hi" }, { 2, "Bye" } }), 
      Map::to_vector( items ) );
}

TEST( Map, keys )
{
    using Mapped = std::map<int, std::string>;
    using Keys = std::set<int>;
    
    Mapped items{ { 1, "Hi" }, { 2, "Bye" } };
    
    EXPECT_EQ( (Keys{ 1, 2 }), 
      Map::keys( items ) );
}

TEST( Map, values )
{
    using Mapped = std::map<int, std::string>;
    using Values = std::vector<std::string>;
    
    Mapped items{ { 1, "Hi" }, { 2, "Bye" } };
    
    EXPECT_EQ( (Values{ "Hi", "Bye" }), 
      Map::values( items ) );
}

TEST( Map, filter )
{
    using Mapped = std::map<int, int>;
    
    Mapped items{ {1,1}, {2,4}, {3,9}, {4,16} };
    
    Mapped even_keys{ {2,4},{4,16} };
    auto lambda = []( const auto& keyvaluepair ){ return keyvaluepair.first % 2 == 0; };
    
    // Map::filter will mutate parameter 'items'
    EXPECT_EQ( even_keys, 
      Map::filter( items, lambda ) );
    EXPECT_EQ( 2, items.size() );
}

TEST( Map, filter_copy )
{
    using Mapped = std::map<int, int>;
    
    Mapped items{ {1,1}, {2,4}, {3,9}, {4,16} };
    
    Mapped even_keys{ {2,4},{4,16} };
    auto lambda = []( const auto& keyvaluepair ){ return keyvaluepair.first % 2 == 0; };
    
    // Map::filter_copy creates a copy, so parameter 'items' is untouched
    EXPECT_EQ( even_keys, 
      Map::filter_copy( items, lambda ) );
    EXPECT_EQ( 4, items.size() );
}

TEST( Map, fold )
{
    using Mapped = std::map<int, int>;
    
    Mapped items{ {1,2}, {3,4}, {5,6} };
    
    // Map::fold takes F( A&, pair<K,V> ) -> void
    EXPECT_EQ( 21, 
      Map::fold( items, 
        []( int& acc, const auto& kvp ){ acc += kvp.first + kvp.second; }, 0 ) );
}

TEST( Map, transform_mutable_values_only )
{
    using Transformed = std::map<int, int>;

    // Map::transform over the values, mutating them
    // Takes F(K,V) -> V i.e. the type of the return value must be V
    // "items" must be a named variable because parameter is non-const
    // (we will mutate it)
    Transformed items = { {1,1}, {2,2}, {3,3} };
    
    EXPECT_EQ( (Transformed{ {1,1}, {2,4}, {3,9} }),
      Map::transform( items, 
        []( int _, int v ){ return v*v; } ) );
}

TEST( Map, map_keys_and_values )
{
    using Mapped = std::map<int,double>;
    
    // Map::map over the pairs<key,value>
    // Takes F( pair<K,V> ) -> pair<K',V'> 
    // i.e. both key and value types can change
    auto lambda =
        []( const auto& kvp )
        {
            return std::make_pair( kvp.first + kvp.second,
                                   (double)kvp.second / (double)kvp.first );
        };
    
    // Map::map - new keys and values, not mutating the original collection, 
    // can be passed as unnamed temporary
    EXPECT_EQ( (Mapped{ {2,1}, {6,2} }),
      Map::map( std::map<int,int>{ {1,1}, {2,4} }, lambda ) );
}

TEST( Vector, filter )
{
    std::vector<int> items{ 1,2,3,4,5,4,3,2,1 };

    // Vector::filter will mutate the input collection
    EXPECT_EQ( (std::vector<int>{1,2,2,1}),
      Vector::filter( items,
        [](int i){ return 0 <= i && i <= 2; } ) );
    EXPECT_EQ( 4, items.size() );
}

TEST( Vector, filter_copy )
{
    std::vector<int> items{ 1,2,3,4,5,4,3,2,1 };
    auto untouched_size = items.size();
    
    // Vector::filter_copy creates a copy, so parameter 'items' is untouched
    EXPECT_EQ( (std::vector<int>{1,2,2,1}),
      Vector::filter_copy( items,
        [](int i){ return 0 <= i && i <= 2; } ) );
    EXPECT_EQ( untouched_size, items.size() );
}

TEST( Vector, fold )
{
    std::vector<int> items{ 1, 2, 3, 4, 5, -4, -6, -2, -1 };
    // Vector::fold takes F( A&, T ) -> void
    auto accumulate_squares = []( std::set<int>& acc, int i ){ acc.insert(i*i); };
    std::set<int> expected{1, 4, 9, 16, 25, 36};
    EXPECT_EQ( expected, 
      Vector::fold( items, accumulate_squares, std::set<int>{} ) );
}

TEST( Vector, sort )
{
    // Vector::sort mutates the input, hence input is non-const reference
    std::vector<int> items{1,2,1,3};
    EXPECT_EQ( (std::vector<int>{1,1,2,3}), 
      Vector::sort( items ) );
}

TEST( Vector, sort_copy )
{
    // Vector::sort_copy copies the input collection,
    // so collection parameter is const& (and can be an unnamed temporary)
    EXPECT_EQ( (std::vector<int>{1,1,2,3}), 
      Vector::sort_copy( std::vector<int>{1,2,1,3} ) );
}

TEST( Vector, map )
{
    // Vector::map takes F(T) -> U
    // Input collection is const and a new collection is returned
    EXPECT_EQ( (std::vector<double>{ 1.1, 2.1, 3.1 }),
      Vector::map( std::vector<int>{ 1,2,3 },
        []( int i ){ return (double)i + 0.1; } ) );
}

TEST( Vector, transform )
{
    // Vector::transform takes F(T) -> T
    // Input collection is mutated
    std::vector<int> items{ 1, 2, 3 };
    EXPECT_EQ( (std::vector<int>{ 2, 4, 6 }),
      Vector::transform( items, []( int i ){ return 2*i; } ) );
}

int main(int argc, char* argv[]) 
{    
    InitGoogleMock( &argc, argv );
    return RUN_ALL_TESTS();   
}

Notice that we can bypass immutability in C++, so whereas in F# Map::filter would always create a copy, it could be preferable in C++ to filter in-place. With that in mind, I’ve written both filter and filter_copy variations. There’s a similar dilemma for map operations – if you want free rein over the output types, then use Map::map or Vector::map. But if you want to transform the data in place (sticking to the existing types), use Map::transform or Vector::transform.

That covers the most popular functions for just Map and Vector, but it would be straight-forward to extend the library to cover List, Set and others. Similarly, I’d like to extend it to include higher-order functions like Choose, but I’ll need C++17’s std::optional for that.

2 Comments

Filed under C++, C++ Code, Programming

How to calculate Fibonacci sequence for large numbers

Firstly, here’s a basic implementation of Fibonacci using recursion:

    using Numbers = std::vector<unsigned long long>;
    Numbers cache( 5000, 0 );
    
    unsigned long long fibonacciRecursiveImpl( size_t n, Numbers& cache )
    {
        if ( n == 0 ) return 0ULL;
        if ( n == 1 ) return 1ULL;
        
        if ( cache[n] != 0 ) return cache[n];
        
        auto num = fibonacciRecursiveImpl( n-1, cache ) + fibonacciRecursiveImpl( n-2, cache );
        cache[n] = num;
        return num;
    }
    
    unsigned long long fibonacciRecursive( size_t n )
    {
        return fibonacciRecursiveImpl( n, cache );
    }

This works fine for small numbers (e.g. up to 20). I was interested to know where it would fail first, data overflow (due to the size of the numbers involved) or stack overflow (due to the recursive approach)? I implemented this on a MacBook using Apple LLVM 7.0 and the C++14 flag, but no other special switches set. It turns out that the overflow problems kick in for n > 93, but there was no sign of stack overflows, even up for n ~ 2000.

Even if there was a stack overflow, you could still use recursion, but change to a tail recursive solution (many C++ compilers will eliminate tail calls, essentially turning this into an iterative solution):

    // Implement recursive approach with possibility for "Tail call elimination",
    // avoids any concerns about stack overflow
    unsigned long long fibonacciTailRecursiveImpl( size_t n, unsigned long long prevprev, unsigned long long prev )
    {
        if ( n == 0 ) return prevprev;
        if ( n == 1 ) return prev;
        
        return fibonacciTailRecursiveImpl( n-1, prev, prevprev + prev );
    }
    
    unsigned long long fibonacciTailRecursive( size_t n )
    {
        return fibonacciTailRecursiveImpl( n, 0, 1 );
    }

So how to avoid the data overflow for Fibonacci over n = 93? At that point, you need to introduce a large number type with its own internal representation of large integers and suitable operator+ implementation. I used one such, BigNum, for this HackerRank challenge. The implementation stores the large integer as a string, with each character taking values in the range [0,100]. This reduces the number of digit-by-digit additions by a factor of 10 compared to a straight decimal representation.

I replaced unsigned long long with BigNum in the recursive solution above and verified that it returns the correct answers for n ~ 2000, with no stack overflows. Here, I’ll show it in an iterative solution (if you don’t want to keep a cache around, this is highly memory efficient, because you only need to store the previous two values):

    BigNum fibonacciBigNumIterativeImpl( size_t n )
    {
        BigNum prevprev(0);
        BigNum prev(1);
        
        if ( n == 0 ) return prevprev;
        if ( n == 1 ) return prev;
        
        BigNum num(0);
        for ( size_t i = 2; i <= n; ++i )
        {
            num = prevprev + prev;
            prevprev = prev;
            prev = num;
        }
        return num;
    }
    
    std::string fibonacciBigNumIterative( size_t n )
    {
        auto result = fibonacciBigNumIterativeImpl( n );
        return result.toString();
    }

Leave a comment

Filed under C++, C++ Code, Programming

How to approximate Pi using C++11 random number generators

The other day I learnt a method to approximate Pi if you have a random number generator for the range [0,1]. Consider a unit circle centred on (0,0) in 2D coordinates. Then one quarter of the area of the circle lies in the quadrant formed by x and y in range [0,1]. The area of the quarter circle is Pi * R^2/4, and here, R = 1 (it’s a unit circle). So we can generate a bunch of random 2D points, calculate the ratio between those that fall inside the circle and those in the outer unit square, then multiply by 4 to approximate Pi.
approximatepi
That sounds like a neat test case for the C++11 random number generators, so I thought I’d try it out. It turns out to work pretty well, if you’re prepared to use a sufficiently large number of random values.

    double approxPi( size_t points )
    {
        std::random_device rand_device;
        
        // mersenne_twister_engine is a random number engine 
        // based on Mersenne Twister algorithm.
        std::mt19937 generator( rand_device() );
        
        // We want random values uniformly distributed in [0,1]
        std::uniform_real_distribution<> unif_zero_one(0, 1);
        
        size_t points_inside{0};
        
        for (int i = 0; i < points; ++i )
        {
            auto x = unif_zero_one( generator );
            auto y = unif_zero_one( generator );
            double d = std::sqrt( x*x + y*y );
            
            if ( d <= 1.0 )
                ++points_inside;
        }
                
        return 4.0 * (static_cast<double>(points_inside) / static_cast<double>(points));
    }
}

void testApproximatePi()
{
    SHOULD_BE_APPROX( 3.14159, 0.3, approxPi( 100 ) );
    SHOULD_BE_APPROX( 3.14159, 0.1, approxPi( 1000 ) );
    SHOULD_BE_APPROX( 3.14159, 0.01, approxPi( 100000 ) );
    SHOULD_BE_APPROX( 3.14159, 0.001, approxPi( 10000000 ) );
    
    std::cout << "\n";
}

A typical run gives reasonable approximations once you get over 100,000 points:
screen-shot-2016-09-16-at-13-30-26

Leave a comment

Filed under C++, C++ Code, Programming