From: Serdar GOKHAN on
Regarding with the DAE example ex12, (Manuscript of 'Solving Index-1 DAEs in MATLAB and Simulink')
http://www.mathworks.com/matlabcentral/fileexchange/7481-manuscript-of-solving-index-1-daes-in-matlab-and-simulink

The PDEs are of the form,

0 = D^2(u1)/Dx^2 - exp(-2*u2)
D(u2)/Dt = D^2(u2)/Dx^2 + exp(-2*u2) - exp(-u1)

The boundary conditions are x = 0 are D(u1)/Dx = exp(-u1) and
u2 = log(p+t). The boundary conditions at x = 1 are u1 = log(1+p+t)
and D(u2)/Dx = exp(-u1).

Regarding with the solution, My first question is ,
1. How is the Mass matrix constructed in the form

M = diag([zeros(N+1,1); ones(N-1,1)]);

In my opnion, since D(u1)/Dt =0,
zeros(N+1,1) ----> would be zeros(N,1) since they are DAE

and, since
f2(1) = u2(1) - log(P + t); ODE
f2(N) = (exp(-u1(N)) - (u2(N) - u2(N-1))/h)/h + exp(-2*u2(N)) + exp(-u1(N)); DAE

ones(N-1,1)---->would be [ones(N,1); zeros(1,1)]
This is my confusion point.

2. How are the f(i) functions are formulated for interior mesh points. Because, when formulated

f1(i) = (u1(i+1) - 2*u1(i) + u1(i-1))/h^2 + (why not - minus ) exp(-2*u2(i));
f2(i) = (u2(i+1) - 2*u2(i) + u2(i-1))/h^2 + exp(-2*u2(i)) + ((why not - minus ))exp(-u1(i));

Any kind of help is highly appreciated.
Thank you