Problem about applying the Newton's method to a system

Question

Problem The amount of pressure needed to sink an object in a bed of homogeneous soft soil that lies above a base of hard soil can be estimated using the pressure needed to sink smaller objects into the soft soil. Suppose the pressure needed to sink circular plates of radius $r$ into the soft soil to a depth of $d$ is given by

$$p=c_1e^{c_2r}+c_3r$$

with the hard base located at a distance $D>d$ from the surface. Here, the $c_i$'s are constants and $c_2>0$.

Using Newton's method, approximate the values of the constants if you are given that a plate of radius 1 inch requires a pressure of $10lb/in^2$ to sink it into a depth of a foot, a plate of radius 2 inches requires a pressure of $12lb/in^2$ to sink it into the same depth, and a plate of radius $3$ inches requires a pressure of $15lb/in^2$ to sink it into this depth. You will need to implement the Newton's method.

So far what I did: $\textbf{c}=(c_1,c_2,c_3)^T$

\begin{eqnarray*} \textbf{P}(\textbf{c})=\begin{pmatrix}p_1\\ p_2\\ p_3\end{pmatrix}=\begin{pmatrix}c_1e^{c_2}+c_3-10\\ c_1e^{2c_2}+2c_3-12\\ c_1e^{3c_2}+3c_3-15\end{pmatrix} \end{eqnarray*}

and the Jacobian matrix is

$$J_P(\textbf{c})=\begin{pmatrix}e^{c_2} & c_1e^{c_2} & 1\\ e^{2c_2} & 2c_1e^{2c_2} & 2\\ e^{3c_2} & 3c_1e^{3c_2} & 3\end{pmatrix}$$

Then I tried to do the iteration on MATLAB like this:

c0=0.5.*ones(3,1); %initial iterate

tol=1.0e-15;

sol=[c0];

c=c0-J(c0)\P(c0);

while norm(c-c0,inf)>tol
    c=c0-J(c0)\P(c0);
    sol=[sol c];
    c0=c;
end

sol

function p=P(c)

p=[c(1,1).*exp(c(2,1))+c(3,1)-10;
    c(1,1).*exp(2.*c(2,1))+2.*c(3,1)-12;
    c(1,1).*exp(3.*c(2,1))+3.*c(3,1)-15];

function j=J(c)

j=[exp(c(2,1)) c(1,1).*exp(c(2,1)) 1;
    exp(2.*c(2,1)) 2.*c(1,1).*exp(2.*c(2,1)) 2;
    exp(3.*c(2,1)) 3.*c(1,1).*exp(3.*c(2,1)) 3];

But this is not giving me the desired solution. I tried testing with the "for" loop then it showed some warnings about matrix being possibly not properly scaled or being too sparse. What am I doing wrong here?

Answer 1

Subtract twice the first row from the second, then thrice the first row from the third. Solve the two resulting equations for $c_1$ thus obtaining a quadratic for $e^{c_2}$. I get $c_2=\ln\left(\frac{15\pm\sqrt{33}}{16}\right)$ so $c_1=-\frac8{e^{2c_2}-2e^{c_2}}$ and $c_3=10-c_1e^{c_2}$. If I try either of

c0 = [8.8 0.26 -1.4]'
c0 = [9.7 -0.55 4.4]'

in your program it converges to a solution.