%% Melinda 254 Well Parameters

dt = 0.090;
tubingHeadPressure = 250.0;
fluidGradient = 0.450;
stuffingBoxFriction = 100.0;
pumpDiameter = 1.500;
tubingInnerDiameter = 2.441;
tubingOuterDiameter = 2.875;
tubingAnchorDepth = 9209.0;
structuralRating = 320000;

rodLengths = [4375.0, 3600.0, 400.0, 150.0];
rodDiameters = [1.225, 1.000, 0.875, 1.500];
rodMaterials = ["fiberglass", "steel", "steel", "steel"];
rodDampingFactors = [0.08, 0.08, 0.08, 0.08];

%
%% Run the Well Setup Function from algorithm.m

melinda252Well = wellSetup(dt, rodLengths, rodDiameters, ...
    rodMaterials, rodDampingFactors, tubingHeadPressure, ...
    stuffingBoxFriction, fluidGradient, pumpDiameter, tubingOuterDiameter, ...
    tubingInnerDiameter, tubingAnchorDepth, structuralRating);

% Alternatively, use the Well Class to define the well.

melinda252 = Well(dt, rodLengths, rodDiameters, ...
    rodMaterials, rodDampingFactors, tubingHeadPressure, ...
    stuffingBoxFriction, fluidGradient, pumpDiameter, tubingOuterDiameter, ...
    tubingInnerDiameter, tubingAnchorDepth, structuralRating);

%
%% Import the Barney Well data into a CSV and setup some persistent variables.

melinda252InputData = csvread("welldata/melinda252.csv");
melinda252InputData = [melinda252InputData; melinda252InputData];

melinda252Data = zeros(length(melinda252InputData), 4);
count = zeros(1, melinda252Well.numTapers);
topPosArray = zeros(melinda252Well.numTapers+1, 100);
topLoadArray = zeros(melinda252Well.numTapers+1, 100);

%
%% Loop through the Barney Well data and calculate downhole data

outputIndex = 1;
for ptI = 2:length(melinda252InputData)
   [pumpPosition, pumpLoad, status, topPosArray, topLoadArray, ...
    count] = pocAlgorithm(melinda252InputData(ptI,1), melinda252InputData(ptI-1,1), ...
    melinda252InputData(ptI,2), count, melinda252Well.dt, melinda252Well.a, melinda252Well.rodDampingFactors, ...
    melinda252Well.factorArray, melinda252Well.rodLengths, double(melinda252Well.lagIndex), ...
    melinda252Well.rodYMs * 10.0^6, melinda252Well.area, melinda252Well.lengthRequired, ...
    melinda252Well.centerPoint, melinda252Well.rodWeightFluidTotal, ...
    melinda252Well.stuffingBoxFriction, melinda252Well.force, topPosArray, topLoadArray);

    melinda252.evalSim(melinda252InputData(ptI,1), melinda252InputData(ptI,2));
    
    if status == 1
        melinda252Data(outputIndex, 1) = melinda252InputData(ptI,1);
        melinda252Data(outputIndex, 2) = melinda252InputData(ptI,2);
        melinda252Data(outputIndex, 3) = pumpPosition;
        melinda252Data(outputIndex, 4) = pumpLoad;
        outputIndex = outputIndex + 1;
    end
    
end



%
%% Plot the surface and downhole cards.

figure;
ax1 = subplot(2,1,1);
plot(melinda252Data(:,1), melinda252Data(:,2));
title('Surface Card');

ax2 = subplot(2,1,2);
plot(melinda252Data(:,3), melinda252Data(:,4));
title('Downhole Card');
linkaxes([ax1, ax2], 'x');
melinda252.plotCards();

%

