While loop stops short of meeting the condition
2 visualizzazioni (ultimi 30 giorni)
Mostra commenti meno recenti
Milos Krsmanovic
il 12 Ott 2022
Commentato: Milos Krsmanovic
il 19 Ott 2022
I have a loop function which is supposed to find the x=0 intercept of two trigonometric functions x=x(t). The two functions alter depending on whether the x is greater or lower than 0. Although I'm using the same while loop for both x(t)>0 and x(t)<0 domains, the while loop unexpectedly stops short of reaching x(t)=0 for only x(t)<0 intervals of the function. That is, only one half of the same while loop works properly.
Please forgive me for not uploading the complete code - I have half a dozen .m files, and also there is an issue of intellectual property (I'm not the owner, just a collaborator). I will try to explain the best I can what the problem is.
The subfunction and the while loop are as follows:
function[x] = get_x(f_dn, X, phi, wndn, wddn, zetadn, ttemp)
if f_dn == 'u'
x = X*cos(wndn*ttemp-phi);
elseif f_dn == 'd'
x = X*exp(-zetadn*wndn*ttemp)*cos(wddn*ttemp-phi);
else fprintf('Incorrect f_up/f_dn entry (should be U or D).');
end
end
ttemp = 1.05*t; % t is the lower limit of the interval, where x(t)=0,...
% already provided as an argument
xapproach = 0.001;
if rem(m,2) == 0 % even interval
xup = xapproach;
xdn = -xapproach;
else
xup = xapproach;
xdn = -xapproach;
end
while (abs(x) > xapproach)
if rem(m,2) == 0 % even interval
if x > xup
ttemp = ttemp + abs(x)*ttemp/m;
elseif x < xdn
ttemp = ttemp - abs(x)*ttemp/m;
end
elseif rem(m,2) ~= 0 % odd interval
if x < xup
ttemp = ttemp + abs(x)*ttemp/m;
elseif x > xdn
ttemp = ttemp - abs(x)*ttemp/m;
end
else
disp('Some error.');
end
x = get_x(f_dn, X, phi, wndn, wddn, zetadn, ttemp)
end
Arguments such as X, phi, etc. are inherited from other functions, defined in other .m files.
The issue is visible from the plot:
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1153113/image.jpeg)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1153118/image.jpeg)
As you can tell, each odd (red) interval starting from the third one onwards ends prematurely, which in turn affects the following even (blue) interval's start. At the same time, even (blue) intervals complete just fine, and I am able to find the x(t)=0 upper limits of those intervals just fine. But both of these are the part of the same while loop. I'm really confused as to why is this happening.
Let us look at a typical command window print of x for the odd (red) interval:
x = -0.1280
x = -0.2757
x = -0.5566
x = -0.8112
x = -0.1367
x = 0.0759
x = -0.0463
x = 0.0270
x = -0.0162
x = 0.0096
x = -0.0057
x = 0.0034
x = -0.0020
x = 0.0012
So from what I can tell my function starts right, it sets a point on the odd (red) interval so that x(t)<0, then starts moving the point right along the red line. If it overshoots, it corrects and moves the point to the left by an increasingly smaller increment. However, it ultimately suddenly stops at x = 0.0012, where the condition abs(x) > xapproach is still ans = 1. Increasing the precision of xapproach produces the same result. As noted, the very same while loop finds the upper limit x(t)=0 of the blue interval just fine, as you can alo see from the plot.
Any ideas on what I can investigate in more detail to understand what the issue is? Quite a few similar questions here had an issue of long precision, but I do not see that being an issue here (again, the same while loop, etc.) Thank you in advance.
2 Commenti
Ghazwan
il 12 Ott 2022
Given the provided information, it looks like it is an issue of rounding the double numbers when it is approaching zero. A possible solution is that you equate the outputs when reaching these points. If not, increasing the number of calculating points would reduce the gap.
Risposta accettata
Mathieu NOE
il 12 Ott 2022
hello @Milos Krsmanovic
why not using this code ? as you see the zero (or any other threshold value) points are obtained with linear interpolation, so better accuracy as using nearest data point.
hope it helps
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1153403/image.png)
% dummy data
n = 100;
x= 25*(0:n-1)/n;
y1 = -cos(x).*exp(-0.1*x);
threshold = 0; % "zero" crossing or any other value here
[t0_pos1,s0_pos1,t0_neg1,s0_neg1]= crossing_V7(y1,x,threshold,'linear'); % positive (pos) and negative (neg) slope crossing points
% ind => time index (samples)
% t0 => corresponding time (x) values
% s0 => corresponding function (y) values , obviously they must be equal to "threshold"
figure(1)
plot(x,y1,'b.-',t0_pos1,s0_pos1,'*r',t0_neg1,s0_neg1,'*g','linewidth',2,'markersize',12);grid on
legend('signal','signal positive slope crossing points','signal negative slope crossing points');
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [t0_pos,s0_pos,t0_neg,s0_neg] = crossing_V7(S,t,level,imeth)
% [ind,t0,s0,t0close,s0close] = crossing_V6(S,t,level,imeth,slope_sign) % older format
% CROSSING find the crossings of a given level of a signal
% ind = CROSSING(S) returns an index vector ind, the signal
% S crosses zero at ind or at between ind and ind+1
% [ind,t0] = CROSSING(S,t) additionally returns a time
% vector t0 of the zero crossings of the signal S. The crossing
% times are linearly interpolated between the given times t
% [ind,t0] = CROSSING(S,t,level) returns the crossings of the
% given level instead of the zero crossings
% ind = CROSSING(S,[],level) as above but without time interpolation
% [ind,t0] = CROSSING(S,t,level,par) allows additional parameters
% par = {'none'|'linear'}.
% With interpolation turned off (par = 'none') this function always
% returns the value left of the zero (the data point thats nearest
% to the zero AND smaller than the zero crossing).
%
% check the number of input arguments
error(nargchk(1,4,nargin));
% check the time vector input for consistency
if nargin < 2 | isempty(t)
% if no time vector is given, use the index vector as time
t = 1:length(S);
elseif length(t) ~= length(S)
% if S and t are not of the same length, throw an error
error('t and S must be of identical length!');
end
% check the level input
if nargin < 3
% set standard value 0, if level is not given
level = 0;
end
% check interpolation method input
if nargin < 4
imeth = 'linear';
end
% make row vectors
t = t(:)';
S = S(:)';
% always search for zeros. So if we want the crossing of
% any other threshold value "level", we subtract it from
% the values and search for zeros.
S = S - level;
% first look for exact zeros
ind0 = find( S == 0 );
% then look for zero crossings between data points
S1 = S(1:end-1) .* S(2:end);
ind1 = find( S1 < 0 );
% bring exact zeros and "in-between" zeros together
ind = sort([ind0 ind1]);
% and pick the associated time values
t0 = t(ind);
s0 = S(ind);
if ~isempty(ind)
if strcmp(imeth,'linear')
% linear interpolation of crossing
for ii=1:length(t0)
%if abs(S(ind(ii))) >= eps(S(ind(ii))) % MATLAB V7 et +
if abs(S(ind(ii))) >= eps*abs(S(ind(ii))) % MATLAB V6 et - EPS * ABS(X)
% interpolate only when data point is not already zero
NUM = (t(ind(ii)+1) - t(ind(ii)));
DEN = (S(ind(ii)+1) - S(ind(ii)));
slope = NUM / DEN;
slope_sign(ii) = sign(slope);
t0(ii) = t0(ii) - S(ind(ii)) * slope;
s0(ii) = level;
end
end
end
% extract the positive slope crossing points
ind_pos = find(sign(slope_sign)>0);
t0_pos = t0(ind_pos);
s0_pos = s0(ind_pos);
% extract the negative slope crossing points
ind_neg = find(sign(slope_sign)<0);
t0_neg = t0(ind_neg);
s0_neg = s0(ind_neg);
else
% empty output
ind_pos = [];
t0_pos = [];
s0_pos = [];
ind_neg = [];
t0_neg = [];
s0_neg = [];
end
end
Più risposte (0)
Vedere anche
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!