%clear all;
%close all;
%clc;

%% Solid fuel and air
kappa_H     = 1.42;
kappa_O     = 0.66;
kappa_N     = 79/21;
n_C         = 12.011*10^-3;
n_H         = 1.00797*10^-3;
n_O         = 15.999*10^-3;
n_N         = 14.007*10^-3;
n_H2O       = 2 * n_H + 1 * n_O;
n_CO2       = n_C + 2 * n_O;
n_CO        = n_C + n_O;
eta_cho     = 0.87;
eta_h2o     = 1-eta_cho;
gamma_1     = .12;
gamma_2     = 2.0;
volume2mass = 1.58;

%% Constants
sigma  = 5.67*10^-8;
chi    = 0.8;
m_w    = 0.35;
c_C    = (eta_cho*n_C)/(n_C+kappa_H*n_H+kappa_O*n_O);
c_H    = (eta_cho*kappa_H*n_H)/(n_C+kappa_H*n_H+kappa_O*n_O);
c_sf_O = (eta_cho*kappa_O*n_O)/(n_C+kappa_H*n_H+kappa_O*n_O);
c_N    = (kappa_N*n_N)/(2*kappa_N*n_N+2*n_O);
c_a_O  = (n_O)/(2*kappa_N*n_N+2*n_O);

%% Relation between duty cycle and flow
beta_motor   = 117.07*10^-6;
gamma_motor = 1000*beta_motor;
beta_blower  = 0.000266557533394774;
gamma_blower = 1000*beta_blower;

%% Density
rho_w   = 983.0167;
rho_cp  = 2800;
rho_eg  = 0.9482;
rho_gc1 = 0.8357;
rho_gc2 = 0.7471;
rho_bc  = 0.4564;

%% Heat transfer coefficients + Area:
alpha_A_bcw   = 45;
epsilon_A_bcw = 0.375;
alpha_A_gcw   = 18;
alpha_A_egw   = 7;

%% Masses
M_w     = 89;
M_bc    = ((0.305+0.048+0.03)*0.2*0.344)*rho_bc; 
M_gc    = ((0.305*0.344*(0.218+0.048))+(0.05*0.344*(0.0775+0.145/2))-(0.048*0.344*0.255))*rho_gc2+(0.065*0.344*(0.470+0.003+0.218+0.048))*rho_gc1;
M_eg    = ((0.115*0.344*(0.770-0.573))+(0.065*0.344*0.573))*rho_eg;
M_cp    = (((0.08085*0.03)+(0.0777*0.04))*rho_cp);
M_total = M_bc + M_gc + M_eg;

%% Specifik heat capacities
c_w  = 4185;
c_iw = 4185;
c_bc = 1067.8;
c_a  = 1005.2;
c_gc = 1034.8;
c_eg = 1017.4;
c_cp = 950;

%% Temperatures [C]
T_a  = 20;

%% Heat values [J/kg]
h_c_1 = 9203.2*1000;
h_c_2 = 32766*1000;
h_c_5 = 10103*1000;
h_c_6 = 141500*1000;

%% Operation points
T_w_bar      = 65;
T_bc_bar     = 420;
T_gc_bar     = 190;
T_eg_bar     = 132;
m_sf_bar     = 9*beta_motor;
m_a_bar      = 17*beta_blower;
M_C_char_bar = 0.0019123; 
O2_bar       = 0.0675;%.068175;
T_iw_bar     = 65;

%% Load non-linear simulation data.
load('.\non_linear_mdl_wheat.mat')

M_C_sim = [modelOutput.time modelOutput.signals.values(:,1)];
T_bc_sim = [modelOutput.time modelOutput.signals.values(:,2)];
T_gc_sim = [modelOutput.time modelOutput.signals.values(:,3)];
T_eg_sim = [modelOutput.time modelOutput.signals.values(:,4)];
T_w_sim = [modelOutput.time modelOutput.signals.values(:,5)];
O2_sim = [modelOutput.time modelOutput.signals.values(:,6)];

%% Save data
save constants_wheat_lin_mdl.mat

disp('Constans for the linear model using wheat are now loaded')