updatong calvin.py

msdogan · Jan 11, 2024 · aeddd57 · aeddd57
1 parent 21723a2
commit aeddd57
Show file tree

Hide file tree

Showing 17 changed files with 74,499 additions and 201 deletions.
diff --git a/back-up/concrete_model/calvin/calvin.py b/back-up/concrete_model/calvin/calvin.py
@@ -1,20 +1,50 @@
+import os
+import sys
+import logging
+
 from pyomo.environ import *
 from pyomo.opt import TerminationCondition
 import numpy as np
 import pandas as pd
-import os
 
 class CALVIN():
 
-  def __init__(self, linksfile, ic=None):
+  def __init__(self, linksfile, ic=None, log_name="calvin"):
     """
     Initialize CALVIN model object.
 
     :param linksfile: (string) CSV file containing network link information
     :param ic: (dict) Initial storage conditions for surface reservoirs
                 only used for annual optimization
+    :param log_name: A name for a logger - will be used to keep logs from different model runs separate in files.
+                Defaults to "calvin", which results in a log file in the current working directory named "calvin.log".
+                You can change this each time you instantiate the CALVIN class if you want to output separate logs
+                for different runs. Otherwise, all results will be appended to the log file (not overwritten). If you
+                run multiple copies of CALVIN simultaneously, make sure to change this, or you could get errors writing
+                to the log file.
+
+                Do not provide a full path to a log file here because this value is also used in a way that is *not* a
+                file path. If being able to specify a full path is important for your workflow, please raise a GitHub
+                issue. It could be supported, but there is no need at this moment.
     :returns: CALVIN model object
     """
+
+    # set up logging code
+    self.log = logging.getLogger(log_name)
+    if not self.log.hasHandlers():  # hasHandlers will only be True if someone already called CALVIN with the same log_name in the same session
+      self.log.setLevel("DEBUG")
+      screen_handler = logging.StreamHandler(sys.stdout)
+      screen_handler.setLevel(logging.INFO)
+      screen_formatter = logging.Formatter('%(levelname)s - %(message)s')
+      screen_handler.setFormatter(screen_formatter)
+      self.log.addHandler(screen_handler)
+
+      file_handler = logging.FileHandler("{}.log".format(log_name))
+      file_handler.setLevel(logging.DEBUG)
+      file_formatter = logging.Formatter('%(asctime)s - %(name)s - %(levelname)s - %(message)s')
+      file_handler.setFormatter(file_formatter)
+      self.log.addHandler(file_handler)
+
     df = pd.read_csv(linksfile)
     df['link'] = df.i.map(str) + '_' + df.j.map(str) + '_' + df.k.map(str)
     df.set_index('link', inplace=True)
@@ -111,19 +141,19 @@ def networkcheck(self):
       lb_out[l[0]] += lb
       ub_out[l[0]] += ub
 
-      # if lb > ub:
-      #   raise ValueError('lb > ub for link %s' % (l[0]+'-'+l[1]))
+      if lb > ub:
+        raise ValueError('lb > ub for link %s' % (l[0]+'-'+l[1]))
 
-    # for n in nodes:
-    #   if num_in[n] == 0 and n not in ['SOURCE','SINK']:
-    #     raise ValueError('no incoming link for ' + n)
-    #   if num_out[n] == 0 and n not in ['SOURCE','SINK']:
-    #     raise ValueError('no outgoing link for ' + n)
+    for n in nodes:
+      if num_in[n] == 0 and n not in ['SOURCE','SINK']:
+        raise ValueError('no incoming link for ' + n)
+      if num_out[n] == 0 and n not in ['SOURCE','SINK']:
+        raise ValueError('no outgoing link for ' + n)
 
