|
2329
|
1 *+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ |
|
|
2 |
|
|
3 SUBROUTINE LSGSET( PRBTYP, LINOBJ, SINGLR, UNITGZ, UNITQ, |
|
|
4 $ N, NCLIN, NFREE, |
|
|
5 $ NROWA, NQ, NROWR, NRANK, NZ, NZ1, |
|
|
6 $ ISTATE, KX, |
|
|
7 $ BIGBND, TOLRNK, NUMINF, SUMINF, |
|
|
8 $ BL, BU, A, RES, FEATOL, |
|
|
9 $ GQ, CQ, R, X, WTINF, ZY, WRK ) |
|
|
10 |
|
|
11 IMPLICIT DOUBLE PRECISION(A-H,O-Z) |
|
|
12 CHARACTER*2 PRBTYP |
|
|
13 LOGICAL LINOBJ, SINGLR, UNITGZ, UNITQ |
|
|
14 INTEGER ISTATE(*), KX(N) |
|
|
15 DOUBLE PRECISION BL(*), BU(*), A(NROWA,*), |
|
|
16 $ RES(*), FEATOL(*) |
|
|
17 DOUBLE PRECISION GQ(N), CQ(*), R(NROWR,*), X(N), WTINF(*), |
|
|
18 $ ZY(NQ,*) |
|
|
19 DOUBLE PRECISION WRK(N) |
|
|
20 |
|
|
21 ************************************************************************ |
|
|
22 * LSGSET finds the number and weighted sum of infeasibilities for |
|
|
23 * the bounds and linear constraints. An appropriate transformed |
|
|
24 * gradient vector is returned in GQ. |
|
|
25 * |
|
|
26 * Positive values of ISTATE(j) will not be altered. These mean |
|
|
27 * the following... |
|
|
28 * |
|
|
29 * 1 2 3 |
|
|
30 * a'x = bl a'x = bu bl = bu |
|
|
31 * |
|
|
32 * Other values of ISTATE(j) will be reset as follows... |
|
|
33 * a'x lt bl a'x gt bu a'x free |
|
|
34 * - 2 - 1 0 |
|
|
35 * |
|
|
36 * If x is feasible, LSGSET computes the vector Q(free)'g(free), |
|
|
37 * where g is the gradient of the the sum of squares plus the |
|
|
38 * linear term. The matrix Q is of the form |
|
|
39 * ( Q(free) 0 ), |
|
|
40 * ( 0 I(fixed)) |
|
|
41 * where Q(free) is the orthogonal factor of A(free) and A is |
|
|
42 * the matrix of constraints in the working set. The transformed |
|
|
43 * gradients are stored in GQ. |
|
|
44 * |
|
|
45 * Systems Optimization Laboratory, Stanford University. |
|
|
46 * Original version written 31-October-1984. |
|
|
47 * Level 2 Blas added 11-June-1986. |
|
|
48 * This version of LSGSET dated 24-June-1986. |
|
|
49 ************************************************************************ |
|
|
50 EXTERNAL DDOT , IDRANK |
|
|
51 INTRINSIC ABS , MAX , MIN |
|
|
52 PARAMETER ( ZERO = 0.0D+0, HALF = 0.5D+0, ONE = 1.0D+0 ) |
|
|
53 |
|
|
54 BIGUPP = BIGBND |
|
|
55 BIGLOW = - BIGBND |
|
|
56 |
|
|
57 NUMINF = 0 |
|
|
58 SUMINF = ZERO |
|
|
59 CALL DLOAD ( N, ZERO, GQ, 1 ) |
|
|
60 |
|
|
61 DO 200 J = 1, N+NCLIN |
|
|
62 IF (ISTATE(J) .LE. 0) THEN |
|
|
63 FEASJ = FEATOL(J) |
|
|
64 IF (J .LE. N) THEN |
|
|
65 CTX = X(J) |
|
|
66 ELSE |
|
|
67 K = J - N |
|
|
68 CTX = DDOT ( N, A(K,1), NROWA, X, 1 ) |
|
|
69 END IF |
|
|
70 ISTATE(J) = 0 |
|
|
71 |
|
|
72 * See if the lower bound is violated. |
|
|
73 |
|
|
74 IF (BL(J) .GT. BIGLOW) THEN |
|
|
75 S = BL(J) - CTX |
|
|
76 IF (S .GT. FEASJ ) THEN |
|
|
77 ISTATE(J) = - 2 |
|
|
78 WEIGHT = - WTINF(J) |
|
|
79 GO TO 160 |
|
|
80 END IF |
|
|
81 END IF |
|
|
82 |
|
|
83 * See if the upper bound is violated. |
|
|
84 |
|
|
85 IF (BU(J) .GE. BIGUPP) GO TO 200 |
|
|
86 S = CTX - BU(J) |
|
|
87 IF (S .LE. FEASJ ) GO TO 200 |
|
|
88 ISTATE(J) = - 1 |
|
|
89 WEIGHT = WTINF(J) |
|
|
90 |
|
|
91 * Add the infeasibility. |
|
|
92 |
|
|
93 160 NUMINF = NUMINF + 1 |
|
|
94 SUMINF = SUMINF + ABS( WEIGHT ) * S |
|
|
95 IF (J .LE. N) THEN |
|
|
96 GQ(J) = WEIGHT |
|
|
97 ELSE |
|
|
98 CALL DAXPY ( N, WEIGHT, A(K,1), NROWA, GQ, 1 ) |
|
|
99 END IF |
|
|
100 END IF |
|
|
101 200 CONTINUE |
|
|
102 |
|
|
103 * ------------------------------------------------------------------ |
|
|
104 * Install GQ, the transformed gradient. |
|
|
105 * ------------------------------------------------------------------ |
|
|
106 SINGLR = .FALSE. |
|
|
107 UNITGZ = .TRUE. |
|
|
108 |
|
|
109 IF (NUMINF .GT. 0) THEN |
|
|
110 CALL CMQMUL( 6, N, NZ, NFREE, NQ, UNITQ, KX, GQ, ZY, WRK ) |
|
|
111 ELSE IF (NUMINF .EQ. 0 .AND. PRBTYP .EQ. 'FP') THEN |
|
|
112 CALL DLOAD ( N, ZERO, GQ, 1 ) |
|
|
113 ELSE |
|
|
114 |
|
|
115 * Ready for the Optimality Phase. |
|
|
116 * Set NZ1 so that Rz1 is nonsingular. |
|
|
117 |
|
|
118 IF (NRANK .EQ. 0) THEN |
|
|
119 IF (LINOBJ) THEN |
|
|
120 CALL DCOPY ( N, CQ, 1, GQ, 1 ) |
|
|
121 ELSE |
|
|
122 CALL DLOAD ( N, ZERO, GQ, 1 ) |
|
|
123 END IF |
|
|
124 NZ1 = 0 |
|
|
125 ELSE |
|
|
126 |
|
|
127 * Compute GQ = - R' * (transformed residual) |
|
|
128 |
|
|
129 CALL DCOPY ( NRANK, RES, 1, GQ, 1 ) |
|
|
130 CALL DSCAL ( NRANK, (-ONE), GQ, 1 ) |
|
|
131 CALL DTRMV ( 'U', 'T', 'N', NRANK, R, NROWR, GQ, 1 ) |
|
|
132 IF (NRANK .LT. N) |
|
|
133 $ CALL DGEMV( 'T', NRANK, N-NRANK, -ONE,R(1,NRANK+1),NROWR, |
|
|
134 $ RES, 1, ZERO, GQ(NRANK+1), 1 ) |
|
|
135 |
|
|
136 IF (LINOBJ) CALL DAXPY ( N, ONE, CQ, 1, GQ, 1 ) |
|
|
137 UNITGZ = .FALSE. |
|
|
138 NZ1 = IDRANK( MIN(NRANK, NZ), R, NROWR+1, TOLRNK ) |
|
|
139 END IF |
|
|
140 END IF |
|
|
141 |
|
|
142 RETURN |
|
|
143 |
|
|
144 * End of LSGSET. |
|
|
145 |
|
|
146 END |
