comparison src/det.cc @ 636:fae2bd91c027

[project @ 1994-08-23 18:39:50 by jwe]
author jwe
date Tue, 23 Aug 1994 18:39:50 +0000
parents 8e4e7e5f307e
children 0a81458ef677
comparison
equal deleted inserted replaced
635:5338832d2cf6 636:fae2bd91c027
30 30
31 #include "tree-const.h" 31 #include "tree-const.h"
32 #include "user-prefs.h" 32 #include "user-prefs.h"
33 #include "gripes.h" 33 #include "gripes.h"
34 #include "error.h" 34 #include "error.h"
35 #include "utils.h"
35 #include "help.h" 36 #include "help.h"
36 #include "defun-dld.h" 37 #include "defun-dld.h"
37 38
38 DEFUN_DLD ("det", Fdet, Sdet, 2, 1, 39 DEFUN_DLD ("det", Fdet, Sdet, 2, 1,
39 "det (X): determinant of a square matrix") 40 "det (X): determinant of a square matrix")
46 { 47 {
47 print_usage ("det"); 48 print_usage ("det");
48 return retval; 49 return retval;
49 } 50 }
50 51
51 tree_constant tmp = args(1).make_numeric ();; 52 tree_constant arg = args(1);
52 53
53 int nr = tmp.rows (); 54 int nr = arg.rows ();
54 int nc = tmp.columns (); 55 int nc = arg.columns ();
55 if (nr == 0 || nc == 0) 56
57 if (nr == 0 && nc == 0)
56 { 58 {
57 int flag = user_pref.propagate_empty_matrices; 59 retval = 1.0;
58 if (flag < 0) 60 return retval;
59 gripe_empty_arg ("det", 0);
60 else if (flag == 0)
61 gripe_empty_arg ("det", 1);
62 } 61 }
63 62
64 if (nr == 0 && nc == 0) 63 if (empty_arg ("det", nr, nc) < 0)
65 return 1.0; 64 return retval;
66 65
67 if (tmp.is_real_matrix ()) 66 if (nr != nc)
68 { 67 {
69 Matrix m = tmp.matrix_value (); 68 gripe_square_matrix_required ("det");
70 if (m.rows () == m.columns ()) 69 return retval;
70 }
71
72 if (arg.is_real_type ())
73 {
74 Matrix m = arg.matrix_value ();
75
76 if (! error_state)
71 { 77 {
72 int info; 78 int info;
73 double rcond = 0.0; 79 double rcond = 0.0;
80
74 DET det = m.determinant (info, rcond); 81 DET det = m.determinant (info, rcond);
82
75 double d = 0.0; 83 double d = 0.0;
84
76 if (info == -1) 85 if (info == -1)
77 warning ("det: matrix singular to machine precision, rcond = %g", 86 warning ("det: matrix singular to machine precision, rcond = %g",
78 rcond); 87 rcond);
79 else 88 else
80 d = det.value (); 89 d = det.value ();
81 90
82 retval = d; 91 retval = d;
83 } 92 }
84 else
85 gripe_square_matrix_required ("det");
86 } 93 }
87 else if (tmp.is_complex_matrix ()) 94 else if (arg.is_complex_matrix ())
88 { 95 {
89 ComplexMatrix m = tmp.complex_matrix_value (); 96 ComplexMatrix m = arg.complex_matrix_value ();
90 if (m.rows () == m.columns ()) 97
98 if (! error_state)
91 { 99 {
92 int info; 100 int info;
93 double rcond = 0.0; 101 double rcond = 0.0;
102
94 ComplexDET det = m.determinant (info, rcond); 103 ComplexDET det = m.determinant (info, rcond);
104
95 Complex c = 0.0; 105 Complex c = 0.0;
106
96 if (info == -1) 107 if (info == -1)
97 warning ("det: matrix singular to machine precision, rcond = %g", 108 warning ("det: matrix singular to machine precision, rcond = %g",
98 rcond); 109 rcond);
99 else 110 else
100 c = det.value (); 111 c = det.value ();
101 112
102 retval = c; 113 retval = c;
103 } 114 }
104 else
105 gripe_square_matrix_required ("det");
106 }
107 else if (tmp.is_real_scalar ())
108 {
109 double d = tmp.double_value ();
110 retval = d;
111 }
112 else if (tmp.is_complex_scalar ())
113 {
114 Complex c = tmp.complex_value ();
115 retval = c;
116 } 115 }
117 else 116 else
118 { 117 {
119 gripe_wrong_type_arg ("det", tmp); 118 gripe_wrong_type_arg ("det", arg);
120 } 119 }
121 120
122 return retval; 121 return retval;
123 } 122 }
124 123