H Doug and Wolfgang,

At the time when we introduced the AD framework I took a stab at quickly adding 
support for to our Vector class. To summarise, I hit the same issue and 
couldn’t find an elegant solution past to get past it (so hopefully this is the 
only outstanding problem to add support to Vector & friends). We have all of 
the tooling in place to convert AD numbers to something printable, but the 
issue is that we create circular inclusions in our headers when trying to use 
them, since the exceptions.h is a very fundamental header and number.h and 
almost all other headers require it.

My suggestion would be to not attack the problem in exceptions.h, but rather in 
vector.templates.h and la_parallel_vector.templates.h itself. So each call to 
AssertIsFinite(), for example over here 
<https://github.com/dealii/dealii/blob/master/include/deal.II/lac/vector.templates.h#L288>,
 should have sent into it a number that’s pre-converted to a 
std::complex<double>. You should be able to do this by invoking the conversion 
tool that will extract the value data for all AD numbers (so for free you’d 
include support for all Sacado and ADOL-C types):
AssertIsFinite(internal::NumberType<std::complex<double>>::value(s));

I’m pretty sure that will work. It’s a bit annoying, but I really think that 
might be the best way forward. If we go with this approach then I would suggest 
that we also add a note to the documentation of AssertIsFinite() to highlight 
this solution.

I hope this helps.

Best,
Jean-Paul

