1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
|
#! /usr/bin/env python
import openturns as ot
import openturns.testing as ott
import sys
ot.TESTPREAMBLE()
ot.PlatformInfo.SetNumericalPrecision(3)
ot.Log.Show(ot.Log.ALL)
def progress(percent):
sys.stderr.write("-- progress=" + str(percent) + "%\n")
def stop():
sys.stderr.write("-- stop?\n")
return False
n = 3
m = 10
x = [[0.5 + i] for i in range(m)]
model = ot.SymbolicFunction(["a", "b", "c", "x"], ["a + b * exp(c * x)"])
p_ref = [2.8, 1.2, 0.5] # a, b, c
modelx = ot.ParametricFunction(model, [0, 1, 2], p_ref)
y = modelx(x)
def residualFunction_py(p):
modelx = ot.ParametricFunction(model, [0, 1, 2], p)
return [modelx(x[i])[0] - y[i, 0] for i in range(m)]
residualFunction = ot.PythonFunction(n, m, residualFunction_py)
bounds = ot.Interval([0, 0, 0], [2.5, 8.0, 19])
for bound in [True, False]:
problem = ot.LeastSquaresProblem(residualFunction)
if bound:
problem.setBounds(bounds)
startingPoint = [1.0] * n
algo = ot.CMinpack(problem)
algo.setStartingPoint(startingPoint)
# algo.setProgressCallback(progress)
# algo.setStopCallback(stop)
algo.run()
result = algo.getResult()
# print(result.getInputSample())
# print(result.getOutputSample())
x_star = result.getOptimalPoint()
print(result)
if bound:
assert x_star in bounds, "optimal point not in bounds"
else:
ott.assert_almost_equal(x_star, p_ref)
|
