OptimizationProblems.jl/src/ADNLPProblems/clnlbeam.jl at f60efc25b69ca8de7bdb877bbb15221665379559 · JuliaSmoothOptimizers/OptimizationProblems.jl · GitHub

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
export clnlbeam

function clnlbeam(args...; n::Int = default_nvar, type::Type{T} = Float64, kwargs...) where {T}
  n_orig = n
  N = div(n - 3, 3)
  n = 3 * N + 3
  @adjust_nvar_warn("clnlbeam", n_orig, n)
  h = 1 // N
  alpha = 350
  function f(y; N = N, h = h, alpha = alpha)
    @views t, x, u = y[1:(N + 1)], y[(N + 2):(2 * N + 2)], y[(2 * N + 3):end]
    return sum(
      1 // 2 * h * (u[i + 1]^2 + u[i]^2) + 1 // 2 * alpha * h * (cos(t[i + 1]) + cos(t[i])) for
      i = 1:N
    )
  end
  function c!(cx, y; N = N)
    @views t, x, u = y[1:(N + 1)], y[(N + 2):(2 * N + 2)], y[(2 * N + 3):end]
    for i = 1:N
      cx[i] = x[i + 1] - x[i] - 1 // 2 * h * (sin(t[i + 1]) + sin(t[i]))
    end
    return cx
  end
  lvar = vcat(-ones(T, N + 1), -T(0.05) * ones(T, N + 1), -T(Inf) * ones(T, N + 1))
  uvar = vcat(ones(T, N + 1), T(0.05) * ones(T, N + 1), T(Inf) * ones(T, N + 1))
  xi = zeros(T, 3 * N + 3)

  # [t[i + 1] - t[i] - T(0.5) * h * u[i + 1] - T(0.5) * h * u[i] for i = 1:N]
  clinrows = vcat(
    1:N, # t[i + 1]
    1:N,
    1:N,
    1:N,
  )
  clincols = vcat(
    [i + 1 for i = 1:N], # t[i + 1]
    [i for i = 1:N], # t[i]
    [2 * N + 2 + i + 1 for i = 1:N], # u[i + 1]
    [2 * N + 2 + i for i = 1:N], # u[i]
  )
  clinvals = vcat(
    T[1 for i = 1:N], # t[i + 1]
    T[-1 for i = 1:N], # t[i]
    [-T(0.5) * h for i = 1:N], # u[i + 1]
    [-T(0.5) * h for i = 1:N], # u[i]
  )
  return ADNLPModels.ADNLPModel!(
    f,
    xi,
    lvar,
    uvar,
    clinrows,
    clincols,
    clinvals,
    c!,
    zeros(T, 2 * N),
    zeros(T, 2 * N),
    name = "clnlbeam";
    kwargs...,
  )
end