(************** Content-type: application/mathematica ************** CreatedBy='Mathematica 5.2' Mathematica-Compatible Notebook This notebook can be used with any Mathematica-compatible application, such as Mathematica, MathReader or Publicon. The data for the notebook starts with the line containing stars above. To get the notebook into a Mathematica-compatible application, do one of the following: * Save the data starting with the line of stars above into a file with a name ending in .nb, then open the file inside the application; * Copy the data starting with the line of stars above to the clipboard, then use the Paste menu command inside the application. Data for notebooks contains only printable 7-bit ASCII and can be sent directly in email or through ftp in text mode. Newlines can be CR, LF or CRLF (Unix, Macintosh or MS-DOS style). NOTE: If you modify the data for this notebook not in a Mathematica- compatible application, you must delete the line below containing the word CacheID, otherwise Mathematica-compatible applications may try to use invalid cache data. For more information on notebooks and Mathematica-compatible applications, contact Wolfram Research: web: http://www.wolfram.com email: info@wolfram.com phone: +1-217-398-0700 (U.S.) Notebook reader applications are available free of charge from Wolfram Research. *******************************************************************) (*CacheID: 232*) (*NotebookFileLineBreakTest NotebookFileLineBreakTest*) (*NotebookOptionsPosition[ 13725, 395]*) (*NotebookOutlinePosition[ 14369, 417]*) (* CellTagsIndexPosition[ 14325, 413]*) (*WindowFrame->Normal*) Notebook[{ Cell["\<\ The following workbook will walk you through the steps necessary to take when \ linearizing a system of ODEs. \ \>", "Text"], Cell["\<\ Let's say we have a system that consists of a CSTR with heat exchange (same \ example modeled by the ODE & Excel CSTR Model w/ Heat Exchange wiki). The \ following equations are used to describe the system:\ \>", "Text"], Cell[CellGroupData[{ Cell[BoxData[ \(dCa\ = \ \(f1\ = \ \((m/\((rho*V)\))\)*\((Ca0 - Ca)\) - \(\(k0\)\(*\)\(Ca\)\(*\)\(Exp[\(-E\)/\((R* T)\)]\)\(\ \)\)\)\)], "Input"], Cell[BoxData[ \(\(-Ca\)\ \[ExponentialE]\^\(-\(\[ExponentialE]\/\(R\ T\)\)\)\ k0 + \ \(\((\(-Ca\) + Ca0)\)\ m\)\/\(rho\ V\)\)], "Output"] }, Open ]], Cell[CellGroupData[{ Cell[BoxData[ \(dT\ = \ \(f2\ = \ \((m*Cp*\((T0 - T)\) - V*H*Ca*k0*Exp[\(-E\)/\((R*T)\)] + U*A*\((Tc - T)\))\)/\((V*rho* Cp)\)\)\)], "Input"], Cell[BoxData[ \(\(Cp\ m\ \((\(-T\) + T0)\) + A\ \((\(-T\) + Tc)\)\ U - Ca\ \ \[ExponentialE]\^\(-\(\[ExponentialE]\/\(R\ T\)\)\)\ H\ k0\ V\)\/\(Cp\ rho\ V\ \)\)], "Output"] }, Open ]], Cell[CellGroupData[{ Cell[BoxData[ \(k\ = \ \(g1\ = \ \ k0*Exp[\(-E\)/\((R*T)\)]\)\)], "Input"], Cell[BoxData[ \(\[ExponentialE]\^\(-\(\[ExponentialE]\/\(R\ T\)\)\)\ k0\)], "Output"] }, Open ]], Cell[TextData[{ "The definitions of the variables used above are as follows:\n\tm = mass \ feed rate and product rate\n\tk0 = rate constant\n\tU = overall heat transfer \ coefficient\n\tA = heat transfer area\n\tT0 = feed temperature\n\tTc = \ coolant temperature\n\tV = reactor volume\n\trho = density of reactor feed\n\t\ H = heat of reaction\n\tCa = concentration of A\n\tT = temperature of reactor\ \n\tR = gas constant\n\tE = activation energy\n\tCp = heat capacity of feed\n\ \nThe reaction being modeled is the irreversible A --> B reaction\n\nIn the \ above equations, dCa describes the change in Ca with respect to time, dT \ describes the change in T with respect to time, and the k equation describes \ the dependance of k on k0. For this example, Ca and T will be treated as \ state variables, which means they are used to describe the conditions of the \ system at any given time. k will be treated as an output variable, which \ means that it is not directly related to the system inputs, but it is related \ to the state variables.\n\nAs was discussed earlier, the first step in \ linearizing these equations is to obtain a truncated Taylor series of each \ function. In the system above, all the variables are all constants except \ for Ca, T, and k. Thus, the partial derivatives we are concerned with are \ only the derivatives with respect to Ca and T. These derivatives can be \ found using the D[ ] function in ", StyleBox["Mathematica", FontSlant->"Italic"], ". This function allows one to find the partial derivatives of a function \ or number of functions with respect to any variable. The syntax employed by \ this function is as follows: D[{f1, f2, f3}, x1]. This expression will \ return a vector of the partial derivatives of f1, f2, and f3 with respect to \ x1.\n\nBecause there are multiple equations and multiple variables in this \ system, the partial derivatives of the functions can be arranged into a \ matrix known as the Jacobian matrix. The Jacobian matrix is arranged as \ follows, for state variables:\n \t\t\t\t\t\t\t\t\t\t\t\t\ J=", Cell[BoxData[ FormBox[ RowBox[{"(", GridBox[{ {\(df1\/dCa\), \(df1\/dT\)}, {\(df2\/dCa\), \(df2\/dT\)} }], ")"}], TraditionalForm]]], "\n " }], "Text", TextAlignment->Left], Cell[CellGroupData[{ Cell[BoxData[ \(M = {{f1}, {f2}}\)], "Input"], Cell[BoxData[ \({{\(-Ca\)\ \[ExponentialE]\^\(-\(\[ExponentialE]\/\(R\ T\)\)\)\ k0 + \(\ \((\(-Ca\) + Ca0)\)\ m\)\/\(rho\ V\)}, {\(Cp\ m\ \((\(-T\) + T0)\) + A\ \ \((\(-T\) + Tc)\)\ U - Ca\ \[ExponentialE]\^\(-\(\[ExponentialE]\/\(R\ \ T\)\)\)\ H\ k0\ V\)\/\(Cp\ rho\ V\)}}\)], "Output"] }, Open ]], Cell[BoxData[""], "Input"], Cell[CellGroupData[{ Cell[BoxData[ \(column1\ = \ D[M, Ca]\)], "Input"], Cell[BoxData[ \({{\(-\[ExponentialE]\^\(-\(\[ExponentialE]\/\(R\ T\)\)\)\)\ k0 - m\/\(rho\ V\)}, {\(-\(\(\[ExponentialE]\^\(-\(\[ExponentialE]\/\(R\ \ T\)\)\)\ H\ k0\)\/\(Cp\ rho\)\)\)}}\)], "Output"] }, Open ]], Cell[CellGroupData[{ Cell[BoxData[ \(column2\ = \ D[M, T]\)], "Input"], Cell[BoxData[ \({{\(-\(\(Ca\ \[ExponentialE]\^\(1 - \[ExponentialE]\/\(R\ T\)\)\ k0\)\/\ \(R\ T\^2\)\)\)}, {\(\(-Cp\)\ m - A\ U - \(Ca\ \[ExponentialE]\^\(1 - \ \[ExponentialE]\/\(R\ T\)\)\ H\ k0\ V\)\/\(R\ T\^2\)\)\/\(Cp\ rho\ V\)}}\)], \ "Output"] }, Open ]], Cell[BoxData[ \(<< \ LinearAlgebra`MatrixManipulation`\)], "Input"], Cell[CellGroupData[{ Cell[BoxData[ \(J\ = \ AppendRows[column1, \ column2]\)], "Input"], Cell[BoxData[ \({{\(-\[ExponentialE]\^\(-\(\[ExponentialE]\/\(R\ T\)\)\)\)\ k0 - m\/\(rho\ V\), \(-\(\(Ca\ \[ExponentialE]\^\(1 - \ \[ExponentialE]\/\(R\ T\)\)\ k0\)\/\(R\ T\^2\)\)\)}, \ {\(-\(\(\[ExponentialE]\^\(-\(\[ExponentialE]\/\(R\ T\)\)\)\ H\ k0\)\/\(Cp\ \ rho\)\)\), \(\(-Cp\)\ m - A\ U - \(Ca\ \[ExponentialE]\^\(1 - \[ExponentialE]\ \/\(R\ T\)\)\ H\ k0\ V\)\/\(R\ T\^2\)\)\/\(Cp\ rho\ V\)}}\)], "Output"] }, Open ]], Cell[CellGroupData[{ Cell[BoxData[ \(MatrixForm[J]\)], "Input"], Cell[BoxData[ TagBox[ RowBox[{"(", "\[NoBreak]", GridBox[{ {\(\(-\[ExponentialE]\^\(-\(\[ExponentialE]\/\(R\ T\)\)\)\)\ k0 - m\/\(rho\ V\)\), \(-\(\(Ca\ \[ExponentialE]\^\(1 - \ \[ExponentialE]\/\(R\ T\)\)\ k0\)\/\(R\ T\^2\)\)\)}, {\(-\(\(\[ExponentialE]\^\(-\(\[ExponentialE]\/\(R\ T\)\)\)\ H\ \ k0\)\/\(Cp\ rho\)\)\), \(\(\(-Cp\)\ m - A\ U - \(Ca\ \[ExponentialE]\^\(1 - \[ExponentialE]\/\(R\ \ T\)\)\ H\ k0\ V\)\/\(R\ T\^2\)\)\/\(Cp\ rho\ V\)\)} }, RowSpacings->1, ColumnSpacings->1, ColumnAlignments->{Left}], "\[NoBreak]", ")"}], Function[ BoxForm`e$, MatrixForm[ BoxForm`e$]]]], "Output"] }, Open ]], Cell["\<\ The above matrix, J, is the Jacobian for the state variables in this problem. \ When using this matrix in linearization applications, T and Ca would be \ substituted with their steady state values (Ts and Cas). The new matrix \ corresponds to the variable \"A\" in the simpler example on the wiki. Once \ this has been done, the Jacobian is then multiplied by the deviation matrix, \ called SS and shown below, which quantifies the difference between the actual \ T and Ca, and their steady state values. \ \>", "Text"], Cell[CellGroupData[{ Cell[BoxData[ \(SS = {{Ca - Cas}, {T - Ts}}\)], "Input"], Cell[BoxData[ \({{Ca - Cas}, {T - Ts}}\)], "Output"] }, Open ]], Cell[CellGroupData[{ Cell[BoxData[ \(MatrixForm[SS]\)], "Input"], Cell[BoxData[ TagBox[ RowBox[{"(", "\[NoBreak]", GridBox[{ {\(Ca - Cas\)}, {\(T - Ts\)} }, RowSpacings->1, ColumnSpacings->1, ColumnAlignments->{Left}], "\[NoBreak]", ")"}], Function[ BoxForm`e$, MatrixForm[ BoxForm`e$]]]], "Output"] }, Open ]], Cell["\<\ To show the behavior of the whole system and how it deviates from steady \ state, the Jacobian and deviation matrices are multiplied together, and yield \ the following result:\ \>", "Text"], Cell[BoxData[ RowBox[{"\t\t\t\t\t", RowBox[{ RowBox[{ RowBox[{"(", GridBox[{ {\(\(d \((Ca - Cas)\)\)\/dt\)}, {\(\(d \((T - Ts)\)\)\/dt\)} }], ")"}], "=", RowBox[{ SubscriptBox[ TagBox[ RowBox[{"(", "\[NoBreak]", GridBox[{ {\(\(-\[ExponentialE]\^\(-\(\[ExponentialE]\/\(R\ T\)\)\ \)\)\ k0 - m\/\(rho\ V\)\), \(-\(\(Ca\ \[ExponentialE]\^\(1 - \ \[ExponentialE]\/\(R\ T\)\)\ k0\)\/\(R\ T\^2\)\)\)}, {\(-\(\(\[ExponentialE]\^\(-\(\[ExponentialE]\/\(R\ T\)\ \)\)\ H\ k0\)\/\(Cp\ rho\)\)\), \(\(\(-Cp\)\ m - A\ U - \(Ca\ \[ExponentialE]\^\(1 - \ \[ExponentialE]\/\(R\ T\)\)\ H\ k0\ V\)\/\(R\ T\^2\)\)\/\(Cp\ rho\ V\)\)} }, RowSpacings->1, ColumnSpacings->1, ColumnAlignments->{Left}], "\[NoBreak]", ")"}], Function[ BoxForm`e$, MatrixForm[ BoxForm`e$]]], \(Cas, Ts\)], RowBox[{"(", GridBox[{ {\(Ca - Cas\)}, {\(T - Ts\)} }], ")"}]}]}], "\[IndentingNewLine]", "\[IndentingNewLine]", "\t\t\t\t\t"}]}]], "Input"], Cell[BoxData[ RowBox[{\(For\ the\ output\ variable\), ",", " ", "k", ",", " ", RowBox[{ "the", " ", "Jacobian", " ", "will", " ", "take", " ", "the", " ", "form", " ", RowBox[{ "of", ":", "\[IndentingNewLine]", " \ ", RowBox[{"(", GridBox[{ {\(dg1\/dCa\), \(dg1\/dT\)} }], ")"}]}]}]}]], "Text"], Cell[CellGroupData[{ Cell[BoxData[ \(Jk\ = \ {D[g1, Ca], D[g1, T]}\)], "Input"], Cell[BoxData[ \({0, \(\[ExponentialE]\^\(1 - \[ExponentialE]\/\(R\ T\)\)\ k0\)\/\(R\ \ T\^2\)}\)], "Output"] }, Open ]], Cell[CellGroupData[{ Cell[BoxData[ \(MatrixForm[Jk]\)], "Input"], Cell[BoxData[ TagBox[ RowBox[{"(", "\[NoBreak]", TagBox[GridBox[{ {"0"}, {\(\(\[ExponentialE]\^\(1 - \[ExponentialE]\/\(R\ T\)\)\ k0\)\/\ \(R\ T\^2\)\)} }, RowSpacings->1, ColumnAlignments->{Left}], Column], "\[NoBreak]", ")"}], Function[ BoxForm`e$, MatrixForm[ BoxForm`e$]]]], "Output"] }, Open ]], Cell[CellGroupData[{ Cell[BoxData[ \(SSk\ = \ {k - ks}\)], "Input"], Cell[BoxData[ \({\[ExponentialE]\^\(-\(\[ExponentialE]\/\(R\ T\)\)\)\ k0 - ks}\)], "Output"] }, Open ]], Cell[CellGroupData[{ Cell[BoxData[ \(MatrixForm[SSk]\)], "Input"], Cell[BoxData[ TagBox[ RowBox[{"(", "\[NoBreak]", TagBox[GridBox[{ {\(\[ExponentialE]\^\(-\(\[ExponentialE]\/\(R\ T\)\)\)\ k0 - ks\)} }, RowSpacings->1, ColumnAlignments->{Left}], Column], "\[NoBreak]", ")"}], Function[ BoxForm`e$, MatrixForm[ BoxForm`e$]]]], "Output"] }, Open ]], Cell[BoxData[ RowBox[{\(To\ show\ the\ behavior\ of\ k\ and\ how\ it\ deviates\ from\ \ it' s\ steady\ state\ value\), ",", " ", \(the\ Jacobian\ and\ deviation\ matrices\ are\ multiplied\ \ together\), ",", " ", RowBox[{\(and\ yield\ the\ following\ \(result : \[IndentingNewLine]\t\t\ \t\t\t\t\t\t\((k - ks)\)\)\), " ", "=", " ", RowBox[{ SubscriptBox[ TagBox[ RowBox[{"(", "\[NoBreak]", TagBox[GridBox[{ {"0"}, {\(\(\[ExponentialE]\^\(1 - \[ExponentialE]\/\(R\ T\)\)\ \ k0\)\/\(R\ T\^2\)\)} }, RowSpacings->1, ColumnAlignments->{Left}], Column], "\[NoBreak]", ")"}], Function[ BoxForm`e$, MatrixForm[ BoxForm`e$]]], \(Ts, Cas\)], "*", TagBox[ RowBox[{"(", "\[NoBreak]", TagBox[GridBox[{ {\(\[ExponentialE]\^\(-\(\[ExponentialE]\/\(R\ T\)\)\)\ \ k0 - ks\)} }, RowSpacings->1, ColumnAlignments->{Left}], Column], "\[NoBreak]", ")"}], Function[ BoxForm`e$, MatrixForm[ BoxForm`e$]]]}]}]}]], "Text"] }, FrontEndVersion->"5.2 for Microsoft Windows", ScreenRectangle->{{0, 1280}, {0, 941}}, WindowSize->{1272, 907}, WindowMargins->{{0, Automatic}, {Automatic, 0}} ] (******************************************************************* Cached data follows. If you edit this Notebook file directly, not using Mathematica, you must remove the line containing CacheID at the top of the file. The cache data will then be recreated when you save this file from within Mathematica. *******************************************************************) (*CellTagsOutline CellTagsIndex->{} *) (*CellTagsIndex CellTagsIndex->{} *) (*NotebookFileOutline Notebook[{ Cell[1754, 51, 135, 3, 33, "Text"], Cell[1892, 56, 231, 4, 33, "Text"], Cell[CellGroupData[{ Cell[2148, 64, 188, 3, 30, "Input"], Cell[2339, 69, 141, 2, 44, "Output"] }, Open ]], Cell[CellGroupData[{ Cell[2517, 76, 176, 3, 30, "Input"], Cell[2696, 81, 176, 3, 53, "Output"] }, Open ]], Cell[CellGroupData[{ Cell[2909, 89, 80, 1, 30, "Input"], Cell[2992, 92, 89, 1, 35, "Output"] }, Open ]], Cell[3096, 96, 2353, 39, 576, "Text"], Cell[CellGroupData[{ Cell[5474, 139, 49, 1, 30, "Input"], Cell[5526, 142, 291, 4, 53, "Output"] }, Open ]], Cell[5832, 149, 26, 0, 30, "Input"], Cell[CellGroupData[{ Cell[5883, 153, 55, 1, 30, "Input"], Cell[5941, 156, 215, 3, 53, "Output"] }, Open ]], Cell[CellGroupData[{ Cell[6193, 164, 54, 1, 30, "Input"], Cell[6250, 167, 253, 4, 64, "Output"] }, Open ]], Cell[6518, 174, 71, 1, 30, "Input"], Cell[CellGroupData[{ Cell[6614, 179, 71, 1, 30, "Input"], Cell[6688, 182, 430, 6, 64, "Output"] }, Open ]], Cell[CellGroupData[{ Cell[7155, 193, 46, 1, 30, "Input"], Cell[7204, 196, 723, 15, 109, "Output"] }, Open ]], Cell[7942, 214, 532, 9, 90, "Text"], Cell[CellGroupData[{ Cell[8499, 227, 60, 1, 30, "Input"], Cell[8562, 230, 56, 1, 29, "Output"] }, Open ]], Cell[CellGroupData[{ Cell[8655, 236, 47, 1, 30, "Input"], Cell[8705, 239, 324, 10, 55, "Output"] }, Open ]], Cell[9044, 252, 200, 4, 33, "Text"], Cell[9247, 258, 1328, 30, 146, "Input"], Cell[10578, 290, 483, 11, 56, "Text"], Cell[CellGroupData[{ Cell[11086, 305, 63, 1, 30, "Input"], Cell[11152, 308, 112, 2, 52, "Output"] }, Open ]], Cell[CellGroupData[{ Cell[11301, 315, 47, 1, 30, "Input"], Cell[11351, 318, 398, 12, 75, "Output"] }, Open ]], Cell[CellGroupData[{ Cell[11786, 335, 51, 1, 30, "Input"], Cell[11840, 338, 105, 2, 39, "Output"] }, Open ]], Cell[CellGroupData[{ Cell[11982, 345, 48, 1, 30, "Input"], Cell[12033, 348, 383, 11, 49, "Output"] }, Open ]], Cell[12431, 362, 1290, 31, 88, "Text"] } ] *) (******************************************************************* End of Mathematica Notebook file. *******************************************************************)