more abstract statespaces

[baggepinnen]

No description provided.

[JuliaControlBot]

Something failed when generating plots. See the log at https://github.com/JuliaControl/ControlExamplePlots.jl/runs/875149431?check_suite_focus=true for more details.

[baggepinnen]

I'm a bit puzzled by
https://github.com/JuliaControl/ControlSystems.jl/pull/521/checks?check_run_id=2671828259#step:6:272
@olof3 it seems you touched this line last, what is this supposed to test? It is actually supposed to not work? In that case, @test_broken is the wrong thing to use.

[olof3]

I'm a bit puzzled by https://github.com/JuliaControl/ControlSystems.jl/pull/521/checks?check_run_id=2671828259#step:6:272
@olof3 it seems you touched this line last, what is this supposed to test?

It's supposed to test feedback between a system and a number, when there are direct terms in the system.

I guess it was failing previously, although we want it to pass. That is the use case for @test_broken, right?

The improved, more involved, feedback method that I have added, although not very much integrated, handles this case. The older, simpler feedback method does not handle direct terms in both systems. So it seems like your changes has changed how the dispatch is done.

Good that it works, but perhaps we would like a simpler, more specialized, and faster feedback method for this standard use case.

I have been looking into this, but I think we need to discuss how much it makes sense to specialize.

[baggepinnen]

I guess it was failing previously, although we want it to pass. That is the use case for @test_broken, right?

Correct, but the line below indicates that positiv feedback should provide the answer that we tested as broken on this line, i.e., both these can't pass?

@test_broken feedback(G3, 1) == ss(-1.5, 0.5, 0.5, 0.5) # Old feedback method
@test feedback(G3, 1, pos_feedback=false) == ss(-1.5, 0.5, 0.5, 0.5)

[baggepinnen]

Ah sorry, yes they should both pass. I'll change the test to not broken anymore :)

[baggepinnen]

In fact, the change made all calls to feedback dispatch to the general method as opposed to the old feedback(::StateSpace, ::StateSpace). All tests seem to pass except this one, so I guess the general method you implemented @olof3 is working well. In that case, we can just remove the simple method and let the advanced method do all the work?

[baggepinnen]

Good that it works, but perhaps we would like a simpler, more specialized, and faster feedback method for this standard use case.
I have been looking into this, but I think we need to discuss how much it makes sense to specialize.

The simple feedback specialized on StateSpace. If one actually wants performance, using StaticArrays is going to be required and the simple feedback does not help (unless adopted for HeteroStateSpace with static arrays). I think it's fine to only have the special method. Some timings, where simple_feedback is the old implementation and the systems have 10 poles and 1 y/u

@btime simple_feedback($s1,$s2)
# 3.327 μs (17 allocations: 10.16 KiB)  with normal arrays
# 262.156 ns (0 allocations: 0 bytes) modified to work with static arrays

@btime feedback($s1,$s2) # 6.752 μs (60 allocations: 14.39 KiB) New, more general feedback

[codecov]

Codecov Report

❗ No coverage uploaded for pull request base (dev@490ba4d). Click here to learn what that means.
The diff coverage is n/a.

@@          Coverage Diff           @@
##             dev     #521   +/-   ##
======================================
  Coverage       ?   84.87%           
======================================
  Files          ?       31           
  Lines          ?     3187           
  Branches       ?        0           
======================================
  Hits           ?     2705           
  Misses         ?      482           
  Partials       ?        0           

---

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 490ba4d...2fc8340. Read the comment docs.

[olof3]

If one actually wants performance, using StaticArrays is going to be required and the simple feedback does not help (unless adopted for HeteroStateSpace with static arrays). I think it's fine to only have the special method. Some timings, where simple_feedback is the old implementation and the systems have 10 poles and 1 y/u

I am a little bit worried that the general feedback method is a little bit too involved . Also, 50% faster for free (less for larger system, more for smaller ones), without having to worry about static arrays seems worthwhile (although I think the performance of the general version could be improved, at the cost of readability) .

The alternative would be to update the simple feedback version to handle direct terms (and probably also positive feedback), although the additional code required is a bit off-putting. I am not sure if there is a simple way to get the dispatch to work, but one could consider renaming the general version to feedback_general or something.

All in all, I don't have that much of a preference. What do you think @mfalt ?

[mfalt]

If we could dispatch to a faster method for SS that would be nice, otherwise there is little reason to keep it. If the faster method only works for some cases I guess we could settle for something like this: (psuedo code:)

function feedback(s1::SS, S2::SS)
 if one_is_zero(s1.D,s2.D)
   feedback_simple(s1,s2)
else
 feedback_general(s1,s2)
end
end

Edit: this might have some potential downsides, but since we know all the types it would be type-stable. It might be worth it.

[olof3]

If we could dispatch to a faster method for SS

with SS you also mean HeteroStateSpace? There should be no difference in the implementation. I think it was just the old type signature that was limiting.

one_is_zero(s1.D,s2.D)

This would not be the main difference. The general version checks and exploits this.

The drawback with the generic version is that it allocates a bit more and does unnecessary multiplications with matrices with 0-length dimensions.

[baggepinnen]

I'm usually for more performant implementations, but in this case it does not feel worth it to me. The call to feedback is probably not going to be the bottleneck in any performance sensitive code, and if you are concerned with the 2x performance hit you can actually get a 10x improvement by not using the regular state-space with heap-allocated arrays.

[olof3]

Yes, only performance wise I agree it's not worth it.
Seems like a little bit of an overwhelming method for a first time user though, and in the docs it seems like it could be difficult to get the main idea across for the basic use case while still capturing the advanced capabilities. Always nice to look at the code one is using and that it is easy to understand.

