how to make a loop with an if condition for a system of equations

So I have to make a for loop with an if condition for stoping the loop and I do not know where to put (i) and also in the if statement I do not know if I wrote the code right, any help in this problem would be appreciated.

z1 =[0;0.0335;0.1250;0.2500;0.3750;0.4665;0.5000];
z2 =[-0.049331610813001978382442370903885;-0.049073798515771506144130540704333;...
    -0.045591949992497543227192386552303;-0.034894245861879646501743852954066;...
    -0.018884472445635220242643376899604;-0.0051843520464032920772057925149539;...
    -0.000020337158849323484471603084704967];

z3 =[-0.00000017038751456961245595141592474853;-0.0042702578341115820738294669113829;...
    -0.0086318152310887342247093465740075;-0.013727659494014945417238515704958;...
    -0.01713291306879930987613700012663;-0.018241529138321719994971901104077;...
    -0.020949265945605807885946000891226];

z4 =[-0.049331610813001978382442370903885;-0.049073798515771506144130540704333;...
    -0.045591949992497543227192386552303;-0.034894245861879646501743852954066;...
    -0.018884472445635220242643376899604;-0.0051843520464032920772057925149539;...
    -0.000020337158849323484471603084704967];

w0onX1 =[0.005, 0.0049723456191038223899642948557202, 0.0046193976625564337806409159469839, 0.0035355339059327376220042218105242,...
    0.001913417161825448858642299920152, 0.00052514668824912603824992069930841, 0];

A11 =765403.71102841025915215223734811;
Nu =0.2850;
D11 =0.000000090096847435747048489786220240126;
A55 =225053.1567503304238293379021462;
N =7;

C{1,1} =[-48.6666666666667,59.7128129211020,-16.0000000000000,8.00000000000000,-5.33333333333333,4.28718707889796,-2.00000000000000;
    -14.9282032302755,6.92820323027552,10.9282032302755,-4.61880215351701,2.92820323027551,-2.30940107675850,1.07179676972449;
    4.00000000000000,-10.9282032302755,1.33333333333333,8.00000000000000,-4.00000000000000,2.92820323027551,-1.33333333333333;
    -2.00000000000000,4.61880215351701,-8.00000000000000,2.22044604925031e-15,8.00000000000000,-4.61880215351700,2.00000000000000;
    1.33333333333334,-2.92820323027551,4.00000000000001,-8.00000000000001,-1.33333333333331,10.9282032302755,-4.00000000000000;
    -1.07179676972449,2.30940107675850,-2.92820323027551,4.61880215351701,-10.9282032302755,-6.92820323027556,14.9282032302755;
    2.00000000000000,-4.28718707889796,5.33333333333333,-8.00000000000000,16.0000000000000,-59.7128129211020,48.6666666666667];

C{1,2} =[1381.33333333333,-2246.42709737007,1301.33333333333,-714.666666666666,490.666666666666,-398.906235963265,186.666666666667;
    684.553755053225,-1002.66666666667,390.276877526613,-106.666666666667,57.7231224733878,-42.6666666666667,19.4462449467756;
    -53.3333333333334,209.709376403674,-277.333333333333,149.333333333333,-42.6666666666666,24.9572902629932,-10.6666666666667;
    16.0000000000000,-42.6666666666667,128.000000000000,-202.666666666667,128.000000000000,-42.6666666666666,16.0000000000000;
    -10.6666666666666,24.9572902629931,-42.6666666666666,149.333333333333,-277.333333333333,209.709376403673,-53.3333333333335;
    19.4462449467757,-42.6666666666669,57.7231224733881,-106.666666666667,390.276877526613,-1002.66666666667,684.553755053226;
    186.666666666667,-398.906235963266,490.666666666667,-714.666666666667,1301.33333333333,-2246.42709737007,1381.33333333333];
C_1 = cell2mat(C(1));

syms U [N,1]
syms PHI [N,1]
syms W [N,1]

coef1=(((2*A11*z1.^2).*(C{1,2}))+(((2*A11*z1).*(1+z4)).*(C{1,1}))+(2*A11*((Nu*z4)-1)));
coef2=(((A11*z1).*(((1-Nu)*z3)+((2*z2).*z1)+(z3.*z4))).*(C{1,1}));
coef3=(((D11*z1.^2).*(C{1,2}))+((D11*z1).*(C{1,1}))-(+D11+(A55*z1.^2)));
coef4=(-((A55*z1.^2).*(C{1,1})));
coef5=(((2*A55*z1).*(C{1,2}))+(((2*A55)+((A11*z3).*(((3*z2).*z1)+z3))).*(C{1,1})));
coef6=(((2*A55*z1).*(C{1,1}))+(2*A55));
coef7=((((2*A11*z3).*z1).*(C{1,2}))+(((2*A11)*((z1.*z2)+((1+Nu)*z3))).*(C{1,1}))+...
    (2*Nu*A11*z2));

