HållfProjektSenaste

%EquilibriumSolverV0.1
clear all, close all, clc
%% Värden
v=70;
R=20;
m=150;
A=0.6;
df=800;
db=200;
h=20;
h1=30;
b=1000;
b1=100;
bd=100;
bb=100;
rd=80;
rb=80;
a1=5;
a2=-8;
c=0.25;
dh=120;
ns=1;
sigmas=(200e6)/4;
g=9.81;
rho=1.21;
v=v*1000/3600;
df=df*1e-3;
db=db*1e-3;
h=h*1e-2;
h1=h1*1e-2;
b=b*1e-3;
b1=b1*1e-3;
bd=bd*1e-3;
bb=bb*1e-3;
rd=rd*1e-3;
rb=rb*1e-3;
dh=dh*1e-3;
d1=b1;
d2=b/2+bd;
d3=b-b1;
d4=b;

%% Jämvikt Fall D
FL=0.5*rho*c*A*(0.7*v)^2;
FD=FL;
Nb=(FD*(h+h1)+m*g*df)/(db+df);
Nf=m*g-Nb;
Vi=m*g*0.5-(h*(m*(0.7*v)^2)/(R*b));
Vy=m*g*0.5+(h*(m*(0.7*v)^2)/(R*b));
Hb=(m*(0.7*v)^2)/(R*(df+db))*df;
Vbi=Vi*(1/(1+(Nf/Nb)));
Vby=Vy*(1/(1+(Nf/Nb)));
Hby=Hb*(1/(1+(Vbi/Vby)));
Hbi=Hb-Hby;
R1x=Hbi+Hby;
Fk=FD*dh*0.5*(1/rd);
R1y=-Fk*(b/2-b1-bd)/(b-2*b1)+(FD*(2*b1-b))/(2*(b-2*b1));
R2y=-Fk*(b/2-b1+bd)/(b-2*b1)-(FD*(b-b1))/(2*(b-2*b1))+(FD*b1*0.5)/(b-2*b1);
R2z=(Vbi*b1+Hbi*0.5*dh+Hby*0.5*dh-Vby*(b-b1))/(b-2*b1);
R1z=(Vby*b1-Vbi*(b-b1)-(Hbi+Hby)*0.5*dh)/(b-2*b1);
%% Snittkrafter
MHbi=Hbi*dh/2;
MHby=Hby*dh/2;

X1y=0:0.001:b1;
X2y=b1+0.001:0.001:(0.5*b+bd);
X3y=(0.5*b+bd+0.001):0.001:(b-b1);
X4y=(b-b1+0.001):0.001:b;
Xiy=[X1y X2y X3y X4y];
%http://pastebin.com/5RgCS36F
%% Fall D

%Z-led
Tz1D=-Vbi*[1:length(X1y)].^0; My1D=(FD/2)*X1y; %från 0 till <b1
Tz2D=(-Vbi-R1z)*[1:length(X2y)].^0; My2D=(FD/2)*X2y+R1y*(X2y-b1);  %från b1 till <0.5*b+bd
Tz3D=(-Vbi-R1z)*[1:length(X3y)].^0; My3D=(FD/2)*X3y+R1y*(X3y-b1)+Fk*(X3y-0.5*b-bd); %Från 0.5*b+bd till <b-b1
Tz4D=(-Vbi-R1z-R2z)*[1:length(X4y)].^0; My4D=0.5*FD*X4y+R1y*(X4y-b1)+Fk*(X4y-0.5*b-bd)+R2y*(X4y-b+b1); %från b-b1 till <b

TzD=[Tz1D Tz2D Tz3D Tz4D]; %dessa ska vi plotta mha av
MyD=[My1D My2D My3D My4D]; %same goes for this badboy


%y-led
Ty1D=-FD/2*([1:length(X1y)].^0); Mz1D=-MHbi-Vbi*X1y; %0 till <b1
Ty2D=((-FD/2)-R1y)*([1:length(X2y)].^0);  Mz2D=-MHbi-Vbi*X2y-R1z*(X2y-b1); %b1 till <b-b1
Ty3D=((-FD/2)-R1y-Fk)*([1:length(X3y)].^0); Mz3D=-MHbi-Vbi*X3y-R1z*(X3y-b1);%Momentberäkning bör göras om %b-b1 till <L
Ty4D=(-FD*0.5-R1y-Fk-R2y)*([1:length(X4y)].^0);  Mz4D=-MHbi-Vbi*X4y-R1z*(X4y-b1)-R2z*(X4y-b+b1);%Educational purpose only

