%% *************************************MINIPROJEKT**********************************************
%temat projektu: Ogrzewanie w klimatyzowanym domu z piwnica
%autor: Kewin Galuszka
%numer indeksu: 241624




%% WARUNKI POCZ¥TKOWE, ZADANE WARTOŒCI
clear all;
close all;
hold on;

cp = 1000; %ciep³o w³aœciwe powietrza
rop = 1.2; %gêstoœæ powietrza

%przyjêto Kelwiny jako jednostkê temperatury

TzewN = -20 + 273.15; %temperatura na zewenatrz budynku
TwewN = 20 + 273.15; %temperatura wewn¹trz pokoju
TpN = 10 + 273.15; %temperatura wewn¹trz piwnicy
TkzN = 35 + 273.15; %temperatura nawiewu klimatyzacji
TgrN = (TzewN+(-5+273.15))/2 ; %temperatura gruntu 

qkN = 480000/24; %moc klimatyzatora w warunkach za³ozonych
fpN = qkN/(cp*rop*TkzN); %nawiew

%% WSPOLCZYNNIKI PRZEWODNICTWA CIEPLNEGO

%wspolczynnik przewodzenia ciep³a (dom - zewn¹trz)
k1=(cp*rop*fpN*TkzN - cp*rop*fpN*TwewN)/((TwewN-TzewN)+0.55*(TpN-TgrN));

%wspolczynnik przewodzenia ciep³a (piwnica - grunt)
k2= 0.55 * k1;

%wspolczynnik przewodzenia ciep³a (dom - piwnica)
kp = k2*(TpN-TgrN)/(TwewN-TpN); 

%% STANY USTALONE - charakterystyki statyczne i obliczanie punktów
Tzew_1=-20;

Tkz_1 = 35;

    TkzN_z= [Tkz_1 + 273.15,Tkz_1  + 273.15,Tkz_1  + 273.15]; %sta³a temperatura klimatyzacji
    TzewN_z= [-20 + 273.15,-15 + 273.15,-10 + 273.15]; %zmienna temperatura zewnêtrzna
    
    for (i=1:3)
        
    TwewOB(i)= ((cp*rop*fpN)*TkzN_z(i)*(k2+kp)+TzewN_z(i)*k1*(k2+kp)+k2*TgrN*kp)/((cp*rop*fpN)*(k2+kp)+k1*(k2+kp)+k2*kp); 
    TpOB(i)= ((cp*rop*fpN)*(TkzN_z(i)*kp+k2*TgrN)+TzewN_z(i)*k1*kp+k2*TgrN*(k1+kp))/((cp*rop*fpN)*(k2+kp)+k1*(k2+kp)+k2*kp);
    end
    
    figure(1);
    subplot(221);
    plot(TzewN_z-273.15,TwewOB-273.15,'r');
    grid minor;
    xlabel('Temperatura zewnêtrzna [st. C]');
    ylabel('Temperatura w pokoju [st. C]');
    subplot(222);
    plot(TzewN_z-273.15,TpOB-273.15,'r');
    grid minor;
    xlabel('Temperatura zewnêtrzna [st. C]');
    ylabel('Temperatura w piwnicy [st. C]');
    
     for (i=1:3)
    TkzN_z2= [35 + 273.15,40 + 273.15,45 + 273.15]; %zmienna temperatura klimatyzacji
    TzewN_z2= [Tzew_1 + 273.15,Tzew_1 + 273.15,Tzew_1 + 273.15]; %sta³a temperatura zewnêtrzna
    
    TwewOB2(i)= ((cp*rop*fpN)*TkzN_z2(i)*(k2+kp)+TzewN_z2(i)*k1*(k2+kp)+k2*TgrN*kp)/((cp*rop*fpN)*(k2+kp)+k1*(k2+kp)+k2*kp) ;
    TpOB2(i)= ((cp*rop*fpN)*(TkzN_z2(i)*kp+k2*TgrN)+TzewN_z2(i)*k1*kp+k2*TgrN*(k1+kp))/((cp*rop*fpN)*(k2+kp)+k1*(k2+kp)+k2*kp);
     end
     
    subplot(223);
    plot(TkzN_z2-273.15,TwewOB2-273.15,'r');
    grid minor;
    xlabel('Temperatura nawiewu klimatyzacji [st. C]');
    ylabel('Temperatura w pokoju[st. C]');
    subplot(224);
    plot(TkzN_z2-273.15,TpOB2-273.15,'r');
    grid minor;
    xlabel('Temperatura nawiewu klimatyzacji [st. C]');
    ylabel('Temperatura w piwnicy [st. C]');
    
    
%% ZADAWANE DANE, WYMUSZENIA
czas_skok_w =[10000,20000,30000];
czas_skok= 10000;

Tgr = TgrN;        %temperatura gruntu 

Tzew0 = TzewN;     %temperatura na zewenatrz budynku
Twew0 = TwewN;     %temperatura wewn¹trz pokoju
Tp0 = TpN;         %temperatura wewn¹trz piwnicy
Tkz0 = TkzN;       %temperatura nawiewu klimatyzacji UWAGA TUTAJ WEKTOR, DLATEGO BLEDY DALEJ
fp0 = fpN;


dTzew =10;
dTkz = 0;

dfp =0;

Vw = 60;
Vp = 60;
dTg = 0;

