diff --git a/build/Jamfile.v2 b/build/Jamfile.v2
index 14439ad42..15966862c 100644
--- a/build/Jamfile.v2
+++ b/build/Jamfile.v2
@@ -16,7 +16,7 @@ project
#gcc:-fvisibility=hidden
intel-linux:-fvisibility=hidden
sun:-xldscope=hidden
- [ check-target-builds ../config//has_gcc_visibility : gcc:-fvisibility=hidden : ]
+ [ check-target-builds ../config//has_gcc_visibility "gcc visibility" : gcc:-fvisibility=hidden : ]
;
cpp-pch pch : ../src/tr1/pch.hpp : ../src/tr1 shared:BOOST_MATH_TR1_DYN_LINK=1 ;
@@ -96,7 +96,7 @@ lib boost_math_tr1l : ../src/tr1/$(TR1_SOURCES)l.cpp pch
shared:BOOST_MATH_TR1_DYN_LINK=1
../config//has_long_double_support
../src/tr1
- [ check-target-builds ../config//has_long_double_support : : no ]
+ [ check-target-builds ../config//has_long_double_support "long double support" : : no ]
;
lib boost_math_c99 : ../src/tr1/$(C99_SOURCES).cpp pch
@@ -116,7 +116,7 @@ lib boost_math_c99l : ../src/tr1/$(C99_SOURCES)l.cpp pch
shared:BOOST_MATH_TR1_DYN_LINK=1
../config//has_long_double_support
../src/tr1
- [ check-target-builds ../config//has_long_double_support : : no ]
+ [ check-target-builds ../config//has_long_double_support "long double support" : : no ]
;
boost-install boost_math_c99 boost_math_c99f boost_math_c99l boost_math_tr1 boost_math_tr1f boost_math_tr1l ;
diff --git a/doc/sf_and_dist/constants.qbk b/doc/sf_and_dist/constants.qbk
index f194caadc..27377934a 100644
--- a/doc/sf_and_dist/constants.qbk
+++ b/doc/sf_and_dist/constants.qbk
@@ -1,75 +1,710 @@
-[section:constants Numeric Constants]
+[section:constants Mathematical Constants]
-[h4 Synopsis]
+[section:intro Introduction]
-``
-#include
-``
+Boost.Math provides a collection of mathematical constants.
- namespace boost{ namespace math{ namespace constants{
-
- template T pi();
- template T root_pi();
- template T root_half_pi();
- template T root_two_pi();
- template T root_ln_four();
- template T e();
- template T half();
- template T euler();
- template T root_two();
- template T ln_two();
- template T ln_ln_two();
- template T third();
- template T twothirds();
- template T pi_minus_three();
- template T four_minus_pi();
-
- }}} // namespaces
+[h4 Why use Boost.Math mathematical constants?]
-[h4 Description]
+* Readable. For the very many jobs just using built-in like `double`, you can just write expressions like
+``double area = pi * r * r;``
+(If that's all you want, jump direct to [link math_toolkit.constants.tutorial.non_templ use in non-template code]!)
+* Effortless - avoiding a search of reference sources.
+* Usable with both builtin floating point types, and user-defined, possibly extended precision, types such as
+NTL, MPFR/GMP, mp_float: in the latter case the constants are computed to the necessary precision and then cached.
+* Accurate - ensuring that the values are as accurate as possible for the
+chosen floating-point type
+ * No loss of accuracy from repeated rounding of intermediate computations.
+ * Result is computed with higher precision and only rounded once.
+ * Less risk of inaccurate result from functions pow, trig and log at [@http://en.wikipedia.org/wiki/Corner_case corner cases].
+ * Less risk of [@http://docs.oracle.com/cd/E19957-01/806-3568/ncg_goldberg.html cancellation error].
+* Faster - can avoid (re-)calculation at runtime. This can be significant if:
+ * Functions pow, trig or log are used.
+ * Inside an inner loop.
+ * Using a high-precision UDT.
+ * Compiler optimizations possible with built-in types, especially `double`, are not available.
+* Portable - as possible between different systems using different floating-point precisions:
+see [link math_toolkit.constants.tutorial.templ use in template code].
+* Tested - by comparison with other published sources, or separately computed at long double precision.
-The header `boost/math/constants/constants.hpp` contains some numeric constants
-that we have found useful in the development of this library. New constants are
-added on an ad-hoc basis based on need.
+[endsect] [/section:intro Introduction]
-Usage is like this:
-
- template
- T circumference(T r)
+[section:tutorial Tutorial]
+
+[section:non_templ Use in non-template code]
+
+When using the math constants at your chosen fixed precision in non-template code,
+you can simply add a `using` declaration, for example, `using boost::math::double_constants`,
+to make the constants of the correct precision for your code
+visible in the current scope, and then use each constant ['as a simple variable]:
+
+ #include
+
+ double area(double r)
{
- return 2 * boost::math::constants::pi() * r;
+ using boost::math::double_constants;
+ return pi * r * r;
}
-All the constants are accurate to at least the 34 decimal digits required for 128-bit
-long doubles, and most are accurate to 100 digits or more when used with
-a suitable arbitrary precision type.
+Had our function been written as taking a `float` rather than a `double`,
+we could have written instead:
-The following table summarises the constants we have at present:
+ #include
+
+ float area(float r)
+ {
+ using boost::math::float_constants;
+ return pi * r * r;
+ }
+
+Likewise, constants that are suitable for use at `long double` precision
+are available in the namespace `boost::math::long_double_constants`.
+
+You can see the full list of available constants at [link math_toolkit.constants.constants].
+
+[endsect]
+
+[section:templ Use in template code]
+
+When using the constants inside a function template, we need to ensure that
+we use a constant of the correct precision for our template parameters.
+We can do this by calling the function-template versions, `pi()`, of the constants
+like this:
+
+ #include
+
+ template
+ Real area(Real r)
+ {
+ using namespace boost::math::constants;
+ return pi() * r * r;
+ }
+
+Although this syntax is a little less "cute" than the non-template version,
+the code is no less efficient
+(at least for the built-in types `float`, `double` and `long double`) :
+the function template versions of the constants are simple inline functions that
+return a constant of the correct precision for the type used. In addition, these
+functions are declared `constexp` for those compilers that support this, allowing
+the result to be used in constant-expressions provided the template argument is a literal type.
+
+[tip Keep in mind the difference between the variable version,
+just `pi`, and the template-function version:
+the template-function requires both a <[~floating-point-type]>
+and function call `()` brackets, for example: `pi()`.
+You cannot write `double p = pi<>()`, nor `double p = pi()`.]
+
+[note You can always use [*both] variable and template-function versions
+[*provided calls are fully qualified], for example:
+``
+double my_pi1 = boost::math::constants::pi();
+double my_pi2 = boost::math::double_constants::pi;
+``
+]
+
+[warning It may be tempting to simply define
+``
+using namespace boost::math::double_constants;
+using namespace boost::math::constants;
+``
+but if you do define two namespaces, this will, of course, create ambiguity!
+``
+double my_pi = pi(); // error C2872: 'pi' : ambiguous symbol
+double my_pi2 = pi; // Context does not allow for disambiguation of overloaded function
+``
+Although the mistake above is fairly obvious,
+it is also not too difficult to do this accidentally, or worse, create it in someone elses code.
+
+Therefore is it prudent to avoid this risk by [*localising the scope of such definitions], as shown above.]
+
+[tip Be very careful with the type provided as parameter.
+For example, providing an [*integer] instead of a floating-point type can be disastrous (a C++ feature).
+
+``cout << "Area = " << area(2) << endl; // Area = 12!!!``
+
+You should get a compiler warning
+[pre
+warning : 'return' : conversion from 'double' to 'int', possible loss of data
+] [/pre]
+Failure to heed this warning can lead to very wrong answers!
+
+You can also avoid this by being explicit about the type of `Area`.
+``cout << "Area = " << area(2) << endl; // Area = 12.566371``
+]
+
+[endsect] [/section:templ Use in template code]
+
+[section:user_def Use With User Defined Types]
+
+The syntax for using the function-call constants with user-defined types is the same
+as it is in the template class, which is to say we use:
+
+ #include
+
+ boost::math::constants::pi();
+
+However, since the precision of the user-defined type may be much greater than that
+of the built-in floating pointer types, how the value returned is created is as follows:
+
+* If the precision of the type is known at compile time:
+ * If the precision is less than or equal to that of a `float` and the type is constructable from a `float`
+ then our code returns a `float` literal. If the user-defined type is a literal type
+ then the function call that returns the constant will be a `constexp`.
+ * If the precision is less than or equal to that of a `double` and the type is constructable from a `double`
+ then our code returns a `double` literal. If the user-defined type is a literal type
+ then the function call that returns the constant will be a `constexp`.
+ * If the precision is less than or equal to that of a `long double` and the type is constructable from a `long double`
+ then our code returns a `long double` literal. If the user-defined type is a literal type
+ then the function call that returns the constant will be a `constexp`.
+ * If the precision is less than 100 decimal digits, then the constant will be constructed
+ (just the once, then cached in a thread-safe manner) from a string representation of the constant.
+ * Otherwise the value is computed (just once, then cached in a thread-safe manner).
+* If the precision is unknown at compile time then:
+ * If the runtime precision (obtained from a call to `boost::math::tools::digits()`) is
+ less than 100 decimal digits, then the constant is constructed "on the fly" from the string
+ representation of the constant.
+ * Otherwise the value is constructed "on the fly" by calculating then value of the constant
+ using the current default precision of the type. Note that this can make use of the constants
+ rather expensive.
+
+In addition, it is possible to pass a `Policy` type as a second template argument, and use this to control
+the precision:
+
+ #include
+
+ typedef boost::math::policies::policy > my_policy_type;
+ boost::math::constants::pi();
+
+[note Boost.Math doesn't know how to control the internal precision of `MyType`, the policy
+just controls how the selection process above is carried out, and the calculation precision
+if the result is computed.]
+
+It is also possible to control which method is used to construct the constant by specialising
+the traits class `construction_traits`:
+
+ namespace boost{ namespace math{ namespace constant{
+
+ template
+ struct construction_traits
+ {
+ typedef mpl::int_ type;
+ };
+
+ }}} // namespaces
+
+Where ['N] takes one of the following values:
[table
-[[Constant][Meaning][Value]]
-[[pi][[pi]][3.1415926535897932384...]]
-[[root_pi][[radic][pi]][1.772453850905516027...]]
-[[root_half_pi][[radic]([pi]/2)][1.253314137315500251...]]
-[[root_two_pi][[radic](2*[pi])][2.506628274631000502...]]
-[[root_ln_four][[radic](ln(4))][1.17741002251547469...]]
-[[e][['e]][2.71828182845904523536...]]
-[[half][0.5][0.5]]
-[[euler][Euler's constant][0.577215664901532860606]]
-[[root_two][[radic]2][1.4142135623730950488...]]
-[[ln_two][ln(2)][0.6931471805599453094...]]
-[[ln_ln_two][ln(ln(2))][-0.3665129205816643...]]
-[[third][1/3][0.333333333333333333...]]
-[[twothirds][2/3][0.666666666666666666...]]
-[[pi_minus_three][[pi]-3][0.14159265358979323846...]]
-[[four_minus_pi][4-[pi]][0.85840734641020676153735...]]
+[[['N]][Meaning]]
+[[0][The precision is unavailable at compile time;
+either construct from a decimal digit string or calculate on the fly depending upon the runtime precision.]]
+[[1][Return a float precision constant.]]
+[[2][Return a double precision constant.]]
+[[3][Return a long double precision constant.]]
+[[4][Construct the result from the string representation, and cache the result.]]
+[[Any other value ['N]][Sets the compile time precision to ['N] bits.]]
]
-[endsect][/section Numeric Constants]
+[h5 Custom Specializing a constant]
-[/
- Copyright 2008 John Maddock and Paul A. Bristow.
- Distributed under the Boost Software License, Version 1.0.
- (See accompanying file LICENSE_1_0.txt or copy at
- http://www.boost.org/LICENSE_1_0.txt).
+In addition, for user-defined types that need special handling, it's possible to \[partially-\] specialize
+the internal structure used by each constant. For example, suppose we're using the C++ wrapper around MPFR
+`mpfr_class`: this has its own representation of Pi which we may well wish to use in place of the above
+mechanism. We can achieve this by specialising the class template `boost::math::constants::detail::constant_pi`:
+
+ namespace boost{ namespace math{ namespace constants{ namespace detail{
+
+ template<>
+ struct constant_pi
+ {
+ template
+ inline T get(const mpl::int_&)
+ {
+ // The template param N is one of the values in the table above,
+ // we can either handle all cases in one as is the case here,
+ // or overload "get" for the different options.
+ mpfr_class result;
+ mpfr_const_pi(result.get_mpfr_t(), GMP_RNDN);
+ return result;
+ }
+ };
+
+ }}}} // namespaces
+
+[h5 Diagnosing what meta-programmed code is doing]
+
+Finally, since it can be tricky to diagnose what meta-programmed code is doing, there is a
+diagnostic routine that prints information about how this library will handle a specific type,
+it can be used like this:
+
+ #include
+
+ int main()
+ {
+ boost::math::constants::print_info_on_type();
+ }
+
+If you wish, you can also pass an optional std::ostream argument to the `print_info_on_type` function.
+Typical output for a user-defined type looks like this:
+
+[pre
+Information on the Implementation and Handling of
+Mathematical Constants for Type class boost::math::concepts::real_concept
+
+Checking for std::numeric_limits specialisation: no
+boost::math::policies::precision
+reports that there is no compile type precision available.
+boost::math::tools::digits()
+reports that the current runtime precision is
+53 binary digits.
+No compile time precision is available, the construction method
+will be decided at runtime and results will not be cached
+- this may lead to poor runtime performance.
+Current runtime precision indicates that
+the constant will be constructed from a string on each call.
]
+
+[endsect] [/section:user_def Use With User Defined Types]
+
+[endsect] [/section:tutorial Tutorial]
+
+[section:constants The Mathematical Constants]
+
+This section lists the mathematical constants, their use(s) (and sometimes rationale for their inclusion).
+[table Mathematical Constants
+[[name] [formula] [Value (6 decimals)] [Uses and Rationale]]
+[[[*Rational fractions]] [] [] [] ]
+[[half] [1/2] [0.5] [] ]
+[[third] [1/3] [0.333333] [] ]
+[[two_thirds] [2/3] [0.66667] [] ]
+[[three_quarters] [3/4] [0.75] [] ]
+
+[[[*two and related]] [] [] [] ]
+[[root_two] [[radic]2] [1.41421] [] ]
+[[root_three] [[radic]3] [1.73205] [] ]
+[[half_root_two] [[radic]2 /2] [0.707106] [] ]
+[[ln_two] [ln(2)] [0.693147] [] ]
+[[ln_ten] [ln(10)] [2.30258] [] ]
+[[ln_ln_two] [ln(ln(2))] [-0.366512] [Gumbel distribution median] ]
+[[root_ln_four] [[radic]ln(4)] [1.177410] [] ]
+[[one_div_root_two] [1/[radic]2] [0.707106] [] ]
+
+[[[*[pi] and related]] [] [] [] ]
+[[pi] [pi] [3.14159] [Ubiquitous. Archimedes constant [@http://en.wikipedia.org/wiki/Pi [pi]]]]
+[[half_pi] [[pi]/2] [1.570796] [] ]
+[[third_pi] [[pi]/3] [1.04719] [] ]
+[[sixth_pi] [[pi]/6] [0.523598] [] ]
+[[two_pi] [2[pi]] [6.28318] [Many uses, most simply, circumference of a circle]]
+[[two_thirds_pi] [2/3 [pi]] [2.09439] [[@http://en.wikipedia.org/wiki/Sphere#Volume_of_a_sphere volume of a hemi-sphere] = 4/3 [pi] r[cubed]]]
+[[three_quarters_pi] [3/4 [pi]] [2.35619] [[@http://en.wikipedia.org/wiki/Sphere#Volume_of_a_sphere volume of a hemi-sphere] = 4/3 [pi] r[cubed]]]
+[[four_thirds_pi] [4/3 [pi]] [4.18879] [[@http://en.wikipedia.org/wiki/Sphere#Volume_of_a_sphere volume of a sphere] = 4/3 [pi] r[cubed]]]
+[[one_div_two_pi] [1/(2[pi])] [1.59155] [Widely used]]
+[[root_pi] [[radic][pi]][1.77245] [Widely used]]
+[[root_half_pi] [[radic] [pi]/2] [1.25331] [Widely used]]
+[[root_two_pi][[radic] [pi]*2] [2.50662] [Widely used]]
+[[one_div_root_pi] [1/[radic][pi]] [0.564189] [] ]
+[[one_div_root_two_pi] [1/[radic](2[pi])] [0.398942] [] ]
+[[root_one_div_pi] [[radic](1/[pi]] [0.564189] [] ]
+[[pi_minus_three] [[pi]-3] [1.41593] [] ]
+[[four_minus_pi] [4 -[pi]] [0.858407] [] ]
+[[pow23_four_minus_pi] [4[super 2/3] - [pi]] [0.795316] [] ]
+[[pi_pow_e] [[pi][super e]] [22.4591] [] ]
+
+[[pi_sqr] [[pi][super 2]] [9.86960] [] ]
+[[pi_sqr_div_six] [[pi][super 2]/6] [1.64493] [] ]
+[[pi_cubed] [[pi][super 3]] [31.00627] [] ]
+[[cbrt_pi] [[radic][super 3] [pi]] [1.46459] [] ]
+[[one_div_cbrt_pi] [1/[radic][super 3] [pi]] [0.682784] [] ]
+
+[[[*Euler's e and related]] [] [] [] ]
+[[e] [e] [2.71828] [[@http://en.wikipedia.org/wiki/E_(mathematical_constant) Euler's constant e]] ]
+[[exp_minus_half] [e [super -1/2]] [0.606530] [] ]
+[[e_pow_pi] [e [super [pi]]] [23.14069] [] ]
+[[root_e] [[radic] e] [1.64872] [] ]
+[[log10_e] [log10(e)] [0.434294] [] ]
+[[one_div_log10_e] [1/log10(e)] [2.30258] [] ]
+
+[[[*Trigonometric]] [] [] [] ]
+[[degree] [radians = [pi] / 180] [0.017453] [] ]
+[[radian] [degrees = 180 / [pi]] [57.2957] [] ]
+[[sin_one] [sin(1)] [0.841470] [] ]
+[[cos_one] [cos(1)] [0.54030] [] ]
+[[sinh_one] [sinh(1)] [1.17520] [] ]
+[[cosh_one] [cosh(1)] [1.54308] [] ]
+
+[[[*Phi]] [ Phidias golden ratio] [[@http://en.wikipedia.org/wiki/Golden_ratio Phidias golden ratio]] [] ]
+[[phi] [(1 + [radic]5) /2] [1.61803] [finance] ]
+[[ln_phi] [ln([phi])] [0.48121] [] ]
+[[one_div_ln_phi] [1/ln([phi])] [2.07808] [] ]
+
+[[[*Euler's Gamma]] [] [] [] ]
+[[euler] [euler] [0.577215] [[@http://en.wikipedia.org/wiki/Euler%E2%80%93Mascheroni_constant Euler-Mascheroni gamma constant]] ]
+[[one_div_euler] [1/euler] [1.73245] [] ]
+[[euler_sqr] [euler[super 2]] [0.333177] [] ]
+
+[[[*Misc]] [] [] [] ]
+[[zeta_two] [[zeta](2)] [1.64493] [[@http://en.wikipedia.org/wiki/Riemann_zeta_function Riemann zeta function]] ]
+[[zeta_three] [[zeta](3)] [1.20205] [[@http://en.wikipedia.org/wiki/Riemann_zeta_function Riemann zeta function]] ]
+[[catalan] [['K]] [0.915965] [[@http://mathworld.wolfram.com/CatalansConstant.html Catalan (or Glaisher) combinatorial constant] ]]
+[[glaisher] [['A]] [1.28242] [[@https://oeis.org/A074962/constant Decimal expansion of Glaisher-Kinkelin constant] ]]
+[[khinchin] [['k]] [2.685452] [[@https://oeis.org/A002210/constant Decimal expansion of Khinchin constant] ]]
+
+[[extreme_value_skewness] [12[radic]6 [zeta](3)/ [pi][super 3]] [1.139547] [Extreme value distribution] ]
+[[rayleigh_skewness] [2[radic][pi]([pi]-3)/(4 - [pi])[super 3/2]] [0.631110] [Rayleigh distribution skewness] ]
+[[rayleigh_kurtosis_excess] [-(6[pi][super 2]-24[pi]+16)/(4-[pi])[super 2]] [0.245089] [[@http://en.wikipedia.org/wiki/Rayleigh_distribution Rayleigh distribution kurtosis excess]] ]
+[[rayleigh_kurtosis] [3+(6[pi][super 2]-24[pi]+16)/(4-[pi])[super 2]] [3.245089] [Rayleigh distribution kurtosis] ]
+
+] [/table]
+
+
+[note Integer values are *not included* in this list of math constants, however interesting,
+because they can be so easily and exactly constructed, even for UDT, for example: `static_cast(42)`.]
+
+[tip If you know the approximate value of the constant, you can search for the value to find Boost.Math chosen name in this table.]
+
+[endsect] [/section:constants The constants]
+
+[section:new_const Defining New Constants]
+
+The library provides some helper code to assist in defining new constants;
+the process for defining a constant called `my_constant` goes like this:
+
+1. [*Define a function that calculates the value of the constant].
+This should be a template function, and be placed in `boost/math/constants/calculate_constants.hpp`
+if the constant is to be added to this library,
+or else defined at the top of your source file if not.
+
+The function should look like this:
+
+ namespace boost{ namespace math{ namespace constants{ namespace detail{
+
+ template
+ template
+ Real constant_my_constant::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_))
+ {
+ int required_precision = N ? N : tools::digits();
+ Real result = /* value computed to required_precision bits */ ;
+ return result;
+ }
+
+ }}}} // namespaces
+
+Then define a placeholder for the constant itself:
+
+ namespace boost{ namespace math{ namespace constants{
+
+ BOOST_DEFINE_MATH_CONSTANT(my_constant, 0.0, "0");
+
+ }}}
+
+
+For example, to calculate [pi]/2, add to `boost/math/constants/calculate_constants.hpp`
+
+ template
+ template
+ inline T constant_half_pi::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_))
+ {
+ BOOST_MATH_STD_USING
+ return pi > >() / static_cast(2);
+ }
+
+Then to `boost/math/constants/constants.hpp` add:
+
+ BOOST_DEFINE_MATH_CONSTANT(half_pi, 0.0, "0"); // Actual values are temporary, we'll replace them later.
+
+[note Previously defined constants like pi and e can be used, but by *not simply calling* `pi()`;
+specifying the precision via the policy
+`pi > >()`
+is essential to ensure full accuracy.]
+
+[warning Newly defined constants can only be used once they are included in
+`boost/math/constants/constants.hpp`. So if you add
+`template T constant_my_constant{...}`,
+then you cannot define `constant_my_constant`
+until you add the temporary `BOOST_DEFINE_MATH_CONSTANT(my_constant, 0.0, "0")`.
+Failing to do this will result in surprising compile errors:
+``
+ error C2143: syntax error : missing ';' before '<'
+ error C2433: 'constant_root_two_div_pi' : 'inline' not permitted on data declarations
+ error C2888: 'T constant_root_two_div_pi' : symbol cannot be defined within namespace 'detail'
+ error C2988: unrecognizable template declaration/definition
+``
+]
+
+2. [*You will need an arbitrary precision type to use to calculate the value]. This library
+currently supports either `cpp_float`, `NTL::RR` or `mpfr_class` used via the bindings in `boost/math/bindings`.
+The default is to use `NTL::RR` unless you define an alternate macro, for example,
+`USE_MPFR` or `USE_CPP_FLOAT` at the start of your program.
+
+3. It is necessary to link to the Boost.Regex library,
+and probably to your chosen arbitrary precision type library.
+
+4. The complete program to generate the constant `half_pi` using function `calculate_half_pi` is then:
+
+ #define USE_CPP_FLOAT // If required.
+ #include
+
+ int main()
+ {
+ BOOST_CONSTANTS_GENERATE(half_pi);
+ }
+
+The output from the program is a snippet of C++ code
+(actually a macro call) that can be cut and pasted
+into `boost/math/constants/constants.hpp` or else into your own code, for example:
+
+[pre
+ BOOST_DEFINE_MATH_CONSTANT(half_pi, 1.570796326794896619231321691639751442e+00, "1.57079632679489661923132169163975144209858469968755291048747229615390820314310449931401741267105853399107404326e+00");
+]
+
+This macro BOOST_DEFINE_MATH_CONSTANT inserts a C++ struct code snippet that
+declares the `float`, `double` and `long double` versions of the constant,
+plus a decimal digit string representation correct to 100 decimal
+digits, and all the meta-programming machinery needed to select between them.
+
+The result of an expanded macro for Pi is shown below.
+
+[import ./pp_pi.hpp]
+
+[preprocessed_pi]
+
+
+[endsect] [/section:new_const Defining New Constants]
+
+[section:FAQ FAQs]
+
+[h4 Why are ['these] Constants Chosen?]
+It is, of course, impossible to please everyone with a list like this.
+
+Some of the criteria we have used are:
+
+* Used in Boost.Math.
+* Commonly used.
+* Expensive to compute.
+* Requested by users.
+* [@http://en.wikipedia.org/wiki/Mathematical_constant Used in science and mathematics.]
+* No integer values (because so cheap to construct).[br]
+(You can easily define your own if found convenient, for example: `FPT one =static_cast(42);`).
+
+[h4 How are constants named?]
+* Not macros, so no upper case.
+* All lower case (following C++ standard names).
+* No CamelCase.
+* Underscore as _ delimiter between words.
+* Numbers spelt as words rather than decimal digits (except following pow).
+* Abbreviation conventions:
+ * root for square root.
+ * cbrt for cube root.
+ * pow for pow function using decimal digits like pow23 for n[super 2/3].
+ * div for divided by or operator /.
+ * minus for operator -, plus for operator +.
+ * sqr for squared.
+ * cubed for cubed n[super 3].
+ * words for greek, like [pi], [zeta] and [Gamma].
+ * words like half, third, three_quarters, sixth for fractions. (Digit(s) can get muddled).
+ * log10 for log[sub 10]
+ * ln for log[sub e]
+
+[h4 How are the constants derived?]
+
+The constants have all been calculated using high-precision software working
+with up to 300-bit precision giving about 100 decimal digits.
+(The precision can be arbitrarily chosen and is limited only by compute time).
+
+[h4 How Accurate are the constants?]
+The minimum accuracy chosen (100 decimal digits) exceeds the
+accuracy of reasonably-foreseeable floating-point hardware (256-bit)
+and should meet most high-precision computations.
+
+[h4 How are the constants tested?]
+
+# Comparison using Boost.Test BOOST_CHECK_CLOSE_FRACTION using long double literals,
+with at least 35 decimal digits, enough to be accurate for all long double implementations.
+The tolerance is usually twice `long double epsilon`.
+
+# Comparison with calculation at long double precision.
+This often requires a slightly higher tolerance than two epsilon
+because of computational noise from round-off etc,
+especially when trig and other functions are called.
+
+# Comparison with independent published values,
+for example, using [@http://oeis.org/ The On-Line Encyclopedia of Integer Sequences (OEIS)]
+again using at least 35 decimal digits strings.
+
+# Comparison with independely calculated values using arbitrary precision tools like
+[@http://www.wolfram.com/mathematica/ Mathematica], again using at least 35 decimal digits literal strings.
+
+[warning We have not yet been able to [*check] that
+[*all] constants are accurate at the full arbitrary precision,
+at present 100 decimal digits.
+But certain key values like `e` and `pi` appear to be accurate
+and internal consistencies suggest that others are this accurate too.
+]
+
+[h4 Why is Portability important?]
+
+Code written using math constants is easily portable even when using different
+floating-point types with differing precision.
+
+It is a mistake to expect that results of computations will be [*identical], but
+you can achieve the [*best accuracy possible for the floating-point type in use].
+
+This has no extra cost to the user, but reduces irritating,
+and often confusing and very hard-to-trace effects,
+caused by the intrinsically limited precision of floating-point calculations.
+
+A harmless symptom of this limit is a spurious least-significant digit;
+at worst, slightly inaccurate constants sometimes cause iterating algorithms
+to diverge wildly because internal comparisons just fail.
+
+[h4 What is the Internal Format of the constants, and why?]
+
+See [link math_toolkit.constants.tutorial tutorial] above for normal use,
+but this FAQ explains the internal details used for the constants.
+
+Constants are stored as 100 decimal digit values.
+However, some compilers do not accept decimal digits strings as long as this.
+So the constant is split into two parts, with the first containing at least
+128-bit long double precision (35 decimal digits),
+and for consistency should be in scientific format with a signed exponent.
+
+The second part is the value of the constant expressed as a string literal,
+accurate to at least 100 decimal digits (in practice that means at least 102 digits).
+Again for consistency use scientific format with a signed exponent.
+
+For types with precision greater than a long double,
+then if T is constructible `T `is constructible from a `const char*`
+then it's directly constructed from the string,
+otherwise we fall back on lexical_cast to convert to type `T`.
+(Using a string is necessary because you can't use a numeric constant
+since even a `long double` might not have enough digits).
+
+So, for example, a constant like pi is internally defined as
+
+ BOOST_DEFINE_MATH_CONSTANT(pi, 3.141592653589793238462643383279502884e+00, "3.14159265358979323846264338327950288419716939937510582097494459230781640628620899862803482534211706798214808651e+00");
+
+In this case the significand is 109 decimal digits, ensuring 100 decimal digits are exact, and exponent is zero.
+
+See [link math_toolkit.constants.new_const defining new constants] to calculate new constants.
+
+A macro definition like this can be pasted into user code where convenient,
+or into `boost/math/constants.hpp` if it is to be added to the Boost.Math library.
+
+[h4 What Floating-point Types could I use?]
+
+Apart from the built-in floating-point types `float`, `double`, `long double`,
+there are several arbitrary precision floating-point classes available,
+but most are not licensed for commercial use.
+
+[h5 Boost.Multiprecision by Christopher Kormanyos]
+
+This work is based on an earlier work called e-float:
+Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations,
+in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011.
+[@http://doi.acm.org/10.1145/1916461.1916469]
+[@https://svn.boost.org/svn/boost/sandbox/e_float/ e_float]
+but is now re-factored and available under the Boost license in the Boost-sandbox at
+[@https://svn.boost.org/svn/boost/sandbox/multiprecision/ multiprecision]
+where it is being refined and prepared for review.
+
+[h5 Boost.cpp_float by John Maddock using Expression Templates]
+
+[@https://svn.boost.org/svn/boost/sandbox/big_number/ Big Number]
+which is a reworking of [@https://svn.boost.org/svn/boost/sandbox/e_float/ e_float]
+by Christopher Kormanyos to use expression templates for faster execution.
+
+[h5 NTL class quad_float]
+
+[@http://shoup.net/ntl/ NTL] by Victor Shoup has fixed and arbitrary high precision fixed and floating-point types.
+However none of these are licenced for commercial use.
+
+ #include // quad precision 106-bit, about 32 decimal digits.
+ using NTL::to_quad_float; // Less precise than arbitrary precision NTL::RR.
+
+NTL class `quad_float`, which gives a form of quadruple precision,
+106-bit significand (but without an extended exponent range.)
+With an IEC559/IEEE 754 compatible processor,
+for example Intel X86 family, with 64-bit double, and 53-bit significand,
+using the significands of [*two] 64-bit doubles,
+if `std::numeric_limits::digits10` is 16,
+then we get about twice the precision,
+so `std::numeric_limits::digits10()` should be 32.
+(the default `std::numeric_limits::digits10()` should be about 40).
+(which seems to agree with experiments).
+We output constants (including some noisy bits,
+an approximation to `std::numeric_limits::max_digits10()`)
+by adding 2 extra decimal digits, so using `quad_float::SetOutputPrecision(32 + 2);`
+
+Apple Mac/Darwin uses a similar ['doubledouble] 106-bit for its built-in `long double` type.
+
+[note The precision of all `doubledouble` floating-point types is rather odd and values given are only approximate.]
+
+[h5 NTL class RR]
+
+Arbitrary precision floating point with NTL class RR,
+default is 150 bit (about 50 decimal digits)
+used here with 300 bit to output 100 decimal digits,
+enough for many practical non-'number-theoretic' C++ applications.
+
+NTL is [*not licenced for commercial use].
+
+This class is used in Boost.Math and an option when using big_number projects to calculate new math constants.
+
+[h5 GMP and MPFR]
+
+[@gmplib.org GMP] and [@http://www.mpfr.org/ MPFR] have also been used to compute constants,
+but are licensed under the [@http://www.gnu.org/copyleft/lesser.html Lesser GPL license]
+and are [*not licensed for commercial use].
+
+[h4 What happened to a previous collection of constants proposed for Boost?]
+
+A review concluded that the way in which the constants were presented did not meet many peoples needs.
+None of the methods proposed met many users' essential requirement to allow writing simply `pi` rather than `pi()`.
+Many science and engineering equations look difficult to read when because function call brackets can be confused
+with the many other brackets often needed. All the methods then proposed of avoiding the brackets failed to meet all needs,
+often on grounds of complexity and lack of applicability to various realistic scenarios.
+
+So the simple namespace method, proposed on its own, but rejected at the first review,
+has been added to allow users to have convenient access to float, double and long double values,
+but combined with template struct and functions to allow simultaneous use
+with other non-built-in floating-point types.
+
+[h4 Why do the constants (internally) have a struct rather than a simple function?]
+
+A function mechanism was provided by in previous versions of Boost.Math.
+
+The new mechanism is to permit partial specialization. See Custom Specializing a constant above.
+It should also allow use with other packages like [@http://www.ttmath.org/ ttmath Bignum C++ library.]
+
+[h4 Where can I find other high precision constants?]
+
+# Constants with very high precision and good accuracy (>40 decimal digits)
+from Simon Plouffe's web based collection [@http://pi.lacim.uqam.ca/eng/].
+# [@https://oeis.org/ The On-Line Encyclopedia of Integer Sequences (OEIS)]
+# Checks using printed text optically scanned values and converted from:
+D. E. Knuth, Art of Computer Programming, Appendix A, Table 1, Vol 1, ISBN 0 201 89683 4 (1997)
+# M. Abrahamovitz & I. E. Stegun, National Bureau of Standards, Handbook of Mathematical Functions,
+a reference source for formulae now superceded by
+# Frank W. Olver, Daniel W. Lozier, Ronald F. Boisvert, Charles W. Clark, NIST Handbook of Mathemetical Functions, Cambridge University Press, ISBN 978-0-521-14063-8, 2010.
+# John F Hart, Computer Approximations, Kreiger (1978) ISBN 0 88275 642 7.
+# Some values from Cephes Mathematical Library, Stephen L. Moshier
+and CALC100 100 decimal digit Complex Variable Calculator Program, a DOS utility.
+# Xavier Gourdon, Pascal Sebah, 50 decimal digits constants at [@http://numbers.computation.free.fr/Constants/constants.html Number, constants and computation].
+
+[h4 Where are Physical Constants?]
+
+Not here in this Boost.Math collection, because physical constants:
+
+* Are measurements.
+* Are not truly constant and keeping changing as mensuration technology improves.
+* Have a instrinsic uncertainty.
+* Mathematical constants are stored and represented at varying precision, but should never be inaccurate.
+
+Some physical constants may be available in Boost.Units.
+
+[endsect] [/section:FAQ FAQ]
+
+[endsect]
+
+
+
diff --git a/doc/sf_and_dist/credits.qbk b/doc/sf_and_dist/credits.qbk
index 624f18f45..fa747ca36 100644
--- a/doc/sf_and_dist/credits.qbk
+++ b/doc/sf_and_dist/credits.qbk
@@ -64,6 +64,8 @@ Thanks to Mark Coleman and Georgi Boshnakov for spot test values
from __Mathematica, and of course,
to Eric Weissten for nurturing __Mathworld, an invaluable resource.
+The Skew-normal distribution and Owen's t function were written by Benjamin Sobotta.
+
Plots of the functions and distributions were prepared in
[@http://www.w3.org/ W3C] standard
[@http://www.svg.org/ Scalable Vector Graphic (SVG)] format
diff --git a/doc/sf_and_dist/dist_reference.qbk b/doc/sf_and_dist/dist_reference.qbk
index a73974c6e..3a7265a2c 100644
--- a/doc/sf_and_dist/dist_reference.qbk
+++ b/doc/sf_and_dist/dist_reference.qbk
@@ -30,6 +30,7 @@
[include distributions/pareto.qbk]
[include distributions/poisson.qbk]
[include distributions/rayleigh.qbk]
+[include distributions/skew_normal.qbk]
[include distributions/students_t.qbk]
[include distributions/triangular.qbk]
[include distributions/uniform.qbk]
diff --git a/doc/sf_and_dist/distributions/normal.qbk b/doc/sf_and_dist/distributions/normal.qbk
index cb437c214..98d96d4a7 100644
--- a/doc/sf_and_dist/distributions/normal.qbk
+++ b/doc/sf_and_dist/distributions/normal.qbk
@@ -35,13 +35,22 @@ is known as the ['Standard Normal Distribution].
Given mean [mu][space]and standard deviation [sigma][space]it has the PDF:
-[equation normal_ref1]
+[space] [equation normal_ref1]
The variation the PDF with its parameters is illustrated
in the following graph:
[graph normal_pdf]
+The cumulative distribution function is given by
+
+[space] [equation normal_cdf]
+
+and illustrated by this graph
+
+[graph normal_cdf]
+
+
[h4 Member Functions]
normal_distribution(RealType mean = 0, RealType sd = 1);
@@ -94,6 +103,7 @@ and /s/ is its standard deviation.
[[quantile from the complement][Using the relation: x = m + s * sqrt(2) * __erfc_inv(2*p)]]
[[mean and standard deviation][The same as `dist.mean()` and `dist.standard_deviation()`]]
[[mode][The same as the mean.]]
+[[median][The same as the mean.]]
[[skewness][0]]
[[kurtosis][3]]
[[kurtosis excess][0]]
@@ -102,7 +112,7 @@ and /s/ is its standard deviation.
[endsect][/section:normal_dist Normal]
[/ normal.qbk
- Copyright 2006, 2007 John Maddock and Paul A. Bristow.
+ Copyright 2006, 2007, 2012 John Maddock and Paul A. Bristow.
Distributed under the Boost Software License, Version 1.0.
(See accompanying file LICENSE_1_0.txt or copy at
http://www.boost.org/LICENSE_1_0.txt).
diff --git a/doc/sf_and_dist/distributions/skew_normal.qbk b/doc/sf_and_dist/distributions/skew_normal.qbk
new file mode 100644
index 000000000..fdd281fac
--- /dev/null
+++ b/doc/sf_and_dist/distributions/skew_normal.qbk
@@ -0,0 +1,193 @@
+[section:skew_normal_dist Skew Normal Distribution]
+
+``#include ``
+
+ namespace boost{ namespace math{
+
+ template
+ class skew_normal_distribution;
+
+ typedef skew_normal_distribution<> normal;
+
+ template
+ class skew_normal_distribution
+ {
+ public:
+ typedef RealType value_type;
+ typedef Policy policy_type;
+ // Constructor:
+ skew_normal_distribution(RealType location = 0, RealType scale = 1, RealType shape = 0);
+ // Accessors:
+ RealType location()const; // mean if normal.
+ RealType scale()const; // width, standard deviation if normal.
+ RealType shape()const; // The distribution is right skewed if shape > 0 and is left skewed if shape < 0.
+ // The distribution is normal if shape is zero.
+ };
+
+ }} // namespaces
+
+The skew normal distribution is a variant of the most well known
+Gaussian statistical distribution.
+
+The skew normal distribution with shape zero resembles the
+[@http://en.wikipedia.org/wiki/Normal_distribution Normal Distribution],
+hence the latter can be regarded as a special case of the more generic skew normal distribution.
+
+If the standard (mean = 0, scale = 1) normal distribution probability density function is
+
+[space][space][equation normal01_pdf]
+
+and the cumulative distribution function
+
+[space][space][equation normal01_cdf]
+
+then the [@http://en.wikipedia.org/wiki/Probability_density_function PDF]
+of the [@http://en.wikipedia.org/wiki/Skew_normal_distribution skew normal distribution]
+with shape parameter [alpha], defined by O'Hagan and Leonhard (1976) is
+
+[space][space][equation skew_normal_pdf0]
+
+Given [@http://en.wikipedia.org/wiki/Location_parameter location] [xi],
+[@http://en.wikipedia.org/wiki/Scale_parameter scale] [omega],
+and [@http://en.wikipedia.org/wiki/Shape_parameter shape] [alpha],
+it can be
+[@http://en.wikipedia.org/wiki/Skew_normal_distribution transformed],
+to the form:
+
+[space][space][equation skew_normal_pdf]
+
+and [@http://en.wikipedia.org/wiki/Cumulative_distribution_function CDF]:
+
+[space][space][equation skew_normal_cdf]
+
+where ['T(h,a)] is Owen's T function, and ['[Phi](x)] is the normal distribution.
+
+The variation the PDF and CDF with its parameters is illustrated
+in the following graphs:
+
+[graph skew_normal_pdf]
+[graph skew_normal_cdf]
+
+[h4 Member Functions]
+
+ skew_normal_distribution(RealType location = 0, RealType scale = 1, RealType shape = 0);
+
+Constructs a skew_normal distribution with location [xi],
+scale [omega] and shape [alpha].
+
+Requires scale > 0, otherwise __domain_error is called.
+
+ RealType location()const;
+
+returns the location [xi] of this distribution,
+
+ RealType scale()const;
+
+returns the scale [omega] of this distribution,
+
+ RealType shape()const;
+
+returns the shape [alpha] of this distribution.
+
+(Location and scale function match other similar distributions,
+allowing the functions `find_location` and `find_scale` to be used generically).
+
+[note While the shape parameter may be chosen arbitrarily (finite),
+the resulting [*skewness] of the distribution is in fact limited to about (-1, 1);
+strictly, the interval is (-0.9952717, 0.9952717).
+
+A parameter [delta] is related to the shape [alpha] by
+[delta] = [alpha] / (1 + [alpha][pow2]),
+and used in the expression for skewness
+[equation skew_normal_skewness]
+] [/note]
+
+[h4 References]
+
+* [@http://azzalini.stat.unipd.it/SN/ Skew-Normal Probability Distribution] for many links and bibliography.
+* [@http://azzalini.stat.unipd.it/SN/Intro/intro.html A very brief introduction to the skew-normal distribution]
+by Adelchi Azzalini (2005-11-2).
+* See a [@http://www.tri.org.au/azzalini.html skew-normal function animation].
+
+[h4 Non-member Accessors]
+
+All the [link math_toolkit.dist.dist_ref.nmp usual non-member accessor functions]
+that are generic to all distributions are supported: __usual_accessors.
+
+The domain of the random variable is ['-[max_value], +[min_value]].
+Infinite values are not supported.
+
+There are no [@http://en.wikipedia.org/wiki/Closed-form_expression closed-form expression]
+known for the mode and median, but these are computed for the
+
+* mode - by finding the maximum of the PDF.
+* median - by computing `quantile(1/2)`.
+
+The maximum of the PDF is sought through searching the root of f'(x)=0.
+
+Both involve iterative methods that will have lower accuracy than other estimates.
+
+[h4 Testing]
+
+__R using library(sn) described at
+[@http://azzalini.stat.unipd.it/SN/ Skew-Normal Probability Distribution],
+and at [@http://cran.r-project.org/web/packages/sn/sn.pd R skew-normal(sn) package].
+
+Package sn provides functions related to the skew-normal (SN)
+and the skew-t (ST) probability distributions,
+both for the univariate and for the the multivariate case,
+including regression models.
+
+__Mathematica was also used to generate some more accurate spot test data.
+
+[h4 Accuracy]
+
+The skew_normal distribution with shape = zero is implemented as a special case,
+equivalent to the normal distribution in terms of the
+[link math_toolkit.special.sf_erf.error_function error function],
+and therefore should have excellent accuracy.
+
+The PDF and mean, variance, skewness and kurtosis are also accurately evaluated using
+[@http://en.wikipedia.org/wiki/Analytical_expression analytical expressions].
+The CDF requires [@http://en.wikipedia.org/wiki/Owen%27s_T_function Owen's T function]
+that is evaluated using a Boost C++ __owens_t implementation of the algorithms of
+M. Patefield and D. Tandy, Journal of Statistical Software, 5(5), 1-25 (2000);
+the complicated accuracy of this function is discussed in detail at __owens_t.
+
+The median and mode are calculated by iterative root finding, and both will be less accurate.
+
+[h4 Implementation]
+
+In the following table, [xi] is the location of the distribution,
+and [omega] is its scale, and [alpha] is its shape.
+
+[table
+[[Function][Implementation Notes]]
+[[pdf][Using:[equation skew_normal_pdf] ]]
+[[cdf][Using: [equation skew_normal_cdf]\n
+where ['T(h,a)] is Owen's T function, and ['[Phi](x)] is the normal distribution. ]]
+[[cdf complement][Using: complement of normal distribution + 2 * Owens_t]]
+[[quantile][Maximum of the pdf is sought through searching the root of f'(x)=0]]
+[[quantile from the complement][-quantile(SN(-location [xi], scale [omega], -shape[alpha]), p)]]
+[[location][location [xi]]]
+[[scale][scale [omega]]]
+[[shape][shape [alpha]]]
+[[median][quantile(1/2)]]
+[[mean][[equation skew_normal_mean]]]
+[[mode][Maximum of the pdf is sought through searching the root of f'(x)=0]]
+[[variance][[equation skew_normal_variance] ]]
+[[skewness][[equation skew_normal_skewness] ]]
+[[kurtosis][kurtosis excess-3]]
+[[kurtosis excess] [ [equation skew_normal_kurt_ex] ]]
+] [/table]
+
+[endsect] [/section:skew_normal_dist skew_Normal]
+
+[/ skew_normal.qbk
+ Copyright 2012 Bejamin Sobotta, John Maddock and Paul A. Bristow.
+ Distributed under the Boost Software License, Version 1.0.
+ (See accompanying file LICENSE_1_0.txt or copy at
+ http://www.boost.org/LICENSE_1_0.txt).
+]
+
diff --git a/doc/sf_and_dist/equations/beta1.svg.svg b/doc/sf_and_dist/equations/beta1.svg.svg
new file mode 100644
index 000000000..877ea89ca
--- /dev/null
+++ b/doc/sf_and_dist/equations/beta1.svg.svg
@@ -0,0 +1,2 @@
+
+beta(a,b)=B(a,b)=Γ(a)Γ(b)Γ(a+b)
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/hankel1.mml b/doc/sf_and_dist/equations/hankel1.mml
new file mode 100644
index 000000000..d10ee157f
--- /dev/null
+++ b/doc/sf_and_dist/equations/hankel1.mml
@@ -0,0 +1,57 @@
+
+
diff --git a/doc/sf_and_dist/equations/hankel1.png b/doc/sf_and_dist/equations/hankel1.png
new file mode 100644
index 000000000..7ec65d9b6
Binary files /dev/null and b/doc/sf_and_dist/equations/hankel1.png differ
diff --git a/doc/sf_and_dist/equations/hankel1.svg b/doc/sf_and_dist/equations/hankel1.svg
new file mode 100644
index 000000000..564245a82
--- /dev/null
+++ b/doc/sf_and_dist/equations/hankel1.svg
@@ -0,0 +1,2 @@
+
+Jν(−z)=(−z)νz−νJν(z)
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/hankel2.mml b/doc/sf_and_dist/equations/hankel2.mml
new file mode 100644
index 000000000..1b46c4f6a
--- /dev/null
+++ b/doc/sf_and_dist/equations/hankel2.mml
@@ -0,0 +1,152 @@
+
+
diff --git a/doc/sf_and_dist/equations/hankel2.png b/doc/sf_and_dist/equations/hankel2.png
new file mode 100644
index 000000000..614f7efbd
Binary files /dev/null and b/doc/sf_and_dist/equations/hankel2.png differ
diff --git a/doc/sf_and_dist/equations/hankel2.svg b/doc/sf_and_dist/equations/hankel2.svg
new file mode 100644
index 000000000..ced038128
--- /dev/null
+++ b/doc/sf_and_dist/equations/hankel2.svg
@@ -0,0 +1,2 @@
+
+Yν(−z)=zνYν(z)(−z)−ν+((−z)νz−ν−(−z)−νzν)Jν(z)cot(πν);ν∉ℤ
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/hankel3.mml b/doc/sf_and_dist/equations/hankel3.mml
new file mode 100644
index 000000000..8b23d8f66
--- /dev/null
+++ b/doc/sf_and_dist/equations/hankel3.mml
@@ -0,0 +1,115 @@
+
+
diff --git a/doc/sf_and_dist/equations/hankel3.png b/doc/sf_and_dist/equations/hankel3.png
new file mode 100644
index 000000000..af2605a4b
Binary files /dev/null and b/doc/sf_and_dist/equations/hankel3.png differ
diff --git a/doc/sf_and_dist/equations/hankel3.svg b/doc/sf_and_dist/equations/hankel3.svg
new file mode 100644
index 000000000..21a9bb6ee
--- /dev/null
+++ b/doc/sf_and_dist/equations/hankel3.svg
@@ -0,0 +1,2 @@
+
+Yν(−z)=(−1)ν(Yν(z)−2π(log(z)−log(−z))Jν(z));ν∈ℤ
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/hankel4.mml b/doc/sf_and_dist/equations/hankel4.mml
new file mode 100644
index 000000000..43bb18bde
--- /dev/null
+++ b/doc/sf_and_dist/equations/hankel4.mml
@@ -0,0 +1,69 @@
+
+]>
+
+hankel4
+
+
+
+
+
+
diff --git a/doc/sf_and_dist/equations/hankel4.png b/doc/sf_and_dist/equations/hankel4.png
new file mode 100644
index 000000000..f1ba66f46
Binary files /dev/null and b/doc/sf_and_dist/equations/hankel4.png differ
diff --git a/doc/sf_and_dist/equations/hankel4.svg b/doc/sf_and_dist/equations/hankel4.svg
new file mode 100644
index 000000000..1126b948e
--- /dev/null
+++ b/doc/sf_and_dist/equations/hankel4.svg
@@ -0,0 +1,2 @@
+
+hv(1)(x)=Ï€21xHv+12(1)(x)
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/hankel5.mml b/doc/sf_and_dist/equations/hankel5.mml
new file mode 100644
index 000000000..70935cb0a
--- /dev/null
+++ b/doc/sf_and_dist/equations/hankel5.mml
@@ -0,0 +1,69 @@
+
+]>
+
+hankel5
+
+
+
+
+
+
diff --git a/doc/sf_and_dist/equations/hankel5.png b/doc/sf_and_dist/equations/hankel5.png
new file mode 100644
index 000000000..822dac13e
Binary files /dev/null and b/doc/sf_and_dist/equations/hankel5.png differ
diff --git a/doc/sf_and_dist/equations/hankel5.svg b/doc/sf_and_dist/equations/hankel5.svg
new file mode 100644
index 000000000..3dd64979a
--- /dev/null
+++ b/doc/sf_and_dist/equations/hankel5.svg
@@ -0,0 +1,2 @@
+
+hv(2)(x)=Ï€21xHv+12(2)(x)
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/normal01_cdf.mml b/doc/sf_and_dist/equations/normal01_cdf.mml
new file mode 100644
index 000000000..67951f026
--- /dev/null
+++ b/doc/sf_and_dist/equations/normal01_cdf.mml
@@ -0,0 +1,62 @@
+
+]>
+
+normal_cdf
+
+
+
+
+
+
diff --git a/doc/sf_and_dist/equations/normal01_cdf.png b/doc/sf_and_dist/equations/normal01_cdf.png
new file mode 100644
index 000000000..6dda9c7a7
Binary files /dev/null and b/doc/sf_and_dist/equations/normal01_cdf.png differ
diff --git a/doc/sf_and_dist/equations/normal01_cdf.svg b/doc/sf_and_dist/equations/normal01_cdf.svg
new file mode 100644
index 000000000..42910865a
--- /dev/null
+++ b/doc/sf_and_dist/equations/normal01_cdf.svg
@@ -0,0 +1,2 @@
+
+Φ(x)=∫−∞xφ(t)dt=12[1+erf(x2)]
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/normal01_pdf.mml b/doc/sf_and_dist/equations/normal01_pdf.mml
new file mode 100644
index 000000000..eb145176e
--- /dev/null
+++ b/doc/sf_and_dist/equations/normal01_pdf.mml
@@ -0,0 +1,45 @@
+
+]>
+
+skew_normal_pdf
+
+
+
+
+
+
+
diff --git a/doc/sf_and_dist/equations/normal01_pdf.png b/doc/sf_and_dist/equations/normal01_pdf.png
new file mode 100644
index 000000000..0b7e762c4
Binary files /dev/null and b/doc/sf_and_dist/equations/normal01_pdf.png differ
diff --git a/doc/sf_and_dist/equations/normal01_pdf.svg b/doc/sf_and_dist/equations/normal01_pdf.svg
new file mode 100644
index 000000000..1907d0cfc
--- /dev/null
+++ b/doc/sf_and_dist/equations/normal01_pdf.svg
@@ -0,0 +1,2 @@
+
+φ(x)=1(2π)e−x22
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/normal_cdf.mml b/doc/sf_and_dist/equations/normal_cdf.mml
new file mode 100644
index 000000000..506d94813
--- /dev/null
+++ b/doc/sf_and_dist/equations/normal_cdf.mml
@@ -0,0 +1,50 @@
+
+]>
+
+normal_cdf
+
+
+
+
+
+
diff --git a/doc/sf_and_dist/equations/normal_cdf.png b/doc/sf_and_dist/equations/normal_cdf.png
new file mode 100644
index 000000000..280573726
Binary files /dev/null and b/doc/sf_and_dist/equations/normal_cdf.png differ
diff --git a/doc/sf_and_dist/equations/normal_cdf.svg b/doc/sf_and_dist/equations/normal_cdf.svg
new file mode 100644
index 000000000..494f760ed
--- /dev/null
+++ b/doc/sf_and_dist/equations/normal_cdf.svg
@@ -0,0 +1,2 @@
+
+12{1+erf((x−μ)(2σ2))]
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/owens_t.mml b/doc/sf_and_dist/equations/owens_t.mml
new file mode 100644
index 000000000..7ad848c56
--- /dev/null
+++ b/doc/sf_and_dist/equations/owens_t.mml
@@ -0,0 +1,105 @@
+
+
+]>
+
+owens_t
+
+
+
+
+
+
+
+
+
+
+
+
+
+
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/owens_t.png b/doc/sf_and_dist/equations/owens_t.png
new file mode 100644
index 000000000..980988da7
Binary files /dev/null and b/doc/sf_and_dist/equations/owens_t.png differ
diff --git a/doc/sf_and_dist/equations/owens_t.svg b/doc/sf_and_dist/equations/owens_t.svg
new file mode 100644
index 000000000..c1c51cec2
--- /dev/null
+++ b/doc/sf_and_dist/equations/owens_t.svg
@@ -0,0 +1,2 @@
+
+T(h,a)=1(2π)∫0aexp{−12h2(1+x2)}(1+x2)dx;(−∞<h,a<+∞)
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/skew_normal_cdf.mml b/doc/sf_and_dist/equations/skew_normal_cdf.mml
new file mode 100644
index 000000000..c3306690b
--- /dev/null
+++ b/doc/sf_and_dist/equations/skew_normal_cdf.mml
@@ -0,0 +1,56 @@
+
+]>
+
+skew_normal_cdf
+
+
+
+
+
+
+
+
+
diff --git a/doc/sf_and_dist/equations/skew_normal_cdf.png b/doc/sf_and_dist/equations/skew_normal_cdf.png
new file mode 100644
index 000000000..0578921e8
Binary files /dev/null and b/doc/sf_and_dist/equations/skew_normal_cdf.png differ
diff --git a/doc/sf_and_dist/equations/skew_normal_cdf.svg b/doc/sf_and_dist/equations/skew_normal_cdf.svg
new file mode 100644
index 000000000..85aade534
--- /dev/null
+++ b/doc/sf_and_dist/equations/skew_normal_cdf.svg
@@ -0,0 +1,2 @@
+
+Φ((x−ξ)ω)−2T(((x−ξ)ω),α)
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/skew_normal_kurt_ex.mml b/doc/sf_and_dist/equations/skew_normal_kurt_ex.mml
new file mode 100644
index 000000000..c7871801e
--- /dev/null
+++ b/doc/sf_and_dist/equations/skew_normal_kurt_ex.mml
@@ -0,0 +1,62 @@
+
+]>
+
+skew_normal_kurt_ex
+
+
+
+
+
+
+
diff --git a/doc/sf_and_dist/equations/skew_normal_kurt_ex.png b/doc/sf_and_dist/equations/skew_normal_kurt_ex.png
new file mode 100644
index 000000000..fbe5909b5
Binary files /dev/null and b/doc/sf_and_dist/equations/skew_normal_kurt_ex.png differ
diff --git a/doc/sf_and_dist/equations/skew_normal_kurt_ex.svg b/doc/sf_and_dist/equations/skew_normal_kurt_ex.svg
new file mode 100644
index 000000000..cf5a749a1
--- /dev/null
+++ b/doc/sf_and_dist/equations/skew_normal_kurt_ex.svg
@@ -0,0 +1,2 @@
+
+2(π−3)(δ(2π))4(1−2δ2π)2
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/skew_normal_mean.mml b/doc/sf_and_dist/equations/skew_normal_mean.mml
new file mode 100644
index 000000000..c2bbb9a77
--- /dev/null
+++ b/doc/sf_and_dist/equations/skew_normal_mean.mml
@@ -0,0 +1,53 @@
+
+]>
+
+skew_normal_mean
+
+
+
+
+
+
diff --git a/doc/sf_and_dist/equations/skew_normal_mean.png b/doc/sf_and_dist/equations/skew_normal_mean.png
new file mode 100644
index 000000000..9141e7921
Binary files /dev/null and b/doc/sf_and_dist/equations/skew_normal_mean.png differ
diff --git a/doc/sf_and_dist/equations/skew_normal_mean.svg b/doc/sf_and_dist/equations/skew_normal_mean.svg
new file mode 100644
index 000000000..44d437f91
--- /dev/null
+++ b/doc/sf_and_dist/equations/skew_normal_mean.svg
@@ -0,0 +1,2 @@
+
+ξ+ωδ(2π)whereδ=α(1+α2)
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/skew_normal_pdf.mml b/doc/sf_and_dist/equations/skew_normal_pdf.mml
new file mode 100644
index 000000000..e1bab6d91
--- /dev/null
+++ b/doc/sf_and_dist/equations/skew_normal_pdf.mml
@@ -0,0 +1,89 @@
+
+]>
+
+skew_normal_pdf
+
+
+
+
+
+
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/skew_normal_pdf.png b/doc/sf_and_dist/equations/skew_normal_pdf.png
new file mode 100644
index 000000000..ae399030a
Binary files /dev/null and b/doc/sf_and_dist/equations/skew_normal_pdf.png differ
diff --git a/doc/sf_and_dist/equations/skew_normal_pdf.svg b/doc/sf_and_dist/equations/skew_normal_pdf.svg
new file mode 100644
index 000000000..01d2d5451
--- /dev/null
+++ b/doc/sf_and_dist/equations/skew_normal_pdf.svg
@@ -0,0 +1,2 @@
+
+1(ωπ)e−(x−ξ)2(2ω2)∫−∞(α((x−ξ)ω))e−t22dt
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/skew_normal_pdf0.mml b/doc/sf_and_dist/equations/skew_normal_pdf0.mml
new file mode 100644
index 000000000..c1c9dd917
--- /dev/null
+++ b/doc/sf_and_dist/equations/skew_normal_pdf0.mml
@@ -0,0 +1,36 @@
+
+]>
+
+skew_normal_pdf0
+
+
+
+
+
+
diff --git a/doc/sf_and_dist/equations/skew_normal_pdf0.png b/doc/sf_and_dist/equations/skew_normal_pdf0.png
new file mode 100644
index 000000000..3cb5d7929
Binary files /dev/null and b/doc/sf_and_dist/equations/skew_normal_pdf0.png differ
diff --git a/doc/sf_and_dist/equations/skew_normal_pdf0.svg b/doc/sf_and_dist/equations/skew_normal_pdf0.svg
new file mode 100644
index 000000000..64141ea40
--- /dev/null
+++ b/doc/sf_and_dist/equations/skew_normal_pdf0.svg
@@ -0,0 +1,2 @@
+
+f(x)=2φ(x)Φ(αx)
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/skew_normal_skewness.mml b/doc/sf_and_dist/equations/skew_normal_skewness.mml
new file mode 100644
index 000000000..accdb9321
--- /dev/null
+++ b/doc/sf_and_dist/equations/skew_normal_skewness.mml
@@ -0,0 +1,71 @@
+
+]>
+
+skew_normal_skewness
+
+
+
+
+
+
+
diff --git a/doc/sf_and_dist/equations/skew_normal_skewness.png b/doc/sf_and_dist/equations/skew_normal_skewness.png
new file mode 100644
index 000000000..f942af364
Binary files /dev/null and b/doc/sf_and_dist/equations/skew_normal_skewness.png differ
diff --git a/doc/sf_and_dist/equations/skew_normal_skewness.svg b/doc/sf_and_dist/equations/skew_normal_skewness.svg
new file mode 100644
index 000000000..9d533486a
--- /dev/null
+++ b/doc/sf_and_dist/equations/skew_normal_skewness.svg
@@ -0,0 +1,2 @@
+
+(4−π)2(δ(2π))3(1−2δ2π)(32)
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/skew_normal_variance.mml b/doc/sf_and_dist/equations/skew_normal_variance.mml
new file mode 100644
index 000000000..a9eedebbf
--- /dev/null
+++ b/doc/sf_and_dist/equations/skew_normal_variance.mml
@@ -0,0 +1,33 @@
+
+]>
+
+skew_normal_variance
+
+
+
+
+
+
\ No newline at end of file
diff --git a/doc/sf_and_dist/equations/skew_normal_variance.png b/doc/sf_and_dist/equations/skew_normal_variance.png
new file mode 100644
index 000000000..4531751b7
Binary files /dev/null and b/doc/sf_and_dist/equations/skew_normal_variance.png differ
diff --git a/doc/sf_and_dist/equations/skew_normal_variance.svg b/doc/sf_and_dist/equations/skew_normal_variance.svg
new file mode 100644
index 000000000..a09a91c5a
--- /dev/null
+++ b/doc/sf_and_dist/equations/skew_normal_variance.svg
@@ -0,0 +1,2 @@
+
+ω2(1−2δ2π)
\ No newline at end of file
diff --git a/doc/sf_and_dist/fp_facets.qbk b/doc/sf_and_dist/fp_facets.qbk
index 0f285712c..6d31ebbd5 100644
--- a/doc/sf_and_dist/fp_facets.qbk
+++ b/doc/sf_and_dist/fp_facets.qbk
@@ -338,7 +338,7 @@ and also string representations used by other programming languages.
[h5 signed_zero]
If the `signed_zero` flag is used with `nonfinite_num_put`,
-then the facet will distinguish between positive and negative zero.
+then the facet will always distinguish between positive and negative zero.
It will format positive zero as "0" or "+0" and negative zero as "-0".
The string representation of positive zero can be controlled
with the `showpos` and `noshowpos` manipulators.
@@ -348,6 +348,16 @@ The input facet `nonfinite_num_get` always parses "0" and "+0"
as positive zero and "-0" as negative zero,
as do most implementations of `std::num_get`.
+[note If the `signed_zero` flag is not set (the default), then a negative zero value
+will be displayed on output in whatever way the platform normally handles it.
+For most platforms, this it will format positive zero as "0" or "+0" and negative zero as "-0".
+But setting the `signed_zero` flag may be more portable.]
+
+[tip A negative zero value can be portably produced using the changesign function
+`(changesign)(static_cast(0))` where `ValType` is `float`, `double` or `long double`,
+ or a User-Defined floating-point type (UDT) provided that this UDT has a sign
+and that the changesign function is implemented.]
+
[h5 trap_infinity]
If the `trap_infinity` flag is used with `nonfinite_num_put`,
@@ -361,6 +371,7 @@ If the `trap_infinity` flag is used with `nonfinite_num_get`,
then the facet will set the `fail bit` of the stream when an attempt is made
to parse a string that represents positive or negative infinity.
+
(See Design Rationale below for a discussion of this inconsistency.)
[h5 trap_nan]
@@ -502,6 +513,9 @@ Example [@../../../example/nonfinite_loopback_ok.cpp] shows loopback works OK.
Example [@../../../example/nonfinite_num_facet.cpp] shows output and re-input
of various finite and nonfinite values.
+A simple example of trapping nonfinite output is at
+[@../../../example/nonfinite_num_facet_trap.cpp nonfinite_num_facet_trap.cpp].
+
A very basic example of using Boost.Archive is at
[@../../../example/nonfinite_serialization_archives.cpp].
diff --git a/doc/sf_and_dist/graphs/dist_graphs.cpp b/doc/sf_and_dist/graphs/dist_graphs.cpp
index 4fc4d0946..535108d65 100644
--- a/doc/sf_and_dist/graphs/dist_graphs.cpp
+++ b/doc/sf_and_dist/graphs/dist_graphs.cpp
@@ -10,7 +10,7 @@
\author John Maddock and Paul A. Bristow
*/
// Copyright John Maddock 2008.
-// Copyright Paul A. Bristow 2008, 2009
+// Copyright Paul A. Bristow 2008, 2009, 2012
// Use, modification and distribution are subject to the
// Boost Software License, Version 1.0. (See accompanying file
// LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
@@ -499,6 +499,26 @@ int main()
rayleigh_cdf_plotter.add(boost::math::rayleigh_distribution<>(10), "σ=10");
rayleigh_cdf_plotter.plot("Rayleigh Distribution CDF", "rayleigh_cdf.svg");
+ distribution_plotter >
+ skew_normal_plotter;
+ skew_normal_plotter.add(boost::math::skew_normal_distribution<>(0,1,0), "{0,1,0}");
+ skew_normal_plotter.add(boost::math::skew_normal_distribution<>(0,1,1), "{0,1,1}");
+ skew_normal_plotter.add(boost::math::skew_normal_distribution<>(0,1,4), "{0,1,4}");
+ skew_normal_plotter.add(boost::math::skew_normal_distribution<>(0,1,20), "{0,1,20}");
+ skew_normal_plotter.add(boost::math::skew_normal_distribution<>(0,1,-2), "{0,1,-2}");
+ skew_normal_plotter.add(boost::math::skew_normal_distribution<>(-2,0.5,-1), "{-2,0.5,-1}");
+ skew_normal_plotter.plot("Skew Normal Distribution PDF", "skew_normal_pdf.svg");
+
+ distribution_plotter >
+ skew_normal_cdf_plotter(false);
+ skew_normal_cdf_plotter.add(boost::math::skew_normal_distribution<>(0,1,0), "{0,1,0}");
+ skew_normal_cdf_plotter.add(boost::math::skew_normal_distribution<>(0,1,1), "{0,1,1}");
+ skew_normal_cdf_plotter.add(boost::math::skew_normal_distribution<>(0,1,4), "{0,1,4}");
+ skew_normal_cdf_plotter.add(boost::math::skew_normal_distribution<>(0,1,20), "{0,1,20}");
+ skew_normal_cdf_plotter.add(boost::math::skew_normal_distribution<>(0,1,-2), "{0,1,-2}");
+ skew_normal_cdf_plotter.add(boost::math::skew_normal_distribution<>(-2,0.5,-1), "{-2,0.5,-1}");
+ skew_normal_cdf_plotter.plot("Skew Normal Distribution CDF", "skew_normal_cdf.svg");
+
distribution_plotter >
triangular_plotter;
triangular_plotter.add(boost::math::triangular_distribution<>(-1,0,1), "{-1,0,1}");
@@ -622,4 +642,4 @@ int main()
hypergeometric_plotter2.add(boost::math::hypergeometric_distribution<>(450, 50, 500), "N=500, r=50, n=450");
hypergeometric_plotter2.plot("Hypergeometric Distribution PDF", "hypergeometric_pdf_2.svg");
-}
+} // int main()
diff --git a/doc/sf_and_dist/graphs/owens_integration_area.png b/doc/sf_and_dist/graphs/owens_integration_area.png
new file mode 100644
index 000000000..8ad21cb4a
Binary files /dev/null and b/doc/sf_and_dist/graphs/owens_integration_area.png differ
diff --git a/doc/sf_and_dist/graphs/owens_integration_area.svg b/doc/sf_and_dist/graphs/owens_integration_area.svg
new file mode 100644
index 000000000..d0aef8a89
--- /dev/null
+++ b/doc/sf_and_dist/graphs/owens_integration_area.svg
@@ -0,0 +1,121 @@
+
+
+
+
diff --git a/doc/sf_and_dist/graphs/plot_owens_3d_xyp.png b/doc/sf_and_dist/graphs/plot_owens_3d_xyp.png
new file mode 100644
index 000000000..73d7a03b3
Binary files /dev/null and b/doc/sf_and_dist/graphs/plot_owens_3d_xyp.png differ
diff --git a/doc/sf_and_dist/graphs/plot_owens_t.png b/doc/sf_and_dist/graphs/plot_owens_t.png
new file mode 100644
index 000000000..678a22b4e
Binary files /dev/null and b/doc/sf_and_dist/graphs/plot_owens_t.png differ
diff --git a/doc/sf_and_dist/graphs/skew_normal_cdf.png b/doc/sf_and_dist/graphs/skew_normal_cdf.png
new file mode 100644
index 000000000..73d28d8cb
Binary files /dev/null and b/doc/sf_and_dist/graphs/skew_normal_cdf.png differ
diff --git a/doc/sf_and_dist/graphs/skew_normal_cdf.svg b/doc/sf_and_dist/graphs/skew_normal_cdf.svg
new file mode 100644
index 000000000..ceb9bc869
--- /dev/null
+++ b/doc/sf_and_dist/graphs/skew_normal_cdf.svg
@@ -0,0 +1,81 @@
+
+
diff --git a/doc/sf_and_dist/graphs/skew_normal_pdf.png b/doc/sf_and_dist/graphs/skew_normal_pdf.png
new file mode 100644
index 000000000..09256c776
Binary files /dev/null and b/doc/sf_and_dist/graphs/skew_normal_pdf.png differ
diff --git a/doc/sf_and_dist/graphs/skew_normal_pdf.svg b/doc/sf_and_dist/graphs/skew_normal_pdf.svg
new file mode 100644
index 000000000..71f0056a2
--- /dev/null
+++ b/doc/sf_and_dist/graphs/skew_normal_pdf.svg
@@ -0,0 +1,80 @@
+
+
diff --git a/doc/sf_and_dist/hankel.qbk b/doc/sf_and_dist/hankel.qbk
new file mode 100644
index 000000000..014a615b3
--- /dev/null
+++ b/doc/sf_and_dist/hankel.qbk
@@ -0,0 +1,136 @@
+
+[section:hankel Hankel Functions]
+[section:cyl_hankel Cyclic Hankel Functions]
+
+[h4 Synopsis]
+
+ template
+ std::complex<``__sf_result``> cyl_hankel_1(T1 v, T2 x);
+
+ template
+ std::complex<``__sf_result``> cyl_hankel_1(T1 v, T2 x, const ``__Policy``&);
+
+ template
+ std::complex<``__sf_result``> cyl_hankel_2(T1 v, T2 x);
+
+ template
+ std::complex<``__sf_result``> cyl_hankel_2(T1 v, T2 x, const ``__Policy``&);
+
+
+[h4 Description]
+
+The functions __cyl_hankel_1 and __cyl_hankel_2 return the result of the
+[@http://dlmf.nist.gov/10.2#P3 Hankel functions] of the first and second kind respectively:
+
+[:['cyl_hankel_1(v, x) = H[sub v][super (1)](x) = J[sub v](x) + i Y[sub v](x)]]
+
+[:['cyl_hankel_2(v, x) = H[sub v][super (2)](x) = J[sub v](x) - i Y[sub v](x)]]
+
+where:
+
+['J[sub v](x)] is the Bessel function of the first kind, and ['Y[sub v](x)] is the Bessel function of the second kind.
+
+The return type of these functions is computed using the __arg_pomotion_rules
+when T1 and T2 are different types. The functions are also optimised for the
+relatively common case that T1 is an integer.
+
+[optional_policy]
+
+Note that while the arguments to these functions are real values, the results are complex.
+That means that the functions can only be instantiated on types `float`, `double` and `long double`.
+The functions have also been extended to operate over the whole range of ['v] and ['x]
+(unlike __cyl_bessel_j and __cyl_neumann).
+
+[h4 Performance]
+
+These functions are generally more efficient than two separate calls to the underlying Bessel
+functions as internally Bessel J and Y can be computed simultaneously.
+
+[h4 Testing]
+
+There are just a few spot tests to exercise all the special case handling - the bulk of the testing is done
+on the Bessel functions upon which these are based.
+
+[h4 Accuracy]
+
+Refer to __cyl_bessel_j and __cyl_neumann.
+
+[h4 Implementation]
+
+For ['x < 0] the following reflection formulae are used:
+
+[@http://functions.wolfram.com/Bessel-TypeFunctions/BesselJ/16/01/01/ [equation hankel1]]
+
+[@http://functions.wolfram.com/Bessel-TypeFunctions/BesselY/16/01/01/ [equation hankel2]]
+
+[@http://functions.wolfram.com/Bessel-TypeFunctions/BesselY/16/01/01/ [equation hankel3]]
+
+Otherwise the implementation is trivially in terms of the Bessel J and Y functions.
+
+Note however, that the Hankel functions compute the Bessel J and Y functions simultaneously,
+and therefore a single Hankel function call is more efficient than two Bessel function calls.
+The one exception is when ['v] is a small positive integer, in which case the usual Bessel function
+routines for integer order are used.
+
+[endsect]
+
+
+[section:sph_hankel Spherical Hankel Functions]
+
+[h4 Synopsis]
+
+ template
+ std::complex<``__sf_result``> sph_hankel_1(T1 v, T2 x);
+
+ template
+ std::complex<``__sf_result``> sph_hankel_1(T1 v, T2 x, const ``__Policy``&);
+
+ template
+ std::complex<``__sf_result``> sph_hankel_2(T1 v, T2 x);
+
+ template
+ std::complex<``__sf_result``> sph_hankel_2(T1 v, T2 x, const ``__Policy``&);
+
+
+[h4 Description]
+
+The functions __sph_hankel_1 and __sph_hankel_2 return the result of the
+[@http://dlmf.nist.gov/10.47#P1 spherical Hankel functions] of the first and second kind respectively:
+
+[equation hankel4]
+
+[equation hankel5]
+
+The return type of these functions is computed using the __arg_pomotion_rules
+when T1 and T2 are different types. The functions are also optimised for the
+relatively common case that T1 is an integer.
+
+[optional_policy]
+
+Note that while the arguments to these functions are real values, the results are complex.
+That means that the functions can only be instantiated on types `float`, `double` and `long double`.
+The functions have also been extended to operate over the whole range of ['v] and ['x]
+(unlike __cyl_bessel_j and __cyl_neumann).
+
+[h4 Testing]
+
+There are just a few spot tests to exercise all the special case handling - the bulk of the testing is done
+on the Bessel functions upon which these are based.
+
+[h4 Accuracy]
+
+Refer to __cyl_bessel_j and __cyl_neumann.
+
+[h4 Implementation]
+
+These functions are trivially implemented in terms of __cyl_hankel_1 and __cyl_hankel_2.
+
+[endsect]
+[endsect]
+
+[/
+ Copyright 2012 John Maddock.
+ Distributed under the Boost Software License, Version 1.0.
+ (See accompanying file LICENSE_1_0.txt or copy at
+ http://www.boost.org/LICENSE_1_0.txt).
+]
diff --git a/doc/sf_and_dist/html/index.html b/doc/sf_and_dist/html/index.html
index 0ea313a77..a4af5cdd0 100644
--- a/doc/sf_and_dist/html/index.html
+++ b/doc/sf_and_dist/html/index.html
@@ -48,9 +48,12 @@