eq1 = ((coef1*U)+(coef2*W)) == 0;
eq2 = ((coef3*PHI)+(coef4*W)) == 0;
eq3 = ((coef5*W)+(coef6*PHI)+(coef7*U)) == 0;
vars = [U;PHI;W];
[Omega,b] = equationsToMatrix([eq1,eq2,eq3],vars);

Omega(1,1) = 1; Omega(N+1,N+1) = 1; Omega(3*N,3*N) = 1;
Omega((2*N)+1,((2*N)+1:3*N)) = C_1(1,(1:N));
Omega(2*N,(N+1:(2*N)-1)) = C_1(N,(1:N-1)); Omega(2*N,2*N) = C_1(N,N)+Nu;
Omega(1,(2:3*N)) = 0; Omega(N+1,(N+2:3*N)) = 0; Omega(2*N,((2*N)+1:3*N)) = 0;
Omega((2*N)+1,(1:2*N)) = 0; Omega(3*N,(1:(3*N)-1)) = 0;

b(N) = (2*z1(N)*(-100000));
sol = Omega\b;
solU = vpa(sol(1:N));
solW = vpa(sol(1+(2*N):N+(2*N)));
solPHI = vpa(sol(1+N:N+N));
Wtotal = w0onX1'+solW;
% some kind of loop like this:
% for i=1:50
% solving equations for up to 50 times and:
% in every step this must happen:
% z2 = C{1,2}*solW & z3 = C{1,1}*solW
% if abs((((Wtotal(i+1)-Wtotal(i))/Wtotal(i))*100)<0.01)
% end, end
% so its like for the first step everything is ready but for second step
% and so on solW that got calculated from the equations will be used in
% calculating the z2 and z3 and this has to happen until the if condition
% satisfies and stops the loop

Answer 1

[David Hill]

z2 =[-0.17162382699849727786214078518756;-0.091002558325694825093302376780363;...
    -0.031615696698776074425432150707406;-0.045376435830610389907683099983263;...
    -0.0049082191520506403009254218790972;-0.047113111927557759176851979732869;...
    -0.12231255346758490388052515365015];
z3 =[-0.00000017038751456961245595141592474853;-0.0042702578341115820738294669113829;...
    -0.0086318152310887342247093465740075;-0.013727659494014945417238515704958;...
    -0.01713291306879930987613700012663;-0.018241529138321719994971901104077;...
    -0.020949265945605807885946000891226];
C{1,1} =[-48.6666666666667,59.7128129211020,-16.0000000000000,8.00000000000000,-5.33333333333333,4.28718707889796,-2.00000000000000;
    -14.9282032302755,6.92820323027552,10.9282032302755,-4.61880215351701,2.92820323027551,-2.30940107675850,1.07179676972449;
    4.00000000000000,-10.9282032302755,1.33333333333333,8.00000000000000,-4.00000000000000,2.92820323027551,-1.33333333333333;
    -2.00000000000000,4.61880215351701,-8.00000000000000,2.22044604925031e-15,8.00000000000000,-4.61880215351700,2.00000000000000;
    1.33333333333334,-2.92820323027551,4.00000000000001,-8.00000000000001,-1.33333333333331,10.9282032302755,-4.00000000000000;
    -1.07179676972449,2.30940107675850,-2.92820323027551,4.61880215351701,-10.9282032302755,-6.92820323027556,14.9282032302755;
    2.00000000000000,-4.28718707889796,5.33333333333333,-8.00000000000000,16.0000000000000,-59.7128129211020,48.6666666666667];
