I'm using the realFFT library with OpenModelica to analyze the frequencies in a PWM signal. When I analyze frequencies below 10 kHz everything works fine. But as soon as I set my maximum frequency to more than 10 kHz my simulation either calculates wrong results, crashes or says that it simulated but doesn't show the results.
What I found out so far: I have a the samplePeriod=8.3us and the number of samples ns=6000 for a max frequency of f_max=12kHz and a resolution of f_res=20Hz which results in a stopTime>=(6000-1)*8.3us --> stopTime>=0.05s (according to the library). So I simulated with this stopTime and it worked but when I set stopTime=0.1s it doesn't show any results but says that it simulated fine. When I further increased the stopTime to 0.2s it calculated wrong results.
This doesn't make sense to me, why should it fail with an increased stopTime? Might this be another OpenModelica issue and works fine with Dymola?
Here are my models:
model TestRealFFT
Modelica.Blocks.Sources.Pulse pulse1(amplitude = 1, offset = 0, period = 100e-6, width = 40) annotation(
Placement(visible = true, transformation(origin = {-168, 16}, extent = {{-60, -60}, {60, 60}}, rotation = 0)));
FFTmath fFTmath1(f_max = 12000, f_res = 20) annotation(
Placement(visible = true, transformation(origin = {176, 8}, extent = {{-80, -80}, {80, 80}}, rotation = 0)));
equation
connect(pulse1.y, fFTmath1.u) annotation(
Line(points = {{-102, 16}, {68, 16}, {68, 8}, {80, 8}}, color = {0, 0, 127}));
annotation(
stopTime = 2,
Diagram(coordinateSystem(extent = {{-300, -250}, {300, 250}})),
Icon(coordinateSystem(extent = {{-300, -250}, {300, 250}})),
__OpenModelica_commandLineOptions = "",
experiment(StartTime = 0, StopTime = 0.2, Tolerance = 1e-06, Interval = 8e-06));
end TestRealFFT;
and
model FFTmath
import Modelica.Constants.{pi};
import Modelica.Math.FastFourierTransform.*;
import Modelica.SIunits.*;
Modelica.Blocks.Interfaces.RealInput u annotation(
Placement(visible = true, transformation(origin = {-120, 0}, extent = {{-20, -20}, {20, 20}}, rotation = 0), iconTransformation(origin = {-120, 0}, extent = {{-20, -20}, {20, 20}}, rotation = 0)));
parameter Frequency f_max = 12000 "Maximum frequency of interest";
parameter Frequency f_res = 20 "Frequency resolution";
final parameter Integer ns = realFFTsamplePoints(f_max, f_res, f_max_factor = 5) "Number of samples";
final parameter Integer nf = div(ns, 2) + 1 "Number of frequency points";
final parameter Integer nfi = max(1, min(integer(ceil(f_max / f_res)) + 1, nf));
final parameter Frequency f_i[nfi](each fixed = false) "FFT frequencies of interested frequency points";
parameter Time samplePeriod = 1 / (2 * f_res * div(ns, 2));
output Integer info(start = 0, fixed = true) "Information flag from FFT computation";
Integer iTick(start = 0, fixed = true);
discrete Real Buf[ns](start = zeros(ns), each fixed = true) "Input buffer";
Real A_i[nfi](start = zeros(nfi), each fixed = true) "FFT amplitudes";
Real Phi_i[nfi](start = zeros(nfi), each fixed = true) "FFT phases";
Real y(start = 0, fixed = true);
// "Signal from which FFT is computed";
initial equation
for i in 1:nfi loop
f_i[i] = (i - 1) * f_res;
end for;
equation
algorithm
when sample(0, samplePeriod) then
iTick := iTick + 1;
y := u;
if iTick <= ns then
Buf[iTick] := y;
end if;
if iTick == ns then
(info, A_i, Phi_i) := realFFT(Buf, nfi);
end if;
end when;
//3 * sin(2 * pi * f1 * time) + sin(2 * pi * f2 * time);
annotation(
Documentation(Icon(graphics = {Text(origin = {-42, 62}, extent = {{110, -78}, {-30, 18}}, textString = "FFT"), Rectangle(origin = {0, -79}, fillPattern = FillPattern.Solid, extent = {{-80, -1}, {80, 1}}), Rectangle(origin = {-79, -49}, fillPattern = FillPattern.Solid, extent = {{-1, -29}, {1, 29}}), Polygon(origin = {-79, -15}, fillPattern = FillPattern.Solid, points = {{0, -5}, {-6, -5}, {0, 5}, {6, -5}, {6, -5}, {0, -5}}), Polygon(origin = {85, -79}, rotation = -90, fillPattern = FillPattern.Solid, points = {{0, -5}, {-6, -5}, {0, 5}, {6, -5}, {6, -5}, {0, -5}}), Rectangle(origin = {-59, -65}, fillPattern = FillPattern.Solid, extent = {{-1, 23}, {1, -15}}), Ellipse(origin = {-59, -39}, fillPattern = FillPattern.Solid, extent = {{-3, 3}, {3, -3}}, endAngle = 360), Ellipse(origin = {-49, -61}, fillPattern = FillPattern.Solid, extent = {{-3, 3}, {3, -3}}, endAngle = 360), Ellipse(origin = {-19, -53}, fillPattern = FillPattern.Solid, extent = {{-3, 3}, {3, -3}}, endAngle = 360), Ellipse(origin = {25, -67}, fillPattern = FillPattern.Solid, extent = {{-3, 3}, {3, -3}}, endAngle = 360), Ellipse(origin = {31, -49}, fillPattern = FillPattern.Solid, extent = {{-3, 3}, {3, -3}}, endAngle = 360), Rectangle(origin = {-49, -65}, fillPattern = FillPattern.Solid, extent = {{-1, 1}, {1, -15}}), Rectangle(origin = {31, -65}, fillPattern = FillPattern.Solid, extent = {{-1, 15}, {1, -15}}), Rectangle(origin = {-19, -63}, fillPattern = FillPattern.Solid, extent = {{-1, 7}, {1, -15}}), Rectangle(origin = {25, -63}, fillPattern = FillPattern.Solid, extent = {{-1, -5}, {1, -15}}), Ellipse(origin = {-67, -61}, fillPattern = FillPattern.Solid, extent = {{-3, 3}, {3, -3}}, endAngle = 360), Rectangle(origin = {-67, -63}, fillPattern = FillPattern.Solid, extent = {{-1, 1}, {1, -15}}), Rectangle(origin = {37, -63}, fillPattern = FillPattern.Solid, extent = {{-1, -5}, {1, -15}}), Ellipse(origin = {37, -67}, fillPattern = FillPattern.Solid, extent = {{-3, 3}, {3, -3}}, endAngle = 360), Line(points = {{-100, 100}, {100, 100}, {100, -100}, {-100, -100}, {-100, 100}, {-100, 100}}, thickness = 0.5)}),
Diagram,
__OpenModelica_commandLineOptions = "");
end FFTmath;
I'm very thankful for any help!
Looks like a bug in OpenModelica. Please open a ticket about it on trac.openmodelica.org.
Problem 1 Huge result files are breaking
When simulating with OpenModelica (v1.14.0-dev-26633-gd9901afc5b) for 0.05 sec simulation time the result mat-file is already 1.06 GB. For 0.1 sec you get around 2.1 GB and a corrupt header.
Problem 2 Wrong results
When you simulate interval [0, 0.2] the solution is completely wrong. Values of y are somewhere around 1.4e+31 and time goes up to 5e+218. Could be because of the broken result file.
But even if you reduce tolerance and interval it won't simulate correctly, maybe because of the huge amount of events.
Works in Dymola 2019, but needs a lot of time to open the result files.