annotate src/inv.cc @ 751:c9f3cb77478d

[project @ 1994-09-30 18:46:17 by jwe]
author jwe
date Fri, 30 Sep 1994 18:46:27 +0000
parents e81d3a66725e
children 33ae9f667c6d
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
7
d68036bcad4c [project @ 1993-08-08 03:53:43 by jwe]
jwe
parents: 1
diff changeset
1 // f-inv.cc -*- C++ -*-
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
2 /*
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
3
277
1781618ec883 [project @ 1994-01-07 19:19:54 by jwe]
jwe
parents: 240
diff changeset
4 Copyright (C) 1993, 1994 John W. Eaton
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
5
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
6 This file is part of Octave.
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
7
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
8 Octave is free software; you can redistribute it and/or modify it
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
9 under the terms of the GNU General Public License as published by the
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
10 Free Software Foundation; either version 2, or (at your option) any
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
11 later version.
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
12
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
13 Octave is distributed in the hope that it will be useful, but WITHOUT
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
14 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
15 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
16 for more details.
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
17
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
18 You should have received a copy of the GNU General Public License
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
19 along with Octave; see the file COPYING. If not, write to the Free
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
20 Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
21
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
22 */
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
23
240
a99f28f5e351 [project @ 1993-11-30 20:24:36 by jwe]
jwe
parents: 217
diff changeset
24 #ifdef HAVE_CONFIG_H
a99f28f5e351 [project @ 1993-11-30 20:24:36 by jwe]
jwe
parents: 217
diff changeset
25 #include "config.h"
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
26 #endif
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
27
453
393e95f46b51 [project @ 1994-06-06 00:05:20 by jwe]
jwe
parents: 277
diff changeset
28 #include "dMatrix.h"
393e95f46b51 [project @ 1994-06-06 00:05:20 by jwe]
jwe
parents: 277
diff changeset
29 #include "CMatrix.h"
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
30
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
31 #include "tree-const.h"
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
32 #include "user-prefs.h"
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
33 #include "gripes.h"
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
34 #include "error.h"
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
35 #include "utils.h"
544
20fbad23ae51 [project @ 1994-07-22 05:04:44 by jwe]
jwe
parents: 519
diff changeset
36 #include "help.h"
519
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
37 #include "defun-dld.h"
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
38
701
0a81458ef677 [project @ 1994-09-15 02:23:24 by jwe]
jwe
parents: 636
diff changeset
39 DEFUN_DLD_BUILTIN ("inv", Finv, Sinv, 2, 1,
519
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
40 "inv (X): inverse of a square matrix")
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
41 {
519
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
42 Octave_object retval;
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
43
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
44 int nargin = args.length ();
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
45
712
36ba0576bd1b [project @ 1994-09-19 14:18:15 by jwe]
jwe
parents: 701
diff changeset
46 if (nargin != 1)
519
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
47 {
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
48 print_usage ("inv");
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
49 return retval;
b9284136189a [project @ 1994-07-19 14:40:20 by jwe]
jwe
parents: 516
diff changeset
50 }
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
51
712
36ba0576bd1b [project @ 1994-09-19 14:18:15 by jwe]
jwe
parents: 701
diff changeset
52 tree_constant arg = args(0);
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
53
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
54 int nr = arg.rows ();
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
55 int nc = arg.columns ();
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
56
718
e81d3a66725e [project @ 1994-09-21 14:58:18 by jwe]
jwe
parents: 712
diff changeset
57 int arg_is_empty = empty_arg ("inverse", nr, nc);
e81d3a66725e [project @ 1994-09-21 14:58:18 by jwe]
jwe
parents: 712
diff changeset
58
e81d3a66725e [project @ 1994-09-21 14:58:18 by jwe]
jwe
parents: 712
diff changeset
59 if (arg_is_empty < 0)
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
60 return retval;
718
e81d3a66725e [project @ 1994-09-21 14:58:18 by jwe]
jwe
parents: 712
diff changeset
61 else if (arg_is_empty > 0)
e81d3a66725e [project @ 1994-09-21 14:58:18 by jwe]
jwe
parents: 712
diff changeset
62 return Matrix ();
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
63
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
64 if (nr != nc)
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
65 {
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
66 gripe_square_matrix_required ("inverse");
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
67 return retval;
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
68 }
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
69
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
70 if (arg.is_real_type ())
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
71 {
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
72 Matrix m = arg.matrix_value ();
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
73
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
74 if (! error_state)
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
75 {
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
76 int info;
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
77 double rcond = 0.0;
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
78
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
79 Matrix minv = m.inverse (info, rcond);
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
80
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
81 if (info == -1)
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
82 warning ("inverse: matrix singular to machine precision,\
217
cd2e49fb4f6b [project @ 1993-11-13 09:16:01 by jwe]
jwe
parents: 162
diff changeset
83 rcond = %g", rcond);
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
84 else
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
85 retval = minv;
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
86 }
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
87 }
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
88 else if (arg.is_complex_type ())
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
89 {
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
90 ComplexMatrix m = arg.complex_matrix_value ();
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
91
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
92 if (! error_state)
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
93 {
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
94 int info;
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
95 double rcond = 0.0;
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
96
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
97 ComplexMatrix minv = m.inverse (info, rcond);
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
98
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
99 if (info == -1)
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
100 warning ("inverse: matrix singular to machine precision,\
217
cd2e49fb4f6b [project @ 1993-11-13 09:16:01 by jwe]
jwe
parents: 162
diff changeset
101 rcond = %g", rcond);
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
102 else
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
103 retval = minv;
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
104 }
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
105 }
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
106 else
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
107 {
636
fae2bd91c027 [project @ 1994-08-23 18:39:50 by jwe]
jwe
parents: 620
diff changeset
108 gripe_wrong_type_arg ("inv", arg);
620
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
109 }
8e4e7e5f307e [project @ 1994-08-16 04:36:32 by jwe]
jwe
parents: 544
diff changeset
110
1
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
111 return retval;
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
112 }
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
113
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
114 /*
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
115 ;;; Local Variables: ***
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
116 ;;; mode: C++ ***
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
117 ;;; page-delimiter: "^/\\*" ***
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
118 ;;; End: ***
78fd87e624cb [project @ 1993-08-08 01:13:40 by jwe]
jwe
parents:
diff changeset
119 */