fibonacci.cpp

英特尔&#174 线程构建模块（英特尔&#174 TBB）是一个屡获殊荣的 C++ 运行时库

    //! fill task arguments
    QueueInsertTask( int n, QueueStream &s ) : my_n(n), my_stream(s) { }
    //! executing task
    /*override*/ task* execute() {
        // Execute of parallel pushing of n-1 initial matrices
        parallel_for( blocked_range<int>( 1, my_n, 10 ), parallel_forFibBody(my_stream) ); 
        my_stream.producer_is_done = true; 
        return 0;
    }
};
//! Parallel queue's processing task
struct QueueProcessTask: public task {
    QueueStream &my_stream;
    //! fill task argument
    QueueProcessTask( QueueStream &s ) : my_stream(s) { }
    //! executing task
    /*override*/ task* execute() {
        while( !my_stream.producer_is_done || my_stream.Queue.size()>1 ) {
            parallel_while<parallel_whileFibBody> w; // run while loop in parallel
            w.run( my_stream, parallel_whileFibBody( w, my_stream ) );
        }
        return 0;
    }
};
//! Root function
value ParallelQueueFib(int n)
{
    QueueStream stream;
    stream.producer_is_done = false;
    task_list list;
    list.push_back(*new(task::allocate_root()) QueueInsertTask( n, stream ));
    list.push_back(*new(task::allocate_root()) QueueProcessTask( stream ));
    // last runs first - scheduler is as LIFO
    task::spawn_root_and_wait(list);
    assert(stream.Queue.size() == 1); // it is easy to lose some work
    Matrix2x2 M; stream.Queue.pop( M ); // get last matrix
    return M.v[0][0]; // and result number
}

// *** Queue with pipeline *** //

//! filter to fills queue
class InputFilter: public filter {
    atomic<int> N; //< index of Fibonacci number
public:
    concurrent_queue<Matrix2x2> Queue;
    //! fill filter arguments
    InputFilter( int n ) : filter(false /*is not serial*/) { N = n; }
    //! executing filter
    /*override*/ void* operator()(void*)
    {
        int n = --N;
        if(n <= 0) return 0;
        Queue.push( Matrix1110 );
        return n == 1? 0 : &Queue; // one less multiplications
    }
};
//! filter to process queue
class MultiplyFilter: public filter {
public:
    MultiplyFilter( ) : filter(false /*is not serial*/) { }
    //! executing filter
    /*override*/ void* operator()(void*p)
    {
        concurrent_queue<Matrix2x2> &Queue = *static_cast<concurrent_queue<Matrix2x2> *>(p);
        Matrix2x2 m1, m2;
        Queue.pop( m1 ); Queue.pop( m2 ); // get two elements
        m1 = m1 * m2; // process them
        Queue.push( m1 ); // and push back
        return this; // just nothing
    }
};
//! Root function
value ParallelPipeFib(int n)
{
    InputFilter input( n );
    MultiplyFilter process;
    // Create the pipeline
    pipeline pipeline;
    // add filters
    pipeline.add_filter( input ); // first
    pipeline.add_filter( process ); // second

    // Run the pipeline
    pipeline.run( n ); // must be larger then max threads number
    pipeline.clear(); // do not forget clear the pipeline

    Matrix2x2 M; input.Queue.pop( M ); // get last element
    return M.v[0][0]; // get value
}

// *** parallel_reduce *** //

//! Functor for parallel_reduce
struct parallel_reduceFibBody {
    Matrix2x2 sum;
    int splitted;  //< flag to make one less operation for splitted bodies
    //! Constructor fills sum with initial matrix
    parallel_reduceFibBody() : sum( Matrix1110 ), splitted(0) { }
    //! Splitting constructor
    parallel_reduceFibBody( parallel_reduceFibBody& other, split ) : sum( Matrix1110 ), splitted(1/*note that it is splitted*/) {}
    //! Join point
    void join( parallel_reduceFibBody &s ) {
        sum = sum * s.sum;
    }
    //! Process multiplications
    void operator()( const blocked_range<int> &r ) {
        for( int k = r.begin() + splitted; k < r.end(); ++k )
            sum = sum * Matrix1110;
        splitted = 0; // reset flag, because this method can be reused for next range
    }
};
//! Root function
value parallel_reduceFib(int n)
{
    parallel_reduceFibBody b;
    parallel_reduce(blocked_range<int>(2, n, 3), b); // do parallel reduce on range [2, n) for b
    return b.sum.v[0][0];
}

// *** parallel_scan *** //

//! Functor for parallel_scan
struct parallel_scanFibBody {
    Matrix2x2 sum;
    int first;  // flag to make one less operation for first range
    //! Constructor fills sum with initial matrix
    parallel_scanFibBody() : sum( Matrix1110 ), first(1) {}
    //! Splitting constructor
    parallel_scanFibBody( parallel_scanFibBody &b, split) : sum( Matrix1110 ), first(1) {}
    //! Join point
    void reverse_join( parallel_scanFibBody &a ) {
        sum = sum * a.sum;
    }
    //! Assign point
    void assign( parallel_scanFibBody &b ) {
        sum = b.sum;
    }
    //! Process multiplications. For two tags
    template<typename T>
    void operator()( const blocked_range<int> &r, T) {
        // see tag.is_final_scan() for what tag is used
        for( int k = r.begin() + first; k < r.end(); ++k )
            sum = sum * Matrix1110;
        first = 0; // reset flag, because this method can be reused for next range
    }
};
//! Root function
value parallel_scanFib(int n)
{
    parallel_scanFibBody b;
    parallel_scan(blocked_range<int>(1/*one less, because body skip first*/, n, 3), b);
    return b.sum.v[0][0];
}

