1#  ___________________________________________________________________________
2#
3#  Pyomo: Python Optimization Modeling Objects
4#  Copyright 2017 National Technology and Engineering Solutions of Sandia, LLC
5#  Under the terms of Contract DE-NA0003525 with National Technology and
6#  Engineering Solutions of Sandia, LLC, the U.S. Government retains certain
7#  rights in this software.
8#  This software is distributed under the 3-clause BSD License.
9#  ___________________________________________________________________________
10
11import os.path
12from pyomo.common.fileutils import this_file_dir, import_file
13import pyomo.common.unittest as unittest
14import pyomo.environ as pyo
15from pyomo.common.dependencies import attempt_import
16from pyomo.common.log import LoggingIntercept
17from pyomo.opt import TerminationCondition
18from io import StringIO
19import logging
20
21from pyomo.contrib.pynumero.dependencies import (
22    numpy as np, numpy_available, scipy, scipy_available
23)
24from pyomo.common.dependencies.scipy import sparse as spa
25
26if not (numpy_available and scipy_available):
27    raise unittest.SkipTest("Pynumero needs scipy and numpy to run CyIpopt tests")
28
29pandas, pandas_available = attempt_import(
30    'pandas',
31    'One of the tests below requires a recent version of pandas for'
32    ' comparing with a tolerance.',
33    minimum_version='1.1.0',
34    defer_check=False)
35
36from pyomo.contrib.pynumero.asl import AmplInterface
37if not AmplInterface.available():
38    raise unittest.SkipTest(
39        "Pynumero needs the ASL extension to run CyIpopt tests")
40
41import pyomo.contrib.pynumero.algorithms.solvers.cyipopt_solver as cyipopt_solver
42if not cyipopt_solver.cyipopt_available:
43    raise unittest.SkipTest(
44        "PyNumero needs CyIpopt installed to run CyIpopt tests")
45import cyipopt as cyipopt_core
46
47example_dir = os.path.join(this_file_dir(), '..')
48
49class TestPyomoCyIpoptSolver(unittest.TestCase):
50    def test_status_maps(self):
51        # verify that all status messages from cyipopy can be cleanly
52        # mapped back to a Pyomo TerminationCondition
53        for msg in cyipopt_core.STATUS_MESSAGES.values():
54            self.assertIn(msg, cyipopt_solver._cyipopt_status_enum)
55        for status in cyipopt_solver._cyipopt_status_enum.values():
56            self.assertIn(status, cyipopt_solver._ipopt_term_cond)
57
58
59class TestExamples(unittest.TestCase):
60    def test_external_grey_box_react_example_maximize_cb_outputs(self):
61        ex = import_file(os.path.join(example_dir, 'external_grey_box', 'react_example', 'maximize_cb_outputs.py'))
62        m = ex.maximize_cb_outputs()
63        self.assertAlmostEqual(pyo.value(m.reactor.inputs['sv']), 1.34381, places=3)
64        self.assertAlmostEqual(pyo.value(m.reactor.outputs['cb']), 1072.4372, places=2)
65
66    def test_external_grey_box_react_example_maximize_cb_outputs_scaling(self):
67        ex = import_file(os.path.join(example_dir, 'external_grey_box', 'react_example', 'maximize_cb_ratio_residuals.py'))
68        aoptions={'nlp_scaling_method': 'user-scaling',
69                 'output_file': '_cyipopt-external-greybox-react-scaling.log',
70                 'file_print_level':10}
71        m = ex.maximize_cb_ratio_residuals_with_output_scaling(additional_options=aoptions)
72        self.assertAlmostEqual(pyo.value(m.reactor.inputs['sv']), 1.26541996, places=3)
73        self.assertAlmostEqual(pyo.value(m.reactor.inputs['cb']), 1071.7410089, places=2)
74        self.assertAlmostEqual(pyo.value(m.reactor.outputs['cb_ratio']), 0.15190409266, places=3)
75
76        with open('_cyipopt-external-greybox-react-scaling.log', 'r') as fd:
77            solver_trace = fd.read()
78        os.remove('_cyipopt-external-greybox-react-scaling.log')
79
80        self.assertIn('nlp_scaling_method = user-scaling', solver_trace)
81        self.assertIn('output_file = _cyipopt-external-greybox-react-scaling.log', solver_trace)
82        self.assertIn('objective scaling factor = 1', solver_trace)
83        self.assertIn('x scaling provided', solver_trace)
84        self.assertIn('c scaling provided', solver_trace)
85        self.assertIn('d scaling provided', solver_trace)
86        self.assertIn('DenseVector "x scaling vector" with 7 elements:', solver_trace)
87        self.assertIn('x scaling vector[    2]= 1.2000000000000000e+00', solver_trace)
88        self.assertIn('x scaling vector[    7]= 1.7000000000000000e+00', solver_trace)
89        self.assertIn('x scaling vector[    6]= 1.1000000000000001e+00', solver_trace)
90        self.assertIn('x scaling vector[    1]= 1.3000000000000000e+00', solver_trace)
91        self.assertIn('x scaling vector[    3]= 1.3999999999999999e+00', solver_trace)
92        self.assertIn('x scaling vector[    4]= 1.5000000000000000e+00', solver_trace)
93        self.assertIn('x scaling vector[    5]= 1.6000000000000001e+00', solver_trace)
94        self.assertIn('DenseVector "c scaling vector" with 6 elements:', solver_trace)
95        self.assertIn('c scaling vector[    1]= 4.2000000000000000e+01', solver_trace)
96        self.assertIn('c scaling vector[    2]= 1.0000000000000001e-01', solver_trace)
97        self.assertIn('c scaling vector[    3]= 2.0000000000000001e-01', solver_trace)
98        self.assertIn('c scaling vector[    4]= 2.9999999999999999e-01', solver_trace)
99        self.assertIn('c scaling vector[    5]= 4.0000000000000002e-01', solver_trace)
100        self.assertIn('c scaling vector[    6]= 1.0000000000000000e+01', solver_trace)
101
102    def test_external_grey_box_react_example_maximize_with_output(self):
103        ex = import_file(os.path.join(example_dir, 'external_grey_box', 'react_example', 'maximize_cb_ratio_residuals.py'))
104        m = ex.maximize_cb_ratio_residuals_with_output()
105        self.assertAlmostEqual(pyo.value(m.reactor.inputs['sv']), 1.26541996, places=3)
106        self.assertAlmostEqual(pyo.value(m.reactor.inputs['cb']), 1071.7410089, places=2)
107        self.assertAlmostEqual(pyo.value(m.reactor.outputs['cb_ratio']), 0.15190409266, places=3)
108
109    def test_external_grey_box_react_example_maximize_with_hessian_with_output(self):
110        ex = import_file(os.path.join(example_dir, 'external_grey_box', 'react_example', 'maximize_cb_ratio_residuals.py'))
111        m = ex.maximize_cb_ratio_residuals_with_hessian_with_output()
112        self.assertAlmostEqual(pyo.value(m.reactor.inputs['sv']), 1.26541996, places=3)
113        self.assertAlmostEqual(pyo.value(m.reactor.inputs['cb']), 1071.7410089, places=2)
114        self.assertAlmostEqual(pyo.value(m.reactor.outputs['cb_ratio']), 0.15190409266, places=3)
115
116    def test_external_grey_box_react_example_maximize_with_hessian_with_output_pyomo(self):
117        ex = import_file(os.path.join(example_dir, 'external_grey_box', 'react_example', 'maximize_cb_ratio_residuals.py'))
118        m = ex.maximize_cb_ratio_residuals_with_hessian_with_output_pyomo()
119        self.assertAlmostEqual(pyo.value(m.sv), 1.26541996, places=3)
120        self.assertAlmostEqual(pyo.value(m.cb), 1071.7410089, places=2)
121        self.assertAlmostEqual(pyo.value(m.cb_ratio), 0.15190409266, places=3)
122
123    def test_pyomo_react_example_maximize_with_obj(self):
124        ex = import_file(os.path.join(example_dir, 'external_grey_box', 'react_example', 'maximize_cb_ratio_residuals.py'))
125        m = ex.maximize_cb_ratio_residuals_with_obj()
126        self.assertAlmostEqual(pyo.value(m.reactor.inputs['sv']), 1.26541996, places=3)
127        self.assertAlmostEqual(pyo.value(m.reactor.inputs['cb']), 1071.7410089, places=2)
128        self.assertAlmostEqual(pyo.value(m.obj), 0.15190409266, places=3)
129
130    def test_external_grey_box_react_example_maximize_with_additional_pyomo_variables(self):
131        ex = import_file(os.path.join(example_dir, 'external_grey_box', 'react_example', 'maximize_cb_ratio_residuals.py'))
132        m = ex.maximize_cb_ratio_residuals_with_pyomo_variables()
133        self.assertAlmostEqual(pyo.value(m.reactor.inputs['sv']), 1.26541996, places=3)
134        self.assertAlmostEqual(pyo.value(m.reactor.inputs['cb']), 1071.7410089, places=2)
135        self.assertAlmostEqual(pyo.value(m.cb_ratio), 0.15190409266, places=3)
136
137    @unittest.skipIf(not pandas_available, "Test uses pandas for data")
138    def test_parameter_estimation(self):
139        data_fname = os.path.join(example_dir, 'external_grey_box', 'param_est', 'smalldata.csv')
140        baseline = pandas.read_csv(data_fname)
141
142        # test the data generator
143        ex = import_file(os.path.join(example_dir, 'external_grey_box', 'param_est', 'generate_data.py'))
144        df1 = ex.generate_data(5, 200, 5, 42)
145        df2 = ex.generate_data_external(5, 200, 5, 42)
146        pandas.testing.assert_frame_equal(df1, baseline, atol=1e-3)
147        pandas.testing.assert_frame_equal(df2, baseline, atol=1e-3)
148
149        # test the estimation
150        ex = import_file(os.path.join(example_dir, 'external_grey_box', 'param_est', 'perform_estimation.py'))
151
152        m = ex.perform_estimation_external(data_fname, solver_trace=False)
153        self.assertAlmostEqual(pyo.value(m.UA), 204.43761, places=3)
154
155        m = ex.perform_estimation_pyomo_only(data_fname, solver_trace=False)
156        self.assertAlmostEqual(pyo.value(m.UA), 204.43761, places=3)
157
158    def test_cyipopt_callbacks(self):
159        ex = import_file(os.path.join(example_dir, 'callback', 'cyipopt_callback.py'))
160
161        output = StringIO()
162        with LoggingIntercept(output, 'pyomo', logging.INFO):
163            ex.main()
164
165        self.assertIn("Residuals for iteration 2",
166                          output.getvalue().strip())
167
168    @unittest.skipIf(not pandas_available, "pandas needed to run this example")
169    def test_cyipopt_functor(self):
170        ex = import_file(os.path.join(example_dir, 'callback', 'cyipopt_functor_callback.py'))
171        df = ex.main()
172        self.assertEqual(df.shape, (7, 5))
173        # check one of the residuals
174        s = df['ca_bal']
175        self.assertAlmostEqual(s.iloc[6], 0, places=3)
176
177    def test_cyipopt_callback_halt(self):
178        ex = import_file(os.path.join(example_dir, 'callback', 'cyipopt_callback_halt.py'))
179        status = ex.main()
180        self.assertEqual(status.solver.termination_condition, TerminationCondition.userInterrupt)
181
182