In this particular application Java code is the final output.\u00a0 However, it is trivial to convert this example to output C\/C++ code.<\/p>\n
The Math<\/h2>\n
To demonstrate these techniques a complex real-world problem is used.\u00a0 This problem comes from geometric computer vision and understanding all the math involved is not necessary, but you can find the paper it is based on here [1].\u00a0 Having an understanding of linear algebra will help you understand all the math terminology.\u00a0 Briefly stated the problem is to solve for a 3×3 Essential matrix given 5 point correspondences. \u00a0 If that means nothing to you, don’t worry.<\/p>\n
The first step in solving this problem is computing the null space of a 5 by 9 linear system.\u00a0 The null space is then stored in four column vectors 
where Again if you don’t understand why we are doing the math, don’t panic.\u00a0 Just know that there are some ugly complex equations which need to be solved for in Java.\u00a0 In summary we need to do the following:<\/p>\n This is actually only 2\/3 of the solution.\u00a0 Since the focus of this article is on the process and not about the details of this particular problem we will skip the rest.<\/p>\n All code examples below are written in Python and need to be run inside of Sage’s interactive shell.\u00a0 Do not try to to run the scripts from Python directly, use Sage.\u00a0 All python files need the following imports:<\/p>\n The first step in solving this problem is creating several 3×3 symbolic matrices.\u00a0 Each element in the matrix is unique.\u00a0 I couldn’t find a function built into Sage, so I wrote my own:<\/p>\n The ‘var’ function is used by Sage to define symbolic variables.\u00a0 Inside of the Sage interactive shell you can do the following:<\/p>\n Now that you can create a symbolic matrix you can solve rest of the problem as follows:<\/p>\n The functions det<\/strong> and trace<\/strong> compute the determinant and trace of a matrix respectively, A.T<\/strong> is the transpose of A<\/strong>.\u00a0 Most of that is straight forward algebra, but what about that last bit? Well the printData() function massages those equations to produce usable Java output, which is the subject of the next section.<\/p>\n What the Java code is going to do is solve a linear system. To do that the equations above need to reformatted and converted into matrix format. The variable eqs<\/strong> in the above python code contains equations which have a similar appearance to the equation below on the left, but they need to be expanded out so that they look like the equation on the right.<\/p>\n After the equations have been expanded, the knowns and unknowns are separated and put into matrix format. Consider the system of equations below:<\/p>\n After it has been converted into matrix format it will look like:<\/p>\n At this point the unknowns can be solved for using standard linear algebra tools, e.g. Ax=B.<\/p>\n Part of this process is shown below.\u00a0 First eq1<\/strong> is expanded out and then the coefficients of The first two statements are just there to show you how long and ugly these equations are.\u00a0 Below is the python code for extractVarEq<\/strong>:<\/p>\n The function expandPower<\/strong>(w) converts statements like a^3 into a*a*a which can be parsed by Java.\u00a0 After processing all 10 equations most of the hard work is done,\u00a0\u00a0 Addition processing is done to simplify the equations, but the code is a bit convoluted and isn’t discussed here.\u00a0 The final output it sent into a text file and then pasted into a Java application.<\/p>\n Want to see what the final output looks like? Well here is a small sniplet:<\/p>\n The full source code can be downloaded from:<\/p>\n [1]<\/strong> David Nister “An Efficient Solution to the Five-Point Relative Pose Problem” Pattern Analysis and Machine Intelligence, 2004<\/em><\/p>\n","protected":false},"excerpt":{"rendered":" Abstract The following article demonstrates how Sage\/Python can be used to solve, format, and output auto generated code for complex symbolic equations. \u00a0 The auto-generated code can be targeted for languages such as C\/C++ or Java.\u00a0 Sage is a computer algebra system (CAS) that uses the Python programming language.\u00a0 By using Python it is easy […]<\/p>\n","protected":false},"author":2,"featured_media":0,"comment_status":"open","ping_status":"open","sticky":false,"template":"","format":"standard","meta":{"_monsterinsights_skip_tracking":false,"_monsterinsights_sitenote_active":false,"_monsterinsights_sitenote_note":"","_monsterinsights_sitenote_category":0,"footnotes":""},"categories":[4],"tags":[16,19,14,15],"class_list":["post-107","post","type-post","status-publish","format-standard","hentry","category-programming","tag-math","tag-python","tag-sage","tag-symbolic"],"jetpack_featured_media_url":"","_links":{"self":[{"href":"https:\/\/peterabeles.com\/blog\/index.php?rest_route=\/wp\/v2\/posts\/107","targetHints":{"allow":["GET"]}}],"collection":[{"href":"https:\/\/peterabeles.com\/blog\/index.php?rest_route=\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/peterabeles.com\/blog\/index.php?rest_route=\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/peterabeles.com\/blog\/index.php?rest_route=\/wp\/v2\/users\/2"}],"replies":[{"embeddable":true,"href":"https:\/\/peterabeles.com\/blog\/index.php?rest_route=%2Fwp%2Fv2%2Fcomments&post=107"}],"version-history":[{"count":9,"href":"https:\/\/peterabeles.com\/blog\/index.php?rest_route=\/wp\/v2\/posts\/107\/revisions"}],"predecessor-version":[{"id":155,"href":"https:\/\/peterabeles.com\/blog\/index.php?rest_route=\/wp\/v2\/posts\/107\/revisions\/155"}],"wp:attachment":[{"href":"https:\/\/peterabeles.com\/blog\/index.php?rest_route=%2Fwp%2Fv2%2Fmedia&parent=107"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/peterabeles.com\/blog\/index.php?rest_route=%2Fwp%2Fv2%2Fcategories&post=107"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/peterabeles.com\/blog\/index.php?rest_route=%2Fwp%2Fv2%2Ftags&post=107"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}
\n
Solving with Sage<\/h2>\n
from sage.all import *\r\nfrom numpy.core.fromnumeric import var<\/pre>\n
def symMatrix( numRows , numCols , letter ):\r\n A = matrix(SR,numRows,numCols)\r\n for i in range(0,numRows):\r\n for j in range(0,numCols):\r\n A[i,j] = var('%s%d%d' % (letter,i,j) )\r\nreturn A<\/pre>\nsage: A = symMatrix(3,3,'a')\r\nsage: A\r\n[a00 a01 a02]\r\n[a10 a11 a12]\r\n[a20 a21 a22]<\/pre>\n
x,y,z = var('x','y','z')\r\n\r\nX = symMatrix( 3, 3 , 'X')\r\nY = symMatrix( 3, 3 , 'Y')\r\nZ = symMatrix( 3, 3 , 'Z')\r\nW = symMatrix( 3, 3 , 'W')\r\n\r\nE = x*X + y*Y + z*Z + 1*W\r\nEE = E*E.T\r\n\r\neq1=det(E)\r\neq2=EE*E-0.5*EE.trace()*E\r\n\r\neqs = (eq1,eq2[0,0],eq2[0,1],eq2[0,2],eq2[1,0],eq2[1,1],eq2[1,2],eq2[2,0],eq2[2,1],eq2[2,2])\r\n\r\nkeysA = ('x^3','y^3','x^2*y','x*y^2','x^2*z','x^2','y^2*z','y^2','x*y*z','x*y')\r\nkeysB = ('x*z^2','x*z','x','y*z^2','y*z','y','z^3','z^2','z','')\r\n\r\n# print out machine code for the linear system\r\nprintData('A',eqs,keysA)\r\nprintData('B',eqs,keysB)<\/pre>\nGenerating Code<\/h2>\n
![A=\\left[ \\begin{array}{ccc} 1 & 2 & 3 \\\\ 5 & 6 & 7 \\\\ 9 & 10 & 11 \\\\ \\end{array} \\right] \\;\\; B=\\left[ \\begin{array}{c} 4 \\\\ 8 \\\\ 12 \\\\ \\end{array} \\right]](https:\/\/s0.wp.com\/latex.php?latex=A%3D%5Cleft%5B++%5Cbegin%7Barray%7D%7Bccc%7D++1+%26+2+%26+3+%5C%5C++5+%26+6+%26+7+%5C%5C++9+%26+10+%26+11+%5C%5C++%5Cend%7Barray%7D++%5Cright%5D++%5C%3B%5C%3B++B%3D%5Cleft%5B++%5Cbegin%7Barray%7D%7Bc%7D++4+%5C%5C++8+%5C%5C++12+%5C%5C++%5Cend%7Barray%7D++%5Cright%5D&bg=f7f3ee&fg=333333&s=0 "A=\\left[ \\begin{array}{ccc} 1 & 2 & 3 \\\\ 5 & 6 & 7 \\\\ 9 & 10 & 11 \\\\ \\end{array} \\right] \\;\\; B=\\left[ \\begin{array}{c} 4 \\\\ 8 \\\\ 12 \\\\ \\end{array} \\right]")
sage: len(str(eq1)) \r\n465\r\nsage: len(str(eq1.expand()))\r\n7065\r\nsage: extractVarEq(eq1.expand(),'x^3')\r\n'X00*X11*X22 - X00*X12*X21 - X01*X10*X22 + X01*X12*X20 + X02*X10*X21 - X02*X11*X20'<\/pre>\n
def extractVarEq( expression , key ):\r\n chars = set('xyz')\r\n # expand out and convert into a string\r\n expression = str(expression.expand())\r\n # Make sure negative symbols are not stripped and split into multiplicative blocks\r\n s = expression.replace('- ','-').split(' ')\r\n # Find blocks multiplied by the key and remove the key from the string\r\n if len(key) == 0:\r\n var = [w for w in s if len(w) != 1 and not any( c in chars for c in w )]\r\n else:\r\n var = [w[:-(1+len(key))] for w in s if (w.endswith(key) and not any(c in chars for c in w[:-len(key)])) ]\r\n\r\n # Expand out power\r\n var = [expandPower(w) for w in var] \r\n\r\n # construct a string which can be compiled\r\n ret = var[0]\r\n for w in var[1:]:\r\n if w[0] == '-':\r\n ret += ' - '+w[1:]\r\n else:\r\n ret += ' + '+w\r\n\r\n return ret<\/pre>\nA.data[0] = X20*( X01*X12 - X02*X11 ) + X10*( -X01*X22 + X02*X21 ) + X00*( X11*X22 - X12*X21 );\r\nA.data[1] = Y02*( Y10*Y21 - Y11*Y20 ) + Y00*( Y11*Y22 - Y12*Y21 ) + Y01*( -Y10*Y22 + Y12*Y20 );\r\nA.data[2] = X22*( X00*Y11 - X01*Y10 - X10*Y01 + X11*Y00 ) + X20*( X01*Y12 - X02*Y11 - X11*Y02 + X12*Y01 ) + X21*( -X00*Y12 + X02*Y10 + X10*Y02 - X12*Y00 ) + X01*( -X10*Y22 + X12*Y20 ) + Y21*( -X00*X12 + X02*X10 ) + X11*( X00*Y22 - X02*Y20 );<\/pre>\n
\n