AguaClara / aguaclara

An open-source Python package for designing and performing research on AguaClara water treatment plants.
https://aguaclara.github.io/aguaclara/
MIT License
24 stars 13 forks source link

Wrong units for ConcClay in floc_model.num_clay() #199

Closed monroews closed 5 years ago

monroews commented 5 years ago

Units on num_clay function are incorrect. This points out the dangers of coercing units! The units on concentration are mass per volume, not 1/volume as shown below!

@u.wraps(1/u.m3, [u.kg/u.m3, u.m], False) def num_clay(ConcClay, material): """Return the number of clay particles in suspension.""" return ConcClay / ((material.Density np.pi material.Diameter**3) / 6)

HannahSi commented 5 years ago

Thank you for letting us know about this! We'll make a fix for this right away.

HannahSi commented 5 years ago

@monroews The unit of 1/volume is actually for the output of the function. The required unit for ConcClay is the first term in the brackets, u.kg/u.m**3. Did you get a specific error when calling the function?

HannahSi commented 5 years ago

If you did get an error, it may be that the function is wrongly coercing the material parameter's units to u.m. I think its units should be ignored, since material is an object, not a quantity.

monroews commented 5 years ago

The function returns this value for a clay concentration of 5 NTU. 1.544e+10 / kilogram 1/meter3

It should return 15443606337.758038 1/meter3

The correct function should be

def num_clay(ConcClay, material):
    """Return the number of clay particles in suspension."""
    return (ConcClay / ((material.Density * np.pi * material.Diameter**3) / 6)).to(1/u.m**3)

Note that this function should be generalized to be generic for any material and should be Particle_number_concentration. Then it can work for any material that has a particle size. Okay... to make it completely generic it would need to have a fractal dimension as an input too. (I think we have a floc number concentration function that uses the fractal dimension).

HannahSi commented 5 years ago

Hi Monroe, I'm sorry we've been caught up with preparing the Python analysis MNC and haven't gotten to work on this issue in the past week! We found that the problem was actually with the u.wraps function, which was applying the units of 1/m^3 on top of the units of 1/kg:

>>> fm.num_clay(5*u.NTU, fm.Clay)
<Quantity(1.544e+10 / kilogram, '1 / meter ** 3')>
>>> x = fm.num_clay(5*u.NTU, fm.Clay)
>>> x.units
<Unit('1 / meter ** 3')>

To solve this, I think we will remove the u.wraps wrapper function and rewrite the function as:

def num_clay(ConcClay, material):
    return ConcClay.to(material.Density.units) / ((material.Density * np.pi * material.Diameter**3) / 6)
HannahSi commented 5 years ago

As for generalizing the function to generic materials, would this be along the lines of what you were suggesting?

def particle_number_concentration(ConcMat, material):
    """Return the number of particles in suspension.

    :param ConcMat: Concentration of the material
    :type ConcMat: float
    :param material: The material in solution
    :type material: floc_model.Material
    """
    ...
monroews commented 5 years ago

yes!

HannahSi commented 5 years ago

Great, thank you!

HannahSi commented 5 years ago

The num_clay() function will be changed to particle_number_concentration() and have the right units in version v0.0.25.