用matlab感受数学之美

最新推荐文章于 2024-05-11 07:11:02 发布

manderous

最新推荐文章于 2024-05-11 07:11:02 发布

阅读量1.3k

点赞数 2

文章标签：数学 matlab

本文链接：https://blog.csdn.net/manderous/article/details/53027704

版权

相比于之前的课程，这个学期，应用数学开了繁多的课程，运筹、数组、matlab这几门的实验课更是增添了很多课程设计、课程论文要写。再配上实变、复变、常微分，觉得这个学期简直完美（括弧微笑脸）。
近期，在写matlab小课程设计，打算从分形入手。想到n多年前，在微博上截过一张图，没想到如今派上了用场。

于是，今天我便特意搜索了matrix67的这篇文章《用三段 140 字符以内的代码生成一张 1024×1024 的图片》
这篇文章中有画Mandelbrot集（曼德勃罗集），不过是c语言代码，于是，我将其修改为matlab代码。通过matlab感受数学之美。

1.Levy 曲线

function levy(n)
% levy(16),n为levy曲线迭代次数
%x1,y1,x2,y2为初始线段两端点坐标,nn为迭代次数
n=16;
x1=0;y1=0;
x2=1;y2=0;
%第i-1次迭代时由各条线段产生的新两条线段的三端点横、纵坐标存储在数组X、Y中
[X,Y]=levy1(x1,y1,x2,y2);
for i=1:n 
    for j=1:length(X)/3
        w=levy2(X(1+3*(j-1):3*j),Y(1+3*(j-1):3*j));
        [XX(3*2*(j-1)+1:3*2*(j-1)+3),YY(3*2*(j-1)+1:3*2*(j-1)+3)]=levy1(w(1,1),w(1,2),w(1,3),w(1,4));
        [XX(3*2*(j-1)+3+1:3*2*(j-1)+3+3),YY(3*2*(j-1)+3+1:3*2*(j-1)+3+3)]=levy1(w(2,1),w(2,2),w(2,3),w(2,4));
    end
    X=XX;
    Y=YY;
end
plot(X,Y)
hold on
axis equal
end

%由以(x1,y1),(x2,y2)为端点的线段生成新的中间点坐标并把(x1,y1),(x2,y2)连同新点横、纵坐标依次分别存储在数组X,Y中
function [X,Y]=levy1(x1,y1,x2,y2)
x3=1/2*(x1+x2+y1-y2);
y3=1/2*(-x1+x2+y1+y2);
X=[x1,x3,x2];
Y=[y1,y3,y2];
end

%把由函数levy1生成的三点横、纵坐标X,Y顺次划分为两组,分别对应两条折线段中每条线段两端点的坐标,
%并依次分别存储在2*4阶矩阵w中,w中第i(i=1,2)行数字代表第i条线段两端点的坐标
function w=levy2(X,Y)
a11=X(1);b11=Y(1);
a12=X(2);b12=Y(2);
a21=X(2);b21=Y(2);
a22=X(3);b22=Y(3);
w=[a11,b11,a12,b12;a21,b21,a22,b22];
end

2.二维元胞自动机算法画Sierpinski三角形

function sierpinski_ca2(m,n)
%sierpinski_ca2(400,400)
m=400;n=400;
t=1;w=zeros(2,m*n);
s=zeros(m,n);
s(m/2,n/2)=1;
for i=[m/2:-1:2,m/2:m-1]
    for j=[n/2:-1:2,n/2:n-1]
        if mod(s(i-1,j-1)+s(i,j-1)+s(i+1,j-1)+s(i-1,j)+s(i+1,j)+s(i-1,j+1)+s(i,j+1)+s(i+1,j+1),2)==1
            s(i,j)=1;
            w(1,t)=i;
            w(2,t)=j;
            t=t+1;
        end
    end
end
plot(w(1,:),w(2,:),'.','markersize',0.1)
end

3.Mandelbrot集

