comparison src/xpow.cc @ 9103:10bed8fbec99

optimize scalar .^ range operation
author Jaroslav Hajek <highegg@gmail.com>
date Wed, 08 Apr 2009 10:53:06 +0200
parents a7c00773a089
children f39b98237d5c
comparison
equal deleted inserted replaced
9102:9dc516d36175 9103:10bed8fbec99
39 #include "fDiagMatrix.h" 39 #include "fDiagMatrix.h"
40 #include "dMatrix.h" 40 #include "dMatrix.h"
41 #include "PermMatrix.h" 41 #include "PermMatrix.h"
42 #include "mx-cm-cdm.h" 42 #include "mx-cm-cdm.h"
43 #include "oct-cmplx.h" 43 #include "oct-cmplx.h"
44 #include "Range.h"
44 #include "quit.h" 45 #include "quit.h"
45 46
46 #include "error.h" 47 #include "error.h"
47 #include "oct-obj.h" 48 #include "oct-obj.h"
48 #include "utils.h" 49 #include "utils.h"
689 } 690 }
690 691
691 return result; 692 return result;
692 } 693 }
693 694
695 octave_value
696 elem_xpow (double a, const Range& r)
697 {
698 octave_value retval;
699
700 if (r.nelem () <= 0)
701 retval = Matrix ();
702 else if (a < 0.0 && ! r.all_elements_are_ints ())
703 {
704 ComplexMatrix mat (1, r.nelem ());
705 Complex atmp (a);
706 Complex *pmat = mat.fortran_vec ();
707
708 pmat[0] = std::pow (atmp, r.base ());
709 Complex mul = std::pow (atmp, r.inc ());
710 for (octave_idx_type i = 1; i < r.nelem (); i++)
711 pmat[i] = pmat[i-1] * mul;
712
713 retval = mat;
714 }
715 else
716 {
717 Matrix mat (1, r.nelem ());
718 double *pmat = mat.fortran_vec ();
719
720 double base = std::pow (a, r.base ());
721 pmat[0] = base;
722 double mul = std::pow (a, r.inc ());
723 for (octave_idx_type i = 1; i < r.nelem (); i++)
724 pmat[i] = (base *= mul);
725
726 retval = mat;
727 }
728
729 return retval;
730 }
731
694 // -*- 3 -*- 732 // -*- 3 -*-
695 octave_value 733 octave_value
696 elem_xpow (const Matrix& a, double b) 734 elem_xpow (const Matrix& a, double b)
697 { 735 {
698 octave_value retval; 736 octave_value retval;
884 OCTAVE_QUIT; 922 OCTAVE_QUIT;
885 result (i, j) = std::pow (a, b (i, j)); 923 result (i, j) = std::pow (a, b (i, j));
886 } 924 }
887 925
888 return result; 926 return result;
927 }
928
929 octave_value
930 elem_xpow (const Complex& a, const Range& r)
931 {
932 octave_value retval;
933
934 if (r.nelem () <= 0)
935 retval = Matrix ();
936 else
937 {
938 ComplexMatrix mat (1, r.nelem ());
939 Complex *pmat = mat.fortran_vec ();
940
941 pmat[0] = std::pow (a, r.base ());
942 Complex mul = std::pow (a, r.inc ());
943 for (octave_idx_type i = 1; i < r.nelem (); i++)
944 pmat[i] = pmat[i-1] * mul;
945
946 retval = mat;
947 }
948
949 return retval;
889 } 950 }
890 951
891 // -*- 9 -*- 952 // -*- 9 -*-
892 octave_value 953 octave_value
893 elem_xpow (const ComplexMatrix& a, double b) 954 elem_xpow (const ComplexMatrix& a, double b)