区域内随机圆形的生成

生成的方法比较简单，同样是使用二维填充等高线函数contourf。

1.设置范围及指定圆的数目

范围设置为10x10的正方形。而圆的数量通过其所占的面积指定，若全部的圆占总面积的\(40\%\)，设圆的半径为\(r\)，则圆的数量为\([\frac{40}{\pi r^2}]\)。

% 指定区域范围
x = 0:0.01:10;
y = 0:0.01:10;
[X, Y] = meshgrid(x, y);

% 设置圆半径及空隙分数
% 计算圆的数量
r = 0.5;
porosity = 0.6;
Acir = (1 - porosity) * 100;
cir_num = round(Acir/(r^2*pi));

2.随机选择圆心并判断是否重叠

每产生一个圆心，就遍历之前所有的圆心，计算他们之间的距离\(dis = \sqrt{x_{dis}^2 + y_{dis}^2}\)，如果\(dis < 2r+0.1\)则说明它们重合。加0.1为了保证至少有0.1的空隙。
还要保证圆的范围不能超过边界，即\(0.5<x_{coor}<9.5\ \&\&\ 0.5<y_{coor}<9.5\)。

% 生成第一个圆心
Center(1, 1) = 10 * rand(1);
Center(1, 2) = 10 * rand(1);
counter_next = 2;

% 循环生成剩余圆心
while (counter_next <= cir_num)
    x_coor = 10 * rand(1);
    y_coor = 10 * rand(1);
    if (x_coor <= 0.5 || x_coor >= 9.5 ...
            || y_coor <= 0.5 || y_coor >= 9.5)  % 防止圆越过边界
        continue;
    end
    for (i = counter_next-1:-1:1)
        x_dis = abs(x_coor - Center(i, 1));
        y_dis = abs(y_coor - Center(i, 2));
        dis = sqrt(x_dis^2+y_dis^2);
        if (dis <= 1.1)
            break;
        elseif (i == 1)
            Center(counter_next, 1) = x_coor;
            Center(counter_next, 2) = y_coor;
            counter_next = counter_next + 1;  % 一开始写在while下面出错
        end
    end
end

3. 使用min函数结合所有圆并做出等高线

完整程序如下：

clear all;
clc;
% 指定区域范围
x = 0:0.01:10;
y = 0:0.01:10;
[X, Y] = meshgrid(x, y);

% 设置圆半径及空隙分数
% 计算圆的数量
r = 0.5;
porosity = 0.6;
Acir = (1 - porosity) * 100;
cir_num = round(Acir/(r^2*pi));

% 声明圆心坐标矩阵并初始化
Center = -ones(cir_num, 2);

% 生成第一个圆心
Center(1, 1) = 10 * rand(1);
Center(1, 2) = 10 * rand(1);
counter_next = 2;

% 循环生成剩余圆心
while (counter_next <= cir_num)
    x_coor = 10 * rand(1);
    y_coor = 10 * rand(1);
    if (x_coor <= 0.5 || x_coor >= 9.5 ...
            || y_coor <= 0.5 || y_coor >= 9.5)  % 防止圆越过边界
        continue;
    end
    for (i = counter_next-1:-1:1)
        x_dis = abs(x_coor - Center(i, 1));
        y_dis = abs(y_coor - Center(i, 2));
        dis = sqrt(x_dis^2+y_dis^2);
        if (dis <= 1.1)
            break;
        elseif (i == 1)
            Center(counter_next, 1) = x_coor;
            Center(counter_next, 2) = y_coor;
            counter_next = counter_next + 1;  % 一开始写在while下面出错
        end
    end
end

% 循环使用min函数将圆形组合
Z = (X - Center(1, 1)).^2 + (Y - Center(1, 2)).^2 - 0.25;
for i = 2:1:cir_num
    Z = min(Z, (X - Center(i, 1)).^2 + (Y - Center(i, 2)).^2 - 0.25);
end

% 使用contourf做出函数值为0的等高线
ax = figure;
[M, C] = contourf(X, Y, Z, [-0.25, 0]); % -0.25是为了之后填充
C.LineWidth = 1;
C.ShowText = 'off';

map = [0, 0.6, 1;
        1, 1, 1];
colormap(ax, map);

结果如下：