MATLAB | 一起来感受数学之美，第一届迷你黑客大赛回顾

Hey真的是好久不见，最近确实是比较忙更新频率也下来了，过段时间应该能恢复正常更新速度，之前给大家解说过今年举办的math is beautiful迷你黑客大赛，但这其实是第二届大赛，本期推送带大家回顾一下第一期大赛，同时也是MATHWORKS社区成立二十周年纪念比赛。

另外由于比赛代码有280个字符的限制，因此很多代码会写的晦涩难懂，以下发布的代码将是本人改编且加注释的版本：

第一届比赛大家的创造性很明显没有第二届足，但是依旧萌生了很多非常有创意的作品：

pale blue dot

原作者：Adam Danz
链接：
https://ww2.mathworks.cn/matlabcentral/communitycontests/contests/4/entries/1116

改写加注释代码：

rng(0)
% 创建黑色背景colormap为粉色的坐标区域
axes('Colormap',pink,'View',-[15,24],'CLim',[-.5,8],'Color','k');               
hold on

% 绘制圆形星球
[x,y,z]=sphere(98);
surf(x,y,z,'EdgeColor','none') 

% 绘制一个圆盘
% 让x,y,z随便某组数据缺失(此处x缺失)
% 使得圆盘变成一个个有空隙的圆环
% 使用hypot获取(x,y)坐标到中心点距离映射为颜色
x(randi(99,7),:)=NaN;
surf(x*2,y*2,z*0,hypot(x,y)*6,'EdgeColor','none')
axis equal

% 生成一些白色(`w`)随机散点和一个青色(`c`)点
n=[1,2E3];
scatter3(rand(n)*4,ones(n)+1,rand(n)*7-4,randg(1,n),'w','filled')  
scatter3(2,2,0,20,'c','filled')
camva(2)

Hi 😃

链接：
https://ww2.mathworks.cn/matlabcentral/communitycontests/contests/4/entries/536

[x y]=meshgrid(linspace(-2.5,2.5,70),linspace(-5.5,5.5,70));
z=exp(-x.^2-0.5*y.^2).*cos(4*x) + exp(-3*((x+0.5).^2+0.5*y.^2));
idx=(abs(z)>0.0005 );
z(idx)=0.001*sign(z(idx));
patch(surf2patch(surf(x,y,z)),'FaceColor','interp');
view(35,65);
colormap(cool);
grid,axis off
camlight headlight

这里只给出链接是因为这个函数并不是链接中的人所发现的，事实上这个神奇的函数：

$y)=e^{-x^2-\frac{y^2}{2}} \cos (4 x)+e^{-3\left((x+0.5)^2+\frac{y^2}{2}\right)}$

是工程师Mike Croucher于2007年在Walking Randomly上使用Mathematica发布的：

http://walkingrandomly.com/?p=19

2010年他本人也给出了俩MATLAB版本代码：

http://walkingrandomly.com/?p=2720

[x y] = meshgrid( linspace(-3,3,50), linspace(-5,5,50) );
z = exp(-x.^2-0.5*y.^2).*cos(4*x) + exp(-3*((x+0.5).^2+0.5*y.^2));
idx = ( abs(z)>0.001 );
z(idx) = 0.001 * sign(z(idx)); 

figure('renderer','opengl')
patch(surf2patch(surf(x,y,z)), 'FaceColor','interp');
set(gca, 'Box','on', ...
    'XColor',[.3 .3 .3], 'YColor',[.3 .3 .3], 'ZColor',[.3 .3 .3], 'FontSize',8)
title('$e^{-x^2 - \frac{y^2}{2}}\cos(4x) + e^{-3((x+0.5)^2+\frac{y^2}{2})}$', ...
    'Interpreter','latex', 'FontSize',15,'Color','b') 

view(35,65)
colormap( [flipud(cool);cool] )
camlight headlight, lighting phong

