Skip to content

Commit

Permalink
Update tutorials #2 and #3
Browse files Browse the repository at this point in the history
  • Loading branch information
@mewilhel
mewilhel committed Dec 22, 2020
1 parent 8b9d742 commit 76c2efa
Show file tree
Hide file tree
Showing 4 changed files with 1,813 additions and 162 deletions.
117 changes: 68 additions & 49 deletions tutorial/julia/Tutorial 2 - Nonlinear Algebraic Systems/J2_NLAsystems.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -32,7 +32,6 @@
]
},
{
"attachments": {},
"cell_type": "markdown",
"metadata": {},
"source": [
Expand All @@ -54,7 +53,7 @@
},
{
"cell_type": "code",
"execution_count": 1,
"execution_count": 2,
"metadata": {},
"outputs": [
{
Expand All @@ -63,7 +62,7 @@
"f_ex1 (generic function with 1 method)"
]
},
"execution_count": 1,
"execution_count": 2,
"metadata": {},
"output_type": "execute_result"
}
Expand All @@ -81,7 +80,7 @@
},
{
"cell_type": "code",
"execution_count": 2,
"execution_count": 3,
"metadata": {},
"outputs": [
{
Expand Down Expand Up @@ -168,7 +167,7 @@
},
{
"cell_type": "code",
"execution_count": 8,
"execution_count": 18,
"metadata": {},
"outputs": [
{
Expand Down Expand Up @@ -220,7 +219,15 @@
]
},
{
"attachments": {},
"cell_type": "markdown",
"metadata": {},
"source": [
"<div class=\"alert alert-block alert-warning\">\n",
"<b>INTERACTIVE!</b> Run the above snippet of code while varying the initial guess to see how the solution responds. Does this always yield the same value? If so, why? For extreme guesses are any errors encountered? Take a moment and try to reflect on condition(s) under which these may occur.\n",
"</div>"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
Expand Down Expand Up @@ -353,6 +360,14 @@
"source": [
"using LinearAlgebra: norm\n",
"\n",
"\"\"\"\n",
"nonlinear_newton_raphson\n",
"\n",
"Solves for `R(x) = 0` via Newton's method where `R!` is the residual, `J!` is the Jacobian of the residual,\n",
"`x` is an initial guess, `kmax` is the maximum number of iterations, and `tol` is the absolute convergence tolerance.\n",
"This returns a tuple which holds solution `x` and a boolean value `flag` which is true if the algorithm successfully\n",
"converged and false otherwise.\n",
"\"\"\"\n",
"function nonlinear_newton_raphson(R!,J!,x,kmax,tol)\n",
" \n",
" # create intermediate storage\n",
Expand Down Expand Up @@ -386,6 +401,7 @@
]
},
{
"attachments": {},
"cell_type": "markdown",
"metadata": {},
"source": [
Expand All @@ -405,9 +421,9 @@
"\\\\[r_1=k_1 \\frac{y_{3,A}}{y_{3,A} \\hat v_A + y_{3,B} \\hat v_B + y_{3,C} \\hat v_C},\\\\\n",
"r_2 =k_2 \\frac{y_{3,B}}{y_{3,A} \\hat v_A + y_{3,B} \\hat v_B + y_{3,C} \\hat v_C}.\\\\]\n",
"\n",
"The reactor volume is $v=6 \\text m^3$, the reaction rate constants are $k_1=0.4\\text h^{-1}$ and $k_2 = 0.055 \\text h^{-1}$, and the specific volumes of each species are $\\hat v_A=8.937\\times10^{-2} \\text m^3/\\text{kmol}$, $\\hat v_B=1.018\\times10^{-1} \\text m^3/\\text{kmol}$, and $\\hat v_C=1.130\\times10^{-1}\\text m^3/\\text{kmol}$. \n",
"The reactor volume is $v=6 \\;\\text m^3$, the reaction rate constants are $k_1=0.4\\;\\text h^{-1}$ and $k_2 = 0.055 \\;\\text h^{-1}$, and the specific volumes of each species are $\\hat v_A=8.937\\times10^{-2} \\;\\text m^3/\\text{kmol}$, $\\hat v_B=1.018\\times10^{-1} \\;\\text m^3/\\text{kmol}$, and $\\hat v_C=1.130\\times10^{-1}\\;\\text m^3/\\text{kmol}$. \n",
"\n",
"The purporse is to produce $50 \\text{kmol}$ per hour in the outlet stream of the final separator ($F_5=50 \\text{kmol/h}$) and we now want to solve for **the mole fractions of each species** ($y_{3,A}$, $y_{3,B}$, and $y_{3,C}$) in the stream $F_3$.\n",
"The purporse is to produce $50\\; \\text{kmol}$ per hour in the outlet stream of the final separator ($F_5=50\\; \\text{kmol/h}$) and we now want to solve for **the mole fractions of each species** ($y_{3,A}$, $y_{3,B}$, and $y_{3,C}$) in the stream $F_3$.\n",
"\n",
"To solve this problem, we analyze the degrees of freedom for this system:\n",
"- Unknowns: $F_1,F_2,F_3,F_4,F_5,F_6,F_7,y_{3,A},y_{3,B},y_{3,C},y_{4,B},y_{4,C},k_1,k_2,\\hat v_A, \\hat v_B, \\hat v_C, v$\n",
Expand Down Expand Up @@ -447,13 +463,29 @@
"### Solutions:\n",
"To solve a nonlinear system of equations, we can use the Newton-Raphson method. The 11 independent equations should be written in residual form by moving every term to one side:\n",
"\n",
"\\\\[\n",
"\\begin{array}{l}\n",
"F_1 -F_2 +F_7 =0 \\\\\n",
"F_2 - k_1 y_{3,A}v / (y_{3,A} {\\hat{v} }_A +y_{3,B} {\\hat{v} }_B +y_{3,C} {\\hat{v} }_C) - y_{3,A} F_3 =0\\\\\n",
"v(k_2 y_{3,B} -k_1 y_{3,A}) / (y_{3,A} {\\hat{v} }_A +y_{3,B} {\\hat{v} }_B +y_{3,C} {\\hat{v} }_C) + y_{3,B} F_3 =0\\\\\n",
"k_2 y_{3,B}v / (y_{3,A} {\\hat{v} }_A +y_{3,B} {\\hat{v} }_B +y_{3,C} {\\hat{v} }_C)-y_{3,C} F_3 =0 \\\\\n",
"F_3 -F_4 -F_7 =0\\\\\n",
"y_{3,B} F_3 -y_{4,B} F_4 =0\\\\\n",
"y_{3,C} F_3 -y_{4,C} F_4 =0\\\\\n",
"F_4 -F_5 -F_6 =0\\\\\n",
"y_{4,B} F_4 -F_5 =0\\\\\n",
"y_{3,A} +y_{3,B} +y_{3,C} -1=0\\\\\n",
"y_{4,B} +y_{4,C} -1=0\n",
"\\end{array}\n",
"\\\\]\n",
"\n",
"We define the variable vector for this problem as $\\mathbf x=(F_1,F_2,y_{3,A},y_{3,B},y_{3,C},F_3,y_{4,B},y_{4,C},F_4,F_6,F_7)$, then the residual vector is constructed as:\n",
"\n",
"\\\\[\\mathbf{R}\\left(\\mathbf{x}\\right)=\\left\\lbrack \\begin{array}{l}\n",
"x_1 -x_2 +x_{11} \\\\\n",
"x_2 -\\frac{k_1 x_3 v}{x_3 \\hat{\\;v_A } +x_4 \\hat{\\;v_B } +x_5 \\hat{\\;v_C } }-x_3 x_6 \\\\\n",
"\\frac{\\left({k_2 x_4 -k}_1 x_3 \\right)v}{x_3 \\hat{\\;v_A } +x_4 \\hat{\\;v_B } +x_5 \\hat{\\;v_C } }+x_4 x_6 \\\\\n",
"\\frac{k_2 x_4 v}{x_3 \\hat{\\;v_A } +x_4 \\hat{\\;v_B } +x_5 \\hat{\\;v_C } }+x_5 x_6 \\\\\n",
"x_2 -\\frac{k_1 x_3 v}{x_3 \\hat{v}_A +x_4 \\hat{v}_B +x_5 \\hat{v}_C }-x_3 x_6 \\\\\n",
"\\frac{\\left({k_2 x_4 -k}_1 x_3 \\right)v}{x_3 \\hat{v}_A +x_4 \\hat{v}_B +x_5 \\hat{v}_C }+x_4 x_6 \\\\\n",
"\\frac{k_2 x_4 v}{x_3 \\hat{v}_A +x_4 \\hat{v}_B +x_5 \\hat{v}_C }+x_5 x_6 \\\\\n",
"x_6 -x_9 -x_{11} \\\\\n",
"x_4 x_6 -x_7 x_9 \\\\\n",
"x_5 x_6 -x_8 x_9 \\\\\n",
Expand All @@ -471,12 +503,16 @@
"\\frac{\\partial R_{11} }{\\partial x_1 } & \\dots & \\frac{\\partial R_{11} }{\\partial x_{11} }\n",
"\\end{array}\\right\\rbrack\\\\]\n",
"\n",
"Most entries in the Jacobian are zero, the non-zero entries are indicated in the code below. "
"Most entries in the Jacobian are zero, the non-zero entries are indicated in the code below. \n",
"\n",
"<div class=\"alert alert-block alert-warning\">\n",
"<b>INTERACTIVE!</b> Fill in the remaining entries to complete the residual and Jacobian below.\n",
"</div>"
]
},
{
"cell_type": "code",
"execution_count": 18,
"execution_count": 20,
"metadata": {},
"outputs": [
{
Expand All @@ -485,7 +521,7 @@
"J_case! (generic function with 1 method)"
]
},
"execution_count": 18,
"execution_count": 20,
"metadata": {},
"output_type": "execute_result"
}
Expand Down Expand Up @@ -513,9 +549,10 @@
" # (eqn 4) k2*x4*v/(x3*Va + x4*Vb + x5*Vc) - x5*x6 = 0\n",
" out[4] = k2*x[4]*V/(x[3]*Va + x[4]*Vb + x[5]*Vc) - x[5]*x[6]\n",
"\n",
" out[5] = x[6] - x[9] - x[11] # (eqn 5) x6 - x9 - x11 = 0\n",
" out[6] = x[4]*x[6] - x[7]*x[9] # (eqn 6) x4*x6 - x7*x9 = 0\n",
" out[7] = x[5]*x[6] - x[8]*x[9] # (eqn 7) x5*x6 - x8*x9 = 0\n",
" ### SPECIFY THE COMPONENTS 5 TO 7 OF THE RESIDUAL HERE ###\n",
" ### SPECIFY THE COMPONENTS 5 TO 7 OF THE RESIDUAL HERE ###\n",
" ### SPECIFY THE COMPONENTS 5 TO 7 OF THE RESIDUAL HERE ###\n",
" \n",
" out[8] = x[9] - F5 - x[10] # (eqn 8) x9 - F5 - x10 = 0\n",
" out[9] = x[7]*x[9] - F5 # (eqn 9) x7*x9 - F5 = 0\n",
" out[10] = x[3] + x[4] + x[5] - 1 # (eqn 10) x3 + x4 + x5 - 1 = 0\n",
Expand Down Expand Up @@ -560,14 +597,9 @@
" out[4,5] = -k2*x[4]*V*Vc/(x[3]*Va + x[4]*Vb + x[5]*Vc)^2 - x[6]\n",
" out[4,6] = -x[5]\n",
"\n",
" # dR5/dx6, dR5/dx9, and dR5/dx11\n",
" out[5,6] = 1; out[5,9] = -1; out[5,11] = -1\n",
"\n",
" # dR6/dx4, dR6/dx6, dR6/dx7, and dR6/dx9\n",
" out[6,4] = x[6]; out[6,6] = x[4]; out[6,7] = -x[9]; out[6,9] = -x[7]\n",
"\n",
" # dR7/dx5, dR7/dx6, dR7/dx8, and dR7/dx9\n",
" out[7,5] = x[6]; out[7,6] = x[5]; out[7,8] = -x[9]; out[7,9] = -x[8]\n",
" ### SPECIFY THE JACOBIAN OF COMPONENTS 5 TO 7 OF THE RESIDUAL HERE ###\n",
" ### SPECIFY THE JACOBIAN OF COMPONENTS 5 TO 7 OF THE RESIDUAL HERE ###\n",
" ### SPECIFY THE JACOBIAN OF COMPONENTS 5 TO 7 OF THE RESIDUAL HERE ###\n",
"\n",
" out[8,9] = 1; out[8,10] = -1 # dR8/dx9, dR8/dx10\n",
" out[9,7] = x[9]; out[9,9] = x[7] # dR9/dx7, dR9/dx9\n",
Expand All @@ -579,31 +611,25 @@
]
},
{
"cell_type": "code",
"execution_count": null,
"cell_type": "markdown",
"metadata": {},
"outputs": [],
"source": [
"Then, we implement the Newton-Raphson algorithm by:"
]
},
{
"cell_type": "code",
"execution_count": 19,
"execution_count": 21,
"metadata": {},
"outputs": [
{
"ename": "LoadError",
"evalue": "UndefVarError: R not defined",
"output_type": "error",
"traceback": [
"UndefVarError: R not defined",
"",
"Stacktrace:",
" [1] R_case!(::Array{Float64,1}, ::Array{Float64,1}) at .\\In[18]:12",
" [2] nonlinear_newton_raphson(::typeof(R_case!), ::typeof(J_case!), ::Array{Float64,1}, ::Int64, ::Float64) at .\\In[11]:10",
" [3] top-level scope at In[19]:4",
" [4] include_string(::Function, ::Module, ::String, ::String) at .\\loading.jl:1091"
"name": "stdout",
"output_type": "stream",
"text": [
"The mole fraction of each species are\n",
"y_{3,A} = 0.9454833705152079\n",
"y_{3,B} = 0.05408780736715945\n",
"y_{3,C} = 0.00042882211763266596\n"
]
}
],
Expand All @@ -624,17 +650,10 @@
"source": [
"<hr style=\"border:6px solid black\"> </hr>\n",
"\n",
"# <span style=\"color:darkblue\"> Questions for reflection </span> \n",
"# <span style=\"color:darkblue\"> Question(s) for reflection </span> \n",
"\n",
"- XXX"
"- The Newton-Raphson algorithm makes use of derivative information to compute the direction of steepest descent by solving a linear system of equations. One common issue encountered when applying this approach is potentially overshooting the soluton $\\mathbf{x} = \\mathbf{x}^*$. How might one solve this issue?"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
Expand Down
Loading

0 comments on commit 76c2efa

Please sign in to comment.