C{1,2} =[1381.33333333333,-2246.42709737007,1301.33333333333,-714.666666666666,490.666666666666,-398.906235963265,186.666666666667;
    684.553755053225,-1002.66666666667,390.276877526613,-106.666666666667,57.7231224733878,-42.6666666666667,19.4462449467756;
    -53.3333333333334,209.709376403674,-277.333333333333,149.333333333333,-42.6666666666666,24.9572902629932,-10.6666666666667;
    16.0000000000000,-42.6666666666667,128.000000000000,-202.666666666667,128.000000000000,-42.6666666666666,16.0000000000000;
    -10.6666666666666,24.9572902629931,-42.6666666666666,149.333333333333,-277.333333333333,209.709376403673,-53.3333333333335;
    19.4462449467757,-42.6666666666669,57.7231224733881,-106.666666666667,390.276877526613,-1002.66666666667,684.553755053226;
    186.666666666667,-398.906235963266,490.666666666667,-714.666666666667,1301.33333333333,-2246.42709737007,1381.33333333333];
[Wtotal{1},solW]=calculate_W(z2,z3);
for i=2:50
  z2 = C{1,2}*solW;
  z3 = C{1,1}*solW;
  [Wtotal{i},solW]=calculate_W(z2,z3);
  ref=abs(((Wtotal{i}-Wtotal{i-1})./Wtotal{i-1})*100);
  if all((ref<0.01)|isnan(ref))
      break;
  end
end
solW %print solW

solW = 

function [Wtotal,solW] = calculate_W(z2,z3)
z1 =[0;0.0335;0.1250;0.2500;0.3750;0.4665;0.5000];
z4 =[-0.049331610813001978382442370903885;-0.049073798515771506144130540704333;...
    -0.045591949992497543227192386552303;-0.034894245861879646501743852954066;...
    -0.018884472445635220242643376899604;-0.0051843520464032920772057925149539;...
    -0.000020337158849323484471603084704967];
w0onX1 =[0.005, 0.0049723456191038223899642948557202, 0.0046193976625564337806409159469839, 0.0035355339059327376220042218105242,...
    0.001913417161825448858642299920152, 0.00052514668824912603824992069930841, 0];
A11 =765403.71102841025915215223734811;
Nu =0.2850;
D11 =0.000000090096847435747048489786220240126;
A55 =225053.1567503304238293379021462;
N =7;
C{1,1} =[-48.6666666666667,59.7128129211020,-16.0000000000000,8.00000000000000,-5.33333333333333,4.28718707889796,-2.00000000000000;
    -14.9282032302755,6.92820323027552,10.9282032302755,-4.61880215351701,2.92820323027551,-2.30940107675850,1.07179676972449;
    4.00000000000000,-10.9282032302755,1.33333333333333,8.00000000000000,-4.00000000000000,2.92820323027551,-1.33333333333333;
    -2.00000000000000,4.61880215351701,-8.00000000000000,2.22044604925031e-15,8.00000000000000,-4.61880215351700,2.00000000000000;
    1.33333333333334,-2.92820323027551,4.00000000000001,-8.00000000000001,-1.33333333333331,10.9282032302755,-4.00000000000000;
    -1.07179676972449,2.30940107675850,-2.92820323027551,4.61880215351701,-10.9282032302755,-6.92820323027556,14.9282032302755;
    2.00000000000000,-4.28718707889796,5.33333333333333,-8.00000000000000,16.0000000000000,-59.7128129211020,48.6666666666667];
C{1,2} =[1381.33333333333,-2246.42709737007,1301.33333333333,-714.666666666666,490.666666666666,-398.906235963265,186.666666666667;
    684.553755053225,-1002.66666666667,390.276877526613,-106.666666666667,57.7231224733878,-42.6666666666667,19.4462449467756;
    -53.3333333333334,209.709376403674,-277.333333333333,149.333333333333,-42.6666666666666,24.9572902629932,-10.6666666666667;
    16.0000000000000,-42.6666666666667,128.000000000000,-202.666666666667,128.000000000000,-42.6666666666666,16.0000000000000;
    -10.6666666666666,24.9572902629931,-42.6666666666666,149.333333333333,-277.333333333333,209.709376403673,-53.3333333333335;
    19.4462449467757,-42.6666666666669,57.7231224733881,-106.666666666667,390.276877526613,-1002.66666666667,684.553755053226;
    186.666666666667,-398.906235963266,490.666666666667,-714.666666666667,1301.33333333333,-2246.42709737007,1381.33333333333];
