diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index c0dc91ef..40e0bf07 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -54,18 +54,25 @@ jobs: name: Documentation runs-on: ubuntu-latest steps: - - uses: actions/checkout@v2 + - uses: actions/checkout@v4 - uses: julia-actions/setup-julia@v1 with: - version: '1.6' - #- run: | - # git config --global user.name name - # git config --global user.email email - # git config --global github.user username - - name: Install dependencies - run: julia --project=docs/ -e 'using Pkg; Pkg.develop(PackageSpec(path=pwd())); Pkg.instantiate()' - - name: Build and deploy + version: '1' + - name: Configure doc environment + shell: julia --project=docs --color=yes {0} + run: | + using Pkg + Pkg.develop(PackageSpec(path=pwd())) + Pkg.instantiate() + - uses: julia-actions/julia-buildpkg@v1 + - uses: julia-actions/julia-docdeploy@v1 env: - GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }} # For authentication with GitHub Actions token - DOCUMENTER_KEY: ${{ secrets.DOCUMENTER_KEY }} # For authentication with SSH deploy key - run: julia --project=docs/ docs/make.jl + GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }} + DOCUMENTER_KEY: ${{ secrets.DOCUMENTER_KEY }} + - name: Run doctests + shell: julia --project=docs --color=yes {0} + run: | + using Documenter: DocMeta, doctest + using PowerDynamics + DocMeta.setdocmeta!(PowerDynamics, :DocTestSetup, :(using PowerDynamics); recursive=true) + doctest(PowerDynamics) diff --git a/.gitignore b/.gitignore index 105057f6..276ebba9 100644 --- a/.gitignore +++ b/.gitignore @@ -21,3 +21,4 @@ Manifest.toml target/ .vscode/ +/docs/src/generated/ diff --git a/Project.toml b/Project.toml index ebe67e9d..abcc4fa4 100644 --- a/Project.toml +++ b/Project.toml @@ -5,6 +5,7 @@ version = "3.1.6" [deps] BlockSystems = "4663d367-db1f-4bef-81e7-dc1bd7f7b428" +ControlSystems = "a6e380b2-a6ca-5380-bf3e-84a91bcd477e" DiffEqBase = "2b5f629d-d688-5b77-993f-72d75c75574e" DiffEqCallbacks = "459566f4-90b8-5000-8ac3-15dfb0a30def" Graphs = "86223c79-3864-5bf0-83f7-82e725a168b6" @@ -15,6 +16,7 @@ LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" Logging = "56ddb016-857b-54e1-b83d-db4d58db5568" MacroTools = "1914dd2f-81c6-5fcd-8719-6d5c9610ff09" Markdown = "d6f4376e-aef5-505a-96c1-9c027394607a" +ModelingToolkit = "961ee093-0014-501f-94e3-6117800e7a78" NLsolve = "2774e3e8-f4cf-5e23-947b-6d7e65073b56" NetworkDynamics = "22e9dc34-2a0d-11e9-0de0-8588d035468b" OrderedCollections = "bac558e1-5e72-5ebc-8fee-abe8a469f55d" @@ -28,9 +30,11 @@ SciMLNLSolve = "e9a6253c-8580-4d32-9898-8661bb511710" Setfield = "efcf1570-3423-57d1-acb7-fd33fddbac46" StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" SteadyStateDiffEq = "9672c7b4-1e72-59bd-8a11-6ac3964bc41f" +Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d" [compat] BlockSystems = "0.4.2" +ControlSystems = "1" DiffEqBase = "^6.9.4" DiffEqCallbacks = "^2" Graphs = "1.4" @@ -38,6 +42,7 @@ Ipopt = "0.7, 0.8, 0.9, 1" JSON = "^0.21.0" Lazy = "^0.14.0, 0.15" MacroTools = "^0.5.1" +ModelingToolkit = "8" NLsolve = "^4.1.0" NetworkDynamics = "0.8" OrderedCollections = "1.3" @@ -51,6 +56,7 @@ SciMLNLSolve = "0.1" Setfield = "1" StaticArrays = "1" SteadyStateDiffEq = "1, 2" +Unitful = "1" julia = "1.6" [extras] diff --git a/README.md b/README.md index 072c9bb3..60792cef 100644 --- a/README.md +++ b/README.md @@ -22,3 +22,9 @@ If you use PowerDynamics.jl in your research publications, [please cite our pape publisher={Elsevier} } ``` + +## Funding + +Development of `PowerDynamics.jl` was funded by the *German Federal Ministry for Economic Affairs and Climate Action* + + diff --git a/bmwk_logo_en.svg b/bmwk_logo_en.svg new file mode 100644 index 00000000..5bed8c7c --- /dev/null +++ b/bmwk_logo_en.svg @@ -0,0 +1 @@ + \ No newline at end of file diff --git a/docs/localmake.jl b/docs/localmake.jl index 8e39f39b..88b360ac 100755 --- a/docs/localmake.jl +++ b/docs/localmake.jl @@ -10,12 +10,6 @@ At the end of each run the user is prompted to rerun the make process. Using rev use the updated `*.md` and source files. This way the Julia session keeps alive and the individual builds are much faster. =# -using Revise -using LiveServer -using REPL.TerminalMenus - -port = isempty(ARGS) ? 8000 : parse(Int, ARGS[1]) -@assert 8000 ≤ port ≤ 9000 "port has to be in range 8000..9000!" using Pkg Pkg.activate(@__DIR__) @@ -23,6 +17,13 @@ Pkg.develop(PackageSpec(path=dirname(@__DIR__))) # adds the package this script Pkg.instantiate() Pkg.update() +using Revise +using LiveServer +using REPL.TerminalMenus + +port = isempty(ARGS) ? 8000 : parse(Int, ARGS[1]) +@assert 8000 ≤ port ≤ 9000 "port has to be in range 8000..9000!" + @info "Start server..." @async serve(;dir=joinpath(@__DIR__, "build"), port) diff --git a/docs/make.jl b/docs/make.jl index bc51b875..efeaafc1 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -1,9 +1,12 @@ using Documenter using PowerDynamics using BlockSystems +using Literate DocMeta.setdocmeta!(PowerDynamics, :DocTestSetup, :(using PowerDynamics); recursive=true) +Literate.markdown(joinpath(@__DIR__, "src", "MARiE.jl"), joinpath(@__DIR__, "src", "generated")) + makedocs( modules = [PowerDynamics], authors = "Tim Kittel, Jan Liße, Sabine Auer and further contributors.", @@ -16,6 +19,7 @@ makedocs( "powergrid_model.md", "node_types.md", "custom_node_types.md", + "Modular Inverters" => "generated/MARiE.md", "line_types.md", "states_solutions.md", "simulations.md", diff --git a/docs/src/MARiE.jl b/docs/src/MARiE.jl new file mode 100644 index 00000000..c38ce98f --- /dev/null +++ b/docs/src/MARiE.jl @@ -0,0 +1,145 @@ +#= +```@meta +CurrentModule = PowerDynamics.ModularInverter +``` +# Modular Inverter Framework +This section describes the modular inverter framework developed in the MARiE-Project funded by the +*German Federal Ministry for Economic Affairs and Climate Action*. + +## Structure + +The inverter model is split in two parts: the inner loop and the outer loop. + +The *inner loop* either acts as a voltage source +``` + +---------+ +u_d_ref -->| | +u_q_ref -->| Voltage |--> u_d +i_d_dis -->| Source |--> u_q +i_q_dis -->| | + +---------+ +``` +which gets a `dq` voltage reference and tries to achive this reference despite the `dq` current disturbance. + +Alternatively, the inner loop can mimic a current source, which tries to realsize a given current setpoint despite a disturbance voltage at the output terminals: +``` + +---------+ +i_d_ref -->| | +i_q_ref -->| Current |--> i_d +u_d_dis -->| Source |--> i_q +u_q_dis -->| | + +---------+ +``` + +The *outer loop* provides references to the inner loop blocks. +As an input, outerloops have acces to the measured currents and voltages at the point of common coupling. +They either eject a dq-reference for the inner voltage source or the inner current soruce: + +``` + +---------+ +i_d_meas -->| | +i_q_meas -->| Outer |--> u_d_ref or i_d_ref +u_d_meas -->| Loop |--> u_q_ref or i_q_ref +u_q_meas -->| | + +---------+ +``` + +`PowerDynamics.jl` provides predefined outer loop as well as inner loop models as `IOBlocks` provide by `BlockSystems.jl`, which can be used to create `IONodes` (see [Custom-Nodes-using-BlockSystems.jl](@ref)). + +## Example +In order to construct a modular inverter we need to import the `ModularInverter` module to the namespace: +=# +using PowerDynamics +using BlockSystems +using PowerDynamics.ModularInverter +nothing #hide +#= +First, we need to define the inner loop as an `IOBlock` object. +You can either use your own or make use of the [Predefined Inner Loops](@ref). +=# +BlockSystems.WARN[] = false # hide +inner = CascadedVoltageControl() +#= +As you can see, the input/output structure of that innerloop adheres to the volten source innerloop convention. + +Next you define the outer loop, either create you own based on the interface or pick from the [Predefined Outer Loops](@ref). +=# +outer = DroopControl() +#= +There are still a lot of open parameters, we need to assign specific values to create the `IONode`(@ref) later. +We can do so by using the keyword aguments of the constructor: +=# +outer = DroopControl(; + ω_ref = 2π*50, # frequency reference + V_ref = 1, # voltage reference (pu) + P_ref = 1, # active power reference (pu) + Q_ref = 0, # reactive power reference (pu) + K_P = 0.4, # active droop gain + K_Q = 0.004, # reactive droop gain + τ_P = 0.1, # active power filter constant + τ_Q = 0.1, # reactive power filter constant +) + +#= +The full inverter model is created by combining outer and inner model: +=# +inv = Inverter(inner, outer) +#= +Which satisfies the interface for `IONodes`. To create the node +=# +node = IONode(inv); + +#= +Which is a valid node definition which can be used in any `PowerDynamics.jl` context. +=# + +#= +## Model Reduction + +`PowerDynamics.jl` also implements automated model order reduction of linear models using the balance residualization technique. +Internaly, we make use of the method implementaion in `ControlSystems.jl`. For that, a linear block (such as the cascaded inner loops) can be decomposed in to its `A`, `B`, `C` and `D` matrix for regular LTI representation. + +For example, we could creat a 10th order representation of the previously defined inner loop by calling +=# +reduced_inner = balresid(inner, 10; reconstruct=false) +#= +To obtain the reduced order model we can construct the inverter again: +=# +reduced_inverter = Inverter(reduced_inner, outer) +# ... and the reduced node +reduced_node = IONode(reduced_inverter); + +#= +```@docs +balresid +``` +=# + +#= +## Predefined Inner Loops +```@docs +CascadedVoltageControl +CascadedVoltageControlCS +CascadedCurrentControl +PT1Source +IdealSource +IdealCSource +``` + +## Predefined Outer Loops +```@docs +FixedVoltage +DroopControl +Synchronverter +PLLCurrent +ConstantPower +FiltConstantPower +``` + +## Inverter Construction +```@docs +Inverter +InverterCS +``` + +=# diff --git a/src/IONodes/IOComponents.jl b/src/IONodes/IOComponents.jl index 1712f2d5..43f9f313 100644 --- a/src/IONodes/IOComponents.jl +++ b/src/IONodes/IOComponents.jl @@ -262,4 +262,34 @@ function ImpedanceConstraint(;name=gensym(:rxconstraint), renamings...) return isempty(renamings) ? block : rename_vars(block; renamings...) end +""" + Cart2Polar(;name=:c2p, renamings...) + +(X, Y) ↦ (mag, arg) transformation +""" +function Cart2Polar(;name=:c2p, renamings...) + @variables t arg(t) mag(t) + @parameters x(t) y(t) + block = IOBlock([mag ~ √(x^2 + y^2), + arg ~ atan(y, x)], + [x, y], [mag, arg]; name) + + replace_vars(block; renamings...) +end + +""" + Polar2Cart(;name=:p2c, renamings...) + +(mag, arg) ↦ (X, Y) transformation +""" +function Polar2Cart(;name=:p2c, renamings...) + @variables t x(t) y(t) + @parameters arg(t) mag(t) + block = IOBlock([x ~ mag * cos(arg), + y ~ mag * sin(arg)], + [mag, arg], [x, y]; name) + + replace_vars(block; renamings...) +end + end # module diff --git a/src/IONodes/IONode.jl b/src/IONodes/IONode.jl index fa8fddbd..6df598ce 100644 --- a/src/IONodes/IONode.jl +++ b/src/IONodes/IONode.jl @@ -51,6 +51,7 @@ function IONode(bp::BlockPara) end IONode(blk::IOBlock, para::Dict) = IONode(BlockPara(blk, para)) +IONode(blk::IOBlock) = IONode(BlockPara(blk, Dict())) # extend the necessary functions for the `AbstractNode` interface function construct_vertex(ion::IONode) diff --git a/src/IONodes/LTI.jl b/src/IONodes/LTI.jl new file mode 100644 index 00000000..21749ebd --- /dev/null +++ b/src/IONodes/LTI.jl @@ -0,0 +1,258 @@ +export identify_lti, balresid + +using ControlSystems +using BlockSystems: @check +using ModelingToolkit: value +using ModelingToolkit.SymbolicUtils: Symbolic +using LinearAlgebra + +function _state_matrix(expr, vars) + if eltype(expr) <: Equation + expr = getproperty.(expr, :rhs) + end + + A = Matrix(undef, length(expr), length(vars)) + fill!(A, NaN) + for (ieq, ex) in enumerate(expr) + # ex = simplify(ex, expand=true) + for (istate, state) in enumerate(vars) + coeff = (ex - substitute(ex, state => 0))/state |> simplify + if !isempty(Set(get_variables(coeff)) ∩ Set(vars)) + error("Could not properly extract $state from $(ex), got $coeff") + end + A[ieq, istate] = coeff + end + end + return _narrow_type(A) +end + +""" + A,B,C,D,x,y,u=identify_lti(blk::IOBlock) + +Identify the matrices A,B,C and D from IOBlock. +""" +function identify_lti(blk::IOBlock) + inputs = blk.inputs + outputs = blk.outputs + output_rhs = similar(outputs, Any) + odestates = Symbolic[] + odeidx = Int[] + algstates = Symbolic[] + algidx = Int[] + for (i, eq) in enumerate(equations(blk)) + (type, var) = BlockSystems.eq_type(eq) + if type == :explicit_diffeq + push!(odestates, var) + push!(odeidx, i) + elseif type == :explicit_algebraic + @assert var ∈ Set(outputs) "Explicit algebraic equations musst represent outputs!" + push!(algstates, var) + push!(algidx, i) + else + error("Only pure LTI systems. Can not handle $type-type equations!") + end + #check if variable is an output + outidx = findfirst(isequal(var), outputs) + if !isnothing(outidx) # is output variable + if type == :explicit_diffeq + output_rhs[outidx] = var + elseif type == :explicit_algebraic + output_rhs[outidx] = eq.rhs + else + error() + end + end + end + + for eq in equations(blk) + @assert isempty(Set(get_variables(eq.rhs)) ∩ algstates) "The rhs of the equations must not contain algebraic states!" + end + + A = _state_matrix(equations(blk)[odeidx], odestates) + B = _state_matrix(equations(blk)[odeidx], inputs) + C = _state_matrix(output_rhs, odestates) + D = _state_matrix(output_rhs, inputs) + + Avars = isempty(A) ? Set() : Set(mapreduce(get_variables, vcat, A)) + Bvars = isempty(B) ? Set() : Set(mapreduce(get_variables, vcat, B)) + Cvars = isempty(C) ? Set() : Set(mapreduce(get_variables, vcat, C)) + Dvars = isempty(D) ? Set() : Set(mapreduce(get_variables, vcat, D)) + + @assert Avars ⊆ Set(blk.iparams) "Matrix A contains non-parameters. Thats an error! $Avars" + @assert Bvars ⊆ Set(blk.iparams) "Matrix B contains non-parameters. Thats an error! $Bvars" + @assert Cvars ⊆ Set(blk.iparams) "Matrix C contains non-parameters. Thats an error! $Cvars" + @assert Dvars ⊆ Set(blk.iparams) "Matrix D contains non-parameters. Thats an error! $Dvars" + + if eltype(A) == Any + # @warn "There appear to be `Any` terms in A" + A = fixdiv(A) + end + if eltype(B) == Any + # @warn "There appear to be `Any` terms in B" + B = fixdiv(B) + end + if eltype(C) == Any + # @warn "There appear to be `Any` terms in C" + C = fixdiv(C) + end + if eltype(D) == Any + # @warn "There appear to be `Any` terms in D" + D = fixdiv(D) + end + + return A, B, C, D, odestates, outputs, inputs +end + +function fixdiv(A) + for idx in eachindex(A) + if A[idx] isa SymbolicUtils.Div + try + A[idx] = eval(Meta.parse(repr(A[idx]))) + catch e + if !isa(e, UndefVarError) + rethrow(e) + end + end + end + end + A = _narrow_type(A) +end + +function IOBlock(A::Matrix,B::Matrix,C::Matrix,D::Matrix,x,y,u; name=gensym(:LTI), warn=BlockSystems.WARN[], rem_eqs=Equation[]) + @assert size(A)[1] == size(A)[2] == size(B)[1] == size(C)[2] == length(x) + @assert size(B)[2] == size(D)[2] == length(u) + @assert size(C)[1] == size(D)[1] == length(y) + iv_candidate = unique(vcat(Symbolics.arguments.(value.(x)), + Symbolics.arguments.(value.(y)), + Symbolics.arguments.(value.(u)))) + @assert length(iv_candidate) == 1 && length(iv_candidate[begin]) == 1 + iv = iv_candidate[1][1] + dt = Differential(iv) + eqs = vcat(dt.(x) .~ A*x + B*u, + y .~ C*x + D*u) + @show + filter!(eq -> !isequal(eq.lhs, eq.rhs), eqs) + + IOBlock(eqs, u, y; name, warn, rem_eqs) +end + +""" + balresid(blk::IOBlock, order; warn=BlockSystems.WARN[], verbose=false, reconstruct=false) + +This function creates a model reduced version of `blk` of given order using the balanced residualization method. +Only works for linear systems. Internaly, it uses the `baltrunc` function from `ControlSystems.jl` +to apply the method to an LTI. The LTI system matrices A,B,C and D are determined symbolicially from the +`IOBlock` object. + +If `reconstruct=true` the resulting system will include `removed_eqs` to reconstruct the original states from the projected states `z_i`. + +""" +function balresid(blk::IOBlock, order; warn=BlockSystems.WARN[], verbose=false, reconstruct=false, getT=false) + @check isempty(blk.iparams) "In order to balresid a system it can not have any internal parameters!" + A, B, C, D, x, y, u = identify_lti(blk) + verbose && @info "Block is LTI of order $(length(x))" + + if length(x) <= order + throw(ArgumentError("Cannot reduce a system of order $(length(x)) to order $order.")) + end + + ss = StateSpace(A,B,C,D) + ss_bal, G, _ = ControlSystems.baltrunc(ss; residual=true, n=order) + _, _, T = ControlSystems.balreal(ss) + + t = get_iv(blk) + z = Num[] + for i in 1:length(x) + zs = subscript(:z, i) + append!(z, @variables $zs(t)) + end + + z_r = z[begin:order] + z_t = z[order+1:end] + A_r = ss_bal.A + B_r = ss_bal.B + C_r = ss_bal.C + D_r = ss_bal.D + + ## check results + An = T*A*inv(T) + Bn = T*B + Cn = C*inv(T) + Dn = D + + A11 = An[1:order, 1:order] + A12 = An[1:order, order+1:end] + A21 = An[order+1:end, 1:order] + A22 = An[order+1:end, order+1:end] + @assert [A11 A12; A21 A22] == An + + B1 = Bn[1:order, :] + B2 = Bn[order+1:end, :] + @assert [B1; B2] == Bn + + C1 = Cn[:, 1:order] + C2 = Cn[:, order+1:end] + @assert [C1 C2] == Cn + + # get matrices for reduced model + @assert A_r ≈ A11 - A12*inv(A22)*A21 + @assert B_r ≈ B1 - A12*inv(A22)*B2 + @assert C_r ≈ C1 - C2*inv(A22)*A21 + @assert D_r ≈ Dn - C2*inv(A22)*B2 + ## + + verbose && @info "Truncated $(length(x)-order) states!" + + # in the sigular pertubation we assume dot(z_t)=0 but z_t is given by constraint + # we calculate z_t based on z_r and u + balA = T*A*inv(T) + A21 = balA[order+1:end, 1:order] + A22 = balA[order+1:end, order+1:end] + B2 = (T*B)[order+1:end, :] + z_t = -inv(A22)*A21*z_r -inv(A22)*B2*u + + outidx = findall(x -> x ∈ Set(y), x) + state_estimates = x .~ inv(T)*vcat(z_r, z_t) + deleteat!(state_estimates, outidx) # outputs are still present in the full system + + if reconstruct + substitutions = Dict(eq.lhs => eq.rhs for eq in state_estimates) + old_rem_eqs = map(eq->BlockSystems.eqsubstitute(eq, substitutions), blk.removed_eqs) + rem_eqs = vcat(state_estimates, old_rem_eqs) + else + rem_eqs = Equation[] + end + + name = Symbol(string(blk.name) * "_trunc_$order") + + blk = IOBlock(A_r, B_r, C_r, D_r, z_r, y, u; rem_eqs, name, warn) + return getT ? (blk, T) : blk +end + +function ControlSystems.balreal(blk::IOBlock) + A, B, C, D, x, y, u = identify_lti(blk) + ss = StateSpace(A,B,C,D) + ss_bal, G, T = ControlSystems.balreal(ss) + + return G, T, x +end + +function ControlSystems.StateSpace(iob::IOBlock) + A,B,C,D = identify_lti(iob) + StateSpace(A,B,C,D) +end + +function _narrow_type(A::AbstractArray) + isempty(A) && return A + elt = mapreduce(typeof, promote_type, A) + convert.(elt, A) +end + +function subscript(s, i, add...) + Symbol(s, _subscript(i), add...) +end + +function _subscript(i::Int) + dig = reverse(digits(i)) + String(map(i -> Char(0x02080 + i), dig)) +end diff --git a/src/IONodes/MIComponents.jl b/src/IONodes/MIComponents.jl new file mode 100644 index 00000000..8c9c95ad --- /dev/null +++ b/src/IONodes/MIComponents.jl @@ -0,0 +1,243 @@ +module MIComponents + +using BlockSystems +using LinearAlgebra + +""" + VRefGen(; name=:Vrefgen, renamings...) + +Create dq-reference from (ω, V) reference. +""" +function VRefGen(; name=:Vrefgen, renamings...) + @variables t δ(t) u_ref_r(t) u_ref_i(t) + @parameters V(t) ω(t) + dt = Differential(t) + block = IOBlock([dt(δ) ~ ω, + u_ref_r ~ V*cos(δ), + u_ref_i ~ V*sin(δ)], + [V, ω], + [u_ref_r, u_ref_i]; + name) + replace_vars(block; renamings...) +end + +""" + UIMeas(; name=:uimeas, renamings...) + +Creates a very simple block which is mainly there for renaming. + +TODO: is this necessary? Maybe make i_i and i_r "globalp" inputs +""" +function UIMeas(; name=:ui_meas, renamings...) + @variables t u_meas_r(t) u_meas_i(t) i_meas_r(t) i_meas_i(t) + @parameters u_r(t) u_i(t) i_i(t) i_r(t) + block = IOBlock([u_meas_r ~ u_r, + u_meas_i ~ u_i, + i_meas_r ~ i_r, + i_meas_i ~ i_i], + [u_r, u_i, i_r, i_i], + [u_meas_r, u_meas_i, i_meas_r, i_meas_i]; + name) + + replace_vars(block; renamings...) +end + +function JuanPLL(; name=:pll, renamings...) + @variables t δ_pll(t) ω_pll(t) u_q(t) + @parameters u_r(t) u_i(t) Kp Ki + dt = Differential(t) + block = IOBlock([u_q ~ 1/sqrt(u_r^2+u_i^2)*(u_r*sin(-δ_pll) + u_i*cos(-δ_pll)), + dt(δ_pll) ~ ω_pll + Kp*u_q, + dt(ω_pll) ~ Ki*u_q + ], + [u_r, u_i], + [δ_pll, ω_pll]; + name) + block = substitute_algebraic_states(block) + replace_vars(block; renamings...) +end + +function L(; kwargs...) + @parameters ω0 Rf Rg Lf Lg + @variables t i_f_r(t) i_f_i(t) + @parameters V_I_r(t) V_I_i(t) V_C_r(t) V_C_i(t) + dt = Differential(t) + + W = [0 1; -1 0] + i_f = [i_f_r, i_f_i] + + x = [i_f_r, i_f_i] + u = [V_I_r, V_I_i, V_C_r, V_C_i] + A = W*ω0 - Rf/Lf * I + B = [1/Lf*I(2) -1/Lf*I] + + L = IOBlock(dt.(x) .~ A*x + B*u, + [V_I_r, V_I_i, V_C_r, V_C_i], + [i_f_r, i_f_i], + name = :L) + replace_vars(L; kwargs...) +end + +function LC(; name=:LC, kwargs...) + @parameters ω0 Rf Rg Lf Lg C + @variables t i_f_r(t) i_f_i(t) V_C_r(t) V_C_i(t) + @parameters V_I_r(t) V_I_i(t) i_g_r(t) i_g_i(t) + dt = Differential(t) + + W = [0 1; -1 0] + i_f = [i_f_r, i_f_i] + V_C = [V_C_r, V_C_i] + + x = [i_f_r, i_f_i, V_C_r, V_C_i] + u = [V_I_r, V_I_i, i_g_r, i_g_i] + A = [W*ω0 - Rf/Lf*I -1/Lf*I ; + 1/C*I W*ω0] + B = [1/Lf*I zeros(2,2); + zeros(2,2) -1/C*I] + + LC = IOBlock(dt.(x) .~ A*x + B*u, + [i_g_r, i_g_i, V_I_r, V_I_i], + [i_f_r, i_f_i, V_C_r, V_C_i]; + name) + + replace_vars(LC; kwargs...) +end + +function LCL(; kwargs...) + @parameters ω0 Rf Rg Lf Lg C + @variables t i_f_r(t) i_f_i(t) i_g_r(t) i_g_i(t) V_C_r(t) V_C_i(t) + @parameters V_I_r(t) V_I_i(t) V_g_r(t) V_g_i(t) + dt = Differential(t) + W = [0 1; -1 0] + + i_f = [i_f_r, i_f_i] + V_C = [V_C_r, V_C_i] + i_g = [i_g_r, i_g_i] + + x = [i_f_r, i_f_i, V_C_r, V_C_i, i_g_r, i_g_i] + u = [V_I_r, V_I_i, V_g_r, V_g_i] + A = [W*ω0 - Rf/Lf*I -1/Lf*I zeros(2,2) ; + 1/C*I W*ω0 -1/C*I ; + zeros(2,2) 1/Lg*I W*ω0 - Rg/Lg*I ] + B = [1/Lf*I zeros(2,2); + zeros(2,2) zeros(2,2); + zeros(2,2) -1/Lg*I ] + + LCL = IOBlock(dt.(x) .~ A*x + B*u, + [V_g_r, V_g_i, V_I_r, V_I_i], + [i_f_r, i_f_i, V_C_r, V_C_i, i_g_r, i_g_i], + name = :LCL) + replace_vars(LCL; kwargs...) +end + +function CC1(; kwargs...) + @parameters ω0 Rf Rg Lf Lg C + @variables t γ_r(t) γ_i(t) V_I_r(t) V_I_i(t) + @parameters KP KI i_f_r(t) i_f_i(t) i_f_ref_r(t) i_f_ref_i(t) V_C_r(t) V_C_i(t) F + dt = Differential(t) + W = [0 1; -1 0] + + γ = [γ_r, γ_i] + i_f = [i_f_r, i_f_i] + i_f_ref = [i_f_ref_r, i_f_ref_i] + V_I = [V_I_r, V_I_i] + V_C = [V_C_r, V_C_i] + + CC1 = IOBlock(vcat(dt.(γ) .~ i_f_ref - i_f, + V_I .~ -Lf*ω0*W*i_f + KP*(i_f_ref - i_f) + KI*γ + F*V_C), + [i_f_r, i_f_i, i_f_ref_r, i_f_ref_i, V_C_r, V_C_i], + [V_I_r, V_I_i], + name=:CC1) + + replace_vars(CC1; kwargs...) +end + +function VC(; kwargs...) + @parameters ω0 Rf Rg Lf Lg C + @variables t γ_r(t) γ_i(t) i_f_ref_r(t) i_f_ref_i(t) + @parameters KP KI F V_C_r(t) V_C_i(t) V_C_ref_r(t) V_C_ref_i(t) i_g_r(t) i_g_i(t) ω0 C + dt = Differential(t) + W = [0 1; -1 0] + + γ = [γ_r, γ_i] + V_C = [V_C_r, V_C_i] + i_g = [i_g_r, i_g_i] + V_C_ref = [V_C_ref_r, V_C_ref_i] + i_f_ref = [i_f_ref_r, i_f_ref_i] + + VC = IOBlock(vcat(dt.(γ) .~ V_C_ref - V_C, + i_f_ref .~ - C*ω0*W*V_C + KP*(V_C_ref - V_C) + KI*γ + F*i_g), + [i_g_r, i_g_i, V_C_r, V_C_i, V_C_ref_r, V_C_ref_i], + [i_f_ref_r, i_f_ref_i], + name=:VC) + replace_vars(VC; kwargs...) +end + +function CC2(;kwargs...) + @parameters ω0 Rf Rg Lf Lg C + @variables t γ_r(t) γ_i(t) V_C_ref_r(t) V_C_ref_i(t) + @parameters KP KI i_g_r(t) i_g_i(t) i_g_ref_r(t) i_g_ref_i(t) V_g_r(t) V_g_i(t) F + dt = Differential(t) + W = [0 1; -1 0] + + γ = [γ_r, γ_i] + i_g = [i_g_r, i_g_i] + i_g_ref = [i_g_ref_r, i_g_ref_i] + V_C_ref = [V_C_ref_r, V_C_ref_i] + V_g = [V_g_r, V_g_i] + + CC2 = IOBlock(vcat(dt.(γ) .~ i_g_ref - i_g, + V_C_ref .~ -Lf*ω0*W*i_g + KP*(i_g_ref - i_g) + KI*γ + F*V_g), + [i_g_r, i_g_i, i_g_ref_r, i_g_ref_i, V_g_r, V_g_i], + [V_C_ref_r, V_C_ref_i], + name=:CC2) + replace_vars(CC2; kwargs...) +end + +function CC1_PR() + dq = CC1() + + renamings = (; u_r=:i_f_r, u_ref_r=:i_f_ref_r, u_i=:i_f_i, u_ref_i=:i_f_ref_i) + pr1 = PR(;name=:PR1, renamings...) + pr2 = PR(;name=:PR2, renamings...) + pr3 = PR(;name=:PR3, renamings...) + pr4 = PR(;name=:PR4, renamings...) + + @variables t + @parameters dq_V_I_r(t) dq_V_I_i(t) pr_V_I_r(t) pr_V_I_i(t) + @variables V_I_r(t) V_I_i(t) + + add = IOBlock([V_I_r ~ dq_V_I_r + pr_V_I_r, + V_I_i ~ dq_V_I_i + pr_V_I_i], + [dq_V_I_r, dq_V_I_i, pr_V_I_r, pr_V_I_i], + [V_I_r, V_I_i], + name=:DQPR_add) + + CC1pr = IOSystem([dq.V_I_r => add.dq_V_I_r, + dq.V_I_i => add.dq_V_I_i, + pr1.y_r + pr2.y_r + pr3.y_r + pr4.y_r => add.pr_V_I_r, + pr1.y_i + pr2.y_i + pr3.y_i + pr4.y_i => add.pr_V_I_i], + [add, dq, pr1, pr2, pr3, pr4], + globalp=[:i_f_r, :i_f_i, :i_f_ref_r, :i_f_ref_i], + outputs=[add.V_I_r, add.V_I_i], + namespace_map=[add.V_I_r=>:V_I_r, + add.V_I_i=>:V_I_i], name=:CC1_pr) + connect_system(CC1pr; name=:CC1) +end + +function PR(;name=:PR, kwargs...) + @variables t x₁_r(t) x₂_r(t) x₁_i(t) x₂_i(t) y_r(t) y_i(t) + @parameters K ω ωc u_r(t) u_i(t) u_ref_r(t) u_ref_i(t) + dt = Differential(t) + blk = IOBlock([dt(x₁_r) ~ x₂_r, + dt(x₂_r) ~ -ω^2*x₁_r - 2*ωc*x₂_r + u_ref_r - u_r, + y_r ~ 2*K*ωc*x₂_r, + dt(x₁_i) ~ x₂_i, + dt(x₂_i) ~ -ω^2*x₁_i - 2*ωc*x₂_i + u_ref_r - u_i, + y_i ~ 2*K*ωc*x₂_i], + [u_r, u_i, u_ref_r, u_ref_i], [y_r, y_i]; name) + + replace_vars(blk; kwargs...) +end + +end # module diff --git a/src/IONodes/ModularInverter.jl b/src/IONodes/ModularInverter.jl new file mode 100644 index 00000000..9dc3ab15 --- /dev/null +++ b/src/IONodes/ModularInverter.jl @@ -0,0 +1,508 @@ +using Unitful +using BlockSystems +using .MIComponents +using PowerDynamics.IOComponents + +export MIParameters +export Inverter, InverterCS +export VSwithLoad, DroopControl, Synchronverter, FixedVoltage, FixedCurrent, PLLCurrent +export IdealSource, IdealCSource, PT1Source, CascadedVoltageControl, CascadedCurrentControl, CascadedVoltageControlCS, ConstantPower, FiltConstantPower + +MIBase = let ω0base = 2π*50u"rad/s", + Sbase = 10000u"W", + Vbase = sqrt(3)*230u"V", + Ibase = Sbase/(Vbase) |> u"A", + Cbase = Ibase/Vbase, + Lbase = Vbase/Ibase, + Rbase = (Vbase^2)/Sbase + (;ω0base, Sbase, Vbase, Ibase, Cbase, Lbase, Rbase) +end + +MIParameters = Dict( + # electrical parameters + :Rf => uconvert(NoUnits, 0.01u"Ω"/MIBase.Rbase), + :Rg => uconvert(NoUnits, 0.01u"Ω"/MIBase.Rbase), + :Lf => ustrip(u"s", 350e-6u"H"/MIBase.Lbase), + :Lg => ustrip(u"s", 350e-6u"H"/MIBase.Lbase), + :C => ustrip(u"s", 100e-6u"F"/MIBase.Cbase), + :ω0 => ustrip(u"rad/s", MIBase.ω0base), + # CC1 control parameters + :CC1₊KP => 1 * ustrip(u"A/V", MIBase.Ibase/MIBase.Vbase), + :CC1₊KI => 1000 * ustrip(u"A/V", MIBase.Ibase/MIBase.Vbase), + # PR configs + :CC1₊PR1₊K => 50*ustrip(u"A/V", MIBase.Ibase/MIBase.Vbase), + :CC1₊PR1₊ω => 3*2π*50, + :CC1₊PR1₊ωc => 5, + :CC1₊PR2₊K => 25*ustrip(u"A/V", MIBase.Ibase/MIBase.Vbase), + :CC1₊PR2₊ω => 6*2π*50, + :CC1₊PR2₊ωc => 5, + :CC1₊PR3₊K => 10*ustrip(u"A/V", MIBase.Ibase/MIBase.Vbase), + :CC1₊PR3₊ω => 9*2π*50, + :CC1₊PR3₊ωc => 5, + :CC1₊PR4₊K => 5*ustrip(u"A/V", MIBase.Ibase/MIBase.Vbase), + :CC1₊PR4₊ω => 12*2π*50, + :CC1₊PR4₊ωc => 5, + # VC control parameters + :VC₊KP => 0.06 * ustrip(u"V/A", MIBase.Vbase/MIBase.Ibase), + :VC₊KI => 20 * ustrip(u"V/A", MIBase.Vbase/MIBase.Ibase), + # CC2 control parameters + :CC2₊KP => 3 * ustrip(u"A/V", MIBase.Ibase/MIBase.Vbase), + :CC2₊KI => 10 * ustrip(u"A/V", MIBase.Ibase/MIBase.Vbase) +) + +""" + Inverter(inner, outer; name)) + +Combines an inner and an outer loop to a closed loop model of an voltage-srouce-inverter. + +Both inner and outer loops musst be `IOBlock` objects adhering to the interface + +Inner loop inputs/outputs: +``` + +---------+ +u_r_ref -->| | +u_i_ref -->| Voltage |--> u_r + i_r -->| Source |--> u_i + i_i -->| | + +---------+ +``` +Outer loop inputs/outputs: + +``` + +---------+ +i_r_meas -->| | +i_i_meas -->| Outer |--> u_r_ref +u_r_meas -->| Loop |--> u_i_ref +u_i_meas -->| | + +---------+ +``` +""" +function Inverter(inner, outer; name=Symbol(string(outer.name)*"_"*string(inner.name))) + @assert BlockSpec([:i_i, :i_r, :u_ref_r, :u_ref_i], [:u_r, :u_i])(inner) "Inner ctrl loop :$(inner.name) does not meet expectation." + @assert BlockSpec([], [:u_ref_r, :u_ref_i]; in_strict=false)(outer) "Outer ctrl loop :$(outer.name) does not meet expectation." + + outerinputs = ModelingToolkit.getname.(outer.inputs) + :u_r ∈ outerinputs && @error "Outer shoud use :u_meas_r instead of :u_r" + :u_i ∈ outerinputs && @error "Outer shoud use :u_meas_i instead of :u_i" + :i_r ∈ outerinputs && @error "Outer shoud use :i_meas_r instead of :i_r" + :i_i ∈ outerinputs && @error "Outer shoud use :i_meas_i instead of :i_i" + + meas = MIComponents.UIMeas() + + sys = IOSystem(:autocon, [outer, inner, meas]; + outputs=[inner.u_r, inner.u_i], + globalp=[:i_r, :i_i], + name) + closed = connect_system(sys) + + @assert BlockSpec([:i_i, :i_r], [:u_r, :u_i])(closed) "Closed loop does not match expectation! $closed" + + return closed +end + +""" + InverterCS(inner, outer; name)) + +Combines an inner and an outer loop to a closed loop model of an current-source-inverter. + +Both inner and outer loops musst be `IOBlock` objects adhering to the interface + +Inner loop inputs/outputs: +``` + +---------+ +i_r_ref -->| | +i_i_ref -->| Current |--> i_r + u_r -->| Source |--> i_i + u_i -->| | + +---------+ +``` +Outer loop inputs/outputs: + +``` + +---------+ +i_r_meas -->| | +i_i_meas -->| Outer |--> i_r_ref +u_r_meas -->| Loop |--> i_i_ref +u_i_meas -->| | + +---------+ +``` +""" +function InverterCS(inner, outer; name=Symbol(string(outer.name)*"_"*string(inner.name))) + # @assert BlockSpec([:u_i, :u_r, :i_ref_r, :i_ref_i], [:i_r, :i_i])(inner) "Inner ctrl loop :$(inner.name) does not meet expectation." + # @assert BlockSpec([], [:i_ref_r, :i_ref_i]; in_strict=false)(outer) "Outer ctrl loop :$(outer.name) does not meet expectation." + + outerinputs = ModelingToolkit.getname.(outer.inputs) + :u_r ∈ outerinputs && @error "Outer shoud use :u_meas_r instead of :u_r" + :u_i ∈ outerinputs && @error "Outer shoud use :u_meas_i instead of :u_i" + :i_r ∈ outerinputs && @error "Outer shoud use :i_meas_r instead of :i_r" + :i_i ∈ outerinputs && @error "Outer shoud use :i_meas_i instead of :i_i" + + meas = MIComponents.UIMeas() + sys = IOSystem(:autocon, [outer, inner, meas]; + outputs=[inner.i_r, inner.i_i], + globalp=[:u_r, :u_i], + name) + closed = connect_system(sys) + + @assert BlockSpec([:u_i, :u_r], [:i_r, :i_i])(closed) "Closed loop does not match expectation! $closed" + + return closed +end + +#### +#### Outer Control Loops +#### +""" + DroopControl(; params...) + +Return block for droop control outer control. + + +``` + ω_ref, V_ref, P_ref, Q_ref, K_P, K_Q, τ_P, τ_Q + v + +-----------------------------------------------------+ + | P_meas = u_meas_r*i_meas_r + i_meas_i*u_meas_i | +i_r_meas -->| Q_meas = -u_meas_r*i_meas_i + i_meas_r*u_meas_i | +i_i_meas -->| d/dt P_fil = (P_meas - P_fil) / τ_P |--> u_r_ref +u_r_meas -->| d/dt Q_fil = (Q_meas - Q_fil) / τ_P |--> u_i_ref +u_i_meas -->| d/dt δ = ω_ref - K_P*(-P_ref + P_fil) | + | V = V_ref - K_Q*(-Q_ref + Q_fil) | + +-----------------------------------------------------+ +``` + +""" +function DroopControl(; params...) + @named Pfil = IOComponents.LowPassFilter(; τ=:τ_P, input=:P_meas, output=:P_fil) + @named Pdroop = IOComponents.DroopControl(; x_ref=:P_ref, K=:K_P, u_ref=:ω_ref, u=:ω, x=:P_fil) + Psys = @connect Pfil.P_fil => Pdroop.P_fil outputs=:remaining name=:Psys + + @named Qfil = IOComponents.LowPassFilter(; τ=:τ_Q, input=:Q_meas, output=:Q_fil) + @named Qdroop = IOComponents.DroopControl(; x_ref=:Q_ref, K=:K_Q, u_ref=:V_ref, u=:V, x=:Q_fil) + Qsys = @connect Qfil.Q_fil => Qdroop.Q_fil outputs=:remaining name=:Qsys + + @named PQmeas = IOComponents.Power(; u_r=:u_meas_r, u_i=:u_meas_i, + i_r=:i_meas_r, i_i=:i_meas_i, + P=:P_meas, Q=:Q_meas) + equations(PQmeas) + + + refgen = MIComponents.VRefGen() + + @named droop = IOSystem(:autocon, [PQmeas, Psys, Qsys, refgen], outputs=:remaining) + con = connect_system(droop) + + return replace_vars(con, params) +end + +""" + Synchronverter(; params...) + +Synchronverter outer loop control for a voltage controlled converter. + +``` + ω0, V_ref, P_ref, Q_ref, Dp, Kq, Kv + v + +---------------------------------------------------------------------+ + | V = sqrt(u_meas_r^2 + u_meas_i^2) | + | Q = -sqrt(3/2) * MfIf * (ω0 + ω) | + | * (-sin(θ)*i_meas_r + cos(θ)*i_meas_i) | +i_r_meas -->| Te = sqrt(3/2) * MfIf * (cos(θ)*i_meas_r + sin(θ)*i_meas_i) | +i_i_meas -->| d/dt ω = 1/J * P_ref/ω0 - Dp*ω - Te, |--> u_r_ref +u_r_meas -->| d/dt θ = ω |--> u_i_ref +u_i_meas -->| d/dt MfIf = (Q_ref - Q + Dq*(V_ref - V))/Kv], | + | u_ref_r = sqrt(3/2) * MfIf * (ω0 + ω) * cos(θ) | + | u_ref_i = sqrt(3/2) * MfIf * (ω0 + ω) * sin(θ) | + +---------------------------------------------------------------------+ +``` +""" +function Synchronverter(; params...) + # Frequency Loop + @variables t ΔT(t) ω(t) θ(t) + @parameters Te(t) P_ref Dp ω0 J + dt = Differential(t) + @named floop = IOBlock([ΔT ~ P_ref/ω0 - Dp*ω - Te, + dt(ω) ~ 1/J * ΔT, + dt(θ) ~ ω], + [Te], + [ω, θ]) + + # Voltage Loop + @variables t ΔQ(t) MfIf(t) V(t) + @parameters u_meas_r(t) u_meas_i(t) Q(t) Q_ref Dq V_ref Kv + @named vloop = IOBlock([V ~ sqrt(u_meas_r^2 + u_meas_i^2), + ΔQ ~ Q_ref - Q + Dq*(V_ref - V), + dt(MfIf) ~ ΔQ/Kv], + [Q, u_meas_r, u_meas_i], + [MfIf]) + BlockSystems.WARN[] = true + + # main model + @variables t Te(t) u_ref_r(t) u_ref_i(t) Q(t) + @parameters i_meas_r(t) i_meas_i(t) MfIf(t) ω(t) θ(t) ω0 + @named machine = IOBlock([Te ~ sqrt(3/2) * MfIf * (cos(θ)*i_meas_r + sin(θ)*i_meas_i), + u_ref_r ~ sqrt(3/2) * MfIf * (ω0 + ω) * cos(θ), + u_ref_i ~ sqrt(3/2) * MfIf * (ω0 + ω) * sin(θ), + Q ~ -sqrt(3/2) * MfIf * (ω0 + ω) * (-sin(θ)*i_meas_r + cos(θ)*i_meas_i)], + [ω, θ, i_meas_r, i_meas_i, MfIf], + [u_ref_r, u_ref_i, Q, Te]) + + @named syncvert = IOSystem(:autocon, [floop, vloop, machine], outputs=:remaining, globalp=[:ω0]) + con = connect_system(syncvert) + return replace_vars(con, params) +end + +""" + FixedVoltage(; params...) + +Outer loop control for a voltage source which tries to keep voltage constant (slack like behavior). + + +``` + u_d_fix, u_i_fix + v ++-------------------+ +| | +| u_r_ref = u_d_fix |--> u_r_ref +| u_i_ref = u_i_fix |--> u_i_ref +| | ++-------------------+ +``` +""" +function FixedVoltage(; params...) + @variables t u_ref_r(t) u_ref_i(t) + @parameters u_fix_r u_fix_i + blk = IOBlock([u_ref_r ~ u_fix_r, + u_ref_i ~ u_fix_i], + [], [u_ref_r, u_ref_i], + name=:FixedU) + replace_vars(blk, params) +end + +""" + PLLCurrent(; params...) + +Outer loop control for a current source inverter, which employs a dynamic PLL to estimate +the phase and frequency of the voltage at PCC. The single parameter, `i_mag_ref` is interpeted +as current reference in `d` component (pure active power injection with fixed current). +""" +function PLLCurrent(; params...) + pll = MIComponents.JuanPLL(; Kp=250, Ki=1000, u_r=:u_meas_r, u_i=:u_meas_i) + + pol2cart = IOComponents.Polar2Cart(;arg=:δ_pll, mag=:i_mag_ref, x=:i_ref_r, y=:i_ref_i) + + con = IOSystem([pll.δ_pll => pol2cart.δ_pll], + [pll, pol2cart]; + outputs = [pol2cart.i_ref_r, pol2cart.i_ref_i], + name=:PLL_current) |> connect_system + con = make_iparam(con, :i_mag_ref) +end + +""" + ConstantPower(; params...) + +Outer loop control for a current source inverter, which measures the voltage at PCC +and commands the `i_ref_r` and `i_ref_i` such that the desired power at PCC is achived. +""" +function ConstantPower(; params...) + @variables t i_ref_r(t) i_ref_i(t) + @parameters u_meas_r(t) u_meas_i(t) P_ref Q_ref + blk = IOBlock([i_ref_r ~ ( P_ref*u_meas_r + Q_ref*u_meas_i)/(u_meas_r^2+u_meas_i^2), + i_ref_i ~ ( P_ref*u_meas_i - Q_ref*u_meas_r)/(u_meas_r^2+u_meas_i^2)], + [u_meas_r, u_meas_i], [i_ref_r, i_ref_i], + name = :ConstantPower) +end + +""" + FiltConstantPower(; params...) + +Outer loop control for a current source inverter, which measures the voltage at PCC +and commands the `i_ref_r` and `i_ref_i` such that the desired power at PCC is achived. +In contrast to the `ConstantPower`-outer loop, the current references follow a PT1 behavior +with time constant τ. +""" +function FiltConstantPower(; params...) + lpf_r = IOComponents.LowPassFilter(input=:u_meas_r, output=:u_filt_r) + lpf_i = IOComponents.LowPassFilter(input=:u_meas_i, output=:u_filt_i) + + @variables t i_ref_r(t) i_ref_i(t) + @parameters u_filt_r(t) u_filt_i(t) P_ref Q_ref + blk = IOBlock([i_ref_r ~ ( P_ref*u_filt_r + Q_ref*u_filt_i)/(u_filt_r^2+u_filt_i^2), + i_ref_i ~ ( P_ref*u_filt_i - Q_ref*u_filt_r)/(u_filt_r^2+u_filt_i^2)], + [u_filt_r, u_filt_i], [i_ref_r, i_ref_i], + name = :PowerCalc) + sys = IOSystem(:autocon, [lpf_r, lpf_i, blk]; globalp=[:τ], outputs=:remaining) + con = connect_system(sys) + if !isempty(params) + con = replace_vars(con, params) + end + con +end + +#### +#### Inner Control Loops +#### +""" + PT1Source(; params...) + +Create Voltage source which follows angle directly but +amplitude with a PT1-lag. +""" +function PT1Source(; params...) + @variables t A(t) u_r(t) u_i(t) + @parameters τ u_ref_r(t) u_ref_i(t) i_r(t) i_i(t) + dt = Differential(t) + blk = IOBlock([dt(A) ~ 1/τ*(√(u_ref_r^2 + u_ref_i^2) - A), + u_r ~ A/√(u_ref_r^2 + u_ref_i^2) * u_ref_r, + u_i ~ A/√(u_ref_r^2 + u_ref_i^2) * u_ref_i], + [u_ref_r, u_ref_i, i_r, i_i], + [u_r, u_i], + name=:PT1Src, + warn=false) + + if !isempty(params) + blk = replace_vars(blk, params) + end + return blk +end + +""" + IdealSource() + +Ideal voltage source which follows the reference directly. +""" +function IdealSource(; params...) + @variables t u_r(t) u_i(t) + @parameters u_ref_r(t) u_ref_i(t) i_r(t) i_i(t) + Vsource = IOBlock([u_r ~ u_ref_r, + u_i ~ u_ref_i], + [u_ref_r, u_ref_i, i_r, i_i], + [u_r, u_i], + name=:VSrc, + warn=false) +end + +""" + IdealCSource() + +Ideal current source which follows the reference directly. +""" +function IdealCSource(; params...) + @variables t i_r(t) i_i(t) + @parameters i_ref_r(t) i_ref_i(t) u_r(t) u_i(t) + Vsource = IOBlock([i_r ~ i_ref_r, + i_i ~ i_ref_i], + [i_ref_r, i_ref_i, u_r, u_i], + [i_r, i_i], + name=:CSrc, + warn=false) +end + + +""" + CascadedVoltageControl(;pr=true, ff=true) + +Cascaded voltage controler over a LC filter, which controlls the output +voltage using 2 levels of cascaded dq controllers. +""" +function CascadedVoltageControl(;pr=true, ff=true) + LC = MIComponents.LC() + CC1 = if pr + MIComponents.CC1_PR() + else + MIComponents.CC1() + end + VC = MIComponents.VC() + vcinner = IOSystem(:autocon, + [LC,CC1,VC]; + globalp=[:ω0, :Rf, :Rg, :Lf, :Lg, :C, :i_g_r, :i_g_i], + name=:VCFOM, + outputs=[:LC₊V_C_r, :LC₊V_C_i], + autopromote=false) |> connect_system + + vcinner = replace_vars(vcinner; + :LC₊V_C_r => :u_r, + :LC₊V_C_i => :u_i, + :VC₊V_C_ref_r => :u_ref_r, + :VC₊V_C_ref_i => :u_ref_i, + :i_g_r => :i_r, + :i_g_i => :i_i) + params = copy(MIParameters) + params[:CC1₊F] = false + params[:VC₊F] = ff + vcinner = replace_vars(vcinner, params; warn=false) +end + +""" + CascadedVoltageControlCS(;pr=true, ff=true) + +Cascaded voltage controler over a LCL filter, which controlls the filter +voltage using 2 levels of cascaded dq controllers. + +Acts as a current source because of the second L in the filter. +""" +function CascadedVoltageControlCS(;pr=true, ff=true) + LCL = MIComponents.LCL() + CC1 = if pr + MIComponents.CC1_PR() + else + MIComponents.CC1() + end + VC = MIComponents.VC() + vcinner = IOSystem(:autocon, + [LCL,CC1,VC]; + globalp=[:ω0, :Rf, :Rg, :Lf, :Lg, :C, :V_g_r, :V_g_i], + name=:VCFOM, + outputs=[:LCL₊i_g_r, :LCL₊i_g_i], + autopromote=false) |> connect_system + + vcinner = replace_vars(vcinner; + :V_g_r => :u_r, + :V_g_i => :u_i, + :VC₊V_C_ref_r => :u_ref_r, + :VC₊V_C_ref_i => :u_ref_i, + :LCL₊i_g_r => :i_r, + :LCL₊i_g_i => :i_i) + params = copy(MIParameters) + params[:CC1₊F] = false + params[:VC₊F] = ff + vcinner = replace_vars(vcinner, params; warn=false) +end + +""" + CascadedCurrentControl(;pr=true) + +Cascaded current controler over a LCL filter, which controlls the output +current using 3 levels of cascaded dq controllers. +""" +function CascadedCurrentControl(; pr=true) + LCL = MIComponents.LCL() + CC1 = if pr + MIComponents.CC1_PR() + else + MIComponents.CC1() + end + VC = MIComponents.VC() + CC2 = MIComponents.CC2() + ccinner = IOSystem(:autocon, + [LCL,CC1,VC,CC2]; + globalp=[:ω0, :Rf, :Rg, :Lf, :Lg, :C, :V_g_r, :V_g_i], + name=:CCFOM, + outputs=[:LCL₊i_g_r, :LCL₊i_g_i], + autopromote=false) |> connect_system + + ccinner = replace_vars(ccinner; + :V_g_r => :u_r, + :V_g_i => :u_i, + :CC2₊i_g_ref_r => :i_ref_r, + :CC2₊i_g_ref_i => :i_ref_i, + :LCL₊i_g_r => :i_r, + :LCL₊i_g_i => :i_i) + + ccparams = copy(MIParameters) + ccparams[:CC1₊F] = false + ccparams[:VC₊F] = false + ccparams[:CC2₊F] = true + + replace_vars(ccinner, ccparams; warn=false) +end diff --git a/src/PowerDynamics.jl b/src/PowerDynamics.jl index b6123be5..966240d5 100644 --- a/src/PowerDynamics.jl +++ b/src/PowerDynamics.jl @@ -49,6 +49,12 @@ include("IONodes/IOComponents.jl") include("IONodes/GFI_MTK.jl") include("IONodes/MTK_Load.jl") +module ModularInverter +include("IONodes/MIComponents.jl") +include("IONodes/ModularInverter.jl") +include("IONodes/LTI.jl") +end + # all line types include("lines/AbstractLine.jl")