From 80982cde38774837ee39bfdf740319e4aa094d0b Mon Sep 17 00:00:00 2001 From: Jean-Paul Pelteret Date: Fri, 5 Oct 2018 16:05:03 +0200 Subject: [PATCH] Rework implementation of taping within ADHelperBase. This commit abstracts the notion of taping almost completely away from the ADHelperBase class. The specifics of the taping process are now fully contained within the driver classes (at the expense of having to have some non-static interfaces to the driver classes). However, this should better enable further extension to use other taped AD implementations in the future. --- .../deal.II/differentiation/ad/ad_drivers.h | 706 ++++++++++++++++-- .../deal.II/differentiation/ad/ad_helpers.h | 95 +-- source/differentiation/ad/ad_helpers.cc | 220 ++---- 3 files changed, 722 insertions(+), 299 deletions(-) diff --git a/include/deal.II/differentiation/ad/ad_drivers.h b/include/deal.II/differentiation/ad/ad_drivers.h index cef0c87d65..4d58a94b9f 100644 --- a/include/deal.II/differentiation/ad/ad_drivers.h +++ b/include/deal.II/differentiation/ad/ad_drivers.h @@ -1,6 +1,6 @@ // --------------------------------------------------------------------- // -// Copyright (C) 2017 by the deal.II authors +// Copyright (C) 2017 - 2018 by the deal.II authors // // This file is part of the deal.II library. // @@ -109,6 +109,7 @@ namespace Differentiation using tape_buffer_sizes = unsigned int; } // namespace types + /** * A collection of special numbers used within the context of * auto-differentiable numbers. @@ -142,6 +143,7 @@ namespace Differentiation #endif // DEAL_II_WITH_ADOLC } // namespace numbers + /** * A driver class for taped auto-differentiable numbers. * @@ -172,37 +174,75 @@ namespace Differentiation //@{ /** - * Enable the recording mode for a given tape. + * Returns whether or not this class is tracking calculations performed + * with its marked independent variables. + */ + bool + is_recording() const; + + /** + * Returns the tape number which is currently activated for recording or + * reading. + */ + types::tape_index + active_tape_index() const; + + /** + * Returns whether or not a tape number has already been used + * or registered. + */ + bool + is_registered_tape(const types::tape_index tape_index) const; + + /** + * Returns whether or not the numerical values of all independent + * variables are recorded in the tape buffer. + */ + bool + keep_independent_values() const; + + /** + * Set the buffer sizes for the next active tape. * - * @param[in] tape_index The index of the tape to be written. - * @param[in] keep_independent_values Determines whether the numerical - * values of all independent variables are recorded in the - * tape buffer. + * This function must be called before start_recording_operations() for + * it to have any influence on the memory allocated to the next recorded + * tape. + * + * @note This function only has an effect when using ADOL-C numbers. As + * stated by the ADOL-C manual, it may be desirable to create a file + * ".adolcrc" in the program run directory and set the buffer size + * therein. + * Alternatively, this function can be used to override the settings for + * any given tape, or can be used in the event that no ".adolcrc" file + * exists. + * The default value for each buffer is set at 64 megabytes, a + * heuristically chosen value thought to be appropriate for use within + * the context of finite element analysis when considering coupled + * problems with multiple vector-valued fields discretised by higher + * order shape functions, as well as complex constitutive laws. + * + * @param[in] obufsize ADOL-C operations buffer size + * @param[in] lbufsize ADOL-C locations buffer size + * @param[in] vbufsize ADOL-C value buffer size + * @param[in] tbufsize ADOL-C Taylor buffer size */ - static void - start_taping(const types::tape_index tape_index, - const bool keep_independent_values); + void + set_tape_buffer_sizes(const types::tape_buffer_sizes obufsize = 67108864, + const types::tape_buffer_sizes lbufsize = 67108864, + const types::tape_buffer_sizes vbufsize = 67108864, + const types::tape_buffer_sizes tbufsize = 67108864); /** - * Enable the recording mode for a given tape, using the run-time - * chosen tape buffer sizes. + * Enable the recording mode for a given tape. * * @param[in] tape_index The index of the tape to be written. * @param[in] keep_independent_values Determines whether the numerical * values of all independent variables are recorded in the * tape buffer. - * @param[in] obufsize ADOL-C operations buffer size - * @param[in] lbufsize ADOL-C locations buffer size - * @param[in] vbufsize ADOL-C value buffer size - * @param[in] tbufsize ADOL-C Taylor buffer size */ - static void - start_taping(const types::tape_index tape_index, - const bool keep_independent_values, - const types::tape_buffer_sizes obufsize, - const types::tape_buffer_sizes lbufsize, - const types::tape_buffer_sizes vbufsize, - const types::tape_buffer_sizes tbufsize); + void + start_taping(const types::tape_index tape_index, + const bool keep_independent_values); /** * Disable the recording mode for a given tape. @@ -212,10 +252,46 @@ namespace Differentiation * @param[in] write_tapes_to_file A flag that specified whether the tape * should be written to file or kept in memory. */ - static void + void stop_taping(const types::tape_index active_tape_index, const bool write_tapes_to_file); + /** + * Select a tape to record to or read from. + * + * This function activates a tape, but depending on whether @p read_mode + * is set, the tape is either taken as previously written to (and put + * into read-only mode), or cleared for (re-)taping. + * + * @param[in] tape_index The index of the tape to be written to/read + * from. + * + * @note The chosen tape index must be greater than + * numbers::invalid_tape_index and less than numbers::max_tape_index. + */ + void + activate_tape(const types::tape_index tape_index); + + /** + * Reset the state of the class. + * + * @note This also resets the active tape number to an invalid number, and + * deactivates the recording mode for taped variables. + */ + void + reset(const bool clear_registered_tapes); + + /** + * Prints the status of all queriable data. Exactly what is printed and + * its format depends on the @p ad_type, as is determined by the + * @p ADNumberTypeCode template parameter. + * + * @param[in] stream The output stream to which the values are to be + * written. + */ + void + print(std::ostream &stream) const; + /** * Print the statistics regarding the usage of the tapes. * @@ -389,6 +465,34 @@ namespace Differentiation //@} + /** + * Operation status + */ + //@{ + + /** + * Set a flag that states that we can safely mark dependent variables + * within the current phase of operations. + */ + void + allow_dependent_variable_marking(); + + /** + * Set a flag that states that we cannot safely mark dependent variables + * within the current phase of operations. + */ + void + prevent_dependent_variable_marking(); + + /** + * Query a flag as to whether or not dependent variables can be marked + * within the current phase of operations. + */ + bool + is_dependent_variable_marking_allowed() const; + + //@} + /** * @name Drivers for scalar functions */ @@ -504,20 +608,45 @@ namespace Differentiation template - void - TapedDrivers::start_taping( - const types::tape_index, - const bool) + bool + TapedDrivers::is_recording() const { AssertThrow(false, ExcRequiresADNumberSpecialization()); + return false; + } + + + template + types::tape_index + TapedDrivers::active_tape_index() const + { + AssertThrow(false, ExcRequiresADNumberSpecialization()); + return numbers::invalid_tape_index; + } + + + template + bool + TapedDrivers::is_registered_tape( + const types::tape_index) const + { + AssertThrow(false, ExcRequiresADNumberSpecialization()); + return false; + } + + + template + bool + TapedDrivers::keep_independent_values() const + { + AssertThrow(false, ExcRequiresADNumberSpecialization()); + return false; } template void - TapedDrivers::start_taping( - const types::tape_index, - const bool, + TapedDrivers::set_tape_buffer_sizes( const types::tape_buffer_sizes, const types::tape_buffer_sizes, const types::tape_buffer_sizes, @@ -527,6 +656,16 @@ namespace Differentiation } + template + void + TapedDrivers::start_taping( + const types::tape_index, + const bool) + { + AssertThrow(false, ExcRequiresADNumberSpecialization()); + } + + template void TapedDrivers::stop_taping( @@ -537,6 +676,31 @@ namespace Differentiation } + template + void + TapedDrivers::activate_tape( + const types::tape_index) + { + AssertThrow(false, ExcRequiresADNumberSpecialization()); + } + + + template + void + TapedDrivers::reset(const bool) + { + AssertThrow(false, ExcRequiresADNumberSpecialization()); + } + + + template + void + TapedDrivers::print(std::ostream &) const + { + AssertThrow(false, ExcRequiresADNumberSpecialization()); + } + + template void TapedDrivers::print_tape_stats( @@ -622,32 +786,120 @@ namespace Differentiation { using scalar_type = double; + TapedDrivers() + : active_tape(numbers::invalid_tape_index) + , keep_values(true) + , is_recording_flag(false) + , use_stored_taped_buffer_sizes(false) + , obufsize(0u) + , lbufsize(0u) + , vbufsize(0u) + , tbufsize(0u) + {} + // === Taping === - static void - start_taping(const types::tape_index tape_index, - const bool keep_independent_values) + bool + is_recording() const { - trace_on(tape_index, keep_independent_values); + return is_recording_flag; } - static void - start_taping(const types::tape_index tape_index, - const bool keep_independent_values, - const types::tape_buffer_sizes obufsize, - const types::tape_buffer_sizes lbufsize, - const types::tape_buffer_sizes vbufsize, - const types::tape_buffer_sizes tbufsize) + types::tape_index + active_tape_index() const { - trace_on(tape_index, - keep_independent_values, - obufsize, - lbufsize, - vbufsize, - tbufsize); + return active_tape; } - static void + bool + keep_independent_values() const + { + return keep_values; + } + + bool + is_registered_tape(const types::tape_index tape_index) const + { + // Sigh... This is a mess :-/ + // The most succinct way to get this piece of information, would be to + // use the getTapeInfos() function, but would come at the expense of + // creating an inactive tape data within ADOL-C's global store. For + // hints as to why this is the way it is, see the getTapeInfos() + // function in + // https://gitlab.com/adol-c/adol-c/blob/master/ADOL-C/src/tape_handling.cpp + // An alternative solution would be to manually access the tape data; + // see the removeTape() function in + // https://gitlab.com/adol-c/adol-c/blob/master/ADOL-C/src/tape_handling.cpp + // with the consideration of the #defines in + // https://gitlab.com/adol-c/adol-c/blob/master/ADOL-C/src/taping_p.h + // as to how this would be performed. + // Doing things "manually" (the second way) without creating the + // additional data object would be a lot more work... + // + // Both of the above solutions would be possible IF ADOL-C exposed + // this data object or a method to access it to the outside world. + // But they don't :-( + // Instead, what we'll have to do is get the statistics for this tape + // and make our own determination as to whether or not this tape exists. + // This effectively executes the first solution, with even more + // overhead! If the tape is in existence, we can assume that it should a + // non-zero number of dependent and independent variables. Those are + // stored in the zeroth and first entries of the statistics vector. + // + // But, oh wait... Surprise! This will trigger an error if the tape + // doesn't exist at all! So lets first check their tape info cache to + // see if the tape REALLY exists (i.e. has been touched, even if nothing + // has been written to it) before trying to access it. It'll only take + // an O(n) search, but at this point who really cares about efficiency? + const std::vector registered_tape_indices = + get_registered_tape_indices(); + const auto it = std::find(registered_tape_indices.begin(), + registered_tape_indices.end(), + tape_index); + if (it == registered_tape_indices.end()) + return false; + + std::vector counts(STAT_SIZE); + ::tapestats(tape_index, counts.data()); + Assert(counts.size() >= 2, ExcInternalError()); + return counts.data()[0] > 0 && counts.data()[1] > 0; + } + + void + set_tape_buffer_sizes(const types::tape_buffer_sizes in_obufsize, + const types::tape_buffer_sizes in_lbufsize, + const types::tape_buffer_sizes in_vbufsize, + const types::tape_buffer_sizes in_tbufsize) + { + // When valid for the chosen AD number type, these values will be used + // the next time start_recording_operations() is called. + obufsize = in_obufsize; + lbufsize = in_lbufsize; + vbufsize = in_vbufsize; + tbufsize = in_tbufsize; + use_stored_taped_buffer_sizes = true; + } + + void + start_taping(const types::tape_index tape_index, + const bool keep_independent_values) + { + if (use_stored_taped_buffer_sizes) + trace_on(tape_index, + keep_independent_values, + obufsize, + lbufsize, + vbufsize, + tbufsize); + else + trace_on(tape_index, keep_independent_values); + + // Set some other flags to their indicated / required values + keep_values = keep_independent_values; + is_recording_flag = true; + } + + void stop_taping(const types::tape_index active_tape_index, const bool write_tapes_to_file) { @@ -655,6 +907,55 @@ namespace Differentiation trace_off(active_tape_index); // Slow else trace_off(); // Fast(er) + + // Now that we've turned tracing off, we've definitely + // stopped all tape recording. + is_recording_flag = false; + + // If the keep_values flag is set, then we expect the user to use this + // tape immediately after recording it. There is therefore no need to + // invalidate it. However, there is now also no way to double-check + // that the newly recorded tape is indeed the active tape. + if (keep_independent_values() == false) + active_tape = numbers::invalid_tape_index; + } + + void + activate_tape(const types::tape_index tape_index) + { + active_tape = tape_index; + } + + void + reset(const bool clear_registered_tapes) + { + active_tape = numbers::invalid_tape_index; + is_recording_flag = false; + if (clear_registered_tapes) + { + const std::vector registered_tape_indices = + get_registered_tape_indices(); + for (const auto &tape_index : registered_tape_indices) + removeTape(tape_index, TapeRemovalType::ADOLC_REMOVE_COMPLETELY); + } + } + + void + print(std::ostream &stream) const + { + const std::vector registered_tape_indices = + get_registered_tape_indices(); + stream << "Registered tapes: "; + auto it_registered_tape = registered_tape_indices.begin(); + for (unsigned int i = 0; i < registered_tape_indices.size(); + ++i, ++it_registered_tape) + stream << *it_registered_tape + << (i < (registered_tape_indices.size() - 1) ? "," : ""); + stream << "\n"; + + stream << "Keep values? " << keep_independent_values() << "\n"; + stream << "Use stored tape buffer sizes? " + << use_stored_taped_buffer_sizes << "\n"; } static void @@ -815,6 +1116,76 @@ namespace Differentiation independent_variables.data(), J.data()); } + + protected: + /** + * Index of the tape that is currently in use. It is this tape that will + * be recorded to or read from when performing computations using "taped" + * auto-differentiable numbers. + */ + types::tape_index active_tape; + + /** + * Mark whether we're going to inform taped data structures to retain + * the coefficients ("Taylors" in ADOL-C nomenclature) stored on the + * tape so that they can be evaluated again at a later stage. + */ + bool keep_values; + + /** + * Mark whether we're currently recording a tape. Dependent on the state + * of this flag, only a restricted set of operations are allowable. + */ + bool is_recording_flag; + + /** + * A flag indicating that we should preferentially use the user-defined + * taped buffer sizes as opposed to either the default values selected + * by the AD library (or, in the case of ADOL-C, defined in an + * ".adolcrc" file). + */ + bool use_stored_taped_buffer_sizes; + + /** + * ADOL-C operations buffer size. + */ + types::tape_buffer_sizes obufsize; + + /** + * ADOL-C locations buffer size. + */ + types::tape_buffer_sizes lbufsize; + + /** + * ADOL-C value buffer size. + */ + types::tape_buffer_sizes vbufsize; + + /** + * ADOL-C Taylor buffer size. + */ + types::tape_buffer_sizes tbufsize; + + /** + * Return a list of registered tape indices + */ + std::vector + get_registered_tape_indices() const + { + // We've chosen to use unsigned shorts for the tape + // index type (a safety precaution) so we need to + // perform a conversion betwwen ADOL-C's native tape + // index type and that chosen by us. + std::vector registered_tape_indices_s; + cachedTraceTags(registered_tape_indices_s); + + std::vector registered_tape_indices( + registered_tape_indices_s.size()); + std::copy(registered_tape_indices_s.begin(), + registered_tape_indices_s.end(), + registered_tape_indices.begin()); + return registered_tape_indices; + } }; # else // DEAL_II_WITH_ADOLC @@ -838,29 +1209,73 @@ namespace Differentiation // === Taping === - static void - start_taping(const types::tape_index, const bool) + bool + is_recording() const { AssertThrow(false, ExcRequiresAdolC()); + return false; } - static void - start_taping(const types::tape_index, - const bool, - const types::tape_buffer_sizes, - const types::tape_buffer_sizes, - const types::tape_buffer_sizes, - const types::tape_buffer_sizes) + types::tape_index + active_tape_index() const { AssertThrow(false, ExcRequiresAdolC()); + return numbers::invalid_tape_index; } - static void + bool + keep_independent_values() const + { + AssertThrow(false, ExcRequiresAdolC()); + return false; + } + + bool + is_registered_tape(const types::tape_index tape_index) const + { + AssertThrow(false, ExcRequiresAdolC()); + return false; + } + + void + set_tape_buffer_sizes(const types::tape_buffer_sizes, + const types::tape_buffer_sizes, + const types::tape_buffer_sizes, + const types::tape_buffer_sizes) + { + AssertThrow(false, ExcRequiresAdolC()); + } + + void + start_taping(const types::tape_index, const bool) + { + AssertThrow(false, ExcRequiresAdolC()); + } + + void stop_taping(const types::tape_index, const bool) { AssertThrow(false, ExcRequiresAdolC()); } + void + activate_tape(const types::tape_index) + { + AssertThrow(false, ExcRequiresAdolC()); + } + + void + reset(const bool) + { + AssertThrow(false, ExcRequiresAdolC()); + } + + void + print(std::ostream &stream) const + { + AssertThrow(false, ExcRequiresAdolC()); + } + static void print_tape_stats(const types::tape_index, std::ostream &) { @@ -936,43 +1351,84 @@ namespace Differentiation // === Taping === - static void - start_taping(const types::tape_index tape_index, - const bool keep_independent_values) + bool + is_recording() const { // ADOL-C only supports 'double', not 'float', so we can forward to // the 'double' implementation of this function - TapedDrivers::start_taping( - tape_index, keep_independent_values); + return taped_driver.is_recording(); } - static void - start_taping(const types::tape_index tape_index, - const bool keep_independent_values, - const types::tape_buffer_sizes obufsize, - const types::tape_buffer_sizes lbufsize, - const types::tape_buffer_sizes vbufsize, - const types::tape_buffer_sizes tbufsize) + types::tape_index + active_tape_index() const { // ADOL-C only supports 'double', not 'float', so we can forward to // the 'double' implementation of this function - TapedDrivers::start_taping( - tape_index, - keep_independent_values, - obufsize, - lbufsize, - vbufsize, - tbufsize); + return taped_driver.active_tape_index(); } - static void + bool + keep_independent_values() const + { + return taped_driver.keep_independent_values(); + } + + bool + is_registered_tape(const types::tape_index tape_index) const + { + // ADOL-C only supports 'double', not 'float', so we can forward to + // the 'double' implementation of this function + return taped_driver.is_registered_tape(tape_index); + } + + void + set_tape_buffer_sizes(const types::tape_buffer_sizes obufsize, + const types::tape_buffer_sizes lbufsize, + const types::tape_buffer_sizes vbufsize, + const types::tape_buffer_sizes tbufsize) + { + // ADOL-C only supports 'double', not 'float', so we can forward to + // the 'double' implementation of this function + taped_driver.set_tape_buffer_sizes(obufsize, + lbufsize, + vbufsize, + tbufsize); + } + + void + start_taping(const types::tape_index tape_index, + const bool keep_independent_values) + { + // ADOL-C only supports 'double', not 'float', so we can forward to + // the 'double' implementation of this function + taped_driver.start_taping(tape_index, keep_independent_values); + } + + void stop_taping(const types::tape_index active_tape_index, const bool write_tapes_to_file) { // ADOL-C only supports 'double', not 'float', so we can forward to // the 'double' implementation of this function - TapedDrivers::stop_taping(active_tape_index, - write_tapes_to_file); + taped_driver.stop_taping(active_tape_index, write_tapes_to_file); + } + + void + activate_tape(const types::tape_index tape_index) + { + taped_driver.activate_tape(tape_index); + } + + void + reset(const bool clear_registered_tapes) + { + taped_driver.reset(clear_registered_tapes); + } + + void + print(std::ostream &stream) const + { + taped_driver.print(stream); } static void @@ -1073,6 +1529,11 @@ namespace Differentiation std::copy(in.begin(), in.end(), out.begin()); return out; } + + /** + * The object that actually takes care of the taping + */ + TapedDrivers taped_driver; }; @@ -1087,6 +1548,31 @@ namespace Differentiation AssertThrow(false, ExcRequiresADNumberSpecialization()); } + template + void + TapelessDrivers:: + allow_dependent_variable_marking() + { + AssertThrow(false, ExcRequiresADNumberSpecialization()); + } + + template + void + TapelessDrivers:: + prevent_dependent_variable_marking() + { + AssertThrow(false, ExcRequiresADNumberSpecialization()); + } + + template + bool + TapelessDrivers:: + is_dependent_variable_marking_allowed() const + { + AssertThrow(false, ExcRequiresADNumberSpecialization()); + return false; + } + template ScalarType @@ -1297,6 +1783,10 @@ namespace Differentiation ADNumberTraits::type_code == NumberTypes::sacado_rad_dfad>::type> { + TapelessDrivers() + : dependent_variable_marking_safe(false) + {} + // === Configuration === static void @@ -1306,6 +1796,25 @@ namespace Differentiation n_independent_variables); } + // === Operation status === + void + allow_dependent_variable_marking() + { + dependent_variable_marking_safe = true; + } + + void + prevent_dependent_variable_marking() + { + dependent_variable_marking_safe = false; + } + + bool + is_dependent_variable_marking_allowed() const + { + return dependent_variable_marking_safe; + } + // === Scalar drivers === static ScalarType @@ -1460,6 +1969,13 @@ namespace Differentiation } } } + + private: + /** + * A flag that states whether or not dependent variables can be marked + * within the current phase of operations. + */ + bool dependent_variable_marking_safe; }; @@ -1479,6 +1995,10 @@ namespace Differentiation ADNumberTraits::type_code == NumberTypes::sacado_dfad_dfad>::type> { + TapelessDrivers() + : dependent_variable_marking_safe(false) + {} + // === Configuration === static void @@ -1488,6 +2008,25 @@ namespace Differentiation n_independent_variables); } + // === Operation status === + void + allow_dependent_variable_marking() + { + dependent_variable_marking_safe = true; + } + + void + prevent_dependent_variable_marking() + { + dependent_variable_marking_safe = false; + } + + bool + is_dependent_variable_marking_allowed() const + { + return dependent_variable_marking_safe; + } + // === Scalar drivers === static ScalarType @@ -1613,6 +2152,13 @@ namespace Differentiation ADNumberTraits::get_directional_derivative( dependent_variables[i], j)); } + + private: + /** + * A flag that states whether or not dependent variables can be marked + * within the current phase of operations. + */ + bool dependent_variable_marking_safe; }; } // namespace AD diff --git a/include/deal.II/differentiation/ad/ad_helpers.h b/include/deal.II/differentiation/ad/ad_helpers.h index 4596bd386f..1ec78b6b7e 100644 --- a/include/deal.II/differentiation/ad/ad_helpers.h +++ b/include/deal.II/differentiation/ad/ad_helpers.h @@ -207,7 +207,7 @@ namespace Differentiation /** * Destructor */ - virtual ~ADHelperBase(); + virtual ~ADHelperBase() = default; //@} @@ -347,11 +347,11 @@ namespace Differentiation is_recording() const; /** - * Returns the tape number which is currently activated for recording or + * Returns the tape index which is currently activated for recording or * reading. */ types::tape_index - active_tape() const; + active_tape_index() const; /** * Returns whether or not a tape number has already been used @@ -424,10 +424,10 @@ namespace Differentiation * @param[in] overwrite_tape Express whether tapes are allowed to be * overwritten. If true then any existing tape with a given * @p tape_index will be destroyed and a new tape traced over it. - * @param[in] keep_values Determines whether the numerical values of all - * independent variables are recorded in the tape buffer. If true, then - * the tape can be immediately used to perform computations after - * recording is complete. + * @param[in] keep_independent_values Determines whether the numerical + * values of all independent variables are recorded in the tape buffer. + * If true, then the tape can be immediately used to perform computations + * after recording is complete. * * @note During the recording phase, no value(), gradient(), hessian(), * or jacobian() operations can be performed. @@ -438,7 +438,7 @@ namespace Differentiation bool start_recording_operations(const types::tape_index tape_index, const bool overwrite_tape = false, - const bool keep_values = true); + const bool keep_independent_values = true); /** * Disable recording mode for a given tape. The use of this function is @@ -469,66 +469,25 @@ namespace Differentiation protected: /** - * @name Taping + * @name Drivers and taping */ //@{ /** - * Index of the tape that is currently in use. It is this tape that will - * be recorded to or read from when performing computations using "taped" - * auto-differentiable numbers. - */ - types::tape_index active_tape_index; - - /** - * A collection of tapes that have been recorded to on this process. + * An object used to help manage stored tapes. * - * It is important to keep track of this so that one doesn't accidentally - * record over a tape (unless specifically instructed to) and that one - * doesn't try to use a tape that doesn't exist. - */ - static std::set registered_tapes; - - /** - * Mark whether we're going to inform taped data structures to retain - * the coefficients ("Taylors" in ADOL-C nomenclature) stored on the - * tape so that they can be evaluated again at a later stage. - */ - bool keep_values; - - /** - * Mark whether we're currently recording a tape. Dependent on the state - * of this flag, only a restricted set of operations are allowable. - */ - bool is_recording_flag; - - /** - * A flag indicating that we should preferentially use the user-defined - * taped buffer sizes as opposed to either the default values selected - * by the AD library (or, in the case of ADOL-C, defined in an - * ".adolcrc" file). - */ - bool use_stored_taped_buffer_sizes; - - /** - * ADOL-C operations buffer size. - */ - types::tape_buffer_sizes obufsize; - - /** - * ADOL-C locations buffer size. - */ - types::tape_buffer_sizes lbufsize; - - /** - * ADOL-C value buffer size. + * In the event that the @p ad_type is a tapeless AD type, then the object + * constructed here is, effectively, a dummy one. */ - types::tape_buffer_sizes vbufsize; + TapedDrivers taped_driver; /** - * ADOL-C Taylor buffer size. + * An object used to help manage tapeless data structures. + * + * In the event that the @p ad_type is a taped AD type, then the object + * constructed here is, effectively, a dummy one. */ - types::tape_buffer_sizes tbufsize; + TapelessDrivers tapeless_driver; /** * Select a tape to record to or read from. @@ -715,24 +674,6 @@ namespace Differentiation //@} - private: - /** - * @name Miscellaneous - */ - //@{ - - /** - * A counter keeping track of the number of helpers in existence. - * - * This is only important information for when we use taped number types. - * As the tapes are stored in some global register, they exist independent - * of these helpers. However, it is assumed that when all helpers go - * out of scope then the tapes can be written over. - */ - static unsigned int n_helpers; - - //@} - }; // class ADHelperBase diff --git a/source/differentiation/ad/ad_helpers.cc b/source/differentiation/ad/ad_helpers.cc index ab22f87b82..7437d82fdb 100644 --- a/source/differentiation/ad/ad_helpers.cc +++ b/source/differentiation/ad/ad_helpers.cc @@ -34,38 +34,17 @@ namespace Differentiation - template - std::set - ADHelperBase::registered_tapes; - - - - template - unsigned int ADHelperBase::n_helpers = 0; - - - template ADHelperBase::ADHelperBase( const unsigned int n_independent_variables, const unsigned int n_dependent_variables) - : active_tape_index(numbers::invalid_tape_index) - , keep_values(true) - , is_recording_flag(false) - , use_stored_taped_buffer_sizes(false) - , obufsize(0u) - , lbufsize(0u) - , vbufsize(0u) - , tbufsize(0u) - , independent_variable_values( + : independent_variable_values( n_independent_variables, dealii::internal::NumberType::value(0.0)) , registered_independent_variable_values(n_independent_variables, false) , registered_marked_independent_variables(n_independent_variables, false) , registered_marked_dependent_variables(n_dependent_variables, false) { - ++n_helpers; - // Tapeless mode must be configured before any active live // variables are created. if (AD::is_tapeless_ad_number::value) @@ -84,20 +63,6 @@ namespace Differentiation - template - ADHelperBase::~ADHelperBase() - { - --n_helpers; - - // Clear any static data when there are no more helpers in scope. - // This means that we've effectively marked any of these tapes to - // be safely overwritten. - if (n_helpers == 0) - registered_tapes.clear(); - } - - - template void ADHelperBase::is_taped == true) { - Assert(this->active_tape() != numbers::invalid_tape_index, + Assert(this->active_tape_index() != numbers::invalid_tape_index, ExcMessage("Invalid tape index")); } Assert( @@ -182,7 +147,7 @@ namespace Differentiation if (ADNumberTraits::is_taped == true) { - Assert(active_tape() != numbers::invalid_tape_index, + Assert(active_tape_index() != numbers::invalid_tape_index, ExcMessage("Invalid tape index")); Assert(is_recording() == true, ExcMessage( @@ -232,7 +197,7 @@ namespace Differentiation { if (ADNumberTraits::is_taped == true) { - Assert(active_tape() != numbers::invalid_tape_index, + Assert(active_tape_index() != numbers::invalid_tape_index, ExcMessage("Invalid tape index")); } Assert(is_recording() == false, @@ -295,16 +260,22 @@ namespace Differentiation bool ADHelperBase::is_recording() const { - return is_recording_flag; + if (AD::is_taped_ad_number::value) + return taped_driver.is_recording(); + else + return tapeless_driver.is_dependent_variable_marking_allowed(); } template types::tape_index - ADHelperBase::active_tape() const + ADHelperBase::active_tape_index() const { - return active_tape_index; + if (AD::is_taped_ad_number::value) + return taped_driver.active_tape_index(); + else + return 1; } @@ -314,7 +285,10 @@ namespace Differentiation ADHelperBase::is_registered_tape( const types::tape_index tape_index) const { - return registered_tapes.find(tape_index) != registered_tapes.end(); + if (AD::is_taped_ad_number::value) + return taped_driver.is_registered_tape(tape_index); + else + return true; } @@ -329,20 +303,13 @@ namespace Differentiation // Set stream to print booleans as "true"/"false" stream.setf(std::ios_base::boolalpha); - stream << "Active tape index: " << active_tape_index << "\n"; - stream << "Registered tapes: "; - auto it_registered_tape = registered_tapes.begin(); - for (unsigned int i = 0; i < registered_tapes.size(); - ++i, ++it_registered_tape) - stream << *it_registered_tape - << (i < (registered_tapes.size() - 1) ? "," : ""); - stream << "\n"; + stream << "Active tape index: " << active_tape_index() << "\n"; stream << "Recording? " << is_recording() << "\n"; - stream << "Keep values? " << keep_values << "\n"; - stream << "Use stored tape buffer sizes? " - << use_stored_taped_buffer_sizes << "\n"; stream << std::flush; + if (AD::is_taped_ad_number::value) + taped_driver.print(stream); + stream << "Registered independent variables: " << "\n"; for (unsigned int i = 0; i < n_independent_variables(); i++) @@ -452,11 +419,8 @@ namespace Differentiation configure_tapeless_mode(new_n_independent_variables, false /*ensure_persistent_setting*/); } - - active_tape_index = numbers::invalid_tape_index; - is_recording_flag = false; - if (clear_registered_tapes) - registered_tapes.clear(); + if (AD::is_taped_ad_number::value) + taped_driver.reset(clear_registered_tapes); independent_variable_values = std::vector( new_n_independent_variables, @@ -521,7 +485,7 @@ namespace Differentiation ExcMessage("Invalid tape index")); Assert(tape_index < numbers::max_tape_index, ExcMessage("Tape index exceeds maximum allowable value")); - active_tape_index = tape_index; + taped_driver.activate_tape(tape_index); reset_registered_independent_variables(); // A tape may have been defined by a different ADHelper, so in this @@ -541,7 +505,6 @@ namespace Differentiation { Assert(ADNumberTraits::is_tapeless == true, ExcInternalError()); - active_tape_index = tape_index; // This is, in effect, a dummy value } } @@ -550,18 +513,18 @@ namespace Differentiation template void ADHelperBase::set_tape_buffer_sizes( - const types::tape_buffer_sizes in_obufsize, - const types::tape_buffer_sizes in_lbufsize, - const types::tape_buffer_sizes in_vbufsize, - const types::tape_buffer_sizes in_tbufsize) + const types::tape_buffer_sizes obufsize, + const types::tape_buffer_sizes lbufsize, + const types::tape_buffer_sizes vbufsize, + const types::tape_buffer_sizes tbufsize) { // When valid for the chosen AD number type, these values will be used the // next time start_recording_operations() is called. - obufsize = in_obufsize; - lbufsize = in_lbufsize; - vbufsize = in_vbufsize; - tbufsize = in_tbufsize; - use_stored_taped_buffer_sizes = true; + if (ADNumberTraits::is_taped == true) + taped_driver.set_tape_buffer_sizes(obufsize, + lbufsize, + vbufsize, + tbufsize); } @@ -571,7 +534,7 @@ namespace Differentiation ADHelperBase::start_recording_operations( const types::tape_index tape_index, const bool overwrite_tape, - const bool keep) + const bool keep_independent_values) { // Define this here for clarity when this flag is used later. const bool read_mode = false; @@ -583,62 +546,41 @@ namespace Differentiation Assert(is_recording() == false, ExcMessage("Already recording...")); } + + // Check conditions to enable tracing if (is_registered_tape(tape_index) == false || overwrite_tape == true) { - // Register the tape... - registered_tapes.insert(tape_index); - // ... and setup the data structures for this class in the + // Setup the data structures for this class in the // appropriate manner activate_tape(tape_index, read_mode); // Start taping - if (use_stored_taped_buffer_sizes) - { - Assert(obufsize > 0u, - ExcMessage("Buffer value not initialized.")); - Assert(lbufsize > 0u, - ExcMessage("Buffer value not initialized.")); - Assert(vbufsize > 0u, - ExcMessage("Buffer value not initialized.")); - Assert(tbufsize > 0u, - ExcMessage("Buffer value not initialized.")); - TapedDrivers::start_taping( - active_tape(), - keep, - obufsize, - lbufsize, - vbufsize, - tbufsize); - - // Reset this, as we don't assume that the same tape buffer - // sizes are used for every tape - use_stored_taped_buffer_sizes = false; - obufsize = 0u; - lbufsize = 0u; - vbufsize = 0u; - tbufsize = 0u; - } - else - { - TapedDrivers::start_taping( - active_tape(), keep); - } + taped_driver.start_taping(active_tape_index(), + keep_independent_values); // Clear the flags that state which independent and // dependent variables have been registered reset_registered_independent_variables(); reset_registered_dependent_variables(); + } + else + { + Assert(is_recording() == false, + ExcMessage("Tape recording is unexpectedly still enabled.")); - // Set some other flags to their indicated / required values - keep_values = keep; - is_recording_flag = true; + // Now we activate the pre-recorded tape so that its immediately + // available for use + activate_recorded_tape(tape_index); } } else { Assert(ADNumberTraits::is_tapeless == true, ExcInternalError()); - is_recording_flag = true; + + // Set the flag that states that we can safely mark dependent + // variables within this current phase of operations + tapeless_driver.allow_dependent_variable_marking(); // Dummy call to ensure that the intuitively correct // value for the active tape (whether "valid" or not) @@ -646,7 +588,7 @@ namespace Differentiation activate_tape(tape_index, read_mode); } - return is_recording_flag; + return is_recording(); } @@ -664,20 +606,8 @@ namespace Differentiation if (ADNumberTraits::is_taped == true) { - // Stop taping - TapedDrivers::stop_taping(active_tape(), - write_tapes_to_file); - - // Now that we've turned tracing off, we've definitely - // stopped all tape recording. - is_recording_flag = false; - - // If the keep_values flag is set, then we expect the user to use this - // tape immediately after recording it. There is therefore no need to - // invalidate it. However, there is now also no way to double-check - // that the newly recorded tape is indeed the active tape. - if (keep_values == false) - active_tape_index = numbers::invalid_tape_index; + // Stop tracing + taped_driver.stop_taping(active_tape_index(), write_tapes_to_file); } else { @@ -689,12 +619,11 @@ namespace Differentiation // By changing this flag, we ensure that the we can no longer // legally alter the values of the dependent variables using - // set_independent_variable(). This is important because the value of + // set_dependent_variable(). This is important because the value of // the tapeless independent variables are set and finalized when - // mark_independent_variable() is called. So we cannot allow this to - // be done when not in the "recording" phase - is_recording_flag = false; - active_tape_index = numbers::invalid_tape_index; + // mark_dependent_variable() is called. So we cannot allow this to + // be done when not in the "recording" phase. + tapeless_driver.prevent_dependent_variable_marking(); } } @@ -713,7 +642,7 @@ namespace Differentiation if (ADNumberTraits::is_taped == true) { - Assert(active_tape() != numbers::invalid_tape_index, + Assert(active_tape_index() != numbers::invalid_tape_index, ExcMessage("Invalid tape index")); Assert(is_recording() == true, ExcMessage( @@ -772,7 +701,7 @@ namespace Differentiation { if (ADNumberTraits::is_taped == true) { - Assert(this->active_tape() != numbers::invalid_tape_index, + Assert(this->active_tape_index() != numbers::invalid_tape_index, ExcMessage("Invalid tape index")); } @@ -796,7 +725,7 @@ namespace Differentiation { if (ADNumberTraits::is_taped == true) { - Assert(this->active_tape() != numbers::invalid_tape_index, + Assert(this->active_tape_index() != numbers::invalid_tape_index, ExcMessage("Invalid tape index")); } @@ -818,7 +747,7 @@ namespace Differentiation { if (ADNumberTraits::is_taped == true) { - Assert(this->active_tape() != numbers::invalid_tape_index, + Assert(this->active_tape_index() != numbers::invalid_tape_index, ExcMessage("Invalid tape index")); } Assert(values.size() == this->n_independent_variables(), @@ -862,7 +791,8 @@ namespace Differentiation ADHelperEnergyFunctional::compute_energy() const { - if (this->keep_values == false || + if ((ADNumberTraits::is_taped == true && + this->taped_driver.keep_independent_values() == false) || ADNumberTraits::is_tapeless == true) { Assert( @@ -882,7 +812,7 @@ namespace Differentiation if (ADNumberTraits::is_taped == true) { - Assert(this->active_tape() != numbers::invalid_tape_index, + Assert(this->active_tape_index() != numbers::invalid_tape_index, ExcMessage("Invalid tape index")); Assert(this->is_recording() == false, ExcMessage( @@ -893,7 +823,7 @@ namespace Differentiation this->n_independent_variables())); return TapedDrivers::value( - this->active_tape(), this->independent_variable_values); + this->active_tape_index(), this->independent_variable_values); } else { @@ -915,7 +845,8 @@ namespace Differentiation ADHelperEnergyFunctional::compute_residual( Vector &gradient) const { - if (this->keep_values == false || + if ((ADNumberTraits::is_taped == true && + this->taped_driver.keep_independent_values() == false) || ADNumberTraits::is_tapeless == true) { Assert( @@ -941,7 +872,7 @@ namespace Differentiation if (ADNumberTraits::is_taped == true) { - Assert(this->active_tape() != numbers::invalid_tape_index, + Assert(this->active_tape_index() != numbers::invalid_tape_index, ExcMessage("Invalid tape index")); Assert(this->is_recording() == false, ExcMessage( @@ -952,7 +883,9 @@ namespace Differentiation this->n_independent_variables())); TapedDrivers::gradient( - this->active_tape(), this->independent_variable_values, gradient); + this->active_tape_index(), + this->independent_variable_values, + gradient); } else { @@ -979,7 +912,8 @@ namespace Differentiation "Cannot computed function Hessian: AD number type does" "not support the calculation of second order derivatives.")); - if (this->keep_values == false) + if ((ADNumberTraits::is_taped == true && + this->taped_driver.keep_independent_values() == false)) { Assert( this->n_registered_independent_variables() == @@ -1006,7 +940,7 @@ namespace Differentiation if (ADNumberTraits::is_taped == true) { - Assert(this->active_tape() != numbers::invalid_tape_index, + Assert(this->active_tape_index() != numbers::invalid_tape_index, ExcMessage("Invalid tape index")); Assert(this->is_recording() == false, ExcMessage( @@ -1017,7 +951,9 @@ namespace Differentiation this->n_independent_variables())); TapedDrivers::hessian( - this->active_tape(), this->independent_variable_values, hessian); + this->active_tape_index(), + this->independent_variable_values, + hessian); } else { -- 2.39.5