Calculating the greater circle distance from source point to centre of all grid cells in a domain/reference grid

I am trying to calculate the greater circle distance in metres from a source point (specified as a latitude/longitude) to the centre of all grid cells within a domain (or reference grid - calculated using makerefmat).

Would anyone be able to help me with this?

% Source point location (decimal degrees)
sourcelon = -7.340612; 
sourcelat = 57.005755;

% Reference grid
% Domain size 
DDX = 4800; % Varation in the domain size for X 
DDY = 5760; % Variation in the domain size for Y 

% Define limits of latitude and longitude for the area of interest
latlim = [53.5, 59.5]; 
lonlim = [-9, -4];

% Create vectors of latitude and longitude which are divided by the same 
% interval sizes as the domain
latvec = linspace(latlim(1),latlim(2),DDY); 
lonvec = linspace(lonlim(1),lonlim(2),DDX); 

% To define dy and dx, divide the initial latitude minus the final latitude
% (or longitude) by the number of elements in the domain array i.e.DDX or DDY 
dy = abs(latlim(1)-latlim(2))/DDY; %dy = difference in latitude (y) between pixels in successive rows
% (degrees) 
dx = abs(lonlim(1)-lonlim(2))/DDX; % dx = difference in longitude (x) between pixels in successive columns
% (degrees) 

% Specify the coordinates for the left hand corner of the map domain 
x11 = lonlim(1);
y11 = latlim(1);

% Create the referenced matrix grid (R) 
R = makerefmat(x11, y11, dx, dy);