Defining Python and symbolic functions: a quick start introduction to functions

Abstract

In this example, we show how to define Python and symbolic functions. Such functions can be evaluated by the library and used, for example, to propagate uncertainties. We focus on functions which have a vector input and a vector output.

Introduction

We consider the following example: * three input variables, * two outputs.

Moreover, we assume that the input distribution has independent Gaussian marginals.

The function is defined by the equations:

Y_1 = X_1 + X_2 + X_3

and

Y_2 = X_1 - X_2 X_3

for any X_1,X_2,X_3\mathbb{R}.

The exact expectation and standard deviation of the output random variable are presented in the following table.

Variable

Expectation

Standard deviation

Y_1

0

1.732

Y_2

0

1.415

[1]:
import openturns as ot

We first define the input random vector of the function.

[2]:
X0 = ot.Normal(0.,1.)
X1 = ot.Normal(0.,1.)
X2 = ot.Normal(0.,1.)
inputDistribution = ot.ComposedDistribution((X0,X1,X2))
inputRandomVector = ot.RandomVector(inputDistribution)

The Python function

Based on a Python function defined with the def keyword, the PythonFunction class creates a function.

The simplest constructor of the PythonFunction class is:

PythonFunction ( nbInputs , nbOutputs , myPythonFunc )

where

  • nbInputs: the number of inputs,

  • nbOutputs: the number of outputs,

  • myPythonFunc: a Python function.

The goal of the PythonFunction class are:

  • to easily create a function in Python,

  • use all the power of the Python libraries in order to evaluate the function.

The function mySimulator has the calling sequence y=mySimulator(x) where:

  • x: the input of the function, a vector with nbInputs dimensions,

  • y: the output of the function, a vector with nbOutputs dimensions.

[3]:
def mySimulator(x):
    y0=x[0]+x[1]+x[2]
    y1=x[0]-x[1]*x[2]
    y=[y0,y1]
    return y

We now define the PythonFunction object.

[4]:
myfunction = ot.PythonFunction (3 ,2 , mySimulator )

This function can be evaluated using parentheses. It produces the same outputs as the mySimulator function.

[5]:
myfunction([1.,2.,3.])
[5]:

[6,-5]

However, the newly created myfunction has services that the basic Python function did not have. For example, we can create a CompositeRandomVector on top of it, by associating it to the input random vector.

[6]:
outputVect =  ot.CompositeRandomVector(myfunction, inputRandomVector)

In the following example, we estimate the output mean based on a simple Monte-Carlo experiment using 10000 function evaluations.

[7]:
montecarlosize = 10000
outputSample = outputVect.getSample(montecarlosize)
[8]:
empiricalMean = outputSample.computeMean()
print(empiricalMean)
empiricalSd = outputSample.computeStandardDeviationPerComponent()
print(empiricalSd)
[-0.0166778,-0.0123527]
[1.73359,1.39888]

What types for x and for y ?

Not all types are possible for the inputs and outputs of the mySimulator function. The following table present some of the available types. All in all, any type which can be converted to or from a “vector” can be managed by the PythonFunction class.

Type

Input X

Output Y

list (Python)

tuple (Python)

array (NumPy)

Point (OpenTURNS)

The vectorized Python function

The PythonFunction class has a func_sample option which vectorizes the computation so that all the output values in the sample are computed from a single function call, without any for loop. To make this possible, the input and output is then a Sample instead of a Point. This often boosts performance (but not always).

The calling sequence of a vectorized Python function is:

def mySimulator (x):
    [...]
    return y

myfunction = PythonFunction(nbInputs, nbOutputs, func_sample = mySimulator)

where

  • x: the input of the function, a Sample with size nbExperiments (getSize) and dimension nbInputs (getDimension),

  • y: the output of the function

    • a numpy array with nbExperiments rows and nbOutputs columns

    • or a Sample with size nbExperiments and dimension nbOutputs

In the following, we present an vectorization example based on the numpy module.

