function [NatNeighborBoundaryPairs, Vertices, VoronoiCells, D_Inv, Directions, MidPoints, NatNeighborPairs, Global_Area] = ...
    UpdateDiscretizationMetrics(BoundaryCandidatesIdx, PosX, PosY, InsideOldMask, NumNodesX, NumNodesY, InteriorCandidatesIdx,MeshSize, InteriorCandidatesMask)


% update the irregular boundaries voronoi mesh 
[NatNeighborBoundaryPairs, Vertices, VoronoiCells] =...
    GetVoronoiDecomposition(BoundaryCandidatesIdx, PosX, PosY, InsideOldMask, [], NumNodesX, NumNodesY);

% update boundary CVs mutual inverse distance, midpoints and directions
[D_Inv_Boundary, Directions, MidPoints] = GetBoundaryIndependentMetrics(NatNeighborBoundaryPairs, PosX, PosY);

% update the mesh neighboring info
NatNeighborPairs = [InteriorCandidatesIdx-1 InteriorCandidatesIdx
                            InteriorCandidatesIdx+1 InteriorCandidatesIdx
                            InteriorCandidatesIdx-NumNodesX InteriorCandidatesIdx
                            InteriorCandidatesIdx+NumNodesX InteriorCandidatesIdx
                            NatNeighborBoundaryPairs];
                        
% update inverse distance for neighboring CVs
D_Inv = [1/MeshSize*ones(4*length(InteriorCandidatesIdx),1); D_Inv_Boundary];                        
                

% create new area array and fill it with regular CVs areas
Global_Area = zeros(NumNodesX*NumNodesY,1);                        
Global_Area(InteriorCandidatesMask) = ones(length(InteriorCandidatesIdx),1)*MeshSize*MeshSize;

end
