What am i doing wrong? This expression is violating basic Mathematics.

If you run this code, you will get

SI(1,1) = -1 [SI is the sine term of the main formula mentioned below]

sqroot(1,1) = 0.0026

MI = 0.0026

Hence, according to %Main Formula

pos1(1,1) = 0.0026 + 0.0026*(-1) = 0

But the answer matlab is returning me is 0.002

How is this even possible?

Is there any error in the for loop?

%Defining values of parameters
Fg = 2.1188;
Fi = 2.7834;
Ft = Fg + Fi;
h = 0.002;
k = 1864;
e = 0.7;
m = 0.001;

% Initializing for x
x = zeros(11,1);
x(1,1) = 0;
l = 0.002;

%Formula 1 - determining new positions
for i = 1:10
    con1 = ((Ft*h)-(k*(h^2)*0.5))*(1-(e^2));
    x(i+1,1) = (Ft-(sqrt((Ft^2)-(2*k*((con1)+((e^2)*(Ft*x(i,1)-(k*0.5*(x(i,1)^2)))))))))/k;
end

%Initializing for constants N & S
N = zeros(10,1);
S = zeros(10,1);
W = (2*Ft)/m;

%Calculating constants N & S
for i = 1:10
    S(i) = (((W*l)-((k*(l^2))/m))*((e^2-1))-((e^2)*(x(i,1))*(W+((k*(x(i,1)))/m))));
end

sininv = zeros(10,1);
sqroot = zeros(10,1);
pos1 = zeros(5000,1);
pos2 = zeros(5000,1);
time = zeros(5000,1);
M = 2*Ft; 
time(1,1) = 0;

for i = 1:10
     N(i,1) = ((k*((x(i,1))^2)) - (2*Ft*(x(i,1))));
     sininv(i,1) = asin((x(i,1)-(M/2*k))/(sqrt(((M/2*k)^2)+(N(i,1)/k))));
     sqroot(i,1) = sqrt((((M/(2*k))^2)+(N(i,1)/k)));

    for j = 1:5000
        
        time(j+1,1) = time(j,1)+0.0001;
        
        %Formula 2 - pos1
        SI = (sin(((sqrt(k/m))*time(1,1))+sininv(1,1)));
        MI = M/(2*k);
        
        %MAIN FORMULA
        pos1(j,1) = MI+(sqroot(i,1))*(sin(((sqrt(k/m))*time(j,1))+sininv(i,1)));
    end
end

[Walter Roberson]

        SI = (sin(((sqrt(k/m))*time(1,1))+sininv(1,1)));

Why are you calculating SI there? You do not use the value

  time(j+1,1) = time(j,1)+0.0001;

It would be easier to go before the loops and define

  time = (0:4999).' * 0.0001;

[Sandip Ghatge]

I was just checking for SI Value, whether it gets -1 or not

and regarding time i tried multiple times as below, but still the issue remains unresolved,

%Defining values of parameters
Fg = 2.1188;
Fi = 2.7834;
Ft = Fg + Fi;
h = 0.002;
k = 1864;
e = 0.7;
m = 0.001;

% Initializing for x
x = zeros(11,1);
x(1,1) = 0;
l = 0.002;

%Formula 1 - determining new positions
for i = 1:10
    con1 = ((Ft*h)-(k*(h^2)*0.5))*(1-(e^2));
    x(i+1,1) = (Ft-(sqrt((Ft^2)-(2*k*((con1)+((e^2)*(Ft*x(i,1)-(k*0.5*(x(i,1)^2)))))))))/k;
end

%Initializing for constants N & S
N = zeros(10,1);
S = zeros(10,1);
W = (2*Ft)/m;

%Calculating constants N & S
for i = 1:10
    S(i) = (((W*l)-((k*(l^2))/m))*((e^2-1))-((e^2)*(x(i,1))*(W+((k*(x(i,1)))/m))));
end

sininv = zeros(10,1);
sqroot = zeros(10,1);
% pos1 = zeros(5000,1);
pos2 = zeros(5000,1);
M = 2*Ft; 
time = (0:4999).' * 0.0001;

for i = 1:10
     N(i,1) = ((k*((x(i,1))^2)) - (2*Ft*(x(i,1))));
     sininv(i,1) = asin((x(i,1)-(M/2*k))/(sqrt(((M/2*k)^2)+(N(i,1)/k))));
     sqroot(i,1) = sqrt((((M/(2*k))^2)+(N(i,1)/k)));
     
     MI = M/(2*k);
     
     pos1 = MI+(sqroot(i,1))*(sin(((sqrt(k/m)).*time)+sininv(i,1)));
end

Answer 1

[Stephen23]

Your code has two loops, loop j is nested inside loop i.

You allocate data to pos1 only using index j, which means that you replace the pos1 data ten times (once for each i iteration), and only keep the last value assigned. Here is a small peice of code to place just after the pos1 allocation and you will see the ten values that you allocate and then replace, once for each i iteration:

if j==1
    disp(pos1(1))
end

When I add that code, these are the ten values that you allocate to pos1(1):

0
0.0007348123729260462
0.001229147082940499
0.001551112338708252
0.001750919824749765
0.001867711745054821
0.001932038559298222
0.001965829383562728
0.001983030898594019
0.001991628840290395

The first nine are simply overwritten by the next value. Only the last value is finally stored. That value corresponds to j==1 and i==10, which you did not take into account when you did the calculation by hand. Using the correct i and j values gives exactly the same value as pos1:

>> tmp = (sqroot(1,1))*(sin(((sqrt(k/m))*time(1,1))+sininv(10,1))) % j==10
tmp = -0.002629935049104

[Sandip Ghatge]

Yes, this has solved my problem. Thank you very much

Before accepting can you give me a hint, to plot all these pos1 values against time, because i can still see my 5000*1 array starting with 0.002, when i click on the workspace.

[Stephen23]

Possibly you should use the i index as well:

pos1 = zeros(5000,10);
...
    pos1(j,i) = ...
...
plot(pos1)

But only you know if this is correct for your algorithm.

[Sandip Ghatge]

Yes, i used the i index and was able to plot. Thank you

Answer 2

[Cris LaPierre]

Nothing is wrong with your equation. However, think through your code. The nested for-loops are messing up your logic. The variable pos1 gets overwritten everytime the outerloop increments. You are checking the value of pos1(1,1) once your code has finished executing. That means you should be checking the math using sqroot(i,1) where i=10 (last value of outer loop counter).

SI(1,1) = -1

sqroot(10,1) = 0.0006

MI = 0.0026

pos1(1,1) = 0.0026 + 0.0006*(-1) = 0.002

[Sandip Ghatge]

Thank you very much, i got the point.