libdspl-2.0
Библиотека алгоритмов цифровой обработки сигналов
ellip_asn_cmplx.c
1 /*
2 * Copyright (c) 2015-2024 Sergey Bakhurin
3 * Digital Signal Processing Library [http://dsplib.org]
4 *
5 * This file is part of DSPL.
6 *
7 * is free software: you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation, either version 3 of the License, or
10 * (at your option) any later version.
11 *
12 * DSPL is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with Foobar. If not, see <http://www.gnu.org/licenses/>.
19 */
20 
21 
22 #include <stdio.h>
23 #include <stdlib.h>
24 #include <string.h>
25 #include <math.h>
26 #include "dspl.h"
27 
28 
29 
30 #ifdef DOXYGEN_ENGLISH
31 
64 #endif
65 #ifdef DOXYGEN_RUSSIAN
66 
104 #endif
105 int DSPL_API ellip_asn_cmplx(complex_t* w, int n, double k, complex_t* u)
106 {
107  double lnd[ELLIP_ITER], t;
108  complex_t tmp0, tmp1;
109  int i, m;
110 
111  if(!u || !w)
112  return ERROR_PTR;
113  if(n<1)
114  return ERROR_SIZE;
115  if(k < 0.0 || k>= 1.0)
116  return ERROR_ELLIP_MODULE;
117 
118  ellip_landen(k,ELLIP_ITER, lnd);
119 
120  for(m = 0; m < n; m++)
121  {
122  RE(u[m]) = RE(w[m]);
123  IM(u[m]) = IM(w[m]);
124  for(i = 1; i < ELLIP_ITER; i++)
125  {
126  RE(tmp0) = lnd[i-1]*RE(u[m]);
127  IM(tmp0) = lnd[i-1]*IM(u[m]);
128  RE(tmp1) = 1.0 - CMRE(tmp0, tmp0);
129  IM(tmp1) = - CMIM(tmp0, tmp0);
130 
131  sqrt_cmplx(&tmp1, 1, &tmp0);
132  RE(tmp0) += 1.0;
133 
134  RE(tmp1) = RE(tmp0) * (1.0 + lnd[i]);
135  IM(tmp1) = IM(tmp0) * (1.0 + lnd[i]);
136 
137  t = 2.0 / ABSSQR(tmp1);
138 
139  RE(tmp0) = t * CMCONJRE(u[m], tmp1);
140  IM(tmp0) = t * CMCONJIM(u[m], tmp1);
141 
142  RE(u[m]) = RE(tmp0);
143  IM(u[m]) = IM(tmp0);
144  }
145 
146  asin_cmplx(&tmp0, 1, u+m);
147  t = 2.0 / M_PI;
148  RE(u[m]) *= t;
149  IM(u[m]) *= t;
150  }
151  return RES_OK;
152 }
#define RE(x)
Макрос определяющий реальную часть комплексного числа.
Definition: dspl.h:420
#define ERROR_PTR
Ошибка указателя. Данная ошибка означает, что один из обязательных указателей (память под который дол...
Definition: dspl.h:610
#define ERROR_SIZE
Ошибка при передаче размера массива. Данная ошибка возникает когда помимо указателя на массив входных...
Definition: dspl.h:618
int ellip_asn_cmplx(complex_t *w, int n, double k, complex_t *u)
Обратная эллиптическая функция Якоби комплексного аргумента
double complex_t[2]
Описание комплексного типа данных.
Definition: dspl.h:86
#define RES_OK
Функция завершилась корректно. Ошибки отсутствуют.
Definition: dspl.h:558
#define ABSSQR(x)
Макрос возвращает квадрат модуля комплексного числа x.
Definition: dspl.h:534
int sqrt_cmplx(complex_t *x, int n, complex_t *y)
Квадратный корень из комплексного вектора x (поэлементный).
Definition: sqrt_cmplx.c:136
#define ERROR_ELLIP_MODULE
Модуль эллиптического интеграла Якоби должен быть от 0 до 1. Данная ошибка возникает при расчете элли...
Definition: dspl.h:569
int asin_cmplx(complex_t *x, int n, complex_t *y)
Арксинус комплексного аргумента x.
Definition: asin_cmplx.c:138
#define IM(x)
Макрос определяющий мнимую часть комплексного числа.
Definition: dspl.h:478
int ellip_landen(double k, int n, double *y)
Расчет коэффициентов ряда полного эллиптического интеграла.
Definition: ellip_landen.c:175