Problem with the Jacobian matrix in a BVP problem I have to solve with the finite difference method the following BVP:
$$-\frac{d}{dx}((1+x)\frac{d}{dx}u(x))=1, \quad x \in (0,1),$$ 
where $$u(0)=u(1)=0.$$
I want to write the equation in the form $$F(\vec{u})=0.$$
So I start discretizing the differential equation: I call $D$ the matrix which discretize the first derivate. 
I got: $$F(\vec{u})=-D*((1+x).*D*\vec{u})-1=0$$ (I wrote with MatLab syntax, hope it's right).
For the boundary conditions there's no problem. The real problem is to write the Jacobian matrix $JF$. I thought it was $$JF=-D1*(c.*D1).$$ But I've seen in MatLab that this is not the right solution.
Any help?
 A: MATLAB is an excellent tool, but you are not assembling the correct equations. I will demonstrate an approach which is slower, but also much safer.
Your equation is
$$ -\left[(1+x)u'(x) \right]'= 1$$
By the product rule of differentiation it is equivalent to 
$$ -(1+x)u''(x) - u'(x) = 1$$
A finite difference discretization can now be read off using the standard space central approximations for $u''$ and $u'$. We have
$$ -(1+x_j) \frac{u_{j+1}- 2u_j + u_{j+1}}{h^2} - \frac{u_{j+1}-u_{j-1}}{2h}, \quad j=1,2,\dotsc,N-1.$$
Here I have assumed a uniform stepsize $h>0$ as well as $Nh=1$. 
It is important to recognize that this is a linear equation in the unknowns $u_j$. In particular, there is no need to apply Newton's method. As you have already pointed out there is no discussion as to the value of $u_0$ and $u_N$ which must be zero.
As a first step, I recommend that you assemble the necessary tridiagonal matrix as a dense matrix by looping over the rows. As soon as you are getting sensible numbers, then switch to sparse matrices for greater computational speed.