TyD=[Ty1D Ty2D Ty3D Ty4D];
MzD=[Mz1D Mz2D Mz3D Mz4D];

MtotD=sqrt(MyD.^2+MzD.^2);
%%

% Fall C
%% Jämvikt Fall C
vc=0;
FLC=0.5*rho*c*A*(vc)^2;
FDC=FLC+m*a1;
NbC=(FDC*(h+h1)+m*(g*df-a1*h1))/(db+df);
NfC=m*g-NbC;
FkC=(FDC*dh)/(2*rd);
R2yC=(FkC*(d2-d1)+FDC*(d4/2-d1))/(d3-d1);
R1yC=FkC+FDC-R2yC;
R1zC=NbC*(1-(d4/2-d1)/(d3-d1));
R2zC=R1zC;



%z-led
Tz1C=(-FDC/2)*([1:length(X1y)].^0); Mz1C=(FDC/2)*X1y;
Tz2C= (R1yC-FDC/2)*([1:length(X2y)].^0); Mz2C=(R1yC*d1)-(R1yC-FDC/2)*X2y;
Tz3C=(R1yC-FkC-FDC/2)*([1:length(X3y)].^0); Mz3C=R1yC*d1-FkC*d2-(R1yC-FkC-FDC/2)*X3y;
Tz4C=(R1yC+R2yC-FkC-FDC/2)*([1:length(X4y)].^0); Mz4C=R1yC*d1+R2yC*d3-FkC*d2-(R1yC+R2yC-FkC-FDC/2)*X4y;

%y-led
Ty1C=(-NbC/2)*([1:length(X1y)].^0); My1C=(-NbC/2)*X1y;
Ty2C=(R1zC-NbC/2)*([1:length(X2y)].^0); My2C=(R1zC-NbC/2)*X2y-R1zC*d1;
Ty3C=(R1zC-NbC/2)*([1:length(X3y)].^0); My3C=(R1zC-NbC/2)*X3y-R1zC*d1;
Ty4C=(R1zC+R2zC-NbC/2)*([1:length(X4y)].^0); My4C=(R1zC+R2zC-NbC/2)*X4y-R1zC*d1-R2zC*d3;

TzC=[Tz1C Tz2C Tz3C Tz4C];
MzC=[Mz1C Mz2C Mz3C Mz4C];

TyC=[Ty1C Ty2C Ty3C Ty4C];
MyC=[My1C My2C My3C My4C];

MtotC=sqrt(MyC.^2+MzC.^2);


%%
%% FALL B
a1=0;

FLB=0.5*rho*c*A*(v)^2;
FDB=FLB+m*a1;
NbB=(FDB*(h+h1)+m*(g*df-a1*h1))/(db+df);
NfB=m*g-NbB;
FkB=(FDB*dh)/(2*rd);
R2yB=(FkB*(d2-d1)+FDB*((d4/2)-d1))/(d3-d1);
R1yB=FkB+FDB-R2yB;
R1zB=NbB*(1-((d4/2)-d1)/(d3-d1));
R2zB=R1zB;

%z-led
Tz1B=(-FDB/2)*([1:length(X1y)].^0); Mz1B=(FDB/2)*X1y;
Tz2B= (R1yB-FDB/2)*([1:length(X2y)].^0); Mz2B=(R1yB*d1)-(R1yB-FDB/2)*X2y;
Tz3B=(R1yB-FkB-FDB/2)*([1:length(X3y)].^0); Mz3B=R1yB*d1-FkB*d2-(R1yB-FkB-FDB/2)*X3y;
Tz4B=(R1yB+R2yB-FkB-FDB/2)*([1:length(X4y)].^0); Mz4B=R1yB*d1+R2yB*d3-FkB*d2-(R1yB+R2yB-FkB-FDB/2)*X4y;

%y-led
Ty1B=(-NbB/2)*([1:length(X1y)].^0); My1B=(-NbB/2)*X1y;
Ty2B=(R1zB-NbB/2)*([1:length(X2y)].^0); My2B=(R1zB-NbB/2)*X2y-R1zB*d1;
Ty3B=(R1zB-NbB/2)*([1:length(X3y)].^0); My3B=(R1zB-NbB/2)*X3y-R1zB*d1;
Ty4B=(R1zB+R2zB-NbB/2)*([1:length(X4y)].^0); My4B=(R1zB+R2zB-NbB/2)*X4y-R1zB*d1-R2zB*d3;

