1%% Copyright (C) 2014-2017, 2019 Colin B. Macdonald 2%% 3%% This file is part of OctSymPy. 4%% 5%% OctSymPy is free software; you can redistribute it and/or modify 6%% it under the terms of the GNU General Public License as published 7%% by the Free Software Foundation; either version 3 of the License, 8%% or (at your option) any later version. 9%% 10%% This software is distributed in the hope that it will be useful, 11%% but WITHOUT ANY WARRANTY; without even the implied warranty 12%% of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See 13%% the GNU General Public License for more details. 14%% 15%% You should have received a copy of the GNU General Public 16%% License along with this software; see the file COPYING. 17%% If not, see <http://www.gnu.org/licenses/>. 18 19%% -*- texinfo -*- 20%% @documentencoding UTF-8 21%% @defmethod @@sym double (@var{x}) 22%% Convert symbolic to doubles. 23%% 24%% Example: 25%% @example 26%% @group 27%% x = sym(1) / 3 28%% @result{} x = (sym) 1/3 29%% @c doctest: +SKIP_IF(compare_versions (OCTAVE_VERSION(), '6.0.0', '<')) 30%% double (x) 31%% @result{} ans = 0.3333 32%% @end group 33%% @end example 34%% 35%% Despite the name, this is one way to convert a complex sym to 36%% floating point: 37%% @example 38%% @group 39%% z = sym(4i) - 3; 40%% double (z) 41%% @result{} ans = -3 + 4i 42%% @end group 43%% @end example 44%% 45%% If conversion fails, you get an error: 46%% @example 47%% @group 48%% syms x 49%% double (x) 50%% @print{} ??? ... can't convert expression ... 51%% @end group 52%% @end example 53%% 54%% @seealso{sym, vpa} 55%% @end defmethod 56 57 58function y = double(x) 59 60 % FIXME: port to uniop? 61 62 if ~(isscalar(x)) 63 % sympy N() works fine on matrices but it gives objects like "Matrix([[1.0,2.0]])" 64 y = zeros(size(x)); 65 for j = 1:numel(x) 66 % temp = x(j) (Issue #17) 67 idx.type = '()'; 68 idx.subs = {j}; 69 temp = double(subsref(x,idx)); 70 if (isempty(temp)) 71 y = []; 72 return 73 end 74 y(j) = temp; 75 end 76 return 77 end 78 79 cmd = { '(x,) = _ins' 80 'if x == zoo:' % zoo -> Inf + Infi 81 ' return (float(sp.oo), float(sp.oo))' 82 'if x == nan:' 83 ' return (float(nan), 0.0)' 84 'x = complex(x)' 85 'return (x.real, x.imag)' 86 }; 87 88 [A, B] = pycall_sympy__ (cmd, x); 89 90 %y = A + B*i; % not quite the same for Inf + InFi 91 if (B == 0.0) 92 y = A; 93 else 94 y = complex(A, B); 95 end 96 97end 98 99 100%!test 101%! % numeric scalar 102%! a = double(sym(10)); 103%! assert (a == 10) 104%! assert (isa (a, 'double')) 105 106%!test 107%! % numeric vectors 108%! a = double(sym([10 12])); 109%! assert (isequal (a, [10 12])) 110%! assert (isa (a, 'double')) 111 112%!test 113%! % complex 114%! a = 3 + 4i; 115%! b = sym(a); 116%! assert (isequal (double (b), a)) 117 118%!xtest 119%! % unexpected, precisely same floating point 120%! a = 3 + 4i; 121%! b = sym(a); 122%! assert (isequal (double (b/pi), a/pi)) 123 124%!test 125%! % floating point 126%! x = sqrt(sym(2)); 127%! assert( abs(double(x) - sqrt(2)) < 2*eps) 128%! x = sym(pi); 129%! assert( abs(double(x) - pi) < 2*eps) 130 131%!test 132%! oo = sym(inf); 133%! assert( double(oo) == inf ) 134%! assert( double(-oo) == -inf ) 135%! assert( isnan(double(0*oo)) ) 136 137%!test 138%! zoo = sym('zoo'); 139%! assert (double(zoo) == complex(inf, inf)) 140 141%!test 142%! zoo = sym('zoo'); 143%! assert (double(-zoo) == double(zoo) ) 144%! assert( isnan(double(0*zoo)) ) 145 146%!test 147%! % nan 148%! snan = sym(nan); 149%! assert( isnan(double(snan))) 150 151%!test 152%! % don't want NaN+NaNi 153%! snan = sym(nan); 154%! assert (isreal (double (snan))) 155 156%!test 157%! % arrays 158%! a = [1 2; 3 4]; 159%! assert( isequal( double(sym(a)), a )) 160%! assert( isequal( double(sym(a)), a )) 161 162%! % should fail with error for non-double 163%!error <TypeError> syms x; double(x) 164%!error <TypeError> syms x; double([1 2 x]) 165