Breaking out of a for-while loop

Managed to fix my mistake earlier. The code seems to be running in an infinite loop for the variable "pkpitch". May I ask, what should I do or take note of when doing for-while loops? Code attached.

clear all
close all
clc

Airfoildensity = 2000;
c = 0.5;
b = 0.25;
xf = 0.45*c;
t = 0.02;
l = 2;
e = xf - 0.25*c;
omegah = 2*2*pi;
omegaalpha = 6*2*pi;
Mass = 2000*t*c*l;
S = -(Mass)*(xf - c/2);
Ialpha = ((1/12)*(Mass)*(c)*(c)) + (Mass)*((xf - c/2))^2;
Density = 1.225;
Kh = Mass*omegah*omegah;
Kalpha = Ialpha*omegaalpha*omegaalpha;
r = xf - c/2;
s = (0.75*c) - xf;
a = xf - b;

%U = 1:1:80;
pkplunge = zeros(1,80);
pkpitch = zeros(1,80);
pkDRplunge = zeros(1,80);
pkDRpitch = zeros(1,80);
UpdatedPitchNatFreq = 0;
UpdatedPlungeNatFreq = 0;


j = 1;
for U = 1:1:80;
    Mmatrix = [Mass S; S Ialpha];
    Kmatrix = [Kh 0; 0 Kalpha];
    pkplunge(1,1) = 2*2*pi;
    pkpitch(1,1) = 6*2*pi;
    PlungeChangeRatio = abs(pkplunge(1,1))/(pkplunge(1,1));
    PitchChangeRatio = abs(pkpitch(1,1))/(pkpitch(1,1));
    UpdatingPlungeNatFreq = zeros(1,80);
    UpdatingPitchNatFreq = zeros(1,80);
    UpdatingPlungeDampRatio = zeros(1,80);
    UpdatingPitchDampRatio = zeros(1,80);
    
    ii = 1;
    while PitchChangeRatio > 0.001
        kpitch = (pkpitch(1,ii))*b/U;
        CtheoPitch = 1-(0.165/(1-(0.0455i/kpitch)))-(0.335/(1-(0.3i/kpitch)));
        DkmatrixPitch = Density*pi*b^2*[-1 a; a -(a^2+0.125*b^2)];
        Ematrix11Pitch = (-2*pi*Density*U*CtheoPitch);
        Ematrix12Pitch = (-Density*pi*b*b*U) + (2*pi*Density*b*U*CtheoPitch)*(a - 0.5*b);
        Ematrix21Pitch = (Density*pi*b*b*U) - (Density*pi*b*U)*(b - (2*a+b)*CtheoPitch);
        Ematrix22Pitch = (-Density*pi*b*b*U*(c/4)) + (Density*pi*b*b*U*(c/4)) + (Density*pi*b*U)*((b - (2*a+b))*CtheoPitch*(a-b/2));
        EmatrixPitch = [Ematrix11Pitch Ematrix12Pitch; Ematrix21Pitch Ematrix22Pitch];
        Fkmatrix12Pitch = -2*pi*Density*b*U*U*CtheoPitch;
        Fkmatrix22Pitch = (Density*pi*b*b*U*U) - (Density*pi*b*U*U)*(b - (2*a-b)*CtheoPitch);
        FkmatrixPitch = [0 Fkmatrix12Pitch; 0 Fkmatrix22Pitch];
        
        Amatrix11Pitch = zeros(2);
        Amatrix12Pitch = eye(2);
        Amatrix21Pitch = (inv(Mmatrix - DkmatrixPitch))*(FkmatrixPitch - Kmatrix);
        Amatrix22Pitch = (inv(Mmatrix - DkmatrixPitch))*EmatrixPitch;
        AmatrixPitch = [Amatrix11Pitch Amatrix12Pitch; Amatrix21Pitch Amatrix22Pitch];
        Amatrix_pitch = AmatrixPitch;
        [V,D] = eig(Amatrix_pitch);
        diagonals = diag(D);
        PitchEValues = sort(imag(diagonals), 'ascend');
        [~,index] = sort((PitchEValues),"ascend")
        UpdatedPitching = abs(PitchEValues(2));
        UpdatedPitchNatFreq = abs(PitchEValues(index(2)));
        
        PitchChangeRatio = abs(UpdatedPitchNatFreq - pkpitch(1,ii))/abs(pkpitch(1,ii));
        pkpitch(1,(ii+1)) = PitchEValues(2);
        
        ii = ii + 1;
    end
    UpdatingPitchNatFreq(j) = UpdatedPitchNatFreq/(2*pi);
    UpdatingPitchDampRatio(j) = -real(PitchEValues(index(2)))/UpdatedPitchNatFreq;
    ii = 1;
    while PlungeChangeRatio > 0.001
        kplunge = ((PlungeChangeRatio(1,ii)*b)/U); 
        Ctheoplunge = 1 - (0.165/(1-(0.0455/kplunge)*i)) - (0.335/(1-(0.3/kplunge)*i));
        Dkmatrixplunge = Density*pi*b^2*[-1 a; a -(a^2+0.125*b^2)];
        Ematrix11plunge = (-2*pi*Density*U*Ctheoplunge);
        Ematrix12plunge = (-Density*pi*b*b*U) + (2*pi*Density*b*U*Ctheoplunge)*(a - 0.5*b);
        Ematrix21plunge = (Density*pi*b*b*U) - (Density*pi*b*U)*(b - (2*a+b)*Ctheoplunge);
        Ematrix22plunge = (-Density*pi*b*b*U*(c/4)) + (Density*pi*b*b*U*(c/4)) + (Density*pi*b*U)*((b - (2*a+b))*Ctheoplunge*(a-b/2));
        Ematrixplunge = [Ematrix11plunge Ematrix12plunge; Ematrix21plunge Ematrix22plunge];
        Fkmatrix12plunge = -2*pi*Density*b*U*U*Ctheoplunge;
        Fkmatrix22plunge = (Density*pi*b*b*U*U) - (Density*pi*b*U*U)*(b - (2*a-b)*Ctheoplunge);
        Fkmatrixplunge = [0 Fkmatrix12plunge; 0 Fkmatrix22plunge];
        
        Amatrix11plunge = zeros(2);
        Amatrix12plunge = eye(2);
        Amatrix21plunge = (inv(Mmatrix - Dkmatrixplunge))*(Fkmatrixplunge - Kmatrix);
        Amatrix22plunge = (inv(Mmatrix - Dkmatrixplunge))*Ematrixplunge;
        Amatrixplunge = [Amatrix11plunge Amatrix12plunge; Amatrix21plunge Amatrix22plunge];
        Amatrix_plunge = Amatrixplunge;
        [V,D] = eig(Amatrix_plunge);
        diagonals = diag(D);
        PlungeEValues = sort(imag(diagonals), 'ascend');
        [~,index] = sort((PlungeEValues),"ascend")
        UpdatedPlunging = abs(PlungeEValues(1));
        UpdatedPlungeNatFreq = abs(PlungeEValues(index(1)));
        
        PlungeChangeRatio = abs(UpdatedPlungeNatFreq - pkplunge(1,ii))/abs(pkplunge(1,ii));
        pkplunge(1,(ii+1)) = PlungeEValues(1);
        ii = ii+1;
        
    end
    UpdatingPlungeNatFreq(j) = UpdatedPlungeNatFreq/(2*pi);
    UpdatingPlungeDampRatio(j) = -real(PitchEValues(index(1)))/UpdatedPlungeNatFreq;        
