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
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
|
#! /usr/bin/env python
import openturns as ot
import persalys
code = """
import openturns as ot
def _exec(x1, x2, x3, x4, x5, x6, x7, x8, x9, x10, x11, x12, x13, x14, x15, x16, x17, x18, x19, x20):
ot.RandomGenerator.SetSeed(123)
b0 = ot.DistFunc.rNormal()
alpha = ot.DistFunc.rNormal(10)
beta = ot.DistFunc.rNormal(6*14)
gamma = ot.DistFunc.rNormal(20*14)
# b1
b1 = [20] * 10 + list(alpha)
# b2
b2 = [[0] *20] * 20
for i in range(6):
for j in range(6):
b2[i][j] = -15.0
# Take into account beta
k = 0
for i in range(6):
for j in range(14):
b2[i][j + 6] = beta[k]
k = k + 1
# Take into account gamma
k = 0
for i in range(6, 20):
for j in range(20):
b2[i][j] = gamma[k]
# b3
b3 = [[[0]*20]*20]*20
for i in range(5):
for j in range(5):
for k in range(5):
b3[i][j][k] = -10.0
# b4
b4 = [[[[0]*20]*20]*20]*20
for i in range(4):
for j in range(4):
for k in range(4):
for l in range(4):
b4[i][j][k][l] = 5
# X is a list, transform it into array
X = ot.Point([x1, x2, x3, x4, x5, x6, x7, x8, x9, x10, x11, x12, x13, x14, x15, x16, x17, x18, x19, x20])
w = (X - [0.5]*20)*2
for k in [2,4,6]:
w[k] = 2.0 * (1.1 * X[k] / (X[k] + 0.1) - 0.5)
y = w.dot(b1)
# Morris function
for i in range(19):
for j in range(i + 1, 20):
y += b2[i][j] * w[i] * w[j]
for i in range(18):
for j in range(i + 1, 19):
for k in range(j + 1, 20):
y += b3[i][j][k] * w[i] * w[j] * w[k]
for i in range(17):
for j in range(i + 1, 18):
for k in range(j + 1, 20):
for l in range(k + 1, 20):
y += b4[i][j][k][l] * w[i] * w[j] * w[k] * w[l]
y = y + b0
y_fake = 0.
return y
"""
myStudy = persalys.Study("myStudy")
inputs = []
for i in range(20):
inputs.append(persalys.Input("x" + str(i + 1), ot.Uniform(0, 1)))
y = persalys.Output("y")
y_fake = persalys.Output("y_fake")
print(inputs)
model = persalys.PythonPhysicalModel("MorrisModel", inputs, [y, y_fake], code)
model.setProcessNumber(1)
myStudy.add(model)
# Morris ##
analysis = persalys.MorrisAnalysis("aMorris", model)
analysis.setInterestVariables(["y"])
analysis.setLevel(4)
analysis.setTrajectoriesNumber(4)
analysis.setSeed(2)
myStudy.add(analysis)
print(analysis)
analysis.run()
print("error message=", analysis.getErrorMessage())
print("result=", analysis.getResult())
print("meanAbsEE =", analysis.getResult().getMeanAbsoluteElementaryEffects(0))
# script
script = myStudy.getPythonScript()
persalys.Study.Add(myStudy)
print(script)
exec(script)
|
