Geology Reference
In-Depth Information
GO TO 50
49 E(I)=1.D0
C Store e(r)/e(0) IN EI.
50 J=NK(M)+I
EI(J)=E(I)
48 CONTINUE
39 CONTINUE
C Find surface value e(d).
C Generate new coefficient of quadratic.
BB=-(14.D0+7.D0*ETAI(NT)+6.D0*FRD)/4.D0
ED=-(BB+DSQRT(BB*BB-4.D0*CC))/2.D0
C Find e(0).
E0=ED/EI(NT)
C Find final e(r).
DO 51 I=1,NT
EI(I)=E0*EI(I)
51 CONTINUE
C Begin calculation of kappa by Runge-Kutta integration.
DO 52 M=1,4
MP1=M+1
C Fix number of model points for region M.
N1=NM(M)
C Restore interpolation matrix for region M.
DO 53 I=1,N1
DO 54 J=1,N1
C(I,J)=D(I,J,M)
54 CONTINUE
53 CONTINUE
C Put radius, density, mean density, ETA and e in
C active locations.
DO 55 I=1,N1
J=NK(M)+I
R(I)=RI(J)
RHO(I)=RHOI(J)
RHOB(I)=RHOBI(J)
ETA(I)=ETAI(J)
E(I)=EI(J)
55 CONTINUE
C Begin calculation of the complementary function YYC.
C Set IC switch to unity for complementary function.
IC=1
C Set initial value of YYC vector at beginning of region.
IF(M.EQ.1)GO TO 56
YYC(1,1)=INITYC(1,M)
YYC(2,1)=INITYC(2,M)
GO TO 57
56 YYC(1,1)=0.D0
YYC(2,1)=0.D0
57 CONTINUE
C Do Runge-Kutta extension of the complementary function through region.
CALL RK(YYC,IC,RHS,RHO,RHOB,E,NIS,N1,C,R,M1)
C Store starting value for next region.
INITYC(1,MP1)=YYC(1,N1)
INITYC(2,MP1)=YYC(2,N1)
C Store YYC in YYCI.
DO 58 I=1,N1
J=NK(M)+I
YYCI(1,J)=YYC(1,I)
YYCI(2,J)=YYC(2,I)
58 CONTINUE
C Begin calculation of the particular integral YYP.
Search WWH ::
Custom Search