function out = sysDAE(w,u,d)

dwdt = -0.25*(2*u+w);                                   % differential state
y    = u + w + d;                                         % algebraic state
out  = [dwdt y];                                              % DAE system