-
Notifications
You must be signed in to change notification settings - Fork 0
/
fixed-point.jl
77 lines (54 loc) · 2.11 KB
/
fixed-point.jl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
# Stopping can also be used for fixed point methods
# Example here concerns the AlternatingDirections Algorithm to find
# a feasible point in the intersection of 2 convex sets A and B.
# This algorithm relies on a fixed point argument, hence it stopped if it finds
# a fixed point.
# Example:
# A={ (x,y) | x=y} and B = {(x,y) | y=0}
# Clearly the unique intersection point is (0,0)
# Note that in this case the projection on A and the projection on B are trivial
# Takeaway: the 2nd scenario illustrates a situation where the algorithm stalls
# as it reached a personal success, i.e. `optimal_sub_pb` is true.
using ADNLPModels, LinearAlgebra, NLPModels, Stopping, Test
# Main algorithm
function AlternatingDirections(stp)
xk = stp.current_state.x
OK = update_and_start!(stp, cx = cons(stp.pb, x0))
while !OK
#First projection
xk1 = 0.5 * (xk[1] + xk[2]) * ones(2)
#Second projection
xk2 = [xk1[1],0.0]
#check if we have a fixed point
Fix = dot(xk-xk2,xk-xk2)
if Fix <= min(eps(Float64),stp.meta.atol) stp.meta.suboptimal = true end
#call the stopping
OK = update_and_stop!(stp, x = xk2, cx = cons(stp.pb, xk2))
xk = xk2
end
return stp
end
# We model the problem using the NLPModels without objective function
# Formulate the problem with NLPModels
c(x) = [x[1] - x[2], x[2]]
lcon = [0.0, 0.0]
ucon = [0.0, 0.0]
nlp = ADNLPModel(x->0.0, zeros(2), c, lcon, ucon)
# 1st scenario: we solve the problem
printstyled("1st scenario:\n")
#Prepare the Stopping
x0 = [0.0, 5.0]
state = NLPAtX(x0)
# Recall that for the optimality_check function x is the pb and y is the state
# Here we take the infinite norm of the residual.
stop = NLPStopping(nlp, state, optimality_check = (x,y) -> norm(y.cx, Inf))
AlternatingDirections(stop)
@test status(stop) == :Optimal
# 2nd scenario: the user gives an irrealistic optimality condition
printstyled("2nd scenario:\n")
reinit!(stop, rstate = true, x = x0)
stop.meta.optimality_check = (x, y) -> norm(y.cx, Inf) + 0.5
AlternatingDirections(stop)
# In this scenario, the algorithm stops because it attains a fixed point
# Hence, status is :SubOptimal.
@test status(stop) == :SubOptimal