C_1 = cell2mat(C(1));
syms U [N,1]
syms PHI [N,1]
syms W [N,1]
coef1=(((2*A11*z1.^2).*(C{1,2}))+(((2*A11*z1).*(1+z4)).*(C{1,1}))+(2*A11*((Nu*z4)-1)));
coef2=(((A11*z1).*(((1-Nu)*z3)+((2*z2).*z1)+(z3.*z4))).*(C{1,1}));
coef3=(((D11*z1.^2).*(C{1,2}))+((D11*z1).*(C{1,1}))-(+D11+(A55*z1.^2)));
coef4=(-((A55*z1.^2).*(C{1,1})));
coef5=(((2*A55*z1).*(C{1,2}))+(((2*A55)+((A11*z3).*(((3*z2).*z1)+z3))).*(C{1,1})));
coef6=(((2*A55*z1).*(C{1,1}))+(2*A55));
coef7=((((2*A11*z3).*z1).*(C{1,2}))+(((2*A11)*((z1.*z2)+((1+Nu)*z3))).*(C{1,1}))+...
    (2*Nu*A11*z2));
eq1 = ((coef1*U)+(coef2*W)) == 0;
eq2 = ((coef3*PHI)+(coef4*W)) == 0;
eq3 = ((coef5*W)+(coef6*PHI)+(coef7*U)) == 0;
vars = [U;PHI;W];
[Omega,b] = equationsToMatrix([eq1,eq2,eq3],vars);
Omega(1,1) = 1; Omega(N+1,N+1) = 1; Omega(3*N,3*N) = 1;
Omega((2*N)+1,((2*N)+1:3*N)) = C_1(1,(1:N));
Omega(2*N,(N+1:(2*N)-1)) = C_1(N,(1:N-1)); Omega(2*N,2*N) = C_1(N,N)+Nu;
Omega(1,(2:3*N)) = 0; Omega(N+1,(N+2:3*N)) = 0; Omega(2*N,((2*N)+1:3*N)) = 0;
Omega((2*N)+1,(1:2*N)) = 0; Omega(3*N,(1:(3*N)-1)) = 0;
b(N) = (2*z1(N)*(-100000));
sol = Omega\b;
solU = vpa(sol(1:N));
solW = vpa(sol(1+(2*N):N+(2*N)));
solPHI = vpa(sol(1+N:N+N));
Wtotal = w0onX1'+solW;
end

[mahdi]

@David Hill

Thank you for answering, the code you wrote seems to be working and it stops in like iteration 5 because if statement gets satisfied but could you also please explain your method used?

[David Hill]

I placed most of your code into a function -- calculate_W. The only thing changing in the calculations are z2 and z3 which I provided as input. Outputs need to be at a minimum Wtotal and solW because they are used in future calculations of z2 and z3. Then the script (after providing the initial inputs of z2 and z3 and constant C) just performs the loop until meeting the condition. I had to modify your condition because the last Wtotal is zero causing division by zero to become nan.

[Wtotal{1},solW]=calculate_W(z2,z3);%perform initial calculation
for i=2:50
  z2 = C{1,2}*solW;%update z2 with new solW
  z3 = C{1,1}*solW;%update z3 with new solW
  [Wtotal{i},solW]=calculate_W(z2,z3);%reperform calculation
  ref=abs(((Wtotal{i}-Wtotal{i-1})./Wtotal{i-1})*100);
  if all((ref<0.01)|isnan(ref))%check conditional statement to break
      break;
  end
end

[mahdi]

@David Hill

I got it now, thank you.

Answer 2

[Pratheek]

Hi Mahdi,

I understand that you want to end the loop based on an if condition, you can leverage break statement to achieve this.

I'm attaching a sample code to do the same

for i = 1:50
    % Solving equations
    % ...

    % Calculate z2 and z3 using solW
    z2 = C{1,2} * solW;
    z3 = C{1,1} * solW;

    % Check the condition and stop the loop if it's satisfied
    if abs(((Wtotal(i+1) - Wtotal(i)) / Wtotal(i)) * 100) < 0.01
        break;
    end
end

1. The for loop iterates from i = 1 to i = 50.
2. Inside the loop, you can solve the equations and perform any necessary calculations.
3. After solving the equations, calculate z2 and z3 using solW.
4. Check the condition using the if statement. If the condition is satisfied, the loop will be stopped using the break statement.
5. If the condition is not met, the loop will continue to the next iteration until i = 50 or the condition is satisfied.

Hope this solves your issue!

[mahdi]

Thanks for answering but the if statement does not work and gives an error like this:

Index exceeds the number of array elements (7).

Error in sym/subsref (line 902)

R_tilde = builtin('subsref',L_tilde,Idx);

Error in forloophelp (line 76)

if abs(((Wtotal(i+1) - Wtotal(i)) / Wtotal(i)) * 100) < 0.01