牛!Mathematica还能这样用!自己制作马赛克拼图

    Wolfram的Blog上更新了一段非常牛的Mathematica代码,真的让我大开眼界。只需要三行代码,你就可以自己做一个马赛克拼图。
imagePool = Map[With[{i = Import[#]}, {i, Mean[Flatten[N[i[[1, 1]]], 1]]}] &, FileNames["Pool/*.jpg"]];
closeMatch[c_] := RandomChoice[Take[SortBy[imagePool, Norm[c - #[[2]]] &], 20]][[1]];
Grid[Reverse[Map[closeMatch, Import["MasterImage.tif"][[1, 1]], {2}]], Spacings -> {0, 0}]

    其中,”Pool/*.jpg”是你的图库,我估计最少也得有几百张吧。我用Photoshop把我的collection全部处理成了35×35的小图;为了让最终效果更佳,我特地把它们全部处理成单色的,并且减小了对比度。”MasterImage.tif”是你的目标图片,Mathematica会把这个图片中的每一个像素用图库中一个合适的图来代替。我把我的照片剪裁了一下,然后压成19×22的大小。Mathematica首先把所有照片以及每个照片的RGB值的中位数存成一个list,函数closeMatch将图片按照RGB值的均方根排序,然后随机从头20个中选出一个。第三行用Grid函数输出我们所要的马赛克拼图。最后我们就得到了——由众MM图所组成的Matrix67的肖像画!!如果你还看不出来的话,站远点儿眯着眼睛就能看出来了。

Read more…

Geek饰物DIY:粘土工艺之Sierpinski三角形

    一提到水果,人们首先想到的往往是苹果;一提到AV女优,最先想到的总是武腾兰;同样地,一提到分形图形,大多数人都会首先想起Sierpinski三角形。Sierpinski三角形可能是最具有代表性的分形图形了,随身佩戴一个Sierpinski三角形绝对够酷。回想Sierpinski三角形的构造方法,将三个同样的三角形的边长缩小一半,再与一个空白的倒三角相拼即可得到一个更高阶的Sierpinski三角形。这种构造方法非常简单,它是在现实生活中最容易构造的分形图形之一,你所需要的仅仅是一种可以拉伸变形的材料。

  
1. 准备好两种颜色的软陶泥(比如蓝色和白色);
2. 捏出四个三角形的长条,三个蓝色的,一个白色的;
3. 把这四个长条拼成一个大三角形

  
4. 把这个长条拉长到原来的四倍(因此横截面积缩小到原来的1/4)

  
5. 切下三段一样长的长条,再捏一个同样大小的白色三角形长条
6. 重复步骤3到5

  
    第7次迭代后,我们得到了2187个三角形,很多细节已经看不清了,此时你可以把它近似地看作一个Sierpinski三角形。接下来你要做的,就是用金属把它串起来,最后烘烤成形即可。钥匙链、手机链、耳环、项链……想拿它干啥就干啥吧。

查看更多:http://www.evilmadscientist.com/article.php/fimofractals

Linux下的数学工具Maxima 简明教程(下)

三角运算
(%i1) trigexpand(sin(10*x+y));
(%o1)                 cos(10 x) sin(y) + sin(10 x) cos(y)
(%i2) trigexpand(sin(2*x));
(%o2)                           2 cos(x) sin(x)
(%i3) trigsimp(2*cos(x)^2+sin(x)^2);
                                     2
(%o3)                             cos (x) + 1
(%i4) trigreduce(-sin(x)^2+3*cos(x)^2+x);
                      cos(2 x)      cos(2 x)   1        1
(%o4)                 -------- + 3 (-------- + -) + x - -
                         2             2       2        2

代数推理
(%i1) assume(x>0,y<-1,z>=0);
(%o1)                      [x > 0, y < - 1, z >= 0]
(%i2) assume(a<b and b<c);
(%o2)                           [b > a, c > b]
(%i3) facts();
(%o3)               [x > 0, - 1 > y, z >= 0, b > a, c > b]
(%i4) is(a>c);
(%o4)                                false
(%i5) is(z-y>0);
(%o5)                                true
(%i6) is(z-x>0);

Maxima was unable to evaluate the predicate:
z - x > 0
-- an error.  Quitting.  To debug this try debugmode(true);
(%i7) prederror:false;
(%o7)                                false
(%i8) is(z-x>0);
(%o8)                               unknown
(%i9) forget(a<b);
(%o9)                               [b > a]
(%i10) is(a>c);
(%o10)                              unknown

级数计算
(%i1) sum(i,i,1,5);
(%o1)                                 15
(%i2) sum(i^2,i,1,5);
(%o2)                                 55
(%i3) sum(1/2^i,i,1,inf);
                                   inf
                                   ====
                                        1
(%o3)                               >    --
                                   /      i
                                   ====  2
                                   i = 1
(%i4) sum(1/2^i,i,1,inf),simpsum;
(%o4)                                  1
(%i5) sum(1/i^2,i,1,inf),simpsum;
                                        2
                                     %pi
(%o5)                                ----
                                      6
(%i6) sum(1/i,i,1,inf),simpsum;
(%o6)                                 inf

微积分
(%i1) limit(1/x,x,inf);
(%o1)                                  0
(%i2) limit(sin(x)/x,x,0);
(%o2)                                  1
(%i3) limit(sin(x),x,inf);
(%o3)                            &n
bsp;    ind
(%i4) diff(3*x^2+x+5/x,x);
                                       5
(%o4)                            6 x - -- + 1
                                        2
                                       x
(%i5) diff(sin(x)*tan(x),x);
                                           2
(%o5)                   cos(x) tan(x) + sec (x) sin(x)
(%i6) diff(%e^(a*x),x);
                                        a x
(%o6)                               a %e
(%i7) integrate(sin(x)^3,x);
                                  3
                               cos (x)
(%o7)                          ------- - cos(x)
                                  3
(%i8) integrate(x^3,x,1,3);
(%o8)                                 20
(%i9) taylor(%e^x,x,0,3);
                                     2    3
                                    x    x
(%o9)/T/                    1 + x + -- + -- + . . .
                                    2    6
(%i10) taylor(sin(x),x,0,5);
                                  3    5
                                 x    x
(%o10)/T/                    x - -- + --- + . . .
                                 6    120
(%i11) taylor(sqrt(x+1),x,1,3);
                                                     2                  3
                    sqrt(2) (x - 1)   sqrt(2) (x - 1)    sqrt(2) (x - 1)
(%o11)/T/ sqrt(2) + --------------- - ---------------- + ----------------
                           4                 32                128
                                                                        + . . .
(%i12) ratsimp(%);
                      3              2
             sqrt(2) x  - 7 sqrt(2) x  + 43 sqrt(2) x + 91 sqrt(2)
(%o12)       -----------------------------------------------------
                                      128

矩阵运算
(%i1) f[i,j]:=i+j;
(%o1)                           f     := i + j
                                 i, j
(%i2) genmatrix(f,3,3);
                                  [ 2  3  4 ]
                                  [         ]
(%o2)                       &nbs
p;     [ 3  4  5 ]
                                  [         ]
                                  [ 4  5  6 ]
(%i3) g[i,j]:=i-2^j;
                                              j
(%o3)                           g     := i - 2
                                 i, j
(%i4) genmatrix(g,3,3);
                               [ - 1  - 3  - 7 ]
                               [               ]
(%o4)                          [  0   - 2  - 6 ]
                               [               ]
                               [  1   - 1  - 5 ]
(%i5) %o2+%o4;
                                 [ 1  0  - 3 ]
                                 [           ]
(%o5)                            [ 3  2  - 1 ]
                                 [           ]
                                 [ 5  4   1  ]
(%i6) %o2.%o4;
                               [ 2  - 16  - 52 ]
                               [               ]
(%o6)                          [ 2  - 22  - 70 ]
                               [               ]
                               [ 2  - 28  - 88 ]
(%i7) %o2^^3;
                               [ 360  474  588 ]
                               [               ]
(%o7)                          [ 474  624  774 ]
                               [               ]
                               [ 588  774  960 ]
(%i8) x:matrix([17, 3],[-8, 11]);
                                  [ 17   3  ]
(%o8)                             [         ]
                                  [ - 8  11 ]
(%i9) x^^-1;
                                [ 11      3  ]
                                [ ---  - --- ]
                                [ 211    211 ]
(%o9)                           [            ]
                                [  8    17   ]
            &n
bsp;                   [ ---   ---  ]
                                [ 211   211  ]

想了解更多请阅读官方文档:
http://maxima.sourceforge.net/docs/manual/en/maxima.html

做人要厚道
转贴请注明出处

Linux下的数学工具Maxima 简明教程(上)

    这个Blog里曾经多次提到过超强数学软件Mathematica,但目前为止我还没发现它的Linux版,Wine似乎也没有用。其实,在Linux下也有很多类似于Mathematica的数学软件,其中Maxima是我用的最多的一个。这里简单介绍一下Maxima的各个函数供大家参考,也方便我自己今后查询。

安装:sudo apt-get install maxima maxima-share
运行:maxima
退出:quit();

基本运算
(%i1) 2+3;
(%o1)                                  5
(%i2) 5*6;
(%o2)                                 30
(%i3) %+2;
(%o3)                                 32
(%i4) %o1*%o3;
(%o4)                                 160
(%i5) 4/7+3/4;
                                      37
(%o5)                                 --
                                      28
(%i6) float(%);
(%o6)                          1.321428571428571
(%i7) 2^32;
(%o7)                             4294967296
(%i8) 30!;
(%o8)                  265252859812191058636308480000000
(%i9) float(sqrt(2));
(%o9)                          1.414213562373095

三角函数和对数函数
(%i1) float(sin(1));
(%o1)                           0.8414709848079
(%i2) sin(%pi/2);
(%o2)                                  1
(%i3) sin(%pi/2)+cos(%pi/3);
                                       3
(%o3)                                  -
                                       2
(%i4) float(sec(%pi/3)+csc(%pi/3));
(%o4)                          3.154700538379252
(%i5) log(1);
(%o5)                                  0
(%i6) float(log(10));
(%o6)                          2.302585092994046
(%i7) log(%e);
(%o7)                                  1
(%i8) log(2^a);
(%o8)                              log(2) a
(%i9) %e^log(2);
(%o9)                                 2

变量操作
(%i1) a^2-b^2;
                                     2    2
(%o1)                               a  - b
(%i2) a:3;
(%o2)                                  3
(%i3) a^2-b^2;
                                         2
(%o3)                               9 - b
(%i4) b:2;
(%o4)                                  2
(%i5) a^2-b^2;
(%o5)                                  5
(%i6) kill(a);
(%o6)                                done
(%i7) kill(b);
(%o7)                                done
(%i8) a^2-b^2;
                                     2    2
(%o8)                      &nb
sp;        a  - b

函数操作
(%i1) f(x):=x^2-1;
                                         2
(%o1)                           f(x) := x  - 1
(%i2) f(2);
(%o2)                                  3
(%i3) f(100);
(%o3)                                9999
(%i4) float(f(2/3));
(%o4)                         - 0.55555555555556
(%i5) a:4/5;
                                       4
(%o5)                                  -
                                       5
(%i6) f(a);
                                       9
(%o6)                                - --
                                       25

多项式运算(展开、合并、化简和消元)
(%i1) expand((a+b)^3);
                            3        2      2      3
(%o1)                      b  + 3 a b  + 3 a  b + a
(%i2) factor(a^2-b^2);
(%o2)                          - (b - a) (b + a)
(%i3) ratsimp((x^2-1)/(x+1));
(%o3)                                x - 1
(%i4) eliminate([x^2+x*y+z=0,3*x+5*y+z=0,x-y-2*z^2=1],[y,z]);
                             4      3       2
(%o4)               [- x (8 x  - 2 x  + 19 x  - 50 x + 25)]

解方程
(%i1) solve(x^2-3*x+4/x=5,x);
                         sqrt(5) + 1      sqrt(5) - 1
(%o1)             [x = - -----------, x = -----------, x = 4]
                              2                2
(%i2) funcsolve(f(n)*(n+1)+2*n=1-f(n)/n,f(n));
                                      n (2 n - 1)
(%o2)                        f(n) = - -----------
                                       2
                                      n  + n + 1
(%i3) solve([x+3*y=10,1/x+x*y=4],[x,y]);
                              sqrt(69) - 9      4 sqrt(3) sqrt(23) - 34
(%o3) [[x = 1, y = 3], [x = - ------------, y = -----------------------],
                                   2            9 sqrt(3) sqrt(23) - 75
                                    sqrt(69) + 9      4 sqrt(3) sqrt(23) + 34
                               [x = ------------, y = -----------------------]]
                                         2            9 sqrt(3) sqrt(23) + 75
(%i4) solve(x^2+b*x+c=0,x);
                           2                       2
                     sqrt(b  - 4 c) + b      sqrt(b  - 4 c) - b
(%o4)         [x = - ------------------, x = ------------------]
                      &nbs
p;      2                       2
(%i5) find_root(x^x=2,x,1,2);
(%o5)                          1.559610469462369
(%i6) find_root(sin(x)=x/2,x,0.1,%pi);
(%o6)                          1.895494267033981

数论相关
(%i1) mod(100,7);
(%o1)                                  2
(%i2) primep(3214567);
(%o2)                                true
(%i3) next_prime(200);
(%o3)                                 211
(%i4) factor(1001);
(%o4)                               7 11 13
(%i5) factor(30!);
                        26  14  7  4   2   2
(%o5)                  2   3   5  7  11  13  17 19 23 29
(%i6) gcd(200,780);
(%o6)                                 20
(%i7) binomial(7,4);
(%o7)                                 35
(%i8) fib(7);
(%o8)                                 13

画函数图像
(%i1) plot2d(x^3+2*x^2-3,[x,-2,2]);
*** X11 output driver not found, switching to dumb terminal!
*** If you want to use the X11 output, please install the gnuplot-x11 package

  14 ++-------+--------+--------+--------+-------+--------+--------+-------++
     +        +        +        +        +       +       x^3+2*x^2-3 $$$$$$ $
  12 ++                                                                    $+
     |                                                                    $ |
  10 ++                                                                  $ ++
     |                                                                  $   |
     |                                                                  $   |
   8 ++                                                                $   ++
     |                                                                $     |
   6 ++                                                             $$     ++
     |                                                             $$       |
   4 ++                                                          $$        ++
     |                                                          $$          |
   2 ++   
                                                     $$          ++
     |                                                      $$$             |
     |                                                     $$               |
   0 ++                                                 $$$                ++
     |                                               $$$$                   |
  -2 ++$$$$$$$$$$$$$$$$$$$$$$$$$$               $$$$$                      ++
     $$       +        +        $$$$$$$$$$$$$$$$ +        +        +        +
  -4 ++-------+--------+--------+--------+-------+--------+--------+-------++
    -2      -1.5      -1      -0.5       0      0.5       1       1.5       2

(%o1)

你可以通过安装gnuplot-x11让maxima在X上画图,安装方法是:
sudo apt-get install gnuplot-x11
maxima也可以画3D图像,比如执行下面代码可以画出sin(x)cos(y)的图像,我就不贴图了,大家自己试试。
plot3d(sin(x)*cos(y),[x,-2,2],[y,-2,2]);

做人要厚道
转贴请注明出处

C语言速成手册(六):其它问题、后记

预处理指令
    以一个井号开头的行都叫做预处理指令。除了#include指令外,我们还经常用到#define指令。#define指令可以告诉编译器,编译时把代码中出现的特定标识当作什么来处理。例如,我们可以这样写:
#define NAME_OF_MY_POTENTIAL_GF "ZPR"
    这样,编译器会在编译前把代码中出现NAME_OF_MY_POTENTIAL_GF的地方全部替换成"ZPR"。这种替换是无条件的,但是有一个例外:当指定的标识属于某个字符串(被引号引起来)时替换不会发生。例如,下面两行代码会输出NAME_OF_MY_POTENTIAL_GF defined as: ZPR
printf("NAME_OF_MY_POTENTIAL_GF defined as: ");
printf(NAME_OF_MY_POTENTIAL_GF);

    其中,后面那个NAME_OF_MY_POTENTIAL_GF被自动替换为"ZPR"。如果哪一天ZPR不要我了,我就可以非常方便地让整个程序适用于另一个MM。

    C语言中通常会用#define代替const。例如,下面的代码假设了输入数据n<=2000。
#include <stdio.h>
#define MAX 2000

int main()
{
   int f[MAX][MAX];
   int i,j,n;
   scanf("%d",&n);
   for ( i=0; i<n; i++ )
   {
      for ( j=0; j<=i; j++ )
      {
          f[i][j] = j ? (f[i-1][j] + f[i-1][j-1]) % 10000 : 1;
          printf( "%5d" , f[i][j] );
      }
      printf("n");
   }

   return 0;
}

    下面的这些指令也是合法的:
#define begin {
#define end }
#define and &&
#define or ||

    #define定义的指令允许带参数。例如,下面的定义也是合法的:
#define sqr(x) x*x
    观察下面的这个程序:
#include <stdio.h>
#define begin {
#define end }
#define writeln(num) printf("%dn",num)
#define sqr(x) x*x

int main()
begin
   writeln(sqr(100));
   writeln(sqr(10+2));
end

    程序输出:
10000
32

    为什么第二个输出的数是32不是144?不要忘了sqr中的x不是一个变量,编译器仅仅是把x替换为10+2,因此sqr(10+2)的结果是10+2*10+2,当然是32咯。为了避免这种情况,这样写就没问题了:
#define sqr(x) ( (x) * (x) )
    下面这个定义很常用:
#define MAX(a,b) ( ((a) > (b)) ? (a) : (b) )

    如果你想写出一个很有个性的C代码,反复使用#define是一个不错的选择。例如,这段代码就极具个性,一个光棍的形象跃然于屏幕上。然而,真正把#define发挥得淋漓尽致的,还是要数这段代码

static声明
    在函数中的变量声明前加一个static可以使这个变量具有“记忆性”。观察下面的程序:
#include <stdio.h>
void printNum()
{
   int a=1;
   static int b=1;
   printf("%d %dn", a++, b++);
}
int main()
{
   int i;
   for ( i=1; i<=5; i++ )
       printNum();
   return 0;
}

    程序输出:
1 1
1 2
1 3
1 4
1 5

short类型和int类型的范围
    最初我们列出的C语言类型和Pascal类型的对比只能提供一个参考。事实上不同的编译器中short和int的范围可能不同。你可以查一下前面说过的limits.h来确定这些类型的实际范围。通常short是16位整数,long是32位整数。在Windows下Dev-C++中int类型是32位。

对64位整型的处理
    和Free Pascal一样,对64位整数类型的处理总是比较麻烦。
    首先,对long long赋值很可能会发生错误,你可以在常数后添加一个LL表明这是long long类型。其次,C语言中有些函数是要区分数据类型的,你需要根据数据类型选用恰当的函数。最后,long long类型的输出很可能也有问题,此时你可以用"%lld"来替换"%d",表明输出的是一个long long类型。在Windows下总要装点怪,我在Windows编译时非要用"%I64d"才行。
    下面的程序代码在Windows下Dev-C++中一切正常。
#include <stdio.h>
#include <stdlib.h>

int main()
{
    long long a;
    a = -5841314520LL;
    a = llabs(a);
    a = a + 1;
    printf("%I64d",a);
    return 0;
}

查漏补缺
    这个系列到这里就结束了。还有我没有说到的语法点吗?欢迎大家补充。

后记
    C语言速成手册到这里就结束了。这很可能是网上现有的原创C语言教材中讲解最快,篇幅最短的,因为它只适合已经学过其它语言,了解程序设计基础知识的人。这一系列的文章略过了大量的概念讲解、示例代码和习题,你可以自己在网上阅读一些C语言程序作为补充。以后我可能还会写一些类似的文章介绍其它语言。下一步我计划写C与C++的区别,对象和类的介绍以及C++的新特性。再以后我可能会向Java或者Ruby的方向发展。
    祝各位努力转C的OIer暑假愉快。

Matrix67原创
转贴请注明出处