40 void poly(
const vec &r, vec &p)
44 p.set_size(n + 1,
false);
48 for (
int i = 0; i < n; i++)
49 p.set_subvector(1, p(1, i + 1) - r(i)*p(0, i));
52 void poly(
const cvec &r, cvec &p)
56 p.set_size(n + 1,
false);
60 for (
int i = 0; i < n; i++)
61 p.set_subvector(1, p(1, i + 1) - r(i)*p(0, i));
66 void roots(
const vec &p, cvec &r)
68 int n = p.size(), m, l;
69 ivec f =
find(p != 0.0);
75 v = v(f(0), f(m - 1));
81 A.set_row(0, -v(1, l - 1) / v(0));
91 r.set_size(n - f(m - 1) - 1,
false);
99 void roots(
const cvec &p, cvec &r)
101 int n = p.size(), m, l;
105 for (
int i = 0; i < n; i++)
114 if (m > 0 && n > 1) {
115 v = v(f(0), f(m - 1));
120 A.set_row(0, -v(1, l - 1) / v(0));
126 r.set_size(n - f(m - 1) - 1,
false);
131 r.set_size(0,
false);
137 it_error_if(p.size() == 0,
"polyval: size of polynomial is zero");
138 it_error_if(x.size() == 0,
"polyval: size of input value vector is zero");
144 for (
int i = 1; i < p.size(); i++)
152 it_error_if(p.size() == 0,
"polyval: size of polynomial is zero");
153 it_error_if(x.size() == 0,
"polyval: size of input value vector is zero");
159 for (
int i = 1; i < p.size(); i++)
160 out = std::complex<double>(p(i)) +
elem_mult(x, out);
167 it_error_if(p.size() == 0,
"polyval: size of polynomial is zero");
168 it_error_if(x.size() == 0,
"polyval: size of input value vector is zero");
174 for (
int i = 1; i < p.size(); i++)
182 it_error_if(p.size() == 0,
"polyval: size of polynomial is zero");
183 it_error_if(x.size() == 0,
"polyval: size of input value vector is zero");
189 for (
int i = 1; i < p.size(); i++)
197 it_assert((n >= 0),
"cheb(): need a non-negative order n!");
199 if (x < 1.0 && x > -1.0) {
214 for (
int i = 0; i < x.size(); ++i) {
215 out(i) =
cheb(n, x(i));
222 it_assert_debug((x.rows() > 0) && (x.cols() > 0),
"cheb(): empty matrix");
224 mat out(x.rows(), x.cols());
225 for (
int i = 0; i < x.rows(); ++i) {
226 for (
int j = 0; j < x.cols(); ++j) {
227 out(i, j) =
cheb(n, x(i, j));
#define it_error_if(t, s)
Abort if t is true.
Various functions on vectors and matrices - header file.
vec cosh(const vec &x)
Cosine hyperbolic function.
vec acos(const vec &x)
Inverse cosine function.
Definitions of eigenvalue decomposition functions.
bool is_even(int x)
Return true if x is an even integer.
ivec find(const bvec &invector)
Return a integer vector with indicies where bvec == 1.
Mat< Num_T > elem_mult(const Mat< Num_T > &m1, const Mat< Num_T > &m2)
Element wise multiplication of two matrices.
#define it_assert(t, s)
Abort if t is not true.
#define it_assert_debug(t, s)
Abort if t is not true and NDEBUG is not defined.
double cheb(int n, double x)
Chebyshev polynomial of the first kindChebyshev polynomials of the first kind can be defined as follo...
ITPP_EXPORT cvec zeros_c(int size)
A Double Complex vector of zeros.
Definitions of converters between different vector and matrix types.
vec acosh(const vec &x)
Inverse cosine hyperbolic function.
vec polyval(const vec &p, const vec &x)
Evaluate polynomialEvaluate the polynomial p (of length at the points x The output is given by ...
ITPP_EXPORT vec ones(int size)
A float vector of ones.
IT++ compatibility types and functions.
bool eig(const mat &A, cvec &d, cmat &V)
Calculates the eigenvalues and eigenvectors of a real non-symmetric matrix.
Definitions of special vectors and matrices.
void roots(const vec &p, cvec &r)
Calculate the roots of the polynomialCalculate the roots r of the polynomial p.
cvec to_cvec(const Vec< T > &v)
Converts a Vec<T> to cvec.
Mat< T > diag(const Vec< T > &v, const int K=0)
Create a diagonal matrix using vector v as its diagonal.
vec cos(const vec &x)
Cosine function.
ITPP_EXPORT cvec ones_c(int size)
A float Complex vector of ones.
void poly(const vec &r, vec &p)
Create a polynomial of the given rootsCreate a polynomial p with roots r.
const Array< T > concat(const Array< T > &a, const T &e)
Append element e to the end of the Array a.