% Musterlösung
%
% Numerische Mathematik, SS 2013
%
% Simulation eines Wasserkreislaufes
%
% Ingolf Busch, Juni 2013

% Input: t_end  Anzahl der Zeitschritte
%        p      Parameter p
%        c      Parameter [c_1,c_2,...,c_5]
%        S      Startwerte [U_0,V_0,W_0,R_0]

% Output: t  Zeitschritte
%         U  Füllmengen in U
%         V  Füllmengen in V
%         W  Füllmengen in W
%         R  Füllmengen in R

function [t,U,V,W,R] = my_wasser(t_end,p,c,S)

% Initialisierung
t = 0:t_end;
U = zeros(t_end+1,1);
V = zeros(t_end+1,1);
W = zeros(t_end+1,1);
R = zeros(t_end+1,1);

%Anfangswerte
U(1) = S(1); 
V(1) = S(2); 
W(1) = S(3); 
R(1) = S(4);

% Zeitschritte

for n=1:t_end
  U(n+1) = U(n) - (c(1)+c(2))*sqrt(U(n)) + p;
  V(n+1) = V(n) + c(1)*sqrt(U(n)) -(c(3)+c(4))*sqrt(V(n));
  W(n+1) = W(n) + c(2)*sqrt(U(n)) + c(4)*sqrt(V(n)) - c(5)*sqrt(W(n));
  R(n+1) = R(n) + c(3)*sqrt(V(n)) + c(5)*sqrt(W(n)) - p;
end