Recurring Error in Simulation

So, I have a fairly complicated Simulink Model which run again and again, to generate different outputs. This has worked fine for several weeks, and now suddenly I have this Error. I didn't change anything in the model, only a random variable definition in the matlab script which initiates the simulation. The variable is not used in the model.

I was able to get the Model to run once by closing the window, reopening it and reloading some variables, but upon a second unaltered simulation it crashed again. The Error message I get is the following:

Error:An error occurred while running the simulation and the simulation was terminated
Caused by:
['My Model/Solver Configuration3']: Transient initialization at time 3.116414621485137e-05, solving for consistent states and modes, failed to converge.
Nonlinear solver: failed to converge, residual norm too large.
Here is the set of components with unconverged equations:

'My Model/Transmission Line3' % an ideal transmission line
Equation location is:
'C:\Program Files\MATLAB\R2023b\toolbox\physmod\elec\library\m\+ee\+passive\+lines\ideal_transmission_line.sscp'(no line number info)

'My Model/Simscape Component1'  % this is a component I made myself by copying the Capacitor component and allowing negative capacitances
Equation location is:
'\\My Path\+NegativeCapacitor\NegativeCapacitor.ssc' (line 29)
Equation location is:
'\\My Path\+NegativeCapacitor\NegativeCapacitor.ssc' (line 29)

'My Model/Simscape Component' % same negative Capacitor
Equation location is:
'\\My Path\+NegativeCapacitor\NegativeCapacitor.ssc' (line 29)
Equation location is:
'\\My Path\+NegativeCapacitor\NegativeCapacitor.ssc' (line 29)

'My Model/Resistor3' % random resistor 1 Ohm
Equation location is:
'C:\Program Files\MATLAB\R2023b\toolbox\physmod\simscape\library\m\+foundation\+electrical\branch.ssc' (line 23)

'My Model/Resistor5'    % random resistor ca 120 Ohm
Equation location is:
'C:\Program Files\MATLAB\R2023b\toolbox\physmod\simscape\library\m\+foundation\+electrical\branch.ssc' (line 23)

I've encountered this error several times before, the last time several months ago, but could so far avoid it by reloading previous versions of the file, or closing and opening again. My usual strategies aren't working this time and I'd like to know why this error occurs so I can fix it once and for all.

A similar question suggested checking for algebraic loops, but the Code

Simulink.BlockDiagram.getAlgebraicLoops('My Model')

did not find any algebraic loops.

I'm using Automatic Solver Selection with Variable Step. My model is pretty idealized, so abrupt jumps occur, however this has not been an issue in the past and I don't understand why it would be now, if this is why the error occurs (this has been suggested in another answer to a similar question). I have found that after a few microseconds the current in my circuit jumps to about 10^300, which is where the simulation stops, understandably. I just don't understand why the current would suddenly jump so high and why at that moment.

If anyone has any ideas on what is causing this error and how to fix it, I would greatly appreciate the help.

Answer 1

[Umar]

Hi Annika,

To troubleshoot, check the component equations for inconsistencies, ensure proper initialization values, and verify the model's stability. Review the custom Simscape components for any unconventional behavior causing abrupt changes. Consider refining the solver settings or adjusting component parameters to prevent sudden spikes in values. Debugging the model step by step and monitoring key variables can help pinpoint the root cause of the convergence problem and facilitate a solution.

Hope this will help resolve your issue.

[Umar]

To debug this issue, we can start by checking the calculations of the material parameters. We can print out the values of these parameters for each selected material to ensure they are correct. Additionally, we can check if any of these parameters are being used in further calculations that could lead to the sudden jump in current.

% InitFunc and all Parameter calculations for % Galliere_neg_Kond_V27_240105_opt_WIP

%% Geometry definition %piezo: material = 'P52_datasheet'; r = 0.003; d = 0.005; %transducer transducer = 'F8V8'; %'TEW20', 'F7G3_TEW21', 'F10G4', 'F8V8' FL_mat = 'water';

global A Am; A = r^2*pi; Am = r*pi*d;

