PLUTO
4.4-patch2
|
Cooling main header file. More...
Go to the source code of this file.
Functions | |
double | CompEquil (double, double, double *) |
void | CoolingSource (const Data *, double, timeStep *, Grid *) |
double | GetMaxRate (double *, double *, double) |
void | Jacobian (double *, double *, double **) |
void | Numerical_Jacobian (double *, double **) |
void | PowerLawCooling (Data_Arr, double, timeStep *, Grid *) |
void | Radiat (double *, double *) |
double | SolveODE_CK45 (double *, double *, double *, double, double, intList *) |
double | SolveODE_RKF23 (double *, double *, double *, double, intList *) |
double | SolveODE_RKF12 (double *, double *, double *, double, intList *) |
double | SolveODE_RK4 (double *, double *, double *, double, intList *) |
void | H2RateTables (double, double *) |
Contains function prototypes and basic definitions for the cooling modules.
double CompEquil | ( | double | N, |
double | T, | ||
double * | v0 | ||
) |
[in] | n | the particle number density (not needed, but kept for compatibility) |
[in] | T | the temperature (in K) for which equilibrium must be found. |
[in,out] | v | an array of primitive variables. On input, only density needs to be defined. On output, fractions will be updated to the equilibrium values. |
Compute the equilibrium ionization balance for (rho,T)
double GetMaxRate | ( | double * | v0, |
double * | k1, | ||
double | T0 | ||
) |
Return an estimate of the maximum rate (dimension 1/time) in the chemical network. This will serve as a "stiffness" detector in the main ode integrator.
For integration to be carried explicitly all the time, return a small value (1.e-12).
void H2RateTables | ( | double | T, |
double * | krvals | ||
) |
On first call, compute and store rate coefficient arrays needed for the H2_COOL module. On subsequent calls, use linear interpolation between adjacent tabulated values to obtain the coefficients at the desired temperature T.
void Jacobian | ( | double * | v, |
double * | rhs, | ||
double ** | dfdy | ||
) |
Compute the jacobian J(k,l) = dfdy
k = row index l = column index
J(0,0) J(0,1) ... J(0, n-1) J(1,0) J(1,1) .... J(1, n-1) . . . . . . . . . J(n-1,0) .... J(n-1, n-1)
or,
+--------------------—+
void Numerical_Jacobian | ( | double * | v, |
double ** | J | ||
) |
Compute Jacobian matrix numerically and get eigenvector and eigenvalue decomposition
The purpose of this function is to detect stiffness.
Note: This function is EXTREMELY time consuming.
void Radiat | ( | double * | v, |
double * | rhs | ||
) |
Cooling for optically thin plasma up to about 200,000 K Plasma composition: H, HeI-II, CI-V, NI-V, OI-V, NeI-V, SI-V Assumed abundances in elem_ab Uses S : Array = Variables vector x line points rhs : output for the system of ODE ibeg, iend : begin and end points of the current line
Cooling for neutral or singly ionized gas: good up to about 35,000 K in equilibrium or shocks in neutral gas up to about 80 km/s. Assumed abundances in ab Uses t : Kelvin dene : electron density cm*-3 fneut : hydrogen neutral fraction (adimensionale) ci,cr : H ionization and recombination rate coefficients
em(1) = TOTAL EMISSIVITY : (ergs cm**3 s**-1) em(2) = Ly alpha + two photon continuum: Aggarwal MNRAS 202, 10**4.3 K em(3) = H alpha: Aggarwal, Case B em(4) = He I 584 + two photon + 623 (all n=2 excitations): Berrington &Kingston,JPB 20:
em(5) = C I 9850 + 9823: Mendoza, IAU 103, 5000 K em(6) = C II, 156 micron: Mendoza, 10,000 K em(7) = C II] 2325 A: Mendoza, 15,000 K em(8) = N I 5200 A: Mendoza, 7500 K em(9) = N II 6584 + 6548 A: Mendoza em(10) = O I 63 micron: Mendoza,2500 K em(11) = O I 6300 A + 6363 A: Mendoza, 7500 K em(12) = O II 3727: Mendoza em(13) = Mg II 2800: Mendoza em(14) = Si II 35 micron: Dufton&Kingston, MNRAS 248 em(15) = S II 6717+6727: Mendoza em(16) = Fe II 25 micron: Nussbaumer&Storey em(17) = Fe II 1.6 micron em(18) = thermal energy lost by ionization em(19) = thermal energy lost by recombination (2/3 kT per recombination. The ionization energy lost is not included here.
Provide r.h.s. for tabulated cooling.
double SolveODE_CK45 | ( | double * | v0, |
double * | k1, | ||
double * | v5th, | ||
double | dt, | ||
double | tol, | ||
intList * | vars | ||
) |
Explicit 5th order adaptive step size Cash-Karp integrator. Attempt to integrate the cooling network with a single step if the error is within the given tolerance. Otherwise, reduce the time step and repeat the integration as many times as necessary.
[in] | v_in | initial condition array |
[in] | k_in | right hand side, R(v_in) |
[out] | v5th | final solution array |
[in] | dt | the time step we want to achieve |
[in] | tol | the relative tolerance |
[in] | *vars | the list of dependent variables |
double SolveODE_RK4 | ( | double * | v0, |
double * | k1, | ||
double * | v4th, | ||
double | dt, | ||
intList * | var_list | ||
) |
Solve the system of ODE with a standard RK4 integrator (no adaptive step).
double SolveODE_RKF12 | ( | double * | v0, |
double * | k1, | ||
double * | v2nd, | ||
double | dt, | ||
intList * | vars | ||
) |
Explicit 2nd order integrator with 1st-order embedded solution. Attempt to integrate the cooling network with a single step if the error is within the given tolerance. The error is computed from the difference between 1st and 2nd order solutions.
[in] | v0 | initial condition array |
[in] | k1 | right hand side, R(v0) |
[out] | v2nd | final solution array |
[in] | dt | the time step we want to achieve |
[in] | *vars | the list of dependent variables |
double SolveODE_RKF23 | ( | double * | v0, |
double * | k1, | ||
double * | v3rd, | ||
double | dt, | ||
intList * | vars | ||
) |
Explicit 3rd order integrator with 2nd-order embedded solution. Attempt to integrate the cooling network with a single step if the error is within the given tolerance.
The error is computed from the difference between 1st and 2nd order solutions.
[in] | v0 | initial condition array |
[in] | k1 | right hand side, R(v0) |
[out] | v2nd | final solution array |
[in] | dt | the time step we want to achieve |
[in] | *vars | the list of dependent variables |