TzB=[Tz1B Tz2B Tz3B Tz4B];
MzB=[Mz1B Mz2B Mz3B Mz4B];

TyB=[Ty1B Ty2B Ty3B Ty4B];
MyB=[My1B My2B My3B My4B];

MtotB=sqrt(MyB.^2+MzB.^2);



% Fall A
%% Jämvikt Fall A
v=0;
a1=0;
FLA=0.5*rho*c*A*(v)^2;
FDA=FLA+m*a1;
NbA=(FDA*(h+h1)+m*(g*df-a1*h1))/(db+df);
NfA=m*g-NbA;
FkA=(FDA*dh)/(2*rd);
R2yA=(FkA*(d2-d1)+FDA*(d4/2-d1))/(d3-d1);
R1yA=FkA+FDA-R2yA;
R1zA=NbA*(1-(d4/2-d1)/(d3-d1));
R2zA=R1zA;

%z-led
Tz1A=(-FDA/2)*([1:length(X1y)].^0); Mz1A=(FDA/2)*X1y;
Tz2A= (R1yA-FDA/2)*([1:length(X2y)].^0); Mz2A=(R1yA*d1)-(R1yA-FDA/2)*X2y;
Tz3A=(R1yA-FkA-FDA/2)*([1:length(X3y)].^0); Mz3A=R1yA*d1-FkA*d2-(R1yA-FkA-FDA/2)*X3y;
Tz4A=(R1yA+R2yA-FkA-FDA/2)*([1:length(X4y)].^0); Mz4A=R1yA*d1+R2yA*d3-FkA*d2-(R1yA+R2yA-FkA-FDA/2)*X4y;

%y-led
Ty1A=(-NbA/2)*([1:length(X1y)].^0); My1A=(-NbA/2)*X1y;
Ty2A=(R1zA-NbA/2)*([1:length(X2y)].^0); My2A=(R1zA-NbA/2)*X2y-R1zA*d1;
Ty3A=(R1zA-NbA/2)*([1:length(X3y)].^0); My3A=(R1zA-NbA/2)*X3y-R1zA*d1;
Ty4A=(R1zA+R2zA-NbA/2)*([1:length(X4y)].^0); My4A=(R1zA+R2zA-NbA/2)*X4y-R1zA*d1-R2zA*d3;

TzA=[Tz1A Tz2A Tz3A Tz4A];
MzA=[Mz1A Mz2A Mz3A Mz4A];

TyA=[Ty1A Ty2A Ty3A Ty4A];
MyA=[My1A My2A My3A My4A];

MtotA=sqrt(MyA.^2+MzA.^2);


% %% Figure 2 PLOTS
% h2=figure('name','Mtot');
% plot(Xiy,MtotA,'b')
% hold on
% plot(Xiy,MtotB,'r')
% hold on
% plot(Xiy,MtotC,'g')
% hold on
% plot(Xiy,MtotD,'k')
%
% legend('Total moment Fall A', 'Total moment Fall B', 'Total moment Fall C', 'Total moment Fall D')
% %% Figure 1 %PLOTS
% h=figure('name','T och M-plots');
% %Plot 1
% subplot(2,2,1)
% plot(Xiy,TzA,'b')
% hold on
% plot(Xiy,TzB,'r')
% hold on
% plot(Xiy,TzC,'g')
% hold on
% plot(Xiy,TyD,'k')
% legend('Tvärkraft, z-rikt. Fall A', 'Tvärkraft, z-rikt. Fall B', 'Tvärkraft, z-rikt. Fall C', 'Tvärkraft, z-rikt. Fall D')
% %--------------------
% %Plot 2
% subplot(2,2,2)
% plot(Xiy,TyA,'b')
% hold on
% plot(Xiy,TyB,'r')
% hold on
% plot(Xiy,TyC,'g')
% hold on
% plot(Xiy,TzD,'k')
% legend('Tvärkraft, y-rikt. Fall A', 'Tvärkraft, y-rikt. Fall B', 'Tvärkraft, y-rikt. Fall C', 'Tvärkraft, y-rikt. Fall D')
% %--------------------
% %Plot 3
% subplot(2,2,3)
% plot(Xiy,MyA,'b')
% hold on
% plot(Xiy,MyB,'r')
% hold on
% plot(Xiy,MyC,'g')
% hold on
% plot(Xiy,MzD,'k')
% legend('Moment, y-rikt. Fall A', 'Moment, y-rikt. Fall B', 'Moment, y-rikt. Fall C', 'Moment, y-rikt. Fall D')
% %--------------------
% %Plot 4
% subplot(2,2,4)
% plot(Xiy,MzA,'b')
% hold on
% plot(Xiy,MzB,'r')
% hold on
% plot(Xiy,MzC,'g')
% hold on
% plot(Xiy,MyD,'k')
% legend('Moment, z-rikt. Fall A', 'Moment, z-rikt. Fall B', 'Moment, z-rikt. Fall C', 'Moment, z-rikt. Fall D')

