1************************************************************************
2* SUBROUTINE PLIPU              ALL SYSTEMS                   97/01/22
3* PURPOSE :
4* EASY TO USE SUBROUTINE FOR LARGE-SCALE UNCONSTRAINED MINIMIZATION.
5*
6* PARAMETERS :
7*  II  NF  NUMBER OF VARIABLES.
8*  RI  X(NF)  VECTOR OF VARIABLES.
9*  II  IPAR(7)  INTEGER PAREMETERS:
10*      IPAR(1)  MAXIMUM NUMBER OF ITERATIONS.
11*      IPAR(2)  MAXIMUM NUMBER OF FUNCTION EVALUATIONS.
12*      IPAR(3)  THIS PARAMETER IS NOT USED IN THE SUBROUTINE PLIP.
13*      IPAR(4)  ESTIMATION INDICATOR. IPAR(4)=0-MINIMUM IS NOT
14*         ESTIMATED. IPAR(4)=1-MINIMUM IS ESTIMATED BY THE VALUE
15*         RPAR(6).
16*      IPAR(5)  METHOD USED. IPAR(5)=1-RANK-ONE METHOD.
17*         IPAR(5)=2-RANK-TWO METHOD.
18*      IPAR(6)  THIS PARAMETER IS NOT USED IN THE SUBROUTINE PLIP.
19*      IPAR(7)  MAXIMUM NUMBER OF VARIABLE METRIC UPDATES.
20*  RI  RPAR(9)  REAL PARAMETERS:
21*      RPAR(1)  MAXIMUM STEPSIZE.
22*      RPAR(2)  TOLERANCE FOR THE CHANGE OF VARIABLES.
23*      RPAR(3)  TOLERANCE FOR THE CHANGE OF FUNCTION VALUES.
24*      RPAR(4)  TOLERANCE FOR THE FUNCTION FALUE.
25*      RPAR(5)  TOLERANCE FOR THE GRADIENT NORM.
26*      RPAR(6)  ESTIMATION OF THE MINIMUM FUNCTION VALUE.
27*      RPAR(7)  THIS PARAMETER IS NOT USED IN THE SUBROUTINE PLIP.
28*      RPAR(8)  THIS PARAMETER IS NOT USED IN THE SUBROUTINE PLIP.
29*      RPAR(9)  THIS PARAMETER IS NOT USED IN THE SUBROUTINE PLIP.
30*  RO  F  VALUE OF THE OBJECTIVE FUNCTION.
31*  RO  GMAX  MAXIMUM PARTIAL DERIVATIVE.
32*  II  IPRNT  PRINT SPECIFICATION. IPRNT=0-NO PRINT.
33*         ABS(IPRNT)=1-PRINT OF FINAL RESULTS.
34*         ABS(IPRNT)=2-PRINT OF FINAL RESULTS AND ITERATIONS.
35*         IPRNT>0-BASIC FINAL RESULTS. IPRNT<0-EXTENDED FINAL
36*         RESULTS.
37*  IO  ITERM  VARIABLE THAT INDICATES THE CAUSE OF TERMINATION.
38*         ITERM=1-IF ABS(X-XO) WAS LESS THAN OR EQUAL TO TOLX IN
39*                   MTESX (USUALLY TWO) SUBSEQUENT ITERATIONS.
40*         ITERM=2-IF ABS(F-FO) WAS LESS THAN OR EQUAL TO TOLF IN
41*                   MTESF (USUALLY TWO) SUBSEQUENT ITERATIONS.
42*         ITERM=3-IF F IS LESS THAN OR EQUAL TO TOLB.
43*         ITERM=4-IF GMAX IS LESS THAN OR EQUAL TO TOLG.
44*         ITERM=6-IF THE TERMINATION CRITERION WAS NOT SATISFIED,
45*                   BUT THE SOLUTION OBTAINED IS PROBABLY ACCEPTABLE.
46*         ITERM=11-IF NIT EXCEEDED MIT. ITERM=12-IF NFV EXCEEDED MFV.
47*         ITERM=13-IF NFG EXCEEDED MFG. ITERM<0-IF THE METHOD FAILED.
48*
49* VARIABLES IN COMMON /STAT/ (STATISTICS) :
50*  IO  NRES  NUMBER OF RESTARTS.
51*  IO  NDEC  NUMBER OF MATRIX DECOMPOSITIONS.
52*  IO  NIN  NUMBER OF INNER ITERATIONS.
53*  IO  NIT  NUMBER OF ITERATIONS.
54*  IO  NFV  NUMBER OF FUNCTION EVALUATIONS.
55*  IO  NFG  NUMBER OF GRADIENT EVALUATIONS.
56*  IO  NFH  NUMBER OF HESSIAN EVALUATIONS.
57*
58* SUBPROGRAMS USED :
59*  S   PLIP  LIMITED MEMORY SHIFTED VARIABLE METRIC METHOD IN THE
60*            PRODUCT FORM.
61*
62* EXTERNAL SUBROUTINES :
63*  SE  OBJ  COMPUTATION OF THE VALUE OF THE OBJECTIVE FUNCTION.
64*         CALLING SEQUENCE: CALL OBJ(NF,X,FF) WHERE NF IS THE NUMBER
65*         OF VARIABLES, X(NF) IS THE VECTOR OF VARIABLES AND FF IS
66*         THE VALUE OF THE OBJECTIVE FUNCTION.
67*  SE  DOBJ  COMPUTATION OF THE GRADIENT OF THE OBJECTIVE FUNCTION.
68*         CALLING SEQUENCE: CALL DOBJ(NF,X,GF) WHERE NF IS THE NUMBER
69*         OF VARIABLES, X(NF) IS THE VECTOR OF VARIABLES AND GF(NF)
70*         IS THE GRADIENT OF THE OBJECTIVE FUNCTION.
71*
72      SUBROUTINE PLIPU(NF,X,IPAR,RPAR,F,GMAX,IPRNT,ITERM)
73      INTEGER NF,IPAR(7),IPRNT,ITERM
74      DOUBLE PRECISION X(*),RPAR(9),F,GMAX
75      INTEGER MF,NB,LGF,LS,LXO,LGO,LSO,LXM,LXR,LGR
76      INTEGER NRES,NDEC,NIN,NIT,NFV,NFG,NFH
77      COMMON /STAT/ NRES,NDEC,NIN,NIT,NFV,NFG,NFH
78      DOUBLE PRECISION RA(:)
79      ALLOCATABLE RA
80      MF=IPAR(7)
81      IF (MF.LE.0) MF=10
82      ALLOCATE (RA(5*NF+NF*MF+2*MF))
83      NB=0
84*
85*     POINTERS FOR AUXILIARY ARRAYS
86*
87      LGF=1
88      LS=LGF+NF
89      LXO=LS+NF
90      LGO=LXO+NF
91      LSO=LGO+NF
92      LXM=LSO+NF
93      LXR=LXM+NF*MF
94      LGR=LXR+MF
95      CALL PLIP(NF,NB,X,IPAR,RA,RA,RA(LGF),RA(LS),RA(LXO),RA(LGO),
96     & RA(LSO),RA(LXM),RA(LXR),RA(LGR),RPAR(1),RPAR(2),RPAR(3),RPAR(4),
97     & RPAR(5),RPAR(6),GMAX,F,IPAR(1),IPAR(2),IPAR(4),IPAR(5),MF,IPRNT,
98     & ITERM)
99      DEALLOCATE (RA)
100      RETURN
101      END
102************************************************************************
103* SUBROUTINE PLIPS              ALL SYSTEMS                   97/01/22
104* PURPOSE :
105* EASY TO USE SUBROUTINE FOR LARGE-SCALE BOX CONSTRAINED MINIMIZATION.
106*
107* PARAMETERS :
108*  II  NF  NUMBER OF VARIABLES.
109*  RI  X(NF)  VECTOR OF VARIABLES.
110*  II  IX(NF)  VECTOR CONTAINING TYPES OF BOUNDS. IX(I)=0-VARIABLE
111*         X(I) IS UNBOUNDED. IX(I)=1-LOWER BOUND XL(I).LE.X(I).
112*         IX(I)=2-UPPER BOUND X(I).LE.XU(I). IX(I)=3-TWO SIDE BOUND
113*         XL(I).LE.X(I).LE.XU(I). IX(I)=5-VARIABLE X(I) IS FIXED.
114*  RI  XL(NF)  VECTOR CONTAINING LOWER BOUNDS FOR VARIABLES.
115*  RI  XU(NF)  VECTOR CONTAINING UPPER BOUNDS FOR VARIABLES.
116*  II  IPAR(7)  INTEGER PAREMETERS:
117*      IPAR(1)  MAXIMUM NUMBER OF ITERATIONS.
118*      IPAR(2)  MAXIMUM NUMBER OF FUNCTION EVALUATIONS.
119*      IPAR(3)  THIS PARAMETER IS NOT USED IN THE SUBROUTINE PLIP.
120*      IPAR(4)  ESTIMATION INDICATOR. IPAR(4)=0-MINIMUM IS NOT
121*         ESTIMATED. IPAR(4)=1-MINIMUM IS ESTIMATED BY THE VALUE
122*         RPAR(6).
123*      IPAR(5)  METHOD USED. IPAR(5)=1-RANK-ONE METHOD.
124*         IPAR(5)=2-RANK-TWO METHOD.
125*      IPAR(6)  THIS PARAMETER IS NOT USED IN THE SUBROUTINE PLIP.
126*      IPAR(7)  MAXIMUM NUMBER OF VARIABLE METRIC UPDATES.
127*  RI  RPAR(9)  REAL PARAMETERS:
128*      RPAR(1)  MAXIMUM STEPSIZE.
129*      RPAR(2)  TOLERANCE FOR THE CHANGE OF VARIABLES.
130*      RPAR(3)  TOLERANCE FOR THE CHANGE OF FUNCTION VALUES.
131*      RPAR(4)  TOLERANCE FOR THE FUNCTION FALUE.
132*      RPAR(5)  TOLERANCE FOR THE GRADIENT NORM.
133*      RPAR(6)  ESTIMATION OF THE MINIMUM FUNCTION VALUE.
134*      RPAR(7)  THIS PARAMETER IS NOT USED IN THE SUBROUTINE PLIP.
135*      RPAR(8)  THIS PARAMETER IS NOT USED IN THE SUBROUTINE PLIP.
136*      RPAR(9)  THIS PARAMETER IS NOT USED IN THE SUBROUTINE PLIP.
137*  RO  F  VALUE OF THE OBJECTIVE FUNCTION.
138*  RO  GMAX  MAXIMUM PARTIAL DERIVATIVE.
139*  II  IPRNT  PRINT SPECIFICATION. IPRNT=0-NO PRINT.
140*         ABS(IPRNT)=1-PRINT OF FINAL RESULTS.
141*         ABS(IPRNT)=2-PRINT OF FINAL RESULTS AND ITERATIONS.
142*         IPRNT>0-BASIC FINAL RESULTS. IPRNT<0-EXTENDED FINAL
143*         RESULTS.
144*  IO  ITERM  VARIABLE THAT INDICATES THE CAUSE OF TERMINATION.
145*         ITERM=1-IF ABS(X-XO) WAS LESS THAN OR EQUAL TO TOLX IN
146*                   MTESX (USUALLY TWO) SUBSEQUENT ITERATIONS.
147*         ITERM=2-IF ABS(F-FO) WAS LESS THAN OR EQUAL TO TOLF IN
148*                   MTESF (USUALLY TWO) SUBSEQUENT ITERATIONS.
149*         ITERM=3-IF F IS LESS THAN OR EQUAL TO TOLB.
150*         ITERM=4-IF GMAX IS LESS THAN OR EQUAL TO TOLG.
151*         ITERM=6-IF THE TERMINATION CRITERION WAS NOT SATISFIED,
152*                   BUT THE SOLUTION OBTAINED IS PROBABLY ACCEPTABLE.
153*         ITERM=11-IF NIT EXCEEDED MIT. ITERM=12-IF NFV EXCEEDED MFV.
154*         ITERM=13-IF NFG EXCEEDED MFG. ITERM<0-IF THE METHOD FAILED.
155*
156* VARIABLES IN COMMON /STAT/ (STATISTICS) :
157*  IO  NRES  NUMBER OF RESTARTS.
158*  IO  NDEC  NUMBER OF MATRIX DECOMPOSITIONS.
159*  IO  NIN  NUMBER OF INNER ITERATIONS.
160*  IO  NIT  NUMBER OF ITERATIONS.
161*  IO  NFV  NUMBER OF FUNCTION EVALUATIONS.
162*  IO  NFG  NUMBER OF GRADIENT EVALUATIONS.
163*  IO  NFH  NUMBER OF HESSIAN EVALUATIONS.
164*
165* SUBPROGRAMS USED :
166*  S   PLIP  LIMITED MEMORY SHIFTED VARIABLE METRIC METHOD IN THE
167*            PRODUCT FORM.
168*
169* EXTERNAL SUBROUTINES :
170*  SE  OBJ  COMPUTATION OF THE VALUE OF THE OBJECTIVE FUNCTION.
171*         CALLING SEQUENCE: CALL OBJ(NF,X,FF) WHERE NF IS THE NUMBER
172*         OF VARIABLES, X(NF) IS THE VECTOR OF VARIABLES AND FF IS
173*         THE VALUE OF THE OBJECTIVE FUNCTION.
174*  SE  DOBJ  COMPUTATION OF THE GRADIENT OF THE OBJECTIVE FUNCTION.
175*         CALLING SEQUENCE: CALL DOBJ(NF,X,GF) WHERE NF IS THE NUMBER
176*         OF VARIABLES, X(NF) IS THE VECTOR OF VARIABLES AND GF(NF)
177*         IS THE GRADIENT OF THE OBJECTIVE FUNCTION.
178*
179      SUBROUTINE PLIPS(NF,X,IX,XL,XU,IPAR,RPAR,F,GMAX,IPRNT,ITERM)
180      INTEGER NF,IX(*),IPAR(7),IPRNT,ITERM
181      DOUBLE PRECISION X(*),XL(*),XU(*),RPAR(9),F,GMAX
182      INTEGER MF,NB,LGF,LS,LXO,LGO,LSO,LXM,LXR,LGR
183      INTEGER NRES,NDEC,NIN,NIT,NFV,NFG,NFH
184      COMMON /STAT/ NRES,NDEC,NIN,NIT,NFV,NFG,NFH
185      DOUBLE PRECISION RA(:)
186      ALLOCATABLE RA
187      MF=IPAR(7)
188      IF (MF.LE.0) MF=10
189      ALLOCATE (RA(5*NF+NF*MF+2*MF))
190      NB=1
191*
192*     POINTERS FOR AUXILIARY ARRAYS
193*
194      LGF=1
195      LS=LGF+NF
196      LXO=LS+NF
197      LGO=LXO+NF
198      LSO=LGO+NF
199      LXM=LSO+NF
200      LXR=LXM+NF*MF
201      LGR=LXR+MF
202      CALL PLIP(NF,NB,X,IX,XL,XU,RA(LGF),RA(LS),RA(LXO),RA(LGO),
203     & RA(LSO),RA(LXM),RA(LXR),RA(LGR),RPAR(1),RPAR(2),RPAR(3),RPAR(4),
204     & RPAR(5),RPAR(6),GMAX,F,IPAR(1),IPAR(2),IPAR(4),IPAR(5),MF,IPRNT,
205     & ITERM)
206      DEALLOCATE (RA)
207      RETURN
208      END
209************************************************************************
210* SUBROUTINE PLIP               ALL SYSTEMS                   01/09/22
211* PURPOSE :
212* GENERAL SUBROUTINE FOR LARGE-SCALE BOX CONSTRAINED MINIMIZATION THAT
213* USE THE SHIFTED LIMITED MEMORY VARIABLE METRIC METHOD BASED ON THE
214* PRODUCT FORM UPDATES.
215*
216* PARAMETERS :
217*  II  NF  NUMBER OF VARIABLES.
218*  II  NB  CHOICE OF SIMPLE BOUNDS. NB=0-SIMPLE BOUNDS SUPPRESSED.
219*         NB>0-SIMPLE BOUNDS ACCEPTED.
220*  RI  X(NF)  VECTOR OF VARIABLES.
221*  II  IX(NF)  VECTOR CONTAINING TYPES OF BOUNDS. IX(I)=0-VARIABLE
222*         X(I) IS UNBOUNDED. IX(I)=1-LOVER BOUND XL(I).LE.X(I).
223*         IX(I)=2-UPPER BOUND X(I).LE.XU(I). IX(I)=3-TWO SIDE BOUND
224*         XL(I).LE.X(I).LE.XU(I). IX(I)=5-VARIABLE X(I) IS FIXED.
225*  RI  XL(NF)  VECTOR CONTAINING LOWER BOUNDS FOR VARIABLES.
226*  RI  XU(NF)  VECTOR CONTAINING UPPER BOUNDS FOR VARIABLES.
227*  RA  GF(NF)  GRADIENT OF THE OBJECTIVE FUNCTION.
228*  RO  S(NF)  DIRECTION VECTOR.
229*  RU  XO(NF)  VECTORS OF VARIABLES DIFFERENCE.
230*  RI  GO(NF)  GRADIENTS DIFFERENCE.
231*  RA  SO(NF)  AUXILIARY VECTOR.
232*  RA  XM(NF*MF)  AUXILIARY VECTOR.
233*  RA  XR(MF)  AUXILIARY VECTOR.
234*  RA  GR(MF)  AUXILIARY VECTOR.
235*  RI  XMAX  MAXIMUM STEPSIZE.
236*  RI  TOLX  TOLERANCE FOR CHANGE OF VARIABLES.
237*  RI  TOLF  TOLERANCE FOR CHANGE OF FUNCTION VALUES.
238*  RI  TOLB  TOLERANCE FOR THE FUNCTION VALUE.
239*  RI  TOLG  TOLERANCE FOR THE GRADIENT NORM.
240*  RI  FMIN  ESTIMATION OF THE MINIMUM FUNCTION VALUE.
241*  RO  GMAX  MAXIMUM PARTIAL DERIVATIVE.
242*  RO  F  VALUE OF THE OBJECTIVE FUNCTION.
243*  II  MIT  MAXIMUM NUMBER OF ITERATIONS.
244*  II  MFV  MAXIMUM NUMBER OF FUNCTION EVALUATIONS.
245*  II  IEST  ESTIMATION INDICATOR. IEST=0-MINIMUM IS NOT ESTIMATED.
246*         IEST=1-MINIMUM IS ESTIMATED BY THE VALUE FMIN.
247*  II  MET  METHOD USED. MET=1-RANK-ONE METHOD. MET=2-RANK-TWO
248*         METHOD.
249*  II  MF  NUMBER OF LIMITED MEMORY STEPS.
250*  II  IPRNT  PRINT SPECIFICATION. IPRNT=0-NO PRINT.
251*         ABS(IPRNT)=1-PRINT OF FINAL RESULTS.
252*         ABS(IPRNT)=2-PRINT OF FINAL RESULTS AND ITERATIONS.
253*         IPRNT>0-BASIC FINAL RESULTS. IPRNT<0-EXTENDED FINAL
254*         RESULTS.
255*  IO  ITERM  VARIABLE THAT INDICATES THE CAUSE OF TERMINATION.
256*         ITERM=1-IF ABS(X-XO) WAS LESS THAN OR EQUAL TO TOLX IN
257*                   MTESX (USUALLY TWO) SUBSEQUEBT ITERATIONS.
258*         ITERM=2-IF ABS(F-FO) WAS LESS THAN OR EQUAL TO TOLF IN
259*                   MTESF (USUALLY TWO) SUBSEQUEBT ITERATIONS.
260*         ITERM=3-IF F IS LESS THAN OR EQUAL TO TOLB.
261*         ITERM=4-IF GMAX IS LESS THAN OR EQUAL TO TOLG.
262*         ITERM=6-IF THE TERMINATION CRITERION WAS NOT SATISFIED,
263*                   BUT THE SOLUTION OBTAINED IS PROBABLY ACCEPTABLE.
264*         ITERM=11-IF NIT EXCEEDED MIT. ITERM=12-IF NFV EXCEEDED MFV.
265*         ITERM=13-IF NFG EXCEEDED MFG. ITERM<0-IF THE METHOD FAILED.
266*
267* VARIABLES IN COMMON /STAT/ (STATISTICS) :
268*  IO  NRES  NUMBER OF RESTARTS.
269*  IO  NDEC  NUMBER OF MATRIX DECOMPOSITION.
270*  IO  NIN  NUMBER OF INNER ITERATIONS.
271*  IO  NIT  NUMBER OF ITERATIONS.
272*  IO  NFV  NUMBER OF FUNCTION EVALUATIONS.
273*  IO  NFG  NUMBER OF GRADIENT EVALUATIONS.
274*  IO  NFH  NUMBER OF HESSIAN EVALUATIONS.
275*
276* SUBPROGRAMS USED :
277*  S   PCBS04  ELIMINATION OF BOX CONSTRAINT VIOLATIONS.
278*  S   PS1L01  STEPSIZE SELECTION USING LINE SEARCH.
279*  S   PULSP3  SHIFTED VARIABLE METRIC UPDATE.
280*  S   PULVP3  SHIFTED LIMITED-MEMORY VARIABLE METRIC UPDATE.
281*  S   PYADC0  ADDITION OF A BOX CONSTRAINT.
282*  S   PYFUT1  TEST ON TERMINATION.
283*  S   PYRMC0  DELETION OF A BOX CONSTRAINT.
284*  S   PYTRCD  COMPUTATION OF PROJECTED DIFFERENCES FOR THE VARIABLE METRIC
285*         UPDATE.
286*  S   PYTRCG  COMPUTATION OF THE PROJECTED GRADIENT.
287*  S   PYTRCS  COMPUTATION OF THE PROJECTED DIRECTION VECTOR.
288*  S   MXDRMM  MULTIPLICATION OF A ROWWISE STORED DENSE RECTANGULAR
289*         MATRIX A BY A VECTOR X.
290*  S   MXDCMD  MULTIPLICATION OF A COLUMNWISE STORED DENSE RECTANGULAR
291*         MATRIX A BY A VECTOR X AND ADDITION OF THE SCALED VECTOR
292*         ALF*Y.
293*  S   MXUCOP  COPYING OF A VECTOR.
294*  S   MXUDIR  VECTOR AUGMENTED BY THE SCALED VECTOR.
295*  RF  MXUDOT  DOT PRODUCT OF TWO VECTORS.
296*  S   MXUNEG  COPYING OF A VECTOR WITH CHANGE OF THE SIGN.
297*  S   MXUZER  VECTOR ELEMENTS CORRESPONDING TO ACTIVE BOUNDS ARE SET
298*         TO ZERO.
299*  S   MXVCOP  COPYING OF A VECTOR.
300*
301* EXTERNAL SUBROUTINES :
302*  SE  OBJ  COMPUTATION OF THE VALUE OF THE OBJECTIVE FUNCTION.
303*         CALLING SEQUENCE: CALL OBJ(NF,X,FF) WHERE NF IS THE NUMBER
304*         OF VARIABLES, X(NF) IS THE VECTOR OF VARIABLES AND FF IS
305*         THE VALUE OF THE OBJECTIVE FUNCTION.
306*  SE  DOBJ  COMPUTATION OF THE GRADIENT OF THE OBJECTIVE FUNCTION.
307*         CALLING SEQUENCE: CALL DOBJ(NF,X,GF) WHERE NF IS THE NUMBER
308*         OF VARIABLES, X(NF) IS THE VECTOR OF VARIABLES AND GF(NF)
309*         IS THE GRADIENT OF THE OBJECTIVE FUNCTION.
310*
311* METHOD :
312* HYBRID METHOD WITH SPARSE MARWIL UPDATES FOR SPARSE LEAST SQUARES
313* PROBLEMS.
314*
315      SUBROUTINE PLIP(NF,NB,X,IX,XL,XU,GF,S,XO,GO,SO,XM,XR,GR,XMAX,
316     & TOLX,TOLF,TOLB,TOLG,FMIN,GMAX,F,MIT,MFV,IEST,MET,MF,IPRNT,ITERM)
317      INTEGER NF,NB,IX(*),MIT,MFV,IEST,MET,MF,IPRNT,ITERM
318      DOUBLE PRECISION X(*),XL(*),XU(*),GF(*),S(*),XO(*),GO(*),SO(*),
319     & XM(*),XR(*),GR(*),XMAX,TOLX,TOLF,TOLG,TOLB,FMIN,GMAX,F
320      INTEGER ITERD,ITERS,ITERH,KD,LD,NTESX,NTESF,MTESX,MTESF,MRED,KIT,
321     & IREST,KBF,MEC,MES,MES1,MES2,MES3,MAXST,ISYS,ITES,INITS,KTERS,
322     & IRES1,IRES2,NRED,INEW,IOLD,I,NN,N,MFG,META,MET3
323      DOUBLE PRECISION R,RO,RP,FO,FP,P,PO,PP,GNORM,SNORM,RMIN,RMAX,
324     & UMAX,FMAX,DMAX,ETA0,ETA9,EPS8,EPS9,ALF1,ALF2,PAR1,PAR2,PAR,TOLD,
325     & TOLS,TOLP
326      DOUBLE PRECISION MXUDOT
327      INTEGER NRES,NDEC,NIN,NIT,NFV,NFG,NFH
328      COMMON /STAT/ NRES,NDEC,NIN,NIT,NFV,NFG,NFH
329      IF (ABS(IPRNT).GT.1) WRITE(6,'(1X,''ENTRY TO PLIP :'')')
330*
331*     INITIATION
332*
333      KBF=0
334      IF (NB.GT.0) KBF=2
335      NRES=0
336      NDEC=0
337      NIN=0
338      NIT=0
339      NFV=0
340      NFG=0
341      NFH=0
342      ISYS=0
343      ITES=1
344      MTESX=2
345      MTESF=2
346      INITS=2
347      ITERM=0
348      ITERD=0
349      ITERS=2
350      ITERH=0
351      KTERS=3
352      IREST=0
353      IRES1=999
354      IRES2=0
355      MRED=10
356      META=1
357      MET3=4
358      MEC=4
359      MES=4
360      MES1=2
361      MES2=2
362      MES3=2
363      ETA0=1.0D-15
364      ETA9=1.0D 120
365      EPS8=1.00D 0
366      EPS9=1.00D-8
367      ALF1=1.0D-10
368      ALF2=1.0D 10
369      RMAX=ETA9
370      DMAX=ETA9
371      FMAX=1.0D 20
372      IF (IEST.LE.0) FMIN=-1.0D 60
373      IF (IEST.GT.0) IEST=1
374      IF (XMAX.LE.0.0D 0) XMAX=1.0D 16
375      IF (TOLX.LE.0.0D 0) TOLX=1.0D-16
376      IF (TOLF.LE.0.0D 0) TOLF=1.0D-14
377      IF (TOLG.LE.0.0D 0) TOLG=1.0D-6
378      IF (TOLB.LE.0.0D 0) TOLB=FMIN+1.0D-16
379      TOLD=1.0D-4
380      TOLS=1.0D-4
381      TOLP=0.9D 0
382      IF (MET.LE.0) MET=2
383      IF (MIT.LE.0) MIT=9000
384      IF (MFV.LE.0) MFV=9000
385      MFG=MFV
386      KD= 1
387      LD=-1
388      KIT=-(IRES1*NF+IRES2)
389      FO=FMIN
390*
391*     INITIAL OPERATIONS WITH SIMPLE BOUNDS
392*
393      IF (KBF.GT.0) THEN
394      DO 2 I = 1,NF
395      IF ((IX(I).EQ.3.OR.IX(I).EQ.4) .AND. XU(I).LE.XL(I)) THEN
396      XU(I) = XL(I)
397      IX(I) = 5
398      ELSE IF (IX(I).EQ.5 .OR. IX(I).EQ.6) THEN
399      XL(I) = X(I)
400      XU(I) = X(I)
401      IX(I) = 5
402      END IF
403    2 CONTINUE
404      CALL PCBS04(NF,X,IX,XL,XU,EPS9,KBF)
405      CALL PYADC0(NF,N,X,IX,XL,XU,INEW)
406      END IF
407      IF (ITERM.NE.0) GO TO 11190
408      CALL OBJ(NF,X,F)
409      NFV=NFV+1
410      CALL DOBJ(NF,X,GF)
411      NFG=NFG+1
41211120 CONTINUE
413      CALL PYTRCG(NF,NF,IX,GF,UMAX,GMAX,KBF,IOLD)
414      IF (ABS(IPRNT).GT.1)
415     & WRITE (6,'(1X,''NIT='',I5,2X,''NFV='',I5,2X,''NFG='',I5,2X,
416     & ''F='', G16.9,2X,''G='',E10.3)') NIT,NFV,NFG,F,GMAX
417      CALL PYFUT1(NF,F,FO,UMAX,GMAX,DMAX,TOLX,TOLF,TOLB,TOLG,KD,
418     & NIT,KIT,MIT,NFV,MFV,NFG,MFG,NTESX,MTESX,NTESF,MTESF,ITES,
419     & IRES1,IRES2,IREST,ITERS,ITERM)
420      IF (ITERM.NE.0) GO TO 11190
421      IF (KBF.GT.0.AND.RMAX.GT.0.0D 0) THEN
422      CALL PYRMC0(NF,N,IX,GF,EPS8,UMAX,GMAX,RMAX,IOLD,IREST)
423      END IF
42411130 CONTINUE
425      IF (IREST.GT.0) THEN
426      NN=0
427      PAR=1.0D 0
428      LD=MIN(LD,1)
429      IF (KIT.LT.NIT) THEN
430        NRES=NRES+1
431        KIT = NIT
432      ELSE
433        ITERM=-10
434        IF (ITERS.LT.0) ITERM=ITERS-5
435      END IF
436      END IF
437      IF (ITERM.NE.0) GO TO 11190
438*
439*     DIRECTION DETERMINATION
440*
441      GNORM=SQRT(MXUDOT(NF,GF,GF,IX,KBF))
442*
443*     NEWTON LIKE STEP
444*
445      CALL MXUNEG(NF,GF,S,IX,KBF)
446      CALL MXDRMM(NF,NN,XM,S,GR)
447      CALL MXDCMD(NF,NN,XM,GR,PAR,S,S)
448      CALL MXUZER(NF,S,IX,KBF)
449      ITERD=1
450      SNORM=SQRT(MXUDOT(NF,S,S,IX,KBF))
451*
452*     TEST ON DESCENT DIRECTION AND PREPARATION OF LINE SEARCH
453*
454      IF (KD.GT.0) P=MXUDOT(NF,GF,S,IX,KBF)
455      IF (ITERD.LT.0) THEN
456        ITERM=ITERD
457      ELSE
458*
459*     TEST ON DESCENT DIRECTION
460*
461      IF (SNORM.LE.0.0D 0) THEN
462        IREST=MAX(IREST,1)
463      ELSE IF (P+TOLD*GNORM*SNORM.LE.0.0D 0) THEN
464        IREST=0
465      ELSE
466*
467*     UNIFORM DESCENT CRITERION
468*
469      IREST=MAX(IREST,1)
470      END IF
471      IF (IREST.EQ.0) THEN
472*
473*     PREPARATION OF LINE SEARCH
474*
475        NRED = 0
476        RMIN=ALF1*GNORM/SNORM
477        RMAX=MIN(ALF2*GNORM/SNORM,XMAX/SNORM)
478      END IF
479      END IF
480      IF (ITERM.NE.0) GO TO 11190
481      IF (IREST.NE.0) GO TO 11130
482      CALL PYTRCS(NF,X,IX,XO,XL,XU,GF,GO,S,RO,FP,FO,F,PO,P,RMAX,ETA9,
483     & KBF)
484      IF (RMAX.EQ.0.0D 0) GO TO 11175
48511170 CONTINUE
486      CALL PS1L01(R,RP,F,FO,FP,P,PO,PP,FMIN,FMAX,RMIN,RMAX,
487     & TOLS,TOLP,PAR1,PAR2,KD,LD,NIT,KIT,NRED,MRED,MAXST,IEST,
488     & INITS,ITERS,KTERS,MES,ISYS)
489      IF (ISYS.EQ.0) GO TO 11174
490      CALL MXUDIR(NF,R,S,XO,X,IX,KBF)
491      CALL PCBS04(NF,X,IX,XL,XU,EPS9,KBF)
492      CALL OBJ(NF,X,F)
493      NFV=NFV+1
494      CALL DOBJ(NF,X,GF)
495      NFG=NFG+1
496      P=MXUDOT(NF,GF,S,IX,KBF)
497      GO TO 11170
49811174 CONTINUE
499      IF (ITERS.LE.0) THEN
500      R=0.0D 0
501      F=FO
502      P=PO
503      CALL MXVCOP(NF,XO,X)
504      CALL MXVCOP(NF,GO,GF)
505      IREST=MAX(IREST,1)
506      LD=KD
507      GO TO 11130
508      END IF
509      CALL MXUNEG(NF,GO,S,IX,KBF)
510      CALL PYTRCD(NF,X,IX,XO,GF,GO,R,F,FO,P,PO,DMAX,KBF,KD,LD,ITERS)
511      CALL MXUCOP(NF,GF,SO,IX,KBF)
512      IF (NN.LT.MF) THEN
513      CALL PULSP3(NF,NN,MF,XM,GR,XO,GO,R,PO,PAR,ITERH,MET3)
514      ELSE
515      CALL PULVP3(NF,NN,XM,XR,GR,S,SO,XO,GO,R,PO,PAR,ITERH,MEC,MET3,
516     & MET)
517      END IF
51811175 CONTINUE
519      IF (ITERH.NE.0) IREST=MAX(IREST,1)
520      IF (KBF.GT.0) CALL PYADC0(NF,N,X,IX,XL,XU,INEW)
521      GO TO 11120
52211190 CONTINUE
523      IF (IPRNT.GT.1.OR.IPRNT.LT.0)
524     & WRITE(6,'(1X,''EXIT FROM PLIP :'')')
525      IF (IPRNT.NE.0)
526     & WRITE (6,'(1X,''NIT='',I5,2X,''NFV='',I5,2X,''NFG='',I5,2X,
527     & ''F='', G16.9,2X,''G='',E10.3,2X,''ITERM='',I3)') NIT,NFV,NFG,
528     & F,GMAX,ITERM
529      IF (IPRNT.LT.0)
530     & WRITE (6,'(1X,''X='',5(G14.7,1X):/(3X,5(G14.7,1X)))')
531     & (X(I),I=1,NF)
532      RETURN
533      END
534