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
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
|
// SPDX-License-Identifier: EPL-2.0 OR GPL-2.0-or-later
// SPDX-FileCopyrightText: Bradley M. Bell <bradbell@seanet.com>
// SPDX-FileContributor: 2003-22 Bradley M. Bell
// ----------------------------------------------------------------------------
# include <cppad/cppad.hpp>
# include <omp.h>
namespace {
typedef CPPAD_TESTVECTOR(double) d_vector;
typedef CPPAD_TESTVECTOR( CppAD::AD<double> ) ad_vector;
// algorithm that we are checkpoingint
const size_t length_of_sum_ = 5000;
void long_sum_algo(const ad_vector& ax, ad_vector& ay)
{ ay[0] = 0.0;
for(size_t i = 0; i < length_of_sum_; ++i)
ay[0] += ax[0];
return;
}
// inform CppAD if we are in parallel mode
bool in_parallel(void)
{ return omp_in_parallel() != 0; }
//
// inform CppAD of the current thread number
size_t thread_num(void)
{ return static_cast<size_t>( omp_get_thread_num() ); }
}
// multi_thread_checkpoint
bool multi_chkpoint_one(void)
{ bool ok = true;
// OpenMP setup
size_t num_threads = 4; // number of threads
omp_set_dynamic(0); // turn off dynamic thread adjustment
omp_set_num_threads( int(num_threads) ); // set number of OMP threads
// check that multi-threading is possible on this machine
if( omp_get_max_threads() < 2 )
{ std::cout << "This machine does not support multi-threading: ";
}
// create checkpoint version of algorithm
size_t n(1), m(1);
ad_vector ax(n), ay(m);
ax[0] = 2.0;
CppAD::atomic_base<double>::option_enum sparsity =
CppAD::atomic_base<double>::set_sparsity_enum;
bool optimize = false;
CppAD::checkpoint<double> atom_fun(
"long_sum", long_sum_algo, ax, ay, sparsity, optimize
);
// setup for using CppAD in parallel mode
CppAD::thread_alloc::parallel_setup(num_threads, in_parallel, thread_num);
CppAD::thread_alloc::hold_memory(true);
CppAD::parallel_ad<double>();
// place to hold result for each thread
d_vector y(num_threads);
for(size_t thread = 0; thread < num_threads; thread++)
y[thread] = 0.0;
# pragma omp parallel for
for(int thread = 0; thread < int(num_threads); thread++)
{ ad_vector au(n), av(m);
au[0] = 1.0;
CppAD::Independent(au);
atom_fun(au, av);
CppAD::ADFun<double> f(au, av);
//
d_vector x(n), v(m);
x[0] = double( thread + 1 );
v = f.Forward(0, x);
//
// this assignment has false sharing; i.e., will case cache resets
// (conversion avoids boost vector conversion warning)
y[size_t(thread)] = v[0];
}
// check the results
for(size_t thread = 0; thread < num_threads; thread++)
{ double check = double( length_of_sum_ * (thread + 1) );
ok &= check == y[thread];
}
return ok;
}
|
