Solving a system of linear equations with computer

[cb951303]

Hello, I'm a senior mechanical engineering student.

I'm trying to write an application that plots the system response for a given block diagram in frequency domain just for fun. Take the system shown on the attachment for example.

What I have in mind is to

• find the Y/U transfer function with block diagram reduction (simplification)
• find it's differential equation
• separate it to first order ODEs
• solve the ODEs and plot Y(t)

I'm aware that this is not the right forum for this but bear with me for a second.
While thinking how to make a computer do block diagram reduction I've noticed that a block diagram is nothing but a system of linear equations.
For example, the block diagram in the attachment can be represented with 3 linear equations with 3 unknowns such as this:

$$C=\frac{B}{E}$$
$$G=\frac{Y}{B}$$
$$E=U - Y$$

The problem is, that the terms of these equations are not numbers but transfer functions such as:

$$G = \frac{1}{2s + 1}$$

Since numerical libraries like LAPACK and NumPy only accepts numeric values for the coefficient matrix how can I make a computer solve this system?

Thank you in advance.

 

[Antiphon]

The computer solves the problem for a single value of complex s.

The easy way is to turn it into a circuit and use spice. That will integrate the matrix in time and plot the system response to arbitrary excitation.

 

[cb951303]

The computer solves the problem for a single value of complex s.

can you elaborate what you mean.

The easy way is to turn it into a circuit and use spice. That will integrate the matrix in time and plot the system response to arbitrary excitation.

My aim is to write a program that does this not use simulations like SPICE otherwise I could plot it in 5 seconds with Matlab/Simulink :)

 

[Filip Larsen]

Perhaps you can get something useful by searching for "transfer function to state space", which amongst other hits give:

https://www.physicsforums.com/showthread.php?t=304447
https://www.physicsforums.com/showthread.php?t=236737

Alternatively, you could grab a textbook on control theory that also describe how to model in state space (most new textbook do, I believe).

 

[cb951303]

Perhaps you can get something useful by searching for "transfer function to state space", which amongst other hits give:

https://www.physicsforums.com/showthread.php?t=304447
https://www.physicsforums.com/showthread.php?t=236737

Alternatively, you could grab a textbook on control theory that also describe how to model in state space (most new textbook do, I believe).

thanks for your input. I think I didn't make it clear.

you see the problem is not about converting transfer function to state space. I already know how to do it. The problem is to make a computer to find the transfer function of a system from its block diagram. it's more like a "finding the appropriate numerical method for the job" problem. and the best method I came up with is to present the block diagram as a system of linear equations and make the computer solve it. the solution of the system will be the transfer function. the problem is: libraries like LAPACK and NumPy only solves linear equation systems with a coefficiant matrix of numeric values. Since my coefficient matrix consists of fractions with s term and not numeric values, I can't use them. I'm asking how can I overcome this?

Thank you.

 

[Filip Larsen]

OK, I think I can see where you are driving at.

I'm not sure you can formulate a linear system that will allow you to combine all the block transfer functions using only numerical linear systems, but even if it can be done I doubt it would be worth the trouble.

I'm not expert on this, but if you ask me, I can see two fairly easy routes you can take.

One option is to make (or find) code that allow you to express and solve a linear system of rational polynomials. Seems you will need a library for general polynomial arithmetic and some code for doing, say, Gauss-Jordan, on such elements.

The other option would be to make recursive structural reductions of your block diagram, which also would require you to express and combine rational polynomials, but now only using a small number of reduction rules that can be hard-coded so you don't need a full library for general polynomial arithmetic.

Maybe others here have a better option?

 

[cb951303]

that's exactly the kind of feedback I'm looking for.

block diagram reduction seems a bit messy for a computer to do. maybe it can be achieved with a directed graph data structure to hold the block diagram but other than that I fear it would be very slow.

but I'll see what can I do with rational polynomials thank you.

meanwhile if anyone else have any feedback I would appreciate it.

 

[Lord Crc]

Not my field at all, so feel free to ignore :)

Are you saying that the solution you're after is or is not a linear combination of the "blocks variables" U, C and G? Ie something on the form
$$Y = \frac{CGU}{CG+1}$$

 

[trambolin]

use "sysic" command of MATLAB Robust Control Toolbox for arbitrary interconnections. Also see "append" and related commands.

 

[cb951303]

Not my field at all, so feel free to ignore :)

Are you saying that the solution you're after is or is not a linear combination of the "blocks variables" U, C and G? Ie something on the form
$$Y = \frac{CGU}{CG+1}$$

I'm saying that the solution (Y)is a linear combination of the blocks(C, G) and the input(U). So yes, your solution is true but by the form of it I think you found it by applying reduction? I'm trying to do the same with solving linear equations which is I think easier for computers.

