headers
Charles Karney | 5 Feb 16:06
Favicon

Re: Exact sampling from the normal distribution

I've implemented MPFR algorithms for sampling from the unit exponential
and unit normal distributions.  The code is in

 
http://randomlib.git.sourceforge.net/git/gitweb.cgi?p=randomlib/randomlib;a=blob;f=examples/RandomMPFR.hpp

This is a header-only implementation via C++ classes.
MPFR{Test,Check}.cpp exercise the implementation.  The three files
RandomMPFR.hpp and MPFR{Test,Check}.cpp are self-contained.  You can
download these and compile and link MPFR{Test,Check}.cpp against MPFR
(version 3.x?).  Converting this to C code should be relatively easy.
RandomMPFR.hpp contains some implementation details at the top.  Here I
raise a few points/questions:

* I explicitly handle rounding, overflow, and underflow.  Partly this is
   because I am confused by the MPFR documentation on the rounding.  It
   is claimed that the ternary value with MPFR_RNDU should always be 0 or
   1.  However, surely that's not the case if there's underflow (from a
   positive number to 0) or overflow (to -inf)?

* The implementation of normal sampling is a *lot* faster than
   mpfr_grandom (MPFR version 3.1.0).  Here's a table showing the times
   (microsec per sample) obtained on an Intel Xeon, x86_42, 2.66GHz
   computer with g++ version 4.6.1 (-O3).  The columns are

    1. RandomLib::NormalDistribution (using the ratio method)
    2. RandomLib::ExactNormal (same algorithm as RandomMPFR.hpp)
    3. RandomMPFR.hpp
    4. mpfr_grandom (MPFR version 3.1.0)
(Continue reading)

Permalink | Reply | 4 comments |
headers
Damien Stehlé | 8 Feb 09:26
Picon

Re: Exact sampling from the normal distribution

Hello Charles,

OK, thank you for the precisions. This is very valuable.
I'll think some more over it, and may come back to you later.

Regards,
Damien

2012/2/5 Charles Karney <charles.karney@...>:
> I've implemented MPFR algorithms for sampling from the unit exponential
> and unit normal distributions.  The code is in
>
>
> http://randomlib.git.sourceforge.net/git/gitweb.cgi?p=randomlib/randomlib;a=blob;f=examples/RandomMPFR.hpp
>
> This is a header-only implementation via C++ classes.
> MPFR{Test,Check}.cpp exercise the implementation.  The three files
> RandomMPFR.hpp and MPFR{Test,Check}.cpp are self-contained.  You can
> download these and compile and link MPFR{Test,Check}.cpp against MPFR
> (version 3.x?).  Converting this to C code should be relatively easy.
> RandomMPFR.hpp contains some implementation details at the top.  Here I
> raise a few points/questions:
>
> * I explicitly handle rounding, overflow, and underflow.  Partly this is
>  because I am confused by the MPFR documentation on the rounding.  It
>  is claimed that the ternary value with MPFR_RNDU should always be 0 or
>  1.  However, surely that's not the case if there's underflow (from a
>  positive number to 0) or overflow (to -inf)?
>
> * The implementation of normal sampling is a *lot* faster than
(Continue reading)

Permalink | Reply | 3 comments |
headers
Charles Karney | 8 Feb 11:10
Favicon

Re: Exact sampling from the normal distribution

On 02/08/12 03:26, Damien Stehlé wrote:
> Hello Charles,
>
> OK, thank you for the precisions. This is very valuable.
> I'll think some more over it, and may come back to you later.

Damien,

A trivial change in the code (adding all the bits needed to get to the
target precision in one shot) gives a further improvement in the times.

               time(us)
     prec     old    new
       8      1.1    1.1
      24      1.1    1.1
      53      1.1    1.1
      64      1.1    1.1
     2^8      1.4    1.2
     2^10     2.2    1.4
     2^12     6.2    2.2
     2^14      32    5.4
     2^16     290     18
     2^18             76
     2^20            290

