jueves, 23 de julio de 2020

Calculation of adiabatic saturation temperature

We will review how to calculate the adiabatic saturation temperature. The definition of this quantity was given in a previous post.

Problem definition


Given a gas-vapor mixture, if the inlet temperature and humidity are known, find the adiabatic saturation temperature.  The following equation is used:

From the inlet conditions we know $T_{\text{air}}$, Y and C. The main problem is that the humidity at saturation $Y_{\text{sat}}$ and the enthalpy of vaporization $\lambda$ depend on $T_{\text{as}}$. So an iterative calculation needs to be employed. Fortunately it can be easily implemented in Excel.

Example

This example comes from Seader's book. Air with inlet temperature 140 [F] and 12.5% relative humidity of water enters an process. Find the adiabatic saturation temperature of this current.
To solve this we also need enthalpy of vaporization and vapor pressure of water, and specific heat of both air and water, all as functions of temperature. All the required correlations can be found on Perry's book.

After a few trials by hand, the result is found using solver, and is equal to 87 [F]. The detailed calculations are on this file.




lunes, 13 de julio de 2020

The wet bulb temperature

I started to study drying of solids, so is a good opportunity to revisit this important concept. The following results were adapted mainly from Treybal's book.

The humidity

Absolute humidity (Y) is simply the ratio of water mass as vapor and the air mass, expressed here in terms of partial and total pressures. Based  on the equality of the molar and pressure ratio:

Where M are the respective molecular weights.

Measurement of wet bulb temperature

Imagine a wet cloth that is exposed to an air stream. If the air is not saturated, water will evaporate from the cloth. The energy required for evaporation comes initially from the rest of the water, decreasing its temperature. Now due to the difference in air and water temperatures, heat will flow from the air to the remaining water.
As more water evaporates the remaining water will keep cooling, increasing the heat transfer rate from the air, until an equilibrium is reached. At this point the heat flow from the air is just enough to sustain the water vaporization. The following relation is satisfied:


 

 Adiabatic cooling

Now, imagine another process, this time air flows over a liquid surface. Similar to the previous case, water evaporates, taking energy from the air in the form of heat.  If this process is carried without external heat flow, it is termed adiabatic.

Before, air temperature and humidity were assumed constant, as the cloth is small compared to the air flow. This time however, both the temperature and humidity of the air change appreciably. The air humidity increases until the air becomes saturated.

A heat balance for the process gives:


Where C is the heat capacity of the inlet humid air. This equation says that the heat lost by the air together with its starting humidity is equal to the heat needed to add the additional humidity. 

Relation between wet bulb temperature and adiabatic cooling


Now lets see the equations for both processes side by side:

It turns out that the following ratio, called the Lewis relation is approximately one for an air-water mixture. So the wet bulb temperature (Teq) is approximately equal to the adiabatic saturation temperature (Tas).

sábado, 6 de junio de 2020

Solving a flow distribution problem

Over the past few weeks I had to study fluid mechanics again, so I thought that for this week post I could share an example of that.

Problem description

Given a network of pipes, where the length and diameter of each section is known, calculate the flow of water in the network (the Bernoulli boundary values must also be known, details later). Lets solve for the network and data that follows (the construction material is PVC).
Section Length[m] Diameter[m]
1 2800 0.25
2 1300 0.18
3 800 0.18
4 1300 0.13

Theory

The Bernoulli equation says that the energy difference between two nodes is due to the friction losses by fluid flow:




The friction losses (h) are proportional to the length of the section (L) and the adimensional head loss (J). Several models are available for J, for water the Hazen-Williams correlation is:

Where the coefficient C is specific for the material. There are four unknown flow rates (Q), we can write three Bernoulli differences between pairs of boundary nodes, and the last equation is the flow rate balance at the central node.


This kind of problems are solved easily using Excel's solver, giving as result Q1=92[l/s], Q2=15[l/s], Q3=55[l/s] and Q4=23[l/s]. The Excel sheet is available on this link.

martes, 17 de diciembre de 2019

Bubble point calculation

The bubble point for a mixture of compounds is a useful parameter for designing distillation systems.

At a given pressure P, the bubble point is the temperature at which the first bubble of vapor forms when a liquid mixture is heated. If a gaseous mixture satisfies Dalton's law the total pressure is equal to the sum of the partial pressures of the components of the mixture:

Where \( P_i \) is the partial pressure of component i. Additionally, if Raoult's law is satisfied, the partial pressure can be expressed in terms of the liquid composition and the vapor pressure of the pure compound:

