开始学习使用emacs

 　　GNU Emacs是Linux里鼎鼎大名的文本编辑器，其实它的历史比Linux长得多，诞生于1975年，和微软公司都差不多年纪。其作者就是“GNU教主”--Richard Stallman。我以前好几次想好好学习这个工具的使用，可是每次都感觉入门困难，vi还好一点(其实vi的强大功能我都没用熟，有点好高鹜远，呵呵）。这两个编辑器与Windows和Linux其它的编辑器差别都很大，上手都感觉不轻松。

　　不过从网上那些很多用过Emacs，除了那些知难而退的，都对它的功能非常满意。学校也有一个老师会用Emacs，他对Emacs评价也是说非常好用，我也一直很眼馋他们。现在才开学不久，有些闲，所以趁此机会学习一下。

　　不过Emacs作为一个古老的程序，它的界面的字体配置就是一个麻烦事。我是一个非常注重字体的人，不管是什么程序，网页，乃至生活中的书，字体对于我来说都很重要的。对于要长久面对的程序，字体的漂亮当然是很重要的。目前在Linux里，作为代码和终端的字体，我唯一满意的只有Courier 10-Pitch和Consolas字体，其中前者是Magic Linux里自带的，后者是从网上下载的，其它的字体显示起来非常不好看。不过这两个字体也是只有放在Qt程序里才好看，在GTK+程序里，这两个字体总是显示得不好看。

　　我从网上查了不少资料，可能Emacs的图形界面程序是基于GTK+的原因，它显示Courier 10-Pitch和Consolas也是奇丑无比。没法，现在唯一能看就是文泉驿的字体了，文泉驿字体显示小号的汉字非常好看，但是英文相对而言就差一些，不过现在一时也找不到更好的选择了。

　　除了字体，还有一个就是配色，我本人还是喜欢那种深色调的。与终端类似的风格，我从FreeMind找到一个链接，它告诉一个Emacs配色扩展，名为color-theme，里面有很多主题，不过其中比较我的喜爱就是是arjen配色，目前就是这个配色为主。

　　FreeMind的的博主pluskid还在Google上弄了一个扩展:yasnippet，主要可以用来为各种编程语言添加代码片断，安装和使用都很方便。所以我也把它加到了我的.emacs文件里。待我学习的深入，.emacs文件也会越来越长，有时间的话，我还可以自己编写模式：scons-mode，可以用scons作为工程管理的工具，进行编译。

Emacs emacs Comments(2) 2008年9月12日 15:11

使用scons软件构建工具

　　scons是一个比较有意思的软件构建工具，其本意是想替代著名的autotools系列，我没有用过autotools，不过看到网上的autotools使用教程就头晕，后来听说了scons，就仔细下载了它的文档看了一下，觉得很符合我的口味。scons与autotools和cmake等工具不同的是，它并不是从一个类似工程配置的文件中生成Makefile，再调用make进行编译，而是直接控制编译器，从这个意义上来说，它其实并不是替代autotools，而是make！下面是一个简单的例子，如果你想编译著名的hello.c文件:

  1. #include <stdio.h>
  2.  
  3. int main(){
  4.     printf("hello,world\n");
  5.     return 0;
  6. }

　　只需在hello.c的当前目录下建立一个叫SConstruct的文件：

  1. Program('hello,'hello.c')

　　然后运行scons就可以了，scons会自动调用gcc（不是make）来编译生成可执行程序hello。如果hello.c中引用了foo.h这个头文件，而foo.h声明的函数定义在foo.c，则SConstrut可以写成：

  1. Program('hello',['hello.c','foo.c'])

　　Program的第一个参数是目标，后面可以是单个源文件名或源文件名数组。注意，scons是很智能的，它知道hello.c依赖于foo.h，所以当foo.h发生改变时，就会重新编译hello.c。这个依赖关系不需要我们额外写代码。如果源文件太多了，全写在那个数组里太长了，而源文件都是在同一个目录以.c为后辍时，我们可以简写如下：

  1. Program('hello',Glob("*.c"))

　　如果程序引入了其它的库，库名为libbar.a，库的目录为/usr/local/lib:

  1. Program('hello',Glob("*.c"),LIBS='bar',LIBPATH='/usr/local/lib')

　　如果使用了pkg-config来管理库，可以用环境(Enviroment)如下，例如我使用了GTK+ 2.0：

  1. env=Enviroment()
  2. env.ParseConfig("pkg-config --libs --cflags gtk+-2.0")
  3. env.Program('hello',Glob("*.c"))

　　scons除了可以编译C/C++程序以外，它还可以编译LaTeX，如

  1. DVI('paper.dvi','paper.tex')

　　就可以用latex命令来编译paper.dvi文件，scons还可以编译Java等其它类型的程序，具体内容见scons的官方文档。
 　　scons本身有一些设计得很好的地方，SConstruct的本身是一个合法的Python程序，它并没有像c，make那样使用自有的配置文件格式，可以使用所有Python的实用库来完成更高级的运用。还有它判断文件更改是默认情况下不是根据时间戳，而是MD5码！所以上面foo.h被touch了，并不会导致重新编译。scons的缺点是其效率太低，速度太慢。
　　scons还可以与现有的一些IDE整合，如Eclipse，KDevelop，Visual Studio等。见这里。

 

Linux scons Comments(1) 2008年9月12日 10:18

线性代数学习六，通过初等行变换对矩阵求逆

　　这是另一种矩阵求逆的方法，其基本思想是若对A矩阵求逆，则在其右边附加一个同型的单位矩阵，然后对这个合并了的矩阵进行初等行变换至行最简形，此时再分割这个矩阵，矩阵左边的应该是一个单位矩阵，而右边的则是一就是原来A矩阵的逆。下面是相关的Lua代码，与以前代码相同的部分省略。

  1. -- Merge matrices horizontally
  2. MatrixMeta.merge=function(this,other)
  3.         if this.row~=other.row or this.column ~= other.column then
  4.                 return nil
  5.         end
  6.  
  7.         local data={}
  8.         for i=1,this.row do
  9.                 for j=1,this.column do
  10.                         table.insert(data,this:at(i,j))
  11.                 end
  12.  
  13.                 for j=1,this.column do
  14.                         table.insert(data,other:at(i,j))
  15.                 end
  16.         end
  17.  
  18.         return Matrix(data,this.row,this.row*2)
  19. end
  20.  
  21. -- Divide the matrix into two parts
  22. MatrixMeta.divide=function(this)
  23.         if this.column~=2*this.row then
  24.                 return nil
  25.         end
  26.  
  27.         local a={}
  28.         local b={}
  29.         local remainder
  30.         for i,v in ipairs(this.data) do
  31.                 remainder=i%this.column
  32.                 if remainder~=0 and remainder <= this.column/2 then
  33.                         table.insert(a,v)
  34.                 else
  35.                         table.insert(b,v)
  36.                 end
  37.         end
  38.        
  39.         a=Matrix(a,this.row,this.row)
  40.         b=Matrix(b,this.row,this.row)
  41.         return a,b
  42. end
  43.  
  44. -- Another way to get inversed matrix
  45. MatrixMeta.inverse2=function(this)
  46.         if this.row ~= this.column or Determinant(this.data):eval()==0 then
  47.                 return nil
  48.         end
  49.  
  50.         local e={}
  51.         for i=1,this.row-1 do      
  52.                 table.insert(e,1)              
  53.                 for j=1,this.row do
  54.                         table.insert(e,0)
  55.                 end
  56.         end
  57.         table.insert(e,1)
  58.        
  59.         e=Matrix(e,this.row,this.column)       
  60.         local merged=this:merge(e)
  61.         merged:reduce()
  62.         local a,b=merged:divide()
  63.         return b
  64. end
  65.  
  66.  
  67.  

Algorithm inversed matrix linear algebra Comments(0) 2008年8月27日 12:18

无法重装Windows

　　Windows使用久了，速度就会越来越慢。而且程序装多了，速度也会明显下降很多，真是麻烦。于是我决定重装Windows。
把Windows C盘里的东西转到D盘后，周围的同学都无一不使用Windows，可是有XP盘的居然一个也找不出来，我以前买过很多次，后来都丢掉了，反正留着有用，还是去外面再买一张Windows XP的重装盘算了。

　　 先买了一张番茄花园的，结果出了问题，安装时选手动安装，结果上面出现一行“Setup is inspecting your computer hardware configuration”，然后屏幕就黑了，无论等多久都不亮了，难道是盘的问题，于是我又换了张，结果仍然不行，出现的问题根第一张一模一样，于是我只好跑到乐购附近的地中海数码中心买了张，也不能用，问题和先前的盘一样。

　　我不禁怀疑是不是我的电脑问题。我又在上面装了ArchLinux，这个以简洁号称的发行版不知道怎么回事，总是无法启动，GRUB的启动参数改了很多次都没用，于是后来又在上面装了一个FreeBSD，这次没问题，我以前也装过，只不过由于不能上网，只好又删了，现在由于对于网络的配置更熟了，所以解决了这个问题，可是在FreeBSD什么程序启动的程序实在是太慢了，只好又删了。

　　不过不管怎么说其它的安装是没问题的，我后来又找隔壁刘钰同学弄了张，不行，然后我从别的同学那里弄了张ghost的盘，ghost也出了问题，毛志盛后来也给我下载刻了张，都是同样的问题，真惨。我觉得屏幕黑是不是我的显示器的原因，于是借了台显示器，还是一样的问题。算了，还是等Destiny回来让他给我用Ghost恢复吧。

　　想起来Windows系统毛病还真不少，Linux无论多少装多少程序运行的速度都不会下降，而且不会因为用久了速度下降，Linux的浏览器虽然兼容性差点，有些网站有问题，但是多开页面也不会卡，Windows就不行了。

Life windows Comments(3) 2008年8月25日 09:56

线性代数学习五，向量的运算与用施密特正交化方法从线性无关组导出正交向量

这里充分使用了Lua的元表来实现Lua风格的面向对象程序设计，用Lua编程的感觉还是不错的。

  1. #!/usr/bin/env lua
  2.  
  3. VectorMeta={}
  4.  
  5. VectorMeta.__add=function(a,b)
  6.         local size=table.maxn(a)
  7.         assert(size==table.maxn(b))
  8.  
  9.         local data={}
  10.         for i=1,size do
  11.                 table.insert(data,a[i]+b[i])
  12.         end
  13.  
  14.         return Vector(data)
  15. end
  16.  
  17. VectorMeta.__sub=function(a,b)
  18.         local size=table.maxn(a)
  19.         assert(size==table.maxn(b))
  20.  
  21.         local data={}
  22.         for i=1,size do
  23.                 table.insert(data,a[i]-b[i])
  24.         end
  25.  
  26.         return Vector(data)
  27. end
  28.  
  29. VectorMeta.__mul=function(a,b)
  30.         if type(a)=="number" then
  31.                 local data={}
  32.                 for _,v in ipairs(b) do
  33.                         table.insert(data,a*v)
  34.                 end
  35.                 return Vector(data)
  36.         else
  37.                 local size=table.maxn(a)
  38.                 assert(size==table.maxn(b))
  39.                
  40.                 local sum=0
  41.                 for i=1,size do
  42.                         sum=sum+a[i]*b[i]
  43.                 end
  44.  
  45.                 return sum
  46.         end
  47. end
  48.  
  49. VectorMeta.__index=function(data,index)
  50.         if type(index)=="number" then
  51.                 return data[index]
  52.         else
  53.                 return VectorMeta[index]
  54.         end
  55. end
  56.  
  57. VectorMeta.length=function(data)
  58.         local sum=0
  59.         for _,v in ipairs(data) do
  60.                 sum=sum+v*v
  61.         end
  62.         return math.sqrt(sum)
  63. end
  64.  
  65. VectorMeta.print=function(this)
  66.         io.write('[')
  67.         for i=1,table.maxn(this)-1 do
  68.                 io.write(this[i],',')
  69.         end
  70.         io.write(this[table.maxn(this)],"]\n")
  71. end
  72.  
  73. -- Use Schmitt method to orthogonalize vectors
  74. VectorMeta.orthogonalize=function(...)
  75.         local alphas=arg
  76.         local betas={arg[1]}
  77.         for i=2,table.maxn(alphas) do
  78.                 local beta=alphas[i]
  79.                 for j=1,i-1 do
  80.                         beta=beta-((alphas[i]*betas[j])/(betas[j]*betas[j]))*betas[j]
  81.                 end
  82.                 assert(type(beta)=="table")
  83.                 table.insert(betas,beta)
  84.         end
  85.         return betas
  86. end
  87.  
  88. Vector=function(data)
  89.         setmetatable(data,VectorMeta)
  90.         return data
  91. end
  92.  

Algorithm Schmitt Orthogonalization Comments(0) 2008年8月25日 09:17

无福消受的KDE4

　　KDE4我以前也不是没有用过，主要是通过Kubuntu这个发行版来感觉。不过感觉都不好，字体等等各方面都不如意。今天一时心血来潮，就从MagicLinux的网站上通过apt-get下载了一下。
KDE4的RPM包是nihui编译的，其它的基于KDE4的发行版有很多配置上的问题，不过对于MagicLinux，我就放心多了。KDE4还依赖于Qt4，这是必然的，不过我的机器上已经自行编译好了Qt4.4，可apt-get依然给我下载了Qt4的编译包。
KDE4+Qt4比我预想的要小，只有200多M，我装好后，退回KDM，进入KDE4，启动画面真不错，黑色的基调，鼠标也是水晶黑，出来的图标也有模糊淡入效果，非常漂亮，只不过速度太慢了些。进入桌面，桌面的背景也非常好看，至于Oxygen主题，也是很不错，可是本来桌面上有一些原来的桌面图标，不过字体是文泉驿，因为都是英文，所以显示起来不好看。那些图标全放在一个框里，我本想把那个框拉大一点，结果出来什么Plasma崩溃什么的，一下子那个装图标的框就没了。
接着，我又去点KDE4的开始菜单，结果一点整个屏幕就白了，然后弹出一个Plasma的崩溃提示信息，连点几次也是这样。我于是右击把它改成经典菜单，结果就可以了，大部分的软件都在，只不过有些菜单项的图标是问号，还有几个有乱码。忍了。
然后开始启动Konsole，这个程序是我进KDE必启动的程序之一。它默认启动的是KDE4的那个Konsole，字体也是文泉驿，由于终端里主要的信息都是英文，所以显示起来也非常不好看，我于是改了一下设置，想改成Courier 10-Pitch，改是改成了，可是字体根本没变，真该死！字体不好看是一点，输入光标总与实际位置要往右偏，这也是一个问题，真麻烦。
然后我又启动了我写那个Qucifer播放器，结果它启动的是Bespin风格，这是我的Qt4程序风格设置不正确时的默认选择，真麻烦，可能是两个Qt4有冲突吧。
真麻烦，这个KDE4我是无福消受了，于是我就把花了半天心思建立起来的KDE4又给卸载了

Linux KDE4 Comments(3) 2008年8月24日 04:01

线性代数学习四，矩阵的运算和求逆

矩阵的运算是通过Lua的metatable来完成的，用起来和C++的操作符重载一样方便，而矩阵的求逆就是通过求伴随矩阵再除以其行列式的值而得。即

A^{-1}=\frac{1}{|A|}A^{*}

在实现过程中遇到了一些麻烦，有一个是我的对于伴随矩阵的理解有问题，原来伴随矩阵不仅仅要每个元素由其代数余子式代替，还得转置一下。

 

  1. #!/usr/bin/env lua
  2.  
  3. require 'determinant'
  4.  
  5. function Matrix(data,row,column)
  6.         local size=table.maxn(data)
  7.         assert(size==row*column)
  8.  
  9.         local matrix={}
  10.         matrix.data=data
  11.         matrix.row=row
  12.         matrix.column=column
  13.  
  14.         matrix.at=function(this,row,column)
  15.                 return this.data[(row-1)*this.column+column]
  16.         end
  17.  
  18.         matrix.set=function(this,row,column,num)
  19.                 this.data[(row-1)*this.column+column]=num
  20.         end
  21.  
  22.         matrix.swaprow=function(this,r1,r2)
  23.                 local temp
  24.                 for i=1,this.column do
  25.                         temp=this:at(r1,i)
  26.                         this:set(r1,i,this:at(r2,i))
  27.                         this:set(r2,i,temp)
  28.                 end
  29.         end
  30.  
  31.         matrix.simplifycolumn=function(this,column,layer)
  32.                 if this:at(layer,column) == 0 then
  33.                         for i=layer+1,this.row do
  34.                                 if this:at(i,column) ~= 0 then
  35.                                         this:swaprow(i,layer)
  36.                                         break
  37.                                 end
  38.                         end
  39.                 end
  40.  
  41.                 if this:at(layer,column) == 0 then
  42.                         return layer
  43.                 end
  44.  
  45.                 for i=1,this.row do
  46.                         if i~= layer then
  47.                                 local tomultiply=-this:at(i,column)/this:at(layer,column)
  48.                                 for j=1,this.column do
  49.                                         this:set(i,j,this:at(layer,j)*tomultiply+this:at(i,j))
  50.                                 end
  51.                         end
  52.                 end
  53.  
  54.                 return layer+1
  55.         end
  56.  
  57.         matrix.swapcolumn=function(this,c1,c2)
  58.                 local temp
  59.                 for i=1,this.row do
  60.                         temp=this:at(i,c1)
  61.                         this:set(i,c1,this:at(i,c2))
  62.                         this:set(i,c2,temp)
  63.                 end
  64.         end
  65.  
  66.         matrix.simplifyrow=function(this,row,layer)
  67.                 if this:at(row,layer) == 0 then
  68.                         for i=layer+1,this.column do
  69.                                 if this:at(row,i) ~= 0 then
  70.                                         this:swapcolumn(i,layer)
  71.                                         break
  72.                                 end
  73.                         end
  74.                 end
  75.                
  76.                 if this:at(row,layer) == 0 then
  77.                         return layer
  78.                 end
  79.  
  80.                 for i=1,this.column do
  81.                         if i ~= layer then
  82.                                 local tomultiply=-this:at(row,i)/this:at(row,layer)
  83.                                 for j=1,this.row do
  84.                                         this:set(j,i,this:at(j,layer)*tomultiply+this:at(j,i))
  85.                                 end
  86.                         end
  87.                 end
  88.  
  89.                 return layer+1
  90.         end
  91.  
  92.         -- Make the matrix to reduced form
  93.         matrix.reduce=function(this)
  94.                 local layer=1
  95.                 for i=1,this.column do
  96.                         layer=this:simplifycolumn(i,layer)
  97.                         if layer>this.row then
  98.                                 break
  99.                         end
  100.                 end
  101.                
  102.                
  103.                 for i=1,this.row do
  104.                         local tomultiply=nil
  105.                         for j=1,this.column do
  106.                                 if tomultiply then
  107.                                         this:set(i,j,this:at(i,j)*tomultiply)
  108.                                 elseif this:at(i,j)~= 0 then
  109.                                         tomultiply=1/this:at(i,j)
  110.                                         this:set(i,j,1)
  111.                                 end
  112.                         end
  113.                 end     
  114.         end
  115.  
  116.         -- Make the matrix to canonical form
  117.         matrix.canonicalize=function(this)
  118.                 this:reduce()
  119.                 local layer=1
  120.                 for i=1,this.row do
  121.                         layer=this:simplifyrow(i,layer)
  122.                         if layer>this.column then
  123.                                 break
  124.                         end
  125.                 end     
  126.         end
  127.  
  128.         -- Get the rank of the matrix
  129.         -- Matrix should be reduced before.
  130.         matrix.rank=function(this)
  131.                 local rank=0
  132.                 for i=1,this.row do
  133.                         for j=1,this.column do
  134.                                 if this:at(i,j) ~= 0 then
  135.                                         rank=rank+1
  136.                                         break
  137.                                 end
  138.                         end
  139.                 end
  140.  
  141.                 return rank
  142.         end
  143.  
  144.         matrix.print=function(this)
  145.                 for i=1,this.row do
  146.                         for j=1,this.column do
  147.                                 io.write(this:at(i,j),'\t')
  148.                         end
  149.                         io.write('\n')
  150.                 end
  151.         end
  152.  
  153.         -- Get the transposed matrix
  154.         matrix.transpose=function(this)
  155.                 local data={}
  156.  
  157.                 for i=1,this.column do
  158.                         for j=1,this.row do
  159.                                 table.insert(data,this:at(j,i))
  160.                         end
  161.                 end
  162.  
  163.                 return Matrix(data,this.column,this.row)
  164.         end
  165.  
  166.         matrix.cofactor=function(this,row,column)
  167.                 local data={}
  168.  
  169.                 for i=1,this.row do
  170.                         for j=1,this.column do
  171.                                 if i~=row and j~=column then
  172.                                         table.insert(data,this:at(i,j))
  173.                                 end
  174.                         end
  175.                 end
  176.  
  177.                 if (row+column)%2==0 then
  178.                         return Determinant(data):eval()
  179.                 else
  180.                         return -Determinant(data):eval()
  181.                 end
  182.         end
  183.  
  184.         -- Get the adjoint matrix of this matrix
  185.         matrix.adjoint=function(this)
  186.                 local data={}
  187.                 for i=1,this.row do
  188.                         for j=1,this.column do
  189.                                 table.insert(data,this:cofactor(i,j))
  190.                         end
  191.                 end
  192.  
  193.                 return Matrix(data,this.row,this.column):transpose()
  194.         end
  195.  
  196.         -- Get the inverse matrix of this matrix
  197.         -- If this matrix is not invertible, nil is returned
  198.         matrix.inverse=function(this)
  199.                 -- Only array can be invertible
  200.                 if this.row ~= this.column then
  201.                         return nil
  202.                 end
  203.  
  204.                 local det=Determinant(this.data):eval();
  205.                 if det==0 then
  206.                         return nil
  207.                 end
  208.                
  209.                 return (1/det)*this:adjoint()
  210.         end
  211.  
  212.         if MatrixCalc then
  213.                 setmetatable(matrix,MatrixCalc)
  214.         end
  215.  
  216.         return matrix
  217. end
  218.  
  219. -- The matrix calculation metatable
  220. MatrixCalc={}
  221.  
  222. MatrixCalc.__add=function(a,b)
  223.         assert(a.row==b.row and a.column==b.column)
  224.         local data={}
  225.  
  226.         for i=1,a.row*a.column do
  227.                 table.insert(data,a.data[i]+b.data[i])
  228.         end
  229.  
  230.         return Matrix(data,a.row,a.column)
  231. end
  232.  
  233. MatrixCalc.__sub=function(a,b)
  234.         assert(a.row==b.row and a.column==b.column)
  235.         local data={}
  236.  
  237.         for i=1,a.row*a.column do
  238.                 table.insert(data,a.data[i]-b.data[i])
  239.         end
  240.  
  241.         return Matrix(data,a.row,a.column)
  242. end
  243.  
  244. MatrixCalc.__unm=function(a)
  245.         local data={}
  246.  
  247.         for i=1,a.row*a.column do
  248.                 table.insert(data,-a.data[i])
  249.         end
  250.  
  251.         return Matrix(data,a.row,a.column)
  252. end
  253.  
  254. MatrixCalc.__mul=function(a,b)
  255.         if tonumber(a) then
  256.                 -- a is a number
  257.                 local data={}
  258.                 for i=1,b.row*b.column do
  259.                         table.insert(data,a*b.data[i])
  260.                 end
  261.                 return Matrix(data,b.row,b.column)
  262.         else
  263.                 -- a is a matrix
  264.                 assert(a.column==b.row)
  265.                
  266.                 local data={}
  267.                 for i=1,a.row do
  268.                         for j=1,b.column do
  269.                                 local value=0
  270.                                 for m=1,a.column do
  271.                                         value=value+a:at(i,m)*b:at(m,j)
  272.                                 end
  273.                                 table.insert(data,value)
  274.                         end
  275.                 end
  276.  
  277.                 return Matrix(data,a.row,b.column)
  278.         end
  279. end
  280.  
  281. MatrixCalc.__eq=function(a,b)
  282.         if a.row~=b.row or a.row~=b.column then
  283.                 return false
  284.         end
  285.  
  286.         for i=1,a.row*a.column do
  287.                 if a.data[i]~=b.data[i] then
  288.                         return false
  289.                 end
  290.         end
  291.  
  292.         return true
  293. end

Algorithm Matrix operation inverse Comments(0) 2008年8月23日 07:47

线性代数学习三，将矩阵变为行最简形，标准形以及计算矩阵的秩

  1. #!/usr/bin/env lua
  2.  
  3. function Matrix(data,row,column)
  4.         local size=table.maxn(data)
  5.         assert(size==row*column)
  6.  
  7.         local matrix={}
  8.         matrix.data=data
  9.         matrix.row=row
  10.         matrix.column=column
  11.  
  12.         matrix.at=function(this,row,column)
  13.                 return this.data[(row-1)*this.column+column]
  14.         end
  15.  
  16.         matrix.set=function(this,row,column,num)
  17.                 this.data[(row-1)*this.column+column]=num
  18.         end
  19.  
  20.         matrix.swaprow=function(this,r1,r2)
  21.                 local temp
  22.                 for i=1,this.column do
  23.                         temp=this:at(r1,i)
  24.                         this:set(r1,i,this:at(r2,i))
  25.                         this:set(r2,i,temp)
  26.                 end
  27.         end
  28.  
  29.         matrix.simplifycolumn=function(this,column,layer)
  30.                 if this:at(layer,column) == 0 then
  31.                         for i=layer+1,this.row do
  32.                                 if this:at(i,column) ~= 0 then
  33.                                         this:swaprow(i,layer)
  34.                                         break
  35.                                 end
  36.                         end
  37.                 end
  38.  
  39.                 if this:at(layer,column) == 0 then
  40.                         return layer
  41.                 end
  42.  
  43.                 for i=1,this.row do
  44.                         if i~= layer then
  45.                                 local tomultiply=-this:at(i,column)/this:at(layer,column)
  46.                                 for j=1,this.column do
  47.                                         this:set(i,j,this:at(layer,j)*tomultiply+this:at(i,j))
  48.                                 end
  49.                         end
  50.                 end
  51.  
  52.                 return layer+1
  53.         end
  54.  
  55.         matrix.swapcolumn=function(this,c1,c2)
  56.                 local temp
  57.                 for i=1,this.row do
  58.                         temp=this:at(i,c1)
  59.                         this:set(i,c1,this:at(i,c2))
  60.                         this:set(i,c2,temp)
  61.                 end
  62.         end
  63.  
  64.         matrix.simplifyrow=function(this,row,layer)
  65.                 if this:at(row,layer) == 0 then
  66.                         for i=layer+1,this.column do
  67.                                 if this:at(row,i) ~= 0 then
  68.                                         this:swapcolumn(i,layer)
  69.                                         break
  70.                                 end
  71.                         end
  72.                 end
  73.                
  74.                 if this:at(row,layer) == 0 then
  75.                         return layer
  76.                 end
  77.  
  78.                 for i=1,this.column do
  79.                         if i ~= layer then
  80.                                 local tomultiply=-this:at(row,i)/this:at(row,layer)
  81.                                 for j=1,this.row do
  82.                                         this:set(j,i,this:at(j,layer)*tomultiply+this:at(j,i))
  83.                                 end
  84.                         end
  85.                 end
  86.  
  87.                 return layer+1
  88.         end
  89.  
  90.         -- Make the matrix to reduced form
  91.         matrix.reduce=function(this)
  92.                 local layer=1
  93.                 for i=1,this.column do
  94.                         layer=this:simplifycolumn(i,layer)
  95.                         if layer>this.row then
  96.                                 break
  97.                         end
  98.                 end
  99.                
  100.                
  101.                 for i=1,this.row do
  102.                         local tomultiply=nil
  103.                         for j=1,this.column do
  104.                                 if tomultiply then
  105.                                         this:set(i,j,this:at(i,j)*tomultiply)
  106.                                 elseif this:at(i,j)~= 0 then
  107.                                         tomultiply=1/this:at(i,j)
  108.                                         this:set(i,j,1)
  109.                                 end
  110.                         end
  111.                 end     
  112.         end
  113.  
  114.         -- Make the matrix to canonical form
  115.         matrix.canonicalize=function(this)
  116.                 this:reduce()
  117.                 local layer=1
  118.                 for i=1,this.row do
  119.                         layer=this:simplifyrow(i,layer)
  120.                         if layer>this.column then
  121.                                 break
  122.                         end
  123.                 end     
  124.         end
  125.  
  126.         -- Get the rank of the matrix
  127.         -- Matrix should be reduced before.
  128.         matrix.rank=function(this)
  129.                 local rank=0
  130.                 for i=1,this.row do
  131.                         for j=1,this.column do
  132.                                 if this:at(i,j) ~= 0 then
  133.                                         rank=rank+1
  134.                                         break
  135.                                 end
  136.                         end
  137.                 end
  138.  
  139.                 return rank
  140.         end
  141.  
  142.         matrix.print=function(this)
  143.                 for i=1,this.row do
  144.                         for j=1,this.column do
  145.                                 io.write(this:at(i,j),'\t')
  146.                         end
  147.                         io.write('\n')
  148.                 end
  149.         end
  150.  
  151.         return matrix
  152. end

Algorithm matrix reduced form canonical form rank Comments(0) 2008年8月23日 05:00

线性代数学习二，用克拉默法则求解线性方程组

　　这是书上１.４的内容，进度得加快

  1. #!/usr/bin/env lua
  2.  
  3. require 'determinant'
  4.  
  5. function Cramer(n)
  6.         assert(n>0)
  7.  
  8.         local cramer={}
  9.  
  10.         -- the coefficient determinant
  11.         cramer.D={}
  12.  
  13.         -- this D1,D2,D3...
  14.         cramer.Ds={}
  15.         for i=1,n do
  16.                 table.insert(cramer.Ds,{})
  17.         end
  18.  
  19.         -- this number of unknown numbers
  20.         cramer.n=n
  21.  
  22.         -- add an equation
  23.         -- for example:
  24.         -- x1+2x2-x3+3x4=2
  25.         -- use cramer:add(1,2,-1,3,2)
  26.         cramer.add=function(this,...)
  27.                 local toadd={...}
  28.                 assert(table.maxn(toadd)==this.n+1)
  29.  
  30.                 for i=1,this.n do
  31.                         table.insert(this.D,toadd[i])
  32.                         for j=1,this.n do
  33.                                 table.insert(this.Ds[j],toadd[i])
  34.                         end
  35.                 end
  36.  
  37.                 local last=toadd[this.n+1];
  38.                 local times=table.maxn(this.D)/this.n
  39.                 for i=1,this.n do
  40.                         this.Ds[i][(times-1)*this.n+i]=last
  41.                 end
  42.         end
  43.  
  44.         -- Use Cramer Rules to solve linear equations
  45.         -- When solution does not exist, nil is returned,
  46.         -- and solution table returned otherwise.
  47.         cramer.solve=function(this)
  48.                 assert(table.maxn(this.D)==this.n*this.n)
  49.  
  50.                 local D=Determinant(this.D):eval();
  51.                 if D==0 then
  52.                         return nil
  53.                 end
  54.  
  55.                 local solutions={}
  56.                 local Dn;
  57.                 for i=1,this.n do
  58.                         Dn=Determinant(this.Ds[i]):eval();
  59.                         table.insert(solutions,Dn/D)
  60.                 end
  61.  
  62.                 return solutions;
  63.         end
  64.  
  65.         return cramer
  66. end

下面是一个例子：

 

  1. cramer=Cramer(4)
  2. cramer:add(1,2,-1,3,2)
  3. cramer:add(2,-1,3,-2,7)
  4. cramer:add(0,3,-1,1,6)
  5. cramer:add(1,-1,1,4,-4)
  6.  
  7. solutions=cramer:solve()
  8. for _,v in ipairs(solutions) do
  9.         print(v)
  10. end

输出是１，３，２，-１，符合书上答案

Algorithm Cramer Rules Linear Equations Comments(0) 2008年8月22日 11:00

线性代数学习一，计算行列式

复习线性代数，写了一个计算行列式的Lua脚本：

  1. #!/usr/bin/env lua
  2.  
  3. -- This is the constructor of a determinant
  4. -- It get an array of numbers as argument
  5. -- And add methods for it
  6. function Determinant(data)
  7.         local det={}
  8.         det.row=math.sqrt(table.maxn(data))
  9.  
  10.         -- If the number is not proper, raise error
  11.         if math.floor (det.row) ~= det.row then
  12.                 error("Determinant's data number must be n^2!")
  13.         end
  14.  
  15.         det.data=data
  16.        
  17.         -- index function
  18.         det.at=function(this,row,column)
  19.                 return this.data[(row-1)*this.row+column];
  20.         end
  21.        
  22.         det.cofactor=function(this,row,column)
  23.                 local cofactor={}
  24.                 local r=1
  25.                 local c=0
  26.                 for _,v in ipairs(this.data) do
  27.                         c=c+1;
  28.                        
  29.                         if c==this.row+1 then
  30.                                 r=r+1
  31.                                 c=1
  32.                         end
  33.  
  34.                         if not (r==row or c==column) then
  35.                                 table.insert(cofactor,v)
  36.                         end
  37.                 end
  38.                
  39.                 return cofactor
  40.         end     
  41.  
  42.         -- evaluation function
  43.         det.eval=function(this)
  44.                 local a=this.data
  45.                 -- this is the simplest case
  46.                 if this.row==2 then
  47.                         return a[1]*a[4]-a[2]*a[3]
  48.                
  49.                 -- for accerleration
  50.                 elseif this.row==3 then
  51.                         return a[1]*a[5]*a[9]+a[3]*a[4]*a[8]+a[2]*a[6]*a[7]-a[3]*a[5]*a[7]-a[1]*a[6]*a[8]-a[2]*a[4]*a[9]
  52.                 else
  53.                         local value=0
  54.                         local cofactor=0
  55.                         local add=true
  56.  
  57.                         for i=1,this.row do
  58.                                 cofactor=Determinant(this:cofactor(1,i))
  59.                                 cofactor=this.data[i]*cofactor:eval()
  60.                                
  61.                                 if add then
  62.                                         value=value+cofactor
  63.                                 else
  64.                                         value=value-cofactor
  65.                                 end
  66.                                 add = not add
  67.                         end
  68.  
  69.                         return value
  70.                 end
  71.         end                         
  72.        
  73.         return det
  74. end
  75.  
  76.  

 

Algorithm Determinant Comments(0) 2008年8月22日 09:15