unit IHCP_Unit;

interface

uses
  Windows, Messages, SysUtils, Variants, Classes, Graphics, Controls, Forms,
  Dialogs, StdCtrls, Spin, ComCtrls, ShellCtrls, DBCtrls, FileCtrl,
  Math, HCArrays, HCConduct, HCMaterial, HCModel, UPohMath, UPohStrConv,
     BasicUnit, DirectUnit, InverseUnit, dblookup, Grids, Outline, DirOutln;

type
  TFIHCP = class(TForm)
    MODEL_label: TLabel;
    DATA: TLabel;
    Label1: TLabel;
    Label2: TLabel;
    BtnBeck: TButton;
    BtnFull: TButton;
    Label3: TLabel;
    SelectFile: TOpenDialog;
    BtnData: TButton;
    BtnResults: TButton;
    Label4: TLabel;
    Material: TComboBox;
    SampleL: TEdit;
    termoclanek: TEdit;
    Freq: TEdit;
    EData: TEdit;
    EResults: TEdit;
    BtnFullReg: TButton;
    Label5: TLabel;
    Ewindow: TEdit;
    Label6: TLabel;
    BtnWindow: TButton;
      
    procedure BtnBeckClick(Sender: TObject);
    procedure BtnDataClick(Sender: TObject);
    procedure BtnResultsClick(Sender: TObject);
    procedure BtnFullClick(Sender: TObject);
    procedure BtnFullRegClick(Sender: TObject);
    procedure BtnWindowClick(Sender: TObject);
  private
    { Private declarations }
  public
    { Public declarations }
  end;

var
  FIHCP: TFIHCP;
  StartTime: Dword;
  Time:double;

  MatDtb: tMatDB;
  model: tHCModel;
  MatID: integer;
  str: string;
  DataFile, ResultFile:TextFile;

  L, dt, hloubka, InitTemp: double;
  pom_temp: array [0..10000] of double;
  q, T, T_comp, T_surf: array of double;
  n, i: integer;

const
  mm=0.001;
  Cel_Kel=273.15;
  max_volL = 0.05*mm;

implementation

{$R *.dfm}

procedure TFIHCP.BtnDataClick(Sender: TObject);
begin
  SelectFile:= TOpenDialog.Create(Self);
  SelectFile.InitialDir:=GetCurrentDir;
  if SelectFile.Execute then begin
     EData.Text:=SelectFile.FileName;
     AssignFile(DataFile, SelectFile.FileName);
     end
     else if EData.Text='' then EData.Text:='Soubor nebyl vybrán!';
  SelectFile.Free;
end;

procedure TFIHCP.BtnResultsClick(Sender: TObject);
begin
  SelectFile:=TOpenDialog.Create(Self);
  SelectFile.InitialDir:=GetCurrentDir;
  if SelectFile.Execute then begin
     EResults.Text:=SelectFile.FileName;
     AssignFile(ResultFile,SelectFile.FileName);
     end
     else if EResults.Text='' then EResults.Text:='Soubor nebyl vybrán!';
  SelectFile.Free;
end;

procedure TFIHCP.BtnBeckClick(Sender: TObject);
var
  FTsteps: integer;

begin
  StartTime:=GetTickCount();
  BtnBeck.Cursor:=crHourGlass;

// Nacteni dialogoveho okna
  str:=Material.Text;
  L := StrToFloat(SampleL.Text)*mm;           // mm
  hloubka:=StrToFloat(termoclanek.Text)*mm;   // mm
  dt:= 1/ StrToFloat(Freq.Text);              //s

// Nacteni datoveho souboru a souboru pro ukladani vysledku
  Reset(DataFile);
  Rewrite(ResultFile);

// Nacteni databaze a materialu
  MatID:=InitMat(MatDtb, str);
  if (MatID < 0) then raise Exception.Create('Material "' + str + '" was not found.');

// Nacteni teplot ze souboru
  ReadData(DataFile, pom_temp, n);
  n:=n+1;
  SetLength(T,n);
  for i:=0 to (n-1) do begin
    T[i]:=pom_temp[i];
  end;
  InitTemp:=T[0];

// Volba poctu doprednych kroku + IHCP Beck
  FTsteps:=FindFTsteps(n, InitTemp, L, MatID, MatDtb, max_volL, dt, hloubka);
  SetLength(q,(n-FTsteps)); SetLength(T_comp,(n-FTsteps)); SetLength(T_surf,(n-FTsteps));
  q[0]:=0;
  IHCP_Beck(MatDtb, L, InitTemp, MatID, max_volL, n, FTsteps, T, q, T_comp, T_surf, dt, hloubka);

