Skip to content

Instantly share code, notes, and snippets.

@svycka

svycka/matlab.m

Last active August 29, 2015 14:07
Show Gist options

  • Select an option

    Learn more about clone URLs
  • Save svycka/5a76d6b66ea0c41b08fe to your computer and use it in GitHub Desktop.

Select an option

Learn more about clone URLs
Save svycka/5a76d6b66ea0c41b08fe to your computer and use it in GitHub Desktop.
Download ZIP
lab1
Raw
matlab.m
function pvz_2
clc, close all, clear all
%KONSTRUKCIJOS DUOMENYS
% MAZGU (MASIU) KONSTANTOS
m1 = 1; m2 = 3;
% ELEMENTU (SPYRUOKLIU) KONSTANTOS
k1 = 2000;
% VEIKIANCIU JEGOS (N)
f1=800; f2=400;
% JEGU VEIKIMO PRADZIA (s)
tf1 = 0.2; tf2 = 0.1;
% MAZGU ITVIRTINIMO MASYVAS
IS = [ 1, 0, 0, 0, 0, 0];
% 2 MAZG0 KINEMATINIIS POVEIKIIS (m)
% u2_deltaU = 0.5;
% % 2 MAZGO KINEMATINIO POVEIKIO PRADZIA (s)
% u2_t_start = 0.2;
% % 2 MAZGO KINEMATINIO POVEIKIO TRUKME (s)
% u2_deltaT = 0.1;
% MAZGAI (MASES) - MATRICU PAVIDALU
% MAZGU MASES (kg)
m= [ 1, m1, m1, m1, m2, m1];
% MAZGUS VEIKIAN?IOS J?GOS (N)
F = [ 0, -f2, 0, 0, 0, -f1];
% JEGU POVEIKIO PRADZIA (s)
tf = [ 0, tf2, 0, 0, 0, tf1];
% ELEMENTAI (SPYRUOKLES)
% STANDUMO KOEFICIENTAI
k=[k1, k1, k1, k1, k1];
% KLAMPIOS TRINTIES KOEFICIENTAI
% c=[10, 10, 10, 10, 10];
c=[ 1, m1, m1, m1, m2, m1];
% MAZGAI (MASES), KURIUOS JUNGIA ELEMENTAI (SPYRUOKLES)
% 3 PAREMTRAS SKIRTAS SPYRUOKLIU ATVAIZDAVIMO AUKSCIUI
ind=[1, 3, 3;
2, 4, 1;
3, 5, 3;
4, 5, 1;
5, 6, 2];
% ATVAIZDAVIMO DUOMENYS
% MAZGU X KORDINATES
x=[2, 2, 4, 4, 6, 8];
% Mazgu Y KOORDINATES
y=[3, 1, 3, 1, 2, 2];
% MASIU DIAMETRAI
rad1=0.2; rad2=1.5;
rads = [rad1, rad1;
rad1, rad1;
rad1, rad1;
rad1, rad1;
rad1, rad2;
rad1, rad1];
nmz=length(m); % MAZGU SKAICIUS (PAVYZDYJE = 6)
nel=length(k); % ELEMENTU SKAICIUS (PAVYZDYJE = 5)
%PRADINIAI GREICIAI IR POSLINKIAI LYGUS 0
U=zeros(nmz,1); % PRADINIAI POSLINKIAI = 0
DU=zeros(nmz,1); % PRADINIAI GREICIAI = 0
% U(2) = u2_deltaU;
% atvaizduojami pradiniai duomenys
% vaizdavimas(x, y, ind, U, rads, 0);
K=zeros(nmz,nmz);
for iii=1:nel
i = ind(iii,1); j = ind(iii,2);
Kel = [k1, -k1;
-k1, k1];
K([i,j],[i,j]) = K([i,j],[i,j])+Kel;
end
F = [ 0, 0, 0, 0, 0, -f1];
K00 = K(find(~IS), find(~IS));
U0 = U(find(~IS));
F0 = F(find(~IS));
K01 = K(find(~IS), find(IS));
U1 = U(find(IS));
% nezinomi poslinkiai
U0 = K00 \ (F0' - K01*U1);
U(find(~IS)) = U0;
U(find(IS)) = U1;
'poslinkiai atlikus skaitin analize'
disp(U');
vaizdavimas(x, y, ind, U, rads, 1);
F = [ 0, -f2, 0, 0, 0, -f1];
U = zeros(nmz,1);
%SKAITINIS INTEGRAVIMAS
Urez=zeros(nmz,1); % rezultatu masyvas, bus pleciamas dinamiskai
TT = 20; %Integravimo laikas
dt=0.0005; % skaitinio integravimo zingsnis
for t=0:dt:TT % skaitinio integravimo ciklas
% PAGREICIU APSKAICIAVIMAS
DDU=pagreitis(m,k,c,ind,F,tf, U,DU,t,@F_laiko_funkcija);
% GREICIU ATNAUJINIMAS
DU=DU+dt*DDU;
% PRISKIRIAMOS GREICIU KRASTINES SALYGOS
DU(find(IS))=0;
% DU(2) = du_laiko_funkcija(u2_t_start, u2_deltaT, u2_deltaU, t);
% POSLINKIU ATNAUJINIMAS
U=U+dt*DU;
% ISSAUGOMI TARPINIAI REZULTATAI
jj=round(t/dt)+1;
% POSLINKIU REIKSMIU ISSAUGOJIMAS
Urez(:,jj)=U;
% ATVAIZDUOJAMA KONSTRUKCIJOS BUSENA
% figure(1); vaizdavimas(x, y, ind, U, rads,t);
% pause(0.01);
end
'poslinkiai atlikus i?reik?tin? analize'
disp(U');
% ATVAIZDUOJAMAS POSLINKIU NUO LAIKO GRAFIKAS
figure(2);hold on; grid on;
axis([0 ,TT,-2 ,0.5]);
spalva={'b-';'r-';'g-';'m-';'c-';'k-';};
for i=1:nmz
plot([0:dt:TT],Urez(i,:),spalva{i});
end
% plot([0.05, 0.05], [-0.2, 0.2], 'g--');
% plot([0.1, 0.1], [-0.2, 0.2], 'k--');
% plot([0.2, 0.2], [-0.2, 0.2], 'r--');
legend('1 mazgas','2 mazgas', '3 mazgas', '4 mazgas', '5 mazgas', '6 mazgas');
xlabel('Laikas (s)');
ylabel('Poslinkiai (m)');
end
% JEGU LAIKO FUNKCIJA
function Ft=F_laiko_funkcija(F, tf, t)
nmz=length(F);
Ft = zeros(nmz, 1);
for i=1:nmz,
if(tf(i) < t)
Ft(i) = F(i);
end
end
if t>0.1 && t<0.2
Ft(2) = F(2);
else
Ft(2)=0;
end
end
% POSLINKU LAIKO FUNKCIJA
function du=du_laiko_funkcija(t_start, deltaT, deltaU, t)
if t >= t_start % JEI POSLINKIS JAU PRASIDEJO
if t <= t_start + deltaT, %JEI POSLINKIS ATLIEKAMAS
% POSLINKIS VYKSTA, APSKAICIUOJAMAS MAZGO GREITIS
omega=(pi/2)/deltaT;
du = deltaU * omega * cos(omega *(t - t_start));
else,
du = 0; % POSLINKIS JAU IVYKO, MAZGO GREITIS = 0
end
else
du = 0; % POSLINKIS DAR NEPRASIDEJO, MAZGO GREITIS = 0
end
return
end
% *********************************************
function DDU=pagreitis(m,k,c,ind,F,tf,U,DU,t,Ft)
nel=length(k);nmz=length(m);
DDU=zeros(nmz,1);
DDU= Ft(F, tf, t); % SUTEIKIAMI JEGU POVEIKIAI
for iii=1:nel % PRIDEDAMOS ELEMENTU VEIKIAMOS JEGOS
i=ind(iii,1);j=ind(iii,2);
T=(U(j)-U(i))*k(iii)+(DU(j)-DU(i))*c(iii);
DDU(i)=DDU(i)+T;
DDU(j)=DDU(j)-T;
end
a = 3;
for iii=1:nmz
DDU(i) = DDU(i)+ a * DU(iii) * m(i);
end
DDU=DDU./m';
return
end
% *********************************************
function vaizdavimas(x,y,ind,U,rads,t)
%PARUOSIAMA FIGURA ATVAIZDAVIMUI
clf; hold on; grid on;axis equal;
axis([min(x)-1 ,max(x)+1,min(y)-1 ,max(y)+1]);
title(['t =', num2str(t), '(s)']);
xlabel('X');
ylabel('Y');
%MAZGU IR ELEMENTU SKAICIUS
nmz=length(x);nel=size(ind,1);
% ATVAIZDUOJAMI ELEMENTAI
for i=1:nel
line([x(ind(i,1))+U(ind(i,1)), ...
x(ind(i,2))+U(ind(i,2)) ], ...
[ind(i,3), ind(i,3)],...
'Color','blue','Linewidth',2);
end
%ATVAIZDUOJAMI MAZGAI
for i=1:nmz
rectangle('Position',[x(i)+U(i)-rads(i,1), ...
y(i)-rads(i,2),rads(i,1)*2,rads(i,2)*2]....
,'Curvature',[1,1],'FaceColor',[0.4 0.6 1]);
end
return
end
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment