SciML / ModelingToolkit.jl

An acausal modeling framework for automatically parallelized scientific machine learning (SciML) in Julia. A computer algebra system for integrated symbolics for physics-informed machine learning and automated transformations of differential equations
https://mtk.sciml.ai/dev/
Other
1.42k stars 207 forks source link

Transform NonlinearSystem to OptimizationSystem? #1077

Open hzgzh opened 3 years ago

hzgzh commented 3 years ago

I have a quation about transform NonlinearSystem and OptimizationSystem. I construct a system based on NonlinearSystem,a source and two turbine and a sink, change source.p and sink.p,we can get the flow through the system. if I want fix the flow to some const value, modify the source.p or sink.p, then need construct the optimization problem. so I want to know can I transform the nonlinearsystem to optimizationsystem, the I provide the loss function like loss = abs(turbine1.a.m - 700.), and solve the optimization problem in ModelingToolkit.

using XSteam 
using ModelingToolkit
using NonlinearSolve
using GalacticOptim

@register XSteam.phs(p, h)
@register XSteam.psh(p, s)
@register XSteam.pht(p, h)
@register XSteam.phs(p, h)
@register XSteam.phv(p, h)
function Port(;name)
    @variables p m h 
     NonlinearSystem(Equation[], [p, m, h], [], name=name, defaults=Dict(p=>1.0, m=>0.0, h=>1000))
end

function Source(;name, p, h)
    val_p = p
    val_h = h
    @named a = Port()
    @parameters p h
    eqs = [
        a.p ~ p
        a.h ~ h
    ]
    NonlinearSystem(eqs, [], [p, h], systems=[a], defaults=Dict(p => val_p, h => val_h), name=name)
end

function Sink(;name, p)
    val_p = p
    @named a = Port()
    @parameters p
    eqs = [
        a.p ~ p
    ]
    NonlinearSystem(eqs, [], [p], systems=[a], defaults=Dict(p => val_p), name=name)
end

function connect_ports(ps...)
    eqs = [
           0 ~ sum(p->p.m, ps) # sum(mdot) = 0
          ]
     for i in 2:length(ps)
        push!(eqs, ps[1].p ~ ps[i].p)
        push!(eqs, ps[1].h ~ ps[i].h)
    end

    return eqs
end

function Turbine(;name, p_in_ref, p_out_ref, h_in_ref, η_ref, m_ref)
    val_p₀ = p_in_ref
    val_p₁ = p_out_ref
    val_h₀ = h_in_ref
    val_η  = η_ref
    val_m₀ = m_ref
    @named a = Port()
    @named b = Port()
    @variables power
    @parameters p_in_ref,p_out_ref,h_in_ref,η_ref,m_ref

    eqs = [
           0 ~ a.m - m_ref*sqrt((a.p^2 - b.p^2)/(p_in_ref^2 - p_out_ref^2)) * sqrt(p_in_ref * XSteam.phv(a.p, a.h) / a.p / XSteam.phv(p_in_ref, h_in_ref))
           η  ~ (a.h - b.h) / (a.h - XSteam.psh(b.p, XSteam.phs(a.p, a.h)))
           0 ~ a.m + b.m
           0 ~ power - a.m * (a.h - b.h)
          ]
    NonlinearSystem(eqs,[power], [p_in_ref, p_out_ref, h_in_ref, η_ref, m_ref], systems=[a, b]
    , defaults=Dict(power => 0.0, p_in_ref => val_p₀, p_out_ref => val_p₁, h_in_ref => val_h₀, η_ref => val_η, m_ref => val_m₀), name=name)
end

const p₀ = 270
const h₀ = 3475.1
const p₁ = 79.09
const m₀ = 800.253
const η = 0.8831

@named source = Source(p = 200.0, h = h₀)
@named turbine1 = Turbine(p_in_ref = p₀, p_out_ref = p₁, h_in_ref = h₀, η_ref = η, m_ref = m₀)
@named turbine2 = Turbine(p_in_ref = 79.09, p_out_ref = 60.0, h_in_ref = 3069.3,  η_ref = 0.918, m_ref = m₀ -41.0)
@named sink = Sink(p = 60.0)

cnnt_eqs = [
            connect_ports(source.a, turbine1.a)
            connect_ports(turbine1.b, turbine2.a)
            connect_ports(turbine2.b, sink.a)
            ]

@named tes_model = NonlinearSystem(cnnt_eqs, [], [], systems=[source, turbine1, turbine2, sink])

sys = structural_simplify(tes_model)

u0 = [
     turbine1.b.p => 65.0
     turbine1.b.h => 3106.0
     turbine1.b.h => 3039.0
]
# these some way to transfrom the sys to a OptimizationSystem?

prob = NonlinearProblem(sys, u0)

sol = solve(prob, NewtonRaphson())

Did OptimizationSystem have signature like OptimizatonSystem(sys,loss,u,p)?

ChrisRackauckas commented 3 years ago

What would this system do? Have the nonlinear systems be the constraints?

hzgzh commented 3 years ago

the system calculate the flow throught the turbine according the input and output pressure. yes, the nonlinear system be the optimization system constraint, can implement in ModelingToolkit?

ChrisRackauckas commented 3 years ago

Not right now.