OpenWAM
TCompresorDep.cpp
1 /* --------------------------------------------------------------------------------*\
2 ==========================|
3  \\ /\ /\ // O pen | OpenWAM: The Open Source 1D Gas-Dynamic Code
4  \\ | X | // W ave |
5  \\ \/_\/ // A ction | CMT-Motores Termicos / Universidad Politecnica Valencia
6  \\/ \// M odel |
7  ----------------------------------------------------------------------------------
8  License
9 
10  This file is part of OpenWAM.
11 
12  OpenWAM is free software: you can redistribute it and/or modify
13  it under the terms of the GNU General Public License as published by
14  the Free Software Foundation, either version 3 of the License, or
15  (at your option) any later version.
16 
17  OpenWAM is distributed in the hope that it will be useful,
18  but WITHOUT ANY WARRANTY; without even the implied warranty of
19  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20  GNU General Public License for more details.
21 
22  You should have received a copy of the GNU General Public License
23  along with OpenWAM. If not, see <http://www.gnu.org/licenses/>.
24 
25 
26  \*-------------------------------------------------------------------------------- */
27 
28 // ---------------------------------------------------------------------------
29 #pragma hdrstop
30 #include <iostream>
31 #ifdef __BORLANDC__
32 #include <vcl.h>
33 #endif
34 //#include <cmath>
35 #include "TCompresorDep.h"
36 
37 #include "TDeposito.h"
38 // ---------------------------------------------------------------------------
39 // ---------------------------------------------------------------------------
40 
41 TCompresorDep::TCompresorDep(int i, nmTipoCalculoEspecies SpeciesModel, int numeroespecies,
42  nmCalculoGamma GammaCalculation, bool ThereIsEGR) :
43  TCompresor(i, SpeciesModel, numeroespecies, GammaCalculation, ThereIsEGR) {
44 
45  FModeloCompresor = nmCompPlenums;
46  Mapa = new TMapaComp(FNumeroCompresor);
47 }
48 
49 // ---------------------------------------------------------------------------
50 // ---------------------------------------------------------------------------
51 
52 TCompresorDep::~TCompresorDep() {
53 }
54 
55 // ---------------------------------------------------------------------------
56 // ---------------------------------------------------------------------------
57 
58 void TCompresorDep::LeeCompresor(const char *FileWAM, fpos_t &filepos) {
59  int modelo = 0;
60  try {
61  FILE *fich = fopen(FileWAM, "r");
62  fsetpos(fich, &filepos);
63 
64  Mapa->LeeMapa(fich);
65 
66  fgetpos(fich, &filepos);
67  fclose(fich);
68  } catch(exception & N) {
69  std::cout << "ERROR: TCompresorDep::LeeCompresor en el compresor: " << FNumeroCompresor << std::endl;
70  std::cout << "Tipo de error: " << N.what() << std::endl;
71  throw Exception("ERROR: LeeCompresor en el compresor: " + std::to_string(FNumeroCompresor) + N.what());
72  }
73 }
74 
75 // ---------------------------------------------------------------------------
76 // ---------------------------------------------------------------------------
77 
78 void TCompresorDep::RelacionDepositoCompresor(TDeposito *DepositoRot, TDeposito *DepositoEst) {
79  try {
80 
81  FDepositoRot = DepositoRot;
82  FDepositoEst = DepositoEst;
83 
84  FDepositoRot->AsignaCompresor(this, -1);
85  FDepositoEst->AsignaCompresor(this, 1);
86 
87  } catch(exception & N) {
88  std::cout << "ERROR: TCompTubos::RelacionTubos en el compresor: " << FNumeroCompresor << std::endl;
89  std::cout << "Tipo de error: " << N.what() << std::endl;
90  throw Exception("ERROR: LeeCompresor en el compresor: " + std::to_string(FNumeroCompresor) + N.what());
91  }
92 }
93 
94 void TCompresorDep::Initialize() {
95 
96  FTemperatura10 = __units::degCToK(FDepositoRot->getTemperature());
97 
98 }
99 
100 // ---------------------------------------------------------------------------
101 // ---------------------------------------------------------------------------
102 
103 // ---------------------------------------------------------------------------
104 // ----------------- CalculaGasto -----------------------------------------//
105 // Calculo del flujo que transcurre a travas de los 2 depositos del //
106 // compresor //
107 // ---------------------------------------------------------------------------
108 
109 void TCompresorDep::CalculaGasto(double TrabajoInsTurbina, double TiempoActual) {
110  int NumGastos = 0;
111  double Aux0 = 0.;
112  double Aux1 = 0.;
113  double Masa0 = 0.;
114  double Masa1 = 0.;
115  double MasaX = 0.;
116  double RC0 = 0.;
117  // double RC1 = 0.;
118  double RCX = 0.;
119  double ErrorMasa = 0.;
120  double ErrorRC = 0.;
121  double *RelComp;
122  double *Gastos;
123  double *Trab;
124  double *Rnd;
125  double *DeltaWork;
126  double Min = 0.;
127  int Contador = 0;
128  int *PuntosDeCruce;
129  int PuntoMin = 0;
130  double DescorrigeGasto = 0.;
131  try {
132 
133  // NO ESTAN PUESTAS LAS ESPECIES !!!!!!! SOLO ESTA HECHO PARA QUE COMPILE.
134  FGamma = FDepositoRot->getGamma(); // Se ha de hacer una media de entrada y salida
135  FRMezcla = FDepositoRot->getR();
136  FGamma4 = __Gamma::G4(FGamma);
137  FGamma5 = __Gamma::G5(FGamma);
138 
139  FPresion10 = FDepositoRot->getPressure();
140  FPresion20 = FDepositoEst->getPressure();
141  FTemperatura10 = pow2(FDepositoRot->getSpeedsound() * __cons::ARef) / FRMezcla / FGamma;
142  FRelacionCompresion = FPresion20 / FPresion10;
143  FDeltaTiempo = TiempoActual - FTiempo0;
144  FTiempo0 = TiempoActual;
145  DescorrigeGasto = __units::BarToPa(FPresion10) / Mapa->getPresionRef() / sqrt(FTemperatura10 / Mapa->getTempRef());
146 
147  if(FRelacionCompresion <= 1.) {
148  FGastoCorregido = Mapa->getGastoRelComp1();
149  FGastoCompresor = Mapa->getGastoRelComp1() * DescorrigeGasto;
150  FRendimiento = Mapa->EvaluaRendSplines(Mapa->getGastoRelComp1());
151  FTrabajo = FGastoCompresor * FRMezcla * FGamma4 / 2 * FTemperatura10 * (pow(FRelacionCompresion,
152  2. * FGamma5) - 1.) / FRendimiento;
153  FTempGasto = FTemperatura10
154  /* +(pow(FRelacionCompresion,Gamma5*2)-1)*FTemperatura10/FRendimiento */;
155  } else {
156  PuntosDeCruce = new int[Mapa->getNumPuntos()];
157  Aux0 = Mapa->GetRelCompInt(0) - FRelacionCompresion;
158  NumGastos = 0;
159  for(int i = 0; i < Mapa->getNumPuntos(); i++) {
160  Aux1 = Mapa->GetRelCompInt(i + 1) - FRelacionCompresion;
161  if(Aux0 * Aux1 < 0) {
162  PuntosDeCruce[NumGastos] = i;
163  ++NumGastos;
164  }
165  Aux0 = Aux1;
166  }
167  if(NumGastos == 0) {
168  FRendimiento = Mapa->EvaluaRendSplines(Mapa->getGastoBombeo());
169  FTrabajo = Mapa->getGastoBombeo() * DescorrigeGasto * FRMezcla * FGamma4 / 2 * FTemperatura10 * (pow(
170  Mapa->getRelCompBombeo(), 2. * FGamma5) - 1.) / FRendimiento;
171  FTempGasto = FTemperatura10 + (pow(FRelacionCompresion, FGamma5 * 2) - 1) * FTemperatura10 / FRendimiento;
172  FGastoCompresor = 0.;
173  } else {
174  RelComp = new double[NumGastos];
175  Gastos = new double[NumGastos];
176  Trab = new double[NumGastos];
177  Rnd = new double[NumGastos];
178  for(int i = 0; i < NumGastos; i++) {
179  ErrorMasa = 1.;
180  ErrorRC = 1.;
181  Contador = 0;
182  Masa0 = Mapa->GetGastoInt(PuntosDeCruce[0]);
183  RC0 = Mapa->GetRelCompInt(PuntosDeCruce[0]);
184  Masa1 = Mapa->GetGastoInt(PuntosDeCruce[0] + 1);
185  // RC1=Mapa->GetRelCompInt(PuntosDeCruce[0]+1);
186  while(ErrorMasa > 1e-5 && ErrorRC > 1e-5 && Contador < 100) {
187  MasaX = (Masa0 + Masa1) / 2.;
188  RCX = Mapa->EvaluaRCHermite(MasaX);
189  if((RC0 - FRelacionCompresion) * (RCX - FRelacionCompresion) < 0) {
190  Masa1 = MasaX;
191  // RC1=RCX;
192  ErrorMasa = fabs(Masa0 - MasaX);
193  } else {
194  ErrorMasa = fabs(Masa0 - MasaX);
195  Masa0 = MasaX;
196  RC0 = RCX;
197  }
198  ErrorRC = fabs(RCX - FRelacionCompresion);
199  ++Contador;
200  }
201  RelComp[i] = RCX;
202  Gastos[i] = MasaX;
203  Rnd[i] = Mapa->EvaluaRendSplines(Gastos[i]);
204  Trab[i] = Gastos[i] * DescorrigeGasto * FRMezcla * FGamma4 / 2. * FTemperatura10 * (pow(RelComp[i],
205  FGamma5 * 2.) - 1.) / Rnd[i];
206  }
207  if(NumGastos > 1) {
208  DeltaWork = new double[NumGastos];
209  for(int i = 0; i < NumGastos; i++) {
210  DeltaWork[i] = fabs(Trab[i] - TrabajoInsTurbina);
211  }
212  PuntoMin = 0;
213  Min = DeltaWork[0];
214  for(int j = 1; j < NumGastos; j++) {
215  if(DeltaWork[j] < Min) {
216  Min = DeltaWork[j];
217  PuntoMin = j;
218  }
219  }
220  FGastoCorregido = Gastos[PuntoMin];
221  FGastoCompresor = Gastos[PuntoMin] * DescorrigeGasto;
222  FTrabajo = Trab[PuntoMin];
223  FTempGasto = FTemperatura10 + (pow(RelComp[PuntoMin], FGamma5 * 2) - 1) * FTemperatura10 / Rnd[PuntoMin];
224  FRendimiento = Rnd[PuntoMin];
225  } else {
226  FGastoCorregido = Gastos[0];
227  FGastoCompresor = Gastos[0] * DescorrigeGasto;
228  FTrabajo = Trab[0];
229  FTempGasto = FTemperatura10 + (pow(RelComp[0], FGamma5 * 2) - 1) * FTemperatura10 / Rnd[0];
230  FRendimiento = Rnd[0];
231  }
232  }
233  }
234  FASonidoSalida = sqrt(FGamma * FRMezcla * FTempGasto);
235  FTrabajo *= FDeltaTiempo;
236  FPotencia = FTrabajo / FDeltaTiempo;
237  FTrabajoPaso += FTrabajo;
238  FDeltaTPaso += FDeltaTiempo;
239  FRegimenCorregido = Mapa->getRegimenCorregido();
240 
241  } catch(exception & N) {
242  std::cout << "ERROR: TCompresorDep::CalculaGasto en el compresor: " << FNumeroCompresor << std::endl;
243  std::cout << "Tipo de error: " << N.what() << std::endl;
244  throw Exception("ERROR: CalculaGasto en el compresor: " + std::to_string(FNumeroCompresor) + N.what());
245  }
246 }
247 
248 // ---------------------------------------------------------------------------
249 // ---------------------------------------------------------------------------
250 
251 #pragma package(smart_init)
TMapaComp
Definition: TMapaComp.h:40
TCompresor
Definition: TCompresor.h:47
TDeposito
Definition: TDeposito.h:44
Exception
Custom exception class.
Definition: Exception.hpp:39
pow2
T pow2(T x)
Returns x to the power of 2.
Definition: Math_wam.h:88