Reference documentation for deal.II version 9.2.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
symengine_number_visitor_internal.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2020 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE at
12 // the top level of the deal.II distribution.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_differentiation_sd_symengine_number_visitor_internal_h
17 #define dealii_differentiation_sd_symengine_number_visitor_internal_h
18 
19 #include <deal.II/base/config.h>
20 
21 #ifdef DEAL_II_WITH_SYMENGINE
22 
24 // Low level
25 # include <symengine/basic.h>
26 # include <symengine/dict.h>
27 # include <symengine/symengine_exception.h>
28 # include <symengine/symengine_rcp.h>
29 
30 // Visitor
31 # include <symengine/visitor.h>
32 
34 
35 # include <deal.II/base/exceptions.h>
36 # include <deal.II/base/numbers.h>
37 
40 
41 # include <boost/serialization/split_member.hpp>
42 
43 
45 
46 
47 namespace Differentiation
48 {
49  namespace SD
50  {
51  namespace internal
52  {
62  template <typename ReturnType, typename ExpressionType>
64  {
65  using symbol_vector_pair =
66  std::vector<std::pair<SD::Expression, SD::Expression>>;
67 
68  public:
69  /*
70  * Constructor.
71  */
72  CSEDictionaryVisitor() = default;
73 
74  /*
75  * Destructor.
76  */
77  virtual ~CSEDictionaryVisitor() = default;
78 
91  void
92  init(const types::symbol_vector &dependent_functions);
93 
131  void
132  call(ReturnType * output_values,
133  const types::symbol_vector &independent_symbols,
134  const ReturnType * substitution_values);
135 
140  template <class Archive>
141  void
142  save(Archive &archive, const unsigned int version) const;
143 
148  template <class Archive>
149  void
150  load(Archive &archive, const unsigned int version);
151 
152 # ifdef DOXYGEN
153 
157  template <class Archive>
158  void
159  serialize(Archive &archive, const unsigned int version);
160 # else
161  // This macro defines the serialize() method that is compatible with
162  // the templated save() and load() method that have been implemented.
163  BOOST_SERIALIZATION_SPLIT_MEMBER()
164 # endif
165 
171  template <typename StreamType>
172  void
173  print(StreamType &stream) const;
174 
178  bool
179  executed() const;
180 
185  unsigned int
187 
191  unsigned int
192  n_reduced_expressions() const;
193 
194  protected:
206  void
207  init(const SymEngine::vec_basic &dependent_functions);
208 
233  void
234  call(ReturnType * output_values,
235  const SymEngine::vec_basic &independent_symbols,
236  const ReturnType * substitution_values);
237 
238  private:
239  // Note: It would be more efficient to store this data in native
240  // SymEngine types, as it would prevent some copying of the data
241  // structures. However, this makes serialization more difficult,
242  // so we use our own serializable types instead, and lose a bit
243  // of efficiency.
244 
249 
254  };
255 
256 
257 
269  template <typename ReturnType, typename ExpressionType>
271  : public SymEngine::BaseVisitor<
272  DictionarySubstitutionVisitor<ReturnType, ExpressionType>>
273  {
274  public:
275  /*
276  * Constructor.
277  */
278  DictionarySubstitutionVisitor() = default;
279 
280  /*
281  * Destructor.
282  */
283  virtual ~DictionarySubstitutionVisitor() override = default;
284 
308  void
310  const Expression & dependent_function,
311  const bool use_cse = false);
312 
335  // The following definition is required due to base class CRTP.
336  void
337  init(const SymEngine::vec_basic &independent_symbols,
338  const SymEngine::Basic & dependent_function,
339  const bool use_cse = false);
340 
361  void
364  const bool use_cse = false);
365 
366 
388  // The following definition is required due to base class CRTP.
389  void
390  init(const SymEngine::vec_basic &independent_symbols,
391  const SymEngine::vec_basic &dependent_functions,
392  const bool use_cse = false);
393 
426  void
427  call(ReturnType *output_values, const ReturnType *substitution_values);
428 
445  // The following definition is required due to base class CRTP.
446  ReturnType
447  call(const std::vector<ReturnType> &substitution_values);
448 
453  template <class Archive>
454  void
455  save(Archive &archive, const unsigned int version) const;
456 
461  template <class Archive>
462  void
463  load(Archive &archive, const unsigned int version);
464 
465 # ifdef DOXYGEN
466 
470  template <class Archive>
471  void
472  serialize(Archive &archive, const unsigned int version);
473 # else
474  // This macro defines the serialize() method that is compatible with
475  // the templated save() and load() method that have been implemented.
476  BOOST_SERIALIZATION_SPLIT_MEMBER()
477 # endif
478 
492  template <typename StreamType>
493  void
494  print(StreamType &stream,
495  const bool print_independent_symbols = false,
496  const bool print_dependent_functions = false,
497  const bool print_cse_reductions = false) const;
498 
499 # ifndef DOXYGEN
500  // The following definitions are required due to base class CRTP.
501  // Since these are not used, and therefore not important to
502  // understand, we'll define them in the most concise manner possible.
503  // We also won't bother to document their existence, since they cannot
504  // be used.
505 # define IMPLEMENT_DSV_BVISIT(Argument) \
506  void bvisit(const Argument &) \
507  { \
508  AssertThrow(false, ExcNotImplemented()); \
509  }
510 
511  IMPLEMENT_DSV_BVISIT(SymEngine::Basic)
512  IMPLEMENT_DSV_BVISIT(SymEngine::Symbol)
513  IMPLEMENT_DSV_BVISIT(SymEngine::Constant)
514  IMPLEMENT_DSV_BVISIT(SymEngine::Integer)
515  IMPLEMENT_DSV_BVISIT(SymEngine::Rational)
516  IMPLEMENT_DSV_BVISIT(SymEngine::RealDouble)
517  IMPLEMENT_DSV_BVISIT(SymEngine::ComplexDouble)
518  IMPLEMENT_DSV_BVISIT(SymEngine::Add)
519  IMPLEMENT_DSV_BVISIT(SymEngine::Mul)
520  IMPLEMENT_DSV_BVISIT(SymEngine::Pow)
521  IMPLEMENT_DSV_BVISIT(SymEngine::Log)
522  IMPLEMENT_DSV_BVISIT(SymEngine::Sin)
523  IMPLEMENT_DSV_BVISIT(SymEngine::Cos)
524  IMPLEMENT_DSV_BVISIT(SymEngine::Tan)
525  IMPLEMENT_DSV_BVISIT(SymEngine::Csc)
526  IMPLEMENT_DSV_BVISIT(SymEngine::Sec)
527  IMPLEMENT_DSV_BVISIT(SymEngine::Cot)
528  IMPLEMENT_DSV_BVISIT(SymEngine::ASin)
529  IMPLEMENT_DSV_BVISIT(SymEngine::ACos)
530  IMPLEMENT_DSV_BVISIT(SymEngine::ATan)
531  IMPLEMENT_DSV_BVISIT(SymEngine::ATan2)
532  IMPLEMENT_DSV_BVISIT(SymEngine::ACsc)
533  IMPLEMENT_DSV_BVISIT(SymEngine::ASec)
534  IMPLEMENT_DSV_BVISIT(SymEngine::ACot)
535  IMPLEMENT_DSV_BVISIT(SymEngine::Sinh)
536  IMPLEMENT_DSV_BVISIT(SymEngine::Cosh)
537  IMPLEMENT_DSV_BVISIT(SymEngine::Tanh)
538  IMPLEMENT_DSV_BVISIT(SymEngine::Csch)
539  IMPLEMENT_DSV_BVISIT(SymEngine::Sech)
540  IMPLEMENT_DSV_BVISIT(SymEngine::Coth)
541  IMPLEMENT_DSV_BVISIT(SymEngine::ASinh)
542  IMPLEMENT_DSV_BVISIT(SymEngine::ACosh)
543  IMPLEMENT_DSV_BVISIT(SymEngine::ATanh)
544  IMPLEMENT_DSV_BVISIT(SymEngine::ACsch)
545  IMPLEMENT_DSV_BVISIT(SymEngine::ACoth)
546  IMPLEMENT_DSV_BVISIT(SymEngine::ASech)
547  IMPLEMENT_DSV_BVISIT(SymEngine::Abs)
548  IMPLEMENT_DSV_BVISIT(SymEngine::Gamma)
549  IMPLEMENT_DSV_BVISIT(SymEngine::LogGamma)
550  IMPLEMENT_DSV_BVISIT(SymEngine::Erf)
551  IMPLEMENT_DSV_BVISIT(SymEngine::Erfc)
552  IMPLEMENT_DSV_BVISIT(SymEngine::Max)
553  IMPLEMENT_DSV_BVISIT(SymEngine::Min)
554 
555 # undef IMPLEMENT_DSV_BVISIT
556 # endif // DOXYGEN
557 
558  private:
559  // Note: It would be more efficient to store this data in native
560  // SymEngine types, as it would prevent some copying of the data
561  // structures. However, this makes serialization more difficult,
562  // so we use our own serializable types instead, and lose a bit
563  // of efficiency.
564 
569 
574 
581  };
582 
583 
584 
585  /* ------------------ inline and template functions ------------------ */
586 
587 
588 # ifndef DOXYGEN
589 
590  /* -------------- CommonSubexpressionEliminationVisitor -------------- */
591 
592 
593  template <typename ReturnType, typename ExpressionType>
594  void
596  const SD::types::symbol_vector &dependent_functions)
597  {
599  dependent_functions));
600  }
601 
602 
603 
604  template <typename ReturnType, typename ExpressionType>
605  void
607  const SymEngine::vec_basic &dependent_functions)
608  {
609  // After the next call, the data stored in replacements is structured
610  // as follows:
611  //
612  // replacements[i] := [f, f(x)]
613  // replacements[i].first = intermediate function label "f"
614  // replacements[i].second = intermediate function definition "f(x)"
615  //
616  // It is to be evaluated top down (i.e. index 0 to
617  // replacements.size()), with the results going back into the
618  // substitution map for the next levels. So for each "i", "x" are the
619  // superset of the input values and the previously evaluated [f_0(x),
620  // f_1(x), ..., f_{i-1}(x)].
621  //
622  // The final result is a set of reduced expressions
623  // that must be computed after the replacement
624  // values have been computed.
625  SymEngine::vec_pair se_replacements;
626  SymEngine::vec_basic se_reduced_exprs;
627  SymEngine::cse(se_replacements, se_reduced_exprs, dependent_functions);
628 
629  intermediate_symbols_exprs =
631  se_replacements);
633  se_reduced_exprs);
634  }
635 
636 
637 
638  template <typename ReturnType, typename ExpressionType>
639  void
641  ReturnType * output_values,
642  const SD::types::symbol_vector &independent_symbols,
643  const ReturnType * substitution_values)
644  {
645  call(output_values,
647  independent_symbols),
648  substitution_values);
649  }
650 
651 
652 
653  template <typename ReturnType, typename ExpressionType>
654  void
656  ReturnType * output_values,
657  const SymEngine::vec_basic &independent_symbols,
658  const ReturnType * substitution_values)
659  {
660  Assert(n_reduced_expressions() > 0, ExcInternalError());
661 
662  // First we add the input values into the substitution map...
663  SymEngine::map_basic_basic substitution_value_map;
664  for (unsigned i = 0; i < independent_symbols.size(); ++i)
665  substitution_value_map[independent_symbols[i]] =
666  static_cast<const SymEngine::RCP<const SymEngine::Basic> &>(
667  ExpressionType(substitution_values[i]));
668 
669  // ... followed by any intermediate evaluations due to the application
670  // of CSE. These are fed directly back into the substitution map...
671  for (unsigned i = 0; i < intermediate_symbols_exprs.size(); ++i)
672  {
673  const SymEngine::RCP<const SymEngine::Basic> &cse_symbol =
674  intermediate_symbols_exprs[i].first;
675  const SymEngine::RCP<const SymEngine::Basic> &cse_expr =
676  intermediate_symbols_exprs[i].second;
677  Assert(substitution_value_map.find(cse_symbol) ==
678  substitution_value_map.end(),
679  ExcMessage(
680  "Reduced symbol already appears in substitution map. "
681  "Is there a clash between the reduced symbol name and "
682  "the symbol used for an independent variable?"));
683  substitution_value_map[cse_symbol] =
684  static_cast<const SymEngine::RCP<const SymEngine::Basic> &>(
685  ExpressionType(ExpressionType(cse_expr)
686  .template substitute_and_evaluate<ReturnType>(
687  substitution_value_map)));
688  }
689 
690  // ... followed by the final reduced expressions
691  for (unsigned i = 0; i < reduced_exprs.size(); ++i)
692  output_values[i] = ExpressionType(reduced_exprs[i])
693  .template substitute_and_evaluate<ReturnType>(
694  substitution_value_map);
695  }
696 
697 
698 
699  template <typename ReturnType, typename ExpressionType>
700  template <class Archive>
701  void
703  Archive &ar,
704  const unsigned int /*version*/) const
705  {
706  // The reduced expressions depend on the intermediate expressions,
707  // so we serialize the latter before the former.
708  ar &intermediate_symbols_exprs;
709  ar &reduced_exprs;
710  }
711 
712 
713 
714  template <typename ReturnType, typename ExpressionType>
715  template <class Archive>
716  void
718  Archive &ar,
719  const unsigned int /*version*/)
720  {
721  Assert(intermediate_symbols_exprs.empty(), ExcInternalError());
722  Assert(reduced_exprs.empty(), ExcInternalError());
723 
724  // The reduced expressions depend on the intermediate expressions,
725  // so we deserialize the latter before the former.
726  ar &intermediate_symbols_exprs;
727  ar &reduced_exprs;
728  }
729 
730 
731 
732  template <typename ReturnType, typename ExpressionType>
733  template <typename StreamType>
734  void
736  StreamType &stream) const
737  {
738  stream << "Common subexpression elimination: \n";
739  stream << " Intermediate reduced expressions: \n";
740  for (unsigned i = 0; i < intermediate_symbols_exprs.size(); ++i)
741  {
742  const SymEngine::RCP<const SymEngine::Basic> &cse_symbol =
743  intermediate_symbols_exprs[i].first;
744  const SymEngine::RCP<const SymEngine::Basic> &cse_expr =
745  intermediate_symbols_exprs[i].second;
746  stream << " " << i << ": " << cse_symbol << " = " << cse_expr
747  << "\n";
748  }
749 
750  stream << " Final reduced expressions for dependent variables: \n";
751  for (unsigned i = 0; i < reduced_exprs.size(); ++i)
752  stream << " " << i << ": " << reduced_exprs[i] << "\n";
753 
754  stream << std::flush;
755  }
756 
757 
758 
759  template <typename ReturnType, typename ExpressionType>
760  bool
762  {
763  // For dictionary substitution, the CSE algorithm moves
764  // ownership of the dependent function expression definition
765  // to the entries in reduced_exprs. So its size thus determines
766  // whether CSE has been executed or not.
767  return (n_reduced_expressions() > 0) ||
768  (n_intermediate_expressions() > 0);
769  }
770 
771 
772 
773  template <typename ReturnType, typename ExpressionType>
774  unsigned int
775  CSEDictionaryVisitor<ReturnType,
776  ExpressionType>::n_intermediate_expressions() const
777  {
778  return intermediate_symbols_exprs.size();
779  }
780 
781 
782 
783  template <typename ReturnType, typename ExpressionType>
784  unsigned int
786  const
787  {
788  return reduced_exprs.size();
789  }
790 
791 
792 
793  /* ------------------ DictionarySubstitutionVisitor ------------------ */
794 
795 
796  template <typename ReturnType, typename ExpressionType>
797  void
799  const types::symbol_vector &inputs,
800  const SD::Expression & output,
801  const bool use_cse)
802  {
803  init(inputs, types::symbol_vector{output}, use_cse);
804  }
805 
806 
807 
808  template <typename ReturnType, typename ExpressionType>
809  void
811  const SymEngine::vec_basic &inputs,
812  const SymEngine::Basic & output,
813  const bool use_cse)
814  {
816  SD::Expression(output.rcp_from_this()),
817  use_cse);
818  }
819 
820 
821 
822  template <typename ReturnType, typename ExpressionType>
823  void
825  const SymEngine::vec_basic &inputs,
826  const SymEngine::vec_basic &outputs,
827  const bool use_cse)
828  {
831  use_cse);
832  }
833 
834 
835 
836  template <typename ReturnType, typename ExpressionType>
837  void
839  const types::symbol_vector &inputs,
840  const types::symbol_vector &outputs,
841  const bool use_cse)
842  {
843  independent_symbols.clear();
844  dependent_functions.clear();
845 
846  independent_symbols = inputs;
847 
848  // Perform common subexpression elimination if requested
849  // Note: After this is done, the results produced by
850  // dependent_functions and cse.reduced_exprs should be
851  // the same. We could keep the former so that we can print
852  // out the original expressions if we wish to do so.
853  if (use_cse == false)
854  dependent_functions = outputs;
855  else
856  {
857  cse.init(outputs);
858  }
859  }
860 
861 
862 
863  template <typename ReturnType, typename ExpressionType>
864  ReturnType
866  const std::vector<ReturnType> &substitution_values)
867  {
868  Assert(
869  dependent_functions.size() == 1,
870  ExcMessage(
871  "Cannot use this call function when more than one symbolic expression is to be evaluated."));
872  Assert(
873  substitution_values.size() == independent_symbols.size(),
874  ExcMessage(
875  "Input substitution vector does not match size of symbol vector."));
876 
877  ReturnType out = ::internal::NumberType<ReturnType>::value(0.0);
878  call(&out, substitution_values.data());
879  return out;
880  }
881 
882 
883 
884  template <typename ReturnType, typename ExpressionType>
885  void
887  ReturnType * output_values,
888  const ReturnType *substitution_values)
889  {
890  // Check to see if CSE has been performed
891  if (cse.executed())
892  {
893  cse.call(output_values, independent_symbols, substitution_values);
894  }
895  else
896  {
897  // Build a substitution map.
898  SymEngine::map_basic_basic substitution_value_map;
899  for (unsigned i = 0; i < independent_symbols.size(); ++i)
900  substitution_value_map[independent_symbols[i]] =
901  static_cast<const SymEngine::RCP<const SymEngine::Basic> &>(
902  ExpressionType(substitution_values[i]));
903 
904  // Since we don't know how to definitively evaluate the
905  // input number type, we create a generic Expression
906  // with the given symbolic expression and ask it to perform
907  // substitution and evaluation for us.
908  Assert(dependent_functions.size() > 0, ExcInternalError());
909  for (unsigned i = 0; i < dependent_functions.size(); ++i)
910  output_values[i] =
911  ExpressionType(dependent_functions[i])
912  .template substitute_and_evaluate<ReturnType>(
913  substitution_value_map);
914  }
915  }
916 
917 
918 
919  template <typename ReturnType, typename ExpressionType>
920  template <class Archive>
921  void
923  Archive & ar,
924  const unsigned int version) const
925  {
926  // Add some dynamic information to determine if CSE has been used,
927  // without relying on the CSE class when deserializing.
928  // const bool used_cse = cse.executed();
929  // ar &used_cse;
930 
931  // CSE and dependent variables both require the independent
932  // symbols, so we serialize them first. The dependent variables
933  // might depend on the outcome of CSE, so we have to serialize
934  // them last.
935  ar &independent_symbols;
936  cse.save(ar, version);
937  ar &dependent_functions;
938  }
939 
940 
941 
942  template <typename ReturnType, typename ExpressionType>
943  template <class Archive>
944  void
946  Archive & ar,
947  const unsigned int version)
948  {
949  Assert(cse.executed() == false, ExcInternalError());
950  Assert(cse.n_intermediate_expressions() == 0, ExcInternalError());
951  Assert(cse.n_reduced_expressions() == 0, ExcInternalError());
952 
953  // CSE and dependent variables both require the independent
954  // symbols, so we deserialize them first. The dependent variables
955  // might depend on the outcome of CSE, so we have to deserialize
956  // them last.
957  ar &independent_symbols;
958  cse.load(ar, version);
959  ar &dependent_functions;
960  }
961 
962 
963 
964  template <typename ReturnType, typename ExpressionType>
965  template <typename StreamType>
966  void
968  StreamType &stream,
969  const bool print_independent_symbols,
970  const bool print_dependent_functions,
971  const bool print_cse_reductions) const
972  {
973  if (print_independent_symbols)
974  {
975  stream << "Independent variables: \n";
976  for (unsigned i = 0; i < independent_symbols.size(); ++i)
977  stream << " " << i << ": " << independent_symbols[i] << "\n";
978 
979  stream << std::flush;
980  }
981 
982  // Check to see if CSE has been performed
983  if (print_cse_reductions && cse.executed())
984  {
985  cse.print(stream);
986  }
987  else
988  {
989  Assert(dependent_functions.size() > 0, ExcInternalError());
990 
991  if (print_dependent_functions)
992  {
993  stream << "Dependent variables: \n";
994  for (unsigned i = 0; i < dependent_functions.size(); ++i)
995  stream << " " << i << dependent_functions[i] << "\n";
996 
997  stream << std::flush;
998  }
999  }
1000  }
1001 
1002 # endif // DOXYGEN
1003 
1004  } // namespace internal
1005  } // namespace SD
1006 } // namespace Differentiation
1007 
1008 
1010 
1011 #endif // DEAL_II_WITH_SYMENGINE
1012 
1013 #endif // dealii_differentiation_sd_symengine_number_visitor_internal_h
Differentiation::SD::internal::CSEDictionaryVisitor::init
void init(const types::symbol_vector &dependent_functions)
Differentiation::SD::internal::DictionarySubstitutionVisitor::cse
CSEDictionaryVisitor< ReturnType, ExpressionType > cse
Definition: symengine_number_visitor_internal.h:580
Differentiation::SD::internal::CSEDictionaryVisitor::print
void print(StreamType &stream) const
Differentiation::SD::internal::CSEDictionaryVisitor::intermediate_symbols_exprs
symbol_vector_pair intermediate_symbols_exprs
Definition: symengine_number_visitor_internal.h:248
Differentiation::SD::Utilities::convert_basic_pair_vector_to_expression_pair_vector
std::vector< std::pair< Expression, Expression > > convert_basic_pair_vector_to_expression_pair_vector(const SymEngine::vec_pair &symbol_value_vector)
Differentiation::SD::Utilities::convert_expression_vector_to_basic_vector
SymEngine::vec_basic convert_expression_vector_to_basic_vector(const SD::types::symbol_vector &symbol_vector)
Differentiation::SD::internal::DictionarySubstitutionVisitor::dependent_functions
SD::types::symbol_vector dependent_functions
Definition: symengine_number_visitor_internal.h:573
internal::NumberType::value
static constexpr const T & value(const T &t)
Definition: numbers.h:703
Differentiation::SD::internal::DictionarySubstitutionVisitor::print
void print(StreamType &stream, const bool print_independent_symbols=false, const bool print_dependent_functions=false, const bool print_cse_reductions=false) const
Differentiation::SD::types::symbol_vector
std::vector< SD::Expression > symbol_vector
Definition: symengine_types.h:71
Differentiation::SD::internal::CSEDictionaryVisitor::n_reduced_expressions
unsigned int n_reduced_expressions() const
Differentiation::SD::internal::CSEDictionaryVisitor::n_intermediate_expressions
unsigned int n_intermediate_expressions() const
Differentiation::SD::internal::DictionarySubstitutionVisitor
Definition: symengine_number_visitor_internal.h:270
symengine_number_types.h
Differentiation::SD::internal::DictionarySubstitutionVisitor::independent_symbols
SD::types::symbol_vector independent_symbols
Definition: symengine_number_visitor_internal.h:568
Differentiation::SD::internal::CSEDictionaryVisitor
Definition: symengine_number_visitor_internal.h:63
Differentiation::SD::internal::CSEDictionaryVisitor::CSEDictionaryVisitor
CSEDictionaryVisitor()=default
DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
#define DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
Definition: config.h:409
Differentiation::SD::internal::DictionarySubstitutionVisitor::call
void call(ReturnType *output_values, const ReturnType *substitution_values)
Differentiation::SD::internal::CSEDictionaryVisitor::serialize
void serialize(Archive &archive, const unsigned int version)
Differentiation::SD::internal::DictionarySubstitutionVisitor::save
void save(Archive &archive, const unsigned int version) const
Differentiation::SD::internal::DictionarySubstitutionVisitor::load
void load(Archive &archive, const unsigned int version)
StandardExceptions::ExcMessage
static ::ExceptionBase & ExcMessage(std::string arg1)
Differentiation::SD::internal::DictionarySubstitutionVisitor::init
void init(const types::symbol_vector &independent_symbols, const Expression &dependent_function, const bool use_cse=false)
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:358
Differentiation
Definition: numbers.h:645
Differentiation::SD::internal::CSEDictionaryVisitor::~CSEDictionaryVisitor
virtual ~CSEDictionaryVisitor()=default
Threads::internal::call
void call(const std::function< RT()> &function, internal::return_value< RT > &ret_val)
Definition: thread_management.h:607
Differentiation::SD::internal::CSEDictionaryVisitor::reduced_exprs
types::symbol_vector reduced_exprs
Definition: symengine_number_visitor_internal.h:253
Differentiation::SD::internal::CSEDictionaryVisitor::executed
bool executed() const
exceptions.h
Differentiation::SD::Expression
Definition: symengine_number_types.h:179
Differentiation::SD::internal::CSEDictionaryVisitor::call
void call(ReturnType *output_values, const types::symbol_vector &independent_symbols, const ReturnType *substitution_values)
StandardExceptions::ExcInternalError
static ::ExceptionBase & ExcInternalError()
Assert
#define Assert(cond, exc)
Definition: exceptions.h:1419
Differentiation::SD::Utilities::convert_basic_vector_to_expression_vector
SD::types::symbol_vector convert_basic_vector_to_expression_vector(const SymEngine::vec_basic &symbol_vector)
Differentiation::SD::internal::CSEDictionaryVisitor::load
void load(Archive &archive, const unsigned int version)
Differentiation::SD::internal::DictionarySubstitutionVisitor::DictionarySubstitutionVisitor
DictionarySubstitutionVisitor()=default
config.h
internal
Definition: aligned_vector.h:369
Differentiation::SD::internal::CSEDictionaryVisitor::save
void save(Archive &archive, const unsigned int version) const
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:359
Differentiation::SD::internal::DictionarySubstitutionVisitor::serialize
void serialize(Archive &archive, const unsigned int version)
Differentiation::SD::internal::DictionarySubstitutionVisitor::~DictionarySubstitutionVisitor
virtual ~DictionarySubstitutionVisitor() override=default
Differentiation::SD::internal::CSEDictionaryVisitor::symbol_vector_pair
std::vector< std::pair< SD::Expression, SD::Expression > > symbol_vector_pair
Definition: symengine_number_visitor_internal.h:66
numbers.h
symengine_utilities.h
DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
#define DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
Definition: config.h:372