1 # ifndef CPPAD_CORE_CHECKPOINT_HPP
2 # define CPPAD_CORE_CHECKPOINT_HPP
246 # define CPPAD_TMP_MULTI_THREAD 0
248 # ifdef CPPAD_FOR_TMB
249 # undef CPPAD_TMP_MULTI_THREAD
250 # define CPPAD_TMP_MULTI_THREAD 1
255 # if ! CPPAD_TMP_MULTI_THREAD
258 template <
class Base>
282 bool transpose =
false;
283 bool dependency =
true;
284 size_t n =
f_.Domain();
285 size_t m =
f_.Range();
291 for(
size_t j = 0; j < n; j++)
295 f_.ForSparseJacCheckpoint(
298 f_.size_forward_set(0);
303 for(
size_t i = 0; i < m; i++)
307 f_.RevSparseJacCheckpoint(
317 bool transpose =
false;
318 bool dependency =
true;
319 size_t n =
f_.Domain();
320 size_t m =
f_.Range();
325 for(
size_t j = 0; j < n; j++)
326 {
for(
size_t i = 0; i < n; i++)
327 identity[ i * n + j ] = (i == j);
330 n, identity, transpose, dependency
332 f_.size_forward_bool(0);
336 for(
size_t j = 0; j < m; j++)
337 {
for(
size_t i = 0; i < m; i++)
338 identity[ i * m + j ] = (i == j);
341 m, identity, transpose, dependency
351 size_t n =
f_.Domain();
352 size_t m =
f_.Range();
356 for(
size_t i = 0; i < m; i++)
362 for(
size_t j = 0; j < n; j++)
368 bool transpose =
false;
369 bool dependency =
false;
370 f_.ForSparseJacCheckpoint(
373 f_.RevSparseHesCheckpoint(
380 f_.size_forward_set(0);
385 size_t n =
f_.Domain();
386 size_t m =
f_.Range();
390 for(
size_t i = 0; i < m; i++)
395 for(
size_t j = 0; j < n; j++)
396 {
for(
size_t i = 0; i < n; i++)
397 identity[ i * n + j ] = (i == j);
401 bool transpose =
false;
402 bool dependency =
false;
403 f_.ForSparseJac(n, identity, transpose, dependency);
408 f_.size_forward_bool(0);
418 template <
class sparsity_type>
421 const sparsity_type& r ,
425 size_t m =
f_.Range();
426 size_t n =
f_.Domain();
437 for(
size_t i = 0; i < m; i++)
438 {
for(
size_t k = 0; k < q; k++)
439 s[i * q + k] =
false;
442 for(
size_t i = 0; i < m; i++)
443 {
for(
size_t k = 0; k < q; k++)
447 for(
size_t j = 0; j < n; j++)
449 bool R_jk = r[j * q + k ];
450 s_ik |= ( J_ij & R_jk );
463 template <
class sparsity_type>
466 const sparsity_type& rt ,
470 size_t m =
f_.Range();
471 size_t n =
f_.Domain();
483 for(
size_t i = 0; i < q; i++)
484 {
for(
size_t j = 0; j < n; j++)
488 for(
size_t k = 0; k < m; k++)
490 bool R_ik = rt[ k * q + i ];
492 s_ij |= (R_ik & J_kj);
495 st[ j * q + i ] = s_ij;
505 template <
class sparsity_type>
511 const sparsity_type& r ,
512 const sparsity_type& u ,
515 {
size_t n =
f_.Domain();
517 size_t m =
f_.Range();
538 t =
f_.RevSparseJac(1, s);
540 for(
size_t j = 0; j < n; j++)
549 bool transpose =
true;
550 sparsity_type a(n * q);
551 a =
f_.RevSparseJac(q, u, transpose);
556 for(
size_t i = 0; i < n; i++)
557 {
for(
size_t k = 0; k < q; k++)
561 for(
size_t j = 0; j < n; j++)
563 bool r_jk = r[j * q + k];
564 h_ik |= ( f_ij & r_jk );
572 for(
size_t i = 0; i < n; i++)
573 {
for(
size_t k = 0; k < q; k++)
575 v[ i * q + k ] =
bool(v[ i * q + k]) | bool(a[ i * q + k]);
607 template <
class Algo,
class ADVector>
617 { CheckSimpleVector< CppAD::AD<Base> , ADVector>();
630 f_.check_for_nan(
false);
637 f_.compare_change_count(0);
644 {
return f_.size_var(); }
663 template <
class ADVector>
664 void operator()(
const ADVector& ax, ADVector& ay,
size_t id = 0)
667 "checkpoint: id is non-zero in atom_fun(ax, ay, id)"
684 {
size_t n =
f_.Domain();
685 size_t m =
f_.Range();
715 for(
size_t i = 0; i < m; i++)
724 vy[i] |=
static_cast<bool>( vx[j] );
737 for(
size_t i = 0; i < m; i++)
739 for(
size_t j = 0; j < n; j++)
743 vy[i] |=
static_cast<bool>( vx[j] );
750 ty =
f_.Forward(q, tx);
757 f_.capacity_order(c, r);
774 size_t n =
f_.Domain();
775 size_t m =
f_.Range();
795 for(i = 0; i < m; i++)
796 {
for(k = 0; k <= q; k++)
803 px =
f_.Reverse(q+1, py);
809 f_.capacity_order(c, r);
823 {
return for_sparse_jac< vectorBool >(q, r, s, x);
835 {
return for_sparse_jac< vector<bool> >(q, r, s, x);
844 const vector< std::set<size_t> >& r ,
845 vector< std::set<size_t> >& s ,
848 size_t m =
f_.Range();
849 size_t n =
f_.Domain();
861 for(
size_t i = 0; i < m; i++)
865 for(
size_t i = 0; i < m; i++)
872 { std::set<size_t>::const_iterator itr_j;
873 const std::set<size_t>& r_j( r[j] );
874 for(itr_j = r_j.begin(); itr_j != r_j.end(); itr_j++)
896 {
return rev_sparse_jac< vectorBool >(q, rt, st, x);
908 {
return rev_sparse_jac< vector<bool> >(q, rt, st, x);
917 const vector< std::set<size_t> >& rt ,
918 vector< std::set<size_t> >& st ,
921 size_t m =
f_.Range();
922 size_t n =
f_.Domain();
935 for(
size_t j = 0; j < n; j++)
940 for(
size_t i = 0; i < m; i++)
942 std::set<size_t>::const_iterator itr_i;
943 const std::set<size_t>& r_i( rt[i] );
944 for(itr_i = r_i.begin(); itr_i != r_i.end(); itr_i++)
979 {
return rev_sparse_hes< vectorBool >(vx, s, t, q, r, u, v, x);
995 {
return rev_sparse_hes< vector<bool> >(vx, s, t, q, r, u, v, x);
1007 const vector< std::set<size_t> >& r ,
1008 const vector< std::set<size_t> >& u ,
1009 vector< std::set<size_t> >& v ,
1011 {
size_t n =
f_.Domain();
1013 size_t m =
f_.Range();
1034 t =
f_.RevSparseJac(1, s);
1036 for(
size_t j = 0; j < n; j++)
1046 bool transpose =
true;
1048 a =
f_.RevSparseJac(q, u, transpose);
1053 for(
size_t i = 0; i < n; i++)
1058 size_t j = *set_itr;
1060 { std::set<size_t>::const_iterator itr_j;
1061 const std::set<size_t>& r_j( r[j] );
1062 for(itr_j = r_j.begin(); itr_j != r_j.end(); itr_j++)
1063 {
size_t k = *itr_j;
1070 std::set<size_t>::const_iterator itr;
1071 for(
size_t i = 0; i < n; i++)
1072 {
for(itr = a[i].begin(); itr != a[i].end(); itr++)
1085 # define THREAD omp_get_thread_num()
1086 template <
class Base>
1087 class checkpoint :
public atomic_base<Base> {
1091 typedef typename atomic_base<Base>::option_enum
option_enum;
1094 vector< ADFun<Base> >
f_;
1105 {
return static_cast< atomic_base<Base>*
>(
this)->
sparsity(); }
1110 bool transpose =
false;
1111 bool dependency =
true;
1112 size_t n =
f_[THREAD].Domain();
1113 size_t m =
f_[THREAD].Range();
1117 { local::sparse_list identity;
1118 identity.resize(n, n);
1119 for(
size_t j = 0; j < n; j++)
1121 identity.add_element(j, j);
1123 f_[THREAD].ForSparseJacCheckpoint(
1126 f_[THREAD].size_forward_set(0);
1129 { local::sparse_list identity;
1130 identity.resize(m, m);
1131 for(
size_t i = 0; i < m; i++)
1133 identity.add_element(i, i);
1135 f_[THREAD].RevSparseJacCheckpoint(
1145 bool transpose =
false;
1146 bool dependency =
true;
1147 size_t n =
f_[THREAD].Domain();
1148 size_t m =
f_[THREAD].Range();
1152 { vectorBool identity(n * n);
1153 for(
size_t j = 0; j < n; j++)
1154 {
for(
size_t i = 0; i < n; i++)
1155 identity[ i * n + j ] = (i == j);
1158 n, identity, transpose, dependency
1160 f_[THREAD].size_forward_bool(0);
1163 { vectorBool identity(m * m);
1164 for(
size_t j = 0; j < m; j++)
1165 {
for(
size_t i = 0; i < m; i++)
1166 identity[ i * m + j ] = (i == j);
1169 m, identity, transpose, dependency
1179 size_t n =
f_[THREAD].Domain();
1180 size_t m =
f_[THREAD].Range();
1183 vector<bool> all_one(m);
1184 for(
size_t i = 0; i < m; i++)
1188 local::sparse_list identity;
1189 identity.resize(n, n);
1190 for(
size_t j = 0; j < n; j++)
1192 identity.add_element(j, j);
1196 bool transpose =
false;
1197 bool dependency =
false;
1198 f_[THREAD].ForSparseJacCheckpoint(
1201 f_[THREAD].RevSparseHesCheckpoint(
1208 f_[THREAD].size_forward_set(0);
1213 size_t n =
f_[THREAD].Domain();
1214 size_t m =
f_[THREAD].Range();
1217 vectorBool all_one(m);
1218 for(
size_t i = 0; i < m; i++)
1222 vectorBool identity(n * n);
1223 for(
size_t j = 0; j < n; j++)
1224 {
for(
size_t i = 0; i < n; i++)
1225 identity[ i * n + j ] = (i == j);
1229 bool transpose =
false;
1230 bool dependency =
false;
1231 f_[THREAD].ForSparseJac(n, identity, transpose, dependency);
1236 f_[THREAD].size_forward_bool(0);
1246 template <
class sparsity_type>
1249 const sparsity_type& r ,
1251 const vector<Base>& x )
1253 size_t m =
f_[THREAD].Range();
1254 size_t n =
f_[THREAD].Domain();
1265 for(
size_t i = 0; i < m; i++)
1266 {
for(
size_t k = 0; k < q; k++)
1267 s[i * q + k] =
false;
1270 for(
size_t i = 0; i < m; i++)
1271 {
for(
size_t k = 0; k < q; k++)
1275 for(
size_t j = 0; j < n; j++)
1277 bool R_jk = r[j * q + k ];
1278 s_ik |= ( J_ij & R_jk );
1280 s[i * q + k] = s_ik;
1291 template <
class sparsity_type>
1294 const sparsity_type& rt ,
1296 const vector<Base>& x )
1298 size_t m =
f_[THREAD].Range();
1299 size_t n =
f_[THREAD].Domain();
1311 for(
size_t i = 0; i < q; i++)
1312 {
for(
size_t j = 0; j < n; j++)
1316 for(
size_t k = 0; k < m; k++)
1318 bool R_ik = rt[ k * q + i ];
1320 s_ij |= (R_ik & J_kj);
1323 st[ j * q + i ] = s_ij;
1333 template <
class sparsity_type>
1335 const vector<bool>& vx ,
1336 const vector<bool>& s ,
1339 const sparsity_type& r ,
1340 const sparsity_type& u ,
1342 const vector<Base>& x )
1343 {
size_t n =
f_[THREAD].Domain();
1345 size_t m =
f_[THREAD].Range();
1366 t =
f_[THREAD].RevSparseJac(1, s);
1368 for(
size_t j = 0; j < n; j++)
1377 bool transpose =
true;
1378 sparsity_type a(n * q);
1379 a =
f_[THREAD].RevSparseJac(q, u, transpose);
1384 for(
size_t i = 0; i < n; i++)
1385 {
for(
size_t k = 0; k < q; k++)
1389 for(
size_t j = 0; j < n; j++)
1391 bool r_jk = r[j * q + k];
1392 h_ik |= ( f_ij & r_jk );
1395 v[i * q + k] = h_ik;
1400 for(
size_t i = 0; i < n; i++)
1401 {
for(
size_t k = 0; k < q; k++)
1403 v[ i * q + k ] =
bool(v[ i * q + k]) | bool(a[ i * q + k]);
1435 template <
class Algo,
class ADVector>
1439 const ADVector& ax ,
1442 atomic_base<Base>::pack_sparsity_enum ,
1443 bool optimize =
true
1446 f_( omp_get_max_threads() ) ,
1452 CheckSimpleVector< CppAD::AD<Base> , ADVector>();
1461 f_[0].Dependent(ay);
1465 f_[0].check_for_nan(
false);
1472 f_[0].compare_change_count(0);
1475 for(
int i = 1; i < omp_get_max_threads(); ++i)
1483 {
return f_[THREAD].size_var(); }
1502 template <
class ADVector>
1503 void operator()(
const ADVector& ax, ADVector& ay,
size_t id = 0)
1506 "checkpoint: id is non-zero in atom_fun(ax, ay, id)"
1519 const vector<bool>& vx ,
1521 const vector<Base>& tx ,
1523 {
size_t n =
f_[THREAD].Domain();
1524 size_t m =
f_[THREAD].Range();
1533 if( vx.size() == 0 )
1548 if(
sparsity() == atomic_base<Base>::set_sparsity_enum )
1554 for(
size_t i = 0; i < m; i++)
1559 size_t j = *set_itr;
1563 vy[i] |=
static_cast<bool>( vx[j] );
1576 for(
size_t i = 0; i < m; i++)
1578 for(
size_t j = 0; j < n; j++)
1582 vy[i] |=
static_cast<bool>( vx[j] );
1589 ty =
f_[THREAD].Forward(q, tx);
1596 f_[THREAD].capacity_order(c, r);
1607 const vector<Base>& tx ,
1608 const vector<Base>& ty ,
1610 const vector<Base>& py )
1613 size_t n =
f_[THREAD].Domain();
1614 size_t m =
f_[THREAD].Range();
1626 f_[THREAD].Forward(q, tx);
1633 vector<Base> check_ty =
f_[THREAD].Forward(q, tx);
1634 for(i = 0; i < m; i++)
1635 {
for(k = 0; k <= q; k++)
1636 { j = i * (q+1) + k;
1642 px =
f_[THREAD].Reverse(q+1, py);
1648 f_[THREAD].capacity_order(c, r);
1659 const vectorBool& r ,
1661 const vector<Base>& x )
1662 {
return for_sparse_jac< vectorBool >(q, r, s, x);
1671 const vector<bool>& r ,
1673 const vector<Base>& x )
1674 {
return for_sparse_jac< vector<bool> >(q, r, s, x);
1683 const vector< std::set<size_t> >& r ,
1684 vector< std::set<size_t> >& s ,
1685 const vector<Base>& x )
1687 size_t m =
f_[THREAD].Range();
1688 size_t n =
f_[THREAD].Domain();
1700 for(
size_t i = 0; i < m; i++)
1704 for(
size_t i = 0; i < m; i++)
1709 size_t j = *set_itr;
1711 { std::set<size_t>::const_iterator itr_j;
1712 const std::set<size_t>& r_j( r[j] );
1713 for(itr_j = r_j.begin(); itr_j != r_j.end(); itr_j++)
1714 {
size_t k = *itr_j;
1732 const vectorBool& rt ,
1734 const vector<Base>& x )
1735 {
return rev_sparse_jac< vectorBool >(q, rt, st, x);
1744 const vector<bool>& rt ,
1746 const vector<Base>& x )
1747 {
return rev_sparse_jac< vector<bool> >(q, rt, st, x);
1756 const vector< std::set<size_t> >& rt ,
1757 vector< std::set<size_t> >& st ,
1758 const vector<Base>& x )
1760 size_t m =
f_[THREAD].Range();
1761 size_t n =
f_[THREAD].Domain();
1774 for(
size_t j = 0; j < n; j++)
1779 for(
size_t i = 0; i < m; i++)
1781 std::set<size_t>::const_iterator itr_i;
1782 const std::set<size_t>& r_i( rt[i] );
1783 for(itr_i = r_i.begin(); itr_i != r_i.end(); itr_i++)
1792 size_t j = *set_itr;
1810 const vector<bool>& vx ,
1811 const vector<bool>& s ,
1814 const vectorBool& r ,
1815 const vectorBool& u ,
1817 const vector<Base>& x )
1818 {
return rev_sparse_hes< vectorBool >(vx, s, t, q, r, u, v, x);
1826 const vector<bool>& vx ,
1827 const vector<bool>& s ,
1830 const vector<bool>& r ,
1831 const vector<bool>& u ,
1833 const vector<Base>& x )
1834 {
return rev_sparse_hes< vector<bool> >(vx, s, t, q, r, u, v, x);
1842 const vector<bool>& vx ,
1843 const vector<bool>& s ,
1846 const vector< std::set<size_t> >& r ,
1847 const vector< std::set<size_t> >& u ,
1848 vector< std::set<size_t> >& v ,
1849 const vector<Base>& x )
1850 {
size_t n =
f_[THREAD].Domain();
1852 size_t m =
f_[THREAD].Range();
1873 t =
f_[THREAD].RevSparseJac(1, s);
1875 for(
size_t j = 0; j < n; j++)
1885 bool transpose =
true;
1886 vector< std::set<size_t> > a(n);
1887 a =
f_[THREAD].RevSparseJac(q, u, transpose);
1892 for(
size_t i = 0; i < n; i++)
1897 size_t j = *set_itr;
1899 { std::set<size_t>::const_iterator itr_j;
1900 const std::set<size_t>& r_j( r[j] );
1901 for(itr_j = r_j.begin(); itr_j != r_j.end(); itr_j++)
1902 {
size_t k = *itr_j;
1909 std::set<size_t>::const_iterator itr;
1910 for(
size_t i = 0; i < n; i++)
1911 {
for(itr = a[i].begin(); itr != a[i].end(); itr++)
1925 # undef CPPAD_TMP_MULTI_THREAD
atomic_base< Base >::option_enum option_enum
same as option_enum in base class
ADFun< Base > f_
AD function corresponding to this checkpoint object.
#define CPPAD_ASSERT_KNOWN(exp, msg)
Check that exp is true, if not print msg and terminate execution.
size_t end(void) const
Fetch end for this vector of sets object.
virtual bool rev_sparse_hes(const vector< bool > &vx, const vector< bool > &s, vector< bool > &t, size_t q, const vector< bool > &r, const vector< bool > &u, vector< bool > &v, const vector< Base > &x)
Link from user_atomic to reverse sparse Hessian bool.
Vector of sets of positive integers, each set stored as a singly linked list.
Class used to hold function objects.
void add_element(size_t i, size_t element)
Add one element to a set.
Vector of sets of positive integers stored as singly linked lists with the element values strictly in...
virtual bool rev_sparse_jac(size_t q, const vector< bool > &rt, vector< bool > &st, const vector< Base > &x)
Link from user_atomic to reverse sparse Jacobian bool.
vectorBool hes_sparse_bool_
atomic_base(void)
Base class for atomic_user functions.
virtual bool for_sparse_jac(size_t q, const vector< std::set< size_t > > &r, vector< std::set< size_t > > &s, const vector< Base > &x)
Link from user_atomic to forward sparse Jacobian sets.
virtual bool for_sparse_jac(size_t q, const vectorBool &r, vectorBool &s, const vector< Base > &x)
Link from user_atomic to forward sparse Jacobian pack.
Vector of sets of positive integers stored as a packed array of bools.
void clear(void)
free memory and set number of elements to zero
void set_jac_sparse_set(void)
set jac_sparse_set_
void Independent(VectorAD &x, size_t abort_op_index)
Declaration of independent variables.
virtual bool for_sparse_jac(size_t q, const vector< bool > &r, vector< bool > &s, const vector< Base > &x)
Link from user_atomic to forward sparse Jacobian bool.
bool for_sparse_jac(size_t q, const sparsity_type &r, sparsity_type &s, const vector< Base > &x)
Link from user_atomic to forward sparse Jacobian pack and bool.
size_t size(void) const
number of elements currently in this vector.
void set_jac_sparse_bool(void)
set jac_sparse_bool_
void operator()(const ADVector &ax, ADVector &ay, size_t id=0)
Implement the user call to afun(ax, ay) and old_atomic call to afun(ax, ay, id).
bool rev_sparse_hes(const vector< bool > &vx, const vector< bool > &s, vector< bool > &t, size_t q, const sparsity_type &r, const sparsity_type &u, sparsity_type &v, const vector< Base > &x)
Link from user_atomic to reverse sparse Hessian bools.
virtual bool reverse(size_t q, const vector< Base > &tx, const vector< Base > &ty, vector< Base > &px, const vector< Base > &py)
Link from user_atomic to reverse mode.
virtual bool rev_sparse_jac(size_t q, const vector< std::set< size_t > > &rt, vector< std::set< size_t > > &st, const vector< Base > &x)
Link from user_atomic to reverse Jacobian sets.
vectorBool jac_sparse_bool_
void set_hes_sparse_bool(void)
set hes_sparse_bool_
#define CPPAD_ASSERT_UNKNOWN(exp)
Check that exp is true, if not terminate execution.
virtual bool forward(size_t p, size_t q, const vector< bool > &vx, vector< bool > &vy, const vector< Base > &tx, vector< Base > &ty)
Link from user_atomic to forward mode.
size_t size_var(void)
Implement the user call to atom_fun.size_var().
checkpoint(const char *name, Algo &algo, const ADVector &ax, ADVector &ay, option_enum sparsity=atomic_base< Base >::pack_sparsity_enum, bool optimize=true)
Constructor of a checkpoint object.
void resize(size_t n_set, size_t end)
Start a new vector of sets.
local::sparse_list jac_sparse_set_
sparsity for entire Jacobian f(x)^{(1)} does not change so can cache it
static void clear(void)
Free all thread_alloc static memory held by atomic_base (avoids reallocations).
virtual bool rev_sparse_hes(const vector< bool > &vx, const vector< bool > &s, vector< bool > &t, size_t q, const vectorBool &r, const vectorBool &u, vectorBool &v, const vector< Base > &x)
Link from user_atomic to reverse sparse Hessian pack.
void set_hes_sparse_set(void)
set hes_sparse_set_
size_t n_set(void) const
Fetch n_set for vector of sets object.
bool rev_sparse_jac(size_t q, const sparsity_type &rt, sparsity_type &st, const vector< Base > &x)
Link from user_atomic to reverse sparse Jacobian pack and bool.
local::sparse_list hes_sparse_set_
sparsity for sum_i f_i(x)^{(2)} does not change so can cache it
virtual bool rev_sparse_hes(const vector< bool > &vx, const vector< bool > &s, vector< bool > &t, size_t q, const vector< std::set< size_t > > &r, const vector< std::set< size_t > > &u, vector< std::set< size_t > > &v, const vector< Base > &x)
Link from user_atomic to reverse sparse Hessian sets.
sparse_list_const_iterator const_iterator
declare a const iterator
cons_iterator for one set of positive integers in a sparse_list object.
void operator()(const ADVector &ax, ADVector &ay, size_t id=0)
Implement the user call to atom_fun(ax, ay).
virtual bool rev_sparse_jac(size_t q, const vectorBool &rt, vectorBool &st, const vector< Base > &x)
Link from user_atomic to reverse sparse Jacobian pack.
size_t size(void) const
number of elements in this vector
option_enum sparsity(void)