Mercurial > hg > octave-lojdl
changeset 14956:466cb8673653
maint: periodic merge of default to jit
author | Max Brister <max@2bass.com> |
---|---|
date | Mon, 18 Jun 2012 16:45:10 -0500 |
parents | 609dcc297db5 (current diff) 1230d5d58d2d (diff) |
children | 9d60ee38c6e2 |
files | build-aux/mkinstalldirs configure.ac src/ov-base.h |
diffstat | 31 files changed, 300 insertions(+), 83 deletions(-) [+] |
line wrap: on
line diff
--- a/.hgsubstate +++ b/.hgsubstate @@ -1,1 +1,1 @@ -f9813bce2c06a6130a68db4478d1b16ddadaf276 gnulib +33f823397dbb0edb57503f2f6dad2362456bc6a9 gnulib
--- a/NEWS +++ b/NEWS @@ -74,11 +74,10 @@ ** Other new functions added in 3.8.0: - colorcube splinefit - lines tetramesh - rgbplot shrinkfaces - findfigs - + colorcube lines splinefit + erfcinv rgbplot tetramesh + findfigs shrinkfaces + ** Deprecated functions. The following functions were deprecated in Octave 3.4 and have been
--- a/build-aux/mkinstalldirs +++ b/build-aux/mkinstalldirs @@ -81,9 +81,9 @@ echo "mkdir -p -- $*" exec mkdir -p -- "$@" else - # On NextStep and OpenStep, the 'mkdir' command does not + # On NextStep and OpenStep, the `mkdir' command does not # recognize any option. It will interpret all options as - # directories to create, and then abort because '.' already + # directories to create, and then abort because `.' already # exists. test -d ./-p && rmdir ./-p test -d ./--version && rmdir ./--version
--- a/configure.ac +++ b/configure.ac @@ -654,6 +654,8 @@ OCTAVE_CXX_COMPLEX_SETTERS OCTAVE_CXX_COMPLEX_REFERENCE_ACCESSORS +OCTAVE_CARBON_CGDISPLAYBITSPERPIXEL + ### Check for the QHull library OCTAVE_CHECK_LIBRARY(qhull, QHull,
--- a/doc/interpreter/arith.txi +++ b/doc/interpreter/arith.txi @@ -271,6 +271,8 @@ @DOCSTRING(erfinv) +@DOCSTRING(erfcinv) + @DOCSTRING(gamma) @DOCSTRING(gammainc)
--- a/doc/interpreter/contributors.in +++ b/doc/interpreter/contributors.in @@ -102,6 +102,7 @@ Kurt Hornik Christopher Hulbert Cyril Humbert +John Hunt Teemu Ikonen Alan W. Irwin Geoff Jacobsen @@ -157,6 +158,7 @@ Rob Mahurin Ricardo Marranita Orestes Mas +Axel Mathéi Makoto Matsumoto Tatsuro Matsuoka Laurent Mazet @@ -217,6 +219,7 @@ Petter Risholm Matthew W. Roberts Andrew Ross +Fabio Rossi Mark van Rossum Joe Rothweiler Kevin Ruland
--- a/libcruft/misc/lo-error.c +++ b/libcruft/misc/lo-error.c @@ -110,8 +110,8 @@ void liboctave_fatal_with_id (const char *id, const char *fmt, ...) { + va_list args; (void) id; /*unused*/ - va_list args; va_start (args, fmt); verror ("fatal", fmt, args); va_end (args); @@ -131,8 +131,8 @@ void liboctave_warning_with_id (const char *id, const char *fmt, ...) { + va_list args; (void) id; /*unused*/ - va_list args; va_start (args, fmt); verror ("warning", fmt, args); va_end (args);
--- a/libcruft/slatec-fn/derf.f +++ b/libcruft/slatec-fn/derf.f @@ -69,7 +69,7 @@ C C ERF(X) = 1.0 - ERFC(X) FOR -1.0 .LE. X .LE. 1.0 C - IF (Y.LE.SQEPS) DERF = 2.0D0*X*X/SQRTPI + IF (Y.LE.SQEPS) DERF = 2.0D0*X/SQRTPI IF (Y.GT.SQEPS) DERF = X*(1.0D0 + DCSEVL (2.D0*X*X-1.D0, 1 ERFCS, NTERF)) RETURN
--- a/liboctave/lo-specfun.cc +++ b/liboctave/lo-specfun.cc @@ -3221,6 +3221,76 @@ return do_erfinv (x, false); } +static double do_erfcinv (double x, bool refine) +{ + // Coefficients of rational approximation. + static const double a[] = + { -2.806989788730439e+01, 1.562324844726888e+02, + -1.951109208597547e+02, 9.783370457507161e+01, + -2.168328665628878e+01, 1.772453852905383e+00 }; + static const double b[] = + { -5.447609879822406e+01, 1.615858368580409e+02, + -1.556989798598866e+02, 6.680131188771972e+01, + -1.328068155288572e+01 }; + static const double c[] = + { -5.504751339936943e-03, -2.279687217114118e-01, + -1.697592457770869e+00, -1.802933168781950e+00, + 3.093354679843505e+00, 2.077595676404383e+00 }; + static const double d[] = + { 7.784695709041462e-03, 3.224671290700398e-01, + 2.445134137142996e+00, 3.754408661907416e+00 }; + + static const double spi2 = 8.862269254527579e-01; // sqrt(pi)/2. + static const double pi = 3.14159265358979323846; + static const double pbreak = 0.95150; + double y; + + // Select case. + if (x <= 1+pbreak && x >= 1-pbreak) + { + // Middle region. + const double q = 0.5*(1-x), r = q*q; + const double yn = (((((a[0]*r + a[1])*r + a[2])*r + a[3])*r + a[4])*r + a[5])*q; + const double yd = ((((b[0]*r + b[1])*r + b[2])*r + b[3])*r + b[4])*r + 1.0; + y = yn / yd; + } + else if (x < 2.0 && x > 0.0) + { + // Tail region. + const double q = x < 1 ? sqrt (-2*log (0.5*x)) : sqrt (-2*log (0.5*(2-x))); + const double yn = ((((c[0]*q + c[1])*q + c[2])*q + c[3])*q + c[4])*q + c[5]; + const double yd = (((d[0]*q + d[1])*q + d[2])*q + d[3])*q + 1.0; + y = yn / yd; + if (x < 1-pbreak) + y *= -1; + } + else if (x == 0.0) + return octave_Inf; + else if (x == 2.0) + return -octave_Inf; + else + return octave_NaN; + + if (refine) + { + // One iteration of Halley's method gives full precision. + double u = (erf (y) - (1-x)) * spi2 * exp (y*y); + y -= u / (1 + y*u); + } + + return y; +} + +double erfcinv (double x) +{ + return do_erfcinv (x, true); +} + +float erfcinv (float x) +{ + return do_erfcinv (x, false); +} + // Implementation based on the Fortran code by W.J.Cody // see http://www.netlib.org/specfun/erf. // Templatized and simplified workflow.
--- a/liboctave/lo-specfun.h +++ b/liboctave/lo-specfun.h @@ -593,6 +593,9 @@ extern OCTAVE_API double erfinv (double x); extern OCTAVE_API float erfinv (float x); +extern OCTAVE_API double erfcinv (double x); +extern OCTAVE_API float erfcinv (float x); + extern OCTAVE_API double erfcx (double x); extern OCTAVE_API float erfcx (float x);
--- a/m4/acinclude.m4 +++ b/m4/acinclude.m4 @@ -158,6 +158,25 @@ AC_LANG_POP(C++) ]) dnl +dnl See if the Carbon Framework defines CGDisplayBitsPerPixel. +dnl +AC_DEFUN([OCTAVE_CARBON_CGDISPLAYBITSPERPIXEL], +[AC_CACHE_CHECK([whether CGDisplayBitsPerPixel is defined in the Carbon Framework], +octave_cv_carbon_cgdisplaybitsperpixel, +[AC_LANG_PUSH(C++) +AC_COMPILE_IFELSE([AC_LANG_PROGRAM([[ +#include <Carbon/Carbon.h> +]], [[ +CGDirectDisplayID display = CGMainDisplayID (); +size_t depth = CGDisplayBitsPerPixel (display); +]])], +octave_cv_carbon_cgdisplaybitsperpixel=yes, octave_cv_carbon_cgdisplaybitsperpixel=no)]) +if test $octave_cv_carbon_cgdisplaybitsperpixel = yes; then +AC_DEFINE(HAVE_CARBON_CGDISPLAYBITSPERPIXEL,1,[Define if Carbon Framework has CGDisplayBitsPerPixel]) +fi +AC_LANG_POP(C++) +]) +dnl dnl The following test is from Karl Berry's Kpathseach library. I'm dnl including it here in case we someday want to make the use of dnl kpathsea optional.
--- a/scripts/help/get_first_help_sentence.m +++ b/scripts/help/get_first_help_sentence.m @@ -50,7 +50,7 @@ error ("get_first_help_sentence: NAME must be a string"); endif - if (!isnumeric (max_len) || max_len <= 0 || max_len != fix (max_len)) + if (! isnumeric (max_len) || max_len <= 0 || max_len != fix (max_len)) error ("get_first_help_sentence: MAX_LEN must be positive integer"); endif @@ -80,10 +80,13 @@ ## This function extracts the first sentence from a plain text help text function [text, status] = first_sentence_plain_text (help_text, max_len) - ## Extract first line by searching for a period or a double line-end. - period_idx = find (help_text == '.', 1); - line_end_idx = strfind (help_text, "\n\n"); - text = help_text (1:min ([period_idx(:); line_end_idx(:); max_len; length(help_text)])); + ## Extract first line by searching for a period followed by a space class + ## character (to support periods in numbers or words) ... + period_idx = regexp (help_text, '\.\s', "once"); + ## ... or a double end-of-line (we subtract 1 because we are not interested + ## in capturing the first newline). + line_end_idx = regexp (help_text, "\n\n", "once") - 1; + text = help_text (1:min ([period_idx; line_end_idx; max_len; length(help_text)])); status = 0; endfunction @@ -91,27 +94,27 @@ ## The function works by removing @def* from the texinfo text. After this, we ## render the text to plain text using makeinfo, and then extract the first line. function [text, status] = first_sentence_texinfo (help_text, max_len) - ## Lines ending with "@\n" are continuation lines, so they should be concatenated - ## with the following line. + ## Lines ending with "@\n" are continuation lines, so they should be + ## concatenated with the following line. help_text = strrep (help_text, "@\n", " "); ## Find, and remove, lines that start with @def. This should remove things ## such as @deftypefn, @deftypefnx, @defvar, etc. keep = true (size (help_text)); def_idx = strfind (help_text, "@def"); - if (!isempty (def_idx)) + if (! isempty (def_idx)) endl_idx = find (help_text == "\n"); for k = 1:length (def_idx) - endl = endl_idx (find (endl_idx > def_idx (k), 1)); + endl = endl_idx (find (endl_idx > def_idx(k), 1)); if (isempty (endl)) - keep (def_idx (k):end) = false; + keep (def_idx(k) : end) = false; else - keep (def_idx (k):endl) = false; + keep (def_idx(k) : endl) = false; endif endfor ## Remove the @end ... that corresponds to the @def we removed above - def1 = def_idx (1); + def1 = def_idx(1); space_idx = find (help_text == " "); space_idx = space_idx (find (space_idx > def1, 1)); bracket_idx = find (help_text == "{" | help_text == "}"); @@ -120,20 +123,20 @@ error ("get_first_help_sentence: couldn't parse texinfo"); endif sep_idx = min (space_idx, bracket_idx); - def_type = help_text (def1+1:sep_idx-1); + def_type = help_text(def1+1:sep_idx-1); end_idx = strfind (help_text, sprintf ("@end %s", def_type)); if (isempty (end_idx)) error ("get_first_help_sentence: couldn't parse texinfo"); endif - endl = endl_idx (find (endl_idx > end_idx, 1)); + endl = endl_idx(find (endl_idx > end_idx, 1)); if (isempty (endl)) - keep (end_idx:end) = false; + keep(end_idx:end) = false; else - keep (end_idx:endl) = false; + keep(end_idx:endl) = false; endif - help_text = help_text (keep); + help_text = help_text(keep); endif ## Run makeinfo to generate plain text @@ -159,8 +162,8 @@ %% Test input validation %!error get_first_help_sentence () %!error get_first_help_sentence (1, 2, 3) -%!error get_first_help_sentence (1) -%!error get_first_help_sentence ('ls', 'a') -%!error get_first_help_sentence ('ls', 0) -%!error get_first_help_sentence ('ls', 80.1) +%!error <NAME must be a string> get_first_help_sentence (1) +%!error <MAX_LEN must be positive integer> get_first_help_sentence ("ls", "a") +%!error <MAX_LEN must be positive integer> get_first_help_sentence ("ls", 0) +%!error <MAX_LEN must be positive integer> get_first_help_sentence ("ls", 80.1)
--- a/scripts/help/unimplemented.m +++ b/scripts/help/unimplemented.m @@ -169,7 +169,6 @@ "echodemo", "ellipj", "ellipke", - "erfcinv", "errordlg", "evalc", "exifread",
--- a/scripts/plot/colorbar.m +++ b/scripts/plot/colorbar.m @@ -162,7 +162,7 @@ "yliminclude", "off", "zliminclude", "off", "deletefcn", {@deletecolorbar, cax, obj}); - set (cax, "deletefcn", {@resetaxis, obj}); + set (cax, "deletefcn", {@resetaxis, ax, obj}); addlistener (ax, "clim", {@update_colorbar_clim, hi, vertical}); addlistener (ax, "plotboxaspectratio", {@update_colorbar_axis, cax, obj}); @@ -195,13 +195,15 @@ endif endfunction -function resetaxis (h, d, orig_props) - if (ishandle (h) && strcmp (get (h, "type"), "axes") - && (isempty (gcbf()) || strcmp (get (gcbf(), "beingdeleted"),"off")) - && ishandle (get (h, "axes"))) - set (get (h, "axes"), "position", orig_props.position, ... - "outerposition", orig_props.outerposition, ... - "activepositionproperty", orig_props.activepositionproperty); +function resetaxis (cax, d, ax, orig_props) + if (ishandle (ax) && strcmp (get (ax, "type"), "axes")) + dellistener (ax, "position") + units = get (ax, "units"); + set (ax, "units", orig_props.units) + set (ax, "position", orig_props.position, ... + "outerposition", orig_props.outerposition, ... + "activepositionproperty", orig_props.activepositionproperty); + set (ax, "units", units) endif endfunction
--- a/scripts/plot/print.m +++ b/scripts/plot/print.m @@ -74,8 +74,18 @@ ## orientation specified. This options is equivalent to changing ## the figure's "paperorientation" property. ## +## @item -TextAlphaBits=@var{n} +## @itemx -GraphicsAlphaBits=@var{n} +## Octave is able to produce output for various printers, bitmaps, and +## vector formats by using ghostscript. +## For bitmap and printer output antialiasing is applied using +## ghostscript's TextAlphaBits and GraphicsAlphaBits options. +## The default number of bits for each is 4. +## Allowed values, for @var{N}, are 1, 2, or 4. +## ## @item -d@var{device} -## Output device, where @var{device} is one of: +## The available output format is specified by the option @var{device}, +## and is one of: ## ## @table @code ## @item ps
--- a/scripts/plot/private/__ghostscript__.m +++ b/scripts/plot/private/__ghostscript__.m @@ -32,6 +32,8 @@ opts.device = ""; opts.epscrop = false; opts.antialiasing = false; + opts.antialiasing_textalphabits = 4;, + opts.antialiasing_graphicsalphabits = 4; opts.resolution = 150; opts.papersize = ""; opts.pageoffset = [0 0]; @@ -70,7 +72,9 @@ if (opts.antialiasing && isempty (strfind (opts.device, "write"))) ## Apply anti-aliasing to all bitmap formats/devices - gs_opts = sprintf ("%s -dTextAlphaBits=4 -dGraphicsAlphaBits=4", gs_opts); + gs_opts = sprintf ("%s -dTextAlphaBits=%d -dGraphicsAlphaBits=%d", + gs_opts, opts.antialiasing_textalphabits, + opts.antialiasing_graphicsalphabits); gs_opts = sprintf ("%s -r%dx%d", gs_opts, [1, 1] * opts.resolution); elseif (any (strcmp (opts.device, {"pswrite", "ps2write", "pdfwrite"}))) gs_opts = sprintf ("%s -dEmbedAllFonts=true", gs_opts); @@ -135,7 +139,7 @@ elseif (isempty (opts.output)) cmd = sprintf ("%s %s", opts.binary, gs_opts); else - cmd = sprintf ("%s %s -sOutputFile=%s", opts.binary, gs_opts, opts.output); + cmd = sprintf ("%s %s -sOutputFile=\"%s\"", opts.binary, gs_opts, opts.output); endif if (! isempty (opts.prepend) && any (strcmpi (opts.device, {"pswrite", "ps2write", "pdfwrite"})))
--- a/scripts/plot/private/__print_parse_opts__.m +++ b/scripts/plot/private/__print_parse_opts__.m @@ -49,6 +49,8 @@ arg_st.ghostscript.pageoffset = []; arg_st.ghostscript.resolution = 150; arg_st.ghostscript.antialiasing = false; + arg_st.ghostscript.antialiasing_textalphabits = 4; + arg_st.ghostscript.antialiasing_graphicsalphabits = 4; arg_st.loose = false; arg_st.lpr_binary = __quote_path__ (__find_binary__ ("lpr")); arg_st.name = ""; @@ -118,6 +120,20 @@ arg_st.fig2dev_binary = arg{10:end}; elseif (strncmp (arg, "-PSTOEDIT:", 9)) arg_st.pstoedit_binary = arg{10:end}; + elseif (strncmpi (arg, "-textalphabits=", 15)) + n = find (arg == "="); + if (! isempty (n) && n == numel (arg) - 1 && ismember (arg(end), "124")) + arg_st.ghostscript.antialiasing_textalphabits = str2num (arg(end)); + else + error ("print: improper syntax, or value, for TextAlphaBits") + endif + elseif (strncmpi (arg, "-graphicsalphabits=", 19)) + n = find (arg == "="); + if (! isempty (n) && n == numel (arg) - 1 && ismember (arg(end), "124")) + arg_st.ghostscript.antialiasing_graphicsalphabits = str2num (arg(end)); + else + error ("print: improper syntax, or value, for GraphicsAlphaBits") + endif elseif ((length (arg) > 2) && arg(1:2) == "-G") arg_st.ghostscript.binary = file_in_path (getenv ("PATH"), arg(3:end)); if (isempty (arg_st.ghostscript.binary))
--- a/scripts/sparse/gmres.m +++ b/scripts/sparse/gmres.m @@ -32,7 +32,7 @@ ## @code{min (10, numel (b) / restart)} is used. ## ## @item @var{x0} is the initial guess, -## if not given or set to [] the default value @code{zeros(size (b))} is used. +## if not given or set to [] the default value @code{zeros (size (b))} is used. ## ## @item @var{m} is the restart parameter, ## if not given or set to [] the default value @code{numel (b)} is used. @@ -57,7 +57,7 @@ ## ## @item 2 : unused, but skipped for compatibility ## -## @item 3 : algorithm reached stagnation +## @item 3 : algorithm reached stagnation (no change between iterations) ## @end table ## ## @item @var{relres} is the final value of the relative residual. @@ -72,7 +72,7 @@ ## @seealso{bicg, bicgstab, cgs, pcg} ## @end deftypefn -function [x, flag, presn, it, resids] = gmres (A, b, restart, rtol, maxit, M1, M2, x0) +function [x, flag, relres, it, resvec] = gmres (A, b, restart, rtol, maxit, M1, M2, x0) if (nargin < 2 || nargin > 8) print_usage (); @@ -142,10 +142,10 @@ ## begin loop iter = 1; restart_it = restart + 1; - resids = zeros (maxit, 1); - resids(1) = presn; + resvec = zeros (maxit, 1); + resvec(1) = presn; prec_b_norm = norm (Pm1x (b), 2); - flag = 1; + flag = 1; # Default flag is maximum # of iterations exceeded while (iter <= maxit * restart && presn > rtol * prec_b_norm) @@ -174,26 +174,36 @@ x = x_old + V(:, 1:restart_it) * Y(1:restart_it); - resids(iter) = presn; + resvec(iter) = presn; if (norm (x - x_old, inf) <= eps) - flag = 3; - break + flag = 3; # Stagnation: no change between iterations + break; endif restart_it++ ; iter++; endwhile - if (presn > rtol * prec_b_norm) - flag = 0; + if (nargout > 1) + ## Calculate extra outputs as requested + relres = presn / prec_b_norm; + if (relres <= rtol) + flag = 0; # Converged to solution within tolerance + endif + + resvec = resvec(1:iter-1); + it = [ceil(iter / restart), rem(iter, restart)]; endif - resids = resids(1:iter-1); - it = [ceil(iter / restart), rem(iter, restart)]; - endfunction +%!demo +%! dim = 20; +%! A = spdiags ([-ones(dim,1) 2*ones(dim,1) ones(dim,1)], [-1:1], dim, dim); +%! b = ones (dim, 1); +%! [x, flag, relres, iter, resvec] = gmres (A, b, 10, 1e-10, dim, @(x) x ./ diag (A), [], b) + %!shared A, b, dim %! dim = 100; %!test @@ -210,7 +220,7 @@ %! A = spdiags ([[1./(2:2:2*(dim-1)) 0]; 1./(1:2:2*dim-1); [0 1./(2:2:2*(dim-1))]]', -1:1, dim, dim); %! A = A'*A; %! b = rand (dim, 1); -%! [x, resids] = gmres (@(x) A*x, b, dim, 1e-10, dim, @(x) x./diag (A), [], []); +%! [x, resvec] = gmres (@(x) A*x, b, dim, 1e-10, dim, @(x) x./diag (A), [], []); %! assert (x, A\b, 1e-9*norm (x, Inf)); %! x = gmres (@(x) A*x, b, dim, 1e-10, 1e6, @(x) diag (diag (A)) \ x, [], []); %! assert (x, A\b, 1e-9*norm (x, Inf)); @@ -218,3 +228,10 @@ %! x = gmres (@(x) A*x, b, dim, 1e-10, 1e6, @(x) x./diag(A), [], []); %! assert (x, A\b, 1e-7*norm (x, Inf)); + +%!error gmres (1) +%!error gmres (1,2,3,4,5,6,7,8,9) +%!error <A must be> gmres ({1},2) +%!error <A must be a function or matrix> gmres ({1},2) +%!error <M1 must be a function or matrix> gmres (1,2,3,4,5,{6}) +%!error <M2 must be a function or matrix> gmres (1,2,3,4,5,6,{7})
--- a/scripts/statistics/distributions/logncdf.m +++ b/scripts/statistics/distributions/logncdf.m @@ -74,15 +74,15 @@ %!shared x,y %! x = [-1 0 1 e Inf]; %! y = [0, 0, 0.5, 1/2+1/2*erf(1/2), 1]; -%!assert (logncdf (x, zeros (1,5), sqrt(2)*ones (1,5)), y) -%!assert (logncdf (x, 0, sqrt(2)*ones (1,5)), y) -%!assert (logncdf (x, zeros (1,5), sqrt(2)), y) -%!assert (logncdf (x, [0 1 NaN 0 1], sqrt(2)), [0 0 NaN y(4:5)]) -%!assert (logncdf (x, 0, sqrt(2)*[0 NaN Inf 1 1]), [NaN NaN NaN y(4:5)]) -%!assert (logncdf ([x(1:3) NaN x(5)], 0, sqrt(2)), [y(1:3) NaN y(5)]) +%!assert (logncdf (x, zeros (1,5), sqrt(2)*ones (1,5)), y, eps) +%!assert (logncdf (x, 0, sqrt(2)*ones (1,5)), y, eps) +%!assert (logncdf (x, zeros (1,5), sqrt(2)), y, eps) +%!assert (logncdf (x, [0 1 NaN 0 1], sqrt(2)), [0 0 NaN y(4:5)], eps) +%!assert (logncdf (x, 0, sqrt(2)*[0 NaN Inf 1 1]), [NaN NaN NaN y(4:5)], eps) +%!assert (logncdf ([x(1:3) NaN x(5)], 0, sqrt(2)), [y(1:3) NaN y(5)], eps) %% Test class of input preserved -%!assert (logncdf ([x, NaN], 0, sqrt(2)), [y, NaN]) +%!assert (logncdf ([x, NaN], 0, sqrt(2)), [y, NaN], eps) %!assert (logncdf (single ([x, NaN]), 0, sqrt(2)), single ([y, NaN]), eps ("single")) %!assert (logncdf ([x, NaN], single (0), sqrt(2)), single ([y, NaN]), eps ("single")) %!assert (logncdf ([x, NaN], 0, single (sqrt(2))), single ([y, NaN]), eps ("single"))
--- a/scripts/statistics/distributions/stdnormal_inv.m +++ b/scripts/statistics/distributions/stdnormal_inv.m @@ -37,7 +37,7 @@ error ("stdnormal_inv: X must not be complex"); endif - inv = sqrt (2) * erfinv (2 * x - 1); + inv = - sqrt (2) * erfcinv (2 * x); endfunction
--- a/src/display.cc +++ b/src/display.cc @@ -41,6 +41,23 @@ display_info *display_info::instance = 0; +#if defined (HAVE_FRAMEWORK_CARBON) && ! defined (HAVE_CARBON_CGDISPLAYBITSPERPIXEL) +// FIXME - This will only work for MacOS > 10.5. For earlier versions +// this code is not needed (use CGDisplayBitsPerPixel instead). +size_t DisplayBitsPerPixel (CGDirectDisplayID display) +{ + CGDisplayModeRef mode = CGDisplayCopyDisplayMode (display); + CFStringRef pixelEncoding = CGDisplayModeCopyPixelEncoding (mode); + + if (CFStringCompare (pixelEncoding, CFSTR (IO32BitDirectPixels), 0) == 0) + return 32; + else if (CFStringCompare (pixelEncoding, CFSTR (IO16BitDirectPixels), 0) == 0) + return 16; + else + return 8; +} +#endif + void display_info::init (bool query) { @@ -72,16 +89,21 @@ if (display) { +# if defined (HAVE_CARBON_CGDISPLAYBITSPERPIXEL) + // For MacOS < 10.7 use the line below dp = CGDisplayBitsPerPixel (display); +# else + // For MacOS > 10.5 use the line below + dp = DisplayBitsPerPixel (display); +# endif ht = CGDisplayPixelsHigh (display); wd = CGDisplayPixelsWide (display); CGSize sz_mm = CGDisplayScreenSize (display); - - // On modern Mac systems (>= 10.5) CGSize is a struct keeping 2 - // CGFloat values, but the CGFloat typedef is not present on - // older systems, so use double instead. + // For MacOS >= 10.6, CGSize is a struct keeping 2 CGFloat values, + // but the CGFloat typedef is not present on older systems, + // so use double instead. double ht_mm = sz_mm.height; double wd_mm = sz_mm.width;
--- a/src/help.cc +++ b/src/help.cc @@ -144,14 +144,14 @@ "-*- texinfo -*-\n\ @deftypefn {Operator} {} #\n\ Begin comment character.\n\ -@seealso{%, #@{}\n\ +@seealso{%, #@\\{}\n\ @end deftypefn"), pair_type ("%", "-*- texinfo -*-\n\ @deftypefn {Operator} {} %\n\ Begin comment character.\n\ -@seealso{#, %@{}\n\ +@seealso{#, %@\\{}\n\ @end deftypefn"), pair_type ("#{", @@ -160,7 +160,7 @@ Begin block comment. There must be nothing else, other than\n\ whitespace, in the line both before and after @code{#@{}.\n\ It is possible to nest block comments.\n\ -@seealso{%@{, #@}, #}\n\ +@seealso{%@\\{, #@\\}, #}\n\ @end deftypefn"), pair_type ("%{", @@ -169,7 +169,7 @@ Begin block comment. There must be nothing else, other than\n\ whitespace, in the line both before and after @code{%@{}.\n\ It is possible to nest block comments.\n\ -@seealso{#@{, %@}, %}\n\ +@seealso{#@\\{, %@\\}, %}\n\ @end deftypefn"), pair_type ("#}", @@ -178,7 +178,7 @@ Close block comment. There must be nothing else, other than\n\ whitespace, in the line both before and after @code{#@}}.\n\ It is possible to nest block comments.\n\ -@seealso{%@}, #@{, #}\n\ +@seealso{%@\\}, #@\\{, #}\n\ @end deftypefn"), pair_type ("%}", @@ -187,7 +187,7 @@ Close block comment. There must be nothing else, other than\n\ whitespace, in the line both before and after @code{%@}}.\n\ It is possible to nest block comments.\n\ -@seealso{#@}, %@{, %}\n\ +@seealso{#@\\}, %@\\{, %}\n\ @end deftypefn"), pair_type ("...",
--- a/src/mappers.cc +++ b/src/mappers.cc @@ -547,7 +547,7 @@ @end example\n\ \n\ @end ifnottex\n\ -@seealso{erfc, erfcx, erfinv}\n\ +@seealso{erfc, erfcx, erfinv, erfcinv}\n\ @end deftypefn") { octave_value retval; @@ -596,7 +596,7 @@ @example\n\ erf (@var{y}) == @var{x}\n\ @end example\n\ -@seealso{erf, erfc, erfcx}\n\ +@seealso{erf, erfc, erfcx, erfcinv}\n\ @end deftypefn") { octave_value retval; @@ -611,10 +611,10 @@ /* ## middle region %!assert (erf (erfinv ([-0.9 -0.3 0 0.4 0.8])), [-0.9 -0.3 0 0.4 0.8], eps) -%!assert (erf (erfinv (single ([-0.9 -0.3 0 0.4 0.8]))), single ([-0.9 -0.3 0 0.4 0.8]), 1e-8) +%!assert (erf (erfinv (single ([-0.9 -0.3 0 0.4 0.8]))), single ([-0.9 -0.3 0 0.4 0.8]), eps ("single")) ## tail region %!assert (erf (erfinv ([-0.999 -0.99 0.9999 0.99999])), [-0.999 -0.99 0.9999 0.99999], eps) -%!assert (erf (erfinv (single ([-0.999 -0.99 0.9999 0.99999]))), single ([-0.999 -0.99 0.9999 0.99999]), 1e-8) +%!assert (erf (erfinv (single ([-0.999 -0.99 0.9999 0.99999]))), single ([-0.999 -0.99 0.9999 0.99999]), eps ("single")) ## backward - loss of accuracy %!assert (erfinv (erf ([-3 -1 -0.4 0.7 1.3 2.8])), [-3 -1 -0.4 0.7 1.3 2.8], -1e-12) %!assert (erfinv (erf (single ([-3 -1 -0.4 0.7 1.3 2.8]))), single ([-3 -1 -0.4 0.7 1.3 2.8]), -1e-4) @@ -626,6 +626,44 @@ %!error erfinv (1, 2) */ +DEFUN (erfcinv, args, , + "-*- texinfo -*-\n\ +@deftypefn {Mapping Function} {} erfcinv (@var{x})\n\ +Compute the inverse complementary error function, i.e., @var{y} such that\n\ +\n\ +@example\n\ +erfc (@var{y}) == @var{x}\n\ +@end example\n\ +@seealso{erfc, erf, erfcx, erfinv}\n\ +@end deftypefn") +{ + octave_value retval; + if (args.length () == 1) + retval = args(0).erfcinv (); + else + print_usage (); + + return retval; +} + +/* +## middle region +%!assert (erfc (erfcinv ([1.9 1.3 1 0.6 0.2])), [1.9 1.3 1 0.6 0.2], eps) +%!assert (erfc (erfcinv (single ([1.9 1.3 1 0.6 0.2]))), single ([1.9 1.3 1 0.6 0.2]), eps ("single")) +## tail region +%!assert (erfc (erfcinv ([0.001 0.01 1.9999 1.99999])), [0.001 0.01 1.9999 1.99999], eps) +%!assert (erfc (erfcinv (single ([0.001 0.01 1.9999 1.99999]))), single ([0.001 0.01 1.9999 1.99999]), eps ("single")) +## backward - loss of accuracy +%!assert (erfcinv (erfc ([-3 -1 -0.4 0.7 1.3 2.8])), [-3 -1 -0.4 0.7 1.3 2.8], -1e-12) +%!assert (erfcinv (erfc (single ([-3 -1 -0.4 0.7 1.3 2.8]))), single ([-3 -1 -0.4 0.7 1.3 2.8]), -1e-4) +## exceptional +%!assert (erfcinv ([2, 0, -0.1, 2.1]), [-Inf, Inf, NaN, NaN]) +%!error erfcinv (1+2i) + +%!error erfcinv () +%!error erfcinv (1, 2) +*/ + DEFUN (erfc, args, , "-*- texinfo -*-\n\ @deftypefn {Mapping Function} {} erfc (@var{z})\n\ @@ -636,7 +674,7 @@ @ifnottex\n\ @w{@code{1 - erf (@var{z})}}.\n\ @end ifnottex\n\ -@seealso{erfcx, erf, erfinv}\n\ +@seealso{erfcinv, erfcx, erf, erfinv}\n\ @end deftypefn") { octave_value retval; @@ -673,7 +711,7 @@ @end example\n\ \n\ @end ifnottex\n\ -@seealso{erfc, erf, erfinv}\n\ +@seealso{erfc, erf, erfinv, erfcinv}\n\ @end deftypefn") { octave_value retval;
--- a/src/ov-base.cc +++ b/src/ov-base.cc @@ -1191,6 +1191,7 @@ "cosh", "erf", "erfinv", + "erfcinv", "erfc", "exp", "expm1",
--- a/src/ov-base.h +++ b/src/ov-base.h @@ -688,6 +688,7 @@ umap_cosh, umap_erf, umap_erfinv, + umap_erfcinv, umap_erfc, umap_erfcx, umap_exp,
--- a/src/ov-float.cc +++ b/src/ov-float.cc @@ -291,6 +291,7 @@ SCALAR_MAPPER (atanh, rc_atanh); SCALAR_MAPPER (erf, ::erff); SCALAR_MAPPER (erfinv, ::erfinv); + SCALAR_MAPPER (erfcinv, ::erfcinv); SCALAR_MAPPER (erfc, ::erfcf); SCALAR_MAPPER (erfcx, ::erfcx); SCALAR_MAPPER (gamma, xgamma);
--- a/src/ov-flt-re-mat.cc +++ b/src/ov-flt-re-mat.cc @@ -784,6 +784,7 @@ RC_ARRAY_MAPPER (atanh, FloatComplex, rc_atanh); ARRAY_MAPPER (erf, float, ::erff); ARRAY_MAPPER (erfinv, float, ::erfinv); + ARRAY_MAPPER (erfcinv, float, ::erfcinv); ARRAY_MAPPER (erfc, float, ::erfcf); ARRAY_MAPPER (erfcx, float, ::erfcx); ARRAY_MAPPER (gamma, float, xgamma);
--- a/src/ov-re-mat.cc +++ b/src/ov-re-mat.cc @@ -907,6 +907,7 @@ RC_ARRAY_MAPPER (atanh, Complex, rc_atanh); ARRAY_MAPPER (erf, double, ::erf); ARRAY_MAPPER (erfinv, double, ::erfinv); + ARRAY_MAPPER (erfcinv, double, ::erfcinv); ARRAY_MAPPER (erfc, double, ::erfc); ARRAY_MAPPER (erfcx, double, ::erfcx); ARRAY_MAPPER (gamma, double, xgamma);
--- a/src/ov-re-sparse.cc +++ b/src/ov-re-sparse.cc @@ -910,6 +910,7 @@ ARRAY_MAPPER (atanh, Complex, rc_atanh); ARRAY_MAPPER (erf, double, ::erf); ARRAY_MAPPER (erfinv, double, ::erfinv); + ARRAY_MAPPER (erfcinv, double, ::erfcinv); ARRAY_MAPPER (erfc, double, ::erfc); ARRAY_MAPPER (gamma, double, xgamma); ARRAY_MAPPER (lgamma, Complex, rc_lgamma);