// NOTE: compile with g++ filename.cpp -std=c++11

#include <iostream>
#include <cmath>
#include <cstdlib>
#define DIM 1024
#define DM1 (DIM-1)
#define _sq(x) ((x)*(x)) // square
#define _cb(x) abs((x)*(x)*(x)) // absolute value of cube
#define _cr(x) (unsigned char)(pow((x),1.0/3.0)) // cube root

unsigned char GR(int,int);
unsigned char BL(int,int);

unsigned char RD(int i,int j){
float x=0,y=0;int k;for(k=0;k++<256;){float a=x*x-y*y+(i-768.0)/512;y=2*x*y+(j-512.0)/512;x=a;if(x*x+y*y>4)break;}return log(k)*47;
}

unsigned char GR(int i,int j){
float x=0,y=0;int k;for(k=0;k++<256;){float a=x*x-y*y+(i-768.0)/512;y=2*x*y+(j-512.0)/512;x=a;if(x*x+y*y>4)break;}return log(k)*47;
}

unsigned char BL(int i,int j){
float x=0,y=0;int k;for(k=0;k++<256;){float a=x*x-y*y+(i-768.0)/512;y=2*x*y+(j-512.0)/512;x=a;if(x*x+y*y>4)break;}return 128-log(k)*23;
}

void pixel_write(int,int);
FILE *fp;
int main(){
    fp = fopen("MathPic.ppm","wb");
    fprintf(fp, "P6\n%d %d\n255\n", DIM, DIM);
    for(int j=0;j<DIM;j++)
        for(int i=0;i<DIM;i++)
            pixel_write(i,j);
    fclose(fp);
    return 0;
}
void pixel_write(int i, int j){
    static unsigned char color[3];
    color[0] = RD(i,j)&255;
    color[1] = GR(i,j)&255;
    color[2] = BL(i,j)&255;
    fwrite(color, 1, 3, fp);
}

4.Mandelbrot集

// NOTE: compile with g++ filename.cpp -std=c++11

#include <iostream>
#include <cmath>
#include <cstdlib>
#define DIM 1024
#define DM1 (DIM-1)
#define _sq(x) ((x)*(x)) // square
#define _cb(x) abs((x)*(x)*(x)) // absolute value of cube
#define _cr(x) (unsigned char)(pow((x),1.0/3.0)) // cube root

unsigned char GR(int,int);
unsigned char BL(int,int);

unsigned char RD(int i,int j){
double a=0,b=0,c,d,n=0;
while((c=a*a)+(d=b*b)<4&&n++<880)
{b=2*a*b+j*8e-9-.645411;a=c-d+i*8e-9+.356888;}
return 255*pow((n-80)/800,3.);
}

unsigned char GR(int i,int j){
double a=0,b=0,c,d,n=0;
while((c=a*a)+(d=b*b)<4&&n++<880)
{b=2*a*b+j*8e-9-.645411;a=c-d+i*8e-9+.356888;}
return 255*pow((n-80)/800,.7);
}

unsigned char BL(int i,int j){
double a=0,b=0,c,d,n=0;
while((c=a*a)+(d=b*b)<4&&n++<880)
{b=2*a*b+j*8e-9-.645411;a=c-d+i*8e-9+.356888;}
return 255*pow((n-80)/800,.5);
}

void pixel_write(int,int);
FILE *fp;
int main(){
    fp = fopen("MathPic.ppm","wb");
    fprintf(fp, "P6\n%d %d\n255\n", DIM, DIM);
    for(int j=0;j<DIM;j++)
        for(int i=0;i<DIM;i++)
            pixel_write(i,j);
    fclose(fp);
    return 0;
}
void pixel_write(int i, int j){
    static unsigned char color[3];
    color[0] = RD(i,j)&255;
    color[1] = GR(i,j)&255;
    color[2] = BL(i,j)&255;
    fwrite(color, 1, 3, fp);
}