diff --git a/epyt/epanet.py b/epyt/epanet.py index aa1f688..4aac8b3 100644 --- a/epyt/epanet.py +++ b/epyt/epanet.py @@ -311,6 +311,25 @@ def __init__(self): EN_R_IS_CLOSED = 2 EN_R_IS_ACTIVE = 3 + # MSX Constants + MSX_NODE = 0 + MSX_LINK = 1 + MSX_TANK = 2 + MSX_SPECIES = 3 + MSX_TERM = 4 + MSX_PARAMETER = 5 + MSX_CONSTANT = 6 + MSX_PATTERN = 7 + + MSX_BULK = 0 + MSX_WALL = 1 + + MSX_NOSOURCE = -1 + MSX_CONCEN = 0 + MSX_MASS = 1 + MSX_SETPOINT = 2 + MSX_FLOWPACED = 3 + class EpytValues: @@ -463,7 +482,7 @@ def isList(var): class epanet: """ EPyt main functions class """ - def __init__(self, *argv, version=2.2, loadfile=False): + def __init__(self, *argv, version=2.2, loadfile=False, msx=False): # Constants # Demand model types. DDA #0 Demand driven analysis, # PDA #1 Pressure driven analysis. @@ -523,6 +542,7 @@ def __init__(self, *argv, version=2.2, loadfile=False): # Initial attributes self.classversion = __version__ + self.api = epanetapi(version, msx=msx) print(f'EPANET version {self.getVersion()} ' f'loaded (EPyT version {self.classversion}).') @@ -592,6 +612,10 @@ def __init__(self, *argv, version=2.2, loadfile=False): plt.rcParams['figure.constrained_layout.use'] = True plt.rcParams['figure.max_open_warning'] = 30 + if msx: + self.msx = epanetmsxapi() + MSX_SPECIES = 3 + def addControls(self, control, *argv): """ Adds a new simple control. @@ -4810,7 +4834,6 @@ def getNodeSourceTypeIndex(self, *argv): else: return value - def getNodeTankBulkReactionCoeff(self, *argv): """ Retrieves the tank bulk rate coefficient. @@ -11115,6 +11138,806 @@ def __setNodeDemandPattern(self, fun, propertyCode, value, *argv): eval('self.api.' + fun + '(i, categ, param[j])') j += 1 + """MSX Funtions""" + + def loadMSXfile(self, msxname): + self.realmsx = msxname + self.MSXPythonSetup(msxname) + + def unloadMSX(self): + self.msx.MSXclose() + + def getMSXSpeciesCount(self): + MSX_SPECIES = self.ToolkitConstants.MSX_SPECIES + return self.msx.MSXgetcount(MSX_SPECIES) + + def getMSXConstantsCount(self): + MSX_CONSTANT = self.ToolkitConstants.MSX_CONSTANT + return self.msx.MSXgetcount(MSX_CONSTANT) + + def getMSXParametersCount(self): + MSX_PARAMETER = self.ToolkitConstants.MSX_PARAMETER + return self.msx.MSXgetcount(MSX_PARAMETER) + + def getMSXPatternsCount(self): + MSX_PATTERN = self.ToolkitConstants.MSX_PATTERN + return self.msx.MSXgetcount(MSX_PATTERN) + + def saveMSXFile(self, msxname): + self.msx.MSXsavemsxfile(msxname) + + def saveMSXQualityFile(self, outfname): + self.msx.MSXsaveoutfile(outfname) + + def solveMSXCompleteHydraulics(self): + self.msx.MSXsolveH() + + def solveMSXCompleteQuality(self): + self.msx.MSXsolveQ() + + def writeMSXReport(self): + self.msx.MSXreport() + + def useMSXHydraulicFile(self, hydname): + self.msx.MSXusehydfile(hydname) + + def getMSXPatternValue(self, patternIndex, patternStep): + return self.msx.MSXgetpatternvalue(patternIndex, patternStep) + + def initializeMSXQualityAnalysis(self, flag): + self.msx.MSXinit(flag) + + def stepMSXQualityAnalysisTimeLeft(self): + t, tleft = self.msx.MSXstep() + return t, tleft + + def getMSXError(self, code): + self.msx.MSXgeterror(code) + + def getMSXOptions(self, param="", getall=False): + msxname=self.msxname + options = {} + options["AREA_UNITS"] = "FT2" + options["RATE_UNITS"] = "HR" + options["SOLVER"] = "EUL" + options["TIMESTEP"] = 300 + options["ATOL"] = 0.01 + options["RTOL"] = 0.001 + options["COUPLING"] = "NONE" + options["COMPILER"] = "NONE" + + try: + with open(msxname, 'r') as f: + sect = 0 + for line in f: + if line.startswith("[OPTIONS]"): + sect = 1 + continue + elif line.startswith("[END") or line.startswith("[REPORTS]"): + break + elif sect == 1: + if not line.strip(): + continue + if line.startswith("["): + return options + key, value = line.split(None, 1) + if key == param or not param: + if key == "TIMESTEP": + options["TIMESTEP"] = float(value) + elif key == "AREA_UNITS": + options["AREA_UNITS"] = value.strip() + elif key == "RATE_UNITS": + options["RATE_UNITS"] = value.strip() + elif key == "SOLVER": + options["SOLVER"] = value.strip() + elif key == "RTOL": + options["RTOL"] = float(value) + elif key == "ATOL": + options["ATOL"] = float(value) + elif key == "COUPLING": + options["COUPLING"] = value.strip() + elif key == "COMPILER": + options["COMPILER"] = value.strip() + else: + options[key] = value + + if not getall and param: + return options[param] + + except FileNotFoundError: + warnings.warn("Please load MSX File.") + return {} + return options + + def getMSXTimeStep(self): + return self.getMSXOptions("TIMESTEP") + + def getMSXRateUnits(self): + return self.getMSXOptions("RATE_UNITS") + + def getMSXAreaUnits(self): + return self.getMSXOptions("AREA_UNITS") + + def getMSXCompiler(self): + return self.getMSXOptions("COMPILER") + + def getMSXCoupling(self): + return self.getMSXOptions("COUPLING") + + def getMSXSolver(self): + return self.getMSXOptions("SOLVER") + + def getMSXAtol(self): + return self.getMSXOptions("ATOL") + + def getMSXRtol(self): + return self.getMSXOptions("RTOL") + + def getMSXConstantsNameID(self, varagin=None): + x = self.getMSXConstantsCount() + value = {} + if varagin is None: + varagin = {} + for i in range(1, x + 1): + varagin[i] = i + 1 + MSX_CONSTANT = self.ToolkitConstants.MSX_CONSTANT + if x > 0: + for i in varagin: + len = self.msx.MSXgetIDlen(MSX_CONSTANT, i) + value[i] = self.msx.MSXgetID(MSX_CONSTANT, i, len) + return value + + def getMSXParametersNameID(self, varagin=None): + x = self.getMSXParametersCount() + value = {} + if varagin is None: + varagin = {} + for i in range(1, x + 1): + varagin[i] = i + 1 + MSX_PARAMETER = self.ToolkitConstants.MSX_PARAMETER + if x > 0: + for i in varagin: + len = self.msx.MSXgetIDlen(MSX_PARAMETER, i) + value[i] = self.msx.MSXgetID(MSX_PARAMETER, i, len) + return value + + def getMSXPatternsNameID(self, varagin=None): + x = self.getMSXPatternsCount() + value = {} + if varagin is None: + varagin = {} + for i in range(1, x + 1): + varagin[i] = i + 1 + MSX_PATTERN = self.ToolkitConstants.MSX_PATTERN + if x > 0: + for i in varagin: + len = self.msx.MSXgetIDlen(MSX_PATTERN, i) + value[i] = self.msx.MSXgetID(MSX_PATTERN, i, len) + return value + + def getMSXSpeciesNameID(self, varagin=None): + x = self.getMSXSpeciesCount() + value = {} + if varagin is None: + varagin = {} + for i in range(1, x + 1): + varagin[i] = i + 1 + MSX_SPECIES = self.ToolkitConstants.MSX_SPECIES + if x > 0: + for i in varagin: + len = self.msx.MSXgetIDlen(MSX_SPECIES, i) + value[i] = self.msx.MSXgetID(MSX_SPECIES, i, len) + return value + + def getMSXParametersIndex(self, varagin=None): + x = self.getMSXParametersCount() + value = {} + if varagin is None: + varagin = {} + varagin = self.getMSXParametersNameID() + y = [value for value in varagin.values()] + else: + y = varagin + MSX_PARAMETER = self.ToolkitConstants.MSX_PARAMETER + if x > 0: + for i in y: + value[i] = self.msx.MSXgetindex(MSX_PARAMETER, i) + return value + + def getMSXSpeciesIndex(self, varagin=None): + x = self.getMSXSpeciesCount() + value = {} + if varagin is None: + varagin = {} + varagin = self.getMSXSpeciesNameID() + y = [value for value in varagin.values()] + else: + y = varagin + MSX_SPECIES = self.ToolkitConstants.MSX_SPECIES + if x > 0: + for i in y: + value[i] = self.msx.MSXgetindex(MSX_SPECIES, i) + return value + + def getMSXPatternsIndex(self, varagin=None): + x = self.getMSXSpeciesCount() + value = {} + if varagin is None: + varagin = {} + varagin = self.getMSXPatternsNameID() + y = [value for value in varagin.values()] + else: + y = varagin + MSX_PATTERN = self.ToolkitConstants.MSX_PATTERN + if x > 0: + for i in y: + value[i] = self.msx.MSXgetindex(MSX_PATTERN, i) + return value + + def getMSXConstantsIndex(self, varagin=None): + x = self.getMSXConstantsCount() + value = {} + if varagin is None: + varagin = {} + varagin = self.getMSXConstantsNameID() + y = [value for value in varagin.values()] + else: + y = varagin + MSX_CONSTANT = self.ToolkitConstants.MSX_CONSTANT + if x > 0: + for i in y: + value[i] = self.msx.MSXgetindex(MSX_CONSTANT, i) + return value + + def getMSXConstantsValue(self, varagin=None): + x = self.getMSXConstantsCount() + value = {} + if varagin is None: + varagin = {} + for i in range(1, x + 1): + varagin[i] = i + 1 + if x > 0: + for i in varagin: + value[i] = self.msx.MSXgetconstant(i) + return value + + def getMSXParametersPipesValue(self): + x = self.getLinkPipeCount() + y = self.getMSXParametersCount() + value = [] + for i in range(1, x + 1): + value_row = [] + for j in range(1, y + 1): + param = self.msx.MSXgetparameter(1, i, j) + value_row.append(param) + value.append(value_row) + return value + + def getMSXParametersTanksValue(self): + x = self.getNodeTankIndex() + y = self.getMSXParametersCount() + value = {} + for i in x: + value[i] = [] + for j in range(1, y + 1): + param = self.msx.MSXgetparameter(0, i, j) + value[i].append(param) + return value + + def getMSXPatternsLengths(self, varagin=None): + x = self.getMSXPatternsCount() + value = {} + if varagin is None: + for i in range(1, x + 1): + value[i] = self.msx.MSXgetpatternlen(i) + else: + if x > 0: + for i in varagin: + value[i] = self.msx.MSXgetpatternlen(i) + return value + + def getMSXPattern(self): + z = self.getMSXPatternsCount() + if z == 0 : + return + val = self.getMSXPatternsLengths() + y = [value for value in val.values()] + tmpmaxlen = max(y) + value = [[0]*tmpmaxlen for _ in range(self.getMSXPatternsCount())] + for i in range(1, self.getMSXPatternsCount() + 1): + z = self.getMSXPatternsLengths([i]) + tmplength = [value for value in z.values()] + for j in range(1, tmplength[0] +1): + value[i-1][j-1] = self.msx.MSXgetpatternvalue(i, j) + if tmplength[0] < tmpmaxlen: + for j in range(tmplength + 1, tmpmaxlen + 1): + value[i - 1][j - 1] = value[i - 1][j - tmplength - 1] + + return value + + def getMSXSpeciesType(self, varagin=None): + + x = self.getMSXSpeciesCount() + value = [] + + if varagin is None: + varagin = {} + for i in range(1, x + 1): + varagin[i] = i + 1 + if x > 0: + for i in varagin: + y = {} + y = self.msx.MSXgetspecies(i) + value.append(y[0]) + return value + + def getMSXSpeciesUnits(self, varagin=None): + x = self.getMSXSpeciesCount() + value = [] + + if varagin is None: + varagin = {} + for i in range(1, x + 1): + varagin[i] = i + 1 + if x > 0: + for i in varagin: + y = {} + y = self.msx.MSXgetspecies(i) + value.append(y[1]) + return value + + def getEquations(self): + msxname = self.msxname + Terms = {} + Pipes = {} + Tanks = {} + with open(msxname, 'r') as f: + + sect = 0 + i = 1 + t = 1 + k = 1 + while True: + tline = f.readline() + if not tline: + break + tline = tline.strip() + if not tline: + continue + tok = tline.split()[0] + + if not tok: + continue + if tok[0] == ';': + continue + + if tok[0] == '[': + + if tok[1:6].upper() == 'TERMS': + sect = 1 + continue + elif tok[1:6].upper() == 'PIPE]': + sect = 2 + continue + elif tok[1:6].upper() == 'TANK]': + sect = 3 + continue + elif tok[1:6].upper() == '[END': + break + else: + sect = 0 + continue + if sect == 0: + continue + + elif sect == 1: + Terms[i] = tline + i = i + 1 + elif sect == 2: + Pipes[t] = tline + t = t + 1 + elif sect == 3: + Tanks[k] = tline + k = k + 1 + return Terms, Pipes, Tanks + + def getMSXEquationsTerms(self): + x, y, z = self.getEquations() + x = list(x.values()) + return x + + def getMSXEquationsPipes(self): + x, y, z = self.getEquations() + y = list(y.values()) + return y + + def getMSXEquationsTanks(self): + x, y, z = self.getEquations() + z = list(z.values()) + return z + + def getMSXSources(self): + value = [] + for i in range(1, self.getNodeCount() + 1): + value_row = [] + for j in range(1, self.getMSXSpeciesCount() + 1): + y = self.msx.MSXgetsource(i, j) + value_row.append(y) + value.append(value_row) + return value + + def getMSXSourceType(self, varagin=None): + value = [] + if varagin == None: + for i in range(1, self.getNodeCount() + 1): + value_row = [] + for j in range(1, self.getMSXSpeciesCount() + 1): + y = self.msx.MSXgetsource(i, j) + value_row.append(y) + value.append(value_row) + else: + for i in varagin: + value_row = [] + for j in range(1, self.getMSXSpeciesCount() + 1): + y = self.msx.MSXgetsource(i, j) + value_row.append(y) + value.append(value_row) + source = [] + for i in value: + source.append([item[0] for item in i]) + return source + + def getMSXSourceLevel(self, varagin=None): + value = [] + if varagin == None: + for i in range(1, self.getNodeCount() + 1): + value_row = [] + for j in range(1, self.getMSXSpeciesCount() + 1): + y = self.msx.MSXgetsource(i, j) + value_row.append(y) + value.append(value_row) + else: + for i in varagin: + value_row = [] + for j in range(1, self.getMSXSpeciesCount() + 1): + y = self.msx.MSXgetsource(i, j) + value_row.append(y) + value.append(value_row) + sourcelevel = [] + for i in value: + sourcelevel.append([item[1] for item in i]) + return sourcelevel + + def getMSXSourcePatternIndex(self, varagin=None): + value = [] + if varagin == None: + for i in range(1, self.getNodeCount() + 1): + value_row = [] + for j in range(1, self.getMSXSpeciesCount() + 1): + y = self.msx.MSXgetsource(i, j) + value_row.append(y) + value.append(value_row) + else: + for i in varagin: + value_row = [] + for j in range(1, self.getMSXSpeciesCount() + 1): + y = self.msx.MSXgetsource(i, j) + value_row.append(y) + value.append(value_row) + sourcepatternindex = [] + for i in value: + sourcepatternindex.append([item[2] for item in i]) + return sourcepatternindex + + def getMSXLinkInitqualValue(self, varagin=None): + value = [] + if varagin == None: + for i in range(1, self.getLinkCount() + 1): + value_row = [] + for j in range(1, self.getMSXSpeciesCount() + 1): + y = self.msx.MSXgetinitqual(1, i, j) + value_row.append(y) + value.append(value_row) + else: + for i in varagin: + value_row = [] + for j in range(1, self.getMSXSpeciesCount() + 1): + y = self.msx.MSXgetinitqual(1, i, j) + value_row.append(y) + value.append(value_row) + return value + + def getMSXNodeInitqualValue(self, varagin=None): + value = [] + if varagin == None: + for i in range(1, self.getNodeCount() + 1): + value_row = [] + for j in range(1, self.getMSXSpeciesCount() + 1): + y = self.msx.MSXgetinitqual(0, i, j) + value_row.append(y) + value.append(value_row) + else: + for i in varagin: + value_row = [] + for j in range(1, self.getMSXSpeciesCount() + 1): + y = self.msx.MSXgetinitqual(0, i, j) + value_row.append(y) + value.append(value_row) + return value + + def getMSXSpeciesATOL(self): + value = [] + for i in range(1, self.getMSXSpeciesCount()): + Atol = [] + value.append(self.msx.MSXgetspecies(i)) + Atol.append([item[2] for item in value]) + return Atol + + def getMSXSpeciesRTOL(self): + value = [] + for i in range(1, self.getMSXSpeciesCount()): + Rtol = [] + value.append(self.msx.MSXgetspecies(i)) + Rtol.append([item[3] for item in value]) + return Rtol + + def getMSXSpeciesConcentration(self, type, index, species): + return self.msx.MSXgetqual(type, index, species) + + def getMSXSourceNodeNameID(self): + nodes = [] + for i in range(1, self.getNodeCount() + 1): + source = [] + value = [] + flag = 0 + value_row = [] + for j in range(1, self.getMSXSpeciesCount() + 1): + y = self.msx.MSXgetsource(i, j) + value_row.append(y) + value.append(value_row) + for k in value: + source.append([item[0] for item in k]) + for sublist in source: + for item in sublist: + if item!='NOSOURCE': + flag = 1 + break + if flag == 1: + nodes.append(i) + return nodes + + def MSXPythonSetup(self, msxname): + + self.msxname = msxname[:-4] + '_temp.msx' + copyfile(msxname,self.msxname) + + self.msx.MSXopen(self.msxname) + + self.MSXEquationsTerms = self.getMSXEquationsTerms() + self.MSXEquationsPipes = self.getMSXEquationsPipes() + self.MSXEquationsTanks = self.getMSXEquationsTanks() + self.MSXSpeciesCount = self.getMSXSpeciesCount() + self.MSXConstantsCount = self.getMSXConstantsCount() + self.MSXParametersCount = self.getMSXParametersCount() + self.MSXPatternsCount = self.getMSXPatternsCount() + self.MSXSpeciesIndex = self.getMSXSpeciesIndex() + self.MSXSpeciesNameID = self.getMSXSpeciesNameID() + self.MSXSpeciesType = self.getMSXSpeciesType() + self.MSXSpeciesUnits = self.getMSXSpeciesUnits() + self.MSXSpeciesATOL = self.getMSXSpeciesATOL() + self.MSXSpeciesRTOL = self.getMSXSpeciesRTOL() + self.MSXConstantsNameID = self.getMSXConstantsNameID() + self.MSXConstantsValue = self.getMSXConstantsValue() + self.MSXConstantsIndex = self.getMSXConstantsIndex() + self.MSXParametersNameID = self.getMSXParametersNameID() + self.MSXParametersIndex = self.getMSXParametersIndex() + self.MSXParametersTanksValue = self.getMSXParametersTanksValue() + self.MSXParametersPipesValue = self.getMSXParametersPipesValue() + self.MSXPatternsNameID = self.getMSXPatternsNameID() + self.MSXPatternsIndex = self.getMSXPatternsIndex() + self.MSXPatternsLengths = self.getMSXPatternsLengths() + self.MSXNodeInitqualValue = self.getMSXNodeInitqualValue() + self.MSXLinkInitqualValue = self.getMSXLinkInitqualValue() + self.MSXSources = self.getMSXSources() + self.MSXSourceType = self.getMSXSourceType() + self.MSXSourceLevel = self.getMSXSourceLevel() + self.MSXSourcePatternIndex = self.getMSXSourcePatternIndex() + self.MSXSourceNodeNameID = self.getMSXSourceNodeNameID() + self.MSXPattern = self.getMSXPattern() + + #options + self.solver = self.getMSXSolver() + self.areaunits = self.getMSXAreaUnits() + self.rateunits = self.getMSXRateUnits() + self.rtol = self.getMSXRtol() + self.atol = self.getMSXAtol() + self.timestep = self.getMSXTimeStep() + self.coupling = self.getMSXCoupling() + self.compiler = self.getMSXCompiler() + + def setMSXOptions(self, *args): + + + for i in range(len(args) // 2): + argument = args[2 * i].lower() + if argument == 'areaunits': + self.areaunits = args[2 * i + 1] + self.changeMSXOptions("AREA_UNITS",args[2 * i + 1]) + elif argument == 'rateunits': + self.rateunits = args[2 * i + 1] + self.changeMSXOptions("RATE_UNITS", args[2 * i + 1]) + elif argument == 'solver': + self.solver = args[2 * i + 1] + self.changeMSXOptions("SOLVER", args[2 * i + 1]) + elif argument == 'timestep': + self.timestep = args[2 * i + 1] + self.changeMSXOptions("TIMESTEP", args[2 * i + 1]) + elif argument == 'atol': + self.atol = args[2 * i + 1] + self.changeMSXOptions("ATOL", args[2 * i + 1]) + elif argument == 'rtol': + self.rtol = args[2 * i + 1] + self.changeMSXOptions("RTOL", args[2 * i + 1]) + elif argument == 'coupling': + self.coupling = args[2 * i + 1] + self.changeMSXOptions("COUPLING", args[2 * i + 1]) + elif argument == 'compiler': + self.compiler = args[2 * i + 1] + self.changeMSXOptions("COMPILER", args[2 * i + 1]) + else: + print('Invalid property found.') + return + + def changeMSXOptions(self, param, change): + msxname=self.msxname + f = open(msxname, 'r+') + lines = f.readlines() + for i, line in enumerate(lines): + + if line.startswith(param): + + lines[i] = param + "\t" + str(change) + "\n" + + f.seek(0) + f.writelines(lines) + f.close() + + def setMSXAreaUnitsCM2(self): + self.changeMSXOptions("AREA_UNITS","CM2") + def setMSXAreaUnitsFT2(self): + self.changeMSXOptions("AREA_UNITS", "FT2") + + def setMSXAreaUnitsM2(self): + self.changeMSXOptions("AREA_UNITS", "M2") + + def setMSXAtol(self, value): + self.changeMSXOptions("ATOL", value) + + def setMSXRtol(self, value): + self.changeMSXOptions("RTOL", value) + + def setMSXCompilerGC(self): + self.changeMSXOptions("COMPILER", "GC") + + def setMSXCompilerVC(self): + self.changeMSXOptions("COMPILER", "VC") + + def setMSXCompilerNONE(self): + self.changeMSXOptions("COMPILER", "NONE") + + def setMSXCouplingFULL(self): + self.changeMSXOptions("COUPLING", "FULL") + + def setMSXCouplingNONE(self): + self.changeMSXOptions("COUPLING", "NONE") + + def setMSXRateUnitsDAY(self): + self.changeMSXOptions("RATE_UNITS", "DAY") + + def setMSXRateUnitsHR(self): + self.changeMSXOptions("RATE_UNITS", "HR") + + def setMSXRateUnitsMIN(self): + self.changeMSXOptions("RATE_UNITS", "MIN") + + def setMSXRateUnitsSEC(self): + self.changeMSXOptions("RATE_UNITS", "SEC") + + def setMSXSolverEUL(self): + self.changeMSXOptions("SOLVER", "EUL") + + def setMSXSolverRK5(self): + self.changeMSXOptions("SOLVER", "RK5") + + def setMSXSolverROS2(self): + self.changeMSXOptions("SOLVER", "ROS2") + + def setMSXTimeStep(self, value): + self.changeMSXOptions("TIMESTEP", value) + + def setMSXPatternValue(self, index, patternTimeStep, patternFactor): + self.msx.MSXsetpatternvalue(index, patternTimeStep, patternFactor) + + def setMSXPattern(self, index, patternVector): + nfactors = len(patternVector) + self.msx.MSXsetpattern(index, patternVector, nfactors) + + def setMSXParametersTanksValue(self, NodeTankIndex, paramindex, value): + self.msx.MSXsetparameter(0, NodeTankIndex, paramindex, value) + + def setMSXParametersPipesValue(self, pipeIndex, value): + for i in range(len(value)): + self.msx.MSXsetparameter(1, pipeIndex, i+1, value[i]) + + def setMSXConstantsValue(self, value): + for i in range(len(value)): + self.msx.MSXsetconstant(i+1, value[i]) + + def addMSXpattern(self, *args): + index = -1 + MSX_PATTERN = self.ToolkitConstants.MSX_PATTERN + if len(args) == 1: + self.msx.MSXaddpattern(args[0]) + index = self.msx.MSXgetindex(MSX_PATTERN, args[0]) + elif len(args) == 2: + self.msx.MSXaddpattern(args[0]) + index = self.msx.MSXgetindex(MSX_PATTERN, args[0]) + self.setMSXPattern(index, args[1]) + return index + + def getMSXComputedQualitySpecie(self, *args): + if self.getMSXSpeciesCount() == 0: + return 0 + if not args: + specie = self.getMSXspeciesNameID() + else: + specie = args[0] + + link_indices = range(1, self.getLinkCount() + 1) + node_indices = range(1, self.getNodeCount() + 1) + speciename = self.getMSXSpeciesIndex(specie) + + node_quality = np.empty((1, len(node_indices), len(speciename))) + link_quality = np.empty((1, len(node_indices), len(speciename))) + + self.solveMSXCompleteHydraulics() + self.initializeMSXQualityAnalysis(0) + + k = 1 + tleft = 1 + t = 0 + time = [0] + if node_indices[-1] < link_indices[-1]: + for i in range(len(speciename)): + for lnk in link_indices: + print(lnk) + link_quality[k, lnk - 1, i] = self.getMSXLinkInitqualValue([lnk])[(speciename)[1]] + if lnk < node_indices[-1] + 1: + node_quality[k, lnk - 1, i] = self.getMSXNodeInitqualValue([lnk])[speciename[i]] + else: + for i in range(len(speciename)): + for lnk in node_indices: + node_quality[k, lnk - 1, i] = self.MSXNodeInitqualValue([lnk])[speciename[i]] + if lnk < link_indices[-1] + 1: + link_quality[k, lnk - 1, i]= self.getMSXLinkInitqualValue([lnk])[speciename[i]] + time_sim = self.getTimeSimulationDuration() + while tleft > 0 and time_sim != t: + k = k + 1 + t, tleft = self.stepMSXQualityAnalysisTimeLeft() + for i in range(len(speciename)): + if node_indices[-1] < link_indices[-1]: + for lnk in link_indices: + link_quality[k, lnk -1, i] = self.getMSXSpeciesConcentration(1, lnk, speciename[i]) + if lnk < node_indices[-1] + 1: + node_quality[k, lnk - 1, i] = self.getMSXSpeciesConcentration(0, lnk, speciename[i]) + else: + for lnk in node_indices: + node_quality[k, lnk - 1, i] = self.getMSXSpeciesConcentration(0, lnk, speciename[i]) + if lnk < link_indices[-1] + 1: + link_quality[k, lnk - 1, i] = self.getMSXSpeciesConcentration(1, lnk, speciename[i]) + time.append(t) + return {'NodeQuality': node_quality, 'LinkQuality': link_quality, 'Time' : np.array(time)} + class epanetapi: """ @@ -13967,3 +14790,630 @@ def ENwriteline(self, line): self.errcode = self._lib.ENwriteline(line.encode("utf-8")) self.ENgeterror() + + +class epanetmsxapi: + """example msx = epanetmsxapi()""" + + def __init__(self, filename=None): + ops = platform.system().lower() + if ops in ["windows"]: + dll_path1 = resource_filename("epyt", os.path.join("libraries", "win", 'epanet2_2', '64bit', + f"epanetmsx.dll")) + elif ops in ["darwin"]: + dll_path1 = resource_filename("epyt", os.path.join("libraries", "mac", 'epanet2_2', '64bit', + f"epanetmsx.dll")) + else: + dll_path1 = resource_filename("epyt", os.path.join("libraries", "glnx", 'epanet2_2', '64bit', + f"epanetmsx.dll")) + + self.msx_lib = cdll.LoadLibrary(dll_path1) + # msx opens starts here + if filename is not None: + """ Open .msx file + msx.MSXopen(filename) + msx.MSXopen(Arsenite.msx)""" + """ filename example : Arsinite.msx or use full path """ + print("Opening MSX file:", filename) + if not os.path.exists(filename): + raise FileNotFoundError(f"File not found: {filename}") + + filename = c_char_p(filename.encode('utf-8')) + err = self.msx_lib.MSXopen(filename) + if err != 0: + self.MSXerror(err) + if err == 503: + print("Error 503 may indicate a problem with the MSX file or the MSX library.") + else: + print("MSX file opened successfully.") + # msx open ends here + + # Error ~ function + self.msx_error = self.msx_lib.MSXgeterror + self.msx_error.argtypes = [c_int, c_char_p, c_int] + + def MSXopen(self, filename): + """ Open .msx file + msx.MSXopen(filename) + msx.MSXopen(Arsenite.msx)""" + """ filename example : Arsinite.msx or use full path """ + print("Opening MSX file:", filename) + if not os.path.exists(filename): + raise FileNotFoundError(f"File not found: {filename}") + + filename = c_char_p(filename.encode('utf-8')) + err = self.msx_lib.MSXopen(filename) + if err != 0: + self.MSXerror(err) + if err == 503: + print("Error 503 may indicate a problem with the MSX file or the MSX library.") + else: + print("MSX file opened successfully.") + # msx open ends here + + def MSXclose(self): + """ Close .msx file + example : msx.MSXclose()""" + err = self.msx_lib.MSXclose() + if err != 0: + self.MSXerror(err) + return err + + def MSXerror(self, err_code): + """ Function that every other function uses in case of an error """ + errmsg = create_string_buffer(256) + self.msx_error(err_code, errmsg, 256) + print(errmsg.value.decode()) + + def MSXgetindex(self, obj_type, obj_id): + """ Retrieves the number of objects of a specific type + MSXgetcount(obj_type, obj_id) + + Parameters: + obj_type: code type of object being sought and must be one of the following + pre-defined constants: + MSX_SPECIES (for a chemical species) the number 3 + MSX_CONSTANT (for a reaction constant) the number 6 + MSX_PARAMETER (for a reaction parameter) the number 5 + MSX_PATTERN (for a time pattern) the number 7 + + obj_id: string containing the object's ID name + Returns: + The index number (starting from 1) of object of that type with that specific name.""" + obj_type = c_int(obj_type) + # obj_id=c_char_p(obj_id) + index = c_int() + err = self.msx_lib.MSXgetindex(obj_type, obj_id.encode("utf-8"), byref(index)) + if err != 0: + Warning(self.MSXerror(err)) + return index.value + + def MSXgetID(self, obj_type, index, id_len=80): + """ Retrieves the ID name of an object given its internal + index number + msx.MSXgetID(obj_type, index, id_len) + print(msx.MSXgetID(3,1,8)) + + Parameters: + obj_type: type of object being sought and must be on of the + following pre-defined constants: + MSX_SPECIES (for chemical species) + MSX_CONSTANT(for reaction constant) + MSX_PARAMETER(for a reaction parameter) + MSX_PATTERN (for a time pattern) + + index: the sequence number of the object (starting from 1 + as listed in the MSX input file) + + id_len: the maximum number of characters that id can hold + + Returns: + id object's ID name""" + + obj_id = create_string_buffer(id_len + 1) + err = self.msx_lib.MSXgetID(obj_type, index, obj_id, id_len) + if err != 0: + Warning(self.MSXerror(err)) + return obj_id.value.decode() + + def MSXgetIDlen(self, obj_type, index): + """Retrieves the number of characters in the ID name of an MSX + object given its internal index number + msx.MSXgetIDlen(obj_type, index) + print(msx.MSXgetIDlen(3,3)) + Parameters: + obj_type: type of object being sought and must be on of the + following pre-defined constants: + MSX_SPECIES (for chemical species) + MSX_CONSTANT(for reaction constant) + MSX_PARAMETER(for a reaction parameter) + MSX_PATTERN (for a time pattern) + + index: the sequence number of the object (starting from 1 + as listed in the MSX input file) + + Returns : the number of characters in the ID name of MSX object + + """ + len = c_int() + err = self.msx_lib.MSXgetIDlen(obj_type, index, byref(len)) + if err: + Warning(self.MSXerror(err)) + return len.value + + def MSXgetspecies(self, index): + """ Retrieves the attributes of a chemical species given its + internal index number + msx.MSXgetspecies(index) + msx.MSXgetspecies(1) + Parameters: + index : integer -> sequence number of the species + + Returns: + type : is returned with one of the following pre-defined constants: + MSX_BULK (defined as 0) for a bulk water species , or + MSX_WALL (defined as 1) for a pipe wall surface species + units: mass units that were defined for the species in question + atol : the absolute concentration tolerance defined for the species. + rtol : the relative concentration tolerance defined for the species. """ + type = c_int() + units = create_string_buffer(16) + atol = c_double() + rtol = c_double() + + err = self.msx_lib.MSXgetspecies( + index, byref(type), units, byref(atol), byref(rtol)) + + if type.value == 0: + type = 'BULK' + elif type.value == 1: + type = 'WALL' + + if err: + Warning(self.MSXerror(err)) + return type, units.value.decode("utf-8"), atol.value, rtol.value + + def MSXgetcount(self, code): + """ Retrieves the number of objects of a specific type + MSXgetcount(code) + + Parameters: + code type of object being sought and must be one of the following + pre-defined constants: + MSX_SPECIES (for a chemical species) the number 3 + MSX_CONSTANT (for a reaction constant) the number 6 + MSX_PARAMETER (for a reaction parameter) the number 5 + MSX_PATTERN (for a time pattern) the number 7 + Returns: + The count number of object of that type. + """ + count = c_int() + err = self.msx_lib.MSXgetcount(code, byref(count)) + if err: + Warning(self.MSXerror(err)) + return count.value + + def MSXgetconstant(self, index): + """ Retrieves the value of a particular rection constant """ + """msx.MSXgetconstant(index) + msx.MSXgetconstant(1)""" + """" Parameters: + index : integer is the sequence number of the reaction + constant ( starting from 1 ) as it + appeared in the MSX input file + + Returns: value -> the value assigned to the constant. """ + value = c_double() + err = self.msx_lib.MSXgetconstant(index, byref(value)) + if err: + Warning(self.MSXerror(err)) + return value.value + + def MSXgetparameter(self, obj_type, index, param): + """Retrieves the value of a particular reaction parameter for a given + pipe + msx.MSXgetparameter(obj_type, index, param) + msx.MSXgetparameter(1,1,1) + Parameters: + obj_type: is type of object being queried and must be either: + MSX_NODE (defined as 0) for a node or + MSX_LINK(defined as 1) for alink + + index: is the internal sequence number (starting from 1) + assigned to the node or link + + param: the sequence number of the parameter (starting from 1 + as listed in the MSX input file) + + Returns: + value : the value assigned to the parameter for the node or link + of interest. """ + value = c_double() + err = self.msx_lib.MSXgetparameter(obj_type, index, param, byref(value)) + if err: + Warning(self.MSXerror(err)) + return value.value + + def MSXgetpatternlen(self, pattern_index): + """Retrieves the number of time periods within a source time pattern + + MSXgetpatternlen(pattern_index) + + Parameters: + pattern_index: the internal sequence number (starting from 1) + of the pattern as it appears in the MSX input file. + + Returns: + len: the number of time periods (and therefore number of multipliers) + that appear in the pattern.""" + len = c_int() + err = self.msx_lib.MSXgetpatternlen(pattern_index, byref(len)) + if err: + Warning(self.MSXerror(err)) + return len.value + + def MSXgetpatternvalue(self, pattern_index, period): + """ Retrieves the multiplier at a specific time period for a + given source time pattern + msx.MSXgetpatternvalue(pattern_index, period) + msx.MSXgetpatternvalue(1,1) + Parameters: + pattern_index: the internal sequence number(starting from 1) + of the pattern as it appears in the MSX input file + + period: the index of the time period (starting from 1) whose + multiplier is being sought """ + value = c_double() + err = self.msx_lib.MSXgetpatternvalue(pattern_index, period, byref(value)) + if err: + Warning(self.MSXerror(err)) + return value.value + + def MSXgetinitqual(self, obj_type, index, species): + """ Retrieves the intial concetration of a particular chemical species + assigned to a specific node or link of the pipe network + msx.MSXgetinitqual(obj_type, index) + msx.MSXgetinitqual(1,1,1) + Parameters: + + type : type of object being queeried and must be either: + MSX_NODE (defined as 0) for a node or , + MSX_LINK (defined as 1) for a link + + index : the internal sequence number (starting from 1) assigned + to the node or link + + species: the sequence number of the species (starting from 1) + + Returns: + value: the initial concetration of the species at the node or + link of interest.""" + # obj_type = c_int(obj_type) + value = c_double() + # species = c_int(species) + # index = c_int(index) + err = self.msx_lib.MSXgetinitqual(obj_type, index, species, byref(value)) + if err: + Warning(self.MSXerror(err)) + return value.value + + def MSXgetsource(self, node_index, species_index): + """ Retrieves information on any external source of a particular + chemical species assigned to a specific node or link of the pipe + network. + msx.MSXgetsource(node_index, species_index) + msx.MSXgetsource(1,1) + + Parameters: + node_index: the internal sequence number (starting from 1) + assigned to the node of interest. + + species_index: the sequence number of the species of interest + (starting from 1 as listed in MSX input file) + Returns: + + type: the type of external source to be utilized and will be one of + the following predefined constants: + MSX_NOSOURCE (defined as -1) for no source + MSX_CONCEN (defined as 0) for a concetration sourc + MSX_MASS (defined as 1) for a mass booster source + MSX_SETPOINT (defined as 2) for a setpoint source + MSX_FLOWPACE (defined as 3) for a flow paced source + + level: the baseline concentration ( or mass flow rate) of the source) + + pat : the index of the time pattern used to add variability to the + the source's baseline level (and will be 0 if no pattern + was defined for the source) + """ + type = c_int() + level = c_double() + pattern = c_int() + node_index = c_int(node_index) + err = self.msx_lib.MSXgetsource(node_index, species_index, + byref(type), byref(level), byref(pattern)) + + if type.value == -1: + type = 'NOSOURCE' + elif type.value == 0: + type = 'CONCEN' + elif type.value == 1: + type = 'MASS' + elif type.value == 2: + type = 'SETPOINT' + elif type.value == 3: + type = 'FLOWPACED' + + if err: + Warning(self.MSXerror(err)) + + return type, level.value, pattern.value + + def MSXsaveoutfile(self, filename): + """ Saves water quality results computed for each node, link + and reporting time period to a named binary file. + msx.MSXsaveoutfile(filename) + msx.MSXsaveoufile(Arsenite.msx) + + Parameters: + filename: name of the permanent output results file""" + err = self.msx_lib.MSXsaveoutfile(filename.encode()) + if err: + Warning(self.MSXerror(err)) + + def MSXsavemsxfile(self, filename): + """ Saves the data associated with the current MSX project into a new + MSX input file + msx.MSXsavemsxfile(filename) + msx.MSXsavemsxfile(Arsenite.msx) + + Parameters: + filename: name of the file to which data are saved""" + err = self.msx_lib.MSXsavemsxfile(filename.encode()) + if err: + Warning(self.MSXerror(err)) + + def MSXsetconstant(self, index, value): + """ Assigns a new value to a specific reaction constant + msx.MSXsetconstant(index, value) + msx.MSXsetconstant(1,10)""" + """" Parameters + index : integer -> is the sequence number of the reaction + constant ( starting from 1 ) as it appeared in the MSX + input file + + Value: float -> the new value to be assigned to the constant.""" + + value = c_double(value) + err = self.msx_lib.MSXsetconstant(index, value) + if err: + Warning(self.MSXerror(err)) + + def MSXsetparameter(self, obj_type, index, param, value): + """ Assigns a value to a particular reaction parameter for a given pipe + or tank within the pipe network + msx.MSXsetparameter(obj_type, index, param, value) + msx.MSXsetparameter(1,1,1,15) + Parameters: + obj_type: is type of object being queried and must be either: + MSX_NODE (defined as 0) for a node or + MSX_LINK (defined as 1) for a link + + index: is the internal sequence number (starting from 1) + assigned to the node or link + + param: the sequence number of the parameter (starting from 1 + as listed in the MSX input file) + + value: the value to be assigned to the parameter for the node or + link of interest. """ + value = c_double(value) + err = self.msx_lib.MSXsetparameter(obj_type, index, param, value) + if err: + Warning(self.MSXerror(err)) + + def MSXsetinitqual(self, obj_type, index, species, value): + """ Assigns an initial concetration of a particular chemical species + node or link of the pipe network + msx.MSXsetinitqual(obj_type, index, species, value) + msx.MSXsetinitqual(1,1,1,15) + Parameters: + type: type of object being queried and must be either : + MSX_NODE(defined as 0) for a node or + MSX_LINK(defined as 1) for a link + index: integer -> the internal sequence number (starting from 1) + assigned to the node or link + + species: the sequence number of the species (starting from 1 as listed in + MASx input file) + + value: float -> the initial concetration of the species to be applied at the node or link + of interest. + """ + + value = c_double(value) + err = self.msx_lib.MSXsetinitqual(obj_type, index, species, value) + if err: + Warning(self.MSXerror(err)) + + def MSXsetpattern(self, index, factors, nfactors): + """Assigns a new set of multipliers to a given MSX source time pattern + MSXsetpattern(index,factors,nfactors) + + Parameters: + index: the internal sequence number (starting from 1) + of the pattern as it appers in the MSX input file + factors: an array of multiplier values to replace those previously used by + the pattern + nfactors: the number of entries in the multiplier array/ vector factors""" + index = c_int(index) + nfactors = c_int(nfactors) + DoubleArray = c_double * len(factors) + mult_array = DoubleArray(*factors) + err = self.msx_lib.MSXsetpattern(index, mult_array, nfactors) + if err: + Warning(self.MSXerror(err)) + + def MSXsetpatternvalue(self, pattern, period, value): + """Assigns a new value to the multiplier for a specific time period + in a given MSX source time pattern. + msx.MSXsetpatternvalue(pattern, period, value) + msx.MSXsetpatternvalue(1,1,10) + Parameters: + pattern: the internal sequence number (starting from 1) of the + pattern as it appears in the MSX input file. + + period: the time period (starting from 1) in the pattern to be replaced + value: the new multiplier value to use for that time period.""" + value = c_double(value) + err = self.msx_lib.MSXsetpatternvalue(pattern, period, value) + if err: + Warning(self.MSXerror(err)) + + def MSXsolveQ(self): + """ Solves for water quality over the entire simulation period + and saves the results to an internal scratch file + msx.MSXsolveQ()""" + err = self.msx_lib.MSXsolveQ() + if err: + Warning(self.MSXerror(err)) + + def MSXsolveH(self): + """ Solves for system hydraulics over the entire simulation period + saving results to an internal scratch file + msx.MSXsolveH() """ + err = self.msx_lib.MSXsolveH() + if err: + Warning(self.MSXerror(err)) + + def MSXaddpattern(self, pattern_id): + """Adds a newm empty MSX source time pattern to an MSX project + MSXaddpattern(pattern_id) + Parameters: + pattern_id: the name of the new pattern """ + err = self.msx_lib.MSXaddpattern(pattern_id.encode("utf-8")) + if err: + Warning(self.MSXerror(err)) + + def MSXusehydfile(self, filename): + """ """ + err = self.msx_lib.MSXusehydfile(filename.encode()) + if err: + Warning(self.MSXerror(err)) + + def MSXstep(self): + """Advances the water quality solution through a single water quality time + step when performing a step-wise simulation + + t, tleft = MSXstep() + Returns: + t : current simulation time at the end of the step(in secconds) + tleft: time left in the simulation (in secconds) + """ + t = c_int() + tleft = c_int() + err = self.msx_lib.MSXstep(byref(t), byref(tleft)) + + if err: + Warning(self.MSXerror(err)) + + return t.value, tleft.value + + def MSXinit(self, flag): + """Initialize the MSX system before solving for water quality results + in the step-wise fashion + + MSXinit(flag) + + Parameters: + flag: Set the flag to 1 if the water quality results should be saved + to a scratch binary file, or 0 if not + + """ + err = self.msx_lib.MSXinit(flag) + if err: + Warning(self.MSXerror(err)) + + def MSXreport(self): + """ Writes water quality simulations results as instructed by + MSX input file to a text file. + msx.MSXreport()""" + err = self.msx_lib.MSXreport() + if err: + Warning(self.MSXerror(err)) + + def MSXgetqual(self, type, index, species): + """Retrieves a chemical species concentration at a given node + or the average concentration along a link at the current sumulation + time step. + + MSXgetqual(type, index, species) + + Parameters: + type: type of object being queried and must be either: + MSX_NODE ( defined as 0) for a node, + MSX_LINK (defined as 1) for a link + index: then internal sequence number (starting from 1) + assigned to the node or link. + species is the sequence number of the species (starting from 1 + as listed in the MSX input file) + + Returns: + The value of the computed concentration of the species at the current + time period. + """ + + value = 0 + value = c_double(value) + err = self.msx_lib.MSXgetqual(type, index, species, byref(value)) + if err: + Warning(self.MSXerror(err)) + return value.value + + def MSXsetsource(self, node, species, type, level, pat): + """"Sets the attributes of an external source of particular chemical + species to specific node of the pipe network + msx.setsource(node, species, type, level, pat) + msx.MSXsetsource(1,1,3,10.565,1) + Parameters: + node: the internal sequence number (starting from1) assigned + to the node of interest. + + species: the sequence number of the species of interest (starting + from 1 as listed in the MSX input file) + + type: the type of external source to be utilized and will be one of + the following predefined constants: + MSX_NOSOURCE (defined as -1) for no source + MSX_CONCEN (defined as 0) for a concetration source + MSX_MASS (defined as 1) for a mass booster source + MSX_SETPOINT (defined as 2) for a setpoint source + MSX_FLOWPACE (defined as 3) for a flow paced source + + level: the baseline concetration (or mass flow rate) of the source + + pat: the index of the time pattern used to add variability to the + source's baseline level ( use 0 if the source has a constant strength) """ + level = c_double(level) + pat = c_int(pat) + type = c_int(type) + err = self.msx_lib.MSXsetsource(node, species, type, level, pat) + if err: + Warning(self.MSXerror(err)) + + def MSXgeterror(self, err): + """Returns the text for an error message given its error code. + msx.MSXgeterror(err) + msx.MSXgeterror(516) + Parameters: + err: the code number of an error condition generated by EPANET-MSX + + Returns: + errmsg: the text of the error message corresponding to the error code""" + errmsg = create_string_buffer(80) + e = self.msx_lib.MSXgeterror(err, errmsg, 80) + + if e: + Warning(errmsg.value.decode()) + + return errmsg.value.decode() diff --git a/epyt/examples/python/Toolkit_api_EX2_using_MSX_functions.py b/epyt/examples/python/Toolkit_api_EX2_using_MSX_functions.py new file mode 100644 index 0000000..529a887 --- /dev/null +++ b/epyt/examples/python/Toolkit_api_EX2_using_MSX_functions.py @@ -0,0 +1,66 @@ +from epyt.epanet import epanet, epanetmsxapi +from epyt import networks +import numpy as np +import os + +# Create EPANET object using INP file and MSX object using MSX file +dirname = os.path.dirname(networks.__file__) +inpname = os.path.join(dirname, 'msx-examples', 'net2-cl2.inp') +msxname = os.path.join(dirname, 'msx-examples', 'net2-cl2.msx') + +d = epanet(inpname, msx=True) +msx = epanetmsxapi(msxname) +MSX_SPECIES = 3 +ss = list(range(1, d.LinkCount + 1)) +uu = list(range(1, msx.MSXgetcount(MSX_SPECIES) + 1)) + +# Initialized quality and time +link_count = d.getLinkCount() +msx_time_step = 300 +time_steps = int(d.getTimeSimulationDuration()/msx_time_step) +quality = [np.zeros((time_steps, 1)) for _ in range(link_count)] +time = np.zeros((time_steps, 1)) +data = { + 'Quality': quality, + 'Time': time +} + +# Obtain a hydraulic solution +msx.MSXsolveH() + +# Run a step-wise water quality analysis without saving results to file +msx.MSXinit(0) + +# Retrieve species concentration at node +for i, nl in enumerate(ss, start=1): + for j_idx, j in enumerate(uu, start=1): + data['Quality'][i - 1][j_idx - 1] = msx.MSXgetinitqual(1, i, j) + +k = 0 +tleft = 1 +# Initialized data time with 0 +data['Time'][k] = [0] +while tleft > 0: + t, tleft = msx.MSXstep() + if t > msx_time_step: + for i, nl in enumerate(ss, start=1): + for g, j in enumerate(uu, start=1): + x = msx.MSXgetqual(1, nl, j) + data['Quality'][i - 1][k, g - 1] = x + data['Time'][k] = t + k += 1 + +# Plot quality over time (in hrs) for link 1 +hrs_time = data['Time'] / 3600 +d.plot_ts(X=hrs_time, Y=data['Quality'][0], + title=f'Quality vs Time Link 1', legend_location='best', marker=None, + xlabel='Time (hrs)', ylabel=f'CL2 Concentration (ppm)', figure_size=[4, 3]) + +# Plot quality over time (in hrs) for link 36 +d.plot_ts(X=hrs_time, Y=data['Quality'][35], + title=f'Quality vs Time Link 36', legend_location='best', marker=None, + xlabel='Time (hrs)', ylabel=f'CL2 Concentration (ppm)', figure_size=[4, 3]) + +# Unload MSX library and EN library. +msx.MSXclose() +d.unload() diff --git a/epyt/libraries/win/epanet2_2/64bit/epanetmsx.dll b/epyt/libraries/win/epanet2_2/64bit/epanetmsx.dll new file mode 100644 index 0000000..e48ed3d Binary files /dev/null and b/epyt/libraries/win/epanet2_2/64bit/epanetmsx.dll differ diff --git a/epyt/libraries/win/epanet2_2/64bit/epanetmsx.exe b/epyt/libraries/win/epanet2_2/64bit/epanetmsx.exe new file mode 100644 index 0000000..1a64ac9 Binary files /dev/null and b/epyt/libraries/win/epanet2_2/64bit/epanetmsx.exe differ diff --git a/epyt/networks/msx-examples/Net3-NH2CL.inp b/epyt/networks/msx-examples/Net3-NH2CL.inp new file mode 100644 index 0000000..a522419 --- /dev/null +++ b/epyt/networks/msx-examples/Net3-NH2CL.inp @@ -0,0 +1,472 @@ +[TITLE] +EPANET Example Network 3 + +[JUNCTIONS] +;ID Elev Demand Pattern + 10 147 0 ; + 15 32 1 3 ; + 20 129 0 ; + 35 12.5 1 4 ; + 40 131.9 0 ; + 50 116.5 0 ; + 60 0 0 ; + 601 0 0 ; + 61 0 0 ; + 101 42 189.95 ; + 103 43 133.2 ; + 105 28.5 135.37 ; + 107 22 54.64 ; + 109 20.3 231.4 ; + 111 10 141.94 ; + 113 2 20.01 ; + 115 14 52.1 ; + 117 13.6 117.71 ; + 119 2 176.13 ; + 120 0 0 ; + 121 -2 41.63 ; + 123 11 1 2 ; + 125 11 45.6 ; + 127 56 17.66 ; + 129 51 0 ; + 131 6 42.75 ; + 139 31 5.89 ; + 141 4 9.85 ; + 143 -4.5 6.2 ; + 145 1 27.63 ; + 147 18.5 8.55 ; + 149 16 27.07 ; + 151 33.5 144.48 ; + 153 66.2 44.17 ; + 157 13.1 51.79 ; + 159 6 41.32 ; + 161 4 15.8 ; + 163 5 9.42 ; + 164 5 0 ; + 166 -2 2.6 ; + 167 -5 14.56 ; + 169 -5 0 ; + 171 -4 39.34 ; + 173 -4 0 ; + 177 8 58.17 ; + 179 8 0 ; + 181 8 0 ; + 183 11 0 ; + 184 16 0 ; + 185 16 25.65 ; + 187 12.5 0 ; + 189 4 107.92 ; + 191 25 81.9 ; + 193 18 71.31 ; + 195 15.5 0 ; + 197 23 17.04 ; + 199 -2 119.32 ; + 201 0.1 44.61 ; + 203 2 1 5 ; + 204 21 0 ; + 205 21 65.36 ; + 206 1 0 ; + 207 9 69.39 ; + 208 16 0 ; + 209 -2 0.87 ; + 211 7 8.67 ; + 213 7 13.94 ; + 215 7 92.19 ; + 217 6 24.22 ; + 219 4 41.32 ; + 225 8 22.8 ; + 229 10.5 64.18 ; + 231 5 16.48 ; + 237 14 15.61 ; + 239 13 44.61 ; + 241 13 0 ; + 243 14 4.34 ; + 247 18 70.38 ; + 249 18 0 ; + 251 30 24.16 ; + 253 36 54.52 ; + 255 27 40.39 ; + 257 17 0 ; + 259 25 0 ; + 261 0 0 ; + 263 0 0 ; + 265 0 0 ; + 267 21 0 ; + 269 0 0 ; + 271 6 0 ; + 273 8 0 ; + 275 10 0 ; + +[RESERVOIRS] +;ID Head Pattern + 4 220.0 ; + 5 167.0 ; + +[TANKS] +;ID Elevation InitLevel MinLevel MaxLevel Diameter MinVol VolCurve + 1 131.9 13.1 .1 32.1 85 0 ; + 2 116.5 23.5 6.5 40.3 50 0 ; + 3 129.0 29.0 4.0 35.5 164 0 ; + +[PIPES] +;ID Node1 Node2 Length Diameter Roughness MinorLoss Status + 20 3 20 99 99 199 0 Open ; + 40 1 40 99 99 199 0 Open ; + 50 2 50 99 99 199 0 Open ; + 60 4 60 1231 24 140 0 Open ; + 101 10 101 14200 18 110 0 Open ; + 103 101 103 1350 16 130 0 Open ; + 105 101 105 2540 12 130 0 Open ; + 107 105 107 1470 12 130 0 Open ; + 109 103 109 3940 16 130 0 Open ; + 111 109 111 2000 12 130 0 Open ; + 112 115 111 1160 12 130 0 Open ; + 113 111 113 1680 12 130 0 Open ; + 114 115 113 2000 8 130 0 Open ; + 115 107 115 1950 8 130 0 Open ; + 116 113 193 1660 12 130 0 Open ; + 117 263 105 2725 12 130 0 Open ; + 119 115 117 2180 12 130 0 Open ; + 120 119 120 730 12 130 0 Open ; + 121 120 117 1870 12 130 0 Open ; + 122 121 120 2050 8 130 0 Open ; + 123 121 119 2000 30 141 0 Open ; + 125 123 121 1500 30 141 0 Open ; + 129 121 125 930 24 130 0 Open ; + 131 125 127 3240 24 130 0 Open ; + 133 20 127 785 20 130 0 Open ; + 135 127 129 900 24 130 0 Open ; + 137 129 131 6480 16 130 0 Open ; + 145 129 139 2750 8 130 0 Open ; + 147 139 141 2050 8 130 0 Open ; + 149 143 141 1400 8 130 0 Open ; + 151 15 143 1650 8 130 0 Open ; + 153 145 141 3510 12 130 0 Open ; + 155 147 145 2200 12 130 0 Open ; + 159 147 149 880 12 130 0 Open ; + 161 149 151 1020 8 130 0 Open ; + 163 151 153 1170 12 130 0 Open ; + 169 125 153 4560 8 130 0 Open ; + 171 119 151 3460 12 130 0 Open ; + 173 119 157 2080 30 141 0 Open ; + 175 157 159 2910 30 141 0 Open ; + 177 159 161 2000 30 141 0 Open ; + 179 161 163 430 30 141 0 Open ; + 180 163 164 150 14 130 0 Open ; + 181 164 166 490 14 130 0 Open ; + 183 265 169 590 30 141 0 Open ; + 185 167 169 60 8 130 0 Open ; + 186 187 204 99.9 8 130 0 Open ; + 187 169 171 1270 30 141 0 Open ; + 189 171 173 50 30 141 0 Open ; + 191 271 171 760 24 130 0 Open ; + 193 35 181 30 24 130 0 Open ; + 195 181 177 30 12 130 0 Open ; + 197 177 179 30 12 130 0 Open ; + 199 179 183 210 12 130 0 Open ; + 201 40 179 1190 12 130 0 Open ; + 202 185 184 99.9 8 130 0 Open ; + 203 183 185 510 8 130 0 Open ; + 204 184 205 4530. 12 130 0 Open ; + 205 204 185 1325. 12 130 0 Open ; + 207 189 183 1350 12 130 0 Open ; + 209 189 187 500 8 130 0 Open ; + 211 169 269 646 12 130 0 Open ; + 213 191 187 2560 12 130 0 Open ; + 215 267 189 1230 12 130 0 Open ; + 217 191 193 520 12 130 0 Open ; + 219 193 195 360 12 130 0 Open ; + 221 161 195 2300 8 130 0 Open ; + 223 197 191 1150 12 130 0 Open ; + 225 111 197 2790 12 130 0 Open ; + 229 173 199 4000 24 141 0 Open ; + 231 199 201 630 24 141 0 Open ; + 233 201 203 120 24 130 0 Open ; + 235 199 273 725 12 130 0 Open ; + 237 205 207 1200 12 130 0 Open ; + 238 207 206 450 12 130 0 Open ; + 239 275 207 1430 12 130 0 Open ; + 240 206 208 510 12 130 0 Open ; + 241 208 209 885 12 130 0 Open ; + 243 209 211 1210 16 130 0 Open ; + 245 211 213 990 16 130 0 Open ; + 247 213 215 4285 16 130 0 Open ; + 249 215 217 1660 16 130 0 Open ; + 251 217 219 2050 14 130 0 Open ; + 257 217 225 1560 12 130 0 Open ; + 261 213 229 2200 8 130 0 Open ; + 263 229 231 1960 12 130 0 Open ; + 269 211 237 2080 12 130 0 Open ; + 271 237 229 790 8 130 0 Open ; + 273 237 239 510 12 130 0 Open ; + 275 239 241 35 12 130 0 Open ; + 277 241 243 2200 12 130 0 Open ; + 281 241 247 445 10 130 0 Open ; + 283 239 249 430 12 130 0 Open ; + 285 247 249 10 12 130 0 Open ; + 287 247 255 1390 10 130 0 Open ; + 289 50 255 925 10 130 0 Open ; + 291 255 253 1100 10 130 0 Open ; + 293 255 251 1100 8 130 0 Open ; + 295 249 251 1450 12 130 0 Open ; + 297 120 257 645 8 130 0 Open ; + 299 257 259 350 8 130 0 Open ; + 301 259 263 1400 8 130 0 Open ; + 303 257 261 1400 8 130 0 Open ; + 305 117 261 645 12 130 0 Open ; + 307 261 263 350 12 130 0 Open ; + 309 265 267 1580 8 130 0 Open ; + 311 193 267 1170 12 130 0 Open ; + 313 269 189 646 12 130 0 Open ; + 315 181 271 260 24 130 0 Open ; + 317 273 275 2230 8 130 0 Open ; + 319 273 205 645 12 130 0 Open ; + 321 163 265 1200 30 141 0 Open ; + 323 201 275 300 12 130 0 Open ; + 325 269 271 1290 8 130 0 Open ; + 329 61 123 45500 30 140 0 Open ; + 330 60 601 1 30 140 0 Closed ; + 333 601 61 1 30 140 0 Open ; + +[PUMPS] +;ID Node1 Node2 Parameters + 10 5 10 HEAD 1 ; + 335 60 61 HEAD 2 ; + +[VALVES] +;ID Node1 Node2 Diameter Type Setting MinorLoss + +[TAGS] + +[DEMANDS] +;Junction Demand Pattern Category + +[STATUS] +;ID Status/Setting + 10 Closed + +[PATTERNS] +;ID Multipliers +; + 1 1.34 1.94 1.46 1.44 .76 .92 + 1 .85 1.07 .96 1.1 1.08 1.19 + 1 1.16 1.08 .96 .83 .79 .74 + 1 .64 .64 .85 .96 1.24 1.67 +; + 2 0 0 0 0 0 1219 + 2 0 0 0 1866 1836 1818 + 2 1818 1822 1822 1817 1824 1816 + 2 1833 1817 1830 1814 1840 1859 +; + 3 620 620 620 620 620 360 + 3 360 0 0 0 0 360 + 3 360 360 360 360 0 0 + 3 0 0 0 0 360 360 +; + 4 1637 1706 1719 1719 1791 1819 + 4 1777 1842 1815 1825 1856 1801 + 4 1819 1733 1664 1620 1613 1620 + 4 1616 1647 1627 1627 1671 1668 +; + 5 4439 4531 4511 4582 4531 4582 + 5 4572 4613 4643 4643 4592 4613 + 5 4531 4521 4449 4439 4449 4460 + 5 4439 4419 4368 4399 4470 4480 + +[CURVES] +;ID X-Value Y-Value +;PUMP: Pump curve for Pump 10 + 1 0 104. + 1 2000. 92. + 1 4000. 63. +;PUMP: Pump curve for Pump 335 + 2 0 200. + 2 8000. 138. + 2 14000. 86. + +[CONTROLS] +Link 10 OPEN AT CLOCKTIME 6 am +Link 10 CLOSED AT CLOCKTIME 8 pm +Link 335 OPEN IF Node 1 BELOW 17.1 +Link 335 CLOSED IF Node 1 ABOVE 19.1 +Link 330 CLOSED IF Node 1 BELOW 17.1 +Link 330 OPEN IF Node 1 ABOVE 19.1 + + +[RULES] + +[ENERGY] + Global Efficiency 75 + Global Price 0 + Demand Charge 0 + +[EMITTERS] +;Junction Coefficient + +[QUALITY] +;Node InitQual + 5 1.0 + +[SOURCES] +;Node Type Quality Pattern + +[REACTIONS] +;Type Pipe/Tank Coefficient + + +[REACTIONS] + Order Bulk 1 + Order Tank 1 + Order Wall 1 + Global Bulk 0 + Global Wall 0 + Limiting Potential 0 + Roughness Correlation 0 + +[MIXING] +;Tank Model + +[TIMES] + Duration 360 + Hydraulic Timestep 1:00 + Quality Timestep 0:05 + Pattern Timestep 1:00 + Pattern Start 0:00 + Report Timestep 1:00 + Report Start 312 + Start ClockTime 6 AM + Statistic Average + +[REPORT] + Status No + Summary No + Page 0 + +[OPTIONS] + Units GPM + Headloss H-W + Specific Gravity 1 + Viscosity 1 + Trials 40 + Accuracy 0.001 + Unbalanced Continue 10 + Pattern 1 + Demand Multiplier 1.0 + Emitter Exponent 0.5 + Quality S5 mg/L + Diffusivity 1 + Tolerance 0.01 + +[COORDINATES] +;Node X-Coord Y-Coord + 10 9.00 27.85 + 15 38.68 23.76 + 20 29.44 26.91 + 35 25.46 10.52 + 40 27.02 9.81 + 50 33.01 3.01 + 60 23.90 29.94 + 601 23.00 29.49 + 61 23.71 29.03 + 101 13.81 22.94 + 103 12.96 21.31 + 105 16.97 21.28 + 107 18.45 20.46 + 109 17.64 18.92 + 111 20.21 17.53 + 113 22.04 16.61 + 115 20.98 19.18 + 117 21.69 21.28 + 119 23.70 22.76 + 120 22.08 23.10 + 121 23.54 25.50 + 123 23.37 27.31 + 125 24.59 25.64 + 127 29.29 26.40 + 129 30.32 26.39 + 131 37.89 29.55 + 139 33.28 24.54 + 141 35.68 23.08 + 143 37.47 21.97 + 145 33.02 19.29 + 147 30.24 20.38 + 149 29.62 20.74 + 151 28.29 21.39 + 153 28.13 22.63 + 157 24.85 20.16 + 159 23.12 17.50 + 161 25.10 15.28 + 163 25.39 14.98 + 164 25.98 15.14 + 166 26.48 15.13 + 167 25.88 12.98 + 169 25.68 12.74 + 171 26.65 11.80 + 173 26.87 11.59 + 179 25.71 10.40 + 181 25.72 10.74 + 183 25.45 10.18 + 184 25.15 9.52 + 185 25.01 9.67 + 187 23.64 11.04 + 189 24.15 11.37 + 191 22.10 14.07 + 193 22.88 14.35 + 195 23.18 14.72 + 197 20.97 15.18 + 199 29.42 8.44 + 201 30.89 8.57 + 203 31.14 8.89 + 204 23.80 10.90 + 205 29.20 6.46 + 206 31.66 6.64 + 207 31.00 6.61 + 208 32.54 6.81 + 209 33.76 6.59 + 211 34.20 5.54 + 213 35.26 6.16 + 215 39.95 8.73 + 217 42.11 8.67 + 219 44.86 9.32 + 225 43.53 7.38 + 229 36.16 3.49 + 231 38.38 2.54 + 237 35.37 3.08 + 239 35.76 2.31 + 241 35.87 2.11 + 243 37.04 0.00 + 247 35.02 2.05 + 249 35.02 1.81 + 251 34.15 1.10 + 253 32.16 2.06 + 255 33.51 2.45 + 257 21.17 23.32 + 259 20.80 23.40 + 261 20.79 21.45 + 263 20.32 21.57 + 265 25.39 13.60 + 267 23.38 12.95 + 269 25.03 12.14 + 271 25.97 11.00 + 273 29.16 7.38 + 275 31.07 8.29 + 4 24.15 31.06 + 5 8.00 27.53 + 1 27.46 9.84 + 2 32.99 3.45 + 3 29.41 27.27 + +[VERTICES] +;Link X-Coord Y-Coord + +[LABELS] +;X-Coord Y-Coord Label & Anchor Node + 8.00 30.00 "LAKE" 5 + 25.00 31.10 "RIVER" 4 + +[BACKDROP] + DIMENSIONS 6.16 -1.55 46.71 32.61 + UNITS None + FILE + OFFSET 0.00 0.00 + +[END] diff --git a/epyt/networks/msx-examples/Net3-NH2CL.msx b/epyt/networks/msx-examples/Net3-NH2CL.msx new file mode 100644 index 0000000..4aa86f1 --- /dev/null +++ b/epyt/networks/msx-examples/Net3-NH2CL.msx @@ -0,0 +1,109 @@ +[TITLE] +NET3 Chloramine Decay Example + +[OPTIONS] +AREA_UNITS FT2 +RATE_UNITS HR +SOLVER ROS2 +COUPLING NONE +COMPILER NONE +TIMESTEP 300 +RTOL 0.0001 +ATOL 1.0e-8 + +[SPECIES] +BULK HOCL MMOL +BULK NH3 MMOL +BULK NH2CL MMOL +BULK NHCL2 MMOL +BULK I MMOL +BULK OCL MMOL +BULK NH4 MMOL +BULK ALK MMOL ;total alkalinity +BULK H MMOL ;hydrogen ion +BULK OH MMOL ;hydroxide ion +BULK CO3 MMOL ;carbonate ion +BULK HCO3 MMOL ;bicarbonate ion +BULK H2CO3 MMOL ;carbonic acid +BULK TOC MG ;total organic carbon +BULK cNH3 MG ;ammonia as mg/L +BULK cNH2CL MG ;monochloramine as mg/L + + +[COEFFICIENTS] +PARAMETER k1 1.5e10 +PARAMETER k2 7.6e-2 +PARAMETER k3 1.0e6 +PARAMETER k4 2.3e-3 +PARAMETER k6 2.2e8 +PARAMETER k7 4.0e5 +PARAMETER k8 1.0e8 +PARAMETER k9 3.0e7 +PARAMETER k10 55.0 +PARAMETER kDOC1 3.0E4 +PARAMETER kDOC2 6.5E5 +CONSTANT S1 0.02 +CONSTANT S2 0.50 + +[TERMS] + k5 (2.5e7*H) + (4.0e4*H2CO3) + (800*HCO3) + a1 k1*HOCL*NH3 + a2 k2*NH2CL + a3 k3*HOCL*NH2CL + a4 k4*NHCL2 + a5 k5*NH2CL*NH2CL + a6 k6*NHCL2*NH3*H + a7 k7*NHCL2*OH + a8 k8*I*NHCL2 + a9 k9*I*NH2CL + a10 k10*NH2CL*NHCL2 + a11 kDOC1*S1*TOC*NH2CL/12000 + a12 kDOC2*S2*TOC*HOCL/12000 + +[PIPE] +RATE HOCL -a1 + a2 - a3 + a4 + a8 - a12 +RATE NH3 -a1 + a2 + a5 - a6 + a11 +RATE NH2CL a1 - a2 - a3 + a4 - a5 + a6 - a9 - a10 - a11 +RATE NHCL2 a3 - a4 + a5 - a6 - a7 - a8 - a10 +RATE I a7 - a8 - a9 +RATE H 0 +RATE ALK 0 +RATE TOC 0 +EQUIL OCL H*OCL - 3.16E-8*HOCL +EQUIL NH4 H*NH3 - 5.01E-10*NH4 +EQUIL CO3 H*CO3 - 5.01E-11*HCO3 +EQUIL H2CO3 H*HCO3 - 5.01E-7*H2CO3 +EQUIL HCO3 ALK - HCO3 - 2*CO3 - OH + H +EQUIL OH H*OH - 1.0E-14 +FORMULA cNH3 17000*NH3 +FORMULA cNH2CL 51500*NH2CL + + +[SOURCES] +;CONC/MASS/FLOW/SETPOINT () +SETPOINT 10 TOC 2 +CONC 15 HOCL 0.8 +MASS 20 NH3 0.5 +FLOW 35 NH2CL 4.5 + +[QUALITY] +NODE 4 NH2CL 0.05E-3 +NODE 5 NH2CL 0.05E-3 +NODE 4 TOC 1.5 +NODE 5 TOC 3.0 +GLOBAL ALK 0.004 +GLOBAL H 2.818E-8 +GLOBAL OH 3.55E-7 + +[REPORT] +NODES ALL +LINKS ALL +SPECIE cNH3 YES 4 +SPECIE cNH2CL YES 4 + +[PATTERNS] +PAT1 190 +PAT2 20 + + + diff --git a/epyt/networks/msx-examples/net2-cl2.inp b/epyt/networks/msx-examples/net2-cl2.inp new file mode 100644 index 0000000..8d02ce4 --- /dev/null +++ b/epyt/networks/msx-examples/net2-cl2.inp @@ -0,0 +1,399 @@ +[TITLE] +EPANET Example Network 2 + +[JUNCTIONS] + 1 50 -694.4 2 ; + 2 100 8 ; + 3 60 14 ; + 4 60 8 ; + 5 100 8 ; + 6 125 5 ; + 7 160 4 ; + 8 110 9 ; + 9 180 14 ; + 10 130 5 ; + 11 185 34.78 ; + 12 210 16 ; + 13 210 2 ; + 14 200 2 ; + 15 190 2 ; + 16 150 20 ; + 17 180 20 ; + 18 100 20 ; + 19 150 5 ; + 20 170 19 ; + 21 150 16 ; + 22 200 10 ; + 23 230 8 ; + 24 190 11 ; + 25 230 6 ; + 27 130 8 ; + 28 110 0 ; + 29 110 7 ; + 30 130 3 ; + 31 190 17 ; + 32 110 17 ; + 33 180 1.5 ; + 34 190 1.5 ; + 35 110 0 ; + 36 110 1 ; + +[RESERVOIRS] + +[TANKS] + 26 235 56.7 50 70 50 0 ; + +[PIPES] + 1 1 2 2400 12 100 0 Open ; + 2 2 5 800 12 100 0 Open ; + 3 2 3 1300 8 100 0 Open ; + 4 3 4 1200 8 100 0 Open ; + 5 4 5 1000 12 100 0 Open ; + 6 5 6 1200 12 100 0 Open ; + 7 6 7 2700 12 100 0 Open ; + 8 7 8 1200 12 140 0 Open ; + 9 7 9 400 12 100 0 Open ; + 10 8 10 1000 8 140 0 Open ; + 11 9 11 700 12 100 0 Open ; + 12 11 12 1900 12 100 0 Open ; + 13 12 13 600 12 100 0 Open ; + 14 13 14 400 12 100 0 Open ; + 15 14 15 300 12 100 0 Open ; + 16 13 16 1500 8 100 0 Open ; + 17 15 17 1500 8 100 0 Open ; + 18 16 17 600 8 100 0 Open ; + 19 17 18 700 12 100 0 Open ; + 20 18 32 350 12 100 0 Open ; + 21 16 19 1400 8 100 0 Open ; + 22 14 20 1100 12 100 0 Open ; + 23 20 21 1300 8 100 0 Open ; + 24 21 22 1300 8 100 0 Open ; + 25 20 22 1300 8 100 0 Open ; + 26 24 23 600 12 100 0 Open ; + 27 15 24 250 12 100 0 Open ; + 28 23 25 300 12 100 0 Open ; + 29 25 26 200 12 100 0 Open ; + 30 25 31 600 12 100 0 Open ; + 31 31 27 400 8 100 0 Open ; + 32 27 29 400 8 100 0 Open ; + 34 29 28 700 8 100 0 Open ; + 35 22 33 1000 8 100 0 Open ; + 36 33 34 400 8 100 0 Open ; + 37 32 19 500 8 100 0 Open ; + 38 29 35 500 8 100 0 Open ; + 39 35 30 1000 8 100 0 Open ; + 40 28 35 700 8 100 0 Open ; + 41 28 36 300 8 100 0 Open ; + +[PUMPS] + +[VALVES] + +[DEMANDS] + +[STATUS] + +[PATTERNS] +; + 1 1.26 + 1 1.04 + 1 .97 + 1 .97 + 1 .89 + 1 1.19 + 1 1.28 + 1 .67 + 1 .67 + 1 1.34 + 1 2.46 + 1 .97 + 1 .92 + 1 .68 + 1 1.43 + 1 .61 + 1 .31 + 1 .78 + 1 .37 + 1 .67 + 1 1.26 + 1 1.56 + 1 1.19 + 1 1.26 + 1 .6 + 1 1.1 + 1 1.03 + 1 .73 + 1 .88 + 1 1.06 + 1 .99 + 1 1.72 + 1 1.12 + 1 1.34 + 1 1.12 + 1 .97 + 1 1.04 + 1 1.15 + 1 .91 + 1 .61 + 1 .68 + 1 .46 + 1 .51 + 1 .74 + 1 1.12 + 1 1.34 + 1 1.26 + 1 .97 + 1 .82 + 1 1.37 + 1 1.03 + 1 .81 + 1 .88 + 1 .81 + 1 .81 + +; + 2 .96 + 2 .96 + 2 .96 + 2 .96 + 2 .96 + 2 .96 + 2 .62 + 2 0 + 2 0 + 2 0 + 2 0 + 2 0 + 2 .8 + 2 1 + 2 1 + 2 1 + 2 1 + 2 .15 + 2 0 + 2 0 + 2 0 + 2 0 + 2 0 + 2 0 + 2 .55 + 2 .92 + 2 .92 + 2 .92 + 2 .92 + 2 .9 + 2 .9 + 2 .45 + 2 0 + 2 0 + 2 0 + 2 0 + 2 0 + 2 .7 + 2 1 + 2 1 + 2 1 + 2 1 + 2 .2 + 2 0 + 2 0 + 2 0 + 2 0 + 2 0 + 2 0 + 2 .74 + 2 .92 + 2 .92 + 2 .92 + 2 .92 + 2 .92 + +; + 3 .98 + 3 1.02 + 3 1.05 + 3 .99 + 3 .64 + 3 .46 + 3 .35 + 3 .35 + 3 .35 + 3 .35 + 3 .35 + 3 .35 + 3 .17 + 3 .17 + 3 .13 + 3 .13 + 3 .13 + 3 .15 + 3 .15 + 3 .15 + 3 .15 + 3 .15 + 3 .15 + 3 .15 + 3 .15 + 3 .12 + 3 .1 + 3 .08 + 3 .11 + 3 .09 + 3 .09 + 3 .08 + 3 .08 + 3 .08 + 3 .08 + 3 .08 + 3 .08 + 3 .09 + 3 .07 + 3 .07 + 3 .09 + 3 .09 + 3 .09 + 3 .09 + 3 .09 + 3 .09 + 3 .09 + 3 .09 + 3 .09 + 3 .08 + 3 .35 + 3 .72 + 3 .82 + 3 .92 + 3 1 + + +[CURVES] + +[CONTROLS] + +[RULES] + +[ENERGY] + Global Efficiency 85 + Global Price 0.0 + Demand Charge 0.0 + +[EMITTERS] + +[QUALITY] + 1 0.5 + 2 0.5 + 3 0.5 + 4 0.5 + 5 0.5 + 6 0.5 + 7 0.5 + 8 0.5 + 9 0.5 + 10 0.5 + 11 0.5 + 12 0.5 + 13 0.5 + 14 0.5 + 15 0.5 + 16 0.5 + 17 0.5 + 18 0.5 + 19 0.5 + 20 0.5 + 21 0.5 + 22 0.5 + 23 0.5 + 24 0.5 + 25 0.5 + 27 0.5 + 28 0.5 + 29 0.5 + 30 0.5 + 31 0.5 + 32 0.5 + 33 0.5 + 34 0.5 + 35 0.5 + 36 0.5 + 26 0.1 + +[SOURCES] + 1 CONCEN 0.8 + +[REACTIONS] + Order Bulk 1 + Order Wall 1 + Global Bulk -0.3 + Global Wall -1.0 + Limiting Potential 0.0 + Roughness Correlation 0.0 + +[MIXING] + +[TIMES] + Duration 55.00 hours + Hydraulic Timestep 1.00 hours + Quality Timestep 5.00 min + Pattern Timestep 1.00 hours + Report Timestep 1.00 hours + Report Start 0.00 hours + Start ClockTime 12 am + + +[OPTIONS] + Units GPM + Headloss H-W + Specific Gravity 1.0 + Viscosity 1.0 + Trials 40 + Accuracy 0.001 + Unbalanced Stop + Quality Chlorine mg/L + Diffusivity 1.0 + Segments 1000 + Tolerance 0.01 + +[COORDINATES] + 1 21.00 4.00 + 2 19.00 20.00 + 3 11.00 21.00 + 4 14.00 28.00 + 5 19.00 25.00 + 6 28.00 23.00 + 7 36.00 39.00 + 8 38.00 30.00 + 9 36.00 42.00 + 10 37.00 23.00 + 11 37.00 49.00 + 12 39.00 60.00 + 13 38.00 64.00 + 14 38.00 66.00 + 15 37.00 69.00 + 16 27.00 65.00 + 17 27.00 69.00 + 18 23.00 68.00 + 19 21.00 59.00 + 20 45.00 68.00 + 21 51.00 62.00 + 22 54.00 69.00 + 23 35.00 74.00 + 24 37.00 71.00 + 25 35.00 76.00 + 27 39.00 87.00 + 28 49.00 85.00 + 29 42.00 86.00 + 30 47.00 80.00 + 31 37.00 80.00 + 32 23.00 64.00 + 33 56.00 73.00 + 34 56.00 77.00 + 35 43.00 81.00 + 36 53.00 87.00 + 26 33.00 76.00 + +[LABELS] + 24.00 7.00 "Pump" 24.00 7.00 + 24.00 4.00 "Station" 24.00 4.00 + 23.00 77.00 "Tank" 23.00 77.00 + +[END] diff --git a/epyt/networks/msx-examples/net2-cl2.msx b/epyt/networks/msx-examples/net2-cl2.msx new file mode 100644 index 0000000..ce21d11 --- /dev/null +++ b/epyt/networks/msx-examples/net2-cl2.msx @@ -0,0 +1,49 @@ +[TITLE] +NET2 Chlorine Example + + +[OPTIONS] +AREA_UNITS FT2 +RATE_UNITS DAY +SOLVER EUL +TIMESTEP 300 +RTOL 0.001 +ATOL 0.001 + + +[SPECIES] +BULK CL2 MG 0.01 0.001 + + +[COEFFICIENTS] +PARAMETER Kb 0.3 +PARAMETER Kw 1.0 + + +[TERMS] +Kf 1.5826e-4 * RE^0.88 / D + +[PIPE] +RATE CL2 -Kb*CL2 - (4/D)*Kw*Kf/(Kw+Kf)*CL2 + + +[TANK] +RATE CL2 -Kb*CL2 + + +[SOURCES] +;CONC/MASS/FLOW/SETPOINT () +CONC 1 CL2 0.8 + + +[QUALITY] +GLOBAL CL2 0.5 +NODE 26 CL2 0.1 + +[PARAMETERS] +;PIPE +;TANK + +[REPORT] +NODES 2 20 23 26 +SPECIE CL2 YES diff --git a/epyt/tests/testMSXunit.py b/epyt/tests/testMSXunit.py new file mode 100644 index 0000000..c8a36e7 --- /dev/null +++ b/epyt/tests/testMSXunit.py @@ -0,0 +1,291 @@ +from epyt.epanet import epanet, epanetmsxapi +import numpy as np +import unittest +import os + + +class MSXtest(unittest.TestCase): + + def setUp(self): + """Call before every test case.""" + # Create EPANET object using the INP file + inpname = os.path.join(os.getcwd(), 'epyt', 'networks', 'Net3-NH2CL.inp') + self.epanetClass = epanet(inpname, msx=True) + file_path = os.path.join(os.getcwd(), 'epyt', 'Net3-NH2CL.msx') + self.msxClass = epanetmsxapi(file_path) + + #self.msxClass.MSXopen(file_path) + + def tearDown(self): + """Call after every test case.""" + self.msxClass.MSXclose() + self.epanetClass.unload() + + """ ------------------------------------------------------------------------- """ + + def test_MSXgetsource(self): + self.assertEqual(self.msxClass.MSXgetsource(1, 14), + ('SETPOINT', 2.0, 0), 'Wrong source comment output') + self.assertEqual(self.msxClass.MSXgetsource(2, 1), + ('CONCEN', 0.8, 0), 'Wrong source comment output') + self.assertEqual(self.msxClass.MSXgetsource(3, 2), + ('MASS', 0.5, 0), 'Wrong source comment output') + self.assertEqual(self.msxClass.MSXgetsource(4, 3), + ('FLOWPACED', 4.5, 0), 'Wrong source comment output') + + def test_MSXsetsource(self): + self.msxClass.MSXsetsource(1, 1, -1, 10.565, 1) + self.assertEqual(self.msxClass.MSXgetsource(1, 1), + ('NOSOURCE', 10.565, 1), 'Wrong source comment output') + self.msxClass.MSXsetsource(1, 1, 0, 10.565, 1) + self.assertEqual(self.msxClass.MSXgetsource(1, 1), + ('CONCEN', 10.565, 1), 'Wrong source comment output') + self.msxClass.MSXsetsource(1, 1, 1, 10.565, 1) + self.assertEqual(self.msxClass.MSXgetsource(1, 1), + ('MASS', 10.565, 1), 'Wrong source comment output') + self.msxClass.MSXsetsource(1, 1, 2, 10.565, 1) + self.assertEqual(self.msxClass.MSXgetsource(1, 1), + ('SETPOINT', 10.565, 1), 'Wrong source comment output') + self.msxClass.MSXsetsource(1, 1, 3, 10.565, 1) + self.assertEqual(self.msxClass.MSXgetsource(1, 1), + ('FLOWPACED', 10.565, 1), 'Wrong source comment output') + # set an integer value for level test_MSXsetsource + self.msxClass.MSXsetsource(1, 1, 3, 10, 1) + self.assertEqual(self.msxClass.MSXgetsource(1, 1), + ('FLOWPACED', 10, 1), 'Wrong source comment output') + def test_MSXgetspecies(self): + self.assertEqual(self.msxClass.MSXgetspecies(1), ('BULK', 'MMOL', 1.0e-8, 0.0001), + 'Wrong species comment output') + self.assertEqual(self.msxClass.MSXgetspecies(2), ('BULK', 'MMOL', 1.0e-8, 0.0001), + 'Wrong species comment output') + self.assertEqual(self.msxClass.MSXgetspecies(3), ('BULK', 'MMOL', 1.0e-8, 0.0001), + 'Wrong species comment output') + self.assertEqual(self.msxClass.MSXgetspecies(4), ('BULK', 'MMOL', 1.0e-8, 0.0001), + 'Wrong species comment output') + self.assertEqual(self.msxClass.MSXgetspecies(5), ('BULK', 'MMOL', 1.0e-8, 0.0001), + 'Wrong species comment output') + + + def test_MSXgetconstant(self): + + self.assertEqual(self.msxClass.MSXgetconstant(1), 0.02, + 'Wrong constant comment output') + self.assertEqual(self.msxClass.MSXgetconstant(2), 0.50, + 'Wrong constant comment output') + + + def test_MSXsetconstant(self): + + # set an integer value + self.msxClass.MSXsetconstant(1, 69) + self.assertEqual(self.msxClass.MSXgetconstant(1), 69, + 'Wrong set/get constant comment output') + # set a float value + self.msxClass.MSXsetconstant(1, 69.420) + self.assertEqual(self.msxClass.MSXgetconstant(1), 69.420, + 'Wrong set/get constant comment output') + + def testMSXgetinitqual(self): + self.assertEqual(self.msxClass.MSXgetinitqual(1, 1, 9), 2.818e-08, + 'Wrong get init qual comment output') + self.assertEqual(self.msxClass.MSXgetinitqual(1, 1, 10), 3.55e-7, + 'Wrong get init qual comment output') + self.assertEqual(self.msxClass.MSXgetinitqual(1, 1, 8), 0.004, + 'Wrong get init qual comment output') + self.assertEqual(self.msxClass.MSXgetinitqual(0, 1, 9), 2.818e-08, + 'Wrong get init qual comment output') + self.assertEqual(self.msxClass.MSXgetinitqual(0, 1, 10), 3.55e-7, + 'Wrong get init qual comment output') + self.assertEqual(self.msxClass.MSXgetinitqual(0, 1, 8), 0.004, + 'Wrong get init qual comment output') + + def testMSXsetinitqual(self): + # set value as integer testMSXsetinitqual + self.msxClass.MSXsetinitqual(1, 1, 1, 69) + self.assertEqual(self.msxClass.MSXgetinitqual(1, 1, 1), 69, + 'Wrong set/get init qual comment output') + # set value as float testMSXsetinitqual + self.msxClass.MSXsetinitqual(1, 1, 1, 69.420) + self.assertEqual(self.msxClass.MSXgetinitqual(1, 1, 1), 69.420, + 'Wrong set/get init qual comment output') + + def test_MSXsetinitqual(self): + #set value as integer test_MSXsetinitqual + self.msxClass.MSXsetinitqual(1,1,1,69) + self.assertEqual(self.msxClass.MSXgetinitqual(1, 1, 1), 69, + 'Wrong set/get init qual comment output') + #set value as float test_MSXsetinitqual + self.msxClass.MSXsetinitqual(1, 1, 1, 69.420) + self.assertEqual(self.msxClass.MSXgetinitqual(1, 1, 1), 69.420, + 'Wrong set/get init qual comment output') + + def test_MSXsetpatternvalue(self): + + #set value as integer test_MSXsetpatternvalue + self.msxClass.MSXsetpatternvalue(1,1,69) + self.assertEqual(self.msxClass.MSXgetpatternvalue(1, 1), 69, + 'Wrong set/get patternvalue comment output') + #set value as float test_MSXsetpatternvalue + self.msxClass.MSXsetpatternvalue(1, 1, 69.420) + self.assertEqual(self.msxClass.MSXgetpatternvalue(1, 1), 69.420, + 'Wrong set/get init patternvalue comment output') + + + def test_MSXgetIDlen(self): + + self.assertEqual(self.msxClass.MSXgetIDlen(3, 1), 4, + 'Wrong get ID len comment output') + self.assertEqual(self.msxClass.MSXgetIDlen(3, 2), 3, + 'Wrong get ID len comment output') + self.assertEqual(self.msxClass.MSXgetIDlen(3, 3), 5, + 'Wrong get ID len comment output') + + def test_MSXgetID(self): + + self.assertEqual(self.msxClass.MSXgetID(3, 1,4), 'HOCL', + 'Wrong get ID comment output') + self.assertEqual(self.msxClass.MSXgetID(3, 2,3), 'NH3', + 'Wrong get ID comment output') + self.assertEqual(self.msxClass.MSXgetID(3, 3,5), 'NH2CL', + 'Wrong get ID comment output') + + + def test_MSXgeterror(self): + + self.assertEqual(self.msxClass.MSXgeterror(516), 'Error 516 - reference made to an illegal object index.', + 'Wrong error comment output') + self.assertEqual(self.msxClass.MSXgeterror(505), 'Error 505 - could not read hydraulic results file.', + 'Wrong error comment output') + self.assertEqual(self.msxClass.MSXgeterror(503), 'Error 503 - could not open MSX input file.', + 'Wrong error comment output') + + + def test_MSXgetparameter(self): + self.assertEqual(self.msxClass.MSXgetparameter(1,1,2), 0.076, + 'Wrong error comment output') + self.assertEqual(self.msxClass.MSXgetparameter(1, 1,4), 2.3e-3, + 'Wrong error comment output') + + def test_MSXsetparameter(self): + #set value as integer test_MSXsetparameter + self.msxClass.MSXsetparameter(1,1,1,69) + self.assertEqual(self.msxClass.MSXgetparameter(1, 1, 1), 69, + 'Wrong error comment output') + #set value as float test_MSXsetparameter + self.msxClass.MSXsetparameter(1, 1, 1, 69.420) + self.assertEqual(self.msxClass.MSXgetparameter(1, 1, 1), 69.420, + 'Wrong error comment output') + + def test_MSXgetcount(self): + #for species + self.assertEqual(self.msxClass.MSXgetcount(3), 16, + 'Wrong error get count output') + #for parameters + self.assertEqual(self.msxClass.MSXgetcount(5), 11, + 'Wrong error get count output') + #for constants + self.assertEqual(self.msxClass.MSXgetcount(6), 2, + 'Wrong error get count output') + #for patterns + self.assertEqual(self.msxClass.MSXgetcount(7), 2, + 'Wrong error get count output') + + def test_MSXgetindex(self): + #testing 4 species (number 3) first , last, 1 char & more than 1 char and middle + self.assertEqual(self.msxClass.MSXgetindex(3,"HOCL"), 1, + 'Wrong error get count output') + self.assertEqual(self.msxClass.MSXgetindex(3,"cNH2CL"), 16, + 'Wrong error get count output') + self.assertEqual(self.msxClass.MSXgetindex(3,"H"), 9, + 'Wrong error get count output') + self.assertEqual(self.msxClass.MSXgetindex(3,"OH"), 10, + 'Wrong error get count output') + #testing parameters (number 5) first , last and one middle case + self.assertEqual(self.msxClass.MSXgetindex(5, "k1"), 1, + 'Wrong error get count output') + self.assertEqual(self.msxClass.MSXgetindex(5, "kDOC2"), 11, + 'Wrong error get count output') + self.assertEqual(self.msxClass.MSXgetindex(5, "k6"), 5, + 'Wrong error get count output') + #testing constants (number 6) + self.assertEqual(self.msxClass.MSXgetindex(6, "S1"), 1, + 'Wrong error get count output') + self.assertEqual(self.msxClass.MSXgetindex(6, "S2"), 2, + 'Wrong error get count output') + #testing paterns (number 7) + self.assertEqual(self.msxClass.MSXgetindex(7, "PAT1"), 1, + 'Wrong error get count output') + self.assertEqual(self.msxClass.MSXgetindex(7, "PAT2"), 2, + 'Wrong error get count output') + + def test_MSXaddpattern(self): + + y=self.msxClass.MSXgetcount(7) + self.msxClass.MSXaddpattern("Johnnys") + self.assertEqual(self.msxClass.MSXgetcount(7), y+1, + 'Wrong error add patter output') + + def test_MSXsetpatter(self): + self.msxClass.MSXaddpattern("JohnLegend") + x = self.msxClass.MSXgetindex(7, "JohnLegend") + mult = [0.5, 0.8, 1.2, 1.0, 0.7, 0.3] + self.msxClass.MSXsetpattern(x, mult, 6) + + self.assertEqual(self.msxClass.MSXgetpatternvalue(x, 1), 0.5, + 'Wrong set/get patternvalue comment output') + self.assertEqual(self.msxClass.MSXgetpatternvalue(x, 2), 0.8, + 'Wrong set/get patternvalue comment output') + self.assertEqual(self.msxClass.MSXgetpatternvalue(x, 3), 1.2, + 'Wrong set/get patternvalue comment output') + self.assertEqual(self.msxClass.MSXgetpatternvalue(x, 4), 1.0, + 'Wrong set/get patternvalue comment output') + self.assertEqual(self.msxClass.MSXgetpatternvalue(x, 5), 0.7, + 'Wrong set/get patternvalue comment output') + self.assertEqual(self.msxClass.MSXgetpatternvalue(x, 6), 0.3, + 'Wrong set/get patternvalue comment output') + + def test_MSXsavesmsxfile(self): + filename = "JohntheLegend.msx" + self.msxClass.MSXsavemsxfile(filename) + full_path = os.path.join(os.getcwd(), filename) + if os.path.exists(full_path): + print(f"The file {filename} exists in the current directory.") + else: + print(f"The file {filename} does not exist in the current directory.") + + + def test_MSXgetqual(self): + + self.msxClass.MSXclose() + self.epanetClass.unload() + inpname = os.path.join(os.getcwd(), 'epyt', 'networks', 'asce-tf-wdst', 'net2-cl2.inp') + self.epanetClass = epanet(inpname, msx=True) + file_path = os.path.join(os.getcwd(), 'epyt', 'net2-cl2.msx') + self.msxClass = epanetmsxapi(file_path) + + self.msxClass.MSXsolveH() + self.msxClass.MSXsolveQ() + t = 0 + tleft = 0 + self.msxClass.MSXinit(0) + c = 0 + while True: + t, tleft = self.msxClass.MSXstep() + c = c + 1 + if c == 1: + self.assertEqual(self.msxClass.MSXgetqual(0, 1,1), 0.8000000188349043, + 'Wrong get qual comment output') + if c == 85: + self.assertEqual(self.msxClass.MSXgetqual(0,1,1), 0.7991662288393907, + 'Wrong get qual comment output') + if c == 660: + self.assertEqual(self.msxClass.MSXgetqual(0,1,1), 0.7999999830262526, + 'Wrong get qual comment output') + if tleft <= 0: + break + + + + +if __name__ == "__main__": + unittest.main() # run all tests