I don't know if it's valid for all kinds of diagrams though. Maybe it won't work on complicated ones, or on multi input, multi output systems.

@trambolin

that's great, I think it's exactly what I'm trying to do *but* unless it's open source it won't help. my aim is to write it myself not to use already written ones.

thanks everyone

 

[Okefenokee]

You need to transform the transfer functions into ODE's, right?

Take this one for instance.

$$G = \frac{1}{2s + 1}$$

we can do a little manipulation to get a regular old fashioned diff e

$$G = \frac{Y}{B} = \frac{1}{2s + 1}$$

$$(Y \cdot 2s + Y) = B$$

$$(2(y'(t) + y(0)) + y(t)) = b(t)$$

You could use that to make state space blocks, by the way. You could make b(t) into a state and start building the matrices.

I'm not quite sure why you don't want to use state space equations. If you think it over, I bet that state space equations will do exactly whatever you want in a very robust way.

http://en.wikibooks.org/wiki/Control_Systems/State-Space_Equations" I use it for reference all the time. Maybe you could take state space equations from individual blocks, combine them, and then transform the result into the canonical form to get the final transfer function, or a diff-e if that's what you want.

 

[cb951303]

You need to transform the transfer functions into ODE's, right?

Take this one for instance.

$$G = \frac{1}{2s + 1}$$

we can do a little manipulation to get a regular old fashioned diff e

$$G = \frac{Y}{B} = \frac{1}{2s + 1}$$

$$(Y \cdot 2s + Y) = B$$

$$(2(y'(t) + y(0)) + y(t)) = b(t)$$

You could use that to make state space blocks, by the way. You could make b(t) into a state and start building the matrices.

I'm not quite sure why you don't want to use state space equations. If you think it over, I bet that state space equations will do exactly whatever you want in a very robust way.

http://en.wikibooks.org/wiki/Control_Systems/State-Space_Equations" I use it for reference all the time. Maybe you could take state space equations from individual blocks, combine them, and then transform the result into the canonical form to get the final transfer function, or a diff-e if that's what you want.

you misunderstood, am I really that unclear :) ??

the problem is not converting transfer function to ODEs it's to make the computer find the transfer function of a system from its block diagram.

state space and ODEs comes after I solve this.

PS: maybe I'm missing something if so can you pls elaborate?

 

[Okefenokee]

the problem is not converting transfer function to ODEs it's to make the computer find the transfer function of a system from its block diagram.

What I'm thinking is that you should look into doing it the other way around.

(1) Convert your individual block functions into state space equations-------> (2) Combine all the individual state space equations into a complete matrix for the whole system (I'll briefly show you how to do it in a bit) ---------> (3) Use linear algebra to convert the matrix into the canonical controller form (the form is described in that wiki-book)-------> (4) Read the total transfer function directly from the canonical form.

By the way, it sounds like you want to get T(s) so that you can simulate an output. You can do that directly from the state space equation without getting the transfer function.

Let's look at step one. That easy enough. Just use http://en.wikibooks.org/wiki/Control_Systems/State-Space_Equations#From_Transfer_Functions".

Now for step two. If you want to combine two state space blocks in series, http://en.wikibooks.org/wiki/Control_Systems/Block_Diagrams#Series_State_Space" If you want to add two state space blocks in parallel (an adder) then you simply stack the matrices like this:

stack the A matrices vertically (second one on the bottom)
stack the B matrices vertically (second one on the bottom)
combine the C matrices horizontally (second one on the right)
combine the D matrices horizontally (second one on the right)

Step 3. Use something similar to the steps of Gauss-Jordan elimination. Except, don't create a reduced-row-echelon form. Instead, get the matrices into canonical form. Make sure that you take any row manipulations of the A matrix and copy it as column manipulations in the C matrix.

Step 4. You should have a canonical form now. Just read the transfer function from it.Do you have any interest in Java programming? I've been working on a math library so that I can use it without worrying about licensing issues. It's not even in an alpha state right now but it's coming along.

 

[cb951303]

Now that makes sense!
I must say I had no idea about the step 2 (combining state space matrices) so I didn't even think about state space form. it should work much faster and I can finally use a library like CLAPACK for the matrix operations without having to write them myself as I would have to do with solving a linear system with a coefficient matrix of polynomial fractions.

thank you!

I only program in C but have you seen LAPACK? it's a well known math library written in FORTRAN. Even Matlab uses it. It also has C version (CLAPACK) converted from FORTRAN with f2c. So it should be pretty easy to create a JAVA extension with http://en.wikipedia.org/wiki/Java_Native_Interface" . Best of all LAPACK is licensed under BSD. Basically, you can do anything you want with it.

 

[trambolin]

I really recommend balanced realizations instead of the Gauss elimination for Step 3 in Okefenokee's post. It is numerically much better for reduced headache.