Alfoa/rom to pyomo

dependencies.xml

@@ -64,8 +64,14 @@ Note all install methods after "main" take
     <line_profiler optional='True'/>
     <!-- <ete3 optional='True'/> -->
     <pywavelets optional='True'>1.1</pywavelets>
+    <numdifftools source="pip">0.9.39</numdifftools>
     <fmpy optional='True'/>
     <xmlschema source="pip"/>
+    <pyomo optional='True'>6.4</pyomo>
+    <glpk skip_check='True' optional='True'/>
+    <ipopt skip_check='True' optional='True'/>
+    <cyipopt optional='True'/>
+    <pyomo-extensions source="pyomo" skip_check='True' optional='True'/>
     <setuptools/>
   </main>
   <alternate name="pip">

doc/user_manual/rom.tex

@@ -191,3 +191,121 @@ \subsection{ROM}
 %%%%% ROM Model - TensorFlow-Keras Interface  %%%%%%%
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 \input{kerasROM.tex}
+
+\subsubsection{SerializePyomo}
+\label{subsubsec:serializepyomo}
+For the purpose of exporting RAVEN trained ROMs in an universal-accessible format for Pyomo to solve, the platform to prepare
+the python syntax including the process to setup concrete models and constraints is created. The python syntax allows to
+to retrieve pickled (serialized) ROMs and solve the defined concrete model via Pynumero's GreyBoxModel, an extension of Pyomo.
+Details on how GreyBoxModel works can be found in the following site:
+
+\begin{lstlisting}
+  https://pyomo.readthedocs.io/en/stable/contributed_packages/pynumero/index.html
+\end{lstlisting}
+
+Two files are required to run the optimization: the pickled rom and the generated python file. The workflow to create both files is
+provided in the example below located at the end of this subsection. The following command line is the defined format to interface both files:
+
+\begin{lstlisting}
+  python3 printed_python_file.py -r rom_pickled_file.pk -f location_of_raven_framework -order order_of_derivative
+\end{lstlisting}
+
+The file names `printed\_python\_file', `rom\_pickled\_file', and `order\_of\_derivative' are arbitrary names, and are not requirements.
+The `location\_of\_raven\_framework' is the relative path to where `raven\_framework' is located. The created python file includes modules from
+RAVEN and is required dealing with derivatives for GreyBoxModel to solve. The rest of the command line syntax `-r', `-f', `-order' are place holders for
+the python file to locate where the pickled ROMs and the python file itself are. Caution, make sure the `order\_of\_derivative' is an integer. If not provided,
+the default value is 1. When running the command a text file named `GreyModelOutput\_cyipopt.txt' will be created showing the results of the optimization.
+
+Since the code developed in ROMs is not a subtype, the code activation for python interface file print is triggered by adding the input subnode with type `Pyomo':
+
+\begin{lstlisting}[style=XML,morekeywords={name,subType}]
+  <Input name="rom_out.py" type="Pyomo">rom_out.py</Input>
+\end{lstlisting}
+
+To verify the output files in the given test file in:
+
+\begin{lstlisting}[style=XML,morekeywords={name,subType}]
+  raven/tests/framework/Models/External/serialize_pyomo.xml
+\end{lstlisting}
+
+and the associated example optimization case:
+
+\begin{lstlisting}[style=XML,morekeywords={name,subType}]
+  python3 rom_out.py -r rom_pickle.pk -f ../../../..
+\end{lstlisting}
+
+are runnable with GreyBoxModel, several optional RAVEN dependencies are required. The following are the required
+optional libraries: `pyomo 6.4' or over, `cmake', `glpk', `ipopt', `cyipopt', and `pyomo-extensions'. Although almost all required
+libraries are standard packages (available on conda or pipy install), the `pyomo-extensions' is a custom library created by RAVEN developers to install `Pynumero'.
+It functions to run the following commands if `pyomo-extensions' is part of the RAVEN library installation:
+
+\begin{lstlisting}
+  pyomo download-extensions
+  pyomo build-extensions
+\end{lstlisting}
+
+For this reason, libraries `Pyomo' and `cmake' are required for successful `Pynumero' installation. Note, although all RAVEN dependency installs will
+be placed in the RAVEN library directory, the current `Pyomo' version saves extensions on the local python local directory.
+Other libraries are optimization solvers for `Pyomo' and `Pynumero'. Make sure the RAVEN conda library is activated when testing the example cases.
+The method to include optional RAVEN libraries is by adding the `optional' flag to the RAVEN library installation command:
+
+\begin{lstlisting}
+  cd raven
+  ./scripts/establish_conda_env.sh --install --optional='pyomo cmake ipopt cyipopt pyomo-extensions'
+\end{lstlisting}
+
+If the user plans to run the optimization case outside of RAVEN libraries and the pickled ROM has already been generated, instead install `pyomo', `cmake',
+`glpk', `ipopt', and `cyipopt' to your system. After this is done, on a terminal, run the `Pyomo' download and extension command introduced earlier.
+This will install `Pynumero' to your system.
+
+If the user plans to use the RAVEN library, it is recommended to activate RAVEN libraries using the following command
+on the terminal:
+
+\begin{lstlisting}
+  source scripts/establish_conda_env.sh --load
+\end{lstlisting}
+
+Make sure that the optional libraries are installed before running the pickled ROM and printed python file. Regardless of whether `Pynumero' is
+installed inside or outside RAVEN libraries, depending on the OS, C++ libraries may be outdated for `Pynumero' to run. It has been reported that
+the following OS has failed to run `Pynumero':
+
+\begin{itemize}
+  \item CentOS 7, 8
+  \item Ubuntu 16, 18
+\end{itemize}
+
+\textbf{Example:} For this example the external model is trained and further loaded as `out\_rom'. The loaded trained model is separately
+processed to `rom\_out.py' and `rom\_pickle.pk'.
+\begin{lstlisting}[style=XML,morekeywords={name,subType}]
+  <Simulation>
+    ...
+    <Files>
+      <Input name="rom_out.py" type="Pyomo">rom_out.py</Input>
+      <Input name="rom_pickle.pk" type="">rom_pickle.pk</Input>
+    </Files>
+    ...
+    <Steps>
+      ...
+      <MultiRun name="sample">
+        ...
+        <Model class="Models" type="ExternalModel">attenuate</Model>
+        ...
+        <Output class="DataObjects" type="PointSet">samples</Output>
+      </MultiRun>
+      <RomTrainer name="train">
+        <Input class="DataObjects" type="PointSet">samples</Input>
+        <Output class="Models" type="ROM">out_rom</Output>
+      </RomTrainer>
+      <IOStep name="serialize">
+        <Input class="Models" type="ROM">out_rom</Input>
+        <Output class="Files" type="">out_rom.py</Output>
+      </IOStep>
+      <IOStep name="pickle">
+        <Input class="Models" type="ROM">out_rom</Input>
+        <Output class="Files" type="">out_pickle.pk</Output>
+      </IOStep>
+      ...
+    </Steps>
+    ...
+  </Simulation>
+\end{lstlisting}

