GSoC July Recap

Hello everyone, a month has passed since the last update.

Here are the various things I have done.

Stability conditions

The PR #28021(merged) has finally been merged. Now, it is possible for SymPy to calculate stability conditions for a continuous-time system (i.e., compute conditions to ensure that every root of a polynomial has a negative real part). The great thing is that these conditions are processed without division, which is very important for the symbolic case, since we don’t have to deal with degenerate cases.

In[1]: p = Poly(symbols("a:f"), x)
In[2]: for c in negative_real_part_conditions(p, x): print(c)
a*b > 0
b**2*(-a*d + b*c) > 0
(-a*d + b*c)**2*(-b**2*(-a*f + b*e) + b*d*(-a*d + b*c)) > 0
(-b**2*(-a*f + b*e) + b*d*(-a*d + b*c))**2*(-b*f*(-a*d + b*c)**3 + (-a*d + b*c)*(-a*f + b*e)*(-b**2*(-a*f + b*e) + b*d*(-a*d + b*c))) > 0
b*f*(-a*d + b*c)*(-b**2*(-a*f + b*e) + b*d*(-a*d + b*c))*(-b*f*(-a*d + b*c)**3 + (-a*d + b*c)*(-a*f + b*e)*(-b**2*(-a*f + b*e) + b*d*(-a*d + b*c)))**3 > 0

In my last update, I also cited two issues.

The first one, #28176(open), is about the simplification of the conditions, since there are a lot of repeating factors. I opened a PR with possible ideas for solving that issue #28265(draft).

Note

I wrote improvements for the cases in which the polynomial is represented by the class Expr.
@oscarbenjamin helped me a lot writing the version for Poly.
Expr and Poly are going to have different purposes.
The first one should be used for cases in which there are a lot of symbols and the computation speed is important.
The second one is very useful because it can reach the best simplification of the inequalities.

In this example we can easily compare the two and see how much better it is than the first implementation:

In [1]: p = Poly(symbols("a:f"), x)
In [2]: for c in negative_real_part_conditions(p, x): print(c) # algorithm for Poly
a*b > 0
a*d - b*c > 0
a*b**2*f - a*b*d**2 - b**3*e + b**2*c*d > 0
a**2*f**2 - 2*a*b*e*f - a*c*d*f + a*d**2*e + b**2*e**2 + b*c**2*f - b*c*d*e > 0
b*f > 0
In [3]: for c in negative_real_part_conditions(p, x, domain = EXRAW): print(c) # algorithm for expressions
a*b > 0
-a*d + b*c > 0
b*(b*(a*f - b*e) - d*(a*d - b*c)) > 0
b*(a*d - b*c)*(f*(a*d - b*c)**2 + (a*f - b*e)*(b*(a*f - b*e) - d*(a*d - b*c))) > 0
b*f > 0

The other issue I opened is #28180(open), which is about performance issues for systems with a lot of symbols. As I said in the last update, the problem was in the characteristic polynomial calculation of the state-space A matrix. I found that calculating the characteristic polynomial using the domain EXRAW, which emulates Expr, fixes the problem. However, before writing the best interface, I need to finish the previous PR.

Discrete-time systems

Other minor changes have been made to the PR #28115(open) about discrete-time transfer functions and I added a usage example to the tutorials.

Locally, I finished writing the discrete counterpart of StateSpace. I’m just waiting for the PR about transfer functions to be merged to upload these changes.

Other work

Some time ago, I opened PR #27986(open), which is about the decomposition of state spaces in the observability and controllability form and any other generic transformation using the method apply_similarity_transform().
It also fixes the way in which the observable subspace of a system is calculated.

Working in that PR, I also found that the function rank, for symbolic matrices, doesn’t count every case of symbol, so, for specific combinations of symbols, the output will be incorrect.

After some discussions, I was convinced that this output is acceptable, since it covers almost all possibilities, but there is a need for a method to find the degenerate cases.

I opened a PR for that purpose #28232(open). Finding a way to solve this issue will be very important, because it will be possible to make the functions is_observable and is_controllable work in all cases.