[9]:
import numpy as np
[10]:
def mySimulatorVect (x):
    # Convert NumericalSample > Array Numpy
    x = np.array (x)
    x0 = x[: ,0] # Extract column 0
    x1 = x[: ,1]
    x2 = x[: ,2]
    y0 = x0 + x1 + x2
    y1 = x0 - x1 * x2
    # Stack the two rows
    y = np.vstack ((y0 ,y1 ))
    y = y.transpose ()
    return y
[11]:
myfunctionVect = ot.PythonFunction (3, 2, func_sample = mySimulatorVect )
[12]:
outputVect =  ot.CompositeRandomVector(myfunctionVect, inputRandomVector)
[13]:
montecarlosize = 10000
outputSample = outputVect.getSample(montecarlosize)

empiricalMean = outputSample.computeMean()
print(empiricalMean)
empiricalSd = outputSample.computeStandardDeviationPerComponent()
print(empiricalSd)
[-0.00077008,-0.0165189]
[1.74231,1.40859]

How to manage the history

The MemoizeFunction class defines a history system to store the calls to the function.

Methods

Description

enableHistory()

enables the history (it is enabled by default)

disableHistory()

disables the history

clearHistory()

deletes the content of the history

getHistoryInput()

a Sample, the history of inputs X

getHistoryOutput()

a Sample, the history of outputs Y

[14]:
myfunction = ot.PythonFunction (3 ,2 , mySimulator )
myfunction = ot.MemoizeFunction(myfunction)
[15]:
outputVariableOfInterest =  ot.CompositeRandomVector(myfunction, inputRandomVector)
montecarlosize = 10
outputSample = outputVariableOfInterest.getSample(montecarlosize)

Get the history of input points.

[16]:
inputs = myfunction.getInputHistory()
inputs
[16]:
v0v1v2
0-0.3426049407409974-0.6896369668677911.2430871864340118
1-2.1386534433938884-1.46352279717364841.1026255045622082
2-0.62127165276042030.05683584819327807-1.0104143995758796
30.140945823728705450.55327862464260220.6904336896063326
40.3056622172571174-1.13650495256445931.692226509358007
5-1.51879905690062840.9695327429694076-1.7656395366104742
61.86032637836427121.21487963389033160.812872250567669
70.32770564349706716-0.6888171650020032-2.011791373087136
8-0.395283255649119871.841280936503913-1.3685050979414375
90.45610178231733167-0.13137430685484950.572819632121423

Symbolic functions

The SymbolicFunction class can create symbolic functions whenever the function is simple enough to be expressed as a string.

This has at least two significant advantages.

  • It generally improves the performance.

  • This allows to automatically evaluate the exact gradient and Hessian matrix.

In practice, some functions cannot be expressed as a symbolic function: in this case, the Python function cannot be avoided.

The calling sequence is the following:

myfunction = SymbolicFunction( list_of_inputs, list_of_formulas)

where

  • list_of_inputs: a list of nbInputs strings, the names of the input variables,

  • list_of_formulas: a list of nbOutputs strings, the equations.

[17]:
myfunction = ot.SymbolicFunction(("x0","x1","x2"),("x0 + x1 + x2","x0 - x1 * x2"))

A SymbolicFunction, like any other function, can also have a history.

[18]:
myfunction = ot.MemoizeFunction(myfunction)
[19]:
outputVect =  ot.CompositeRandomVector(myfunction, inputRandomVector)
[20]:
montecarlosize = 10000
outputSample = outputVect.getSample(montecarlosize)
empiricalMean = outputSample.computeMean()
print(empiricalMean)
[-0.000337898,-0.00609942]

Since the history is enabled, we can get the history of outputs of the function.

[21]:
outputs = myfunction.getOutputHistory()
outputs[1:10,:]
[21]:
v0v1
00.8350684490856450.6082917234676658
1-1.5331063675958376-0.3603107009663977
20.3294989655386499-0.5302618334408122
31.69777518925271350.09262297601997455
41.0671706780737071-0.7204318270150004
5-0.3294656311476475-0.043150705383642796
6-1.4948567574549438-0.45129933759593865
71.0942955644067511.9061369435653874
80.326106500781734030.14681166047125288