ravenframework/BaseClasses/BaseEntity.py

@@ -246,7 +246,6 @@ def _validateSolutionExportVariables(self, solutionExport):
     # get acceptable names
     fromSolnExport = set(self.getSolutionExportVariableNames())
     acceptable = set(self._formatSolutionExportVariableNames(fromSolnExport))
-
     # remove registered solution export names first
     remaining = requested - acceptable
     # anything remaining is unknown!

ravenframework/Models/Model.py

@@ -13,13 +13,17 @@
 # limitations under the License.
 """
 Module where the base class and the specialization of different type of Model are
+
+@author crisrab, alfoa
+
 """
 #External Modules------------------------------------------------------------------------------------
 import copy
 import numpy as np
 import abc
 import sys
 import importlib
+import pickle
 #External Modules End--------------------------------------------------------------------------------
 
 #Internal Modules------------------------------------------------------------------------------------
@@ -43,7 +47,6 @@ def loadFromPlugins(cls):
     """
     cls.plugins = importlib.import_module(".ModelPlugInFactory","ravenframework.Models")
 
-
   @classmethod
   def getInputSpecification(cls):
     """
@@ -381,6 +384,23 @@ def initialize(self,runInfo,inputs,initDict=None):
     """
     pass
 
+  def serialize(self,fileObjIn,**kwargs):
+    """
+      This method is the base class method that is aimed to serialize the model (and derived) instances.
+      @ In, fileObjIn, str or File object, the filename of the output serialized binary file or the RAVEN File instance
+      @ In, kwargs, dict, dictionary of options that the derived class might require
+      @ Out, None
+    """
+    import cloudpickle
+    if isinstance(fileObjIn,str):
+      fileObj = open(filename, mode='wb+')
+    else:
+      fileObj = fileObjIn # if issues occur add 'isintance(fileObjIn,Files)'.
+      fileObj.open(mode='wb+')
+    cloudpickle.dump(self,fileObj, protocol=pickle.HIGHEST_PROTOCOL)
+    fileObj.flush()
+    fileObj.close()
+
   @abc.abstractmethod
   def createNewInput(self,myInput,samplerType,**kwargs):
     """

ravenframework/Models/ROM.py

