Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 0 additions & 3 deletions Project.toml
Original file line number Diff line number Diff line change
Expand Up @@ -14,9 +14,6 @@ ExaModels = "1037b233-b668-4ce9-9b63-f9f681f55dd2"
COPSBenchmarkExaModels = "ExaModels"
COPSBenchmarkJuMP = "JuMP"

[sources]
ExaModels = {url = "https://github.com/exanauts/ExaModels.jl", rev = "main"}

[compat]
JuMP = "^1.19"
julia = "1.11"
42 changes: 21 additions & 21 deletions ext/COPSBenchmarkExaModels/catmix.jl
Original file line number Diff line number Diff line change
Expand Up @@ -7,65 +7,67 @@
ne = 2
nc = 3

tf = 1
h = tf / nh # Final time
h = T(1) / T(nh) # Final time / nh (was Int/Int → Float64)

rho = [
rho = T[
0.11270166537926,
0.50000000000000,
0.88729833462074,
]
bc = [1.0, 0.0] # Boundary conditions for x
alpha = 0.0 # Smoothing parameter
bc = T[1, 0] # Boundary conditions for x
alpha = T(0) # Smoothing parameter
neg_one = -T(1)
ten_T = T(10)
one_T = T(1)
rho_index = [(i, rho[i]) for i in 1:nc]

c = ExaModels.ExaCore(T; backend = backend, concrete = Val(true))
ExaModels.@add_var(c, u, nh, nc; lvar = zeros(nh, nc), uvar = ones(nh, nc), start = zeros(nh, nc))
ExaModels.@add_var(c, v, nh, ne; start = [mod(j, ne) for i in 1:nh, j in 1:ne])
ExaModels.@add_var(c, w, nh, nc, ne; start = zeros(nh, nc, ne))
ExaModels.@add_var(c, pp, nh, nc, ne; start = [mod(k, ne) for i in 1:nh, j in 1:nc, k in 1:ne])
ExaModels.@add_var(c, Dpp, nh, nc, ne; start = zeros(nh, nc, ne))
ExaModels.@add_var(c, ppf, ne; start = [mod(i,ne) for i in 1:ne])

ExaModels.@add_obj(c, -1.0 + ppf[1] + ppf[2])
ExaModels.@add_var(c, u, nh, nc; lvar = zeros(T, nh, nc), uvar = ones(T, nh, nc), start = zeros(T, nh, nc))
ExaModels.@add_var(c, v, nh, ne; start = T[mod(j, ne) for i in 1:nh, j in 1:ne])
ExaModels.@add_var(c, w, nh, nc, ne; start = zeros(T, nh, nc, ne))
ExaModels.@add_var(c, pp, nh, nc, ne; start = T[mod(k, ne) for i in 1:nh, j in 1:nc, k in 1:ne])
ExaModels.@add_var(c, Dpp, nh, nc, ne; start = zeros(T, nh, nc, ne))
ExaModels.@add_var(c, ppf, ne; start = T[mod(i,ne) for i in 1:ne])

ExaModels.@add_obj(c, neg_one + ppf[1] + ppf[2])
ExaModels.@add_obj(c, alpha/h*(u[i+1, j] - u[i, j])^2 for i in 1:nh-1, j in 1:nc)

ExaModels.@add_con(
c,
c1,
pp[i, k, s] - v[i, s] - h*sum(w[i, j, s]*(rho^j/factorial(j)) for j in 1:nc) for i=1:nh, (k, rho) in rho_index, s=1:ne
pp[i, k, s] - v[i, s] - h*sum(w[i, j, s]*(rho^j/T(factorial(j))) for j in 1:nc) for i=1:nh, (k, rho) in rho_index, s=1:ne
)

ExaModels.@add_con(
c,
c2,
Dpp[i, k, s] - sum(w[i, j, s]*(rho^(j-1)/factorial(j-1)) for j in 1:nc) for i=1:nh, (k, rho) in rho_index, s=1:ne
Dpp[i, k, s] - sum(w[i, j, s]*(rho^(j-1)/T(factorial(j-1))) for j in 1:nc) for i=1:nh, (k, rho) in rho_index, s=1:ne
)

ExaModels.@add_con(
c,
c3,
ppf[s] - v[nh, s] - h * sum(w[nh, j, s] / factorial(j) for j in 1:nc) for s in 1:ne
ppf[s] - v[nh, s] - h * sum(w[nh, j, s] / T(factorial(j)) for j in 1:nc) for s in 1:ne
)

ExaModels.@add_con(
c,
c4,
v[i, s] + sum(w[i, j, s] * h / factorial(j) for j in 1:nc) - v[i+1, s] for i in 1:nh-1, s in 1:ne
v[i, s] + sum(w[i, j, s] * h / T(factorial(j)) for j in 1:nc) - v[i+1, s] for i in 1:nh-1, s in 1:ne
)



ExaModels.@add_con(
c,
c5,
Dpp[i,j,1] - u[i,j] * (10.0*pp[i,j,2] - pp[i,j,1]) for i=1:nh, j=1:nc
Dpp[i,j,1] - u[i,j] * (ten_T*pp[i,j,2] - pp[i,j,1]) for i=1:nh, j=1:nc
)

ExaModels.@add_con(
c,
c6,
Dpp[i,j,2] - u[i,j] * (pp[i,j,1] - 10.0*pp[i,j,2]) + (1 - u[i,j])*pp[i,j,2] for i=1:nh, j=1:nc
Dpp[i,j,2] - u[i,j] * (pp[i,j,1] - ten_T*pp[i,j,2]) + (one_T - u[i,j])*pp[i,j,2] for i=1:nh, j=1:nc
)


Expand All @@ -77,5 +79,3 @@

return ExaModels.ExaModel(c; kwargs...)
end


30 changes: 15 additions & 15 deletions ext/COPSBenchmarkExaModels/gasoil.jl
Original file line number Diff line number Diff line change
Expand Up @@ -12,18 +12,19 @@
nm = 21 # number of measurements

# roots of k-th degree Legendre polynomial
rho = [0.06943184420297, 0.33000947820757, 0.66999052179243, 0.93056815579703]
rho = T[0.06943184420297, 0.33000947820757, 0.66999052179243, 0.93056815579703]
# ODE initial conditions
bc = [1, 1, 2, 0]
# times at which observations made
tau = [0.0, 0.025, 0.05, 0.075, 0.10, 0.125, 0.150, 0.175, 0.20, 0.225, 0.250, 0.30, 0.35, 0.40, 0.45, 0.50, 0.55, 0.65, 0.75, 0.85, 0.95]
tau = T[0.0, 0.025, 0.05, 0.075, 0.10, 0.125, 0.150, 0.175, 0.20, 0.225, 0.250, 0.30, 0.35, 0.40, 0.45, 0.50, 0.55, 0.65, 0.75, 0.85, 0.95]
# ODEs defined in [0,tf]
tf = tau[nm]
# uniform interval length
h = tf / nh
t = [(i-1)*h for i in 1:nh+1]
h = tf / T(nh)
t = T[T(i-1)*h for i in 1:nh+1]
zero_T = T(0)

itau = Int[min(nh, floor(tau[i]/h)+1) for i in 1:nm]
itau = Int[min(nh, Int(floor(tau[i]/h))+1) for i in 1:nm]

# Concentrations
z = reshape(T[
Expand All @@ -50,10 +51,10 @@
0.0690, 0.0100,
], ne, nm)'

v0 = zeros(nh, ne)
v0 = zeros(T, nh, ne)
# Starting-value
for i in 1:itau[1], s in 1:ne
v0[i, s] = bc[s]
v0[i, s] = T(bc[s])
end
for j in 2:nm, i =itau[j-1]+1:itau[j], s in 1:ne
v0[i, s] = z[j, s]
Expand All @@ -65,31 +66,31 @@
core = ExaModels.ExaCore(T; backend= backend, concrete = Val(true))

# ODE parameters
ExaModels.@add_var(core, theta, 1:np; lvar = 0.0, start=0.0)
ExaModels.@add_var(core, theta, 1:np; lvar = zero_T, start=zero_T)
# The collocation approximation u is defined by the parameters v and w.
# uc and Duc are, respectively, u and u' evaluated at the collocation points.
ExaModels.@add_var(core, v, 1:nh, 1:ne; start=[v0[i, s] for i =1:nh, s = 1:ne])
ExaModels.@add_var(core, w, 1:nh, 1:nc, 1:ne; start=0.0)
ExaModels.@add_var(core, w, 1:nh, 1:nc, 1:ne; start=zero_T)
ExaModels.@add_var(core, uc, 1:nh, 1:nc, 1:ne; start=[v0[i, s] for i =1:nh, j=1:nc, s = 1:ne])
ExaModels.@add_var(core, Duc, 1:nh, 1:nc, 1:ne; start=0.0)
ExaModels.@add_var(core, Duc, 1:nh, 1:nc, 1:ne; start=zero_T)

itr = [(j, s, itau[j], tau[j], t[itau[j]], z[j,s]) for j=1:nm, s in 1:ne]
itr2 = [(j,rho[j]) for j=1:nc]

# L2 error
ExaModels.@add_obj(core, (v[itauj,s] + sum(w[itauj,k,s]*(tauj-tj)^k/(factorial(k)*h^(k-1)) for k in 1:nc) - zjs)^2 for (j,s,itauj, tauj, tj, zjs) in itr)
ExaModels.@add_obj(core, (v[itauj,s] + sum(w[itauj,k,s]*(tauj-tj)^k/(T(factorial(k))*h^(k-1)) for k in 1:nc) - zjs)^2 for (j,s,itauj, tauj, tj, zjs) in itr)

# Collocation model
ExaModels.@add_con(
core,
c1,
- uc[i, j, s] + v[i,s] + h*sum(w[i,k,s]*(rhoj^k/factorial(k)) for k in 1:nc)
- uc[i, j, s] + v[i,s] + h*sum(w[i,k,s]*(rhoj^k/T(factorial(k))) for k in 1:nc)
for i=1:nh, (j,rhoj) in itr2, s=1:ne
)
ExaModels.@add_con(
core,
c2,
- Duc[i, j, s] + sum(w[i,k,s]*(rhoj^(k-1)/factorial(k-1)) for k in 1:nc)
- Duc[i, j, s] + sum(w[i,k,s]*(rhoj^(k-1)/T(factorial(k-1))) for k in 1:nc)
for i=1:nh, (j,rhoj) in itr2, s=1:ne
)

Expand All @@ -100,7 +101,7 @@
ExaModels.@add_con(
core,
c4,
v[i, s] + sum(w[i, j, s]*h/factorial(j) for j in 1:nc) - v[i+1, s]
v[i, s] + sum(w[i, j, s]*h/T(factorial(j)) for j in 1:nc) - v[i+1, s]
for i=1:nh-1, s=1:ne
)
ExaModels.@add_con(
Expand All @@ -117,4 +118,3 @@
)
return ExaModels.ExaModel(core; kwargs...)
end

66 changes: 35 additions & 31 deletions ext/COPSBenchmarkExaModels/glider.jl
Original file line number Diff line number Diff line change
Expand Up @@ -6,52 +6,56 @@
# COPS 3.1 - March 2004

@inline function COPSBenchmark.glider_model(::ExaModelsBackend, nh; T = Float64, backend = nothing, kwargs...)
# Design parameters
x_0 = 0.0
y_0 = 1000.0
y_f = 900.0
vx_0 = 13.23
vx_f = 13.23
vy_0 = -1.288
vy_f = -1.288
u_c = 2.5
r_0 = 100.0
m = 100.0
g = 9.81
cd0 = 0.034
cd1 = 0.069662
S = 14.0
rho = 1.13
cL_min = 0.0
cL_max = 1.4
cL0 = cL_max / 2
# Design parameters (T-typed)
x_0 = T(0)
y_0 = T(1000)
y_f = T(900)
vx_0 = T(13.23)
vx_f = T(13.23)
vy_0 = T(-1.288)
vy_f = T(-1.288)
u_c = T(2.5)
r_0 = T(100)
m = T(100)
g = T(9.81)
cd0 = T(0.034)
cd1 = T(0.069662)
S = T(14)
rho = T(1.13)
cL_min = T(0)
cL_max = T(1.4)
cL0 = cL_max / T(2)
half = T(0.5)
inv_nh = T(1) / T(nh)
r_offset = T(2.5)
one_T = T(1)

c = ExaModels.ExaCore(T; backend = backend, concrete = Val(true))

ExaModels.@add_var(c, t_f, 1; lvar = 0.0, start = 1.0)
ExaModels.@add_var(c, x, nh+1; lvar = zeros(nh+1), start = [x_0 + vx_0*(k/nh) for k in 0:nh])
ExaModels.@add_var(c, y, nh+1; start = [y_0 + (k/nh)*(y_f - y_0) for k in 0:nh])
ExaModels.@add_var(c, vx, nh+1; lvar = zeros(nh+1), start = fill(vx_0, nh+1))
ExaModels.@add_var(c, t_f, 1; lvar = T(0), start = T(1))
ExaModels.@add_var(c, x, nh+1; lvar = zeros(T, nh+1), start = [x_0 + vx_0*(T(k)*inv_nh) for k in 0:nh])
ExaModels.@add_var(c, y, nh+1; start = [y_0 + (T(k)*inv_nh)*(y_f - y_0) for k in 0:nh])
ExaModels.@add_var(c, vx, nh+1; lvar = zeros(T, nh+1), start = fill(vx_0, nh+1))
ExaModels.@add_var(c, vy, nh+1; start = fill(vy_0, nh+1))
ExaModels.@add_var(c, cL, nh+1; lvar = fill(cL_min,nh+1), uvar = fill(cL_max, nh+1), start = fill(cL0, nh+1))

# Expressions (matching JuMP @expressions)
ExaModels.@add_expr(c, r, (x[i]/r_0 - 2.5)^2 for i in 1:nh+1)
ExaModels.@add_expr(c, u, u_c*(1 - r[i])*exp(-r[i]) for i in 1:nh+1)
ExaModels.@add_expr(c, r, (x[i]/r_0 - r_offset)^2 for i in 1:nh+1)
ExaModels.@add_expr(c, u, u_c*(one_T - r[i])*exp(-r[i]) for i in 1:nh+1)
ExaModels.@add_expr(c, w, vy[i] - u[i] for i in 1:nh+1)
ExaModels.@add_expr(c, v, sqrt(vx[i]^2 + w[i]^2) for i in 1:nh+1)
ExaModels.@add_expr(c, D, 0.5*(cd0+cd1*cL[i]^2)*rho*S*v[i]^2 for i in 1:nh+1)
ExaModels.@add_expr(c, L, 0.5*cL[i]*rho*S*v[i]^2 for i in 1:nh+1)
ExaModels.@add_expr(c, D, half*(cd0+cd1*cL[i]^2)*rho*S*v[i]^2 for i in 1:nh+1)
ExaModels.@add_expr(c, L, half*cL[i]*rho*S*v[i]^2 for i in 1:nh+1)
ExaModels.@add_expr(c, vx_dot, (-L[i]*(w[i]/v[i]) - D[i]*(vx[i]/v[i]))/m for i in 1:nh+1)
ExaModels.@add_expr(c, vy_dot, (L[i]*(vx[i]/v[i]) - D[i]*(w[i]/v[i]))/m - g for i in 1:nh+1)

ExaModels.@add_obj(c, -x[nh+1])

# Dynamics
ExaModels.@add_con(c, c1, x[j] - (x[j-1] + 0.5 * t_f[1]/nh * (vx[j] + vx[j-1])) for j in 2:nh+1)
ExaModels.@add_con(c, c2, y[j] - (y[j-1] + 0.5 * t_f[1]/nh * (vy[j] + vy[j-1])) for j in 2:nh+1)
ExaModels.@add_con(c, c3, vx[j] - (vx[j-1] + 0.5 * t_f[1]/nh * (vx_dot[j] + vx_dot[j-1])) for j in 2:nh+1)
ExaModels.@add_con(c, c4, vy[j] - (vy[j-1] + 0.5 * t_f[1]/nh * (vy_dot[j] + vy_dot[j-1])) for j in 2:nh+1)
ExaModels.@add_con(c, c1, x[j] - (x[j-1] + half * t_f[1]*inv_nh * (vx[j] + vx[j-1])) for j in 2:nh+1)
ExaModels.@add_con(c, c2, y[j] - (y[j-1] + half * t_f[1]*inv_nh * (vy[j] + vy[j-1])) for j in 2:nh+1)
ExaModels.@add_con(c, c3, vx[j] - (vx[j-1] + half * t_f[1]*inv_nh * (vx_dot[j] + vx_dot[j-1])) for j in 2:nh+1)
ExaModels.@add_con(c, c4, vy[j] - (vy[j-1] + half * t_f[1]*inv_nh * (vy_dot[j] + vy_dot[j-1])) for j in 2:nh+1)

# Boundary constraints
ExaModels.@add_con(c, c5, x[1] - x_0)
Expand Down
33 changes: 16 additions & 17 deletions ext/COPSBenchmarkExaModels/marine.jl
Original file line number Diff line number Diff line change
Expand Up @@ -10,14 +10,15 @@
ne = 8 # number of differential equations
nm = 21 # number of measurements

rho = [0.5] # roots of k-th degree Legendre polynomial
tau = collect(range(0.0, 10.0, 21)) # times at which observations made
tf = tau[nm] # ODEs defined in [0,tf]
h = tf / nh # uniform interval length
t = [(i-1)*h for i in 1:nh+1] # partition
rho = T[0.5] # roots of k-th degree Legendre polynomial
tau = collect(T, range(T(0), T(10), 21)) # times at which observations made
tf = tau[nm] # ODEs defined in [0,tf]
h = tf / T(nh) # uniform interval length
t = T[T(i-1)*h for i in 1:nh+1] # partition
zero_T = T(0)

# itau[i] is the largest integer k with t[k] <= tau[i]
itau = Int[min(nh, floor(tau[i]/h)+1) for i in 1:nm]
itau = Int[min(nh, Int(floor(tau[i]/h))+1) for i in 1:nm]

# Observation
z = reshape(T[
Expand All @@ -44,7 +45,7 @@
1.0, 12.0, 198.0, 707.0, 2562.0, 3163.0, 3232.0, 5566.0,
], ne, nm)'

v0 = zeros(nh, ne)
v0 = zeros(T, nh, ne)
# Starting-value
for i in 1:itau[1], s in 1:ne
v0[i, s] = z[1, s]
Expand All @@ -59,15 +60,15 @@
core = ExaModels.ExaCore(T; backend= backend, concrete = Val(true))

# Growth rates
ExaModels.@add_var(core, g, 1:ne-1; lvar = 0.0)
ExaModels.@add_var(core, g, 1:ne-1; lvar = zero_T)
# Mortality rates
ExaModels.@add_var(core, m, 1:ne; lvar = 0.0)
ExaModels.@add_var(core, m, 1:ne; lvar = zero_T)
# The collocation approximation u is defined by the parameters v and w.
# uc and Duc are, respectively, u and u' evaluated at the collocation points.
ExaModels.@add_var(core, v, 1:nh, 1:ne; start=[v0[i, s] for i=1:nh, s=1:ne])
ExaModels.@add_var(core, w, 1:nh, 1:nc, 1:ne; start=0.0)
ExaModels.@add_var(core, w, 1:nh, 1:nc, 1:ne; start=zero_T)
ExaModels.@add_var(core, uc, 1:nh, 1:nc, 1:ne; start=[v0[i, s] for i=1:nh, s=1:ne])
ExaModels.@add_var(core, Duc, 1:nh, 1:nc, 1:ne; start=0.0)
ExaModels.@add_var(core, Duc, 1:nh, 1:nc, 1:ne; start=zero_T)

# error

Expand All @@ -76,26 +77,26 @@
itr = [(j, s, itau[j], tau[j], t[itau[j]], z[j,s]) for j=1:nm, s in 1:ne]
itr2 = [(j,rho[j]) for j=1:nc]

ExaModels.@add_obj(core, (v[itauj,s] + sum(w[itauj,k,s]*(tauj-tj)^k/(factorial(k)*h^(k-1)) for k in 1:nc) - zjs)^2 for (j,s,itauj,tauj, tj, zjs) in itr)
ExaModels.@add_obj(core, (v[itauj,s] + sum(w[itauj,k,s]*(tauj-tj)^k/(T(factorial(k))*h^(k-1)) for k in 1:nc) - zjs)^2 for (j,s,itauj,tauj, tj, zjs) in itr)

# Collocation model
ExaModels.@add_con(
core,
c1,
- uc[i, j, s] + v[i,s] + h*sum(w[i,k,s]*(rhoj^k/factorial(k)) for k in 1:nc)
- uc[i, j, s] + v[i,s] + h*sum(w[i,k,s]*(rhoj^k/T(factorial(k))) for k in 1:nc)
for i=1:nh, (j,rhoj) in itr2, s=1:ne
)
ExaModels.@add_con(
core,
c2,
- Duc[i, j, s] + sum(w[i,k,s]*(rhoj^(k-1)/factorial(k-1)) for k in 1:nc)
- Duc[i, j, s] + sum(w[i,k,s]*(rhoj^(k-1)/T(factorial(k-1))) for k in 1:nc)
for i=1:nh, (j,rhoj) in itr2, s=1:ne
)
# Continuity
ExaModels.@add_con(
core,
c3,
v[i, s] + sum(w[i, j, s]*h/factorial(j) for j in 1:nc) - v[i+1, s]
v[i, s] + sum(w[i, j, s]*h/T(factorial(j)) for j in 1:nc) - v[i+1, s]
for i=1:nh-1, s=1:ne
)
# Boundary conditions
Expand All @@ -121,5 +122,3 @@

return ExaModels.ExaModel(core; kwargs...)
end


Loading
Loading