You are absolutely right, and that is a fair point. What we have learned throughout this process about numerical simulation has already been very valuable in itself (things we would never have had to think about using commercial tools). And from a research group perspective, the right attitude is exactly that: contributing to making the tools better rather than just being end users. That is precisely why we intend to keep working on this and contribute back to the ecosystem wherever we can.
Some progress: I simplified the coupled system, and now all works fine:
This is with a two-phase rectifier. I still have the warning:
┌ Warning: Did not converge after `maxiters = 100` substitutions. Either there \
│ is a cycle in the rules or `maxiters` needs to be higher.
└ @ Symbolics ~/.julia/packages/Symbolics/oZEAe/src/variable.jl:451
But this warning does not affect the result. I also created a flattened version of the same system (without components). With that version there is no warning.
I currently use this solver call:
sol = solve(prob, Rodas5P(), abstol = 1e-9, reltol = 1e-12, dense = false, saveat = 0.0001, maxiters = 10000000)
Another solver might work better, but this already works.
Update:
The original 3-phase example now works, too:
I had to switch to the RadauIIA5 solver:
sol = solve(prob, RadauIIA5(), abstol = 1e-8, reltol = 1e-8, saveat = 0.0001)
You can find the code at: GitHub - ufechner7/MTKCircuits.jl · GitHub
Update:
It still does not work for the original value of aero = Constant(k = -97). It does work for k= -9.7, though.
Great results! The waveforms look really clean, especially the 3-phase case.
Did you manage to pinpoint what was actually causing the issue? Was it something related to the topology of the system (e.g. how the diode switching interacts with the connections), or does it come down entirely to solver choice and configuration? I’m asking because the fact that RadauIIA5 works while other solvers get stuck suggests the problem might be more about stiffness handling during the switching events rather than a structural issue with the model itself, but curious if you dug deeper into it.
Yes, stiffness handling during switching events is the main problem. So the problem is not solved yet. It still fails if the generator torque and thus the voltage become too high.
UPDATE:
The simulation now also works with the original parameters.
For higher DC voltages (here 928 V), it is necessary to add a small (here 650 pF) capacitor in parallel to the diodes, and also a very small resistor (here 1 mOhm) in series.
Example:
include("examples/coupled_system2.jl)
The following solver parameters were used:
sol = solve(prob, RadauIIA5(κ = 0.005), abstol = 1e-7, reltol = 1e-8, saveat = 0.000025, maxiters = 5e6)
Update II:
Some of the examples did not work on Windows. Now they do. We had to add initial guesses. Furthermore, there is now a bin/install script and a script to create a system image.
Thanks! I was travelling at the time, but I’m keeping tabs on the topic. Since it has been marked resolved, I’ll chime in with my experience with audio circuits.
Some general tips:
- Successful simulation: ramp sources from zero, use known starting solution from SPICE, tune tolerances, switch solvers
- Sometimes helps: serial/parallel R for diodes, clamping exponentials
- Initialization: I haven’t seen ramping fail for any audio circuit yet
- Agents: use them for mechanical workflows which are not quite search-and-replace, I used to do this manually for years
- Switching: use the DiffEq event system, use known good runs for
affect!()state projections, change Newton’s method to Halley’s method as a get-out-of-jail-free card (it behaves better away from the root, which is where switching sends you)
I think PSpice/Simulink are less rigorous than the DAE-based MTK simulation. This can be an advantage if you are interested in a lax-but-close-enough result that always outputs something.
For non-pathological cases (e.g. floating ideal transformers would fail!), MTK can be made to work, but it’s going to annoy you with its rigor in plenty of large, practical circuits.
The DAE-based equivalent is not Simulink but Simscape.
Thank for the tips! They are really useful and give us a clearer picture of what to expect from MTK in this kind of application.
With all the contributions in this thread we have managed to get the diode-based circuits working reasonably well, especially rectifiers, although we are still exploring how well they scale to more complex topologies. On top of that we are starting to test IGBT-based systems, which add another layer of complexity since EMT simulations of power electronics typically involve carrier signals in the kHz range, we are curious to see how MTK handles that.
As @ufechner7 has set up a GitHub repository to collect these models, I plan to keep uploading more implementations there as we make progress, so that anyone facing similar problems has something to work from. The next few weeks will be quite busy with other commitments and some travel, but I will update the repo whenever I get a free moment.
Nice! Making better circuit simulations in MTK will help a lot of applications. My work on solar PV will also benefit from this, since it is essentially another diode-based circuit.
If you want to simulate switching-mode power converters, you should use a linearized state-space representation as PLECS does. See this very old paper: PLECS-piece-wise linear electrical circuit simulation for Simulink | IEEE Conference Publication | IEEE Xplore
Thanks for the reference! I will have a look at it. Once I have new results I will let you know.