Avoid many for-loops when calculating contact lengths between elements

I built a thermal model composed of N elements (N~10^5-10^6), grouped in a smaller number of groups (nGroups~10^2-10^3). The 1xN vector elementToGroup defines to which group each element belongs. The position of the first element of each group is zero.

The elements are disposed along one direction and have non-uniform length. The 1xN vector ds defines their length.

Some of the groups are in thermal contact. The Mx2 matrix groupPairs defines pairs of groups in thermal contact (the number M can vary).

My goal is to find:

1. elementPairs: the Qx2 matrix defining pairs of elements that are in thermal contact, i.e. that in groups that are in thermal contact and are encompassing the same length at the same distance from zero (the number Q can vary).
2. pairContactLength: the Qx1 vector defining the length for which the two elements of each pair are touching.

I have a prototype code that does what it is supposed to do on a small number of elements (N=20). However, the code is not optimized and quite inelegant. I believe it'll scale up very badly.

Do you have suggestions for improving the code and make it more efficient and readable? Some hints would be really appreciated. Thank you!

clear all

% % % Inputs
N=20;
ds = [12 6 6 8 16 20 17 2 7 5 6 20 11 7 17 3 18 2 15 15];
% ds=round(rand(1,N)*20);

sStartPositions=[0, cumsum(ds(1,1:end-1))];
sEndPositions=cumsum( ds(1,:) );
sAvePositions=(sStartPositions+sEndPositions)/2;

elementToGroup( 1: 3)=1;
elementToGroup(4: 10)=2;
elementToGroup(11: 14)=3;
elementToGroup(15: 20)=4;

groupPairs = [ 1,2; 2,3; 3,4; 1,3;];


% % % Visuailze the elements in each group
nGroups=max(elementToGroup);
sStartPositions_perElement=cell(nGroups,1);
for g =1:nGroups
    idxGroupLast=find(elementToGroup==g,1,'last');
    idxGroup=elementToGroup==g;
    if g==1
        sStartPositions_perElement{g} = [sStartPositions(idxGroup), sEndPositions(idxGroupLast)];
    elseif g==nGroups
        idxGroupPrevious=find(elementToGroup==g-1,1,'last');
        sStartPositions_perElement{g} = [sStartPositions(idxGroup), sEndPositions(N)]-sEndPositions(idxGroupPrevious);
    else
        idxGroupPrevious=find(elementToGroup==g-1,1,'last');
        sStartPositions_perElement{g} = [sStartPositions(idxGroup), sEndPositions(idxGroupLast)]-sEndPositions(idxGroupPrevious);
    end
end

% % % Plot element positions
figure
hold on, grid on, grid minor
for g =1:max(elementToGroup)
    plot(sStartPositions_perElement{g},g*ones(1,length(sStartPositions_perElement{g})),'o-')
end
ylim([0 nGroups+1])


% Inizialize with a high number of elements (matrix and vector will be truncated afterwards
elementPairs = zeros(N*2,2);
pairContactLength = zeros(N*2,1);
current_pair_of_elements = 1;
for current_pair_of_groups = 1:length(groupPairs)
    % Find elements in the current two groups to pair
    element_from = find(elementToGroup==groupPairs(current_pair_of_groups,1));
    element_to   = find(elementToGroup==groupPairs(current_pair_of_groups,2));

    % Calculate the position the elements of the current groups to pair, starting from zero
    position_from_zero_from = cumsum(ds(element_from));
    position_from_zero_to   = cumsum(ds(element_to));
          
    current_element_from = 1;
    current_element_to = 1;
    current_position = 0;
    while current_element_from<=length(element_from) && current_element_to<=length(element_to) && (position_from_zero_from(current_element_from)-ds(element_from(current_element_from)) <= position_from_zero_to(end))
        % Pair the current elements from the two groups
        elementPairs(current_pair_of_elements,:) = [element_from(current_element_from), element_to(current_element_to)];

        % Calculate the length for which the two elements are touching
        if position_from_zero_from(current_element_from) < position_from_zero_to(current_element_to)
            pairContactLength(current_pair_of_elements) = position_from_zero_from(current_element_from) - current_position;
            current_position = position_from_zero_from(current_element_from);
            current_element_from = current_element_from+1;
            current_pair_of_elements = current_pair_of_elements+1;
        elseif position_from_zero_from(current_element_from) > position_from_zero_to(current_element_to)
            pairContactLength(current_pair_of_elements) = position_from_zero_to(current_element_to) - current_position;
            current_position = position_from_zero_to(current_element_to);
            current_element_to = current_element_to+1;
            current_pair_of_elements = current_pair_of_elements+1;
        elseif position_from_zero_from(current_element_from) == position_from_zero_to(current_element_to)
            current_element_from = current_element_from+1;
        end
    end
    % now all elements of the current group have been paired
end


% % % Display output
clc
elementPairs = elementPairs(1:current_pair_of_elements-1,:)
pairContactLength = pairContactLength(1:current_pair_of_elements-1)