added function ones and polyroots.

Changes to be committed:
modified:   dspl/doc/html/formula.repository
modified:   dspl/dox/ru/error_list.dox
modified:   dspl/src/array.c
modified:   dspl/src/polyval.c
new file:   examples/src/polyroots_test.c
modified:   include/dspl.c
modified:   include/dspl.h
pull/6/merge
Dsplib 2020-05-18 20:44:38 +03:00
rodzic 74b1b0cece
commit 66ab45251c
7 zmienionych plików z 222 dodań i 102 usunięć

Wyświetl plik

@ -63,83 +63,83 @@
\form#62:\[ k_i = \left( \frac{k_{i-1}} { 1+\sqrt{1-k_{i-1}^2} } \right)^2 \]
\form#63:$ k<1 $
\form#64:$ y = \textrm{sn}(u K(k), k)$
\form#65:\[ Y(k) = \frac{1}{N} \sum_{m = 0}^{n-1} x(m) \exp \left( j \frac{2\pi}{n} m k \right), \]
\form#66:$n = n_0 \times n_1 \times n_2 \times n_3 \times \ldots \times n_p \times m$
\form#67:$n = n_0 \times n_1 \times n_2 \ldots \times n_p \times m$
\form#68:$ n = 725760 $
\form#69:$725760 = 2 \cdot 3 \cdot 4 \cdot 5 \cdot 6 \cdot 7 \cdot 9 \cdot 16 $
\form#70:$ n = 172804 = 43201 \cdot 4 $
\form#71:$ n = 13 \cdot 17 \cdot 23 \cdot 13 = 66079 $
\form#72:$\sqrt{2^{31}} = 46340.95$
\form#73:$x_0$
\form#74:$x_1$
\form#75:$x(k) = x_0 + k \cdot dx$
\form#76:$dx = \frac{x_1 - x_0}{n-1}$
\form#77:$k = 0 \ldots n-1.$
\form#78:$dx = \frac{x_1 - x_0}{n}$
\form#79:$10^{x_0}$
\form#80:$10^{x_1}$
\form#81:$x(k) = 10^{x_0} \cdot dx^k$
\form#82:$dx = \sqrt[n-1]{10^{x_1 - x_0}}$
\form#83:$dx = \sqrt[n]{10^{x_1 - x_0}}$
\form#84:$ H(j \omega) $
\form#85:$ H(j \omega)$
\form#86:$ H(s) $
\form#87:\[ H(s) = \frac {\sum_{k = 0}^{N} b_k s^k} {\sum_{m = 0}^{N} a_m s^m}, \]
\form#88:$ N $
\form#89:$ s = j \omega $
\form#90:$ \omega $
\form#91:$H(s)$
\form#92:$ H \left(\mathrm{e}^{j\omega} \right) $
\form#93:$ 2\pi $
\form#94:$ \pi $
\form#95:$ -\pi $
\form#96:$ H \left(e^{j \omega} \right)$
\form#97:\[ H(z) = \frac {\sum_{k = 0}^{N} b_k z^{-k}} {\sum_{m = 0}^{N} a_m z^{-m}}, \]
\form#98:$N$
\form#99:$z = e^{j \omega} $
\form#100:$\omega$
\form#101:$ 2 \pi-$
\form#102:$2 \pi$
\form#103:$-\pi$
\form#104:$ \pi$
\form#105:$ H \left(e^{j \omega} \right) = H^* \left(e^{-j \omega} \right)$
\form#106:$\pi$
\form#107:$ -R_p $
\form#108:$ H(s)$
\form#109:$-R_p$
\form#110:$ R_p $
\form#111:$-R_s$
\form#112:$H(j\cdot 1) = -R_s$
\form#113:\[ H(s) = \frac{\sum_{n = 0}^{N_z} b_n \cdot s^n}{\sum_{m = 0}^{N_p} a_m \cdot s^m} = \frac{\prod_{n = 0}^{N_z}(s-z_n)}{\prod_{m = 0}^{N_p} (s-p_m)} \]
\form#114:\[ H(z) = \sum_{n = 0}^{ord} h_n z^{-n} \]
\form#115:$ F(s) $
\form#116:$F(s)$
\form#117:$Y(s) = (H \circ F)(s) = H(F(s))$
\form#118:\[ H(s) = \frac{\sum\limits_{m = 0}^{n} b_m s^m} {\sum\limits_{k = 0}^{n} a_k s^k}, \quad F(s) = \frac{\sum\limits_{m = 0}^{p} d_m s^m} {\sum\limits_{k = 0}^{p} c_k s^k}, \quad Y(s) = \frac{\sum\limits_{m = 0}^{n p} \beta_m s^m} {\sum\limits_{k = 0}^{n p} \alpha_k s^k} \]
\form#119:$Y(s) = (H \circ F)(s)$
\form#120:\[ s \leftarrow \frac{1 - z^{-1}}{1 - z^{-1}}. \]
\form#121:$\Omega$
\form#122:\[ \Omega = \tan(\omega / 2). \]
\form#123:\[ s(t) = \sum\limits_{n = 0}^{n_{\omega}-1} S(\omega_n) \exp(j\omega_n t) \]
\form#124:$\omega_n$
\form#125:$S(\omega_n)$
\form#126:$ I_0(x)$
\form#127:$ x $
\form#128:$[0 \ 3]$
\form#129:$ \textrm{sinc}(x,a) = \frac{\sin(ax)}{ax}$
\form#130:\[ \textrm{Si}(x) = \int_{0}^{x} \frac{\sin(x)}{x} \, dx\]
\form#131:$[-6\pi \ 6\pi]$
\form#132:$P_N(x)$
\form#133:$N-$
\form#134:\[ P_N(x) = a_0 + a_1 \cdot x + a_2 \cdot x^2 + a_3 \cdot x^3 + ... a_N \cdot x^N. \]
\form#135:\[ P_N(x) = a_0 + x \cdot (a_1 + x \cdot (a_2 + \cdot ( \ldots x \cdot (a_{N-1} + x\cdot a_N) \ldots ))) \]
\form#136:$\mu$
\form#137:$\sigma$
\form#138:$\sigma^2$
\form#139:$\mu = 0$
\form#140:$\sigma=1$
\form#141:$\frac{P}{Q}$
\form#142:$P$
\form#143:$Q$
\form#144:$1/F_{\textrm{s}}$
\form#65:$n = n_0 \times n_1 \times n_2 \ldots \times n_p \times m$
\form#66:$\frac{P}{Q}$
\form#67:$P$
\form#68:$Q$
\form#69:$1/F_{\textrm{s}}$
\form#70:\[ Y(k) = \frac{1}{N} \sum_{m = 0}^{n-1} x(m) \exp \left( j \frac{2\pi}{n} m k \right), \]
\form#71:$n = n_0 \times n_1 \times n_2 \times n_3 \times \ldots \times n_p \times m$
\form#72:$ n = 725760 $
\form#73:$725760 = 2 \cdot 3 \cdot 4 \cdot 5 \cdot 6 \cdot 7 \cdot 9 \cdot 16 $
\form#74:$ n = 172804 = 43201 \cdot 4 $
\form#75:$ n = 13 \cdot 17 \cdot 23 \cdot 13 = 66079 $
\form#76:$\sqrt{2^{31}} = 46340.95$
\form#77:$x_0$
\form#78:$x_1$
\form#79:$x(k) = x_0 + k \cdot dx$
\form#80:$dx = \frac{x_1 - x_0}{n-1}$
\form#81:$k = 0 \ldots n-1.$
\form#82:$dx = \frac{x_1 - x_0}{n}$
\form#83:$10^{x_0}$
\form#84:$10^{x_1}$
\form#85:$x(k) = 10^{x_0} \cdot dx^k$
\form#86:$dx = \sqrt[n-1]{10^{x_1 - x_0}}$
\form#87:$dx = \sqrt[n]{10^{x_1 - x_0}}$
\form#88:$ H(j \omega) $
\form#89:$ H(j \omega)$
\form#90:$ H(s) $
\form#91:\[ H(s) = \frac {\sum_{k = 0}^{N} b_k s^k} {\sum_{m = 0}^{N} a_m s^m}, \]
\form#92:$ N $
\form#93:$ s = j \omega $
\form#94:$ \omega $
\form#95:$H(s)$
\form#96:$ H \left(\mathrm{e}^{j\omega} \right) $
\form#97:$ 2\pi $
\form#98:$ \pi $
\form#99:$ -\pi $
\form#100:$ H \left(e^{j \omega} \right)$
\form#101:\[ H(z) = \frac {\sum_{k = 0}^{N} b_k z^{-k}} {\sum_{m = 0}^{N} a_m z^{-m}}, \]
\form#102:$N$
\form#103:$z = e^{j \omega} $
\form#104:$\omega$
\form#105:$ 2 \pi-$
\form#106:$2 \pi$
\form#107:$-\pi$
\form#108:$ \pi$
\form#109:$ H \left(e^{j \omega} \right) = H^* \left(e^{-j \omega} \right)$
\form#110:$\pi$
\form#111:$ -R_p $
\form#112:$ H(s)$
\form#113:$-R_p$
\form#114:$ R_p $
\form#115:$-R_s$
\form#116:$H(j\cdot 1) = -R_s$
\form#117:\[ H(s) = \frac{\sum_{n = 0}^{N_z} b_n \cdot s^n}{\sum_{m = 0}^{N_p} a_m \cdot s^m} = \frac{\prod_{n = 0}^{N_z}(s-z_n)}{\prod_{m = 0}^{N_p} (s-p_m)} \]
\form#118:\[ H(z) = \sum_{n = 0}^{ord} h_n z^{-n} \]
\form#119:$ F(s) $
\form#120:$F(s)$
\form#121:$Y(s) = (H \circ F)(s) = H(F(s))$
\form#122:\[ H(s) = \frac{\sum\limits_{m = 0}^{n} b_m s^m} {\sum\limits_{k = 0}^{n} a_k s^k}, \quad F(s) = \frac{\sum\limits_{m = 0}^{p} d_m s^m} {\sum\limits_{k = 0}^{p} c_k s^k}, \quad Y(s) = \frac{\sum\limits_{m = 0}^{n p} \beta_m s^m} {\sum\limits_{k = 0}^{n p} \alpha_k s^k} \]
\form#123:$Y(s) = (H \circ F)(s)$
\form#124:\[ s \leftarrow \frac{1 - z^{-1}}{1 - z^{-1}}. \]
\form#125:$\Omega$
\form#126:\[ \Omega = \tan(\omega / 2). \]
\form#127:\[ s(t) = \sum\limits_{n = 0}^{n_{\omega}-1} S(\omega_n) \exp(j\omega_n t) \]
\form#128:$\omega_n$
\form#129:$S(\omega_n)$
\form#130:$ I_0(x)$
\form#131:$ x $
\form#132:$[0 \ 3]$
\form#133:$ \textrm{sinc}(x,a) = \frac{\sin(ax)}{ax}$
\form#134:\[ \textrm{Si}(x) = \int_{0}^{x} \frac{\sin(x)}{x} \, dx\]
\form#135:$[-6\pi \ 6\pi]$
\form#136:$P_N(x)$
\form#137:$N-$
\form#138:\[ P_N(x) = a_0 + a_1 \cdot x + a_2 \cdot x^2 + a_3 \cdot x^3 + ... a_N \cdot x^N. \]
\form#139:\[ P_N(x) = a_0 + x \cdot (a_1 + x \cdot (a_2 + \cdot ( \ldots x \cdot (a_{N-1} + x\cdot a_N) \ldots ))) \]
\form#140:$\mu$
\form#141:$\sigma$
\form#142:$\sigma^2$
\form#143:$\mu = 0$
\form#144:$\sigma=1$

