testingAlgorithm1 (trafficMatrix, phys, parameters)

function [exitMsg exitFlag netState] = testingAlgorithm1 (trafficMatrix, phys, parameters)

    % Initialization 
    exitMsg='';
    exitFlag=0;
    netState = struct('lightpathTable',[],'lightpathRoutingMatrix',[],'flowTable',[],'flowRoutingMatrix',[],'numberOfOccupiedTxs',zeros(1,phys.N),'numberOfOccupiedRxs',zeros(1,phys.N),'numberOfOccupiedTWCs',zeros(1,phys.N));
    
    % Begin
    netState.flowTable = zeros (0,7);
    serialNumber = 1;
    for ingressNode = 1:phys.N
        for egressNode = 1:phys.N
            if (ingressNode == egressNode), continue; end
            if (trafficMatrix (ingressNode , egressNode) > 0)
                netState.flowTable = [netState.flowTable ; serialNumber ingressNode egressNode trafficMatrix(ingressNode,egressNode) -1 -1 -1];            
                serialNumber = serialNumber + 1;
            end    
        end
    end
    netState.flowRoutingMatrix = zeros (size(netState.flowTable,1),0);
    netState.lightpathTable = zeros (0,3);
    netState.lightpathRoutingMatrix = zeros (0,phys.M);

    residualFlowTable = netState.flowTable; % I remove from here the carried traffic

    while (1)
        initialResidualTraffic = sum(residualFlowTable(:,4));
        initialNumberOfLps = size(netState.lightpathTable,1);

        [residualFlowTable netState] = routeTrafficInVT (residualFlowTable , netState , phys);

        [maximumResidualTrafficValue , flowId] = max(residualFlowTable(:,4));
        
        if (maximumResidualTrafficValue <= 0), break; end % everything routed, great!!

        % try a lightpath between any two pairs
        [residualFlowTable netState] = createNewLp (residualFlowTable , netState , phys);

        if (initialResidualTraffic == sum(residualFlowTable(:,4)) && (initialNumberOfLps == size(netState.lightpathTable,1)))
            break;
        end
    end
end

function [residualFlowTable netState] = routeTrafficInVT (residualFlowTable , netState , phys)

    freeCapacityOfLightpaths = ones(1,size(netState.flowRoutingMatrix,2))*phys.lightpathCapacity - sum(netState.flowRoutingMatrix,1);
    [aux orderedListOfFlowIdsPerResidualTraffic] = sort(residualFlowTable(:,4), 'descend');

    for flowId = reshape(orderedListOfFlowIdsPerResidualTraffic,1,numel(orderedListOfFlowIdsPerResidualTraffic))
        % 1. take input output pair with more residual traffic
        ingressNode = residualFlowTable(flowId,2);
        egressNode = residualFlowTable(flowId,3);
        residualTraffic = residualFlowTable(flowId,4);
        if residualTraffic<=0, break; end

        % 2. try to route the max of traffic in the residual capacity, for
        % a given ingress-egress pair
        while (1)
            if residualFlowTable(flowId,4) <= 0, break; end
            notFullLightpaths = find (freeCapacityOfLightpaths > 0);
            lpTable_notFullLps = netState.lightpathTable(notFullLightpaths,:);
            lpTable_notFullLps = [lpTable_notFullLps ones(size(lpTable_notFullLps,1),1)];
            [sequenceOfSPFLpIds,sequenceOfSPFNodeIds , totalCost] = libraryGraph_shortestPath...
                (lpTable_notFullLps(:,2:4) ,ingressNode , egressNode);
            
            if (numel(sequenceOfSPFLpIds) == 0), break; end
            
            indexesOfTraversedLps = lpTable_notFullLps(sequenceOfSPFLpIds,1);
            minimumFreeCapacityInSP = min (freeCapacityOfLightpaths(indexesOfTraversedLps));
            carriedTraffic = min (minimumFreeCapacityInSP , residualFlowTable(flowId,4));
            residualFlowTable(flowId,4) = residualFlowTable(flowId,4) - carriedTraffic;

            for lpId=reshape(indexesOfTraversedLps,1,numel(indexesOfTraversedLps))
                netState.flowRoutingMatrix (flowId,lpId) = netState.flowRoutingMatrix (flowId,lpId) + carriedTraffic;
                freeCapacityOfLightpaths(lpId) = freeCapacityOfLightpaths(lpId) - carriedTraffic;
            end
        end
    end
end



function [residualFlowTable netState] = createNewLp (residualFlowTable , netState , phys)
    
    [aux orderedListOfFlowIdsPerResidualTraffic] = sort(residualFlowTable(:,4), 'descend');

    for flowId = reshape(orderedListOfFlowIdsPerResidualTraffic,1,numel(orderedListOfFlowIdsPerResidualTraffic))
        ingressNode = residualFlowTable(flowId,2);
        egressNode = residualFlowTable(flowId,3);
        residualTraffic = residualFlowTable(flowId,4);

        if (residualTraffic <= 0), break; end
        
        if(netState.numberOfOccupiedTxs(ingressNode) >= phys.numberTxPerNode(ingressNode) ||...
                netState.numberOfOccupiedRxs(egressNode) >= phys.numberRxPerNode(egressNode))
            continue;
        end
        
        [exitFlag sequenceOfSPFLinkIds sequenceOfSPFNodeIds sequenceOfWavelengths occupiedTWCsPerNode newLightpathRoutingEntry] = ...
        libraryRWA_kShortestRWA(1 , ingressNode, egressNode, phys, netState);

        if (exitFlag == 0)
            
            availableTWCsPerNode = phys.numberTWCPerNode - netState.numberOfOccupiedTWCs;
            if sum(availableTWCsPerNode >= occupiedTWCsPerNode') == phys.N

                netState.lightpathRoutingMatrix(end+1,:) = newLightpathRoutingEntry;
                serialNumber = max(netState.lightpathTable(:,1))+1;
                if isempty(serialNumber), serialNumber = 1; end
                netState.lightpathTable(end+1,:) = [serialNumber ingressNode egressNode];
                netState.flowRoutingMatrix(:,end+1) = zeros(size(netState.flowRoutingMatrix,1),1);
                netState.numberOfOccupiedTxs(ingressNode) = netState.numberOfOccupiedTxs(ingressNode) + 1;
                netState.numberOfOccupiedRxs(egressNode) = netState.numberOfOccupiedRxs(egressNode) + 1;
                netState.numberOfOccupiedTWCs = netState.numberOfOccupiedTWCs + occupiedTWCsPerNode';

                % try to route the residual traffic in this same input output
                % pair
                carriedTraffic = min (phys.lightpathCapacity , residualFlowTable (flowId,4));
                netState.flowRoutingMatrix (flowId,end) = carriedTraffic;
                residualFlowTable (flowId,4) = residualFlowTable (flowId,4) - carriedTraffic;
                break;
            end
        end
    end
end