text:
#using Ipopt: optimize!

using Base: Int64
using XLSX: length
using Pkg
#=
Pkg.add("Tables")
Pkg.add("TableView")
Pkg.add("DataFrames")
Pkg.add("XLSX")
Pkg.add("ExcelReaders")
Pkg.add("JuMP")
Pkg.add("Ipopt")
Pkg.add("NamedArrays")
Pkg.add("DelimitedFiles")
=#
using XLSX, ExcelReaders, DataFrames, Tables, JuMP, Ipopt, NamedArrays, DelimitedFiles;

#filepath cross system compatability code
if Sys.KERNEL === :NT || Sys.KERNEL === :Windows
    pathmark = "\\";
else
    pathmark = "/";
end

#functions for ease of display

function tableprint(c)
    withTitles = [["Titles" permutedims(names(c)[2])];[names(c)[1] c]]
end

function tablecsv(c::Any, d::AbstractString)
    withTitles = [[d permutedims(names(c)[2])];[names(c)[1] c]]
    writedlm("tableCSV"*pathmark*d*".csv", withTitles, ',');
end

IOSource = ExcelReaders.readxlsheet("."*pathmark*"data"*pathmark*"IOData"*pathmark*"5209055001DO001_201819.xls", "Table 5");

#indexing vectors for initial data import groups
intermediaryTotalsCol = findall(x -> occursin("T4", x), string.(IOSource[3,:]));
intermediaryTotalsRow = findall(x -> occursin("T1", x), string.(IOSource[:,1]));
finalTotalsCol = findall(x -> occursin("T6", x), string.(IOSource[3,:]));
finalTotalsRow = findall(x -> occursin("Australian Production", x), string.(IOSource[:,2]));
finalDemandCol = findall(x -> occursin('Q', x), string.(IOSource[3,:]));
factorRow = findall(x -> occursin('P', x), string.(IOSource[:,1]));
IOSourceCol = vcat(intermediaryTotalsCol, finalDemandCol, finalTotalsCol);
IOSourceRow = vcat(intermediaryTotalsRow, factorRow[1:length(factorRow)-1], finalTotalsRow);
#initialising IO
IO = zeros(length(IOSourceRow), length(IOSourceCol));
#import numerical data into IO
IO[1:length(IOSourceRow), 1:length(IOSourceCol)] = IOSource[IOSourceRow, IOSourceCol];
#creating vectors of titles for IO
IOCodeRow = string.(IOSource[IOSourceRow, 1]);
IOCodeRow[length(IOSourceRow)] = "T2";
IOCodeCol = string.(IOSource[3, IOSourceCol]);
IONameRow = string.(IOSource[IOSourceRow, 2]);
IONameCol = string.(IOSource[2, IOSourceCol]);
IO = NamedArray(IO);
setnames!(IO, string.(IOCodeRow), 1);
setnames!(IO, string.(IOCodeCol), 2);
#code to sum public and private entities into one column
investment = findall(x -> occursin("Capital Formation", x), IONameCol);
IO[:, investment[1]]=sum(eachcol(IO[:, investment[1:2]]));
IO = IO[:,Not(investment[2])];
#alter title vectors accordingly (include Q in total investment column in IOcode)
IOCodeCol[investment[1]] = "Q3+Q4";
IOCodeCol = string.(IOCodeCol[Not(investment[2])]);
IONameCol[investment[1]] = "Private and Public Gross Fixed Capital Formation";
IONameCol = IONameCol[Not(investment[2])];
#creating a dictionary for the index of each column and row in IO by IOCode
IOColDict = Dict(IOCodeCol .=> [1:1:8;]);
IORowDict = Dict(IOCodeRow .=> [1:1:8;]);
IOCapForm = findall(x -> occursin("Capital Formation", x), IONameCol);
IOChangeInv = findall(x -> occursin("Changes in Inventories", x), IONameCol);

#importing relevant ASNA data for table 5
ASNAHouseCap = ExcelReaders.readxlsheet("."*pathmark*"data"*pathmark*"ASNAData"*pathmark*"5204039_Household_Capital_Account.xls", "Data1");
ASNANonFinCap = ExcelReaders.readxlsheet("."*pathmark*"data"*pathmark*"ASNAData"*pathmark*"5204018_NonFin_Corp_Capital_Account.xls", "Data1");
ASNAFinCap = ExcelReaders.readxlsheet("."*pathmark*"data"*pathmark*"ASNAData"*pathmark*"5204026_Fin_Corp_Capital_Account.xls", "Data1");
ASNAGovCap = ExcelReaders.readxlsheet("."*pathmark*"data"*pathmark*"ASNAData"*pathmark*"5204032_GenGov_Capital_Account.xls", "Data1");
ASNAYearRow = findall(x -> occursin("2019", x), string.(ASNAHouseCap[:,1]));

#table 5
#creating table 5a - allocation of investment expenditure (broken into subsections for dict referencing purposes)
#subsection a is fixed capital expenditure
table5aNameCol = ["Households", "Non-Financial Corporations", "Financial Corporations", "General Government", "Total"];
table5aNameRow = ["Domestic Commodities", "Imported Commodities, complementary", "Imported Commodities, competing",
"Taxes less subsidies on products", "Other taxes less subsidies on investment", "Total indirect taxes",
"Total fixed capital expenditure"];
table5aRowDict = Dict(table5aNameRow .=> [1:1:length(table5aNameRow);]);
table5aColDict = Dict(table5aNameCol .=> [1:1:length(table5aNameCol);]);
table5a = zeros(length(table5aNameRow), length(table5aNameCol));
table5a = NamedArray(table5a);
setnames!(table5a, string.(table5aNameRow), 1);
setnames!(table5a, string.(table5aNameCol), 2);

#filling in totals column from corresponding IO data
table5a[table5aRowDict["Domestic Commodities"], table5aColDict["Total"]] = sum(IO[IORowDict["T1"],IOCapForm]);
table5a[table5aRowDict["Imported Commodities, complementary"], table5aColDict["Total"]] = sum(IO[IORowDict["P5"],IOCapForm]);
table5a[table5aRowDict["Imported Commodities, competing"], table5aColDict["Total"]] = sum(IO[IORowDict["P6"],IOCapForm]);
table5a[table5aRowDict["Taxes less subsidies on products"], table5aColDict["Total"]] = sum(IO[IORowDict["P3"],IOCapForm]);
table5a[table5aRowDict["Other taxes less subsidies on investment"], table5aColDict["Total"]] = sum(IO[IORowDict["P4"],IOCapForm]);
table5aTaxes = findall(x -> occursin("taxes", lowercase(x)), table5aNameRow);
table5aTaxes = table5aTaxes[Not(3)];
table5a[table5aRowDict["Total indirect taxes"], table5aColDict["Total"]] = sum(table5a[table5aTaxes,table5aColDict["Total"]]);
table5a[table5aRowDict["Total fixed capital expenditure"], table5aColDict["Total"]] = sum(table5a[Not(table5aRowDict["Total indirect taxes"]),table5aColDict["Total"]]);

#creating index variables for the measurements that we want
ASNAHouseCapTotCapForm = findall(x -> occursin("Gross fixed capital formation ;", x), string.(ASNAHouseCap[1,:]));
ASNANonFinCapTotCapForm = findall(x -> occursin("Gross fixed capital formation ;", x), string.(ASNANonFinCap[1,:]));
ASNAFinCapTotCapForm = findall(x -> occursin("Gross fixed capital formation ;", x), string.(ASNAFinCap[1,:]));
ASNAGenGovCapTotCapForm = findall(x -> occursin("General government ;  Gross fixed capital formation ;", x), string.(ASNAGovCap[1,:]));

#filling in totals row from ASNA Data
table5a[table5aRowDict["Total fixed capital expenditure"], table5aColDict["Households"]] = first(ASNAHouseCap[ASNAYearRow, ASNAHouseCapTotCapForm]);
table5a[table5aRowDict["Total fixed capital expenditure"], table5aColDict["Non-Financial Corporations"]] = first(ASNANonFinCap[ASNAYearRow, ASNANonFinCapTotCapForm]);
table5a[table5aRowDict["Total fixed capital expenditure"], table5aColDict["Financial Corporations"]] = first(ASNAFinCap[ASNAYearRow, ASNAFinCapTotCapForm]);
table5a[table5aRowDict["Total fixed capital expenditure"], table5aColDict["General Government"]] = first(ASNAGovCap[ASNAYearRow, ASNAGenGovCapTotCapForm]);

#filling in non-total values note- i think its including the taxes subtotal in the distribution
for i in [1:1:length(table5aRowDict)-1;];
    for ring in [1:1:length(table5aColDict)-1;];
        table5a[i,ring] = table5a[table5aRowDict["Total fixed capital expenditure"],ring]*table5a[i,table5aColDict["Total"]] / table5a[length(table5aRowDict),length(table5aColDict)];
        table5a[table5aRowDict["Total indirect taxes"], ring] = sum(table5a[table5aTaxes, ring]);
    end
end

#creating table 5b - allocation of investment expenditure (broken into subsections for dict referencing purposes)
#subsection b is fixed capital expenditure
table5bNameCol = ["Households", "Non-Financial Corporations", "Financial Corporations", "General Government", "Total"];
table5bNameRow = ["Domestic Commodities", "Imported Commodities, complementary", "Imported Commodities, competing",
"Taxes less subsidies on products", "Total change in inventories"];
table5bRowDict = Dict(table5bNameRow .=> [1:1:length(table5bNameRow);]);
table5bColDict = Dict(table5bNameCol .=> [1:1:length(table5bNameCol);]);
table5b = zeros(length(table5bNameRow), length(table5bNameCol));
table5b = NamedArray(table5b);
setnames!(table5b, string.(table5bNameRow), 1);
setnames!(table5b, string.(table5bNameCol), 2);

#filling in totals column from corresponding IO data
table5b[table5bRowDict["Domestic Commodities"], table5bColDict["Total"]] = sum(IO[IORowDict["T1"],IOChangeInv]);
table5b[table5bRowDict["Imported Commodities, complementary"], table5bColDict["Total"]] = sum(IO[IORowDict["P5"],IOChangeInv]);
table5b[table5bRowDict["Imported Commodities, competing"], table5bColDict["Total"]] = sum(IO[IORowDict["P6"],IOChangeInv]);
table5b[table5bRowDict["Taxes less subsidies on products"], table5bColDict["Total"]] = sum(IO[IORowDict["P3"],IOChangeInv]);
table5b[table5bRowDict["Total change in inventories"], table5bColDict["Total"]] = sum(table5b[:,table5bColDict["Total"]]);

#creating index variables for the measurements that we want
ASNAHouseCapChangeInv = findall(x -> occursin("Changes in inventories ;", x), string.(ASNAHouseCap[1,:]));
ASNANonFinCapChangeInv = findall(x -> occursin("Changes in inventories ;", x), string.(ASNANonFinCap[1,:]));
ASNAFinCapChangeInv = findall(x -> occursin("Changes in inventories ;", x), string.(ASNAFinCap[1,:]));
ASNAGenGovCapChangeInv = findall(x -> occursin("General government ;  Changes in inventories ;", x), string.(ASNAGovCap[1,:]));

#filling in totals row from ASNA Data
table5b[table5bRowDict["Total change in inventories"], table5bColDict["Households"]] = first(ASNAHouseCap[ASNAYearRow, ASNAHouseCapChangeInv]);
table5b[table5bRowDict["Total change in inventories"], table5bColDict["Non-Financial Corporations"]] = first(ASNANonFinCap[ASNAYearRow, ASNANonFinCapChangeInv]);
table5b[table5bRowDict["Total change in inventories"], table5bColDict["Financial Corporations"]] = first(ASNAFinCap[ASNAYearRow, ASNAFinCapChangeInv]);
table5b[table5bRowDict["Total change in inventories"], table5bColDict["General Government"]] = first(ASNAGovCap[ASNAYearRow, ASNAGenGovCapChangeInv]);

#calculate non-total values with lagrangian optimisation
table5bScalingFact = abs(minimum(table5b)) * 2;
mod5b = Model(Ipopt.Optimizer);
@variable(mod5b, x[1:(length(table5bNameRow)-1), 1:(length(table5bNameCol)-1)]);
@NLobjective(mod5b, Min, sum((x[i,j] - table5bScalingFact) ^ 2 for i in 1:(length(table5bNameRow)-1), j in 1:(length(table5bNameCol)-1)));
for i in 1:(length(table5bNameRow)-1);
    @constraint(mod5b, sum(x[:,i]) == table5b[table5bRowDict["Total change in inventories"],i]+table5bScalingFact);
end;
for i in 1:(length(table5bNameCol)-1);
    @constraint(mod5b, sum(x[i,:]) == table5b[i,table5bColDict["Total"]]+table5bScalingFact);
end;
optimize!(mod5b);

#plug back into table 5b
table5b[1:(length(table5bNameRow)-1),1:(length(table5bNameCol)-1)]=value.(x).-table5bScalingFact/4;

#creating table 5c - allocation of investment expenditure (broken into subsections for dict referencing purposes)
#subsection c is totals
table5cNameCol = ["Households", "Non-Financial Corporations", "Financial Corporations", "General Government", "Total"];
table5cNameRow = ["Domestic Commodities", "Imported Commodities", "Taxes less subsidies on products", "Other taxes less subsidies on investment", "Total investment expenditure"];
table5cRowDict = Dict(table5cNameRow .=> [1:1:length(table5cNameRow);]);
table5cColDict = Dict(table5cNameCol .=> [1:1:length(table5cNameCol);]);
table5c = zeros(length(table5cNameRow), length(table5cNameCol));
table5c = NamedArray(table5c);
setnames!(table5c, string.(table5cNameRow), 1);
setnames!(table5c, string.(table5cNameCol), 2);

#do totals calcuations to get all values in 5c
table5c[table5cRowDict["Domestic Commodities"],:] = (table5a[table5aRowDict["Domestic Commodities"],:] +
table5b[table5bRowDict["Domestic Commodities"],:]);
table5c[table5cRowDict["Imported Commodities"],:] = sum(eachrow(table5a[[table5aRowDict["Imported Commodities, competing"],table5aRowDict["Imported Commodities, complementary"]],:] +
table5b[[table5bRowDict["Imported Commodities, competing"],table5bRowDict["Imported Commodities, complementary"]],:]));
table5c[table5cRowDict["Taxes less subsidies on products"],:] = table5a[table5aRowDict["Taxes less subsidies on products"],:];
table5c[table5cRowDict["Other taxes less subsidies on investment"],:] = (table5a[table5aRowDict["Other taxes less subsidies on investment"],:] +
table5b[table5bRowDict["Taxes less subsidies on products"],:]);
table5c[table5cRowDict["Total investment expenditure"],:] = (table5a[table5aRowDict["Total fixed capital expenditure"],:] +
table5b[table5bRowDict["Total change in inventories"],:]);

