Mercurial > octave
view libinterp/corefcn/__luinc__.cc @ 22898:9baa19102908
refactor display and disp functions (bug #49794)
* pr-output.cc (Fdisp, Ffdisp): Tag with dispatch classes.
(Fdisplay): New function.
* ov-class.cc (octave_class::print_with_name):
Simply call octave_base_value::print_with_name.
* ov-classdef.cc (octave_classdef::print): Simply call print_raw.
(octave_classdef::print_with_name):
Simply call octave_base_value::print_with_name.
* variables.cc (bind_ans): Call display function to print result.
* pt-assign.cc (tree_simple_assignment::rvalue1,
tree_multi_assignment::rvalue): Likewise.
* pt-id.cc (tree_identifier::rvalue): Likewise.
* display.m: Delete.
* scripts/general/module.mk: Update.
author | John W. Eaton <jwe@octave.org> |
---|---|
date | Fri, 16 Dec 2016 00:10:27 -0500 |
parents | 3a2b891d0b33 |
children | ef4d915df748 |
line wrap: on
line source
/* Copyright (C) 2005-2016 David Bateman This file is part of Octave. Octave is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 3 of the License, or (at your option) any later version. Octave is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with Octave; see the file COPYING. If not, see <http://www.gnu.org/licenses/>. */ #if defined (HAVE_CONFIG_H) # include "config.h" #endif #include "defun.h" #include "error.h" #include "errwarn.h" #include "ovl.h" #include "utils.h" #include "oct-map.h" #include "MatrixType.h" #include "sparse-lu.h" #include "ov-re-sparse.h" #include "ov-cx-sparse.h" // FIXME: Deprecated in 4.0 and should be removed in 4.4. DEFUN (__luinc__, args, nargout, doc: /* -*- texinfo -*- @deftypefn {} {[@var{L}, @var{U}, @var{P}, @var{Q}] =} __luinc__ (@var{A}, '0') @deftypefnx {} {[@var{L}, @var{U}, @var{P}, @var{Q}] =} __luinc__ (@var{A}, @var{droptol}) @deftypefnx {} {[@var{L}, @var{U}, @var{P}, @var{Q}] =} __luinc__ (@var{A}, @var{opts}) Internal implementation of @code{luinc}. See documentation for @code{luinc}. @seealso{luinc} @end deftypefn */) { int nargin = args.length (); if (nargin < 2 || nargin > 3) print_usage (); if (! args(0).is_sparse_type ()) error ("luinc: matrix A must be sparse"); bool zero_level = false; bool milu = false; bool udiag = false; Matrix thresh; double droptol = -1.0; bool vecout = false; if (args(1).is_string ()) { if (args(1).string_value () == "0") zero_level = true; else error ("luinc: unrecognized string argument"); } else if (args(1).is_map ()) { octave_scalar_map map = args(1).xscalar_map_value ("luinc: OPTS must be a scalar structure"); octave_value tmp; tmp = map.getfield ("droptol"); if (tmp.is_defined ()) droptol = tmp.double_value (); tmp = map.getfield ("milu"); if (tmp.is_defined ()) { double val = tmp.double_value (); milu = (val == 0.0 ? false : true); } tmp = map.getfield ("udiag"); if (tmp.is_defined ()) { double val = tmp.double_value (); udiag = (val == 0.0 ? false : true); } tmp = map.getfield ("thresh"); if (tmp.is_defined ()) { thresh = tmp.matrix_value (); if (thresh.numel () == 1) { thresh.resize (1, 2); thresh(1) = thresh(0); } else if (thresh.numel () != 2) error ("luinc: THRESH must be a 1 or 2-element vector"); } } else droptol = args(1).double_value (); if (nargin == 3) { std::string tmp = args(2).string_value (); if (tmp == "vector") vecout = true; else error ("luinc: unrecognized string argument"); } // FIXME: Add code for zero-level factorization if (zero_level) error ("luinc: zero-level factorization not implemented"); octave_value_list retval; if (args(0).is_real_type ()) { SparseMatrix sm = args(0).sparse_matrix_value (); octave_idx_type sm_nr = sm.rows (); octave_idx_type sm_nc = sm.cols (); ColumnVector Qinit (sm_nc); for (octave_idx_type i = 0; i < sm_nc; i++) Qinit(i) = i; switch (nargout) { case 0: case 1: case 2: { octave::math::sparse_lu<SparseMatrix> fact (sm, Qinit, thresh, false, true, droptol, milu, udiag); SparseMatrix P = fact.Pr (); SparseMatrix L = P.transpose () * fact.L (); retval(1) = octave_value (fact.U (), MatrixType (MatrixType::Upper)); retval(0) = octave_value (L, MatrixType (MatrixType::Permuted_Lower, sm_nr, fact.row_perm ())); } break; case 3: { octave::math::sparse_lu<SparseMatrix> fact (sm, Qinit, thresh, false, true, droptol, milu, udiag); if (vecout) retval(2) = fact.Pr_vec (); else retval(2) = fact.Pr_mat (); retval(1) = octave_value (fact.U (), MatrixType (MatrixType::Upper)); retval(0) = octave_value (fact.L (), MatrixType (MatrixType::Lower)); } break; case 4: default: { octave::math::sparse_lu<SparseMatrix> fact (sm, Qinit, thresh, false, false, droptol, milu, udiag); if (vecout) { retval(3) = fact.Pc_vec (); retval(2) = fact.Pr_vec (); } else { retval(3) = fact.Pc_mat (); retval(2) = fact.Pr_mat (); } retval(1) = octave_value (fact.U (), MatrixType (MatrixType::Upper)); retval(0) = octave_value (fact.L (), MatrixType (MatrixType::Lower)); } break; } } else { SparseComplexMatrix sm = args(0).sparse_complex_matrix_value (); octave_idx_type sm_nr = sm.rows (); octave_idx_type sm_nc = sm.cols (); ColumnVector Qinit (sm_nc); for (octave_idx_type i = 0; i < sm_nc; i++) Qinit(i) = i; switch (nargout) { case 0: case 1: case 2: { octave::math::sparse_lu<SparseComplexMatrix> fact (sm, Qinit, thresh, false, true, droptol, milu, udiag); SparseMatrix P = fact.Pr (); SparseComplexMatrix L = P.transpose () * fact.L (); retval(1) = octave_value (fact.U (), MatrixType (MatrixType::Upper)); retval(0) = octave_value (L, MatrixType (MatrixType::Permuted_Lower, sm_nr, fact.row_perm ())); } break; case 3: { octave::math::sparse_lu<SparseComplexMatrix> fact (sm, Qinit, thresh, false, true, droptol, milu, udiag); if (vecout) retval(2) = fact.Pr_vec (); else retval(2) = fact.Pr_mat (); retval(1) = octave_value (fact.U (), MatrixType (MatrixType::Upper)); retval(0) = octave_value (fact.L (), MatrixType (MatrixType::Lower)); } break; case 4: default: { octave::math::sparse_lu<SparseComplexMatrix> fact (sm, Qinit, thresh, false, false, droptol, milu, udiag); if (vecout) { retval(3) = fact.Pc_vec (); retval(2) = fact.Pr_vec (); } else { retval(3) = fact.Pc_mat (); retval(2) = fact.Pr_mat (); } retval(1) = octave_value (fact.U (), MatrixType (MatrixType::Upper)); retval(0) = octave_value (fact.L (), MatrixType (MatrixType::Lower)); } break; } } return retval; }