As last ingredient we need a model for the vapor pressure as a function of temperature, in this example we will use an exponential:



Example

A mixture of 50% toluene and 50% benzene is kept at a pressure of 760 [mmHg]. Determine the bubble point of the mixture using the following parameter values for the vapor pressure (Where \( T \) is in Celsius degrees and the vapor pressure value has units of [mmHg]).

Compound A B
Benzene -324.01 10.46
Toluene -364.43 9.95

One key point of this calculation is to note that at bubble point, the first single bubble is a negligible fraction of the total mass of the system, so the total composition \( z_i \) is equal to the liquid composition \( x_i \). We must find then the value of temperature that satisfies:



For this we use a short python code:

from scipy import optimize
import math

def vapor_press(T, A, B):
    return math.exp(A/T + B)

def system_press(T):
    return 760 - (0.5 * vapor_press(T, -324.01, 10.46) + 0.5 * vapor_press(T, -364.43, 9.95))

T_b = optimize.bisect(system_press, 85, 110)
print('Bubble point temperature: {:.2f} C'.format(T_b))
 
 
The result is 93.57 [C]. Due to the problem being nonlinear, we used the bisection algorithm to find a solution, it looks for a root in a specified interval. The bubble point must be located between the boiling point of both pure compounds, which can be obtained from:


Solving for each compound separately you obtain the approximate limits 85 [C] for benzene and 110 [C] for toluene.

sábado, 8 de julio de 2017

A recursive finite difference implementation


In this entry I will show how to program recursive functions on Python, using as an example the determination of finite difference expressions.

 Theory

The book Computational fluid dynamics, Fourth edition, by Klaus A. Hoffmann and Steve T. Chiang provides recursive definitions for finite difference approximations of derivatives. First they define the forward difference and backward difference operators:



Higher order operators are defined recursively:


Finally they define higher order derivatives around a node i in terms of these higher order operators.

Forward difference


Backward difference


Central difference (n even)



Central difference (n odd)



