diff --git a/src/Engine Updates.txt b/src/Engine Updates.txt index a7d7ad4..47c987b 100644 --- a/src/Engine Updates.txt +++ b/src/Engine Updates.txt @@ -1,76 +1,57 @@ ----------------------- -Build 2.00.11 (8/15/07) +Build 2.00.12 (2/25/08) ----------------------- Computational Engine Changes (epanet2.dll and epanet2d.exe): =============================================================================== CODE MODULES CHANGES =============================================================================== -All All variables previously declared as floats were re-defined as - doubles except for those used to write results to binary output - files and those used as arguments in the toolkit functions. +EPANET.C Temporary files are now written to the users's current working +INPUT1.C directory or to the designated Temporary Directory when running +VARS.H the Windows GUI, thus avoiding file access problems for users +FUNCS.H who do not have administrative privileges on their machine. ------------------------------------------------------------------------------- -EPANET.C The symbols CLE (for command line executable) and SOL (for shared - object library) were introduced in addition to the existing - symbol DLL (for dynamic linked library) to support conditional - compliation for shared object libraries on Unix/Linux. +EPANET.C The following tank parameters were made available for retrieval +TOOLKIT.H using the Toolkit's ENgetnodevalue function: tank diameter, +EPANET2.H minimum volume, index of the tank's volume curve, the initial, + minimum, and maximum water levels, the fraction of a 2-compart- + ment tank devoted to the mixing zone, and the tank's bulk + reaction rate coefficient. These same parameters (with the + exception of the volume curve) were made modifiable using the + ENsetnodevalue function, as was the choice of mixing model. + See the Toolkit Help file for details. ------------------------------------------------------------------------------- -EPANET.C EN_INITVOLUME, EN_MIXMODEL and EN_MIXZONEVOL were added to the -EPANET2.H list of tank parameters that can be retrieved with -TOOLKIT.H ENgetnodevalue. +EPANET.C A new Toolkit function, ENaddpattern, was added that allows +TOOLKIT.C programmers to add a new time pattern to the network in any +EPANET2.H custom code that they write. ------------------------------------------------------------------------------- -EPANET.C Missing code was added to the ENgetnodevalue function to return - a value for EN_SOURCEPAT when requested. +HYDRAUL.C A DAMPLIMIT parameter was added to control at what point during +INPUT1.C the hydraulic solution iterations status checks on PRVs and PSVs +INPUT3.C begin and subsequent flow changes are dampened. +REPORT.C +VARS.H +TEXT.H ------------------------------------------------------------------------------- -EPANET.H The function declarations in these "include" files were modified -TOOLKIT.H to support conditional compliation for shared object libraries - on Unix/Linux. +HYDRAUL.C Demands for tanks (net inflows/outflows) were not always being + set to zero when the tank was full or empty which was causing + problems in the water quality routing solutions. ------------------------------------------------------------------------------- -INPUT3.C The keyword "HEADLOSS" is no longer confused with "HEAD" when - parsing reporting variable names in the [REPORT] section of the - input file. +QUALITY.C The water quality functions for all of the tank mixing models + were modified so as to produce more robust results for cases + where there is a significant amount of volume exchange over a + water quality time step. ------------------------------------------------------------------------------- -RULES.C The input values for a tank's FILLTIME or DRAINTIME in a rule - premise are now correctly converted to seconds in the - newpremise() function. +EPANET.C A problem with the system not recognizing that water quality +QUALITY.C reactions could begin after some initial period of time when +VARS.H the Toolkit was used to modify reaction rate coefficients was + fixed. ------------------------------------------------------------------------------- -HYDRAUL.C During hydraulic balancing, status checks on control valves are - now made only when the convergence error drops to less than a - factor of 10 of the accuracy limit and a damping factor of 0.6 - is applied to flow changes computed at each iteration beyond - this point. +INPFILE.C A problem with extraneous lines being written to the [REPORT] + section of an EPANET input file created with the ENsaveinpfile + function was fixed. Also the number of decimal places for some + items written to the saved file, such as nodal demands, was + increased. ------------------------------------------------------------------------------- -HYDRAUL.C The logic for determining the status of PRVs and PSVs was changed - to produce more robust solutions, particularly for the case of - multi-stage regulator stations where there are two or more PRVs - connected in parallel. -------------------------------------------------------------------------------- -HYDRAUL.C The matrix coefficients for fully open control valves are now set - directly rather than by assuming the valve has a certain length - and friction factor. -------------------------------------------------------------------------------- -HYDRAUL.C Changes were made to how the P-coefficient for FCVs and the Y- - coefficient for GPVs are calculated. -------------------------------------------------------------------------------- -HYDRAUL.C An extraneous "if" statement was removed from the resistance() - function; the Y-coefficient value in gpvcoeff() was corrected; - the P-coefficient value in fcvcoeff() was corrected. -------------------------------------------------------------------------------- -QUALITY.C The memory pool used for water quality routing segments was given - the name SegPool and declared as a static global variable. -------------------------------------------------------------------------------- -QUALITY.C In the release() function, the upstream node quality is now mixed - together with that of the upstream pipe segment when the quality - difference between the two is less than the CTOL tolerance. -------------------------------------------------------------------------------- -INPFILE.C The ENsaveinpfile toolkit function now writes disabled reporting - variables to the [REPORT] section of the generated input file. -------------------------------------------------------------------------------- -INPFILE.C The appearance of an extraneous character at the end of the .INP - file produced by the toolkit function ENsaveinpfile was fixed. -------------------------------------------------------------------------------- -TYPES.H MAXID was increased to allow ID names to contain up to 31 -HASH.C characters. -------------------------------------------------------------------------------- -TYPES.H The code version was changed to 20011. +TYPES.H The code version was changed to 20012 and INT4 was explicitly + defined as an "int" data type. =============================================================================== diff --git a/src/Readme.txt b/src/Readme.txt index 2717685..ba91c01 100644 --- a/src/Readme.txt +++ b/src/Readme.txt @@ -48,7 +48,7 @@ Also included are the following header files: The comments at the top of each file lists the date when the latest update was made, and these updates can be located in the code by searching for comments with the phrase "/*** Updated" or with the -release number (e.g., 2.00.11) in them. +release number (e.g., 2.00.12) in them. Other useful documentation that can be consulted includes the EPANET Programmers Toolkit Help file and the EPANET Version 2 Users Manual. diff --git a/src/epanet.c b/src/epanet.c index ad56353..267a0de 100644 --- a/src/epanet.c +++ b/src/epanet.c @@ -11,6 +11,7 @@ DATE: 5/30/00 11/19/01 6/24/02 8/15/07 (2.00.11) + 2/14/08 (2.00.12) AUTHOR: L. Rossman US EPA - NRMRL @@ -114,18 +115,26 @@ execute function x and set the error code equal to its return value. /*** Following lines are deprecated ***/ //(2.00.11 - LR) //#ifdef DLL //#include - -/*** Updated 9/7/00 ***/ -#include - +//#include //#endif +/*** Need to define WINDOWS to use the getTmpName function ***/ //(2.00.12 - LR) +// --- define WINDOWS +#undef WINDOWS +#ifdef _WIN32 + #define WINDOWS +#endif +#ifdef __WIN32__ + #define WINDOWS +#endif +/************************************************************/ #include #include #include #include #include +#include //(2.00.12 - LR) #include "hash.h" #include "text.h" #include "types.h" @@ -135,7 +144,6 @@ execute function x and set the error code equal to its return value. #include "vars.h" #include "toolkit.h" - void (* viewprog) (char *); /* Pointer to progress viewing function */ @@ -349,10 +357,21 @@ int DLLEXPORT ENclose() { if (Openflag) writetime(FMT105); freedata(); + + if (TmpOutFile != OutFile) //(2.00.12 - LR) + { //(2.00.12 - LR) + if (TmpOutFile != NULL) fclose(TmpOutFile); //(2.00.12 - LR) + remove(TmpFname); //(2.00.12 - LR) + } //(2.00.12 - LR) + if (InFile != NULL) fclose(InFile); if (RptFile != NULL) fclose(RptFile); if (HydFile != NULL) fclose(HydFile); if (OutFile != NULL) fclose(OutFile); + + if (Hydflag == SCRATCH) remove(HydFname); //(2.00.12 - LR) + if (Outflag == SCRATCH) remove(OutFname); //(2.00.12 - LR) + Openflag = FALSE; OpenHflag = FALSE; SaveHflag = FALSE; @@ -1403,6 +1422,59 @@ int DLLEXPORT ENgetnodevalue(int index, int code, float *value) v = C[index]*Ucf[QUALITY]; break; +/*** New parameters added for retrieval begins here ***/ //(2.00.12 - LR) +/*** (Thanks to Nicolas Basile of Ecole Polytechnique ***/ +/*** de Montreal for suggesting some of these.) ***/ + + case EN_TANKDIAM: + v = 0.0; + if ( index > Njuncs ) + { + v = 4.0/PI*sqrt(Tank[index-Njuncs].A)*Ucf[ELEV]; + } + break; + + case EN_MINVOLUME: + v = 0.0; + if ( index > Njuncs ) v = Tank[index-Njuncs].Vmin * Ucf[VOLUME]; + break; + + case EN_VOLCURVE: + v = 0.0; + if ( index > Njuncs ) v = Tank[index-Njuncs].Vcurve; + break; + + case EN_MINLEVEL: + v = 0.0; + if ( index > Njuncs ) + { + v = (Tank[index-Njuncs].Hmin - Node[index].El) * Ucf[ELEV]; + } + break; + + case EN_MAXLEVEL: + v = 0.0; + if ( index > Njuncs ) + { + v = (Tank[index-Njuncs].Hmax - Node[index].El) * Ucf[ELEV]; + } + break; + + case EN_MIXFRACTION: + v = 1.0; + if ( index > Njuncs && Tank[index-Njuncs].Vmax > 0.0) + { + v = Tank[index-Njuncs].V1max / Tank[index-Njuncs].Vmax; + } + break; + + case EN_TANK_KBULK: + v = 0.0; + if (index > Njuncs) v = Tank[index-Njuncs].Kb * SECperDAY; + break; + +/*** New parameter additions ends here. ***/ //(2.00.12 - LR) + default: return(251); } *value = (float)v; @@ -1849,6 +1921,85 @@ int DLLEXPORT ENsetnodevalue(int index, int code, float v) } break; +/*** New parameters added for retrieval begins here ***/ //(2.00.12 - LR) +/*** (Thanks to Nicolas Basile of Ecole Polytechnique ***/ +/*** de Montreal for suggesting some of these.) ***/ + + case EN_TANKDIAM: + if (value <= 0.0) return(202); + j = index - Njuncs; + if (j > 0 && Tank[j].A > 0.0) + { + value /= Ucf[ELEV]; + Tank[j].A = PI*SQR(value)/4.0; + Tank[j].Vmin = tankvolume(j, Tank[j].Hmin); + Tank[j].V0 = tankvolume(j, Tank[j].H0); + Tank[j].Vmax = tankvolume(j, Tank[j].Hmax); + } + break; + + case EN_MINVOLUME: + if (value < 0.0) return(202); + j = index - Njuncs; + if (j > 0 && Tank[j].A > 0.0) + { + Tank[j].Vmin = value/Ucf[VOLUME]; + Tank[j].V0 = tankvolume(j, Tank[j].H0); + Tank[j].Vmax = tankvolume(j, Tank[j].Hmax); + } + break; + + case EN_MINLEVEL: + if (value < 0.0) return(202); + j = index - Njuncs; + if (j > 0 && Tank[j].A > 0.0) + { + if (Tank[j].Vcurve > 0) return(202); + Tank[j].Hmin = value/Ucf[ELEV] + Node[index].El; + Tank[j].Vmin = tankvolume(j, Tank[j].Hmin); + } + break; + + case EN_MAXLEVEL: + if (value < 0.0) return(202); + j = index - Njuncs; + if (j > 0 && Tank[j].A > 0.0) + { + if (Tank[j].Vcurve > 0) return(202); + Tank[j].Hmax = value/Ucf[ELEV] + Node[index].El; + Tank[j].Vmax = tankvolume(j, Tank[j].Hmax); + } + break; + + case EN_MIXMODEL: + j = ROUND(value); + if (j < MIX1 || j > LIFO) return(202); + if (index > Njuncs && Tank[index-Njuncs].A > 0.0) + { + Tank[index-Njuncs].MixModel = (char)j; + } + break; + + case EN_MIXFRACTION: + if (value < 0.0 || value > 1.0) return(202); + j = index - Njuncs; + if (j > 0 && Tank[j].A > 0.0) + { + Tank[j].V1max = value*Tank[j].Vmax; + } + break; + + case EN_TANK_KBULK: + j = index - Njuncs; + if (j > 0 && Tank[j].A > 0.0) + { + Tank[j].Kb = value/SECperDAY; + Reactflag = 1; + } + break; + +/*** New parameter additions ends here. ***/ //(2.00.12 - LR) + default: return(251); } return(0); @@ -1953,11 +2104,19 @@ int DLLEXPORT ENsetlinkvalue(int index, int code, float v) break; case EN_KBULK: - if (Link[index].Type <= PIPE) Link[index].Kb = value/SECperDAY; + if (Link[index].Type <= PIPE) + { + Link[index].Kb = value/SECperDAY; + Reactflag = 1; //(2.00.12 - LR) + } break; case EN_KWALL: - if (Link[index].Type <= PIPE) Link[index].Kw = value/SECperDAY; + if (Link[index].Type <= PIPE) + { + Link[index].Kw = value/SECperDAY; + Reactflag = 1; //(2.00.12 - LR) + } break; default: return(251); @@ -1965,6 +2124,74 @@ int DLLEXPORT ENsetlinkvalue(int index, int code, float v) return(0); } + +int DLLEXPORT ENaddpattern(char *id) +/*---------------------------------------------------------------- +** Input: id = ID name of the new pattern +** Output: none +** Returns: error code +** Purpose: adds a new time pattern appended to the end of the +** existing patterns. +**---------------------------------------------------------------- +*/ +{ + int i, j, n, err = 0; + Spattern *tmpPat; + +/* Check if a pattern with same id already exists */ + + if ( !Openflag ) return(102); + if ( ENgetpatternindex(id, &i) == 0 ) return(215); + +/* Check that id name is not too long */ + + if (strlen(id) > MAXID) return(250); + +/* Allocate memory for a new array of patterns */ + + n = Npats + 1; + tmpPat = (Spattern *) calloc(n+1, sizeof(Spattern)); + if ( tmpPat == NULL ) return(101); + +/* Copy contents of old pattern array to new one */ + + for (i=0; i<=Npats; i++) + { + strcpy(tmpPat[i].ID, Pattern[i].ID); + tmpPat[i].Length = Pattern[i].Length; + tmpPat[i].F = (double *) calloc(Pattern[i].Length, sizeof(double)); + if (tmpPat[i].F == NULL) err = 1; + else for (j=0; j 0) Outflag = SAVE; //(2.00.12 - LR) + else Outflag = SCRATCH; //(2.00.12 - LR) /* Check that file names are not identical */ if (strcomp(f1,f2) || strcomp(f1,f3) || strcomp(f2,f3)) @@ -2271,7 +2500,8 @@ int openhydfile() HydFile = NULL; switch(Hydflag) { - case SCRATCH: HydFile = tmpfile(); + case SCRATCH: getTmpName(HydFname); //(2.00.12 - LR) + HydFile = fopen(HydFname, "w+b"); //(2.00.12 - LR) break; case SAVE: HydFile = fopen(HydFname,"w+b"); break; @@ -2337,9 +2567,13 @@ int openoutfile() if (TmpOutFile != NULL) fclose(TmpOutFile); TmpOutFile = NULL; + if (Outflag == SCRATCH) remove(OutFname); //(2.00.12 - LR) + remove(TmpFname); //(2.00.12 - LR) + /* If output file name was supplied, then attempt to */ /* open it. Otherwise open a temporary output file. */ - if (strlen(OutFname) != 0) + //if (strlen(OutFname) != 0) //(2.00.12 - LR) + if (Outflag == SAVE) //(2.00.12 - LR) { if ( (OutFile = fopen(OutFname,"w+b")) == NULL) { @@ -2347,10 +2581,15 @@ int openoutfile() errcode = 304; } } - else if ( (OutFile = tmpfile()) == NULL) + //else if ( (OutFile = tmpfile()) == NULL) //(2.00.12 - LR) + else //(2.00.12 - LR) { - writecon(FMT08); - errcode = 304; + getTmpName(OutFname); //(2.00.12 - LR) + if ( (OutFile = fopen(OutFname,"w+b")) == NULL) //(2.00.12 - LR) + { + writecon(FMT08); + errcode = 304; + } } /* Save basic network data & energy usage results */ @@ -2364,7 +2603,10 @@ int openoutfile() { if (Tstatflag != SERIES) { - if ( (TmpOutFile = tmpfile()) == NULL) errcode = 304; + //if ( (TmpOutFile = tmpfile()) == NULL) errcode = 304; //(2.00.12 - LR) + getTmpName(TmpFname); //(2.00.12 - LR) + TmpOutFile = fopen(TmpFname, "w+b"); //(2.00.12 - LR) + if (TmpOutFile == NULL) errcode = 304; //(2.00.12 - LR) } else TmpOutFile = OutFile; } @@ -2638,6 +2880,50 @@ void freedata() ---------------------------------------------------------------- */ +/*** New function for 2.00.12 ***/ //(2.00.12 - LR) +char* getTmpName(char* fname) +// +// Input: fname = file name string +// Output: returns pointer to file name +// Purpose: creates a temporary file name with path prepended to it. +// +{ + char name[MAXFNAME+1]; + int n; + + // --- for Windows systems: + #ifdef WINDOWS + // --- use system function tmpnam() to create a temporary file name + tmpnam(name); + + // --- if user supplied the name of a temporary directory, + // then make it be the prefix of the full file name + n = strlen(TmpDir); + if ( n > 0 ) + { + strcpy(fname, TmpDir); + if ( fname[n-1] != '\\' ) strcat(fname, "\\"); + } + + // --- otherwise, use the relative path notation as the file name + // prefix so that the file will be placed in the current directory + else + { + strcpy(fname, ".\\"); + } + + // --- now add the prefix to the file name + strcat(fname, name); + + // --- for non-Windows systems: + #else + // --- use system function mkstemp() to create a temporary file name + strcpy(fname, "enXXXXXX"); + mkstemp(fname); + #endif + return fname; +} + int strcomp(char *s1, char *s2) /*--------------------------------------------------------------- diff --git a/src/epanet2.def b/src/epanet2.def index 797fc34..7ee86fa 100644 --- a/src/epanet2.def +++ b/src/epanet2.def @@ -1,57 +1,58 @@ LIBRARY EPANET2.DLL EXPORTS - ENclose @4 ; ENclose - ENcloseH @11 ; ENcloseH - ENcloseQ @20 ; ENcloseQ - ENepanet @1 ; ENepanet - ENgetcontrol @26 ; ENgetcontrol - ENgetcount @27 ; ENgetcount - ENgeterror @36 ; ENgeterror - ENgetflowunits @30 ; ENgetflowunits - ENgetlinkid @42 ; ENgetlinkid - ENgetlinkindex @41 ; ENgetlinkindex - ENgetlinknodes @44 ; ENgetlinknodes - ENgetlinktype @43 ; ENgetlinktype - ENgetlinkvalue @45 ; ENgetlinkvalue - ENgetnodeid @38 ; ENgetnodeid - ENgetnodeindex @37 ; ENgetnodeindex - ENgetnodetype @39 ; ENgetnodetype - ENgetnodevalue @40 ; ENgetnodevalue - ENgetoption @28 ; ENgetoption - ENgetpatternid @32 ; ENgetpatternid - ENgetpatternindex @31 ; ENgetpatternindex - ENgetpatternlen @33 ; ENgetpatternlen - ENgetpatternvalue @34 ; ENgetpatternvalue - ENgetqualtype @35 ; ENgetqualtype - ENgettimeparam @29 ; ENgettimeparam - ENgetversion @25 ; ENgetversion - ENinitH @8 ; ENinitH - ENinitQ @16 ; ENinitQ - ENnextH @10 ; ENnextH - ENnextQ @18 ; ENnextQ - ENopen @2 ; ENopen - ENopenH @7 ; ENopenH - ENopenQ @15 ; ENopenQ - ENreport @22 ; ENreport - ENresetreport @23 ; ENresetreport - ENrunH @9 ; ENrunH - ENrunQ @17 ; ENrunQ - ENsaveH @6 ; ENsaveH - ENsavehydfile @12 ; ENsavehydfile - ENsaveinpfile @3 ; ENsaveinpfile - ENsetcontrol @46 ; ENsetcontrol - ENsetlinkvalue @48 ; ENsetlinkvalue - ENsetnodevalue @47 ; ENsetnodevalue - ENsetoption @52 ; ENsetoption - ENsetpattern @49 ; ENsetpattern - ENsetpatternvalue @50 ; ENsetpatternvalue - ENsetqualtype @54 ; ENsetqualtype - ENsetreport @24 ; ENsetreport - ENsetstatusreport @53 ; ENsetstatusreport - ENsettimeparam @51 ; ENsettimeparam - ENsolveH @5 ; ENsolveH - ENsolveQ @14 ; ENsolveQ - ENstepQ @19 ; ENstepQ - ENusehydfile @13 ; ENusehydfile - ENwriteline @21 ; ENwriteline + ENaddpattern = _ENaddpattern@4 + ENclose = _ENclose@0 + ENcloseH = _ENcloseH@0 + ENcloseQ = _ENcloseQ@0 + ENepanet = _ENepanet@16 + ENgetcontrol = _ENgetcontrol@24 + ENgetcount = _ENgetcount@8 + ENgeterror = _ENgeterror@12 + ENgetflowunits = _ENgetflowunits@4 + ENgetlinkid = _ENgetlinkid@8 + ENgetlinkindex = _ENgetlinkindex@8 + ENgetlinknodes = _ENgetlinknodes@12 + ENgetlinktype = _ENgetlinktype@8 + ENgetlinkvalue = _ENgetlinkvalue@12 + ENgetnodeid = _ENgetnodeid@8 + ENgetnodeindex = _ENgetnodeindex@8 + ENgetnodetype = _ENgetnodetype@8 + ENgetnodevalue = _ENgetnodevalue@12 + ENgetoption = _ENgetoption@8 + ENgetpatternid = _ENgetpatternid@8 + ENgetpatternindex = _ENgetpatternindex@8 + ENgetpatternlen = _ENgetpatternlen@8 + ENgetpatternvalue = _ENgetpatternvalue@12 + ENgetqualtype = _ENgetqualtype@8 + ENgettimeparam = _ENgettimeparam@8 + ENgetversion = _ENgetversion@4 + ENinitH = _ENinitH@4 + ENinitQ = _ENinitQ@4 + ENnextH = _ENnextH@4 + ENnextQ = _ENnextQ@4 + ENopen = _ENopen@12 + ENopenH = _ENopenH@0 + ENopenQ = _ENopenQ@0 + ENreport = _ENreport@0 + ENresetreport = _ENresetreport@0 + ENrunH = _ENrunH@4 + ENrunQ = _ENrunQ@4 + ENsaveH = _ENsaveH@0 + ENsavehydfile = _ENsavehydfile@4 + ENsaveinpfile = _ENsaveinpfile@4 + ENsetcontrol = _ENsetcontrol@24 + ENsetlinkvalue = _ENsetlinkvalue@12 + ENsetnodevalue = _ENsetnodevalue@12 + ENsetoption = _ENsetoption@8 + ENsetpattern = _ENsetpattern@12 + ENsetpatternvalue = _ENsetpatternvalue@12 + ENsetqualtype = _ENsetqualtype@16 + ENsetreport = _ENsetreport@4 + ENsetstatusreport = _ENsetstatusreport@4 + ENsettimeparam = _ENsettimeparam@8 + ENsolveH = _ENsolveH@0 + ENsolveQ = _ENsolveQ@0 + ENstepQ = _ENstepQ@4 + ENusehydfile = _ENusehydfile@4 + ENwriteline = _ENwriteline@4 diff --git a/src/epanet2.h b/src/epanet2.h index 2abefd7..8d3447b 100644 --- a/src/epanet2.h +++ b/src/epanet2.h @@ -3,7 +3,7 @@ ** ** C/C++ header file for EPANET Programmers Toolkit ** -** Last updated on 8/15/07 (2.00.11) +** Last updated on 2/14/08 (2.00.12) */ #ifndef EPANET2_H @@ -29,6 +29,14 @@ #define EN_MIXMODEL 15 #define EN_MIXZONEVOL 16 +#define EN_TANKDIAM 17 +#define EN_MINVOLUME 18 +#define EN_VOLCURVE 19 +#define EN_MINLEVEL 20 +#define EN_MAXLEVEL 21 +#define EN_MIXFRACTION 22 +#define EN_TANK_KBULK 23 + #define EN_DIAMETER 0 /* Link parameters */ #define EN_LENGTH 1 #define EN_ROUGHNESS 2 @@ -113,6 +121,11 @@ #define EN_MAXIMUM 3 #define EN_RANGE 4 +#define EN_MIX1 0 /* Tank mixing models */ +#define EN_MIX2 1 +#define EN_FIFO 2 +#define EN_LIFO 3 + #define EN_NOSAVE 0 /* Save-results-to-file flag */ #define EN_SAVE 1 #define EN_INITFLOW 10 /* Re-initialize flow flag */ @@ -205,6 +218,7 @@ int DLLEXPORT ENsetcontrol(int, int, int, float, int, float); int DLLEXPORT ENsetnodevalue(int, int, float); int DLLEXPORT ENsetlinkvalue(int, int, float); + int DLLEXPORT ENaddpattern(char *); int DLLEXPORT ENsetpattern(int, float *, int); int DLLEXPORT ENsetpatternvalue(int, int, float); int DLLEXPORT ENsettimeparam(int, long); diff --git a/src/funcs.h b/src/funcs.h index 20881c8..31ac700 100644 --- a/src/funcs.h +++ b/src/funcs.h @@ -9,6 +9,7 @@ DATE: 5/8/00 10/25/00 12/29/00 3/1/01 + 2/14/08 (2.00.12) AUTHOR: L. Rossman US EPA - NRMRL @@ -33,7 +34,8 @@ int openfiles(char *,char *,char *); /* Opens input & report files */ int openhydfile(void); /* Opens hydraulics file */ int openoutfile(void); /* Opens binary output file */ int strcomp(char *, char *); /* Compares two strings */ -double interp(int, double *, /* Interpolates a data curve */ +char* getTmpName(char* fname); /* Gets temporary file name */ //(2.00.12 - LR) +double interp(int, double *, /* Interpolates a data curve */ double *, double); int findnode(char *); /* Finds node's index from ID */ int findlink(char *); /* Finds link's index from ID */ diff --git a/src/hydraul.c b/src/hydraul.c index dc55247..722b40c 100644 --- a/src/hydraul.c +++ b/src/hydraul.c @@ -12,6 +12,7 @@ DATE: 6/5/00 11/19/01 6/24/02 8/15/07 (2.00.11) + 2/14/08 (2.00.12) AUTHOR: L. Rossman US EPA - NRMRL @@ -1114,12 +1115,14 @@ int netsolve(int *iter, double *relerr) ** using Todini's Gradient algorithm ** *** Updated 9/7/00 *** -*** Updated 2.00.11 *** +*** Updated 2.00.11 *** +*** Updated 2.00.12 *** ** Notes: Status checks on CVs, pumps and pipes to tanks are made ** every CheckFreq iteration, up until MaxCheck iterations ** are reached. Status checks on control valves are made -** when within a factor of 10 of convergence. Also at this -** point, future computed flow changes are only 60% of +** every iteration if DampLimit = 0 or only when the +** convergence error is at or below DampLimit. If DampLimit +** is > 0 then future computed flow changes are only 60% of ** their full value. A complete status check on all links ** is made when convergence is achieved. If convergence is ** not achieved in MaxIter trials and ExtraIter > 0 then @@ -1135,7 +1138,6 @@ int netsolve(int *iter, double *relerr) int nextcheck; /* Next status check trial */ int maxtrials; /* Max. trials for convergence */ double newerr; /* New convergence error */ - double hacc10 = 10.0*Hacc; /* Relaxed convergence criterion */ int valveChange; /* Valve status change flag */ int statChange; @@ -1179,18 +1181,18 @@ int netsolve(int *iter, double *relerr) /* Write convergence error to status report if called for */ if (Statflag == FULL) writerelerr(*iter,*relerr); - /* Apply relaxation & check valve status if close to convergence */ - - if ( *relerr <= hacc10) + /* Apply solution damping & check for change in valve status */ + RelaxFactor = 1.0; + valveChange = FALSE; + if ( DampLimit > 0.0 ) { - RelaxFactor = 0.6; - valveChange = valvestatus(); - } - else - { - RelaxFactor = 1.0; - valveChange = FALSE; + if( *relerr <= DampLimit ) + { + RelaxFactor = 0.6; + valveChange = valvestatus(); + } } + else valveChange = valvestatus(); /* Check for convergence */ if (*relerr <= Hacc) @@ -1283,7 +1285,7 @@ int valvestatus() ** Input: none ** Output: returns 1 if any pressure or flow control valve //(2.00.11 - LR) ** changes status, 0 otherwise //(2.00.11 - LR) -** Purpose: updates status for PRVs, PSVs & FCVs whose status //(2.00.11 - LR) +** Purpose: updates status for PRVs & PSVs whose status //(2.00.12 - LR) ** is not fixed to OPEN/CLOSED **----------------------------------------------------------------- */ @@ -1314,8 +1316,9 @@ int valvestatus() S[k] = psvstatus(k,s,hset,H[n1],H[n2]); break; - case FCV: S[k] = fcvstatus(k,s,H[n1],H[n2]); //(2.00.11 - LR) - break; //(2.00.11 - LR) +//// FCV status checks moved back into the linkstatus() function //// //(2.00.12 - LR) +// case FCV: S[k] = fcvstatus(k,s,H[n1],H[n2]); //(2.00.12 - LR) +// break; //(2.00.12 - LR) default: continue; } @@ -1341,7 +1344,7 @@ int linkstatus() **-------------------------------------------------------------- ** Input: none ** Output: returns 1 if any link changes status, 0 otherwise -** Purpose: determines new status for pumps, CVs, & pipes //(2.00.11 - LR) +** Purpose: determines new status for pumps, CVs, FCVs & pipes //(2.00.12 - LR) ** to tanks. **-------------------------------------------------------------- */ @@ -1370,8 +1373,8 @@ int linkstatus() S[k] = pumpstatus(k,-dh); /* Check for status changes in non-fixed FCVs */ - //if (Link[k].Type == FCV && K[k] != MISSING) //(2.00.11 - LR) - // S[k] = fcvstatus(k,status,H[n1],H[n2]); //(2.00.11 - LR) + if (Link[k].Type == FCV && K[k] != MISSING) //(2.00.12 - LR)// + S[k] = fcvstatus(k,status,H[n1],H[n2]); //(2.00.12 - LR)// /* Check for flow into (out of) full (empty) tanks */ if (n1 > Njuncs || n2 > Njuncs) tankstatus(k,n1,n2); @@ -1381,7 +1384,7 @@ int linkstatus() { change = TRUE; if (Statflag == FULL) writestatchange(k,status,S[k]); - + //if (S[k] <= CLOSED) Q[k] = QZERO; //(2.00.11 - LR) //else setlinkflow(k, dh); //(2.00.11 - LR) } @@ -1763,8 +1766,11 @@ double newflows() dqsum += ABS(dq); /* Update net flows to tanks */ - if (n1 > Njuncs) D[n1] -= Q[k]; - if (n2 > Njuncs) D[n2] += Q[k]; + if ( S[k] > CLOSED ) //(2.00.12 - LR) + { + if (n1 > Njuncs) D[n1] -= Q[k]; + if (n2 > Njuncs) D[n2] += Q[k]; + } } diff --git a/src/inpfile.c b/src/inpfile.c index 7692d11..204d61f 100644 --- a/src/inpfile.c +++ b/src/inpfile.c @@ -9,6 +9,7 @@ DATE: 5/8/00 11/19/01 6/24/02 8/15/07 (2.00.11) + 2/14/08 (2.00.12) AUTHOR: L. Rossman US EPA - NRMRL @@ -43,8 +44,12 @@ extern char *RptFlagTxt[]; extern char *SectTxt[]; -/*** Updated 6/24/02 ***/ -void savemapinfo(FILE *ftmp) +void saveauxdata(FILE *f) //(2.00.12 - LR) +/* +------------------------------------------------------------ + Writes auxilary data from original input file to new file. +------------------------------------------------------------ +*/ { int sect,newsect; char *tok; @@ -55,14 +60,10 @@ void savemapinfo(FILE *ftmp) rewind(InFile); while (fgets(line,MAXLINE,InFile) != NULL) { - /* Skip blank lines & those beginning with a comment */ - strcpy(s,line); - tok = strtok(line,SEPSTR); - if (tok == NULL) continue; - if (*tok == ';') continue; - /* Check if line begins with a new section heading */ - if (*tok == '[') + strcpy(s,line); + tok = strtok(s,SEPSTR); + if (tok != NULL && *tok == '[') { newsect = findmatch(tok,SectTxt); if (newsect >= 0) @@ -76,7 +77,7 @@ void savemapinfo(FILE *ftmp) case _VERTICES: case _LABELS: case _BACKDROP: - case _TAGS: fwrite(s, sizeof(char), strlen(s), ftmp); + case _TAGS: fprintf(f, "%s", line); //(2.00.12 - LR) } continue; } @@ -91,13 +92,14 @@ void savemapinfo(FILE *ftmp) case _VERTICES: case _LABELS: case _BACKDROP: - case _TAGS: fwrite(s, sizeof(char), strlen(s), ftmp); + case _TAGS: fprintf(f, "%s", line); //(2.00.12 - LR) } } } -/*** End of update ***/ +//// This function was heavily modified. //// //(2.00.12 - LR) + int saveinpfile(char *fname) /* ------------------------------------------------- @@ -105,38 +107,16 @@ int saveinpfile(char *fname) ------------------------------------------------- */ { - int i,j,n; - double d,kc,ke,km,ucf; - -/*** Updated 6/24/02 ***/ - char s[MAXLINE+1], s1[MAXLINE+1], s2[MAXLINE+1]; - + int i,j,n; + double d,kc,ke,km,ucf; + char s[MAXLINE+1], s1[MAXLINE+1], s2[MAXLINE+1]; Pdemand demand; Psource source; - FILE *f; + FILE *f; -/*** Updated 6/24/02 ***/ - FILE *ftmp; +/* Open the new text file */ -/* Copy [RULES], [COORDS], [VERTICES], [LABELS], [BACKDROP] & [TAGS] */ -/* sections from original input file to new input file */ - - ftmp = NULL; - if (InFile) - { - ftmp = tmpfile(); - if (ftmp) savemapinfo(ftmp); - } - -/* Open text file */ - - if ((f = fopen(fname,"wt")) == NULL) - { - if (ftmp) fclose(ftmp); - return(102); - } - -/*** End of update ***/ + if ((f = fopen(fname,"wt")) == NULL) return(308); /* Write [TITLE] section */ @@ -151,7 +131,7 @@ int saveinpfile(char *fname) fprintf(f,"\n\n[JUNCTIONS]"); for (i=1; i<=Njuncs; i++) - fprintf(f,"\n %-15s %12.2f", Node[i].ID, Node[i].El*Ucf[ELEV]); + fprintf(f,"\n %-31s %12.4f", Node[i].ID, Node[i].El*Ucf[ELEV]); /* Write [RESERVOIRS] section */ @@ -161,9 +141,9 @@ int saveinpfile(char *fname) if (Tank[i].A == 0.0) { n = Tank[i].Node; - sprintf(s," %-15s %12.2f",Node[n].ID, Node[n].El*Ucf[ELEV]); + sprintf(s," %-31s %12.4f",Node[n].ID, Node[n].El*Ucf[ELEV]); if ((j = Tank[i].Pat) > 0) - sprintf(s1," %-15s",Pattern[j].ID); + sprintf(s1," %-31s",Pattern[j].ID); else strcpy(s1,""); fprintf(f, "\n%s %s", s,s1); @@ -178,7 +158,7 @@ int saveinpfile(char *fname) if (Tank[i].A > 0.0) { n = Tank[i].Node; - sprintf(s," %-15s %12.2f %12.2f %12.2f %12.2f %12.2f %12.2f", + sprintf(s," %-31s %12.4f %12.4f %12.4f %12.4f %12.4f %12.4f", Node[n].ID, Node[n].El*Ucf[ELEV], (Tank[i].H0 - Node[n].El)*Ucf[ELEV], @@ -187,7 +167,7 @@ int saveinpfile(char *fname) sqrt(4.0*Tank[i].A/PI)*Ucf[ELEV], Tank[i].Vmin*SQR(Ucf[ELEV])*Ucf[ELEV]); if ((j = Tank[i].Vcurve) > 0) - sprintf(s1,"%-15s",Curve[j].ID); + sprintf(s1,"%-31s",Curve[j].ID); else strcpy(s1,""); fprintf(f, "\n%s %s", s,s1); @@ -205,21 +185,18 @@ int saveinpfile(char *fname) kc = Link[i].Kc; if (Formflag == DW) kc = kc*Ucf[ELEV]*1000.0; km = Link[i].Km*SQR(d)*SQR(d)/0.02517; - -/*** Updated 6/24/02 ***/ - sprintf(s," %-15s %-15s %-15s %12.2f %12.2f", + sprintf(s," %-31s %-31s %-31s %12.4f %12.4f", Link[i].ID, Node[Link[i].N1].ID, Node[Link[i].N2].ID, Link[i].Len*Ucf[LENGTH], d*Ucf[DIAM]); if (Formflag == DW) sprintf(s1, "%12.4f %12.4f", kc, km); - else sprintf(s1, "%12.2f %12.4f", kc, km); + else sprintf(s1, "%12.4f %12.4f", kc, km); if (Link[i].Type == CV) sprintf(s2,"CV"); else if (Link[i].Stat == CLOSED) sprintf(s2,"CLOSED"); else strcpy(s2,""); fprintf(f,"\n%s %s %s",s,s1,s2); -/*** End of update ***/ } } @@ -229,14 +206,14 @@ int saveinpfile(char *fname) for (i=1; i<=Npumps; i++) { n = Pump[i].Link; - sprintf(s," %-15s %-15s %-15s", + sprintf(s," %-31s %-31s %-31s", Link[n].ID, Node[Link[n].N1].ID, Node[Link[n].N2].ID); /* Pump has constant power */ if (Pump[i].Ptype == CONST_HP) - sprintf(s1, " POWER %.2f", Link[n].Km); + sprintf(s1, " POWER %.4f", Link[n].Km); /* Pump has a head curve */ else if ((j = Pump[i].Hcurve) > 0) @@ -245,7 +222,7 @@ int saveinpfile(char *fname) /* Old format used for pump curve */ else { - fprintf(f, "\n%s %12.2f %12.2f %12.2f 0.0 %12.2f",s, + fprintf(f, "\n%s %12.4f %12.4f %12.4f 0.0 %12.4f",s, -Pump[i].H0*Ucf[HEAD], (-Pump[i].H0 - Pump[i].R*pow(Pump[i].Q0,Pump[i].N))*Ucf[HEAD], Pump[i].Q0*Ucf[FLOW], @@ -260,7 +237,7 @@ int saveinpfile(char *fname) strcat(s,s1); if (Link[n].Kc != 1.0) - sprintf(s1, " SPEED %.2f", Link[n].Kc); + sprintf(s1, " SPEED %.4f", Link[n].Kc); else strcpy(s1,""); strcat(s,s1); @@ -285,7 +262,7 @@ int saveinpfile(char *fname) } km = Link[n].Km*SQR(d)*SQR(d)/0.02517; - sprintf(s," %-15s %-15s %-15s %12.2f %5s", + sprintf(s," %-31s %-31s %-31s %12.4f %5s", Link[n].ID, Node[Link[n].N1].ID, Node[Link[n].N2].ID, @@ -293,8 +270,8 @@ int saveinpfile(char *fname) LinkTxt[Link[n].Type]); if (Link[n].Type == GPV && (j = ROUND(Link[n].Kc)) > 0) - sprintf(s1,"%-15s %12.2f", Curve[j].ID, km); - else sprintf(s1,"%12.2f %12.2f",kc,km); + sprintf(s1,"%-31s %12.4f", Curve[j].ID, km); + else sprintf(s1,"%12.4f %12.4f",kc,km); fprintf(f, "\n%s %s", s,s1); } @@ -302,15 +279,12 @@ int saveinpfile(char *fname) /* Write [DEMANDS] section */ fprintf(f, "\n\n[DEMANDS]"); - -/*** Updated 11/19/01 ***/ ucf = Ucf[DEMAND]; - for (i=1; i<=Njuncs; i++) { for (demand = Node[i].D; demand != NULL; demand = demand->next) { - sprintf(s," %-15s %12.2f",Node[i].ID,ucf*demand->Base); + sprintf(s," %-31s %14.6f",Node[i].ID,ucf*demand->Base); if ((j = demand->Pat) > 0) sprintf(s1," %s",Pattern[j].ID); else strcpy(s1,""); fprintf(f,"\n%s %s",s,s1); @@ -324,7 +298,7 @@ int saveinpfile(char *fname) { if (Node[i].Ke == 0.0) continue; ke = Ucf[FLOW]/pow(Ucf[PRESSURE]*Node[i].Ke,(1.0/Qexp)); - fprintf(f,"\n %-15s %12.2f",Node[i].ID,ke); + fprintf(f,"\n %-31s %14.6f",Node[i].ID,ke); } /* Write [STATUS] section */ @@ -335,7 +309,7 @@ int saveinpfile(char *fname) if (Link[i].Type <= PUMP) { if (Link[i].Stat == CLOSED) - fprintf(f, "\n %-15s %s",Link[i].ID,StatTxt[CLOSED]); + fprintf(f, "\n %-31s %s",Link[i].ID,StatTxt[CLOSED]); /* Write pump speed here for pumps with old-style pump curve input */ else if (Link[i].Type == PUMP) @@ -346,7 +320,7 @@ int saveinpfile(char *fname) Pump[n].Ptype != CONST_HP && Link[i].Kc != 1.0 ) - fprintf(f, "\n %-15s %-.2f",Link[i].ID, Link[i].Kc); + fprintf(f, "\n %-31s %-.4f",Link[i].ID, Link[i].Kc); } } @@ -354,9 +328,9 @@ int saveinpfile(char *fname) else if (Link[i].Kc == MISSING) { if (Link[i].Stat == OPEN) - fprintf(f, "\n %-15s %s",Link[i].ID,StatTxt[OPEN]); + fprintf(f, "\n %-31s %s",Link[i].ID,StatTxt[OPEN]); if (Link[i].Stat == CLOSED) - fprintf(f, "\n%-15s %s",Link[i].ID,StatTxt[CLOSED]); + fprintf(f, "\n%-31s %s",Link[i].ID,StatTxt[CLOSED]); } } @@ -368,7 +342,7 @@ int saveinpfile(char *fname) { for (j=0; jType], source->C0); @@ -467,7 +441,7 @@ int saveinpfile(char *fname) for (i=1; i<=Ntanks; i++) { if (Tank[i].A == 0.0) continue; - fprintf(f, "\n %-15s %-8s %12.4f", + fprintf(f, "\n %-31s %-8s %12.4f", Node[Tank[i].Node].ID, MixTxt[Tank[i].MixModel], (Tank[i].V1max/Tank[i].Vmax)); @@ -479,25 +453,25 @@ int saveinpfile(char *fname) fprintf(f, "\n ORDER BULK %-.2f", BulkOrder); fprintf(f, "\n ORDER WALL %-.0f", WallOrder); fprintf(f, "\n ORDER TANK %-.2f", TankOrder); - fprintf(f, "\n GLOBAL BULK %-.4f", Kbulk*SECperDAY); - fprintf(f, "\n GLOBAL WALL %-.4f", Kwall*SECperDAY); + fprintf(f, "\n GLOBAL BULK %-.6f", Kbulk*SECperDAY); + fprintf(f, "\n GLOBAL WALL %-.6f", Kwall*SECperDAY); if (Climit > 0.0) - fprintf(f, "\n LIMITING POTENTIAL %-.4f", Climit); + fprintf(f, "\n LIMITING POTENTIAL %-.6f", Climit); if (Rfactor != MISSING && Rfactor != 0.0) - fprintf(f, "\n ROUGHNESS CORRELATION %-.4f",Rfactor); + fprintf(f, "\n ROUGHNESS CORRELATION %-.6f",Rfactor); for (i=1; i<=Nlinks; i++) { if (Link[i].Type > PIPE) continue; if (Link[i].Kb != Kbulk) - fprintf(f, "\n BULK %-15s %-.4f",Link[i].ID,Link[i].Kb*SECperDAY); + fprintf(f, "\n BULK %-31s %-.6f",Link[i].ID,Link[i].Kb*SECperDAY); if (Link[i].Kw != Kwall) - fprintf(f, "\n WALL %-15s %-.4f",Link[i].ID,Link[i].Kw*SECperDAY); + fprintf(f, "\n WALL %-31s %-.6f",Link[i].ID,Link[i].Kw*SECperDAY); } for (i=1; i<=Ntanks; i++) { if (Tank[i].A == 0.0) continue; if (Tank[i].Kb != Kbulk) - fprintf(f, "\n TANK %-15s %-.4f",Node[Tank[i].Node].ID, + fprintf(f, "\n TANK %-31s %-.6f",Node[Tank[i].Node].ID, Tank[i].Kb*SECperDAY); } @@ -508,20 +482,18 @@ int saveinpfile(char *fname) fprintf(f, "\n GLOBAL PRICE %-.4f", Ecost); if (Epat != 0) fprintf(f, "\n GLOBAL PATTERN %s", Pattern[Epat].ID); - fprintf(f, "\n GLOBAL EFFIC %-.2f", Epump); + fprintf(f, "\n GLOBAL EFFIC %-.4f", Epump); fprintf(f, "\n DEMAND CHARGE %-.4f", Dcost); for (i=1; i<=Npumps; i++) { if (Pump[i].Ecost > 0.0) - fprintf(f, "\n PUMP %-15s PRICE %-.4f", + fprintf(f, "\n PUMP %-31s PRICE %-.4f", Link[Pump[i].Link].ID,Pump[i].Ecost); if (Pump[i].Epat > 0.0) - fprintf(f, "\n PUMP %-15s PATTERN %s", + fprintf(f, "\n PUMP %-31s PATTERN %s", Link[Pump[i].Link].ID,Pattern[Pump[i].Epat].ID); - -/*** Updated 3/1/01 ***/ if (Pump[i].Ecurve > 0.0) - fprintf(f, "\n PUMP %-15s EFFIC %s", + fprintf(f, "\n PUMP %-31s EFFIC %s", Link[Pump[i].Link].ID,Curve[Pump[i].Ecurve].ID); } @@ -558,22 +530,22 @@ int saveinpfile(char *fname) if (Qualflag == CHEM) fprintf(f, "\n QUALITY %s %s", ChemName, ChemUnits); if (Qualflag == TRACE) - fprintf(f, "\n QUALITY TRACE %-15s", Node[TraceNode].ID); + fprintf(f, "\n QUALITY TRACE %-31s", Node[TraceNode].ID); if (Qualflag == AGE) fprintf(f, "\n QUALITY AGE"); if (Qualflag == NONE) fprintf(f, "\n QUALITY NONE"); - fprintf(f, "\n DEMAND MULTIPLIER %-.2f", Dmult); - -/*** Updated 11/19/01 ***/ - fprintf(f, "\n EMITTER EXPONENT %-.2f", 1.0/Qexp); - - fprintf(f, "\n VISCOSITY %-.4f", Viscos/VISCOS); - fprintf(f, "\n DIFFUSIVITY %-.4f", Diffus/DIFFUS); - fprintf(f, "\n SPECIFIC GRAVITY %-.4f", SpGrav); + fprintf(f, "\n DEMAND MULTIPLIER %-.4f", Dmult); + fprintf(f, "\n EMITTER EXPONENT %-.4f", 1.0/Qexp); + fprintf(f, "\n VISCOSITY %-.6f", Viscos/VISCOS); + fprintf(f, "\n DIFFUSIVITY %-.6f", Diffus/DIFFUS); + fprintf(f, "\n SPECIFIC GRAVITY %-.6f", SpGrav); fprintf(f, "\n TRIALS %-d", MaxIter); fprintf(f, "\n ACCURACY %-.8f", Hacc); fprintf(f, "\n TOLERANCE %-.8f", Ctol*Ucf[QUALITY]); + fprintf(f, "\n CHECKFREQ %-d", CheckFreq); + fprintf(f, "\n MAXCHECK %-d", MaxCheck); + fprintf(f, "\n DAMPLIMIT %-.8f", DampLimit); /* Write [REPORT] section */ @@ -622,33 +594,27 @@ int saveinpfile(char *fname) } } } - for (i=0; i -BIG) - fprintf(f, "\n %-20sABOVE %.4f", Field[i].Name, Field[i].RptLim[HI]); + fprintf(f, "\n %-20sABOVE %.6f", Field[i].Name, Field[i].RptLim[HI]); } - else fprintf(f, "\n %-20sNO", Field[i].Name); -/********************************************************************/ + else fprintf(f, "\n %-20sNO", Field[i].Name); } fprintf(f, "\n"); -/*** Updated *****************************************/ //(2.00.11 - LR) -/* Copy data from scratch file to new input file */ - if (ftmp != NULL) - { - fseek(ftmp, 0, SEEK_SET); - while ( (j = fgetc(ftmp)) != EOF ) fputc(j, f); - fclose(ftmp); - } -/*****************************************************/ +/* Save auxilary data to new input file */ + + saveauxdata(f); - fprintf(f, "\n\n[END]"); +/* Close the new input file */ + + fprintf(f, "\n[END]"); fclose(f); return(0); } diff --git a/src/input1.c b/src/input1.c index eace922..4ef941c 100644 --- a/src/input1.c +++ b/src/input1.c @@ -8,6 +8,7 @@ DATE: 5/30/00 9/7/00 11/19/01 6/24/02 + 2/14/08 (2.00.12) AUTHOR: L. Rossman US EPA - NRMRL @@ -36,7 +37,7 @@ AUTHOR: L. Rossman --------------------- Module Global Variables ---------------------- */ -#define MAXITER 40 /* Default max. # hydraulic iterations */ +#define MAXITER 200 /* Default max. # hydraulic iterations */ //(2.00.12 - LR) #define HACC 0.001 /* Default hydraulics convergence ratio */ #define HTOL 0.0005 /* Default hydraulic head tolerance (ft) */ @@ -57,6 +58,7 @@ AUTHOR: L. Rossman #define RQTOL 1E-7 /* Default low flow resistance tolerance */ #define CHECKFREQ 2 /* Default status check frequency */ #define MAXCHECK 10 /* Default # iterations for status checks */ +#define DAMPLIMIT 0 /* Default damping threshold */ //(2.00.12 - LR) extern char *Fldname[]; /* Defined in enumstxt.h in EPANET.C */ extern char *RptFlowUnitsTxt[]; @@ -96,6 +98,8 @@ void setdefaults() strncpy(Title[0],"",MAXMSG); strncpy(Title[1],"",MAXMSG); strncpy(Title[2],"",MAXMSG); + strncpy(TmpDir,"",MAXFNAME); //(2.00.12 - LR) + strncpy(TmpFname,"",MAXFNAME); //(2.00.12 - LR) strncpy(HydFname,"",MAXFNAME); strncpy(MapFname,"",MAXFNAME); strncpy(ChemName,t_CHEMICAL,MAXID); @@ -146,6 +150,7 @@ void setdefaults() RQtol = RQTOL; /* Default hydraulics parameters */ CheckFreq = CHECKFREQ; MaxCheck = MAXCHECK; + DampLimit = DAMPLIMIT; //(2.00.12 - LR) } /* End of setdefaults */ diff --git a/src/input3.c b/src/input3.c index 7a2771d..a8fdfe0 100644 --- a/src/input3.c +++ b/src/input3.c @@ -10,6 +10,7 @@ DATE: 5/30/00 3/1/01 6/24/02 8/15/07 (2.00.11) + 2/14/08 (2.00.12) AUTHOR: L. Rossman US EPA - NRMRL @@ -1599,7 +1600,8 @@ int optionvalue(int n) ** QTOL value ** RQTOL value ** CHECKFREQ value -** MAXCHECK value +** MAXCHECK value +** DAMPLIMIT value //(2.00.12 - LR) **-------------------------------------------------------------- */ { @@ -1633,6 +1635,13 @@ int optionvalue(int n) return(0); } +/* Check for Damping Limit option */ //(2.00.12 - LR) + if (match(Tok[0],w_DAMPLIMIT)) + { + DampLimit = y; + return(0); + } + /* All other options must be > 0 */ if (y <= 0.0) return(213); diff --git a/src/quality.c b/src/quality.c index 15be2a1..28469de 100644 --- a/src/quality.c +++ b/src/quality.c @@ -8,6 +8,7 @@ DATE: 5/29/00 9/7/00 10/25/00 8/15/07 (2.00.11) + 2/14/08 (2.00.12) AUTHOR: L. Rossman US EPA - NRMRL @@ -75,7 +76,10 @@ double *MassIn; /* Total mass inflow to node */ double Sc; /* Schmidt Number */ double Bucf; /* Bulk reaction units conversion factor */ double Tucf; /* Tank reaction units conversion factor */ -char Reactflag; /* Reaction indicator */ + +/*** Moved to vars.h ***/ //(2.00.12 - LR) +//char Reactflag; /* Reaction indicator */ + char OutOfMemory; /* Out of memory indicator */ static alloc_handle_t *SegPool; // Memory pool for water quality segments //(2.00.11 - LR) @@ -688,10 +692,12 @@ void accumulate(long dt) i = UP_NODE(k); /* Upstream node */ j = DOWN_NODE(k); /* Downstream node */ v = ABS(Q[k])*dt; /* Flow volume */ + +//// Start of deprecated code segment //// //(2.00.12 - LR) /* If link volume < flow volume, then transport upstream */ /* quality to downstream node and remove all link segments. */ - if (LINKVOL(k) < v) +/* if (LINKVOL(k) < v) { VolIn[j] += v; seg = FirstSeg[k]; @@ -700,10 +706,14 @@ void accumulate(long dt) MassIn[j] += v*cseg; removesegs(k); } - +*/ /* Otherwise remove flow volume from leading segments */ /* and accumulate flow mass at downstream node */ - else while (v > 0.0) + //else + +//// End of deprecated code segment. //// //(2.00.12 - LR) + + while (v > 0.0) //(2.00.12 - LR) { /* Identify leading segment in pipe */ seg = FirstSeg[k]; @@ -1003,6 +1013,37 @@ void updatetanks(long dt) } +//// Deprecated version of tankmix1 //// //(2.00.12 - LR) +//void tankmix1(int i, long dt) +/* +**--------------------------------------------- +** Input: i = tank index +** dt = current WQ time step +** Output: none +** Purpose: complete mix tank model +**--------------------------------------------- +*/ +//{ +// int n; +// double cin; + +// /* Blend inflow with contents */ +// n = Tank[i].Node; +// if (VolIn[n] > 0.0) cin = MassIn[n]/VolIn[n]; +// else cin = 0.0; +// if (Tank[i].V > 0.0) +// Tank[i].C = tankreact(Tank[i].C,Tank[i].V,Tank[i].Kb,dt) + +// (cin - Tank[i].C)*VolIn[n]/Tank[i].V; +// else Tank[i].C = cin; +// Tank[i].C = MAX(0.0, Tank[i].C); + +// /* Update tank volume & nodal quality */ +// Tank[i].V += D[n]*dt; +// C[n] = Tank[i].C; +//} + + +//// New version of tankmix1 //// //(2.00.12 - LR) void tankmix1(int i, long dt) /* **--------------------------------------------- @@ -1015,23 +1056,32 @@ void tankmix1(int i, long dt) { int n; double cin; + double c, cmax, vold, vin; - /* Blend inflow with contents */ + /* React contents of tank */ + c = tankreact(Tank[i].C,Tank[i].V,Tank[i].Kb,dt); + + /* Determine tank & volumes */ + vold = Tank[i].V; n = Tank[i].Node; - if (VolIn[n] > 0.0) cin = MassIn[n]/VolIn[n]; - else cin = 0.0; - if (Tank[i].V > 0.0) - Tank[i].C = tankreact(Tank[i].C,Tank[i].V,Tank[i].Kb,dt) + - (cin - Tank[i].C)*VolIn[n]/Tank[i].V; - else Tank[i].C = cin; - Tank[i].C = MAX(0.0, Tank[i].C); - - /* Update tank volume & nodal quality */ Tank[i].V += D[n]*dt; + vin = VolIn[n]; + + /* Compute inflow concen. */ + if (vin > 0.0) cin = MassIn[n]/vin; + else cin = 0.0; + cmax = MAX(c, cin); + + /* Mix inflow with tank contents */ + if (vin > 0.0) c = (c*vold + cin*vin)/(vold + vin); + c = MIN(c, cmax); + c = MAX(c, 0.0); + Tank[i].C = c; C[n] = Tank[i].C; } /*** Updated 10/25/00 ***/ +//// New version of tankmix2 //// //(2.00.12 - LR) void tankmix2(int i, long dt) /* **------------------------------------------------ @@ -1045,11 +1095,11 @@ void tankmix2(int i, long dt) */ { int k,n; - long tstep, tstar; - double cin, /* Inflow quality */ - qin, /* Inflow rate */ - qout, /* Outflow rate */ - qnet; /* Net flow rate */ + double cin, /* Inflow quality */ + vin, /* Inflow volume */ + vt, /* Transferred volume */ + vnet, /* Net volume change */ + v1max; /* Full mixing zone volume */ Pseg seg1,seg2; /* Compartment segments */ /* Identify segments for each compartment */ @@ -1058,88 +1108,64 @@ void tankmix2(int i, long dt) seg2 = FirstSeg[k]; if (seg1 == NULL || seg2 == NULL) return; + /* React contents of each compartment */ + seg1->c = tankreact(seg1->c,seg1->v,Tank[i].Kb,dt); + seg2->c = tankreact(seg2->c,seg2->v,Tank[i].Kb,dt); + /* Find inflows & outflows */ n = Tank[i].Node; - qnet = D[n]; - qin = VolIn[n]/(double)dt; - qout = qin - qnet; - if (qin > 0.0) cin = MassIn[n]/VolIn[n]; + vnet = D[n]*dt; + vin = VolIn[n]; + if (vin > 0.0) cin = MassIn[n]/vin; else cin = 0.0; - Tank[i].V += qnet*dt; + v1max = Tank[i].V1max; - /* Case of no net volume change */ - if (ABS(qnet) < TINY) + /* Tank is filling */ + vt = 0.0; + if (vnet > 0.0) { - seg1->c = tankreact(seg1->c,seg1->v,Tank[i].Kb,dt); - seg2->c = tankreact(seg2->c,seg2->v,Tank[i].Kb,dt); - } - - /* Case of net filling (qnet > 0) */ - else if (qnet > 0.0) - { - - /* Case where ambient zone empty & mixing zone filling */ - if (seg2->v <= 0.0) + vt = MAX(0.0, (seg1->v + vnet - v1max)); + if (vin > 0.0) { - tstar = (long) ((Tank[i].V1max-(seg1->v))/qnet); - tstep = MIN(dt, tstar); - if (seg1->v > 0.0) - seg1->c = tankreact(seg1->c,seg1->v,Tank[i].Kb,tstep) + - qin*tstep*(cin-(seg1->c))/(seg1->v); - else seg1->c = cin; - seg1->c = MAX(0.0, seg1->c); - seg1->v += qnet*tstep; - seg2->c = 0.0; - dt -= tstep; + seg1->c = ((seg1->c)*(seg1->v) + cin*vin) / (seg1->v + vin); } - - /* Case where mixing zone full & ambient zone filling */ - if (dt > 1) + if (vt > 0.0) { - seg1->c = tankreact(seg1->c,seg1->v,Tank[i].Kb,dt) + - qin*dt*(cin-(seg1->c))/(seg1->v); - seg1->c = MAX(0.0, seg1->c); - if (seg2->v <= 0.0) seg2->c = seg1->c; - else - seg2->c = tankreact(seg2->c,seg2->v,Tank[i].Kb,dt) + - qnet*dt*((seg1->c)-(seg2->c))/(seg2->v); - seg2->c = MAX(0.0, seg2->c); - seg2->v += qnet*dt; - } - } - - /* Case of net emptying (qnet < 0) */ - else if (qnet < 0.0 && seg1->v > 0.0) - { - - /* Case where mixing zone full & ambient zone draining */ - if ((seg2->v) > 0.0) - { - tstar = (long)(seg2->v/-qnet); - tstep = MIN(dt, tstar); - seg1->c = tankreact(seg1->c,seg1->v,Tank[i].Kb,tstep) + - (qin*cin - qnet*(seg2->c) - qout*(seg1->c))* - tstep/(seg1->v); - seg2->c = tankreact(seg2->c,seg2->v,Tank[i].Kb,tstep); - seg2->v += qnet*tstep; - seg1->c = MAX(0.0, seg1->c); - seg2->c = MAX(0.0, seg2->c); - seg2->v = MAX(0.0, seg2->v); - dt -= tstep; - } - - /* Case where ambient zone empty & mixing zone draining */ - if (dt > 1) - { - seg1->c = tankreact(seg1->c,seg1->v,Tank[i].Kb,dt) + - qin*dt*(cin-(seg1->c))/(seg1->v); - seg1->v += qnet*dt; - seg1->c = MAX(0.0, seg1->c); - seg1->v = MAX(0.0, seg1->v); - seg2->c = 0.0; + seg2->c = ((seg2->c)*(seg2->v) + (seg1->c)*vt) / (seg2->v + vt); } } + /* Tank is emptying */ + if (vnet < 0.0) + { + if (seg2->v > 0.0) + { + vt = MIN(seg2->v, (-vnet)); + } + if (vin + vt > 0.0) + { + seg1->c = ((seg1->c)*(seg1->v) + cin*vin + (seg2->c)*vt) / + (seg1->v + vin + vt); + } + } + + /* Update segment volumes */ + if (vt > 0.0) + { + seg1->v = v1max; + if (vnet > 0.0) seg2->v += vt; + else seg2->v = MAX(0.0, ((seg2->v)-vt)); + } + else + { + seg1->v += vnet; + seg1->v = MIN(seg1->v, v1max); + seg1->v = MAX(0.0, seg1->v); + seg2->v = 0.0; + } + Tank[i].V += vnet; + Tank[i].V = MAX(0.0, Tank[i].V); + /* Use quality of mixed compartment (seg1) to */ /* represent quality of tank since this is where */ /* outflow begins to flow from */ @@ -1185,6 +1211,7 @@ void tankmix3(int i, long dt) if (vin > 0.0) cin = MassIn[n]/VolIn[n]; else cin = 0.0; Tank[i].V += vnet; + Tank[i].V = MAX(0.0, Tank[i].V); //(2.00.12 - LR) /* Withdraw flow from first segment */ vsum = 0.0; @@ -1201,10 +1228,13 @@ void tankmix3(int i, long dt) vout -= vseg; /* Remaining flow volume */ if (vout >= 0.0 && vseg >= seg->v) /* Seg used up */ { - FirstSeg[k] = seg->prev; - if (FirstSeg[k] == NULL) LastSeg[k] = NULL; - seg->prev = FreeSeg; - FreeSeg = seg; + if (seg->prev) //(2.00.12 - LR) + { //(2.00.12 - LR) + FirstSeg[k] = seg->prev; + //if (FirstSeg[k] == NULL) LastSeg[k] = NULL; //(2.00.12 - LR) + seg->prev = FreeSeg; + FreeSeg = seg; + } //(2.00.12 - LR) } else /* Remaining volume in segment */ { @@ -1271,6 +1301,7 @@ void tankmix4(int i, long dt) if (vin > 0.0) cin = MassIn[n]/VolIn[n]; else cin = 0.0; Tank[i].V += vnet; + Tank[i].V = MAX(0.0, Tank[i].V); //(2.00.12 - LR) Tank[i].C = LastSeg[k]->c; /* If tank filling, then create new last seg */ @@ -1317,10 +1348,13 @@ void tankmix4(int i, long dt) vnet -= vseg; if (vnet >= 0.0 && vseg >= seg->v) /* Seg used up */ { - LastSeg[k] = seg->prev; - if (LastSeg[k] == NULL) FirstSeg[k] = NULL; - seg->prev = FreeSeg; - FreeSeg = seg; + if (seg->prev) //(2.00.12 - LR) + { //(2.00.12 - LR) + LastSeg[k] = seg->prev; + //if (LastSeg[k] == NULL) FirstSeg[k] = NULL; //(2.00.12 - LR) + seg->prev = FreeSeg; + FreeSeg = seg; + } //(2.00.12 - LR) } else /* Remaining volume in segment */ { diff --git a/src/report.c b/src/report.c index 895c844..ad22934 100644 --- a/src/report.c +++ b/src/report.c @@ -7,6 +7,7 @@ VERSION: 2.00 DATE: 5/30/00 6/24/02 8/15/07 (2.00.11) + 2/14/08 (2.00.12) AUTHOR: L. Rossman US EPA - NRMRL @@ -200,6 +201,14 @@ void writesummary() writeline(s); sprintf(s,FMT27,Hacc); writeline(s); + + sprintf(s,FMT27a,CheckFreq); //(2.00.12 - LR) + writeline(s); //(2.00.12 - LR) + sprintf(s,FMT27b,MaxCheck); //(2.00.12 - LR) + writeline(s); //(2.00.12 - LR) + sprintf(s,FMT27c,DampLimit); //(2.00.12 - LR) + writeline(s); //(2.00.12 - LR) + sprintf(s,FMT28,MaxIter); writeline(s); if (Qualflag == NONE || Dur == 0.0) @@ -1210,5 +1219,4 @@ int getnodetype(int i) return(2); } - /********************* END OF REPORT.C ********************/ diff --git a/src/text.h b/src/text.h index 9815b64..dfdee67 100644 --- a/src/text.h +++ b/src/text.h @@ -7,6 +7,7 @@ VERSION: 2.00 DATE: 5/8/00 10/25/00 8/15/07 (2.00.11) + 2/14/08 (2.00.12) AUTHOR: L. Rossman US EPA - NRMRL @@ -137,6 +138,7 @@ AUTHOR: L. Rossman #define w_RQTOL "RQTOL" #define w_CHECKFREQ "CHECKFREQ" #define w_MAXCHECK "MAXCHECK" +#define w_DAMPLIMIT "DAMPLIMIT" //(2.00.12 - LR) #define w_SECONDS "SEC" #define w_MINUTES "MIN" @@ -328,7 +330,7 @@ AUTHOR: L. Rossman #define LOGO4 \ "* Analysis for Pipe Networks *" #define LOGO5 \ -"* Version 2.00.11 *" //(2.00.11 - LR) +"* Version 2.00.12 *" //(2.00.12 - LR) #define LOGO6 \ "******************************************************************" #define FMT01 "\n... EPANET Version 2.0\n" @@ -357,6 +359,11 @@ AUTHOR: L. Rossman #define FMT25 " Headloss Formula .................. %s" #define FMT26 " Hydraulic Timestep ................ %-.2f %s" #define FMT27 " Hydraulic Accuracy ................ %-.6f" + +#define FMT27a " Status Check Frequency ............ %-d" //(2.00.12 - LR) +#define FMT27b " Maximum Trials Checked ............ %-d" //(2.00.12 - LR) +#define FMT27c " Damping Limit Threshold ........... %-.6f" //(2.00.12 - LR) + #define FMT28 " Maximum Trials .................... %-d" #define FMT29 " Quality Analysis .................. None" #define FMT30 " Quality Analysis .................. %s" @@ -393,6 +400,10 @@ AUTHOR: L. Rossman #define FMT57 " %s %s switched from %s to %s" #define FMT58 "%10s: Balanced after %-d trials" #define FMT59 "%10s: Unbalanced after %-d trials (flow change = %-.6f)" + +#define FMT60a " Max. flow imbalance is %.4f %s at Node %s" //(2.00.12 - LR) +#define FMT60b " Max. head imbalance is %.4f %s at Link %s" //(2.00.12 - LR) + #define FMT61 "%10s: Valve %s caused ill-conditioning" #define FMT62 "%10s: System ill-conditioned at node %s" #define FMT63 "%10s: %s %s changed by rule %s" diff --git a/src/toolkit.h b/src/toolkit.h index 580de2c..f1a6f5d 100644 --- a/src/toolkit.h +++ b/src/toolkit.h @@ -8,6 +8,7 @@ DATE: 5/8/00 10/25/00 3/1/01 8/15/07 (2.00.11) + 2/14/08 (2.00.12) AUTHOR: L. Rossman US EPA - NRMRL @@ -50,6 +51,14 @@ AUTHOR: L. Rossman #define EN_MIXMODEL 15 #define EN_MIXZONEVOL 16 +#define EN_TANKDIAM 17 +#define EN_MINVOLUME 18 +#define EN_VOLCURVE 19 +#define EN_MINLEVEL 20 +#define EN_MAXLEVEL 21 +#define EN_MIXFRACTION 22 +#define EN_TANK_KBULK 23 + #define EN_DIAMETER 0 /* Link parameters */ #define EN_LENGTH 1 #define EN_ROUGHNESS 2 @@ -134,6 +143,11 @@ AUTHOR: L. Rossman #define EN_MAXIMUM 3 #define EN_RANGE 4 +#define EN_MIX1 0 /* Tank mixing models */ +#define EN_MIX2 1 +#define EN_FIFO 2 +#define EN_LIFO 3 + #define EN_NOSAVE 0 /* Save-results-to-file flag */ #define EN_SAVE 1 @@ -200,6 +214,7 @@ AUTHOR: L. Rossman int DLLEXPORT ENsetcontrol(int, int, int, float, int, float); int DLLEXPORT ENsetnodevalue(int, int, float); int DLLEXPORT ENsetlinkvalue(int, int, float); + int DLLEXPORT ENaddpattern(char *); int DLLEXPORT ENsetpattern(int, float *, int); int DLLEXPORT ENsetpatternvalue(int, int, float); int DLLEXPORT ENsettimeparam(int, long); diff --git a/src/types.h b/src/types.h index 43f16f8..37030f4 100644 --- a/src/types.h +++ b/src/types.h @@ -11,6 +11,7 @@ DATE: 5/8/00 12/6/01 6/24/02 8/15/07 (2.00.11) + 2/14/08 (2.00.12) AUTHOR: L. Rossman US EPA - NRMRL @@ -20,6 +21,13 @@ AUTHOR: L. Rossman /*********************************************************/ /* All floats have been re-declared as doubles (7/3/07). */ /*********************************************************/ +/* +------------------------------------------- + Definition of 4-byte integers & reals +------------------------------------------- +*/ +typedef float REAL4; //(2.00.11 - LR) +typedef int INT4; //(2.00.12 - LR) /* ----------------------------- @@ -27,7 +35,7 @@ AUTHOR: L. Rossman ----------------------------- */ /*** Updated ***/ -#define CODEVERSION 20011 //(2.00.11 - LR) +#define CODEVERSION 20012 //(2.00.12 - LR) #define MAGICNUMBER 516114521 #define VERSION 200 #define EOFMARK 0x1A /* Use 0x04 for UNIX systems */ @@ -118,8 +126,6 @@ AUTHOR: L. Rossman Global Data Structures ------------------------------------------------------ */ -typedef float REAL4; //(2.00.11 - LR) -typedef long INT4; //(2.00.11 - LR) struct IDstring /* Holds component ID labels */ { diff --git a/src/vars.h b/src/vars.h index 5740f2c..1c92fcd 100644 --- a/src/vars.h +++ b/src/vars.h @@ -5,6 +5,7 @@ VERSION: 2.00 DATE: 5/8/00 6/24/02 + 2/14/08 (2.00.12) AUTHOR: L. Rossman US EPA - NRMRL @@ -25,6 +26,8 @@ EXTERN char Msg[MAXMSG+1], /* Text of output message */ HydFname[MAXFNAME+1], /* Hydraulics file name */ OutFname[MAXFNAME+1], /* Binary output file name */ MapFname[MAXFNAME+1], /* Map file name */ + TmpFname[MAXFNAME+1], /* Temporary file name */ //(2.00.12 - LR) + TmpDir[MAXFNAME+1], /* Temporary directory name */ //(2.00.12 - LR) Title[MAXTITLE][MAXMSG+1], /* Problem title */ ChemName[MAXID+1], /* Name of chemical */ ChemUnits[MAXID+1], /* Units of chemical */ @@ -33,8 +36,10 @@ EXTERN char Msg[MAXMSG+1], /* Text of output message */ /*** Updated 6/24/02 ***/ Atime[13], /* Clock time (hrs:min:sec) */ + Outflag, /* Output file flag */ //(2.00.12 - LR) Hydflag, /* Hydraulics flag */ Qualflag, /* Water quality flag */ + Reactflag, /* Reaction indicator */ //(2.00.12 - LR) Unitsflag, /* Unit system flag */ Flowflag, /* Flow units flag */ Pressflag, /* Pressure units flag */ @@ -95,6 +100,7 @@ EXTERN double Ucf[MAXVAR], /* Unit conversion factors */ Qexp, /* Exponent in orifice formula */ Dmult, /* Demand multiplier */ Hacc, /* Hydraulics solution accuracy */ + DampLimit, /* Solution damping threshold */ //(2.00.12 - LR) BulkOrder, /* Bulk flow reaction order */ WallOrder, /* Pipe wall reaction order */ TankOrder, /* Tank reaction order */