-    #   if ub_in[n] < lb_out[n]:
-    #     raise ValueError('ub_in < lb_out for %s (%d < %d)' % (n, ub_in[n], lb_out[n]))
-    #   if lb_in[n] > ub_out[n]:
-    #     raise ValueError('lb_in > ub_out for %s (%d > %d)' % (n, lb_in[n], ub_out[n]))
+      if ub_in[n] < lb_out[n]:
+        raise ValueError('ub_in < lb_out for %s (%d < %d)' % (n, ub_in[n], lb_out[n]))
+      if lb_in[n] > ub_out[n]:
+        raise ValueError('lb_in > ub_out for %s (%d > %d)' % (n, lb_in[n], ub_out[n]))
 
   def add_ag_region_sinks(self):
     """
@@ -143,7 +173,7 @@ def add_ag_region_sinks(self):
       links.upper_bound = maxub
       links['link'] = links.i.map(str) + '_' + links.j.map(str) + '_' + links.k.map(str)
       links.set_index('link', inplace=True)
-      self.df = self.df.append(links.drop_duplicates())
+      self.df = self.df._append(links.drop_duplicates())
 
 
   def fix_hydropower_lbs(self):
@@ -199,16 +229,16 @@ def create_pyomo_model(self, debug_mode=False, debug_cost=2e7):
     else:
       df = self.df
 
-    print('Creating Pyomo Model (debug=%s)' % debug_mode)
+    self.log.info('Creating Pyomo Model (debug=%s)' % debug_mode)
 
     model = ConcreteModel()
 
     model.N = Set(initialize=self.nodes)
     model.k = Set(initialize=range(15))
     model.A = Set(within=model.N*model.N*model.k, 
                   initialize=self.links, ordered=True)
-    model.source = Param(initialize='SOURCE')
-    model.sink = Param(initialize='SINK')
+    model.source = Param(initialize='SOURCE', within=Any)
+    model.sink = Param(initialize='SINK', within=Any)
 
     def init_params(p):
       if p == 'cost' and debug_mode:
@@ -289,7 +319,7 @@ def solve_pyomo_model(self, solver='glpk', nproc=1, debug_mode=False, maxiter=10
     from pyomo.opt import SolverFactory
     opt = SolverFactory(solver)
 
-    if nproc > 1 and solver is not 'glpk':
+    if nproc > 1 and solver != 'glpk':
       opt.options['threads'] = nproc
 
     if debug_mode:
@@ -298,25 +328,25 @@ def solve_pyomo_model(self, solver='glpk', nproc=1, debug_mode=False, maxiter=10
       vol_total = 0
 
       while run_again and i < maxiter:
-        print('-----Solving Pyomo Model (debug=%s)' % debug_mode)
+        self.log.info('-----Solving Pyomo Model (debug=%s)' % debug_mode)
         self.results = opt.solve(self.model)
-        print('Finished. Fixing debug flows...')
+        self.log.info('Finished. Fixing debug flows...')
         run_again,vol = self.fix_debug_flows()
         i += 1
         vol_total += vol
 
       if run_again:
-        print(('Warning: Debug mode maximum iterations reached.'
+        self.log.info(('Warning: Debug mode maximum iterations reached.'
                ' Will still try to solve without debug mode.'))
       else:
-        print('All debug flows eliminated (iter=%d, vol=%0.2f)' % (i,vol_total))
+        self.log.info('All debug flows eliminated (iter=%d, vol=%0.2f)' % (i,vol_total))
 
     else:
-      print('-----Solving Pyomo Model (debug=%s)' % debug_mode)
-      self.results = opt.solve(self.model, tee=True)
+      self.log.info('-----Solving Pyomo Model (debug=%s)' % debug_mode)
+      self.results = opt.solve(self.model, tee=False)
 
       if self.results.solver.termination_condition == TerminationCondition.optimal:
-        print('Optimal Solution Found (debug=%s).' % debug_mode)
+        self.log.info('Optimal Solution Found (debug=%s).' % debug_mode)
         self.model.solutions.load_from(self.results)
       else:
         raise RuntimeError('Problem Infeasible. Run again starting from debug mode.')
@@ -357,7 +387,7 @@ def fix_debug_flows(self, tol=1e-7):
             v = model.X[s].value*1.2
             model.u[s2].value += v
             vol_total += v
-            print('%s UB raised by %0.2f (%0.2f%%)' % (l[0]+'_'+l[1], v, v*100/iv))
+            self.log.info('%s UB raised by %0.2f (%0.2f%%)' % (l[0]+'_'+l[1], v, v*100/iv))
             df.loc['_'.join(str(x) for x in l[0:3]), 'upper_bound'] = model.u[s2].value
 
         # if we need to bring in extra water
@@ -367,7 +397,7 @@ def fix_debug_flows(self, tol=1e-7):
 
         if 'DBUGSRC' in dbl[0]:
           vol_to_reduce = max(model.X[s].value*1.2, 0.5)
-          print('Volume to reduce: %.2e' % vol_to_reduce)
+          self.log.info('Volume to reduce: %.2e' % vol_to_reduce)
 
           children = [dbl[1]]
           for i in range(max_depth):
@@ -397,14 +427,14 @@ def fix_debug_flows(self, tol=1e-7):
               model.l[s2].value -= v
               vol_to_reduce -= v
               vol_total += v
-              print('%s LB reduced by %.2e (%0.2f%%). Dual=%.2e' % (l[0]+'_'+l[1], v, v*100/iv, dl))
+              self.log.info('%s LB reduced by %.2e (%0.2f%%). Dual=%.2e' % (l[0]+'_'+l[1], v, v*100/iv, dl))
               df.loc['_'.join(str(x) for x in l[0:3]), 'lower_bound'] = model.l[s2].value
 
             if vol_to_reduce == 0:
               break
 
           if vol_to_reduce > 0:
-            print('Debug -> %s: could not reduce full amount (%.2e left)' % (dbl[1],vol_to_reduce))
+            self.log.info('Debug -> %s: could not reduce full amount (%.2e left)' % (dbl[1],vol_to_reduce))
 
     self.df, self.model = df, model
-    return run_again,vol_total
+    return run_again, vol_total
diff --git a/concrete_model/calvin/calvin.py b/concrete_model/calvin/calvin.py
@@ -1,20 +1,50 @@
+import os
+import sys
+import logging
+
 from pyomo.environ import *
 from pyomo.opt import TerminationCondition
 import numpy as np
 import pandas as pd
-import os
 
 class CALVIN():
 
-  def __init__(self, linksfile, ic=None):
+  def __init__(self, linksfile, ic=None, log_name="calvin"):
     """
     Initialize CALVIN model object.
 
     :param linksfile: (string) CSV file containing network link information
     :param ic: (dict) Initial storage conditions for surface reservoirs
                 only used for annual optimization