#table 6
#importing relevant ASNA data
ASNAHouseInc = ExcelReaders.readxlsheet("."*pathmark*"data"*pathmark*"ASNAData"*pathmark*"5204036_Household_Income_Account.xls", "Data1");
ASNANonFinInc = ExcelReaders.readxlsheet("."*pathmark*"data"*pathmark*"ASNAData"*pathmark*"5204017_NonFin_Corp_Income_Account.xls", "Data1");
ASNAFinInc = ExcelReaders.readxlsheet("."*pathmark*"data"*pathmark*"ASNAData"*pathmark*"5204025_Fin_Corp_Income_Account.xls", "Data1");
ASNAGovInc = ExcelReaders.readxlsheet("."*pathmark*"data"*pathmark*"ASNAData"*pathmark*"5204030_GenGov_Income_Account.xls", "Data1");
ASNAExtInc = ExcelReaders.readxlsheet("."*pathmark*"data"*pathmark*"ASNAData"*pathmark*"5204043_External_Accounts.xls", "Data1");
#initialising table
tableName = ["Households", "Non-Financial Corporations", "Financial Corporations", "General Government", "External", "Total"];
tableDict = Dict(tableName .=> [1:1:length(tableName);]);
table6 = zeros(length(tableName),length(tableName));
table6 = NamedArray(table6);
setnames!(table6, string.(tableName), 1);
setnames!(table6, string.(tableName), 2);
#allocating total column and row data
table6[tableDict["Total"],tableDict["Households"]] = (first(ASNAHouseInc[ASNAYearRow,findall(x -> occursin("receivable - Interest", x), string.(ASNAHouseInc[1,:]))])
+first(ASNAHouseInc[ASNAYearRow,findall(x -> occursin("receivable - Imputed interest", x), string.(ASNAHouseInc[1,:]))]));
table6[tableDict["Total"],tableDict["Non-Financial Corporations"]] = (first(ASNANonFinInc[ASNAYearRow,findall(x -> occursin("receivable - Interest", x), string.(ASNANonFinInc[1,:]))])
+first(ASNANonFinInc[ASNAYearRow,findall(x -> occursin("Property income attributed to insurance policyholders", x), string.(ASNANonFinInc[1,:]))]));
table6[tableDict["Total"],tableDict["Financial Corporations"]] = first(ASNAFinInc[ASNAYearRow,findall(x -> occursin("receivable - Interest", x), string.(ASNAFinInc[1,:]))]);
table6[tableDict["Total"],tableDict["General Government"]] = first(ASNAGovInc[ASNAYearRow,findall(x -> occursin("General government ;  Property income receivable - Interest ;", x), string.(ASNAGovInc[1,:]))]);
table6[tableDict["Total"],tableDict["External"]] = first(ASNAExtInc[ASNAYearRow,findall(x -> occursin("receivable - Interest", x), string.(ASNAExtInc[1,:]))]);

table6[tableDict["Households"],tableDict["Total"]] = sum(ASNAHouseInc[ASNAYearRow,findall(x -> occursin("Property income payable - Interest", x), string.(ASNAHouseInc[1,:]))]);
table6[tableDict["Non-Financial Corporations"],tableDict["Total"]] = first(ASNANonFinInc[ASNAYearRow,findall(x -> occursin("Property income payable - Interest", x), string.(ASNANonFinInc[1,:]))]);
table6[tableDict["Financial Corporations"],tableDict["Total"]] = (first(ASNAFinInc[ASNAYearRow,findall(x -> occursin("Property income payable - Interest", x), string.(ASNAFinInc[1,:]))])
+first(ASNAFinInc[ASNAYearRow,findall(x -> occursin("payable - Property income attributed to insurance policy holders", x), string.(ASNAFinInc[1,:]))]));
table6[tableDict["General Government"],tableDict["Total"]] = first(ASNAGovInc[ASNAYearRow,findall(x -> occursin("General government ;  Property income payable - Total interest ;", x), string.(ASNAGovInc[1,:]))]);
table6[tableDict["External"],tableDict["Total"]] = first(ASNAExtInc[ASNAYearRow,findall(x -> occursin("Property income payable - Interest", x), string.(ASNAExtInc[1,:]))]);

if 0.98*sum(table6[:,length(tableName)])<sum(table6[length(tableName),:])<1.02*sum(table6[:,length(tableName)])
else
    error("Large discrepancy in row and column total sums in table 6")
end

#assuming gov only receives interest from fin corps
table6[tableDict["Financial Corporations"], tableDict["General Government"]] = table6[tableDict["Total"], tableDict["General Government"]];
#from other ASNA data
table6[tableDict["Financial Corporations"], tableDict["Non-Financial Corporations"]] = first(ASNANonFinInc[ASNAYearRow,findall(x -> occursin("attributed", x), string.(ASNANonFinInc[1,:]))]);
table6[tableDict["Financial Corporations"], tableDict["Households"]] = first(ASNAFinInc[ASNAYearRow,findall(x -> occursin("attributed", x), string.(ASNAFinInc[1,:]))])-first(ASNANonFinInc[ASNAYearRow,findall(x -> occursin("attributed", x), string.(ASNANonFinInc[1,:]))]);
table6[tableDict["General Government"], tableDict["Households"]] = first(ASNAGovInc[ASNAYearRow,findall(x -> occursin("General government ;  Property income payable - Interest - On unfunded superannuation liabilities ;", x), string.(ASNAGovInc[1,:]))]);

#solve for missing values with scaling
table6Step3 = zeros(length(tableName),length(tableName));
table6Step3Row = [tableDict["Non-Financial Corporations"],tableDict["Financial Corporations"],tableDict["General Government"]];
table6Step3Col = [tableDict["Households"],tableDict["External"]];
for i in table6Step3Col;
    for ring in table6Step3Row;
        table6Step3[ring,i] = (table6[tableDict["Total"],i]-sum(table6[1:(length(tableName)-1),i]))*(
            table6[ring,tableDict["Total"]]-sum(table6[ring,1:(length(tableName)-1)]))/sum(table6[
            table6Step3Row,tableDict["Total"]]-sum(eachcol(table6[table6Step3Row,1:(length(tableName)-1)])));
    end
end
table6 = table6+table6Step3;

table6Step4 = zeros(length(tableName),length(tableName));
table6Step4Row = [1:1:(length(tableName)-1);];
table6Step4Col = [tableDict["Non-Financial Corporations"],tableDict["Financial Corporations"]];
for i in table6Step4Col;
    for ring in table6Step4Row;
        table6Step4[ring,i] = (table6[tableDict["Total"],i]-sum(table6[1:(length(tableName)-1),i]))*(
            table6[ring,tableDict["Total"]]-sum(table6[ring,1:(length(tableName)-1)]))/sum(table6[
            table6Step4Row,tableDict["Total"]]-sum(eachcol(table6[table6Step4Row,1:(length(tableName)-1)])));
    end
end
table6 = table6+table6Step4;

#table 7
#initialising table
table7 = zeros(length(tableName),length(tableName));
table7 = NamedArray(table7);
setnames!(table7, string.(tableName), 1);
setnames!(table7, string.(tableName), 2);
#allocating total column and row data
table7[tableDict["Total"],tableDict["Households"]] = first(ASNAHouseInc[ASNAYearRow,findall(x -> occursin("receivable - Dividends", x), string.(ASNAHouseInc[1,:]))]);
table7[tableDict["Total"],tableDict["Non-Financial Corporations"]] = first(ASNANonFinInc[ASNAYearRow,findall(x -> occursin("receivable - Dividends", x), string.(ASNANonFinInc[1,:]))]);
table7[tableDict["Total"],tableDict["Financial Corporations"]] = first(ASNAFinInc[ASNAYearRow,findall(x -> occursin("receivable - Dividends", x), string.(ASNAFinInc[1,:]))]);
table7[tableDict["Total"],tableDict["General Government"]] = first(ASNAGovInc[ASNAYearRow,findall(x -> occursin("General government ;  Property income receivable - Dividends ;", x), string.(ASNAGovInc[1,:]))]);
table7[tableDict["Total"],tableDict["External"]] = first(ASNAExtInc[ASNAYearRow,findall(x -> occursin("receivable - Dividends", x), string.(ASNAExtInc[1,:]))]);

table7[tableDict["Households"],tableDict["Total"]] = 0.0
table7[tableDict["Non-Financial Corporations"],tableDict["Total"]] = first(ASNANonFinInc[ASNAYearRow,findall(x -> occursin("payable - Dividends", x), string.(ASNANonFinInc[1,:]))]);
table7[tableDict["Financial Corporations"],tableDict["Total"]] = first(ASNAFinInc[ASNAYearRow,findall(x -> occursin("payable - Dividends", x), string.(ASNAFinInc[1,:]))]);
table7[tableDict["General Government"],tableDict["Total"]] = 0.0;
table7[tableDict["External"],tableDict["Total"]] = first(ASNAExtInc[ASNAYearRow,findall(x -> occursin("payable - Dividends", x), string.(ASNAExtInc[1,:]))]);

if 0.98*sum(table7[:,length(tableName)])<sum(table7[length(tableName),:])<1.02*sum(table7[:,length(tableName)])
else
    error("Large discrepancy in row and column total sums in table 7")
end

#assuming gov only receives dividends from non-fin corps
table7[tableDict["Non-Financial Corporations"], tableDict["General Government"]] = table7[tableDict["Total"], tableDict["General Government"]];

#solve for missing values with scaling
table7Step3 = zeros(length(tableName),length(tableName));
table7Step3Row = [tableDict["Non-Financial Corporations"],tableDict["Financial Corporations"]];
table7Step3Col = [tableDict["External"]];
for i in table7Step3Col;
    for ring in table7Step3Row;
        table7Step3[ring,i] = (table7[tableDict["Total"],i]-sum(table7[1:(length(tableName)-1),i]))*(
            table7[ring,tableDict["Total"]]-sum(table7[ring,1:(length(tableName)-1)]))/sum(table7[
            table7Step3Row,tableDict["Total"]]-sum(eachcol(table7[table7Step3Row,1:(length(tableName)-1)])));
    end
end
table7 = table7+table7Step3;

table7Step4 = zeros(length(tableName),length(tableName));
table7Step4Row = [tableDict["Non-Financial Corporations"],tableDict["Financial Corporations"], tableDict["External"]];
table7Step4Col = [tableDict["Households"],tableDict["Non-Financial Corporations"],tableDict["Financial Corporations"]];
for i in table7Step4Col;
    for ring in table7Step4Row;
        table7Step4[ring,i] = (table7[tableDict["Total"],i]-sum(table7[1:(length(tableName)-1),i]))*(
            table7[ring,tableDict["Total"]]-sum(table7[ring,1:(length(tableName)-1)]))/sum(table7[
            table7Step4Row,tableDict["Total"]]-sum(eachcol(table7[table7Step4Row,1:(length(tableName)-1)])));
    end
end
table7 = table7+table7Step4;

#table 8
#initialising table
table8 = zeros(length(tableName),length(tableName));
table8 = NamedArray(table8);
setnames!(table8, string.(tableName), 1);
setnames!(table8, string.(tableName), 2);
#allocating total column and row data
table8[tableDict["Total"],tableDict["Households"]] = first(ASNAHouseInc[ASNAYearRow,findall(x -> occursin("receivable - Reinvested", x), string.(ASNAHouseInc[1,:]))]);
table8[tableDict["Total"],tableDict["Non-Financial Corporations"]] = first(ASNANonFinInc[ASNAYearRow,findall(x -> occursin("receivable - Reinvested", x), string.(ASNANonFinInc[1,:]))]);
table8[tableDict["Total"],tableDict["Financial Corporations"]] = first(ASNAFinInc[ASNAYearRow,findall(x -> occursin("receivable - Reinvested", x), string.(ASNAFinInc[1,:]))]);
table8[tableDict["Total"],tableDict["General Government"]] = 0.0;
table8[tableDict["Total"],tableDict["External"]] = first(ASNAExtInc[ASNAYearRow,findall(x -> occursin("receivable - Reinvested", x), string.(ASNAExtInc[1,:]))]);

table8[tableDict["Households"],tableDict["Total"]] = 0.0
table8[tableDict["Non-Financial Corporations"],tableDict["Total"]] = first(ASNANonFinInc[ASNAYearRow,findall(x -> occursin("payable - Reinvested", x), string.(ASNANonFinInc[1,:]))]);
table8[tableDict["Financial Corporations"],tableDict["Total"]] = first(ASNAFinInc[ASNAYearRow,findall(x -> occursin("payable - Reinvested", x), string.(ASNAFinInc[1,:]))]);
table8[tableDict["General Government"],tableDict["Total"]] = 0.0;
table8[tableDict["External"],tableDict["Total"]] = first(ASNAExtInc[ASNAYearRow,findall(x -> occursin("payable - Reinvested", x), string.(ASNAExtInc[1,:]))]);

if 0.98*sum(table8[:,length(tableName)])<sum(table8[length(tableName),:])<1.02*sum(table8[:,length(tableName)])
else
    error("Large discrepancy in row and column total sums in table 8")
end
#assuming that values can be negative
#solve missing values with Ipopt
mod8 = Model(Ipopt.Optimizer);
@variable(mod8, x[1:3, 1:4]>=0);
@NLobjective(mod8, Min, sum((x[i,j]) ^ 2 for i in 1:3, j in 1:4));
@constraint(mod8, sum(x[:,1]) == table8[tableDict["Total"],1]);
@constraint(mod8, sum(x[:,2]) == table8[tableDict["Total"],2]);
@constraint(mod8, sum(x[:,3]) == table8[tableDict["Total"],3]);
@constraint(mod8, sum(x[:,4]) == table8[tableDict["Total"],5]);
@constraint(mod8, sum(x[1,:]) == table8[2,tableDict["Total"]]);
@constraint(mod8, sum(x[2,:]) == table8[3,tableDict["Total"]]);
@constraint(mod8, sum(x[3,:]) == table8[5,tableDict["Total"]]);
#@constraint(mod8, x[3,4] == 0);
optimize!(mod8);
table8[[2,3,5], [1,2,3,5]] = table8[[2,3,5], [1,2,3,5]] + value.(x);

