How does an incoherent beam propagate in free space?

I've written a code that takes a circle function, adds some decoherence and propagates it through space:

L = 500;               %side length 
M = 10000;                %number of samples 
delta = L/M;             %src sample interval 
u = -L/2:delta:L/2-delta;  %src coords 
v = (L/2-delta:-delta:-L/2)';

lambda = 0.65;     %wavelength 
srcw = 50;              %source half width in x (um) 
z = 2750;               %propagation dist and focal length (um) 

fu = -1/(2*delta):1/L:1/(2*delta)-(1/L);
fv = fu;
[FU,FV] = meshgrid(fu,fv);
Hincoh = circ(sqrt(FU.^2+FV.^2)/f0);
OTF = ifft2(abs(fft2(fftshift(Hincoh))).^2);
OTF = abs(OTF/OTF(1,1));

u1 = circ(sqrt(u.^2+v.^2)/srcw);
I1 = abs(u1.^2);        %src intensity
 
figure; 
imagesc(u,u,I1); 
axis square; axis xy; 
colormap('gray'); xlabel('x (um)'); ylabel('y (um)'); 
title('z= 0 m');

% Add incoherence
Ig = I1/max(max(I1));
Gg = fft2(fftshift(Ig));
Gi = Gg.*OTF;
Ii = ifftshift(ifft2(Gi));
% Remove residual imag parts, values < 0
Ii = real(Ii); mask = Ii>=0; Ii = mask.*Ii;
ui = sqrt(Ii);
  
clear I1

% Propagate to lens

zprop1 = 2750;

uprop = Propagate2D(ui,zprop1,lambda,Fv,Fv);
Iprop = abs(uprop.^2);

figure; 
imagesc(u,u,Iprop); 
axis square; axis xy; 
colormap('gray'); xlabel('x (um)'); ylabel('y (um)'); 
title('After propagation');

figure; % intensity slice
plot(u,Iprop(M/2+1,:)');
xlabel('x(um)');ylabel('Intensity')

The circ function is

function[out]=circ(r)
% 
% circle function 
% 
% evaluates circ(r) 
% note: returns odd number of samples for diameter 
% 
out=abs(r)<=1; 
end 

and the Propagate2D function is

function u = Propagate2D(E0,Z,Lambda,Fx,Fy)

[FxGrid, FyGrid] = meshgrid(Fx,Fy);

H = exp(((1i*2*pi*Z)/Lambda)*sqrt(1-(Lambda*FxGrid).^2 - (Lambda*FyGrid).^2));
H = fftshift(H)';

clear FxGrid;
clear FyGrid;

E1 = fft2(fftshift(E0));
clear E0;

E2 = H.*E1;

clear H
clear E1;

u = ifftshift(ifft2(E2));
I = abs(u.^2);

end

The beam I'm actually trying to model is divergent. I'm not sure whether or not the propagation I've shown is correct as it doesn't look how I'd expect. Does anyone have any insight into whether or not this is correct and whether it is possible to include divergence using this method?