@@ -108,7 +108,6 @@ def __init__(self):
     """
     super().__init__()
     self.amITrained = False               # boolean flag, is the ROM trained?
-    self.supervisedEngine = None          # dict of ROM instances (== number of targets => keys are the targets)
     self.printTag = 'ROM MODEL'           # label
     self.cvInstanceName = None            # the name of Cross Validation instance
     self.cvInstance = None                # Instance of provided cross validation
@@ -417,6 +416,38 @@ def evaluate(self, request):
       resultsDict[k] = np.atleast_1d(v)
     return resultsDict
 
+  def derivatives(self, request, feats = None, order=1):
+    """
+      This method is aimed to evaluate the derivatives using this ROM
+      The differentiation method depends on the specific ROM. Numerical
+      differentiation is available by default for any ROM. Some ROMs (e.g.
+      neural networks) use automatic differentiation.
+      @ In, request, dict, coordinate of the central point
+      @ In, feats, list, optional, list of features we need to compute the
+                                   derivative for (default all)
+      @ In, order, int, order of the derivative (1 = first, 2 = second, etc.). Max 10
+      @ Out, derivatives, dict, the dict containing the derivatives for each target
+                                ({'d(feature1)/d(target1)':np.array(size 1 or n_ts),
+                                  'd(feature1)/d(target1)':np.array(...)}. For example,
+                                {"x1/y1":np.array(...),"x2/y1":np.array(...),etc.}
+    """
+    if self.pickled:
+      self.raiseAnError(RuntimeError,'ROM "', self.name, '" has not been loaded yet!  Use an IOStep to load it.')
+    if not self.amITrained:
+      self.raiseAnError(RuntimeError, "ROM ", self.name, " has not been trained yet and, consequentially, can not be evaluated!")
+    derivatives = {}
+    if self.segment:
+      derivatives = mathUtils.derivatives(self.supervisedContainer[0].evaluate, request, var=feats, n=order)
+    else:
+      for rom in self.supervisedContainer:
+        sliceEvaluation = mathUtils.derivatives(rom.evaluate, request, var=feats, n=order)
+        if len(list(derivatives.keys())) == 0:
+          derivatives.update(sliceEvaluation)
+        else:
+          for key in derivatives.keys():
+            derivatives[key] = np.append(derivatives[key],sliceEvaluation[key])
+    return derivatives
+
   def _externalRun(self,inRun):
     """
       Method that performs the actual run of the imported external model (separated from run method for parallelization purposes)
@@ -567,3 +598,133 @@ def writeXML(self, what='all'):
       engines[0].writeXMLPreamble(xml)
       engines[0].writeXML(xml)
     return xml
+
+  def writePyomoGreyModel(self):
+    """
+      @ In, None, called by the OutStreamPrint object to cause the ROM to print itself
+      @ Out, xml, xmlUtils.StaticXmlElement, written meta
+    """
+    template = r"""
+# MODEL GENERATED BY RAVEN (raven.inl.gov)
+import pyomo.environ as pyo
+from pyomo.contrib.pynumero.interfaces.external_grey_box import ExternalGreyBoxModel, ExternalGreyBoxBlock
+from pyomo.contrib.pynumero.dependencies import (numpy as np)
+from pyomo.contrib.pynumero.asl import AmplInterface
+from pyomo.contrib.pynumero.algorithms.solvers.cyipopt_solver import CyIpoptSolver, CyIpoptNLP
+import os, sys, pickle
+from contextlib import redirect_stdout
+# RAVEN ROM PYOMO GREY MODEL CLASS
+
+class ravenROM(ExternalGreyBoxModel):
+
+  def __init__(self, **kwargs):
+    self._rom_file = kwargs.get("rom_file")
+    self._raven_framework = kwargs.get("raven_framework")
+    self._order = kwargs.get("order")
+    self._raven_framework = os.path.abspath(self._raven_framework)
+    if not os.path.exists(self._raven_framework):
+      raise IOError('The RAVEN framework directory does not exist in location "' + str(self._raven_framework)+'" !')
+    if os.path.isdir(os.path.join(self._raven_framework,"ravenframework")):
+      sys.path.append(self._raven_framework)
+    else:
+      raise IOError('The RAVEN framework directory does not exist in location "' + str(self._raven_framework)+'" !')
+    from ravenframework.CustomDrivers import DriverUtils as dutils
+    dutils.doSetup()
+    # de-serialize the ROM
+    self._rom_file = os.path.abspath(self._rom_file)
+    if not os.path.exists(self._rom_file):
+      raise IOError('The serialized (binary) file has not been found in location "' + str(self._rom_file)+'" !')
+    self.rom = pickle.load(open(self._rom_file, mode='rb'))
+    #
+    self.settings = self.rom.getInitParams()
+    # get input names
+    self._input_names = self.settings.get('Features')
+    # n_inputs
+    self._n_inputs = len(self._input_names)
+    # get output names
+    self._output_names = self.settings.get('Target')
+    # n_outputs
+    self._n_outputs = len(self._output_names)
+    # input values storage
+    self._input_values = np.zeros(self._n_inputs, dtype=np.float64)
+
+  def return_train_values(self, feat):
+    return self.rom.trainingSet.get(feat)
+
+  def input_names(self):
+    return self._input_names
+
+  def output_names(self):
+    return self._output_names
+
+  def set_input_values(self, input_values):
+    assert len(input_values) == self._n_inputs
+    np.copyto(self._input_values, input_values)
+
+  def evaluate_equality_constraints(self):
+      raise NotImplementedError('This method should not be called for this model.')
+
+  def evaluate_outputs(self):
+    request = {k:np.asarray(v) for k,v in zip(self._input_names,self._input_values)}
+    outs = self.rom.evaluate(request)
+    eval_outputs = np.asarray([outs[k].flatten() for k in self._output_names], dtype=np.float64)
+    return eval_outputs.flatten()
+
+  def evaluate_jacobian_outputs(self):
+    request = {k:np.asarray(v) for k,v in zip(self._input_names,self._input_values)}
+    derivatives = self.rom.derivatives(request, order=self._order)
+    jac = np.zeros((self._n_outputs, self._n_inputs))
+    for tc, target in enumerate(self._output_names):
+      for fc, feature in enumerate(self._input_names):
+        jac[tc,fc] = derivatives['d{}|d{}'.format(target, feature)]
+    return jac
+
+def pyomoModel(ex_model):
+  m = pyo.ConcreteModel()
+  m.egb = ExternalGreyBoxBlock()
+  m.egb.set_external_model(ex_model)
+  for inp in ex_model.input_names():
+    m.egb.inputs[inp].value = np.mean(ex_model.return_train_values(inp))
+    m.egb.inputs[inp].setlb(np.min(ex_model.return_train_values(inp)))
+    m.egb.inputs[inp].setub(np.max(ex_model.return_train_values(inp)))
+  for out in ex_model.output_names():
+    m.egb.outputs[out].value = np.mean(ex_model.return_train_values(out))
+    m.egb.outputs[out].setlb(np.min(ex_model.return_train_values(out)))
+    m.egb.outputs[out].setub(np.max(ex_model.return_train_values(out)))
+  m.obj = pyo.Objective(expr=m.egb.outputs[out])
+
+  return m
+
+if __name__ == '__main__':
+  for cnt, item in enumerate(sys.argv):
+    if item.lower() == "-r":
+      rom_file = sys.argv[cnt+1]
+    if item.lower() == "-f":
+      raven_framework = sys.argv[cnt+1]
+    if item.lower() == "-order":
+      order = int(sys.argv[cnt+1])
+    else:
+      order = 1
+  ext_model = ravenROM(**{'rom_file':rom_file,'raven_framework':raven_framework,'order':order})
+  concreteModel = pyomoModel(ext_model)
+  ### here you should implement the optimization problem
+  ###
+  solver = pyo.SolverFactory('cyipopt')
+  solver.config.options['hessian_approximation'] = 'limited-memory'
+  results = solver.solve(concreteModel)
+  print(results)
+  count = 0
+  lines = str(results).split("\n")
+  with open ("GreyModelOutput_cyipopt.txt", "w") as textfile:
+    with redirect_stdout(textfile):
+      print("---------------------< Output Summary >---------------------" + "\n")
+      for element in lines:
+        textfile.write(element + "\n")
+      print("---------------------< Optimization Results >---------------------" + "\n")
+      concreteModel.pprint()
+"""
+
+    return template
+
+
+

ravenframework/Steps/IOStep.py

@@ -111,6 +111,8 @@ def _localInitializeStep(self,inDictionary):
       elif isinstance(inDictionary['Input'][i], (Models.ROM, Models.ExternalModel)):
         # ... file
         if isinstance(outputs[i],Files.File):
+          if 'PYOMO' == outputs[i].getType().upper():
+            self.actionType.append('MODEL-PYOMO')
           if 'FMU' == outputs[i].getType().upper():
             self.actionType.append('MODEL-FMU')
           else:
@@ -206,6 +208,12 @@ def _localTakeAstepRun(self, inDictionary):
         cloudpickle.dump(inDictionary['Input'][i], fileobj, protocol=pickle.HIGHEST_PROTOCOL)
         fileobj.flush()
         fileobj.close()
+
+      elif self.actionType[i] == 'MODEL-PYOMO':
+        outfile = open(outputs[i].getAbsFile(),"w")
+        outfile.write(inDictionary['Input'][i].writePyomoGreyModel())
+        outfile.close()
+
       elif self.actionType[i] == 'MODEL-FMU':
         #check the ROM is trained first (if ExternalModel no check it is performed)
         if isinstance(inDictionary['Input'][i],Models.ROM) and not inDictionary['Input'][i].amITrained: