This commit was generated by cvs2svn to compensate for changes in r2,
[u/mrichter/AliRoot.git] / MINICERN / mathlib / gen / d / minsq.F
1 *
2 * $Id$
3 *
4 * $Log$
5 * Revision 1.1.1.1  1996/04/01 15:02:20  mclareni
6 * Mathlib gen
7 *
8 *
9 #include "gen/pilot.h"
10       SUBROUTINE MINSQ (M,N,F,X,E,IPRINT,NFUN ,NW,W,COV,XSTEP)
11       DIMENSION F(M),X(N),E(N),W(NW)
12       IFLAG=1
13       ESCALE=1.
14       MAXFUN=NFUN
15       MPLUSN=M+N
16       KST=N+MPLUSN
17       NPLUS=N+1
18       KINV=NPLUS*(MPLUSN+1)
19       KSTORE=KINV-MPLUSN-1
20       CALL FCN1(M,N,F,X,IFLAG,NW,W,KEND,KEND2N)
21       IFLAG=4
22       NN=N+N
23       K=NN
24       DO 1 I=1,M
25       K=K+1
26       W(K)=F(I)
27     1 CONTINUE
28       IINV=2
29       K=KST
30       I=1
31     2 X(I)=X(I)+E(I)
32       CALL FCN1(M,N,F,X,IFLAG,NW,W,KEND,KEND2N)
33       X(I)=X(I)-E(I)
34       DO 3 J=1,N
35       K=K+1
36       W(K)=0.
37       W(J)=0.
38     3 CONTINUE
39       SUM=0.
40       KK=NN
41       DO 4 J=1,M
42       KK=KK+1
43       F(J)=F(J)-W(KK)
44       SUM=SUM+F(J)*F(J)
45     4 CONTINUE
46       IF (SUM) 5,5,6
47     5 WRITE(6,7)I
48     7 FORMAT(5X,'E(',I3,') UNREASONABLY SMALL')
49       RETURN
50     6 SUM=1./SQRT (SUM)
51       J=K-N+I
52       W(J)=E(I)*SUM
53       DO 9 J=1,M
54       K=K+1
55       W(K)=F(J)*SUM
56       KK=NN+J
57       DO 11 II=1,I
58       KK=KK+MPLUSN
59       W(II)=W(II)+W(KK)*W(K)
60    11 CONTINUE
61     9 CONTINUE
62       ILESS=I-1
63       IGAMAX=N+I-1
64       INCINV=N-ILESS
65       INCINP=INCINV+1
66       IF (ILESS) 13,13,14
67    13 W(KINV)=1.
68       GO TO 15
69    14 B=1.
70       DO 16 J=NPLUS,IGAMAX
71       W(J)=0.
72    16 CONTINUE
73       KK=KINV
74       DO 17 II=1,ILESS
75       IIP=II+N
76       W(IIP)=W(IIP)+W(KK)*W(II)
77       JL=II+1
78       IF (JL-ILESS) 18,18,19
79    18 DO 20 JJ=JL,ILESS
80       KK=KK+1
81       JJP=JJ+N
82       W(IIP)=W(IIP)+W(KK)*W(JJ)
83       W(JJP)=W(JJP)+W(KK)*W(II)
84    20 CONTINUE
85    19 B=B-W(II)*W(IIP)
86       KK=KK+INCINP
87    17 CONTINUE
88       B=1./B
89       KK=KINV
90       DO 21 II=NPLUS,IGAMAX
91       BB=-B*W(II)
92       DO 22 JJ=II,IGAMAX
93       W(KK)=W(KK)-BB*W(JJ)
94       KK=KK+1
95    22 CONTINUE
96       W(KK)=BB
97       KK=KK+INCINV
98    21 CONTINUE
99       W(KK)=B
100    15 GO TO (27,24),IINV
101    24 I=I+1
102       IF (I-N) 2,2,25
103    25 IINV=1
104       FF=0.
105       KL=NN
106       DO 26 I=1,M
107       KL=KL+1
108       F(I)=W(KL)
109       FF=FF+F(I)*F(I)
110    26 CONTINUE
111       ICONT = 1
112       ISS=1
113       MC=N+1
114       IPP=IPRINT*(IPRINT-1)
115       ITC=0
116       IPS=1
117       IPC=0
118    27 IPC=IPC-IPRINT
119       IF (IPC) 28,29,29
120    28 WRITE(6,30)ITC,MC,FF
121    30 FORMAT(//5X,'ITERATION',I4,I9,' CALLS OF FCN   ',5X,'F=',E24.14)
122       WRITE(6,31)(X(I),I=1,N)
123    31 FORMAT(5X,'VARIABLES',/(5E24.14))
124       WRITE(6,32)(F(I),I=1,M)
125    32 FORMAT(5X,'FUNCTIONS',/(5E24.14))
126       IPC=IPP
127       GO TO (29,33),IPS
128    29 GO TO (34,35),ICONT
129    35 IF (CHANGE-1.) 10,10,36
130    10 CONTINUE
131 C*UL   37 WRITE(6,38)
132       WRITE(6,38)
133    38 FORMAT(//5X,'      FINAL VALUES OF FUNCTIONS AND VARIABLES')
134       IPS=2
135       GO TO 28
136    33 IFLAG=3
137       CALL FCN1(M,N,F,X,IFLAG,NW,W,KEND,KEND2N)
138       IF(COV.EQ.0.) RETURN
139       IWC1 = 2*N+N*(N+M)
140       IWC=IWC1
141       DO 91 I=1,N
142       DO 91 J=1,I
143       IWC=IWC+1
144       W(IWC)=0.
145       DO 90 MR=1,N
146       JCOV=N+J+MR*(N+M)
147       ICOV=N+I+MR*(N+M)
148    90 W(IWC)=W(IWC)+W(ICOV)*W(JCOV)
149    91 CONTINUE
150       WRITE(6,1000)
151  1000 FORMAT(///40X,'VARIANCE-COVARIANCE MATRIX'///)
152       IB=IWC1
153       DO 92 I=1,N
154       IA=IB+1
155       IB=IB+I
156    92 WRITE(6,1001)(W(J),J=IA,IB)
157  1001 FORMAT(10X,5G20.6//)
158       RETURN
159    36 ICONT=1
160    34 ITC=ITC+1
161       K=N
162       KK=KST
163       DO 39 I=1,N
164       K=K+1
165       W(K)=0.
166       KK=KK+N
167       W(I)=0.
168       DO 40 J=1,M
169       KK=KK+1
170       W(I)=W(I)+W(KK)*F(J)
171    40 CONTINUE
172    39 CONTINUE
173       DM=0.
174       K=KINV
175       DO 41 II=1,N
176       IIP=II+N
177       W(IIP)=W(IIP)+W(K)*W(II)
178       JL=II+1
179       IF (JL-N) 42,42,43
180    42 DO 44 JJ=JL,N
181       JJP=JJ+N
182       K=K+1
183       W(IIP)=W(IIP)+W(K)*W(JJ)
184       W(JJP)=W(JJP)+W(K)*W(II)
185    44 CONTINUE
186       K=K+1
187    43 IF (DM-ABS (W(II)*W(IIP))) 45,41,41
188    45 DM=ABS (W(II)*W(IIP))
189       KL=II
190    41 CONTINUE
191       II=N+MPLUSN*KL
192       CHANGE=0.
193       DO 46 I=1,N
194       JL=N+I
195       W(I)=0.
196       DO 47 J=NPLUS,NN
197       JL=JL+MPLUSN
198       W(I)=W(I)+W(J)*W(JL)
199    47 CONTINUE
200       II=II+1
201       W(II)=W(JL)
202       W(JL)=X(I)
203       IF (ABS (E(I)*CHANGE)-ABS (W(I))) 48,48,46
204    48 CHANGE=ABS (W(I)/E(I))
205    46 CONTINUE
206       DO 49 I=1,M
207       II=II+1
208       JL=JL+1
209       W(II)=W(JL)
210       W(JL)=F(I)
211    49 CONTINUE
212       FC=FF
213       ACC=0.1/CHANGE
214       IT=3
215       XC=0.
216       XL=0.
217       IS=3
218       IF (CHANGE-1.) 50,50,51
219    50 ICONT=2
220    51 CALL VD01A (IT,XC,FC,20,ACC,0.1,XSTEP)
221       GO TO (52,53,53,53),IT
222    52 MC=MC+1
223       IF (MC-MAXFUN) 54,54,55
224    55 WRITE(6,56)MAXFUN
225    56 FORMAT(5X,I6,' CALLS OF FCN')
226       ISS=2
227       GO TO 53
228    54 XL=XC-XL
229       DO 57 J=1,N
230       X(J)=X(J)+XL*W(J)
231    57 CONTINUE
232       XL=XC
233       CALL FCN1(M,N,F,X,IFLAG,NW,W,KEND,KEND2N)
234       FC=0.
235       DO 58 J=1,M
236       FC=FC+F(J)*F(J)
237    58 CONTINUE
238       GO TO (59,59,60),IS
239    60 K=N
240       IF (FC-FF) 61,51,62
241    61 IS=2
242       FMIN=FC
243       FSEC=FF
244       GO TO 63
245    62 IS=1
246       FMIN=FF
247       FSEC=FC
248       GO TO 63
249    59 IF (FC-FSEC) 64,51,51
250    64 K=KSTORE
251       GO TO (75,74),IS
252    75 K=N
253    74 IF (FC-FMIN) 65,51,66
254    66 FSEC=FC
255       GO TO 63
256    65 IS=3-IS
257       FSEC=FMIN
258       FMIN=FC
259    63 DO 67 J=1,N
260       K=K+1
261       W(K)=X(J)
262    67 CONTINUE
263       DO 68 J=1,M
264       K=K+1
265       W(K)=F(J)
266    68 CONTINUE
267       GO TO 51
268    53 K=KSTORE
269       KK=N
270       GO TO (69,70,69),IS
271    70 K=N
272       KK=KSTORE
273    69 SUM=0.
274       DM=0.
275       JJ=KSTORE
276       DO 71 J=1,N
277       K=K+1
278       KK=KK+1
279       JJ=JJ+1
280       X(J)=W(K)
281       W(JJ)=W(K)-W(KK)
282    71 CONTINUE
283       DO 72 J=1,M
284       K=K+1
285       KK=KK+1
286       JJ=JJ+1
287       F(J)=W(K)
288       W(JJ)=W(K)-W(KK)
289       SUM=SUM+W(JJ)*W(JJ)
290       DM=DM+F(J)*W(JJ)
291    72 CONTINUE
292       GO TO (73,10),ISS
293    73 J=KINV
294       KK=NPLUS-KL
295       DO 76 I=1,KL
296       K=J+KL-I
297       J=K+KK
298       W(I)=W(K)
299       W(K)=W(J-1)
300    76 CONTINUE
301       IF (KL-N) 77,78,78
302    77 KL=KL+1
303       JJ=K
304       DO 79 I=KL,N
305       K=K+1
306       J=J+NPLUS-I
307       W(I)=W(K)
308       W(K)=W(J-1)
309    79 CONTINUE
310       W(JJ)=W(K)
311       B=1./W(KL-1)
312       W(KL-1)=W(N)
313       GO TO 88
314    78 B=1./W(N)
315    88 K=KINV
316       DO 80 I=1,ILESS
317       BB=B*W(I)
318       DO 81 J=I,ILESS
319       W(K)=W(K)-BB*W(J)
320       K=K+1
321    81 CONTINUE
322       K=K+1
323    80 CONTINUE
324       IF (FMIN-FF) 82,83,83
325    83 CHANGE=0.
326       GO TO 84
327    82 FF=FMIN
328       CHANGE=ABS (XC)*CHANGE
329    84 XL=-DM/FMIN
330       SUM=1./SQRT (SUM+DM*XL)
331       K=KSTORE
332       DO 85 I=1,N
333       K=K+1
334       W(K)=SUM*W(K)
335       W(I)=0.
336    85 CONTINUE
337       DO 86 I=1,M
338       K=K+1
339       W(K)=SUM*(W(K)+XL*F(I))
340       KK=NN+I
341       DO 87 J=1,N
342       KK=KK+MPLUSN
343       W(J)=W(J)+W(KK)*W(K)
344    87 CONTINUE
345    86 CONTINUE
346       GO TO 14
347       END