switch transducer case 'TEW20' d_back = 0.003; d_front = 0.0028; BL_mat = 'air'; front_mat = 'epoxy'; back_mat = 'alu'; back_mat2 = 'air'; d_back2 = 0; case 'F7G3_TEW21' d_back = 0.003; d_front = 0.0028; BL_mat = 'air'; front_mat = 'epoxy'; back_mat = 'steel'; back_mat2 = 'air'; d_back2 = 0; case 'F10G4' d_back = 0.002; d_front = 0.003; front_mat = 'epoxy'; BL_mat = 'air'; back_mat = 'steel'; back_mat2 = 'air'; d_back2 = 0; case 'F8V8' d_back = 0.003; d_front = 0.0005; d_back2 = 0.001; back_mat2 = 'air'; BL_mat = 'air'; front_mat = 'alu'; back_mat = 'epoxy'; back_mat2 = 'alu'; end

define_params('BL', BL_mat, false); global Z_BL c_BL rho_BL ZR; define_params('back', back_mat, false); global Z_back c_back rho_back; define_params('front', front_mat, true); global Z_front c_front rho_front; define_params('FL', FL_mat, false); global Z_FL c_FL rho_FL;

global Z_back2 c_back2 rho_back2; define_params('back2', back_mat2, false);

%% Material Parameters switch material case 'P52' g33 = 0.0142188; g31 = -0.0333414; c33D = 1.362e11; c11D = 7.268e10; h33 = 0.0142188*1.362e11; h31 = -0.0333414*7.268e10; eps33S = 9.249e-9; rho = 14672.908; case 'P5' c33D = 1.295e11; c11D = 8.587e10; h33 = 1.302e9; h31 = -1.684e9; eps33S = 8.711e-9; rho = 17435.538; case 'P52_datasheet' g33 = 23e-3; g31 = -10e-30; c11D = 7.1428e10; c33D = 1e11; h33 = g33*c33D; h31 = g31*c11D; eps33S = 11.15627e-9; rho = 7500; end

k = 1.5*d/r; C0 = A*eps33S/d; Z33 = A*sqrt(c33D/rho)*rho; t33 = d/sqrt(c33D/rho); Z11 = 2*pi*d*r*sqrt(c11D/rho)*rho; t11 = 2*r/sqrt(c11D/rho);

t_back = d_back/c_back; t_back2 = d_back2/c_back2; t_front = d_front/c_front;

function define_params(name, mat, mantel) global A Am geomat = mat; switch geomat case 'alu' geo_c = 5037; geo_rho = 2700; geo_Z = geo_c * geo_rho *A; geo_ZR = geo_c * geo_rho *Am; case 'steel' geo_c = 5752; geo_rho = 7800; geo_Z = geo_c * geo_rho *A; geo_ZR = geo_c * geo_rho *Am; case 'epoxy' geo_c = 2450; geo_rho = 1200; geo_Z = geo_c * geo_rho *A; geo_ZR = geo_c * geo_rho *Am; case 'air' geo_c = 330; geo_rho = 1.2; geo_Z = geo_c * geo_rho *A; geo_ZR = geo_c * geo_rho *Am; case 'water' geo_c = 1500; geo_rho = 1000; geo_Z = geo_c * geo_rho *A; geo_ZR = geo_c * geo_rho *Am; end eval(append('global Z_', name, ' c_', name, ' rho_', name, ' ZR ;')); eval(append('Z_', name, ' = ', sprintf('%.10f',geo_Z), ';')); eval(append('c_', name, ' = ', sprintf('%.10f',geo_c), ';')); eval(append('rho_', name, ' = ', sprintf('%.10f',geo_rho), ';')); if mantel == true eval(append('ZR = ', sprintf('%.10f',geo_ZR), ';')); end end

% Print material parameters fprintf('Material Parameters:\n'); fprintf('g33: %.10f\n', g33); fprintf('g31: %.10f\n', g31); fprintf('c33D: %.10f\n', c33D); fprintf('c11D: %.10f\n', c11D); fprintf('h33: %.10f\n', h33); fprintf('h31: %.10f\n', h31); fprintf('eps33S: %.10f\n', eps33S); fprintf('rho: %.10f\n', rho);

[Umar]

Also, another issue lies in the line where the capacitance parameter is defined as C = { 1, 'F' };. The value of the capacitance is set to 1 Farad, which might lead to numerical instability during the simulation, especially when the frequency of operation is high. To address the bug and improve the stability of the simulation, it is recommended to adjust the initial capacitance value to a more reasonable value. You can modify the capacitance parameter initialization as follows:

C = { 1e-6, 'F' }; % Set the initial capacitance to 1 microfarad