zL = 0.001;
x = linspace(-3,3,50);
y = linspace(-5,5,50);
[X,Y]=meshgrid(x,y);
Z = exp(-X.^2-Y.^2/2).*cos(4*X) + exp(-3*((X+0.5).^2+Y.^2/2));
Z(Z>0.001) = zL;
Z(Z<-0.001) = -zL;
surf(X,Y,Z);
title({'Secret Messages Hidden Inside Equations',...
       '$e^{-x^2 - \frac{y^2}{2}}\cos(4x) + e^{-3((x+0.5)^2+\frac{y^2}{2})}$'},...
       'Interpreter','latex', 'FontSize',15,'Color','b') 
colormap(flipud(cool))
view([1 -1.5 2])

Cumulus

原作者：Jenny Bosten
链接：
https://ww2.mathworks.cn/matlabcentral/communitycontests/contests/4/entries/4211

改写加注释代码：

rng(333)
X=-164:.646:164;
[t,r]=cart2pol(X,X'-90);
s=(285-r)/285;
% 二维傅里叶逆变换生成噪声
image(cat(3,s*.3+.12,s*.5+.21,s*.36+.57),'AlphaData',5E4*abs(ifft2(r.^-1.7.*cos(6*rand(508)))));
B=getframe(gcf,[78 52 426 250]);

% 将每一行像素随机左右平移生成水波
m=B.cdata;
d=flip(m);
for k=1:250
    m=[m;circshift(d(k,:,:),round(15*randn),2)];
end
image(m)
axis off

首先就是通过二维傅里叶逆变换生成噪声：

之后将图形下半部分每一行都随机左右移动模拟水波:

A Time-Lapse of the Night Sky

原作者：Adam Danz
链接：
https://ww2.mathworks.cn/matlabcentral/communitycontests/contests/4/entries/4841

改写加注释代码：

% 使用q/q的方式生成单位列向量与q相乘
% 从而生成从上到下数值逐渐减小的矩阵用作背景
q=50:-.5:1;
hold
imagesc(q'.*(q./q),[20,80])

% 固定每一个星轨角度变化
% 随机生成一些星轨初始角度和初始位置半径
% 并转化为x,y坐标
t=0:2.5E-5:.1;
n=4001;
[x,y]=pol2cart(6*rand(1,n)+t',rand(1,n).*(t*0+142));
m=99;
plot(x+m,y+m,'w')
% 拼凑出房子和草地
x=[1,1:.1:m,m]; 
y=[1,cos(0:.1/98:1)*30+randg(1,1,981),1]; 
fill(x,y,'k')
fill([4 4 3 6 9 8 8]*5,[5 7 7 8 7 7 5]*5,'k')
axis([1 m 1 m])
camva(4) 
colormap gray

Above the clouds

原作者：Tim
链接：
https://ww2.mathworks.cn/matlabcentral/communitycontests/contests/4/entries/8273

改写加注释代码：

a=200;
b=(.5:a)'/a;
c=(-cos(b*2*pi)+1).^.2;
d=ones(a);
f=b-.5;
r=f'.^2+f.^2;
% 生成不对称柏林噪声
m=50;
surf(b,b',abs(ifftn(exp(7i*rand(a))./r.^.9)).*(c*c')*30)
l=(m:-1:1)/m;
% 绘制云彩
hold on
for n = 1:m
    surf(b,b',d*n,d+cat(3,1,1,1),'EdgeAlpha',0,'FaceAlpha',max(.2,l(n)));
end
zlim([-a/2,a])
shading flat
colormap(flip([b,b,b]))
camva(5)
axis off

感觉是本次比赛最酷的一个了，非常像是水墨画，同时此代码还有很多不同数值和配色的Remix Tree：

Ghost Pentagon Flower

原作者：Daniel Pereira
链接：
https://ww2.mathworks.cn/matlabcentral/communitycontests/contests/4/entries/124

代码：

theta = linspace(0,360,6);
c=bone(20);
for i=20:-1:1
    patch(i*sind(theta+18*i),i*cosd(theta+18*i),c(i,:),'edgecolor','none');
end
axis equal off; set(gcf,'color','w');

White Hole

原作者：Michal Halon
链接：
https://ww2.mathworks.cn/matlabcentral/communitycontests/contests/4/entries/5828

代码：

n=2000;
x=rand(4,n);
i=1;
for j=1:n
    x(4,i)=0.5-sqrt((x(1,i)-0.5)^2+(x(2,i)-0.5)^2);
    if x(4,i)>0.5
        x(:,i)=[];
    else
        if x(4,i)<0.08
            x(4,i)=0;
        end
        i=i+1;
    end
end
colormap bone
bubblechart(x(1,:),x(2,:),x(3,:),x(4,:),'MarkerFaceAlpha',0.35,'MarkerEdgeAlpha',0);
set(gcf,'Color','k')
axis square
axis off

代码想法贼简单，就是生成一些点，如果在圆外面就删掉，然后在每个点处生成随机大小半透明圆点，结束。
不过代码属实是写麻烦了，以下给出我的改写版：

n=2000;
x=rand(4,n);
x(4,:)=0.5-vecnorm(x(1:2,:)-.5);
x(:,x(4,:)>.5)=[];
x(4,x(4,:)<0.08)=0;
colormap bone
bubblechart(x(1,:),x(2,:),x(3,:),x(4,:),'MarkerFaceAlpha',0.35,'MarkerEdgeAlpha',0);
set(gcf,'Color','k')
axis square
axis off

Universe in a Nutshell

原作者：Sebastian Kraemer on
链接：
https://ww2.mathworks.cn/matlabcentral/communitycontests/contests/4/entries/5036

这真的是本次比赛注释最全，代码最美观的选手了，注释比较全就不翻译了，仅仅对代码做了微调：

% set background
figure('Color','k','InvertHardcopy','off')
set(gca,'Color','k')
hold
% function shortcuts
m=@(a,c)max(min(a,c),-c);
s=@vecnorm;
n=@(p)p./s(p);
% random uniformly distributed values on the unit sphere
P=n(randn(3,1e3));
% iterate gravitational force exchange
A=copper(320);
for k=1:320
    Q=P;
    % distance of all points to each other
    D=reshape(P,3,1,[])-P;
    % calculate gravitational exchange force
    U=sum(m(D./s(D).^2,10),3)/1e4-1e-3;
    % add bounded, cumulative displacement and project result to unit sphere
    P=n(P+n(U).*m(s(U),.01));
    % plot imitating more points than actually calculated
    h=scatter(P(1,:),P(2,:),m(s(P-Q).^-3/1e8,80),A(k,:),'filled');
    alpha(h,.3)
end
axis off equal

Unknown pleasures - album cover

原作者： Alberto Cuadra Lara
链接：
https://ww2.mathworks.cn/matlabcentral/communitycontests/contests/4/entries/6323

这张仿图真的太太经典了，此图曾经被作为Joy Division乐队的专辑封面

% A. Cuadra - 20 Oct 2021
figure('Color','k');
hold on;
set(gca,'Color','k','XColor','none','YColor','none')
x=linspace(-9,13);
y=4e3;% y0 - initial position
for i=1:80% 80 number of lines
s=randi([1,2]);
n=.04*rand(1,1e2)+.4/s*exp(-(x-randi([1,4])).^2/2/s^2);%noise vector + normal pdf
patch(x,y+n*9e2,'k')% fill the gap in black
plot(x,y+n*9e2,'w')% plot line
y=y-5e1;% move y axis
end

The Lighthouse

原作者：Jr
链接：
https://ww2.mathworks.cn/matlabcentral/communitycontests/contests/4/entries/3626

此代码直接用了一些特殊符号绘制灯塔，这点需要比较新的版本才能做到：

c=1200;
n=800;
hold;
% 生成一系列正弦函数填充图形
for i=18:-1:1
    y=i/(8+3)*n;
    u=c*(1-(i-1)/20);
    x=linspace(0,c,n);
    p=y/n*12;
    t=x*6/u;
    q=sin(p+t)+sin(p+.3*t);
    r=y+q*160*.8^(i-1);
    v=[x' r';c 0;0 0];
    fill(v(:,1),v(:,2),i/18*[.1 .7 1],'EdgeColor','n');
end
% 特殊符号绘制灯塔
text(9,440,'♜','Fontsi',53);
text(111,438,'◍','color','y');
axis off;
ylim([0,n])