unit uArcos_ConPerdidasCuadraticas;
interface

uses
  Classes,
  SysUtils, xmatdefs, usimplex, uglobs,
  uNodos,
  uActorNodal,
  ucosa, ufechas, uconstantesSimSEE, uranddispos,
  uFichasLPD, Math;

type

  { TFicha_ArcoPerdidasCuadraticas }

  TFicha_ArcoPerdidasCuadraticas = class(TFichaLPD)
  public
    rendimiento: NReal;
    peaje: NReal; //[USD/MWh] costo variable por transito
    fd: NReal; // factor de disponibilidad
    tRepHoras: NReal;
    PMax: NReal; //[MW] potencia máxima admisible

    constructor Create(fecha: TFecha; periodicidad: TPeriodicidad;
      rendimiento: NReal; peaje: NReal;
    //[USD/MWh] costo variable por transito
      fd: NReal; // factor de disponibilidad
      tRepHoras: NReal; PMax: NReal //[MW] potencia máxima admisible
      );
    constructor Create_ReadFromText(f: TArchiTexto); override;
    procedure WriteToText(f: TArchiTexto); override;
    function infoAd_: string; override;
    procedure Free; override;
  end;

  { TArcoPerdidasCuadraticas }

  TArcoPerdidasCuadraticas = class(TActorBiNodal)
  public
    NLineasDisponibles: integer;
    PMaxDisponible: NReal;
    // resultados del despacho
    P_Entrante: TDAOfNReal; //[MW] Potencia Entrante

    (**************************************************************************)
    (*               A T R I B U T O S   P E R S I S T E N T E S              *)
    (**************************************************************************)

    //Lo dejo fuera del iterar para no perder la informacion si se carga una sala
    //en un editor compilado sin ITERAR
    aplicarPerdidasCuadraticas: boolean;
    (**************************************************************************)
    k: NReal;
    P_PasoAnterior: TDAofNReal; //Guarda las potencias despachadas en el paso
    //anterior. Las usa para iterar y aproximar las
    //pérdidas cuadráticas

    //Costo por poste de la transmisión
    //costo[iposte]:= P[iposte] * pa.peaje * globs.durpos[iposte];
    costo: TDAofNReal;
    // valor por poste del multiplicador de Lagrange de la restrición de PMax
    costo_congestion: TDAOfNReal; //[USD/MW]

    pa: TFicha_ArcoPerdidasCuadraticas;

    constructor Create(nombre: string; nacimiento, muerte: TFecha;
      lpdUnidades, lpd: TFichasLPD; Entrada, Salida: TNodo;
      aplicarPerdidasCuadraticas: boolean);

    constructor Create_ReadFromText(f: TArchiTexto); override;
    procedure WriteToText(f: TArchiTexto); override;

    function Rec: TCosa_RecLnk; override;
    procedure BeforeRead(version, id_hilo: integer); override;
    procedure AfterRead(f:TArchiTexto); override;

    procedure PrepararMemoria(globs: TGlobs); override;
    procedure RegistrarParametrosDinamicos; override;
    function InfoAd_: string; override;
    class function DescClase: string; override;

    procedure Sim_Cronica_Inicio; override;
    procedure SorteosDelPaso(sortear: boolean); override;

    procedure opt_nvers(var ivar, ivae, ires: integer); override;
    procedure opt_cargue(s: TSimplex); override;
    procedure opt_fijarRestriccionesDeCaja(s: TSimplex); override;
    procedure opt_leerSolucion(s: TSimplex); override;

    function opt_necesitoIterar(kIteracion: integer;
      var errRelativo: NReal): boolean;
      override;

    function getNombreVar(ivar: integer; var nombre: string): boolean; override;
    function getNombreRes(ires: integer; var nombre: string): boolean; override;

    procedure PosicionarseEnEstrellita; override;
    procedure AcumAux1(peso: NReal); override;
    procedure SetAux1; override;

    procedure optx_nvxs(var ixr, ixd, iauxr, iauxd: integer); override;

    procedure dump_Variables(var f: TextFile; charIndentacion: char); override;

    //    function infoAd : String; virtual;
    procedure CambioFichaPD; override;
    procedure PubliVars; override;
    procedure Free; override;


  end;

procedure TArcoPerdidasCuadraticas_cambioFichaPD(Actor: TCosa);

procedure AlInicio;
procedure AlFinal;

implementation

uses uActores;


//-------------------
// Métodos de TArcoPerdidasCuadraticas
//===================

constructor TArcoPerdidasCuadraticas.Create(nombre: string; nacimiento, muerte: TFecha;
  lpdUnidades, lpd: TFichasLPD; Entrada, Salida: TNodo;
  aplicarPerdidasCuadraticas: boolean);