end

Answer 1

[Walter Roberson]

        [V,D] = eig(Amatrix_pitch);
        diagonals = diag(D);

The eigenvalues are complex, and have positive and negative imaginary components

        PitchEValues = sort(imag(diagonals), 'ascend');

Positive and negative imaginary components are extracted and sorted. The output PitchEValues is signed.

        [~,index] = sort((PitchEValues),"ascend")

You are re-sorting something that you just sorted, so the index order is going to be the 1 to the length

        UpdatedPitching = abs(PitchEValues(2));

The code explicitly acknowledges knowing that the PitchEValues are signed, by taking absolute value here, which is not something that it would do if the person writing the code had (mistakenly) believed that they were all positive.

        UpdatedPitchNatFreq = abs(PitchEValues(index(2)));

Again explicitly acknowledges that the values are signed. Also, since the sort that produced the index was applied to something already sorted, index(2) is always going to be 2, so UpdatedPitchNatFreq is always going to equal UpdatedPitching, and both are positive.

        PitchChangeRatio = abs(UpdatedPitchNatFreq - pkpitch(1,ii))/abs(pkpitch(1,ii));

Is pkpitch signed or unsigned? We will hold off on judging that for a second -- but when the code divides by the abs() then the code implicitly believes it might be negative. We know that UpdatedPitchNatFreq is positive. Hypothesize that pkpitch is negative: then positive minus negative gives a more positive number, and you can be sure that the output of the ratio would be greater than 1. If pkpitch is positive, then potentially PitchChangeRatio could decrease.

       pkpitch(1,(ii+1)) = PitchEValues(2);

We discussed above that PitchEValues is signed, so the stored pkpitch certainly has the potential to be signed.

In particular, after a few iterations, what happens in practice is that the calculation settles down and pkpitch does not change... but it is negative, and the updated PitchEValues is its positive counterpart, so you have a case of abs((positive X)-(negative X))/abs(negative X) which is going to be 2*(positive X) / (positive X) which is going to be 2. So the PitchChangeRatio locks in at 2, and your loop never ends.