1#! /usr/bin/env python 2 3from __future__ import print_function 4from openturns import * 5# from math import * 6 7TESTPREAMBLE() 8RandomGenerator.SetSeed(0) 9 10try: 11 # 12 # Physical model 13 # 14 15 EtatLimite = Function( 16 ['X1', 'X2'], ['G'], ["-0.5*(X1-X2)*(X1-X2) - (X1+X2)/(sqrt(2)) + 3"]) 17 18 dim = EtatLimite.getInputDimension() 19 print(dim) 20 21 # 22 # probabilistic model 23 # 24 25 mean = [0.] * dim 26 sigma = [1.] * dim 27 R = IdentityMatrix(dim) 28 29 myDistribution = Normal(mean, sigma, R) 30 31 start = myDistribution.getMean() 32 Covariance = myDistribution.getCovariance() 33 34 # 35 # limit state 36 # 37 38 vect = RandomVector(myDistribution) 39 40 output = RandomVector(EtatLimite, vect) 41 42 myEvent = ThresholdEvent(output, Less(), 0.0) 43 44 # 45 # Calculation 46 # 47 # 48 # FORM/SORM Cobyla 49 myCobyla = Cobyla() 50 myCobyla.setMaximumEvaluationNumber(100 * dim) 51 myCobyla.setMaximumAbsoluteError(1.0e-10) 52 myCobyla.setMaximumRelativeError(1.0e-10) 53 myCobyla.setMaximumResidualError(1.0e-10) 54 myCobyla.setMaximumConstraintError(1.0e-10) 55 56 myAlgoC = FORM(myCobyla, myEvent, start) 57 myAlgoC2 = SORM(myCobyla, myEvent, start) 58 59 myAlgoC.run() 60 myAlgoC2.run() 61 62 resultC = FORMResult(myAlgoC.getResult()) 63 resultC2 = SORMResult(myAlgoC2.getResult()) 64 65 # 66 # FORM/SORM Abdo Rackwitz 67 myAbdoRackwitz = AbdoRackwitz() 68 myAbdoRackwitz.setMaximumIterationNumber(100 * dim) 69 myAbdoRackwitz.setMaximumAbsoluteError(1.0e-10) 70 myAbdoRackwitz.setMaximumRelativeError(1.0e-10) 71 myAbdoRackwitz.setMaximumResidualError(1.0e-10) 72 myAbdoRackwitz.setMaximumConstraintError(1.0e-10) 73 74 myAlgoAR = FORM(myAbdoRackwitz, myEvent, start) 75 myAlgoAR2 = SORM(myAbdoRackwitz, myEvent, start) 76 77 myAlgoAR.run() 78 # myAlgoAR2.run() 79 80 resultAR = FORMResult(myAlgoAR.getResult()) 81 # resultAR2 = SORMResult(myAlgoAR2.getResult()) 82 83 # 84 # Monte Carlo 85 CoV_MC = 0.5 86 myMC = MonteCarlo(myEvent) 87 myMC.setMaximumOuterSampling(100000) 88 myMC.setBlockSize(1) 89 myMC.setMaximumCoefficientOfVariation(CoV_MC) 90 myMC.run() 91 92 # 93 # LHS 94 CoV_LHS = 0.1 95 myLHS = LHS(myEvent) 96 myLHS.setMaximumOuterSampling(100000) 97 myLHS.setBlockSize(1) 98 myLHS.setMaximumCoefficientOfVariation(CoV_LHS) 99 myLHS.run() 100 101 # 102 # Outputs 103 # 104 print("") 105 print("") 106 print( 107 "************************************************************************************************") 108 print( 109 "***************************************** FORM COBYLA *****************************************") 110 print( 111 "************************************************************************************************") 112 print("event probability = %.5e" % PfC) 113 print("generalized reliability index = %.5f" % Beta_generalizedC) 114 print( 115 "************************************************************************************************") 116 for i in range(u_starC.getDimension()): 117 print("standard space design point = %.5f" % u_starC[i]) 118 print( 119 "************************************************************************************************") 120 for i in range(x_starC.getDimension()): 121 print("physical space design point = %.5f" % x_starC[i]) 122 print( 123 "************************************************************************************************") 124 print("is standard point origin in failure space? ", PtC) 125 print( 126 "************************************************************************************************") 127 for i in range(gammaC.getDimension()): 128 print("importance factors = %.5f" % gammaC[i]) 129 print( 130 "************************************************************************************************") 131 print("Hasofer reliability index = %.5f" % beta_hasoferC) 132 print( 133 "************************************************************************************************") 134 for i in range(SensitivityC.getSize()): 135 for j in range(SensitivityC[i].getDimension()): 136 print("Pf sensitivity = %.5f" % SensitivityC[i][j]) 137 print( 138 "************************************************************************************************") 139 print("") 140 print( 141 "************************************************************************************************") 142 print( 143 "************************************** FORM ABDO RACKWITZ **************************************") 144 print( 145 "************************************************************************************************") 146 print("event probability = %.5e" % PfAR) 147 print("generalized reliability index = %.5f" % Beta_generalizedAR) 148 print( 149 "************************************************************************************************") 150 for i in range(u_starAR.getDimension()): 151 print("standard space design point = %.5f" % u_starAR[i]) 152 print( 153 "************************************************************************************************") 154 for i in range(x_starAR.getDimension()): 155 print("physical space design point = %.5f" % x_starAR[i]) 156 print( 157 "************************************************************************************************") 158 print("is standard point origin in failure space? ", PtAR) 159 print( 160 "************************************************************************************************") 161 for i in range(gammaAR.getDimension()): 162 print("importance factors = %.5f" % gammaAR[i]) 163 print( 164 "************************************************************************************************") 165 print("Hasofer reliability index = %.5f" % beta_hasoferAR) 166 print( 167 "************************************************************************************************") 168 for i in range(SensitivityAR.getSize()): 169 for j in range(SensitivityAR[i].getDimension()): 170 print("Pf sensitivity = %.5f" % SensitivityAR[i][j]) 171 print( 172 "************************************************************************************************") 173 print("") 174 print( 175 "************************************************************************************************") 176 print( 177 "***************************************** SORM COBYLA *****************************************") 178 print( 179 "************************************************************************************************") 180 print("Breitung event probability = %.5e" % PFBreitC2) 181 print("Breitung generalized reliability index = %.5f" % BetaBreitC2) 182 print("Hohenbichler event probability = %.5e" % PFHBC2) 183 print("Hohenbichler generalized reliability index = %.5f" % BetaHBC2) 184 print("Tvedt event probability = %.5e" % PFTvedtC2) 185 print("Tvedt generalized reliability index = %.5f" % BetaTvedtC2) 186 print( 187 "************************************************************************************************") 188 for i in range(CurvC2.getDimension()): 189 print("sorted curvatures = %.5f" % cleanScalar(CurvC2[i])) 190 print( 191 "************************************************************************************************") 192 for i in range(u_starC2.getDimension()): 193 print("standard space design point = %.5f" % u_starC2[i]) 194 print( 195 "************************************************************************************************") 196 for i in range(x_starC2.getDimension()): 197 print("physical space design point = %.5f" % x_starC2[i]) 198 print( 199 "************************************************************************************************") 200 print( 201 "************************************************************************************************") 202 print("is standard point origin in failure space? ", PtC2) 203 print( 204 "************************************************************************************************") 205 for i in range(gammaC2.getDimension()): 206 print("importance factors = %.5f" % gammaC2[i]) 207 print( 208 "************************************************************************************************") 209 print("Hasofer reliability index = %.5f" % beta_hasoferC2) 210 print( 211 "************************************************************************************************") 212 print("") 213 print( 214 "************************************************************************************************") 215 print( 216 "************************************** SORM ABDO RACKWITZ **************************************") 217 print( 218 "************************************************************************************************") 219 print("Breitung event probability = %.5e" % PFBreitAR2) 220 print("Breitung generalized reliability index = %.5f" % BetaBreitAR2) 221 print("Hohenbichler event probability = %.5e" % PFHBAR2) 222 print("Hohenbichler generalized reliability index = %.5f" % BetaHBAR2) 223 print("Tvedt event probability = %.5e" % PFTvedtAR2) 224 print("Tvedt generalized reliability index = %.5f" % BetaTvedtAR2) 225 print( 226 "************************************************************************************************") 227 for i in range(CurvAR2.getDimension()): 228 print("sorted curvatures = %.5f" % cleanScalar(CurvAR2[i])) 229 print( 230 "************************************************************************************************") 231 for i in range(u_starAR2.getDimension()): 232 print("standard space design point = %.5f" % u_starAR2[i]) 233 print( 234 "************************************************************************************************") 235 for i in range(x_starAR2.getDimension()): 236 print("physical space design point = %.5f" % x_starAR2[i]) 237 print( 238 "************************************************************************************************") 239 print( 240 "************************************************************************************************") 241 print("is standard point origin in failure space? ", PtAR2) 242 print( 243 "************************************************************************************************") 244 for i in range(gammaAR2.getDimension()): 245 print("importance factors = %.5f" % gammaAR2[i]) 246 print( 247 "************************************************************************************************") 248 print("Hasofer reliability index = %.5f" % beta_hasoferAR2) 249 print( 250 "************************************************************************************************") 251 print("") 252 print( 253 "************************************************************************************************") 254 print( 255 "**************************************** MONTE CARLO *******************************************") 256 print( 257 "************************************************************************************************") 258 print("Pf estimation = %.5e" % PFMC) 259 print("Pf Variance estimation = %.5e" % Variance_PF_MC) 260 print("CoV = %.5f" % CVMC) 261 print("90% Confidence Interval =", "%.5e" % length90MC) 262 print("CI at 90% =[", "%.5e" % (PFMC - 0.5 * length90MC), 263 "; %.5e" % (PFMC + 0.5 * length90MC), "]") 264 print( 265 "************************************************************************************************") 266 print("") 267 print( 268 "************************************************************************************************") 269 print( 270 "******************************************* L H S **********************************************") 271 print( 272 "************************************************************************************************") 273 print("Pf estimation = %.5e" % PFLHS) 274 print("Pf Variance estimation = %.5e" % Variance_PF_LHS) 275 print("CoV = %.5f" % CVLHS) 276 print("90% Confidence Interval =", "%.5e" % length90LHS) 277 print("CI at 90% =[", "%.5e" % (PFLHS - 0.5 * length90LHS), 278 "; %.5e" % (PFLHS + 0.5 * length90LHS), "]") 279 print( 280 "************************************************************************************************") 281 print("") 282 283except: 284 import sys 285 print("t_Waarts_concave.py", sys.exc_info()[0], sys.exc_info()[1]) 286