1## Copyright 2014-2016 Oliver Heimlich 2## 3## This program is free software; you can redistribute it and/or modify 4## it under the terms of the GNU General Public License as published by 5## the Free Software Foundation; either version 3 of the License, or 6## (at your option) any later version. 7## 8## This program is distributed in the hope that it will be useful, 9## but WITHOUT ANY WARRANTY; without even the implied warranty of 10## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 11## GNU General Public License for more details. 12## 13## You should have received a copy of the GNU General Public License 14## along with this program; if not, see <http://www.gnu.org/licenses/>. 15 16## -*- texinfo -*- 17## @documentencoding UTF-8 18## @defmethod {@@infsup} cos (@var{X}) 19## 20## Compute the cosine in radians. 21## 22## Accuracy: The result is a tight enclosure. 23## 24## @example 25## @group 26## cos (infsup (1)) 27## @result{} ans ⊂ [0.5403, 0.54031] 28## @end group 29## @end example 30## @seealso{@@infsup/acos, @@infsup/sec, @@infsup/cosh} 31## @end defmethod 32 33## Author: Oliver Heimlich 34## Keywords: interval 35## Created: 2014-10-05 36 37function x = cos (x) 38 39 if (nargin ~= 1) 40 print_usage (); 41 return 42 endif 43 44 l = u = sinsignl = sinsignu = zeros (size (x.inf)); 45 46 ## Check, if wid (x) is certainly greater than 2*pi. This can save the 47 ## computation if some cosine values. 48 width = mpfr_function_d ('minus', -inf, x.sup, x.inf); 49 persistent pi = infsup ("pi"); 50 persistent twopi = 2 .* pi; 51 certainlyfullperiod = width >= sup (twopi); 52 l(certainlyfullperiod) = -1; 53 u(certainlyfullperiod) = 1; 54 55 possiblynotfullperiod = not (certainlyfullperiod); 56 if (__check_crlibm__ ()) 57 l(possiblynotfullperiod) = min (crlibm_function ('cos', -inf, ... 58 x.inf(possiblynotfullperiod)), ... 59 crlibm_function ('cos', -inf, ... 60 x.sup(possiblynotfullperiod))); 61 u(possiblynotfullperiod) = max (crlibm_function ('cos', inf, ... 62 x.inf(possiblynotfullperiod)), ... 63 crlibm_function ('cos', inf, ... 64 x.sup(possiblynotfullperiod))); 65 66 ## We use sign (-sin) to know the gradient at the boundaries. 67 sinsignl(possiblynotfullperiod) = sign (-crlibm_function ('sin', .5, ... 68 x.inf(possiblynotfullperiod))); 69 sinsignu(possiblynotfullperiod) = sign (-crlibm_function ('sin', .5, ... 70 x.sup(possiblynotfullperiod))); 71 else 72 l(possiblynotfullperiod) = min (mpfr_function_d ('cos', -inf, ... 73 x.inf(possiblynotfullperiod)), ... 74 mpfr_function_d ('cos', -inf, ... 75 x.sup(possiblynotfullperiod))); 76 u(possiblynotfullperiod) = max (mpfr_function_d ('cos', inf, ... 77 x.inf(possiblynotfullperiod)), ... 78 mpfr_function_d ('cos', inf, ... 79 x.sup(possiblynotfullperiod))); 80 81 ## We use sign (-sin) to know the gradient at the boundaries. 82 sinsignl(possiblynotfullperiod) = sign (-mpfr_function_d ('sin', .5, ... 83 x.inf(possiblynotfullperiod))); 84 sinsignu(possiblynotfullperiod) = sign (-mpfr_function_d ('sin', .5, ... 85 x.sup(possiblynotfullperiod))); 86 endif 87 88 ## In case of sign (-sin) == 0, we conservatively use sign (-sin) of nextout. 89 sinsignl(sinsignl == 0) = (-1) .* sign (l(sinsignl == 0)); 90 sinsignu(sinsignu == 0) = sign (u(sinsignu == 0)); 91 92 containsinf = possiblynotfullperiod & ((sinsignl == -1 & sinsignu == 1) | ... 93 (sinsignl == sinsignu & ... 94 width >= sup (pi))) ... 95 & ne (0, x); 96 l(containsinf) = -1; 97 98 containssup = possiblynotfullperiod & ((sinsignl == 1 & sinsignu == -1) | ... 99 (sinsignl == sinsignu & ... 100 width >= sup (pi))); 101 u(containssup) = 1; 102 103 emptyresult = isempty (x); 104 l(emptyresult) = inf; 105 u(emptyresult) = -inf; 106 107 l(l == 0) = -0; 108 109 x.inf = l; 110 x.sup = u; 111 112endfunction 113 114%!# from the documentation string 115%!assert (cos (infsup (1)) == "[0x1.14A280FB5068Bp-1, 0x1.14A280FB5068Cp-1]"); 116 117%!shared testdata 118%! # Load compiled test data (from src/test/*.itl) 119%! testdata = load (file_in_loadpath ("test/itl.mat")); 120 121%!test 122%! # Scalar evaluation 123%! testcases = testdata.NoSignal.infsup.cos; 124%! for testcase = [testcases]' 125%! assert (isequaln (... 126%! cos (testcase.in{1}), ... 127%! testcase.out)); 128%! endfor 129 130%!test 131%! # Vector evaluation 132%! testcases = testdata.NoSignal.infsup.cos; 133%! in1 = vertcat (vertcat (testcases.in){:, 1}); 134%! out = vertcat (testcases.out); 135%! assert (isequaln (cos (in1), out)); 136 137%!test 138%! # N-dimensional array evaluation 139%! testcases = testdata.NoSignal.infsup.cos; 140%! in1 = vertcat (vertcat (testcases.in){:, 1}); 141%! out = vertcat (testcases.out); 142%! # Reshape data 143%! i = -1; 144%! do 145%! i = i + 1; 146%! testsize = factor (numel (in1) + i); 147%! until (numel (testsize) > 2) 148%! in1 = reshape ([in1; in1(1:i)], testsize); 149%! out = reshape ([out; out(1:i)], testsize); 150%! assert (isequaln (cos (in1), out)); 151