Implementation

 A recursive definition is a perfect example to be programmed. A Python module finiteDifference.py was implemented (https://drive.google.com/drive/folders/0B7_nEWx0hv56UE9SLWVQUmQyTjQ?usp=sharing). The module contains the basic class definitions needed to obtain higher order finite difference approximations. The basic class is GridPoint, which allows to define a function evaluation at some specified grid location. Here it is __init__ method:

class GridPoint(object):
       def __init__(self, coeff = 1, ind = 0):
....

This allows you to define a gridPoint instance and inspect it using print:

>>> import finiteDifference as fd
>>> a=fd.GridPoint(coeff=2,ind=1)
>>> print(a)
  2.0*f(i+1)

The next class is Den, which is the denominator of a finite difference approximation, it is a very simple class:

class Den(object):
 def __init__(self, power = 1):
...
 

>>> b=fd.Den(power=2)
>>> print(b)
(dh)^2

Next is the Num class,  which is the numerator of a finite difference approximation. It is constructed from a list of instances of the GridPoint class, as the other classes, it can be inspected by using print on it:


class Num(object):
 def __init__(self, exp = None):
...

>>> c=fd.GridPoint(coeff=2,ind=2)
>>> d=fd.GridPoint(coeff=1.0,ind=-1)
>>> e=fd.Num([c,d])
>>> print(e)
  2.0*f(i+2)+  1.0*f(i-1)

The Num by default groups GridPoint instances terms with equal ind parameter:

>>> f=fd.GridPoint(coeff=-1,ind=3)
>>> print(f)
 -1.0*f(i+3)
>>> g=fd.GridPoint(coeff=2,ind=3)
>>> print(g)
  2.0*f(i+3)
>>> h=fd.Num([f,g])
>>> print(h)
  1.0*f(i+3)

Finally, the class FinDiff is the class used to obtain the finite difference approximation. FindDiff makes internal use of the rest of the classes defined in the module (parameter method can be one of backward, forward and central):


class FinDiff(object):
 def __init__(self, method = 'central', order = 1, point = GridPoint()):
...

FindDiff class has as attributes numer and denom that are instances of classes Num and Den, respectively, which must be printed separately:

>>> i=fd.FinDiff(method='forward',order=1,point=fd.GridPoint())
>>> print(i.numer)
  1.0*f(i+1) -1.0*f(i)
>>> print(i.denom)
(dh)^1

Results

Hoffmann and Chiang offer the following tables with the coefficients of the different grid point evaluations involved  in the finite difference representations of the derivatives, for every method, with the order ranging from 1 to 4.

Forward difference

Order$f_i$$f_{i+1}$$f_{i+2}$$f_{i+3}$$f_{i+4}$
1-11
21-21
3-13-31
41-46-41

Lets see the results using FinDiff:

>>> fd1=fd.FinDiff(method='forward',order=1,point=fd.GridPoint())
>>> fd2=fd.FinDiff(method='forward',order=2,point=fd.GridPoint())
>>> fd3=fd.FinDiff(method='forward',order=3,point=fd.GridPoint())
>>> fd4=fd.FinDiff(method='forward',order=4,point=fd.GridPoint())
>>> fd5=fd.FinDiff(method='forward',order=5,point=fd.GridPoint())
>>> print(fd1.numer)
  1.0*f(i+1) -1.0*f(i)
>>> print(fd2.numer)
  1.0*f(i+2) -2.0*f(i+1)+  1.0*f(i)
>>> print(fd3.numer)
  1.0*f(i+3) -3.0*f(i+2)+  3.0*f(i+1) -1.0*f(i)
>>> print(fd4.numer)
  1.0*f(i+4) -4.0*f(i+3)+  6.0*f(i+2) -4.0*f(i+1)+  1.0*f(i)
>>> print(fd5.numer)
  1.0*f(i+5) -5.0*f(i+4)+ 10.0*f(i+3)-10.0*f(i+2)+  5.0*f(i+1) -1.0*f(i)

The results agree up to fourth order and additionally the result for order 5 is shown.

 Backwards difference

Order$f_{i-4}$$f_{i-3}$$f_{i-2}$$f_{i-1}$$f_{i}$
1-11
21-21
3-13-31
41-46-41

>>> bd1=fd.FinDiff(method='backward',order=1,point=fd.GridPoint())
>>> bd2=fd.FinDiff(method='backward',order=2,point=fd.GridPoint())
>>> bd3=fd.FinDiff(method='backward',order=3,point=fd.GridPoint())
>>> bd4=fd.FinDiff(method='backward',order=4,point=fd.GridPoint())
>>> bd5=fd.FinDiff(method='backward',order=5,point=fd.GridPoint())
>>> print(bd1.numer)
  1.0*f(i) -1.0*f(i-1)
>>> print(bd2.numer)
  1.0*f(i) -2.0*f(i-1)+  1.0*f(i-2)
>>> print(bd3.numer)
  1.0*f(i) -3.0*f(i-1)+  3.0*f(i-2) -1.0*f(i-3)
>>> print(bd4.numer)
  1.0*f(i) -4.0*f(i-1)+  6.0*f(i-2) -4.0*f(i-3)+  1.0*f(i-4)
>>> print(bd5.numer)
  1.0*f(i) -5.0*f(i-1)+ 10.0*f(i-2)-10.0*f(i-3)+  5.0*f(i-4) -1.0*f(i-5)

Central difference

Order$f_{i-2}$$f_{i-1}$$f_{i}$$f_{i+1}$$f_{i+2}$
1-0.500.5
21-21
3-0.510-10.5
41-46-41

>>> cd1=fd.FinDiff(method='central',order=1,point=fd.GridPoint())
>>> cd2=fd.FinDiff(method='central',order=2,point=fd.GridPoint())
>>> cd3=fd.FinDiff(method='central',order=3,point=fd.GridPoint())
>>> cd4=fd.FinDiff(method='central',order=4,point=fd.GridPoint())
>>> cd5=fd.FinDiff(method='central',order=5,point=fd.GridPoint())
>>> print(cd1.numer)
  0.5*f(i+1)+  0.0*f(i) -0.5*f(i-1)
>>> print(cd2.numer)
  1.0*f(i+1) -2.0*f(i)+  1.0*f(i-1)
>>> print(cd3.numer)
  0.5*f(i+2) -1.0*f(i+1)+  0.0*f(i)+  1.0*f(i-1) -0.5*f(i-2)
>>> print(cd4.numer)
  1.0*f(i+2) -4.0*f(i+1)+  6.0*f(i) -4.0*f(i-1)+  1.0*f(i-2)
>>> print(cd5.numer)
  0.5*f(i+3) -2.0*f(i+2)+  2.5*f(i+1)+  0.0*f(i) -2.5*f(i-1)+  2.0*f(i-2) -0.5*f(i-3)
 

Starting with IDAES: Steady state CSTR