By setting the initial capacitance to a smaller value (1 microfarad in this case), the simulation is less likely to encounter numerical issues, especially at higher frequencies.

Step by Step guidance

Step by Step:

Save the text from NegativeCapacitor.txt as a .ssc file.

Create a folder named +NegativeCapacitor and place the .ssc file and the corresponding .jpg file in it.

Add the folder to the Matlab path.

Run InitFunc.m to initialize the necessary functions and variables.

Define the frequency variable by executing freqHz = 200000; in the Matlab console.

Experiment with different frequency values between 1e5 and 1e6 to observe the system behavior.

Run the simulation to analyze the performance of the NegativeCapacitor component at the specified frequency.

By following these steps and adjusting the initial capacitance value, you can enhance the stability and accuracy of the simulation process in Matlab.

[Annika Michel]

First of all, thanks for your answer! I'm going a little crazy with this and I'm very glad your trying to help.

You wrote above that I could try adjusting the solver settings. What settings would you recommend? I have tried lowering the max step time, however this didn't fix the issue and the simulation got of course much slower. I don't know the Solvers well enough to know how they differ, so if you can recommend one that is better than the auto selection, I would be glad to know.

As for the variables in the block parameters: Those all check out. I've checked them, and again, this exact simulation ran like this for months with these exact parameters. I don't know where you see a capacitance defined as 1F, in my model all capacitances are of the order 1e-8F or lower. My custom component is of course a bit of a wildcard, but since I literally changed only "C>0" to "abs(C>0)" in the capacitor for this I don't think it is the source of the issue. Again, I have been using this simulation in this form for months and never had an issue.

I may try upgrading to 2024a and see if that helps, if you have any specific advice beyond "check everything again" I'd be glad to hear it.

Thanks again for your answer!

Answer 2

[Umar]

Hi Annika,

No problem, so I will try to tackle your questions by providing solutions to them to best of my knowledge. So, your first question was

What settings would you recommend? I have tried lowering the max step time, however this didn't fix the issue and the simulation got of course much slower. I don't know the Solvers well enough to know how they differ, so if you can recommend one that is better than the auto selection, I would be glad to know.

In my opinion, default solver selection is usually reliable, but for specific cases, manual selection might be beneficial. Here are some common solvers and their characteristics:

ode45: A good general-purpose solver for non-stiff problems. ode15s: Ideal for stiff systems where fast transients are present. ode23: Balances efficiency and accuracy for medium-sized problems. You can manually select a solver using the odeset function. For instance:

options = odeset('RelTol',1e-5,'AbsTol',1e-8); [t,y] = ode45(@odefun,tspan,y0,options);

As for the variables in the block parameters: Those all check out. I've checked them, and again, this exact simulation ran like this for months with these exact parameters. I don't know where you see a capacitance defined as 1F, in my model all capacitances are of the order 1e-8F or lower. My custom component is of course a bit of a wildcard, but since I literally changed only "C>0" to "abs(C>0)" in the capacitor for this I don't think it is the source of the issue. Again, I have been using this simulation in this form for months and never had an issue.

If you suspect a custom component might be causing the issue, it's essential to scrutinize its implementation. Verify the changes made, especially the alteration from "C>0" to "abs(C>0)" in the capacitor. Ensure that this modification aligns with the expected behavior and doesn't introduce unintended consequences.

I may try upgrading to 2024a and see if that helps, if you have any specific advice beyond "check everything again" I'd be glad to hear it.

That is a valid option because newer versions often come with bug fixes and performance enhancements. Before upgrading, ensure compatibility with your existing codebase and any specific features you rely on.

While you have validated the parameters and settings, a systematic troubleshooting approach can help pinpoint the root cause. Check for any recent changes, unexpected interactions, or boundary conditions that might have been overlooked.

Let me know if I can assist further.

[Annika Michel]

Thanks for your help! I saved all the workspace variables in a working iteration and then with the error and compared them one by one. And I found my mistake. I was accidentally using a variable name twice, so I overwrote a gain value with a much higher number. This explains why I suddenly had an incredibly high current, and also why it worked when I ran the model outside of my iterating code.

So yeah, my mistake and thanks for your help :D I am an idiot.

Answer 3

[Umar]

Hi Annika,

Glad to help out. People make mistakes but when we work as a team together, problems get resolved. I am glad that you are happy and problem is resolved. Good luck with your future endeavors.