13.5.9. Macro “sensitivitySobolFunctionFlowrate.C”
13.5.9.1. Objective
The objective of this macro is to perform Sobol sensitivity analysis on a set of eight parameters used
in the flowrateModel model described in Presentation of the problem.
13.5.9.2. Macro Uranie
// Define the DataServer
TDataServer *tds = new TDataServer("tdsflowreate", "DataBase flowreate");
gROOT->LoadMacro("UserFunctions.C");
tds->addAttribute( new TUniformDistribution("rw", 0.05, 0.15));
tds->addAttribute( new TUniformDistribution("r", 100.0, 50000.0));
tds->addAttribute( new TUniformDistribution("tu", 63070.0, 115600.0));
tds->addAttribute( new TUniformDistribution("tl", 63.1, 116.0));
tds->addAttribute( new TUniformDistribution("hu", 990.0, 1110.0));
tds->addAttribute( new TUniformDistribution("hl", 700.0, 820.0));
tds->addAttribute( new TUniformDistribution("l", 1120.0, 1680.0));
tds->addAttribute( new TUniformDistribution("kw", 9855.0, 12045.0));
Int_t ns = 100000;
TSobol * tsobol = new TSobol(tds, "flowrateModel", ns, "rw:r:tu:tl:hu:hl:l:kw", "flowrateModel", "pouet");
tsobol->setDrawProgressBar(kFALSE);
tsobol->computeIndexes();
tsobol->getResultTuple()->Scan("*","Algo==\"--first--\" || Algo==\"--total--\"");
TCanvas *cc = new TCanvas("c1", "histgramme",5,64,1270,667);
TPad *pad = new 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");
gSystem->Rename("_sobol_launching_.dat","ref_sobol_launching_.dat");
tds->exportData("_onlyMandN_sobol_launching_.dat","rw:r:tu:tl:hu:hl:l:kw:flowrateModel","sobol__n__iter__tdsflowreate < 100");
The function flowrateModel is loaded from the macro UserFunctions.C (the file can be found in
${URANIESYS}/share/uranie/macros)
gROOT->LoadMacro("UserFunctions.C");
Each parameter is related to the TDataServer as a TAttribute and obeys an uniform law on specific interval:
TDataServer *tds = new TDataServer("tdsflowreate", "DataBase flowreate");
tds->addAttribute( new TUniformDistribution("rw", 0.05, 0.15));
tds->addAttribute( new TUniformDistribution("r", 100.0, 50000.0));
tds->addAttribute( new TUniformDistribution("tu", 63070.0, 115600.0));
tds->addAttribute( new TUniformDistribution("tl", 63.1, 116.0));
tds->addAttribute( new TUniformDistribution("hu", 990.0, 1110.0));
tds->addAttribute( new TUniformDistribution("hl", 700.0, 820.0));
tds->addAttribute( new TUniformDistribution("l", 1120.0, 1680.0));
tds->addAttribute( new TUniformDistribution("kw", 9855.0, 12045.0));
To instantiate the TSobol, one uses the TDataServer, the name of the function and the number of
samplings needed to perform sensitivity analysis (here ns=600):
TSobol * tsobol = new TSobol(tds, "flowrateModel", ns, "rw:r:tu:tl:hu:hl:l:kw", "flowrateModel", "pouet");
Computation of the sensitivity indexes:
tsobol->computeIndexes();
The automatic backup of data (the file _sobol_launching_.dat) is renamed so that it can be used in
other macros (see Macro “sensitivitySobolRe-estimation.C” and
Macro “sensitivitySobolLoadFile.C”) while the tds contains is exported (only the \(M\)
and \(N\) matrices content) to be also used in another macro (Macro “sensitivitySobolWithData.C”) :
gSystem->Rename("_sobol_launching_.dat","ref_sobol_launching_.dat");
tds->exportData("_onlyMandN_sobol_launching_.dat","rw:r:tu:tl:hu:hl:l:kw:flowrateModel","sobol__n__iter__tdsflowreate < 100");
13.5.9.3. Graph
Figure 13.25 Graph of the macro “sensitivitySobolFunctionFlowrate.C”
13.5.9.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[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 10000
** Input att [rw] First [0.830033] Total Order[0.865762]
** Input att [r] First [0] Total Order[0.000102212]
** Input att [tu] First [0] Total Order[0.0001]
** Input att [tl] First [0] Total Order[0.000110756]
** Input att [hu] First [0.0417298] Total Order[0.0554922]
** Input att [hl] First [0.0367345] Total Order[0.0526188]
** Input att [l] First [0.0384214] Total Order[0.0535728]
** Input att [kw] First [0.00669831] Total Order[0.0132316]
************************************************************************************************************
* Row * Out.Out * Inp.Inp * Order.Ord * Method.Me * Algo.Algo * Value.Val * CILower.C * CIUpper.C *
************************************************************************************************************
* 0 * flowrateM * rw * First * Sobol * --first-- * 0.8300331 * 0.8238356 * 0.8360321 *
* 4 * flowrateM * rw * Total * Sobol * --total-- * 0.8657619 * 0.8465652 * 0.8850599 *
* 8 * flowrateM * r * First * Sobol * --first-- * 0 * 0 * 0.0196004 *
* 12 * flowrateM * r * Total * Sobol * --total-- * 0.0001022 * 9.828e-05 * 0.0001062 *
* 16 * flowrateM * tu * First * Sobol * --first-- * 0 * 0 * 0.0196004 *
* 20 * flowrateM * tu * Total * Sobol * --total-- * 0.0001000 * 9.615e-05 * 0.0001039 *
* 24 * flowrateM * tl * First * Sobol * --first-- * 0 * 0 * 0.0196004 *
* 28 * flowrateM * tl * Total * Sobol * --total-- * 0.0001107 * 0.0001064 * 0.0001151 *
* 32 * flowrateM * hu * First * Sobol * --first-- * 0.0417297 * 0.0221474 * 0.0612800 *
* 36 * flowrateM * hu * Total * Sobol * --total-- * 0.0554921 * 0.0534156 * 0.0576470 *
* 40 * flowrateM * hl * First * Sobol * --first-- * 0.0367345 * 0.0171464 * 0.0562944 *
* 44 * flowrateM * hl * Total * Sobol * --total-- * 0.0526187 * 0.0506469 * 0.0546652 *
* 48 * flowrateM * l * First * Sobol * --first-- * 0.0384214 * 0.0188352 * 0.0579782 *
* 52 * flowrateM * l * Total * Sobol * --total-- * 0.0535727 * 0.0515661 * 0.0556552 *
* 56 * flowrateM * kw * First * Sobol * --first-- * 0.0066983 * 0 * 0.0262952 *
* 60 * flowrateM * kw * Total * Sobol * --total-- * 0.0132315 * 0.0127261 * 0.0137570 *
* 64 * flowrateM * __sum__ * First * Sobol * --first-- * 0.9536171 * -1 * -1 *
* 65 * flowrateM * __sum__ * Total * Sobol * --total-- * 1.0409902 * -1 * -1 *
************************************************************************************************************
==> 18 selected entries