van_der_waals.cpp

气体热力性质计算程序

#include "..\\include\\van_der_Waals.h"

namespace XZGas
{
	//the derivative of v above T with p as a constant
	inline double vdWGasModel::dvdTp(double pp,double vv)
    {
		double v3=vv*vv*vv;
        
		return Rg*v3/(pp*v3-a*vv+a*b*2);
    }

    vdWGasModel::vdWGasModel(const GasData &gas):
    GasModel(gas)
    {
        a=Rg*Rg*gas.Tc*gas.Tc/gas.pc*27/64;
        b=Rg*gas.Tc/8/gas.pc;
    }
    
    vdWGasModel::~vdWGasModel()
    {
    }
    
    inline double vdWGasModel::Cp(double TT,double pp)
    {
		//from the relationship between Cp and Cv: Cp-Cv=T*dvdTp*dpdTv
		double vv=v(TT,pp);
		double cv=Cv(TT,vv);
		double cp_cv=TT*dvdTp(pp,vv)*Rg/(vv-b);

	    return cv+cp_cv;
    }

	inline double vdWGasModel::Cv(double TT,double vv)
	{		
		//Cv=Cv(T)+T*SdpdTv2*dv
		double cv=_Cp(TT);
		cv-=Rg;
		
		return cv;
	}

	inline double vdWGasModel::Gama(double TT,double pp)
	{
		//Gama=Cp/Cv

		double vv=v(TT,pp);
		double cv=Cv(TT,pp);
		return (cv+TT*dvdTp(pp,vv)*Rg/(vv-b))/cv;
	}

	inline double vdWGasModel::h(double TT,double pp)
	{	
		//from definition: h=u+pv	
		
		double vv=v(TT,pp),v0=v(T0,p0);
		double uu=u(TT,vv);
		
		return uu+pp*vv-p0*v0+h0;
	}

	inline double vdWGasModel::u(double TT,double vv)
	{	
		//from the formula: du=CvdT+(T*dpdTv-p)dv

		//the first part of calculus above temperature
		double u1=_ICp(TT)-_ICp(T0)-Rg*(TT-T0);

		//the second part of calculus above special volume
		double v0=v(T0,p0);
		double u2=a*(1/v0-1/vv);
		
		return u1+u2+u0;
	}

	inline double vdWGasModel::s(double TT,double pp)
	{		
		double vv=v(TT,pp);
        		
		double s1=_ICpT(TT)-_ICpT(T0)-Rg*log(TT/T0);
		
		double v0=v(T0,p0);
		double s2=Rg*log((vv-b)/(v0-b));

		return s1+s2+s0;
	}

    //calculate the pressure with 
    inline double vdWGasModel::p(double TT,double vv)
    {
        return (Rg*TT/(vv-b)-a/vv/vv);
    }    
    
	//calculate the value of v with iteration
    inline double vdWGasModel::v(double TT,double pp)
    {
        double _v1=1.0,_v2=1.0,d;   
        short n=0;
        //this iteration converges as fast as usually for less than 5 times
        do
        {
            _v1=b+Rg*TT/(pp+a/_v1/_v1);
            d=_v1-_v2;
            ++n;
            /*if this iteration runs more than 7 times, it may run out of 
            the precision of double varibles*/
            if(n>7)
                break;
			_v2=_v1;
        }while(d>1e-10||d<-1e-10);
        return _v1;
    }
    
	//to compute the temperature under certain pressure and volume
    inline double vdWGasModel::T(double pp,double vv)
    {
        return (pp+a/vv/vv)*(vv-b)/Rg;    
    }

	inline double vdWGasModel::Alphav(double pp,double vv)
	{
		return dvdTp(pp,vv)/vv;
	}

	inline double vdWGasModel::KT(double TT,double vv)
	{
		//variable Reciprocal is the reciprocal of the return value,
		//which has a much simpler form
		double Reciprocal=vv*Rg*TT/(vv-b)/(vv-b)-2*a/vv/vv;

		return 1.0/Reciprocal;
	}

	inline double vdWGasModel::Beta(double pp,double vv)
	{
		return Rg/(vv-b);
	}

	inline void vdWGasModel::ChangeGas(const GasData &newgas)
	{
		//the same kind of gas
		if(strcmp(GasName,newgas.Name) == 0)
			return ;

		GasModel::ChangeGas(newgas);

		a=Rg*Rg*newgas.Tc*newgas.Tc/newgas.pc*27/64;
        b=Rg*newgas.Tc/8/newgas.pc;
	}
		
}