Everglades Landscape Model (ELM) Home Page
Main Page | Data Structures | Directories | File List | Data Fields | Globals

watmgmt.h File Reference

Header file for the Water Managment calculations. More...

#include "globals.h"

Include dependency graph for watmgmt.h:

Include dependency graph

This graph shows which files directly or indirectly include this file:

Included by dependency graph

Go to the source code of this file.

Data Structures

struct  Chan_reach
struct  Chan
struct  Cells
struct  Structure
 Water control structure attributes. More...
struct  Schedule
struct  Points
struct  HistData

Defines

#define STEP(x)   ( ((x)>=0) ? 1.0 : 0.0 )
#define MAX_H_STRUCT   180
#define LEFT   1
#define RIGHT   -1
#define NONE   0
#define EAST   1
#define SOUTH   2
#define ALL   3

Functions

void CanalReInit ()
 Re-initialize the depths and concentrations in canals under the Positional Analysis mode.
void Canal_Network_Init (float baseDatum, float *elev)
 Initialize the water managment network topology.
void Run_Canal_Network (float *SWH, float *ElevMap, float *MC, float *GWV, float *poros, float *GWcond, double *NA, double *PA, double *SA, double *GNA, double *GPA, double *GSA, float *Unsat, float *sp_yield)
 Runs the water management network.
void ReadStructures (char *name, float BASE_DATUM)
 Read attribute data on water control structures.
void Channel_configure (float *ElevMap, struct Chan *channel_first)
 Configure attributes of grid cells interacting with canal vectors.
void getCanalElev (int chan_n)
 Assign elevation to canal using slope & distance from start point.
void FluxChannel (int chan_n, float *SWH, float *ElevMap, float *MC, float *GWH, float *poros, float *GWcond, double *NA, double *PA, double *SA, double *GNA, double *GPA, double *GSA, float *Unsat, float *sp_yield)
 Runs the iterative algorithm over the network of canals.
void Flows_in_Structures (float *SWH, float *Elev, float *MC, double *NA, double *PA, double *SA)
 Determine flows through water control structures.
float f_Manning (float delta, float SWater, float SW_coef)
 Surface water exchange between cell and canal.
float f_Ground (float dh, float height, float GW_coef, float l_Length)
 Subsurface ground water exchange between cell and canal.
float GetGraph (struct Schedule *graph, float x)
 Graph to read the time dependent schedules.
ChanReadChanStruct (char *filename)
 Reads the information about canals from data file.
CellsMarkCell (int x, int y, int index, float length, struct Cells *cell, int ch_number, int levee, int xm, int ym, int c_num, int *marked, float distTot)
 Mark cell function.
ScheduleRead_schedule (char *sched_name, char *filename, float BASE_DATUM)
 Read target stage schedule as a recurring time-series graph.
int Wrdcmp (char *s, char *t)
 Compare two words.
float UTM2kmx (double UTM)
 Converter from meters in UTM to grid cell row location (zero origin).
float UTM2kmy (double UTM)
 Converter from meters in UTM to grid cell column location (zero origin).
double julday (int mon, int day, int year, int h, int mi, double se)
 Determine the Julian calendar day from a Gregorian date.
void init_pvar (VOIDP Map, UCHAR *mask, unsigned char Mtype, float iv)
 Initialize a variable to a value.
VOIDP nalloc (unsigned mem_size, const char var_name[])
 Allocate memory for a variable.
char * Scip (char *s, char SYM)
 Skip ahead in a string until next field.
int isInteger (char *target_str)
 Determine if an integer is present in string.
int isFloat (char *target_str)
 Determine if a float is present in string.
float FMOD (float x, float y)
 Modulus of a pair of (double) arguments.
int Round (float x)
 truncate (not rounding) to an integer
float Flux_SWcells (int i0, int i1, int j0, int j1, float *SWater, float *Elevation, float *MC)
 Surface water flux calculations.
void Flux_SWstuff (int i0, int i1, int j0, int j1, float Flux, float *SURFACE_WAT, double *STUF1, double *STUF2, double *STUF3)
 Flux of surface water constituents.

Variables

FILE * schedFile
FILE * ChanInFile
FILE * CanalOutFile
FILE * CanalOutFile_P
FILE * CanalOutFile_S
FILE * WstructOutFile
FILE * WstructOutFile_P
FILE * WstructOutFile_S
FILE * canDebugFile
FILE * CanalCellInteractDebug
FILE * F_struct_wat
FILE * F_struct_TP
FILE * F_struct_TS
int CHo
int CHANNEL_MAX_ITER
int UTM_EOrigin
int UTM_NOrigin
int printchan = 0
int N_c_iter = 0
int num_cell = 0
int C_Mark = 0
float C_F
float F_ERROR
float MINDEPTH
int num_chan
int num_struct_hist
int numTPser
int numTSser
float * MCopen
char modelFileName [300]
Chan ** Chan_list
Cellscell_last
Cellscell
Structurestruct_first
Structurestructs
HistData Hist_data [MAX_H_STRUCT]
char * ModelPath
char * OutputPath
char * ProjName
char modelName [20]
char modelVers [10]
char SimAlt [20]
char SimModif [20]
double Jdate_init
double Jdate_end
int PORnumday
float can_Intvl
int canPrint
unsigned char * HAB
int * basn
float GP_DetentZ
float GP_MinCheck
float GP_SLRise
float GP_mannDepthPow
float GP_mannHeadPow
float GP_calibGWat
int numBasn
basnDef ** basn_list
basnDefbasins
double * TOT_VOL_CAN
double * VOL_IN_STR
double * VOL_IN_OVL
double * VOL_IN_SPG
double * VOL_IN_GW
double * VOL_OUT_STR
double * VOL_OUT_OVL
double * VOL_OUT_SPG
double * VOL_OUT_GW
double * TOT_P_CAN
double * P_IN_STR
double * P_IN_OVL
double * P_IN_SPG
double * P_IN_GW
double * P_OUT_STR
double * P_OUT_OVL
double * P_OUT_SPG
double * P_OUT_GW
double * TOT_S_CAN
double * S_IN_STR
double * S_IN_OVL
double * S_IN_SPG
double * S_IN_GW
double * S_OUT_STR
double * S_OUT_OVL
double * S_OUT_SPG
double * S_OUT_GW


Detailed Description

Header file for the Water Managment calculations.

This defines or declares variables & functions that are global to WatMgmt.c.

Note: documented with Doxygen, which expects specific syntax within special comments.

The Everglades Landscape Model (ELM).
last updated: Jan 2005

Definition in file watmgmt.h.


Define Documentation

#define STEP  )     ( ((x)>=0) ? 1.0 : 0.0 )
 

Definition at line 23 of file watmgmt.h.

Referenced by Channel_configure().

#define MAX_H_STRUCT   180
 

maximum number of water control structures that have historical or other data as input

Definition at line 24 of file watmgmt.h.

#define LEFT   1
 

Definition at line 26 of file watmgmt.h.

Referenced by Channel_configure().

#define RIGHT   -1
 

Definition at line 27 of file watmgmt.h.

Referenced by Channel_configure().

#define NONE   0
 

Definition at line 29 of file watmgmt.h.

Referenced by Channel_configure().

#define EAST   1
 

Definition at line 30 of file watmgmt.h.

Referenced by Channel_configure().

#define SOUTH   2
 

Definition at line 31 of file watmgmt.h.

Referenced by Channel_configure().

#define ALL   3
 

Definition at line 32 of file watmgmt.h.

Referenced by Channel_configure().


Function Documentation

void CanalReInit  ) 
 

Re-initialize the depths and concentrations in canals under the Positional Analysis mode.

Definition at line 268 of file WatMgmt.c.

00269 {
00270     int i;
00271     for( i = 0; i < num_chan; i++ )
00272 
00273     {
00274         Chan_list[i]->wat_depth = Chan_list[i]->ic_depth;
00275         Chan_list[i]->P = Chan_list[i]->ic_P_con * Chan_list[i]->area * Chan_list[i]->ic_depth; /* kg/m3 * m2 * m */
00276         Chan_list[i]->N = Chan_list[i]->ic_N_con * Chan_list[i]->area * Chan_list[i]->ic_depth; /* kg/m3 * m2 * m */
00277         Chan_list[i]->S = Chan_list[i]->ic_S_con * Chan_list[i]->area * Chan_list[i]->ic_depth; /* kg/m3 * m2 * m */
00278          
00279     }
00280 }

void Canal_Network_Init float  baseDatum,
float *  elev
 

Initialize the water managment network topology.

Initialize/configure the network of canals/levees and and water control structures. The canal information stored in the global array Chan_list

Parameters:
baseDatum GP_DATUM_DISTANCE
elev SED_ELEV

Definition at line 72 of file WatMgmt.c.

00074 { struct Structure *structs;
00075   int i,j;
00076   char filename[30];
00077 
00078   sprintf(filename,"CanalData");
00079 
00080 /* Input the canals geometry file, return pointer to first vector               */
00081 /* Configure the canals, marking the cells that interact with them              */ 
00082 /* and storing the pointers to the cells for each of the canals in Chan_list array */
00083   Channel_configure(elev, ReadChanStruct (filename));  
00084 
00085   for( i = 0; i < num_chan; i++ )
00086   {
00087       Chan_list[i]->P = Chan_list[i]->ic_P_con * Chan_list[i]->area * Chan_list[i]->ic_depth; /* kg/m3 * m2 * m */
00088       Chan_list[i]->N = Chan_list[i]->ic_N_con * Chan_list[i]->area * Chan_list[i]->ic_depth; /* kg/m3 * m2 * m */
00089       Chan_list[i]->S = Chan_list[i]->ic_S_con * Chan_list[i]->area * Chan_list[i]->ic_depth; /* kg/m3 * m2 * m */
00090   }
00091 
00092   usrErr ("\tCanal geometry configured OK");
00093   
00094 /* Read data on the water control structures, return pointer to the first structure */
00095   ReadStructures  (filename, baseDatum);
00096 
00097   ON_MAP[T(2,2)] = 0;
00098   
00099   MCopen = (float *) nalloc(sizeof(float)*(s0+2)*(s1+2),"MCopen"); /* open-water manning's coef, used in structure cell-cell bridge flow */
00100   init_pvar(MCopen,NULL,'f',0.05); /* this is a hard-wired parameter - but do not do anything much with it now for the bridge flows */
00101   
00102  
00103  if ( debug > 0 )       /* Output the canal/levee modifications to ON_MAP */
00104  {
00105         sprintf( modelFileName, "%s/%s/Output/Debug/ON_MAP_CANAL.txt", OutputPath, ProjName );
00106    
00107 
00108   if ( ( ChanInFile = fopen( modelFileName, "w" ) ) ==  NULL )
00109   {
00110      printf( "Can't open the %s canal/levee debug file!\n", modelFileName );
00111      exit(-1) ;
00112   }
00113   fprintf ( ChanInFile, "ROWS=%d\nCOLUMNS=%d\nFORMAT=text\nINFO=\"Debug data: Overland flow restrictions due to levees. \"\nMISSING=0\n", s0, s1 );
00114 
00115   for ( i = 1; i <= s0 ; i++ ) 
00116     { for ( j = 1 ; j <= s1 ; j++ )
00117          fprintf( ChanInFile, "%4d\t", ON_MAP[T(i,j)]) ;
00118       fprintf ( ChanInFile, "\n" );
00119     }   
00120   fclose( ChanInFile ) ;
00121  }
00122 
00123 
00124 /* Open file for writing structure flows */
00125     { if(H_OPSYS==UNIX) 
00126                 sprintf( modelFileName, "%s/%s/Output/Canal/structsOut", OutputPath, ProjName );
00127       else sprintf( modelFileName, "%s%s:Output:structsOut", OutputPath, ProjName );
00128     }
00129 
00130         /* Open file to print structure flows to */
00131     if ( ( WstructOutFile = fopen( modelFileName, "w" ) ) ==  NULL )
00132         {
00133        printf( "Can't open %s file!\n", modelFileName );
00134        exit(-1) ;
00135     }
00136 /* Print The Header Line For This File */
00137     fprintf ( WstructOutFile,
00138               "%s %s %s %s scenario: Sum of water flows at structures (ft^3/sec)\n     Date\t", &modelName, &modelVers, &SimAlt, &SimModif );
00139     structs = struct_first;
00140     while ( structs != NULL )
00141     {   fprintf( WstructOutFile, "%7s\t", structs->S_nam);       
00142             structs = structs->next_in_list;
00143     }
00144 
00145 /* Open file for writing structure P concs */
00146     { if(H_OPSYS==UNIX) 
00147                 sprintf( modelFileName, "%s/%s/Output/Canal/structsOut_P", OutputPath, ProjName );
00148       else sprintf( modelFileName, "%s%s:Output:structsOut_P", OutputPath, ProjName );
00149     }
00150 
00151         /*  Open file to print structure flows to */
00152     if ( ( WstructOutFile_P = fopen( modelFileName, "w" ) ) ==  NULL )
00153         {
00154        printf( "Can't open %s file!\n", modelFileName );
00155        exit(-1) ;
00156     }
00157 /* Print The Header Line For This File */
00158     fprintf ( WstructOutFile_P,
00159               "%s %s %s %s scenario: Flow weighted mean TP conc at structures ( mg/L) \n     Date\t", &modelName, &modelVers, &SimAlt, &SimModif );
00160     structs = struct_first;
00161     while ( structs != NULL ) 
00162     {   fprintf( WstructOutFile_P, "%7s\t", structs->S_nam);      
00163             structs = structs->next_in_list; 
00164     }
00165 
00166 /* Open file for writing structure Salin concs */
00167     { if(H_OPSYS==UNIX) 
00168                 sprintf( modelFileName, "%s/%s/Output/Canal/structsOut_S", OutputPath, ProjName );
00169       else sprintf( modelFileName, "%s%s:Output:structsOut_S", OutputPath, ProjName );
00170     }
00171 
00172         /*  Open file to print structure flows to */
00173     if ( ( WstructOutFile_S = fopen( modelFileName, "w" ) ) ==  NULL )
00174         {
00175        printf( "Can't open %s file!\n", modelFileName );
00176        exit(-1) ;
00177     }
00178 /* Print The Header Line For This File */
00179     fprintf ( WstructOutFile_S,
00180               "%s %s %s %s scenario: Flow weighted mean tracer concentration at structures ( g/L) \n     Date\t", &modelName, &modelVers, &SimAlt, &SimModif );
00181     structs = struct_first;
00182     while ( structs != NULL ) 
00183     {   fprintf( WstructOutFile_S, "%7s\t", structs->S_nam);      
00184             structs = structs->next_in_list; 
00185     }
00186 
00187 /****/
00188     
00189         /* Open file to print canal stage info to */
00190     { if(H_OPSYS==UNIX) 
00191                 sprintf( modelFileName, "%s/%s/Output/Canal/CanalOut", OutputPath, ProjName );
00192       else sprintf( modelFileName, "%s%s:Output:CanalOut", OutputPath, ProjName );
00193     }
00194 
00195     if ( ( CanalOutFile = fopen( modelFileName, "w" ) ) ==  NULL )
00196         {
00197        printf( "Can't open %s file!\n", modelFileName );
00198        exit(-1) ;
00199     }
00200 
00201     fprintf ( CanalOutFile, "%s %s %s %s scenario: Instantaneous water depth (meters, from bottom of canal) in canal Reaches\n    Date\t", &modelName, &modelVers, &SimAlt, &SimModif ); /* print the header line for this file */
00202     /*  channels with negative widths are reserved for new canals */
00203     /* v2.2 put a "R_" prefix in front of the canal ID# */
00204     for( i = 0; i < num_chan; i++ ) {if (Chan_list[i]->width > 0) fprintf ( CanalOutFile, "R_%d\t",Chan_list[i]->number);}
00205     
00206 
00207         /* Open file to print canal phosph info to */
00208     { if(H_OPSYS==UNIX) 
00209                 sprintf( modelFileName, "%s/%s/Output/Canal/CanalOut_P", OutputPath, ProjName );
00210       else sprintf( modelFileName, "%s%s:Output:CanalOut_P", OutputPath, ProjName );
00211     }
00212 
00213     if ( ( CanalOutFile_P = fopen( modelFileName, "w" ) ) ==  NULL )
00214         {
00215        printf( "Can't open %s file!\n", modelFileName );
00216        exit(-1) ;
00217     }
00218 
00219     fprintf ( CanalOutFile_P, "%s %s %s %s scenario: Instantaneous TP conc (mg/L) in canal Reaches\n    Date\t", &modelName, &modelVers, &SimAlt, &SimModif ); /* print the header line for this file */
00220     /* v2.2 put a "R_" prefix in front of the canal ID# */
00221     for( i = 0; i < num_chan; i++ ) {if (Chan_list[i]->width > 0) fprintf ( CanalOutFile_P, "R_%d\t",Chan_list[i]->number);}
00222     
00223         /* Open file to print canal salinity/tracer info to */
00224     { if(H_OPSYS==UNIX) 
00225                 sprintf( modelFileName, "%s/%s/Output/Canal/CanalOut_S", OutputPath, ProjName );
00226       else sprintf( modelFileName, "%s%s:Output:CanalOut_S", OutputPath, ProjName );
00227     }
00228 
00229     if ( ( CanalOutFile_S = fopen( modelFileName, "w" ) ) ==  NULL )
00230         {
00231        printf( "Can't open %s file!\n", modelFileName );
00232        exit(-1) ;
00233     }
00234 
00235     fprintf ( CanalOutFile_S, "%s %s %s %s scenario: Instantaneous tracer conc (g/L) in canal Reaches\n    Date\t", &modelName, &modelVers, &SimAlt, &SimModif ); /* print the header line for this file */
00236     /* v2.2 put a "R_" prefix in front of the canal ID# */
00237     for( i = 0; i < num_chan; i++ ) {if (Chan_list[i]->width > 0) fprintf ( CanalOutFile_S, "R_%d\t",Chan_list[i]->number);}
00238     
00239  
00240 
00241 /* Open file for canal-cell flux debugging purposes */
00242   if( debug > 0 ) 
00243   { 
00244      if(H_OPSYS==UNIX) 
00245                 sprintf( modelFileName, "%s/%s/Output/Debug/Can_debug", OutputPath, ProjName );
00246       else sprintf( modelFileName, "%s%s:Output:Can_debug", OutputPath, ProjName );
00247     
00248 
00249         /* Open file for debug */
00250     if ( ( canDebugFile = fopen( modelFileName, "w" ) ) ==  NULL )
00251         {
00252        printf( "Can't open %s file!\n", modelFileName );
00253        exit(-1) ;
00254     }
00255   
00256      fprintf ( canDebugFile, "Depending on level of the debug parameter, data on canal-cell flows printed here. \n" ); /* v2.5 made use of this file pointer */
00257      fflush(canDebugFile); 
00258     
00259   }
00260   
00261   return;
00262 
00263 }

void Run_Canal_Network float *  SWH,
float *  Elevation,
float *  MC,
float *  GW,
float *  poros,
float *  GWcond,
double *  NA,
double *  PA,
double *  SA,
double *  GNA,
double *  GPA,
double *  GSA,
float *  Unsat,
float *  sp_yield
 

Runs the water management network.

Invoke the calls to water control structure flows and the canal-cell interactions.

Parameters:
SWH SURFACE_WAT
Elevation SED_ELEV
MC HYD_MANNINGS_N
GW SAT_WATER
poros HP_HYD_POROSITY
GWcond HYD_RCCONDUCT
NA DINdummy
PA TP_SF_WT
SA SALT_SURF_WT
GNA DINdummy
GPA TP_SED_WT
GSA SALT_SED_WT
Unsat UNSAT_WATER
sp_yield HP_HYD_SPEC_YIELD

Definition at line 306 of file WatMgmt.c.

00310 { int i;
00311  float CH_vol;
00312 
00313  if ( canPrint && (( N_c_iter++ % hyd_iter ) == 0.0) ) {
00314         if (SensiOn) { 
00315             printchan = (SimTime.IsSimulationEnd) ? (1) : (0); /* flag to indicate to print canal/struct data */
00316          }
00317          else printchan = 1;
00318  }
00319  else  printchan = 0;
00320 
00321  if (printchan == 1) {
00322      
00323      fprintf( WstructOutFile, "\n%d/%d/%d\t",SimTime.yr[0],SimTime.mo[0],SimTime.da[0] );
00324      fprintf( WstructOutFile_P, "\n%d/%d/%d\t",SimTime.yr[0],SimTime.mo[0],SimTime.da[0] );
00325 /*      fprintf( WstructOutFile_N, "\n%d/%d/%d\t",SimTime.yr[0],SimTime.mo[0],SimTime.da[0] ); */
00326      fprintf( WstructOutFile_S, "\n%d/%d/%d\t",SimTime.yr[0],SimTime.mo[0],SimTime.da[0] );
00327 
00328      fprintf (CanalOutFile, "\n%d/%d/%d\t",SimTime.yr[0],SimTime.mo[0],SimTime.da[0] );
00329      fprintf (CanalOutFile_P, "\n%d/%d/%d\t",SimTime.yr[0],SimTime.mo[0],SimTime.da[0] );
00330 /*      fprintf (CanalOutFile_N, "\n%d/%d/%d\t",SimTime.yr[0],SimTime.mo[0],SimTime.da[0] ); */
00331      fprintf (CanalOutFile_S, "\n%d/%d/%d\t",SimTime.yr[0],SimTime.mo[0],SimTime.da[0] );
00332  }
00333     
00334   
00335  for ( i = 0; i < num_chan; i++ )
00336  {   /* set the sum of historical/other-data flows to/from a canal per iteration to 0 */        
00337      Chan_list[i]->sumHistIn = Chan_list[i]->sumHistOut = 0.0;
00338      Chan_list[i]->sumRuleIn = Chan_list[i]->sumRuleOut = 0.0;
00339  }
00340 
00341 /* calculate flows through water control structures */
00342 Flows_in_Structures ( SWH, Elevation, MC, NA, PA, SA ); 
00343   
00344 /* output  header for canal-cell flux output canal debug file, at relatively high debug levels (writes a LOT of stuff) */
00345  if( debug > 3 ) 
00346  {
00347      fprintf( canDebugFile, "N%3d  WatDepth     Pconc     Sconc        flux_S        flux_G        flux_L           Qin          Qout #_iter      Area    error \n", N_c_iter );
00348       
00349  }
00350 
00351 /* Flux the water along a canal using the relaxation procedure */
00352   
00353  for( i = 0; i < num_chan; i++ )
00354      if (Chan_list[i]->width > 0) /* a neg width canal is skipped */
00355 
00356      { 
00357          FluxChannel ( i, SWH, Elevation, MC, GW, poros, GWcond, NA, PA, SA, GNA, GPA, GSA, Unsat, sp_yield );
00358              /* calculate total canal volume after canal fluxes */
00359          CH_vol = Chan_list[i]->area*Chan_list[i]->wat_depth; 
00360 
00361              /* sum the storages on iteration that performs budget checks */
00362          if ( !( FMOD(N_c_iter,(1/canstep) )) && (SimTime.IsBudgEnd ) ) { 
00363              TOT_VOL_CAN[Chan_list[i]->basin] += CH_vol;
00364              TOT_VOL_CAN[0] += CH_vol;
00365              TOT_P_CAN[Chan_list[i]->basin] += Chan_list[i]->P;
00366              TOT_P_CAN[0] += Chan_list[i]->P;
00367              TOT_S_CAN[Chan_list[i]->basin] += Chan_list[i]->S;
00368              TOT_S_CAN[0] += Chan_list[i]->S;
00369          }
00370                                                                                                
00371      
00372          if ( printchan )
00373          {
00374                  /* report N&P concentration in mg/L (kg/m3==g/L * 1000 = mg/L) */
00375              Chan_list[i]->P_con = (CH_vol > 0) ? (Chan_list[i]->P/CH_vol*1000.0) : 0.0;
00376 /*              Chan_list[i]->N_con = (CH_vol > 0) ? (Chan_list[i]->N/CH_vol*1000.0) : 0.0; */
00377              Chan_list[i]->S_con = (CH_vol > 0) ? (Chan_list[i]->S/CH_vol       ) : 0.0;
00378  
00379              /* v2.2 increased decimal places in output (from 6.3f to 7.4f) */
00380              fprintf( CanalOutFile, "%6.2f\t", Chan_list[i]->wat_depth );
00381              fprintf( CanalOutFile_P, "%7.4f\t", Chan_list[i]->P_con );
00382 /*          fprintf( CanalOutFile_N, "%7.4f\t", Chan_list[i]->N_con ); */
00383              fprintf( CanalOutFile_S, "%7.4f\t", Chan_list[i]->S_con ); 
00384          }
00385      } 
00386  if ( printchan )
00387  {
00388      fflush(CanalOutFile);
00389      fflush(CanalOutFile_P);
00390 /*       fflush(CanalOutFile_N); */
00391      fflush(CanalOutFile_S);
00392  }
00393   
00394  return;
00395 }

void ReadStructures char *  filename,
float  BASE_DATUM
 

Read attribute data on water control structures.

Parameters:
filename Water control structures attributes file name
BASE_DATUM GP_DATUM_DISTANCE

Variables local to function
disAggNum Total # of disaggregated water control structures
IDhist ID of historal data set
lincnt Flow data input file record # starting at the point where the model reads values
structName Water control structure name in time series file
histTP Input data field that either has TP conc (ug/L=ppb) that is constant across time or has (any) string to indicate the use of TimeSeries data from a datafile
histTS Input data field that either has TS (salt) conc (g/L=ppt) that is constant across time or has (any) string to indicate the use of TimeSeries data from a datafile

Definition at line 582 of file WatMgmt.c.

References HistData::aggCnt, HistData::aggID, HistData::arrayN, HistData::arrayP, HistData::arrayS, HistData::arrayWat, Structure::canal_fr, Structure::canal_to, Structure::cell_i_fr, Structure::cell_i_HW, Structure::cell_i_to, Structure::cell_i_TW, Structure::cell_j_fr, Structure::cell_j_HW, Structure::cell_j_to, Structure::cell_j_TW, Chan_list, ChanInFile, CHo, F_struct_TP, F_struct_TS, F_struct_wat, Structure::flag, HistData::flag, H_OPSYS, Hist_data, Structure::histID, Structure::HW_graph, Structure::HW_stage, isFloat(), isInteger(), julday(), modelFileName, modelName, ModelPath, msgStr, nalloc(), Structure::next_in_list, num_chan, num_struct_hist, numTPser, numTSser, ON_MAP, PORnumday, ProjName, Read_schedule(), s0, Structure::S_nam, Scip(), SimAlt, SimModif, Structure::str_cell_i, Structure::str_cell_j, struct_first, structs, Structure::Sum_N, Structure::Sum_P, Structure::Sum_S, Structure::SumFlow, T, TAB, Structure::TN, Structure::TNser, Structure::TP, Structure::TPser, Structure::TS, Structure::TSser, Structure::TW_graph, Structure::TW_stage, usrErr(), usrErr0(), Structure::w_coef, Chan::width, and WriteMsg().

Referenced by Canal_Network_Init().