#=spread the external receivable totals between fin and non-fin Corporations
table8Step3 = zeros(length(tableName),length(tableName));
table8Step3Row = [tableDict["Non-Financial Corporations"],tableDict["Financial Corporations"]];
table8Step3Col = [tableDict["External"]];
for i in table8Step3Col;
    for ring in table8Step3Row;
        table8Step3[ring,i] = (table8[tableDict["Total"],i]-sum(table8[1:(length(tableName)-1),i]))*(
            table8[ring,tableDict["Total"]]-sum(table8[ring,1:(length(tableName)-1)]))/sum(table8[
            table8Step3Row,tableDict["Total"]]-sum(eachcol(table8[table8Step3Row,1:(length(tableName)-1)])));
    end
end
table8 = table8+table8Step3;
=#

#table 9
#initialising table
table9 = zeros(length(tableName),length(tableName));
table9 = NamedArray(table9);
setnames!(table9, string.(tableName), 1);
setnames!(table9, string.(tableName), 2);
#allocating total column and row data
table9[tableDict["Total"],tableDict["Households"]] = first(ASNAHouseInc[ASNAYearRow,findall(x -> occursin("receivable - Rent on natural", x), string.(ASNAHouseInc[1,:]))]);
table9[tableDict["Total"],tableDict["Non-Financial Corporations"]] = first(ASNANonFinInc[ASNAYearRow,findall(x -> occursin("receivable - Rent on natural", x), string.(ASNANonFinInc[1,:]))]);
table9[tableDict["Total"],tableDict["Financial Corporations"]] = first(ASNAFinInc[ASNAYearRow,findall(x -> occursin("receivable - Rent on natural", x), string.(ASNAFinInc[1,:]))]);
table9[tableDict["Total"],tableDict["General Government"]] = first(ASNAGovInc[ASNAYearRow,findall(x -> occursin("General government ;  Property income receivable - Rent on natural assets ;", x), string.(ASNAGovInc[1,:]))]);
table9[tableDict["Total"],tableDict["External"]] = 0.0;

table9[tableDict["Households"],tableDict["Total"]] =  first(ASNAHouseInc[ASNAYearRow,findall(x -> occursin("payable - Rent on natural", x), string.(ASNAHouseInc[1,:]))]);
table9[tableDict["Non-Financial Corporations"],tableDict["Total"]] = first(ASNANonFinInc[ASNAYearRow,findall(x -> occursin("payable - Rent on natural", x), string.(ASNANonFinInc[1,:]))]);
table9[tableDict["Financial Corporations"],tableDict["Total"]] = first(ASNAFinInc[ASNAYearRow,findall(x -> occursin("payable - Rent on natural", x), string.(ASNAFinInc[1,:]))]);
table9[tableDict["General Government"],tableDict["Total"]] =  0.0;
table9[tableDict["External"],tableDict["Total"]] = 0.0;

if 0.98*sum(table9[:,length(tableName)])<sum(table9[length(tableName),:])<1.02*sum(table9[:,length(tableName)])
else
    error("Large discrepancy in row and column total sums in table 9")
end

#solve missing values with Ipopt
mod9 = Model(Ipopt.Optimizer);
@variable(mod9, x[1:(length(tableName)-1), 1:(length(tableName)-1)]>=0);
@NLobjective(mod9, Min, sum((x[i,j]) ^ 2 for i in 1:(length(tableName)-1), j in 1:(length(tableName)-1)));
for i in 1:(length(tableName)-1);
    @constraint(mod9, sum(x[:,i]) == table9[tableDict["Total"],i]-sum(table9[1:(length(tableName)-1),i]));
    @constraint(mod9, sum(x[i,:]) == table9[i,tableDict["Total"]]-sum(table9[i,1:(length(tableName)-1)]));
end;
@constraint(mod9, x[tableDict["Households"],tableDict["Households"]] == 0);
@constraint(mod9, x[tableDict["Households"],tableDict["Non-Financial Corporations"]] == 0);
optimize!(mod9);
table9[1:(length(tableName)-1), 1:(length(tableName)-1)] = table9[1:(length(tableName)-1), 1:(length(tableName)-1)] + value.(x);

#table 10
#initialising table
table10 = zeros(length(tableName),length(tableName));
table10 = NamedArray(table10);
setnames!(table10, string.(tableName), 1);
setnames!(table10, string.(tableName), 2);
#allocating total column and row data
table10[tableDict["Total"],tableDict["Households"]] = first(ASNAHouseInc[ASNAYearRow,findall(x -> occursin("Social assistance benefits", x), string.(ASNAHouseInc[1,:]))]);
table10[tableDict["General Government"],tableDict["Total"]] =  table10[tableDict["Total"],tableDict["Households"]];
#filling in missing values
table10[tableDict["General Government"],tableDict["Households"]] =  table10[tableDict["Total"],tableDict["Households"]];

#table 11
#initialising table
table11 = zeros(length(tableName),length(tableName));
table11 = NamedArray(table11);
setnames!(table11, string.(tableName), 1);
setnames!(table11, string.(tableName), 2);
#allocating total column and row data
#assuming that the large discrepancy in the paid and received claims is a result of households not reporting small claims
#table11[tableDict["Total"],tableDict["Households"]] = first(ASNAHouseInc[ASNAYearRow,findall(x -> occursin("receivable - Non-life", x), string.(ASNAHouseInc[1,:]))]);
table11[tableDict["Total"],tableDict["Non-Financial Corporations"]] = first(ASNANonFinInc[ASNAYearRow,findall(x -> occursin("receivable - Non-life", x), string.(ASNANonFinInc[1,:]))]);
table11[tableDict["Financial Corporations"],tableDict["Total"]] = first(ASNAFinInc[ASNAYearRow,findall(x -> occursin("payable - Non-life", x), string.(ASNAFinInc[1,:]))]);
table11[tableDict["Total"],tableDict["Households"]]= table11[tableDict["Financial Corporations"],tableDict["Total"]] - table11[tableDict["Total"],tableDict["Non-Financial Corporations"]];
if 0.98*sum(table11[:,length(tableName)])<sum(table11[length(tableName),:])<1.02*sum(table11[:,length(tableName)])
else
    error("Large discrepancy in row and column total sums in table 11")
end
#filling in missing values
table11[tableDict["Financial Corporations"],tableDict["Households"]]=table11[tableDict["Total"],tableDict["Households"]];
table11[tableDict["Financial Corporations"],tableDict["Non-Financial Corporations"]] = table11[tableDict["Total"],tableDict["Non-Financial Corporations"]];

#table 12
#initialising table
table12 = zeros(length(tableName),length(tableName));
table12 = NamedArray(table12);
setnames!(table12, string.(tableName), 1);
setnames!(table12, string.(tableName), 2);
#allocating total column and row data
#also assuming in the same manner as in table 11 that the missing data is from the household side
table12[tableDict["Total"],tableDict["Households"]] = 0.0;
table12[tableDict["Total"],tableDict["Non-Financial Corporations"]] = 0.0;
table12[tableDict["Total"],tableDict["Financial Corporations"]] = first(ASNAFinInc[ASNAYearRow,findall(x -> occursin("receivable - Net non-life insurance premiums", x), string.(ASNAFinInc[1,:]))]);
table12[tableDict["Total"],tableDict["General Government"]] = 0.0;
table12[tableDict["Total"],tableDict["External"]] = first(ASNAExtInc[ASNAYearRow,findall(x -> occursin("receivable - Non-life insurance transfers", x), string.(ASNAExtInc[1,:]))]);

#table12[tableDict["Households"],tableDict["Total"]] =  first(ASNAHouseInc[ASNAYearRow,findall(x -> occursin("payable - Net non-life insurance premiums", x), string.(ASNAHouseInc[1,:]))]);
table12[tableDict["Non-Financial Corporations"],tableDict["Total"]] = first(ASNANonFinInc[ASNAYearRow,findall(x -> occursin("payable - Net non-life insurance premiums", x), string.(ASNANonFinInc[1,:]))]);
table12[tableDict["Financial Corporations"],tableDict["Total"]] = 0.0;
table12[tableDict["General Government"],tableDict["Total"]] =  0.0;
table12[tableDict["External"],tableDict["Total"]] = first(ASNAExtInc[ASNAYearRow,findall(x -> occursin("payable - Non-life insurance transfers", x), string.(ASNAExtInc[1,:]))]);
table12[tableDict["Households"],tableDict["Total"]] = sum(table12[length(tableName),:])-sum(table12[:,length(tableName)]);

if 0.98*sum(table12[:,length(tableName)])<sum(table12[length(tableName),:])<1.02*sum(table12[:,length(tableName)])
else
    error("Large discrepancy in row and column total sums in table 12")
end

#filling in empty values
table12[tableDict["External"],tableDict["Financial Corporations"]] = table12[tableDict["External"],tableDict["Total"]];
#=
mod12 = Model(Ipopt.Optimizer);
@variable(mod12, x[1:(length(tableName)-1), 1:(length(tableName)-1)]>=0);
@NLobjective(mod12, Min, sum((x[i,j])^2 for i in 1:(length(tableName)-1), j in 1:(length(tableName)-1)));
for i in 1:(length(tableName)-1);
    @constraint(mod12, sum(x[:,i]) == table12[tableDict["Total"],i]-sum(table12[1:(length(tableName)-1),i]));
    @constraint(mod12, sum(x[i,:]) == table12[i,tableDict["Total"]]-sum(table12[i,1:(length(tableName)-1)]));
end;
optimize!(mod12);
table12[1:(length(tableName)-1), 1:(length(tableName)-1)] = table12[1:(length(tableName)-1), 1:(length(tableName)-1)] + value.(x);
=#
#spread the external receivable totals between fin and non-fin Corporations
table12Step3 = zeros(length(tableName),length(tableName));
table12Step3Row = [tableDict["Non-Financial Corporations"],tableDict["Households"]];
table12Step3Col = [tableDict["External"],tableDict["Financial Corporations"]];
for i in table12Step3Col;
    for ring in table12Step3Row;
        table12Step3[ring,i] = (table12[tableDict["Total"],i]-sum(table12[1:(length(tableName)-1),i]))*(
            table12[ring,tableDict["Total"]]-sum(table12[ring,1:(length(tableName)-1)]))/sum(table12[
            table12Step3Row,tableDict["Total"]]-sum(eachcol(table12[table12Step3Row,1:(length(tableName)-1)])));
    end
end
table12 = table12+table12Step3;

#table 13
#initialising table
table13 = zeros(length(tableName),length(tableName));
table13 = NamedArray(table13);
setnames!(table13, string.(tableName), 1);
setnames!(table13, string.(tableName), 2);
#allocating total column and row data
table13[tableDict["Total"],tableDict["Households"]] = (first(ASNAHouseInc[ASNAYearRow,findall(x -> occursin("receivable - Total other current transfers ;", x), string.(ASNAHouseInc[1,:]))])
+first(ASNAHouseInc[ASNAYearRow,findall(x -> occursin("receivable - Current transfers to ", x), string.(ASNAHouseInc[1,:]))]));
table13[tableDict["Total"],tableDict["Non-Financial Corporations"]] = first(ASNANonFinInc[ASNAYearRow,findall(x -> occursin("receivable - Other current transfers ;", x), string.(ASNANonFinInc[1,:]))]);
table13[tableDict["Total"],tableDict["Financial Corporations"]] = first(ASNAFinInc[ASNAYearRow,findall(x -> occursin("receivable - Other current transfers ;", x), string.(ASNAFinInc[1,:]))]);
table13[tableDict["Total"],tableDict["General Government"]] = first(ASNAGovInc[ASNAYearRow,findall(x -> occursin("General government ;  Secondary income receivable - Other current transfers ;", x), string.(ASNAGovInc[1,:]))]);
table13[tableDict["Total"],tableDict["External"]] = (first(ASNAExtInc[ASNAYearRow,findall(x -> occursin("receivable - Other current transfers ;", x), string.(ASNAExtInc[1,:]))])
+first(ASNAExtInc[ASNAYearRow,findall(x -> occursin("receivable - Current international cooperation ;", x), string.(ASNAExtInc[1,:]))]));

table13[tableDict["Households"],tableDict["Total"]] =  first(ASNAHouseInc[ASNAYearRow,findall(x -> occursin("payable - Total other current transfers ;", x), string.(ASNAHouseInc[1,:]))]);
table13[tableDict["Non-Financial Corporations"],tableDict["Total"]] = (first(ASNANonFinInc[ASNAYearRow,findall(x -> occursin("payable - Current transfers to non-profit institutions ;", x), string.(ASNANonFinInc[1,:]))])
+ first(ASNANonFinInc[ASNAYearRow,findall(x -> occursin("payable - Other current transfers ;", x), string.(ASNANonFinInc[1,:]))]));
table13[tableDict["Financial Corporations"],tableDict["Total"]] = first(ASNAFinInc[ASNAYearRow,findall(x -> occursin("payable - Other current transfers ;", x), string.(ASNAFinInc[1,:]))]);
table13[tableDict["General Government"],tableDict["Total"]] =  sum(ASNAGovInc[ASNAYearRow,findall(x -> occursin("General government ;  Secondary income payable - Other current transfers -", x), string.(ASNAGovInc[1,:]))]);
table13[tableDict["External"],tableDict["Total"]] = first(ASNAExtInc[ASNAYearRow,findall(x -> occursin("payable - Other current transfers", x), string.(ASNAExtInc[1,:]))]);

if 0.98*sum(table13[:,length(tableName)])<sum(table13[length(tableName),:])<1.02*sum(table13[:,length(tableName)])
else
    error("Large discrepancy in row and column total sums in table 13")
end

#filling in empty values
table13Step3 = zeros(length(tableName),length(tableName));
table13Step3Row = [tableDict["Non-Financial Corporations"],tableDict["Financial Corporations"],tableDict["General Government"],tableDict["External"]];
table13Step3Col = [tableDict["General Government"]];
for i in table13Step3Col;
    for ring in table13Step3Row;
        table13Step3[ring,i] = (table13[tableDict["Total"],i]-sum(table13[1:(length(tableName)-1),i]))*(
            table13[ring,tableDict["Total"]]-sum(table13[ring,1:(length(tableName)-1)]))/sum(table13[
            table13Step3Row,tableDict["Total"]]-sum(eachcol(table13[table13Step3Row,1:(length(tableName)-1)])));
    end
end
table13 = table13+table13Step3;

table13Step3 = zeros(length(tableName),length(tableName));
table13Step3Row = [tableDict["Households"],tableDict["Non-Financial Corporations"],tableDict["Financial Corporations"],tableDict["General Government"],tableDict["External"]];
table13Step3Col = [tableDict["Households"],tableDict["Non-Financial Corporations"],tableDict["Financial Corporations"],tableDict["External"]];
for i in table13Step3Col;
    for ring in table13Step3Row;
        table13Step3[ring,i] = (table13[tableDict["Total"],i]-sum(table13[1:(length(tableName)-1),i]))*(
            table13[ring,tableDict["Total"]]-sum(table13[ring,1:(length(tableName)-1)]))/sum(table13[
            table13Step3Row,tableDict["Total"]]-sum(eachcol(table13[table13Step3Row,1:(length(tableName)-1)])));
    end
