diff options
Diffstat (limited to 'Master/texmf-dist/tex/context/base/m-matrix.mkiv')
-rw-r--r-- | Master/texmf-dist/tex/context/base/m-matrix.mkiv | 495 |
1 files changed, 495 insertions, 0 deletions
diff --git a/Master/texmf-dist/tex/context/base/m-matrix.mkiv b/Master/texmf-dist/tex/context/base/m-matrix.mkiv new file mode 100644 index 00000000000..ccb376e39c7 --- /dev/null +++ b/Master/texmf-dist/tex/context/base/m-matrix.mkiv @@ -0,0 +1,495 @@ +%D \module +%D [ file=m-matrix, +%D version=2014.11.04, % already a year older +%D title=\CONTEXT\ Extra Modules, +%D subtitle=Matrices, +%D author={Jeong Dalyoung \& Hans Hagen}, +%D date=\currentdate, +%D copyright={PRAGMA ADE \& \CONTEXT\ Development Team}] +%C +%C This module is part of the \CONTEXT\ macro||package and is +%C therefore copyrighted by \PRAGMA. See mreadme.pdf for +%C details. + +%D This code is based on a post by Dalyoung on the context list. After that +%D we turned it into a module and improved the code a bit. Feel free to ask +%D us for more. Once we're satisfied, a more general helper l-matrix could +%D be made. Dalyoung does the clever bits, and Hans only cleanes up and +%D optimizes a bit. + +% \registerctxluafile{l-matrix}{1.001} % not yet + +\startmodule[matrix] + +\startluacode + +local settings_to_hash = utilities.parsers.settings_to_hash +local formatters = string.formatters +local copy = table.copy +local insert = table.insert +local remove = table.remove + +local matrix = { } +moduledata.matrix = matrix + +local f_matrix_slot = formatters["%s_{%s%s}"] + +function matrix.symbolic(sym, x, y, nx ,ny) -- symMatrix("a", "m", "n") + local nx = nx or 2 + local ny = ny or nx + local function filled(i,y) + local mrow = { } + for j=1,nx do + mrow[#mrow+1] = f_matrix_slot(sym,i,j) + end + mrow[#mrow+1] = "\\cdots" + mrow[#mrow+1] = f_matrix_slot(sym,i,y) + return mrow + end + local function dummy() + local mrow = { } + for j=1,nx do + mrow[#mrow+1] = "\\vdots" + end + mrow[#mrow+1] = "\\ddots" + mrow[#mrow+1] = "\\vdots" + return mrow + end + -- + local mm = { } + for i=1,ny do + mm[i] = filled(i,y) + end + mm[#mm+1] = dummy() + mm[#mm+1] = filled(x,y) + return mm +end + +-- todo: define a matrix at the tex end so that we have more control + +local fences_p = { + left = "\\left(\\,", + right = "\\,\\right)", +} + +local fences_b = { + left = "\\left[\\,", + right = "\\,\\right]", +} + +function matrix.typeset(m,options) + local options = settings_to_hash(options or "") + context.startmatrix(options.determinant and fences_b or fences_p) + for i=1, #m do + local mi = m[i] + for j=1,#mi do + context.NC(mi[j]) + end + context.NR() + end + context.stopmatrix() +end + +-- interchange two rows (i-th, j-th) + +function matrix.swap(t,i,j) + t[i], t[j] = t[j], t[i] +end + +-- replace i-th row with factor * (i-th row) + +function matrix.multiply(m,i,factor) + local mi = m[i] + for k=1,#mi do + mi[k] = factor * mi[k] + end + return m +end + +-- scalar product "factor * m" + +function matrix.scalar(m, factor) + for i=1,#m do + local mi = m[i] + for j=1,#mi do + mi[j] = factor * mi[j] + end + end + return m +end + +-- replace i-th row with i-th row + factor * (j-th row) + +function matrix.sumrow(m,i,j,factor) + local mi = m[i] + local mj = m[j] + for k=1,#mi do + mi[k] = mi[k] + factor * mj[k] + end +end + +-- transpose of a matrix + +function matrix.transpose(m) + local t = { } + for j=1,#m[1] do + local r = { } + for i=1,#m do + r[i] = m[i][j] + end + t[j] = r + end + return t +end + +-- inner product of two vectors + +function matrix.inner(u,v) + local nu = #u + if nu == 0 then + return 0 + end + local nv = #v + if nv ~= nu then + return 0 + end + local result = 0 + for i=1,nu do + result = result + u[i] * v[i] + end + return result +end + +-- product of two matrices + +function matrix.product(m1,m2) + local product = { } + if #m1[1] == #m2 then + for i=1,#m1 do + local m1i = m1[i] + local mrow = { } + for j=1,#m2[1] do + local temp = 0 + for k=1,#m1[1] do + temp = temp + m1i[k] * m2[k][j] + end + mrow[j] = temp + end + product[i] = mrow + end + end + return product +end + +local function uppertri(m,sign) + local temp = copy(m) + for i=1,#temp-1 do + local pivot = temp[i][i] + if pivot == 0 then + local pRow = i +1 + while temp[pRow][i] == 0 do + pRow = pRow + 1 + if pRow > #temp then -- if there is no nonzero number + return temp + end + end + temp[i], temp[pRow] = temp[pRow], temp[i] + if sign then + sign = -sign + end + end + local mi = temp[i] + for k=i+1, #temp do + local factor = -temp[k][i]/mi[i] + local mk = temp[k] + for l=i,#mk do + mk[l] = mk[l] + factor * mi[l] + end + end + end + if sign then + return temp, sign + else + return temp + end +end + +matrix.uppertri = uppertri + +function matrix.determinant(m) + if #m == #m[1] then + local d = 1 + local t, s = uppertri(m,1) + for i=1,#t do + d = d * t[i][i] + end + return s*d + else + return 0 + end +end + +local function rowechelon(m,r) + local temp = copy(m) + local pRow = 1 + local pCol = 1 + while pRow <= #temp do + local pivot = temp[pRow][pCol] + if pivot == 0 then + local i = pRow + local n = #temp + while temp[i][pCol] == 0 do + i = i + 1 + if i > n then + -- no nonzero number in a column + pCol = pCol + 1 + if pCol > #temp[pRow] then + -- there is no nonzero number in a row + return temp + end + i = pRow + end + end + temp[pRow], temp[i] = temp[i], temp[pRow] + end + local row = temp[pRow] + pivot = row[pCol] + for l=pCol,#row do + row[l] = row[l]/pivot + end + + if r == 1 then + -- make the "reduced row echelon form" + local row = temp[pRow] + for k=1,pRow-1 do + local current = temp[k] + local factor = -current[pCol] + local mk = current + for l=pCol,#mk do + mk[l] = mk[l] + factor * row[l] + end + end + end + -- just make the row echelon form + local row = temp[pRow] + for k=pRow+1, #temp do + local current = temp[k] + local factor = -current[pCol] + local mk = current + for l=pCol,#mk do + mk[l] = mk[l] + factor * row[l] + end + end + pRow = pRow + 1 + pCol = pCol + 1 + + if pRow > #temp or pCol > #temp[1] then + pRow = #temp + 1 + end + end + return temp +end + +matrix.rowechelon = rowechelon +matrix.rowEchelon = rowechelon + +-- solve the linear equation m X = c + +local function solve(m,c) + local n = #m + if n ~= #c then + return copy(m) + end + local newm = copy(m) + local temp = copy(c) + for i=1,n do + insert(newm[i],temp[i]) + end + return rowechelon(newm,1) +end + +matrix.solve = solve + +-- find the inverse matrix of m + +local function inverse(m) + local n = #m + local temp = copy(m) + if n ~= #m[1] then + return temp + end + for i=1,n do + for j=1,n do + insert(temp[i],j == i and 1 or 0) + end + end + temp = rowechelon(temp,1) + for i=1,n do + for j=1,n do + remove(temp[i], 1) + end + end + return temp +end + +matrix.inverse = inverse + +\stopluacode + +\stopmodule + +\unexpanded\def\ctxmodulematrix#1{\ctxlua{moduledata.matrix.#1}} + +\continueifinputfile{m-matrix.mkiv} + +\starttext + +\startluacode +document.DemoMatrixA = { + { 0, 2, 4, -4, 1 }, + { 0, 0, 2, 3, 4 }, + { 2, 2, -6, 2, 4 }, + { 2, 0, -6, 9, 7 }, + { 2, 3, 4, 5, 6 }, + { 6, 6, -6, 6, 6 }, +} + +document.DemoMatrixB = { + { 0, 2, 4, -4, 1 }, + { 0, 0, 2, 3, 4 }, + { 2, 2, -6, 2, 4 }, + { 2, 0, -6, 9, 7 }, + { 2, 2, -6, 2, 4 }, + { 2, 2, -6, 2, 4 }, +} +\stopluacode + +\startsubject[title={A symbolic matrix}] + +\ctxmodulematrix{typeset(moduledata.matrix.symbolic("a", "m", "n"))} +\ctxmodulematrix{typeset(moduledata.matrix.symbolic("a", "m", "n", 4, 8))} + +\stopsubject + +\startsubject[title={Swap two rows (2 and 4)}] + +\startluacode +moduledata.matrix.typeset(document.DemoMatrixA) +context.blank() +moduledata.matrix.swap(document.DemoMatrixA, 2, 4) +context.blank() +moduledata.matrix.typeset(document.DemoMatrixA) +\stopluacode + +\stopsubject + +\startsubject[title={Multiply $3 \times r_2$}] + +\startluacode +moduledata.matrix.typeset(document.DemoMatrixA) +context.blank() +moduledata.matrix.typeset(moduledata.matrix.multiply(document.DemoMatrixA, 2, 3)) +\stopluacode + +\stopsubject + +\startsubject[title={Row 2 + $3 \times r_4$}] + +\startluacode +moduledata.matrix.typeset(document.DemoMatrixA) +context.blank() +moduledata.matrix.sumrow(document.DemoMatrixA, 2, 3, 4) +context.blank() +moduledata.matrix.typeset(document.DemoMatrixA) +\stopluacode + +\stopsubject + +\startsubject[title={Transpose a matrix}] + +\startluacode +moduledata.matrix.typeset(document.DemoMatrixA) +context.blank() +moduledata.matrix.typeset(moduledata.matrix.transpose(document.DemoMatrixA)) +\stopluacode + +\stopsubject + +\startsubject[title={The inner product of two vectors}] + +\startluacode +context(moduledata.matrix.inner({ 1, 2, 3 }, { 3, 1, 2 })) +context.blank() +context(moduledata.matrix.inner({ 1, 2, 3 }, { 3, 1, 2, 4 })) +\stopluacode + +\startsubject[title={The product of two matrices}] + +\startluacode +moduledata.matrix.typeset(document.DemoMatrixA) +context.blank() +moduledata.matrix.typeset(moduledata.matrix.product(document.DemoMatrixA,document.DemoMatrixA)) +\stopluacode + +\stopsubject + +\startsubject[title={An Upper Triangular Matrix}] + +\ctxmodulematrix{typeset(moduledata.matrix.uppertri(document.DemoMatrixB))} + +\startsubject[title={A determinant}] + +\startluacode +local m = { + { 1, 2, 4 }, + { 0, 0, 2 }, + { 2, 2, -6 }, +} +context(moduledata.matrix.determinant(m)) +\stopluacode + +\stopsubject + +\startsubject[title={Row echelon form}] + +\startluacode +local m = { + { 1, 3, -2, 0, 2, 0, 0 }, + { 2, 6, -5, -2, 4, -3, -1 }, + { 0, 0, 5, 10, 0, 15, 5 }, + { 2, 6, 0, 8, 4, 18, 6 }, +} + +moduledata.matrix.typeset(m) +moduledata.matrix.typeset(moduledata.matrix.rowechelon(m,1)) +\stopluacode + +\stopsubject + +\startsubject[title={Solving linear equation}] + +\startluacode +local m = { + { 1, 3, -2, 0 }, + { 2, 0, 1, 2 }, + { 6, -5, -2, 4 }, + { -3, -1, 5, 10 }, +} + +local c = { 5, 2, 6, 8 } + +moduledata.matrix.typeset(moduledata.matrix.solve(m,c)) +\stopluacode + +\stopsubject + +\startsubject[title={Inverse matrix}] + +\startcombination[2*1] + {\ctxlua{moduledata.matrix.typeset { { 1, 1, 1 }, { 0, 2, 3 }, { 3, 2, 1 } }}} {} + {\ctxlua{moduledata.matrix.typeset(moduledata.matrix.inverse { { 1, 1, 1 }, { 0, 2, 3 }, { 3, 2, 1 } })}} {} +\stopcombination + +\stopsubject + +\stoptext |