> On 18 Sep 2019, at 03:40, Doug <doug.shid...@gmail.com> wrote:
> 
> so 'a' is of type Sacado::Fad::DFad<double>. It then needs to call 
> 
>    numbers::is_finite (Sacado::Fad::DFad<double>), which doesn't exist. 
> It probably just tries to go through all of the overloads of 
> numbers::is_finite() and wants to see whether it can convert the 
> argument Sacado::Fad::DFad<double> to the argument type of these 
> overloads. The error message you show then would just explain why this 
> one possibility (namely, numbers::is_finite(std::complex<long double>&)) 
> does not work. But that doesn't mean that this is the right overload 
> anyway -- I suspect that your compiler produces similar error messages 
> above or below the one you show for all of the other overloads, right? 
> 
> 
> True, the error message gets long pretty quickly, but the undefined is_finite 
> was another issue. Even if is_finite exists, the complex constructor is still 
> an issue.
>  
> I *think* that the solution is to simply provide an overload for 
>    numbers::is_finite (const Sacado::Fad::DFad<double> &x) 
> Can you try this? You could declare it in your own .cc file before you 
> #include <deal.II/lac/la_parallel_vector.templates.h> 
> 
>  Attached is a tiny .cc file that simply instantiates the vector. Copy pasted 
> here here convenience.
> 
> #include <Sacado.hpp>
> namespace dealii{
> namespace numbers{
>     bool is_finite(const Sacado::Fad::DFad<double> &x) {
>     ¦   (void) x;
>     ¦   return true;
>     }   
> }}
> #include <deal.II/lac/vector.h>
> #include <deal.II/lac/la_parallel_vector.templates.h>
> int main (int /*argc*/, char * /*argv*/[]){
>     using namespace dealii;
>     dealii::LinearAlgebra::distributed::Vector<double> vector_double;
>     using ADtype = Sacado::Fad::DFad<double>;
>     dealii::LinearAlgebra::distributed::Vector<ADtype> vector_ad;
>     return 0;
> }
> 
> I also provide an function for numbers::is_finite (const 
> Sacado::Fad::DFad<double> &x) to avoid the first set of errors. However, I 
> still get the error below
> 
> In file included from 
> /home/ddong/Libraries/dealii/install/include/deal.II/base/aligned_vector.h:22:0,
>                  from 
> /home/ddong/Libraries/dealii/install/include/deal.II/lac/vector.h:22,
>                  from 
> /home/ddong/Codes/PHiLiP/src/instantiate_vector_ad.cpp:13:
> /home/ddong/Libraries/dealii/install/include/deal.II/lac/la_parallel_vector.templates.h:
>  In instantiation of ‘dealii::LinearAlgebra::distributed::Vector<Number, 
> MemorySpace>& dealii::LinearAlgebra::distributed::Vector<Number, 
> MemorySpace>::operator*=(Number) [with Number = Sacado::Fad::DFad<double>; 
> MemorySpace = dealii::MemorySpace::Host]’:
> /home/ddong/Codes/PHiLiP/src/instantiate_vector_ad.cpp:26:1:   required from 
> here
> /home/ddong/Libraries/dealii/install/include/deal.II/lac/la_parallel_vector.templates.h:1651:7:
>  error: no matching function for call to ‘std::complex<double>::complex(const 
> Sacado::Fad::DFad<double>&)’
>        AssertIsFinite(factor);
> In file included from /usr/include/trilinos/Teuchos_ConfigDefs.hpp:94:0,
>                  from /usr/include/trilinos/Teuchos_PromotionTraits.hpp:45,
>                  from 
> /usr/include/trilinos/Sacado_Fad_Exp_GeneralFadTraits.hpp:139,
>                  from /usr/include/trilinos/Sacado.hpp:52,
>                  from 
> /home/ddong/Codes/PHiLiP/src/instantiate_vector_ad.cpp:1:
> /usr/include/c++/7/complex:1512:3: note: candidate: constexpr 
> std::complex<double>::complex(const std::complex<long double>&)
>    complex<double>::complex(const complex<long double>& __z)
>    ^~~~~~~~~~~~~~~
> /usr/include/c++/7/complex:1512:3: note:   no known conversion for argument 1 
> from ‘const Sacado::Fad::DFad<double>’ to ‘const std::complex<long double>&’
> /usr/include/c++/7/complex:1219:26: note: candidate: constexpr 
> std::complex<double>::complex(const std::complex<float>&)
>        _GLIBCXX_CONSTEXPR complex(const complex<float>& __z)
>                           ^~~~~~~
> /usr/include/c++/7/complex:1219:26: note:   no known conversion for argument 
> 1 from ‘const Sacado::Fad::DFad<double>’ to ‘const std::complex<float>&’
> /usr/include/c++/7/complex:1209:26: note: candidate: constexpr 
> std::complex<double>::complex(double, double)
>        _GLIBCXX_CONSTEXPR complex(double __r = 0.0, double __i = 0.0)
>                           ^~~~~~~
> /usr/include/c++/7/complex:1209:26: note:   no known conversion for argument 
> 1 from ‘const Sacado::Fad::DFad<double>’ to ‘double’
> /usr/include/c++/7/complex:1207:26: note: candidate: constexpr 
> std::complex<double>::complex(std::complex<double>::_ComplexT)
>        _GLIBCXX_CONSTEXPR complex(_ComplexT __z) : _M_value(__z) { }
>                           ^~~~~~~
> /usr/include/c++/7/complex:1207:26: note:   no known conversion for argument 
> 1 from ‘const Sacado::Fad::DFad<double>’ to ‘std::complex<double>::_ComplexT 
> {aka __complex__ double}’
> /usr/include/c++/7/complex:1202:12: note: candidate: constexpr 
> std::complex<double>::complex(const std::complex<double>&)
>      struct complex<double>
>             ^~~~~~~~~~~~~~~
> /usr/include/c++/7/complex:1202:12: note:   no known conversion for argument 
> 1 from ‘const Sacado::Fad::DFad<double>’ to ‘const std::complex<double>&’
> /usr/include/c++/7/complex:1202:12: note: candidate: constexpr 
> std::complex<double>::complex(std::complex<double>&&)
> /usr/include/c++/7/complex:1202:12: note:   no known conversion for argument 
> 1 from ‘const Sacado::Fad::DFad<double>’ to ‘std::complex<double>&&’
> In file included from 
> /home/ddong/Libraries/dealii/install/include/deal.II/base/aligned_vector.h:22:0,
>                  from 
> /home/ddong/Libraries/dealii/install/include/deal.II/lac/vector.h:22,
>                  from 
> /home/ddong/Codes/PHiLiP/src/instantiate_vector_ad.cpp:13:
> 
> 
> So, even if the is_finite function exists, the exception to be thrown still 
> requires the complex constructor std::complex<double>::complex(const 
> Sacado::Fad::DFad<double>&) to exist.
> 
> -- 
> The deal.II project is located at http://www.dealii.org/ 
> <http://www.dealii.org/>
> For mailing list/forum options, see 
> https://groups.google.com/d/forum/dealii?hl=en 
> <https://groups.google.com/d/forum/dealii?hl=en>
> --- 
> You received this message because you are subscribed to the Google Groups 
> "deal.II User Group" group.
> To unsubscribe from this group and stop receiving emails from it, send an 
> email to dealii+unsubscr...@googlegroups.com 
> <mailto:dealii+unsubscr...@googlegroups.com>.
> To view this discussion on the web visit 
> https://groups.google.com/d/msgid/dealii/7c3d83a7-154c-4bd2-a736-6c5c83f6f829%40googlegroups.com
>  
> <https://groups.google.com/d/msgid/dealii/7c3d83a7-154c-4bd2-a736-6c5c83f6f829%40googlegroups.com?utm_medium=email&utm_source=footer>.
> <instantiate_vector_ad.cpp>

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/F1F1A9D9-48F0-4B7D-9EDF-A81786D58E59%40gmail.com.

Reply via email to