Getting Data from a Case

The Case object contains all the information about a specific Case recording whether it was recorded by a problem, driver, system, or solver. Case objects have a number methods for accessing variables and their values.

Example of Getting Variable Data from Case Recording of a Driver

Here are the methods typically used when retrieving data from the recording of a Driver.

Case.get_objectives(scaled=True, use_indices=True)[source]

Get the values of the objectives, as seen by the driver, for this case.

Parameters
scaledbool

If True, then return scaled values.

use_indicesbool

If True, apply indices.

Returns
PromAbsDict

Map of variables to their values.

Case.get_constraints(scaled=True, use_indices=True)[source]

Get the values of the constraints, as seen by the driver, for this case.

Parameters
scaledbool

If True, then return scaled values.

use_indicesbool

If True, apply indices.

Returns
PromAbsDict

Map of variables to their values.

Case.get_design_vars(scaled=True, use_indices=True)[source]

Get the values of the design variables, as seen by the driver, for this case.

Parameters
scaledbool

If True, then return scaled values.

use_indicesbool

If True, apply indices.

Returns
PromAbsDict

Map of variables to their values.

Case.get_responses(scaled=True, use_indices=True)[source]

Get the values of the responses, as seen by the driver, for this case.

Parameters
scaledbool

If True, then return scaled values.

use_indicesbool

If True, apply indices.

Returns
PromAbsDict

Map of variables to their values.

The following example shows how to use these methods to easily check the variables of interest from the first driver iteration.

import openmdao.api as om
from openmdao.test_suite.components.sellar import SellarDerivatives

import numpy as np

prob = om.Problem(model=SellarDerivatives())

model = prob.model
model.add_design_var('z', lower=np.array([-10.0, 0.0]),
                          upper=np.array([10.0, 10.0]))
model.add_design_var('x', lower=0.0, upper=10.0)
model.add_objective('obj')
model.add_constraint('con1', upper=0.0)
model.add_constraint('con2', upper=0.0)

driver = prob.driver = om.ScipyOptimizeDriver(optimizer='SLSQP', tol=1e-9, disp=False)

driver.add_recorder(om.SqliteRecorder("cases.sql"))

driver.recording_options['includes'] = []
driver.recording_options['record_objectives'] = True
driver.recording_options['record_constraints'] = True
driver.recording_options['record_desvars'] = True

prob.setup()
prob.set_solver_print(0)
prob.run_driver()
prob.cleanup()

cr = om.CaseReader("cases.sql")

driver_cases = cr.list_cases('driver')
case = cr.get_case(driver_cases[0])

print(sorted(case.outputs.keys()))
driver
    rank0:ScipyOptimize_SLSQP|0
driver
    rank0:ScipyOptimize_SLSQP|1
driver
    rank0:ScipyOptimize_SLSQP|2
driver
    rank0:ScipyOptimize_SLSQP|3
driver
    rank0:ScipyOptimize_SLSQP|4
driver
    rank0:ScipyOptimize_SLSQP|5
driver
    rank0:ScipyOptimize_SLSQP|6
['con1', 'con2', 'obj', 'x', 'z']
objs = case.get_objectives()
cons = case.get_constraints()
dvs = case.get_design_vars()
rsps = case.get_responses()

# keys() will give you the promoted variable names
print((sorted(objs.keys()), sorted(cons.keys()), sorted(dvs.keys())))
(['obj'], ['con1', 'con2'], ['x', 'z'])
# alternatively, you can get the absolute names
print((sorted(objs.absolute_names()), sorted(cons.absolute_names()), sorted(dvs.absolute_names())))
(['obj_cmp.obj'], ['con_cmp1.con1', 'con_cmp2.con2'], ['x', 'z'])
# you can access variable values using either the promoted or the absolute name
print((objs['obj'], objs['obj_cmp.obj']))
(array([28.58830817]), array([28.58830817]))
print((dvs['x'], dvs['_auto_ivc.v1']))
(array([1.]), array([1.]))
print((rsps['obj'], rsps['obj_cmp.obj']))
(array([28.58830817]), array([28.58830817]))
# you can also access the variables directly from the case object
print((case['obj'], case['obj_cmp.obj']))
(array([28.58830817]), array([28.58830817]))
print((case['x'], case['_auto_ivc.v1']))
(array([1.]), array([1.]))

Note

Note that you can use either the promoted or absolute names when accessing variables.

Getting Variable Data from Case Recording of a Problem

Here are the methods typically used when retrieving data from the recording of a Problem.

Case.list_inputs(values=True, prom_name=False, units=False, shape=False, desc=False, hierarchical=True, print_arrays=False, tags=None, includes=None, excludes=None, out_stream=DEFAULT_OUT_STREAM)[source]

