18#if defined(DEAL_II_WITH_ADOLC) || defined(DEAL_II_TRILINOS_WITH_SACADO)
23# include <type_traits>
37 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
39 const unsigned int n_independent_variables,
40 const unsigned int n_dependent_variables)
41 : independent_variable_values(
42 n_independent_variables,
44 , registered_independent_variable_values(n_independent_variables, false)
45 , registered_marked_independent_variables(n_independent_variables, false)
46 , registered_marked_dependent_variables(n_dependent_variables, false)
55 "Floating point/arithmetic numbers have no derivatives."));
59 "The AD number type does not support the calculation of any derivatives."));
79 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
82 ScalarType>::reset_registered_independent_variables()
84 std::fill(registered_independent_variable_values.begin(),
85 registered_independent_variable_values.end(),
91 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
96 std::fill(registered_marked_dependent_variables.begin(),
97 registered_marked_dependent_variables.end(),
103 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
106 const unsigned int index,
114 if (this->is_recording() ==
false)
115 start_recording_operations(1 );
117 Assert(this->is_recording() ==
true,
119 "Cannot change the value of an independent variable "
120 "of the tapeless variety while this class is not set "
121 "in recording operations."));
125 Assert(this->active_tape_index() !=
130 index < n_independent_variables(),
132 "Trying to set the value of a non-existent independent variable."));
134 independent_variable_values[index] = value;
135 registered_independent_variable_values[index] =
true;
140 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
143 const unsigned int index,
147 Assert(registered_independent_variable_values[index] ==
true,
153 registered_marked_independent_variables[index - 1] ==
true,
155 "Need to extract sensitivities in the order they're created."));
162 Assert(is_recording() ==
true,
164 "The marking of independent variables is only valid "
165 "during recording."));
169 independent_variable_values[index],
171 this->n_independent_variables(),
173 registered_marked_independent_variables[index] =
true;
178 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
181 ScalarType>::finalize_sensitive_independent_variables()
const
184 Assert(n_registered_independent_variables() == n_independent_variables(),
185 ExcMessage(
"Not all values of sensitivities have been recorded!"));
188 if (this->independent_variables.size() == 0)
190 this->independent_variables.resize(
191 this->n_independent_variables(),
195 for (
unsigned int i = 0; i < this->n_independent_variables(); ++i)
196 this->mark_independent_variable(i, this->independent_variables[i]);
202 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
213 Assert(is_recording() ==
false,
215 "The initialization of non-sensitive independent variables is "
216 "only valid outside of recording operations."));
219 Assert(registered_independent_variable_values[index] ==
true,
222 out = independent_variable_values[index];
227 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
230 ScalarType>::n_registered_independent_variables()
const
232 return std::count(registered_independent_variable_values.begin(),
233 registered_independent_variable_values.end(),
239 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
243 return independent_variable_values.size();
248 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
253 return std::count(registered_marked_dependent_variables.begin(),
254 registered_marked_dependent_variables.end(),
260 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
264 return dependent_variables.size();
269 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
274 return taped_driver.is_recording();
276 return tapeless_driver.is_dependent_variable_marking_allowed();
281 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
287 return taped_driver.active_tape_index();
294 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
300 return taped_driver.is_registered_tape(tape_index);
307 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
312 const std::ios_base::fmtflags stream_flags(stream.flags());
314 stream.setf(std::ios_base::boolalpha);
316 stream <<
"Active tape index: " << active_tape_index() <<
"\n";
317 stream <<
"Recording? " << is_recording() <<
"\n";
318 stream << std::flush;
321 taped_driver.print(stream);
323 stream <<
"Registered independent variables: "
325 for (
unsigned int i = 0; i < n_independent_variables(); i++)
326 stream << registered_independent_variable_values[i]
327 << (i < (n_independent_variables() - 1) ?
"," :
"");
330 stream <<
"Independent variable values: "
332 print_values(stream);
334 stream <<
"Registered marked independent variables: "
336 for (
unsigned int i = 0; i < n_independent_variables(); i++)
337 stream << registered_marked_independent_variables[i]
338 << (i < (n_independent_variables() - 1) ?
"," :
"")
342 stream <<
"Dependent variable values: "
344 for (
unsigned int i = 0; i < n_dependent_variables(); i++)
345 stream << dependent_variables[i]
346 << (i < (n_dependent_variables() - 1) ?
"," :
"");
349 stream <<
"Registered dependent variables: "
351 for (
unsigned int i = 0; i < n_dependent_variables(); i++)
352 stream << registered_marked_dependent_variables[i]
353 << (i < (n_dependent_variables() - 1) ?
"," :
"");
357 stream.flags(stream_flags);
362 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
365 std::ostream &stream)
const
367 for (
unsigned int i = 0; i < n_independent_variables(); i++)
368 stream << independent_variable_values[i]
369 << (i < (n_independent_variables() - 1) ?
"," :
"")
376 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
380 std::ostream & stream)
const
385 Assert(is_registered_tape(tape_index),
388 this->taped_driver.print_tape_stats(tape_index, stream);
393 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
396 const unsigned int n_independent_variables,
397 const unsigned int n_dependent_variables,
398 const bool clear_registered_tapes)
400 const unsigned int new_n_independent_variables =
402 n_independent_variables :
403 this->n_independent_variables());
404 const unsigned int new_n_dependent_variables =
406 n_dependent_variables :
407 this->n_dependent_variables());
421 std::vector<ad_type>().swap(independent_variables);
422 std::vector<ad_type>().swap(dependent_variables);
433 taped_driver.reset(clear_registered_tapes);
435 independent_variable_values = std::vector<scalar_type>(
436 new_n_independent_variables,
438 registered_independent_variable_values =
439 std::vector<bool>(new_n_independent_variables,
false);
440 registered_marked_independent_variables =
441 std::vector<bool>(new_n_independent_variables,
false);
442 dependent_variables =
443 std::vector<ad_type>(new_n_dependent_variables,
445 registered_marked_dependent_variables =
446 std::vector<bool>(new_n_dependent_variables,
false);
451 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
454 const unsigned int n_independent_variables,
455 const bool ensure_persistent_setting)
462 n_independent_variables);
464 if (ensure_persistent_setting ==
true)
477 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
482 activate_tape(tape_index,
true );
487 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
495 return taped_driver.requires_retaping(tape_index);
500 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
508 return taped_driver.last_action_requires_retaping();
513 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
520 taped_driver.remove_tape(taped_driver.active_tape_index());
525 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
529 const bool read_mode)
536 ExcMessage(
"Tape index exceeds maximum allowable value"));
537 taped_driver.activate_tape(tape_index);
538 reset_registered_independent_variables();
543 if (read_mode ==
true)
545 Assert(is_registered_tape(tape_index),
547 reset_registered_dependent_variables(
true);
548 Assert(n_registered_dependent_variables() ==
549 n_dependent_variables(),
550 ExcMessage(
"Not all dependent variables have been set!"));
562 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
573 taped_driver.set_tape_buffer_sizes(obufsize,
581 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
585 const bool overwrite_tape,
586 const bool keep_independent_values)
589 const bool read_mode =
false;
593 if (overwrite_tape !=
true)
595 Assert(is_recording() ==
false,
600 if (is_registered_tape(tape_index) ==
false || overwrite_tape ==
true)
604 activate_tape(tape_index, read_mode);
607 taped_driver.start_taping(active_tape_index(),
608 keep_independent_values);
612 reset_registered_independent_variables();
613 reset_registered_dependent_variables();
617 Assert(is_recording() ==
false,
619 "Tape recording is unexpectedly still enabled."));
623 activate_recorded_tape(tape_index);
628 Assert(ADNumberTraits<ad_type>::is_tapeless ==
true,
633 tapeless_driver.allow_dependent_variable_marking();
638 activate_tape(tape_index, read_mode);
641 return is_recording();
646 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
649 const bool write_tapes_to_file)
654 Assert(n_registered_independent_variables() == n_independent_variables(),
655 ExcMessage(
"Not all values of sensitivities have been recorded!"));
660 taped_driver.stop_taping(active_tape_index(), write_tapes_to_file);
667 Assert(n_registered_dependent_variables() == n_dependent_variables(),
668 ExcMessage(
"Not all dependent variables have been set!"));
676 tapeless_driver.prevent_dependent_variable_marking();
682 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
685 const unsigned int index,
689 Assert(registered_marked_dependent_variables[index] ==
false,
691 "This dependent variable has already been registered."));
697 Assert(is_recording() ==
true,
699 "Must be recording when registering dependent variables."));
705 registered_marked_dependent_variables[index] =
true;
714 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
716 const unsigned int n_independent_variables,
717 const unsigned int n_dependent_variables)
718 :
HelperBase<ADNumberTypeCode, ScalarType>(n_independent_variables,
719 n_dependent_variables)
724 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
727 const std::vector<scalar_type> &dof_values)
733 Assert(dof_values.size() == this->n_independent_variables(),
735 "Vector size does not match number of independent variables"));
736 for (
unsigned int i = 0; i < this->n_independent_variables(); ++i)
738 Assert(this->registered_independent_variable_values[i] ==
false,
739 ExcMessage(
"Independent variable value already registered."));
741 set_dof_values(dof_values);
746 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
754 Assert(this->active_tape_index() !=
761 this->finalize_sensitive_independent_variables();
762 Assert(this->independent_variables.size() ==
763 this->n_independent_variables(),
765 this->n_independent_variables()));
767 return this->independent_variables;
772 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
775 const std::vector<scalar_type> &
values)
779 Assert(this->active_tape_index() !=
783 Assert(
values.size() == this->n_independent_variables(),
785 "Vector size does not match number of independent variables"));
786 for (
unsigned int i = 0; i < this->n_independent_variables(); ++i)
797 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
799 const unsigned int n_independent_variables)
800 :
CellLevelBase<ADNumberTypeCode, ScalarType>(n_independent_variables, 1)
805 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
817 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
822 this->taped_driver.keep_independent_values() ==
false) ||
826 this->n_registered_independent_variables() ==
827 this->n_independent_variables(),
829 "Not all values of sensitivities have been registered or subsequently set!"));
831 Assert(this->n_registered_dependent_variables() ==
832 this->n_dependent_variables(),
833 ExcMessage(
"Not all dependent variables have been registered."));
836 this->n_dependent_variables() == 1,
838 "The EnergyFunctional class expects there to be only one dependent variable."));
842 Assert(this->active_tape_index() !=
845 Assert(this->is_recording() ==
false,
847 "Cannot compute value while tape is being recorded."));
848 Assert(this->independent_variable_values.size() ==
849 this->n_independent_variables(),
851 this->n_independent_variables()));
853 return this->taped_driver.value(this->active_tape_index(),
854 this->independent_variable_values);
860 Assert(this->independent_variables.size() ==
861 this->n_independent_variables(),
863 this->n_independent_variables()));
865 return this->tapeless_driver.value(this->dependent_variables);
871 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
877 this->taped_driver.keep_independent_values() ==
false) ||
881 this->n_registered_independent_variables() ==
882 this->n_independent_variables(),
884 "Not all values of sensitivities have been registered or subsequently set!"));
886 Assert(this->n_registered_dependent_variables() ==
887 this->n_dependent_variables(),
888 ExcMessage(
"Not all dependent variables have been registered."));
891 this->n_dependent_variables() == 1,
893 "The EnergyFunctional class expects there to be only one dependent variable."));
898 if (gradient.size() != this->n_independent_variables())
899 gradient.reinit(this->n_independent_variables(),
904 Assert(this->active_tape_index() !=
907 Assert(this->is_recording() ==
false,
909 "Cannot compute gradient while tape is being recorded."));
910 Assert(this->independent_variable_values.size() ==
911 this->n_independent_variables(),
913 this->n_independent_variables()));
915 this->taped_driver.gradient(this->active_tape_index(),
916 this->independent_variable_values,
923 Assert(this->independent_variables.size() ==
924 this->n_independent_variables(),
926 this->n_independent_variables()));
928 this->tapeless_driver.gradient(this->independent_variables,
929 this->dependent_variables,
936 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
943 "Cannot computed function Hessian: AD number type does "
944 "not support the calculation of second order derivatives."));
947 this->taped_driver.keep_independent_values() ==
false))
950 this->n_registered_independent_variables() ==
951 this->n_independent_variables(),
953 "Not all values of sensitivities have been registered or subsequently set!"));
955 Assert(this->n_registered_dependent_variables() ==
956 this->n_dependent_variables(),
957 ExcMessage(
"Not all dependent variables have been registered."));
960 this->n_dependent_variables() == 1,
962 "The EnergyFunctional class expects there to be only one dependent variable."));
967 if (hessian.m() != this->n_independent_variables() ||
968 hessian.n() != this->n_independent_variables())
969 hessian.reinit({this->n_independent_variables(),
970 this->n_independent_variables()},
973 if (ADNumberTraits<ad_type>::is_taped ==
true)
975 Assert(this->active_tape_index() !=
978 Assert(this->is_recording() ==
false,
980 "Cannot compute hessian while tape is being recorded."));
981 Assert(this->independent_variable_values.size() ==
982 this->n_independent_variables(),
984 this->n_independent_variables()));
986 this->taped_driver.hessian(this->active_tape_index(),
987 this->independent_variable_values,
994 Assert(this->independent_variables.size() ==
995 this->n_independent_variables(),
997 this->n_independent_variables()));
999 this->tapeless_driver.hessian(this->independent_variables,
1000 this->dependent_variables,
1010 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
1012 const unsigned int n_independent_variables,
1013 const unsigned int n_dependent_variables)
1014 :
CellLevelBase<ADNumberTypeCode, ScalarType>(n_independent_variables,
1015 n_dependent_variables)
1020 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
1025 Assert(residual.size() == this->n_dependent_variables(),
1027 "Vector size does not match number of dependent variables"));
1028 for (
unsigned int i = 0; i < this->n_dependent_variables(); ++i)
1035 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
1041 this->taped_driver.keep_independent_values() ==
false) ||
1045 this->n_registered_independent_variables() ==
1046 this->n_independent_variables(),
1048 "Not all values of sensitivities have been registered or subsequently set!"));
1050 Assert(this->n_registered_dependent_variables() ==
1051 this->n_dependent_variables(),
1052 ExcMessage(
"Not all dependent variables have been registered."));
1057 if (
values.size() != this->n_dependent_variables())
1058 values.reinit(this->n_dependent_variables(),
1063 Assert(this->active_tape_index() !=
1066 Assert(this->is_recording() ==
false,
1068 "Cannot compute values while tape is being recorded."));
1069 Assert(this->independent_variable_values.size() ==
1070 this->n_independent_variables(),
1072 this->n_independent_variables()));
1074 this->taped_driver.values(this->active_tape_index(),
1075 this->n_dependent_variables(),
1076 this->independent_variable_values,
1083 this->tapeless_driver.values(this->dependent_variables,
values);
1089 template <enum AD::NumberTypes ADNumberTypeCode,
typename ScalarType>
1095 this->taped_driver.keep_independent_values() ==
false) ||
1099 this->n_registered_independent_variables() ==
1100 this->n_independent_variables(),
1102 "Not all values of sensitivities have been registered or subsequently set!"));
1104 Assert(this->n_registered_dependent_variables() ==
1105 this->n_dependent_variables(),
1106 ExcMessage(
"Not all dependent variables have been registered."));
1111 if (jacobian.
m() != this->n_dependent_variables() ||
1112 jacobian.
n() != this->n_independent_variables())
1113 jacobian.reinit({this->n_dependent_variables(),
1114 this->n_independent_variables()},
1119 Assert(this->active_tape_index() !=
1122 Assert(this->is_recording() ==
false,
1124 "Cannot compute hessian while tape is being recorded."));
1125 Assert(this->independent_variable_values.size() ==
1126 this->n_independent_variables(),
1128 this->n_independent_variables()));
1130 this->taped_driver.jacobian(this->active_tape_index(),
1131 this->n_dependent_variables(),
1132 this->independent_variable_values,
1139 Assert(this->independent_variables.size() ==
1140 this->n_independent_variables(),
1142 this->n_independent_variables()));
1144 this->tapeless_driver.jacobian(this->independent_variables,
1145 this->dependent_variables,
1158 typename ScalarType>
1161 const unsigned int n_dependent_variables)
1162 :
HelperBase<ADNumberTypeCode, ScalarType>(n_independent_variables,
1163 n_dependent_variables)
1164 , symmetric_independent_variables(n_independent_variables, false)
1171 typename ScalarType>
1174 const unsigned int n_independent_variables,
1175 const unsigned int n_dependent_variables,
1176 const bool clear_registered_tapes)
1179 n_dependent_variables,
1180 clear_registered_tapes);
1182 const unsigned int new_n_independent_variables =
1184 n_independent_variables :
1185 this->n_independent_variables());
1186 symmetric_independent_variables =
1187 std::vector<bool>(new_n_independent_variables,
false);
1194 typename ScalarType>
1199 Assert(index < symmetric_independent_variables.size(),
1201 return symmetric_independent_variables[index];
1208 typename ScalarType>
1213 return std::count(symmetric_independent_variables.begin(),
1214 symmetric_independent_variables.end(),
1222 typename ScalarType>
1231 Assert(
values.size() == this->n_independent_variables(),
1233 "Vector size does not match number of independent variables"));
1234 for (
unsigned int i = 0; i < this->n_independent_variables(); ++i)
1236 Assert(this->registered_independent_variable_values[i] ==
false,
1237 ExcMessage(
"Independent variable value already registered."));
1239 set_independent_variables(
values);
1246 typename ScalarType>
1249 ScalarType>::ad_type> &
1255 Assert(this->active_tape_index() !=
1263 this->finalize_sensitive_independent_variables();
1264 Assert(this->independent_variables.size() ==
1265 this->n_independent_variables(),
1267 this->n_independent_variables()));
1269 return this->independent_variables;
1276 typename ScalarType>
1280 const bool symmetric_component,
1286 index < this->n_independent_variables(),
1288 "Trying to set the symmetry flag of a non-existent independent variable."));
1289 Assert(index < symmetric_independent_variables.size(),
1291 symmetric_independent_variables[index] = symmetric_component;
1298 typename ScalarType>
1305 Assert(this->active_tape_index() !=
1309 Assert(
values.size() == this->n_independent_variables(),
1311 "Vector size does not match number of independent variables"));
1312 for (
unsigned int i = 0; i < this->n_independent_variables(); ++i)
1325 typename ScalarType>
1327 const unsigned int n_independent_variables)
1329 n_independent_variables,
1337 typename ScalarType>
1351 typename ScalarType>
1356 this->taped_driver.keep_independent_values() ==
false) ||
1360 this->n_registered_independent_variables() ==
1361 this->n_independent_variables(),
1363 "Not all values of sensitivities have been registered or subsequently set!"));
1365 Assert(this->n_registered_dependent_variables() ==
1366 this->n_dependent_variables(),
1367 ExcMessage(
"Not all dependent variables have been registered."));
1370 this->n_dependent_variables() == 1,
1372 "The ScalarFunction class expects there to be only one dependent variable."));
1376 Assert(this->active_tape_index() !=
1379 Assert(this->is_recording() ==
false,
1381 "Cannot compute values while tape is being recorded."));
1382 Assert(this->independent_variable_values.size() ==
1383 this->n_independent_variables(),
1385 this->n_independent_variables()));
1387 return this->taped_driver.value(this->active_tape_index(),
1388 this->independent_variable_values);
1394 return this->tapeless_driver.value(this->dependent_variables);
1401 typename ScalarType>
1407 this->taped_driver.keep_independent_values() ==
false) ||
1411 this->n_registered_independent_variables() ==
1412 this->n_independent_variables(),
1414 "Not all values of sensitivities have been registered or subsequently set!"));
1416 Assert(this->n_registered_dependent_variables() ==
1417 this->n_dependent_variables(),
1418 ExcMessage(
"Not all dependent variables have been registered."));
1421 this->n_dependent_variables() == 1,
1423 "The ScalarFunction class expects there to be only one dependent variable."));
1428 if (gradient.size() != this->n_independent_variables())
1429 gradient.reinit(this->n_independent_variables(),
1434 Assert(this->active_tape_index() !=
1437 Assert(this->is_recording() ==
false,
1439 "Cannot compute gradient while tape is being recorded."));
1440 Assert(this->independent_variable_values.size() ==
1441 this->n_independent_variables(),
1443 this->n_independent_variables()));
1445 this->taped_driver.gradient(this->active_tape_index(),
1446 this->independent_variable_values,
1453 Assert(this->independent_variables.size() ==
1454 this->n_independent_variables(),
1456 this->n_independent_variables()));
1458 this->tapeless_driver.gradient(this->independent_variables,
1459 this->dependent_variables,
1464 for (
unsigned int i = 0; i < this->n_independent_variables(); i++)
1466 if (this->is_symmetric_independent_variable(i) ==
true)
1475 typename ScalarType>
1482 "Cannot computed function Hessian: AD number type does "
1483 "not support the calculation of second order derivatives."));
1486 this->taped_driver.keep_independent_values() ==
false))
1489 this->n_registered_independent_variables() ==
1490 this->n_independent_variables(),
1492 "Not all values of sensitivities have been registered or subsequently set!"));
1494 Assert(this->n_registered_dependent_variables() ==
1495 this->n_dependent_variables(),
1496 ExcMessage(
"Not all dependent variables have been registered."));
1499 this->n_dependent_variables() == 1,
1501 "The ScalarFunction class expects there to be only one dependent variable."));
1506 if (hessian.m() != this->n_independent_variables() ||
1507 hessian.n() != this->n_independent_variables())
1508 hessian.reinit({this->n_independent_variables(),
1509 this->n_independent_variables()},
1514 Assert(this->active_tape_index() !=
1517 Assert(this->is_recording() ==
false,
1519 "Cannot compute Hessian while tape is being recorded."));
1520 Assert(this->independent_variable_values.size() ==
1521 this->n_independent_variables(),
1523 this->n_independent_variables()));
1525 this->taped_driver.hessian(this->active_tape_index(),
1526 this->independent_variable_values,
1533 Assert(this->independent_variables.size() ==
1534 this->n_independent_variables(),
1536 this->n_independent_variables()));
1538 this->tapeless_driver.hessian(this->independent_variables,
1539 this->dependent_variables,
1544 for (
unsigned int i = 0; i < this->n_independent_variables(); i++)
1545 for (
unsigned int j = 0; j < i + 1; j++)
1547 if (this->is_symmetric_independent_variable(i) ==
true &&
1548 this->is_symmetric_independent_variable(j) ==
true)
1550 hessian[i][j] *= 0.25;
1552 hessian[j][i] *= 0.25;
1554 else if ((this->is_symmetric_independent_variable(i) ==
true &&
1555 this->is_symmetric_independent_variable(j) ==
false) ||
1556 (this->is_symmetric_independent_variable(j) ==
true &&
1557 this->is_symmetric_independent_variable(i) ==
false))
1559 hessian[i][j] *= 0.5;
1561 hessian[j][i] *= 0.5;
1570 typename ScalarType>
1589 const std::vector<unsigned int> row_index_set(
1590 internal::extract_field_component_indices<dim>(extractor_row));
1591 const std::vector<unsigned int> col_index_set(
1592 internal::extract_field_component_indices<dim>(extractor_col));
1598 hessian[row_index_set[0]][col_index_set[0]]);
1607 typename ScalarType>
1627 const std::vector<unsigned int> row_index_set(
1628 internal::extract_field_component_indices<dim>(extractor_row));
1629 const std::vector<unsigned int> col_index_set(
1630 internal::extract_field_component_indices<dim>(extractor_col));
1632 for (
unsigned int r = 0; r < row_index_set.size(); ++r)
1633 for (
unsigned int c = 0; c < col_index_set.size(); ++c)
1636 out, r, c, hessian[row_index_set[r]][col_index_set[c]]);
1650 typename ScalarType>
1652 const unsigned int n_independent_variables,
1653 const unsigned int n_dependent_variables)
1655 n_independent_variables,
1656 n_dependent_variables)
1663 typename ScalarType>
1668 Assert(funcs.size() == this->n_dependent_variables(),
1670 "Vector size does not match number of dependent variables"));
1671 for (
unsigned int i = 0; i < this->n_dependent_variables(); ++i)
1680 typename ScalarType>
1686 this->taped_driver.keep_independent_values() ==
false) ||
1690 this->n_registered_independent_variables() ==
1691 this->n_independent_variables(),
1693 "Not all values of sensitivities have been registered or subsequently set!"));
1695 Assert(this->n_registered_dependent_variables() ==
1696 this->n_dependent_variables(),
1697 ExcMessage(
"Not all dependent variables have been registered."));
1702 if (
values.size() != this->n_dependent_variables())
1703 values.reinit(this->n_dependent_variables(),
1708 Assert(this->active_tape_index() !=
1711 Assert(this->is_recording() ==
false,
1713 "Cannot compute values while tape is being recorded."));
1714 Assert(this->independent_variable_values.size() ==
1715 this->n_independent_variables(),
1717 this->n_independent_variables()));
1719 this->taped_driver.values(this->active_tape_index(),
1720 this->n_dependent_variables(),
1721 this->independent_variable_values,
1728 this->tapeless_driver.values(this->dependent_variables,
values);
1736 typename ScalarType>
1742 this->taped_driver.keep_independent_values() ==
false) ||
1746 this->n_registered_independent_variables() ==
1747 this->n_independent_variables(),
1749 "Not all values of sensitivities have been registered or subsequently set!"));
1751 Assert(this->n_registered_dependent_variables() ==
1752 this->n_dependent_variables(),
1753 ExcMessage(
"Not all dependent variables have been registered."));
1758 if (jacobian.
m() != this->n_dependent_variables() ||
1759 jacobian.
n() != this->n_independent_variables())
1760 jacobian.reinit({this->n_dependent_variables(),
1761 this->n_independent_variables()},
1766 Assert(this->active_tape_index() !=
1769 Assert(this->is_recording() ==
false,
1771 "Cannot compute Jacobian while tape is being recorded."));
1772 Assert(this->independent_variable_values.size() ==
1773 this->n_independent_variables(),
1775 this->n_independent_variables()));
1777 this->taped_driver.jacobian(this->active_tape_index(),
1778 this->n_dependent_variables(),
1779 this->independent_variable_values,
1786 Assert(this->independent_variables.size() ==
1787 this->n_independent_variables(),
1789 this->n_independent_variables()));
1791 this->tapeless_driver.jacobian(this->independent_variables,
1792 this->dependent_variables,
1796 for (
unsigned int j = 0; j < this->n_independent_variables(); j++)
1801 if (this->is_symmetric_independent_variable(j) ==
true)
1802 for (
unsigned int i = 0; i < this->n_dependent_variables(); i++)
1803 jacobian[i][j] *= 0.5;
1811 typename ScalarType>
1831 const std::vector<unsigned int> row_index_set(
1832 internal::extract_field_component_indices<dim>(extractor_row));
1833 const std::vector<unsigned int> col_index_set(
1834 internal::extract_field_component_indices<dim>(extractor_col));
1840 jacobian[row_index_set[0]][col_index_set[0]]);
1849 typename ScalarType>
1869 const std::vector<unsigned int> row_index_set(
1870 internal::extract_field_component_indices<dim>(extractor_row));
1871 const std::vector<unsigned int> col_index_set(
1872 internal::extract_field_component_indices<dim>(extractor_col));
1874 for (
unsigned int r = 0; r < row_index_set.size(); ++r)
1875 for (
unsigned int c = 0; c < col_index_set.size(); ++c)
1878 out, r, c, jacobian[row_index_set[r]][col_index_set[c]]);
1890# include "ad_helpers.inst"
1892# ifdef DEAL_II_WITH_ADOLC
1893# include "ad_helpers.inst1"
1895# ifdef DEAL_II_TRILINOS_WITH_SACADO
1896# include "ad_helpers.inst2"
const std::vector< ad_type > & get_sensitive_dof_values() const
typename HelperBase< ADNumberTypeCode, ScalarType >::ad_type ad_type
void set_dof_values(const std::vector< scalar_type > &dof_values)
void register_dof_values(const std::vector< scalar_type > &dof_values)
CellLevelBase(const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
void register_energy_functional(const ad_type &energy)
EnergyFunctional(const unsigned int n_independent_variables)
typename HelperBase< ADNumberTypeCode, ScalarType >::ad_type ad_type
void compute_residual(Vector< scalar_type > &residual) const override
virtual void compute_linearization(FullMatrix< scalar_type > &linearization) const override
scalar_type compute_energy() const
typename HelperBase< ADNumberTypeCode, ScalarType >::scalar_type scalar_type
unsigned int n_registered_dependent_variables() const
void stop_recording_operations(const bool write_tapes_to_file=false)
typename AD::NumberTraits< ScalarType, ADNumberTypeCode >::scalar_type scalar_type
Types< ad_type >::tape_index active_tape_index() const
void print_tape_stats(const typename Types< ad_type >::tape_index tape_index, std::ostream &stream) const
std::size_t n_independent_variables() const
bool start_recording_operations(const typename Types< ad_type >::tape_index tape_index, const bool overwrite_tape=false, const bool keep_independent_values=true)
void mark_independent_variable(const unsigned int index, ad_type &out) const
bool is_registered_tape(const typename Types< ad_type >::tape_index tape_index) const
void reset_registered_dependent_variables(const bool flag=false)
void activate_tape(const typename Types< ad_type >::tape_index tape_index, const bool read_mode)
void set_tape_buffer_sizes(const typename Types< ad_type >::tape_buffer_sizes obufsize=64 *1024 *1024, const typename Types< ad_type >::tape_buffer_sizes lbufsize=64 *1024 *1024, const typename Types< ad_type >::tape_buffer_sizes vbufsize=64 *1024 *1024, const typename Types< ad_type >::tape_buffer_sizes tbufsize=64 *1024 *1024)
void print(std::ostream &stream) const
void set_sensitivity_value(const unsigned int index, const scalar_type &value)
bool active_tape_requires_retaping() const
HelperBase(const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
virtual void reset(const unsigned int n_independent_variables=::numbers::invalid_unsigned_int, const unsigned int n_dependent_variables=::numbers::invalid_unsigned_int, const bool clear_registered_tapes=true)
bool recorded_tape_requires_retaping(const typename Types< ad_type >::tape_index tape_index) const
std::vector< ad_type > dependent_variables
bool is_recording() const
std::size_t n_dependent_variables() const
void print_values(std::ostream &stream) const
void register_dependent_variable(const unsigned int index, const ad_type &func)
static void configure_tapeless_mode(const unsigned int n_independent_variables, const bool ensure_persistent_setting=true)
void initialize_non_sensitive_independent_variable(const unsigned int index, ad_type &out) const
typename AD::NumberTraits< ScalarType, ADNumberTypeCode >::ad_type ad_type
void activate_recorded_tape(const typename Types< ad_type >::tape_index tape_index)
bool is_symmetric_independent_variable(const unsigned int index) const
PointLevelFunctionsBase(const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
void set_independent_variables(const std::vector< scalar_type > &values)
typename HelperBase< ADNumberTypeCode, ScalarType >::scalar_type scalar_type
unsigned int n_symmetric_independent_variables() const
void register_independent_variables(const std::vector< scalar_type > &values)
void set_sensitivity_value(const unsigned int index, const bool symmetric_component, const scalar_type &value)
const std::vector< ad_type > & get_sensitive_variables() const
virtual void reset(const unsigned int n_independent_variables=::numbers::invalid_unsigned_int, const unsigned int n_dependent_variables=::numbers::invalid_unsigned_int, const bool clear_registered_tapes=true) override
ResidualLinearization(const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
virtual void compute_residual(Vector< scalar_type > &residual) const override
virtual void compute_linearization(FullMatrix< scalar_type > &linearization) const override
void register_residual_vector(const std::vector< ad_type > &residual)
static internal::ScalarFieldHessian< dim, scalar_type, ExtractorType_Row, ExtractorType_Col >::type extract_hessian_component(const FullMatrix< scalar_type > &hessian, const ExtractorType_Row &extractor_row, const ExtractorType_Col &extractor_col)
typename HelperBase< ADNumberTypeCode, ScalarType >::ad_type ad_type
void register_dependent_variable(const ad_type &func)
ScalarFunction(const unsigned int n_independent_variables)
void compute_gradient(Vector< scalar_type > &gradient) const
scalar_type compute_value() const
void compute_hessian(FullMatrix< scalar_type > &hessian) const
typename HelperBase< ADNumberTypeCode, ScalarType >::scalar_type scalar_type
void compute_values(Vector< scalar_type > &values) const
void register_dependent_variables(const std::vector< ad_type > &funcs)
void compute_jacobian(FullMatrix< scalar_type > &jacobian) const
typename HelperBase< ADNumberTypeCode, ScalarType >::scalar_type scalar_type
static internal::VectorFieldJacobian< dim, scalar_type, ExtractorType_Row, ExtractorType_Col >::type extract_jacobian_component(const FullMatrix< scalar_type > &jacobian, const ExtractorType_Row &extractor_row, const ExtractorType_Col &extractor_col)
VectorFunction(const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define Assert(cond, exc)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcMessage(std::string arg1)
std::enable_if<!(ADNumberTraits< ADNumberType >::type_code==NumberTypes::adolc_tapeless)>::type configure_tapeless_mode(const unsigned int)
void set_tensor_entry(TensorType &t, const unsigned int unrolled_index, const NumberType &value)
static const unsigned int invalid_unsigned_int
static const Types< ADNumberType >::tape_index invalid_tape_index
static void initialize_global_environment(const unsigned int n_independent_variables)