]>
Commit | Line | Data |
---|---|---|
fe4da5cc | 1 | * |
2 | * $Id$ | |
3 | * | |
4 | * $Log$ | |
5 | * Revision 1.1.1.1 1996/04/01 15:02:28 mclareni | |
6 | * Mathlib gen | |
7 | * | |
8 | * | |
9 | #include "gen/pilot.h" | |
10 | #if defined(CERNLIB_DOUBLE) | |
11 | SUBROUTINE DCHEBN(M,N,A,MDIM,B,TOL,RELERR,X,RESMAX,IRK,ITER,IOCD) | |
12 | ||
13 | IMPLICIT DOUBLE PRECISION (A-H,O-Z) | |
14 | ||
15 | DIMENSION A(MDIM,*),B(*),X(*) | |
16 | ||
17 | PARAMETER (R0 = 0, R1 = 1, R2 = 2) | |
18 | ||
19 | DATA BIG /1E+37/ | |
20 | ||
21 | IRK=N | |
22 | NR=1 | |
23 | RELTMP=RELERR | |
24 | RELERR=0 | |
25 | ||
26 | CALL DVSET(M,R1,A(1,N+1),A(2,N+1)) | |
27 | CALL DVSCL(M,-R1,B(1),B(2),A(1,N+2),A(2,N+2)) | |
28 | DO 10 I = 1,M | |
29 | A(I,N+3)=N+I | |
30 | 10 CONTINUE | |
31 | DO 20 J = 1,N | |
32 | A(M+1,J)=J | |
33 | 20 CONTINUE | |
34 | A(M+1,N+1)=0 | |
35 | ITER=0 | |
36 | IOCD=1 | |
37 | CALL DVSET(N,R0,X(1),X(2)) | |
38 | ||
39 | LEV=1 | |
40 | ||
41 | K=0 | |
42 | NK=N+1 | |
43 | 30 K=K+1 | |
44 | NK=NK-1 | |
45 | MODE=0 | |
46 | CALL DVSET(M-K+1,R1,B(1),B(2)) | |
47 | ||
48 | 50 D=-BIG | |
49 | DO 60 I = K,M | |
50 | IF(B(I) .NE. 0) THEN | |
51 | DD=ABS(A(I,N+2)) | |
52 | IF(DD .GT. D) THEN | |
53 | IQ=I | |
54 | D=DD | |
55 | ENDIF | |
56 | ENDIF | |
57 | 60 CONTINUE | |
58 | IF(K .LE. 1 .AND. D .LE. TOL) THEN | |
59 | RESMAX=0 | |
60 | MODE=2 | |
61 | GOTO 380 | |
62 | ENDIF | |
63 | D=TOL | |
64 | DO 80 J = 1,NK | |
65 | DD=ABS(A(IQ,J)) | |
66 | IF(DD .GT. D) THEN | |
67 | IP=J | |
68 | D=DD | |
69 | ENDIF | |
70 | 80 CONTINUE | |
71 | IF(D .GT. TOL) GOTO 330 | |
72 | ||
73 | B(IQ)=0 | |
74 | IF(MODE .EQ. 1) GOTO 50 | |
75 | DO 100 I = K,M | |
76 | IF(B(I) .NE. 0) THEN | |
77 | DO 90 J = 1,NK | |
78 | IF(ABS(A(I,J)) .GT. TOL) THEN | |
79 | MODE=1 | |
80 | GOTO 50 | |
81 | ENDIF | |
82 | 90 CONTINUE | |
83 | ENDIF | |
84 | 100 CONTINUE | |
85 | IRK=K-1 | |
86 | NR=N+1-IRK | |
87 | IOCD=0 | |
88 | GOTO 160 | |
89 | ||
90 | 110 CALL DVXCH(N+3,A(IQ,1),A(IQ,2),A(K,1),A(K,2)) | |
91 | CALL DVXCH(M+1,A(1,IP),A(2,IP),A(1,NK),A(2,NK)) | |
92 | IF(K .LT. N) GOTO 30 | |
93 | 160 IF(IRK .EQ. M) GOTO 380 | |
94 | ||
95 | LEV=2 | |
96 | ||
97 | D=TOL | |
98 | DO 170 I = IRK+1,M | |
99 | DD=ABS(A(I,N+2)) | |
100 | IF(DD .GT. D) THEN | |
101 | IQ=I | |
102 | D=DD | |
103 | ENDIF | |
104 | 170 CONTINUE | |
105 | ||
106 | IF(D .LE. TOL) THEN | |
107 | RESMAX=0 | |
108 | MODE=3 | |
109 | GOTO 380 | |
110 | ENDIF | |
111 | IF(A(IQ,N+2) .GE. -TOL) THEN | |
112 | A(IQ,N+1)=2-A(IQ,N+1) | |
113 | CALL DVSCL(N+4-NR,-R1,A(IQ,1),A(IQ,2),A(IQ,1),A(IQ,2)) | |
114 | A(IQ,N+1)=-A(IQ,N+1) | |
115 | ENDIF | |
116 | DO 220 J = NR,N | |
117 | IF(A(IQ,J) .GE. TOL) THEN | |
118 | CALL DVSCA(M,R2, | |
119 | 1 A(1,J),A(2,J),A(1,N+1),A(2,N+1),A(1,N+1),A(2,N+1)) | |
120 | CALL DVSCL(M,-R1,A(1,J),A(2,J),A(1,J),A(2,J)) | |
121 | A(M+1,J)=-A(M+1,J) | |
122 | ENDIF | |
123 | 220 CONTINUE | |
124 | IP=N+1 | |
125 | GOTO 330 | |
126 | ||
127 | 230 IF(IRK+1 .EQ. M) GO TO 380 | |
128 | CALL DVXCH(IRK+3,A(IQ,1),A(IQ,2),A(M,1),A(M,2)) | |
129 | ||
130 | LEV=3 | |
131 | ||
132 | 260 D=-TOL | |
133 | H=2*A(M,N+2) | |
134 | DO 280 I = IRK+1,M-1 | |
135 | IF(A(I,N+2) .LT. D) THEN | |
136 | IQ=I | |
137 | D=A(I,N+2) | |
138 | MODE=0 | |
139 | ELSE | |
140 | DD=H-A(I,N+2) | |
141 | IF(DD .LT. D) THEN | |
142 | IQ=I | |
143 | D=DD | |
144 | MODE=1 | |
145 | ENDIF | |
146 | ENDIF | |
147 | 280 CONTINUE | |
148 | IF(D .GE. -TOL) GOTO 380 | |
149 | DD=-D/A(M,N+2) | |
150 | IF(DD .LT. RELTMP) THEN | |
151 | RELERR=DD | |
152 | MODE=4 | |
153 | GOTO 380 | |
154 | ENDIF | |
155 | IF(MODE .NE. 0) THEN | |
156 | CALL DVSCS(IRK+1,R2, | |
157 | 1 A(M,1),A(M,2),A(IQ,1),A(IQ,2),A(IQ,1),A(IQ,2)) | |
158 | A(IQ,N+2)=D | |
159 | A(IQ,N+3)=-A(IQ,N+3) | |
160 | ENDIF | |
161 | D=BIG | |
162 | DO 320 J = NR,N+1 | |
163 | IF(A(IQ,J) .GT. TOL) THEN | |
164 | DD=A(M,J)/A(IQ,J) | |
165 | IF(DD .LT. D) THEN | |
166 | IP=J | |
167 | D=DD | |
168 | ENDIF | |
169 | ENDIF | |
170 | 320 CONTINUE | |
171 | IF(D .LT. BIG) GO TO 330 | |
172 | IOCD=2 | |
173 | GOTO 380 | |
174 | ||
175 | 330 RPVT=1/A(IQ,IP) | |
176 | CALL DVSCL(M,RPVT,A(1,IP),A(2,IP),A(1,IP),A(2,IP)) | |
177 | DO 360 I = 1,M | |
178 | IF(I .NE. IQ) THEN | |
179 | D=A(I,IP) | |
180 | CALL DVSCA(N+3-NR,-D, | |
181 | 1 A(IQ,NR),A(IQ,NR+1),A(I,NR),A(I,NR+1),A(I,NR),A(I,NR+1)) | |
182 | A(I,IP)=D | |
183 | ENDIF | |
184 | 360 CONTINUE | |
185 | CALL DVSCL(IRK+2,-RPVT,A(IQ,NR),A(IQ,NR+1),A(IQ,NR),A(IQ,NR+1)) | |
186 | A(IQ,IP)=RPVT | |
187 | D=A(M+1,IP) | |
188 | A(M+1,IP)=A(IQ,N+3) | |
189 | A(IQ,N+3)=D | |
190 | ITER=ITER+1 | |
191 | GOTO (110,230,260), LEV | |
192 | ||
193 | 380 CALL DVSET(M,R0,B(1),B(2)) | |
194 | IF(MODE .EQ. 2) GOTO 450 | |
195 | DO 400 I = 1,IRK | |
196 | X(INT(A(I,N+3)))=A(I,N+2) | |
197 | 400 CONTINUE | |
198 | IF(MODE .EQ. 3 .OR. IRK .EQ. M) GOTO 450 | |
199 | DO 410 J = NR,N+1 | |
200 | B(INT(ABS(A(M+1,J)))-N)=A(M,N+2)*SIGN(R1,A(M+1,J)) | |
201 | 410 CONTINUE | |
202 | DO 420 I = IRK+1,M-1 | |
203 | B(INT(ABS(A(I,N+3)))-N)=(A(M,N+2)-A(I,N+2))*SIGN(R1,A(I,N+3)) | |
204 | 420 CONTINUE | |
205 | 430 DO 440 J = NR,N+1 | |
206 | IF(ABS(A(M,J)) .LE. TOL) THEN | |
207 | IOCD=0 | |
208 | GOTO 450 | |
209 | ENDIF | |
210 | 440 CONTINUE | |
211 | 450 IF(MODE .NE. 2 .AND. MODE .NE. 3) RESMAX=A(M,N+2) | |
212 | IF(IRK .EQ. M) RESMAX=0 | |
213 | IF(MODE .EQ. 4) RESMAX=RESMAX-D | |
214 | RETURN | |
215 | END | |
216 | #endif |