Now the cost is very close to offset-linear -- approximately the time is
(1.1 + 2.75e-4 * prec) us per sample.  Presumably the time could be
improved for low precisions by working in units of 32 random bits
(instead of 64).  But that would complicate the code somewhat.
(Continue reading)

Permalink | Reply | 2 comments |
headers
keith.briggs | 16 Feb 18:21
Favicon

RE: MPFI-1.5.1 is released

Am I confused or does MPFI really have a contradictory definition of the comparison functions?

MPFR: Compare op1 and op2. Return a positive value if op1 > op2, zero if op1 = op2, and a negative value if op1 <
op2. (http://www.mpfr.org/mpfr-current/mpfr.html )

MPFI: Compares op1 and op2. Return a positive value if op1 < op2, zero if op1 overlaps or contains op2, and a
negative value if op1 > op2. (http://perso.ens-lyon.fr/nathalie.revol/mpfi.html )

Keith

-----Original Message-----
From: philipe theveny [mailto:philippe.theveny@...] 
Sent: 19 January 2012 10:10
To: mpfr@...
Subject: [MPFR] MPFI-1.5.1 is released

Hi,

MPFR users may be interested by the new MPFI version which is mostly a bug-fix release. You can read the
announce here: 
http://lists.gforge.inria.fr/pipermail/mpfi-users/2012-January/000046.html

The MPFI library is based on MPFR and provides interval arithmetic and elementary functions in arbitrary precision.

Best Regards,
The MPFI development team
Permalink | Reply | 6 comments |
headers
Sylvain CHEVILLARD | 16 Feb 18:30
Picon
Picon
Favicon

Re: MPFI-1.5.1 is released

According to a test that I just performed, I think that this is a typo 
in the documentation of MPFI.
Both MPFR and MPFI seem to return a positive value if op1 > op2.

@Nathalie: could you please confirm and correct the documentation?

Sylvain

Le 16/02/2012 18:21, keith.briggs@... a écrit :
> Am I confused or does MPFI really have a contradictory definition of the comparison functions?
>
> MPFR: Compare op1 and op2. Return a positive value if op1>  op2, zero if op1 = op2, and a negative value if op1< 
op2. (http://www.mpfr.org/mpfr-current/mpfr.html )
>
> MPFI: Compares op1 and op2. Return a positive value if op1<  op2, zero if op1 overlaps or contains op2, and a
negative value if op1>  op2. (http://perso.ens-lyon.fr/nathalie.revol/mpfi.html )
>
> Keith
>
>
> -----Original Message-----
> From: philipe theveny [mailto:philippe.theveny@...]
> Sent: 19 January 2012 10:10
> To: mpfr@...
> Subject: [MPFR] MPFI-1.5.1 is released
>
> Hi,
>
> MPFR users may be interested by the new MPFI version which is mostly a bug-fix release. You can read the
announce here:
(Continue reading)

Permalink | Reply | 1 comment |
headers
Vincent Lefevre | 16 Feb 18:40

Re: MPFI-1.5.1 is released

On 2012-02-16 18:30:01 +0100, Sylvain CHEVILLARD wrote:
> According to a test that I just performed, I think that this is a typo in
> the documentation of MPFI.
> Both MPFR and MPFI seem to return a positive value if op1 > op2.
> 
> @Nathalie: could you please confirm and correct the documentation?

Note that the MPFI manual (in the Debian package) is correct.

--

-- 
Vincent Lefèvre <vincent@...> - Web: <http://www.vinc17.net/>
100% accessible validated (X)HTML - Blog: <http://www.vinc17.net/blog/>
Work: CR INRIA - computer arithmetic / AriC project (LIP, ENS-Lyon)
Permalink | Reply | 0 comments |
headers
Nathalie Revol | 16 Feb 19:06
Picon
Picon

Re: MPFI-1.5.1 is released


Dear Keith

ouch, the extremely old version of MPFI (with bugged documentation
and probably bugged implementation) is still available at the URL you
indicate.

I just removed it, please access to MPFI here:
https://gforge.inria.fr/frs/?group_id=157

Sorry for this.

And yes, MPFR and MPFI have compatible definition of the comparison,
which is the one indicated for MPFR in your message below.

Again, sorry for the trouble.
	Nathalie

On Feb 16, 2012, at 6:21 PM, <keith.briggs@...>
<keith.briggs@...> wrote:

> Am I confused or does MPFI really have a contradictory definition of the comparison functions?
> 
> MPFR: Compare op1 and op2. Return a positive value if op1 > op2, zero if op1 = op2, and a negative value if op1 <
op2. (http://www.mpfr.org/mpfr-current/mpfr.html )
> 
> MPFI: Compares op1 and op2. Return a positive value if op1 < op2, zero if op1 overlaps or contains op2, and a
negative value if op1 > op2. (http://perso.ens-lyon.fr/nathalie.revol/mpfi.html )
> 
> Keith
(Continue reading)

Permalink | Reply | 3 comments |
headers
Sylvain CHEVILLARD | 16 Feb 19:14
Picon
Picon
Favicon

Re: MPFI-1.5.1 is released

Dear Nathalie,

is there somewhere a HTML version of the documentation of the current 
version of MPFI? The one pointed by Keith, though very old, was very 
useful for a quick search over the doc while coding with MPFI.

Thanks,

Sylvain

Le 16/02/2012 19:06, Nathalie Revol a écrit :
>
> Dear Keith
>
> ouch, the extremely old version of MPFI (with bugged documentation
> and probably bugged implementation) is still available at the URL you
> indicate.
>
> I just removed it, please access to MPFI here:
> https://gforge.inria.fr/frs/?group_id=157
>
> Sorry for this.
>
> And yes, MPFR and MPFI have compatible definition of the comparison,
> which is the one indicated for MPFR in your message below.
>
> Again, sorry for the trouble.
> 	Nathalie
>
> On Feb 16, 2012, at 6:21 PM,<keith.briggs@...>
(Continue reading)

Permalink | Reply | 0 comments |
headers
keith.briggs | 16 Feb 20:14
Favicon

RE: MPFI-1.5.1 is released

Dear Nathalie and others,

Thank you for the replies.    I just sent the link to the old documentation for simplicity, but I had already
checked that the current documentation is also wrong.  It says
(file:///mpfi-1.5.1/doc/mpfi.html/Interval-Comparison.html ):

The default comparison functions return a positive value if the first interval has all its elements
strictly lower than all elements of the second one, a negative value if the first interval has all its
elements strictly greater than all elements of the second one and 0 otherwise, i.e. if they overlap or if
one is contained in the other. 

Keith
________________________________________
From: Nathalie Revol [Nathalie.Revol@...]
Sent: 16 February 2012 18:06
To: Briggs,KM,Keith,DUB8 R
Cc: mpfr@...
Subject: Re: [MPFR] MPFI-1.5.1 is released

Dear Keith

ouch, the extremely old version of MPFI (with bugged documentation
and probably bugged implementation) is still available at the URL you
indicate.

I just removed it, please access to MPFI here:
https://gforge.inria.fr/frs/?group_id=157

Sorry for this.
(Continue reading)

Permalink | Reply | 1 comment |
headers
philipe theveny | 20 Feb 11:03
Picon
Picon

Re: MPFI-1.5.1 is released

Dear Keith,

> Thank you for the replies.    I just sent the link to the old documentation for simplicity, but I had already
checked that the current documentation is also wrong.  It says
(file:///mpfi-1.5.1/doc/mpfi.html/Interval-Comparison.html ):
>
> The default comparison functions return a positive value if the first interval has all its elements
strictly lower than all elements of the second one, a negative value if the first interval has all its
elements strictly greater than all elements of the second one and 0 otherwise, i.e. if they overlap or if
one is contained in the other.

The MPFI documentation is clearly wrong here. I fixed it in the svn 
repository with the revision r654. The introduction of the section now 
matches the documentation of mpfi_cmp.

Thanks you for your bug report.
Philippe

MPFI users, please notice that there is also a mailing list dedicated to 
mpfi: mpfi-users@...
Permalink | Reply | 0 comments |

Gmane