begin
  inherited Create(nombre, nacimiento, muerte, lpdUnidades, Entrada, Salida);
  self.lpd := lpd;
  Self.aplicarPerdidasCuadraticas := aplicarPerdidasCuadraticas;
end;

function TArcoPerdidasCuadraticas.Rec: TCosa_RecLnk;
begin
  Result:=inherited Rec;
  Result.addCampoDef('lpd', TCosa(lpd));
  Result.addCampoDef('aplicarPerdidasCuadraticas', aplicarPerdidasCuadraticas);
end;

procedure TArcoPerdidasCuadraticas.PrepararMemoria(globs: TGlobs);
begin
  inherited prepararMemoria(globs);
  setlength(P_Entrante, globs.NPostes);
  SetLength(costo, globs.NPostes);
  setlength(costo_congestion, globs.NPostes);

  if aplicarPerdidasCuadraticas then
    SetLength(P_PasoAnterior, globs.NPostes)
  else
    P_PasoAnterior := nil;
end;

procedure TArcoPerdidasCuadraticas.RegistrarParametrosDinamicos;
begin
  inherited registrarParametrosDinamicos;
  lpd.expandirFichas(globs);
  lpd.RegistrarFichasAActualizar(Self, globs.ActualizadorLPD, @pA, nil,
    TArcoPerdidasCuadraticas_cambioFichaPD);
end;

function TArcoPerdidasCuadraticas.InfoAd_: string;
begin
  Result := '';
end;

class function TArcoPerdidasCuadraticas.DescClase: string;
begin
  Result := 'Arco';
end;

procedure TArcoPerdidasCuadraticas.Sim_Cronica_Inicio;
begin
  inherited Sim_Cronica_Inicio;
  //Para que de igual si se ejecuta mas de una vez seguida
  if aplicarPerdidasCuadraticas then
    vclear(P_PasoAnterior);
end;

procedure TArcoPerdidasCuadraticas.SorteosDelPaso(sortear: boolean);
begin
  if globs.ObligarDisponibilidad_1_ then
  begin
    NLineasDisponibles := paUnidades.nUnidades[0];
    PMaxDisponible := pa.PMax * NLineasDisponibles;
  end
  else
  if sortear then
  begin
    NLineasDisponibles := Sorteos_RepRotUnidades;
    PMaxDisponible := pa.PMax * NLineasDisponibles;
  end
  else
  begin
    NLineasDisponibles := paUnidades.nUnidades[0];
    PMaxDisponible := pa.PMax * pa.fd * NLineasDisponibles;
  end;
end;

procedure TArcoPerdidasCuadraticas.opt_nvers(var ivar, ivae, ires: integer);
begin
  Self.ivar := ivar;
  if NLineasDisponibles > 0 then
    ivar := ivar + globs.NPostes;
end;

procedure TArcoPerdidasCuadraticas.opt_cargue(s: TSimplex);
var
  ibaseResSal, ibaseResEnt: integer;
  iposte: integer;
begin
  if NLineasDisponibles > 0 then
  begin

{$IFDEF SPXCONLOG}
    for iposte := 0 to globs.NPostes - 1 do
      s.set_NombreVar(ivar + iposte, Nombre + '_P[MW]' + IntToStr(iposte + 1));
{$ENDIF}

    // Tenemos que aportar a las restricciones de
    // demanda de los nodos de
    ibaseResEnt := NodoA.ires;
    ibaseResSal := NodoB.ires;
    if aplicarPerdidasCuadraticas then
    begin
      for iposte := 0 to globs.NPostes - 1 do
      begin
        s.pon_e(ibaseResEnt + iposte, ivar + iposte, -1);
        s.pon_e(ibaseResSal + iposte, ivar + iposte, 1 - 2 * k * P_Entrante[iposte]);
        s.acum_e(ibaseResSal + iposte, s.nc, k * P_Entrante[iposte] *
          P_Entrante[iposte]);
        // Ahora agregamos a la función de UTILIDAD (-costo)
        s.pon_e(s.nf, ivar + iposte, -pa.peaje * globs.durpos[iposte]);
      end;
    end
    else
      for iposte := 0 to globs.NPostes - 1 do
      begin
        s.pon_e(ibaseResEnt + iposte, ivar + iposte, -1);
        s.pon_e(ibaseResSal + iposte, ivar + iposte, pa.rendimiento);
        // Ahora agregamos a la función de UTILIDAD (-costo)
        s.pon_e(s.nf, ivar + iposte, -pa.peaje * globs.durpos[iposte]);
      end;

  end;
end;

procedure TArcoPerdidasCuadraticas.opt_fijarRestriccionesDeCaja(s: TSimplex);
var
  iposte: integer;