// *** Raw tasks *** //

//! task class which computes Fibonacci numbers by Lucas formula
struct FibTask: public task {
    const int n;
    value& sum;
    value x, y;
    bool second_phase; //< flag of continuation
    // task arguments
    FibTask( int n_, value& sum_ ) : 
        n(n_), sum(sum_), second_phase(false)
    {}
    //! Execute task
    /*override*/ task* execute() {
        // Using Lucas' formula here
        if( second_phase ) { // childs finished
            sum = n&1 ? x*x + y*y : x*x - y*y;
            return NULL;
        }
        if( n <= 2 ) {
            sum = n!=0;
            return NULL;
        } else {
            recycle_as_continuation();  // repeat this task when childs finish
            second_phase = true; // mark second phase
            FibTask& a = *new( allocate_child() ) FibTask( n/2 + 1, x );
            FibTask& b = *new( allocate_child() ) FibTask( n/2 - 1 + (n&1), y );
            set_ref_count(2);
            spawn( a );
            return &b;
        }
    }
};
//! Root function
value ParallelTaskFib(int n) { 
    value sum;
    FibTask& a = *new(task::allocate_root()) FibTask(n, sum);
    task::spawn_root_and_wait(a);
    return sum;
}

/////////////////////////// Main ////////////////////////////////////////////////////

//! A closed range of int.
struct IntRange {
    int low;
    int high;
    void set_from_string( const char* s );
    IntRange( int low_, int high_ ) : low(low_), high(high_) {}
};

void IntRange::set_from_string( const char* s ) {
    char* end;
    high = low = strtol(s,&end,0);
    switch( *end ) {
    case ':': 
        high = strtol(end+1,0,0); 
        break;
    case '\0':
        break;
    default:
        printf("unexpected character = %c\n",*end);
    }
}

//! Tick count for start
static tick_count t0;

typedef value (*MeasureFunc)(int);
//! Measure ticks count in loop [2..n]
value Measure(const char *name, MeasureFunc func, int n)
{
    value result;
    printf(name); t0 = tick_count::now();
    for(int number = 2; number <= n; number++)
        result = func(number);
    printf("\t- in %f msec\n", (tick_count::now() - t0).seconds()*1000);
    return result;
}

//! program entry
int main(int argc, char* argv[])
{
    int NumbersCount = argc>1 ? strtol(argv[1],0,0) : 100;
    IntRange NThread(1,4);// Number of threads to use.
    if(argc>2) NThread.set_from_string(argv[2]);
    unsigned long ntrial = argc>3? (unsigned long)strtoul(argv[3],0,0) : 1;
    value result, sum;

    printf("Fibonacci numbers example. Generating %d numbers..\n",  NumbersCount);

    result = Measure("Serial loop", SerialFib, NumbersCount);
    sum = Measure("Serial matrix", SerialMatrixFib, NumbersCount); assert(result == sum);
    sum = Measure("Serial vector", SerialVectorFib, NumbersCount); assert(result == sum);
    sum = Measure("Serial queue", SerialQueueFib, NumbersCount); assert(result == sum);
    // now in parallel
    for( unsigned long i=0; i<ntrial; ++i ) {
        for(int threads = NThread.low; threads <= NThread.high; threads++)
        {
            task_scheduler_init scheduler_init(threads); printf("Threads number is %d\n", threads);

            sum = Measure("Shared serial (mutex)\t", SharedSerialFib<mutex>, NumbersCount); assert(result == sum);
            sum = Measure("Shared serial (spin_mutex)", SharedSerialFib<spin_mutex>, NumbersCount); assert(result == sum);
            sum = Measure("Shared serial (queuing_mutex)", SharedSerialFib<queuing_mutex>, NumbersCount); assert(result == sum);
            sum = Measure("Shared serial (Conc.HashTable)", ConcurrentHashSerialFib, NumbersCount); assert(result == sum);
            sum = Measure("Parallel while+for/queue", ParallelQueueFib, NumbersCount); assert(result == sum);
            sum = Measure("Parallel pipe/queue\t", ParallelPipeFib, NumbersCount); assert(result == sum);
            sum = Measure("Parallel reduce\t\t", parallel_reduceFib, NumbersCount); assert(result == sum);
            sum = Measure("Parallel scan\t\t", parallel_scanFib, NumbersCount); assert(result == sum);
            sum = Measure("Parallel tasks\t\t", ParallelTaskFib, NumbersCount); assert(result == sum);
        }

    #ifdef __GNUC__
        printf("Fibonacci number #%d modulo 2^64 is %lld\n\n", NumbersCount, result);
    #else
        printf("Fibonacci number #%d modulo 2^64 is %I64d\n\n", NumbersCount, result);
    #endif
    }
    return 0;
}

// Utils

void Matrix2x2Multiply(const value a[2][2], const value b[2][2], value c[2][2])
{
    for( int i = 0; i <= 1; i++)
        for( int j = 0; j <= 1; j++)
            c[i][j] = a[i][0]*b[0][j] + a[i][1]*b[1][j];
}

Matrix2x2 Matrix2x2::operator *(const Matrix2x2 &to) const
{
    Matrix2x2 result;
    Matrix2x2Multiply(v, to.v, result.v);
    return result;
}