13.6.15. Macro “sensitivitySobolLoadFile.py”
13.6.15.1. Objective
The objective of this macro is to perform a full Sobol analysis using the existing file created when the Sobol class is allowed to perform the design-of-experiments and the estimations by itself when this one is not considered accurate enough (from the statistical point of view). The idea is to use anyway all the computations already done and generate some more to increase the statistical precision. In order to do that there are few things to keep in mind:
The problem should obviously be exactly the same: same input and output (name, statistical laws, parameter values…).
The input file (here
ref_sobol_launching_.dat) should contains the internal variable in order to figure out from what matrices every configuration is taken out of (generally an iterator whose name should look likesobol__n__iter__plus the dataserver name).If the pseudo-random generator seed is set to a given value, be sure that you are not be using the same seed than the once used to generate the input file, which would lead to twice the same events.
One can find another discussion for this objective in the third item in the tip box in Computation of Sobol’s sensitivity indices.
Warning
The ref_sobol_launching_.dat file used as input is not provided in the usual sub-directory
“/share/uranie/macros” of the installation folder of Uranie ($URANIESYS) but can be
generated by the user by running the macro discussed in Macro “sensitivitySobolFunctionFlowrate.py”.
13.6.15.2. Macro Uranie
"""
Example of sobol estimation once a file is loaded
"""
from URANIE import DataServer, Sensitivity
import ROOT
ROOT.gROOT.LoadMacro("UserFunctions.C")
# Define the DataServer
tds = DataServer.TDataServer("tdsflowreate", "DataBase flowreate")
tds.addAttribute(DataServer.TUniformDistribution("rw", 0.05, 0.15))
tds.addAttribute(DataServer.TUniformDistribution("r", 100.0, 50000.0))
tds.addAttribute(DataServer.TUniformDistribution("tu", 63070.0, 115600.0))
tds.addAttribute(DataServer.TUniformDistribution("tl", 63.1, 116.0))
tds.addAttribute(DataServer.TUniformDistribution("hu", 990.0, 1110.0))
tds.addAttribute(DataServer.TUniformDistribution("hl", 700.0, 820.0))
tds.addAttribute(DataServer.TUniformDistribution("l", 1120.0, 1680.0))
tds.addAttribute(DataServer.TUniformDistribution("kw", 9855.0, 12045.0))
ns = 100000
tsobol = Sensitivity.TSobol(tds, "flowrateModel", ns,
"rw:r:tu:tl:hu:hl:l:kw",
"flowrateModel", "DummyPython")
tsobol.loadOtherSobolFile("ref_sobol_launching_.dat")
tsobol.setDrawProgressBar(False)
tsobol.computeIndexes()
myfilter = "Algo==\"--first--\" || Algo==\"--total--\""
tsobol.getResultTuple().Scan("*", myfilter)
cc = ROOT.TCanvas("c1", "histgramme", 5, 64, 1270, 667)
pad = ROOT.TPad("pad", "pad", 0, 0.03, 1, 1)
pad.Draw()
pad.Divide(2, 1)
pad.cd(1)
tsobol.drawIndexes("Flowrate", "", "nonewcanv, hist, all")
pad.cd(2)
tsobol.drawIndexes("Flowrate", "", "nonewcanv, pie, first")
As for Macro “sensitivitySobolFunctionFlowrate.py”, the UserFunctions.C file is loaded and all the input
stochastic distribution are chosen wisely. The number of new estimations is set (to 10000, doubling
the statistic as this is also the number used in Macro “sensitivitySobolFunctionFlowrate.py” which provided
the input file that would be used to complete this estimation).
All the code lines are exactly those taken out of Macro Uranie but the one that is used to load a previous estimation which is shown below:
tsobol.loadOtherSobolFile("ref_sobol_launching_.dat")
Once done, the computeIndexes() method is called and few lines are shown to display the results in
the classical plot form, leading to the figure below (it is provided to illustrate how to represent
results). The numerical results are shown in the console below and the improvement in terms of
statistical precision can be seen by comparing the 95 percent confidence intervals going from
Console to Console.
13.6.15.3. Graph
The results of the previous macro is shown in Figure 13.31, where the evolution of the sobol coefficient is shown for all inputs with the uncertainty band, for the first order coefficient and the total one, respectively on the top and bottom panel.
Figure 13.31 Graph of the macro “sensitivitySobolLoadFile.py”
13.6.15.4. Console
--- Uranie v4.11/0 --- Developed with ROOT (6.36.06)
Copyright (C) 2013-2026 CEA/DES
Contact: support-uranie@cea.fr
Date: Thu Feb 12, 2026
<URANIE::WARNING>
<URANIE::WARNING> *** URANIE WARNING ***
<URANIE::WARNING> *** File[${SOURCEDIR}/dataSERVER/souRCE/TDataServer.cxx] Line[760]
<URANIE::WARNING> TDataServer::fileDataRead: Expected iterator tdsflowreate_1__n__iter__ not found but tdsflowreate__n__iter__ looks like an URANIE iterator => Will be used as so.
<URANIE::WARNING> *** END of URANIE WARNING ***
<URANIE::WARNING>
<URANIE::WARNING>
<URANIE::WARNING> *** URANIE WARNING ***
<URANIE::WARNING> *** File[${SOURCEDIR}/dataSERVER/souRCE/TDataServer.cxx] Line[8531]
<URANIE::WARNING> TDataServer::getTuple Error : There is no tree!
<URANIE::WARNING> *** END of URANIE WARNING ***
<URANIE::WARNING>
<URANIE::INFO>
<URANIE::INFO> *** URANIE INFORMATION ***
<URANIE::INFO> *** File[${SOURCEDIR}/meTIER/sampler/souRCE/TSamplerStochastic.cxx] Line[66]
<URANIE::INFO> TSamplerStochastic::init: the TDS [tdsflowreate] contains data: we need to empty it !
<URANIE::INFO> *** END of URANIE INFORMATION ***
<URANIE::INFO>
** Case of Output atty [flowrateModel] nSimPerIndex 20000
** Input att [rw] First [0.823516] Total Order[0.873849]
** Input att [r] First [0] Total Order[5.22438e-05]
** Input att [tu] First [0] Total Order[5e-05]
** Input att [tl] First [0] Total Order[6.08391e-05]
** Input att [hu] First [0.0319889] Total Order[0.0547724]
** Input att [hl] First [0.0304453] Total Order[0.0540884]
** Input att [l] First [0.0315935] Total Order[0.0527791]
** Input att [kw] First [0] Total Order[0.0129065]
************************************************************************************************************
* Row * Out.Out * Inp.Inp * Order.Ord * Method.Me * Algo.Algo * Value.Val * CILower.C * CIUpper.C *
************************************************************************************************************
* 0 * flowrateM * rw * First * Sobol * --first-- * 0.8235162 * 0.8190045 * 0.8279262 *
* 4 * flowrateM * rw * Total * Sobol * --total-- * 0.8738492 * 0.8602341 * 0.8875120 *
* 8 * flowrateM * r * First * Sobol * --first-- * 0 * 0 * 0.0138594 *
* 12 * flowrateM * r * Total * Sobol * --total-- * 5.224e-05 * 5.081e-05 * 5.371e-05 *
* 16 * flowrateM * tu * First * Sobol * --first-- * 0 * 0 * 0.0138594 *
* 20 * flowrateM * tu * Total * Sobol * --total-- * 5.000e-05 * 4.863e-05 * 5.140e-05 *
* 24 * flowrateM * tl * First * Sobol * --first-- * 0 * 0 * 0.0138594 *
* 28 * flowrateM * tl * Total * Sobol * --total-- * 6.083e-05 * 5.917e-05 * 6.254e-05 *
* 32 * flowrateM * hu * First * Sobol * --first-- * 0.0319888 * 0.0181374 * 0.0458280 *
* 36 * flowrateM * hu * Total * Sobol * --total-- * 0.0547723 * 0.0533147 * 0.0562686 *
* 40 * flowrateM * hl * First * Sobol * --first-- * 0.0304453 * 0.0165929 * 0.0442861 *
* 44 * flowrateM * hl * Total * Sobol * --total-- * 0.0540884 * 0.0526485 * 0.0555665 *
* 48 * flowrateM * l * First * Sobol * --first-- * 0.0315935 * 0.0177418 * 0.0454330 *
* 52 * flowrateM * l * Total * Sobol * --total-- * 0.0527791 * 0.0513732 * 0.0542224 *
* 56 * flowrateM * kw * First * Sobol * --first-- * 0 * 0 * 0.0138594 *
* 60 * flowrateM * kw * Total * Sobol * --total-- * 0.0129065 * 0.0125558 * 0.0132669 *
* 64 * flowrateM * __sum__ * First * Sobol * --first-- * 0.9175440 * -1 * -1 *
* 65 * flowrateM * __sum__ * Total * Sobol * --total-- * 1.0485587 * -1 * -1 *
************************************************************************************************************
==> 18 selected entries