Mercurial > octave-nkf

comparison src/DLD-FUNCTIONS/qr.cc @ 10154:40dfc0c99116

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
DLD-FUNCTIONS/*.cc: untabify
author John W. Eaton <jwe@octave.org>
date Wed, 20 Jan 2010 17:33:41 -0500
parents 09da0bd91412
children d0ce5e973937
comparison
equal deleted inserted replaced
10153:2c28f9d0360f 10154:40dfc0c99116
214 bool economy = false; 214 bool economy = false;
215 bool is_cmplx = false; 215 bool is_cmplx = false;
216 int have_b = 0; 216 int have_b = 0;
217 217
218 if (arg.is_complex_type ()) 218 if (arg.is_complex_type ())
219 is_cmplx = true; 219 is_cmplx = true;
220 if (nargin > 1) 220 if (nargin > 1)
221 { 221 {
222 have_b = 1; 222 have_b = 1;
223 if (args(nargin-1).is_scalar_type ()) 223 if (args(nargin-1).is_scalar_type ())
224 { 224 {
225 int val = args(nargin-1).int_value (); 225 int val = args(nargin-1).int_value ();
226 if (val == 0) 226 if (val == 0)
227 { 227 {
228 economy = true; 228 economy = true;
229 have_b = (nargin > 2 ? 2 : 0); 229 have_b = (nargin > 2 ? 2 : 0);
230 } 230 }
231 } 231 }
232 if (have_b > 0 && args(have_b).is_complex_type ()) 232 if (have_b > 0 && args(have_b).is_complex_type ())
233 is_cmplx = true; 233 is_cmplx = true;
234 } 234 }
235 235
236 if (!error_state) 236 if (!error_state)
237 { 237 {
238 if (have_b && nargout < 2) 238 if (have_b && nargout < 2)
239 error ("qr: incorrect number of output arguments"); 239 error ("qr: incorrect number of output arguments");
240 else if (is_cmplx) 240 else if (is_cmplx)
241 { 241 {
242 SparseComplexQR q (arg.sparse_complex_matrix_value ()); 242 SparseComplexQR q (arg.sparse_complex_matrix_value ());
243 if (!error_state) 243 if (!error_state)
244 { 244 {
245 if (have_b > 0) 245 if (have_b > 0)
246 { 246 {
247 retval(1) = q.R (economy); 247 retval(1) = q.R (economy);
248 retval(0) = q.C (args(have_b).complex_matrix_value ()); 248 retval(0) = q.C (args(have_b).complex_matrix_value ());
249 if (arg.rows() < arg.columns()) 249 if (arg.rows() < arg.columns())
250 warning ("qr: non minimum norm solution for under-determined problem"); 250 warning ("qr: non minimum norm solution for under-determined problem");
251 } 251 }
252 else if (nargout > 1) 252 else if (nargout > 1)
253 { 253 {
254 retval(1) = q.R (economy); 254 retval(1) = q.R (economy);
255 retval(0) = q.Q (); 255 retval(0) = q.Q ();
256 } 256 }
257 else 257 else
258 retval(0) = q.R (economy); 258 retval(0) = q.R (economy);
259 } 259 }
260 } 260 }
261 else 261 else
262 { 262 {
263 SparseQR q (arg.sparse_matrix_value ()); 263 SparseQR q (arg.sparse_matrix_value ());
264 if (!error_state) 264 if (!error_state)
265 { 265 {
266 if (have_b > 0) 266 if (have_b > 0)
267 { 267 {
268 retval(1) = q.R (economy); 268 retval(1) = q.R (economy);
269 retval(0) = q.C (args(have_b).matrix_value ()); 269 retval(0) = q.C (args(have_b).matrix_value ());
270 if (args(0).rows() < args(0).columns()) 270 if (args(0).rows() < args(0).columns())
271 warning ("qr: non minimum norm solution for under-determined problem"); 271 warning ("qr: non minimum norm solution for under-determined problem");
272 } 272 }
273 else if (nargout > 1) 273 else if (nargout > 1)
274 { 274 {
275 retval(1) = q.R (economy); 275 retval(1) = q.R (economy);
276 retval(0) = q.Q (); 276 retval(0) = q.Q ();
277 } 277 }
278 else 278 else
279 retval(0) = q.R (economy); 279 retval(0) = q.R (economy);
280 } 280 }
281 } 281 }
282 } 282 }
283 } 283 }
284 else 284 else
285 { 285 {
286 QR::type type = (nargout == 0 || nargout == 1) ? QR::raw 286 QR::type type = (nargout == 0 || nargout == 1) ? QR::raw
287 : (nargin == 2 ? QR::economy : QR::std); 287 : (nargin == 2 ? QR::economy : QR::std);
288 288
289 if (arg.is_single_type ()) 289 if (arg.is_single_type ())
290 { 290 {
291 if (arg.is_real_type ()) 291 if (arg.is_real_type ())
292 { 292 {
293 FloatMatrix m = arg.float_matrix_value (); 293 FloatMatrix m = arg.float_matrix_value ();
294 294
295 if (! error_state) 295 if (! error_state)
296 { 296 {
297 switch (nargout) 297 switch (nargout)
298 { 298 {
299 case 0: 299 case 0:
300 case 1: 300 case 1:
301 { 301 {
302 FloatQR fact (m, type); 302 FloatQR fact (m, type);
303 retval(0) = fact.R (); 303 retval(0) = fact.R ();
304 } 304 }
305 break; 305 break;
306 306
307 case 2: 307 case 2:
308 { 308 {
309 FloatQR fact (m, type); 309 FloatQR fact (m, type);
310 retval(1) = get_qr_r (fact); 310 retval(1) = get_qr_r (fact);
311 retval(0) = fact.Q (); 311 retval(0) = fact.Q ();
312 } 312 }
313 break; 313 break;
314 314
315 default: 315 default:
316 { 316 {
317 FloatQRP fact (m, type); 317 FloatQRP fact (m, type);
318 if (type == QR::economy) 318 if (type == QR::economy)
319 retval(2) = fact.Pvec (); 319 retval(2) = fact.Pvec ();
320 else 320 else
321 retval(2) = fact.P (); 321 retval(2) = fact.P ();
322 retval(1) = get_qr_r (fact); 322 retval(1) = get_qr_r (fact);
323 retval(0) = fact.Q (); 323 retval(0) = fact.Q ();
324 } 324 }
325 break; 325 break;
326 } 326 }
327 } 327 }
328 } 328 }
329 else if (arg.is_complex_type ()) 329 else if (arg.is_complex_type ())
330 { 330 {
331 FloatComplexMatrix m = arg.float_complex_matrix_value (); 331 FloatComplexMatrix m = arg.float_complex_matrix_value ();
332 332
333 if (! error_state) 333 if (! error_state)
334 { 334 {
335 switch (nargout) 335 switch (nargout)
336 { 336 {
337 case 0: 337 case 0:
338 case 1: 338 case 1:
339 { 339 {
340 FloatComplexQR fact (m, type); 340 FloatComplexQR fact (m, type);
341 retval(0) = fact.R (); 341 retval(0) = fact.R ();
342 } 342 }
343 break; 343 break;
344 344
345 case 2: 345 case 2:
346 { 346 {
347 FloatComplexQR fact (m, type); 347 FloatComplexQR fact (m, type);
348 retval(1) = get_qr_r (fact); 348 retval(1) = get_qr_r (fact);
349 retval(0) = fact.Q (); 349 retval(0) = fact.Q ();
350 } 350 }
351 break; 351 break;
352 352
353 default: 353 default:
354 { 354 {
355 FloatComplexQRP fact (m, type); 355 FloatComplexQRP fact (m, type);
356 if (type == QR::economy) 356 if (type == QR::economy)
357 retval(2) = fact.Pvec (); 357 retval(2) = fact.Pvec ();
358 else 358 else
359 retval(2) = fact.P (); 359 retval(2) = fact.P ();
360 retval(1) = get_qr_r (fact); 360 retval(1) = get_qr_r (fact);
361 retval(0) = fact.Q (); 361 retval(0) = fact.Q ();
362 } 362 }
363 break; 363 break;
364 } 364 }
365 } 365 }
366 } 366 }
367 } 367 }
368 else 368 else
369 { 369 {
370 if (arg.is_real_type ()) 370 if (arg.is_real_type ())
371 { 371 {
372 Matrix m = arg.matrix_value (); 372 Matrix m = arg.matrix_value ();
373 373
374 if (! error_state) 374 if (! error_state)
375 { 375 {
376 switch (nargout) 376 switch (nargout)
377 { 377 {
378 case 0: 378 case 0:
379 case 1: 379 case 1:
380 { 380 {
381 QR fact (m, type); 381 QR fact (m, type);
382 retval(0) = fact.R (); 382 retval(0) = fact.R ();
383 } 383 }
384 break; 384 break;
385 385
386 case 2: 386 case 2:
387 { 387 {
388 QR fact (m, type); 388 QR fact (m, type);
389 retval(1) = get_qr_r (fact); 389 retval(1) = get_qr_r (fact);
390 retval(0) = fact.Q (); 390 retval(0) = fact.Q ();
391 } 391 }
392 break; 392 break;
393 393
394 default: 394 default:
395 { 395 {
396 QRP fact (m, type); 396 QRP fact (m, type);
397 if (type == QR::economy) 397 if (type == QR::economy)
398 retval(2) = fact.Pvec (); 398 retval(2) = fact.Pvec ();
399 else 399 else
400 retval(2) = fact.P (); 400 retval(2) = fact.P ();
401 retval(1) = get_qr_r (fact); 401 retval(1) = get_qr_r (fact);
402 retval(0) = fact.Q (); 402 retval(0) = fact.Q ();
403 } 403 }
404 break; 404 break;
405 } 405 }
406 } 406 }
407 } 407 }
408 else if (arg.is_complex_type ()) 408 else if (arg.is_complex_type ())
409 { 409 {
410 ComplexMatrix m = arg.complex_matrix_value (); 410 ComplexMatrix m = arg.complex_matrix_value ();
411 411
412 if (! error_state) 412 if (! error_state)
413 { 413 {
414 switch (nargout) 414 switch (nargout)
415 { 415 {
416 case 0: 416 case 0:
417 case 1: 417 case 1:
418 { 418 {
419 ComplexQR fact (m, type); 419 ComplexQR fact (m, type);
420 retval(0) = fact.R (); 420 retval(0) = fact.R ();
421 } 421 }
422 break; 422 break;
423 423
424 case 2: 424 case 2:
425 { 425 {
426 ComplexQR fact (m, type); 426 ComplexQR fact (m, type);
427 retval(1) = get_qr_r (fact); 427 retval(1) = get_qr_r (fact);
428 retval(0) = fact.Q (); 428 retval(0) = fact.Q ();
429 } 429 }
430 break; 430 break;
431 431
432 default: 432 default:
433 { 433 {
434 ComplexQRP fact (m, type); 434 ComplexQRP fact (m, type);
435 if (type == QR::economy) 435 if (type == QR::economy)
436 retval(2) = fact.Pvec (); 436 retval(2) = fact.Pvec ();
437 else 437 else
438 retval(2) = fact.P (); 438 retval(2) = fact.P ();
439 retval(1) = get_qr_r (fact); 439 retval(1) = get_qr_r (fact);
440 retval(0) = fact.Q (); 440 retval(0) = fact.Q ();
441 } 441 }
442 break; 442 break;
443 } 443 }
444 } 444 }
445 } 445 }
446 else 446 else
447 gripe_wrong_type_arg ("qr", arg); 447 gripe_wrong_type_arg ("qr", arg);
448 } 448 }
449 } 449 }
450 450
451 return retval; 451 return retval;
452 } 452 }
453 453
800 && argu.is_numeric_type () && argv.is_numeric_type ()) 800 && argu.is_numeric_type () && argv.is_numeric_type ())
801 { 801 {
802 if (check_qr_dims (argq, argr, true)) 802 if (check_qr_dims (argq, argr, true))
803 { 803 {
804 if (argq.is_real_type () 804 if (argq.is_real_type ()
805 && argr.is_real_type () 805 && argr.is_real_type ()
806 && argu.is_real_type () 806 && argu.is_real_type ()
807 && argv.is_real_type ()) 807 && argv.is_real_type ())
808 { 808 {
809 // all real case 809 // all real case
810 if (argq.is_single_type () 810 if (argq.is_single_type ()
811 || argr.is_single_type () 811 || argr.is_single_type ()
812 || argu.is_single_type () 812 || argu.is_single_type ()
813 || argv.is_single_type ()) 813 || argv.is_single_type ())
814 { 814 {
815 FloatMatrix Q = argq.float_matrix_value (); 815 FloatMatrix Q = argq.float_matrix_value ();
816 FloatMatrix R = argr.float_matrix_value (); 816 FloatMatrix R = argr.float_matrix_value ();
817 FloatMatrix u = argu.float_matrix_value (); 817 FloatMatrix u = argu.float_matrix_value ();
818 FloatMatrix v = argv.float_matrix_value (); 818 FloatMatrix v = argv.float_matrix_value ();
819 819
820 FloatQR fact (Q, R); 820 FloatQR fact (Q, R);
821 fact.update (u, v); 821 fact.update (u, v);
822 822
823 retval(1) = get_qr_r (fact); 823 retval(1) = get_qr_r (fact);
824 retval(0) = fact.Q (); 824 retval(0) = fact.Q ();
825 } 825 }
826 else 826 else
827 { 827 {
828 Matrix Q = argq.matrix_value (); 828 Matrix Q = argq.matrix_value ();
829 Matrix R = argr.matrix_value (); 829 Matrix R = argr.matrix_value ();
830 Matrix u = argu.matrix_value (); 830 Matrix u = argu.matrix_value ();
831 Matrix v = argv.matrix_value (); 831 Matrix v = argv.matrix_value ();
832 832
833 QR fact (Q, R); 833 QR fact (Q, R);
834 fact.update (u, v); 834 fact.update (u, v);
835 835
836 retval(1) = get_qr_r (fact); 836 retval(1) = get_qr_r (fact);
837 retval(0) = fact.Q (); 837 retval(0) = fact.Q ();
838 } 838 }
839 } 839 }
840 else 840 else
841 { 841 {
842 // complex case 842 // complex case
843 if (argq.is_single_type () 843 if (argq.is_single_type ()
844 || argr.is_single_type () 844 || argr.is_single_type ()
845 || argu.is_single_type () 845 || argu.is_single_type ()
846 || argv.is_single_type ()) 846 || argv.is_single_type ())
847 { 847 {
848 FloatComplexMatrix Q = argq.float_complex_matrix_value (); 848 FloatComplexMatrix Q = argq.float_complex_matrix_value ();
849 FloatComplexMatrix R = argr.float_complex_matrix_value (); 849 FloatComplexMatrix R = argr.float_complex_matrix_value ();
850 FloatComplexMatrix u = argu.float_complex_matrix_value (); 850 FloatComplexMatrix u = argu.float_complex_matrix_value ();
851 FloatComplexMatrix v = argv.float_complex_matrix_value (); 851 FloatComplexMatrix v = argv.float_complex_matrix_value ();
852 852
853 FloatComplexQR fact (Q, R); 853 FloatComplexQR fact (Q, R);
854 fact.update (u, v); 854 fact.update (u, v);
855 855
856 retval(1) = get_qr_r (fact); 856 retval(1) = get_qr_r (fact);
857 retval(0) = fact.Q (); 857 retval(0) = fact.Q ();
858 } 858 }
859 else 859 else
860 { 860 {
861 ComplexMatrix Q = argq.complex_matrix_value (); 861 ComplexMatrix Q = argq.complex_matrix_value ();
862 ComplexMatrix R = argr.complex_matrix_value (); 862 ComplexMatrix R = argr.complex_matrix_value ();
863 ComplexMatrix u = argu.complex_matrix_value (); 863 ComplexMatrix u = argu.complex_matrix_value ();
864 ComplexMatrix v = argv.complex_matrix_value (); 864 ComplexMatrix v = argv.complex_matrix_value ();
865 865
866 ComplexQR fact (Q, R); 866 ComplexQR fact (Q, R);
867 fact.update (u, v); 867 fact.update (u, v);
868 868
869 retval(1) = get_qr_r (fact); 869 retval(1) = get_qr_r (fact);
870 retval(0) = fact.Q (); 870 retval(0) = fact.Q ();
871 } 871 }
872 } 872 }
873 } 873 }
874 else 874 else
875 error ("qrupdate: dimensions mismatch"); 875 error ("qrupdate: dimensions mismatch");
876 } 876 }
877 else 877 else
878 error ("qrupdate: expecting numeric arguments"); 878 error ("qrupdate: expecting numeric arguments");
879 879
880 return retval; 880 return retval;
999 && (col || argx.rows () == 1)) 999 && (col || argx.rows () == 1))
1000 { 1000 {
1001 if (check_index (argj, col)) 1001 if (check_index (argj, col))
1002 { 1002 {
1003 MArray<octave_idx_type> j 1003 MArray<octave_idx_type> j
1004 = argj.octave_idx_type_vector_value (); 1004 = argj.octave_idx_type_vector_value ();
1005 1005
1006 if (argq.is_real_type () 1006 if (argq.is_real_type ()
1007 && argr.is_real_type () 1007 && argr.is_real_type ()
1008 && argx.is_real_type ()) 1008 && argx.is_real_type ())
1009 { 1009 {
1010 // real case 1010 // real case
1011 if (argq.is_single_type () 1011 if (argq.is_single_type ()
1012 || argr.is_single_type () 1012 || argr.is_single_type ()
1013 || argx.is_single_type ()) 1013 || argx.is_single_type ())
1014 { 1014 {
1015 FloatMatrix Q = argq.float_matrix_value (); 1015 FloatMatrix Q = argq.float_matrix_value ();
1016 FloatMatrix R = argr.float_matrix_value (); 1016 FloatMatrix R = argr.float_matrix_value ();
1017 FloatMatrix x = argx.float_matrix_value (); 1017 FloatMatrix x = argx.float_matrix_value ();
1018 1018
1019 FloatQR fact (Q, R); 1019 FloatQR fact (Q, R);
1020 1020
1021 if (col) 1021 if (col)
1022 fact.insert_col (x, j-1); 1022 fact.insert_col (x, j-1);
1023 else 1023 else
1024 fact.insert_row (x.row (0), j(0)-1); 1024 fact.insert_row (x.row (0), j(0)-1);
1025 1025
1026 retval(1) = get_qr_r (fact); 1026 retval(1) = get_qr_r (fact);
1027 retval(0) = fact.Q (); 1027 retval(0) = fact.Q ();
1028 1028
1029 } 1029 }
1030 else 1030 else
1031 { 1031 {
1032 Matrix Q = argq.matrix_value (); 1032 Matrix Q = argq.matrix_value ();
1033 Matrix R = argr.matrix_value (); 1033 Matrix R = argr.matrix_value ();
1034 Matrix x = argx.matrix_value (); 1034 Matrix x = argx.matrix_value ();
1035 1035
1036 QR fact (Q, R); 1036 QR fact (Q, R);
1037 1037
1038 if (col) 1038 if (col)
1039 fact.insert_col (x, j-1); 1039 fact.insert_col (x, j-1);
1040 else 1040 else
1041 fact.insert_row (x.row (0), j(0)-1); 1041 fact.insert_row (x.row (0), j(0)-1);
1042 1042
1043 retval(1) = get_qr_r (fact); 1043 retval(1) = get_qr_r (fact);
1044 retval(0) = fact.Q (); 1044 retval(0) = fact.Q ();
1045 1045
1046 } 1046 }
1047 } 1047 }
1048 else 1048 else
1049 { 1049 {
1050 // complex case 1050 // complex case
1051 if (argq.is_single_type () 1051 if (argq.is_single_type ()
1052 || argr.is_single_type () 1052 || argr.is_single_type ()
1053 || argx.is_single_type ()) 1053 || argx.is_single_type ())
1054 { 1054 {
1055 FloatComplexMatrix Q = argq.float_complex_matrix_value (); 1055 FloatComplexMatrix Q = argq.float_complex_matrix_value ();
1056 FloatComplexMatrix R = argr.float_complex_matrix_value (); 1056 FloatComplexMatrix R = argr.float_complex_matrix_value ();
1057 FloatComplexMatrix x = argx.float_complex_matrix_value (); 1057 FloatComplexMatrix x = argx.float_complex_matrix_value ();
1058 1058
1059 FloatComplexQR fact (Q, R); 1059 FloatComplexQR fact (Q, R);
1060 1060
1061 if (col) 1061 if (col)
1062 fact.insert_col (x, j-1); 1062 fact.insert_col (x, j-1);
1063 else 1063 else
1064 fact.insert_row (x.row (0), j(0)-1); 1064 fact.insert_row (x.row (0), j(0)-1);
1065 1065
1066 retval(1) = get_qr_r (fact); 1066 retval(1) = get_qr_r (fact);
1067 retval(0) = fact.Q (); 1067 retval(0) = fact.Q ();
1068 } 1068 }
1069 else 1069 else
1070 { 1070 {
1071 ComplexMatrix Q = argq.complex_matrix_value (); 1071 ComplexMatrix Q = argq.complex_matrix_value ();
1072 ComplexMatrix R = argr.complex_matrix_value (); 1072 ComplexMatrix R = argr.complex_matrix_value ();
1073 ComplexMatrix x = argx.complex_matrix_value (); 1073 ComplexMatrix x = argx.complex_matrix_value ();
1074 1074
1075 ComplexQR fact (Q, R); 1075 ComplexQR fact (Q, R);
1076 1076
1077 if (col) 1077 if (col)
1078 fact.insert_col (x, j-1); 1078 fact.insert_col (x, j-1);
1079 else 1079 else
1080 fact.insert_row (x.row (0), j(0)-1); 1080 fact.insert_row (x.row (0), j(0)-1);
1081 1081
1082 retval(1) = get_qr_r (fact); 1082 retval(1) = get_qr_r (fact);
1083 retval(0) = fact.Q (); 1083 retval(0) = fact.Q ();
1084 } 1084 }
1085 } 1085 }
1086 1086
1087 } 1087 }
1088 else 1088 else
1089 error ("qrinsert: invalid index"); 1089 error ("qrinsert: invalid index");
1213 if (check_qr_dims (argq, argr, col)) 1213 if (check_qr_dims (argq, argr, col))
1214 { 1214 {
1215 if (check_index (argj, col)) 1215 if (check_index (argj, col))
1216 { 1216 {
1217 MArray<octave_idx_type> j 1217 MArray<octave_idx_type> j
1218 = argj.octave_idx_type_vector_value (); 1218 = argj.octave_idx_type_vector_value ();
1219 1219
1220 if (argq.is_real_type () 1220 if (argq.is_real_type ()
1221 && argr.is_real_type ()) 1221 && argr.is_real_type ())
1222 { 1222 {
1223 // real case 1223 // real case
1224 if (argq.is_single_type () 1224 if (argq.is_single_type ()
1225 || argr.is_single_type ()) 1225 || argr.is_single_type ())
1226 { 1226 {
1227 FloatMatrix Q = argq.float_matrix_value (); 1227 FloatMatrix Q = argq.float_matrix_value ();
1228 FloatMatrix R = argr.float_matrix_value (); 1228 FloatMatrix R = argr.float_matrix_value ();
1229 1229
1230 FloatQR fact (Q, R); 1230 FloatQR fact (Q, R);
1231 1231
1232 if (col) 1232 if (col)
1233 fact.delete_col (j-1); 1233 fact.delete_col (j-1);
1234 else 1234 else
1235 fact.delete_row (j(0)-1); 1235 fact.delete_row (j(0)-1);
1236 1236
1237 retval(1) = get_qr_r (fact); 1237 retval(1) = get_qr_r (fact);
1238 retval(0) = fact.Q (); 1238 retval(0) = fact.Q ();
1239 } 1239 }
1240 else 1240 else
1241 { 1241 {
1242 Matrix Q = argq.matrix_value (); 1242 Matrix Q = argq.matrix_value ();
1243 Matrix R = argr.matrix_value (); 1243 Matrix R = argr.matrix_value ();
1244 1244
1245 QR fact (Q, R); 1245 QR fact (Q, R);
1246 1246
1247 if (col) 1247 if (col)
1248 fact.delete_col (j-1); 1248 fact.delete_col (j-1);
1249 else 1249 else
1250 fact.delete_row (j(0)-1); 1250 fact.delete_row (j(0)-1);
1251 1251
1252 retval(1) = get_qr_r (fact); 1252 retval(1) = get_qr_r (fact);
1253 retval(0) = fact.Q (); 1253 retval(0) = fact.Q ();
1254 } 1254 }
1255 } 1255 }
1256 else 1256 else
1257 { 1257 {
1258 // complex case 1258 // complex case
1259 if (argq.is_single_type () 1259 if (argq.is_single_type ()
1260 || argr.is_single_type ()) 1260 || argr.is_single_type ())
1261 { 1261 {
1262 FloatComplexMatrix Q = argq.float_complex_matrix_value (); 1262 FloatComplexMatrix Q = argq.float_complex_matrix_value ();
1263 FloatComplexMatrix R = argr.float_complex_matrix_value (); 1263 FloatComplexMatrix R = argr.float_complex_matrix_value ();
1264 1264
1265 FloatComplexQR fact (Q, R); 1265 FloatComplexQR fact (Q, R);
1266 1266
1267 if (col) 1267 if (col)
1268 fact.delete_col (j-1); 1268 fact.delete_col (j-1);
1269 else 1269 else
1270 fact.delete_row (j(0)-1); 1270 fact.delete_row (j(0)-1);
1271 1271
1272 retval(1) = get_qr_r (fact); 1272 retval(1) = get_qr_r (fact);
1273 retval(0) = fact.Q (); 1273 retval(0) = fact.Q ();
1274 } 1274 }
1275 else 1275 else
1276 { 1276 {
1277 ComplexMatrix Q = argq.complex_matrix_value (); 1277 ComplexMatrix Q = argq.complex_matrix_value ();
1278 ComplexMatrix R = argr.complex_matrix_value (); 1278 ComplexMatrix R = argr.complex_matrix_value ();
1279 1279
1280 ComplexQR fact (Q, R); 1280 ComplexQR fact (Q, R);
1281 1281
1282 if (col) 1282 if (col)
1283 fact.delete_col (j-1); 1283 fact.delete_col (j-1);
1284 else 1284 else
1285 fact.delete_row (j(0)-1); 1285 fact.delete_row (j(0)-1);
1286 1286
1287 retval(1) = get_qr_r (fact); 1287 retval(1) = get_qr_r (fact);
1288 retval(0) = fact.Q (); 1288 retval(0) = fact.Q ();
1289 } 1289 }
1290 } 1290 }
1291 } 1291 }
1292 else 1292 else
1293 error ("qrdelete: invalid index"); 1293 error ("qrdelete: invalid index");
1294 } 1294 }
1449 1449
1450 if (argq.is_real_type () 1450 if (argq.is_real_type ()
1451 && argr.is_real_type ()) 1451 && argr.is_real_type ())
1452 { 1452 {
1453 // all real case 1453 // all real case
1454 if (argq.is_single_type () 1454 if (argq.is_single_type ()
1455 && argr.is_single_type ()) 1455 && argr.is_single_type ())
1456 { 1456 {
1457 FloatMatrix Q = argq.float_matrix_value (); 1457 FloatMatrix Q = argq.float_matrix_value ();
1458 FloatMatrix R = argr.float_matrix_value (); 1458 FloatMatrix R = argr.float_matrix_value ();
1459 1459
1460 FloatQR fact (Q, R); 1460 FloatQR fact (Q, R);
1461 fact.shift_cols (i-1, j-1); 1461 fact.shift_cols (i-1, j-1);
1462 1462
1463 retval(1) = get_qr_r (fact); 1463 retval(1) = get_qr_r (fact);
1464 retval(0) = fact.Q (); 1464 retval(0) = fact.Q ();
1465 } 1465 }
1466 else 1466 else
1467 { 1467 {
1468 Matrix Q = argq.matrix_value (); 1468 Matrix Q = argq.matrix_value ();
1469 Matrix R = argr.matrix_value (); 1469 Matrix R = argr.matrix_value ();
1470 1470
1471 QR fact (Q, R); 1471 QR fact (Q, R);
1472 fact.shift_cols (i-1, j-1); 1472 fact.shift_cols (i-1, j-1);
1473 1473
1474 retval(1) = get_qr_r (fact); 1474 retval(1) = get_qr_r (fact);
1475 retval(0) = fact.Q (); 1475 retval(0) = fact.Q ();
1476 } 1476 }
1477 } 1477 }
1478 else 1478 else
1479 { 1479 {
1480 // complex case 1480 // complex case
1481 if (argq.is_single_type () 1481 if (argq.is_single_type ()
1482 && argr.is_single_type ()) 1482 && argr.is_single_type ())
1483 { 1483 {
1484 FloatComplexMatrix Q = argq.float_complex_matrix_value (); 1484 FloatComplexMatrix Q = argq.float_complex_matrix_value ();
1485 FloatComplexMatrix R = argr.float_complex_matrix_value (); 1485 FloatComplexMatrix R = argr.float_complex_matrix_value ();
1486 1486
1487 FloatComplexQR fact (Q, R); 1487 FloatComplexQR fact (Q, R);
1488 fact.shift_cols (i-1, j-1); 1488 fact.shift_cols (i-1, j-1);
1489 1489
1490 retval(1) = get_qr_r (fact); 1490 retval(1) = get_qr_r (fact);
1491 retval(0) = fact.Q (); 1491 retval(0) = fact.Q ();
1492 } 1492 }
1493 else 1493 else
1494 { 1494 {
1495 ComplexMatrix Q = argq.complex_matrix_value (); 1495 ComplexMatrix Q = argq.complex_matrix_value ();
1496 ComplexMatrix R = argr.complex_matrix_value (); 1496 ComplexMatrix R = argr.complex_matrix_value ();
1497 1497
1498 ComplexQR fact (Q, R); 1498 ComplexQR fact (Q, R);
1499 fact.shift_cols (i-1, j-1); 1499 fact.shift_cols (i-1, j-1);
1500 1500
1501 retval(1) = get_qr_r (fact); 1501 retval(1) = get_qr_r (fact);
1502 retval(0) = fact.Q (); 1502 retval(0) = fact.Q ();
1503 } 1503 }
1504 } 1504 }
1505 } 1505 }
1506 else 1506 else
1507 error ("qrshift: invalid index"); 1507 error ("qrshift: invalid index");
1508 } 1508 }
1509 else 1509 else
1510 error ("qrshift: dimensions mismatch"); 1510 error ("qrshift: dimensions mismatch");
1511 } 1511 }
1512 else 1512 else
1513 error ("qrshift: expecting numeric arguments"); 1513 error ("qrshift: expecting numeric arguments");
1514 1514
1515 return retval; 1515 return retval;