13.6.8. Macro “modelerFlowratePolynChaosIntegration.C”

13.6.8.1. Objective

The objective of this macro is to build a polynomial chaos expansion in order to get a surrogate model along with a global sensitivity interpretation for the flowrate function, whose purpose and behaviour have been already introduced in Presentation of the problem. The method used here is the regression one, as discussed below.

13.6.8.2. Macro Uranie

    // Define the DataServer
    TDataServer *tds = new TDataServer("tdsFlowrate", "Design of experiments for Flowrate");

    gROOT->LoadMacro("UserFunctions.C");

    // Add the eight attributes of the study with uniform law
    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));

    // Define of TNisp object
    TNisp *nisp = new TNisp(tds);
    nisp->generateSample("Petras", 5); //State that there is a sample ...

    // Compute the output variable
    TLauncherFunction * tlf = new TLauncherFunction(tds, "flowrateModel","*","ymod");
    tlf->setDrawProgressBar(kFALSE);
    tlf->run();

    // build a chaos polynomial
    TPolynomialChaos *pc = new TPolynomialChaos(tds,nisp);

    // compute the pc coefficients using the "Integration" method
    Int_t degree = 4;
    pc->setDegree(degree);
    pc->computeChaosExpansion("Integration");

    // Uncertainty and sensitivity analysis
    cout << "Variable ymod ================" << endl;
    cout << "Mean     = " << pc->getMean("ymod") << endl;
    cout << "Variance = " << pc->getVariance("ymod") << endl;
    cout << "First Order Indices ================" << endl;
    cout << "Indice First Order[1] = " <<  pc->getIndexFirstOrder(0,0) << endl;
    cout << "Indice First Order[2] = " <<  pc->getIndexFirstOrder(1,0) << endl;
    cout << "Indice First Order[3] = " <<  pc->getIndexFirstOrder(2,0) << endl;
    cout << "Indice First Order[4] = " <<  pc->getIndexFirstOrder(3,0) << endl;
    cout << "Indice First Order[5] = " <<  pc->getIndexFirstOrder(4,0) << endl;
    cout << "Indice First Order[6] = " <<  pc->getIndexFirstOrder(5,0) << endl;
    cout << "Indice First Order[7] = " <<  pc->getIndexFirstOrder(6,0) << endl;
    cout << "Indice First Order[8] = " <<  pc->getIndexFirstOrder(7,0) << endl;
    cout << "Total Order Indices ================" << endl;
    cout << "Indice Total Order[1] = " << pc->getIndexTotalOrder("rw","ymod") << endl;
    cout << "Indice Total Order[2] = " << pc->getIndexTotalOrder("r","ymod") << endl;
    cout << "Indice Total Order[3] = " << pc->getIndexTotalOrder("tu","ymod") << endl;
    cout << "Indice Total Order[4] = " << pc->getIndexTotalOrder("tl","ymod") << endl;
    cout << "Indice Total Order[5] = " << pc->getIndexTotalOrder("hu","ymod") << endl;
    cout << "Indice Total Order[6] = " << pc->getIndexTotalOrder("hl","ymod") << endl;
    cout << "Indice Total Order[7] = " << pc->getIndexTotalOrder("l","ymod") << endl;
    cout << "Indice Total Order[8] = " << pc->getIndexTotalOrder("kw","ymod") << endl;

    // Dump main factors up to a certain threshold
    Double_t seuil = 0.99;
    cout<<"Ordered functionnal ANOVA"<<endl;
    pc->getAnovaOrdered(seuil,0);

    cout << "Number of experiments = " << tds->getNPatterns() << endl;

    //save the pv in a program (C langage)
    pc->exportFunction("NispFlowrate","NispFlowrate");

The first part is just creating a TDataServer and providing the attributes needed to define the problem. From there, a TNisp object is created, providing the dataserver that specifies the inputs. This class is used to generate the sample (Petras being a design-of-experiments dedicated to integration problem).

    // Define of TNisp object
    TNisp *nisp = new TNisp(tds);
    nisp->generateSample("Petras", 5); //State that there is a sample ...

The function is launched through a TLauncherFunction instance in order to get the output values that will be needed to train the surrogate model.

    TLauncherFunction * tlf = new TLauncherFunction(tds, "flowrateModel","*","ymod");
    tlf->run();

Finally, a TPolynomialChaos instance is created and the computation of the coefficients is performed by requesting a truncature on the resulting degree of the polynomial expansion (set to 4) and the use of a regression method.

    // build a chaos polynomial
    TPolynomialChaos *pc = new TPolynomialChaos(tds,nisp);
    // compute the pc coefficients using the "Integration" method
    Int_t degree = 4;
    pc->setDegree(degree);
    pc->computeChaosExpansion("Integration");

The rest of the code is showing how to access the resulting sensitivity indices either one-by-one, or ordered up to a chosen threshold of the output variance.

13.6.8.3. 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

Variable ymod ================
Mean     = 77.6511
Variance = 2078.91
First Order Indices ================
Indice First Order[1] = 0.828922
Indice First Order[2] = 1.04477e-06
Indice First Order[3] = 8.94692e-12
Indice First Order[4] = 5.30299e-06
Indice First Order[5] = 0.0413852
Indice First Order[6] = 0.0413852
Indice First Order[7] = 0.0393419
Indice First Order[8] = 0.00952157
Total Order Indices ================
Indice Total Order[1] = 0.866834
Indice Total Order[2] = 2.16178e-06
Indice Total Order[3] = 2.28296e-09
Indice Total Order[4] = 1.07386e-05
Indice Total Order[5] = 0.0541095
Indice Total Order[6] = 0.0541095
Indice Total Order[7] = 0.0520767
Indice Total Order[8] = 0.0127306