clc;
clear;
close all;


%% define parameters

%Q=14.6358;
%S=0.0001;
%n=0.02;
%epsilon=0.01;

Q=str2double(answer1);
S=str2double(answer2);
n=str2double(answer3);
epsilon=str2double(answer4);
N=str2double(answer5);
Z=str2double(answer6);
Y=str2double(answer7);
B=str2double(answer8);

q=Q/B,
yc=(q*q/9.81)^(1/3),
sc=((n*n*Q*Q)/((B*B*yc*yc)*(((B*yc)/(B+2*yc))^(4/3)))),

%mohasebeh omgh normal;
yn=0.01
A=((1/n)*(B*yn)*(((B*yn)/(b+2*yn))^(2/3))*S);
while( abs(A-Q) > 0.01)
    yn=yn+0.0001
    A=((1/n)*(b*yn)*(((b*yn)/(b+2*yn))^(2/3))*S);
    
end

if(S < sc )
    %monhani M ast;
   if (Y< yc)
    %monhani M3 ast; 
       go to line  77; 
     
   elseif(Y<yn)
    %monhani M2 ast; 
     go to line  108;
   else
    %monhani M1 ast;    
       go to line  77 ;
   end
   
    
end
if(S >sc)
    %monhani S ast;
    if (Y< yn)
    %monhani S3 ast; 
       go to line  77; 
     
   elseif(Y<yc)
    %monhani S2 ast; 
     go to line  ;
   else
    %monhani S1 ast;    
       go to line  61 ;
   end
    
end
if(S<0 &S=0)
    if(Y > yc)
      %monhani H2 ya A2 ast;    
       go to line  61 ;      
        
    end
     %monhani H3 ya A3 ast; 
      go to line  61 ; 
end


i=1;
while (x(i)< 15)
    
    y(i+1)=y(i)+0.1,
    z(i+1)=z(i)+(X*S),
    H(i)=z(i)+y(i)+(1.15*(Q*Q)/(2*g*(B*B*y(i)*y(i)))),
    H(i+1)=z(i+1)+y(i+1)+(1.15*(Q*Q)/(2*g*(B*B*y(i+1)*y(i+1)))),
    Sf(i)=((n*n*Q*Q)/((B*B*y(i)*y(i))*(((B*y(i))/(B+2*y(i)))^(4/3)))),
    Sf(i+1)=((n*n*Q*Q)/((B*B*y(i+1)*y(i+1))*(((B*y(i+1))/B+2*y(i+1))^(4/3)))),
    hf(i)=((Sf(i)+Sf(i+1))/2)*X,
    G(i+1)=H(i)+hf(i),
    while (abs(G(i+1)-H(i+1))> epsilon)
        if((G(i+1)-H(i+1))<0)
         y(i+1)=y(i+1)-0.01 ,  
        else
          y(i+1)=y(i+1)+0.01 , 
        end
        H(i+1)=z(i+1)+y(i+1)+(1.15*(Q*Q)/(2*g*(B*B*y(i+1)*y(i+1)))),
        Sf(i+1)=((n*n*Q*Q)/((B*B*y(i+1)*y(i+1))*(((B*y(i+1))/B+2*y(i+1))^(4/3)))),
        hf(i)=((Sf(i)+Sf(i+1))/2)*X,
        G(i+1)=H(i)+hf(i),
    end
   

 i=i+1;
 x(i)= X+x(i-1);
    
end


i=1;
while (x(i)< 15)
    
    y(i+1)=y(i)-0.1,
    z(i+1)=z(i)+(X*S),
    H(i)=z(i)+y(i)+(1.15*(Q*Q)/(2*g*(B*B*y(i)*y(i)))),
    H(i+1)=z(i+1)+y(i+1)+(1.15*(Q*Q)/(2*g*(B*B*y(i+1)*y(i+1)))),
    Sf(i)=((n*n*Q*Q)/((B*B*y(i)*y(i))*(((B*y(i))/(B+2*y(i)))^(4/3)))),
    Sf(i+1)=((n*n*Q*Q)/((B*B*y(i+1)*y(i+1))*(((B*y(i+1))/B+2*y(i+1))^(4/3)))),
    hf(i)=((Sf(i)+Sf(i+1))/2)*X,
    G(i+1)=H(i)+hf(i),
    while (abs(G(i+1)-H(i+1))> epsilon)
        if((G(i+1)-H(i+1))<0)
         y(i+1)=y(i+1)-0.01 ,  
        else
          y(i+1)=y(i+1)+0.01 , 
        end
        H(i+1)=z(i+1)+y(i+1)+(1.15*(Q*Q)/(2*g*(B*B*y(i+1)*y(i+1)))),
        Sf(i+1)=((n*n*Q*Q)/((B*B*y(i+1)*y(i+1))*(((B*y(i+1))/B+2*y(i+1))^(4/3)))),
        hf(i)=((Sf(i)+Sf(i+1))/2)*X,
        G(i+1)=H(i)+hf(i),
    end
   

 i=i+1;
 x(i)= X+x(i-1);
    
end
   
 figure;
 
plot(-x,y);
xlabel('x');
ylabel('y');
%semilogy(x,y,'LineWidth',2);
grid on;