% ----- POLITECHNIKA RZESZOWSKA 2001 -----
% Metody komputerowe
% Metoda różnic skończonych
% ----------------------------------------
close all
clear all
% A B C D
% 1 2 3 4 5 6 7 8 9 10 11
% #.___.___.___.___.___.___.___.___.___.___.
% #| /\ /\
% ==== ====
%
% A - utwierdzenie, B i C - podpory przegubowe, D - swobodny koniec
% liczba fragmentów ---> 10 == długo�ć belki 10[m]
% WARTO�CI WEJ�CIOWE: -----------------------------------------------------
n = 50; % liczba przedziałów, na które jest dzielone każdy fragment
L = 1; % długo�ć jednego fragmentu [m]
h = L/n; % długo�ć jednego przedziału [m]
E = 205e9; % moduł Younga [kN/m2]
% NUMEROWANIE: ------------------------------------------------------------
wym = 10*n+1; % poczštkowy wymiar macierzy współczynników równy
% liczbie wszystkich węzłów
numery = [1:wym]; % pomocniczy wektor okre�lajšcy w których punktach
% sš liczone ugięcia
Nu = [1 4 8]; % numery węzłów w których sš podpory
Nu_rz =(Nu-1)*n+1; % numery węzłów podporowych z uwzględnieniem
% liczby podziałów
% OBCIĽŻĘNIE i SZTYWNO�Ć:--------------------------------------------------
x = 1:wym; % wartosci X
I = ((cos(x/n)+2).*0.005).^4; % funkcja opisujšca zależno�ć moment
% bezwładno�ci od wsp. x
q = sin(x/n)+2; % funkcja opisujšca natężenie obcišżenia w
% zależno�ci od x
P = ((h^4)/E)*q./I; % obliczenie wektora prawych stron
% warto�ci (q(x)*h^4)/E*I(x)
% konstruowanie poczštkowej macierzy współczynników: ----------------------
Ws = zeros(wym,wym);
D4 = [1 -4 6 -4 1]; %operator różnicowy d4/dx4
B = ones(wym,1) *D4;
M_D4 = spdiags( B,[-2:2],wym,wym); %macierz pasmowa
Ws = Ws + M_D4;
% Uwzględnianie warunków brzegowych: -------------------------------------
% D) prawy skrajny koniec belki: swobodny. Warunki brzegowe:
% 1. d2w/dx2=0 (moment zg.= 0),stšd ugięcie w pierwszym punkcie (wym+1)
% w(wym+1)=-w(wym-1)+2w(wym)
% 2. d3w/dx3=0 (siła poprz. = 0),stšd ugięcie w drugim punkcie(wym+2)
% w(wym+2)=w(wym-2)-4w(wym-1)+4w(wym)
%------------------------------------------------------------------------
wolny_P = [0 -1 2; % macierz modyfikujšca fragment
1 0 -4]; % macierzy współczynników
Ws(wym-1:wym,wym-2:wym) = Ws(wym-1:wym,wym-2:wym) + wolny_P;
% A) lewy skrajny koniec belki - warunki brzegowe -----------------------
% 1. dw/dx=0 , stšd ugięcie w m punkcie W(0)
% w(0)=w(2);
% -----------------------------------------------------------------------
Ws(1,2) = Ws(1,2)-4;
% Ws(1,3) = Ws(1,3)+1; dla utwierdzenie z przesuwem pionowym
Ws(2,2) = Ws(2,2)+1;
% A), C) i B) podpory , ugięcie równe 0;---------------------------------
Ws(Nu_rz,:)=[]; % odrzucam zbędne wiersze
numery(Nu_rz)=[];
Ws(:,Nu_rz)=[]; % odrzucam zbędne kolumny
P(Nu_rz)=[]; % uwzględnienie podpór w prawej stronie
% ROZWIĽZANIE UKŁADU ROWNAŃ: --------------------------------------------
% Ws*W=Q co odpowiada równaniu linii ugięcia belki: d4/dx4(w)= q/EI
% W = inv(Ws)*P; %
W=-Ws\P'; % Wektor przemieszczeń
EXTremum=[min(W) max(W)] % obliczenie ekstremów
Ug = zeros(wym,1);
Ug(numery) = W;
figure(1)
plot(Ug);
xlim([1 wym]);
hold on;
plot([1 wym],[0 0],':k'); % linia zerowego ugięcia podpory
plot([Nu_rz(:),Nu_rz(:)],[-max(abs(Ug))*0.1 max(abs(Ug))*0.1],':k');
% plot([I*10e3],'-o','LineWidth',1,'MarkerEdgeColor',[1 0 0])
% plot([q/10e2],'-v','LineWidth',1,'MarkerFaceColor',[1 0 1])