samer@0
|
1 /*
|
samer@0
|
2 * Copyright (c) 2000, Samer Abdallah, King's College London.
|
samer@0
|
3 * All rights reserved.
|
samer@0
|
4 *
|
samer@0
|
5 * This software is provided AS iS and WITHOUT ANY WARRANTY;
|
samer@0
|
6 * without even the implied warranty of MERCHANTABILITY or
|
samer@0
|
7 * FITNESS FOR A PARTICULAR PURPOSE.
|
samer@0
|
8 */
|
samer@0
|
9
|
samer@0
|
10 package samer.maths.opt;
|
samer@0
|
11 import samer.maths.*;
|
samer@0
|
12 import samer.core.*;
|
samer@0
|
13
|
samer@0
|
14 public class ConstrainedGillMurray extends GillMurray
|
samer@0
|
15 {
|
samer@0
|
16 State S;
|
samer@0
|
17 Constraints C;
|
samer@0
|
18
|
samer@0
|
19 public ConstrainedGillMurray(State s, Constraints c) {
|
samer@0
|
20 super(s); S=s; C=c;
|
samer@0
|
21 }
|
samer@0
|
22
|
samer@0
|
23 public void computeHg(double [] g)
|
samer@0
|
24 {
|
samer@0
|
25 C.zeroInactive(S.h);
|
samer@0
|
26 C.mul(S.h,H,g);
|
samer@0
|
27 C.negate(S.h);
|
samer@0
|
28 S.normh = Util.maxabs(S.h);
|
samer@0
|
29 }
|
samer@0
|
30
|
samer@0
|
31 public void steepest() {
|
samer@0
|
32 C.negate(S.P1.g,S.h);
|
samer@0
|
33 C.zeroInactive(S.h);
|
samer@0
|
34 S.normh = Util.maxabs(S.h);
|
samer@0
|
35 }
|
samer@0
|
36
|
samer@0
|
37 public void resetDimension(int k) {
|
samer@0
|
38 double [] hk=H0.getArray()[k];
|
samer@0
|
39 for (int i=0; i<S.n; i++) H[i][k]=H[k][i]=hk[i];
|
samer@0
|
40 }
|
samer@0
|
41
|
samer@0
|
42 /* calculate new hessian and search direction */
|
samer@0
|
43 public void update()
|
samer@0
|
44 {
|
samer@0
|
45 // sparse subtraction and dot product
|
samer@0
|
46 C.sub(u,S.P2.x,S.P1.x);
|
samer@0
|
47 C.sub(v,S.P2.g,S.P1.g);
|
samer@0
|
48 double uv = C.dot(u,v);
|
samer@0
|
49
|
samer@0
|
50 if (uvCheck(uv,C.dot(u,u),C.dot(v,v))) {
|
samer@0
|
51 C.mul(Hv,H,v);
|
samer@0
|
52
|
samer@0
|
53 double a = 1/uv;
|
samer@0
|
54 double b = a*(1+a*C.dot(v,Hv));
|
samer@0
|
55
|
samer@0
|
56 for (int k=0; k<C.m; k++) {
|
samer@0
|
57 int i=C.active[k];
|
samer@0
|
58 for (int l=0; l<k; l++) {
|
samer@0
|
59 int j=C.active[l];
|
samer@0
|
60 H[i][j] = H[i][j] - a*(Hv[i]*u[j]+u[i]*Hv[j]) + b*u[i]*u[j];
|
samer@0
|
61 H[j][i] = H[i][j];
|
samer@0
|
62 }
|
samer@0
|
63 H[i][i] = H[i][i] - a*(2*Hv[i]*u[i]) + b*u[i]*u[i];
|
samer@0
|
64 }
|
samer@0
|
65 }
|
samer@0
|
66 computeHg(S.P2.g);
|
samer@0
|
67 S.setSlope2(C.dot(S.P2.g,S.h));
|
samer@0
|
68 }
|
samer@0
|
69 }
|