First, to clarity, you are working with a defuzzification algorithm described in the paper "Defuzzification of spatial fuzzy sets by feature distance minimization". The term Ck
is initially built for k = 0
as:
Core(f) = {x is a member of a set X | m(x) >= m(y) for all y in X}
where f
is some discrete grayscale 2D image, m(x)
is defined as the membership function with range [0, 1]
, and X
is formed by the pairs of points that belong to the domain of f
. Then, Supp(f)
is formed by those points where m(x) > 0
. For a grayscale image f
where 0 is black, this means anything that is not black is in the Supp(f)
.
So, here is an input f
(which also represents Supp(f)
), Core(f)
, and Supp(f)\Core(f)
respectively. In Matlab terms, Core(f)
is given by: core = f; core(f ~= max(max(f))) = 0;

Now we can solve the problems in the question. First: "... among all pixels p that belong to Supp(F)\Ck being 4-neighbours of Ck ..." can be translated as:
allp = imdilate(core, strel('diamond', 1)) - core;
allp_in_f = (allp/255) .* f;
supposing the input is of type uint8
. The second statement discards points that are not in Supp(f)
by setting them to 0, but in this example there are no such points. The result of this is step is shown in the following image:

Now: "... select the pixel p that minimizes the distance d(F,Ck union {p}) ...". At this moment, any pixel p
in the last image minimizes it, supposing the metric is the L-infinity (Chebyshev, chessboard, ...) but even if we picked L-1 norm, all the pixels except the diagonal ones would satisfy it. Or, another interpretation (if we completely ignore the paper, since it wasn't mentioned in the question): the distance considers the intensities in f
, and the method wants the one that minimizes the weighted (intensities) distance. In this case, and also for future iterations (the method is an iterative one, described in the paper), the simplest way is to perform a flood-fill from the boundary obtained in the last image shown. This flood-fill would be done in a breadth-first manner, so the first time a pixel touches f
you would backtrack and determine the pixel p
that minimized the distance. This method is also known as wave-front propagation.