Solving a Nonlinear BVP using Finite Difference Method

Question

I am trying to write a code to solve a nonlinear BVP using the Finite Difference Method. The BVP is:
$(T^2)\frac{\partial^2 T}{\partial x^2} + T \left( \frac{\partial T}{\partial x}\right)^2 + Q = 0$
The boundary conditions are Dirichlet at $x = 0$ and Robin at $x = L$. The parameter $Q$ is a constant.

My (probably naive) approach is as follows:

1. Use second-order accurate approximations for the first and second derivatives.
2. Use the discretized form to generate a system of coupled equations for the spatial points.
3. Construct a linear system of the form $\mathbf{A}\mathbf{x} = \mathbf{b}$, where $\mathbf{x}$ is a vector of $T$ at each spatial point $x_0, \ldots, x_L$, $\mathbf{A}$ is a matrix of the coefficients of $T$ arising from the discretization, and $\mathbf{b}$ is a solution vector including the source term $Q$.
4. Solve the linear system for $\mathbf{x}$ somehow, I have not gotten this far.

I really don't know what to do about the nonlinearity. My first problem arises with the squared first derivative. If I apply the finite difference approximation to the first derivative, I get:
$T_i\left(\frac{T_{i+1} - T_{i-1}}{2\Delta x} \right)^2$
Developing the above term yields:
$T_i \left(\frac{T_{i+1}^2-2T_{i+1}T_{i-1}+T_{i-1}^2}{4 \Delta x^2} \right)$

I don't know how to handle this term. I need to linearize something at some point (before or after discretization?) but I am not sure how. And say I figure out a way to handle the squared first derivative, I don't know how to deal with the $T_i$ in front of it.

My question is: how do you discretize and then linearize (or vice-versa) these nonlinear terms?

Thank you for your time. Please let me know if I should elucidate anything.

Answer 1

Thanks to the help in the comments I solved the problem. I divided the spatial domain into $N$ grid points and discretized the PDE (without any sort of linearization), which resulted in $N-1$ nonlinear algebraic equations (as Geoff suggested) for the temperature at each grid point.

Then I employed Newton's method (again, as suggested by Geoff):

1. Implement boundary conditions on the $N-1$ equations
2. Develop a physically reasonable initial guess $F(x_0)$
3. Compute the Jacobian of the $N-1$ equations, resulting in an $N-1 \times N-1$ tri-diagonal matrix
4. Perform Newton iteration to solve for $\Delta T$

The method converged quickly and I am left with the correct temperature distribution.

Thanks for the help with this problem.