The subseismic equation and boundary conditions - Earth Dynamics: Deformations and Oscillations of the Rotating Earth - page 439

C Calculate potential condition parameter AND.
AND=-AN/SUM
GO TO 65
C Construct coefficient matrix for N=0.
57 CONTINUE
CMC1(1,1)=GZEROC*Y(1,3)-Y(2,3)/RHOAC-Y(5,3)
CMC1(1,2)=1.D0
CMC1(2,1)=4.D0*PI*G*RHOAC*Y(1,3)
CMC1(2,2)=0.D0
C Solve for linear combination coefficients and AND.
BVC1(1)=0.D0
BVC2(1)=-1.D0/CMC1(2,1)
AND=0.D0
C Calculate unscaled Love numbers.
SUM1=BVC1(1)*Y(1,3)
SUM2=BVC2(1)*Y(1,3)
SUM3=0.D0
65 CONTINUE
C Scale Love numbers to dimensionless values.
BVC1(1)=GZEROC*SUM1
BVC2(1)=GZEROC*SUM2/RIOB
BVC3(1)=GZEROC*SUM3
C Test if N=1.
IF(N.NE.1) GO TO 66
IF(ICOUNT.EQ.3) GO TO 67
IF(ICOUNT.EQ.2) GO TO 68
GAMA=1.D0/BVC3(1)
COR=CORS
ICOUNT=2
GO TO 28
68 CONTINUE
BETA=(1.D0/BVC3(1)-GAMA)/SIG
ANGS=ANGSS
ICOUNT=3
GO TO 28
67 CONTINUE
ALPHA=(1.D0/BVC3(1)-BETA*SIG-GAMA)/SIGS
C Write out internal Love numbers.
66 CONTINUE
WRITE(6,69)N,BVC1(1),BVC2(1)
69 FORMAT(/I15,11X,F16.8,F28.8)
WRITE(2,70)BVC3(1),AND
70 FORMAT(2D23.15)
C End loop over degree for inner core.
26 CONTINUE
C Begin loop over degree for the shell (mantle + crust).
WRITE(6,25)
DO 71 NP1=1,NN
N=NP1-1
C Set up interpolation for the mantle.
N1=NM(3)
N2=2*N1-2
N3=N1-2
C Put mantle values in active locations.
DO 72 I=1,N1
J=NK(3)+I
R(I)=RI(J)
RHO(I)=RHOI(J)
MU(I)=MUI(J)
LAMBDA(I)=LAMBDAI(J)
GZERO(I)=GZEROI(J)