00583 {
00593 struct Structure *structs, *struct_next, *struct_last;
00594 
00595  char ss[322], *s;
00596  char sched_name[20];
00597  int i, j, k;
00598  int disAggNum=0; 
00599  int IDhist[MAX_H_STRUCT]; 
00600  int lincnt; 
00601  double Jdate_read;
00602  char lline[2001], *line;
00603  int yyyy, mm, dd;
00604  char structName[20]; 
00605  char scenario[20], s_modname[20];
00606  char histTP[5]; 
00607  char histTS[5]; 
00608  
00609  num_struct_hist = 0;
00610  numTPser = numTSser = 0;
00611  
00612  { if(H_OPSYS==UNIX) 
00613      sprintf( modelFileName, "%s/%s/Data/%s.struct", ModelPath, ProjName, filename );
00614  else sprintf( modelFileName, "%s%s:Data:%s.struct", ModelPath, ProjName, filename );
00615  }
00616 
00617 /* Open file with structures data */
00618  if ( ( ChanInFile = fopen( modelFileName, "r" ) ) ==  NULL )
00619  {
00620      printf( "Can't open %s file!\n", modelFileName ) ;
00621      exit(-1) ;
00622  }
00623     
00624 /* Allocate memory for structures */
00625  if ( (structs = (struct Structure *) malloc( (size_t) sizeof( struct Structure ))) == NULL )
00626  {
00627      printf( "No memory for first structure\n " ) ;
00628      exit( -2 ) ;
00629  }
00630  struct_first = structs;
00631 
00632  fgets( ss, 320, ChanInFile );
00633  sscanf(ss, "%s\t%s\t%s\t%s",&scenario,&scenario,&s_modname,&scenario); /*  2 junk strings followed by the model name and actual scenario ID */
00634  if (strcmp(scenario,SimAlt) != 0) {
00635          sprintf(msgStr, "The scenario (%s) found in the %s file doesn't match the one (%s) you asked for in Driver.parm!",
00636                  scenario, modelFileName, SimAlt); usrErr(msgStr);
00637      exit(-1);
00638   }
00639  if (strcmp(s_modname,modelName) != 0) {
00640          sprintf(msgStr, "The model name (%s) found in the %s file doesn't match the one (%s) in the canal data file!",
00641                  s_modname, modelFileName, modelName); usrErr(msgStr);
00642      exit(-1);
00643   }
00644  
00645 /* Skip 2'nd header line */
00646  fgets( ss, 320, ChanInFile );
00647 
00648 
00649 /* Read structure descriptors until none left */   
00650  while ( fgets( ss, 320, ChanInFile ) != NULL && !feof( ChanInFile ) )
00651  { 
00652      s = ss; 
00653          
00654          
00655      sscanf( s, "%d", &structs->flag ); 
00656 /* structure flag: <0 - skip structure info, structure not operating; */
00657 /*                 >0 - historical/other-data data array, flag=1 normally in this case, but is >1 for aggregated SFWMM structures; */
00658 /*                  0 - rule-based structure   */
00659      if ( structs->flag < 0 ) continue;
00660      if ( structs->flag > 0 && structs->flag <10 ) {
00661          structs->histID = ++num_struct_hist;
00662          Hist_data[structs->histID-1].flag = structs->flag;
00663      }
00664      
00665          
00666 /*    assign integer to the historical data ID, and increment the counter for structures with historical/other-data data to be read from datafile. */
00667 /*    flags >= 20 represent structures (such as the S-11s) that were aggregated in SFWMM output */
00668 /*    and need to be assigned disaggregated flows for separate ELM structures */
00669 /*    Ex: SFWMM structure name S11 may have flag = 2 in CanalData.struct.  That S11 flow holds the sum of S-11A-C flows. */
00670 /*    The flags on S-11A, S-11B, and S-11C would all have been assigned 20, 20, and 20 in the CanalData.struct file */
00671 
00672 /*           NO LONGER NECESSARY: as a temporary measure until the SFWMM output .dss file is fixed, the structures that have calculated "other" partitions */
00673 /*             (S140,S7,S8) are assigned historical flags of 9 (they are in the dss file read by a SAS routine, the flag */
00674 /*             of 9 allows them to remain in the file, but the structure is turned off because it has flag >1)  */
00675          
00676      
00677      s = Scip( s, TAB );
00678          
00679      if ( *s != TAB ) {
00680              sscanf ( s, "%s", &structs->S_nam );
00681              if ( structs->flag >0 && structs->flag <10) 
00682                  sscanf ( s, "%s", &Hist_data[structs->histID-1].S_nam );
00683              
00684                  /* for solute concentration of flows thru structures, */
00685                  /* data use flag XXser: -1== simulated, 0==static historical, 1= time-varying historical */
00686              s = Scip( s, TAB );
00687              if ( *s == TAB ) structs->TPser = -1; 
00688              else {
00689                  sscanf ( s, "%s", &histTP ); /* read the field that either has TP conc (ug/L=ppb) that is constant across time
00690                                                    or has (any) string to indicate the use of TimeSeries data from a datafile (read later)*/
00691                  if (isInteger(histTP) ) {
00692                       structs->TPser = 0; 
00693                       structs->TP = atof(histTP)*1.0e-6;  /* conc (ug/L=>g/L) that is constant across time */
00694                  }
00695                  else { structs->TPser = 1; numTPser++;} /* otherwise, we will be reading time series data later */
00696              }
00697              
00698              s = Scip( s, TAB );
00699              if ( *s == TAB ) structs->TNser = -1; /* value to indicate use of simulated, not historical, data */
00700              else {
00701                  sscanf ( s, "%f", &structs->TN );
00702                  structs->TNser = 0;
00703                  structs->TN *= 1.0e-6;/* read the TN conc (ug/L=>g/L) that is constant across time */
00704              }
00705              
00706 
00707              s = Scip( s, TAB );
00708              if ( *s == TAB ) structs->TSser = -1; /* value to indicate use of simulated, not historical, data */
00709              /* v2.2 added code to allow use of TimeSeries input data in addition to old method of a constant concentration */
00710              else {
00711                   sscanf ( s, "%s", &histTS ); /* read the field that either has TS (salt) conc (g/L=ppt) that is constant across time
00712                                                    or has (any) string to indicate the use of TimeSeries data from a datafile (read later)*/
00713                  if (isFloat(histTS) ) {
00714                       structs->TSser = 0; 
00715                       structs->TS = atof(histTS)*1.0;  /* conc (no conversion, g/L) that is constant across time */
00716                  }
00717                  else { structs->TSser = 1; numTSser++;} /* otherwise, we will be reading time series data later */
00718 
00719              }
00720              
00721              s = Scip( s, TAB );
00722              if ( *s == TAB ) structs->str_cell_i = 0;  /* non 0 if elevation at the structure location is desired */
00723              else {
00724                  sscanf( s, "%d", &structs->str_cell_i );
00725                  structs->str_cell_i++; /* this converts the data map coords to model coords (map+1) */
00726              }
00727              
00728              s = Scip( s, TAB );
00729 
00730              if ( *s == TAB ) structs->str_cell_j = 0;  /* non 0 if elevation at the structure location is desired */
00731              else {
00732                  sscanf( s, "%d", &structs->str_cell_j );
00733                  structs->str_cell_j++;
00734              }
00735 
00736              if ( structs->str_cell_i != 0 && !ON_MAP[T(structs->str_cell_i, structs->str_cell_j)] )
00737              {  printf ( "Location of water control structure %s is off-map.\n", structs->S_nam );
00738              exit (-2);                 /*check that cells obtained are within on-map boundaries */
00739              }
00740              
00741              
00742      }
00743          
00744      s = Scip( s, TAB );
00745 
00746      if ( *s == TAB ) structs->canal_fr  = 0;   /* non 0 if the structure involves a donor canal */
00747      else sscanf( s, "%d", &structs->canal_fr  );
00748      s = Scip( s, TAB );
00749      if ( *s == TAB) structs->canal_to = 0;             /* non 0 if the structure involves a receiving canal */
00750      else sscanf( s, "%d", &structs->canal_to );
00751      s = Scip( s, TAB );
00752 
00753      if ( structs->canal_fr  - CHo >= num_chan || structs->canal_to - CHo >= num_chan )
00754      {  printf ( "Canal from/to doesn't exist for water control structure %s\n", structs->S_nam );
00755      exit (-2);                 /* check that canals read are consistent with canals data */
00756      }           
00757 /*  check to ensure canals associated with structures are valid */
00758      if ( structs->canal_fr != 0 && Chan_list[structs->canal_fr- CHo]->width <= 0.1 )
00759      { printf ( "The donor canal %d has been turned off for water control structure %s\n", structs->canal_fr, structs->S_nam );
00760      exit (-2);                 
00761      }           
00762      if ( structs->canal_to != 0 && Chan_list[structs->canal_to- CHo]->width <= 0.1 )
00763      {  printf ( "The recipient canal %d has been turned off for water control structure %s\n", structs->canal_to, structs->S_nam );
00764      exit (-2);                 
00765      }           
00766          
00767      if (*s == TAB) structs->cell_i_fr = 0;  /* non 0 if the structure involves a donor cell */
00768      else 
00769      { sscanf( s, "%d", &structs->cell_i_fr );
00770      structs->cell_i_fr++;
00771      }   
00772      s = Scip( s, TAB );
00773      if (*s == TAB) structs->cell_j_fr = 0;
00774      else 
00775      { sscanf( s, "%d", &structs->cell_j_fr );
00776      structs->cell_j_fr++;
00777      }
00778      s = Scip( s, TAB );
00779          
00780 
00781      if ( structs->cell_i_fr > 2 && !ON_MAP[T(structs->cell_i_fr, structs->cell_j_fr)] )
00782      {  printf ( "Donor cell for water control structure %s is off-map.\n", structs->S_nam );
00783      exit (-2);                 /*check that cells obtained are within on-map boundaries */
00784      }
00785 
00786          /* if valid donor cell, mark it on ON_MAP */
00787      if ( structs->cell_j_fr ) ON_MAP[T(structs->cell_i_fr, structs->cell_j_fr)] += 100;
00788          
00789      if (*s == TAB) structs->cell_i_to = 0;  /* non 0 if the structure involves a recipient cell */
00790      else 
00791      { sscanf( s, "%d", &structs->cell_i_to );
00792      structs->cell_i_to++;
00793      }   
00794          
00795      s = Scip( s, TAB );
00796      if (*s == TAB) structs->cell_j_to = 0;
00797      else
00798      { sscanf( s, "%d", &structs->cell_j_to );
00799      structs->cell_j_to++;
00800      }   
00801 
00802 
00803      if ( structs->cell_i_to > 2 && !ON_MAP[T(structs->cell_i_to, structs->cell_j_to)] )
00804      {  printf ( "Recipient cell for water control structure %s is off-map.\n", structs->S_nam );
00805      exit (-2);                 /*check that cells obtained are within on-map boundaries */
00806      }
00807 
00808          /* if valid recieving cell, mark it on ON_MAP */
00809          /* this used to check for cell_j_to > 2, not putting recipient LEC cells ON_MAP */
00810      if ( structs->cell_j_to ) ON_MAP[T(structs->cell_i_to, structs->cell_j_to)] += 100;
00811          
00812      s = Scip( s, TAB );                /* read the targeted water levels */
00813      if (*s == TAB) structs->HW_stage = 0;
00814      else if ( isalpha (*s) ) {  /* if there is a graph for the schedule */
00815          sscanf( s, "%s", sched_name );
00816          structs->HW_graph = Read_schedule ( sched_name, filename, BASE_DATUM );
00817          structs->HW_stage = -99;  
00818      } 
00819      else   { sscanf( s, "%f", &structs->HW_stage );
00820      structs->HW_stage += BASE_DATUM;
00821      }
00822                         
00823      s = Scip( s, TAB );
00824      if (*s == TAB) structs->TW_stage = 0;
00825      else if ( isalpha (*s) ) {  /* if there is a graph for the schedule */
00826          sscanf( s, "%s", sched_name );
00827          structs->TW_graph = Read_schedule ( sched_name, filename, BASE_DATUM );
00828          structs->TW_stage = -99;
00829      } 
00830      else       { sscanf( s, "%f", &structs->TW_stage );
00831      structs->TW_stage -= BASE_DATUM; 
00832      }
00833      s = Scip( s, TAB );
00834          
00835      if (*s == TAB) structs->cell_i_HW = 0;
00836      else               /* read coordinates of the reference cell for head water */
00837      { sscanf( s, "%d", &structs->cell_i_HW );
00838      structs->cell_i_HW++;
00839      }   
00840      s = Scip( s, TAB );
00841      if (*s == TAB) structs->cell_j_HW = 0;
00842      else 
00843      { sscanf( s, "%d", &structs->cell_j_HW );
00844      structs->cell_j_HW++;
00845      }   
00846      s = Scip( s, TAB );
00847 
00848      if ( structs->cell_i_HW > s0 || structs->cell_j_HW > s1 )
00849      {  printf ( "Error in definition of water control structure %s", structs->S_nam );
00850      exit (-2); /*check that cell obtained is within array boundaries */
00851      }
00852 
00853      if (*s == TAB) structs->cell_i_TW = 0;
00854      else               /* read coordinates of the reference cell for tail water */
00855      { sscanf( s, "%d", &structs->cell_i_TW );
00856      structs->cell_i_TW++;
00857      }   
00858      s = Scip( s, TAB );
00859      if (*s == TAB) structs->cell_j_TW = 0;
00860      else 
00861      { sscanf( s, "%d", &structs->cell_j_TW );
00862      structs->cell_j_TW++;
00863      }   
00864      s = Scip( s, TAB );
00865 
00866      if ( structs->cell_i_TW > s0 || structs->cell_j_TW > s1 )
00867      {  printf ( "Error in definition of water control structure %s", structs->S_nam );
00868      exit (-2); /*check that cell obtained is within array boundaries */
00869      }
00870 
00871      if (*s == TAB || *s == '?') structs->w_coef = 0;
00872      else sscanf( s, "%f", &structs->w_coef );
00873          /* the next field after the coef is the date the file was exported, ignored here */
00874         
00875      struct_last = structs;
00876           
00877          /* Allocate memory for next structure */
00878      if ( (struct_next = (struct Structure *) malloc( (size_t) sizeof( struct Structure ))) == NULL )
00879      {
00880          printf( "Failed to allocate memory for next structure (%s)\n ", structs->S_nam ) ;
00881          exit( -2 ) ;
00882      };
00883      
00884      /* initialize the sum of structure  flows and constit solute in flows */
00885      structs->SumFlow = 0.0; 
00886      structs->Sum_P = 0.0; 
00887      structs->Sum_N = 0.0; 
00888      structs->Sum_S = 0.0; 
00889 
00890      structs->next_in_list = struct_next;
00891      structs = struct_next;
00892  }
00893    
00894  struct_last->next_in_list = NULL;
00895  free ((char *)structs);
00896  fclose (ChanInFile);
00897    
00898 
00899  usrErr( "\tWater control structure locs/attributes read OK" );
00900  usrErr0 ( "\tControl structures' water flow time series...");
00901  
00902    
00903 /**********************/
00904 /* read structure flow input file(s) */
00905      { if(H_OPSYS==UNIX) 
00906          sprintf( modelFileName, "%s/%s/Data/%s.struct_wat", ModelPath, ProjName, filename );
00907      else sprintf( modelFileName, "%s%s:Data:%s.struct_wat", ModelPath, ProjName, filename );
00908      }
00909 
00910 /* Open file with structure flow data */
00911      if ( ( F_struct_wat = fopen( modelFileName, "r" ) ) ==  NULL )
00912      {
00913          printf( "Can't open %s file!\n", modelFileName ) ;
00914          exit(-1) ;
00915      }
00916 
00917      fgets(lline, 2000, F_struct_wat); /* first header line with the alternative's name */
00918      line=lline;
00919      sscanf(line, "%s",&scenario);
00920      if (strcmp(scenario,SimAlt) != 0) {
00921          sprintf(msgStr, "The scenario (%s) found in the %s file doesn't match the one (%s) you asked for in Driver.parm!", scenario, modelFileName, SimAlt); usrErr(msgStr);
00922          exit(-1);
00923      }
00924 
00925      fgets(lline, 2000, F_struct_wat); /* read 2'nd header line with column names */
00926      line=lline;
00927      line = Scip( line, TAB );
00928      line = Scip( line, TAB );
00929      line = Scip( line, TAB ); /*read past the three date columns, ready to read the variable names */
00930 
00931 
00932      /* loop through all the control structure columns */
00933      for ( i = 0; i < num_struct_hist; i++ ) { 
00934 
00935          Hist_data[i].arrayWat = (float*) nalloc( ( (PORnumday+2)*sizeof(float) ), "timeseries temp" );
00936          Hist_data[i].arrayP =   (float*) nalloc( ( (PORnumday+2)*sizeof(float) ), "timeseries temp" );
00937           Hist_data[i].arrayN =  (float*) nalloc( ( (PORnumday+2)*sizeof(float) ), "timeseries temp" ); 
00938           Hist_data[i].arrayS =  (float*) nalloc( ( (PORnumday+2)*sizeof(float) ), "timeseries temp" ); 
00939 
00940 /*  allocate memory for dissaggregated structures and provide a couple of attributes */
00941          if ( Hist_data[i].flag > 1 ) {
00942              structs = struct_first;  /* go through all structures to pull out the disaggregated ones */
00943              while ( structs != NULL ) 
00944              {
00945                  /* Ex: the SFWMM aggregated S10's flag is 2, S10AB&C's flags are 20 */
00946                  if ( structs->flag == 10 * Hist_data[i].flag  )  
00947                  { 
00948                         /* disAggNum counts the total number of disaggregated structures for simulation */
00949                      disAggNum++;
00950                      Hist_data[num_struct_hist+disAggNum-1].arrayWat = 
00951                              (float*) nalloc( ( (PORnumday+2)*sizeof(float) ), "timeseries temp" );
00952                      Hist_data[num_struct_hist+disAggNum-1].arrayP    = 
00953                              (float*) nalloc( ( (PORnumday+2)*sizeof(float) ), "timeseries temp" );
00954                      Hist_data[num_struct_hist+disAggNum-1].arrayN    = 
00955                              (float*) nalloc( ( (PORnumday+2)*sizeof(float) ), "timeseries temp" ); 
00956                      Hist_data[num_struct_hist+disAggNum-1].arrayS    = 
00957                              (float*) nalloc( ( (PORnumday+2)*sizeof(float) ), "timeseries temp" ); 
00958 
00959                      Hist_data[num_struct_hist+disAggNum-1].aggID = i; /* the aggregated structure ID is attached to the dissagregated one */
00960                      strcpy( Hist_data[num_struct_hist+disAggNum-1].S_nam, structs->S_nam ); /* the disaggregated structure's name is put into the HistData structure */
00961                      Hist_data[i].aggCnt++; /* increment a count of the number of disagg structs associated w/ each agg struct */
00962                      structs->flag = 1; /* set the disagg struct's processing flag to 1, the only flag used for simulation structure flows (new input style) */
00963                      structs->histID = num_struct_hist+disAggNum; /* store the historical ID value for the disagg structures */
00964                      
00965                  } 
00966    
00967                  structs = structs->next_in_list;
00968              }             
00969          }
00970          
00971          sscanf ( line,"%s\t",&structName); /* read the structure names in the header */
00972          if (strcmp(Hist_data[i].S_nam , structName) != 0) {
00973              printf( "variable name  %s in input file does not match order of CanalData.structs.\n", structName);
00974              exit(-1); 
00975              }
00976          line = Scip( line, TAB );
00977      }
00978 
00979 
00980      lincnt = 0; /*  historical/other-data data input file record #  starting at the point where the model reads values */
00981      while ( fgets(lline, 2000, F_struct_wat) != NULL && !feof( F_struct_wat ) )
00982      {
00983          line = lline;
00984          sscanf( line, "%d %d %d",&yyyy,&mm,&dd); /* read the date of the current record */
00985          
00986         /* julian day, returns a double (args = mon, day, year, hours, minutes, seconds) */
00987          Jdate_read = julday( mm,dd,yyyy,0,0,0.0 ); 
00988          if ((lincnt == 0) && (Jdate_read > Jdate_init) ) {
00989              printf (" \n***Error\nStructure flow data (%s) init date > simulation start date\n", modelFileName); exit (-1);
00990          }
00991          if ((Jdate_read >= Jdate_init) && (Jdate_read <= Jdate_end))
00992          {
00993              line = Scip( line, TAB );
00994              line = Scip( line, TAB );
00995 
00996              for ( i = 0; i < num_struct_hist; i++ ) { /* loop through all the control structure columns in the input file */
00997                  line = Scip( line, TAB );
00998                  if (*line == TAB) Hist_data[i].arrayWat[((int)(lincnt))] = 0;
00999                  else sscanf( line, "%f", &Hist_data[i].arrayWat[((int)(lincnt))] );
01000                  Hist_data[i].arrayWat[((int)(lincnt))] *= 2446.848; /* cfs => m^3/d (0.02832*60*60*24)  */
01001              }
01002          
01003          for (i = num_struct_hist; i < num_struct_hist+disAggNum; i++) { /* dissaggregate all agg structs */
01004              k = Hist_data[i].aggID; /* use the aggregated historical data WCstructure */
01005              Hist_data[i].arrayWat[((int)(lincnt))] =
01006                  Hist_data[k].arrayWat[((int)(lincnt))] / (Hist_data[k].aggCnt) ;
01007 
01008          }
01009          
01010               
01011          lincnt++; /* increment the number of daily data records */
01012      }
01013 
01014 } /* end of reading structure flow records */
01015      
01016 if (lincnt == 0) {
01017     printf (" \n***Error\nStructure flow data (%s) has date problem (0 records read)\n", modelFileName); exit (-1);
01018 }
01019 else if (lincnt != PORnumday ) {
01020     printf (" \n***Error\nStructure flow data (%s) has problem (%d records read, %d expected)\n", modelFileName, lincnt, PORnumday); exit (-1);
01021 }
01022      
01023 sprintf (msgStr, "done, found %d records.",lincnt); usrErr(msgStr);
01024 fclose (F_struct_wat);
01025 
01026 /**********************/
01027 /* read time series structure P concentration input file(s) */
01028 /* there can be situations when we have no time series TP data (numTPser==0) */
01029  if (numTPser>0) {
01030      usrErr0 ( "\tControl structures' TP concen. time series...");
01031 
01032      { if(H_OPSYS==UNIX) 
01033          sprintf( modelFileName, "%s/%s/Data/%s.struct_TP", ModelPath, ProjName, filename );
01034      else sprintf( modelFileName, "%s%s:Data:%s.struct_TP", ModelPath, ProjName, filename );
01035      }
01036 
01037 /* Open file with structure concentration data */
01038      if ( ( F_struct_TP = fopen( modelFileName, "r" ) ) ==  NULL )
01039      {
01040          printf( "Can't open %s file!\n", modelFileName ) ;
01041          exit(-1) ;
01042      }
01043 
01044      fgets(lline, 2000, F_struct_TP); /* first header line with the alternative's name */
01045      line=lline;
01046      sscanf(line, "%s",&scenario);
01047      if (strcmp(scenario,SimAlt) != 0) {
01048          sprintf(msgStr, "The scenario (%s) found in the %s file doesn't match the one (%s) you asked for in Driver.parm!", scenario, modelFileName, SimAlt); usrErr(msgStr);
01049          exit(-1);
01050      }
01051      fgets(lline, 2000, F_struct_TP); /* read 2'nd header line with column names */
01052      line=lline;
01053      line = Scip( line, TAB );
01054      line = Scip( line, TAB );
01055      line = Scip( line, TAB );
01056 
01057 /* now we're past the three date columns, ready to read the variable names */
01058      for ( j = 0; j < numTPser; j++ ) { /* loop through all the data columns */
01059          IDhist[j] = -1;
01060          strcpy(structName,"NULL"); /* a string structName read previously is retained - this NULL assignment will take care of case
01061                                    where nothing was found/read in the time series file, and thus cause an error */
01062          sscanf ( line,"%s\t",&structName); /* read the structure names in the header */
01063          for ( i = 0; i < num_struct_hist+disAggNum; i++ ) { 
01064              if (strcmp(Hist_data[i].S_nam , structName) == 0) {
01065                  IDhist[j] = i;
01066                  line = Scip( line, TAB );
01067              }
01068          }
01069          if (IDhist[j] == -1 ) { printf( "variable name %s in %s not found in CanalData.structs. (if 'NULL', a variable in CanalData.structs defined as having a time series column was not found in %s)\n", structName, modelFileName, modelFileName);
01070          exit(-1); 
01071          }
01072      }
01073 
01074      lincnt = 0; /*  input file record #  starting at the point where the model reads values */
01075 
01076      while ( fgets(lline, 2000, F_struct_TP) != NULL && !feof( F_struct_TP ) )
01077      {
01078          line = lline;
01079          sscanf( line, "%d %d %d",&yyyy,&mm,&dd); /* read the date of the current record */
01080          
01081         /* julian day, returns a double (args = mon, day, year, hours, minutes, seconds) */
01082          Jdate_read = julday( mm,dd,yyyy,0,0,0.0 ); 
01083          if ((lincnt == 0) && (Jdate_read > Jdate_init) ) {
01084              printf (" \n***Error\nStructure TP data (%s) init date > simulation start date\n", modelFileName); exit (-1);
01085          }
01086          if ((Jdate_read >= Jdate_init) && (Jdate_read <= Jdate_end))
01087          {
01088              line = Scip( line, TAB );
01089              line = Scip( line, TAB );
01090              for ( j = 0; j < numTPser; j++ ) { /* loop through all the control structure columns in the input file */
01091                  line = Scip( line, TAB );
01092                  if (*line == TAB) Hist_data[IDhist[j]].arrayP[((int)(lincnt))] = 0;
01093                  else sscanf( line, "%f", &Hist_data[IDhist[j]].arrayP[((int)(lincnt))] );
01094                  Hist_data[IDhist[j]].arrayP[((int)(lincnt))] *= 1.0E-6; /* data read in ug/L, convert to g/L=kg/m3 */
01095                  
01096              }
01097               
01098              lincnt++; /* increment the number of daily data records */
01099          }
01100 
01101      } /* end of reading TP conc data records */
01102      
01103      if (lincnt == 0 ) {
01104              printf (" \n***Error\nTP Time series data (%s) has date problem (0 records read)\n", modelFileName); exit (-1);
01105      }
01106      else if (lincnt != PORnumday ) {
01107              printf (" \n***Error\nTP Time series data (%s) has problem (%d records read, %d expected)\n", modelFileName, lincnt, PORnumday); exit (-1);
01108      }
01109      sprintf (msgStr, "done, found %d records.",lincnt); usrErr(msgStr);
01110      
01111      fclose (F_struct_TP);
01112  } /* end of reading TP time series data (if needed) */
01113  
01114  
01115  
01116 /**********************/
01117 /* read time series structure tracer (salt) concentration input file(s) */
01118 /* there can be situations when we have no time series TS (salt) data (numTSser==0) */
01119  if (numTSser>0) {
01120      usrErr0 ( "\tControl structures' Tracer concen. time series...");
01121 
01122      { if(H_OPSYS==UNIX) 
01123          sprintf( modelFileName, "%s/%s/Data/%s.struct_TS", ModelPath, ProjName, filename );
01124      else sprintf( modelFileName, "%s%s:Data:%s.struct_TS", ModelPath, ProjName, filename );
01125      }
01126 
01127 /* Open file with structure concentration data */
01128      if ( ( F_struct_TS = fopen( modelFileName, "r" ) ) ==  NULL )
01129      {
01130          printf( "Can't open %s file!\n", modelFileName ) ;
01131          exit(-1) ;
01132      }
01133 
01134      fgets(lline, 2000, F_struct_TS); /* first header line with the alternative's name */
01135      line=lline;
01136      sscanf(line, "%s",&scenario);
01137      if (strcmp(scenario,SimAlt) != 0) {
01138          sprintf(msgStr, "The scenario (%s) found in the %s file doesn't match the one (%s) you asked for in Driver.parm!", scenario, modelFileName, SimAlt); usrErr(msgStr);
01139          exit(-1);
01140      }
01141      fgets(lline, 2000, F_struct_TS); /* read 2'nd header line with column names */
01142      line=lline;
01143      line = Scip( line, TAB );
01144      line = Scip( line, TAB );
01145      line = Scip( line, TAB );
01146 
01147 /* now we're past the three date columns, ready to read the variable names */
01148      for ( j = 0; j < numTSser; j++ ) { /* loop through all the data columns */
01149          IDhist[j] = -1;
01150          strcpy(structName,"NULL"); /* a string structName read previously is retained - this NULL assignment will take care of case
01151                                    where nothing was found/read in the time series file, and thus cause an error */
01152          sscanf ( line,"%s\t",&structName); /* read the structure names in the header */
01153          for ( i = 0; i < num_struct_hist+disAggNum; i++ ) { 
01154              if (strcmp(Hist_data[i].S_nam , structName) == 0) {
01155                  IDhist[j] = i;
01156                  line = Scip( line, TAB );
01157              }
01158          }
01159          if (IDhist[j] == -1 ) { printf( "variable name %s in %s not found in CanalData.structs. (if 'NULL', a variable in CanalData.structs defined as having a time series column was not found in %s)\n", structName, modelFileName, modelFileName);
01160          exit(-1); 
01161          }
01162      }
01163 
01164      lincnt = 0; /*  input file record #  starting at the point where the model reads values */
01165 
01166      while ( fgets(lline, 2000, F_struct_TS) != NULL && !feof( F_struct_TS ) )
01167      {
01168          line = lline;
01169          sscanf( line, "%d %d %d",&yyyy,&mm,&dd); /* read the date of the current record */
01170          
01171         /* julian day, returns a double (args = mon, day, year, hours, minutes, seconds) */
01172          Jdate_read = julday( mm,dd,yyyy,0,0,0.0 ); 
01173          if ((lincnt == 0) && (Jdate_read > Jdate_init) ) {
01174              printf (" \n***Error\nStructure Tracer data (%s) init date > simulation start date\n", modelFileName); exit (-1);
01175          }
01176          if ((Jdate_read >= Jdate_init) && (Jdate_read <= Jdate_end))
01177          {
01178              line = Scip( line, TAB );
01179              line = Scip( line, TAB );
01180              for ( j = 0; j < numTSser; j++ ) { /* loop through all the control structure columns in the input file */
01181                  line = Scip( line, TAB );
01182                  if (*line == TAB) Hist_data[IDhist[j]].arrayS[((int)(lincnt))] = 0;
01183                  else sscanf( line, "%f", &Hist_data[IDhist[j]].arrayS[((int)(lincnt))] );
01184                  Hist_data[IDhist[j]].arrayS[((int)(lincnt))] *= 1.0; /* data read in g/L=kg/m3, no conversion*/
01185                  
01186              }
01187               
01188              lincnt++; /* increment the number of daily data records */
01189          }
01190 
01191      } /* end of reading tracer (salt) conc data records */
01192      
01193      if (lincnt == 0 ) {
01194              printf (" \n***Error\nTracer Time series data (%s) has date problem (0 records read)\n", modelFileName); exit (-1);
01195      }
01196      else if (lincnt != PORnumday ) {
01197              printf (" \n***Error\nTracer Time series data (%s) has problem (%d records read, %d expected)\n", modelFileName, lincnt, PORnumday); exit (-1);
01198      }
01199      sprintf (msgStr, "done, found %d records.",lincnt); usrErr(msgStr);
01200      
01201      fclose (F_struct_TS);
01202  } /* end of reading trace (salt) time series data (if needed) */
01203  
01204 
01205  
01206  sprintf(msgStr, "\t**** Evaluate scenario: %s %s ****", &SimAlt, &SimModif ); usrErr(msgStr); WriteMsg(msgStr,1);
01207 
01208  return;
01209 }

Here is the call graph for this function:

void Channel_configure float *  ElevMap,
struct Chan channel_first
 

Configure attributes of grid cells interacting with canal vectors.

Identifies the cells that are crossed by the canal, calculates the length of canal within each cell and marks the cells directly interacting with the canal. Store the pointers to the adjacent cell arrays.

Parameters:
ElevMap SED_ELEV
channel_first Pointer to first canal data structure

Definition at line 1317 of file WatMgmt.c.

References Abs, ALL, Chan::area, C_Mark, CanalCellInteractDebug, cell, cell_last, Chan::cells, Chan_list, CHo, Chan::cond, EAST, Chan::elev_drop, Chan::elev_end, Chan::elev_slope, Chan::elev_start, getCanalElev(), init_pvar(), LEFT, Chan::length, Chan::levee, MarkCell(), Chan::minVol, modelFileName, nalloc(), Cells::next_cell, Chan::next_in_list, Chan_reach::next_reach, NONE, num_cell, num_chan, Chan::num_of_cells, Chan::number, ON_MAP, OutputPath, ProjName, Chan::reaches, RIGHT, Chan::roil, Chan::roir, Round(), s0, s1, sec_per_day, Chan::seg_len, SOUTH, Chan::SPG_flow_coef, STEP, Chan::SW_flow_coef, T, Chan::width, Chan_reach::x0, Chan_reach::x1, Chan_reach::y0, and Chan_reach::y1.

Referenced by Canal_Network_Init().

