|
Текст подпрограммы и версий qsjac_c.zip |
Тексты тестовых примеров tqsjac_c.zip |
Вычисление узлов и весов квадратурной формулы Гаусса - Якоби.
Подпрограмма qsjac_c вычисляeт yзлы xi и веса Wi квадратурной формулы Гаусса - Якоби
b
∫ (1 - x)α (1 + x)β f(x) dx ≈
a
N
≈ ∑ Wi f(xi ) по методу Ньютона .
i =1
int qsjac_c (real *x, real *w, real *alf, real *bta,
integer *n, real *a, real *b, real *c)
Параметры
| x - | вещественный одномерный массив размерности 512 вычисленных узлов формулы Гаусса - Якоби; |
| w - | вещественный одномерный массив размерности 512 вычисленных весов формулы Гаусса - Якоби; |
|
alf - bta | заданные параметры α и β , соответственно, в весовой функции (1 - x)α (1 + x) β (тип: вещественный); |
| n - | заданное число узлов (и весов) интегрирования в формуле Гаусса - Якоби (2 ≤ n ≤ 512, тип: целый); |
| a, b, c - | вещественные одномерные массивы размерности 512, заданных коэффициентов рекуррентного соотношения для полиномов Якоби. |
Версии: нет
Вызываемые подпрограммы
| qsjac1_c - | служебная подпрограмма; |
| qsjac2_c - | вещественная подпрограмма - функция, вычисляющая значение логарифма от гамма - функции ln ( Г (x) ). |
Замечания по использованию
|
Подпрограмма qsjac_c вычисляет узлы и веса квадратурной формулы Гаусса - Якоби n
∑ wi f(xi )
i =1
с абсолютной погрешностью 10 - 10 и для значений 2 ≤ n ≤ 512. Перед обращением к подпрограмме qsjac_c должны быть вычислены коэффициенты Anα , β, Bnα , β, Cnα , β рекуррентного соотношения для полиномов Якоби Jnα , β (x) = ( Anα , β x + Bnα , β ) Jn -1α , β (x) - Cnα , β jn -2α , β (x)
по формулам:
( α + β + 2n ) ( α + β + 2n - 1 )
Anα , β = ------------------------------------------- ,
2n ( α + β + n )
( α + β ) ( α - β ) ( α + β + 2n - 1 )
Bnα , β = ------------------------------------------------ ,
2n ( α + β + n ) ( α + β + 2n - 2 )
( α + β + 2n ) ( α + 2n - 1 ) ( β + 2n - 1 )
Cnα , β = ------------------------------------------------------ .
n ( α + β + n ) ( α + β + 2n - 2 )
|
int main(void)
{
/* System generated locals */
int i__1;
/* Local variables */
static float a[512], b[512], c__[512];
static int k, n;
static float p, w[512], x[512];
extern int qsjac_c(float *, float *, float *, float *,
int *, float *, float *, float *);
static float q1, q2, q3, q4, q5, fn;
static int nn;
static float fn1, dba, alf, sab, bta;
alf = 1.f;
bta = 1.f;
nn = 8;
sab = alf + bta;
dba = alf - bta;
p = sab * dba;
a[0] = 0.f;
b[0] = 0.f;
c__[0] = 0.f;
for (n = 2; n <= 512; ++n) {
fn = (float) n;
fn1 = (float) (n - 1);
q1 = fn * 2.f + sab;
q2 = q1 - 1.f;
q3 = sab + fn;
q4 = fn * q3;
q5 = q4 * (q1 - 2.f);
a[n - 1] = q1 * q2 / (q4 * 2.f);
b[n - 1] = p * q2 / (q5 * 2.f);
c__[n - 1] = q1 * (alf + fn1) * (bta + fn1) / q5;
/* l1: */
}
qsjac_c(x, w, &alf, &bta, &nn, a, b, c__);
i__1 = nn;
for (k = 1; k <= i__1; ++k) {
/* l2: */
printf("\n %16.7f %16.7f \n",x[k - 1],w[k - 1]);
}
return 0;
} /* main */
Результаты:
x1 = 0.9195339081 , w1 = 0.0205900956 ,
x2 = 0.7387738651 , w2 = 0.1021477023 ,
x3 = 0.4779249498 , w3 = 0.2253365549 ,
x4 = 0.1652789576 , w4 = 0.3185923136 ,
x5 = - 0.1652789576 , w5 = 0.3185923136 ,
x6 = - 0.4779249498 , w6 = 0.2253365549 ,
x7 = - 0.7387738651 , w7 = 0.1021477023 ,
x8 = - 0.9195339081 , w8 = 0.0205900956 .