.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_reliability_sensitivity/reliability/plot_form_explained.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_reliability_sensitivity_reliability_plot_form_explained.py: An illustrated example of a FORM probability estimate ===================================================== .. GENERATED FROM PYTHON SOURCE LINES 7-14 Abstract -------- In this example we illustrate the different steps of a FORM/SORM analysis on a simple example. We focus on the different steps and compare them with an analytic computation whenever possible. .. GENERATED FROM PYTHON SOURCE LINES 14-18 .. code-block:: default import openturns as ot import openturns.viewer as otv import numpy as np .. GENERATED FROM PYTHON SOURCE LINES 19-29 Position of the problem ----------------------- We consider a bivariate random vector :math:`X = (X_1, X_2)` with the following independent marginals : - an exponential distribution with parameter :math:`\lambda=1`, :math:`X_1 \sim \mathcal{E}(1.0)` ; - a standard unit gaussian :math:`X_2 \sim \mathcal{N}(0,1)`. The support of the input vector is :math:`[0, +\infty[ \times \mathbb{R}` .. GENERATED FROM PYTHON SOURCE LINES 31-36 .. code-block:: default distX1 = ot.Exponential(1.0) distX2 = ot.Normal() distX = ot.ComposedDistribution([distX1, distX2]) .. GENERATED FROM PYTHON SOURCE LINES 37-38 We can draw the bidimensional PDF of the distribution `distX` over :math:`[0,-10] \times [10,10]` : .. GENERATED FROM PYTHON SOURCE LINES 38-47 .. code-block:: default ot.ResourceMap.SetAsUnsignedInteger("Contour-DefaultLevelsNumber", 8) graphPDF = distX.drawPDF([0, -10], [10, 10]) graphPDF.setTitle(r"2D-PDF of the input variables $(X_1, X_2)$") graphPDF.setXTitle(r"$x_1$") graphPDF.setYTitle(r"$x_2$") graphPDF.setLegendPosition("bottomright") view = otv.View(graphPDF) .. image-sg:: /auto_reliability_sensitivity/reliability/images/sphx_glr_plot_form_explained_001.png :alt: 2D-PDF of the input variables $(X_1, X_2)$ :srcset: /auto_reliability_sensitivity/reliability/images/sphx_glr_plot_form_explained_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 48-50 We consider the model :math:`f : (x_1, x_2) \mapsto x_1 x_2` which maps the random input vector :math:`X` to the output variable :math:`Y=f(X) \in \mathbb{R}`. We also draw the isolines of the model `f`. .. GENERATED FROM PYTHON SOURCE LINES 50-58 .. code-block:: default f = ot.SymbolicFunction(["x1", "x2"], ["x1 * x2"]) graphModel = f.draw([0.0, -10.0], [10.0, 10.0]) graphModel.setXTitle(r"$x_1$") graphModel.setXTitle(r"$x_2$") graphModel.setTitle(r"Isolines of the model : $Y = f(X)$") view = otv.View(graphModel) .. image-sg:: /auto_reliability_sensitivity/reliability/images/sphx_glr_plot_form_explained_002.png :alt: Isolines of the model : $Y = f(X)$ :srcset: /auto_reliability_sensitivity/reliability/images/sphx_glr_plot_form_explained_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 59-67 We want to estimate the probability :math:`P_f` of the output variable to be greater than a prescribed threshold :math:`s=10` : this is the failure event. This probability is simply expressed as an integral : .. math:: P_f = \int_{\mathcal{D}} \mathbf{1}_{\mathcal{D}}(x) df_{X_1,X_2}(x) where :math:`\mathcal{D} = \{ (x_1, x_2) \in [0,+\infty[ \times \mathbb{R} / x_1 x_2 \geq s \}` is the failure domain. In the general case the probability density function :math:`f_{X_1,X_2}` and the domain of integration :math:`\mathcal{D}` are difficult to handle. .. GENERATED FROM PYTHON SOURCE LINES 69-70 We first define RandomVector objects and the failure event associated to the output random variable. .. GENERATED FROM PYTHON SOURCE LINES 70-76 .. code-block:: default vectorX = ot.RandomVector(distX) vectorY = ot.CompositeRandomVector(f, vectorX) s = 10.0 event = ot.ThresholdEvent(vectorY, ot.Greater(), s) .. GENERATED FROM PYTHON SOURCE LINES 77-84 This event can easily be represented with a 1D curve as it is a branch of an hyperbole. If :math:`y = x_1 x_2 = 10.0`, then the boundary of the domain of failure is the curve : .. math:: h : x_1 \mapsto \frac{10.0}{x_1} .. GENERATED FROM PYTHON SOURCE LINES 86-87 We shall represent this curve using a :class:`~openturns.Contour` object. .. GENERATED FROM PYTHON SOURCE LINES 87-96 .. code-block:: default nx, ny = 15, 15 xx = ot.Box([nx], ot.Interval([0.0], [10.0])).generate() yy = ot.Box([ny], ot.Interval([-10.0], [10.0])).generate() inputData = ot.Box([nx, ny], ot.Interval([0.0, -10.0], [10.0, 10.0])).generate() outputData = f(inputData) mycontour = ot.Contour(xx, yy, outputData, [10.0], ["10.0"]) myGraph = ot.Graph("Representation of the failure domain", r"$X_1$", r"$X_2$", True, "") myGraph.add(mycontour) .. GENERATED FROM PYTHON SOURCE LINES 97-103 .. code-block:: default texts = [r" Event : $\mathcal{D} = \{Y \geq 10.0\}$"] myText = ot.Text([[4.0, 4.0]], texts) myText.setTextSize(1) myGraph.add(myText) view = otv.View(myGraph) .. image-sg:: /auto_reliability_sensitivity/reliability/images/sphx_glr_plot_form_explained_003.png :alt: Representation of the failure domain :srcset: /auto_reliability_sensitivity/reliability/images/sphx_glr_plot_form_explained_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 104-106 We can superimpose the event boundary with the 2D-PDF ot the input variables : .. GENERATED FROM PYTHON SOURCE LINES 106-112 .. code-block:: default mycontour.setColor("black") mycontour.setLabels(["event"]) graphPDF.add(mycontour) graphPDF.setLegendPosition("bottomright") view = otv.View(graphPDF) .. image-sg:: /auto_reliability_sensitivity/reliability/images/sphx_glr_plot_form_explained_004.png :alt: 2D-PDF of the input variables $(X_1, X_2)$ :srcset: /auto_reliability_sensitivity/reliability/images/sphx_glr_plot_form_explained_004.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 113-115 From the previous figure we observe that in the failure domain the PDF takes small (and even very small) values. Consequently the probability of the failure, the integral :math:`P_f` is also expected to be small. The FORM/SORM methods estimate this kind of integral. .. GENERATED FROM PYTHON SOURCE LINES 117-133 The FORM approximation ---------------------- The basic steps of the FORM (or SORM) algorithm are : - an isoprobabilistic transform ; - finding the design point : that is the nearest point wrt the origin in the standard space ; - estimating the probability integral. As mentioned, both the density function and the domain of integration are complex in general. The first step of the FORM method makes the density easier to work with and the second step tackles the domain of integration problem. Variable transform ^^^^^^^^^^^^^^^^^^ .. GENERATED FROM PYTHON SOURCE LINES 135-138 OpenTURNS has several isoprobabilistic transforms and the FORM/SORM classes implement the Generalized Nataf and Rosenblatt transforms. In this case the `distX` distribution is not elliptical so the default method is the Rosenblatt transform. .. GENERATED FROM PYTHON SOURCE LINES 138-140 .. code-block:: default print("Is Elliptical ? ", distX.isElliptical()) .. rst-class:: sphx-glr-script-out .. code-block:: none Is Elliptical ? False .. GENERATED FROM PYTHON SOURCE LINES 141-149 We seek an isoprobabilistic transform :math:`T` such as .. math:: T : X \mapsto Z where each component of :math:`Z` is a standard unit gaussian. The isoprobabilistic transform and its inverse are methods of the distribution `distX` : .. GENERATED FROM PYTHON SOURCE LINES 149-152 .. code-block:: default transformation = distX.getIsoProbabilisticTransformation() inverseTransformation = distX.getInverseIsoProbabilisticTransformation() .. GENERATED FROM PYTHON SOURCE LINES 153-157 The main goal of this step is to work with a simpler probability density function of the input variables as they will be standard gaussian unit and uncorrelated. The domain of integration will still be complicated but will be handled with a well chosen approximate. .. GENERATED FROM PYTHON SOURCE LINES 159-164 We detail the Rosenblatt transform in this simple case. In this example we consider independent variables so the transform is simpler, we only have to perform the transformation on each variable. For the second one is already a standard unit gaussian we transform the first variable only. .. GENERATED FROM PYTHON SOURCE LINES 166-168 We draw a realization of the random input vector. This point is said to be in the physical space. We shall focus on the first component. .. GENERATED FROM PYTHON SOURCE LINES 168-170 .. code-block:: default xi = vectorX.getRealization() .. GENERATED FROM PYTHON SOURCE LINES 171-175 The first step of the Rosenblatt transform is to build a random variable :math:`u` with a uniform law in ]0,1[. This is done through an evaluation of the CDF of `distX1` at the given point in the physical space. Once again, please note that the second component is left unchanged. .. GENERATED FROM PYTHON SOURCE LINES 175-177 .. code-block:: default ui = [distX1.computeCDF(xi[0]), xi[1]] .. GENERATED FROM PYTHON SOURCE LINES 178-180 The second step is to build a standard unit gaussian from a uniform variable. This is done by a simple call to the probit function. The point `zi` is said to be in the standard space. .. GENERATED FROM PYTHON SOURCE LINES 180-182 .. code-block:: default zi = [-ot.Normal().computeInverseSurvivalFunction(ui[0])[0], ui[1]] .. GENERATED FROM PYTHON SOURCE LINES 183-184 The sought transform then maps a point in the physical space to the standard space : .. GENERATED FROM PYTHON SOURCE LINES 184-187 .. code-block:: default print(xi, "->", ui, "->", zi) .. rst-class:: sphx-glr-script-out .. code-block:: none [0.481072,-0.267264] -> [0.3818794784091921, -0.2672639149611159] -> [-0.30054830438164254, -0.2672639149611159] .. GENERATED FROM PYTHON SOURCE LINES 188-190 We also build the isoprobabilistic transform :math:`T_1` and its inverse :math:`T_1^{-1}` for the first marginal : .. GENERATED FROM PYTHON SOURCE LINES 190-193 .. code-block:: default transformX1 = distX1.getIsoProbabilisticTransformation() inverseTransformX1 = distX1.getInverseIsoProbabilisticTransformation() .. GENERATED FROM PYTHON SOURCE LINES 194-200 We can check the result of our experiment against : - the 2D-transform :math:`T` ; - the 1D-transform :math:`T_1` and the second component unchanged ; and observe the results are the same. .. GENERATED FROM PYTHON SOURCE LINES 200-207 .. code-block:: default zi1D = [transformX1([xi[0]])[0], xi[1]] zi2D = transformation(xi) print("zi = ", zi) print("zi1D = ", zi1D) print("zi2D = ", zi2D) .. rst-class:: sphx-glr-script-out .. code-block:: none zi = [-0.30054830438164254, -0.2672639149611159] zi1D = [-0.30054830438164254, -0.2672639149611159] zi2D = [-0.300548,-0.267264] .. GENERATED FROM PYTHON SOURCE LINES 208-211 We can represent the boundary of the event in the standard space : that is a composition of the hyperbole :math:`h : x \mapsto 10/x` and the inverse transform :math:`T_1^{-1}` defined by :math:`inverseTransformX1`. .. GENERATED FROM PYTHON SOURCE LINES 211-227 .. code-block:: default failureBoundaryPhysicalSpace = ot.SymbolicFunction(["x"], ["10.0 / x"]) failureBoundaryStandardSpace = ot.ComposedFunction( failureBoundaryPhysicalSpace, inverseTransformX1 ) x = np.linspace(1.1, 5.0, 100) cx = np.array([failureBoundaryStandardSpace([xi])[0] for xi in x]) graphStandardSpace = ot.Graph( "Failure event in the standard space", r"$u_1$", r"$u_2$", True, "" ) curveCX = ot.Curve(x, cx, r"Boundary of the event $\partial \mathcal{D}$") curveCX.setLineStyle("solid") curveCX.setColor("blue") graphStandardSpace.add(curveCX) .. GENERATED FROM PYTHON SOURCE LINES 228-229 We add the origin to the previous graph. .. GENERATED FROM PYTHON SOURCE LINES 229-244 .. code-block:: default cloud = ot.Cloud([0.0], [0.0]) cloud.setColor("black") cloud.setPointStyle("fcircle") cloud.setLegend("origin") graphStandardSpace.add(cloud) graphStandardSpace.setGrid(True) graphStandardSpace.setLegendPosition("bottomright") # Some annotation texts = [r"Event : $\mathcal{D} = \{Y \geq 10.0\}$"] myText = ot.Text([[3.0, 4.0]], texts) myText.setTextSize(1) graphStandardSpace.add(myText) view = otv.View(graphStandardSpace) .. image-sg:: /auto_reliability_sensitivity/reliability/images/sphx_glr_plot_form_explained_005.png :alt: Failure event in the standard space :srcset: /auto_reliability_sensitivity/reliability/images/sphx_glr_plot_form_explained_005.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 245-253 The design point ^^^^^^^^^^^^^^^^ The FORM and SORM methods assume that the failure probability integral has its support in the vicinity of the closest point of the domain to the origin. Then the second step of the method is to find this point, *the design point*, through a minimization problem under constraints. .. GENERATED FROM PYTHON SOURCE LINES 255-256 We configure the Cobyla solver that we use for the optimization : .. GENERATED FROM PYTHON SOURCE LINES 256-263 .. code-block:: default solver = ot.Cobyla() solver.setMaximumIterationNumber(10000) solver.setMaximumAbsoluteError(1.0e-3) solver.setMaximumRelativeError(1.0e-3) solver.setMaximumResidualError(1.0e-3) solver.setMaximumConstraintError(1.0e-3) .. GENERATED FROM PYTHON SOURCE LINES 264-266 We build the FORM algorithm with its basic constructor. The starting point for the optimization algorithm is the mean of the input variables. .. GENERATED FROM PYTHON SOURCE LINES 266-268 .. code-block:: default algoFORM = ot.FORM(solver, event, distX.getMean()) .. GENERATED FROM PYTHON SOURCE LINES 269-270 We are ready to run the algorithm and store the result : .. GENERATED FROM PYTHON SOURCE LINES 270-273 .. code-block:: default algoFORM.run() result = algoFORM.getResult() .. GENERATED FROM PYTHON SOURCE LINES 274-276 The design point can be retrieved in both physical and standard space with respectively the `getPhysicalSpaceDesignPoint` and `getStandardSpaceDesignPoint` methods. .. GENERATED FROM PYTHON SOURCE LINES 276-282 .. code-block:: default designPointPhysicalSpace = result.getPhysicalSpaceDesignPoint() designPointStandardSpace = result.getStandardSpaceDesignPoint() print("Design point in physical space : ", designPointPhysicalSpace) print("Design point in standard space : ", designPointStandardSpace) .. rst-class:: sphx-glr-script-out .. code-block:: none Design point in physical space : [4.84183,2.06513] Design point in standard space : [2.41384,2.06513] .. GENERATED FROM PYTHON SOURCE LINES 283-284 We can get the Hasofer index with the `getHasoferReliabilityIndex` method which is the distance of the design point to the origin : .. GENERATED FROM PYTHON SOURCE LINES 284-288 .. code-block:: default betaHL = result.getHasoferReliabilityIndex() print("Hasofer index : ", betaHL) .. rst-class:: sphx-glr-script-out .. code-block:: none Hasofer index : 3.176696193499833 .. GENERATED FROM PYTHON SOURCE LINES 289-290 We visualize it on the previous graph. .. GENERATED FROM PYTHON SOURCE LINES 290-307 .. code-block:: default cloud = ot.Cloud([designPointStandardSpace[0]], [designPointStandardSpace[1]]) cloud.setColor("red") cloud.setPointStyle("fcircle") cloud.setLegend("design point") graphStandardSpace.add(cloud) graphStandardSpace.setGrid(True) graphStandardSpace.setLegendPosition("bottomright") cc = ot.Curve( [0.0, designPointStandardSpace[0]], [0.0, designPointStandardSpace[1]], r"$\beta_{HL}$ distance", ) cc.setLineStyle("dashed") cc.setColor("black") graphStandardSpace.add(cc) view = otv.View(graphStandardSpace) .. image-sg:: /auto_reliability_sensitivity/reliability/images/sphx_glr_plot_form_explained_006.png :alt: Failure event in the standard space :srcset: /auto_reliability_sensitivity/reliability/images/sphx_glr_plot_form_explained_006.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 308-319 Estimating the failure probability integral ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ The last step of the FORM algorithm is to replace the domain of integration by the half-space at the design point. In this simple example the half-space is delimited by the tangent at the design point in the standard space. The expression of the failure domain boundary in the standard space is the composition of the hyperbole :math:`h:x \mapsto 10/x` and the inverse transform on the first variable. We can compute the gradient (here the first derivative of a 1D function :math:`h(u_0)` ) at any given point with the getGradient method : .. GENERATED FROM PYTHON SOURCE LINES 321-328 .. code-block:: default u0 = [designPointStandardSpace[0]] du0 = failureBoundaryStandardSpace.getGradient().gradient(u0) print("abscissa of the design point u0 = ", u0[0]) print("value of the failure boundary at u0 = ", failureBoundaryStandardSpace(u0)[0]) print("value of the gradient of the failure boundary at u0 = ", du0[0, 0]) .. rst-class:: sphx-glr-script-out .. code-block:: none abscissa of the design point u0 = 2.413844253971037 value of the failure boundary at u0 = 2.065335164471684 value of the gradient of the failure boundary at u0 = -1.1706609709179403 .. GENERATED FROM PYTHON SOURCE LINES 329-335 In the standard space the equation of the tangent :math:`\Pi_{u_0}(x)` is given by .. math:: \Pi_{u_0}(x) = (h \circ T^{-1}) (u_0) + \frac{d}{dx} (h \circ T^{-1}) (u_0) (x-u_0) .. GENERATED FROM PYTHON SOURCE LINES 335-343 .. code-block:: default x = np.linspace(1.1, 5.0, 100) hyperplane = failureBoundaryStandardSpace(u0)[0] + du0[0, 0] * (x - u0) curveHyperplane = ot.Curve(x, hyperplane, r"$\Pi_{u_0}$ (FORM)") curveHyperplane.setLineStyle("dashed") curveHyperplane.setColor("green") graphStandardSpace.add(curveHyperplane) view = otv.View(graphStandardSpace) .. image-sg:: /auto_reliability_sensitivity/reliability/images/sphx_glr_plot_form_explained_007.png :alt: Failure event in the standard space :srcset: /auto_reliability_sensitivity/reliability/images/sphx_glr_plot_form_explained_007.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 344-352 In the standard space the PDF of the input variables is rotationally invariant so .. math:: P_f \approx E(\beta_{HL}), where :math:`E(.)` is the survival function of the standard unit gaussian. .. GENERATED FROM PYTHON SOURCE LINES 352-355 .. code-block:: default pf = ot.Normal().computeSurvivalFunction(betaHL) print("FORM : Pf = ", pf) .. rst-class:: sphx-glr-script-out .. code-block:: none FORM : Pf = 0.0007448149708283228 .. GENERATED FROM PYTHON SOURCE LINES 356-357 This proability of failure is the one computed in the FORMResult and obtained with the `getEventProbability` method : .. GENERATED FROM PYTHON SOURCE LINES 357-361 .. code-block:: default pf = result.getEventProbability() print("Probability of failure (FORM) Pf = ", pf) .. rst-class:: sphx-glr-script-out .. code-block:: none Probability of failure (FORM) Pf = 0.0007448149708283228 .. GENERATED FROM PYTHON SOURCE LINES 362-370 The SORM approximation ---------------------- The SORM approximate uses an osculating paraboloid instead of the half-space delimited by the tangent at the design point. In this case it is a simple parabola we can obtain through Taylor expansion at the design point. However, in the general case one has to manipulate the gradient and the hessian in the standard space which is cumbersome. .. GENERATED FROM PYTHON SOURCE LINES 372-374 We need the value of the second derivative of the failure boundary function at the design point in the standard space : .. GENERATED FROM PYTHON SOURCE LINES 374-380 .. code-block:: default u0 = [designPointStandardSpace[0]] d2u0 = failureBoundaryStandardSpace.getHessian().hessian(u0) print("abscissa of the design point u0 = ", u0[0]) print("value of the hessian of the failure boundary at u0 = ", d2u0[0, 0, 0]) .. rst-class:: sphx-glr-script-out .. code-block:: none abscissa of the design point u0 = 2.413844253971037 value of the hessian of the failure boundary at u0 = 0.9401058369722994 .. GENERATED FROM PYTHON SOURCE LINES 381-387 In the standard space the equation of the osculating parabola :math:`\mathcal{P}_{u_0}(x)` at :math:`u_0` is given by .. math:: \mathcal{P}_{u_0}(x) = h \circ T^{-1} (u_0) + \frac{d}{dx} (h \circ T^{-1})(u_0) (x-u_0) + \frac{1}{2} \frac{d^2}{dx^2} (h \circ T^{-1})(u_0) (x-u_0)^2 .. GENERATED FROM PYTHON SOURCE LINES 387-398 .. code-block:: default x = np.linspace(1.1, 5.0, 100) parabola = ( failureBoundaryStandardSpace(u0)[0] + du0[0, 0] * (x - u0) + 0.5 * d2u0[0, 0, 0] * (x - u0) ** 2 ) curveParabola = ot.Curve(x, parabola, r"$\mathcal{P}_{u_0}$ (SORM)") curveParabola.setLineStyle("dashed") curveParabola.setColor("orange") graphStandardSpace.add(curveParabola) view = otv.View(graphStandardSpace) .. image-sg:: /auto_reliability_sensitivity/reliability/images/sphx_glr_plot_form_explained_008.png :alt: Failure event in the standard space :srcset: /auto_reliability_sensitivity/reliability/images/sphx_glr_plot_form_explained_008.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 399-409 The next step is to estimate the principal curvatures of the osculating paraboloid. For any regular function :math:`g` the curvature :math:`\kappa(x_0)` at the point :math:`x_0` in cartesian coordinates reads as .. math:: \kappa(x_0) = \frac{g''(x_0)}{(1+[g'(x_0)]^2)^{3/2}}. For the oscilating parabola of concern we use the gradient and hessian previously computed : .. GENERATED FROM PYTHON SOURCE LINES 409-413 .. code-block:: default curvature = (d2u0[0, 0, 0]) / (1 + (du0[0, 0]) ** 2) ** (3 / 2) print("Curvature (analytic formula) = ", curvature) .. rst-class:: sphx-glr-script-out .. code-block:: none Curvature (analytic formula) = 0.2575913913877353 .. GENERATED FROM PYTHON SOURCE LINES 414-415 We build the SORM algorithm and run it : .. GENERATED FROM PYTHON SOURCE LINES 415-418 .. code-block:: default algoSORM = ot.SORM(solver, event, distX.getMean()) algoSORM.run() .. GENERATED FROM PYTHON SOURCE LINES 419-420 The SORM result is obtained with the `getResult` method : .. GENERATED FROM PYTHON SOURCE LINES 420-422 .. code-block:: default resultSORM = algoSORM.getResult() .. GENERATED FROM PYTHON SOURCE LINES 423-425 The principal curvatures of the osculating paraboloid at the design point is obtained by the `getSortedCurvatures` method : .. GENERATED FROM PYTHON SOURCE LINES 425-428 .. code-block:: default print("Curvature (estimated) = ", resultSORM.getSortedCurvatures()[1]) .. rst-class:: sphx-glr-script-out .. code-block:: none Curvature (estimated) = 0.25761034541451805 .. GENERATED FROM PYTHON SOURCE LINES 429-437 Once the curvature is obtained there are several ways of approximating the failure probability :math:`P_f`. OpenTURNS implements the Breitung, Hohenbichler and Tvedt estimates. For instance, the Breitung approximation gives .. math:: P_f \approx E(\beta_{HL}) \frac{1}{\sqrt{1+\beta_{HL}\kappa}} .. GENERATED FROM PYTHON SOURCE LINES 437-441 .. code-block:: default coeff = (1.0 + betaHL * curvature) ** (-0.5) pf = (1.0 - ot.Normal().computeCDF(betaHL)) * coeff print("SORM : Pf = ", pf) .. rst-class:: sphx-glr-script-out .. code-block:: none SORM : Pf = 0.0005523531956154747 .. GENERATED FROM PYTHON SOURCE LINES 442-443 We can compare with the different estimators implemented in OpenTURNS : .. GENERATED FROM PYTHON SOURCE LINES 443-451 .. code-block:: default pfBreitung = resultSORM.getEventProbabilityBreitung() pfHohenbichler = resultSORM.getEventProbabilityHohenbichler() pfTvedt = resultSORM.getEventProbabilityTvedt() print("Probability of failure (SORM Breintung) Pf = ", pfBreitung) print("Probability of failure (SORM Hohenbichler) Pf = ", pfHohenbichler) print("Probability of failure (SORM Tvedt) Pf = ", pfTvedt) .. rst-class:: sphx-glr-script-out .. code-block:: none Probability of failure (SORM Breintung) Pf = 0.0005523440504786004 Probability of failure (SORM Hohenbichler) Pf = 0.0005420328660300074 Probability of failure (SORM Tvedt) Pf = 0.0005381057564255441 .. GENERATED FROM PYTHON SOURCE LINES 452-453 Display all figures .. GENERATED FROM PYTHON SOURCE LINES 453-454 .. code-block:: default otv.View.ShowAll() .. rst-class:: sphx-glr-timing **Total running time of the script:** ( 0 minutes 0.843 seconds) .. _sphx_glr_download_auto_reliability_sensitivity_reliability_plot_form_explained.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_form_explained.py ` .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_form_explained.ipynb `