OpenWAM
TCCCilindro.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 // ---------------------------------------------------------------------------
30 // En la Tesis de Jose Miguel Corberan, paginas 104-139
31 // Resumido en el articulo "Solucion a la condicion de contorno de la union
32 // cilindro-conducto de los MCIA" (Esta en el COMETT PROGRAMME 1995)
33 // ---------------------------------------------------------------------------
34 // ---------------------------------------------------------------------------
35 #pragma hdrstop
36 
37 #include "TCCCilindro.h"
38 //#include <cmath>
39 #ifdef __BORLANDC__
40 #include <vcl.h>
41 #endif
42 #include "TTubo.h"
43 
44 #include "TValvula4T.h"
45 #include "TLamina.h"
46 #include "TLumbrera.h"
47 
48 #include "TBloqueMotor.h"
49 #include "TCilindro4T.h"
50 #include "Globales.h"
51 
52 // ---------------------------------------------------------------------------
53 // ---------------------------------------------------------------------------
54 
55 TCCCilindro::TCCCilindro(nmTypeBC TipoCC, int numCC, nmTipoCalculoEspecies SpeciesModel, int numeroespecies,
56  nmCalculoGamma GammaCalculation, bool ThereIsEGR) :
57  TCondicionContorno(TipoCC, numCC, SpeciesModel, numeroespecies, GammaCalculation, ThereIsEGR) {
58 
59  if(TipoCC == nmIntakeValve)
60  FTipoValv = nmValvAdmision;
61  else if(TipoCC == nmExhaustValve)
62  FTipoValv = nmValvEscape;
63  else
64  printf("ERROR en tipo de valvula TCCCilindro en la condicion de contorno: %d\n", FNumeroCC);
65 
66  FTuboExtremo = NULL;
67  FValvula = NULL;
68 
69  FGasto = 0.;
70  FVelocity = 0.;
71  FSonido = 1.;
72  FMomento = 0.;
73  FVelocidadGarganta = 0.;
74  FMachGarganta = 1.;
75  FGastoGarganta = 0.;
76  FRelacionPresionGarganta = 1.;
77 
78  FTime0 = 0.;
79  FTime1 = 0.;
80  FAnguloAnterior = 0.;
81  FAnguloActual = 0.;
82 
83 }
84 
85 // ---------------------------------------------------------------------------
86 // ---------------------------------------------------------------------------
87 
88 TCCCilindro::~TCCCilindro() {
89 
90  delete[] FTuboExtremo;
91 
92  if(FValvula != NULL)
93  delete FValvula;
94 
95  FValvula = NULL;
96 
97 }
98 
99 // ---------------------------------------------------------------------------
100 // ---------------------------------------------------------------------------
101 
102 void TCCCilindro::AsignaTipoValvula(TTipoValvula **Origen, int Valv, int i) {
103  try {
104 
105  switch(Origen[Valv - 1]->getTypeOfValve()) {
106  case nmValvula4T:
107  FValvula = new TValvula4T(dynamic_cast<TValvula4T*>(Origen[Valv - 1]), i);
108  break;
109  case nmLamina:
110  FValvula = new TLamina(dynamic_cast<TLamina*>(Origen[Valv - 1]), i);
111  break;
112  case nmLumbrera2T:
113  FValvula = new TLumbrera(dynamic_cast<TLumbrera*>(Origen[Valv - 1]), i);
114  break;
115  }
116 
117  // FValvula->putDiametroTubo(FTuboExtremo[0].Pipe->GetDiametro(FNodoFin));
118 
119  FValvula->PutPipe(FTuboExtremo[0].Pipe, FNodoFin);
120 
121  FSeccionValvula = __geom::Circle_area(FValvula->getDiametro());
122  FSeccionTubo = __geom::Circle_area(FTuboExtremo[0].Pipe->GetDiametro(FNodoFin));
123 
124  } catch(exception & N) {
125  std::cout << "ERROR: TCCCilindro::AsignaTipoValvula en la condicion de contorno: " << FNumeroCC << std::endl;
126  std::cout << "Tipo de error: " << N.what() << std::endl;
127  throw Exception(N.what());
128  }
129 }
130 
131 // ---------------------------------------------------------------------------
132 // ---------------------------------------------------------------------------
133 
134 void TCCCilindro::ReadBoundaryData(const char *FileWAM, fpos_t &filepos, int NumberOfPipes, TTubo **Pipe, int nDPF,
135  TDPF **DPF) {
136  try {
137  int i = 0;
138  int numid = 0; // Variable necesaria para WAMer.
139 
140  FTuboExtremo = new stTuboExtremo[1];
141  FTuboExtremo[0].Pipe = NULL;
142 
143  FPref = 1;
144 
145  while(FNumeroTubosCC < 1 && i < NumberOfPipes) {
146  if(Pipe[i]->getNodoIzq() == FNumeroCC) {
147  FTuboExtremo[FNumeroTubosCC].Pipe = Pipe[i];
148  FTuboExtremo[FNumeroTubosCC].TipoExtremo = nmLeft;
149  FNodoFin = 0;
150  FIndiceCC = 0;
151  FCC = &(FTuboExtremo[FNumeroTubosCC].Beta);
152  FCD = &(FTuboExtremo[FNumeroTubosCC].Landa);
153  FNumeroTubosCC++;
154  }
155  if(Pipe[i]->getNodoDer() == FNumeroCC) {
156  FTuboExtremo[FNumeroTubosCC].Pipe = Pipe[i];
157  FTuboExtremo[FNumeroTubosCC].TipoExtremo = nmRight;
158  FNodoFin = Pipe[i]->getNin() - 1;
159  FIndiceCC = 1;
160  FCC = &(FTuboExtremo[FNumeroTubosCC].Landa);
161  FCD = &(FTuboExtremo[FNumeroTubosCC].Beta);
162  FNumeroTubosCC++;
163  }
164  i++;
165  }
166  FILE *fich = fopen(FileWAM, "r");
167  fsetpos(fich, &filepos);
168 
169  fscanf(fich, "%d ",
170  &numid); // Esto es un dato que necesita el WAMer. Los usuarios de WAM hacemos la vista gorda hasta que se arregle.
171  fscanf(fich, "%d ", &FNumeroCilindro);
172 
173  fgetpos(fich, &filepos);
174  fclose(fich);
175 
176  // Inicializacion del transporte de especies quimicas
177  FFraccionMasicaEspecie = new double[FNumeroEspecies - FIntEGR];
178  for(int i = 0; i < FNumeroEspecies - FIntEGR; i++) {
179  FFraccionMasicaEspecie[i] = FTuboExtremo[0].Pipe->GetFraccionMasicaInicial(i);
180  }
181 
182  }
183 
184  catch(exception & N) {
185  std::cout << "ERROR: TCCCilindro::LeeCCCilindro en la condicion de contorno: " << FNumeroCC << std::endl;
186  std::cout << "Tipo de error: " << N.what() << std::endl;
187  throw Exception(N.what());
188  }
189 }
190 
191 // ---------------------------------------------------------------------------
192 // ---------------------------------------------------------------------------
193 
194 void TCCCilindro::AsignaCilindro(TBloqueMotor *EngineBlock) {
195  try {
196 
197  FMotor = EngineBlock;
198  FCilindro = EngineBlock->GetCilindro(FNumeroCilindro - 1);
199  FAnguloActual = FCilindro->getAnguloActual();
200  FValvula->PutCylider(FCilindro);
201 
202  } catch(exception & N) {
203  std::cout << "ERROR: TCCCilindro::AsignaCilindro en la condicion de contorno: " << FNumeroCC << std::endl;
204  std::cout << "Tipo de error: " << N.what() << std::endl;
205  throw Exception(N.what());
206  }
207 }
208 
209 // ---------------------------------------------------------------------------
210 // ---------------------------------------------------------------------------
211 
212 void TCCCilindro::CalculaCondicionContorno(double Time) {
213  try {
214  double rel_CCon_Entropia, coef, FraccionMasicaAcum = 0.;
215 
216  FTime0 = FTime1;
217  FTime1 = Time;
218 
219  FGamma = FTuboExtremo[0].Pipe->GetGamma(FNodoFin);
220  FGamma1 = __Gamma::G1(FGamma);
221  FGamma2 = __Gamma::G2(FGamma);
222  FGamma3 = __Gamma::G3(FGamma);
223  FGamma4 = __Gamma::G4(FGamma);
224  FGamma5 = __Gamma::G5(FGamma);
225  FGamma6 = __Gamma::G6(FGamma);
226 
227  // FSeccionValvula = Pi * pow(FValvula->getDiametro(),2.)/ 4.;
228  // FSeccionTubo=Pi*pow(FTuboExtremo[0].Pipe->GetDiametro(FNodoFin),2.)/4.;
229 
230  FAd = pow(FCilindro->getPressure() / FPref, 1. / FGamma4);
231  rel_CCon_Entropia = *FCC / FTuboExtremo[0].Entropia;
232  if(rel_CCon_Entropia / FAd > 1.000005) { // Flujo entrante al cilindro
233  FSentidoFlujo = nmEntrante;
234  FValvula->GetCDin(FTime1);
235  FCDEntrada = FValvula->getCDTubVol();
236  if(FCDEntrada > 0.0001) { /* Abierto */
237  FSeccionEficaz = FCDEntrada * FSeccionValvula;
238  FlujoEntranteCilindro();
239  /* CALCULO DEL MOMENTO ANGULAR ENTRANTE L */
240  if(FGasto < -1e-5) {
241  FCTorbellino = FValvula->getCTorb();
242  if(FTipoValv == nmValvEscape) {
243  coef = FCTorbellino / 4.;
244  } else {
245  coef = FCTorbellino;
246  }
247  FMomento = coef * pow2(__cons::ARef * FGasto) / (200000 * FMotor->getGeometria().Carrera * FGamma *
248  FCilindro->getPressure()) * (pow2(*FCC + *FCD) / 4. + FGamma3 * pow2(
249  ((*FCD - *FCC) / FGamma1)));
250  }
251  // Transporte de especies quimicas.
252  for(int j = 0; j < FNumeroEspecies - 2; j++) {
253  FFraccionMasicaEspecie[j] = FTuboExtremo[0].Pipe->GetFraccionMasicaCC(FIndiceCC, j);
254  FraccionMasicaAcum += FFraccionMasicaEspecie[j];
255  }
256  FFraccionMasicaEspecie[FNumeroEspecies - 2] = 1. - FraccionMasicaAcum;
257  if(FHayEGR)
258  FFraccionMasicaEspecie[FNumeroEspecies - 1] = FTuboExtremo[0].Pipe->GetFraccionMasicaCC(FIndiceCC, FNumeroEspecies - 1);
259  } else { /* Cerrado */
260  FMomento = 0.;
261  FGasto = 0.;
262  *FCD = *FCC;
263  FVelocity = 0.;
264  FSonido = *FCD;
265  FVelocidadGarganta = 0.;
266  FMachGarganta = 0.;
267  FGastoGarganta = 0.;
268  FRelacionPresionGarganta = 0.;
269  FSeccionEficaz = 0.;
270  // La composicion se mantiene, al estar el cilindro cerrado.
271  }
272  } else if(rel_CCon_Entropia / FAd < .999995) { // Flujo saliente del cilindro
273  FSentidoFlujo = nmSaliente;
274  FValvula->GetCDout(FTime1);
275  FCDSalida = FValvula->getCDVolTub();
276  if(FCDSalida > 0.0001) { /* Abierto */
277  FSeccionEficaz = FCDSalida * FSeccionValvula;
278  FlujoSalienteCilindro();
279  /* CALCULO DEL MOMENTO ANGULAR SALIENTE */
280  if(FGasto > 1e-5) {
281  FMomento = -FCilindro->getMomentoAngular() * FGasto / FCilindro->getMasa();
282  }
283 
284  // Transporte de especies quimicas.
285  for(int j = 0; j < FNumeroEspecies - 2; j++) {
286  FFraccionMasicaEspecie[j] = FCilindro->GetComposicionSaliente(j);
287  FraccionMasicaAcum += FFraccionMasicaEspecie[j];
288  }
289  FFraccionMasicaEspecie[FNumeroEspecies - 2] = 1. - FraccionMasicaAcum;
290  if(FHayEGR)
291  FFraccionMasicaEspecie[FNumeroEspecies - 1] = FCilindro->GetComposicionSaliente(FNumeroEspecies - 1);
292  } else { /* Cerrado */
293  FMomento = 0.;
294  FGasto = 0.;
295  *FCD = *FCC;
296  FVelocity = 0.;
297  FSonido = *FCD;
298  FVelocidadGarganta = 0.;
299  FMachGarganta = 0.;
300  FGastoGarganta = 0.;
301  FRelacionPresionGarganta = 0.;
302  FSeccionEficaz = 0.;
303  // La composicion se mantiene, al estar el cilindro cerrado.
304  }
305  } else { // Flujo Parado
306  FSentidoFlujo = nmParado;
307  *FCD = *FCC;
308  FGasto = 0.;
309  FVelocity = 0.;
310  FSonido = *FCD;
311  FVelocidadGarganta = 0.;
312  FMachGarganta = 0.;
313  FGastoGarganta = 0.;
314  FRelacionPresionGarganta = 0.;
315  FSeccionEficaz = 0.;
316  // La composicion se mantiene, al estar el flujo parado.
317  }
318  FValvula->AcumulaCDMedio(Time);
319  } catch(exception & N) {
320  std::cout << "ERROR: TCCCilindro::CalculaCondicionContorno en la condicion de contorno: " << FNumeroCC << std::endl;
321  std::cout << "Tipo de error: " << N.what() << std::endl;
322  throw Exception(N.what());
323  }
324 }
325 
326 // ---------------------------------------------------------------------------
327 // ---------------------------------------------------------------------------
328 
329 void TCCCilindro::FlujoEntranteCilindro() {
330  try {
331  double vel_son_garganta = 0., velocidad_garganta = 0., Mach = 0., xaa2 = 0., ycal = 0., d1 = 0.;
332 
333  Fk = FSeccionTubo / FSeccionEficaz;
334  if(Fk < 1)
335  Fk = 1;
336  vel_son_garganta = FTuboExtremo[0].Entropia * FAd; // Velocity del sonido en la garganta. Adimensional.
337 
338  // Calculo de la velocidad en la garganta.Caso de salto subcritico.
339  FCaso = nmFlujoEntranteSaltoSubcritico;
340  if(Fk == 1) {
341  FSonido = FTuboExtremo[0].Entropia * FAd;
342  FVelocity = (*FCC - FSonido) / FGamma3;
343  } else
344  Resolucion(vel_son_garganta, *FCC, FCaso, &FVelocity, &FSonido);
345 
346  // Ecuacion de la energia
347  velocidad_garganta = sqrt(2. * FGamma6 * (pow2(FSonido) + FGamma3 * pow2(FVelocity) - pow2(vel_son_garganta)));
348  // Se ha calculado la velocidad en la garganta en valor absoluto.
349 
350  // Calculo de la velocidad en la garganta en el caso de salto supercritico
351  if(velocidad_garganta > vel_son_garganta) {
352  FCaso = nmFlujoEntranteSaltoSupercritico;
353  Resolucion(0.0, 1.0, FCaso, &ycal, &Mach);
354  FVelocity = *FCC / (1 / Mach + FGamma3);
355  FSonido = FVelocity / Mach;
356 
357  d1 = Fk * FVelocity * pow(FSonido, 1. / FGamma3);
358  vel_son_garganta = pow(d1, FGamma1 / FGamma2);
359  velocidad_garganta = vel_son_garganta;
360  }
361  double Ga3U = FGamma3 * FVelocity;
362  // Fin caso de salto supercritico
363 
364  xaa2 = pow(FTuboExtremo[0].Entropia, FGamma4);
365  FGasto = __units::BarToPa(-FGamma * FSeccionTubo * pow(FSonido,
366  2 * FGamma6) * FVelocity) / (__cons::ARef * xaa2); // Massflow entrante al cilindro negativo
367  *FCD = FSonido - Ga3U;
368  *FCC = FSonido + Ga3U;
369  FRelacionPresionGarganta = pow(FSonido / (FTuboExtremo[0].Entropia * FAd), FGamma4);
370  FGastoGarganta = FGasto / (FCDEntrada * FSeccionValvula);
371  FMachGarganta = -velocidad_garganta / vel_son_garganta; // Negativo por ser flujo entrante
372  FVelocidadGarganta = velocidad_garganta;
373  }
374 
375  catch(exception & N) {
376  std::cout << "ERROR: TCCCilindro::FlujoEntranteCilindro en la condicion de contorno: " << FNumeroCC << std::endl;
377  std::cout << "Tipo de error: " << N.what() << std::endl;
378  throw Exception(N.what());
379  }
380 }
381 
382 // ---------------------------------------------------------------------------
383 // ---------------------------------------------------------------------------
384 
385 void TCCCilindro::FlujoSalienteCilindro() {
386  try {
387 
388  double a1 = 0., xx = 0., yy = 0., d1 = 0., Ga3U = 0.;
389  double a2cr = 0., val1 = 0., val2 = 0., u2cr = 0., ycal = 0., error = 0., valde = 0., miembro2 = 0.;
390  // Variables para resolver la onda de choque.
391  double relacion_velocidades_son, Mach_tras_ondachoque, Mach, temp_antes_ondachoque, temp_tras_ondachoque;
392  double root_a = 0.;
393 
394  Fk = FSeccionTubo / FSeccionEficaz;
395  if(Fk < 1)
396  Fk = 1.0;
397 
398  double sqrtGa2 = sqrt(2. / FGamma2);
399 
400  /* Calculo del valor de la velocidad del sonido en el extremo del tubo para
401  el cual el salto es critico. */
402  u2cr = FCilindro->getSpeedsound() / __cons::ARef * sqrtGa2 * (sqrt(pow2(Fk) + FGamma1 * FGamma2) - Fk) / FGamma1;
403  a2cr = sqrt(pow2(FCilindro->getSpeedsound() / __cons::ARef) - FGamma3 * pow2(u2cr));
404  // Ecuacion de la energia. Garganta-Cylinder.
405 
406  /* A partir de a2cr se determina el error en el calculo de A2 al suponer salto
407  subcritico. Si es negativo, el salto es supercritico. Si es positivo, el salto
408  es subcritico. */
409  // FSSubcritico(a2cr,&error,&miembro2);
410  stFSSub FSA2(FTuboExtremo[0].Entropia, FAd, FGamma, Fk, *FCC, FCilindro->getSpeedsound() / __cons::ARef);
411 
412  error = FSA2(a2cr);
413 
414  if(error < 0.) { // Salto de presiones supercritico.
415 
416  /* Determinacion del intervalo de iteracion. Para ello se supone que
417  en el extremo del tubo se dan las condiciones criticas. Explicado en
418  los apuntes de Pedro. */
419  a1 = sqrtGa2 * FCilindro->getSpeedsound() / __cons::ARef;
420  FVelocidadGarganta = a1;
421  xx = pow(FAd / (FCilindro->getSpeedsound() / __cons::ARef), FGamma4);
422  yy = pow(a1, 2. / FGamma1);
423  Fcc = FVelocidadGarganta * yy * xx / Fk;
424  // FSSupercritico(FVelocidadGarganta,&val1,&val2);
425 
426  stFSSup FU2(FTuboExtremo[0].Entropia, Fcc, FGamma, Fk, *FCC, FCilindro->getSpeedsound() / __cons::ARef);
427  val1 = FU2(FVelocidadGarganta);
428 
429  if(val1 < 0.)
430  valde = FVelocidadGarganta;
431  else
432  valde = (FCilindro->getSpeedsound() / __cons::ARef) / sqrt(FGamma3);
433 
434  /* Una vez conocido el intervalo de iteracion, se pasa a la resolucion
435  del caso flujo saliente salto supercritico. */
436  FCaso = nmFlujoSalienteSaltoSupercritico;
437  Resolucion(0.0, valde, FCaso, &FVelocity, &FSonido);
438  Ga3U = FVelocity * FGamma3;
439  // Calcula del massflow. Como es saliente del cilindro, siempre es positivo.
440  xx = pow(sqrtGa2, (FGamma2 / FGamma1));
441  yy = pow(FAd, FGamma4);
442  FGasto = __units::BarToPa(FCDSalida * FSeccionValvula * FGamma * xx * yy) / (FCilindro->getSpeedsound());
443 
444  /* Reduccion a flujo subsonico mediante onda de choque plana en el caso
445  de que se hayan obtenido condiciones supersonicas en el extremo del
446  tubo. Explicado en la tesis Corberan (pagina de la 47 a la 52
447  (punto 2.5) y de la 122 a la 129 (lo importante a partir de la 127) */
448  Mach = FVelocity / FSonido;
449  xx = *FCC + Ga3U;
450  FTuboExtremo[0].Entropia = FTuboExtremo[0].Entropia * FSonido / xx;
451  // Ecuacion de la caracteristica incidente.
452  if(Mach > 1.) {
453 
454  /* Las ecuaciones siguientes corresponden a la resolucion de la onda
455  de choque plana. Se pueden encontrar en el punto 2.5 de la tesis
456  de Corberan. */
457  xx = FGamma4 * pow2(Mach) - 1.;
458  Mach_tras_ondachoque = sqrt((pow2(Mach) + 2. / FGamma1) / xx);
459  temp_tras_ondachoque = FGamma3 * pow2(Mach) + 1.;
460  temp_antes_ondachoque = FGamma3 * pow2(Mach_tras_ondachoque) + 1.;
461  relacion_velocidades_son = sqrt(temp_tras_ondachoque / temp_antes_ondachoque);
462  FSonido = FSonido * relacion_velocidades_son;
463  FVelocity = FSonido * Mach_tras_ondachoque;
464  Ga3U = FVelocity * FGamma3;
465  d1 = xx * FGamma1 / FGamma2;
466  FTuboExtremo[0].Entropia = FTuboExtremo[0].Entropia * relacion_velocidades_son / pow(d1, FGamma5);
467  }
468  } else { // Salto de presiones subcritico.
469 
470  // Resolucion del caso de flujo saliente salto subcritico.
471  FCaso = nmFlujoSalienteSaltoSubcritico;
472  Resolucion(a2cr, FCilindro->getSpeedsound() / __cons::ARef, FCaso, &ycal, &FSonido);
473  // Aplicando la Ecuacion de la Energia entre el cilindro y la garganta:
474  root_a = pow2(FCilindro->getSpeedsound() / __cons::ARef) - pow2(FSonido);
475  if(root_a > 0) {
476  FVelocity = sqrt((pow2(FCilindro->getSpeedsound() / __cons::ARef) - pow2(FSonido)) / FGamma3);
477  } else if(root_a > -1e12) {
478  FVelocity = 0;
479  } else {
480  FVelocity = 0.;
481  printf("ERROR: Calculating outflow in boundary %d", FNumeroCC);
482  }
483  Ga3U = FVelocity * FGamma3;
484  // Calculo del massflow. Como es saliente del cilindro, siempre es positivo.
485  xx = *FCC + Ga3U;
486  a1 = FCilindro->getSpeedsound() / __cons::ARef * xx / (FTuboExtremo[0].Entropia * FAd);
487  FVelocidadGarganta = Fk * pow2(a1) * FVelocity / pow2(FSonido);
488  FGasto = __units::BarToPa(FCDSalida * FSeccionValvula * FGamma * pow(FAd / (FCilindro->getSpeedsound() / __cons::ARef),
489  FGamma4) * FVelocidadGarganta * pow(a1, 2. / FGamma1)) / __cons::ARef;
490 
491  FTuboExtremo[0].Entropia = FTuboExtremo[0].Entropia * FSonido / xx;
492  // Ecuacion de la caracteristica incidente.
493  }
494  *FCD = FSonido + Ga3U;
495  *FCC = FSonido - Ga3U;
496  d1 = FSonido / (FTuboExtremo[0].Entropia * FAd);
497  FRelacionPresionGarganta = pow(d1, FGamma4);
498  FMachGarganta = FVelocidadGarganta / a1;
499  // FMachGarganta = FVelocidadGarganta/a1; // Positivo por ser flujo saliente
500  FGastoGarganta = FGasto / (FCDSalida * FSeccionValvula);
501  }
502 
503  catch(exception & N) {
504  std::cout << "ERROR: TCCCilindro::FlujoSalienteCilindro en la condicion de contorno: " << FNumeroCC << std::endl;
505  std::cout << "Tipo de error: " << N.what() << std::endl;
506  throw Exception(N.what());
507  }
508 }
509 
510 // ---------------------------------------------------------------------------
511 // ---------------------------------------------------------------------------
512 
513 void TCCCilindro::Resolucion(double ext1, double ext2, nmCaso Caso, double *u2t, double *a2t) {
514  try {
515  if(Caso == nmFlujoEntranteSaltoSubcritico) {
516  stFESub FEA2(FTuboExtremo[0].Entropia, FAd, FGamma, Fk, *FCC);
517  *a2t = FindRoot(FEA2, ext1, ext2);
518  *u2t = FEA2.U2;
519  } else if(Caso == nmFlujoEntranteSaltoSupercritico) {
520  stFESup FMatch(FGamma, Fk);
521  *a2t = FindRoot(FMatch, ext1, ext2);
522  *u2t = 0.;
523  } else if(Caso == nmFlujoSalienteSaltoSubcritico) {
524  stFSSub FSA2(FTuboExtremo[0].Entropia, FAd, FGamma, Fk, *FCC, FCilindro->getSpeedsound() / __cons::ARef);
525  *a2t = FindRoot(FSA2, ext1, ext2);
526  *u2t = FSA2.U2;
527  } else if(Caso == nmFlujoSalienteSaltoSupercritico) {
528  stFSSup FU2(FTuboExtremo[0].Entropia, Fcc, FGamma, Fk, *FCC, FCilindro->getSpeedsound() / __cons::ARef);
529  *u2t = FindRoot(FU2, ext1, ext2);
530  *a2t = FU2.A2;
531  } else {
532  printf("Error en la definicion del flujo TCCDeposito::Resolucion en la condicion de contorno: %d\n", FNumeroCC);
533  throw Exception("");
534  }
535  }
536 
537  catch(exception & N) {
538  std::cout << "ERROR: TCCCilindro::Resolucion en la condicion de contorno: " << FNumeroCC << std::endl;
539  std::cout << "Tipo de error: " << N.what() << std::endl;
540  throw Exception(N.what());
541  }
542 }
543 
544 // ---------------------------------------------------------------------------
545 // ---------------------------------------------------------------------------
546 
547 // void TCCCilindro::FESubcritico(double vel_son_supuesta,double *u2_1,double *u2_2)
548 // {
549 // try
550 // {
551 //
552 // double xx,yy,u2;
553 //
555 // contorno de la union cilindro-conducto de los MCIA". Ecuacion 4.30 en la
556 // tesis de Corberan */
557 //
558 // xx=vel_son_supuesta/(FTuboExtremo[0].Entropia*FAd);
559 // yy=pow(xx,4.*FGamma6);
560 // yy=pow(Fk,2.)*yy-1.;
561 // *u2_2=FTuboExtremo[0].Entropia*FAd*sqrt(2.*FGamma6*(pow(xx,2.)-1.)/yy); // Valor absoluto
562 //
564 //
565 // *u2_1=(*FCC-vel_son_supuesta)/FGamma3; // En valor absoluto
566 //
567 // }
568 // catch(Exception &N)
569 // {
570 // std::cout << "ERROR: TCCCilindro::FESubcritico en la condicion de contorno: " << FNumeroCC << std::endl;
571 // std::cout << "Tipo de error: " << N.what() << std::endl;
572 // throw Exception(N.what());
573 // }
574 // }
575 //
578 //
579 // void TCCCilindro::FESupercritico(double mach_supuesto,double *miembro1,
580 // double *miembro2)
581 // {
582 // try
583 // {
584 //
585 // double xx,yy;
586 //
588 // contorno de la union cilindro-conducto de los MCIA". Ecuacion (4.31) de
589 // la tesis de Corberan */
590 //
591 // yy=(FGamma2/2.)*pow(Fk*mach_supuesto,2.*FGamma1/FGamma2);
592 // xx=FGamma3*pow(mach_supuesto,2);
593 // *miembro1=xx-yy+1.; // Miembro 1 de la ecuacion (21)
594 // *miembro2=0; // Miembro 2 de la ecuacion (21)
595 //
596 // }
597 // catch(Exception &N)
598 // {
599 // std::cout << "ERROR: TCCCilindro::FESupercritico en la condicion de contorno: " << FNumeroCC << std::endl;
600 // std::cout << "Tipo de error: " << N.what() << std::endl;
601 // throw Exception(N.what());
602 // }
603 // }
604 //
607 //
608 // void TCCCilindro::FSSubcritico(double vel_son_supuesta,double *error,double *miembro2)
609 // {
610 // try
611 // {
612 // double a1,u1,u2;
613 //
614 // *miembro2=0;
615 //
617 // de Corberan. */
618 //
619 // u2 = sqrt((pow(FCilindro->getSpeedsound()/__cons::ARef,2)-pow(vel_son_supuesta,2))/FGamma3);
620 // a1 = FCilindro->getSpeedsound()/__cons::ARef*(*FCC+FGamma3*u2)/(FTuboExtremo[0].Entropia*FAd);
621 // u1 = Fk*u2*pow(a1,2)/pow(vel_son_supuesta,2);
622 // *error=pow(a1,2)+FGamma3*pow(u1,2)-pow(FCilindro->getSpeedsound()/__cons::ARef,2);
623 //
624 // }
625 // catch(Exception &N)
626 // {
627 // std::cout << "ERROR: TCCCilindro::FSSubcritico en la condicion de contorno: " << FNumeroCC << std::endl;
628 // std::cout << "Tipo de error: " << N.what() << std::endl;
629 // throw Exception(N.what());
630 // }
631 // }
632 //
635 //
636 // void TCCCilindro::FSSupercritico(double vel_supuesta,double *a2_1,double *a2_2)
637 // {
638 // try
639 // {
640 //
642 // *a2_1 = sqrt(pow(FCilindro->getSpeedsound()/__cons::ARef,2)-FGamma3*pow(vel_supuesta,2));
643 //
645 // *a2_2 = sqrt(vel_supuesta*pow((*FCC+FGamma3*vel_supuesta)/
646 // FTuboExtremo[0].Entropia,FGamma4)/Fcc);
647 //
648 // }
649 // catch(Exception &N)
650 // {
651 // std::cout << "ERROR: TCCCilindro::FSSupercritico en la condicion de contorno: " << FNumeroCC << std::endl;
652 // std::cout << "Tipo de error: " << N.what() << std::endl;
653 // throw Exception(N.what());
654 // }
655 // }
656 
657 // ---------------------------------------------------------------------------
658 // ---------------------------------------------------------------------------
659 
660 void TCCCilindro::ActualizaAnguloValvula(double TiempoActual, double Regimen) {
661  try {
662 
663  FTime0 = FTime1;
664  FTime1 = TiempoActual;
665  FDeltaT = FTime1 - FTime0;
666  FDeltaAngulo = 360. * Regimen / 60. * FDeltaT;
667  FAnguloAnterior = FAnguloActual;
668  FAnguloActual = FAnguloAnterior + FDeltaAngulo;
669  if(FAnguloActual > 720.) {
670  FAnguloActual -= 720.;
671  }
672 
673  } catch(exception & N) {
674  std::cout << "ERROR: TCCCilindro::ActualizaAnguloValvula en la condicion de contorno: " << FNumeroCC << std::endl;
675  std::cout << "Tipo de error: " << N.what() << std::endl;
676  throw Exception(N.what());
677  }
678 }
679 
680 #pragma package(smart_init)
TTubo
a Finite differences pipe.
Definition: TTubo.h:116
TTubo::GetGamma
double GetGamma(int i) const
Gets the specific heat capacities ratio at a given cell.
Definition: TTubo.cpp:5444
FindRoot
double FindRoot(T &func, const double x1, const double x2)
Finds the root of a function.
Definition: Math_wam.h:459
TCCCilindro::ActualizaAnguloValvula
void ActualizaAnguloValvula(double TiempoActual, double Regimen)
Definition: TCCCilindro.cpp:660
stFESub
Definition: BoundaryFunctions.h:37
stTuboExtremo
Definition: Globales.h:730
TTipoValvula
Definition: TTipoValvula.h:53
TDPF
Definition: TDPF.h:45
TTubo::GetFraccionMasicaCC
double GetFraccionMasicaCC(int j, int i)
Definition: TTubo.h:953
TLamina
Definition: TLamina.h:50
TTubo::GetFraccionMasicaInicial
double GetFraccionMasicaInicial(int i) const
Gets the initial mass fraction of species i.
Definition: TTubo.cpp:5440
TCondicionContorno
Definition: TCondicionContorno.h:54
stFESup
Definition: BoundaryFunctions.h:70
TValvula4T
Definition: TValvula4T.h:45
stFSSub
Definition: BoundaryFunctions.h:93
Exception
Custom exception class.
Definition: Exception.hpp:39
TLumbrera
Definition: TLumbrera.h:45
TTubo.h
TBloqueMotor
Definition: TBloqueMotor.h:43
pow2
T pow2(T x)
Returns x to the power of 2.
Definition: Math_wam.h:88
stFSSup
Definition: BoundaryFunctions.h:125