end
table13 = table13+table13Step3;

#table 14
#initialising table
table14 = zeros(length(tableName),length(tableName));
table14 = NamedArray(table14);
setnames!(table14, string.(tableName), 1);
setnames!(table14, string.(tableName), 2);
#filling in data from ASNA
table14[tableDict["General Government"],tableDict["Households"]] = first(ASNAHouseCap[ASNAYearRow,findall(x -> occursin("Capital transfers, receivable from general government ;", x), string.(ASNAHouseCap[1,:]))]);
table14[tableDict["Total"],tableDict["Households"]] = first(ASNAHouseCap[ASNAYearRow,findall(x -> occursin("Other capital transfers, receivable ;", x), string.(ASNAHouseCap[1,:]))])+table14[tableDict["General Government"],tableDict["Households"]];
table14[tableDict["General Government"],tableDict["Non-Financial Corporations"]] = first(ASNANonFinCap[ASNAYearRow,findall(x -> occursin("Capital transfers, receivable from general government ;", x), string.(ASNANonFinCap[1,:]))]);
table14[tableDict["Total"],tableDict["Non-Financial Corporations"]] = first(ASNANonFinCap[ASNAYearRow,findall(x -> occursin("Other capital transfers, receivable ;", x), string.(ASNANonFinCap[1,:]))])+table14[tableDict["General Government"],tableDict["Non-Financial Corporations"]];
table14[tableDict["General Government"],tableDict["Financial Corporations"]] = first(ASNAFinCap[ASNAYearRow,findall(x -> occursin("Capital transfers, receivable from general government ;", x), string.(ASNAFinCap[1,:]))]);
table14[tableDict["Total"],tableDict["Financial Corporations"]] = first(ASNAFinCap[ASNAYearRow,findall(x -> occursin("Other capital transfers, receivable ;", x), string.(ASNAFinCap[1,:]))])+table14[tableDict["General Government"],tableDict["Financial Corporations"]];
table14[tableDict["Households"],tableDict["General Government"]] = first(ASNAHouseCap[ASNAYearRow,findall(x -> occursin("Capital transfers, payable to general government ;", x), string.(ASNAHouseCap[1,:]))]);
table14[tableDict["Non-Financial Corporations"],tableDict["General Government"]] = first(ASNANonFinCap[ASNAYearRow,findall(x -> occursin("Capital transfers, payable to general government ;", x), string.(ASNANonFinCap[1,:]))]);
table14[tableDict["Financial Corporations"],tableDict["General Government"]] = first(ASNAFinCap[ASNAYearRow,findall(x -> occursin("Capital transfers, payable to general government ;", x), string.(ASNAFinCap[1,:]))]);
table14[tableDict["Total"],tableDict["General Government"]] = first(ASNAGovCap[ASNAYearRow,findall(x -> occursin("General government ;  Capital transfers, receivable ;", x), string.(ASNAGovCap[1,:]))]);
table14[tableDict["Total"],tableDict["External"]] = first(ASNAExtInc[ASNAYearRow,findall(x -> occursin("Capital transfers, receivable ;", x), string.(ASNAExtInc[1,:]))]);
table14[tableDict["Households"],tableDict["Total"]] = first(ASNAHouseCap[ASNAYearRow,findall(x -> occursin("Other capital transfers, payable ;", x), string.(ASNAHouseCap[1,:]))])+table14[tableDict["Households"],tableDict["General Government"]];
table14[tableDict["Non-Financial Corporations"],tableDict["Total"]] = first(ASNANonFinCap[ASNAYearRow,findall(x -> occursin("Other capital transfers, payable ;", x), string.(ASNANonFinCap[1,:]))])+table14[tableDict["Non-Financial Corporations"],tableDict["General Government"]];
table14[tableDict["Financial Corporations"],tableDict["Total"]] = first(ASNAFinCap[ASNAYearRow,findall(x -> occursin("Other capital transfers, payable ;", x), string.(ASNAFinCap[1,:]))])+table14[tableDict["Financial Corporations"],tableDict["General Government"]];
table14[tableDict["General Government"],tableDict["Total"]] = first(ASNAGovCap[ASNAYearRow,findall(x -> occursin("General government ;  Capital transfers, payable ;", x), string.(ASNAGovCap[1,:]))]);
table14[tableDict["External"],tableDict["Total"]] = first(ASNAExtInc[ASNAYearRow,findall(x -> occursin("Capital transfers, payable ;", x), string.(ASNAExtInc[1,:]))]);
if 0.98*sum(table14[:,length(tableName)])<sum(table14[length(tableName),:])<1.02*sum(table14[:,length(tableName)])
else
    error("Large discrepancy in row and column total sums in table 14")
end

#solve missing values with Ipopt
mod14 = Model(Ipopt.Optimizer);
@variable(mod14, x[1:(length(tableName)-1), 1:(length(tableName)-1)]>=0);
@NLobjective(mod14, Min, sum((x[i,j]) ^ 2 for i in 1:(length(tableName)-1), j in 1:(length(tableName)-1)));
for i in 1:(length(tableName)-1);
    @constraint(mod14, sum(x[:,i]) == table14[tableDict["Total"],i]-sum(table14[1:(length(tableName)-1),i]));
    @constraint(mod14, sum(x[i,:]) == table14[i,tableDict["Total"]]-sum(table14[i,1:(length(tableName)-1)]));
    @constraint(mod14, x[i,i] == 0);
end;

optimize!(mod14);
table14[1:(length(tableName)-1), 1:(length(tableName)-1)] = table14[1:(length(tableName)-1), 1:(length(tableName)-1)] + value.(x);

#table 15
#initialising table
table15 = zeros(length(tableName),length(tableName));
table15 = NamedArray(table15);
setnames!(table15, string.(tableName), 1);
setnames!(table15, string.(tableName), 2);
#allocating total column and row data
if first(ASNAHouseCap[ASNAYearRow,findall(x -> occursin("Acquisitions less disposals of non-produced non-financial assets ;", x), string.(ASNAHouseCap[1,:]))])<=0
    table15[tableDict["Households"],tableDict["Total"]] = abs(first(ASNAHouseCap[ASNAYearRow,findall(x -> occursin("Acquisitions less disposals of non-produced non-financial assets ;", x), string.(ASNAHouseCap[1,:]))]));
else
    table15[tableDict["Total"],tableDict["Households"]] = abs(first(ASNAHouseCap[ASNAYearRow,findall(x -> occursin("Acquisitions less disposals of non-produced non-financial assets ;", x), string.(ASNAHouseCap[1,:]))]));
end
if first(ASNANonFinCap[ASNAYearRow,findall(x -> occursin("Acquisitions less disposals of non-produced non-financial assets ;", x), string.(ASNANonFinCap[1,:]))])<=0
    table15[tableDict["Non-Financial Corporations"],tableDict["Total"]] = abs(first(ASNANonFinCap[ASNAYearRow,findall(x -> occursin("Acquisitions less disposals of non-produced non-financial assets ;", x), string.(ASNANonFinCap[1,:]))]));
else
    table15[tableDict["Total"],tableDict["Non-Financial Corporations"]] = abs(first(ASNANonFinCap[ASNAYearRow,findall(x -> occursin("Acquisitions less disposals of non-produced non-financial assets ;", x), string.(ASNANonFinCap[1,:]))]));
end
if first(ASNAFinCap[ASNAYearRow,findall(x -> occursin("Acquisitions less disposals of non-produced non-financial assets ;", x), string.(ASNAFinCap[1,:]))])<=0
    table15[tableDict["Financial Corporations"],tableDict["Total"]] = abs(first(ASNAFinCap[ASNAYearRow,findall(x -> occursin("Acquisitions less disposals of non-produced non-financial assets ;", x), string.(ASNAFinCap[1,:]))]));
else
    table15[tableDict["Total"],tableDict["Financial Corporations"]] = abs(first(ASNAFinCap[ASNAYearRow,findall(x -> occursin("Acquisitions less disposals of non-produced non-financial assets ;", x), string.(ASNAFinCap[1,:]))]));
end
if first(ASNAGovCap[ASNAYearRow,findall(x -> occursin("General government ;  Acquisitions less disposals of non-produced non-financial assets ;", x), string.(ASNAGovCap[1,:]))])<=0
    table15[tableDict["General Government"],tableDict["Total"]] = abs(first(ASNAGovCap[ASNAYearRow,findall(x -> occursin("General government ;  Acquisitions less disposals of non-produced non-financial assets ;", x), string.(ASNAGovCap[1,:]))]));
else
    table15[tableDict["Total"],tableDict["General Government"]] = abs(first(ASNAGovCap[ASNAYearRow,findall(x -> occursin("General government ;  Acquisitions less disposals of non-produced non-financial assets ;", x), string.(ASNAGovCap[1,:]))]));
end
if first(ASNAExtInc[ASNAYearRow,findall(x -> occursin("Acquisitions less disposals of non-produced non-financial assets ;", x), string.(ASNAExtInc[1,:]))])<=0
    table15[tableDict["External"],tableDict["Total"]] = abs(first(ASNAExtInc[ASNAYearRow,findall(x -> occursin("Acquisitions less disposals of non-produced non-financial assets ;", x), string.(ASNAExtInc[1,:]))]));
else
    table15[tableDict["Total"],tableDict["External"]] = abs(first(ASNAExtInc[ASNAYearRow,findall(x -> occursin("Acquisitions less disposals of non-produced non-financial assets ;", x), string.(ASNAExtInc[1,:]))]));
end

if 0.98*sum(table15[:,length(tableName)])<sum(table15[length(tableName),:])<1.02*sum(table15[:,length(tableName)])
else
    error("Large discrepancy in row and column total sums in table 15")
end

#solve missing values with Ipopt
mod15 = Model(Ipopt.Optimizer);
@variable(mod15, x[1:(length(tableName)-1), 1:(length(tableName)-1)]>=0);
@NLobjective(mod15, Min, sum((x[i,j]) ^ 2 for i in 1:(length(tableName)-1), j in 1:(length(tableName)-1)));
for i in 1:(length(tableName)-1);
    @constraint(mod15, sum(x[:,i]) == table15[tableDict["Total"],i]-sum(table15[1:(length(tableName)-1),i]));
    @constraint(mod15, sum(x[i,:]) == table15[i,tableDict["Total"]]-sum(table15[i,1:(length(tableName)-1)]));
    @constraint(mod15, x[i,i] == 0);
end;

optimize!(mod15);
table15[1:(length(tableName)-1), 1:(length(tableName)-1)] = table15[1:(length(tableName)-1), 1:(length(tableName)-1)] + value.(x);

#table 16a
#initialising table
table16aNameRow = ["Netsaving;", "Totalnetcapitaltransfers;", "Grossfixedcapitalformation;", "Netlending(+)/netborrowing(-);"];
table16aDataCol = [ASNAHouseCap, ASNANonFinCap, ASNAFinCap, ASNAGovCap, ASNAExtInc];
table16aRowDict = Dict(table16aNameRow .=> [1:1:length(table16aNameRow);]);
table16aDataDict = Dict(tableName[1:length(tableName)-1] .=> table16aDataCol);
table16a = zeros(length(table16aNameRow),length(tableName));
table16a = NamedArray(table16a);
setnames!(table16a, string.(table16aNameRow), 1);
setnames!(table16a, string.(tableName), 2);
#filling in values
for i in tableName
    if i == "Total"
        for ring in table16aNameRow
            table16a[table16aRowDict[ring],tableDict[i]]=sum(table16a[table16aRowDict[ring],:]);
        end
    elseif i == "General Government"
        for ring in table16aNameRow
            table16a[table16aRowDict[ring],tableDict[i]]=first(table16aDataDict[i][ASNAYearRow,findall(x -> occursin("Generalgovernment;"*ring, x), filter.(x -> !isspace(x), string.(table16aDataDict[i][1,:])))]);
        end
    else
        for ring in table16aNameRow
            if isempty(findall(x -> occursin(ring, x), filter.(x -> !isspace(x),string.(table16aDataDict[i][1,:]))))
                table16a[table16aRowDict[ring],tableDict[i]]=0.0;
            else
                table16a[table16aRowDict[ring],tableDict[i]]=first(table16aDataDict[i][ASNAYearRow,findall(x -> occursin(ring, x), filter.(x -> !isspace(x), string.(table16aDataDict[i][1,:])))]);
            end
        end
    end
end

