RECQR
Computes recurrence coefficients for monic polynomials given a quadrature rule.
Required Arguments
QX — Array of length N containing the quadrature points. (Input)
QW — Array of length N containing the quadrature weights. (Input)
B — Array of length NTERM containing recurrence coefficients. (Output)
C — Array of length NTERM containing recurrence coefficients. (Output)
Optional Arguments
N — Number of quadrature points. (Input)
Default: N = size (QX,1).
NTERM — Number of recurrence coefficients. (Input)
NTERM must be less than or equal to N.
Default: NTERM = size (B,1).
FORTRAN 90 Interface
Generic: CALL RECQR (QX, QW, B, C [, …])
Specific: The specific interface names are S_RECQR and D_RECQR.
FORTRAN 77 Interface
Single: CALL RECQR (N, QX, QW, NTERM, B, C)
Double: The double precision name is DRECQR.
Description
The routine RECQR produces the recurrence coefficients for the orthogonal polynomials given the points and weights for the Gauss quadrature formula. It is assumed that the orthogonal polynomials are monic; hence the three-term recursion may be written
where p0 = 1 and p-1 = 0. It is obvious from this representation that the degree of pi is i and that pi is monic. In order for the recurrence to give rise to a sequence of orthogonal polynomials (with respect to a nonnegative measure), it is necessary and sufficient that ci > 0.
This routine is an inverse routine to GQRCF. Given the recurrence coefficients, the routine GQRCF produces the corresponding Gauss quadrature formula, whereas the routine RECQR produces the recurrence coefficients given the quadrature formula.
Comments
1. Workspace may be explicitly provided, if desired, by use of R2CQR/DR2CQR. The reference is:
CALL R2CQR (N, QX, QW, NTERM, B, C, WK)
The additional argument is:
WK — Work array of length 2 * N.
2. The recurrence coefficients B(I) and C(I) define the monic polynomials via the relation P(I) = (X − B(I + 1)) * P(I − 1) − C(I + 1) * P(I − 2). The zero-th moment
of the weight function is returned in C(1).
Example
To illustrate the use of RECQR, we will input a simple choice of recurrence coefficients, call GQRCF for the quadrature formula, put this information into RECQR, and recover the recurrence coefficients.
USE RECQR_INT
USE UMACH_INT
USE GQRCF_INT
IMPLICIT NONE
INTEGER N
PARAMETER (N=5)
INTEGER I, J, NFIX, NOUT, NTERM
REAL B(N), C(N), FLOAT, QW(N), QX(N), QXFIX(2)
INTRINSIC FLOAT
! Get output unit number
CALL UMACH (2, NOUT)
NFIX = 0
! Set arrays B and C of recurrence
! coefficients
DO 10 J=1, N
B(J) = FLOAT(J)
C(J) = FLOAT(J)/2.0
10 CONTINUE
WRITE (NOUT,99995)
99995 FORMAT (1X, 'Original recurrence coefficients')
WRITE (NOUT,99996) (I,B(I),I,C(I),I=1,N)
99996 FORMAT (5(6X,'B(',I1,') = ',F8.4,7X,'C(',I1,') = ',F8.5,/))
!
! The call to GQRCF will compute the
! quadrature rule from the recurrence
! coefficients given above.
!
CALL GQRCF (N, B, C, QX, QW)
WRITE (NOUT,99997)
99997 FORMAT (/, 1X, 'Quadrature rule from the recurrence coefficients' &
)
WRITE (NOUT,99998) (I,QX(I),I,QW(I),I=1,N)
99998 FORMAT (5(6X,'QX(',I1,') = ',F8.4,7X,'QW(',I1,') = ',F8.5,/))
!
! Call RECQR to recover the original
! recurrence coefficients
NTERM = N
CALL RECQR (QX, QW, B, C)
WRITE (NOUT,99999)
99999 FORMAT (/, 1X, 'Recurrence coefficients determined by RECQR')
WRITE (NOUT,99996) (I,B(I),I,C(I),I=1,N)
!
END
Output
Original recurrence coefficients
B(1) = 1.0000 C(1) = 0.50000
B(2) = 2.0000 C(2) = 1.00000
B(3) = 3.0000 C(3) = 1.50000
B(4) = 4.0000 C(4) = 2.00000
B(5) = 5.0000 C(5) = 2.50000
Quadrature rule from the recurrence coefficients
QX(1) = 0.1525 QW(1) = 0.25328
QX(2) = 1.4237 QW(2) = 0.17172
QX(3) = 2.7211 QW(3) = 0.06698
QX(4) = 4.2856 QW(4) = 0.00790
QX(5) = 6.4171 QW(5) = 0.00012
Recurrence coefficients determined by RECQR
B(1) = 1.0000 C(1) = 0.50000
B(2) = 2.0000 C(2) = 1.00000
B(3) = 3.0000 C(3) = 1.50000
B(4) = 4.0000 C(4) = 2.00000
B(5) = 5.0000 C(5) = 2.50000