Wyświetl plik

@ -218,6 +218,15 @@
*/
/*!
\ingroup ERROR_CODE_GROUP
\def ERROR_MALLOC
\brief Ошибка динамического выделения памяти.
Данная ошибка означает, что при динамическом выделении памяти произошла ошбика.
В результате функция `malloc` в теле вызваемой функции вернула `NULL` указатель.
Дальнейшая обработка функцией невозможна.
*/
/*!
\ingroup ERROR_CODE_GROUP
@ -226,7 +235,6 @@
*/
/*!
\ingroup ERROR_CODE_GROUP
\def ERROR_MIN_MAX
@ -234,7 +242,6 @@
*/
/*!
\ingroup ERROR_CODE_GROUP
\def ERROR_NEGATIVE
@ -245,6 +252,16 @@
/*!
\ingroup ERROR_CODE_GROUP
\def ERROR_POLY_AN
\brief Неверно задан старший коэффициент полинома.
Например при вычислении кореней полинома степени \f$N\f$
\f[ P_N(x) = a_0 + a_1 x + a_2 x^2 + \ldots a_N x^N \f]
старший коэффициентом \f$a_N\f$ не может быть равным нулю.
*/
/*!
\ingroup ERROR_CODE_GROUP

Wyświetl plik

@ -26,6 +26,36 @@
#include "blas.h"
/******************************************************************************
Vector linear transformation
*******************************************************************************/
int array_scale_lin(double* x, int n,
double xmin, double xmax, double dx,
double h, double* y)
{
double kx;
int k;
if(!x)
return ERROR_PTR;
if(n<1)
return ERROR_SIZE;
if(h<0.0)
return ERROR_NEGATIVE;
if(xmin >= xmax)
return ERROR_MIN_MAX;
kx = h / (xmax - xmin);
for(k = 0; k < n; k++)
y[k] = (x[k] - xmin) * kx + dx;
return RES_OK;
}
/******************************************************************************
\fn int concat(void* a, size_t na, void* b, size_t nb, void* c)
\brief
@ -489,32 +519,19 @@ int DSPL_API logspace(double x0, double x1, int n, int type, double* x)
}
/******************************************************************************
Vector linear transformation
/*******************************************************************************
Oned double array
*******************************************************************************/
int array_scale_lin(double* x, int n,
double xmin, double xmax, double dx,
double h, double* y)
int DSPL_API ones(double* x, int n)
{
double kx;
int k;
int i;
if(!x)
return ERROR_PTR;
if(n<1)
return ERROR_SIZE;
if(h<0.0)
return ERROR_NEGATIVE;
if(xmin >= xmax)
return ERROR_MIN_MAX;
kx = h / (xmax - xmin);
for(k = 0; k < n; k++)
y[k] = (x[k] - xmin) * kx + dx;
return RES_OK;
for(i = 0; i < n; i++)
x[i] = 1.0;
return RES_OK;
}