#table 16b
#importing relevant data from ASNA
#importing relevant ASNA data
ASNAHouseFin = ExcelReaders.readxlsheet("."*pathmark*"data"*pathmark*"ASNAData"*pathmark*"5204040_Household_Financial_Account.xls", "Data1");
ASNANonFinFin = ExcelReaders.readxlsheet("."*pathmark*"data"*pathmark*"ASNAData"*pathmark*"5204019_NonFin_Corp_Financial_Account.xls", "Data1");
ASNAFinFin = ExcelReaders.readxlsheet("."*pathmark*"data"*pathmark*"ASNAData"*pathmark*"5204027_Fin_Corp_Financial_Account.xls", "Data1");
ASNAGovFin = ExcelReaders.readxlsheet("."*pathmark*"data"*pathmark*"ASNAData"*pathmark*"5204033_GenGov_Financial_Account.xls", "Data1");
ASNAExtFin = ExcelReaders.readxlsheet("."*pathmark*"data"*pathmark*"ASNAData"*pathmark*"5204044_External_Financial_Account.xls", "Data1");
#initialising table
#make no spaces
table16bNameRow = ["Acquisitionoffinancialassets-MonetarygoldandSDRs;","Acquisitionoffinancialassets-Currencyanddeposits;",
"Acquisitionoffinancialassets-Billsofexchange;","Acquisitionoffinancialassets-Onenamepaper–issuedinAustralia;",
"Acquisitionoffinancialassets-Onenamepaper–issuedoffshore;","Acquisitionoffinancialassets-Bondsetc.–issuedinAustralia;",
"Acquisitionoffinancialassets-Bondsetc.–issuedoffshore;","Acquisitionoffinancialassets-Derivatives;",
"Acquisitionoffinancialassets-Loansandplacements;","Acquisitionoffinancialassets-Sharesandotherequity;",
"Acquisitionoffinancialassets-Totalinsurancetechnicalreserves;","Acquisitionoffinancialassets-Otheraccounts",
"Incurrenceofliabilities-MonetarygoldandSDRs;","Incurrenceofliabilities-Currencyanddeposits;",
"Incurrenceofliabilities-Billsofexchange;","Incurrenceofliabilities-Onenamepaper–issuedinAustralia;",
"Incurrenceofliabilities-Onenamepaper–issuedoffshore;","Incurrenceofliabilities-Bondsetc.–issuedinAustralia;",
"Incurrenceofliabilities-Bondsetc.–issuedoffshore;","Incurrenceofliabilities-Derivatives;",
"Incurrenceofliabilities-Loansandplacements;","Incurrenceofliabilities-Sharesandotherequity;",
"Incurrenceofliabilities-Totalinsurancetechnicalreserves;","Incurrenceofliabilities-Otheraccounts"];
table16bDataCol = [ASNAHouseFin, ASNANonFinFin, ASNAFinFin, ASNAGovFin, ASNAExtFin];
table16bRowDict = Dict(table16bNameRow .=> [1:1:length(table16bNameRow);]);
table16bDataDict = Dict(tableName[1:length(tableName)-1] .=> table16bDataCol);
table16b = zeros(length(table16bNameRow),length(tableName));
table16b = NamedArray(table16b);
setnames!(table16b, string.(table16bNameRow), 1);
setnames!(table16b, string.(tableName), 2);
#filling in values
for i in tableName
    if i == "Total"
        for ring in table16bNameRow
            table16b[table16bRowDict[ring],tableDict[i]]=sum(table16b[table16bRowDict[ring],:]);
        end
    elseif i == "General Government"
        for ring in table16bNameRow
            if isempty(findall(x -> occursin("Generalgovernment-"*ring, x), filter.(x -> !isspace(x), string.(table16bDataDict[i][1,:]))))
                table16b[table16bRowDict[ring],tableDict[i]]=0.0;
            else
                table16b[table16bRowDict[ring],tableDict[i]]=1000*first(table16bDataDict[i][findall(x -> occursin("2019", x), string.(table16bDataDict[i][:,1])),findall(x -> occursin("Generalgovernment-"*ring, x), filter.(x -> !isspace(x), string.(table16bDataDict[i][1,:])))]);
            end
        end
    elseif i == "External"
        for ring in table16bNameRow
            key = ring
            if occursin("liabilities",ring) && occursin("issuedoffshore",ring)
                ring = replace(ring,"–issuedoffshore" => "")
            else
            end
            if isempty(findall(x -> occursin(ring, x), filter.(x -> !isspace(x), string.(table16bDataDict[i][1,:]))))
                table16b[table16bRowDict[key],tableDict[i]]=0.0;
            else
                table16b[table16bRowDict[key],tableDict[i]]=1000*first(table16bDataDict[i][findall(x -> occursin("2019", x), string.(table16bDataDict[i][:,1])),findall(x -> occursin(ring, x), filter.(x -> !isspace(x), string.(table16bDataDict[i][1,:])))]);
            end
        end
    else
        for ring in table16bNameRow
            if isempty(findall(x -> occursin(ring, x), filter.(x -> !isspace(x), string.(table16bDataDict[i][1,:]))))
                table16b[table16bRowDict[ring],tableDict[i]]=0.0;
            else
                table16b[table16bRowDict[ring],tableDict[i]]=1000*first(table16bDataDict[i][findall(x -> occursin("2019", x), string.(table16bDataDict[i][:,1])),findall(x -> occursin(ring, x), filter.(x -> !isspace(x), string.(table16bDataDict[i][1,:])))]);
            end
        end
    end
end

table16bAqcRow = findall(x -> occursin("Acquisition", x), string.(table16bNameRow));
table16bLiaRow = findall(x -> occursin("liabilities", x), string.(table16bNameRow));
table16bAqcTotal = sum(eachrow(table16b[table16bAqcRow,:]));
table16bLiaTotal = sum(eachrow(table16b[table16bLiaRow,:]));
table16SurplusOrDeficit = table16bAqcTotal - table16bLiaTotal;

#initialising table
#make no spaces
eAndONameRow = ["Neterrorsandomissions;"];
eAndODataCol = [ASNAHouseFin, ASNANonFinFin, ASNAFinFin, ASNAGovFin, ASNAExtFin];
eAndORowDict = Dict(eAndONameRow .=> [1:1:length(eAndONameRow);]);
eAndODataDict = Dict(tableName[1:length(tableName)-1] .=> eAndODataCol);
eAndO = zeros(length(eAndONameRow),length(tableName));
#filling in values
for i in tableName
    if i == "Total"
        for ring in eAndONameRow
            eAndO[eAndORowDict[ring],tableDict[i]]=sum(eAndO[eAndORowDict[ring],:]);
        end
    elseif i == "General Government"
        for ring in eAndONameRow
            if isempty(findall(x -> occursin("Generalgovernment-"*ring, x), filter.(x -> !isspace(x), string.(eAndODataDict[i][1,:]))))
                eAndO[eAndORowDict[ring],tableDict[i]]=0.0;
            else
                eAndO[eAndORowDict[ring],tableDict[i]]=1000*first(eAndODataDict[i][findall(x -> occursin("2019", x), string.(eAndODataDict[i][:,1])),findall(x -> occursin("Generalgovernment-"*ring, x), filter.(x -> !isspace(x), string.(eAndODataDict[i][1,:])))]);
            end
        end
    else
        for ring in eAndONameRow
            if isempty(findall(x -> occursin(ring, x), filter.(x -> !isspace(x), string.(eAndODataDict[i][1,:]))))
                eAndO[eAndORowDict[ring],tableDict[i]]=0.0;
            else
                eAndO[eAndORowDict[ring],tableDict[i]]=1000*first(eAndODataDict[i][findall(x -> occursin("2019", x), string.(eAndODataDict[i][:,1])),findall(x -> occursin(ring, x), filter.(x -> !isspace(x), string.(eAndODataDict[i][1,:])))]);
            end
        end
    end
end

table17TableNames = ["MonetarygoldandSDRs;","Currencyanddeposits;","Billsofexchange;",
"Onenamepaper–issuedinAustralia;","Onenamepaper–issuedoffshore;","Bondsetc.–issuedinAustralia;",
"Bondsetc.–issuedoffshore;","Derivatives;","Loansandplacements;","Sharesandotherequity;",
"Totalinsurancetechnicalreserves;","Otheraccounts"];

function myfind(c)
    a = similar(c, Int)
    count = 1
    @inbounds for i in eachindex(c)
        a[count] = i
        count += (c[i] != zero(eltype(c)))
    end
    return resize!(a, count-1)
end
table17=zeros(length(tableName), length(tableName))
for ring in [1:1:length(table17TableNames);]
    println("this round is number ", ring)
    table17=zeros(length(tableName), length(tableName));
    table17[length(tableName),:] = table16b[table16bRowDict["Acquisitionoffinancialassets-"*table17TableNames[ring]],:];
    table17[:,length(tableName)] = table16b[table16bRowDict["Incurrenceofliabilities-"*table17TableNames[ring]],:];
    if sum(table17[length(tableName),1:(length(tableName)-1)]) == sum(table17[1:(length(tableName)-1),length(tableName)])
    else
        table17RowTot = sum(table17[length(tableName),1:(length(tableName)-1)]);
        table17ColTot = sum(table17[1:(length(tableName)-1),length(tableName)]);
        table17TotAv = (table17RowTot + table17ColTot)/2;
        if 0.98*abs(table17RowTot)<abs(table17TotAv)<1.04*abs(table17RowTot)
        else
            error("Large discrepancy in row and column total sums in table 17."*string(ring))
        end
        table17[length(tableName),length(tableName)] = table17TotAv;
        for i in 1:(length(tableName)-1);
            table17[length(tableName),i] = table17[length(tableName),i]*table17TotAv/table17RowTot;
            table17[i,length(tableName)] = table17[i,length(tableName)]*table17TotAv/table17ColTot;
        end
    end
    nonZerosRow = myfind(table17[:,length(tableName)]);
    nonZerosCol = myfind(table17[length(tableName),:]);
    if length(nonZerosRow)-1 == 1;
        table17[nonZerosRow[1],nonZerosCol] = table17[tableDict["Total"],nonZerosCol];
    elseif length(nonZerosCol)-1 == 1;
        table17[nonZerosRow,nonZerosCol[1]] = table17[nonZerosRow,tableDict["Total"]];
    else
        mod17 = Model(Ipopt.Optimizer);
        @variable(mod17, x[1:(length(nonZerosRow)-1), 1:(length(nonZerosCol)-1)]);
        @NLobjective(mod17, Min, sum((x[i,j]) ^ 2 for i in 1:(length(nonZerosRow)-1), j in 1:(length(nonZerosCol)-1)));
        for i in 1:length(nonZerosRow)-1;
            @constraint(mod17, sum(x[i,:]) == table17[nonZerosRow[i],tableDict["Total"]]);
        end
        for i in 1:(length(nonZerosCol)-1);
            @constraint(mod17, sum(x[:,i]) == table17[tableDict["Total"],nonZerosCol[i]]);
        end
        if tableDict["External"] in nonZerosRow && tableDict["External"] in nonZerosCol
            @constraint(mod17, x[length(nonZerosRow)-1,length(nonZerosCol)-1] == 0);
        else
        end
        optimize!(mod17);
        table17[nonZerosRow[1:length(nonZerosRow)-1], nonZerosCol[1:length(nonZerosCol)-1]] = value.(x);
    end
    @eval $(Symbol(:table17_, ring)) = NamedArray(table17);
    setnames!((@eval $(Symbol(:table17_, ring))), string.(tableName), 1);
    setnames!((@eval $(Symbol(:table17_, ring))), string.(tableName), 2);
end
table19 = zeros(length(tableName),length(tableName))
for ring in [1:1:length(table17TableNames);]
    table19 = table19 + (@eval $(Symbol(:table17_, ring)));
end
table19 = NamedArray(table19);
setnames!(table19, string.(tableName), 1);
setnames!(table19, string.(tableName), 2);
table19TotRow = table19[:,tableDict["Total"]]+table16a[table16aRowDict["Netlending(+)/netborrowing(-);"],:];
table19EAndO = table19[tableDict["Total"],:]-table19TotRow;

#final bits and pieces to fill in gaps
ASNAIncomeFromGDP = ExcelReaders.readxlsheet("."*pathmark*"data"*pathmark*"ASNAData"*pathmark*"5204006_Income_from_GDP.xls", "Data1");
productionToHouseholds = (first(ASNAIncomeFromGDP[findall(x -> occursin("2019", x), string.(ASNAIncomeFromGDP[:,1])),findall(x -> occursin("Compensationofemployees;", x), filter.(x -> !isspace(x), string.(ASNAIncomeFromGDP[1,:])))])
-first(ASNAExtInc[findall(x -> occursin("2019", x), string.(ASNAExtInc[:,1])),findall(x -> occursin("receivable-Compensationofemployees;", x), filter.(x -> !isspace(x), string.(ASNAExtInc[1,:])))])
+first(ASNAIncomeFromGDP[findall(x -> occursin("2019", x), string.(ASNAIncomeFromGDP[:,1])),findall(x -> occursin("Dwellingsownedbypersons;Grossoperatingsurplus;", x), filter.(x -> !isspace(x), string.(ASNAIncomeFromGDP[1,:])))])
+first(ASNAIncomeFromGDP[findall(x -> occursin("2019", x), string.(ASNAIncomeFromGDP[:,1])),findall(x -> occursin("Grossmixedincome;", x), filter.(x -> !isspace(x), string.(ASNAIncomeFromGDP[1,:])))]));

#current to current transfers
#from households
currentToCurrent = transpose(table6 + table7 + table8 + table9 + table10 + table11 + table12 + table13);
currentToCurrent[tableDict["General Government"],tableDict["Households"]] = (currentToCurrent[tableDict["General Government"],tableDict["Households"]] + first(ASNAHouseInc[findall(x -> occursin("2019", x), string.(ASNAHouseInc[:,1])),findall(x -> occursin("payable-Incometax;", x), filter.(x -> !isspace(x), string.(ASNAHouseInc[1,:])))])
+ first(ASNAHouseInc[findall(x -> occursin("2019", x), string.(ASNAHouseInc[:,1])),findall(x -> occursin("payable-Othercurrenttaxes", x), filter.(x -> !isspace(x), string.(ASNAHouseInc[1,:])))]) + IO[IORowDict["P3"], IOColDict["Q1"]]);
currentToCurrent[tableDict["External"],tableDict["Households"]] = currentToCurrent[tableDict["External"],tableDict["Households"]] + sum(IO[IORowDict["P5"]:IORowDict["P6"],IOColDict["Q1"]]);
#from non-financial
currentToCurrent[tableDict["General Government"],tableDict["Non-Financial Corporations"]] = currentToCurrent[tableDict["General Government"],tableDict["Non-Financial Corporations"]] + first(ASNANonFinInc[findall(x -> occursin("2019", x), string.(ASNANonFinInc[:,1])),findall(x -> occursin("Secondaryincomepayable-Currenttaxesonincome,wealth,etc-Total;", x), filter.(x -> !isspace(x), string.(ASNANonFinInc[1,:])))]);
#from financial
currentToCurrent[tableDict["General Government"],tableDict["Financial Corporations"]] = currentToCurrent[tableDict["General Government"],tableDict["Financial Corporations"]] + first(ASNAFinInc[findall(x -> occursin("2019", x), string.(ASNAFinInc[:,1])),findall(x -> occursin("Secondaryincomepayable-Currenttaxesonincome,wealth,etc-Total;", x), filter.(x -> !isspace(x), string.(ASNAFinInc[1,:])))]);
#from general Government
currentToCurrent[tableDict["External"],tableDict["General Government"]] = currentToCurrent[tableDict["External"],tableDict["General Government"]] + sum(IO[IORowDict["P5"]:IORowDict["P6"],IOColDict["Q2"]]);
#from external
currentToCurrent[tableDict["General Government"],tableDict["External"]] = currentToCurrent[tableDict["General Government"],tableDict["External"]] + sum(IO[IORowDict["P3"],IOColDict["Q7"]]) + first(ASNAExtInc[findall(x -> occursin("2019", x), string.(ASNAExtInc[:,1])),findall(x -> occursin("payable-Currenttaxesonincome,wealth,etc;", x), filter.(x -> !isspace(x), string.(ASNAExtInc[1,:])))]);
currentToCurrent[tableDict["Households"],tableDict["External"]] = currentToCurrent[tableDict["Households"],tableDict["External"]] + first(ASNAExtInc[findall(x -> occursin("2019", x), string.(ASNAExtInc[:,1])),findall(x -> occursin("payable-Compensationofemployees;", x), filter.(x -> !isspace(x), string.(ASNAExtInc[1,:])))]);
currentToCurrent[tableDict["External"],tableDict["External"]] = currentToCurrent[tableDict["External"],tableDict["External"]] + sum(IO[IORowDict["P6"],IOColDict["Q7"]]);