Return and optionally log a list of input names and other optional information.

Parameters
valuesbool, optional

When True, display/return input values. Default is True.

prom_namebool, optional

When True, display/return the promoted name of the variable. Default is False.

unitsbool, optional

When True, display/return units. Default is False.

shapebool, optional

When True, display/return the shape of the value. Default is False.

descbool, optional

When True, display/return description. Default is False.

hierarchicalbool, optional

When True, human readable output shows variables in hierarchical format.

print_arraysbool, optional

When False, in the columnar display, just display norm of any ndarrays with size > 1. The norm is surrounded by vertical bars to indicate that it is a norm. When True, also display full values of the ndarray below the row. Format is affected by the values set with numpy.set_printoptions Default is False.

tagsstr or list of strs

User defined tags that can be used to filter what gets listed. Only inputs with the given tags will be listed. Default is None, which means there will be no filtering based on tags.

includesiter of str or None

Glob patterns for pathnames to include in the check. Default is None, which includes all.

excludesiter of str or None

Glob patterns for pathnames to exclude from the check. Default is None, which excludes nothing.

out_streamfile-like object

Where to send human readable output. Default is sys.stdout. Set to None to suppress.

Returns
list

list of input names and other optional information about those inputs

Case.list_outputs(explicit=True, implicit=True, values=True, prom_name=False, residuals=False, residuals_tol=None, units=False, shape=False, bounds=False, scaling=False, desc=False, hierarchical=True, print_arrays=False, tags=None, includes=None, excludes=None, list_autoivcs=False, out_stream=DEFAULT_OUT_STREAM)[source]

Return and optionally log a list of output names and other optional information.

Parameters
explicitbool, optional

include outputs from explicit components. Default is True.

implicitbool, optional

include outputs from implicit components. Default is True.

valuesbool, optional

When True, display/return output values. Default is True.

prom_namebool, optional

When True, display/return the promoted name of the variable. Default is False.

residualsbool, optional

When True, display/return residual values. Default is False.

residuals_tolfloat, optional

If set, limits the output of list_outputs to only variables where the norm of the resids array is greater than the given ‘residuals_tol’. Default is None.

unitsbool, optional

When True, display/return units. Default is False.

shapebool, optional

When True, display/return the shape of the value. Default is False.

boundsbool, optional

When True, display/return bounds (lower and upper). Default is False.

scalingbool, optional

When True, display/return scaling (ref, ref0, and res_ref). Default is False.

descbool, optional

When True, display/return description. Default is False.

hierarchicalbool, optional

When True, human readable output shows variables in hierarchical format.

print_arraysbool, optional

When False, in the columnar display, just display norm of any ndarrays with size > 1. The norm is surrounded by vertical bars to indicate that it is a norm. When True, also display full values of the ndarray below the row. Format is affected by the values set with numpy.set_printoptions Default is False.

tagsstr or list of strs

User defined tags that can be used to filter what gets listed. Only outputs with the given tags will be listed. Default is None, which means there will be no filtering based on tags.

includesiter of str or None

Glob patterns for pathnames to include in the check. Default is None, which includes all.

excludesiter of str or None

Glob patterns for pathnames to exclude from the check. Default is None, which excludes nothing.

list_autoivcsbool

If True, include auto_ivc outputs in the listing. Defaults to False.

out_streamfile-like

Where to send human readable output. Default is sys.stdout. Set to None to suppress.

Returns
list

list of output names and other optional information about those outputs

The following example shows how to use these methods.

import openmdao.api as om
from openmdao.test_suite.components.sellar import SellarProblem

prob = SellarProblem()

recorder = om.SqliteRecorder("cases.sql")
prob.model.add_recorder(recorder)
prob.model.recording_options['record_residuals'] = True

prob.setup()

d1 = prob.model.d1
d1.nonlinear_solver = om.NonlinearBlockGS(maxiter=5)
d1.add_recorder(recorder)

prob.run_driver()
prob.cleanup()

cr = om.CaseReader("cases.sql")

system_cases = cr.list_cases('root.d1')

case = cr.get_case(system_cases[1])

# list_inputs will print a report to the screen
case_inputs = sorted(case.list_inputs())
system
    rank0:Driver|0|root._solve_nonlinear|0|d1._solve_nonlinear|0
system
    rank0:Driver|0|root._solve_nonlinear|0|NonlinearBlockGS|1|d1._solve_nonlinear|1
system
    rank0:Driver|0|root._solve_nonlinear|0|NonlinearBlockGS|2|d1._solve_nonlinear|2
system
    rank0:Driver|0|root._solve_nonlinear|0|NonlinearBlockGS|3|d1._solve_nonlinear|3