%--------------------
%% Vridmoment och normalkrafter
%FALL D
N1D=Hbi*([1:length(X1y)].^0);
N2D=(Hbi-R1x)*([1:length(X2y)].^0);
N3D=(Hby)*([1:length(X3y)].^0);
N4D=(Hby)*([1:length(X4y)].^0);
ND=[N1D N2D N3D N4D];
%Nmax=[N1D N2D N3D N4D];
%Nmax=max(Nmax);
Mv1D=((-FD*dh)/4)*([1:length(X1y)].^0);
Mv2D=(Fk*rd)*([1:length(X2y)].^0);
Mv3D=((Fk*rd)-(-FD*b)/4)*([1:length(X3y)].^0);
Mv4D=((Fk*rd)-(-FD*b)/4)*([1:length(X4y)].^0);
MvTotD=[Mv1D Mv2D Mv3D Mv4D];

%FALL C
Mv1C=((FDC*dh)/4)*([1:length(X1y)].^0);
Mv2C=((FDC*dh)/4)*([1:length(X2y)].^0);
Mv3C=(((FDC*dh)/4)-(FkC*rd))*([1:length(X3y)].^0);
Mv4C=(((FDC*dh)/4)-(FkC*rd))*([1:length(X4y)].^0);
MvTotC=[Mv1C Mv2C Mv3C Mv4C];
TmaxC=sqrt(TzC.^2+TyC.^2);
%skit i detta tillsvidare

%Dimensionering

axelnsDelar=[0.6 1 1 0.6];
diameterD=[];
slutDiameterD=[];
diameterC=[];
slutDiameterC=[];
TmaxD=sqrt(TzD.^2+TyD.^2);
NC=0*([1:length(Xiy)].^0);

ghettoMatris=[101 500 300 100];

%ta fram max värden för de olika grejerna
counterMin=0;
counterMax=0;
for i = 1:4
    rotterD=[];
    rotterC=[];
    counterMin=counterMax+1;
    counterMax=counterMax+ghettoMatris(i);

    for k = counterMin:counterMax %axelns alla delar (aka x axeln)
        %Fall D
        rD=roots([1 0 0 0 -1/(pi^2*sigmas^2)*(16*ND(k).^2+48*TmaxD(k).^2)/(axelnsDelar(i)^4) -(256*MtotD(k).*ND(k)+384*MvTotD(k).*TmaxD(k))/(pi^2*sigmas^2)/(axelnsDelar(i)^5) -(1024*MtotD(k).^2+768*MvTotD(k).^2)/(pi^2*sigmas^2)/(axelnsDelar(i)^6)]);
        rD=rD(imag(rD)==0);
        rotterD=[rotterD rD'];

        %Fall C
        rC=roots([1 0 0 0 -1/(pi^2*sigmas^2)*(16*NC(k).^2+48*TmaxC(k).^2)/(axelnsDelar(i)^4) -(256*MtotC(k).*NC(k)+384*MvTotC(k).*TmaxC(k))/(pi^2*sigmas^2)/(axelnsDelar(i)^5) -(1024*MtotC(k).^2+768*MvTotC(k).^2)/(pi^2*sigmas^2)/(axelnsDelar(i)^6)]);
        rC=rC(imag(rC)==0);
        rotterC=[rotterC rC'];
    end
    %Fall D
    diameterD=max(abs(rotterD));
    diameterD=max(diameterD);
    slutDiameterD=[slutDiameterD diameterD]

    %Fall C
    diameterC=max(abs(rotterC));
    diameterC=max(diameterC);
    slutDiameterC=[slutDiameterC diameterC]
end


%--------------------