#capital to capital transfers
capitalToCapital = transpose(table14) + transpose(table15) + table19

ASNAConsFixCap = ExcelReaders.readxlsheet("."*pathmark*"data"*pathmark*"ASNAData"*pathmark*"5204047_Cons_Fixed_Capital_By_Industry.xls", "Data1");

samName = ["Production Activities", "Factors of Production", "Households Current Account", "Non-Financial Corporations Current Account",
"Financial Corporations Current Account", "Government Current Account", "Foreigners Current Account",
"Households Capital Account", "Non-Financial Corporations Capital Account", "Financial Corporations Capital Account",
"Government Capital Account", "Foreigners Capital Account", "Net Errors and Ommisions", "Total"];
samDict = Dict(samName .=> [1:1:length(samName);]);
sam = zeros(length(samName),length(samName));
sam[samDict["Production Activities"],samDict["Production Activities"]] = IO[IORowDict["T1"],IOColDict["T4"]];
sam[samDict["Production Activities"],samDict["Households Current Account"]] = IO[IORowDict["T1"],IOColDict["Q1"]];
sam[samDict["Production Activities"],samDict["Government Current Account"]] = IO[IORowDict["T1"],IOColDict["Q2"]];
sam[samDict["Production Activities"],samDict["Foreigners Current Account"]] = IO[IORowDict["T1"],IOColDict["Q7"]];
sam[samDict["Production Activities"],samDict["Households Capital Account"]:samDict["Government Capital Account"]] = table5c[table5cRowDict["Domestic Commodities"],table5cColDict["Households"]:table5cColDict["General Government"]];
sam[samDict["Factors of Production"], samDict["Production Activities"]] = sum(IO[IORowDict["P1"]:IORowDict["P2"],IOColDict["T4"]]);
sam[samDict["Households Current Account"], samDict["Factors of Production"]] = productionToHouseholds;
sam[samDict["Non-Financial Corporations Current Account"], samDict["Factors of Production"]] = first(ASNANonFinInc[findall(x -> occursin("2019", x), string.(ASNANonFinInc[:,1])),findall(x -> occursin("Grossoperatingsurplus;", x), filter.(x -> !isspace(x), string.(ASNANonFinInc[1,:])))]);
sam[samDict["Financial Corporations Current Account"], samDict["Factors of Production"]] = first(ASNAFinInc[findall(x -> occursin("2019", x), string.(ASNAFinInc[:,1])),findall(x -> occursin("Grossoperatingsurplus;", x), filter.(x -> !isspace(x), string.(ASNAFinInc[1,:])))]);
sam[samDict["Government Current Account"], samDict["Factors of Production"]] = first(ASNAGovInc[findall(x -> occursin("2019", x), string.(ASNAGovInc[:,1])),findall(x -> occursin("Generalgovernment;Consumptionoffixedcapital;", x), filter.(x -> !isspace(x), string.(ASNAGovInc[1,:])))]);
sam[samDict["Foreigners Current Account"], samDict["Factors of Production"]] = first(ASNAExtInc[findall(x -> occursin("2019", x), string.(ASNAExtInc[:,1])),findall(x -> occursin("receivable-Compensationofemployees;", x), filter.(x -> !isspace(x), string.(ASNAExtInc[1,:])))]);
sam[samDict["Government Current Account"], samDict["Production Activities"]] = sum(IO[IORowDict["P3"]:IORowDict["P4"],IOColDict["T4"]]);
sam[samDict["Foreigners Current Account"], samDict["Production Activities"]] = sum(IO[IORowDict["P5"]:IORowDict["P6"],IOColDict["T4"]]);
sam[samDict["Households Current Account"]:samDict["Foreigners Current Account"],samDict["Households Current Account"]:samDict["Foreigners Current Account"]]=currentToCurrent[1:length(tableName)-1,1:length(tableName)-1];
sam[samDict["Households Capital Account"], samDict["Households Current Account"]] = (first(ASNAHouseInc[findall(x -> occursin("2019", x), string.(ASNAHouseInc[:,1])),findall(x -> occursin("Netsaving;", x), filter.(x -> !isspace(x), string.(ASNAHouseInc[1,:])))])
+first(ASNAConsFixCap[findall(x -> occursin("2019", x), string.(ASNAConsFixCap[:,1])),findall(x -> occursin("Households;Dwellingsownedbypersons;", x), filter.(x -> !isspace(x), string.(ASNAConsFixCap[1,:])))])
+first(ASNAConsFixCap[findall(x -> occursin("2019", x), string.(ASNAConsFixCap[:,1])),findall(x -> occursin("Households;Other;", x), filter.(x -> !isspace(x), string.(ASNAConsFixCap[1,:])))]));
sam[samDict["Non-Financial Corporations Capital Account"], samDict["Non-Financial Corporations Current Account"]] = first(ASNANonFinInc[findall(x -> occursin("2019", x), string.(ASNANonFinInc[:,1])),findall(x -> occursin("Consumptionoffixedcapital;", x), filter.(x -> !isspace(x), string.(ASNANonFinInc[1,:])))])+first(ASNANonFinInc[findall(x -> occursin("2019", x), string.(ASNANonFinInc[:,1])),findall(x -> occursin("Netsaving;", x), filter.(x -> !isspace(x), string.(ASNANonFinInc[1,:])))]);
sam[samDict["Financial Corporations Capital Account"], samDict["Financial Corporations Current Account"]] = first(ASNAFinInc[findall(x -> occursin("2019", x), string.(ASNAFinInc[:,1])),findall(x -> occursin("Consumptionoffixedcapital;", x), filter.(x -> !isspace(x), string.(ASNAFinInc[1,:])))])+first(ASNAFinInc[findall(x -> occursin("2019", x), string.(ASNAFinInc[:,1])),findall(x -> occursin("Netsaving;", x), filter.(x -> !isspace(x), string.(ASNAFinInc[1,:])))]);
sam[samDict["Government Capital Account"], samDict["Government Current Account"]] = first(ASNAGovInc[findall(x -> occursin("2019", x), string.(ASNAGovInc[:,1])),findall(x -> occursin("Generalgovernment;Consumptionoffixedcapital;", x), filter.(x -> !isspace(x), string.(ASNAGovInc[1,:])))])+first(ASNAGovInc[findall(x -> occursin("2019", x), string.(ASNAGovInc[:,1])),findall(x -> occursin("Generalgovernment;Netsaving;", x), filter.(x -> !isspace(x), string.(ASNAGovInc[1,:])))]);
sam[samDict["Foreigners Capital Account"], samDict["Foreigners Current Account"]] = first(ASNAExtInc[findall(x -> occursin("2019", x), string.(ASNAExtInc[:,1])),findall(x -> occursin("Balanceonexternalincomeaccount;", x), filter.(x -> !isspace(x), string.(ASNAExtInc[1,:])))]);

sam[samDict["Households Capital Account"]:samDict["Foreigners Capital Account"],samDict["Households Capital Account"]:samDict["Foreigners Capital Account"]]= capitalToCapital[1:length(tableDict)-1,1:length(tableDict)-1];
sam[samDict["Government Current Account"],samDict["Households Capital Account"]:samDict["Government Capital Account"]] = table5c[table5cRowDict["Taxes less subsidies on products"],1:length(table5cNameCol)-1] + table5c[table5cRowDict["Other taxes less subsidies on investment"],1:length(table5cNameCol)-1];
sam[samDict["Foreigners Current Account"],samDict["Households Capital Account"]:samDict["Government Capital Account"]] = table5c[table5cRowDict["Imported Commodities"],1:length(table5cNameCol)-1];
sam[samDict["Net Errors and Ommisions"],samDict["Households Capital Account"]:samDict["Foreigners Capital Account"]] = -1 .* table19EAndO[1:length(tableName)-1] + sum(table19EAndO[1:length(tableName)-1]).*abs.(table19EAndO[1:length(tableName)-1])./sum(abs.(table19EAndO[1:length(tableName)-1]));

for i in 1:length(samName)-1
    sam[i,samDict["Total"]]=sum(sam[i,1:length(samName)-1]);
    sam[samDict["Total"],i]=sum(sam[1:length(samName)-1,i]);
end

