For completeness, we need some examples that include
- integration of the PDE variables over the entire domain or individual elements/material domains, for $m = 0, 1, 2$
- calculation of cumulative integrals, at the knots of the Chebyshev basis (see this question: https://scicomp.stackexchange.com/q/23398/37438)
The coefficients in the Chebyshev expansion can be obtained from the work array wkres
as follows (TOMS Algorithm 660, pg. 195):
C
NPTL = NPOLY + 1
DO 100 I = 1,NEL
C ITH ELEMENT
C IU IS THE COMPONENT OF SOLUTION VECTOR AT LHS
C OF ELEMENT I
IU = (I - 1)*(NPOLY)*NPDE + 1
DO 80 IS = 1,NPTL
DO 80 JK = 1,NPDE
COEFF(JK,IS,I) = 0.0D0
80 CONTINUE
DO 90 IS = 1,NPTL
DO 90 JS = 1,NPTL
DO 90 JK = 2,NPDE
COEFF(JK,IS,I) = COEFF(JK,IS,I) +
1 WKRES(IS + (JS-1)*NPTL)*Y(IU + JS + JK - 2)
90 CONTINUE
100 CONTINUE
C
Or using free-form Fortran:
nptl = npoly + 1
do i = 1, nel
! ith element
! iu is the component of solution vector at LHS of element i
iu = (i - 1)*npoly*npde + 1
do is = 1, nptl
do jk = 1, npde
coeff(jk,is,i) = 0.0d0
end do
end do
do is = 1, nptl
do js = 1, nptl
do jk = 2, npde
coeff(jk,is,i) = coeff(jk,is,i) + &
wkres(is + (js - 1)*nptl) * y(iu + js + jk - 2)
end do
end do
end do
end do
The coefficient $a_{j,n}(t)$ of the $k$-th PDE is stored in COEFF(K,J,N)
, where $j$ runs along elements in the mesh, and $n$ counts the Chebyshev basis polynomials $T_n(x)$.