Jet/cooling.c at master · chomps/Jet · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63

#include "paul.h"

//static double r0;

double get_dV( double * , double * );

void setCoolingParams( struct domain * theDomain ){

//   r0   = theDomain->theParList.Nozzle_r0;

}

void cool_src( double * prim , double * cons , double dVdt ){

   double T0 = 3e-3; //1e5 K
   double Tmin = 0.05*T0;

   double Q = 0.0;
   double Pp  = prim[PPP];
   double rho = prim[RHO];

   double T = Pp/rho;
   if( T > Tmin ){
      Q = -5.0*rho*rho*(T-Tmin)/T0*pow(T/T0,-1.7);
   }

   cons[TAU] += Q*dVdt;

}

void add_cooling( struct domain * theDomain , double dt ){

   struct cell ** theCells = theDomain->theCells;
   int Nt = theDomain->Nt;
   int Np = theDomain->Np;
   int * Nr = theDomain->Nr;
   double * t_jph = theDomain->t_jph;
   double * p_kph = theDomain->p_kph;
   double t = theDomain->t;

   int i,j,k;
   for( j=0 ; j<Nt ; ++j ){
      double thp = t_jph[j];
      double thm = t_jph[j-1];
      double th = .5*(thp+thm);
      for( k=0 ; k<Np ; ++k ){
         double php = p_kph[k];
         double phm = p_kph[k-1];
         int jk = j+Nt*k;
         for( i=0 ; i<Nr[jk] ; ++i ){
            struct cell * c = &(theCells[jk][i]);
            double rp = c->riph;
            double rm = rp - c->dr;
            double xp[3] = {rp,thp,php};
            double xm[3] = {rm,thm,phm};
            double dV = get_dV(xp,xm);
            cool_src( c->prim , c->cons , dV*dt );
         }
      }
   }

}