\n",
"```\n",
"\n",
"```{raw} latex\n",
"%%start demo%%\n",
"```\n",
"\n",
"::::{panels}\n",
":column: col-7 left-side\n",
":card: border-0 shadow-none\n",
"```{raw} latex\n",
"\\begin{minipage}[t]{0.5\\textwidth}\n",
"```\n",
"A system of nonlinear equations is defined by its residual and Jacobian.\n",
"\n",
"```{raw} latex\n",
"\\end{minipage}\\hfill\n",
"```\n",
"---\n",
":column: col-5 right-side\n",
":card: shadow-none comment\n",
"```{raw} latex\n",
"\\begin{minipage}[t]{0.4\\textwidth}\\begin{mdframed}[default]\\small\n",
"```\n",
"Be careful when coding a Jacobian all in one statement. Spaces separate columns, so `x[3]-1` is not the same as `x[3] - 1`.\n",
"```{raw} latex\n",
"\\end{mdframed}\\end{minipage}\n",
"```\n",
"::::"
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "f8be8083",
"metadata": {},
"outputs": [],
"source": [
"function func(x)\n",
" [ exp(x[2]-x[1]) - 2,\n",
" x[1]*x[2] + x[3],\n",
" x[2]*x[3] + x[1]^2 - x[2]\n",
" ];\n",
"end;\n",
" \n",
"function jac(x)\n",
" [ \n",
" -exp(x[2]-x[1]) exp(x[2]-x[1]) 0\n",
" x[2] x[1] 1\n",
" 2*x[1] x[3]-1 x[2]\n",
" ];\n",
"end;"
]
},
{
"cell_type": "markdown",
"id": "4d675973",
"metadata": {},
"source": [
"We will use a `BigFloat` starting value, and commensurately small stopping tolerances, in order to get a sequence long enough to measure convergence."
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "f2387c9b",
"metadata": {},
"outputs": [],
"source": [
"x₁ = BigFloat.([0,0,0])\n",
"ϵ = eps(BigFloat)\n",
"x = FNC.newtonsys(func,jac,x₁,xtol=ϵ,ftol=ϵ);"
]
},
{
"cell_type": "markdown",
"id": "9259e349",
"metadata": {},
"source": [
"Let's compute the residual of the last result in order to check the quality."
]
},
{
"cell_type": "code",
"execution_count": 4,
"id": "c835ba83",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"residual = norm(func(r)) = 0.0\n"
]
}
],
"source": [
"r = x[end]\n",
"@show residual = norm(func(r));"
]
},
{
"cell_type": "markdown",
"id": "5bcc02e5",
"metadata": {},
"source": [
"We take the sequence of norms of errors, applying the log so that we can look at the exponents."
]
},
{
"cell_type": "code",
"execution_count": 5,
"id": "6c3784ab",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"7-element Vector{Float64}:\n",
" 0.7937993447128696\n",
" 3.6959808854483027\n",
" 2.4326597889977153\n",
" 2.3110932374368063\n",
" 2.1325411149310054\n",
" 2.029767250340732\n",
" 2.0340010945857525"
]
},
"execution_count": 5,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"logerr = [ Float64( log(norm(r-x[k])) ) for k in 1:length(x)-1 ]\n",
"[ logerr[k+1]/logerr[k] for k in 1:length(logerr)-1 ]"
]
},
{
"cell_type": "markdown",
"id": "33cc73e7",
"metadata": {},
"source": [
"The ratio is neatly converging toward 2, which is expected for quadratic convergence.\n",
"\n",
"```{raw} html\n",
"

\n",
"```\n",
"\n",
"```{raw} latex\n",
"%%end demo%%\n",
"```\n",
"\n",
"## Exercises\n",
"\n",
"(problem-newtonsys-byhand)=\n",
"1. ✍ Suppose that\n",
" \n",
" ```{math}\n",
" \\mathbf{f}(\\mathbf{x}) =\n",
" \\begin{bmatrix}\n",
" x_1x_2+x_2^2-1 \\\\[1mm] x_1x_2^3 + x_1^2x_2^2 + 1\n",
" \\end{bmatrix}.\n",
" ```\n",
"\n",
" Let $\\mathbf{x}_1=[-2,1]^T$. Use Newton's method to find $\\mathbf{x}_2$.\n",
"\n",
"2. ✍ Suppose that $\\mathbf{f}(\\mathbf{x}) = \\mathbf{A}\\mathbf{x} - \\mathbf{b}$ for a constant $n\\times n$ matrix $\\mathbf{A}$ and constant $n\\times 1$ vector $\\mathbf{b}$. Show that Newton's method converges to the exact root in one iteration.\n",
"\n",
" (problem-newtonsys-spherepotential)=\n",
"3. Two curves in the $(u,v)$ plane are defined implicitly by the equations $u\\log u + v \\log v = -0.3$ and $u^4 + v^2 = 1$.\n",
" \n",
" **(a)** ✍ Write the intersection of these curves in the form $\\mathbf{f}(\\mathbf{x}) = \\boldsymbol{0}$ for two-dimensional $\\mathbf{f}$ and $\\mathbf{x}$.\n",
"\n",
" **(b)** ✍ Find the Jacobian matrix of $\\mathbf{f}$.\n",
"\n",
" **(c)** ⌨ Use {numref}`Function {number}