system
    rank0:Driver|0|root._solve_nonlinear|0|NonlinearBlockGS|4|d1._solve_nonlinear|4
system
    rank0:Driver|0|root._solve_nonlinear|0|NonlinearBlockGS|5|d1._solve_nonlinear|5
system
    rank0:Driver|0|root._solve_nonlinear|0|NonlinearBlockGS|6|d1._solve_nonlinear|6
system
    rank0:Driver|0|root._solve_nonlinear|0|NonlinearBlockGS|7|d1._solve_nonlinear|7
3 Input(s) in 'd1'
------------------

varname  value              
-------  -------------------
z        |5.385164807134504|
x        [1.]               
y2       [12.27257053]
print(case_inputs[0][1]['value']) # d1.x
[1.]
print(case_inputs[1][1]['value']) # d1.y2
[12.27257053]
print(case_inputs[2][1]['value']) # d1.z
[5. 2.]
case_outputs = case.list_outputs(prom_name=True)
1 Explicit Output(s) in 'd1'
----------------------------

varname  value          prom_name
-------  -------------  ---------
y1       [25.54548589]  y1       


0 Implicit Output(s) in 'd1'
----------------------------
print(case_outputs[0][1]['value']) # d1.y1
[25.54548589]

The Case.list_inputs and Case.list_outputs methods have optional arguments that let you filter based on variable names what gets listed. This is shown in these examples.

import openmdao.api as om
from openmdao.core.tests.test_expl_comp import RectangleComp

model = om.Group()
prob = om.Problem(model)
model.add_recorder(om.SqliteRecorder('cases.sql'))

model.add_subsystem('rect', RectangleComp(), promotes=['length', 'width', 'area'])

prob.setup(check=False)

prob.set_val('length', 100.)
prob.set_val('width', 60.0)

prob.run_model()

prob.cleanup()

cr = om.CaseReader("cases.sql")

cases = cr.get_cases()
case = cases[0]

# Inputs with includes
inputs = case.list_inputs(includes=['*length'], out_stream=None)
print(sorted([inp[0] for inp in inputs]))
['rect.length']
# Inputs with multiple tags
inputs = case.list_inputs(excludes=['*length'], out_stream=None)
print(sorted([inp[0] for inp in inputs]))
['rect.width']
# Outputs with includes
outputs = case.list_outputs(includes=['*area'], out_stream=None)
print(sorted([outp[0] for outp in outputs]))
['rect.area']
# Inputs with excludes
inputs = case.list_inputs(excludes=['*length'], out_stream=None)
print(sorted(['rect.width']))
['rect.width']

Finally, you can also make use of the variable tagging feature when getting values from cases. This example shows how to do that.

import openmdao.api as om

class RectangleCompWithTags(om.ExplicitComponent):
    """
    A simple Explicit Component that also has input and output with tags.
    """

    def setup(self):
        self.add_input('length', val=1., tags=["tag1", "tag2"])
        self.add_input('width', val=1., tags=["tag2"])
        self.add_output('area', val=1., tags="tag1")

    def setup_partials(self):
        self.declare_partials('*', '*')

    def compute(self, inputs, outputs):
        outputs['area'] = inputs['length'] * inputs['width']

model = om.Group()
prob = om.Problem(model)
model.add_recorder(om.SqliteRecorder('cases.sql'))

model.add_subsystem('rect', RectangleCompWithTags(), promotes=['length', 'width', 'area'])

prob.setup(check=False)

prob.set_val('length', 100.0)
prob.set_val('width', 60.0)

prob.run_model()

prob.cleanup()

cr = om.CaseReader("cases.sql")

cases = cr.get_cases()
case = cases[0]

# Inputs with tag that matches
inputs = case.list_inputs(out_stream=None, tags="tag1")
print(sorted([inp[0] for inp in inputs]))
['rect.length']
# Inputs with multiple tags
inputs = case.list_inputs(out_stream=None, tags=["tag1", "tag2"])
print(sorted([inp[0] for inp in inputs]))
['rect.length', 'rect.width']
# Outputs with tag that does match
outputs = case.list_outputs(tags="tag1")
1 Explicit Output(s) in 'model'
-------------------------------

varname  value  
-------  -------
rect
  area   [6000.]


0 Implicit Output(s) in 'model'
-------------------------------
print(sorted([outp[0] for outp in outputs]))
['rect.area']

Getting Variable Data from Case By Specifying Variable Name and Units Desired

You can also get variable values from a Case like you would from a Problem using dictionary-like access or, if you want the value in different units, using the get_val method.

Case.get_val(name, units=None, indices=None)[source]

Get an output/input variable.

