Introduction and theoretical background - Earth Dynamics: Deformations and Oscillations of the Rotating Earth - page 79

DO 23 J=1,N2
C(K,J)=C(K,J)-C(IR(I),J)*ELL
23 CONTINUE
22 CONTINUE
21 CONTINUE
13 CONTINUE
C Unscramble equations and record row interchanges.
DO 24 I=1,N
IDET(IR(I))=IC(I)
DO 24 J=1,N2
HOLD(IC(I),J)=C(IR(I),J)
24 CONTINUE
C Put equations back in original matrix.
DO 25 I=1,N
DO 25 J=1,N2
C(I,J)=HOLD(I,J)
25 CONTINUE
C Determine sign of determinant.
ICH=0
DO 26 I=1,NM1
IP1=I+1
DO 26 J=IP1,N
IF(IDET(J).GE.IDET(I))GO TO 26
IHOLD=IDET(J)
IDET(J)=IDET(I)
IDET(I)=IHOLD
ICH=ICH+1
26 CONTINUE
IF((ICH/2)*2.NE.ICH)DET=-DET
C Put inverse in A and solution in B.
DO 27 I=1,N
B(I)=C(I,NP1)
DO 27 J=1,N
JJ=J+NP1
A(I,J)=C(I,JJ)

27 CONTINUE

RETURN
END

1.6 Interpolation and approximation
In the study of Earth's dynamics, we are more often concerned with functions
represented numerically than by simple analytical expressions. Thus, interpolation
and approximation are of central importance.
Commonly, we have a numerical representation of a function, f ( x ), as a table of
function values ( f 1 ,..., f N )at N discrete locations, or nodes, ( x 1 ,..., x N ) along the
x -axis. The question is then how best to represent the function between its tabulated
values. Let us consider its representation s ( x ) on the interval ( x i , x i + 1 ).
The crudest approximation is to take it to be a constant equal to one of the values
at the end points, either f i or f i + 1 , or the mean value, giving

f i + f i + 1
2 .

s ( x ) = (1.542)

The approximating function then appears as a histogram.