01318 { int i, j;                                                       /*coord. of current cell crossed by the canal */
01319  int  cellLoc, cellLoc_i, cellLoc_j;
01320  int HV, k, L_range, R_range;
01321  float T_length;                                                /* length of a straight segment of a canal reach (m) */
01322  float C_length;                                          /*  total length of all of the segments of a canal reach (m) */
01323  float distTot;                                       /*  cumulative (along grid cells) distance along a canal reach, from the starting point (m) */
01324  float A1, B1, A2, B2;                                                   /* parameters of the canal equation */
01325  float L, p_middle, m, i45, x, y, x1, y1, length, RoIL, RoIR;
01326  float cell_step = 1.0;
01327    
01328  float a0, b0, a1, b1;
01329  int c_num, z;
01330  struct Cells  *cell; /* addresses for adjacent cell structures */
01331  struct Chan *channel;
01332  struct Chan_reach *Reach;
01333  
01334  int numChek; /* placeholder for debugging cell elev along canals */
01335    
01336 
01337  int *marked; /* a "temporary" fix used to determine if a cell belonging to a canal has already been marked */
01338  marked = (int *) nalloc(sizeof(int)*(s0+2)*(s1+2),"marked");
01339  init_pvar(marked,NULL,'d',0);
01340 
01341 /* Open file to print canal-cell interaction info to */
01342  sprintf( modelFileName, "%s/%s/Output/Debug/CanalCells_interaction.txt", OutputPath, ProjName );
01343  if ( ( CanalCellInteractDebug = fopen( modelFileName, "w" ) ) ==  NULL )
01344  {
01345      printf( "Can't open %s file!\n", modelFileName ) ;
01346      exit(-1) ;
01347  }
01348 
01349      /* allocate memory for array of pointers to Chan structs */
01350  if ( (Chan_list = 
01351        (struct Chan **) malloc( (size_t) sizeof( struct Chan *) * num_chan)) == NULL )
01352  {
01353      printf( "Failed to allocate memory for next canal (%d)\n ", channel->number ) ;
01354      exit( -2 ) ;
01355  };
01356 
01357      /* allocate memory for first cell structure */
01358  if ( (cell = (struct Cells *) malloc( (size_t) sizeof( struct Cells ))) == NULL )
01359  { 
01360      printf( "Failed to allocate memory for cell structure\n " ) ;
01361      exit( -2 ) ;
01362  }
01363 
01364      /* store the pointer to first cell in the array */
01365  channel = channel_first;
01366 
01367      /* loop over all the canal reaches (== channel, defined as an individual canal reach with an ID number) */
01368  while (channel != NULL)
01369  {
01370      Reach = channel->reaches;
01371      c_num = channel->number;
01372      z  = channel->levee;
01373      RoIL = channel->roil;
01374      RoIR = channel->roir;
01375         
01376      Chan_list[c_num - CHo] = channel;
01377         
01378      channel->cells = cell;
01379      distTot = 0;
01380      C_length = 0;
01381      C_Mark = 0;
01382     
01383          /* a neg width canal is skipped, to SKIPCHAN near end of this cycle along canals */
01384      if (channel->width < 0)   goto SKIPCHAN; 
01385     
01386      fprintf ( CanalCellInteractDebug, "N = %d  levee = %d\n", c_num, z );
01387      
01388 /* get the elevation at the starting reach segment (may be up or down stream) of the canal */
01389      cellLoc_i=(int)(Reach->x0+1);
01390      cellLoc_j=(int)(Reach->y0+1);
01391      cellLoc = T(cellLoc_i,cellLoc_j );
01392      Chan_list[c_num - CHo]->elev_start = ElevMap[cellLoc]; 
01393      if ( !ON_MAP[cellLoc] ) {
01394         printf( "Starting elevation of Reach %d is out of the active domain.  Fix the reach location.\n",c_num);
01395         exit(-2);
01396      }
01397      fprintf ( CanalCellInteractDebug, " %d, %d, StartElev = %f\n", cellLoc_i,cellLoc_j,  Chan_list[c_num - CHo]->elev_start );
01398 
01399                 
01400      while (Reach != NULL)  /* loop through canal reach segments */
01401      {
01402          a0 = Reach->x0;
01403          b0 = Reach->y0;
01404          a1 = Reach->x1;
01405          b1 = Reach->y1;
01406 
01407          x1 = 0;  y1 = 0;
01408 
01409          T_length = sqrt( (b1 - b0)*(b1 - b0) + (a1 - a0)*(a1 - a0) );  /* total reach segment length */
01410          C_length += T_length;                                                                              /* total canal reach length */
01411 
01412          if ( a1 != a0 && b1 != b0 ) 
01413          {
01414              A1 = (b1 - b0)/(a1 - a0);                      /* define the equation: y = A1x + B1 */
01415              B1 = b0 - A1*a0;
01416              A2 = 1/A1;  B2 = - B1/A1;                      /* define the equation: x = A2y + B2 */
01417          } 
01418     
01419          x = a0; y = b0;
01420          i = Round (x);    j = Round (y);   
01421               
01422              /* loop along cells of the canal reach segment */ 
01423          while ( x1 != a1 || y1 != b1 )
01424          {
01425              if ( a1 == a0 )                                     /* if goes straight horizontally */
01426              {  x1 = a0; 
01427              y1 = (j == Round (b1) ) ? b1 : j + STEP(b1 - b0)*cell_step; 
01428              length = y1 - y;
01429              }
01430        
01431              else if ( b1 == b0 )                                   /* if goes straight vertically */
01432              { x1 = (i == Round (a1) ) ? a1 : i + STEP(a1 - a0)*cell_step; 
01433              y1 = b0; 
01434              length = x1 - x;
01435              }
01436 
01437              else if ( i == Round (a1) && j == Round (b1) )
01438              { x1 = a1;  y1 = b1;                                               /* if last cell */ 
01439              length = T_length * (y1 - y)/(b1 - b0);
01440              C_Mark = 0;
01441              }
01442           
01443              else
01444              { y1 = A1*(i + STEP(a1 - a0)*cell_step) + B1;  length = T_length * (y1 - y)/(b1 - b0);
01445              x1 = A2*(j + STEP(b1 - b0)*cell_step) + B2;  L = T_length * (x1 - x)/(a1 - a0);
01446     
01447              if ( length < L )                       /* define the end coordinates and length */
01448                  x1 = i + STEP(a1 - a0)*cell_step;  /* (the next coordinate is the closest)*/
01449              else
01450              {  y1 = j + STEP(b1 - b0)*cell_step;  length = L;
01451              }
01452              }
01453 
01454              if ( Abs (a1 - a0) < Abs (b1 - b0) )
01455                      /* define the direction of canal reach and  */
01456              { p_middle = (x1 + x)/2;      /* calculate the middle point on the canal reach in cell */
01457              m = i + cell_step/2;                                                               /* m - the center of the cell */
01458              HV = 1;  /* to indicate that the canal reach goes closer to the horizontal direction */
01459              }
01460              else
01461              { p_middle = (y1 + y)/2;
01462              m = j + cell_step/2;
01463              HV = 0;
01464              }
01465              
01466                  distTot = distTot + length*celWid;  /*  cumulative (along grid cells) distance along a canal reach, from the starting point (m) */
01467          
01468                  /* mark cells depending on where the center of canal reach is rel. to the cell center*/
01469 
01470              if ( Round (a1) == Round (a0) ) /* if the canal reach goes within the range of one */
01471                      /*cell horizontally, mark the two vertical cells */
01472              { if (b0 < b1)     /* if going from left to right */
01473              {          /* if the middle point of the canal reach is lower than the */
01474                      /* cell center, mark the cell as LEFT else as RIGHT */
01475                  if ( p_middle > m ) 
01476                  { cell = MarkCell ( i, j, LEFT, length, cell, EAST, z, i, j, c_num, marked, distTot );
01477                  L_range = (int)(RoIL - 1); R_range = (int)RoIR;
01478                  }  
01479                  else 
01480                  { cell = MarkCell ( i, j, RIGHT, length, cell, EAST, z, i - 1, j, c_num, marked, distTot );
01481                  L_range = (int) RoIL;  R_range = (int)(RoIR - 1);
01482                  }  
01483                  for ( k = 0; k < L_range; k++ )
01484                      cell = MarkCell ( i - k - 1, j, LEFT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01485                  for ( k = 0; k < R_range; k++ )
01486                      cell = MarkCell ( i + k + 1, j, RIGHT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01487              }
01488              else  /* vice versa if we go from right to left */
01489              {          /* if the middle point of the canal reach is lower than the */
01490                      /* cell center, mark the cell as RIGHT else as LEFT */
01491                  if ( p_middle > m ) 
01492                  { cell = MarkCell ( i, j, RIGHT, length, cell, EAST, z, i, j, c_num, marked, distTot );
01493                  L_range = (int)RoIL; R_range = (int)(RoIR - 1);
01494                  } 
01495                  else 
01496                  { cell = MarkCell ( i, j, LEFT, length, cell, EAST, z, i - 1, j, c_num, marked, distTot );
01497                  L_range = (int)(RoIL - 1); R_range = (int)RoIR;
01498                  }
01499                  for ( k = 0; k < L_range; k++ )
01500                      cell = MarkCell ( i + k + 1, j, LEFT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01501                  for ( k = 0; k < R_range; k++ )
01502                      cell = MarkCell ( i - k - 1, j, RIGHT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01503              }
01504              }
01505              else if ( Round (b1) == Round (b0) )       /* if the canal reach goes vertically */
01506              { if (a0 < a1) /* if going from top to bottom */
01507              {      /* if the middle point of the canal reach is further to the right than the */
01508                      /* cell center, mark the cell as RIGHT else as LEFT */
01509                  if ( p_middle > m ) 
01510                  { cell = MarkCell ( i, j, RIGHT, length, cell, SOUTH, z, i, j, c_num, marked, distTot );
01511                  L_range = (int)RoIL; R_range = (int)(RoIR -1);
01512                  }
01513                  else 
01514                  { cell = MarkCell ( i, j, LEFT, length, cell, SOUTH, z, i, j - 1, c_num, marked, distTot );
01515                  L_range = (int)(RoIL - 1); R_range = (int)RoIR;
01516                  }
01517                  for ( k = 0; k < L_range; k++ )
01518                      cell = MarkCell ( i, j + k + 1, LEFT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01519                  for ( k = 0; k < R_range; k++ )
01520                      cell = MarkCell ( i, j - k - 1, RIGHT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01521              }
01522              else  /* vice versa if we go from bottom to top */
01523              {      /* if the middle point of the canal reach is further to the left than the */
01524                      /* cell center, mark the cell as LEFT else as RIGHT */
01525                  if ( p_middle > m ) 
01526                  { cell = MarkCell ( i, j, LEFT, length, cell, SOUTH, z, i, j, c_num, marked, distTot );
01527                  L_range = (int)(RoIL - 1); R_range = (int)RoIR;
01528                  }
01529                  else 
01530                  { cell = MarkCell ( i, j, RIGHT, length, cell, SOUTH, z, i, j - 1, c_num, marked, distTot );
01531                  L_range = (int)RoIL; R_range = (int)(RoIR - 1);
01532                  }
01533                  for ( k = 0; k < L_range; k++ )
01534                      cell = MarkCell ( i, j - k - 1, LEFT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01535                  for ( k = 0; k < R_range; k++ )
01536                      cell = MarkCell ( i, j + k + 1, RIGHT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01537              }
01538              }
01539              else
01540                      /* if the canal reach crosses more than one cell (in both dimensions), mark diagonal cells */
01541              { if ( A1>0 )                 /* this is an extremely confusing check, but in this */
01542                  i45 = ( HV ) ? p_middle - m : m - p_middle;  /*form it seems to be working OK */
01543              else 
01544                  i45 = p_middle - m;
01545                  /*  b0>b1, a0>a1 \/ b0<b1, a0>a1 */
01546                  /*  b0>b1, a0<a1 /\ b0<b1, a0<a1 */
01547              if ( b0 < b1 )
01548              { if ( a0 < a1 )     /* when a0<a1 the canal reach goes down in the \ direction */
01549              { if ( i45 > 0 ) 
01550              { cell = MarkCell ( i, j, LEFT, length, cell, EAST, z, i, j, c_num, marked, distTot );
01551              cell = MarkCell ( i + 1, j - 1, RIGHT, length, cell, SOUTH, z, i + 1, j - 1, c_num, marked, distTot );
01552              }
01553              else 
01554              { cell = MarkCell ( i, j, RIGHT, length, cell, SOUTH, z, i, j, c_num, marked, distTot );
01555              cell = MarkCell ( i - 1, j + 1, LEFT, length, cell, EAST, z, i - 1, j + 1, c_num, marked, distTot );
01556              }
01557              for ( k = 0; k < RoIR-1; k++ )
01558                  cell = MarkCell ( i + k + 1, j - k - 1, RIGHT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01559              for ( k = 0; k < RoIL-1; k++ )
01560                  cell = MarkCell ( i - k - 1, j + k + 1, LEFT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01561              }
01562              else                                 /* in this case the canal reach goes up in the / direction */
01563              { if ( i45 > 0 ) 
01564              { cell = MarkCell ( i, j, LEFT, length, cell, NONE, z, i, j, c_num, marked, distTot );
01565              cell = MarkCell ( i + 1, j + 1, RIGHT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01566              }
01567              else 
01568              { cell = MarkCell ( i, j, RIGHT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01569              cell = MarkCell ( i - 1, j - 1, LEFT, length, cell, NONE, z, i - 1, j - 1, c_num, marked, distTot );
01570              }
01571              for ( k = 0; k < RoIR-1; k++ )
01572                  cell = MarkCell ( i + k + 1, j + k + 1, RIGHT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01573              for ( k = 0; k < RoIL-1; k++ )
01574                  cell = MarkCell ( i - k - 1, j - k - 1, LEFT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01575              }
01576              }
01577                  /* if going in the opposite direction change left and right */            
01578              else
01579              { if ( a0 > a1 )       /* when a0<a1 the canal reach goes up in the \ direction */
01580              { if ( i45 < 0 ) 
01581              { cell = MarkCell ( i, j, LEFT, length, cell, SOUTH, z, i, j, c_num, marked, distTot );
01582              cell = MarkCell ( i - 1, j + 1, RIGHT, length, cell, EAST, z, i - 1, j + 1, c_num, marked, distTot );
01583              }
01584              else 
01585              { cell = MarkCell ( i, j, RIGHT, length, cell, EAST, z, i, j, c_num, marked, distTot );
01586              cell = MarkCell ( i + 1, j - 1, LEFT, length, cell, SOUTH, z, i + 1, j - 1, c_num, marked, distTot );
01587              }
01588 
01589              for ( k = 0; k < RoIL-1; k++ )
01590                  cell = MarkCell ( i + k + 1, j - k - 1, LEFT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01591              for ( k = 0; k < RoIR-1; k++ )
01592                  cell = MarkCell ( i - k - 1, j + k + 1, RIGHT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01593              }
01594              else                               /* in this case the canal reach goes down in the / direction */
01595              { if ( i45 < 0 ) 
01596              { cell = MarkCell ( i, j, LEFT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01597              cell = MarkCell ( i - 1, j - 1, RIGHT, length, cell, NONE, z, i - 1, j - 1, c_num, marked, distTot );
01598              }
01599              else 
01600              { cell = MarkCell ( i, j, RIGHT, length, cell, NONE, z, i, j, c_num, marked, distTot );
01601              cell = MarkCell ( i + 1, j + 1, LEFT, length, cell, ALL, z, i + 1, j + 1, c_num, marked, distTot );
01602              }
01603 
01604              for ( k = 0; k < (int)(RoIL-1); k++ )
01605                  cell = MarkCell ( i + k + 1, j + k + 1, LEFT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01606              for ( k = 0; k < (int)(RoIR-1); k++ )
01607                  cell = MarkCell ( i - k - 1, j - k - 1, RIGHT, length, cell, ALL, z, i, j, c_num, marked, distTot );
01608              }
01609              }
01610              }
01611          
01612              num_cell++;        /*at this point we're counting the cells only once thinking that all 
01613                                   cells have the same length of interaction irrelevant of their distance from the canal 
01614                                   and that this length will be equal to the canal length / number of cells it crosses */
01615 
01616                  /* prepare for the next step in cycle */
01617              x = x1;  y = y1;
01618              i = Round (x);     if ( a1 < a0 && i == x ) i--;         /* if we're going upwards, we */
01619                  /* need to change the rule of Rounding so that we move ahead along cells */
01620              j = Round (y);     if ( b1 < b0 && j == y ) j--;
01621              C_Mark = 1;
01622     
01623 
01624          } /* end loop along cells of the canal reach segment  */
01625          
01626          Reach = Reach->next_reach;
01627      } /* end loop over all reach segments of a canal */
01628 
01629 /* get the elevation at the ending reach segment (may be up or down stream) of the canal */
01630      cellLoc_i=(int)(a1+1);
01631      cellLoc_j=(int)(b1+1);
01632      cellLoc = T(cellLoc_i,cellLoc_j );
01633      Chan_list[c_num - CHo]->elev_end = ElevMap[cellLoc];  
01634 
01635      Chan_list[c_num - CHo]->num_of_cells = num_cell;
01636      num_cell = 0;
01637          /* store the canal attributes of:
01638             length of entire canal, 
01639             area of entire canal, 
01640             minimum volume allowed in canal,
01641             avg length of cell along reach segment,
01642             overland flow coefficient (not incl. dynamic manning's n, C_F is for sensitivity experiments only (normally=1) ),
01643             seepage flow coefficient (including GP_calibGWat groundwater flow calibration coef),
01644             land surface elevation difference between starting and ending points of canal 
01645             slope of the elevation gradient from start to end points */
01646      Chan_list[c_num - CHo]->length = C_length*celWid;
01647      Chan_list[c_num - CHo]->area = Chan_list[c_num - CHo]->width*Chan_list[c_num - CHo]->length;
01648      Chan_list[c_num - CHo]->minVol = Chan_list[c_num - CHo]->area*MINDEPTH;
01649      Chan_list[c_num - CHo]->seg_len = Chan_list[c_num - CHo]->length/Chan_list[c_num - CHo]->num_of_cells;
01650      Chan_list[c_num - CHo]->SW_flow_coef = sqrt(Chan_list[c_num - CHo]->seg_len) * sec_per_day * C_F;
01651      Chan_list[c_num - CHo]->SPG_flow_coef = Chan_list[c_num - CHo]->cond * GP_calibGWat;
01652      Chan_list[c_num - CHo]->elev_drop = Chan_list[c_num - CHo]->elev_start - Chan_list[c_num - CHo]->elev_end;
01653      Chan_list[c_num - CHo]->elev_slope = Chan_list[c_num - CHo]->elev_drop / Chan_list[c_num - CHo]->length ; /* v2.5 calculated slope */
01654      
01655      if ( !ON_MAP[cellLoc] ) {
01656         printf( "Ending elevation of Reach %d is out of the active domain.  Fix the reach location.\n",c_num);
01657         exit(-2);
01658      }
01659      fprintf ( CanalCellInteractDebug, " %d, %d, EndingElev = %f, ElevSlope= %f, ReachLength= %f\n", cellLoc_i,cellLoc_j,  
01660                                         Chan_list[c_num - CHo]->elev_end, Chan_list[c_num - CHo]->elev_slope, Chan_list[c_num - CHo]->length );
01661      
01662      getCanalElev(c_num - CHo); /* v2.5 */
01663      
01664    SKIPCHAN: /* a single goto comes here, skipping a canal with negative width */
01665      channel = channel->next_in_list;
01666 
01667      cell_last->next_cell = NULL;                       /* close the cell list */
01668 
01669  } /* end loop loop over all the canal reaches */
01670  
01671  fclose ( CanalCellInteractDebug );
01672  
01673  return; 
01674 }

Here is the call graph for this function:

void getCanalElev int  chan_n  ) 
 

Assign elevation to canal using slope & distance from start point.

This function does one thing:
1. For each grid cell along a canal, we determine the elevation of the top of the canal bank, using the slope of the canal and the distance from the beginning point of the canal. (function new in v2.5)

Parameters:
chan_n Canal reach number (adjusted for the beginning ID number)
Returns:
void

Definition at line 1687 of file WatMgmt.c.

References CanalCellInteractDebug, cell, Chan::cells, Chan_list, CHo, Chan::elev_slope, Chan::elev_start, Cells::next_cell, Cells::reachDistDown, Cells::reachElev, T, Cells::x, and Cells::y.

Referenced by Channel_configure().

01688 {
01689      int i;
01690          /* cycle along canal cells */
01691      cell = Chan_list[chan_n]->cells;
01692      i = T(cell->x, cell->y);
01693      
01694      while ( cell != NULL )
01695      { 
01696         cell->reachElev = Chan_list[chan_n]->elev_start - cell->reachDistDown * Chan_list[chan_n]->elev_slope;
01697         fprintf ( CanalCellInteractDebug, " Reach= %d, Cell (%d, %d), ReachElev= %f\n", chan_n+CHo, cell->x, cell->y, cell->reachElev );
01698         cell = cell->next_cell;
01699      }
01700 }

void FluxChannel int  chan_n,
float *  SWH,
float *  ElevMap,
float *  MC,
float *  GWH,
float *  poros,
float *  GWcond,
double *  NA,
double *  PA,
double *  SA,
double *  GNA,
double *  GPA,
double *  GSA,
float *  Unsat,
float *  sp_yield
 

Runs the iterative algorithm over the network of canals.

The iterative relaxation routine calculates the new head in the canal and the water exchange with the adjacent cells, including constituent fluxes.

Parameters:
chan_n Canal ID number
SWH SURFACE_WAT
ElevMap SED_ELEV
MC HYD_MANNINGS_N
GWH SAT_WATER
poros HP_HYD_POROSITY
GWcond HYD_RCCONDUCT
NA DINdummy
PA TP_SF_WT
SA SALT_SURF_WT
GNA DINdummy
GPA TP_SED_WT
GSA SALT_SED_WT
Unsat UNSAT_WATER
sp_yield HP_HYD_SPEC_YIELD

Variables local to function
converged The convergence flag
error The error calculated as the diff between the incoming volume and the absorbed volume
prev_error The error calculated in the previous iteration
iter The iteration number in the iterative relaxation procedure
max_iter The max number of iterations allowed
init_factor The depth increment for the initial guess (now same as "factor")
factor The depth increment for the new guess (now same as "init_factor")
attempt The counter of number of attempts at convergence
CanWatDep The water level height in canal
T_flux_S, T_flux_G, T_flux_L The total net fluxes among canal and cells (Surface, Groundwater, Levee seepage)
fluxS, fluxG, fluxL The flux among canal and cells (Surface, Groundwater, Levee seepage)
fluxO, fluxCk Flow-checks to ensure the canal is not drained below its minimum level
SW_coef, SPG_coef, GW_coef Aggregated flow coefficients (Surface, Groundwater, Levee seepage)
N_fl, P_fl, S_fl Flux of Nitrogen (UNUSED), Phosphorus, and Salt/tracer among canal and cells
I_Length Average length of cell along reach segment
seg_area Average area of reach segment
MIN_VOL Minimum canal volume allowed
GW_head Groundwater head
h_GWflow, h_SPflow Water heights for groundwater and seepage flows
dh Difference between hydraulic heads
Qin, Qout, Qout1 Flow volumes
SWater SURFACE_WAT
tot_head Grid cell water head
CH_vol Total volume of canal at the estimated water depth
CH_vol_R Total canal volume prior to relaxation procedure
CH_bottElev Elevation of bottom of canal at cell location
H_rad_ch, H_rad_cell Hydraulic radii, canal (channel) and cell
can_cell_M, can_cell_A Temp vars used in reducing excess flows between canal<->cell
errorConv Numerical error after convergence

Definition at line 1805 of file WatMgmt.c.

References Abs, Chan::area, Chan::basin, basn, basn_list, Structure::canal_fr, Structure::canal_to, canDebugFile, cell, CELL_SIZE, Chan::cells, Chan_list, CHANNEL_MAX_ITER, CHo, debug, Chan::depth, dynERRORnum, Chan::edgeMann, F_ERROR, f_Ground(), f_Manning(), Chan::family, basndef::family, Structure::flow, HAB, Cells::ind, Chan::levee, Max, Min, MINDEPTH, Chan::minVol, msgStr, Chan::N, Cells::next_cell, Structure::next_in_list, Chan::number, ON_MAP, Chan::P, Chan::P_con, P_IN_GW, P_IN_OVL, P_IN_SPG, P_OUT_GW, P_OUT_OVL, P_OUT_SPG, Chan::parent, basndef::parent, ramp, Cells::reachElev, Chan::S, Chan::S_con, S_IN_GW, S_IN_OVL, S_IN_SPG, S_OUT_GW, S_OUT_OVL, S_OUT_SPG, Chan::seg_len, sgn, SimTime, Chan::SPG_flow_coef, structs, Chan::SW_flow_coef, T, simTime::TIME, True, VOL_IN_GW, VOL_IN_OVL, VOL_IN_SPG, VOL_OUT_GW, VOL_OUT_OVL, VOL_OUT_SPG, Chan::wat_depth, Chan::width, WriteMsg(), Cells::x, and Cells::y.

Referenced by Run_Canal_Network().

01809 { 
01810                 
01842 int converged = 0;                                              
01843  float error = 1.0;                                     
01844  float prev_error = 1.0;                                        
01845  int iter = 0;                                                  
01846  int max_iter = CHANNEL_MAX_ITER;               
01847  float init_factor = 0.5;                               
01848  float factor = init_factor;                    
01849 
01850  int i;
01851  int attempt=0;                                                 
01852  double CanWatDep = Chan_list[chan_n]->wat_depth;
01853  double T_flux_S, T_flux_G, T_flux_L;           
01854  double fluxO, fluxCk;
01855  double fluxS, fluxG, fluxL; 
01856  double SW_coef, SPG_coef, GW_coef;     
01857  double N_fl, P_fl, S_fl;                               
01858  /* double SNfl, GNfl, SPfl, GPfl, SSfl, GSfl; UNUSED */
01859  double I_Length;                                               
01860  double seg_area;                                               
01861  double MIN_VOL;                                                
01862  double GW_head;                                                
01863  double h_GWflow, h_SPflow;                             
01864  double dh;                                                     
01865  double Qin, Qout, Qout1;                               
01866  double SWater;                                                 
01867  double tot_head;                                               
01868  double CH_vol;                                         
01869  double CH_vol_R;                                       
01870  double CH_bottElev;                                    
01871  double H_rad_ch, H_rad_cell;                   
01872  double can_cell_M, can_cell_A;                 
01873  float errorConv;                                               
01874 
01875 /* UNUSED below, unimplemented gwater/sfwater integration vars */ 
01876     int equilib, revers=1;
01877     float H_pot_don, H_pot_rec;
01878     
01879     float fluxTOunsat_don, fluxHoriz; /* vertical and horiz components of the calculated total groundwater Flux */
01880     float UnsatZ_don, UnsatCap_don, UnsatPot_don, Sat_vs_unsat; /* unsaturated zone capacities, potentials, in the donor cell */
01881     float UnsatZ_rec, UnsatCap_rec, UnsatPot_rec; /* unsaturated zone capacities, potentials, in the recipient cell */
01882     float sfTOsat_don, unsatTOsat_don, unsatTOsat_rec, satTOsf_rec, horizTOsat_rec; /* fluxes for vertical integration */
01883     
01884     float sedwat_don, sedwat_rec;
01885     float sed_stuf1fl_rec, sed_stuf2fl_rec, sed_stuf3fl_rec; /* vertical flux of solutes */
01886     float sf_stuf1fl_don, sf_stuf2fl_don, sf_stuf3fl_don; /* vertical flux of solutes */
01887     float sed_stuf1fl_horz, sed_stuf2fl_horz, sed_stuf3fl_horz; /* horiz flux of solutes */
01888 /* UNUSED above, unimplemented gwater/sfwater integration vars */ 
01889  
01890  
01891 /* look at all the structures connected to the canal and add/subtract flow from/to them */
01892  Qin = 0; Qout = 0;
01893  structs = struct_first;
01894  while ( structs != NULL )
01895  {
01896      
01897      if ( structs->canal_fr  == chan_n+CHo )  
01898      {
01899          if ( structs->flow > 0 ) Qout += structs->flow;
01900          else Qin -= structs->flow;
01901      } 
01902      if ( structs->canal_to == chan_n+CHo ) 
01903      {
01904          if ( structs->flow > 0 ) Qin += structs->flow;
01905          else Qout -= structs->flow;
01906      } 
01907     
01908      structs = structs->next_in_list;
01909  }
01910   
01911  CH_vol_R = Chan_list[chan_n]->area * CanWatDep;  /* total canal volume prior to relaxation procedure */
01912  MIN_VOL =  Chan_list[chan_n]->minVol;                          /* minimum canal volume allowed */
01913  SPG_coef = Chan_list[chan_n]->SPG_flow_coef;   /* GP_calibGWat already in this */
01914 
01915  I_Length =   Chan_list[chan_n]->seg_len;          /* avg length of cell along reach segment */
01916  seg_area = I_Length * Chan_list[chan_n]->width;  /* avg area of reach segment */
01917  can_cell_M = Chan_list[chan_n]->area*CELL_SIZE;  /* used in reducing excess flows between canal<->cell */
01918  can_cell_A = Chan_list[chan_n]->area+CELL_SIZE;  /* used in reducing excess flows between canal<->cell */
01919  
01920 
01921  
01922 /* start relaxation procedure */  
01923  do
01924  {
01925      if   ( Abs( error) < F_ERROR ) 
01926      {
01927              converged = 1; /* when we've converged, we'll go through the cell list one more time w/o modifying the depth estimate */
01928                         /* occasionally, we get a modified total flux value, and diff error, a second time with the same depth estimate! */
01929              errorConv=error;
01930              factor = 0;
01931      }
01932      else if  ( iter >= (max_iter-1) )
01933      {
01934          if (attempt == 4) {
01935              converged = 1; /* haven't actually converged, but stop now */
01936              sprintf(msgStr,"Day %6.1f: ERROR - couldn't converge on Canal %d, HALTING the iterative relaxation after 4 tries, error = %f.", 
01937                         SimTime.TIME, Chan_list[chan_n]->number, error); 
01938              WriteMsg( msgStr,True ); dynERRORnum++;
01939          }
01940          else {
01941              if (debug > 3) { sprintf(msgStr,"Day %6.1f: Warning - couldn't converge on Canal %d, redoing the iterative relaxation, error = %f.", 
01942                         SimTime.TIME, Chan_list[chan_n]->number, error); 
01943              WriteMsg( msgStr,True ); }
01944 
01945              attempt++; /* the attempt counter increment */
01946              iter = 0; /* reset iterations to 0, and double factor & max # iterations on first try */ 
01947              if (attempt == 1) { factor = init_factor*2.0;  max_iter = CHANNEL_MAX_ITER * 2; }
01948              else if (attempt == 2) { factor = init_factor*0.5; max_iter = CHANNEL_MAX_ITER * 2; }
01949              else if (attempt == 3) { factor = init_factor*4.0; max_iter = CHANNEL_MAX_ITER * 4; }
01950              else if (attempt == 4) { factor = init_factor*0.05; max_iter = CHANNEL_MAX_ITER * 4; }
01951          }
01952      }
01953 
01954      else
01955      {
01956              /* if error changed sign, we have overshot the target, therefore reverse direction */
01957          if (error * prev_error < 0 && iter != 1) 
01958              factor = - sgn(error) *  Abs( factor ) * 0.5;      /* and reduce the step */
01959          
01960              /* if error contiues to grow after first iter.(#0), reverse direction */   
01961          if (iter == 1 && error > 0) factor = - factor;
01962          
01963      }
01964      
01965      if (iter != 0 && !converged ) 
01966      {
01967          CanWatDep += factor; 
01968          prev_error = error;
01969          if ( CanWatDep < MINDEPTH ) 
01970          { 
01971              if (debug > 3) { sprintf(msgStr,"Day %6.1f: Warning - estimate is below minimum depth, revising estimate. (Canal %d, depth %f)", 
01972                         SimTime.TIME, Chan_list[chan_n]->number,CanWatDep); 
01973              WriteMsg( msgStr,True ); }            
01974 
01975              CanWatDep = MINDEPTH; 
01976             /* converged = 1; */
01977             factor = -factor*0.1; /* reverse the direction and decrease for next estimate*/
01978          }      
01979      }     
01980     
01981      T_flux_G = T_flux_S = T_flux_L = 0.0;           /* at each convergence iter, set the sum flows to 0 */
01982      CH_vol = Chan_list[chan_n]->area * CanWatDep;  /* volume of canal at the estimated water depth */
01983     
01984          /* cycle along canal cells */
01985      cell = Chan_list[chan_n]->cells;
01986      
01987      while ( cell != NULL )
01988      { 
01989          i = T(cell->x, cell->y);
01990          if ( !ON_MAP[i] )       /* check that cell is ON_MAP, else we do nothing */
01991                     /* return to beginning of do loop                */
01992          {  cell = cell->next_cell;  continue;
01993          }
01994          if (basn[i] == 0) {
01995              sprintf(msgStr,"ERROR - Cell (%d, %d) is out-of-system.", cell->x, cell->y); 
01996              WriteMsg( msgStr,True );  dynERRORnum++; }
01997 
01998          fluxS = 0.;  fluxG = 0.; fluxL = 0.; /* overland, subsurface, and seepage cell-canal fluxes set to 0 */
01999            /* v2.5 In the case of some canal reaches, there is a lip or berm, and often associated vegetation, that is 
02000            along the edge of a canal, decreasing flow rates between canal & marsh.  If user provides a non-zero roughness
02001            parameter in the CanalData.chan input file, the the mean roughness of this berm and the cell is used for 
02002            surface exchanges */
02003          SW_coef = ( MC[i] == 0.0 ) ? 0 : Chan_list[chan_n]->SW_flow_coef / ( (Chan_list[chan_n]->edgeMann > 0.0) ?  (MC[i] + Chan_list[chan_n]->edgeMann)/2 : MC[i] );
02004          GW_coef = GWcond[i] * poros[HAB[i]];   
02005          SWater = SWH[i];
02006 
02007          if (debug>0 && (SWater < -GP_MinCheck) ) {
02008          sprintf(msgStr,"Day %6.1f: capacityERR - neg sfwater along canal %d, candepth %f, cell (%d, %d), celldepth %f", 
02009                  SimTime.TIME,Chan_list[chan_n]->number,CanWatDep, cell->x, cell->y, SWater); 
02010          WriteMsg( msgStr,True );  dynERRORnum++; }
02011          
02012          GW_head = GWH[i]/poros[HAB[i]];
02013          tot_head = SWater + ElevMap[i]; 
02014          CH_bottElev = cell->reachElev - Chan_list[chan_n]->depth;  /* elev of bottom of canal at cell location (v2.5 is from fixed vector slope instead of actual grid cell elev) */
02015              
02016          dh = ( CH_bottElev + CanWatDep ) - tot_head;
02017 
02018              /* hydraulic radius of canal for overland flow out of canal */
02019          H_rad_ch =   Max ( ( seg_area * ramp(CanWatDep - Chan_list[chan_n]->depth) +
02020                         SWater * (CELL_SIZE-seg_area) ) / CELL_SIZE, 0.0);   /* v2.5 contrain to >= 0 */
02021             /* hydraulic radius of cell for overland flow into canal (same eqn right now) */
02022          H_rad_cell = Max ( ( seg_area * ramp(CanWatDep - Chan_list[chan_n]->depth) +
02023                         SWater * (CELL_SIZE-seg_area) ) / CELL_SIZE, 0.0);     /* v2.5 contrain to >= 0 */
02024 
02025                 /* determine the heights of the water cross sections associated with the subsurf and seep flows */
02026          if (dh > 0.0) { /* from canal */
02027                 h_GWflow = Min(Chan_list[chan_n]->depth, CanWatDep);
02028                 h_SPflow = Max(CH_bottElev + CanWatDep - ElevMap[i], 0.0);
02029          }
02030          else if (dh < 0.0) { /* from cell */
02031                 h_GWflow = Max(GW_head-CH_bottElev, 0.0);
02032                 h_SPflow = Max(tot_head-ElevMap[i], 0.0);
02033          }
02034 
02035                  
02036          switch( Chan_list[chan_n]->levee )
02037          {      
02038                 
02039              case 2:    /* if levees are on both sides of the canal */
02040                          /*seepage and subsurface flow on both sides */
02041                  fluxL =  (h_SPflow > 0.0) ? (f_Ground ( dh, h_SPflow, SPG_coef, I_Length) ) : (0.0);
02042                  fluxG =  (h_GWflow > 0.0) ? (f_Ground ( dh, h_GWflow, GW_coef, I_Length) ) : (0.0); 
02043                  break;
02044              case 0: /* if there are no levees */
02045                      /*overland flow either side */
02046                  if ( dh > 0 ) {
02047                          /* flux from canal to cell  */
02048                      fluxS = f_Manning ( dh, H_rad_ch, SW_coef);
02049                  }
02050                  else if (SWater>GP_DetentZ) {
02051                      fluxS =  f_Manning ( dh, H_rad_cell, SW_coef) ; 
02052                     /* fluxS =  ( (-fluxS) > (SWater-GP_DetentZ) *CELL_SIZE )  ? ( -(SWater-GP_DetentZ)*CELL_SIZE ) : (fluxS);*/
02053                     if (-fluxS > (SWater-GP_DetentZ) *CELL_SIZE )
02054                         fluxS = -(SWater-GP_DetentZ)*CELL_SIZE;
02055                  }
02056           
02057                                 /* subsurface groundwater flow both sides */
02058                  fluxG =  (h_GWflow > 0.0) ? (f_Ground ( dh, h_GWflow, GW_coef, I_Length) ) : (0.0); 
02059                  break;
02060              case 1:  /* if levee is on the left */ 
02061                  if ( cell->ind < 0  )
02062                  { /*overland flow */
02063                      if ( dh > 0 ) {
02064                              /* flux from canal to cell  */
02065                          fluxS = f_Manning ( dh, H_rad_ch, SW_coef);
02066                      }
02067                      else if (SWater>GP_DetentZ) {
02068                          fluxS = f_Manning ( dh, H_rad_cell, SW_coef) ; 
02069                         /* fluxS =  ( (-fluxS) > (SWater-GP_DetentZ) *CELL_SIZE )  ? ( -(SWater-GP_DetentZ)*CELL_SIZE ) : (fluxS);*/
02070                     if (-fluxS > (SWater-GP_DetentZ) *CELL_SIZE )
02071                         fluxS = -(SWater-GP_DetentZ)*CELL_SIZE;
02072                      }
02073           
02074                  }      
02075                  else 
02076                          /*seepage on the levee side */
02077                      fluxL =  (h_SPflow > 0.0) ? (f_Ground ( dh, h_SPflow, SPG_coef, I_Length) ) : (0.0);
02078                      
02079                                 /* subsurface groundwater flow both sides */
02080                  fluxG =  (h_GWflow > 0.0) ? (f_Ground ( dh, h_GWflow, GW_coef, I_Length) ) : (0.0); 
02081                  break;
02082              case -1: /* if levee is on the right */
02083                  if ( cell->ind > 0 )
02084                  { /*overland flow */
02085                      if ( dh > 0 ) {
02086                              /* flux from canal to cell */
02087                          fluxS = f_Manning( dh, H_rad_ch, SW_coef);
02088                      }
02089                      else if (SWater>GP_DetentZ) {
02090                          fluxS = f_Manning ( dh, H_rad_cell, SW_coef); 
02091                          /*fluxS =  ( (-fluxS) > (SWater-GP_DetentZ) *CELL_SIZE )  ? ( -(SWater-GP_DetentZ)*CELL_SIZE ) : (fluxS);*/
02092                     if (-fluxS > (SWater-GP_DetentZ) *CELL_SIZE )
02093                         fluxS = -(SWater-GP_DetentZ)*CELL_SIZE;
02094                      }
02095                  }      
02096                  else
02097                          /*seepage on the levee side */
02098                      fluxL =  (h_SPflow > 0.0) ? (f_Ground ( dh, h_SPflow, SPG_coef, I_Length) ) : (0.0);
02099 
02100                                 /* subsurface groundwater flow both sides */
02101                  fluxG =  (h_GWflow > 0.0) ? (f_Ground ( dh, h_GWflow, GW_coef, I_Length) ) : (0.0); 
02102                  break;  
02103              default:
02104                  printf ("Error in levee location: Channel N %d", chan_n+CHo );
02105                  exit (-1);
02106                  break;
02107          } /*end switch */
02108           
02109 
02110 
02111              /* can_cell_M and can_cell_A are the areas of the canal and the cell multiplied (M) and added (A) together, respectively */
02112              /* the flux eqn is the same as the virtual weir eqn, equilibrating the heads across two unequal area regions */
02113 
02114              /* Surface flows: making sure that the cell head is not higher than the canal head */
02115          if (fluxS>0 && ( (fluxCk = tot_head + fluxS/CELL_SIZE - (CH_bottElev + CanWatDep - fluxS/Chan_list[chan_n]->area) ) > 0.0 ) ) {
02116              fluxS = (can_cell_M*(CH_bottElev + CanWatDep)/can_cell_A - can_cell_M*tot_head/can_cell_A );
02117              fluxCk = tot_head + fluxS/CELL_SIZE - (CH_bottElev + CanWatDep - fluxS/Chan_list[chan_n]->area);
02118              if (fluxCk > F_ERROR && debug>3) {
02119                  sprintf(msgStr,"Day %6.1f: Warning - excessive head in cell after flux from Canal %d (%f m diff)", 
02120                          SimTime.TIME, Chan_list[chan_n]->number,fluxCk); 
02121                  WriteMsg( msgStr,True );  
02122              }
02123          }
02124          
02125              /* Surface flows: making sure that the canal water level is not flooded much higher than the water level in the cell */
02126          if ( fluxS<0 && ( (fluxCk = (CH_bottElev + CanWatDep - fluxS/Chan_list[chan_n]->area) - (tot_head + fluxS/CELL_SIZE) ) > 0.0 ) ) {
02127              fluxS = (can_cell_M*(CH_bottElev + CanWatDep)/can_cell_A - can_cell_M*tot_head/can_cell_A );
02128              fluxCk =  (CH_bottElev + CanWatDep - fluxS/Chan_list[chan_n]->area) - (tot_head + fluxS/CELL_SIZE);
02129              if (fluxCk > F_ERROR && debug>3) {
02130                  sprintf(msgStr,"Day %6.1f: Warning - excessive head in Canal %d after flux from cell (%f m diff)", 
02131                          SimTime.TIME, Chan_list[chan_n]->number,fluxCk); 
02132                  WriteMsg( msgStr,True );  
02133              }
02134          }
02135 
02136              /* now ensure the canal is not drained below its minimum level */ 
02137          if ( fluxS > 0 || fluxG > 0 || fluxL > 0) /* check gw, seep, and sfwat */
02138          {
02139              fluxO  = (fluxS > 0) ? fluxS : 0;
02140              fluxO += (fluxG > 0) ? fluxG : 0;
02141              fluxO += (fluxL > 0) ? fluxL : 0;
02142     
02143                  /* using vol associated with new estimate (before subtracting losses) */
02144              if ( (fluxCk = CH_vol - MIN_VOL - fluxS - fluxG - fluxL) < 0.0 )
02145              { fluxCk =  (fluxO>0) ? ( ( fluxO + fluxCk )/fluxO) : (0.0);       /* v2.5 removed a divide by zero, but shouldn't be needed */
02146              if (fluxS>0) fluxS *= /* (double) */fluxCk; /* reduce all pos outflows by this proportion */
02147              if (fluxG>0) fluxG *= fluxCk;
02148              if (fluxL>0) fluxL *= fluxCk;
02149              if ( (fluxCk = (CH_vol - MIN_VOL - fluxS - fluxG - fluxL)/Chan_list[chan_n]->area) < -GP_MinCheck ) {
02150                  sprintf(msgStr,"Day %6.1f: capacityERR - excessive draining from Canal %d (%f m below min)", SimTime.TIME,Chan_list[chan_n]->number,fluxCk); 
02151                  WriteMsg( msgStr,True ); dynERRORnum++; }
02152              }
02153          }
02154 
02155 /******** UNUSED *********************/
02156     if (debug>99 && fluxG != 0.0) do {
02157         /* This is a part of the routine used to integrate surface/groundwater
02158          fluxes in the cell-cell Fluxes.c.  It only modifies the allocations if
02159          a cell is a recipient of groundwater flux (does not deal with cells as
02160          donor, which is handled in Fluxes.c).
02161          NOTE:****This routine is NOT fully integrated into the canal iterative solution,
02162          and thus is NOT used for now - get to it some day if it is determined to be needed.  */
02163        fluxHoriz =  fluxG;
02164 
02165 
02166 /**** recipient cell's **pre-flux** capacities */
02167         UnsatZ_rec  =  (fluxG>0) ? ( ElevMap[i] - GWH[i] / poros[HAB[i]] ) : (0.0) ; /* unsat zone depth */
02168         if (debug>2 && (UnsatZ_rec < -0.001 ) ) {
02169             sprintf(msgStr,"Day %6.1f: capacityERR - neg unsat depth (%f m) in recipient cell (%d,%d) in canal-cell gw fluxes!", 
02170                     SimTime.TIME, UnsatZ_rec, cell->x, cell->y ); WriteMsg( msgStr,True ); dynERRORnum++; }
02171             
02172         UnsatCap_rec =  (fluxG>0) ? (UnsatZ_rec * poros[HAB[i]]) : (0.0); /* maximum pore space capacity (UnsatCap)
02173                                                               in the depth (UnsatZ) of new unsaturated zone */
02174         UnsatPot_rec  = (fluxG>0) ? (UnsatCap_rec - Unsat[i] ) : (0.0); /* (height of) the volume of pore space (soil "removed")
02175                                                            that is unoccupied in the unsat zone */
02176      /* sf-unsat-sat fluxes */
02177         horizTOsat_rec = (fluxG>0) ? (fluxHoriz) : (0.0); /* lateral inflow to soil into saturated storage */
02178         satTOsf_rec = (fluxG>0) ? (Max(fluxHoriz - UnsatPot_rec, 0.0) ) : (0.0); /* upwelling beyond soil capacity */
02179         if (fluxG>0) unsatTOsat_rec = (UnsatZ_rec > 0.0) ? /* incorporation of unsat moisture into sat storage with
02180                                                  rising water table due to horiz inflow */
02181             ( ((horizTOsat_rec-satTOsf_rec)/poros[HAB[i]] ) / UnsatZ_rec * Unsat[i] ) :
02182             (0.0);
02183         else
02184             unsatTOsat_rec = 0.0;
02185         
02186 /**** potential new head in recipient cell will be ... */
02187         if (fluxG>0)  H_pot_rec = (GWH[i] + horizTOsat_rec + unsatTOsat_rec - satTOsf_rec)
02188             / poros[HAB[i]] + (SWater + satTOsf_rec);
02189 
02190 
02191         equilib = 1;
02192         
02193     } while  (!equilib); /* end of ***incomplete** routine for integrating groundwater-surface water and preventing head reversals */
02194 /******** UNUSED *********************/
02195 
02196    /*  if (fluxG>0) fluxG=fluxHoriz; */ /* unused, part of unimplemented sf/gwater routine */
02197     
02198          CH_vol -= ( fluxS + fluxG + fluxL); /* subtract flows (pos == loss) from volume associated with new estimate */
02199          
02200          if ( converged )
02201          {      /* the nutrient mass fluxes (pre-flux volumes) (mass=kg, vol=m3)*/
02202                  /* remember pos flux is from canal */
02203                  /* (CH_vol_R is vol prior to in/out fluxes) */
02204         /* surface water canal-cell flux */               
02205              if ( fluxS != 0 ) {
02206              if ( fluxS > 0 ) /* surface water flux from canal to cell */
02207              {  
02208                  /* this is where Disp_Calc comes in */ /* hcf TEST to see about reducing dispersion, halfing the flux of salt (i.e., ca. 500m grid dispersion) */
02209                  N_fl = (CH_vol_R > 0) ? (Chan_list[chan_n]->N/CH_vol_R*fluxS) : (0.0);
02210                  P_fl = (CH_vol_R > 0) ? (Chan_list[chan_n]->P/CH_vol_R*fluxS) : (0.0);
02211                  S_fl = (CH_vol_R > 0) ? (Chan_list[chan_n]->S/CH_vol_R*fluxS) : (0.0);
02212 
02213                  if (Chan_list[chan_n]->basin != basn[i] && debug > -999 ) {
02214                      if ( (Chan_list[chan_n]->family !=  basn_list[basn[i]]->family) && (debug > 0) ) {        /* if the flow is not within the family... */
02215                          sprintf(msgStr,"Day %6.1f: ERROR - no basin-basin canal-cell overland flow! (Reach %d, basn %d, with cell %d,%d, basin %d)", 
02216                                  SimTime.TIME,Chan_list[chan_n]->number, Chan_list[chan_n]->basin, cell->x, cell->y, basn[i]); 
02217                          WriteMsg( msgStr,True );  dynERRORnum++;
02218                      }
02219                      else {    /* if it's flow within a family, just keep
02220                                   track of what the children do among themselves */            
02221                          if  ( !Chan_list[chan_n]->parent  ) { 
02222                          /* then find out about the flow for the family's sake */
02223                              VOL_OUT_OVL[Chan_list[chan_n]->basin] += fluxS;
02224                              P_OUT_OVL[Chan_list[chan_n]->basin] += P_fl;
02225                              S_OUT_OVL[Chan_list[chan_n]->basin] += S_fl;
02226                          }
02227                          if ( !basn_list[basn[i]]->parent ) {                 
02228                                  /* then find out about the flow for the family's sake */
02229                              VOL_IN_OVL[basn[i]] += fluxS; 
02230                              P_IN_OVL[basn[i]] += P_fl; 
02231                              S_IN_OVL[basn[i]] += S_fl; 
02232                          }
02233                      }
02234                  }                
02235              }
02236              else /* surface water flux from cell to canal */
02237              {
02238                  /* this is where Disp_Calc comes in */
02239                  N_fl = ( SWH[i] > 0 ) ? NA[i]*fluxS/CELL_SIZE/SWH[i] : 0.;
02240                  P_fl = ( SWH[i] > 0 ) ? PA[i]*fluxS/CELL_SIZE/SWH[i] : 0.;
02241                  S_fl = ( SWH[i] > 0 ) ? SA[i]*fluxS/CELL_SIZE/SWH[i] : 0.;
02242                  if (Chan_list[chan_n]->basin != basn[i] && debug > -999 ) {
02243 
02244                      if ( (Chan_list[chan_n]->family !=  basn_list[basn[i]]->family) && (debug > 0) ) {        /* if the flow is not within the family... */
02245                          sprintf(msgStr,"Day %6.1f: ERROR - no basin-basin canal-cell overland flow! (Reach %d, basn %d, with cell %d,%d, basin %d)", 
02246                                         SimTime.TIME,Chan_list[chan_n]->number, Chan_list[chan_n]->basin, cell->x, cell->y, basn[i]); 
02247                          WriteMsg( msgStr,True );  dynERRORnum++;
02248                  }
02249                  else {    /* if it's flow within a family, just keep
02250                           track of what the children do among themselves */            
02251                      if  ( !basn_list[basn[i]]->parent  ) { 
02252                          /* then find out about the flow for the family's sake */
02253                          VOL_OUT_OVL[basn[i]] -= fluxS;
02254                          P_OUT_OVL[basn[i]] -= P_fl;
02255                          S_OUT_OVL[basn[i]] -= S_fl;
02256                      }
02257                      if ( !Chan_list[chan_n]->parent ) {                 
02258                          /* then find out about the flow for the family's sake */
02259                          VOL_IN_OVL[Chan_list[chan_n]->basin] -= fluxS; 
02260                          P_IN_OVL[Chan_list[chan_n]->basin] -= P_fl; 
02261                          S_IN_OVL[Chan_list[chan_n]->basin] -= S_fl; 
02262                      }
02263                  }
02264                  }
02265 
02266 
02267              } 
02268                  /* now pass that mass flux to/from the cell/canal */
02269              NA[i] += N_fl;
02270              PA[i] += P_fl;
02271              SA[i] += S_fl;
02272              Chan_list[chan_n]->N -= N_fl;
02273              Chan_list[chan_n]->P -= P_fl;
02274              Chan_list[chan_n]->S -= S_fl;
02275              }
02276 
02277         /* groundwater canal-cell flux */                 
02278              if ( fluxG != 0 ) {
02279              if ( fluxG > 0 ) /* ground water flux from canal to cell */
02280              {  
02281                  N_fl = (CH_vol_R > 0) ? (Chan_list[chan_n]->N/CH_vol_R*fluxG) : (0.0);
02282                  P_fl = (CH_vol_R > 0) ? (Chan_list[chan_n]->P/CH_vol_R*fluxG) : (0.0);
02283                  S_fl = (CH_vol_R > 0) ? (Chan_list[chan_n]->S/CH_vol_R*fluxG) : (0.0);
02284 
02285                  if (Chan_list[chan_n]->basin != basn[i] && debug > -999 ) {
02286                      if (Chan_list[chan_n]->family !=  
02287                          basn_list[basn[i]]->family ) {        /* if the flow is not within the family... */
02288                          if  ( !Chan_list[chan_n]->parent  ) { /* and if the donor or recipient is a child... */
02289                          /* then find out about the flow for the family's sake */
02290                              VOL_OUT_GW[Chan_list[chan_n]->family] += fluxG;
02291                              P_OUT_GW[Chan_list[chan_n]->family] += P_fl;
02292                              S_OUT_GW[Chan_list[chan_n]->family] += S_fl;
02293                          }
02294                          if ( !basn_list[basn[i]]->parent ) {                 
02295                            /* then find out about the flow for the family's sake */
02296                              VOL_IN_GW[basn_list[basn[i]]->family] += fluxG; 
02297                              P_IN_GW[basn_list[basn[i]]->family] += P_fl; 
02298                              S_IN_GW[basn_list[basn[i]]->family] += S_fl; 
02299                          }
02300                      VOL_OUT_GW[Chan_list[chan_n]->basin] += fluxG;
02301                      VOL_IN_GW[basn[i]] += fluxG; 
02302                      P_OUT_GW[Chan_list[chan_n]->basin] += P_fl;
02303                      P_IN_GW[basn[i]] += P_fl; 
02304                      S_OUT_GW[Chan_list[chan_n]->basin] += S_fl;
02305                      S_IN_GW[basn[i]] += S_fl; 
02306                  }
02307                  else {    /* if it's flow within a family, just keep
02308                           track of what the children do among themselves */            
02309                      if  ( !Chan_list[chan_n]->parent  ) { 
02310                          /* then find out about the flow for the family's sake */
02311                          VOL_OUT_GW[Chan_list[chan_n]->basin] += fluxG;
02312                          P_OUT_GW[Chan_list[chan_n]->basin] += P_fl;
02313                          S_OUT_GW[Chan_list[chan_n]->basin] += S_fl;
02314                      }
02315                      if ( !basn_list[basn[i]]->parent ) {                 
02316                          /* then find out about the flow for the family's sake */
02317                          VOL_IN_GW[basn[i]] += fluxG; 
02318                          P_IN_GW[basn[i]] += P_fl; 
02319                          S_IN_GW[basn[i]] += S_fl; 
02320                      }
02321                  }
02322                  }
02323              } /* end of positive flow evaluations */
02324              
02325              
02326              else  /*  ground water flux from cell to canal */
02327              {
02328                  N_fl = ( GWH[i] > 0 ) ? GNA[i]*fluxG/CELL_SIZE/GWH[i] : 0.;
02329                  P_fl = ( GWH[i] > 0 ) ? GPA[i]*fluxG/CELL_SIZE/GWH[i] : 0.;
02330                  S_fl = ( GWH[i] > 0 ) ? GSA[i]*fluxG/CELL_SIZE/GWH[i] : 0.;
02331                  if (Chan_list[chan_n]->basin != basn[i] && debug > -999) {
02332                      if (Chan_list[chan_n]->family !=  
02333                          basn_list[basn[i]]->family ) {        /* if the flow is not within the family... */
02334                          if  ( !basn_list[basn[i]]->parent  ) { /* and if the donor or recipient is a child... */
02335                          /* then find out about the flow for the family's sake */
02336                              VOL_OUT_GW[basn_list[basn[i]]->family] -= fluxG;
02337                              P_OUT_GW[basn_list[basn[i]]->family] -= P_fl;
02338                              S_OUT_GW[basn_list[basn[i]]->family] -= S_fl;
02339                          }
02340                          if ( !Chan_list[chan_n]->parent ) {                 
02341                            /* then find out about the flow for the family's sake */
02342                              VOL_IN_GW[Chan_list[chan_n]->family] -= fluxG; 
02343                              P_IN_GW[Chan_list[chan_n]->family] -= P_fl; 
02344                              S_IN_GW[Chan_list[chan_n]->family] -= S_fl; 
02345                          }
02346                      VOL_IN_GW[Chan_list[chan_n]->basin] -= fluxG;
02347                      VOL_OUT_GW[basn[i]] -= fluxG;
02348                      P_IN_GW[Chan_list[chan_n]->basin] -= P_fl;
02349                      P_OUT_GW[basn[i]] -= P_fl;
02350                      S_IN_GW[Chan_list[chan_n]->basin] -= S_fl;
02351                      S_OUT_GW[basn[i]] -= S_fl;
02352                  }
02353                  else {    /* if it's flow within a family, just keep
02354                           track of what the children do among themselves */            
02355                      if  ( !basn_list[basn[i]]->parent  ) { 
02356                          /* then find out about the flow for the family's sake */
02357                          VOL_OUT_GW[basn[i]] -= fluxG;
02358                          P_OUT_GW[basn[i]] -= P_fl;
02359                          S_OUT_GW[basn[i]] -= S_fl;
02360                      }
02361                      if ( !Chan_list[chan_n]->parent ) {                 
02362                          /* then find out about the flow for the family's sake */
02363                          VOL_IN_GW[Chan_list[chan_n]->basin] -= fluxG; 
02364                          P_IN_GW[Chan_list[chan_n]->basin] -= P_fl; 
02365                          S_IN_GW[Chan_list[chan_n]->basin] -= S_fl; 
02366                      }
02367                  }
02368                  }
02369              } 
02370                  /* now pass that mass flux to/from the cell/canal
02371                   this includes allocating nuts to surface water thru upflow (unimplemented) */
02372              /* unused, part of unimplemented sf/gwater routine */
02373 /*              SNfl = Chan_list[chan_n]->N/CH_vol_R*satTOsf_rec ; */
02374 /*              GNfl = N_fl - SNfl; */
02375 /*              SPfl = Chan_list[chan_n]->P/CH_vol_R*satTOsf_rec ; */
02376 /*              GPfl = P_fl - SPfl; */
02377 /*              SSfl = Chan_list[chan_n]->S/CH_vol_R*satTOsf_rec ; */
02378 /*              GSfl = S_fl - SSfl; */
02379 /*              GNA[i] += GNfl; */
02380 /*              GPA[i] += GPfl; */
02381 /*              GSA[i] += GSfl; */
02382              GNA[i] += N_fl;
02383              GPA[i] += P_fl;
02384              GSA[i] += S_fl;
02385              Chan_list[chan_n]->N -= N_fl;
02386              Chan_list[chan_n]->P -= P_fl;
02387              Chan_list[chan_n]->S -= S_fl;
02388                  /* upflow into surface water nuts (this is actually canal conc., but what the heck) */
02389              /* unused, part of unimplemented sf/gwater routine */
02390 /*              NA[i] += SNfl; */
02391 /*              PA[i] += SPfl; */
02392 /*              SA[i] += SSfl; */
02393              }
02394 
02395         /* seepage canal-cell flux */             
02396              if ( fluxL != 0 ) {
02397              if ( fluxL > 0 ) /* seepage flux from canal to cell */
02398              {  
02399                  N_fl = (CH_vol_R > 0) ? (Chan_list[chan_n]->N/CH_vol_R*fluxL) : (0.0);
02400                  P_fl = (CH_vol_R > 0) ? (Chan_list[chan_n]->P/CH_vol_R*fluxL) : (0.0);
02401                  S_fl = (CH_vol_R > 0) ? (Chan_list[chan_n]->S/CH_vol_R*fluxL) : (0.0);
02402                  if (Chan_list[chan_n]->basin != basn[i] && debug > -999 ) {
02403                      if (Chan_list[chan_n]->family !=  
02404                          basn_list[basn[i]]->family ) {        /* if the flow is not within the family... */
02405                          if  ( !Chan_list[chan_n]->parent  ) { /* and if the donor or recipient is a child... */
02406                          /* then find out about the flow for the family's sake */
02407                              VOL_OUT_SPG[Chan_list[chan_n]->family] += fluxL;
02408                              P_OUT_SPG[Chan_list[chan_n]->family] += P_fl;
02409                              S_OUT_SPG[Chan_list[chan_n]->family] += S_fl;
02410                          }
02411                          if ( !basn_list[basn[i]]->parent ) {                 
02412                            /* then find out about the flow for the family's sake */
02413                              VOL_IN_SPG[basn_list[basn[i]]->family] += fluxL; 
02414                              P_IN_SPG[basn_list[basn[i]]->family] += P_fl; 
02415                              S_IN_SPG[basn_list[basn[i]]->family] += S_fl; 
02416                          }
02417                      VOL_OUT_SPG[Chan_list[chan_n]->basin] += fluxL;
02418                      VOL_IN_SPG[basn[i]] += fluxL; 
02419                      P_OUT_SPG[Chan_list[chan_n]->basin] += P_fl;
02420                      P_IN_SPG[basn[i]] += P_fl; 
02421                      S_OUT_SPG[Chan_list[chan_n]->basin] += S_fl;
02422                      S_IN_SPG[basn[i]] += S_fl; 
02423                  }
02424                  else {    /* if it's flow within a family, just keep
02425                           track of what the children do among themselves */            
02426                      if  ( !Chan_list[chan_n]->parent  ) { 
02427                          /* then find out about the flow for the family's sake */
02428                          VOL_OUT_SPG[Chan_list[chan_n]->basin] += fluxL;
02429                          P_OUT_SPG[Chan_list[chan_n]->basin] += P_fl;
02430                          S_OUT_SPG[Chan_list[chan_n]->basin] += S_fl;
02431                      }
02432                      if ( !basn_list[basn[i]]->parent ) {                 
02433                          /* then find out about the flow for the family's sake */
02434                          VOL_IN_SPG[basn[i]] += fluxL; 
02435                          P_IN_SPG[basn[i]] += P_fl; 
02436                          S_IN_SPG[basn[i]] += S_fl; 
02437                      }
02438                  }
02439                  }
02440                      if ( tot_head > GW_head ) { /* to cell surface water */
02441                                 NA[i] += N_fl;
02442                                 PA[i] += P_fl;
02443                                 SA[i] += S_fl;
02444                                 }
02445                          else { /* to cell groundwater */
02446                                 GNA[i] += N_fl;
02447                                 GPA[i] += P_fl;
02448                                 GSA[i] += S_fl;
02449                         }
02450              }
02451              else if (fluxL < 0) /*  seepage flux from cell to canal */
02452              {
02453                  if ( tot_head > GW_head ) { /* cell's tot head > its groundwater head */
02454                         N_fl = ( SWH[i] > 0 ) ? NA[i]*fluxL/CELL_SIZE/SWH[i] : 0.;
02455                         P_fl = ( SWH[i] > 0 ) ? PA[i]*fluxL/CELL_SIZE/SWH[i] : 0.;
02456                         S_fl = ( SWH[i] > 0 ) ? SA[i]*fluxL/CELL_SIZE/SWH[i] : 0.;
02457                         NA[i] += N_fl;
02458                         PA[i] += P_fl;
02459                         SA[i] += S_fl;
02460                  }
02461                  else {
02462                         N_fl = ( GWH[i] > 0 ) ? GNA[i]*fluxL/CELL_SIZE/GWH[i] : 0.;
02463                         P_fl = ( GWH[i] > 0 ) ? GPA[i]*fluxL/CELL_SIZE/GWH[i] : 0.;
02464                         S_fl = ( GWH[i] > 0 ) ? GSA[i]*fluxL/CELL_SIZE/GWH[i] : 0.;
02465                         GNA[i] += N_fl;
02466                         GPA[i] += P_fl;
02467                         GSA[i] += S_fl;
02468                 }
02469                  if (Chan_list[chan_n]->basin != basn[i] && debug > -999) {
02470                      if (Chan_list[chan_n]->family !=  
02471                          basn_list[basn[i]]->family ) {        /* if the flow is not within the family... */
02472                          if  ( !basn_list[basn[i]]->parent  ) { /* and if the donor or recipient is a child... */
02473                          /* then find out about the flow for the family's sake */
02474                              VOL_OUT_SPG[basn_list[basn[i]]->family] -= fluxL;
02475                              P_OUT_SPG[basn_list[basn[i]]->family] -= P_fl;
02476                              S_OUT_SPG[basn_list[basn[i]]->family] -= S_fl;
02477                          }
02478                          if ( !Chan_list[chan_n]->parent ) {                 
02479                            /* then find out about the flow for the family's sake */
02480                              VOL_IN_SPG[Chan_list[chan_n]->family] -= fluxL; 
02481                              P_IN_SPG[Chan_list[chan_n]->family] -= P_fl; 
02482                              S_IN_SPG[Chan_list[chan_n]->family] -= S_fl; 
02483                          }
02484                      VOL_IN_SPG[Chan_list[chan_n]->basin] -= fluxL;
02485                      VOL_OUT_SPG[basn[i]] -= fluxL;
02486                      P_IN_SPG[Chan_list[chan_n]->basin] -= P_fl;
02487                      P_OUT_SPG[basn[i]] -= P_fl;
02488                      S_IN_SPG[Chan_list[chan_n]->basin] -= S_fl;
02489                      S_OUT_SPG[basn[i]] -= S_fl;
02490                  }
02491                  else {    /* if it's flow within a family, just keep
02492                           track of what the children do among themselves */            
02493                      if  ( !basn_list[basn[i]]->parent  ) { 
02494                          /* then find out about the flow for the family's sake */
02495                          VOL_OUT_SPG[basn[i]] -= fluxL;
02496                          P_OUT_SPG[basn[i]] -= P_fl;
02497                          S_OUT_SPG[basn[i]] -= S_fl;
02498                      }
02499                      if ( !Chan_list[chan_n]->parent ) {                 
02500                          /* then find out about the flow for the family's sake */
02501                          VOL_IN_SPG[Chan_list[chan_n]->basin] -= fluxL; 
02502                          P_IN_SPG[Chan_list[chan_n]->basin] -= P_fl; 
02503                          S_IN_SPG[Chan_list[chan_n]->basin] -= S_fl; 
02504                      }
02505                  }
02506                  }
02507              } 
02508                  /* now pass that mass flux to/from the canal */
02509              Chan_list[chan_n]->N -= N_fl;
02510              Chan_list[chan_n]->P -= P_fl;
02511              Chan_list[chan_n]->S -= S_fl;
02512              }
02513                                 /* now we recalculate the depths of SFWater and Gwater in cell */
02514 /*              SWH[i] += (fluxS+satTOsf_rec)/CELL_SIZE;  */
02515 /*              GWH[i] += (fluxG-satTOsf_rec)/CELL_SIZE;  */
02516              SWH[i] += (fluxS)/CELL_SIZE; 
02517              GWH[i] += (fluxG)/CELL_SIZE; 
02518              if ( tot_head > GW_head ) SWH[i] += fluxL/CELL_SIZE; /* whether to add seepage to surface or */
02519              else GWH[i] += fluxL/CELL_SIZE;                      /* groundwater is not critical, as vertical updates will follow */
02520                                 
02521 
02522          } /* that's all we do when finally converged */
02523                     
02524          T_flux_S += fluxS;
02525          T_flux_G += fluxG;
02526          T_flux_L += fluxL;
02527                   
02528          cell = cell->next_cell;
02529      } /*end cycle along canal cells */
02530 
02531          
02532      error = (CanWatDep - Chan_list[chan_n]->wat_depth) - (Qin - Qout - T_flux_S - T_flux_G - T_flux_L)/ Chan_list[chan_n]->area;
02533      if (converged && error != errorConv ) {
02534          /* error the second time around should be close to that of the first */
02535          if (fabs(error-errorConv)>F_ERROR ) { 
02536              sprintf(msgStr,"Day %6.1f: ERROR - converg error != to error after cell-canal fluxes (Chan %d, est= %5.3f, last depth= %5.3f m)", 
02537                         SimTime.TIME, Chan_list[chan_n]->number, CanWatDep, Chan_list[chan_n]->wat_depth); 
02538              WriteMsg( msgStr,True );  dynERRORnum++;
02539          }
02540          
02541          if (CanWatDep>Chan_list[chan_n]->depth*3.0 || CanWatDep < 0.0) {
02542                 CanWatDep = MINDEPTH; /* lost the mass balance, but let it go for now */
02543                 sprintf(msgStr,"           Reset Chan %d depth to minimum allowed (%5.3f m)", 
02544                                 Chan_list[chan_n]->number, MINDEPTH); 
02545                 WriteMsg( msgStr,True );  }
02546                  
02547                  /*errorConv is the error upon first reaching convergence criteria; */
02548              /* we then recalc the fluxes in the same loop thru cells (w/ same depth est.), */
02549              /* and should have the same fluxes and error as before */
02550          if (debug>4) { /* this here to halt program with high debug */
02551              printf("\nDay %6.1f: ERROR - converg error != to error after cell-canal fluxes (Chan %d, est= %5.3f m)\n", 
02552                         SimTime.TIME, Chan_list[chan_n]->number, CanWatDep);
02553              exit (-1);
02554 
02555          }
02556          
02557      }
02558      
02559      ++iter;
02560 
02561  } while (!converged );  /* end relaxation procedure */
02562 
02563  if ( factor != 0 && Qout > 0) /* if relaxation was abnormally stopped because of hitting the bottom */
02564  { 
02565      sprintf(msgStr,"Day %6.1f: ERROR - hit bottom of Canal %d, stopped the iterative relaxation", SimTime.TIME, Chan_list[chan_n]->number); 
02566      WriteMsg( msgStr,True );     
02567      dynERRORnum++;         
02568             /*#######*******########  this old code should not be reached
02569              the code is questionable - an ERROR message printed
02570              to file will indicate signif problem due to entrance to this */
02571      Qout1 = (Chan_list[chan_n]->wat_depth - CanWatDep)*Chan_list[chan_n]->area + Qin - T_flux_S - T_flux_G - T_flux_L;
02572      structs = struct_first;
02573      while ( structs != NULL )   /* recalculate the outflow in the data for possible use in the next canal inflow */
02574      {  
02575          if ( structs->canal_fr  == chan_n+CHo )   structs->flow *= Qout1/Qout;
02576     
02577          structs = structs->next_in_list;
02578      }
02579      Qout = Qout1;
02580  } 
02581 
02582  Chan_list[chan_n]->wat_depth = CanWatDep;      /* new canal water depth */
02583  CH_vol =  CanWatDep * Chan_list[chan_n]->area; /* new canal volume */
02584 
02585  if( debug > 3 ) 
02586  {
02587          /* concentration in mg/L (kg/m3==g/L, convert P (not S) to mg/L) */
02588      Chan_list[chan_n]->P_con = (CH_vol > 0) ? (Chan_list[chan_n]->P/CH_vol*1000.0) : 0.0;
02589      Chan_list[chan_n]->S_con = (CH_vol > 0) ? (Chan_list[chan_n]->S/CH_vol) : 0.0;
02590      fprintf( canDebugFile, 
02591               "%4d  %8.3f  %8.4f  %8.4f  %12.1f  %12.1f  %12.1f  %12.1f  %12.1f  %5d  %9.1f  %7.5f \n",
02592               Chan_list[chan_n]->number, Chan_list[chan_n]->wat_depth, Chan_list[chan_n]->P_con,  
02593               Chan_list[chan_n]->S_con, T_flux_S, T_flux_G, T_flux_L, Qin, Qout, iter, Chan_list[chan_n]->area, error );
02594  
02595  }
02596 
02597  return;
02598 }

Here is the call graph for this function:

void Flows_in_Structures float *  SWH,
float *  Elev,
float *  MC,
double *  NA,
double *  PA,
double *  SA
 

Determine flows through water control structures.

Parameters:
SWH SURFACE_WAT
Elev SED_ELEV
MC HYD_MANNINGS_N
NA DINdummy
PA TP_SF_WT
SA SALT_SURF_WT
Remarks:
Database (Structs_attr) formatting notes/requirements:
  • 1) No particular sequence is necessary, but the sequence of the structures in the CanalData.struct and CanalData.struct_wat files must be the same (and the prog will catch and whip you if you try otherwise).
  • 2) Due to the cycle that checks for multiple historical/other-data outflows from a canal, be sure that the database (to CanalData.struct) does not erroneously have 2 types (canal and cell) of recipients for a given structure.
  • 3) There need to be 2 header lines (as described in the dbase: run-name and col IDs).
  • 4) Flow is expected to ALWAYS be positive, even though there remain checks for reverse flows in the code. For a two-way structure, designate a separate structure for the reverse flow.
  • 5) As labeled in dbase, the concentration of P and N in structure flows are ug/L (ppb), and S is in g/L (ppt). (The code converts and deals with kg, m^3, and kg/m^3 (==g/L) ).
  • 6) The water flows in the input data are cfs, converted and used here as m^3/d.
  • 7) For rule based structures, don't forget to provide the cell location of the structure itself so that the elevation at a canal-canal flow structure is available.
  • 8) ******* You MAY NOT have multiple historical/other-data demands from a particular cell - only canals are allowed such multiple demands.
  • 9) If you want to apply a concentration to an internal structure that is one of multiple outflows from a canal, you need to apply that conc to all structures (actually the first in the sequence of the database) draining that canal. (This would only be for tracer studies etc., all internal model structure concentrations are calc'd).
Flow to/from the on-map, within-model-domain area, can only occur by a flow to/from a "cell" mapped at (1,1). That's based on our standard (0,0) upper left cell: the code uses a (1,1) origin - so that (1,1) from a map => (2,2) in the code!

Variables local to function
cell_to, cell_fr, can_to, can_fr The cell or canal the flow goes to, or comes from
cell_struct The grid cell location of the water control structure
TWf, HWf The flag to indicate a schedule for the TailWater, HeadWater
HeadH, HeadT, HeadH2, HeadT2 Head in Headwater or Tailwater
CH_vol, cell_vol Canal, cell total water volume
Nconc, Pconc, Sconc Nitrogen (UNUSED), Phosphorus, Salt/tracer concentration
N_fl, P_fl, S_fl Nitrogen (UNUSED), Phosphorus, Salt/tracer flows
chan_over Ratio of canal volume to historical/other-data demands
flow_rev Flow reversal check not invoked in current version
ChanHistOut Sum of the multiple historical/other-data outflows from canal
HeadH_drop, HeadT_drop Elevation difference between begin & end in (Headwater, Tailwater) canals
FFlux Flow (m) calculated between two cells
Note:
v2.3note re-activated used of head and tail-water stage schedules, for tidal boundary conditions

Definition at line 2675 of file WatMgmt.c.

References Abs, Chan::area, HistData::arrayN, HistData::arrayP, HistData::arrayS, HistData::arrayWat, Chan::basin, basn, basn_list, boundcond_depth, Structure::canal_fr, Structure::canal_to, canDebugFile, canstep, Structure::cell_i_fr, Structure::cell_i_to, Structure::cell_i_TW, Structure::cell_j_fr, Structure::cell_j_to, Structure::cell_j_TW, CELL_SIZE, Chan_list, Structure::conc_N, Structure::conc_P, Structure::conc_S, debug, Chan::depth, dynERRORnum, Chan::family, basndef::family, Structure::flag, Structure::flow, Flux_SWcells(), Flux_SWstuff(), FMOD(), GetGraph(), GP_SLRise, Hist_data, Structure::histID, Structure::HW_graph, Structure::HW_stage, Max, MCopen, Min, MINDEPTH, msgStr, Structure::multiOut, Chan::N, Structure::next_in_list, Chan::number, Chan::P, P_IN_STR, P_OUT_STR, ramp, Chan::S, S_IN_STR, Structure::S_nam, S_OUT_STR, SimTime, Structure::str_cell_i, Structure::str_cell_j, structs, Structure::Sum_P, Structure::Sum_S, Structure::SumFlow, Chan::sumHistIn, Chan::sumHistOut, Chan::sumRuleIn, Chan::sumRuleOut, T, simTime::TIME, Structure::TN, Structure::TNser, Structure::TP, Structure::TPser, True, Structure::TS, Structure::TSser, Structure::TW_graph, Structure::TW_stage, VOL_IN_STR, VOL_OUT_STR, Structure::w_coef, Chan::wat_depth, WriteMsg(), WstructOutFile, WstructOutFile_P, and WstructOutFile_S.

Referenced by Run_Canal_Network().

02676 {
02691         struct Structure *structs, *struct_last, *struct_hist_start;
02692     int cell_to, cell_fr, can_to, can_fr;
02693     int cell_struct;
02694     int TWf, HWf;
02695     double HeadH, HeadT, HeadH2, HeadT2;
02696     double CH_vol, cell_vol;
02697     double Nconc, Pconc, Sconc;
02698     double N_fl, P_fl, S_fl; 
02699     double chan_over;
02700     double flow_rev; 
02701     double ChanHistOut;
02702     double HeadH_drop, HeadT_drop;
02703     float FFlux; 
02704 
02705     /* from unitmod.h */
02706     extern float* boundcond_depth; /* depth in external cells */
02707 
02708 
02709 
02710     TWf = 0; HWf = 0; /* 0 when there is no schedule graphs for tail or headwater */
02711     
02712     structs = struct_first; 
02713 
02714     while ( structs != NULL )
02715     {           
02716         if ( structs->flag < 0 || (structs->flag >1) ) {
02717             structs->flow = 0.0; /* used in summing up flows for Qin and Qout later */
02718             structs->conc_P = 0.0;
02719             structs->conc_N = 0.0;
02720             structs->conc_S = 0.0;
02721             goto OUT; 
02722         } /* nothing to do when structure inactive (neg or >1 flag) */
02723         
02725         if ( structs->TW_stage == -99 ) /* this means that there is a graph for TW schedules */
02726         { 
02727         structs->TW_stage = GetGraph ( structs->TW_graph, ( FMOD(SimTime.TIME,365.0) >0.0 ) ? ( FMOD(SimTime.TIME,365.0) ) : ( 0.0) ); 
02728         /* add sea level rise */
02729         structs->TW_stage = structs->TW_stage + GP_SLRise * (SimTime.TIME/365);
02730         TWf = 1; 
02731         }
02732           
02733         if ( structs->HW_stage == -99 ) /* this means that there is a graph for HW schedules */
02734         { 
02735         structs->HW_stage = GetGraph ( structs->HW_graph, ( FMOD(SimTime.TIME,365.0) >0.0 ) ? ( FMOD(SimTime.TIME,365.0) ) : ( 365.0) ); 
02736         /* add sea level rise */
02737         structs->HW_stage = structs->HW_stage + GP_SLRise * (SimTime.TIME/365);
02738         HWf = 1; }
02739         
02740           
02741 /* CASE 1 */
02742 /*check that this is a canal-to-canal structure */
02743 /*  these flows are always internal to model domain */
02744         if ( structs->canal_fr  != 0 && structs->canal_to != 0 ) {
02745             can_fr = structs->canal_fr - CHo;
02746             can_to = structs->canal_to - CHo;
02747             cell_struct = T(structs->str_cell_i,structs->str_cell_j); 
02748                 
02749                 /* first check for historical/sfwmm flow */
02750             if ( structs->flag == 1 ) {
02751                 if (structs->multiOut == 1) goto OUT; /* this structure has already been processed as a multiple outflow from a canal */
02752                     /* FLOW */
02753                     /* don't do much with zero or negative flows */           
02754                 structs->flow = Hist_data[structs->histID-1].arrayWat[((int)(SimTime.TIME))]*canstep;
02755                 if (structs->flow == 0.0) { 
02756                     structs->conc_P = 0.0;
02757                     structs->conc_N = 0.0;
02758                     structs->conc_S = 0.0;
02759                     goto OUT; 
02760                 }   /* zero the reported struct flow conc (for output info only) */
02761                 else if (structs->flow < 0.0) {
02762                     sprintf( msgStr, "Day %6.1f: ERROR IN ASSUMPTION OF NO NEG HIST FLOW from chan %d",
02763                              SimTime.TIME, Chan_list[can_to]->number );
02764                     WriteMsg( msgStr,True );
02765                     dynERRORnum++;
02766                     goto OUT;
02767                 }
02768 
02769                     /* for multiple historical/other-data outflows from a canal, ensure mass balance (sum of outflows !> avail volume) */
02770                     /* all flows from the canal associated with each set must be non-negative, which is currently valid */
02771                     /* the below cycle thru the structures will also process a canal->cell flow structure */
02772 
02773                 struct_hist_start = structs; /* this is the first structure w/ historical info for a particular canal */
02774                 ChanHistOut = 0.0; /* sum of all historical/other-data outflows from a canal */
02775 
02776                 while ( structs != NULL )   
02777                 {  
02778                     if  ( structs->flag == 1 && struct_hist_start->canal_fr  == structs->canal_fr   ) { /* look for outflows from same canal */
02779                         structs->flow = Hist_data[structs->histID-1].arrayWat[((int)(SimTime.TIME))]*canstep;
02780                         ChanHistOut += structs->flow; /* sum the multiple historical/other-data outflows from canal */
02781                     }
02782                         
02783                     structs = structs->next_in_list; /* increment to next structure */
02784                 } /* end cycle thru list of remaining structures in total list */
02785                     
02786                 structs = struct_hist_start; /* now we are pointing back to the first (or only) water control structure draining the canal */
02787                          
02788                     /* concentration in donor canal for structure flow, g/L = kg/m^3 */
02789                     /* this is based on previous depth/volume calc'd in FluxChannel */
02790                     /* or use historical conc */
02791 
02792                     /* CONC */
02793                                 /* canal TOTAL volume before new structure flows */
02794                 CH_vol = Chan_list[can_fr]->area*Chan_list[can_fr]->wat_depth; 
02795                     /* concentrations applied to all structures draining this canal */
02796                 Pconc = (structs->TPser == -1) ? 
02797                     ( (CH_vol>0) ? (Chan_list[can_fr]->P/CH_vol) : (0.0) ) :
02798                     ( (structs->TPser == 0) ? (structs->TP) : ( Hist_data[structs->histID-1].arrayP[((int)(SimTime.TIME))] ) );
02799                 Nconc = (structs->TNser == -1) ? 
02800                     ( (CH_vol>0) ? (Chan_list[can_fr]->N/CH_vol) : (0.0) ) :
02801                     ( (structs->TNser == 0) ? (structs->TN) : ( Hist_data[structs->histID-1].arrayN[((int)(SimTime.TIME))] ) );
02802                 Sconc = (structs->TSser == -1) ? 
02803                     ( (CH_vol>0) ? (Chan_list[can_fr]->S/CH_vol) : (0.0) ) :
02804                     ( (structs->TSser == 0) ? (structs->TS) : ( Hist_data[structs->histID-1].arrayS[((int)(SimTime.TIME))] ) );
02805 
02806 
02807                     /* vol avail for fluxing */
02808                 CH_vol =  Chan_list[can_fr]->area*(ramp(Chan_list[can_fr]->wat_depth-MINDEPTH) );
02809                 chan_over = (ChanHistOut > 0.0 ) ? ( CH_vol/ChanHistOut ) : 1.0;
02810 
02811                 do /* loop again through the remaining sequence of structures to correct any outflows for the canal */
02812                 {   
02813                     if (structs->flag == 1 &&  struct_hist_start->canal_fr  == structs->canal_fr   ) { 
02814 
02815                         structs->multiOut = 1; /* set the flag to indicate this structure has been processed */
02816                             
02817                             /* revised FLOW */
02818                         if ( chan_over < 1.0 ) { /* if hist/other-data demand > canal vol */
02819                             if( debug > 0 ) { 
02820                                 fprintf( canDebugFile, "Day \t%6.1f\t %s:\tTotDemand= \t%7.0f\t; exceeded chan \t%d\t avail vol= \t%7.0f \tby \t%7.0f \tm3 (\t%5.3f \tm deep).\n",
02821                                          SimTime.TIME, structs->S_nam, ChanHistOut, Chan_list[can_fr]->number,
02822                                          CH_vol, (ChanHistOut - CH_vol), Chan_list[can_fr]->wat_depth );
02823                                 
02824                             }  
02825                                 /* decrement the individual structure outflows by it's proportional contrib. */
02826                             structs->flow = (structs->flow) * chan_over ;
02827                         } 
02828                             /* sum of all realized (corrected for excess demand) historical/other-data outflows from a canal during an iteration */
02829                             /* used to determine any allowable rule-based canal flows after the historical/other-data demands */
02830                         Chan_list[can_fr]->sumHistOut += structs->flow;         
02831 
02832                             /* the recipient changes as we cycle thru the structures (may be another canal or a cell) */
02833                         can_to = structs->canal_to - CHo;
02834                         cell_to = T(structs->cell_i_to,structs->cell_j_to); 
02835 
02836                         structs->conc_P = Pconc;/* concentrations applied to all structures draining this canal */
02837                         structs->conc_N = Nconc;
02838                         structs->conc_S = (structs->TSser == -1) ?
02839                             (Sconc) :
02840                             ( (structs->TSser == 0) ? (structs->TS) :
02841                               ( Hist_data[structs->histID-1].arrayS[((int)(SimTime.TIME))] ) ); /* see note below on tracer conc. application */
02842                         
02843                         P_fl = structs->flow * structs->conc_P; 
02844                         N_fl = structs->flow * structs->conc_N;
02845                         S_fl = structs->flow * structs->conc_S;
02846                                 /* change in nut mass in fr canal (kg) */
02847                         Chan_list[can_fr]->N -=  N_fl ;
02848                         Chan_list[can_fr]->P -=  P_fl ;
02849                         Chan_list[can_fr]->S -=  (structs->TSser == -1) ? (S_fl) : (0);
02850                             /* for the salt (S), we sometimes want to apply a tracer to internal canals/cells.
02851                                Thus, in that case, we don't subtract salt from the canal because this is "introduced" salt
02852                                that was not calculated from the available water and salt (i.e., the conc. was from the
02853                                CanalData.struct data file, either fixed or a time series) */
02854                             
02855                         if (can_to > 0) {    /* IF this is a canal recip, change in nut mass in to canal (kg) */
02856                                 /* sum of all realized (corrected for excess demand) historical/other-data outflows into a canal during an iteration */
02857                                 /* used to determine any allowable rule-based canal flows after the historical/other-data demands */
02858                             Chan_list[can_to]->sumHistIn += structs->flow;              
02859 
02860                             Chan_list[can_to]->N +=  N_fl ;
02861                             Chan_list[can_to]->P +=  P_fl ;
02862                             Chan_list[can_to]->S +=  S_fl ;
02863                                 /* mass balance in fr and to canals */
02864                                 /* canals themselves always belong to a parent hydro basin,
02865                                    thus canal-canal struct flows do not keep track of child-parent, child-child flows */
02866                             if (Chan_list[can_fr]->basin != Chan_list[can_to]->basin && debug > -999 ) {
02867                                 VOL_OUT_STR[Chan_list[can_fr]->basin] += structs->flow;
02868                                 VOL_IN_STR[Chan_list[can_to]->basin] += structs->flow; 
02869                                 P_OUT_STR[Chan_list[can_fr]->basin] += P_fl;
02870                                 P_IN_STR[Chan_list[can_to]->basin] += P_fl; 
02871                                 S_OUT_STR[Chan_list[can_fr]->basin] += (structs->TSser == -1) ? (S_fl) : (0); /* see note above on tracer conc. application */
02872                                 S_IN_STR[Chan_list [can_to]->basin] += S_fl; 
02873                             }
02874 
02875 
02876                                              
02877                         }
02878                         else if (cell_to > 0) {
02879                                 /* OR change in nut mass in downstream internal cell (kg) */
02880  
02881                             if (structs->cell_i_to > 2 ) { /* to an on-map cell */
02882                                 SWH[cell_to] += structs->flow/CELL_SIZE;
02883                                 PA[cell_to] += P_fl; 
02884                                 NA[cell_to] += N_fl; 
02885                                 SA[cell_to] += S_fl; 
02886                                     /* mass balance in fr canal and to cell */
02887                                 if (Chan_list[can_fr]->basin != basn[cell_to] && debug > -999 ) {
02888                                     VOL_OUT_STR[Chan_list[can_fr]->basin] += structs->flow;
02889                                     VOL_IN_STR[basn[cell_to]] += structs->flow; 
02890                                     P_OUT_STR[Chan_list[can_fr]->basin] += P_fl;
02891                                     P_IN_STR[basn[cell_to]] += P_fl; 
02892                                     S_OUT_STR[Chan_list[can_fr]->basin] += (structs->TSser == -1) ? (S_fl) : (0); /* see note above on tracer conc. application */
02893                                     S_IN_STR[basn[cell_to]] += S_fl; 
02894 
02895                                         /* all structure flows, EXCEPT canal-canal, must be between parent (family) hydro basins */
02896                                     if ( (Chan_list[can_fr]->family ==  basn_list[basn[cell_to]]->family) && (debug > 0) ) {
02897                                         sprintf( msgStr, "Day %6.1f: ERROR in basin configuration for structure from chan %d to cell (%d,%d): flow must be between diff hydrologic basins",
02898                                                  SimTime.TIME, Chan_list[can_fr]->number,structs->cell_i_to,structs->cell_j_to );
02899                                         WriteMsg( msgStr,True );
02900                                         dynERRORnum++;
02901                                     }
02902                                 }
02903                                 
02904 
02905 
02906                             }
02907                             else if (structs->cell_i_to == 2 && debug > -999 ) { /* to off-map cell, mass balance check */
02908                                 VOL_OUT_STR[Chan_list[can_fr]->basin] += structs->flow; 
02909                                 VOL_OUT_STR[0] += structs->flow; 
02910                                 P_OUT_STR[Chan_list[can_fr]->basin] += P_fl; 
02911                                 P_OUT_STR[0] += P_fl; 
02912                                 S_OUT_STR[Chan_list[can_fr]->basin] += (structs->TSser == -1) ? (S_fl) : (0); /* see note above on tracer conc. application */
02913                                 S_OUT_STR[0] += (structs->TSser == -1) ? (S_fl) : (0); /* see note above on tracer conc. application */
02914                             }
02915                         }
02916                             /* don't do much with zero or negative flows */           
02917                         if (structs->flow == 0.0) {
02918                             structs->conc_P = 0.0;
02919                             structs->conc_N = 0.0;
02920                             structs->conc_S = 0.0;
02921                         } /* for output purposes only, zero the reported conc */
02922                         else if (structs->flow < 0.0) {
02923                             sprintf( msgStr, "Day %6.1f: ERROR IN ASSUMPTION OF NO NEG HIST FLOW from another cell/chan into chan %d",
02924                                      SimTime.TIME,Chan_list[can_fr]->number );
02925                             WriteMsg( msgStr,True );
02926                             dynERRORnum++;
02927                         }
02928                     }
02929 
02930                     structs = structs->next_in_list;
02931                 } while (structs != NULL); /*  cycling to the end of the list again !!!! */
02932 
02933                 structs = struct_hist_start; /* now pointing back to the first historical outflow structure that we dealt with */
02934 
02935                 goto OUT;
02936             } /* finished check for multiple (positive) outflows exceeding canal volume */
02937 
02938             else  { /* Rule- based flow calcs */
02939 
02940                 HeadH_drop = 0.5 * Abs(Chan_list[can_fr]->elev_drop); /* elevation difference between begin&end of canal */
02941                 HeadT_drop = 0.5 * Abs(Chan_list[can_to]->elev_drop);
02942                 
02943                 HeadH = HeadH_drop + /* add portion of the elevation drop along the canal length */
02944                     Elev[cell_struct] - Chan_list[can_fr]->depth + Chan_list[can_fr]->wat_depth 
02945                     + (Chan_list[can_fr]->sumHistIn - Chan_list[can_fr]->sumHistOut)/Chan_list[can_fr]->area; 
02946                                 /* in both head and tail, add net historical/other-data flows before determining hydraulic potential */
02947                     /* in tailwater only, check to see if other virtual struct has added water already, add to head */
02948                 HeadT = -HeadT_drop + /* subtract portion of the elevation drop along the canal length */
02949                     Elev [cell_struct] - Chan_list[can_to]->depth + Chan_list[can_to]->wat_depth
02950                     + ( Chan_list[can_to]->sumRuleIn + Chan_list[can_to]->sumHistIn - Chan_list[can_to]->sumHistOut)/Chan_list[can_to]->area;
02951 
02952                     /* FLOW */
02953                 if ( structs->w_coef == 0.0)  {/* check for a virtual structure, flow only in pos direction */
02954                     if ( HeadH>HeadT ) {         
02955                         structs->flow =
02956                             Chan_list[can_fr]->area*Chan_list[can_to]->area/
02957                             (Chan_list[can_fr]->area+Chan_list[can_to]->area) * (HeadH-HeadT) ;
02958                         if (debug >3) {
02959                             HeadH2 = HeadH - structs->flow/Chan_list[can_fr]->area; 
02960                             HeadT2 = HeadT + structs->flow/Chan_list[can_to]->area;
02961                             sprintf( msgStr, "Day %6.1f:  virtual flow from chan %d (%f m) to chan %d, New headH-headT=%f, Orig HeadH-HeadT=%f",
02962                                      SimTime.TIME,Chan_list[can_fr]->number,HeadH-Elev[cell_struct]+Chan_list[can_fr]->depth,
02963                                      Chan_list[can_to]->number,(HeadH2-HeadT2),(HeadH-HeadT) );
02964                             WriteMsg( msgStr,True ); 
02965                             if (Abs(HeadH2 - HeadT2) > GP_MinCheck) /* virtual flows should exactly equilibrate the two canals */
02966                             { sprintf( msgStr, "Day %6.1f: ERROR virtual flow from chan %d to chan %d, New headH-headT=%f, Orig HeadH-HeadT=%f",
02967                                        SimTime.TIME,Chan_list[can_fr]->number,Chan_list[can_to]->number,(HeadH2-HeadT2),(HeadH-HeadT) );
02968                             WriteMsg( msgStr,True ); dynERRORnum++; }
02969 
02970                         }
02971                         
02972                         CH_vol =  Chan_list[can_fr]->area* ramp((HeadH-HeadH_drop)-Elev[cell_struct]+Chan_list[can_fr]->depth-MINDEPTH); /* avail flow volume */
02973                         if (structs->flow > CH_vol) {
02974                             structs->flow =  CH_vol;
02975                             HeadH2 = HeadH - structs->flow/Chan_list[can_fr]->area; 
02976                             HeadT2 = HeadT + structs->flow/Chan_list[can_to]->area;
02977                             if (debug>3 && Abs(HeadH2 - HeadT2) > GP_MinCheck) {
02978                                 sprintf( msgStr, "Day %6.1f: Warning: stage didn't equilibrate in modified virtual struct flow from chan %d (%f m) to chan %d, New headH-headT=%f, Orig HeadH-HeadT=%f",
02979                                          SimTime.TIME,Chan_list[can_fr]->number,CH_vol/Chan_list[can_fr]->area,Chan_list[can_to]->number,(HeadH2-HeadT2),(HeadH-HeadT) );
02980                                 WriteMsg( msgStr,True ); }
02981                         }
02982                         
02983                             /* revise volume to be TOTAL volume for conc calc */
02984                         CH_vol = CH_vol + MINDEPTH*Chan_list[can_fr]->area;
02985                             /* conc in donor (fr) canal (kg/m3) */
02986                         structs->conc_P = (CH_vol>0) ? (Chan_list[can_fr]->P/CH_vol) : (0.0);
02987                         structs->conc_N = (CH_vol>0) ? (Chan_list[can_fr]->N/CH_vol) : (0.0);
02988                         structs->conc_S = (CH_vol>0) ? (Chan_list[can_fr]->S/CH_vol) : (0.0);
02989                     }
02990                     else  structs->flow = 0.0;
02991                     
02992                 }
02993                 
02994                 else
02995                 {
02996                     printf ("\nVirtual structs are the only rule-based canal-canal flows allowed in this version!\n"); exit(-1);
02997                         /* the code that was here is in older version,
02998                            not fully implemented (ELM driven by historical/other-data flows, only canal-canal virtual structs allowed) */
02999                     
03000                 } /* end of rule based flow calcs */
03001  
03002                     /* don't do much with zero or negative rule-based flows  */           
03003                 if (structs->flow == 0.0) {
03004                     structs->conc_P = 0.0;
03005                     structs->conc_N = 0.0;
03006                     structs->conc_S = 0.0;
03007                     goto OUT;
03008                 } /* for output purposes only, zero the reported conc */
03009                 else if (structs->flow < 0.0 && structs->w_coef != 0) {                 
03010                     sprintf( msgStr, "Day %6.1f: ERROR IN ASSUMPTION OF NO NEG RULE-BASED FLOW from chan %d",
03011                              SimTime.TIME,Chan_list[can_to]->number);
03012                     WriteMsg( msgStr,True );
03013                     goto OUT;
03014                 } 
03015                                                                 
03016 
03017                     /* increment the sum of rule-based inflows/outflows to canals
03018                        for potential multiple virtual inflows into some canals */
03019                 Chan_list[can_fr]->sumRuleOut += structs->flow;         
03020                 Chan_list[can_to]->sumRuleIn += structs->flow;          
03021 
03022                     /* nut flows */
03023                 P_fl = structs->flow * structs->conc_P; 
03024                 N_fl = structs->flow * structs->conc_N;
03025                 S_fl = structs->flow * structs->conc_S;
03026                     /* change in nut mass in fr canal (kg) */
03027                 Chan_list[can_fr]->N -=  N_fl ;
03028                 Chan_list[can_fr]->P -=  P_fl ;
03029                 Chan_list[can_fr]->S -=  S_fl ;
03030                     /* change in nut mass in to canal (kg) */
03031                 Chan_list[can_to]->N +=  N_fl ;
03032                 Chan_list[can_to]->P +=  P_fl ;
03033                 Chan_list[can_to]->S +=  S_fl ;
03034                     /* mass balance in fr and to canals */
03035                     /* virtual structure flows can not (by definition) establish flow among basins */
03036                 if (Chan_list[can_fr]->basin != Chan_list[can_to]->basin && debug > 0 ) {
03037                     sprintf( msgStr, "Day %6.1f: ERROR - no hydrologic basin-basin flows in virtual structure, canals %d - %d",
03038                              SimTime.TIME, Chan_list[can_to]->number, Chan_list[can_fr]->number);
03039                     WriteMsg( msgStr,True ); dynERRORnum++;
03040                 }
03041          
03042                 goto OUT;
03043             }
03044         }
03045 
03046 /* CASE 2 */
03047 /*check that this is a structure between a cell and a cell*/
03048         if ( structs->cell_i_to != 0 && structs->cell_j_to != 0  
03049              && structs->cell_i_fr != 0 && structs->cell_j_fr != 0 ) { 
03050 
03051             cell_to = T(structs->cell_i_to,structs->cell_j_to); 
03052             cell_fr = T(structs->cell_i_fr,structs->cell_j_fr);
03053 
03054                 /* first check for historical (or other data) flow */
03055             if ( structs->flag == 1 ) {
03056                     /* FLOW */
03057                     /* don't do much with zero or negative flows */
03058                 structs->flow = Hist_data[structs->histID-1].arrayWat[((int)(SimTime.TIME))]*canstep;
03059                 if (structs->flow == 0.0) {
03060                     structs->conc_P = 0.0;
03061                     structs->conc_N = 0.0;
03062                     structs->conc_S = 0.0;
03063                     goto OUT;    /* for output purposes only, zero the reported conc */
03064                 }
03065                 else if (structs->flow < 0.0) { /* all flows should be positive, but check */
03066                     sprintf( msgStr, "Day %6.1f: ERROR IN ASSUMPTION OF NO NEG HIST FLOW from cell (%d,%d)",
03067                              SimTime.TIME, structs->cell_i_to, structs->cell_j_to);
03068                     WriteMsg( msgStr,True );
03069                     goto OUT;
03070                 }
03071                     /* for internal cells, check for excessive historical/model demand */
03072                 if (structs->cell_i_fr != 2 ) 
03073                 {
03074 
03075                     cell_vol = SWH[cell_fr] * CELL_SIZE;
03076                     if (structs->flow > cell_vol) { /* if hist/other-data demand > cell vol */
03077                         if( debug > 0 ) { 
03078                             fprintf( canDebugFile, "Day \t%6.1f\t %s:\tTotDemand= \t%7.0f\t; exceeded cell \t(%d,%d)\t vol= \t%7.0f \tby \t%7.0f\t m3.\n",
03079                                      SimTime.TIME, structs->S_nam, structs->flow, structs->cell_i_fr, structs->cell_j_fr, cell_vol, (structs->flow - cell_vol) );
03080                            
03081                         }
03082                             
03083                         structs->flow = cell_vol;
03084                     }}
03085                     /* CONCENTRATION */
03086                     /* now calculate the nutrient conc in the donor cell */
03087                     /* nut conc in flow from outside system using historical data or 0.0 */
03088                 if (structs->cell_i_fr == 2) { 
03089                     structs->conc_P = (structs->TPser == 1)  ?
03090                         Hist_data[structs->histID-1].arrayP[((int)(SimTime.TIME))] :
03091                         ((structs->TPser == 0) ? (structs->TP) : (0.0) );
03092                     structs->conc_N = (structs->TNser == 1)  ?
03093                         Hist_data[structs->histID-1].arrayN[((int)(SimTime.TIME))] :
03094                         ((structs->TNser == 0) ? (structs->TN) : (0.0) );
03095                     structs->conc_S = (structs->TSser == 1)  ?
03096                         Hist_data[structs->histID-1].arrayS[((int)(SimTime.TIME))] :
03097                         ((structs->TSser == 0) ? (structs->TS) : (0.0) );
03098                 }
03099                 else
03100                 { /* nut conc in flow within system - check to see if there are historical concs with the historical flow */
03101                     structs->conc_P = (structs->TPser == -1)  ?
03102                         ( (SWH[cell_fr] > 0) ? PA[cell_fr]/CELL_SIZE/SWH[cell_fr] : 0.0 ) :
03103                         ( (structs->TPser == 0) ? (structs->TP) : (Hist_data[structs->histID-1].arrayP[((int)(SimTime.TIME))] ) );
03104                     structs->conc_N = ( structs->TNser == -1)  ?
03105                         ( (SWH[cell_fr] > 0) ? NA[cell_fr]/CELL_SIZE/SWH[cell_fr] : 0.0 ) :
03106                         ( (structs->TNser == 0) ? (structs->TN) : (Hist_data[structs->histID-1].arrayN[((int)(SimTime.TIME))] ) );
03107                     structs->conc_S = ( structs->TSser == -1)  ?
03108                         ( (SWH[cell_fr] > 0) ? SA[cell_fr]/CELL_SIZE/SWH[cell_fr] : 0.0 ) :
03109                         ( (structs->TSser == 0) ? (structs->TS) : (Hist_data[structs->histID-1].arrayS[((int)(SimTime.TIME))] ) );
03110                 }
03111                     /* the nut and water update done after the rule-based flow calcs */
03112 /* determine mass transfer */
03113                     
03114                 P_fl = structs->flow * structs->conc_P;
03115                 N_fl = structs->flow * structs->conc_N;
03116                 S_fl = structs->flow * structs->conc_S;
03117                 
03118                 
03119                 if (structs->cell_i_to == 2) { /* "to" cell is outside system */
03120                     if (debug > -999) {/* m3 and kg nuts outflow out from the system */
03121                         VOL_OUT_STR[basn[cell_fr]] += structs->flow; /* mass balance */
03122                         VOL_OUT_STR[0] += structs->flow; 
03123                         P_OUT_STR[basn[cell_fr]] += P_fl; 
03124                         P_OUT_STR[0] += P_fl; 
03125                         S_OUT_STR[basn[cell_fr]] += (structs->TSser == -1) ? (S_fl) : (0); /* see below note on tracer application */
03126                         S_OUT_STR[0] += (structs->TSser == -1) ? (S_fl) : (0); /* see below note on tracer application */
03127                     } 
03128                     PA[cell_fr] -= P_fl;
03129                     NA[cell_fr] -= N_fl; 
03130                     SA[cell_fr] -= (structs->TSser == -1) ? (S_fl) : (0);
03131                         /* for the salt (S), we sometimes want to apply a tracer to internal canals/cells.
03132                            Thus, in that case, we don't subtract salt from the canal because this is "introduced" salt
03133                            that was not calculated from the available water and salt (i.e., the conc. was from the
03134                            CanalData.struct data file, either fixed or a time series) */
03135                         /* recalculate the surface water depth in the "from" cell */
03136                     SWH[cell_fr] -= structs->flow/CELL_SIZE;
03137                 }
03138                 else if (structs->cell_i_fr == 2) {  /* "from" cell is outside system */
03139                     if (debug > -999) {/* m3 and kg nuts inflow to the system */
03140                         VOL_IN_STR[basn[cell_to]] += structs->flow; 
03141                         VOL_IN_STR[0] += structs->flow; 
03142                         P_IN_STR[basn[cell_to]] += P_fl; 
03143                         P_IN_STR[0] += P_fl; 
03144                         S_IN_STR[basn[cell_to]] += S_fl; 
03145                         S_IN_STR[0] += S_fl; 
03146                     } 
03147                     PA[cell_to] += P_fl;
03148                     NA[cell_to] += N_fl;
03149                     SA[cell_to] += S_fl;
03150                         /* recalculate the surface water depth in the "to" cell */
03151                     SWH[cell_to] += structs->flow/CELL_SIZE;
03152                 }
03153                 else { /* internal cell-cell m3 and kg nuts flow */
03154                     PA[cell_fr] -= P_fl;
03155                     NA[cell_fr] -= N_fl;        
03156                     SA[cell_fr] -= (structs->TSser == -1) ? (S_fl) : (0);
03157                         /* for the salt (S), we sometimes want to apply a tracer to internal canals/cells.
03158                            Thus, in that case, we don't subtract salt from the cell because this is "introduced" salt
03159                            that was not calculated from the available water and salt (i.e., the conc. was from the
03160                            CanalData.struct data file, either fixed or a time series) */
03161                     
03162                     PA[cell_to] += P_fl;
03163                     NA[cell_to] += N_fl;
03164                     SA[cell_to] += S_fl;
03165                         /* recalculate the surface water depth in the interacting cells */
03166                     SWH[cell_fr] -= structs->flow/CELL_SIZE;
03167                     SWH[cell_to] += structs->flow/CELL_SIZE;
03168                     if (basn[cell_to] != basn[cell_fr] && debug > -999 ) {
03169                         VOL_OUT_STR[basn[cell_fr]] += structs->flow; 
03170                         VOL_IN_STR[basn[cell_to]] += structs->flow; 
03171                         P_OUT_STR[basn[cell_fr]] += P_fl; 
03172                         P_IN_STR[basn[cell_to]] += P_fl; 
03173                         S_OUT_STR[basn[cell_fr]] += (structs->TSser == -1) ? (S_fl) : (0); /* see above note on tracer application */
03174                         S_IN_STR[basn[cell_to]] += S_fl; 
03175 
03176                             /* all structure flows, EXCEPT canal-canal, must be between parent (family) hydro basins */
03177                         if ( (basn_list[basn[cell_fr]]->family ==  basn_list[basn[cell_to]]->family) && (debug > 0) ) {
03178                             sprintf( msgStr, "Day %6.1f: ERROR in basin configuration for structure from cell (%d,%d) to cell (%d,%d): flow must be between diff hydrologic basins",
03179                                      SimTime.TIME, structs->cell_i_fr,structs->cell_j_fr,structs->cell_i_to,structs->cell_j_to );
03180                             WriteMsg( msgStr,True ); dynERRORnum++;
03181                         }
03182                     }
03183                     
03184                 }
03185                 goto OUT;
03186             }
03187             
03188                 /**** Rule- based flows */
03189             else {  
03190                 if ( structs->w_coef != 0.0)  {/* check for a virtual structure  */
03191                     printf ("\nVirtual structs are the only rule-based cell-cell flows in this version!\n"); exit(-1);
03192                 }
03193                 else {
03194                         /* calc'd flow, cell-cell, using the functions defined in Fluxes.c for overland flow, bridge-flow only app right now */
03195                         /* MCopen is a (whole array) of Manning's coefs that are initialized at an open-water value */
03196                     FFlux = Flux_SWcells(structs->cell_i_fr,structs->cell_j_fr,
03197                                       structs->cell_i_to,structs->cell_j_to,SWH,Elev,MCopen);
03198                         /* flow may be negative */
03199                                 /* FFlux units = m */
03200                     structs->flow = FFlux*CELL_SIZE;
03201                     structs->conc_P = (SWH[cell_fr] > 0) ? PA[cell_fr]/CELL_SIZE/SWH[cell_fr] : 0.0 ;
03202 /*                     structs->conc_N = (SWH[cell_fr] > 0) ? NA[cell_fr]/CELL_SIZE/SWH[cell_fr] : 0.0 ; */
03203                     structs->conc_S = (SWH[cell_fr] > 0) ? SA[cell_fr]/CELL_SIZE/SWH[cell_fr] : 0.0 ;
03204 
03205                     if (FFlux != 0.0) Flux_SWstuff ( structs->cell_i_fr,structs->cell_j_fr,
03206                                                   structs->cell_i_to,structs->cell_j_to,FFlux,SWH,SA,NA,PA);
03207                     
03208                     
03209                 }
03210                 goto OUT;
03211             } /* end of rule-based flow calc */
03212             
03213         }
03214 
03215 /* CASE 3 */
03216 /*check that this is a structure between a canal and a cell*/
03217         if ( structs->canal_fr  != 0 && structs->cell_i_to != 0 && structs->cell_j_to != 0 ) { 
03218             can_fr = structs->canal_fr  - CHo;
03219             cell_to = T(structs->cell_i_to,structs->cell_j_to);           
03220             cell_struct = T(structs->str_cell_i,structs->str_cell_j);
03221 
03222                 /* first check for historical (or other data) flow */
03223             if ( structs->flag == 1 ) {
03224 
03225                 if (structs->multiOut == 1) goto OUT; /* this structure has already been processed as a multiple outflow from a canal */
03226                 
03227                 structs->flow = Hist_data[structs->histID-1].arrayWat[((int)(SimTime.TIME))]*canstep;
03228                     /* don't do much with zero or negative flows */           
03229                 if (structs->flow == 0.0) {/* zero the reported struct flow conc (for output info only) */
03230                     structs->conc_P = 0.0;
03231                     structs->conc_N = 0.0;
03232                     structs->conc_S = 0.0;
03233                     goto OUT; 
03234                 }
03235                 else if (structs->flow < 0.0) { 
03236                     sprintf( msgStr, "Day %6.1f: ERROR IN ASSUMPTION OF NO NEG HIST FLOW from cell (%d,%d)",
03237                              SimTime.TIME, structs->cell_i_to, structs->cell_j_to );
03238                     WriteMsg( msgStr,True ); dynERRORnum++;
03239                     goto OUT;
03240                 }
03241                 
03242                     /* for multiple historical/other-data outflows from a canal, ensure mass balance (sum of outflows !> avail volume) */
03243                     /* all flows from the canal associated with each set must be non-negative, which is currently valid */
03244 
03245                 struct_hist_start = structs; /* this is the first structure w/ historical info for a particular canal */
03246                 ChanHistOut = 0.0; /* sum of all historical outflows from a canal */
03247 
03248                 while ( structs != NULL )   
03249                 {  
03250                     if  ( structs->flag == 1 && struct_hist_start->canal_fr  == structs->canal_fr   ) { /* look for outflows from same canal */
03251                         structs->flow = Hist_data[structs->histID-1].arrayWat[((int)(SimTime.TIME))]*canstep;
03252                         ChanHistOut += structs->flow; /* sum the multiple historical outflows from canal */
03253                     }
03254                         
03255                     structs = structs->next_in_list; /* increment to next structure */
03256                 } /* end cycle thru list of remaining structures in total list */
03257                     
03258                 structs = struct_hist_start; /* now we are pointing back to the first (or only) water control structure draining the canal */                 
03259 
03260                     /* concentration in donor canal for structure flow, g/L = kg/m^3 */
03261                     /* this is based on previous depth/volume calc'd in FluxChannel */
03262                     /* or use historical conc */
03263                         /* canal TOTAL volume before new structure flows */
03264                 CH_vol = Chan_list[can_fr]->area*Chan_list[can_fr]->wat_depth; 
03265                     /* concentrations applied to all structures draining this canal */
03266                 Pconc = (structs->TPser == -1) ? 
03267                     ( (CH_vol>0) ? (Chan_list[can_fr]->P/CH_vol) : (0.0) ) :
03268                     ( (structs->TPser == 0) ? (structs->TP) : (Hist_data[structs->histID-1].arrayP[((int)(SimTime.TIME))] ) ) ;
03269                 Nconc = (structs->TNser == -1) ? 
03270                     ( (CH_vol>0) ? (Chan_list[can_fr]->N/CH_vol) : (0.0) ) :
03271                     ( (structs->TNser == 0) ? (structs->TN) : (Hist_data[structs->histID-1].arrayN[((int)(SimTime.TIME))] ) ) ;
03272                 Sconc = (structs->TSser == -1) ? 
03273                     ( (CH_vol>0) ? (Chan_list[can_fr]->S/CH_vol) : (0.0) ) :
03274                     ( (structs->TSser == 0) ? (structs->TS) : (Hist_data[structs->histID-1].arrayS[((int)(SimTime.TIME))] ) ) ;
03275 
03276                    /* vol avail for fluxing */
03277                 CH_vol =  Chan_list[can_fr]->area*(ramp(Chan_list[can_fr]->wat_depth-MINDEPTH) );
03278                 chan_over = (ChanHistOut > 0.0 ) ? ( CH_vol/ChanHistOut ) : 1.0;
03279 
03280                 do /* loop again through the remaining sequence of structures to correct any outflows for the canal */
03281                 {   
03282                     if ( structs->flag == 1 && struct_hist_start->canal_fr  == structs->canal_fr   ) { 
03283 
03284                         structs->multiOut = 1; /* flag to indicate this structure has been processed (below) */
03285 
03286                             /* revised FLOW */
03287                         if ( chan_over < 1.0 ) { /* if hist/other-data demand > canal vol */
03288                             if( debug > 0 ) { 
03289                                 fprintf( canDebugFile, "Day \t%6.1f\t %s:\tTotDemand= \t%7.0f\t; exceeded chan \t%d\t avail vol= \t%7.0f\t by \t%7.0f\t m3 (\t%5.3f\t m deep).\n",
03290                                          SimTime.TIME, structs->S_nam, ChanHistOut, Chan_list[can_fr]->number,
03291                                          CH_vol, (ChanHistOut - CH_vol), Chan_list[can_fr]->wat_depth );
03292                                 
03293                             }  
03294                                 /* decrement the individual structure outflows by it's proportional contrib. */
03295                             structs->flow = (structs->flow) * chan_over ;
03296                         } 
03297                             /* sum of all realized (corrected for excess demand) historical/other-data outflows from a canal during an iteration */
03298                             /* used to determine any allowable rule-based canal flows after the historical/other-data demands */
03299                         Chan_list[can_fr]->sumHistOut += structs->flow;         
03300 
03301                             /* the recipient changes as we cycle thru the structures (may be another cell or a canal) */
03302                         can_to = structs->canal_to - CHo;
03303                         cell_to = T(structs->cell_i_to,structs->cell_j_to); 
03304 
03305                         structs->conc_P = Pconc;/* concentrations applied to all structures draining this canal */
03306                         structs->conc_N = Nconc;
03307                         structs->conc_S = (structs->TSser == -1) ?
03308                             (Sconc) :
03309                             ( (structs->TSser == 0) ? (structs->TS) :
03310                             ( Hist_data[structs->histID-1].arrayS[((int)(SimTime.TIME))] ) ); /* see note below on tracer conc. application */
03311                          
03312                         P_fl = structs->flow * structs->conc_P; 
03313                         N_fl = structs->flow * structs->conc_N;
03314                         S_fl = structs->flow * structs->conc_S;
03315                                 /* change in nut mass in fr canal (kg) */
03316                         Chan_list[can_fr]->N -=  N_fl ;
03317                         Chan_list[can_fr]->P -=  P_fl ;
03318                         Chan_list[can_fr]->S -=  (structs->TSser == -1) ? (S_fl) : (0);
03319                             /* for the salt (S), we sometimes want to apply a tracer to internal canals/cells.
03320                                Thus, in that case, we don't subtract salt from the canal because this is "introduced" salt
03321                                that was not calculated from the available water and salt (i.e., the conc. was from the
03322                                CanalData.struct data file, either fixed or a time series) */
03323   
03324                         if (can_to > 0) {    /* IF this is a canal recip, change in nut mass in to canal (kg) */
03325                             /* sum of all realized (corrected for excess demand) historical/other-data outflows into a canal during an iteration */
03326                             /* used to determine any allowable rule-based canal flows after the historical/other-data demands */
03327                                 Chan_list[can_to]->sumHistIn += structs->flow;          
03328 
03329                             Chan_list[can_to]->N +=  N_fl ;
03330                             Chan_list[can_to]->P +=  P_fl ;
03331                             Chan_list[can_to]->S +=  S_fl ;
03332                                 /* mass balance in fr and to canals */
03333                             if (Chan_list[can_fr]->basin != Chan_list[can_to]->basin && debug > -999 ) {
03334                                 VOL_OUT_STR[Chan_list[can_fr]->basin] += structs->flow;
03335                                 VOL_IN_STR[Chan_list [can_to]->basin] += structs->flow; 
03336                                 P_OUT_STR[Chan_list[can_fr]->basin] += P_fl;
03337                                 P_IN_STR[Chan_list [can_to]->basin] += P_fl; 
03338                                 S_OUT_STR[Chan_list[can_fr]->basin] += (structs->TSser == -1) ? (S_fl) : (0); /* see above note on tracer application */
03339                                 S_IN_STR[Chan_list [can_to]->basin] += S_fl; 
03340                             }
03341                         }
03342                         else if (cell_to > 0) {
03343                                 /* OR change in nut mass in downstream cell (kg) */
03344  
03345                             if (structs->cell_i_to > 2 ) { /* to an on-map cell */
03346                                 SWH[cell_to] += structs->flow/CELL_SIZE;
03347                                 PA[cell_to] += P_fl; 
03348                                 NA[cell_to] += N_fl; 
03349                                 SA[cell_to] += S_fl; 
03350                                     /* mass balance in fr canal and to cell */
03351                                 if (Chan_list[can_fr]->basin != basn[cell_to] && debug > -999 ) {
03352                                     VOL_OUT_STR[Chan_list[can_fr]->basin] += structs->flow;
03353                                     VOL_IN_STR[basn[cell_to]] += structs->flow; 
03354                                     P_OUT_STR[Chan_list[can_fr]->basin] += P_fl;
03355                                     P_IN_STR[basn[cell_to]] += P_fl; 
03356                                     S_OUT_STR[Chan_list[can_fr]->basin] += (structs->TSser == -1) ? (S_fl) : (0); /* see above note on tracer application */
03357                                     S_IN_STR[basn[cell_to]] += S_fl; 
03358 
03359                                         /* all structure flows, EXCEPT canal-canal, must be between parent (family) hydro basins */
03360                                     if ( (Chan_list[can_fr]->family ==  basn_list[basn[cell_to]]->family) && (debug > 0) ) {
03361                                         sprintf( msgStr, "Day %6.1f: ERROR in basin configuration for structure from chan %d to cell (%d,%d): flow must be between diff hydrologic basins",
03362                                                  SimTime.TIME, Chan_list[can_fr]->number,structs->cell_i_to,structs->cell_j_to );
03363                                         WriteMsg( msgStr,True ); dynERRORnum++;
03364                                     }
03365                                 }
03366                             }
03367                             else if (structs->cell_i_to == 2 && debug > -999 ) { /* to off-map cell, mass balance check */
03368                                 VOL_OUT_STR[Chan_list[can_fr]->basin] += structs->flow; 
03369                                 VOL_OUT_STR[0] += structs->flow; 
03370                                 P_OUT_STR[Chan_list[can_fr]->basin] += P_fl; 
03371                                 P_OUT_STR[0] += P_fl; 
03372                                 S_OUT_STR[Chan_list[can_fr]->basin] += (structs->TSser == -1) ? (S_fl) : (0); /* see above note on tracer application */
03373                                 S_OUT_STR[0] += (structs->TSser == -1) ? (S_fl) : (0); /* see above note on tracer application */
03374                             }
03375                         }
03376                             /* don't do much with zero or negative flows */           
03377                         if (structs->flow == 0.0) {
03378                             structs->conc_P = 0.0;
03379                             structs->conc_N = 0.0;
03380                             structs->conc_S = 0.0;
03381                         } /* for output purposes only, zero the reported conc */
03382                         else if (structs->flow < 0.0) {
03383                             sprintf( msgStr, "Day %6.1f: ERROR IN ASSUMPTION OF NO NEG HIST FLOW from another chan/cell into chan %d",
03384                                      SimTime.TIME, Chan_list[can_fr]->number );
03385                             WriteMsg( msgStr,True ); dynERRORnum++;
03386                         }
03387                     }
03388 
03389                     structs = structs->next_in_list;
03390                 } while (structs != NULL); /*  cycling to the end of the list again !!!! */
03391 
03392                 structs = struct_hist_start; /* now pointing back to the first historical outflow structure that we dealt with */
03393     
03394 
03395                 goto OUT;
03396                     /* finished check for multiple (positive) outflows exceeding canal volume */
03397             }
03398                 /* FLOW */
03399                 /**** Rule- based flows */
03400             else { 
03401                 HeadH = Elev [cell_struct] - Chan_list[can_fr]->depth + Chan_list[can_fr]->wat_depth
03402                                 + (Chan_list[can_fr]->sumHistIn - Chan_list[can_fr]->sumHistOut)/Chan_list[can_fr]->area; 
03403           /* in v2.3, the rule-based canal->cell flows use either tide data or SFWMM/other model data for boundary conditions */
03404                 
03405                /* v2.2 and previous eqn: HeadT = ( ( structs->cell_i_to == 2 ) ? (Elev[cell_struct] - 0.1) : (SWH[cell_to]+Elev[cell_to]) ); */
03406                 /* Head at Tailwater location: check that destination is out of model domain (always in v2.3 and prior versions)
03407                         The TWf (TailWater flag) indicates use of annually-recurring tidal data, with Sea Level Rise added.
03408                         If the TWf is false, use relative depth data from the SFWMM or other model output (boundcond_depth data is NOT stage, but depth relative to elev NGVD'29) */
03409                         /* v2.5 use cell_i_TW, cell_j_TW in cell_to location */
03410                 HeadT = ( ( structs->cell_i_to == 2 ) ? 
03411                         ( (TWf) ? (structs->TW_stage) : ( boundcond_depth[T(structs->cell_i_TW,structs->cell_j_TW)]  + Elev[cell_struct]) ): 
03412                         (SWH[cell_to]+Elev[cell_to]) );
03413                                 /* volume available for fluxing - right now, the only structures here drain border seepage canals, leave 30 cm water */
03414                 CH_vol =  Chan_list[can_fr]->area* ramp(HeadH-Elev[cell_struct]+Chan_list[can_fr]->depth-0.3); /* avail flow volume */
03415             
03416                 if (HeadH>HeadT) 
03417                 {
03418                     /* this ****DOES NOT**** work when the HW and TW checks are being made at a remote cell */
03419                     /* ok for now, just fix when driving ELM with rules and checking remote cells */
03420                     structs->flow =  Min(structs->w_coef * ( HeadH - HeadT ) * canstep, CH_vol );                   
03421                         /* all rule based instances of canal->cell flow  go
03422                            to external cells, so this reversal check not invoked */
03423                     flow_rev = ( structs->cell_i_to != 2 ) ? 
03424                         (HeadH - structs->flow/Chan_list[can_fr]->area) - (HeadT + structs->flow/CELL_SIZE) :
03425                         (0.0);
03426                         /* the flux eqn is the same as the (new Jul98) virtual weir eqn, equilibrating the heads across two unequal area regions */
03427                     /* flow is only in positive directioin */
03428                     if (flow_rev < 0.0) {  structs->flow = 
03429                         ramp( Chan_list[can_fr]->area*CELL_SIZE*HeadH/(Chan_list[can_fr]->area+CELL_SIZE) -
03430                         Chan_list[can_fr]->area*CELL_SIZE*HeadT/(Chan_list[can_fr]->area+CELL_SIZE) );
03431                           }
03432                }
03433                 else structs->flow = 0.0; /* undefined flow if heads are equal, thus this stmt added */
03434                 
03435             
03436             
03437                     /* don't do much with zero or negative flows */           
03438                 if (structs->flow == 0.0) {
03439                     structs->conc_P = 0.0;
03440                     structs->conc_N = 0.0;
03441                     structs->conc_S = 0.0;
03442                     goto OUT;
03443                 } /* for output purposes only, zero the reported conc */
03444                 else if (structs->flow < 0) {  
03445                     sprintf( msgStr, "Day %6.1f: ERROR IN ASSUMPTION OF NO NEG FLOW from cell (%d,%d)",
03446                              SimTime.TIME, structs->cell_i_to, structs->cell_j_to );
03447                     WriteMsg( msgStr,True ); dynERRORnum++;
03448                     goto OUT;
03449                 }
03450           
03451                     /* rule-based flows:  calculate the nutrient concentrations */
03452                     /* conc in canal -> cell flow  */
03453                 CH_vol =  Chan_list[can_fr]->area* ramp(HeadH-Elev[cell_struct]+Chan_list[can_fr]->depth); /* TOTAL volume */
03454                     /* conc in donor (fr) canal (kg/m3) */
03455                 structs->conc_P = (CH_vol>0) ? (Chan_list[can_fr]->P/CH_vol) : (0.0);
03456                 structs->conc_N = (CH_vol>0) ? (Chan_list[can_fr]->N/CH_vol) : (0.0);
03457                 structs->conc_S = (CH_vol>0) ? (Chan_list[can_fr]->S/CH_vol) : (0.0);
03458 
03459                 P_fl = structs->flow * structs->conc_P; 
03460                 N_fl = structs->flow * structs->conc_N;
03461                 S_fl = structs->flow * structs->conc_S;
03462                   
03463                     /* increment the sum of rule-based flows to from canals
03464                        for potential multiple virtual flows associated with some canals */
03465                 Chan_list[can_fr]->sumRuleOut += structs->flow;         
03466                     /* change in nut mass in downstream cell (kg) */
03467                 if (structs->cell_i_to > 2 ) { /* to an on-map cell */
03468                     SWH[cell_to] += structs->flow/CELL_SIZE;
03469                     PA[cell_to] += P_fl; 
03470                     NA[cell_to] += N_fl; 
03471                     SA[cell_to] += S_fl; 
03472                         /* mass balance in fr canal and to cell */
03473                     if (Chan_list[can_fr]->basin != basn[cell_to] && debug > -999 ) {
03474                         VOL_OUT_STR[Chan_list[can_fr]->basin] += structs->flow;
03475                         VOL_IN_STR[basn[cell_to]] += structs->flow; 
03476                         P_OUT_STR[Chan_list[can_fr]->basin] += P_fl;
03477                         P_IN_STR[basn[cell_to]] += P_fl; 
03478                         S_OUT_STR[Chan_list[can_fr]->basin] += S_fl;
03479                         S_IN_STR[basn[cell_to]] += S_fl; 
03480 
03481                             /* all structure flows, EXCEPT canal-canal, must be between parent (family) hydro basins */
03482                         if ( (Chan_list[can_fr]->family ==  basn_list[basn[cell_to]]->family) && (debug > 0) ) {
03483                             sprintf( msgStr, "Day %6.1f: ERROR in basin configuration for structure from chan %d to cell (%d,%d): flow must be between diff hydrologic basins",
03484                                      SimTime.TIME, Chan_list[can_fr]->number,structs->cell_i_to,structs->cell_j_to );
03485                             WriteMsg( msgStr,True ); dynERRORnum++;
03486                         }
03487                     }
03488                 }
03489                 else if (structs->cell_i_to == 2 && debug > -999 ) { /* to off-map cell, mass balance check */
03490                     VOL_OUT_STR[Chan_list[can_fr]->basin] += structs->flow; 
03491                     VOL_OUT_STR[0] += structs->flow; 
03492                     P_OUT_STR[Chan_list[can_fr]->basin] += P_fl; 
03493                     P_OUT_STR[0] += P_fl; 
03494                     S_OUT_STR[Chan_list[can_fr]->basin] += S_fl; 
03495                     S_OUT_STR[0] += S_fl; 
03496                 }
03497 
03498             
03499                     /* calc the change in mass in  canal (kg) */
03500                 Chan_list[can_fr]->N -=  N_fl ;
03501                 Chan_list[can_fr]->P -=  P_fl ;
03502                 Chan_list[can_fr]->S -=  S_fl ;
03503 
03504                 goto OUT;
03505             }
03506         }
03507         
03508     
03509 
03510 /* CASE 4 */
03511 /*check that this is a structure between a cell and a canal*/
03512             /* if from 1,1 to canal, always will require historical data */
03513         if ( structs->canal_to != 0 && structs->cell_i_fr != 0 && structs->cell_j_fr != 0 ) { 
03514             can_to = structs->canal_to - CHo;
03515             cell_fr = T(structs->cell_i_fr,structs->cell_j_fr);           
03516             cell_struct = T(structs->str_cell_i,structs->str_cell_j);
03517 
03518                 /* first check for historical (or other data) flow */
03519             if ( structs->flag == 1 ) {
03520   
03521                 if (structs->multiOut == 1) {
03522                     sprintf( msgStr, "Day %6.1f: ERROR IN ASSUMPTION OF NO NEG HIST FLOW from chan %d",
03523                              SimTime.TIME, Chan_list[can_to]->number );
03524                     WriteMsg( msgStr,True ); dynERRORnum++;
03525                     goto OUT; /* this structure has already been processed as a multiple outflow from a canal */
03526                 }
03527                 
03528                 structs->flow = Hist_data[structs->histID-1].arrayWat[((int)(SimTime.TIME))]*canstep;
03529                     /* don't do much with zero or negative flows */           
03530                 if (structs->flow == 0.0) {
03531                     structs->conc_P = 0.0;
03532                     structs->conc_N = 0.0;
03533                     structs->conc_S = 0.0;
03534                     goto OUT; 
03535                 }
03536                 else if (structs->flow <0) { 
03537                     sprintf( msgStr, "Day %6.1f: ERROR IN ASSUMPTION OF NO NEG HIST FLOW from chan %d",
03538                              SimTime.TIME, Chan_list[can_to]->number );
03539                     WriteMsg( msgStr,True ); dynERRORnum++;
03540                     goto OUT;
03541                 } 
03542 
03543                 if (structs->cell_i_fr != 2 ) /* cell (2,2) == (1,1) in maps, not in model domain */
03544                 {
03545 
03546                     cell_vol = SWH[cell_fr] * CELL_SIZE;
03547                     if (structs->flow > cell_vol) {
03548                         if( debug > 0 ) { 
03549                             fprintf( canDebugFile, "Day \t%6.1f\t %s:\tTotDemand= \t%7.0f\t; exceeded cell \t(%d,%d)\t vol= \t%7.0f\t by \t%7.0f\t m3.\n",
03550                                      SimTime.TIME, structs->S_nam, structs->flow, structs->cell_i_fr, structs->cell_j_fr, cell_vol, (structs->flow - cell_vol) );
03551                             
03552                         }
03553                         structs->flow = cell_vol;
03554                     }}
03555                     /* calculate the nutrient conc in the donor cell */
03556                     /* nut conc in flow from outside system using historical data or 0.0 */
03557                 if (structs->cell_i_fr == 2) { 
03558                     structs->conc_P = (structs->TPser == 1)  ?
03559                         Hist_data[structs->histID-1].arrayP[((int)(SimTime.TIME))] :
03560                         ((structs->TPser == 0) ? (structs->TP) : (0.0) );
03561                     structs->conc_N = (structs->TNser == 1)  ?
03562                         Hist_data[structs->histID-1].arrayN[((int)(SimTime.TIME))] :
03563                         ((structs->TNser == 0) ? (structs->TN) : (0.0) );
03564                     structs->conc_S = (structs->TSser == 1)  ?
03565                         Hist_data[structs->histID-1].arrayS[((int)(SimTime.TIME))] :
03566                         ((structs->TSser == 0) ? (structs->TS) : (0.0) );
03567                 }
03568                 else
03569                 { /* nut conc in flow within system - check to see if there are historical concs with the historical flow */
03570                     structs->conc_P = (structs->TPser == -1)  ?
03571                         ( (SWH[cell_fr] > 0) ? PA[cell_fr]/CELL_SIZE/SWH[cell_fr] : 0.0 ) :
03572                         ( (structs->TPser == 0) ? (structs->TP) : (Hist_data[structs->histID-1].arrayP[((int)(SimTime.TIME))] ) );
03573                     structs->conc_N = ( structs->TNser == -1)  ?
03574                         ( (SWH[cell_fr] > 0) ? NA[cell_fr]/CELL_SIZE/SWH[cell_fr] : 0.0 ) :
03575                         ( (structs->TNser == 0) ? (structs->TN) : (Hist_data[structs->histID-1].arrayN[((int)(SimTime.TIME))] ) );
03576                     structs->conc_S = ( structs->TSser == -1)  ?
03577                         ( (SWH[cell_fr] > 0) ? SA[cell_fr]/CELL_SIZE/SWH[cell_fr] : 0.0 ) :
03578                         ( (structs->TSser == 0) ? (structs->TS) : (Hist_data[structs->histID-1].arrayS[((int)(SimTime.TIME))] ) );
03579                 }
03580                     /* now determine the nutrient mass (kg) flow from cell to canal */
03581                     /* 1000 is to convert g/m3 to kg */ 
03582                 P_fl = structs->flow * structs->conc_P;
03583                 N_fl = structs->flow * structs->conc_N;
03584                 S_fl = structs->flow * structs->conc_S;
03585                 if (structs->cell_i_fr == 2 ) {/* "from" cell is outside system */
03586                     if (debug > -999) {  
03587                         VOL_IN_STR[Chan_list[can_to]->basin] += structs->flow; /* mass balance */
03588                         VOL_IN_STR[0] += structs->flow; 
03589                         P_IN_STR[Chan_list[can_to]->basin] += P_fl; 
03590                         P_IN_STR[0] += P_fl; 
03591                         S_IN_STR[Chan_list[can_to]->basin] += S_fl; 
03592                         S_IN_STR[0] += S_fl; 
03593                     } }
03594                 else { /* internal cell-> canal flow */
03595                     SWH[cell_fr] -= structs->flow/CELL_SIZE;
03596                     PA[cell_fr] -= P_fl; 
03597                     NA[cell_fr] -= N_fl; 
03598                     SA[cell_fr] -= (structs->TSser == -1) ? (S_fl) : (0);
03599                             /* for the salt (S), we sometimes want to apply a tracer to internal canals/cells.
03600                                Thus, in that case, we don't subtract salt from the canal because this is "introduced" salt
03601                                that was not calculated from the available water and salt (i.e., the conc. was from the
03602                                CanalData.struct data file, either fixed or a time series) */
03603                         /* mass balance in fr cell and to canal */
03604                     if (basn[cell_fr] != Chan_list[can_to]->basin && debug > -999 ) {
03605                         VOL_OUT_STR[basn[cell_fr]] += structs->flow;
03606                         VOL_IN_STR[Chan_list[can_to]->basin] += structs->flow; 
03607                         P_OUT_STR[basn[cell_fr]] += P_fl;
03608                         P_IN_STR[Chan_list[can_to]->basin] += P_fl; 
03609                         S_OUT_STR[basn[cell_fr]] += (structs->TSser == -1) ? (S_fl) : (0); /* see above note on tracer application */
03610                         S_IN_STR[Chan_list[can_to]->basin] += S_fl; 
03611 
03612                             /* all structure flows, EXCEPT canal-canal, must be between parent (family) hydro basins */
03613                         if ( (Chan_list[can_to]->family ==  basn_list[basn[cell_fr]]->family) && (debug > 0) ) {
03614                             sprintf( msgStr, "Day %6.1f: ERROR in basin configuration for structure from cell (%d,%d) to chan %d: flow must be between diff hydrologic basins",
03615                                      SimTime.TIME, structs->cell_i_fr, structs->cell_j_fr, Chan_list[can_to]->number );
03616                             WriteMsg( msgStr,True ); dynERRORnum++;
03617                         }
03618                     }
03619                 }
03620                 
03621                   /* sum of all realized (corrected for excess demand) historical/other-data outflows into a canal during an iteration */
03622                   /* used to determine any allowable rule-based canal flows after the historical/other-data demands */
03623                 Chan_list[can_to]->sumHistIn += structs->flow;          
03624 
03625                    /* change in nut mass in "to" canal (kg) */
03626                 Chan_list[can_to]->N +=  N_fl ;
03627                 Chan_list[can_to]->P +=  P_fl ;
03628                 Chan_list[can_to]->S +=  S_fl ;
03629 
03630                 goto OUT;
03631             }
03632 
03633             else { /* rule based flow calc */
03634                 if (structs->cell_i_fr != 2  ) {
03635                     sprintf( msgStr, "Day %6.1f: ERROR in flow configuration for structure from cell (%d,%d) to chan %d: this rule-based flow is restricted to flow from outside model domain.",
03636                              SimTime.TIME, structs->cell_i_fr, structs->cell_j_fr, Chan_list[can_to]->number );
03637                     WriteMsg( msgStr,True );   dynERRORnum++;    
03638                 }
03639                 HeadH = (structs->HW_stage) ; 
03640           /* in v2.3, the rule-based cell->canal flows use only tide data for boundary conditions */
03641                 HeadT = Elev [cell_struct] - Chan_list[can_to]->depth + Chan_list[can_to]->wat_depth
03642                     + (Chan_list[can_to]->sumHistIn - Chan_list[can_to]->sumHistOut)/Chan_list[can_to]->area;
03643             
03644                 if (HeadH>HeadT) 
03645                 {
03646                     structs->flow =  Max(structs->w_coef * ( HeadH - HeadT ) * canstep, 0 );                   
03647                 }
03648                 else structs->flow = 0.0; /* undefined flow if heads are equal, thus this stmt added */
03649                 
03650                 /* don't do much with zero or negative flows */           
03651             if (structs->flow == 0.0) {
03652                 structs->conc_P = 0.0;
03653                 structs->conc_N = 0.0;
03654                 structs->conc_S = 0.0;
03655                 goto OUT;
03656             } /* for output purposes only, zero the reported conc */
03657 
03658                 /*  rule-based: calculate the nutrients - this is restricted to calculated tidal boundary condition flows  from outside of system */
03659                 /* concentration in structure flows, g/L = kg/m^3; mass of stock in kg */
03660 
03661                 if (structs->cell_i_fr == 2) { 
03662                     structs->conc_P = (structs->TPser == 1)  ?
03663                         Hist_data[structs->histID-1].arrayP[((int)(SimTime.TIME))] :
03664                         ((structs->TPser == 0) ? (structs->TP) : (0.0) );
03665                     structs->conc_N = (structs->TNser == 1)  ?
03666                         Hist_data[structs->histID-1].arrayN[((int)(SimTime.TIME))] :
03667                         ((structs->TNser == 0) ? (structs->TN) : (0.0) );
03668                     structs->conc_S = (structs->TSser == 1)  ?
03669                         Hist_data[structs->histID-1].arrayS[((int)(SimTime.TIME))] :
03670                         ((structs->TSser == 0) ? (structs->TS) : (0.0) );
03671                 }
03672             P_fl = structs->flow * structs->conc_P; 
03673             N_fl = structs->flow * structs->conc_N;
03674             S_fl = structs->flow * structs->conc_S;
03675              
03676                 /* this calculated flow is restricted to an outside-system cell, which does not need updating */
03677                 /* these are nutrient stocks, not concentrations */  
03678                 /* mass balance in fr cell and to canal */
03679             if (structs->cell_i_fr == 2  ) {
03680                 VOL_IN_STR[Chan_list[can_to]->basin] += structs->flow; 
03681                 VOL_IN_STR[0] += structs->flow; 
03682                 P_IN_STR[Chan_list[can_to]->basin] += P_fl; 
03683                 P_IN_STR[0] += P_fl; 
03684                 S_IN_STR[Chan_list[can_to]->basin] += S_fl; 
03685                 S_IN_STR[0] += S_fl; 
03686             }
03687             else        {  
03688                 sprintf( msgStr, "Day %6.1f: ERROR in flow to chan %d: cannot have calculated cell->canal flows except from outside of model domain. ",
03689                          SimTime.TIME, Chan_list[can_to]->number);
03690                 WriteMsg( msgStr,True ); dynERRORnum++;
03691                 goto OUT; 
03692             }
03693             
03694 
03695                 /* change in nut mass in canal (kg) */
03696             Chan_list[can_to]->N +=  N_fl ;
03697             Chan_list[can_to]->P +=  P_fl ;
03698             Chan_list[can_to]->S +=  S_fl ;
03699             
03700             goto OUT;
03701             }
03702         }
03703         
03704  
03705 /* CASE 5 */
03706 /* if none of the above structure flow conditions */
03707         printf ("\nError in data for structure N %s: impossible condition! \n", structs->S_nam );
03708         exit (-1);
03709         
03710       OUT:
03711         if (TWf) { structs->TW_stage = -99; TWf = 0; }
03712         if (HWf) { structs->HW_stage = -99; HWf = 0; }
03713 
03714             /* sum of flows and solute masses thru structs over the can_Intvl canal summary interval */
03715         structs->SumFlow += structs->flow; 
03716         structs->Sum_P += structs->conc_P * structs->flow; 
03717 /*        structs->Sum_N += structs->conc_N * structs->flow; */
03718         structs->Sum_S += structs->conc_S * structs->flow; 
03719 
03720         structs->multiOut = 0; /* reset the multi-outflow flag to zero after done */
03721         
03722             /* print sum of flows and the flow weighted mean conc thru structures */
03723         if ( printchan ) { 
03724             fprintf( WstructOutFile, "%7.0f\t", (structs->SumFlow)/(2446.848) ); /* convert to cfs (0.02832*60*60*24) */
03725             fprintf( WstructOutFile_P, "%7.4f\t", 
03726                      (structs->SumFlow > 0.0) ? (structs->Sum_P/structs->SumFlow*1000.0) : (0.0) ); /* convert g/L to mg/L */
03727 /*            fprintf( WstructOutFile_N, "%7.4f\t", 
03728               (structs->SumFlow > 0.0) ? (structs->Sum_N/structs->SumFlow*1000.0) : (0.0) );*/ /* convert g/L to mg/L */
03729             fprintf( WstructOutFile_S, "%7.4f\t", 
03730                      (structs->SumFlow > 0.0) ? (structs->Sum_S/structs->SumFlow) : (0.0) ); /* g/L, no conversion */
03731 
03732             structs->SumFlow = 0.0; 
03733             structs->Sum_P = 0.0;
03734 /*            structs->Sum_N = 0.0; */
03735             structs->Sum_S = 0.0;
03736         }
03737         
03738             /* increment to the next structure in the list */
03739         structs = structs->next_in_list;
03740     }
03741 
03742     if (printchan) {
03743         fflush( WstructOutFile );
03744         fflush( WstructOutFile_P );
03745 /*         fflush( WstructOutFile_N ); */
03746         fflush( WstructOutFile_S );
03747     }
03748  
03749     return;
03750     
03751 }

Here is the call graph for this function:

float f_Manning float  delta,
float  Water,
float  SW_coef
 

Surface water exchange between cell and canal.

Parameters:
delta Head difference between cell and canal (m)
Water Hydraulic radius (m)
SW_coef Aggregated flow coefficient (m^0.5 * sec)/m^(1/3)
Returns:
m^3

Definition at line 2609 of file WatMgmt.c.

References Abs, GP_mannDepthPow, and sgn.

Referenced by FluxChannel().

02610 {       float a_delta;
02611     
02612         a_delta = Abs (delta);
02613         return ( sgn( delta ) * SW_coef * pow(Water, GP_mannDepthPow) * sqrt(a_delta) * canstep);
02614 
02615 }

float f_Ground float  dh,
float  height,
float  GW_coef,
float  I_Length
 

Subsurface ground water exchange between cell and canal.

Parameters:
dh Head difference between cell and canal (m)
height Height of water (m)
GW_coef Aggregated flow coefficient (m/d)
I_Length Average length of cell along reach segment
Returns:
m^3

Definition at line 2626 of file WatMgmt.c.

Referenced by FluxChannel().

02627 {       
02628         float  w;
02629 
02630 /* porosity already part of GW_coef */
02631         w = dh * I_Length * GW_coef / (0.5*celWid) * height * canstep; 
02632 
02633         return ( w);
02634 
02635 }

float GetGraph struct Schedule graph,
float  x
 

Graph to read the time dependent schedules.

Parameters:
graph A struct of Schedule
x Time, julian day within recurring year
Returns:
value The target stage (m)

Definition at line 3762 of file WatMgmt.c.

References Schedule::graph_points, Schedule::num_point, Points::time, and Points::value.

Referenced by Flows_in_Structures().

03763 {
03764     float s;
03765     int ig = 0, Np;
03766         
03767     Np = graph->num_point;
03768 
03769     while(1) 
03770     {
03771         if (x <= graph->graph_points[ig].time) 
03772         { if ( ig > 0 ) ig--; 
03773         else return ( graph->graph_points[0].value );
03774         }
03775         else if (x > graph->graph_points[ig].time && x <= graph->graph_points[ig+1].time) 
03776         {
03777             s =         ( graph->graph_points[ig+1].value - graph->graph_points[ig].value )/
03778                 ( graph->graph_points[ig+1].time - graph->graph_points[ig].time );
03779             return ( s * ( x - graph->graph_points[ig].time ) + graph->graph_points[ig].value ); 
03780         }
03781         else if (x > graph->graph_points[ig+1].time) 
03782         { if ( ig < Np ) ig++; 
03783         else return ( graph->graph_points[Np].value );
03784         }
03785     }
03786 }

struct Chan* ReadChanStruct char *  filename  ) 
 

Reads the information about canals from data file.

Parameters:
filename Canal/levee attributes data file.
Returns:
chan_first Pointer to first canal data struct

Definition at line 406 of file WatMgmt.c.

References Chan::basin, basins, basn_list, basndef::basnTxt, C_F, ChanInFile, CHANNEL_MAX_ITER, CHo, Chan::cond, Chan::depth, Chan::edgeMann, F_ERROR, Chan::family, basndef::family, H_OPSYS, Chan::ic_depth, Chan::ic_N_con, Chan::ic_P_con, Chan::ic_S_con, Chan::levee, MINDEPTH, modelFileName, modelName, ModelPath, msgStr, Chan::next_in_list, Chan_reach::next_reach, num_chan, Chan::number, Chan::parent, basndef::parent, ProjName, Chan::reaches, Chan::roil, Chan::roir, SimAlt, usrErr(), UTM2kmx(), UTM2kmy(), UTM_EOrigin, UTM_NOrigin, Chan::wat_depth, Chan::width, Chan_reach::x0, Chan_reach::x1, Chan_reach::y0, and Chan_reach::y1.

Referenced by Canal_Network_Init().

00407 { 
00408  struct Chan *channels, *chan_first, *chan_last;
00409  struct Chan_reach *C_reach, *last_reach, *next_reach;
00410  int i, index, firstReach, firstCanal = 1;
00411  char ss[422], scenario[20], modnam[20], bas_txt[8];
00412  float C_x, C_y, C_x1, C_y1;
00413 
00414  
00415  int ibas, foundBasn;
00416  
00417 
00418  { if(H_OPSYS==UNIX) 
00419      sprintf( modelFileName, "%s/%s/Data/%s.chan", ModelPath, ProjName, filename );
00420  else sprintf( modelFileName, "%s%s:Data:%s.chan", ModelPath, ProjName, filename );
00421  }
00422 
00423 /* Open file with canals data */
00424  if ( ( ChanInFile = fopen( modelFileName, "r" ) ) ==  NULL )
00425  {
00426      sprintf( msgStr,"Can't open %s file! ",modelFileName ) ; usrErr(msgStr);
00427      exit(-1) ;
00428  }
00429   
00430 /* Allocate memory for canal reach structure. A canal reach is the straight canal segment */
00431 /* defined in the canal data file by the coordinates of its firstReaching and ending points.   */
00432 /* Each canal can be made of several canal reaches              */
00433  if ( (C_reach = (struct Chan_reach *) malloc( (size_t) sizeof( struct Chan_reach ))) == NULL )
00434  {
00435      usrErr( "No memory for first canal reach\n " ) ;
00436      exit( -2 ) ;
00437  };
00438 
00439  num_chan = 0;
00440    
00441 /* Read the general canal network information */ 
00442  fgets( ss, 420, ChanInFile ); sscanf( ss, "%d", &UTM_EOrigin );
00443  fgets( ss, 420, ChanInFile ); sscanf( ss, "%d", &UTM_NOrigin );
00444  fgets( ss, 420, ChanInFile );  /* skip comment line */
00445  fgets( ss, 420, ChanInFile ); sscanf( ss,"%s %s %d %d", &modnam, &scenario, &CHo, &CHANNEL_MAX_ITER );
00446  if (strcmp(scenario,SimAlt) != 0) {
00447      sprintf(msgStr, "The scenario (%s) found in the %s file doesn't match the one (%s) you asked for in Driver.parm!",
00448              scenario, modelFileName, SimAlt); usrErr(msgStr);
00449              exit(-1);
00450  }
00451  if (strcmp(modnam,modelName) != 0) {
00452      sprintf(msgStr, "The model name (%s) found in the %s file doesn't match the one (%s) you asked for in Driver.parm!",
00453              modnam, modelFileName, modelName); usrErr(msgStr);
00454              exit(-1);
00455  }
00456   
00457  fgets( ss, 420, ChanInFile );   /* skip comment line */
00458  fgets( ss, 420, ChanInFile ); sscanf ( ss, "%f %f %f ", &F_ERROR, &MINDEPTH, &C_F );
00459      /* the C_F is used only for sensitivity experiments */
00460  fgets( ss, 420, ChanInFile );   /* skip comment line */
00461    
00462 /* Read canal descriptors until none left */   
00463  while ( fgets( ss, 420, ChanInFile ) != NULL && !feof( ChanInFile ) )
00464  { 
00465      if ( *ss == '#' )  /* "#" identifies each new canal in the data file */
00466      { /* Allocate memory for canal structure */
00467          if ( (channels = (struct Chan *) malloc( (size_t) sizeof( struct Chan ))) == NULL )
00468          {
00469              printf( "No memory for first canal\n " ) ;
00470              exit( -2 ) ;
00471          };
00472                 
00473              /* Read canal information:  - canal number, 
00474                 - levee marker (1 - left levee, -1 - right, 0 - none, 2 - both sides),
00475                 - ranges of cell interaction to the left and to the right - functionality removed, should ALWAYS=1,
00476                 - canal depth and width,
00477                 - seepage coefficient through the levee,
00478                 - initial S and P concentrations (doubles) in the canal water,
00479                 - initial water depth in canal (can be >, = or < than canal depth 
00480                 - Manning's n associated w/ edge of canal, to accomodate topographic lip/berm and/or denser veg along canal length
00481                 - hydrologic basin that canal exchanges surface water with */
00482          i = sscanf( ss+1, "%d\t%d\t%d\t%d\t%f\t%f\t%f\t%f\t%f\t%f\t%f\t%s",
00483                      &channels->number, 
00484                      &channels->levee,
00485                      &channels->roil,   &channels->roir, 
00486                      &channels->depth,  &channels->width,  
00487                      &channels->cond,   
00488                      &channels->ic_S_con,  &channels->ic_P_con,  
00489                      &channels->wat_depth, &channels->edgeMann, &bas_txt);
00490 
00491          foundBasn=0;
00492          for (ibas=numBasn;ibas>0;ibas--) {
00493              basins = basn_list[ibas];
00494              if(strcmp(bas_txt,basins->basnTxt)==0) {
00495                  channels->basin = ibas;
00496                  foundBasn=1;
00497                  channels->family = basn_list[ibas]->family;
00498                  channels->parent = basn_list[ibas]->parent;
00499                      /* canals may ONLY belong to parent hydro basins, not to indicator regions */
00500                  if (!channels->parent) {
00501                      printf ("Sorry - %s's canal %d's basin is not a parent hydrologic basin! Please assign to a hydro basin in the CanalData.chan file.\n",
00502                             modelName,channels->number); exit (-1);}     
00503              }
00504          }
00505          if (foundBasn==0) { printf ("Sorry - %s's canal %d's basin does not exist! Please define the basin in the basinIR file.\n",
00506                             modelName,channels->number); exit (-1);}                
00507 
00508          channels->ic_P_con *= 0.001; /* initial input value in mg/L, convert to g/L==kg/m3 */
00509          channels->ic_N_con = 0.00; /* v2.2 "removed" N */ /* initial input value in mg/L, convert to g/L==kg/m3 */
00510          channels->ic_S_con *= 1.0; /* initial input value in g/L, no conversion */
00511              /* store the initial depth for re-initializing under Positional Analysis mode */
00512          channels->ic_depth = channels->wat_depth;
00513          
00514 
00515          channels->reaches = C_reach; 
00516                         
00517          if (i < 10)    /* check for the number of input fields read */
00518          {printf ( " Error in CanalData.chan: reach input %d\n", (num_chan+1) ); exit (-1);}
00519                         
00520                             
00521          num_chan++;            /*count the number of canals (not canal reaches) */
00522          firstReach = 1;
00523                         
00524          if (firstCanal) 
00525          {
00526              firstCanal = 0;  
00527              chan_first = channels;
00528              chan_last = channels;
00529              continue; 
00530          }
00531 
00532          last_reach->next_reach = NULL;
00533          chan_last->next_in_list = channels;
00534          chan_last = channels;     
00535      }
00536      else  
00537      { /* read the pairs of coordinates for the consecutive canal reaches */
00538          i = sscanf ( ss, "%f %f", &C_x1, &C_y1 );
00539          if (i < 2)     /* check for the number of input fields read */
00540          {printf ( " Error in CanalData.chan coords: %d'th reach input \n", (num_chan+1) ); exit (-1);}
00541                 
00542              /* Convert from UTM to km units with 0,0 origin */
00543          C_x = UTM2kmx (C_x1);  C_y = UTM2kmy (C_y1);           
00544          C_reach->x1 = C_x;  C_reach->y1 = C_y;
00545                         
00546          if ( firstReach )                   
00547          {  C_reach->x0 = C_x; C_reach->y0 = C_y; 
00548          firstReach = 0;
00549          }
00550          else 
00551          {  
00552                                 /* Allocate memory for canal reach structure */
00553              if ( (next_reach = (struct Chan_reach *) malloc( (size_t) sizeof( struct Chan_reach ))) == NULL )
00554              {  printf( "No memory for first canal reach\n " ) ;
00555              exit( -2 ) ;
00556              };
00557              C_reach->next_reach = next_reach;
00558              last_reach = C_reach;
00559              C_reach = next_reach;
00560              C_reach->x0 = C_x; C_reach->y0 = C_y;
00561          }                       
00562      }  
00563  }/* end loop for reading */
00564   
00565  chan_last->next_in_list = NULL;
00566  last_reach->next_reach = NULL;
00567  free ((char *) C_reach);
00568  fclose (ChanInFile);
00569    
00570  usrErr( "\tCanal locs/attributes read OK" );
00571    
00572  return (chan_first);
00573 }

Here is the call graph for this function:

struct Cells* MarkCell int  x,
int  y,
int  c_ind,
float  length,
struct Cells cell,
int  direction,
int  levee,
int  xm,
int  ym,
int  c_num,
int *  marked,
float  distTot
 

Mark cell function.

This function does two things:
1. It stores the cell coordinates (x, y) in the interaction array, also identifying the location of this cell with respect to the canal;
2. It marks the cell (xm,ym) on the ON_MAP array to prevent horizontal overland fluxes if there is a levee associated with this canal.

Parameters:
x Grid cell row
y Grid cell column
c_ind Cell index = 1 for left cells and -1 for right cells
length Length of canal associated with this cell
cell A struct of Cells
direction Direction of interaction: 1 - none, 2 -to the East, 3 - to the South, 4 - all directions
levee Levee location attribute
xm Grid cell row marked
ym Grid cell column marked
c_num Canal ID number
marked Denote a marked cell
distTot Cumulative (along grid cells) distance along a canal reach, from the starting point (m)
Returns:
cell Pointer to struct of Cells

Definition at line 1728 of file WatMgmt.c.

References CanalCellInteractDebug, cell, cell_last, Cells::ind, Cells::length, Cells::next_cell, ON_MAP, Cells::reachDistDown, s0, T, Cells::x, and Cells::y.

Referenced by Channel_configure().

01730 { struct Cells *cell_next;
01731   int N, cellLoc;
01732   
01733    if ( x >= s0 || y >= s1 )
01734      {  printf ( "Error in definition of canal in cell x=%d  y=%d ", x, y );
01735         exit (-2);                      /*check that cells obtained are within array boundaries */
01736      }
01737 
01738    N = T(xm+1, ym+1);
01739    if ( !ON_MAP[N] ) return cell;
01740    
01741    if ( levee && direction != ALL )
01742    ON_MAP[N] = ((ON_MAP[N]-1) & (direction+100)) + 1;   /* modifying ON_MAP to take into account the 
01743    allowed direction of fluxing: 00 - none, 01=1 - to the East, 10=2 - to the South, 11=3 - all 
01744    ON_MAP they naturally become 1 - none, 2 -to the East, 3 - to the South, 4 - all directions 
01745    0 is reserved for cells out of the ON_MAP area*/
01746    
01747    if ( !C_Mark ) return cell; /* C_Mark is to make gaps between canals, so that there will be no 
01748                 canals being connected through cells if it turns out that they are linked to the same cells */
01749              
01750        /* a temporary fix used to determine if a cell belonging to a canal has already been marked */
01751        /***********/
01752    cellLoc = T(x+1,y+1);
01753    if (marked[cellLoc] == c_num ) return cell;
01754    marked[cellLoc] = c_num;
01755        /***********/
01756              
01757 
01758    
01759    cell->x = x+1;                                               /* shift coord to match the other maps */
01760    cell->y = y+1;
01761    cell->ind = c_ind;   
01762    cell->length = length*celWid;                /* length of canal associated with this cell * cell width m */
01763    cell->reachDistDown = distTot; /*  cumulative (by grid cells) distance along a canal reach, from the starting point (m) (v2.5) */
01764    
01765    fprintf ( CanalCellInteractDebug, " %d, %d, l/r = %d, distance= %7.1f\n", cell->x, cell->y, c_ind, cell->reachDistDown );
01766    
01767    /* allocate memory for next cell structure */
01768    if ( (cell_next = (struct Cells *) malloc((size_t) sizeof( struct Cells ))) == NULL )
01769        {
01770          printf( "Failed to allocate memory for cell structure\n " ) ;
01771          exit( -2 ) ;
01772        }
01773    cell->next_cell = cell_next;
01774    cell_last = cell;
01775 
01776    return cell_next;
01777 }

struct Schedule* Read_schedule char *  sched_name,
char *  filename,
float  BASE_DATUM
 

Read target stage schedule as a recurring time-series graph.

Parameters:
sched_name The name of the target stage location
filename Input file name for graph
BASE_DATUM GP_DATUM_DISTANCE
Returns:
Sce_Graph struct of the time series schedule graph

Definition at line 1222 of file WatMgmt.c.

References Schedule::graph_points, H_OPSYS, modelFileName, modelName, ModelPath, msgStr, Schedule::num_point, ProjName, schedFile, Scip(), SimAlt, Points::time, usrErr(), Points::value, and Wrdcmp().

Referenced by ReadStructures().

01223 { 
01224 char ss[301], *s, scenario[20], modnam[20];
01225 int i, count = 1, sched_found=0;
01226 struct Schedule *Sce_Graph;
01227 
01228    { if(H_OPSYS==UNIX) 
01229                 sprintf( modelFileName, "%s/%s/Data/%s.graph", ModelPath, ProjName, filename );
01230      else sprintf( modelFileName, "%s%s:Data:%s.graph", ModelPath, ProjName, filename );
01231    }
01232 
01233 /* Open file with schedule graphs */
01234   if ( ( schedFile = fopen( modelFileName, "r" ) ) ==  NULL )
01235         {
01236        printf( "Can't open %s file!\n", modelFileName ) ;
01237        exit(-1) ;
01238     }
01239 
01240  fgets( ss, 300, schedFile );  /* skip comment line */
01241  fgets( ss, 300, schedFile );  /* skip comment line */
01242  fgets( ss, 300, schedFile );  /* skip comment line */
01243  fgets( ss, 300, schedFile ); sscanf( ss,"%s %s", &modnam, &scenario );
01244  if (strcmp(scenario,SimAlt) != 0) {
01245      sprintf(msgStr, "The scenario (%s) found in the %s file doesn't match the one (%s) you asked for in Driver.parm!",
01246              scenario, modelFileName, SimAlt); usrErr(msgStr);
01247              exit(-1);
01248  }
01249  if (strcmp(modnam,modelName) != 0) {
01250      sprintf(msgStr, "The model name (%s) found in the %s file doesn't match the one (%s) you asked for in Driver.parm!",
01251              modnam, modelFileName, modelName); usrErr(msgStr);
01252              exit(-1);
01253  }
01254 
01255 /* Read schedule descriptors until none left */   
01256   while ( fgets( ss, 300, schedFile ) != NULL && !feof( schedFile ) )
01257         { if ( !Wrdcmp ( sched_name, ss ) ) continue;
01258           sched_found = True;
01259           
01260           fgets( ss, 300, schedFile );  
01261           s = ss;
01262           /* count the number of pairs for the graph */
01263           while ( *s++ != EOL ) if ( *s == '(' ) count++;    
01264 
01265           /* Allocate memory for next schedule */
01266       if ( (Sce_Graph = (struct Schedule *) malloc( (size_t) sizeof( struct Schedule ))) == NULL )
01267        {
01268          printf( "Failed to allocate memory for Managed flow Schedule graph in (%s)\n ", sched_name ) ;
01269          exit( -2 ) ;
01270        };
01271      
01272           /* Allocate memory for next schedule graph containing 'count' points*/
01273       if ( (Sce_Graph->graph_points = 
01274                 (struct Points *) malloc( (size_t) sizeof( struct Points ) * count)) == NULL )
01275        {
01276          printf( "Failed to allocate memory for Managed flow Schedule graph in (%s)\n ", sched_name ) ;
01277          exit( -2 ) ;
01278        };
01279           
01280           Sce_Graph->num_point = count;
01281           
01282           /* Read the graph points */
01283           s = ss;
01284           for ( i = 0; i < count && *s!=EOS; i++ )
01285           { 
01286                 s = Scip( s,'(');       
01287                 sscanf( s, "%f", &Sce_Graph->graph_points[i].time ); 
01288                 s = Scip ( s, ',');
01289                 sscanf( s, "%f", &Sce_Graph->graph_points[i].value );
01290                 Sce_Graph->graph_points[i].value += BASE_DATUM;
01291           }
01292         } /* end of while */
01293         /* v2.5 included this error-trap */
01294         if (!sched_found) {
01295              sprintf(msgStr, "The schedule name (%s) was not found in the %s schedule file! (Was named in the CanalData.struct file.)",
01296              sched_name, modelFileName); usrErr(msgStr);
01297              exit(-1);
01298     }
01299 
01300         
01301   fclose (schedFile);
01302   return Sce_Graph;
01303 }

Here is the call graph for this function:

int Wrdcmp char *  s,
char *  t
 

Compare two words.

Parameters:
s A string
t Another string
Returns:
True/false

Definition at line 3815 of file WatMgmt.c.

Referenced by Read_schedule().

03816 { 
03817   for ( ; *s == *t; s++, t++ )
03818          if ( isspace (*(s+1)) || *(s+1) == EOS ) return 1;
03819         
03820   return 0;
03821   
03822 }

float UTM2kmx double  UTM  ) 
 

Converter from meters in UTM to grid cell row location (zero origin).

Parameters:
UTM Geographic coordinate (m, ELM uses UTM zone17, NAD'27, but not important here)
Returns:
Grid cell row number

Definition at line 3804 of file WatMgmt.c.

References Abs, and UTM_NOrigin.

Referenced by ReadChanStruct().

03805 { return ( Abs(UTM - UTM_NOrigin)/celWid );              
03806 }

float UTM2kmy double  UTM  ) 
 

Converter from meters in UTM to grid cell column location (zero origin).

Parameters:
UTM Geographic coordinate (m, ELM uses UTM zone17, NAD'27, but not important here)
Returns:
Grid cell column number

Definition at line 3795 of file WatMgmt.c.

References Abs, and UTM_EOrigin.

Referenced by ReadChanStruct().

03796 { return ( Abs(UTM - UTM_EOrigin)/celWid );  
03797 }

double julday int  mon,
int  day,
int  year,
int  h,
int  mi,
double  se
 

Determine the Julian calendar day from a Gregorian date.

Returns:
julian day + hms as a real number
Takes a date, and returns a Julian day. A Julian day is the number of days since some base date (in the very distant past). Handy for getting date of x number of days after a given Julian date (use jdate to get that from the Gregorian date). Author: Robert G. Tantzen, translator: Nat Howard Translated from the algol original in Collected Algorithms of CACM (This and jdate are algorithm 199). All taken (unmodified) from SFWMD HSM's /vol/hsm/src/libs/xmgr_julday/ directory.

Parameters:
mon Month
day Day
year Year
h Hour
mi Minute
se Second
Returns:
Julian day

Definition at line 1500 of file Driver_Utilities.c.

Referenced by get_parmf(), readSeriesCol(), and ReadStructures().

01501 {
01502   long m = mon, d = day, y = year;
01503   long c, ya, j;
01504   double seconds = h * 3600.0 + mi * 60 + se;
01505 
01506   if (m > 2)
01507     m -= 3;
01508   else {
01509     m += 9;
01510     --y;
01511   }
01512   c = y / 100L;
01513   ya = y - (100L * c);
01514   j = (146097L * c) / 4L + (1461L * ya) / 4L + (153L * m + 2L) / 5L + d + 1721119L;
01515   if (seconds < 12 * 3600.0) {
01516     j--;
01517     seconds += 12.0 * 3600.0;
01518   } 
01519   else {
01520     seconds = seconds - 12.0 * 3600.0;
01521   }
01522   return (j + (seconds / 3600.0) / 24.0);
01523 }

void init_pvar VOIDP  Map,
UCHAR mask,
unsigned char  Mtype,
float  iv
 

Initialize a variable to a value.

Parameters:
Map array of data
mask data mask
Mtype the data type of the map data
iv the value used to initialize variable

Definition at line 925 of file Driver_Utilities.c.

00926 {
00927   int i0, i1;
00928   
00929   switch(Mtype) {
00930   case 'b' :    /* added double for (non-map) basin (b) array budget calcs */
00931     for(i0=0; i0<=numBasn; i0++) {
00932         ((double*)Map)[i0] = iv;
00933       }
00934     break;
00935   case 'l' :    /* added double (l == letter "ell" ) for map arrays */
00936     for(i0=0; i0<=s0+1; i0++) 
00937       for(i1=0; i1<=s1+1; i1++) {
00938         if(mask==NULL) ((double*)Map)[T(i0,i1)] = iv;
00939         else if ( mask[T(i0,i1)] == 0 ) ((double*)Map)[T(i0,i1)] = 0;
00940         else ((double*)Map)[T(i0,i1)] = iv;
00941       }
00942     break;
00943   case 'f' :    
00944     for(i0=0; i0<=s0+1; i0++) 
00945       for(i1=0; i1<=s1+1; i1++) {
00946         if(mask==NULL) ((float*)Map)[T(i0,i1)] = iv;
00947         else if ( mask[T(i0,i1)] == 0 ) ((float*)Map)[T(i0,i1)] = 0;
00948         else ((float*)Map)[T(i0,i1)] = iv;
00949       }
00950     break;
00951   case 'i' :    case 'd' :      
00952     for(i0=0; i0<=s0+1; i0++) 
00953       for(i1=0; i1<=s1+1; i1++) {
00954         if(mask==NULL) ((int*)Map)[T(i0,i1)] = (int)iv;
00955         else if ( mask[T(i0,i1)] == 0 ) ((int*)Map)[T(i0,i1)] = 0;
00956         else ((int*)Map)[T(i0,i1)] = (int)iv;
00957       }
00958     break;
00959   case 'c' :    
00960     for(i0=0; i0<=s0+1; i0++) 
00961       for(i1=0; i1<=s1+1; i1++) {
00962         if(mask==NULL) ((unsigned char*)Map)[T(i0,i1)] = (unsigned char) '\0' + (int) iv;
00963         else if ( mask[T(i0,i1)] == 0 ) ((unsigned char*)Map)[T(i0,i1)] = '\0';
00964         else ((unsigned char*)Map)[T(i0,i1)] = (unsigned char) '\0' + (int) iv;
00965       } 
00966     break;
00967  
00968    default :
00969     printf(" in default case\n");
00970    break;
00971   }
00972 }

VOIDP nalloc unsigned  mem_size,
const char  var_name[]
 

Allocate memory for a variable.

Parameters:
mem_size The size of memory space
var_name The variable's name
Returns:
Pointer to that memory

Definition at line 1774 of file Driver_Utilities.c.

01775 {
01776   VOIDP rp;
01777 
01778   
01779   if(mem_size == 0) return(NULL);
01780   rp = (VOIDP)malloc( mem_size );
01781   total_memory += mem_size;
01782   fasync(stderr);
01783   if( rp == NULL ) {
01784     fprintf(stderr,"Sorry, out of memory(%d): %s\n",mem_size,var_name);
01785     Exit(0);
01786   }
01787   fmulti(stderr);
01788   return(rp);
01789 }

char* Scip char *  s,
char  SYM
 

Skip ahead in a string until next field.

Parameters:
s A character array
SYM The specific symbol being used to skip over
Returns:
Location in string

Definition at line 1705 of file Driver_Utilities.c.

Referenced by BIRinit(), get_hab_parm(), Read_schedule(), and ReadStructures().

01706 {
01707   if(*s == SYM ) return ++s;
01708   while(*s != SYM && *s != EOS ) s++;
01709   if(*s != EOS) return ++s;
01710   else  return s;
01711 }

int isInteger char *  target_str  ) 
 

Determine if an integer is present in string.

Parameters:
target_str Target string
Returns:
true or false (found int or not)

Definition at line 1392 of file Driver_Utilities.c.

Referenced by ReadStructures(), and readViewParms().

01392                                  {
01393 
01394         int i=-1,first_num=0;
01395         char ch;
01396 
01397         while( (ch=target_str[++i]) != '\0' ) { 
01398                         if( isdigit(ch) ) first_num=1;
01399                         if( (ch=='-' || ch=='+') && first_num ) return(0);  
01400                         if( !( isspace(ch) || isdigit(ch) || ch=='-' || ch=='+') ) return(0);  
01401         }
01402         return(1);              
01403 }

int isFloat char *  target_str  ) 
 

Determine if a float is present in string.

Parameters:
target_str Target string
Returns:
true or false (found float or not)

Definition at line 1410 of file Driver_Utilities.c.

Referenced by ReadStructures(), and readViewParms().

01410                                {
01411 
01412         int i=-1,first_num=0;
01413         char ch;
01414 
01415         while( (ch=target_str[++i]) != '\0' ) { 
01416                         if( isdigit(ch) ) first_num=1;
01417                         if( (ch=='-' || ch=='+') && first_num ) return(0);  
01418                         if( !( isspace(ch) || isdigit(ch) || ch=='-' || ch=='+' || ch=='.' || toupper(ch) == 'E') ) return(0);  
01419         }
01420         return(1);              
01421 }

float FMOD float  x,
float  y
 

Modulus of a pair of (double) arguments.

Parameters:
x Numerator
y Denominator
Returns:
modulus (float) result

Definition at line 1592 of file Driver_Utilities.c.

01592                              { 
01593 return (float)fmod((double)x, (double)y); 
01594 } 

int Round float  x  ) 
 

truncate (not rounding) to an integer

Parameters:
x the value being operated on
Returns:
result

Definition at line 1694 of file Driver_Utilities.c.

Referenced by Channel_configure().

01695 { 
01696   int i = (int)x;
01697   return (i);
01698 }

float Flux_SWcells int  i0,
int  i1,
int  j0,
int  j1,
float *  SWater,
float *  Elevation,
float *  MC
 

Surface water flux calculations.

Application of Manning's eqn to calculate flux between two adjacent cells (i0,i1) and (j0,j1) Returns height flux. Flux is positive if flow is from i to j. Checks for available volume, and that flow cannot make the head in recepient cell higher than in the donor one.

Parameters:
i0 Row of cell from (positive flow)
i1 Column of cell from (positive flow)
j0 Row of cell to (positive flow)
j1 Column of cell to (positive flow)
SWater SURFACE_WAT
Elevation SED_ELEV
MC HYD_MANNINGS_N
Returns:
Flux Water flux (m)

Variables local to function
dh, adh The head difference, and absolute value of head difference (m)
MC_cells The mean manning's coefficient for a pair of cells or a cell
Hi, Hj The stage height in the i'th and j'th cells, respectively (m)
cellLoci, cellLocj The array location of the i'th and j'th cells, respectively

Definition at line 217 of file Fluxes.c.

References Abs, BCondFlow, boundcond_depth, GP_DetentZ, GP_mannDepthPow, GP_mannHeadPow, Min, ON_MAP, ramp, and T.

Referenced by Flows_in_Structures(), and Flux_SWater().

00218 {
00225   float dh, adh;
00226   float MC_cells;
00227   float Hi, Hj;
00228   float Flux = 0.;
00229   int cellLoci = T(i0,i1);
00230   int cellLocj = T(j0,j1);
00231 
00232    
00233   MC_cells = (MC[cellLoci] + MC[T(j0,j1)])/2.;
00234 
00235   /* If an on-map cell is marked 3, we are at a model boundary allowing surface water exchange */     
00236   if (!ON_MAP[cellLoci] && BCondFlow[cellLocj] == 3 ) 
00237   {
00238     /* the off-map cell given head 5 cm less than donor */
00239     /* Hi = Elevation[cellLocj] + Max(SWater[cellLocj]-0.05,0.0) ;  */ /* v2.3 not using this Hi */ 
00240 
00241     /* new dynamic boundary condition stage */
00242     Hi = boundcond_depth[cellLoci]; /* ?v2.2 not using this Hi */ 
00243 
00244     MC_cells = MC[cellLocj];       /* the mannings n is not avg, but the value of onmap boundary cell */
00245   }
00246 
00247   else 
00248     Hi = SWater[cellLoci] + Elevation[cellLoci];
00249 
00250   if(BCondFlow[cellLoci] == 3 && !ON_MAP[cellLocj] )   
00251   {  
00252     /* the off-map cell given head 5 cm less than donor  */
00253       /* Hj = Elevation[cellLoci] + Max(SWater[cellLoci]-0.05,0.0); */   /* v2.3 not using this Hj */ 
00254 
00255     Hj = boundcond_depth[cellLocj]; /* ?v2.2 not using this Hj */
00256 
00257     MC_cells = MC[cellLoci];      /* the mannings n is not avg, but the value of onmap boundary cell */
00258   }
00259 
00260   else 
00261     Hj = SWater[cellLocj] + Elevation[cellLocj];
00262 
00263   dh = Hi - Hj;         /* dh is "from --> to" */
00264   adh = Abs (dh);
00265                 
00266   if (dh > 0) 
00267   {
00268     if(SWater[cellLoci] < GP_DetentZ) 
00269       return 0.0; 
00270     /* step_Cell is constant calc'd in Generic_Driver.c at initialization ( m^(-1.5) * sec )*/
00271     Flux = (MC_cells != 0) ? 
00272         (pow(adh,GP_mannHeadPow) / MC_cells  * pow(SWater[cellLoci],GP_mannDepthPow)*step_Cell) : (0.0);
00273 
00274     /* ensure adequate volume avail */
00275     Flux =  ( Flux > ramp(SWater[cellLoci] - GP_DetentZ) ) ? (ramp(SWater[cellLoci] - GP_DetentZ)) : (Flux);
00276 
00277     /* check to ensure no flip/flops associated with depth */
00278     if ( ( Hi - Flux ) < ( Hj + Flux ) )        
00279       Flux = Min ( dh/2.0, ramp(SWater[cellLoci] - GP_DetentZ) );
00280 
00281   } /* end of dh > 0 */
00282 
00283   else
00284   {     
00285     if (SWater[cellLocj] < GP_DetentZ) 
00286       return 0.0; 
00287     /* step_Cell is constant calc'd in Generic_Driver.c at initialization ( m^(-1.5) * sec )*/
00288     /* Flux is negative in this case */
00289     Flux = (MC_cells != 0) ? 
00290            ( - pow(adh,GP_mannHeadPow) / MC_cells  * pow(SWater[cellLocj],GP_mannDepthPow)*step_Cell) : (0.0);
00291 
00292     /* ensure adequate volume avail */
00293     Flux =  ( -Flux > ramp(SWater[cellLocj] - GP_DetentZ) ) ? (-ramp(SWater[cellLocj] - GP_DetentZ)) : (Flux);
00294 
00295     /* check to ensure no flip/flops associated with depth */           
00296     if ( ( Hi - Flux ) > ( Hj + Flux ) )  
00297       Flux = - Min ( adh/2.0, ramp(SWater[cellLocj] - GP_DetentZ) );
00298   }
00299         
00300   return (Flux);        /* returns height flux */
00301 }

void Flux_SWstuff int  i0,
int  i1,
int  j0,
int  j1,
float  Flux,
float *  SURFACE_WAT,
double *  STUF1,
double *  STUF2,
double *  STUF3
 

Flux of surface water constituents.

The constituents (nutrients, salinity, etc), are passed among cells, and surface water variable is updated, along with making budget calcs. These are units of constituent mass being fluxed from i0,i1 to j0,j1 Flux & SURFACE_WAT are in units of height (m)

Parameters:
i0 Row of cell from (positive flow)
i1 Column of cell from (positive flow)
j0 Row of cell to (positive flow)
j1 Column of cell to (positive flow)
Flux Water flux (m) between grid cells
SURFACE_WAT SURFACE_WAT
STUF1 SALT_SURF_WT
STUF2 DINdummy
STUF3 TP_SF_WT

Definition at line 324 of file Fluxes.c.

References basins, basn, basn_list, debug, Disp_Calc(), dynERRORnum, basndef::family, basndef::FLok, Max, Min, msgStr, basndef::numFLok, ON_MAP, P_IN_OVL, P_OUT_OVL, basndef::parent, S_IN_OVL, S_OUT_OVL, sfstep, SimTime, SURFACE_WAT, T, simTime::TIME, True, VOL_IN_OVL, VOL_OUT_OVL, and WriteMsg().

Referenced by Flows_in_Structures(), and Flux_SWater().

00325 {
00326     float m1=0.0, m2=0.0, m3=0.0;
00327     int  ii, flo_chek, cel_i, cel_j;
00328     float disp_num; /* numerical dispersion (m2/d) */
00329     float veloc;    /* velocity (m/d) */
00330     float velocAdj;  /* velocity adjusted for numerical dispersion (m/d) */
00331     float FluxAdj; /* Flux adjusted for numerical dispersioin */
00332     double fl_prop_i, fl_prop_j; /* proportion of surface water constituents that should be fluxed with water (depends on magnitude of dispersion) */
00333     extern float dispParm_scaled; 
00334 
00335     cel_i = T(i0,i1);
00336     cel_j = T(j0,j1);
00337     
00338 /*    veloc =  Abs(Flux) * celWid/( (Flux>0.0)?(SURFACE_WAT[cel_i]):(SURFACE_WAT[cel_j])) / (sfstep) ;  
00339     disp_num = 0.5 * veloc * (celWid - veloc * sfstep)  ; 
00340     velocAdj = (veloc * celWid - disp_num)/celWid; 
00341     FluxAdj = dispParm_scaled * velocAdj * sfstep * ( (Flux>0.0)?(SURFACE_WAT[cel_i]):(SURFACE_WAT[cel_j]))/celWid;
00342 */
00343     FluxAdj = Disp_Calc(Flux, SURFACE_WAT[cel_i], SURFACE_WAT[cel_j], sfstep);
00344     
00345     fl_prop_i = (SURFACE_WAT[cel_i]>0.0) ? (Max(Flux-FluxAdj,0.0)/SURFACE_WAT[cel_i]) : (0.0); /* pos Flux */
00346     fl_prop_j = (SURFACE_WAT[cel_j]>0.0) ? (Min(Flux+FluxAdj,0.0)/SURFACE_WAT[cel_j]) : (0.0); /* neg Flux */
00347         /*  if (i0==41 && i1==12 && veloc>100.0) printf ("\nveloc=%f m/d, disp_num=%f m2/d, FF=%f m; FFadj=%f m \n",veloc, disp_num, Flux, FluxAdj); */
00348     fl_prop_i = Min(fl_prop_i, 1.0);
00349     fl_prop_j = Min(fl_prop_j, 1.0);
00350  
00351  if (Flux >0.0) {
00352      m1 = STUF1[cel_i]*fl_prop_i;
00353      m2 = STUF2[cel_i]*fl_prop_i;
00354      m3 = STUF3[cel_i]*fl_prop_i;
00355  }
00356  else  {
00357      m1 = STUF1[cel_j]*fl_prop_j;
00358      m2 = STUF2[cel_j]*fl_prop_j;
00359      m3 = STUF3[cel_j]*fl_prop_j;
00360  }
00361         
00362  STUF1[cel_j] += m1;  /* add the masses of constituents */
00363  STUF2[cel_j] += m2;
00364  STUF3[cel_j] += m3;
00365  STUF1[cel_i] -= m1;
00366  STUF2[cel_i] -= m2;
00367  STUF3[cel_i] -= m3;
00368  SURFACE_WAT[cel_j] += Flux; /* now update the surfwater depths */
00369  SURFACE_WAT[cel_i] -= Flux;
00370 
00371  if (debug > 2) {
00372      if (STUF3[cel_j] < -GP_MinCheck) {
00373          sprintf(msgStr,"Day %6.1f: capacityERR - neg surface water P (%f kg) in cell (%d,%d) after SWfluxes!", 
00374                  SimTime.TIME, STUF3[cel_j], j0,j1 ); 
00375          WriteMsg( msgStr,True );  dynERRORnum++;}
00376      if (STUF3[cel_i] < -GP_MinCheck) {
00377          sprintf(msgStr,"Day %6.1f: capacityERR - neg surface water P (%f kg) in cell (%d,%d) after SWfluxes!", 
00378                  SimTime.TIME, STUF3[cel_i], i0,i1 ); 
00379          WriteMsg( msgStr,True );  dynERRORnum++; }
00380      if (STUF1[cel_j] < -GP_MinCheck) {
00381          sprintf(msgStr,"Day %6.1f: capacityERR - neg surface water S (%f kg) in cell (%d,%d) after SWfluxes!", 
00382                  SimTime.TIME, STUF1[cel_j], j0,j1 ); 
00383          WriteMsg( msgStr,True );  dynERRORnum++; }
00384      if (STUF1[cel_i] < -GP_MinCheck) {
00385          sprintf(msgStr,"Day %6.1f: capacityERR - neg surface water S (%f kg) in cell (%d,%d) after SWfluxes!", 
00386                  SimTime.TIME, STUF1[cel_i], i0,i1 ); 
00387          WriteMsg( msgStr,True );  dynERRORnum++; }
00388      if (SURFACE_WAT[cel_j] < -GP_MinCheck) {
00389          sprintf(msgStr,"Day %6.1f: capacityERR - negative surface water (%f m) in cell (%d,%d)!", 
00390                  SimTime.TIME, SURFACE_WAT[cel_j], j0,j1 ) ; 
00391          WriteMsg( msgStr,True );  dynERRORnum++;}
00392 
00393      if (SURFACE_WAT[cel_i] < -GP_MinCheck) {
00394          sprintf(msgStr,"Day %6.1f: capacityERR - negative surface water (%f m) in cell (%d,%d)!", 
00395                  SimTime.TIME, SURFACE_WAT[cel_i], i0,i1 ) ; 
00396          WriteMsg( msgStr,True );  dynERRORnum++; }
00397  }
00398         
00399 
00400 /* mass balance sums */
00401  if (basn[cel_i] != basn[cel_j])  { 
00402 
00403         /* first do the normal case where all cells are on-map */
00404      if ( ON_MAP[cel_j] && ON_MAP[cel_i] ) { 
00405          if ( Flux > 0  ) { /* positive fluxes out of basn[cel_i] */
00406              if (basn_list[basn[cel_i]]->family !=  
00407                  basn_list[basn[cel_j]]->family ) {        /* if the flow is not within the family... */
00408                  if  ( !basn_list[basn[cel_i]]->parent  ) { /* and if the donor or recipient is a child... */
00409                      /* then find out about the flow for the family's sake */
00410                      VOL_OUT_OVL[basn_list[basn[cel_i]]->family] += Flux*CELL_SIZE;
00411                      P_OUT_OVL[basn_list[basn[cel_i]]->family] += m3;
00412                      S_OUT_OVL[basn_list[basn[cel_i]]->family] += m1;
00413                  }
00414                  if ( !basn_list[basn[cel_j]]->parent ) {                 
00415                      /* then find out about the flow for the family's sake */
00416                      VOL_IN_OVL[basn_list[basn[cel_j]]->family] += Flux*CELL_SIZE;
00417                      P_IN_OVL[basn_list[basn[cel_j]]->family] += m3;
00418                      S_IN_OVL[basn_list[basn[cel_j]]->family] += m1;
00419                  }
00420                      /* now sum the parents' flows */
00421                  VOL_OUT_OVL[basn[cel_i]] += Flux*CELL_SIZE; 
00422                  P_OUT_OVL[basn[cel_i]] += m3; 
00423                  S_OUT_OVL[basn[cel_i]] += m1; 
00424                  VOL_IN_OVL[basn[cel_j]]  += Flux*CELL_SIZE; 
00425                  P_IN_OVL[basn[cel_j]]  += m3; 
00426                  S_IN_OVL[basn[cel_j]]  += m1; 
00427                  
00428              }
00429              else {    /* if it's flow within a family, just keep
00430                           track of what the children do among themselves */            
00431                  if  ( !basn_list[basn[cel_i]]->parent  ) { 
00432                      VOL_OUT_OVL[basn[cel_i]] += Flux*CELL_SIZE; 
00433                      P_OUT_OVL[basn[cel_i]] += m3; 
00434                      S_OUT_OVL[basn[cel_i]] += m1; 
00435                  }
00436                  if ( !basn_list[basn[cel_j]]->parent ) {                 
00437                      VOL_IN_OVL[basn[cel_j]]  += Flux*CELL_SIZE; 
00438                      P_IN_OVL[basn[cel_j]]  += m3; 
00439                      S_IN_OVL[basn[cel_j]]  += m1; 
00440                  }
00441              }
00442                         
00443              if (debug > 0 && WatMgmtOn) { /* check for basin misconfiguration (allowable basin-basin flows) */
00444                  basins = basn_list[basn[cel_i]];
00445                  flo_chek = 0;
00446                  for (ii=0; ii<basins->numFLok; ii++) { if (basn[cel_j] == basins->FLok[ii] ) flo_chek = 1; }
00447                  if (!flo_chek) {
00448                      sprintf(msgStr,"Day %5.3f: ERROR - no (pos) SW flow from cell (%d,%d) of basin %d into cell (%d,%d) of basin %d!", 
00449                              SimTime.TIME, i0,i1,basn[cel_i], j0,j1, basn[cel_j]); 
00450                      WriteMsg( msgStr,True );  dynERRORnum++; }
00451              }
00452              
00453 
00454          }
00455          else { /* negative fluxes out of basn[cel_j] */
00456              if (basn_list[basn[cel_i]]->family !=  
00457                  basn_list[basn[cel_j]]->family ) {        /* if the flow is not within the family... */
00458                  if  ( !basn_list[basn[cel_j]]->parent  ) { /* and if the donor or recipient is a child... */
00459                      /* then find out about the flow for the family's sake */
00460                      VOL_OUT_OVL[basn_list[basn[cel_j]]->family] -= Flux*CELL_SIZE;
00461                      P_OUT_OVL[basn_list[basn[cel_j]]->family] -= m3;
00462                      S_OUT_OVL[basn_list[basn[cel_j]]->family] -= m1;
00463                  }
00464                  if ( !basn_list[basn[cel_i]]->parent ) {                 
00465                      /* then find out about the flow for the family's sake */
00466                      VOL_IN_OVL[basn_list[basn[cel_i]]->family] -= Flux*CELL_SIZE;
00467                      P_IN_OVL[basn_list[basn[cel_i]]->family] -= m3;
00468                      S_IN_OVL[basn_list[basn[cel_i]]->family] -= m1;
00469                  }
00470                      /* now sum the parents' flows */
00471                  VOL_OUT_OVL[basn[cel_j]] -= Flux*CELL_SIZE; 
00472                  P_OUT_OVL[basn[cel_j]] -= m3; 
00473                  S_OUT_OVL[basn[cel_j]] -= m1; 
00474                  VOL_IN_OVL[basn[cel_i]]  -= Flux*CELL_SIZE; 
00475                  P_IN_OVL[basn[cel_i]]  -= m3; 
00476                  S_IN_OVL[basn[cel_i]]  -= m1; 
00477                  
00478              }
00479              else {    /* if it's flow within a family, just keep
00480                           track of what the children do among themselves */            
00481                  if  ( !basn_list[basn[cel_j]]->parent  ) { 
00482                      VOL_OUT_OVL[basn[cel_j]] -= Flux*CELL_SIZE; 
00483                      P_OUT_OVL[basn[cel_j]] -= m3; 
00484                      S_OUT_OVL[basn[cel_j]] -= m1; 
00485                  }
00486                  if ( !basn_list[basn[cel_i]]->parent ) {                 
00487                      VOL_IN_OVL[basn[cel_i]]  -= Flux*CELL_SIZE; 
00488                      P_IN_OVL[basn[cel_i]]  -= m3; 
00489                      S_IN_OVL[basn[cel_i]]  -= m1; 
00490                  }
00491              }
00492 
00493 
00494              if (debug > 0 && WatMgmtOn) { /* check for basin misconfiguration (allowable basin-basin flows) */
00495                  basins = basn_list[basn[cel_j]];
00496                  flo_chek = 0;
00497                  for (ii=0; ii<basins->numFLok; ii++) { if (basn[cel_i] == basins->FLok[ii] ) flo_chek = 1; }
00498                  if (!flo_chek) {
00499                      sprintf(msgStr,"Day %5.3f: ERROR - no (neg) SW flow from cell (%d,%d) of basin %d into cell (%d,%d) of basin %d!", 
00500                              SimTime.TIME, i0,i1, basn[cel_j], j0,j1, basn[cel_i]); 
00501                      WriteMsg( msgStr,True );  dynERRORnum++; }
00502              }
00503                  
00504              
00505          }
00506      }
00507      else  if ( !ON_MAP[cel_j]) /* so now the j,j cell is off-map, recipient if pos flow */
00508          if (Flux > 0) { 
00509              if ( !basn_list[basn[cel_i]]->parent ) { /* child's play */
00510                  VOL_OUT_OVL[basn_list[basn[cel_i]]->family] += Flux*CELL_SIZE;
00511                  P_OUT_OVL[basn_list[basn[cel_i]]->family] += m3;
00512                  S_OUT_OVL[basn_list[basn[cel_i]]->family] += m1;
00513              }
00514              /* parents' play */
00515              VOL_OUT_OVL[basn[cel_i]] += Flux*CELL_SIZE; 
00516              VOL_OUT_OVL[0]+= Flux*CELL_SIZE;
00517              P_OUT_OVL[basn[cel_i]] += m3; 
00518              P_OUT_OVL[0]+= m3;
00519              S_OUT_OVL[basn[cel_i]] += m1; 
00520              S_OUT_OVL[0]+= m1;
00521          }
00522          else { /* negative flows */
00523              if ( !basn_list[basn[cel_i]]->parent ) {/* child's play */
00524                  VOL_IN_OVL[basn_list[basn[cel_i]]->family] -= Flux*CELL_SIZE;
00525                  P_IN_OVL[basn_list[basn[cel_i]]->family] -= m3;
00526                  S_IN_OVL[basn_list[basn[cel_i]]->family] -= m1;
00527              }
00528              /* parents' play */
00529              VOL_IN_OVL[basn[cel_i]]  -= Flux*CELL_SIZE;
00530              VOL_IN_OVL[0]               -= Flux*CELL_SIZE;
00531              P_IN_OVL[basn[cel_i]]  -= m3;
00532              P_IN_OVL[0]               -= m3;
00533              S_IN_OVL[basn[cel_i]]  -= m1;
00534              S_IN_OVL[0]               -= m1;
00535          }
00536      else  if ( !ON_MAP[cel_i]) /* so now the i,i cell is off-map, donor if pos flow */
00537          if (Flux > 0) { 
00538              if ( !basn_list[basn[cel_j]]->parent ) {/* child's play */
00539                  VOL_IN_OVL[basn_list[basn[cel_j]]->family] += Flux*CELL_SIZE;
00540                  P_IN_OVL[basn_list[basn[cel_j]]->family] += m3;
00541                  S_IN_OVL[basn_list[basn[cel_j]]->family] += m1;
00542              }
00543              /* parents' play */
00544              VOL_IN_OVL[basn[cel_j]] += Flux*CELL_SIZE;
00545              VOL_IN_OVL[0]              += Flux*CELL_SIZE;
00546              P_IN_OVL[basn[cel_j]] += m3;
00547              P_IN_OVL[0]              += m3;
00548              S_IN_OVL[basn[cel_j]] += m1;
00549              S_IN_OVL[0]              += m1;
00550          }
00551          else { /*negative flows */
00552              if ( !basn_list[basn[cel_j]]->parent ) {/* child's play */
00553                  VOL_OUT_OVL[basn_list[basn[cel_j]]->family] -= Flux*CELL_SIZE;
00554                  P_OUT_OVL[basn_list[basn[cel_j]]->family] -= m3;
00555                  S_OUT_OVL[basn_list[basn[cel_j]]->family] -= m1;
00556              }
00557              /* parents' play */
00558              VOL_OUT_OVL[basn[cel_j]]  -= Flux*CELL_SIZE;
00559              VOL_OUT_OVL[0]               -= Flux*CELL_SIZE;
00560              P_OUT_OVL[basn[cel_j]]  -= m3;
00561              P_OUT_OVL[0]               -= m3;
00562              S_OUT_OVL[basn[cel_j]]  -= m1;
00563              S_OUT_OVL[0]               -= m1;
00564          }
00565  }
00566                 
00567  return;
00568         
00569 }

Here is the call graph for this function:


Variable Documentation

FILE* schedFile
 

File pointer to the schedule data file

Definition at line 34 of file watmgmt.h.

Referenced by Read_schedule().

FILE* ChanInFile
 

File pointer to several (opened-closed) input files

Definition at line 35 of file watmgmt.h.

Referenced by Canal_Network_Init(), ReadChanStruct(), and ReadStructures().

FILE* CanalOutFile
 

File pointer to canal hydrologic output file

Definition at line 36 of file watmgmt.h.

Referenced by Canal_Network_Init(), and Run_Canal_Network().

FILE* CanalOutFile_P
 

File pointer to canal phosphorus output file

Definition at line 37 of file watmgmt.h.

Referenced by Canal_Network_Init(), and Run_Canal_Network().

FILE* CanalOutFile_S
 

File pointer to canal salt/tracer output file

Definition at line 38 of file watmgmt.h.

Referenced by Canal_Network_Init(), and Run_Canal_Network().

FILE* WstructOutFile
 

File pointer to water control structure hydrologic output file

Definition at line 39 of file watmgmt.h.

Referenced by Canal_Network_Init(), Flows_in_Structures(), and Run_Canal_Network().

FILE* WstructOutFile_P
 

File pointer to water control structure phosphorus output file

Definition at line 40 of file watmgmt.h.

Referenced by Canal_Network_Init(), Flows_in_Structures(), and Run_Canal_Network().

FILE* WstructOutFile_S
 

File pointer to water control structure salt/tracer output file

Definition at line 41 of file watmgmt.h.

Referenced by Canal_Network_Init(), Flows_in_Structures(), and Run_Canal_Network().

FILE* canDebugFile
 

File pointer to canal debug file

Definition at line 42 of file watmgmt.h.

Referenced by Canal_Network_Init(), Flows_in_Structures(), FluxChannel(), and Run_Canal_Network().

FILE* CanalCellInteractDebug
 

File pointer to canal-cell interaction file

Definition at line 43 of file watmgmt.h.

Referenced by Channel_configure(), getCanalElev(), and MarkCell().

FILE* F_struct_wat
 

File pointer to the water control structure flow data file

Definition at line 44 of file watmgmt.h.

Referenced by ReadStructures().

FILE* F_struct_TP
 

File pointer to the water control structure TP data file

Definition at line 45 of file watmgmt.h.

Referenced by ReadStructures().

FILE* F_struct_TS
 

File pointer to the water control structure TS (salt/tracer) data file

Definition at line 46 of file watmgmt.h.

Referenced by ReadStructures().

int CHo
 

the number of the first canal in Chan.file sequence

Definition at line 48 of file watmgmt.h.

Referenced by Channel_configure(), FluxChannel(), getCanalElev(), ReadChanStruct(), and ReadStructures().

int CHANNEL_MAX_ITER
 

maximal number of iterations in relaxation procedure (input data)

Definition at line 49 of file watmgmt.h.

Referenced by FluxChannel(), and ReadChanStruct().

int UTM_EOrigin
 

Definition at line 50 of file watmgmt.h.

Referenced by ReadChanStruct(), and UTM2kmy().

int UTM_NOrigin
 

Map origin, easting and northing UTM

Definition at line 50 of file watmgmt.h.

Referenced by ReadChanStruct(), and UTM2kmx().

int printchan = 0
 

flag to indicate time to print canal/struct data

Definition at line 51 of file watmgmt.h.

Referenced by Run_Canal_Network().

int N_c_iter = 0
 

cumulative number of iterations through canal code

Definition at line 52 of file watmgmt.h.

Referenced by Run_Canal_Network().

int num_cell = 0
 

counter of number cells along a reach segment

Definition at line 53 of file watmgmt.h.

Referenced by Channel_configure().

int C_Mark = 0
 

indicator that a cell is marked with a canal already

Definition at line 54 of file watmgmt.h.

Referenced by Channel_configure().

float C_F
 

acceleration of overland flow for canal/cell interactions - for sensitivity experiments only

Definition at line 56 of file watmgmt.h.

Referenced by ReadChanStruct().

float F_ERROR
 

convergence error for relaxation procedure in canals

Definition at line 57 of file watmgmt.h.

Referenced by FluxChannel(), and ReadChanStruct().

float MINDEPTH
 

minimal allowed water depth in canal

Definition at line 58 of file watmgmt.h.

Referenced by Flows_in_Structures(), FluxChannel(), and ReadChanStruct().

int num_chan
 

number of channels (canals)

Definition at line 60 of file watmgmt.h.

Referenced by Channel_configure(), ReadChanStruct(), and ReadStructures().

int num_struct_hist
 

number of structures driven by data

Definition at line 61 of file watmgmt.h.

Referenced by ReadStructures().

int numTPser
 

number of structures with hist TP data

Definition at line 62 of file watmgmt.h.

Referenced by ReadStructures().

int numTSser
 

number of structures hist TS (TotalSalt/tracer) data

Definition at line 63 of file watmgmt.h.

Referenced by ReadStructures().

float* MCopen
 

manning's coef for open water, fixed parameter for all domain (used in structure flow for bridges)

Definition at line 65 of file watmgmt.h.

Referenced by Canal_Network_Init(), and Flows_in_Structures().

char modelFileName[300]
 

generic string for input/output file name

Definition at line 66 of file watmgmt.h.

struct Chan** Chan_list
 

Definition at line 121 of file watmgmt.h.

Referenced by Canal_Network_Init(), CanalReInit(), Channel_configure(), Flows_in_Structures(), FluxChannel(), getCanalElev(), ReadStructures(), and Run_Canal_Network().

struct Cells* cell_last
 

Definition at line 135 of file watmgmt.h.

Referenced by Channel_configure(), and MarkCell().

struct Cells* cell
 

Definition at line 136 of file watmgmt.h.

Referenced by Channel_configure(), FluxChannel(), getCanalElev(), HabSwitch(), and MarkCell().

struct Structure* struct_first
 

Definition at line 180 of file watmgmt.h.

Referenced by ReadStructures().

struct Structure* structs
 

Definition at line 181 of file watmgmt.h.

Referenced by Canal_Network_Init(), Flows_in_Structures(), FluxChannel(), and ReadStructures().

struct HistData Hist_data[MAX_H_STRUCT]
 

pointer to array of pointers for the water control structure time series

Definition at line 209 of file watmgmt.h.

Referenced by Flows_in_Structures(), and ReadStructures().

char* ModelPath
 

Environment variables used in model

ModelPath environment variable, base pathname for executable and input data
ProjName environment variable, the name of the model project
DriverPath environment variable, base pathname for source code
OS_TYPE environment variable, the type of operating system being used (informational purpose only, not used in code)

Definition at line 36 of file driver_utilities.h.

char * OutputPath
 

base pathname for all model output (user input)

Definition at line 38 of file driver_utilities.h.

char * ProjName
 

Environment variables used in model

ModelPath environment variable, base pathname for executable and input data
ProjName environment variable, the name of the model project
DriverPath environment variable, base pathname for source code
OS_TYPE environment variable, the type of operating system being used (informational purpose only, not used in code)

Definition at line 36 of file driver_utilities.h.

char modelName[20]
 

Model name/version (user input)

modelName Name given to model implementation (user input)
modelVers Version given to model implementation (e.g., v.2.1) (user input)

Definition at line 44 of file driver_utilities.h.

char modelVers[10]
 

Model name/version (user input)

modelName Name given to model implementation (user input)
modelVers Version given to model implementation (e.g., v.2.1) (user input)

Definition at line 44 of file driver_utilities.h.

char SimAlt[20]
 

simulation scenario/alterative name

Definition at line 35 of file generic_driver.h.

Referenced by BIRoutfiles(), Canal_Network_Init(), get_parmf(), open_point_lists(), Read_schedule(), ReadChanStruct(), and ReadStructures().

char SimModif[20]
 

simulation scenario/alterative modifier name/note

Definition at line 36 of file generic_driver.h.

Referenced by BIRoutfiles(), Canal_Network_Init(), get_parmf(), open_point_lists(), and ReadStructures().

double Jdate_init
 

julian day (since epochs ago) of the initialization of model

Definition at line 42 of file generic_driver.h.

Referenced by get_parmf(), send_point_lists2(), and track_time().

double Jdate_end
 

julian day (since epochs ago) of the ending of model

Definition at line 43 of file generic_driver.h.

Referenced by get_parmf().

int PORnumday
 

number of days of simulation Period Of Record (incl. start day)

Definition at line 58 of file generic_driver.h.

Referenced by get_parmf(), processData(), ReadStructures(), and track_time().

float can_Intvl
 

time (day) interval between canal (water management) summaries

Definition at line 67 of file generic_driver.h.

Referenced by get_parmf(), and track_time().

int canPrint
 

boolean flag to indicate if canal data is to be printed at current time

Definition at line 73 of file generic_driver.h.

Referenced by Run_Canal_Network(), and track_time().

unsigned char* HAB
 

_Units_: dimless; _Brief_: Habitat, or vegetation community type (integer attribute, defining database parameter lookups)

Definition at line 70 of file unitmod_vars.h.

int* basn
 

_Units_: dimless; _Brief_: Map array defining hydrologic basins and Indicator Regions

Definition at line 129 of file unitmod_vars.h.

float GP_DetentZ
 

Global parameter. _Units_: m; _Brief_: ***detention depth in a grid cell, below which surface flows do not occur; _Extended_: scale-dependent relative to topographic heterogeneity

Definition at line 43 of file unitmod_globparms.h.

float GP_MinCheck
 

Global parameter. _Units_: m; _Brief_: ***small threshold number, for relative error-checking (not a multiplier etc); _Extended_: only used in constraining fluxes at extremely mimimal conditions

Definition at line 44 of file unitmod_globparms.h.

float GP_SLRise
 

Global parameter. _Units_: m/yr; _Brief_: ***rate of Sea Level Rise; _Extended_: based on CERP Guidance Memo 016.00

Definition at line 47 of file unitmod_globparms.h.

Referenced by Flows_in_Structures(), and ReadGlobalParms().

float GP_mannDepthPow
 

Global parameter. _Units_: dimless; _Brief_: ***power used in manning's equation water depth; _Extended_: for "true" manning's, use 1.667

Definition at line 33 of file unitmod_globparms.h.

float GP_mannHeadPow
 

Global parameter. _Units_: dimless; _Brief_: ***power used in manning's equation head difference; _Extended_: for "true" manning's, use 2.0

Definition at line 34 of file unitmod_globparms.h.

float GP_calibGWat
 

Global parameter. _Units_: dimless; _Brief_: ***calibration parameter, multiply groundwater cell-cell flow calculation; _Extended_: v2.3=1.05

Definition at line 35 of file unitmod_globparms.h.

Referenced by ReadGlobalParms().

int numBasn
 

number of basins in model domain

Definition at line 32 of file budgstats.h.

basnDef** basn_list
 

Basin/Indicator Region data

Definition at line 30 of file budgstats.h.

basnDef* basins
 

Basin/Indicator Region data

Definition at line 31 of file budgstats.h.

double* TOT_VOL_CAN
 

_Units_: m^3; _Brief_:

Definition at line 36 of file budgstats_birvars.h.

Referenced by alloc_mem_stats(), BIRbudg_reset(), BIRbudg_sumFinal(), and Run_Canal_Network().

double* VOL_IN_STR
 

_Units_: m^3; _Brief_:

Definition at line 47 of file budgstats_birvars.h.

Referenced by alloc_mem_stats(), BIRbudg_print(), BIRbudg_reset(), BIRbudg_sumFinal(), and Flows_in_Structures().

double * VOL_IN_OVL
 

_Units_: m^3; _Brief_:

Definition at line 45 of file budgstats_birvars.h.

double * VOL_IN_SPG
 

_Units_: m^3; _Brief_:

Definition at line 46 of file budgstats_birvars.h.

Referenced by alloc_mem_stats(), BIRbudg_print(), BIRbudg_reset(), BIRbudg_sumFinal(), and FluxChannel().

double * VOL_IN_GW
 

_Units_: m^3; _Brief_:

Definition at line 44 of file budgstats_birvars.h.

double* VOL_OUT_STR
 

_Units_: m^3; _Brief_:

Definition at line 54 of file budgstats_birvars.h.

Referenced by alloc_mem_stats(), BIRbudg_print(), BIRbudg_reset(), BIRbudg_sumFinal(), and Flows_in_Structures().

double * VOL_OUT_OVL
 

_Units_: m^3; _Brief_:

Definition at line 52 of file budgstats_birvars.h.

double * VOL_OUT_SPG
 

_Units_: m^3; _Brief_:

Definition at line 53 of file budgstats_birvars.h.

Referenced by alloc_mem_stats(), BIRbudg_print(), BIRbudg_reset(), BIRbudg_sumFinal(), and FluxChannel().

double * VOL_OUT_GW
 

_Units_: m^3; _Brief_:

Definition at line 51 of file budgstats_birvars.h.

double* TOT_P_CAN
 

_Units_: kgP; _Brief_:

Definition at line 76 of file budgstats_birvars.h.

Referenced by alloc_mem_stats(), BIRbudg_reset(), BIRbudg_sumFinal(), and Run_Canal_Network().

double* P_IN_STR
 

_Units_: kgP; _Brief_:

Definition at line 65 of file budgstats_birvars.h.

Referenced by alloc_mem_stats(), BIRbudg_print(), BIRbudg_reset(), BIRbudg_sumFinal(), and Flows_in_Structures().

double * P_IN_OVL
 

_Units_: kgP; _Brief_:

Definition at line 63 of file budgstats_birvars.h.

double * P_IN_SPG
 

_Units_: kgP; _Brief_:

Definition at line 64 of file budgstats_birvars.h.

Referenced by alloc_mem_stats(), BIRbudg_print(), BIRbudg_reset(), BIRbudg_sumFinal(), and FluxChannel().

double * P_IN_GW
 

_Units_: kgP; _Brief_:

Definition at line 62 of file budgstats_birvars.h.

double* P_OUT_STR
 

_Units_: kgP; _Brief_:

Definition at line 73 of file budgstats_birvars.h.

Referenced by alloc_mem_stats(), BIRbudg_print(), BIRbudg_reset(), BIRbudg_sumFinal(), and Flows_in_Structures().

double * P_OUT_OVL
 

_Units_: kgP; _Brief_:

Definition at line 71 of file budgstats_birvars.h.

double * P_OUT_SPG
 

_Units_: kgP; _Brief_:

Definition at line 72 of file budgstats_birvars.h.

Referenced by alloc_mem_stats(), BIRbudg_print(), BIRbudg_reset(), BIRbudg_sumFinal(), and FluxChannel().

double * P_OUT_GW
 

_Units_: kgP; _Brief_:

Definition at line 70 of file budgstats_birvars.h.

double* TOT_S_CAN
 

_Units_: kgSalt; _Brief_:

Definition at line 146 of file budgstats_birvars.h.

Referenced by alloc_mem_stats(), BIRbudg_reset(), BIRbudg_sumFinal(), and Run_Canal_Network().

double* S_IN_STR
 

_Units_: kgSalt; _Brief_:

Definition at line 136 of file budgstats_birvars.h.

Referenced by alloc_mem_stats(), BIRbudg_print(), BIRbudg_reset(), BIRbudg_sumFinal(), and Flows_in_Structures().

double * S_IN_OVL
 

_Units_: kgSalt; _Brief_:

Definition at line 134 of file budgstats_birvars.h.

double * S_IN_SPG
 

_Units_: kgSalt; _Brief_:

Definition at line 135 of file budgstats_birvars.h.

Referenced by alloc_mem_stats(), BIRbudg_print(), BIRbudg_reset(), BIRbudg_sumFinal(), and FluxChannel().

double * S_IN_GW
 

_Units_: kgSalt; _Brief_:

Definition at line 133 of file budgstats_birvars.h.

double* S_OUT_STR
 

_Units_: kgSalt; _Brief_:

Definition at line 143 of file budgstats_birvars.h.

Referenced by alloc_mem_stats(), BIRbudg_print(), BIRbudg_reset(), BIRbudg_sumFinal(), and Flows_in_Structures().

double * S_OUT_OVL
 

_Units_: kgSalt; _Brief_:

Definition at line 141 of file budgstats_birvars.h.

double * S_OUT_SPG
 

_Units_: kgSalt; _Brief_:

Definition at line 142 of file budgstats_birvars.h.

Referenced by alloc_mem_stats(), BIRbudg_print(), BIRbudg_reset(), BIRbudg_sumFinal(), and FluxChannel().

double * S_OUT_GW
 

_Units_: kgSalt; _Brief_:

Definition at line 140 of file budgstats_birvars.h.


Generated on Thu Jul 6 11:20:51 2006 for ELM source code by  doxygen 1.3.9.1