13.6.7. Macro “modelerFlowratePolynChaosRegression.C”

13.6.7.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.7.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("QmcSobol", 500); //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 "Regression" method
    Int_t degree = 4;
    pc->setDegree(degree);
    pc->computeChaosExpansion("Regression");

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

    // Define of TNisp object
    TNisp *nisp = new TNisp(tds);
    nisp->generateSample("QmcSobol", 500); //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 "Regression" method
    Int_t degree = 4;
    pc->setDegree(degree);
    pc->computeChaosExpansion("Regression");

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.7.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.651
Variance = 2079.68
First Order Indices ================
Indice First Order[1] = 0.828732
Indice First Order[2] = 1.05548e-06
Indice First Order[3] = 9.86805e-08
Indice First Order[4] = 4.77668e-06
Indice First Order[5] = 0.0414015
Indice First Order[6] = 0.0413215
Indice First Order[7] = 0.039401
Indice First Order[8] = 0.00954644
Total Order Indices ================
Indice Total Order[1] = 0.866773
Indice Total Order[2] = 9.27467e-06
Indice Total Order[3] = 7.01768e-06
Indice Total Order[4] = 1.71944e-05
Indice Total Order[5] = 0.0541916
Indice Total Order[6] = 0.0540578
Indice Total Order[7] = 0.0522121
Indice Total Order[8] = 0.0127809