+    :param log_name: A name for a logger - will be used to keep logs from different model runs separate in files.
+                Defaults to "calvin", which results in a log file in the current working directory named "calvin.log".
+                You can change this each time you instantiate the CALVIN class if you want to output separate logs
+                for different runs. Otherwise, all results will be appended to the log file (not overwritten). If you
+                run multiple copies of CALVIN simultaneously, make sure to change this, or you could get errors writing
+                to the log file.
+
+                Do not provide a full path to a log file here because this value is also used in a way that is *not* a
+                file path. If being able to specify a full path is important for your workflow, please raise a GitHub
+                issue. It could be supported, but there is no need at this moment.
     :returns: CALVIN model object
     """
+
+    # set up logging code
+    self.log = logging.getLogger(log_name)
+    if not self.log.hasHandlers():  # hasHandlers will only be True if someone already called CALVIN with the same log_name in the same session
+      self.log.setLevel("DEBUG")
+      screen_handler = logging.StreamHandler(sys.stdout)
+      screen_handler.setLevel(logging.INFO)
+      screen_formatter = logging.Formatter('%(levelname)s - %(message)s')
+      screen_handler.setFormatter(screen_formatter)
+      self.log.addHandler(screen_handler)
+
+      file_handler = logging.FileHandler("{}.log".format(log_name))
+      file_handler.setLevel(logging.DEBUG)
+      file_formatter = logging.Formatter('%(asctime)s - %(name)s - %(levelname)s - %(message)s')
+      file_handler.setFormatter(file_formatter)
+      self.log.addHandler(file_handler)
+
     df = pd.read_csv(linksfile)
     df['link'] = df.i.map(str) + '_' + df.j.map(str) + '_' + df.k.map(str)
     df.set_index('link', inplace=True)
@@ -111,19 +141,19 @@ def networkcheck(self):
       lb_out[l[0]] += lb
       ub_out[l[0]] += ub
 
-      # if lb > ub:
-      #   raise ValueError('lb > ub for link %s' % (l[0]+'-'+l[1]))
+      if lb > ub:
+        raise ValueError('lb > ub for link %s' % (l[0]+'-'+l[1]))
 
-    # for n in nodes:
-    #   if num_in[n] == 0 and n not in ['SOURCE','SINK']:
-    #     raise ValueError('no incoming link for ' + n)
-    #   if num_out[n] == 0 and n not in ['SOURCE','SINK']:
-    #     raise ValueError('no outgoing link for ' + n)
+    for n in nodes:
+      if num_in[n] == 0 and n not in ['SOURCE','SINK']:
+        raise ValueError('no incoming link for ' + n)
+      if num_out[n] == 0 and n not in ['SOURCE','SINK']:
+        raise ValueError('no outgoing link for ' + n)
 
-    #   if ub_in[n] < lb_out[n]:
-    #     raise ValueError('ub_in < lb_out for %s (%d < %d)' % (n, ub_in[n], lb_out[n]))
-    #   if lb_in[n] > ub_out[n]:
-    #     raise ValueError('lb_in > ub_out for %s (%d > %d)' % (n, lb_in[n], ub_out[n]))
+      if ub_in[n] < lb_out[n]:
+        raise ValueError('ub_in < lb_out for %s (%d < %d)' % (n, ub_in[n], lb_out[n]))
+      if lb_in[n] > ub_out[n]:
+        raise ValueError('lb_in > ub_out for %s (%d > %d)' % (n, lb_in[n], ub_out[n]))
 
   def add_ag_region_sinks(self):
     """
@@ -143,7 +173,7 @@ def add_ag_region_sinks(self):
       links.upper_bound = maxub
       links['link'] = links.i.map(str) + '_' + links.j.map(str) + '_' + links.k.map(str)
       links.set_index('link', inplace=True)