Cvw = cp * rop * Vw;
Cvp = cp * rop * Vp; 
Cp0 = cp * rop * fp0;





%% SYMULACJA

for(i_t=1:3)
 sim('miniprojekt_sim_nn');
    figure(2);
    subplot(211);
    plot(t,Tp-273.15,'r');
    hold on;
    grid on;
    xlabel('Czas [s]');
    ylabel('Temperatura [st. C]');
    title('Zmiana temperatury w piwnicy')

    subplot(212);
    plot(t,Twew-273.15,'r');
    hold on;
    grid on;
    xlabel('Czas [s]');
    ylabel('Temperatura [st. C]');
    title('Zmiana temperatury w pokoju')
end
    
   %% Macierz i równania stanu
   fp = fp0;
   for(i_ms=1:3)
       
   Tkz0_ms = TkzN_z(i_ms);
   Tzew0_ms = TzewN_z(i_ms);
   
   X0=[Twew0(i_ms), Tp0(i_ms)]';
   A=[(-cp*rop*fp-k1-kp)/Cvw,kp/Cvw;
      kp/Cvp,(-kp-k2)/Cvp];
   B=[cp*rop*fp/Cvw,k1/Cvw,0;
      0,0,k2/Cvp];
   C = [1,0;0,1];
   D= [0,0,0;0,0,0];
   
   sim('macierz_stanu_nn');
   figure(3);
   subplot(211);
    plot(t_ms,Tp_ms-273.15,'r');
    hold on;
    grid on;
    xlabel('Czas [s]');
    ylabel('Temperatura [st. C]');
    title('Zmiana temperatury w piwnicy (macierz stanu)')

    subplot(212);
    plot(t_ms,Twew_ms-273.15,'r');
    hold on;
    grid on;
    xlabel('Czas [s]');
    ylabel('Temperatura [st. C]');
    title('Zmiana temperatury w pokoju (macierz stanu)')    
   end
    
    %% Transmitancje
    
M = [Cvw * Cvp, Cvw*kp + Cvw*k2 + Cp0*Cvp + Cvp*k1 + Cvp*kp, (kp+k2)*(Cp0) + k1*kp + k1*k2 + kp*k2];

%Pokoj
    L11 = [Cp0*Cvp, Cp0*(kp+k2)];   %Tkz
    L12 = [Cvp*k1, k1*(k2+kp)];     %Tzew
    L13 = [k2*kp];               %Tg

%Piwnica
    L21 = [Cp0*kp];              %Tkz
    L22 = [k1*kp];               %Tzew
    L23 = [Cvw*k2 (k2*Cp0+k1*k2+k2*kp)];%Tg
    
    sim('transmitancje_nn');
   figure(4);
   subplot(211);
   plot(t_tr,Tp_tr-273.15,'r');
    hold on;
    grid on;
    xlabel('Czas [s]');
    ylabel('Temperatura [st. C]');
    title('Zmiana temperatury w piwnicy (transmitancje)')

    subplot(212);
    plot(t_tr,Twew_tr-273.15,'r');
    hold on;
    grid on;
    xlabel('Czas [s]');
    ylabel('Temperatura [st. C]');
    title('Zmiana temperatury w pokoju (transmitancje)')    
    
    %% podsumowanie
    figure(5);
    subplot(321);
    plot(t,Tp-273.15,'r');
    hold on;
    grid on;
    xlabel('Czas [s]');
    ylabel('Temperatura [st. C]');
    title('Zmiana temperatury w piwnicy')

    subplot(322);
    plot(t,Twew-273.15,'r');
    hold on;
    grid on;
    xlabel('Czas [s]');
    ylabel('Temperatura [st. C]');
    title('Zmiana temperatury w pokoju')
    
    for(i_ms=1:3)
       
   Tkz0_ms = TkzN_z(i_ms);
   Tzew0_ms = TzewN_z(i_ms);
   
   X0=[Twew0(i_ms), Tp0(i_ms)]';
   A=[(-cp*rop*fp-k1-kp)/Cvw,kp/Cvw;
      kp/Cvp,(-kp-k2)/Cvp];
   B=[cp*rop*fp/Cvw,k1/Cvw,0;
      0,0,k2/Cvp];
   C = [1,0;0,1];
   D= [0,0,0;0,0,0];
   
   sim('macierz_stanu_nn');
   subplot(323);
    plot(t_ms,Tp_ms-273.15,'r');
    hold on;
    grid on;
    xlabel('Czas [s]');
    ylabel('Temperatura [st. C]');
    title('Zmiana temperatury w piwnicy (macierz stanu)')

    subplot(324);
    plot(t_ms,Twew_ms-273.15,'r');
    hold on;
    grid on;
    xlabel('Czas [s]');
    ylabel('Temperatura [st. C]');
    title('Zmiana temperatury w pokoju (macierz stanu)')    
   end
    
    subplot(325);
   plot(t_tr,Tp_tr-273.15,'r');
    hold on;
    grid on;
    xlabel('Czas [s]');
    ylabel('Temperatura [st. C]');
    title('Zmiana temperatury w piwnicy (transmitancje)')

    subplot(326);
    plot(t_tr,Twew_tr-273.15,'r');
    hold on;
    grid on;
    xlabel('Czas [s]');
    ylabel('Temperatura [st. C]');
    title('Zmiana temperatury w pokoju (transmitancje)')
   