How to change initial conditions and run loop again until R^2 is above 0.8

I am trying to create a combined pendulum model to model the force from a CFD solution and instead of guessing and checking I want to create a script to change conditions such as lengths (L and L2), initial angle(A1 and A2) and damping coefficient (c and c2) and re run the (2 pendulum forces) Fr and Fr2 in my script and compare the FR (combined force pendulum) to the Total Force Variable.

I currently an running a while loop with a nested if statement at the end to chnage the varaible but its not doing anything and the constants remain zero any help would be greatly appreciated.

Thanks, Willie

%%
CFDoutput = readmatrix('Cyl_3T_4_15_Fill_90sec.csv');
% Cut first 5 seconds from CFD
T   = CFDoutput(:,1);
loc = find(T<5);
start = length(loc);
loc2 = find(T==80.0130);
endend = length(loc2);

time = T(start:end);
XForce = CFDoutput(:,8);
xForce = XForce(start:end);
YForce = CFDoutput(:,9);
yForce = YForce(start:end);
ZForce = CFDoutput(:,10);
zForce = ZForce(start:end);

TotalForce = sqrt(xForce.^2 + yForce.^2 + zForce.^2);

%% 2D pendulum motion 
R2 = 0;

 while R2 < 0.8
% System 1 Variables 
mass=67500; % [kg] pendulum bob mass
L=5.9;      % [m] pendulum rod length
g=9.81;     % [m/s^2] acceleration due to gravity
t=90;       % [s] simulation time
c =0;% 0.05;  % Damping coefficient
dt= 2/91;   % Time step
i= t/dt;    % Number of iterations
A1 = 0;     %initial angle deg
% System 2 Variables 
mass2=67500; % [kg] pendulum bob mass
L2=2;      % [m] pendulum rod length
c2 =0;% 0.1;    % Damping coefficient
A2 = 0;     % initial angle deg

% Allocation of Variables  
theta=zeros(round(i),1);
omega=zeros(round(i),1);
alpha=zeros(round(i),1);
Fr=zeros(round(i),1);
theta(1,:)=deg2rad(A1); % initial angular position
omega(1,:)=0;           % initial angular velocity
alpha(1,:)=0;           % initial angular acceleration

% Solve Equations of Motion 
for n=1:i
    theta(n+1,:)= theta(n,:)+omega(n,:)*dt; 
    omega(n+1,:)= omega(n,:)+alpha(n,:)*dt; 
    alpha(n+1,:)=(-g*sin(theta(n+1,:)))/L-c*omega(n+1,:); 
    Fr(n+1,:)= mass*g*cos(theta(n+1,:))+ mass*L*(omega(n+1,:)).^2; 
end

% % 2D pendulum motion 2

% Allocation of variables 2
t=0:i;    
theta2=zeros(round(i),1);
omega2=zeros(round(i),1);
alpha2=zeros(round(i),1);
Fr2=zeros(round(i),1);
theta2(1,:)=deg2rad(A2); % initial angular position
omega2(1,:)=0;           % initial angular velocity
alpha2(1,:)=0;           % initial angular acceleration

% Solve Equations of Motion 
for n=1:i
    theta2(n+1,:)= theta2(n,:)+omega2(n,:)*dt; 
    omega2(n+1,:)= omega2(n,:)+alpha2(n,:)*dt; 
    alpha2(n+1,:)=(-g*sin(theta2(n+1,:)))/L2-c2*omega2(n+1,:); 
    Fr2(n+1,:)= mass2*g*cos(theta2(n+1,:))+ mass2*L2*(omega2(n+1,:)).^2; 
end

% Combined Pendulum Force
FR = (Fr + Fr2)/2;

%R squared Value curve fitting 
Fit = fitlm(TotalForce,FR);
R2 = Fit.Rsquared.Ordinary;
    if R2 < 0.8
        for j=1:i
        c = c + 0.001;
        c2 = c2 + 0.001;
        L = L + 0.001;
        L2 = L2 + 0.001;
        break
    end
end

[Walter Roberson]

Put everything except the initial conditions into a function. Have the function return R2. Now loop around providing new initial conditions and calling the function until R2 meets your criteria.

[Wilfredo Huaman]

I understand almost everything except the part where you say "Now loop around providing new initial conditions" could you elaborate a little more on this?