#!/usr/bin/octave -qf
%
% 2014
%

%
% SetupParameters
%

function SetupParameters()

  global parameters;


  % Arbeitsgas: ARGON
 
  % U_mol = T * S_mol - R * T % molare Innere Energie
  % R_S = R / M_Argon % spezifische Gaskonstante Argon
  % U_S = T * S_S - R_S * T % spezifische Innere Energie Argon
  % S_S = 3,875 % kJ/K pro kg % spezifische Entropie Argon
  parameters.M_S = 0.03995; % kg/mol spezifische Masse für Argon
  parameters.R_S = 8.3145e-3 / parameters.M_S; % kJ / kg * K spezifische Gaskonstante für Argon
  parameters.S_Std = 3.876; % kJ / K spezifische Entropie Argon bei 1 bar und 25°C (siehe D. Joensson - TECHNISCHE MECHANIK, Band 29, Heft 1)
  parameters.T_Std = 25 + 273; % Temperatur für Standard-Entropie
  parameters.p_Std = 1.e5; % Druck für Standard-Entropie in N / m^2

  parameters.C_p = 5 / 2; % isobare Wärmekapazität zweiatomiges Gas
  parameters.C_V = 3 / 2; % isochore Wäremkapazität zweiatomiges Gas (Wasserstoff)
  parameters.kappa = 5 / 3; % C_p / C_V = 7/5 % zweiatomiges Gas

  % Druck = 1e5 / Volumen 1 m^3 / T_Std = 298K
  parameters.Masse = (parameters.p_Std * 1.) / (parameters.R_S * 273) * 1.e-3; % in kg
  
  parameters.T_min = 273 - 200;
  parameters.T_Schritt = 100;
  parameters.T_max = 273 + 800;
  
  parameters.PlotSize = 3;
  
  parameters.Dimensions = "-S1600,1200";

  
  parameters.filename = "Isothermal_and_isentrop_process";

  parameters  
  
endfunction


%
%
% M A I N
%
%

%
% Derived from [[File:Animation Diagramm Carnotprozess.ogv]]
%

begin_time = time();
 
global parameters;

SetupParameters();


Volumen = 0:0.01:parameters.PlotSize;

close(); % Close previous plots
hold("on");

% Plot Isotherm
% Randbedingung ? konstante Temperatur

for temperatur = parameters.T_min:parameters.T_Schritt:parameters.T_max

  % p * V = n * R * T = m * R_S * T
  % p = m * R_S * T / V
  
  Masse_in_g = 1000. * parameters.Masse;
  Druck_isotherm = Masse_in_g  * parameters.R_S * temperatur ./ Volumen;
  
  % "spezifischer" Druck in N / (m^2 kg) ?
  % Menge, Druck, Volumen
  % 1 kg, 1 g ?

  plot(Volumen, Druck_isotherm/1e5, "color", "b");
  
endfor

% Plot Isentrop
% Randbedingung V = 5, Temperatur wie vorher ?

for temperatur = parameters.T_min:parameters.T_Schritt:parameters.T_max

  % Die Temperatur liefert den Bezugspunkt
  Volumen_ref = 1.; % in m^3
  Masse_in_g = 1000. * parameters.Masse;
  Druck_ref = Masse_in_g  * parameters.R_S * temperatur ./ Volumen_ref;
  
  Druck_isentrop = Druck_ref * (Volumen_ref ./ Volumen) .^ parameters.kappa;

  plot(Volumen, Druck_isentrop/1e5, "color", "r");

endfor


% Basis von 1 m^3

for temperatur = parameters.T_min:parameters.T_Schritt:parameters.T_max

  Volumen_ref = 1.; % in m^3
  Masse_in_g = 1000. * parameters.Masse;
  Druck_ref = Masse_in_g  * parameters.R_S * temperatur ./ Volumen_ref;

  plot(Volumen_ref, Druck_ref/1e5, "marker", "o", "markeredgecolor", "black", "markersize", 3.);
  text(Volumen_ref*1.03, Druck_ref/1e5, [num2str(temperatur-273) , " °C"], "color", "black");


endfor


title (["Argon m = ", num2str(parameters.Masse), "kg"]);
ylim([0,6]); % TODO magic number
xlim([0,3]);
xlabel("Volumen in m^3");
ylabel("Druck in 10^5 N / m");
Legend_Text(1,:) = "Isotherm";
Legend_Text(12,:) = "Isentrop";
legend(Legend_Text);

filename = [parameters.filename, ".svg"]
print (filename, "-dsvg"); % , parameters.Dimensions);

printf("Total processing time: %f seconds\n", time()-begin_time);
printf("Finished\n");