Skip to content

Commit

Permalink
Browse files Browse the repository at this point in the history
…nProblems into add-marine
  • Loading branch information
tmigot committed Oct 22, 2021
2 parents 6f6ebab + 795a1d0 commit ec8a501
Show file tree
Hide file tree
Showing 2 changed files with 44 additions and 4 deletions.
46 changes: 43 additions & 3 deletions src/marine.jl
Original file line number Diff line number Diff line change
Expand Up @@ -2,12 +2,28 @@ export marine

# Marine Population Dynamics COPS Problem v.0.3.1
# https://www.mcs.anl.gov/~more//cops/cops3.pdf
<<<<<<< HEAD
function marine(args...; n = 400, T = 10.0, kwargs...)
model = CartesianDiscreteModel((0, T), n)

valuetype = Float64
reffe = ReferenceFE(lagrangian, valuetype, 1)
Vfree = TestFESpace(model, reffe; conformity = :H1)
=======
# n=400, 800, 1600
function marine(args...; n = 400, kwargs...)
T = 10.0
model = CartesianDiscreteModel((0, T), n)

labels = get_face_labeling(model)
add_tag_from_tags!(labels, "diri0", [1]) #initial time condition

valuetype = Float64
reffe = ReferenceFE(lagrangian, valuetype, 1)
Vfree = TestFESpace(model, reffe; conformity = :H1)
VI = TestFESpace(model, reffe; conformity = :H1, labels = labels, dirichlet_tags = ["diri0"])
UI = TrialFESpace(VI, 0.0)
>>>>>>> 795a1d06d1b16edeafe3d23aef2a206e8ce94b7a
Ufree = TrialFESpace(Vfree)
Xpde = MultiFieldFESpace([Vfree, Vfree, Vfree, Vfree, Vfree, Vfree, Vfree, Vfree])
Ypde = MultiFieldFESpace([Ufree, Ufree, Ufree, Ufree, Ufree, Ufree, Ufree, Ufree])
Expand All @@ -21,6 +37,7 @@ function marine(args...; n = 400, T = 10.0, kwargs...)
conv(u, ∇u) = (∇u one(∇u)) u
c(u, v) = conv (v, (u))
function res(y, u, p)
<<<<<<< HEAD
# y1, y2, y3, y4, y5, y6, y7, y8 = y
# m1, m2, m3, m4, m5, m6, m7, m8, g1, g2, g3, g4, g5, g6, g7 = u
# p1, p2, p3, p4, p5, p6, p7, p8 = v
Expand All @@ -29,6 +46,29 @@ function marine(args...; n = 400, T = 10.0, kwargs...)
c(y[j], p[j]) - u[8 + j -1] * y[j - 1] + (u[j] + u[8 + j]) * y[j]
)dΩ for j=2:7) +
( c(y[8], p[8]) - u[15] * y[7] + u[8] * y[8] )dΩ
=======
y1, y2, y3, y4, y5, y6, y7, y8 = y
m1, m2, m3, m4, m5, m6, m7, m8, g1, g2, g3, g4, g5, g6, g7 = u
p1, p2, p3, p4, p5, p6, p7, p8 = p
return ( c(y1, p1) + p1 * (m1 + g1) * y1 +
c(y2, p2) - p2 * g1 * y1 + p2 * (m2 + g2) * y2 +
c(y3, p3) - p3 * g2 * y2 + p3 * (m3 + g3) * y3 +
c(y4, p4) - p4 * g3 * y3 + p4 * (m4 + g4) * y4 +
c(y5, p5) - p5 * g4 * y4 + p5 * (m5 + g5) * y5 +
c(y6, p6) - p6 * g5 * y5 + p6 * (m6 + g6) * y6 +
c(y7, p7) - p7 * g6 * y6 + p7 * (m7 + g7) * y7 +
c(y8, p8) - p8 * g7 * y7 + p8 * m8 * y8 )dΩ
#=
return ∫( c(y[1], p[1]) + (u[1] + u[8 + 1]) * y[1] +
c(y[2], p[2]) - u[8 + 2 -1] * y[2 - 1] + (u[2] + u[8 + 2]) * y[2] +
c(y[3], p[3]) - u[8 + 3 -1] * y[3 - 1] + (u[3] + u[8 + 3]) * y[3] +
c(y[4], p[4]) - u[8 + 4 -1] * y[4 - 1] + (u[4] + u[8 + 4]) * y[4] +
c(y[5], p[5]) - u[8 + 5 -1] * y[5 - 1] + (u[5] + u[8 + 5]) * y[5] +
c(y[6], p[6]) - u[8 + 6 -1] * y[6 - 1] + (u[6] + u[8 + 6]) * y[6] +
c(y[7], p[7]) - u[8 + 7 -1] * y[7 - 1] + (u[7] + u[8 + 7]) * y[7] +
c(y[8], p[8]) - u[15] * y[7] + u[8] * y[8] )dΩ
=#
>>>>>>> 795a1d06d1b16edeafe3d23aef2a206e8ce94b7a
end

trian = Triangulation(model)
Expand Down Expand Up @@ -83,16 +123,16 @@ function marine(args...; n = 400, T = 10.0, kwargs...)
1.0 12.0 198.0 707.0 2562.0 3163.0 3232.0 5566.0
]

objterm = PDEOptimizationProblems.InterpolatedEnergyFETerm(8, 21, zmes, 1, τ, dΩ)
objterm = PDEOptimizationProblems.InterpolatedEnergyFETerm(8, 21, zmes, 1, τ, dΩ, 1 / n)
f = (y, u) -> PDEOptimizationProblems.interpolated_measurement(objterm, y)

ndofs = Gridap.FESpaces.num_free_dofs(Ypde) + Gridap.FESpaces.num_free_dofs(Ycon)
xin = zeros(ndofs)
return GridapPDENLPModel(xin, f, trian, Ypde, Ycon, Xpde, Xcon, op, name = "Marine Population Dynamics n=$n")
end

apinene_meta = Dict(
:name => "apinene",
marine_meta = Dict(
:name => "marine",
:domaindim => UInt8(1),
:pbtype => :yu,
:nθ => 0,
Expand Down
2 changes: 1 addition & 1 deletion test/runtests.jl
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,7 @@ using PDEOptimizationProblems
include("utils.jl")

# Test that every problem can be instantiated.
for prob in setdiff(names(PDEOptimizationProblems), [:PDEOptimizationProblems])
for prob in [:marine] #setdiff(names(PDEOptimizationProblems), [:PDEOptimizationProblems])
@time begin
print(prob)
prob_fn = PDEOptimizationProblems.eval(prob)
Expand Down

0 comments on commit ec8a501

Please # to comment.