Function is used if you want to specify display units.

Parameters
namestr

Promoted or relative variable name in the root system’s namespace.

unitsstr, optional

Units to convert to before upon return.

indicesint or list of ints or tuple of ints or int ndarray or Iterable or None, optional

Indices or slice to return.

Returns
float or ndarray

The requested output/input variable.

This example shows both methods of getting variable data by name.

import openmdao.api as om

model = om.Group()
model.add_recorder(om.SqliteRecorder('cases.sql'))

speed = om.ExecComp('v=x/t', x={'units': 'm'}, t={'units': 's'}, v={'units': 'm/s'})

model.add_subsystem('speed', speed, promotes=['x', 't', 'v'])

prob = om.Problem(model)
prob.setup()

prob.set_val('x', 100., units='m')
prob.set_val('t', 60., units='s')

prob.run_model()
prob.cleanup()

cr = om.CaseReader('cases.sql')
case = cr.get_case(0)

print(case['x'])
[100.]
print(case.get_val('x', units='ft'))
[328.0839895]
print(case['v'])
[1.66666667]
print(case.get_val('v', units='ft/s'))
[5.46806649]

Getting Derivative Data from a Case

A driver has the ability to record derivatives but it is not enabled by default. If you do enable this option, the recorded cases will have a value for the jacobian.

import openmdao.api as om
from openmdao.test_suite.components.sellar_feature import SellarMDA

import numpy as np

prob = om.Problem(model=SellarMDA())

model = prob.model
model.add_design_var('z', lower=np.array([-10.0, 0.0]),
                          upper=np.array([10.0, 10.0]))
model.add_design_var('x', lower=0.0, upper=10.0)
model.add_objective('obj')
model.add_constraint('con1', upper=0.0)
model.add_constraint('con2', upper=0.0)

driver = prob.driver = om.ScipyOptimizeDriver(optimizer='SLSQP', tol=1e-9, disp=False)
driver.recording_options['record_derivatives'] = True

driver.add_recorder(om.SqliteRecorder('cases.sql'))

prob.setup()
prob.set_solver_print(0)
prob.run_driver()
prob.cleanup()

cr = om.CaseReader('cases.sql')

# Get derivatives associated with the last iteration.
derivs = cr.get_case(-1).derivatives

# check that derivatives have been recorded.
print(set(derivs.keys()))
{('con1', 'x'), ('con2', 'z'), ('obj', 'x'), ('con2', 'x'), ('con1', 'z'), ('obj', 'z')}
# Get specific derivative.
print(derivs['obj', 'z'])
[[3.90585345 1.97002049]]

Problem recording can also include recording of the derivatives as this example shows.

import openmdao.api as om
from openmdao.test_suite.components.eggcrate import EggCrate

prob = om.Problem()
model = prob.model

model.add_subsystem('egg_crate', EggCrate(), promotes=['x', 'y', 'f_xy'])
model.add_design_var('x', lower=-50.0, upper=50.0)
model.add_design_var('y', lower=-50.0, upper=50.0)
model.add_objective('f_xy')
prob.driver = om.ScipyOptimizeDriver(disp=False, tol=1e-9)

prob.recording_options['record_derivatives'] = True
recorder = om.SqliteRecorder('cases.sql')
prob.add_recorder(recorder)

prob.setup()

prob.set_solver_print(0)

prob.set_val('x', 2.5)
prob.set_val('y', 2.5)

prob.run_driver()
print(prob.get_val('x'), prob.get_val('y'), prob.get_val('f_xy'))
[3.01960159] [3.01960159] [18.97639468]
case_name_1 = "c1"
prob.record(case_name_1)


prob.set_val('x', 0.1)
prob.set_val('y', -0.1)
prob.run_driver()
print(prob.get_val('x'), prob.get_val('y'), prob.get_val('f_xy'))
[-2.14311975e-08] [2.14312031e-08] [2.388341e-14]
case_name_2 = "c2"
prob.record(case_name_2)
prob.cleanup()

cr = om.CaseReader('cases.sql')

num_problem_cases = len(cr.list_cases('problem'))
print(num_problem_cases)
2
c1 = cr.get_case(case_name_1)
c2 = cr.get_case(case_name_2)

# check that derivatives have been recorded properly.
print(c1.derivatives[('f_xy', 'x')][0])
[-1.50431553e-05]
print(c1.derivatives[('f_xy', 'y')][0])
[-1.50431553e-05]
print(c2.derivatives[('f_xy', 'x')][0])
[-1.11442227e-06]
print(c2.derivatives[('f_xy', 'y')][0])
[1.11442256e-06]

For both Driver and Problem, the recording of the derivatives is not affected by the includes and excludes options.