differencesCol = ["Diff" (sam[14,:]-sam[:,14])' 0.0];
differencesRow = sam[:,14]-sam[14,:];

samShortName = [ "Prod", "Fctr", "HhldCrt", "NfinCrt", "Fin_Crt", "GvntCrt", "ExtlCrt", "HhldCpl", "NFinCpl", "Fin_Cpl", "GvntCpl", "ExtlCpl","ErrOmm", "Total"];
samShortName2 = ["Titles" "Prod" "Fctr" "HhldCrt" "NfinCrt" "Fin_Crt" "GvntCrt" "ExtlCrt" "HhldCpl" "NFinCpl" "Fin_Cpl" "GvntCpl" "ExtlCpl" "ErrOmm" "Total" "Diff"];

sam = NamedArray(sam);
setnames!(sam, samShortName, 1);
setnames!(sam, samShortName, 2);

withTitlesAndDiff=[samShortName2 ; [samShortName sam differencesRow] ; differencesCol];

writedlm("tableCSV"*pathmark*"samWithTitlesAndDiff.csv", withTitlesAndDiff, ',');

#remove the negative values, put them into the correct spot
#create matrix on top to store Changes
samNegs = zeros(length(samName)-2,length(samName)-2);
for i in 1:length(samName)-2
    for j in 1:length(samName)-2
        if sam[i,j] < 0
            samNegs[j,i] = abs(sam[i,j]);
            sam[i,j]=0;
        else
        end
    end
end

sam[1:length(samName)-2,1:length(samName)-2]=sam[1:length(samName)-2,1:length(samName)-2]+samNegs;

function posfind(c)
    a = similar(c, Int)
    count = 1
    @inbounds for i in eachindex(c)
        a[count] = i
        count += (c[i] > 0.0)
    end
    return resize!(a, count-1)
end

#RAS the differences into the sam
rasDiffCol = sam[:,14]-sam[14,:];
rasDiffRow = -1 .* rasDiffCol;
posRow = posfind(rasDiffRow);
posCol = posfind(rasDiffCol);
rasDiffCol = rasDiffCol[posCol];
rasDiffRow = rasDiffRow[posRow];

#initialize matrix
modDiff = Model(Ipopt.Optimizer);
@variable(modDiff, x[1:length(rasDiffRow), 1:length(rasDiffCol)]);
@NLobjective(modDiff, Min, sum((x[i,j]) ^ 2 for i in 1:(length(rasDiffRow)), j in 1:(length(rasDiffCol))));
for i in 1:length(rasDiffRow);
    @constraint(modDiff, sum(x[i,:]) == rasDiffRow[i]);
    for j in 1:length(rasDiffCol);
        if sam[posRow[i],posCol[j]] == 0.0;
            @constraint(modDiff, x[i,j] == 0.0);
        else
        end
    end
end
for i in 1:length(rasDiffCol);
    @constraint(modDiff, sum(x[:,i]) == rasDiffCol[i]);
end

optimize!(modDiff);

sam[posRow,posCol]=sam[posRow,posCol]+value.(x);

for i in 1:length(samName)-1
    sam[i,samDict["Total"]]=sum(sam[i,1:length(samName)-1]);
    sam[samDict["Total"],i]=sum(sam[1:length(samName)-1,i]);
end

sam[length(samName),length(samName)]=sum(sam[length(samName),:]);

sam = round.(sam, digits = 9);

tablecsv(sam, "sam");

#make full SAM with all 114 production activities (splitting up factors of production was quite challenging, that could be returned to)

samFull = zeros(length(sam)+114,length(sam)+114);
samFullName = vcat(IOSource[first(findall(x -> occursin("101", x), string.(IOSource[:,1]))):first(findall(x -> occursin("9502", x), string.(IOSource[:,1]))),1],samName);
samFullDict = Dict(samFullName .=> [1:1:length(samFullName);])
samFull[1:114,1:114] = IOSource[first(findall(x -> occursin("101", x), string.(IOSource[:,1]))):first(findall(x -> occursin("9502", x), string.(IOSource[:,1]))), first(findall(x -> occursin("101", x), string.(IOSource[1,:]))): first(findall(x -> occursin("9502", x), string.(IOSource[1,:])))];
samFull[115,1:114] = IOSource[first(findall(x -> occursin("P1", x), string.(IOSource[:,1]))), first(findall(x -> occursin("101", x), string.(IOSource[1,:]))): first(findall(x -> occursin("9502", x), string.(IOSource[1,:])))]+
IOSource[first(findall(x -> occursin("P2", x), string.(IOSource[:,1]))), first(findall(x -> occursin("101", x), string.(IOSource[1,:]))): first(findall(x -> occursin("9502", x), string.(IOSource[1,:])))];
samFull[samFullDict["Government Current Account"],1:114]=IOSource[first(findall(x -> occursin("P3", x), string.(IOSource[:,1]))), first(findall(x -> occursin("101", x), string.(IOSource[1,:]))): first(findall(x -> occursin("9502", x), string.(IOSource[1,:])))]+
IOSource[first(findall(x -> occursin("P4", x), string.(IOSource[:,1]))), first(findall(x -> occursin("101", x), string.(IOSource[1,:]))): first(findall(x -> occursin("9502", x), string.(IOSource[1,:])))];
samFull[samFullDict["Foreigners Current Account"],1:114]=IOSource[first(findall(x -> occursin("P5", x), string.(IOSource[:,1]))), first(findall(x -> occursin("101", x), string.(IOSource[1,:]))): first(findall(x -> occursin("9502", x), string.(IOSource[1,:])))]+
IOSource[first(findall(x -> occursin("P6", x), string.(IOSource[:,1]))), first(findall(x -> occursin("101", x), string.(IOSource[1,:]))): first(findall(x -> occursin("9502", x), string.(IOSource[1,:])))];
samFull[1:114,samFullDict["Households Current Account"]] = IOSource[first(findall(x -> occursin("101", x), string.(IOSource[:,1]))):first(findall(x -> occursin("9502", x), string.(IOSource[:,1]))), first(findall(x -> occursin("Households ;\n Final Consumption Expenditure", x), string.(IOSource[2,:])))];
samFull[1:114,samFullDict["Government Current Account"]] = IOSource[first(findall(x -> occursin("101", x), string.(IOSource[:,1]))):first(findall(x -> occursin("9502", x), string.(IOSource[:,1]))), first(findall(x -> occursin("General Government ;\n Final Consumption Expenditure", x), string.(IOSource[2,:])))];
samFull[1:114,samFullDict["Foreigners Current Account"]] = IOSource[first(findall(x -> occursin("101", x), string.(IOSource[:,1]))):first(findall(x -> occursin("9502", x), string.(IOSource[:,1]))), first(findall(x -> occursin("Exports of Goods and Services", x), string.(IOSource[2,:])))];

#split up IO table to make a full sam
#do the optimizing as in dual 2by2 example by patrick

#=convert dataframe to dictionary
function increment!( d::Dict{S, T}, k::S, i::T) where {T<:Real, S<:Any}
    if haskey(d, k)
        d[k] += i
    else
        d[k] = i
    end
end
increment!(d::Dict{S, T}, k::S ) where {T<:Real, S<:Any} = increment!( d, k, one(T))

function df2dict( df::DataFrame, key_col::Symbol, val_col::Symbol=:null)
    keytype = typeof(df[1,key_col])
    if val_col == :null
        valtype = Int
    else
        valtype = typeof(df[1,val_col])
    end
    D = Dict{keytype, valtype}()
    for i=1:size(df,1)
        if !ismissing(df[i,key_col])
            if val_col == :null
                increment!( D, df[i,key_col] )
            elseif valtype <: Real
                increment!( D, df[i,key_col], df[i,val_col] )
            else
                if haskey(D, df[i,key_col])
                    @warn("non-unique entry: $(df[i,key_col])")
                else
                    D[df[i,key_col]] = df[i,val_col]
                end
            end
        end
    end
    return D
end
df[!, "IOCode"]=IOcode
insertcols!(df, 2, :name => vector)
D = df2dict(df, :IOCode, :x3)
=#

#notes
#implement loop and more ras-ing - take out spaces in all find functions
#generalise the ASNAYearRow to just search every time
#rename tables from numbers to abbreviations#new filepath test
#fix 101 to 0101 in IOSource IO code collumn and row

max_stars_repo_name: uq-aibe/sam

max_issues_repo_name: uq-aibe/sam

max_forks_repo_name: uq-aibe/sam


text:
using GLMakie,CairoMakie
GLMakie.activate!()

function scatters_in_3D()
    seed!(123)
    xyz = randn(10, 3)
    x, y, z = xyz[:, 1], xyz[:, 2], xyz[:, 3]
    fig = Figure(resolution=(1600, 400))
    ax1 = Axis3(fig[1, 1]; aspect=(1, 1, 1), perspectiveness=0.5)
    ax2 = Axis3(fig[1, 2]; aspect=(1, 1, 1), perspectiveness=0.5)
    ax3 = Axis3(fig[1, 3]; aspect=:data, perspectiveness=0.5)
    scatter!(ax1, x, y, z; markersize=50)
    meshscatter!(ax2, x, y, z; markersize=0.25)
    hm = meshscatter!(ax3, x, y, z; markersize=0.25,
        marker=FRect3D(Vec3f(0), Vec3f(1)), color=1:size(xyz)[2],
        colormap=:plasma, transparency=false)
    Colorbar(fig[1, 4], hm, label="values", height=Relative(0.5))
    fig
end
scatters_in_3D()

ext: jl

max_stars_repo_name: BioMedAI/Pluto-Note

max_issues_repo_name: BioMedAI/Pluto-Note

max_forks_repo_name: BioMedAI/Pluto-Note


text:
using Revise
using JuMP
using Ipopt
using PyPlot
using DataFrames
using StructArrays
using CSV
using LinearAlgebra
using KNITRO
using Gurobi

include("MPPDCs/MPPDC_quantity.jl") #The MPPDC model with quantity bidding
include("MILPs/MILP_quantity.jl") #The MILP model with quantity bidding
include("MPPDCs/NLPs/NLP_dual_MPPDC_quantity.jl") #The NLP model with quantity bidding

###############DECLARE PARTICIPATING AGENTS###########
rel_path = "DATA"
filepath = joinpath(@__DIR__, rel_path)
cd(filepath)
generators_dat = CSV.read("generators_data.csv", header=true, delim= ";")
MC_g = generators_dat[1,2:8]
MC_g_str = generators_dat[1,9]
#CAPACITIES OF GENERATORS
capacity_g = generators_dat[2,2:9]
capacity_g = values(capacity_g) ./ 100
capacity_g_nonstr = values(capacity_g[1:7])
capacity_g_str = values(capacity_g[8])

#RAMPING LIMIT OF THE STRATEGIC GENERATOR:
Rᴳ = 1.0

res_dat = CSV.read("wind_PV.csv", header=true, delim= ";")
eight_res  = res_dat[1:3:end,:] #selecting every 3rd of the res_dat
renewable_nonstr = ["PV"]
renewable_str = ["WIND"]
capacity_re_nonstr = eight_res.PV/100
capacity_re_str =  eight_res.Wind/100

#LOADING THE DEMAND
load_dat = CSV.read("load.csv",header=true, delim= ";")
load = load_dat.load[1:3:end,:]/100

eff_ch = eff_dch = generators_dat.ESS[end]
energy_max = generators_dat.ESS[end-2]/100
power_max = generators_dat.ESS[end-4]/100

MC_re_str = MC_re_nonstr = 0 #Zero marginal cost

str_stor = ["ESS_STR"]
str_gens = ["GEN_STR", "WIND_STR"]
str_agents = vcat(str_stor, str_gens)

generators = names(MC_g)
renewables = [:PV]

non_str_gens = vcat(generators)
non_str_res = vcat(renewables)
non_str_agents = vcat(generators, renewables)

#NOTE:ATM There is no non-str ESS involved in the pool
#MARKET PRICE LIMITS
p_cap = 100
p_floor = 0

function t_interval(lenght, resolution)
    t_steps = [i for i in range(1, stop = lenght, step = resolution)]
    return t_steps
end

t_steps = t_interval(8,1)

function compute_profits()
    lambda = value.(NLP[:λ][:])
    #NON-STRATEGIC GENERATORS
    gen_1 = value.(NLP[:g]["GEN1",:])
    gen_2 = value.(NLP[:g]["GEN2",:])
    gen_3 = value.(NLP[:g]["GEN3",:])
    gen_5 = value.(NLP[:g]["GEN5",:])
    gen_6 = value.(NLP[:g]["GEN6",:])
    gen_7 = value.(NLP[:g]["GEN7",:])
    gen_8 = value.(NLP[:g]["GEN8",:])

    #DEMAND
    d = value.(NLP[:d])
    #NON-STRATEGIC RES
    res_PV = value.(NLP[:w][:PV,:])
    #DISPATCHED QUANTITIES FROMS STR AGENTS
    dch = value.(NLP[:dch_st])
    ch = value.(NLP[:ch_st])
    gen = value.(NLP[:g_st])
    res = value.(NLP[:w_st])
    #CALCULATING PROFITS
    earned_ESS = sum(dch[t] * lambda[t] for t in t_steps)
    paid_ESS = sum(ch[t] * lambda[t] for t in t_steps)
    profit_ESS =  earned_ESS - paid_ESS
    # println("PROFIT_ESS: ", profit_ESS)

    earned_GEN = sum(gen[t] * lambda[t] for t in t_steps)
    paid_GEN = sum(dat_str.GEN_STR.MC_g[t] * gen[t] for t in t_steps)
    profit_GEN =  earned_GEN - paid_GEN
    # println("PROFIT_GEN: ", profit_GEN)

    earned_RES = sum(res[t] * lambda[t] for t in t_steps)
    paid_RES = sum(dat_str.WIND_STR.MC_g[t] * res[t] for t in t_steps)
    profit_RES =  earned_RES - paid_RES
    # println("PROFIT_RES: ", profit_RES)

    return profit_ESS, profit_GEN, profit_RES
end
function compute_profits_MPPDC()
    lambda = value.(MPPDC[:λ][:])
    #NON-STRATEGIC GENERATORS
    gen_1 = value.(MPPDC[:g]["GEN1",:])
    gen_2 = value.(MPPDC[:g]["GEN2",:])
    gen_3 = value.(MPPDC[:g]["GEN3",:])
    gen_5 = value.(MPPDC[:g]["GEN5",:])
    gen_6 = value.(MPPDC[:g]["GEN6",:])
    gen_7 = value.(MPPDC[:g]["GEN7",:])
    gen_8 = value.(MPPDC[:g]["GEN8",:])

    #DEMAND
    d = value.(MPPDC[:d])
    #NON-STRATEGIC RES
    res_PV = value.(MPPDC[:w][:PV,:])
    #DISPATCHED QUANTITIES FROMS STR AGENTS
    dch = value.(MPPDC[:dch_st])
    ch = value.(MPPDC[:ch_st])
    gen = value.(MPPDC[:g_st])
    res = value.(MPPDC[:w_st])
    #CALCULATING PROFITS
    earned_ESS = sum(dch[t] * lambda[t] for t in t_steps)
    paid_ESS = sum(ch[t] * lambda[t] for t in t_steps)
    profit_ESS =  earned_ESS - paid_ESS
    # println("PROFIT_ESS: ", profit_ESS)

    earned_GEN = sum(gen[t] * lambda[t] for t in t_steps)
    paid_GEN = sum(dat_str.GEN_STR.MC_g[t] * gen[t] for t in t_steps)
    profit_GEN =  earned_GEN - paid_GEN
    # println("PROFIT_GEN: ", profit_GEN)

    earned_RES = sum(res[t] * lambda[t] for t in t_steps)
    paid_RES = sum(dat_str.WIND_STR.MC_g[t] * res[t] for t in t_steps)
    profit_RES =  earned_RES - paid_RES
    # println("PROFIT_RES: ", profit_RES)

    return profit_ESS, profit_GEN, profit_RES
end

###########INIT BID, TECH LIMIT DFs#########
function return_data()
    #DATA OF THE STR AGENTS, DIFFERENTIATE BETWEEN BIDS AND TRUE VALUES
    dat_str = DataFrame()
    for i in range(1, stop = length(str_gens)), l in range(1, stop = length(str_stor))
        #This needs to be revisited....
        if i ==  1 #if str agent is GEN
            dat_str[i] =  StructArray(
            P_max = ones(length(t_steps)) * capacity_g_str,
            G = zeros(length(t_steps)),
            c_g = zeros(length(t_steps)),
            MC_g = ones(length(t_steps)) * MC_g_str)
        elseif i == 2 #if str agent is RE
            dat_str[i] =  StructArray(
            P_max = capacity_re_str,
            G = zeros(length(t_steps)),
            c_g = zeros(length(t_steps)),
            MC_g = zeros(length(t_steps)))
        end

        dat_str[i + l] =  StructArray(
        E_max = ones(length(t_steps)) * energy_max,
        P_max = ones(length(t_steps)) * power_max,
        DCH = zeros(length(t_steps)),
        MC_dch = ones(length(t_steps)) * p_cap,
        c_dch = zeros(length(t_steps)),
        CH = zeros(length(t_steps)),
        c_ch = zeros(length(t_steps)),
        MC_ch = ones(length(t_steps)) * p_floor)
    end

    #NOTE: For non-str agents set the quantity bid equal to the maximum capacity
    #NON STR AGENTS THE BIDS ARE EQUAL TO THE TRUE VALUES
    dat_nonstr = DataFrame()
    for i in range(1, stop = length(non_str_gens)), k in range(1, stop = length(non_str_res))
        dat_nonstr[i] =  StructArray(
        G = ones(length(t_steps)) * capacity_g_nonstr[i],
        c_g = ones(length(t_steps)) * MC_g[1+i-1])
        dat_nonstr[i + k] = StructArray(
        G = capacity_re_nonstr,
        c_g = ones(length(t_steps)) * MC_re_nonstr #as zero marginal cost
        )
    end

    #NOW RENAME THEM
    names!(dat_str, vcat(Symbol.(str_gens), Symbol.(str_stor)))
    names!(dat_nonstr, Symbol.(non_str_agents))

    return dat_str, dat_nonstr
end

dat_str, dat_nonstr = return_data()
#
optimizer = (with_optimizer(Ipopt.Optimizer,
                                max_iter = 3000,
                                print_level = 0
                                # tol = 1.0e-16,
                                # acceptable_tol = 1.0e-10
                                ))

#SETTING THE CURRENT AGENT
#NOTE: ACCEPTED STRINGS: "ESS_STR", "GEN_STR", "WIND_STR"
current = "ESS_STR"
bidding = "quantity"
model_type = "MPEC"
model_type = "MPPDC"
model_type = "MILP"
objective_type = "dual" #Dual refers to the linearized objective functon, whereas primal is the original including bi-linear terms
flag_regularize = true

MPPDC = build_MPPDC_quantity(t_steps, non_str_gens,
non_str_res, current, optimizer, objective_type)

println("SOLUTION METHOD:
                ....Agent: $current
                ....Objective_type: $objective_type")

flag_sequanctial_red = true
#NOTE: MIND THE STEP SIZE
results_MPPDC_ESS = DataFrame(ϵ = [], iter = [], profit = [])
results_MPPDC_GEN = DataFrame(ϵ = [], iter = [], profit = [])
results_MPPDC_WIND = DataFrame(ϵ = [], iter = [], profit = [])

#FOR MPPDC THE SEQUANTIAL ALGORITHM
if flag_sequanctial_red == true
    #Rebuild the problem
    MPPDC = build_MPPDC_quantity(t_steps, non_str_gens,
    non_str_res, current, optimizer, objective_type)
    global iter_limit = 100
    global iterations = 0
    global tol = 0.001
    global epsilon_new = 1000
    JuMP.fix(MPPDC[:ϵ], epsilon_new; force = true)

    optimize!(MPPDC)
    objective_value(MPPDC)
    profit_ESS, profit_GEN, profit_RES = compute_profits_MPPDC()
    println("PROFITS:
                    PROFIT_ESS: $profit_ESS,
                    PROFIT_GEN: $profit_GEN,
                    PROFIT_RES: $profit_RES")
    global iterations += 1
    push!(results_MPPDC_ESS, [epsilon_new, iterations, profit_ESS])
    ############THE SMAL SCRIPT SNIPPET FOR UPDATING THE REGULARIZATION TERM
    while epsilon_new >= tol && iterations < iter_limit
        #########SETTING THE WARM START POINT TO THE PREVIOUS SOLUTION#########
        #if terminition status is local optimal, otherwise should not use the solution as starting point
        term_stat = termination_status(MPPDC)
        # if term_stat == MOI.OPTIMAL || term_stat == MOI.LOCALLY_SOLVED
        set_start_value.(all_variables(MPPDC), value.(all_variables(MPPDC)))
        ###################FIXING EPRSILONS TO THE REDUCED PARAMS##############
        global epsilon_new = value.(MPPDC[:ϵ])*0.1
        # else
        #     set_start_value.(all_variables(MPPDC), prev_vals)
        #     global epsilon_new = value.(MPPDC[:ϵ])*0.99
        # end
        fix(MPPDC[:ϵ], epsilon_new; force = true)
        println("\\
        ######### \\
        THE NEW RELAXATION TERMS: \\
                    ϵ: $epsilon_new and \\
                    #ITERATIONS: $iterations \\
                    #TERMINATION: $term_stat
        #########")
        ################OPTIMIZING############
        optimize!(MPPDC)
        profit_ESS, profit_GEN, profit_RES = compute_profits_MPPDC()
        println("PROFITS:
                        PROFIT_ESS: $profit_ESS,
                        PROFIT_GEN: $profit_GEN,
                        PROFIT_RES: $profit_RES")

        ##############CACLCULATING PROFIT############
        # profit = compute_profit()
        # println("PROFIT: ", profit)
        ############PROGRESSING IN TIME##############
        global iterations += 1
        push!(results_MPPDC_ESS, [epsilon_new, iterations, profit_ESS])
        global prev_vals = value.(all_variables(MPPDC))
    end
end

@show res_MPPDC = results_MPPDC_ESS.profit[end]

#FOR NLP
NLP = build_NLP_dual_quantity_mod(t_steps, non_str_gens,
non_str_res, str_agents, current, optimizer)

results_NLP_ESS = DataFrame(ϵ1 = [], ϵ2 = [], iter = [], profit = [])
results_NLP_GEN = DataFrame(ϵ1 = [], ϵ2 = [], iter = [], profit = [])
results_NLP_WIND = DataFrame(ϵ1 = [], ϵ2 = [], iter = [], profit = [])

if flag_sequanctial_red == true
    global iter_limit = 50
    global iterations = 0
    global tol = 0.001
    global epsilon1_new = 100
    global epsilon2_new = 100
    #Rebuild the model, otherwise there is some bias
    NLP = build_NLP_dual_quantity_mod(t_steps, non_str_gens,
    non_str_res, str_agents, current, optimizer)

    JuMP.fix(NLP[:ϵ1], epsilon1_new; force = true)
    JuMP.fix(NLP[:ϵ2], epsilon2_new; force = true)
    optimize!(NLP)
    objective_value(NLP)
    profit_ESS, profit_GEN, profit_RES = compute_profits()
    println("PROFITS:
                    PROFIT_ESS: $profit_ESS,
                    PROFIT_GEN: $profit_GEN,
                    PROFIT_RES: $profit_RES")

    global iterations += 1
    push!(results_NLP_ESS, [epsilon1_new, epsilon2_new, iterations, profit_ESS])
    ############THE SMAL SCRIPT SNIPPET FOR UPDATING THE REGULARIZATION TERM
    while epsilon1_new >= tol || epsilon2_new >= tol
        #After reaching iteration limit exit:
        if iterations >= iter_limit
            break
        end
        #########SETTING THE WARM START POINT TO THE PREVIOUS SOLUTION#########
        #if terminition status is local optimal, otherwise should not use the solution as starting point
        term_stat = termination_status(NLP)
        # if term_stat == MOI.OPTIMAL || term_stat == MOI.LOCALLY_SOLVED
        set_start_value.(all_variables(NLP), value.(all_variables(NLP)))
        ###################FIXING EPRSILONS TO THE REDUCED PARAMS##############
        global epsilon1_new = value.(NLP[:ϵ1])*0.1
        if epsilon2_new >= tol
            global epsilon2_new = value.(NLP[:ϵ2])*0.1
        end
        # else
        #     set_start_value.(all_variables(MPPDC), prev_vals)
        #     global epsilon_new = value.(MPPDC[:ϵ])*0.99
        # end
        fix(NLP[:ϵ1], epsilon1_new; force = true)
        fix(NLP[:ϵ2], epsilon2_new; force = true)

        println("
        #########
        THE NEW RELAXATION TERMS:
                    # ϵ1: $epsilon1_new and
                    ϵ2: $epsilon2_new and
                    #ITERATIONS: $iterations
                    #TERMINATION: $term_stat
        #########")
        ################OPTIMIZING############
        optimize!(NLP)
        ##############CACLCULATING PROFIT############
        profit_ESS, profit_GEN, profit_RES = compute_profits()
        println("PROFITS:
                        PROFIT_ESS: $profit_ESS,
                        PROFIT_GEN: $profit_GEN,
                        PROFIT_RES: $profit_RES")
        global iterations += 1
        push!(results_NLP_ESS, [epsilon1_new, epsilon2_new, iterations, profit_ESS])
        ############PROGRESSING IN TIME##############
        global prev_vals = value.(all_variables(NLP))
    end
end

@show res_NLP = results_NLP_ESS.profit[end]

#OR BUILD WITH KKT CONDITIONS INSTEAD OF MPPDC, NOTE: With the linearized objective, LP solver may be used
optimizer = with_optimizer(Gurobi.Optimizer)
MILP = build_MILP_quantity(t_steps, non_str_gens,
non_str_res, current, optimizer, objective_type)
optimize!(MILP)

@show res_MILP = objective_value(MILP)


text:
function _ossb_build_basis_of_solutions(instance::CombinatorialInstance{T}) where T
  # "Basis": the i-th returned vector must use arm i. The same solution might be returned for
  # several arms.
  d = dimension(instance)
  x = Vector{Vector{T}}(undef, d)
  all_arms = all_arm_indices(instance)
  for arm_to_include in 1:d
    # TODO: rather use DataStructures.jl's DefaultDict, but CombinatorialBandits.jl sometimes uses keys(weights) for a set of arms... Rather use all_arm_indices(instance) in a next refactoring.
    weights = Dict{T, Float64}()
    for (arm, arm_object) in enumerate(all_arms)
      weights[arm_object] = (arm == arm_to_include) ? 1.0 : 0.0
    end
    x[arm_to_include] = solve_linear(instance, weights)
  end
  return x
end

function _ossb_initial_feasible_solution(instance::CombinatorialInstance{T}, θ::Dict{T, Float64}) where T
  # Compute useful results.
  d = dimension(instance)
  m = maximum_solution_length(instance)
  Δmin = estimate_Δmin(instance, θ)
  idx_to_arm = Dict(idx => arm for (idx, arm) in enumerate(all_arm_indices(instance)))
  x = _ossb_build_basis_of_solutions(instance)

  # Ignore repeated solutions in the basis (loses an indexing property, but that's not important).
  x = collect(Set(x))

  # Compute the initial solution.
  v = d * m / Δmin^2
  t = zeros(d)
  for i in 1:d
    for j in 1:length(x)
      t[i] += any(x[j] .== Ref(idx_to_arm[i])) # Ref(.) to avoid broadcasting over the contents of solutions.
    end
  end

  return (t * v, v)
end

function _ossb_confusing_arms(instance::CombinatorialInstance{T}, θ::Dict{T, Float64}; ε::Float64=1.0e-6) where T
  # Precompute useful results.
  d = dimension(instance)
  idx_to_arm = Dict(idx => arm for (idx, arm) in enumerate(all_arm_indices(instance)))
  x_max = solve_linear(copy(instance), θ)
  r_max = sum(θ[arm] for arm in x_max)

  # Compute the set of confusing arms.
  I = Set{T}()
  for i in 1:d
    # Compute a solution that plays arm i *and* maximises the reward.
    # In particular, this is not necessarily a solution computed by _ossb_build_basis_of_solutions, because that
    # basis only ensures that an arm is included in the solution, not that it maximises the reward in any sense.
    weights = copy(θ)
    weights[idx_to_arm[i]] = d^2 * max(maximum(abs.(collect(values(θ)))), 1) # Force i to be part of the solution.
    xi = solve_linear(copy(instance), weights)

    if length(xi) == 0
      # Probably an error when optimising, skip.
      continue
    end

    # If this solution yields a suboptimum reward, then it cannot be part of an optimum solution.
    if sum(θ[arm] for arm in xi) <= r_max * (1.0 - ε)
      push!(I, idx_to_arm[i])
    end
  end

  return I
end

function _ossb_is_solution_feasible(instance::CombinatorialInstance{T}, θ::Dict{T, Float64}, t::Vector{Float64}, solve_all_budgets_at_once::Bool) where T
  # Check whether the constraints ∑i x_i / t_i ≤ Δx² are respected. Procedure: maximise ∑i x_i / t_i - Δx², check the sign.
  Δmax = estimate_Δmax(instance, θ)
  Δmin = estimate_Δmin(instance, θ)
  arm_to_idx = Dict(arm => idx for (idx, arm) in enumerate(all_arm_indices(instance)))
  I = _ossb_confusing_arms(instance, θ)
  # I = all_arm_indices(instance)
  I_all = all_arm_indices(instance) # Only for gap.

  x_max = solve_linear(copy(instance), θ)
  r_max = sum(θ[arm] for arm in x_max)

  # Ensure no component of t is zero.
  t = copy(t)
  t[t .== 0.0] .= minimum(t[t .> 0.0])

  # Actual maximisation.
  # m, _, vars = get_lp_formulation(copy(instance), θ) # A t ≤ b
  # @objective(m, Min, (sum(vars[i] * θ[i] for i in I_all) - r_max) ^2 - sum(vars[i] / t[arm_to_idx[i]] for i in I))
  # # @objective(m, Max, sum(vars[i] / t[arm_to_idx[i]] for i in I) - (sum(vars[i] * θ[i] for i in I) - r_max) ^2)
  # @constraint(m, sum(vars[i] for i in I) >= 1)
  # optimize!(m)
  # best_objective = objective_value(m)
  # @show [value(vars[i]) for i in I]
  # @show best_objective
  # @show (sum(value(vars[i]) * θ[i] for i in I_all) - r_max) ^2

  # println(m)
  rewards = Dict{T, Float64}(i => 1 / t[arm_to_idx[i]] for i in I_all)
  weights = Dict{T, Int}(i => (i in I) ? round(Int, θ[i] * t[arm_to_idx[i]] / Δmin, RoundUp) : 0.0 for i in I_all)
  _, o = _maximise_nonlinear_through_budgeted(instance, rewards, weights, round(Int, Δmax / Δmin),
                                                           (sol, budget) -> begin
                                                              # @show sol
                                                              # @show budget * Δmin
                                                              # @show (budget * Δmin - r_max)^2
                                                              sum(rewards[i] for i in sol) - (budget * Δmin - r_max)^2
                                                            end,
                                                           solve_all_budgets_at_once)

  # @show o

  # best_subobjective = best_objective
  # @show best_subobjective
  # # @show rewards
  # # @show weights
  return o >= 0.0
  # return best_objective >= 0.0
end

function _ossb_compute_feasible_v(instance::CombinatorialInstance{T}, θ::Dict{T, Float64}, t::Vector{Float64}) where T
  # Compute a v such that all constraints are respected, supposing that the input t satisfies the ∑i x_i / t_i ≤ Δx² constraints.
  # Two conditions:
  # - the objective is positive: v r_\max - θ^T t ≥ 0
  # - the point is a convex combination of integer solutions: A t ≤ v b
  arm_to_idx = Dict(arm => idx for (idx, arm) in enumerate(all_arm_indices(instance)))
  I = _ossb_confusing_arms(instance, θ)
  # I = all_arm_indices(instance)
  I_all = all_arm_indices(instance)

  x_max = solve_linear(copy(instance), θ)
  r_max = sum(θ[arm] for arm in x_max)

  # First condition gives a first estimation of v.
  θ⊤t = sum(θ[i] * t[arm_to_idx[i]] for i in I_all)
  v_min = θ⊤t / r_max

  # The second condition requires an explicit model to work with.
  m, _, vars = get_lp_formulation(instance, Dict(arm => 0.0 for arm in all_arm_indices(instance)))
  var_to_arm = Dict(values(vars) .=> keys(vars))

  # Compute the number of constraints.
  n_constraints = 0
  for (F, S) in list_of_constraint_types(m)
    # Ignore the case where this is an integrality/binary/bound constraint, we only need linear constraints.
    if S == MOI.Integer || S == MOI.ZeroOne || F == JuMP.VariableRef
      continue
    end

    n_constraints += num_constraints(m, F, S)
  end

  # Precompute A^T t and b.
  At = zeros(Float64, n_constraints)
  b = zeros(Float64, n_constraints)

  i = 1
  for (F, S) in list_of_constraint_types(m)
    # Ignore the case where this is an integrality/binary/bound constraint, we only need linear constraints.
    if S == MOI.Integer || S == MOI.ZeroOne || F == JuMP.VariableRef
      continue
    end

    for c in all_constraints(m, F, S)
      lhs = value(c, vr -> t[arm_to_idx[var_to_arm[vr]]])
      rhs = MOI.constant(moi_set(constraint_object(c)))

      At[i] = lhs
      b[i] = rhs
      i += 1
    end
  end

  # Check if this value satisfies the second constraint, and increase if need be.
  v = v_min
  if ! all(At .<= v .* b)
    vals = At ./ b
    vals[isinf.(vals)] .= -Inf # In case some b is zero.
    v = maximum(vals)
  end

  return v
end


text:
struct VonNeumann end

capacity(::VonNeumann, d::Int) = 2 * d

##### eachdir - iterating over adjacent sites
@inline function eachdir(::VonNeumann, x::Site{1})
  return ((x[1] - 1,),
          (x[1] + 1,))
end

@inline function eachdir(::VonNeumann, x::Site{2})
  return ((x[1] - 1, x[2]),
          (x[1] + 1, x[2]),
          (x[1], x[2] - 1),
          (x[1], x[2] + 1))
end

@inline function eachdir(::VonNeumann, x::Site{3})
  return ((x[1] - 1, x[2], x[3]),
          (x[1] + 1, x[2], x[3]),
          (x[1], x[2] - 1, x[3]),
          (x[1], x[2] + 1, x[3]),
          (x[1], x[2], x[3] - 1),
          (x[1], x[2], x[3] + 1))
end

##### distance - this could be done with @generated instead
@inline function distance(::VonNeumann, x::Site{1}, y::Site{1})
  return abs(x[1] - y[1])
end

@inline function distance(::VonNeumann, x::Site{2}, y::Site{2})
  return abs(x[1] - y[1]) + abs(x[2] - y[2])
end

@inline function distance(::VonNeumann, x::Site{3}, y::Site{3})
  return abs(x[1] - y[1]) + abs(x[2] - y[2]) + abs(x[3] - y[3])
end