Skip to content

Commit

Permalink
Epa ad fit fix (#187)
Browse files Browse the repository at this point in the history
* add bounds to curve_fit Pe for advective dispersion

* fix import for aeration files and curve fit for advective disperison

* update version
  • Loading branch information
monroews authored Mar 8, 2019
1 parent af9f08c commit 53ac17d
Show file tree
Hide file tree
Showing 2 changed files with 8 additions and 6 deletions.
12 changes: 7 additions & 5 deletions aguaclara/research/environmental_processes_analysis.py
Original file line number Diff line number Diff line change
Expand Up @@ -221,8 +221,8 @@ def aeration_data(DO_column, dirpath):
filepaths = [os.path.join(dirpath, i) for i in filenames]
#DO_data is a list of numpy arrays. Thus each of the numpy data arrays can have different lengths to accommodate short and long experiments
# cycle through all of the files and extract the column of data with oxygen concentrations and the times
DO_data=[Column_of_data(i,0,-1,DO_column,'mg/L') for i in filepaths]
time_data=[(ftime(i,0,-1)).to(u.s) for i in filepaths]
DO_data=[column_of_data(i,0,DO_column,-1,'mg/L') for i in filepaths]
time_data=[(column_of_time(i,0,-1)).to(u.s) for i in filepaths]
aeration_collection = collections.namedtuple('aeration_results','filepaths airflows DO_data time_data')
aeration_results = aeration_collection(filepaths, airflows, DO_data, time_data)
return aeration_results
Expand Down Expand Up @@ -398,7 +398,7 @@ def E_Advective_Dispersion(t, Pe):
# replace any times at zero with a number VERY close to zero to avoid
# divide by zero errors
if isinstance(t, list):
t[t == 0] = 10**(-50)
t[t == 0] = 10**(-10)
return (Pe/(4*np.pi*t))**(0.5)*np.exp((-Pe*((1-t)**2))/(4*t))


Expand Down Expand Up @@ -549,13 +549,15 @@ def Solver_AD_Pe(t_data, C_data, theta_guess, C_bar_guess):
peclet number that best fits the data
"""

#remove time=0 data to eliminate divide by zero error
t_data = t_data[1:-1]
C_data = C_data[1:-1]
C_unitless = C_data.magnitude
C_units = str(C_bar_guess.units)
t_seconds = (t_data.to(u.s)).magnitude
# assume that a guess of 1 reactor in series is close enough to get a solution
p0 = [theta_guess.to(u.s).magnitude, C_bar_guess.magnitude,5]
popt, pcov = curve_fit(Tracer_AD_Pe, t_seconds, C_unitless, p0)
popt, pcov = curve_fit(Tracer_AD_Pe, t_seconds, C_unitless, p0, bounds=(0.01,np.inf))
Solver_theta = popt[0]*u.s
Solver_C_bar = popt[1]*u(C_units)
Solver_Pe = popt[2]
Expand Down
2 changes: 1 addition & 1 deletion setup.py
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
from setuptools import setup, find_packages

setup(name='aguaclara',
version='0.0.20',
version='0.0.21',
description='Open source functions for AguaClara water treatment research and plant design.',
url='https://github.com/AguaClara/aguaclara',
author='AguaClara at Cornell',
Expand Down

0 comments on commit 53ac17d

Please sign in to comment.