-      self.df = self.df.append(links.drop_duplicates())
+      self.df = self.df._append(links.drop_duplicates())
 
 
   def fix_hydropower_lbs(self):
@@ -199,16 +229,16 @@ def create_pyomo_model(self, debug_mode=False, debug_cost=2e7):
     else:
       df = self.df
 
-    print('Creating Pyomo Model (debug=%s)' % debug_mode)
+    self.log.info('Creating Pyomo Model (debug=%s)' % debug_mode)
 
     model = ConcreteModel()
 
     model.N = Set(initialize=self.nodes)
     model.k = Set(initialize=range(15))
     model.A = Set(within=model.N*model.N*model.k, 
                   initialize=self.links, ordered=True)
-    model.source = Param(initialize='SOURCE',within=Any)
-    model.sink = Param(initialize='SINK',within=Any)
+    model.source = Param(initialize='SOURCE', within=Any)
+    model.sink = Param(initialize='SINK', within=Any)
 
     def init_params(p):
       if p == 'cost' and debug_mode:
@@ -298,25 +328,25 @@ def solve_pyomo_model(self, solver='glpk', nproc=1, debug_mode=False, maxiter=10
       vol_total = 0
 
       while run_again and i < maxiter:
-        print('-----Solving Pyomo Model (debug=%s)' % debug_mode)
+        self.log.info('-----Solving Pyomo Model (debug=%s)' % debug_mode)
         self.results = opt.solve(self.model)
-        print('Finished. Fixing debug flows...')
+        self.log.info('Finished. Fixing debug flows...')
         run_again,vol = self.fix_debug_flows()
         i += 1
         vol_total += vol
 
       if run_again:
-        print(('Warning: Debug mode maximum iterations reached.'
+        self.log.info(('Warning: Debug mode maximum iterations reached.'
                ' Will still try to solve without debug mode.'))
       else:
-        print('All debug flows eliminated (iter=%d, vol=%0.2f)' % (i,vol_total))
+        self.log.info('All debug flows eliminated (iter=%d, vol=%0.2f)' % (i,vol_total))
 
     else:
-      print('-----Solving Pyomo Model (debug=%s)' % debug_mode)
-      self.results = opt.solve(self.model, tee=True)
+      self.log.info('-----Solving Pyomo Model (debug=%s)' % debug_mode)
+      self.results = opt.solve(self.model, tee=False)
 
       if self.results.solver.termination_condition == TerminationCondition.optimal:
-        print('Optimal Solution Found (debug=%s).' % debug_mode)
+        self.log.info('Optimal Solution Found (debug=%s).' % debug_mode)
         self.model.solutions.load_from(self.results)
       else:
         raise RuntimeError('Problem Infeasible. Run again starting from debug mode.')
@@ -357,7 +387,7 @@ def fix_debug_flows(self, tol=1e-7):
             v = model.X[s].value*1.2
             model.u[s2].value += v
             vol_total += v
-            print('%s UB raised by %0.2f (%0.2f%%)' % (l[0]+'_'+l[1], v, v*100/iv))
+            self.log.info('%s UB raised by %0.2f (%0.2f%%)' % (l[0]+'_'+l[1], v, v*100/iv))
             df.loc['_'.join(str(x) for x in l[0:3]), 'upper_bound'] = model.u[s2].value
 
         # if we need to bring in extra water
@@ -367,7 +397,7 @@ def fix_debug_flows(self, tol=1e-7):
 
         if 'DBUGSRC' in dbl[0]:
           vol_to_reduce = max(model.X[s].value*1.2, 0.5)
-          print('Volume to reduce: %.2e' % vol_to_reduce)
+          self.log.info('Volume to reduce: %.2e' % vol_to_reduce)
 
           children = [dbl[1]]
           for i in range(max_depth):
@@ -397,14 +427,14 @@ def fix_debug_flows(self, tol=1e-7):
               model.l[s2].value -= v
               vol_to_reduce -= v
               vol_total += v
-              print('%s LB reduced by %.2e (%0.2f%%). Dual=%.2e' % (l[0]+'_'+l[1], v, v*100/iv, dl))
+              self.log.info('%s LB reduced by %.2e (%0.2f%%). Dual=%.2e' % (l[0]+'_'+l[1], v, v*100/iv, dl))
               df.loc['_'.join(str(x) for x in l[0:3]), 'lower_bound'] = model.l[s2].value
 
             if vol_to_reduce == 0:
               break
 
           if vol_to_reduce > 0:
-            print('Debug -> %s: could not reduce full amount (%.2e left)' % (dbl[1],vol_to_reduce))
+            self.log.info('Debug -> %s: could not reduce full amount (%.2e left)' % (dbl[1],vol_to_reduce))
 
     self.df, self.model = df, model
-    return run_again,vol_total
+    return run_again, vol_total