[2563] | 1 | /*************************************************************************
|
---|
| 2 | >>> SOURCE LICENSE >>>
|
---|
| 3 | This program is free software; you can redistribute it and/or modify
|
---|
| 4 | it under the terms of the GNU General Public License as published by
|
---|
| 5 | the Free Software Foundation (www.fsf.org); either version 2 of the
|
---|
| 6 | License, or (at your option) any later version.
|
---|
| 7 |
|
---|
| 8 | This program is distributed in the hope that it will be useful,
|
---|
| 9 | but WITHOUT ANY WARRANTY; without even the implied warranty of
|
---|
| 10 | MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
---|
| 11 | GNU General Public License for more details.
|
---|
| 12 |
|
---|
| 13 | A copy of the GNU General Public License is available at
|
---|
| 14 | http://www.fsf.org/licensing/licenses
|
---|
| 15 |
|
---|
| 16 | >>> END OF LICENSE >>>
|
---|
| 17 | *************************************************************************/
|
---|
| 18 |
|
---|
| 19 | using System;
|
---|
| 20 |
|
---|
| 21 | namespace alglib
|
---|
| 22 | {
|
---|
| 23 | public class legendre
|
---|
| 24 | {
|
---|
| 25 | /*************************************************************************
|
---|
| 26 | Calculation of the value of the Legendre polynomial Pn.
|
---|
| 27 |
|
---|
| 28 | Parameters:
|
---|
| 29 | n - degree, n>=0
|
---|
| 30 | x - argument
|
---|
| 31 |
|
---|
| 32 | Result:
|
---|
| 33 | the value of the Legendre polynomial Pn at x
|
---|
| 34 | *************************************************************************/
|
---|
| 35 | public static double legendrecalculate(int n,
|
---|
| 36 | double x)
|
---|
| 37 | {
|
---|
| 38 | double result = 0;
|
---|
| 39 | double a = 0;
|
---|
| 40 | double b = 0;
|
---|
| 41 | int i = 0;
|
---|
| 42 |
|
---|
| 43 | result = 1;
|
---|
| 44 | a = 1;
|
---|
| 45 | b = x;
|
---|
| 46 | if( n==0 )
|
---|
| 47 | {
|
---|
| 48 | result = a;
|
---|
| 49 | return result;
|
---|
| 50 | }
|
---|
| 51 | if( n==1 )
|
---|
| 52 | {
|
---|
| 53 | result = b;
|
---|
| 54 | return result;
|
---|
| 55 | }
|
---|
| 56 | for(i=2; i<=n; i++)
|
---|
| 57 | {
|
---|
| 58 | result = ((2*i-1)*x*b-(i-1)*a)/i;
|
---|
| 59 | a = b;
|
---|
| 60 | b = result;
|
---|
| 61 | }
|
---|
| 62 | return result;
|
---|
| 63 | }
|
---|
| 64 |
|
---|
| 65 |
|
---|
| 66 | /*************************************************************************
|
---|
| 67 | Summation of Legendre polynomials using Clenshaws recurrence formula.
|
---|
| 68 |
|
---|
| 69 | This routine calculates
|
---|
| 70 | c[0]*P0(x) + c[1]*P1(x) + ... + c[N]*PN(x)
|
---|
| 71 |
|
---|
| 72 | Parameters:
|
---|
| 73 | n - degree, n>=0
|
---|
| 74 | x - argument
|
---|
| 75 |
|
---|
| 76 | Result:
|
---|
| 77 | the value of the Legendre polynomial at x
|
---|
| 78 | *************************************************************************/
|
---|
| 79 | public static double legendresum(ref double[] c,
|
---|
| 80 | int n,
|
---|
| 81 | double x)
|
---|
| 82 | {
|
---|
| 83 | double result = 0;
|
---|
| 84 | double b1 = 0;
|
---|
| 85 | double b2 = 0;
|
---|
| 86 | int i = 0;
|
---|
| 87 |
|
---|
| 88 | b1 = 0;
|
---|
| 89 | b2 = 0;
|
---|
| 90 | for(i=n; i>=0; i--)
|
---|
| 91 | {
|
---|
| 92 | result = (2*i+1)*x*b1/(i+1)-(i+1)*b2/(i+2)+c[i];
|
---|
| 93 | b2 = b1;
|
---|
| 94 | b1 = result;
|
---|
| 95 | }
|
---|
| 96 | return result;
|
---|
| 97 | }
|
---|
| 98 |
|
---|
| 99 |
|
---|
| 100 | /*************************************************************************
|
---|
| 101 | Representation of Pn as C[0] + C[1]*X + ... + C[N]*X^N
|
---|
| 102 |
|
---|
| 103 | Input parameters:
|
---|
| 104 | N - polynomial degree, n>=0
|
---|
| 105 |
|
---|
| 106 | Output parameters:
|
---|
| 107 | C - coefficients
|
---|
| 108 | *************************************************************************/
|
---|
| 109 | public static void legendrecoefficients(int n,
|
---|
| 110 | ref double[] c)
|
---|
| 111 | {
|
---|
| 112 | int i = 0;
|
---|
| 113 |
|
---|
| 114 | c = new double[n+1];
|
---|
| 115 | for(i=0; i<=n; i++)
|
---|
| 116 | {
|
---|
| 117 | c[i] = 0;
|
---|
| 118 | }
|
---|
| 119 | c[n] = 1;
|
---|
| 120 | for(i=1; i<=n; i++)
|
---|
| 121 | {
|
---|
| 122 | c[n] = c[n]*(n+i)/2/i;
|
---|
| 123 | }
|
---|
| 124 | for(i=0; i<=n/2-1; i++)
|
---|
| 125 | {
|
---|
| 126 | c[n-2*(i+1)] = -(c[n-2*i]*(n-2*i)*(n-2*i-1)/2/(i+1)/(2*(n-i)-1));
|
---|
| 127 | }
|
---|
| 128 | }
|
---|
| 129 | }
|
---|
| 130 | }
|
---|