Wyświetl plik

@ -63,6 +63,42 @@ int DSPL_API poly_z2a_cmplx(complex_t* z, int nz, int ord, complex_t* a)
/******************************************************************************
Real polynomial roots calculation
*******************************************************************************/
int DSPL_API polyroots(double* a, int ord, complex_t* r, int* info)
{
complex_t *t = NULL;
int m;
int err;
if(!a || !r)
return ERROR_PTR;
if(ord<0)
return ERROR_POLY_ORD;
if(a[ord] == 0.0)
return ERROR_POLY_AN;
t = (complex_t*)malloc(ord * ord * sizeof(complex_t));
if(!t)
return ERROR_MALLOC;
for(m = 0; m < ord-1; m++)
{
RE(t[m * (ord+1) + 1]) = 1.0;
RE(t[m + ord * (ord - 1)]) = -a[m] / a[ord];
}
RE(t[ord * ord - 1]) = -a[ord-1] / a[ord];
err = matrix_eig_cmplx(t, ord, r, info);
return err;
}
/******************************************************************************
Real polynomial evaluation
*******************************************************************************/

Wyświetl plik

@ -0,0 +1,34 @@
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "dspl.h"
#define N 2
int main(int argc, char* argv[])
{
void* hdspl; /* DSPL handle */
double a[N+1] = {2.0, 2.0, 1.0}; /* P(x) = 0.5 - 2x + x^2 */
complex_t r[N] = {0}; /* roots */
int err, n, info;
hdspl = dspl_load(); /* Load DSPL functions */
if(!hdspl)
{
printf("libdspl loading error!\n");
return -1;
}
/* roots calculation */
err = polyroots(a, N, r, &info);
printf("Error code: %.8x\n");
/* print roots */
for(n = 0; n < N; n++)
printf("r[%d] = % -8.5f% -8.5f j\n", n, RE(r[n]), IM(r[n]));
/* free dspl handle */
dspl_free(hdspl);
return 0;
}