Perhaps a nonstandard doc string would be the solution. With simple use case first, and then separately, below, the current one.

[baggepinnen]

Note that this pr doesn't add nor remove any docstrings, it simply removes one method which was strictly less capable that another existing method. I'm not even sure that maintaining separate implementations for educational purposes is going to be easier for new users, having to figure out which of several implementations of the same function will be called and why is likely much more confusing than just having one generic method.

Docstring can certainly be improved though. The docstring for feedback is probably the one from this package I read the most often and scratch my head about.

[baggepinnen]

I added a block diagram to the docstring in 1658340
see if it makes sense and if it's actually a correct depiction of the signals, I wasn't completely sure :P

[olof3]

I made a small further stab at the docstring and added a diagram for the basic use case as well. Feel free toimprove them further. I guess these methods (and docstrings) should really be for any LTISystem?

[baggepinnen]

LGTM

I guess these methods (and docstrings) should really be for any LTISystem?

The general version currently only works for ss, but the block diagram for the simple usage applies to transfer functions as well. This usecase already has a small block diagram associated with it in the other docstring
https://github.com/JuliaControl/ControlSystems.jl/pull/521/files#diff-40aa26d41b0c99b2f88018c7c8e914ee958ce3cb5e7fb1c11902eef04f06627fR206
so users looking at ?feedback will see both.

[olof3]

The general version currently only works for ss, but the block diagram for the simple usage applies to transfer functions as well. This usecase already has a small block diagram associated with it in the other docstring

Hopefully I'll get around to fix the general version to support DelayLTISystem. Possibly also for transfer functions.

so users looking at ?feedback will see both.

That is not ideal. Not sure what the best way is for getting maintainable non-duplicated docs? Does this happen if the methods for the different types are grouped together?

[baggepinnen]

With the latest commit, ?feedback shows the following

help?> feedback
search: feedback feedback2dof

  feedback(L)
  feedback(P1,P2)

  Returns L/(1+L) or P1/(1+P1*P2)

  ─────────────────────────────────────────────────────────────────────────────────────────────────────────

  feedback(sys)
  feedback(sys1,sys2)

  Forms the negative feedback interconnection

  >-+ sys1 +-->
    |      |
   (-)sys2 +

  If no second system is given, negative identity feedback is assumed

  ─────────────────────────────────────────────────────────────────────────────────────────────────────────

  feedback(sys1::AbstractStateSpace, sys2::AbstractStateSpace)
  feedback(sys1::AbstractStateSpace, sys2::AbstractStateSpace;
           U1=:, Y1=:, U2=:, Y2=:, W1=:, Z1=:, W2=Int[], Z2=Int[],
           Wperm=:, Zperm=:, pos_feedback::Bool=false)

  Basic use feedback(sys1, sys2) forms the feedback interconnection

             ┌──────────────┐
  ◄──────────┤     sys1     │◄──── Σ ◄──────
      │      │              │      │
      │      └──────────────┘      -1
      │                            |
      │      ┌──────────────┐      │
      └─────►│     sys2     ├──────┘
             │              │
             └──────────────┘

  Advanced use feedback also supports more flexible use according to the figure below

                ┌──────────────┐
        z1◄─────┤     sys1     │◄──────w1
   ┌─── y1◄─────┤              │◄──────u1 ◄─┐
   │            └──────────────┘            │
   │                                        α
   │            ┌──────────────┐            │
   └──► u2─────►│     sys2     ├───────►y2──┘
        w2─────►│              ├───────►z2
                └──────────────┘

  U1, W1 specifies the indices of the input signals of sys1 corresponding to u1 and w1 Y1, Z1 specifies the
  indices of the output signals of sys1 corresponding to y1 and z1 U2, W2, Y2, Z2 specifies the
  corresponding signals of sys2

  Specify Wperm and Zperm to reorder the inputs (corresponding to [w1; w2]) and outputs (corresponding to
  [z1; z2]) in the resulting statespace model.

  Negative feedback (α = -1) is the default. Specify pos_feedback=true for positive feedback (α = 1).

  See Zhou, Doyle, Glover (1996) for similar (somewhat less symmetric) formulas.
  ≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡

  ─────────────────────────────────────────────────────────────────────────────────────────────────────────

  feedback(sys)
  feedback(sys1, sys2)

  Forms the negative feedback interconnection

  >-+ sys1 +-->
    |      |
   (-)sys2 +

  If no second system is given, negative identity feedback is assumed

There is clearly some redundancy there.

Splitting the docstring like this

  feedback(sys1::AbstractStateSpace, sys2::AbstractStateSpace)
  feedback(sys1::AbstractStateSpace, sys2::AbstractStateSpace;
           U1=:, Y1=:, U2=:, Y2=:, W1=:, Z1=:, W2=Int[], Z2=Int[],
           Wperm=:, Zperm=:, pos_feedback::Bool=false)

can actually be confusing since the second method would overwrite the first, why is the first row there at all? All keyword arguments have default values so that the first method exists follows from the second signature.

[olof3]

There is clearly some redundancy there.

Yes, it would definitely be good to improve this. I am not sure how to do it in a good way though. I think you know better what works well. Perhaps a method feedback(::LTISystem, ::LTISystem) and putting the docs for the simple case there, removing all the redundant ones? Only keeping the docs for the complicated one for AbstractStateSpace?

why is the first row there at all? All keyword arguments have default values so that the first method exists follows from the second signature.

I agree, that was perhaps not the clearest way of doing it. But I think the basic use be illustrated in the docs somewhere, although perhaps just mentioning, as it is done in the doc text, would be sufficient. Feel free to improve it.