Thank you and sorry about the wrong group post. Do you happen to know, how 
I can crosspost it?

Tobias

On Thursday, January 14, 2016 at 3:35:34 PM UTC, Aaron Meurer wrote:
>
> You may want to crosspost this to the mpmath list. 
>
> Aaron Meurer 
>
> On Thu, Jan 14, 2016 at 10:14 AM, Tobias Hartung 
> <[email protected] <javascript:>> wrote: 
> > Hi everyone, 
> > 
> > I am trying to integrate polynomials on an n-sphere and implemented (see 
> > attached) the algorithm published in Alan Genz "Fully Symmetric 
> > Interpolatory Rules for Multiple Integrals over Hyper-Spherical 
> Surfaces" 
> > Journal of Computational and Applied Mathematics 157: 187-195, 2003. On 
> > double precision, everything is fine; however, I lose precision in some 
> > integrals using more than double prec. 
> > 
> > I am using randomized samples on the sphere which (by construction) 
> should 
> > integrate the polynomials (that I tested) exactly. However, integrating 
> all 
> > polynomials up to degree 6 on a 5-dimensional sphere with points that 
> should 
> > integrate all polynomials up to degree 7 exactly, I obtain the following 
> > test results : 
> > 
> > mp.prec = 1024 (aka mp.dps = 307) 
> > 
> > polynomial  |  degree of polynomial  |  log_10 ( difference of two 
> > randomized integrations )  |  log_10 ( relative error of weights )  | 
> > analytic value of integral 
> > 
> > [0, 0, 0, 0, 0, 0]      0      -inf                     -306.962262 
> > 31.0062767 
> > [1, 0, 0, 0, 0, 0]      1      -308.111506      -306.962262      0.0 
> > [2, 0, 0, 0, 0, 0]      2      -14.519354        -306.962262     
>  5.16771278 
> > [1, 1, 0, 0, 0, 0]      2      -14.8787358      -306.962262      0.0 
> > [3, 0, 0, 0, 0, 0]      3      -308.583548      -306.962262      0.0 
> > [2, 1, 0, 0, 0, 0]      3      -309.157806      -306.962262      0.0 
> > [1, 1, 1, 0, 0, 0]      3      -309.759866      -306.962262      0.0 
> > [4, 0, 0, 0, 0, 0]      4      -14.6442927      -306.962262     
>  1.93789229 
> > [3, 1, 0, 0, 0, 0]      4      -15.3047045      -306.962262      0.0 
> > [2, 2, 0, 0, 0, 0]      4      -15.0673319      -306.962262     
>  0.645964098 
> > [2, 1, 1, 0, 0, 0]      4      -16.0594155      -306.962262      0.0 
> > [1, 1, 1, 1, 0, 0]      4      -31.8146065      -306.962262      0.0 
> > [5, 0, 0, 0, 0, 0]      5      -309.37189        -306.962262      0.0 
> > [4, 1, 0, 0, 0, 0]      5      -309.481766      -306.962262      0.0 
> > [3, 2, 0, 0, 0, 0]      5      -309.739262      -306.962262      0.0 
> > [3, 1, 1, 0, 0, 0]      5      -310.095075      -306.962262      0.0 
> > [2, 2, 1, 0, 0, 0]      5      -310.474598      -306.962262      0.0 
> > [2, 1, 1, 1, 0, 0]      5      -310.714154      -306.962262      0.0 
> > [1, 1, 1, 1, 1, 0]      5      -311.46643        -306.962262      0.0 
> > [6, 0, 0, 0, 0, 0]      6      -14.7692314      -306.962262     
>  0.968946146 
> > [5, 1, 0, 0, 0, 0]      6      -15.6057345      -306.962262      0.0 
> > [4, 2, 0, 0, 0, 0]      6      -15.4314091      -306.962262     
>  0.193789229 
> > [3, 3, 0, 0, 0, 0]      6      -15.8275833      -306.962262      0.0 
> > [4, 1, 1, 0, 0, 0]      6      -16.5822942      -306.962262      0.0 
> > [3, 2, 1, 0, 0, 0]      6      -16.4753651      -306.962262      0.0 
> > [2, 2, 2, 0, 0, 0]      6      -16.0997045      -306.962262 
> > 0.0645964098 
> > [3, 1, 1, 1, 0, 0]      6      -32.3374853      -306.962262      0.0 
> > [2, 2, 1, 1, 0, 0]      6      -17.3552924      -306.962262      0.0 
> > [2, 1, 1, 1, 1, 0]      6      -33.4100353      -306.962262      0.0 
> > [1, 1, 1, 1, 1, 1]      6      -49.1751847      -306.962262      0.0 
> > 
> > It should be noted that the polynomial are written in the following 
> form: a 
> > list p of length n represents the polynomial 
> > 
> > z_0^{p[0]}  z_1^{p[1]}  z_2^{p[2]}  ...  z_{n-1}^{p[n-1]}. 
> > 
> > If we look at the second to last column, we can see that the weights add 
> are 
> > exact with respect to mp.prec=1024 but do not add up to 1.0 exactly. 
> Hence, 
> > the loss in precision is unlikely to be caused by the weights being only 
> on 
> > double prec (and even if it is, then I don't understand why because they 
> are 
> > calculated with mp.prec=1024, as well). 
> > 
> > The really interesting column is the middle one. Here, I took two 
> randomized 
> > set of integration points (both should integrate all polynomials exactly 
> up 
> > to machine error) and printed the decadic logarithm of the absolute 
> value of 
> > their difference. In other words, we are expecting the middle column to 
> be 
> > populated with numbers in the vicinity of -307 (just like the fourth 
> > column). However, we only get these values for odd degree polynomials 
> (and 
> > the constant 1 but that is just the sum of weights, i.e., the 
> randomization 
> > doesn't do anything). For even degree polynomials, we have significantly 
> > reduced precision. 
> > 
> > Does anyone have an idea what might be the reason for such behavior? 
> > 
> > Thank you very much, 
> > Tobias 
> > 
> > -- 
> > You received this message because you are subscribed to the Google 
> Groups 
> > "sympy" group. 
> > To unsubscribe from this group and stop receiving emails from it, send 
> an 
> > email to [email protected] <javascript:>. 
> > To post to this group, send email to [email protected] 
> <javascript:>. 
> > Visit this group at https://groups.google.com/group/sympy. 
> > To view this discussion on the web visit 
> > 
> https://groups.google.com/d/msgid/sympy/cf76f223-96a0-4789-8530-db224622f6bb%40googlegroups.com.
>  
>
> > For more options, visit https://groups.google.com/d/optout. 
>

-- 
You received this message because you are subscribed to the Google Groups 
"sympy" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to [email protected].
To post to this group, send email to [email protected].
Visit this group at https://groups.google.com/group/sympy.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/sympy/015dd6ad-d2e8-495c-82f7-1c157b0c7ab0%40googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Reply via email to