dolfin team mailing list archive

["Garth N. Wells" <gnw20@xxxxxxxxx>]

On 02/06/11 17:33, Benjamin Kehlet wrote:
> On 2 June 2011 17:49, Garth N. Wells <gnw20@xxxxxxxxx> wrote:
>>
>>
>> On 02/06/11 13:41, Anders Logg wrote:
>>> Anyone using or interested in the ODE solvers should take a look
>>> below.
>>>
>>> On Thu, Jun 02, 2011 at 02:17:17PM +0200, Benjamin Kehlet wrote:
>>>> On 2 June 2011 14:02, Anders Logg <logg@xxxxxxxxx> wrote:
>>>>> On Thu, Jun 02, 2011 at 01:10:01PM +0200, Benjamin Kehlet wrote:
>>>>>> On 2 June 2011 11:51, Anders Logg <logg@xxxxxxxxx> wrote:
>>>>>>> On Thu, Jun 02, 2011 at 10:46:29AM +0100, Garth N. Wells wrote:
>>>>>>>>
>>>>>>>>
>>>>>>>> On 02/06/11 10:26, Anders Logg wrote:
>>>>>>>>> On Thu, Jun 02, 2011 at 10:07:59AM +0100, Garth N. Wells wrote:
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> On 01/06/11 23:46, Anders Logg wrote:
>>>>>>>>>>> Have you checked that there is no performance penalty?
>>>>>>>>>>
>>>>>>>>>> I just have - evaluating a Legendgre polynomial 10k times at the same
>>>>>>>>>> point is just noise with both methods (of the order 10^-5 - 10^-4 s).
>>>>>>>>>
>>>>>>>>> It may be noise for some applications, but not for others. I'm not
>>>>>>>>> sure this is a bottle-neck for the ODE code (Benjamin will know) but
>>>>>>>>> we need to evaluate Legendre polynomials of degree > 100 many times
>>>>>>>>> and then it may not be noise.
>>>>>>>>>
>>>>>>>>
>>>>>>>> For very high degree (e.g. 200) Boost is marginally faster.
>>>>>>>
>>>>>>> Sounds promising then.
>>>>>>>
>>>>>>>>>> The Boost code is slightly slower because it doesn't cache the values
>>>>>>>>>> (which is nice not to do), but may be faster if the call is inlined.
>>>>>>>>>> It's not possible to inline it at the moment because of clashes between
>>>>>>>>>> tr1:tuple and boost::tuple (Boost bug, I suspect). Old and new are the
>>>>>>>>>> same when evaluating at different points.
>>>>>>>>>
>>>>>>>>> Let's wait for Benjamin to comment.
>>>>>>>>>
>>>>>>>>
>>>>>>>> The speed is about the same (with scope to improve the speed for Boost)
>>>>>>>> for unique values. The caller should be responsible for caching, if
>>>>>>>> desired, since it can lead to memory blow out.
>>>>>>>>
>>>>>>>> Legendre does not appear in the ode code. It only appears in the
>>>>>>>> computation of quadrature schemes.
>>>>>>>
>>>>>>> True, but the quadrature schemes are used in the ode code.
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>>> Garth
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>>>> Garth
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>>> Benjamin has
>>>>>>>>>>> worked quite hard on optimizing some of the basic math routines (in
>>>>>>>>>>> some cases by many many orders of magnitude).
>>>>>>>>>>>
>>>>>>>>>>> Benjamin, can you take a look that it still works?
>>>>>>
>>>>>> Yes, the performance seems to be about the same, but I'm unable to
>>>>>> compile it with support for GMP.
>>>>>>
>>>>>> /usr/include/boost/math/special_functions/legendre.hpp:178:
>>>>>> instantiated from ‘typename boost::math::tools::promote_args<RT,
>>>>>> float, float, float, float, float>::type boost::math::legendre_p(int,
>>>>>> int, T, const Policy&) [with T = __gmp_expr<__mpf_struct [1],
>>>>>> __mpf_struct [1]>, Policy =
>>>>>> boost::math::policies::policy<boost::math::policies::default_policy,
>>>>>> boost::math::policies::default_policy,
>>>>>> boost::math::policies::default_policy,
>>>>>> boost::math::policies::default_policy,
>>>>>> boost::math::policies::default_policy,
>>>>>> boost::math::policies::default_policy,
>>>>>> boost::math::policies::default_policy,
>>>>>> boost::math::policies::default_policy,
>>>>>> boost::math::policies::default_policy,
>>>>>> boost::math::policies::default_policy,
>>>>>> boost::math::policies::default_policy,
>>>>>> boost::math::policies::default_policy,
>>>>>> boost::math::policies::default_policy>]’
>>>>>> /usr/include/boost/math/special_functions/legendre.hpp:185:
>>>>>> instantiated from ‘typename boost::math::tools::promote_args<RT,
>>>>>> float, float, float, float, float>::type boost::math::legendre_p(int,
>>>>>> int, T) [with T = __gmp_expr<__mpf_struct [1], __mpf_struct [1]>]’
>>>>>> /home/benjamik/fenics/dolfin-wells_gmp/dolfin/math/Legendre.cpp:42:
>>>>>> instantiated from here
>>>>>> /usr/include/boost/math/special_functions/legendre.hpp:167: error: no
>>>>>> matching function for call to ‘pow(__gmp_expr<__mpf_struct [1],
>>>>>> __gmp_binary_expr<long int, __gmp_expr<__mpf_struct [1],
>>>>>> __gmp_binary_expr<__gmp_expr<__mpf_struct [1], __mpf_struct [1]>,
>>>>>> __gmp_expr<__mpf_struct [1], __mpf_struct [1]>,
>>>>>> __gmp_binary_multiplies> >, __gmp_binary_minus> >,
>>>>>> __gmp_expr<__mpf_struct [1], __gmp_binary_expr<__gmp_expr<__mpf_struct
>>>>>> [1], __mpf_struct [1]>, long int, __gmp_binary_divides> >)’
>>>>>> /usr/include/bits/mathcalls.h:154: note: candidates are: double
>>>>>> pow(double, double)
>>>>>> /usr/include/c++/4.4/cmath:358: note:                 float
>>>>>> std::pow(float, float)
>>>>>> /usr/include/c++/4.4/cmath:362: note:                 long double
>>>>>> std::pow(long double, long double)
>>>>>> /usr/include/c++/4.4/cmath:369: note:                 double
>>>>>> std::pow(double, int)
>>>>>> /usr/include/c++/4.4/cmath:373: note:                 float std::pow(float, int)
>>>>>> /usr/include/c++/4.4/cmath:377: note:                 long double
>>>>>> std::pow(long double, int)
>>>>>> [...]
>>>>>>
>>>>>> boost::math::legendre seems to rely on std::pow which is not
>>>>>> templated, only implemented with the most common types.
>>>>>
>>>>> If it's not possible to make it work, we need to revert back.
>>>>
>>>> I don't know of any solution to this. This is the same problem that we
>>>> discussed some months back (then related to Armadillo): Templated
>>>> libraries which rely on non-templated  code (often old and implemented
>>>> i c), so they only support the types which these underlying libraries
>>>> can handle. I think the only solution here is a change in
>>>> boost::math::Legendre.
>>>>
>>>> Of course another solution would be to split the ODE solver from
>>>> Dolfin and let it continue as a separate project, and then import code
>>>> from that when we are going to look at automation/generating code for
>>>> time-dependent problems.
>>>
>>> Yes, perhaps it's time for that. Since it is going to be removed soon
>>> (and replaced by code generation), the best option might be to remove
>>> it before the release of 1.0.
>>>
>>> Are there any objections? Is anyone using the ODE solvers?
>>>
>>
>> No objection, I think that it's a good idea.
>>
>> Once the ODE solvers are out, we can re-design the arbitrary precision
>> interface.
> 
> Good. Just out of curiosity: What is the need of arbitrary precision
> then? 

I don't have anything in particular in mind, but it is a nice feature to
have.

Garth

> I thought the ODE solvers were the only reason for supporting
> that.
> 
> Benjamin
> 
>>
>> Garth
>>
>>> (They will make a comeback later in new form.)
>>>
>>> --
>>> Anders
>>
>>