So to describe the process from the outer layers of the code to the inner layers.
First major change. The way the current loops of the system and their interaction is decribed. A loop will have the beginning and the ending node as the same and will be a continuous list. The problem is if loop manipulations have to be done - that is differences between loops have to be found in terms of branch segments, it is a little difficult.
So all loops and common branch segments are defined in terms of branches that are elements between two nodes. This is the smallest denomination. Also, every one of them has a direction. If a branch is in a loop and is the diagonal element in system_loops, it has the default direction as "forward". If it shows the interaction between loops, it was have the necessary direction.
So before the ODE solver, a pre-processing has been done to break up all elements in system_loops into lists of branches. This is the code (click on "View Raw" below the code box to see it in a new window):
The next major change is in the ODE solver. Three steps here:
- Check for stiff equations - L/R ratio and isolate every stiff branch to a single equation by row operations.
- With each row operation on the matrices A, B and E, perform the same row operations on the system_loops. So remove elements in common and add new elements.
- Since, the objective has been to isolate a stiff loop and approximate it to a static equation without a d/dt, the interactions with other loops has to be eliminated and these interactions will have to be readjusted.
Now in detail. Here is the code (click on "View Raw" below the code box to see it in a new window):
Check if a diagonal element has a small L/R ratio (time constant). Now this judgement is purely a guess and I have chosen it as 0.1*SimulationStep. But may have to be changed. If that is the case, this equation has to be isolated. So, the interaction between this loop and other loops, must be eliminated and first we start with the branch that is causing the problems. So look for off diagonal elements in the same row with small time constants.
Suppose an off diagonal element [c1, c2] has a small time constant. This means the branch between loops c1 and c2 has a small time constant. And this therefore means that loop c2 will also have this branch and will also have a small time constant. So, the way to go is to check which loop has a smaller time constant - c1 or c2 by looking at the diagonal elements c1,c1 and c2,c2. The smaller one remains as it is and the larger one is manipulated to remove this stiff branch. This is by changing row c2 to c2-c1.
Look through the entire matrix this way. The objective is then that only the stiff branches will be in their own isolated loops and the stiff elements removed from all other loops.
Next - loop manipulations by changing the elements.
The reason for doing this is that when the branch currents are measured, they are measured with respect to the branches present in the loops. So with the above row operations, the loops have changed, so the actual loop descriptions must also change. So when we do c2=c2-c1, we need to find out which elements in c1 are present in c2 and delete them from c2 and which elements are present in c1 but not in c2 and add them reversing their direction.
Last step - readjust the matrices. Here is the code (click on "View Raw" below the code box to see it in a new window):
This was actually a bit tricky. The first part is the stiff equation must be connected to an input source because interactions with all other loops has been eliminated. To eliminate the interaction with other loops, the off diagonal elements in matrix A are set to zero in the stiff loop.
Then check if any of the elements in the row of matrix B is nonzero. If yes, leave it as it is. If not, add the stiff loop with another nonstiff loop that is connected to the input. This will cause one of the elements of B matrix in the stiff row as non zero. At the same time, some of the off diagonal elements in the stiff row will now be nonzero as these are the resistances in the nonstiff loop with which the addition occurred. These resistances are essentially the resistances between the source and the stiff branch. So add all those to the diagonal element in A and then get rid of them. So now the diagonal element of A contains the entire resistance as a loop between the source and the stiff branch and is also linked to the source. So the stiff dynamical equation has been approximated to a static isolated equation.
Next, remove all the elements in matrix E - eliminating all d/dt terms and making it a static equation. The last part is to account for the interaction elements that still linger between the stiff loops and non stiff loops. So if Z is the stiff loop and X and Y are nonstiff loops that interact with Z in that they have nonzero elements in the same column as Z, this means that X and Y will interact with each other. This interaction needs to be captured as we want to isolated loop Z. So look for the elements and move them to the offdiagonal entries between X and Y with the appropriate sign.