Generate a pore network model with grain diameters

Hello everyone,

I would like to generate a pore network model with grain diameters having log normal distribution (Grain diameter range- 20 um to 500 um; Mean= 200 um, Std dev, 0.1*Mean). The centre to centre distance for the grains need to be maintained as 550 um (to ensure there is no overlap). The domain size is 10 mm X 10 mm. I am using the following code:

clc;
clear;

% Parameters
domain_size = 10000; % Domain size in micrometers (10 mm = 10000 µm)
c_c_distance = 550; % Center-to-center distance in micrometers (µm)
N = floor(domain_size / c_c_distance); % Grid size (18 x 18 = 324 pores)
num_pores = N * N; % Total number of pores

% Define range for grain diameters
min_diameter = 20; % Minimum diameter in micrometers (µm)
max_diameter = 500; % Maximum diameter in micrometers (µm)

% Log-normal distribution parameters
mu = log(sqrt(min_diameter * max_diameter)); % Mean in log-space
sigma = sqrt(log(max_diameter / min_diameter)) / 3.29; % Standard deviation in log-space

% Generate log-normal distributed diameters
diameters = lognrnd(mu, sigma, [num_pores, 1]);

% Enforce the diameter range [20 µm, 500 µm]
diameters(diameters < min_diameter) = min_diameter;
diameters(diameters > max_diameter) = max_diameter;

% Limit the diameters to avoid overlap (should be less than c_c_distance)
max_valid_diameter = c_c_distance - 10;
diameters(diameters >= max_valid_diameter) = max_valid_diameter;

% Generate grid coordinates
x = (0:N-1) * c_c_distance; % X-coordinates
y = (0:N-1) * c_c_distance; % Y-coordinates
[X, Y] = meshgrid(x, y);

% Flatten the grid and assign diameters
X = X(:);
Y = Y(:);

% Plotting the heterogeneous pore network
figure;
hold on;
axis equal;
xlim([0, domain_size]);
ylim([0, domain_size]);
title('Heterogeneous Pore Network Model (10 mm × 10 mm) with Log-normal Grain Diameters');
xlabel('X Position (µm)');
ylabel('Y Position (µm)');

% Plot pores using circles
for i = 1:num_pores
    % Plot each pore as a circle
    rectangle('Position', [X(i) - diameters(i)/2, Y(i) - diameters(i)/2, diameters(i), diameters(i)], ...
              'Curvature', [1, 1], 'EdgeColor', 'b', 'FaceColor', [0.7, 0.8, 1], 'LineWidth', 0.5);
end

grid on;
hold off;

But seems like there are many small grain diameters than the larger ones. Also, it is seen that the grains are overlapping the axis when I change the size distribution. Can someone help me modify it or give ideas for a better one.

Any help is appreciated

Thanks in advance

Answer 1

[Arnav]

Hi @Mahima Horta,

It seems like the formula that you are using for finding μ and σ is inaccurate. The current formula seems to be finding μ as the logarithm of the geometric mean of the minimum and maximum diameters, but you want the mean to be 200 µm and the standard deviation to be 0.1 times the mean. Note that the parameters μ and σ are not the mean and standard deviation in logarithmic scale. They are simply parameters of the distribution that need to be calculated from the actual mean and standard deviation in linear scale.

The correct formula for μ and σ in terms of (mean) and (standard deviation) for a Log-normal distribution is shown below:

and

You may refer to the following Wikipedia link for the same: https://en.wikipedia.org/wiki/Log-normal_distribution

Therefore, the code for calculating mu and sigma can be changed to the following:

mean_diameter = 200;
std_dev = 0.1 * mean_diameter;
mu = log(mean_diameter^2 / sqrt(std_dev^2 + mean_diameter^2));
sigma = sqrt(log(1 + (std_dev^2 / mean_diameter^2)));

Regarding the issue of overlapping pores with the axis, you can fix this by adjusting the centres by a distance equal to their radius. Here's how you can do that:

X(i) = min(max(X(i) - diameters(i)/2, 0), domain_size - diameters(i));
Y(i) = min(max(Y(i) - diameters(i)/2, 0), domain_size - diameters(i));

The resultant plot is shown below:

Hope it helps!