// Vypis vysledku
  writeln(ResultFile,'Vzorek',#9,'Čas [s]',#9,'T_meas [°C]',#9,'T_comp [°C]',#9,'T_surf [°C]',#9,'q [W/m^2]',#9,'err [°C]');
  for i:=1 to (n-FTsteps-1) do
    writeln(ResultFile,FloatToStr(i),#9,FloatToStr(i*dt),#9,FloatToStr(T[i]),#9,FloatToStr(T_comp[i]),#9,FloatToStr(T_surf[i]),#9,FloatToStr(-q[i]),#9,FloatToStr(T[i]-T_comp[i]));

// free RAM a ukonceni
  fHC_Free(MatDtb);
  CloseFile(ResultFile); CloseFile(DataFile);

  Time:=(GetTickCount()-StartTime)/1000;
  ShowMessage('Výpočet sekvenční metodou byl úspěšně ukončen. Potřebný čas byl: '+FloatToStr(Time)+' sekund.');
  Application.Terminate;
end;


procedure TFIHCP.BtnFullClick(Sender: TObject);
var
  del_theta: array of double;

begin
  StartTime:=GetTickCount();
  BtnFull.Cursor:=crHourGlass;

// Nacteni dialogoveho okna
  str:=Material.Text;
  L := StrToFloat(SampleL.Text)*mm;           // mm
  hloubka:=StrToFloat(termoclanek.Text)*mm;   // mm
  dt:= 1/ StrToFloat(Freq.Text);              //s

// Nacteni datoveho souboru a souboru pro ukladani vysledku
  Reset(DataFile);
  Rewrite(ResultFile);

// Nacteni databaze a materialu
  MatID:=InitMat(MatDtb, str);
  if (MatID < 0) then raise Exception.Create('Material "' + str + '" was not found.');

// Nacteni teplot ze souboru
  ReadData(DataFile, pom_temp, n);
  SetLength(T,n);
  InitTemp:=pom_temp[0];
  for i:=1 to n do begin
    T[i-1]:=pom_temp[i]-InitTemp;
  end;

// Nastaveni koeficientu citlivosti + IHCP FULL DOMAIN
  SetLength(del_theta,n); SetLength(q,n); SetLength(T_comp,n); SetLength(T_surf,n);
  Sensitivity(n, L, MatID, MatDtb, max_volL, dt, del_theta, hloubka);
  IHCP_Full(n, T, q, del_theta);

// VYPOCET PRIME ULOHY PRO ZJISTENI TEPLOT
  // Model
    InitModel(model, L, InitTemp, MatID, max_volL);
    SetBoundCond(model);
    fHC_ClearAccuracy(model); // clear accuracy arrays

  // Prima uloha
    for i:=0 to (n-1) do q[i]:=-q[i];
    MainDirectTask(T_comp, T_surf, MatDtb, model, q, n, dt, hloubka);

// Vypis vysledku
  writeln(ResultFile,'Vzorek',#9,'Čas [s]',#9,'T_meas [°C]',#9,'T_comp [°C]',#9,'T_surf [°C]',#9,'q [W/m^2]',#9,'err [°C]');
  for i:=0 to (n-1) do
    writeln(ResultFile,FloatToStr(i+1),#9,FloatToStr((i+1)*dt),#9,FloatToStr(pom_temp[i+1]),#9,FloatToStr(T_comp[i]),#9,FloatToStr(T_surf[i]),#9,FloatToStr(-q[i]),#9,FloatToStr(pom_temp[i+1]-T_comp[i]));

// free RAM a ukonceni
	fHC_Free(model);
  fHC_Free(MatDtb);
  CloseFile(ResultFile); CloseFile(DataFile);
  
  Time:=(GetTickCount()-StartTime)/1000;
  ShowMessage('Výpočet full-domain metodou byl úspěšně ukončen. Potřebný čas byl: '+FloatToStr(Time)+' sekund.');
  Application.Terminate;
end;


procedure TFIHCP.BtnFullRegClick(Sender: TObject);
var
  del_theta: array of double;

begin
  StartTime:=GetTickCount();
  BtnFullReg.Cursor:=crHourGlass;

// Nacteni dialogoveho okna
  str:=Material.Text;
  L := StrToFloat(SampleL.Text)*mm;           // mm
  hloubka:=StrToFloat(termoclanek.Text)*mm;   // mm
  dt:= 1/ StrToFloat(Freq.Text);              //s

// Nacteni datoveho souboru a souboru pro ukladani vysledku
  Reset(DataFile);
  Rewrite(ResultFile);

// Nacteni databaze a materialu
  MatID:=InitMat(MatDtb, str);
  if (MatID < 0) then raise Exception.Create('Material "' + str + '" was not found.');

// Nacteni teplot ze souboru
  ReadData(DataFile, pom_temp, n);
  SetLength(T,n);
  InitTemp:=pom_temp[0];
  for i:=1 to n do begin
    T[i-1]:=pom_temp[i]-InitTemp;
  end;

// Nastaveni koeficientu citlivosti + IHCP FULL+reg DOMAIN
  SetLength(del_theta,n); SetLength(q,n); SetLength(T_comp,n); SetLength(T_surf,n);
  Sensitivity(n, L, MatID, MatDtb, max_volL, dt, del_theta, hloubka);
  IHCP_Full_reg(n, T, q, del_theta);

// VYPOCET PRIME ULOHY PRO ZJISTENI TEPLOT
  // Model
    InitModel(model, L, InitTemp, MatID, max_volL);
    SetBoundCond(model);
    fHC_ClearAccuracy(model);
  // Direct Task
    for i:=0 to (n-1) do q[i]:=-q[i];
    MainDirectTask(T_comp, T_surf, MatDtb, model, q, n, dt, hloubka);

// Vypis vysledku
  writeln(ResultFile,'Vzorek',#9,'Čas [s]',#9,'T_meas [°C]',#9,'T_comp [°C]',#9,'T_surf [°C]',#9,'q [W/m^2]',#9,'err [°C]');
  for i:=0 to (n-1) do
    writeln(ResultFile,FloatToStr(i+1),#9,FloatToStr((i+1)*dt),#9,FloatToStr(pom_temp[i+1]),#9,FloatToStr(T_comp[i]),#9,FloatToStr(T_surf[i]),#9,FloatToStr(-q[i]),#9,FloatToStr(pom_temp[i+1]-T_comp[i]));

// free RAM a ukonceni
  fHC_Free(MatDtb);
  CloseFile(ResultFile); CloseFile(DataFile);

  Time:=(GetTickCount()-StartTime)/1000;
  ShowMessage('Výpočet regularizovanou full-domain metodou byl úspěšně ukončen. Potřebný čas byl: '+FloatToStr(Time)+' sekund.');
  Application.Terminate;
end;

procedure TFIHCP.BtnWindowClick(Sender: TObject);
var
  vTbck : tFloatArray;
  Max: double;
  Temp, Temp2, T_win, T_swin, T_qnull, T_qnull2, q, q_first, q_sec, q_third, del_theta: array of double;
  window, overlap, j, k: integer;

begin
  StartTime:=GetTickCount();
  BtnWindow.Cursor:=crHourGlass;

// Nacteni dialogoveho okna
  str:=Material.Text;
  L := StrToFloat(SampleL.Text)*mm;           // mm
  hloubka:=StrToFloat(termoclanek.Text)*mm;   // mm
  dt:= 1/ StrToFloat(Freq.Text);              //s
  window := StrToInt(Ewindow.Text);
  Max:=0;

// Nacteni datoveho souboru a souboru pro ukladani vysledku
  Reset(DataFile);
  Rewrite(ResultFile);

// Nacteni databaze a materialu
  MatID:=InitMat(MatDtb, str);
  if (MatID < 0) then raise Exception.Create('Material "' + str + '" was not found.');

// Nacteni teplot ze souboru
  ReadData(DataFile, pom_temp, n);
  SetLength(T,n);
  for i:=1 to n do begin
    T[i-1]:=pom_temp[i];
  end;
  InitTemp:=pom_temp[0];
  if window>n then window:=n;

// Nastaveni koeficientu citlivosti + velikost prekryvu
  SetLength(del_theta, window); SetLength(Temp,window); SetLength(T_qnull,window);
  Sensitivity(window, L, MatID, MatDtb, max_volL, dt, del_theta, hloubka);
  for i:=0 to (window-1) do
    if Abs(del_theta[i])>Max then begin
      Max:=Abs(del_theta[i]);
      overlap:=i;
    end;
  overlap:=250;
  if overlap>window then overlap:=window-2;  // zajisteni, ze prekryv neni delsi nez okno

// Nastaveni poctu oken
  if ((n-overlap+1) mod (window-overlap) = 0) then k:=((n-overlap+1) div (window-overlap))-1
  else k:=((n-overlap+1) div (window-overlap));

// Nacteni teplot pro prvni okno
  for i:=0 to (window-1) do begin
    Temp[i]:=T[i]-InitTemp;
  end;

// IHCP FULL+reg na prvnim okne
  SetLength(q,n); SetLength(T_comp,n); SetLength(T_surf,n);
  SetLength(q_first,window);
  SetLength(q_sec,window-overlap); SetLength(T_win, window-overlap); SetLength(T_swin,window-overlap);
  IHCP_Full_reg(window, Temp, q_first, del_theta);

// zapis q prvniho okna (do pulky prekryvu)
  for i:=0 to (window-overlap-1) do
     q_sec[i]:=-q_first[i];
  for i:=0 to (window-overlap)+(overlap div 2)  do
     q[i]:=q_first[i];

// VYPOCET PRIME ULOHY PRO ZJISTENI TEPLOT v bode zacatku dalsiho okna
  // Model
    InitModel(model, L, InitTemp, MatID, max_volL);
    fHC_Init(vTbck);
    SetBoundCond(model);
    fHC_ClearAccuracy(model);
  // Prima uloha
    MainDirectTask(T_win, T_swin, MatDtb, model, q_sec, (window-overlap), dt, hloubka);
    fHC_Copy(vTbck, model.vT);  // zaloha teplotniho pole

  for j:=1 to (k-1) do begin  // Zacatek cyklu pro dalsi okna
    TempNull(model, vTbck, window, MatDtb, dt, T_qnull, hloubka);  // Spocteni teplot pro q=0 v celem okne
    // Nacteni teplot pro dane okno (bez vlivu z minulosti)
    for i:=0 to (window-1) do
      Temp[i]:=T[j*(window-overlap)+i]-T_qnull[i];

    // IHCP FULL+reg a zapis (od pulkdy prekryvu minuleho okna, do pulky prekryvu nasledujiciho)
    IHCP_Full_reg(window, Temp, q_first, del_theta);
    for i:=0 to (window-overlap-1) do
      q_sec[i]:=-q_first[i];
    for i:= ((overlap div 2)+1) to (window-overlap)+(overlap div 2) do
       q[j*(window-overlap)+i]:=q_first[i];

    // VYPOCET PRIME ULOHY PRO ZJISTENI TEPLOT v bode zacatku dalsiho okna
    // Model
      fHC_Copy(model.vT,vTbck); // vlozeni zalohovaneho teplotniho pole
    // Prima uloha
      MainDirectTask(T_win, T_swin, MatDtb, model, q_sec, (window-overlap), dt, hloubka);
      fHC_Copy(vTbck, model.vT);  // zaloha teplotniho pole
  end;

// Vypocet posledniho (mensiho) okna
  SetLength(Temp2,n-k*(window-overlap)+1); SetLength(T_qnull2,n-k*(window-overlap)+1);
  // Nacteni teplot pro posledni okno (bez vlivu z minulosti)
  TempNull(model, vTbck, n-k*(window-overlap)+1, MatDtb, dt, T_qnull2, hloubka);
  for i:=0 to (n-k*(window-overlap)-1) do
      Temp2[i]:=T[k*(window-overlap)+i]-T_qnull2[i];

  // IHCP FULL+reg a zapis do q
    SetLength(q_third, n-k*(window-overlap));
    IHCP_Full_reg(n-k*(window-overlap), Temp2, q_third, del_theta);  
    for i:=((overlap div 2)+1) to (n-k*(window-overlap)-1) do q[k*(window-overlap)+i]:=q_third[i];

  fHC_Free(model);

// VYPOCET PRIME ULOHY PRO ZJISTENI TEPLOT
  // Model
    InitModel(model, L, InitTemp, MatID, max_volL);
    SetBoundCond(model);
    fHC_ClearAccuracy(model); // clear accuracy arrays
  // Prima uloha
    for i:=0 to (n-1) do q[i]:=-q[i];
    MainDirectTask(T_comp, T_surf, MatDtb, model, q, n, dt, hloubka);

// Vypis vysledku
  writeln(ResultFile,'Vzorek',#9,'Čas [s]',#9,'T_meas [°C]',#9,'T_comp [°C]',#9,'T_surf [°C]',#9,'q [W/m^2]',#9,'err [°C]');
  for i:=0 to (n-1) do
    writeln(ResultFile,FloatToStr(i+1),#9,FloatToStr((i+1)*dt),#9,FloatToStr(pom_temp[i+1]),#9,FloatToStr(T_comp[i]),#9,FloatToStr(T_surf[i]),#9,FloatToStr(-q[i]),#9,FloatToStr(pom_temp[i+1]-T_comp[i]));

// free RAM a ukonceni
	fHC_Free(model);
  fHC_Free(MatDtb);
  CloseFile(ResultFile); CloseFile(DataFile);

  Time:=(GetTickCount()-StartTime)/1000;
  ShowMessage('Výpočet metodou s posuvným oknem byl úspěšně ukončen. Potřebný čas byl: '+FloatToStr(Time)+' sekund.');
  Application.Terminate;
end;

end.
