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

material properties in floc_model are missing their units #184

Closed monroews closed 5 years ago

monroews commented 5 years ago

If I call floc_model.Clay.Density it returns a value of 2650. It should return a value of 2650 kg/m3 The same is true for the other material properties.

If I apply units by fm.Clay.Density = 2650 * u.kg/u.m**3

Then when I call fm.sep_dist_clay(C_clay, fm.Clay) it returns a value of zero which isn't correct.

The material properties must have the correct units.

HannahSi commented 5 years ago

@monroews I apologize for the delay on resolving this issue! Elena and I are currently adding a fix to the "init" function for the different materials and making sure it is tested.

Just to clarify, was the output of your line of code supposed to be nonzero, or 0*u.kg/u.m**3?

monroews commented 5 years ago

Given a real clay concentration the separation distance between clay must be greater than zero. Might there be another error in the equation?

HannahSi commented 5 years ago

That makes sense. What is your C_Clay value? I ran the code below with an arbitrary concentration and obtained

fm.sep_dist_clay(.5*u.kg/u.m**3, fm.Clay)

Output: 9.836853510957468e-05 meter

Currently, the calculation for the separation distance is ( (density/concentration) (pi/6diameter^3) )^(1/3). Unfortunately I'm not familiar with floc modeling. Are you able to confirm whether this is correct?

monroews commented 5 years ago

I believe this is correct. So proceed!

On Sun, 10 Mar 2019, 4:15 p.m. Hannah Si, notifications@github.com wrote:

That makes sense. What is your C_Clay value? I ran the code below with an arbitrary concentration and obtained

fm.sep_dist_clay(.5*u.kg/u.m**3, fm.Clay)

Output: 9.836853510957468e-05 meter

Currently, the calculation for the separation distance is ( (density/concentration) (pi/6diameter^3) )^(1/3). Unfortunately I'm not familiar with floc modeling. Are you able to confirm whether this is correct?

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/AguaClara/aguaclara/issues/184#issuecomment-471340043, or mute the thread https://github.com/notifications/unsubscribe-auth/AbN_65-VHMUJ8rwg1aV9o-Mt7uCZnfmfks5vVWf4gaJpZM4badyg .

elenapertsalis commented 5 years ago

Thank you! This is just a note for us:

We're trying to apply a wrapper function to the init function for Material so that the correct units are automatically correctly applied to the inputs. This is our code:

class Material:
    @u.wraps(None, [None, u.m, u.kg/u.m**3, u.kg/u.mol], False)
    def __init__(self, name, diameter, density, molecWeight):
        self.name = name
        self.Diameter = diameter
        self.Density = density
        self.MolecWeight = molecWeight

When we try running the code:

   Clay = fm.Material('Clay', 7 * 10**-6 * u.m, 2650 * u.kg/u.m**3, None)

We get this error: pint.errors.DimensionalityError: Cannot convert from 'meter' ([length]) to 'kilogram / meter ** 3' ([mass] / [length] ** 3)

HannahSi commented 5 years ago

I think I understand the error now. I'll first rewrite the wrapper function as:

@u.wraps(ret=None, args=[None, u.m, u.kg/u.m**3, u.kg/u.mol], strict=False)

In the __init__ function, diameter is the 3rd argument, so u.wraps tries to convert diameter's units to the 3rd unit in the args list, which is u.kg/u.m3. If we modify args to be `[None, None, u.m, u.kg/u.m3, u.kg/u.mol]`, this error no longer appears.

It turns out though that the wrapper function doesn't actually apply the units to the arguments if they are unitless. It only gives an error if the units are wrong, so for each argument we'll have to check if it has units. If it doesn't, multiply it by the correct units. If it does, we actually still have to convert them to the correct units using the to function.

monroews commented 5 years ago

I don't think it is necessary to convert to the "correct" units. The whole point of using pint is that users can use any units they want and they equations still work. The only reason to use the "to" function is for display purposes if a user wants to see the value with a particular set of units.

HannahSi commented 5 years ago

I see, so to allow users to choose their units if they create their own floc_model.Material, we won't use u.wraps to enforce any units, but simply expect the user to input the correct units. Could you confirm that this all we need to change in floc_model.py? (The initializations for Clay, PACl, Alum, AlOH3, and Humic Acid all have units now.)

class Material:
    def __init__(self, name, diameter, density, molecWeight):
       ...

class Chemical(Material):
    def __init__(self, name, diameter, density, molecWeight, Precipitate,
                 AluminumMPM=None):
       ...
    def define_Precip(self, diameter, density, molecweight, alumMPM):
       ...

Clay = Material('Clay', 7 * 10**-6 * u.m, 2650 * u.kg/u.m**3, None)
PACl = Chemical('PACl', 9 * 10 **-8 * u.m, 1138 * u.kg/u.m**3,
                 1.039 * u.kg/u.mol, 'PACl', AluminumMPM=13)
Alum = Chemical('Alum', 7 * 10 **-8 * u.m, 2420 * u.kg/u.m**3, 
                0.59921 * u.kg/u.mol, 'AlOH3', AluminumMPM=2)
Alum.define_Precip(7 * 10 **-8 * u.m, 2420 * u.kg/u.m**3, 
                0.078 * u.kg/u.mol, 1)
HumicAcid = Chemical('Humic Acid', 72 * 10**-9 * u.m, 1780 * u.kg/u.m**3, None,
                'Humic Acid')
monroews commented 5 years ago

Could you test a couple of the functions that use material and make sure that they work. If they work, then go ahead and merge and publish!

elenapertsalis commented 5 years ago

Hi, for the functions in the floc_model file, do you have recommended test cases that you can give me? The following functions are examples of functions we need test cases for: Thank you! @u.wraps(u.kg/u.m**3, None, False) def dens_alum_nanocluster(coag): def p(C, Cprime): @u.wraps(u.m, [u.dimensionless, u.m, u.dimensionless], False) def diam_fractal(DIM_FRACTAL, DiamInitial, NumCol): @u.wraps(None, u.m, False) def ratio_clay_sphere(RatioHeightDiameter): @u.wraps(u.m, [u.W/u.kg, u.degK], False) def eta_kolmogorov(EnergyDis, Temp):

monroews commented 5 years ago

Create your own test cases. Just make sure that the functions return appropriate units and that the functions that call these functions still work after you've fixed the units.

elenapertsalis commented 5 years ago

For this release, we have not yet tested the functions because we need the textbook to do this. We've tested the constructors thus far.

HannahSi commented 5 years ago

On a different note, the UASB team said that they would like some guidance on how to create materials like Clay, PACl, etc. as we did above or how to change their properties. We should write some example code in team_resources (or the aguaclara wiki when it is working) for doing this.

elenapertsalis commented 5 years ago

I agree, Hannah! We should start working on that after break.

HannahSi commented 5 years ago

The code with the correct units will be released in version 24 (see the v0.0.24 branch). Apologies for the delay!