Wyświetl plik

@ -139,7 +139,10 @@ p_matrix_transpose_cmplx matrix_transpose_cmplx ;
p_matrix_transpose_hermite matrix_transpose_hermite ;
p_minmax minmax ;
p_ones ones ;
p_poly_z2a_cmplx poly_z2a_cmplx ;
p_polyroots polyroots ;
p_polyval polyval ;
p_polyval_cmplx polyval_cmplx ;
@ -337,8 +340,11 @@ void* dspl_load()
LOAD_FUNC(matrix_transpose_cmplx);
LOAD_FUNC(matrix_transpose_hermite);
LOAD_FUNC(minmax);
LOAD_FUNC(ones);
LOAD_FUNC(poly_z2a_cmplx);
LOAD_FUNC(polyroots);
LOAD_FUNC(polyval);
LOAD_FUNC(polyval_cmplx);

Wyświetl plik

@ -130,12 +130,14 @@ typedef struct
/* L 0x12xxxxxx*/
#define ERROR_LAPACK 0x12011601
/* M 0x13xxxxxx*/
#define ERROR_MALLOC 0x13011212
#define ERROR_MATRIX_SIZE 0x13011926
#define ERROR_MIN_MAX 0x13091413
/* N 0x14xxxxxx*/
#define ERROR_NEGATIVE 0x14050701
/* O 0x15xxxxxx*/
/* P 0x16xxxxxx*/
#define ERROR_POLY_AN 0x16150114
#define ERROR_POLY_ORD 0x16151518
#define ERROR_PTR 0x16201800
/* Q 0x17xxxxxx*/
@ -795,11 +797,19 @@ DECLARE_FUNC(int, minmax, double* x
COMMA double* xmin
COMMA double* xmax);
/*----------------------------------------------------------------------------*/
DECLARE_FUNC(int, ones, double* x
COMMA int n);
/*----------------------------------------------------------------------------*/
DECLARE_FUNC(int, poly_z2a_cmplx, complex_t*
COMMA int
COMMA int
COMMA complex_t*);
/*----------------------------------------------------------------------------*/
DECLARE_FUNC(int, polyroots, double* a
COMMA int ord
COMMA complex_t* r
COMMA int* info);
/*----------------------------------------------------------------------------*/
DECLARE_FUNC(int, polyval, double*
COMMA int
COMMA double*