Use finites differences for the BVP 
Use finite differences for the BVP $u(x)''=\tan^{-1}(u(x))+2u(x)+\cos(x)$ with $u(0)=u(1)=0$. 

I don't know what I have to do with $\tan^{-1}(u(x))$. When I applied finite differences, I need to have the term $u(x)$ without depending on the function to find the coefficients of the matrix.
 A: Using a standard second order discretization of $u''(x)$ on a grid $0=x_0<x_1<\cdots<x_{N-1}<x_N=1$ with uniform spacing $h = x_{j+1}-x_j$, the linear system of equations to solve becomes $$\frac{u_{j+1} - 2u_j + u_{j-1}}{h^2} = \arctan(u_j) + 2u_j + \cos(x_j), \qquad j=1,2,\ldots,N-1. \tag{*}$$
Let $\mathbf{u}=[u_1,u_2,\ldots,u_{N-1}]^\top$ be the vector of unknowns.  ($u_0=u_N=0$ aren't included, since their value is determined by the boundary condition... same reasoning that $j=0$ and $j=N$ aren't included in (*).)  We can write (*) in the form $$F(x) = 0, \tag{**}$$ where $F : \mathbb{R}^{N-1} \to \mathbb{R}^{N-1}$, and the $j$th component of $F$ is $$(F(x))_j := \frac{u_{j+1} - 2u_j + u_{j-1}}{h^2} - \arctan(u_j) - 2u_j - \cos(x_j).$$
Now we can use a nonlinear solver!  An approximate solution of (**) (using, for example, Newton's method) will give us an approximation to the solution of (*)!  Continuing with the example of Newton's method, we need an initial guess $\mathbf{u}_0 \in \mathbb{R}^{N-1}$.  Then, the iterative scheme (for $\mathbf{u}_0, \mathbf{u}_1, \mathbf{u}_2, \ldots$) is to solve the linear system of equations $$J(\mathbf{u}_n)(\mathbf{u}_{n+1} - \mathbf{u}_n) = -F(\mathbf{u}_n) \tag{***}$$ for $\mathbf{u}_{n+1}$!  (In other words, we start with initial guess $\mathbf{u}_0$ for the solution of (*), compute $J(\mathbf{u}_0)$ and $F(\mathbf{u}_0)$, and solve the linear system of equations (***) for $\mathbf{u}_1$ -- and then repeat for $\mathbf{u}_2$, $\mathbf{u}_3$, etc.!  Each of these iterates -- $\mathbf{u}_1, \mathbf{u}_2, \mathbf{u}_3, \ldots$ -- is an approximation of the solution $\mathbf{u}$ of (*).) 
In (***), $J$ is the Jacobian of $F$ -- in other words, $$J = \begin{bmatrix} \dfrac{\partial F_1}{\partial u_1} & \cdots & \dfrac{\partial F_1}{\partial u_{N-1}} \\ \vdots & \ddots & \vdots \\ \dfrac{\partial F_{N-1}}{\partial u_1} & \cdots & \dfrac{\partial F_{N-1}}{\partial u_{N-1}} \end{bmatrix},$$ where "$F_j$" is shorthand for $(F(x))_j$.
For example, since $$(F(x))_1 = \frac{u_2 - 2u_1}{h^2} - \arctan(u_1) - 2u_1 - \cos(x_1),$$ then the first two elements of the first row of $J$ will be nonzero, while the remaining elements will be zero... $$J_{1,1} = -\dfrac{2}{h^2} - \dfrac{1}{1 + u_1^2} - 2, \qquad J_{1,2} = \dfrac{1}{h^2}, \qquad J_{1,k} = 0 \text{ for } 3 \leq k \leq N-1.$$