begin
  if NLineasDisponibles > 0 then
    // Le fijamos como cota máxima PMax a la potencia en todos los postes
    for iposte := 0 to globs.NPostes - 1 do
      s.cota_sup_set(ivar + iposte, PMaxDisponible);
end;

procedure TArcoPerdidasCuadraticas.opt_leerSolucion(s: TSimplex);
var
  iposte: integer;
  m: NReal;
  P: NReal;
begin
  CostoDirectoDelPaso := 0;
  vclear(P_Entrante);
  vclear(P_NodoA);
  vclear(P_NodoB);
  vclear(costo_congestion);

  if NLineasDisponibles > 0 then
    // recuperamos los valores de Potencia despachada
    for iposte := 0 to globs.NPostes - 1 do
    begin
      if s.flg_x[ivar + iposte] > 0 then // 1 = TopeActivo 2 = Variable Fijada
      begin
        P := s.xval(ivar + iposte);
        P_Entrante[iposte] := P;
        m := P * pa.peaje * globs.DurPos[iposte];
        costo[iposte] := m;
        CostoDirectoDelPaso := CostoDirectoDelPaso + m;
        costo_congestion[iposte] := -s.xmult(ivar + iposte);
        P_NodoA[iposte] := -P;
        P_NodoB[iposte] := P * pa.rendimiento;
      end;
    end;
end;



function TArcoPerdidasCuadraticas.opt_necesitoIterar(kIteracion: integer; var errRelativo: NReal): boolean;
var
  iposte: integer;
  errRelativoi, pmaxi: NReal;
begin
  if aplicarPerdidasCuadraticas then
  begin
    errRelativo := 0;
    //errRelativo retorna con el máximo error relativo de todos los postes

    for iposte := 0 to globs.NPostes - 1 do
    begin
      pmaxi := max(P_Entrante[iposte], P_PasoAnterior[iposte]);

      if pmaxi > 1 then
      begin
        errRelativoi := abs(P_Entrante[iposte] - P_PasoAnterior[iposte]) / pmaxi;
        if errRelativoi > errRelativo then
          errRelativo := errRelativoi;
      end;
      P_Entrante[iposte] := 0.9 * P_Entrante[iposte] + 0.1 * P_PasoAnterior[iposte];
    end;

    Result := errRelativo > 0.1;
  end
  else
    Result := False;
end;


function TArcoPerdidasCuadraticas.getNombreVar(ivar: integer; var nombre: string): boolean;
begin
  if (ivar >= self.ivar) and (ivar < self.ivar + globs.NPostes) then
  begin
    nombre := self.Nombre + '_P[MW]' + IntToStr(ivar - self.ivar + 1);
    Result := True;
  end
  else
    Result := False;
end;

function TArcoPerdidasCuadraticas.getNombreRes(ires: integer; var nombre: string): boolean;
begin
  Result := False;
end;

procedure TArcoPerdidasCuadraticas.PosicionarseEnEstrellita;
var
  iposte: integer;
begin
  if aplicarPerdidasCuadraticas then
  begin
    for iposte := 0 to globs.NPostes - 1 do
    begin
      self.P_Entrante[iposte] :=
        globs.Auxs_r0[iauxr + iposte][globs.CF.ordinalEstrellaActual];
      //Para que de igual en la optimización monohilo y multihilo
      self.P_PasoAnterior[iposte] := self.P_Entrante[iposte];
    end;
  end;
end;

procedure TArcoPerdidasCuadraticas.AcumAux1(peso: NReal);
var
  iposte: integer;
begin
  if aplicarPerdidasCuadraticas then
  begin
    for iposte := 0 to globs.NPostes - 1 do
    begin
      globs.Auxs_r1[iauxr + iposte][globs.CF.ordinalEstrellaActual] :=
        globs.Auxs_r1[iauxr + iposte][globs.CF.ordinalEstrellaActual] +
        P_Entrante[iposte] * peso;
    end;
  end;
end;

procedure TArcoPerdidasCuadraticas.SetAux1;
var
  iposte: integer;
begin
  if aplicarPerdidasCuadraticas then
  begin
    for iposte := 0 to globs.NPostes - 1 do
    begin
      globs.Auxs_r1[iauxr + iposte][globs.CF.ordinalEstrellaActual] :=
        P_Entrante[iposte];
    end;
  end;
end;

procedure TArcoPerdidasCuadraticas.optx_nvxs(var ixr, ixd, iauxr, iauxd: integer);
begin
  if aplicarPerdidasCuadraticas then
  begin
    Self.iauxr := iauxr;
    iauxr := iauxr + globs.NPostes;
  end;
end;

