Hi Roy

Now it is time to deal with the term:  Nabla . (rho  sum_s (h_s Ds Nabla(c_s) ) ).

I will have to calculate the derivative of h_N and h_N2, so the strategy that I've been using that was just writing h as a function of T cannot be applied now.

In my opinion, we should use the expression for the h_equilibrium (around line 743 of thermodynamic.C):

      HoverRT = -a[0]/T2 + a[1]*lnT/T + a[2] + a[3]*T/2. + a[4]*T2/3. + a[5]*T3/4. + a[6]*T4/5. + a[8]/T;

and multiply it by R*T.

So, I'd have:

h_N = R*T*(a0_N/T^2  + a1_N*lnT/T + a2_N + a3_N*T/2 + a4_N*T^2/3 + a5_N*T^3/4 + a6_N*T^4/5 + a8_N/T)
h_N2 = R*T*(a0_N2/T^2  + a1_N2*lnT/T + a2_N2 + a3_N2*T/2 + a4_N2*T^2/3 + a5_N2*T^3/4 + a6_N2*T^4/5 + a8_N2/T)

Do you agree?

--
Kemelli C. Estacio-Hiroms
Postdoctoral Researcher
Institute for Computational Engineering and Sciences
The University of Texas at Austin
(512) 299-8403
kemelli@ices.utexas.edu