procedure TArcoPerdidasCuadraticas.dump_Variables(var f: TextFile; charIndentacion: char);
begin
  inherited dump_Variables(f, charIndentacion);
  writeln(f, charIndentacion, 'peaje= ', FloatToStrF(pa.peaje, ffFixed, 10, 3));
  writeln(f, charIndentacion, 'PMax= ', FloatToStrF(pa.PMax, ffFixed, 10, 3));
  writeln(f);
end;

procedure TArcoPerdidasCuadraticas.CambioFichaPD;
begin
  ActualizarProbabilidadesReparacionYRotura(pa.fd, pa.tRepHoras);
  if aplicarPerdidasCuadraticas then
    k := (1 - pa.rendimiento) / pa.PMax;
end;

procedure TArcoPerdidasCuadraticas.PubliVars;
begin
  inherited PubliVars;
  PublicarVariableVR('P', '[MW]', 6, 1, P_Entrante, True, True);
  PublicarVariableVR('Costo', '[USD]', 6, 1, costo, True, True);
  PublicarVariableVR('CostoCongestion', '[USD/MWh]', 6, 1, costo_congestion, True, True);
  PublicarVariableNI('NLineasDisponibles', '[u]', NLineasDisponibles, True);
end;

procedure TArcoPerdidasCuadraticas.Free;
begin
  lpd.Free;
  setlength(P_Entrante, 0);
  SetLength(costo, 0);
  setlength(costo_congestion, 0);
  if P_PasoAnterior <> nil then
    SetLength(P_PasoAnterior, 0);
  inherited Free;
end;


//---------------------
//Metodos de TFicha_ArcoPerdidasCuadraticas
//=====================

constructor TFicha_ArcoPerdidasCuadraticas.Create(fecha: TFecha; periodicidad: TPeriodicidad;
  rendimiento: NReal; peaje: NReal;
  //[USD/MWh] costo variable por transito
  fd: NReal; // factor de disponibilidad
  tRepHoras: NReal; PMax: NReal //[MW] potencia máxima admisible
  );
begin
  inherited Create(fecha, periodicidad);
  self.rendimiento := rendimiento;
  self.peaje := peaje;
  self.fd := fd;
  self.tRepHoras := tRepHoras;
  self.PMax := PMax;
end;

constructor TFicha_ArcoPerdidasCuadraticas.Create_ReadFromText(f: TArchiTexto);
begin
  if f.Version < 14 then
  begin
    inherited Create_ReadFromText(f);
    f.IniciarLecturaRetrasada;
    f.rd('rendimiento', rendimiento);
    f.rd('peaje', peaje);
    f.rd('fd', fd);
    f.rd('PMax', PMax);
    f.EjecutarLectura;
    self.tRepHoras := 8;
  end
  else
  begin
    inherited Create_ReadFromText(f);
    f.IniciarLecturaRetrasada;
    f.rd('rendimiento', rendimiento);
    f.rd('peaje', peaje);
    f.rd('fd', fd);
    f.rd('tRepHoras', tRepHoras);
    f.rd('PMax', PMax);
    f.EjecutarLectura;
  end;
end;

procedure TFicha_ArcoPerdidasCuadraticas.WriteToText(f: TArchiTexto);
begin
  inherited WriteToText(f);
  f.wr('rendimiento', rendimiento);
  f.wr('peaje', peaje);
  f.wr('fd', fd);
  f.wr('tRepHoras', tRepHoras);
  f.wr('PMax', PMax);
end;

function TFicha_ArcoPerdidasCuadraticas.infoAd_: string;
begin
  Result := 'PMáx= ' + FloatToStrF(PMax, ffGeneral, 10, 1) + ' MW, ' +
    'fDisp= ' + FloatToStrF(fd, ffGeneral, 10, 2) + ' p.u., ' +
    'tRep= ' + FloatToStrF(tRepHoras, ffGeneral, 10, 1) + 'h, ' +
    'Rendimiento= ' + FloatToStrF(rendimiento, ffGeneral, 10, 2) +
    ' p.u., ' + 'peaje= ' + FloatToStrF(peaje, ffGeneral, 10, 1) + ' USD/MWh';
end;

procedure TFicha_ArcoPerdidasCuadraticas.Free;
begin
  inherited Free;
end;


procedure TArcoPerdidasCuadraticas_cambioFichaPD(Actor: TCosa);
begin
  (Actor as TArcoPerdidasCuadraticas).cambioFichaPD;
end;

procedure AlInicio;
begin
  registrarClaseDeCosa(TArcoPerdidasCuadraticas.ClassName, TArcoPerdidasCuadraticas);
  registrarClaseDeCosa(TFicha_ArcoPerdidasCuadraticas.ClassName, TFicha_ArcoPerdidasCuadraticas);
end;

procedure AlFinal;
begin
end;

end.
