NISP Toolbox Manual - Forge
Contents
1. 0 0 2 0 EN 2 2 2 3 2 4 2 5 2 6 2 7 2 8 2 9 3 0 Figure 5 6 The histogram of a Uniform random variable with 1000 samples 39 Empirical histogram Log Uniform variable 0 18 0 12 7 Figure 5 7 The histogram of a Log Uniform random variable with 1000 samples 40 Chapter 6 The setrandvar class In this chapter we presen the setrandvar class The first section gives a brief outline of the features of this class and the second section present several examples 6 1 Introduction The setrandvar class allows to manage a collection of random variables and to build a Design Of Experiments DOE Several types of DOE are provided e Monte Carlo e Latin Hypercube Sampling e Smolyak Once a DOE is created we can retrieve the information experiment by experiment or the whole matrix of experiments This last feature allows to benefit from the fact that Scilab can natively manage matrices so that we do not have to perform loops to manage the complete DOE Hence good performances can be observed even if the language still i2. Figure 6 10 Latin Hypercube Sampling Uniform random variable the LHS sampling is not a Monte Carlo sampling randvar destroy vu1 randvar destroy vu2 setrandvar destroy srv 6 2 4 Other types of DOEs The following Scilab session allows to generate a Monte Carlo sampling with two uniform variables in the interval 1 1 The figure 6 11 presents this sampling and the figures 6 12 and 6 13 present the histograms of the two uniform random variables vul randvar new Uniforme 1 0 1 0 vu2 randvar new Uniforme 1 0 1 0 Srv setrandvar new setrandvar addrandvar srv vul setrandvar_addrandvar srv vu2 setrandvar_buildsample srv MonteCarlo 1000 sampling setrandvar getsample srv randvar_destroy vul randvar_destroy vu2 setrandvar_destroy srv It is easy to change the type of sampling by modifying the second argument of the setrandvar_buildsa function This way we can create the Petras Quadrature and Sobol sampling presented in figures 6 14 6 15 and 6 16 51 Monte Carlo Sampling Variable 22 Uniforme 1 0 1 0 1 0 0 8 0 6 0 4 0 2 0 0 0 2 0 4 0 6 0 8 1 0 Variable 1 Uniforme 1 0 1 0 Figure 6 11 Monte Carlo Sampling Two uniform variables in the interval 1 1 Variable 1 Uniforme 1 0 1 0 0 84 0 77 0 67
3. scilE 2c Usecase 1 Usecase 2 Randvar 1 Randvar 2 Setrandvar 1 Setrandvar Lhs design Setrandvar Monte Carlo samplin Figure 1 4 Demonstrations provided with the NISP toolbox Chapter 2 Theory In this chapter we extremely briefly present the theory which is used in the library This section is a tutorial introduction to the NISP module 2 1 Sensitivity analysis In this section we present the sensitivity analysis and emphasize the difference between global and local analysis Consider the model Y f X 2 1 where X Dx C R is the input and Y Dy CR is the output of the model The mapping f is presented in figure 2 1 Dx Dy Figure 2 1 Global analysis The assume that the input X is a random variable so that the output variable Y is also a random variable We are interested in measuring the sensitivity of the output depending on the uncertainty of the input More precisely we are interested in knowing e the input variables X which generate the most variability in the output Y e the input variables X which are not significant e a sub space of the input variables where the variability is maximum e if input variables interacts Consider the mapping presented in figure 2 1 The f mapping transforms the domain Dx into the domain Dy If f is sufficiently smooth small perturbations of X generate small perturbations of Y The local sensitivity analysis focu
4. 1 my handle children title text Variable 1 Normale 1 0 0 5 Plot Var 2 my_handle scf clf my handle reset histplot 50 sampling 2 my handle children title text Variable 2 Uniforme 2 0 3 0 Variable 1 Normale 1 0 0 5 Figure 6 9 Latin Hypercube Sampling Normal random variable We can use the mean and variance on each random variable and check that the expected result is computed We insist on the fact that the mean and variance functions are not provided by the NISP library these are pre defined functions which are available in the Scilab library That means that any Scilab function can be now used to process the data generated by the toolbox for ivar 1 2 m mean sampling ivar mprintf Variable d Mean 4f n ivar m v variance sampling ivar mprintf Variable 4d Variance Af n ivar v end The previous script produces the following output Variable 1 Mean 1 000000 Variable 1 Variance 0 249925 Variable 2 Mean 2 500000 Variable 2 Variance 0 083417 Our numerical simulation is now finished but we must destroy the objects so that the memory managed by the toolbox is deleted 50 Variable 2 Uniforme 2 0 3 0
5. Mean 3 500000 Variance 13 842473 First order sensitivity index Variable X1 0 313953 Variable X2 0 442325 Variable X3 0 000000 Total sensitivity index Variable X1 0 557675 62 Variable X2 0 442326 Variable X3 0 243721 We now focus on the variance generated by the variables 1 and 3 We set the group to the empty group with the polychaos setgroupempty function and add variables with the polychaos_setgroupaddvar function groupe 1 3 polychaos setgroupempty pc polychaos setgroupaddvar pc groupe 1 polychaos setgroupaddvar pc groupe 2 mprintf Fraction of the variance of ja group of variablesAn mprintf uuuuGroupeuXiwetuX2u 4fAMn polychaos_getgroupind pc The previous script produces the following output Fraction of the variance of a group of variables Groupe X1 et X2 0 557674 The function polychaos getanova prints the functionnal decomposition of the normalized variance polychaos_getanova pc The previous script produces the following output 00 0 313953 0 442325 1 55229e 009 8 08643e 031 0 243721 7 26213e 031 1 6007e 007 k O HO HO FA k OO kk KA k k ki OO We can compute the density function associated with the output variable of the function In order to compute it we use the polychaos buildsample function and create a Latin Hypercube Sampling with 10000 experiments The polychaos getsample function allows to quiery the polynomial and get th
6. 0 03538 SRC_2 0 12570 SRC_3 0 31817 SRC_4 0 54314 expected 0 expected 0 expected 0 expected 0 03333 13333 30000 53333 SUM 1 00066 expected 1 00000 The previous script also produces the plots 2 2 2 3 2 4 and 2 5 2 5 Sensitivity analysis for nonlinear models Let us focus on one particular input X of the model f with 1 2 particular value say x for example then the variance of the output Y surely decreases because the variable X is not randome anymore We can then measure the conditionnal variance given Xi denoted by V Y X z Since X is a random variable the conditionnal variance V Y X x is a random variable We are interested in the average value of this variance that is we are interested in E V Y X If X has a large weight in the variance V Y then E V Y X is small The theorem of the total variance states that if V Y is finite then V Y V E Y X E V Y X If X has a large weight in the variance V Y then V E Y X is large 14 p If we set X toa Scatter plot for X1 25 Figure 2 2 Scatter plot for an affine model Variable X4 Scatter plot for X2 25 Figure 2 3 Scatter plot for an affine model Variable X5 15 Scatter plot for X3 15 20 15 10 5 0 5 10 15 20 25 Figure 2 4 Scatter plot for an affine model Variable X3 Scatter plot for X4 25 Figure 2 5
7. 0 5 044 0 3 0 2 0 1 0 0 Figure 6 12 Latin Hypercube Sampling First uniform variable in 1 1 52 Variable 2 Uniforme 1 0 1 0 0 8 074 0 67 0 5 1 0 0 8 0 6 0 4 0 2 0 0 0 2 0 4 0 6 0 8 1 0 Figure 6 13 Latin Hypercube Sampling Second uniform variable in 1 1 Petras Sampling 731 371 1 3 TTT 0 7 7 IDE E do H D e o 9 pe H Ho E 5 gt ist i d d 4d 2 8 2 ti eg AG did id i fr a_i ii i L 0 2 0 3 0 4 0 5 0 6 0 7 0 8 0 Variable 1 Uniforme 0 0 1 0 Figure 6 14 Petras sampling Two uniform variables in the interval 1 1 53 Quadrature Sampling 1 0 TT 3 T T T T T T FF rf T Er e e N N e o ee ese o e D e o e e ee 0 9 ho 0 e e e e D e e o o ee kee D D D CJ D D D D D D D D D ee PRES ee ley n K e o e e ee po e o k Ta Se Q9 kes e e e e e o e e ee e S 064 g
8. 44 AA 8 31 5 1 2 Parameters of the Log normal distribution 32 5 1 3 Uniform random number generation 32 a RAR CR AY A 33 Bl CABE ciu o rr AA AS AS A II 33 5 22 The Oriented Object system ss o socs eoc se ORO EN ee A se 33 Do Examples HT 35 Bol Amen session cusa EE ed XO SOR Rw OR we aod BR dot d 36 Soe Variable translormalions lt il RAR 36 6 The setrandvar class A1 Bl ba c x e Le lt cu uud Y O RDA AER de RP ds ADRES RS AER 41 E Tee pce i r a a 41 6 2 1 A Monte Carlo design with 2 variables en 6 2 0 A Monte Carlo design with 2 variables DAS AH UON ao bon Geom at aa atta Dad Otheriypesol DOES o ue e Peed o bed YM Re RU ne La 7 The polychaos class GEI o pie A DE A Dio soc IA A RA ee ee A 7 2 1 Product of two random variables 122 The pene BECA a AAA AAA 8 Thanks Bibliography Chapter 1 Introduction 1 1 The OPUS project The goal of this toolbox is to provide a tool to manage uncertainties in simulated models This toolbox is based on the NISP library where NISP stands for Non Intrusive Spectral Projection This work has been realized in the context of the OPUS project http opus project fr Open Source Platform for Uncertainty treatments in Simulation funded by ANR the french Agence Nationale pour la Recherche http www agence nationale recherche fr The toolbox is released under the Lesser General Public Licenc
9. 5 Monte Carlo Sampling Histogram and exact distribution functions for the second variable function p exppdf x lambda p lambda exp lambda x endfunction The following script creates a Monte Carlo sampling where the first variable is Normal and the second variable is Exponential Then we compare the empirical histogram and the exact distribution function We use the dnorm function defined in the Stixbox module _initseed 0 rvi randvar_new Normale 1 0 0 5 rv2 randvar new Exponentielle 5 Definition d un groupe de variables aleatoires Srv setrandvar new Setrandvar addrandvar srv rvi setrandvar_addrandvar srv rv2 nisp np 1000 setrandvar_buildsample srv MonteCarlo np sampling setrandvar_getsample srv Check sampling of random variable 1 mean sampling 1 variance sampling 1 Check sampling of random variable 2 min sampling 2 max sampling 2 Plot scf histplot 40 sampling 1 x linspace 1 3 1000 p dnorm x 1 0 5 plot x p r 46 xtitle Empirical histogram of X1 X P X legend Empirical Exact scf histplot 40 sampling 2 x linspace 0 2 1000 p exppdf x 5 plot x p r xtitle Empiricaluhistogramiof X2 X P X legend Empirical Exact Clean up Setrandvar destroy srv randvar destroy rv1 randvar destroy rv2 The previous sc
10. Vy expectation mprintf uuuuVariable X2u 14 fu expectationy 4f n S2 S2_expectation re abs S2 S2_expectation S2_expectation mprintf QuuuuuuuRelativeyErrory f n re mprintf Total sensitivity index n STi polychaos_getindextotal pc 1 mprintf j j ji Variable X1 7 4 f Nn ST1 ST2 polychaos_getindextotal pc 2 mprintf j j j Variable X2 7 4 f Nn 8T2 Clean up polychaos destroy pc randvar destroy vul randvar destroy vu2 randvar destroy vx1 randvar destroy vx2 setrandvar_destroy srvu setrandvar_destroy srvx The previous script produces the following output Mean 5 250000 expectation 5 250000 Variance 18 687500 expectation 18 687500 First order sensitivity index Variable X1 0 163880 expectation 0 163880 Relative Error 0 000000 Variable X2 0 752508 expectation Relative Error 0 000000 Total sensitivity index Variable X1 0 247492 Variable X2 0 836120 0 752508 We see that the polynomial chaos performs an exact computation 21 2 7 Sobol decomposition Sobol 4 introduced the sensitivity index based on V E Y X by decomposing the function f as a sum of function with an increasing number of parameters We consider the function f Y PE 2 p 2 69 where X z 25 0 1 P If f can be integrated in 0 1 P then there is a unique decom
11. computed from rv2 into a sample from rv The mean allows to check that the transformed samples have an mean value which corresponds to the random variable rv nisp initseed O a 10 0 b 20 0 rv randvar_new LogUniforme a b rv2 randvar new Uniforme 2 3 nbshots 1000 valuesLou zeros nbshots for i 1 nbshots valuesUni i randvar_getvalue rv2 valuesLou i randvar getvalue rv rv2 valuesUni i end computed mean valuesLou mu b a log b log a expected mu mprintf Expectation 4 5f expected 4 5f n computed expected scf histplot 50 valuesUni xtitle Empirical jhistogram Uniform variable X P X scf histplot 50 valuesLou xtitle Empirical histogram Log Uniform variable X P X randvar_destroy rv randvar destroy rv2 The previous script produces the following output Expectation 14 63075 expected 14 42695 The previous script also produces the figures 5 6 and 5 7 The transformation depends on the mother random variable rvi and on the daughter ran dom variable rv Specific transformations are provided for all many combinations of the two distribution functions T hese transformations will be analysed in the next sections 38 Empirical histogram Uniform variable x 4 amp 08 0 6 7 0 47
12. experiments setrandvar buildsample srv Lhs 1000 At this point the DOE is stored in the memory space of the NISP library but we do not have a direct access to it We now call the setrandvar_getsample function and store that DOE into the sampling matrix sampling setrandvar getsample srv The sampling matrix has 1000 rows corresponding to each experiment and 2 columns cor responding to each input random variable The following script allows to plot the sampling which is is presented in figure 6 8 my handle scf clf my handle reset plot sampling 1 sampling 2 my_handle children children children line_mode off my_handle children children children mark_mode on my_handle children children children mark_size 2 my_handle children title text Latin Hypercube Sampling my handle children x label text Variable 1 Normale 1 0 0 5 my _handle children y_label text Variable 2 Uniforme 2 0 3 0 The following script allows to plot the histogram of the two variables which are presented in figures 6 9 and 6 10 Plot Var 1 48 Latin Hypercube Sampling Variable 2 Uniforme 2 0 3 0 1 0 0 5 0 0 0 5 1 0 1 5 2 0 2 5 3 0 Variable 1 Normale 1 0 0 5 Figure 6 8 Latin Hypercube Sampling The first variable is Normal the second variable is Uniform 49 my handle scf clf my handle reset histplot 50 sampling
13. know the value of the SCIHOME global variable gt SCIHOME SCIHOME C Users baudin AppData Roaming Scilab scilab 5 2 0 beta 1 On my Linux system the Scilab 5 1 startup file is located in home myname Scilab scilab 5 1 scilab On my Windows system the Scilab 5 1 startup file is located in C Users myname AppData Roaming Scilab scilab 5 1 scilab This file is a regular Scilab script which is automatically loaded at Scilab s startup If that file does not already exist create it Copy the following lines into the scilab file and configure the path to the toolboxes stored in the SCILABTBX variable exec C tbxnisp loader sce The following script presents the messages which are generated when the unit tests script of the toolbox is launched exec C tbxnisp runtests sce Tests beginning the 2009 11 18 at 12 47 45 TMPDIR C Users baudin AppData Local Temp SCI_TMP_6372_ 001 004 tbxnisp nisp passed ref created 002 004 tbxnisp polychaosi passed ref created 003 004 tbxnisp randvar i passed ref created 004 004 tbxnisp setrandvari passed ref created Summary tests 4 100 passed 0 0 failed 0 o skipped 0 o length 3 84 sec Tests ending the 2009 11 18 at 12 47 48 end verbatim 29 Chapter 4 Configuration functions In this section we present functions which allow to configure the NISP toolbox The nisp_ functions allo
14. li TE 1 Figure 7 3 Product function Histogram of the output on a LHS design with 10000 experiments b 0 1 si sin x 1 s2 sin x 2 y 1 st a s2 s2 b x 3 x 3 x 3 x 3 s1 endfunction We create 3 uncertain parameters which are uniform in the interval 7 7 and put these random variables into the collection srvu rvul randvar new Uniforme pi pi rvu2 randvar_new Uniforme pi pi rvu3 randvar_new Uniforme pi pi srvu setrandvar_new setrandvar_addrandvar srvu rvul setrandvar_addrandvar srvu rvu2 setrandvar_addrandvar srvu rvu3 The collection of stochastic variables is created with the function setrandvar_new The calling sequence srvx setrandvar_new nx allows to create a collection of nx 3 random variables uniform in the interval 0 1 Then we create a Petras DOE for the stochastic collection srvx and transform it into a DOE for the uncertain parameters srvu nx setrandvar_getdimension srvu Srvx setrandvar new nx degre 9 setrandvar_buildsample srvx Petras degre setrandvar_buildsample srvu srvx 6l We use the polychaos new function and create the new polynomial pc with 3 inputs and 1 output noutput 1 pc polychaos new srvx noutput The next step allows to perform the simulations associated with the DOE prescribed by the collection srvu Here we must perform np 751 experiments np setrandvar
15. ne setrandvar addrandvar srv rv Destructor setrandvar destroy srv Static methods tokenmatrix setrandvar tokens nb setrandvar size srv name np srv k i value i Figure 6 1 Outline of the methods of the setrandvar class 42 variable is associated with a Uniform distribution function The simulation is based on 1000 experiments The function nisp initseed is used to set the value of the seed to zero so that the re sults can be reproduced The setrandvar new function is used to create a new set of ran dom variables Then we create two new random variables with the randvar new function These two variables are added to the set with the setrandvar_addrandvar function The setrandvar buildsample allows to build the design of experiments which can be retrieved as matrix with the setrandvar_getsample function The sampling matrix has np rows and 2 columns one for each input variable nisp_initseed 0 rvul randvar_new Normale 1 3 rvu2 randvar_new Uniforme 2 3 srvu setrandvar_new setrandvar_addrandvar srvu rvul setrandvar_addrandvar srvu rvu2 77 np 5000 setrandvar buildsample srvu MonteCarlo np sampling setrandvar getsample srvu Check sampling of random variable 41 mean sampling 1 Expectation 1 Check sampling of random variable 2 mean sampling 2 Expectation 2 5 ii scf histplot 5
16. normal distribution For the LogNormale law the distribution function is usually defined by the expected value y and the standard deviation c of the underlying Normal random variable But when we create a LogNormale randvar the parameters to pass to the constructor are the expected value of the LogNormal random variable E X and the standard deviation of the underlying Normale random variable o The expected value and the variance of the Log Normal law are given by E X exp 50 5 1 V X exp o 1 exp 2u 07 5 2 It is possible to invert these formulas in the situation where the given parameters are the expected value and the variance of the Log Normal random variable We can invert completely the previous equations and get u In E X In rk 5 3 g EE co l 1 sien e 5 4 In particular the expected value u of with the Normal random variable satisfies the equation p In E X oi 5 5 5 1 3 Uniform random number generation In this section we present the generation of uniform random numbers The goal of this section is to warn users about a current limitation of the library Indeed the random number generator is based on the compiler so that its quality cannot be guaranted The Uniforme law is associated with the parameters a b R with a lt b It produces real values uniform in the interval a b To compute the uniform random number X in the interval a b a uniform rand
17. po o o D D D D D D D D D D D oe Se ERC Yen tas Ue Ce Ce Cu een leg S poo o o D D D D D D D D D D D eee 047 2 B ooo e e o o e e ee S gt 0 30 0 0 e e c H e e e ee eens e e a e 0 oo 027 boss e o e D D e o e o ee ke s 0 e e o o e e ee of ee e D D D D D D D D D D D D D eee o D D e o o o o Mn i i ATI I 0 0 0 1 0 2 0 3 0 4 0 5 0 6 0 7 0 8 0 9 Variable 1 Uniforme 0 0 1 0 Figure 6 15 Quadrature sampling Two uniform variables in the interval 1 1 Sobol Sampling Variable 2 Uniforme 1 0 1 0 Variable 1 Uniforme 1 0 1 0 Figure 6 16 Sobol sampling Two uniform variables in the interval 1 1 54 Chapter 7 The polychaos class 7 1 Introduction The polychaos class allows to manage a polynomial chaos expansion The coefficients of the expansion are computed based on given numerical experiments which creates the association between the inputs and the outputs Once computed the expansion can be used as a regular function The mean standard deviation or quantile can also be directly retrieved The tool allows to get the following results e mean e variance e quantile e correlation etc Moreover we can generate the C source code which computes the output of the polynomial chaos expansion This C source code is stand alone that is it is independent of both the NISP library a
18. quick interactive overview of its features The following script presents the messages which are generated when the builder of the toolbox is launched The builder script performs the following steps e compile the NISP C library e compile the C gateway library the glue between the library and Scilab e generate the Java help files from the xml files e generate the loader script exec C tbxnisp builder sce Building sources Generate a loader file Generate a Makefile Running the Makefile Compilation of utils cpp Compilation of blasi_d cpp Compilation of dcdflib cpp Compilation of faure cpp Compilation of halton cpp Compilation of linpack_d cpp Compilation of niederreiter cpp Compilation of reversehalton cpp Compilation of sobol cpp Building shared library be patient Generate a cleaner file Generate a loader file Generate a Makefile Running the Makefile Compilation of nisp gc cpp Compilation of nisp gva cpp Compilation of nisp ind cpp Compilation of nisp index cpp Compilation of nisp inv cpp Compilation of nisp math cpp Compilation of nisp msg cpp Compilation of nisp_conf cpp Compilation of nisp ort cpp Compilation of nisp pc cpp Compilation of nisp_polyrule cpp 26 Compilation Compilation Compilation Compilation Compilation Compilation of of of of of of nisp qua cpp nisp random cpp nisp smo cpp nisp util cpp nisp_va cpp nisp_smolyak cpp Building shared library be patient Generate a
19. the figure 1 1 as a framework and will follow the steps in order In the following Scilab script we define the function Example which takes a vector of size 2 as input and returns a scalar as output function y Exemple x y 1 x 1 x 2 endfunction We now create a collection of two stochastic normalized random variables Since the ran dom variables are normalized we use the default parameters of the randvar new function The normalized collection is stored in the variable srvx vxl randvar new Normale vx2 randvar new Uniforme srvx setrandvar_new setrandvar_addrandvar srvx vx1 setrandvar_addrandvar srvx vx2 We create a collection of two uncertain parameters We explicitely set the parameters of each random variable that is the first Normal variable is associated with a mean equal to 1 0 and a standard deviation equal to 0 5 while the second Uniform variable is in the interval 1 0 2 5 This collection is stored in the variable srvu vul randvar new Normale 1 0 0 5 vu2 randvar new Uniforme 1 0 2 5 srvu setrandvar new setrandvar_addrandvar srvu vul setrandvar_addrandvar srvu vu2 56 Methods output polychaos gettarget pc np polychaos getsizetarget pc polychaos getsample pc k ovar polychaos getquantile pc k polychaos getsample pc polychaos getquantile pc alpha polychaos getoutput pc polychaos
20. vu2 randvar getvalue vu2 error 4 Undefined variable vu2 It is now impossible to generate values from the variable vu2 Moreover it may be difficult to know exactly what went wrong and what exact variable is lost At any time we can use the randvar_tokens function in order to get the list of current variables Deleting these variables allows to clean the memory space properly without memory loss gt randvar_tokens ans 3 4 gt randvar destroy 3 ans 3 gt randvar destroy 4 ans 4 randvar tokens ans L 5 3 Examples In this section we present to examples of use of the randvar class The first example presents the simulation of a Normal random variable and the generation of 1000 random variables The second example presents the transformation of a Uniform outcome into a LogUniform outcome 35 5 3 1 A sample session We present a sample Scilab session where the randvar class is used to generate samples from the Normale law In the following Scilab session we create a Normale random variable and compute samples from this law The nisp initseed function is used to initialize the seed for the uniform random variable generator Then we use the randvar new function to create a new random variable from the Normale law with mean 1 and standard deviation 0 5 The main loop allows to compute 1000 samples from this law based on calls to the randvar_getvalue function Once the samples are compu
21. 0 sampling 1 xtitle Empirical jhistogram of X1 scf histplot 50 sampling 2 xtitle Empirical histogram of X2 Clean up setrandvar_destroy srvu randvar destroy rvu1 randvar destroy rvu2 The previous script produces the following output gt mean sampling 1 Expectation 1 ans 1 0064346 mean sampling 2 Expectation 2 5 ans 2 5030984 The prevous script also produces the figures 6 2 and 6 3 We may now want to add the exact distribution to these histograms and compare The Normal distribution function is not provided by Scilab but is provided by the Stixbox module Indeed the dnorm function of the Stixbox module computes the Normal probability distribution function 43 Empirical histogram of X1 3 a o o 3 a Figure 6 2 Monte Carlo Sampling Normal random variable Empirical histogram of X2 0 0 2 0 21 2 2 23 2 4 25 2 6 27 2 8 29 3 0 Figure 6 3 Monte Carlo Sampling Uniform random variable 44 In order to install this module we can run the atomsInstall function as in the following script atoms
22. Install stixbox The following script compares the empirical and theoretical distributions scf histplot 50 sampling 1 xtitle Empirical histogram of X1 x linspace 15 15 1000 y dnorm x 1 3 plot x y r legend Empirical Exact The previous script produces the figure 6 4 E Empirical Exact nd lu 0 05 4 a 3 n od o 3 a Figure 6 4 Monte Carlo Sampling Histogram and exact distribution functions for the first variable The following script performs the same comparison for the second variable scf histplot 50 sampling 2 xtitle Empirical histogram of X2 x linspace 2 3 1000 y ones 1000 1 plot x y r The previous script produces the figure 6 5 6 2 2 A Monte Carlo design with 2 variables In this section we create a Monte Carlo design with 2 variables We are going to use the exponential distribution function which is not defined in Scilab The following exppdf function computes the probability distribution function of the exponential distribution function 45 Empirical histogram of X2 Figure 6
23. NISP Toolbox Manual Michael Baudin INRIA Jean Marc Martinez CEA Version 0 2 January 2011 Abstract This document is a brief introduction to the NISP module We present the installation process of the module in binary from ATOMS or from the sources We present the configuration functions and the randvar setrandvar and polychaos classes Contents 1 Introduction 3 LEI The OPUS DOSE 2 2 WE Eod EX Se NN SLT A ene wr 3 l2 The A cun d RG SMS RA ER RD ER RR ou e RA EE RA AU 3 lo The NIOF modules oes t sieto aa ok a a oa kue o kar oih kan aoh E a 4 2 Theory 9 21 PIO so uc EE mie RAR AR kas E A 9 2 2 Standardized regression coefficients of affine models 10 2 3 Link with the linear correlation coefficients 11 ZA Using scatter plots s 4 24 s es du ERR ge ra 12 2 5 Sensitivity analysis for nonlinear models 14 2 6 The effect of the interactions 17 Z1 poboldecomposkion 2 ssc s oe a A a a A 22 205 IRMA MACION uoo e OE E OE R E COE Ea A UR 22 3 Installation 23 ml MCR e er NR x e an Oe Bane be ACE OS bane DS E ace A 23 3 2 Installing the toolbox from ATOMS e sr ssa ou due o X xx x 24 3 3 Installing the toolbox from the sources 25 4 Configuration functions 30 5 The randvar class 31 5 1 The distribution functions 4 4 4 4 a e da d De au 31 SI EDIT le ana EE ary ee d pos Xo 4 X04 Xo 4
24. NISP module The NISP toolbox is available under the following operating systems e Linux 32 bits e Linux 64 bits e Windows 32 bits e Mac OS X The following list presents the features provided by the NISP toolbox 4 e Manage various types of random variables uniform normal exponential log normal Generate random numbers from a given random variable Transform an outcome from a given random variable into another e Manage various Design of Experiments for sets of random variables Monte Carlo Sobol Latin Hypercube Sampling various samplings based on Smolyak designs Manage polynomial chaos expansion and get specific outputs including mean variance quantile correlation etc e Generate the C source code which computes the output of the polynomial chaos expansion This User s Manual completes the online help provided with the toolbox but does not replace it The goal of this document is to provide both a global overview of the toolbox and to give some details about its implementation The detailed calling sequence of each function is provided by the online help and will not be reproduced in this document The inline help is presented in the figure 1 2 For example in order to access to the help associated with the randvar class we type the following statements in the Scilab console help randvar The previous statements opens the Help Browser
25. Scatter plot for an affine model Variable X4 16 The first order sensivity indice of Y to the variable X is defined by s V E Y X5 V 2 30 fori 1 2 p The sensitivity indice measures the part of the variance which is caused by the uncertainty in X We can compute the sensitivity indice when the function f is linear Assume that the output Y depends linearily on the input X Y bo SS PX 2 31 i 1 2 p where 5 ER for 0 1 2 p Then E Y X Bo H E X GiXi 2 32 since the expection of a sum is the sum of expectations Then V E Y X V 6 X 2 33 fgV X 2 34 since the variance of a constant term is zero Therefore the sensitivity index of Y to the variable X 19 BV Xi S BM 2 35 for i 1 2 p Hence if we make the assumption that a model is affine then the empirical linear correlation coefficient can be used to estimate the sensitivity indices 2 6 The effect of the interactions In this example we consider a non linear non additive model made of the product of two inde pendent random variables The goal of this example is to show that in some cases we have to consider the interations between the variables Consider the function Y X1X2 2 36 where X and X gt are two independent normal random variables with mean 44 and uz and variances 2 2 of and 03 Let us compute the expectation of the random variable Y The expectation of Y
26. able is V Xi l V X2 4 2 22 V X3 29 V X 16 2 23 Since the variables are independent the variance of the output Y is V Y 2V X1 V X3 V X3 V X4 30 2 24 The standardized regression coefficient is BV X GRO a 2 25 V Y 30 225 for 1 2 3 4 More specifically we have 1 4 cn 2 26 SRC on SRC ze 2 26 9 16 SRC3 SRC4 2 27 3 30 4 30 We have the following inequalities S HRC gt SRC3 gt SRC gt S RC 2 28 This means that the variable which causes the most variance in the output is X4 while the variable which causes the least variance in the output is X4 The script below performs the analysis with the NISP module The sampling is based on a Latin Hypercube Sampling design with 5000 points 12 function y Exemple x y x 1 x 2 x 3 x 4 endfunction function r lincorrcoef x y Returns the linear correlation coefficient of x and y The variables are expected to be column matrices with the same size x x y y mx mean x my mean y sx sqrt sum x mx 2 sy sqrt sum y my 72 r x mx y my sx sy endfunction Initialisation de la graine aleatoire nisp initseed O Create the random variables rvui randvar new Normale 0 1 rvu2 randvar new Normale 0 2 rvu3 randvar new Normale 0 3 rvu4 randvar_new Normale 0 4 srvu setrandvar_new setrandvar_addran
27. and displays the helps page presented in figure Several demonstration scripts are provided with the toolbox and are presented in the figure 1 4 These demonstrations are available under the question mark in the menu of the Scilab console Finally the unit tests provided with the toolbox cover all the features of the toolbox When we want to know how to use a particular feature and do not find the information we can search in the unit tests which often provide the answer Help Browser NISP Toolbox Table of Contents nisp Functions to configure the NISP library overview An overview of the NISP toolbox polychaos A class to manage a Polynomial Chaos expansion randvar A class to manage a random variable setrandvar A class to manage a set of random variables Figure 1 2 The NISP inline help Help Browser Nom randvar class to manage a random variable SYNOPSIS tokens randvar tokens randvar size randvar new name a b randvar new name a randvar new name randvar destroy rv value randvar getvalue rv value randvar getvalue rv rv2 value2 randvar getlog rv Figure 1 3 The online help of the randvar function GUI Genetic Algorithms Simulated Annealing Graphics Signal Processing CACSD Optimization and Simulation Polynomials Simulation Scicos Metanet Tcl Tk Sound file handling Random Spreadsheet Nisp
28. ar new setrandvar_addrandvar srvu vul setrandvar_addrandvar srvu vu2 5 Create the Design 0f Experiments degre 2 setrandvar buildsample srvx Quadrature degre Setrandvar buildsample srvu srvx 6 Create the polynomial ny 1 pc polychaos new srvx ny np setrandvar getsize srvx mprintf Number jof experiments din np polychaos setsizetarget pc np T Perform the DOE inputdata setrandvar getsample srvu outputdata Exemple inputdata 20 polychaos settarget pc outputdata 8 Compute the coefficients of the polynomial expansion polychaos setdegree pc degre polychaos computeexp pc srvx Integration 9 Get the sensitivity indices average polychaos getmean pc Ey expectation mui mu2 var polychaos getvariance pc Vy expectation mu2 2 sigmai 2 mui 2 sigma2 2 sigmai 2 sigma2 2 mprintf Meany yuu uhf y expectation y f 1n average Ey_expectation mprintf Varianceyyuu sf yu expectation f n var Vy_expectation mprintf First order sensitivity index n S1 polychaos getindexfirst pc 1 S1 expectation mu272 sigmai172 Vy_expectation mprintf uuuVariable X1 f expectation 4f Nn 81 81 expectation re abs S1 S81 expectation S1 expectation mprintf GQuuuuuuuRelativeyErrory 4f n re 2 polychaos getindexfirst pc 2 82 expectation mui 2 sigma2 2
29. bers from a given distribution function 5 1 The distribution functions In this section we present the distribution functions provided by the randvar class We especially present the Log normal distribution function 5 1 1 Overview The table 5 1 gives the list of distribution functions which are available with the randvar class 3 Each distribution functions have zero one or two parameters One random variable can be specified by giving explicitely its parameters or by using default parameters The parameters for all distribution function are presented in the figure 5 2 which also presents the conditions which must be satisfied by the parameters Name DES EX V X Normale VE exp 1 ep LL c Uniforme Co i E i bta Bes Exponentielle SCH FAR y S i LogNormale ox am O P 2 ar p eed exp u 30 exp o 1 exp 2u e 0 ifa lt 0 LogUniforme TI e E i i MORBO soi E x Figure 5 1 Distributions functions of the randvar class The expected value is denoted by E X and the variance is denoted by V X 31 Name Parameter 1 a Parameter 2 b Conditions Normale p 0 o 1 Bion Uniforme a 0 b a b Exponentielle A 1 LogNormale w 0 1 o 1 0 nc LogUniforme a 0 1 b 1 0 a b gt 0 a lt b Figure 5 2 Default parameters for distributions functions 5 1 2 Parameters of the Log
30. cleaner file Building gateway Generate a gateway file Generate a loader file Generate a Makefile Running the Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Makelib makefile of of of of of of of of of of of of of of of of of of of of of of of of of of of of of of of of of of of of nisp_gettoken cpp nisp_gwsupport cpp nisp_PolynomialChaos_map cpp nisp_RandomVariable_map cpp nisp_SetRandomVariable_map cpp sci_nisp_startup cpp sci_nisp_shutdown cpp sci_nisp_verboselevelset cpp Sci nisp verboselevelget cpp Sci nisp initseed cpp Sci randvar new cpp Sci randvar destroy cpp Sci randvar size cpp Sci randvar tokens cpp Sci randvar getlog cpp Sci randvar getvalue cpp Sci setrandvar new cpp Sci setrandvar tokens cpp Sci setrandvar size cpp Sci setrandvar destroy cpp Sci setrandvar freememory cpp Sci setrandvar addrandvar cpp Sci setrandvar getlog cpp sci_setrandvar_getdimension cpp sci_setrandvar_getsize cpp sci_setrandvar_getsample cpp sci_setrandvar_setsample cpp Sci
31. des donn es Bruits et Signaux Master M2 Recherche Astronomie astrophysique 2006 4 I M Sobol Sensitivity estimates for nonlinear mathematical models Mathematical Modelling and Computational Experiments 1 407 414 2003 67
32. dvar srvu rvul setrandvar_addrandvar srvu rvu2 setrandvar_addrandvar srvu rvu3 setrandvar_addrandvar srvu rvu4 Create a sampling by a Latin Hypercube Sampling with size 5000 nbshots 5000 setrandvar_buildsample srvu Lhs nbshots sampling setrandvar_getsample srvu Perform the experiments y Exemple sampling Scatter plots y depending on X_i for k 1 4 scf plot y sampling k rx xistr X string k xtitle Scatteryplotyfory xistr xistr Y end Compute the sample linear correlation coefficients rhol lincorrcoef sampling 1 y D SRC1 rho172 SRCiexpected 1 30 mprintf SRC_1 5f expected 5f n SRC1 SRClexpected 77 rho2 lincorrcoef sampling 2 y 13 SRC2 rho272 SRC2expected 4 30 mprintf SRC_2 5f expected 5f n SRC2 SRC2expected 1 rho3 lincorrcoef sampling 3 y D SRC3 rho372 SRC3expected 9 30 mprintf SRC_3 5f expected 5f n SRC3 SRC3expected rho4 lincorrcoef sampling 4 y SRC4 rho472 SRC4expected 16 30 mprintf SRC_4 5f expected 5f n SRC4 SRC4expected 77 SUM SRC1 SRC2 SRC3 SRC4 SUMexpected 1 mprintf SsUM 5f expected 5f 1n SUM SUMexpected 77 Clean up randvar destroy rvu1 randvar_destroy rvu2 randvar_destroy rvu3 randvar destroy rvu4 setrandvar_destroy srvu The previous script produces the following output SRC_1
33. e LGPL as all components of the OPUS project 1 2 The NISP library The NISP library is based on a set of 3 C classes so that it provides an object oriented framework for uncertainty analysis The Scilab toolbox provides a pseudo object oriented interface to this library so that the two approaches are consistent The NISP library is release under the LGPL licence The NISP library provides three tools which are detailed below e The randvar class allows to manage random variables specified by their distribution law and their parameters Once a random variable is created one can generate random numbers from the associated law e The setrandvar class allows to manage a collection of random variables This collection is associated with a sampling method such as MonteCarlo Sobol Quadrature etc It is possible to build the sample and to get it back so that the experiments can be performed e The polychaos class allows to manage a polynomial representation of the simulated model One such object must be associated with a set of experiments which have been performed This set may be read from a data file The object is linked with a collection of random variables Then the coefficients of the polynomial can be computed by integration quadra ture Once done the mean the variance and the Sobol indices can be directly computed from the coefficients The figure 1 1 presents the NISP methodology The process requires that the us
34. e outputs We plot it with the histplot Scilab graphic function which produces the figure 7 4 polychaos buildsample pc Lhs 10000 0 sample output polychaos getsample pc scf histplot 50 sample_output xtitle Ishigami Histogram We can plot a bar graph of the sensitivity indices as presented in figure 7 5 for i 1 nx indexfirst i polychaos getindexfirst pc i indextotal i polychaos getindextotal pc i end scf bar indextotal 0 2 blue bar indexfirst 0 15 yellow 63 Fonction Ishigami Histogramme normalis 0 12 0 10 0 08 0 06 0 04 0 02 ll 10 5 0 0 00 mJ 10 15 20 Figure 7 4 Ishigami function Histogram of the variable on a LHS design with 10000 experiments 64 legend Total First order pos 1 xtitle Ishigami Sensitivity indices Fonction Ishigami Indice de sensibilit Cc totale CU premier ordre Figure 7 5 Ishigami function Sensitivity indices 65 Chapter 8 Thanks Many thanks to Allan Cornet who helped us many times in the creation of this toolbox 66 Bibliography 1 T Homma and A Saltelli Importance measures in global sensitivity analysis of non linear models Reliability Engineering and System Safety 52 1 17 1996 2 Julien Jacques Contributions l analyse de sensibilit et l analyse discriminante g n ral is e 2005 3 Didier Pelat Bases et m thodes pour le traitement
35. er has a numerical solver which has the form Y f X where X are input uncertain parameters and Y are output random variables The method is based on the following steps e We begin by defining normalized random variables For example we may use a random variables in the interval 0 1 or a Normal random variable with mean 0 and variance 1 This choice allows to define the basis for the polynomial chaos denoted by Yx x gt o Depending on the type of random variable the polynomials Yx x gt o are based on Hermite Legendre or Laguerre polynomials e We can now define a Design Of Experiments DOE and with random variable transforma tions rules we get the physical uncertain parameters X Several types of DOE are available Monte Carlo Latin Hypercube Sampling etc If N experiments are required the DOE define the collection of normalized random variables amp 1 w Transformation rules allows to compute the uncertain parameters X 1 v which are the input of the numerical solver f e We can now perform the simulations that is compute the collection of outputs Yi 1 w where Y f Xi e The variables Y are then projected on the polynomial basis and the coefficients y are computed by integration or regression Random Uncertain Numerical Spectral Variable Parameter gt Solver Projection 6 A Y f X Y 2 y w 5 Figure 1 1 The NISP methodology 13 The
36. ests The current version is based on the NISP Library v2 1 23 3 2 Installing the toolbox from ATOMS The ATOMS component is the Scilab tool which allows to search download install and load toolboxes ATOMS comes with Scilab v5 2 T he Scilab NISP toolbox has been packaged and is provided mainly by the ATOMS component The toolbox is provided in binary form depending on the user s operating system The Scilab NISP toolbox is available for the following platforms e Windows 32 bits e Linux 32 bits 64 bits e Mac OS X The ATOMS component allows to use a toolbox based on compiled source code without having a compiler installed in the system Installing the Scilab NISP toolbox from ATOMS requires the following steps e atomsList prints the list of current toolboxes e atomsShow prints informations about a toolbox e atomsInstall installs a toolbox on the system e atomsLoad loads a toolbox Once installed and loaded the toolbox will be available on the system from session to session so that there is no need to load the toolbox again it will be available right from the start of the session In the following Scilab session we use the atomsList function to print the list of all ATOMS toolboxes atomsList ANN Toolbox ANN Toolbox dde_toolbox Dynamic Data Exchange client for Scilab module lycee Scilab pour les lyc l es NISP Non Intrusive Spectral Projection plotlib Matlab like Plotting librar
37. getmultind pc polychaos_getlog pc polychaos getinvquantile pc threshold polychaos getindextotal pc polychaos getindexfirst pc ny polychaos_getdimoutput pc nx 7 polychaos_getdiminput pc p polychaos getdimexp pc no 7 polychaos_getdegree pc polychaos getcovariance pc polychaos getcorrelation pc polychaos_getanova L pc polychaos generatecode pc filename funname polychaos computeoutput pc polychaos computeexp pc srv method polychaos computeexp pc pc2 invalue varopt polychaos buildsample pc type np order Figure 7 2 More methods from the polychaos class 57 The first design of experiment is build on the stochastic set srvx and based on a Quadrature type of DOE Then this DOE is transformed into a DOE for the uncertain collection of parameters srvu degre 2 setrandvar_buildsample srvx Quadrature degre setrandvar_buildsample srvu srvx The next steps will be to create the polynomial and actually perform the DOE But before doing this we can take a look at the DOE associated with the stochastic and uncertain collection of random variables We can use the setrandvar_getsample twice and get the following output gt setrandvar_getsample srvx ans 1 7320508 0 1127017 1 7320508 0 5 1 7320508 0 8872983 0 0 1127017 0 0 5 0 0 8872983 1 7320508 0 1127017 1 7320508 0 5 1 7320508 0 8872983 gt set
38. getsize srvu polychaos_setsizetarget pc np This is slow for k 1 np inputdata setrandvar getsample srvu k outputdata Exemple inputdata polychaos settarget pc k outputdata end The previous loop works but is slow when np is large Instead the following script is fast because it uses vectorization inputdata setrandvar getsample srvu outputdata Exemple inputdata polychaos settarget pc outputdata We can now compute the polynomial expansion by integration polychaos setdegree pc degre polychaos computeexp pc srvx Integration Everything is now ready so that we can do the sensitivy analysis as in the following script average polychaos getmean pc var polychaos getvariance pc mprintf Meanuuuuuuuu u fAn average mprintf Varianceyyuu sf n var mprintf First order sensitivity index n u Variable X1 f1n polychaos_getindexfirst pc 1 mprintf u uuuyVariable X2 4fAn polychaos_getindexfirst pc 2 mprintf uuuuVariableLX3L u f n polychaos_getindexfirst pc 3 mprintf Total sensitivity index n u Variable X1 4f n polychaos_getindextotal pc 1 uVariable X2 4f n polychaos_getindextotal pc 2 uVariableuX3u u4fAn polychaos_getindextotal pc 3 mprintf mprintf mprintf mprintf y The previous script produces the following output
39. ious inequality by V Y and get po wei _ Ti VW a Therefore the sum of the first order sensitivity indices satisfies the inequality Si 5 lt 1 19 2 55 2 56 2 57 2 58 2 59 2 60 2 61 2 62 2 63 2 64 2 65 2 66 2 67 Hence in this example one part of the variance V Y cannot be explained neither by Xi alone or by X alone because it is caused by the interactions between X and X We define by S12 the sensitivity index associated with the group of variables X X5 as 010 vw S2 1 Si S2 2 68 The following Scilab script performs the sensitivity analysis on the previous example We consider two normal variables where the first variable has mean 1 5 and standard deviation 0 5 while the second variable has mean 3 5 and standard deviation 2 5 function y Exemple x y 1 x 1 x 2 endfunction First variable Normal mul 1 5 sigmal 0 5 Second variable Normal mu2 3 5 sigma2 2 5 1 Two stochastic normalized variables vx1 randvar new Normale vx2 randvar new Normale 2 A collection of stoch variables srvx setrandvar_new setrandvar_addrandvar srvx vxl setrandvar_addrandvar srvx vx2 3 Two uncertain parameters vul randvar_new Normale mul sigmal vu2 randvar_new Normale mu2 sigma2 4 A collection of uncertain parameters srvu setrandv
40. ip 2 4 The linear correlation coefficient between Y and X is Cov Y Xi PY X 2 10 VV Y VV X for i 1 2 p In the particular case of the affine model 2 4 we have CovY Xi Cov fo Xi PB Cov X Xi B2Cov Xa Xi 2 11 BiCov Xi Xi SE gie By Cov Xp Xi 2 12 2 13 Since the random variables X are independent we have Cov X X 0 for any j i There fore Cov Y X BiCov Xi Xi 2 14 BV 2 15 Hence the correlation coefficient can be simplified into x BL 2 16 Pr VOY OO _ VV 2 17 V Y 11 We square the previous equality and get AVX Therefore the square of the linear correlation coefficient is equal to the first order sensitivity index i e 2 18 Pix SRC 2 19 2 4 Using scatter plots In this section we present an example of an affine model where the difference between local and global sensitivity is made clearer by the use of scatter plots Assume four independent random variables X for i 1 2 3 4 We assume that the variables X are normally distributed with zero mean and i variance Let us consider the affine model Notice that the derivative of Y with respect to any of its input is equal to one i e OY 1 2 21 for 1 2 3 4 This means that locally the inputs all have the same effect on the output Let us compute the standardized regression coefficients of this model By hypothesis the variance of each vari
41. is E Y f X4 XoF Xi Xa dx1dx2 2 37 17 where F x x2 is the joint probability distribution function of the variables X and X Since X and X are independent variables we have F x1 2 F X1 F2 X35 2 38 where Fj is the probability distribution function of X and F gt is the probability distribution function of X5 Then we have E Y i X X3 FACX1 F5 X2 dr1dt2 ZS XF An lm D WEEN E X E X3 Therefore The variance of Y is The expectation E Y is EV fS IRE E GII i XEM E Now we have which leads to Therefore E Y V X1 E X1 V X5 E X2 oi ui ea 15 Finally we get V Y ei ui e3 13 pupo 18 2 39 2 40 2 41 2 42 2 43 2 44 2 45 2 46 2 47 We can expand the previous equality and get V Y 0703 0143 1702 pis LST The last two terms of the previous equality can be simplified so that we get V Y 0103 0113 1303 The sensitivity indices can be computed from the definitions z VUE a VECI Es E ANA 40 We have E Y X1 E X2 X u2X1 Similarily F Y X2 11X2 Hence V m2X3 ELSE V X2 S2 VY S We get i pV X1 cU ECT V Y Di V Y Finally the first order sensitivity indices are 3 H207 V Y g pio V Y Since 0202 gt 0 we have 22 2 2 2 2 22 22 1201 109 V Y 0703 01413 H103 We divide the prev
42. l function can be used to install the toolbox for all the users of this computer We finally load the toolbox with the atomsLoad function gt atomsLoad NISP Start NISP Toolbox Load gateways Load help Load demos ans INISP 2 1 D Programs SC3623 1 contrib NISP 2 1 Now that the toolbox is loaded it will be automatically loaded at the next Scilab session 3 3 Installing the toolbox from the sources In this section we present the steps which are reguired in order to install the toolbox from the sources In order to install the toolbox from the sources a compiler is reguired to be installed on the machine This toolbox can be used with Scilab v5 1 and Scilab v5 2 We suppose that the archive has been unpacked in the tbxnisp directory The following is a short list of the steps which are reguired to setup the toolbox 1 build the toolbox run the tbrnisp builder sce script to create the binaries of the library create the binaries for the gateway generate the documentation 2 load the toolbox run the tbrnisp load sce script to load all commands and setup the documentation 25 3 setup the startup configuration file of your Scilab system so that the toolbox is known at startup see below for details 4 run the unit tests run the tbrnisp runtests sce script to perform all unit tests and check that the toolbox is OK 5 run the demos run the tbrnisp rundemos sce script to run all demonstration scripts and get a
43. n transferred to the polynomial pc by the polychaos settarget function for k 1 np inputdata setrandvar getsample srvu k outputdata Exemple inputdata mprintf Experiment d inputy f 4 uoutputyu u fYn k inputdata 1 inputdata 2 outputdata polychaos settarget pc k outputdata end The previous script produces the following output Experiment 1 input 0 133975 1 169052 output 0 156623 Experiment 2 input 0 133975 1 750000 output 0 234456 Experiment 3 input 0 133975 2 330948 output 0 312288 Experiment 4 input 1 000000 1 169052 output 1 169052 Experiment 45 input 1 000000 1 750000 output 1 750000 Experiment 6 input 1 000000 2 330948 output 2 330948 Experiment 7 input 1 866025 1 169052 output 2 181482 Experiment 8 input 1 866025 1 750000 output 3 265544 Experiment 9 input 1 866025 2 330948 output 4 349607 There is actually a much faster way of computing the output Indeed using vectorisation we can compute all the outputs in one single call to the Exemple function inputdata setrandvar getsample srvu outputdata Exemple inputdata polychaos settarget pc outputdata We can compute the polynomial expansion based on numerical integration so that the coeffi cients of the polynomial are determined This is done with the polychaos_computeexp function which stands for compute the expansion polychaos setdegree pc degre polychaos compu
44. nd Scilab It can be used as a meta model The figure 7 1 presents the most commonly used methods available in the polychaos class More methods are presented in figure 7 2 The inline help contains the detailed calling sequence for each function and will not be repeated here More than 50 methods are available and most of them will not be presented here More informations about the Oriented Object system used in this toolbox can be found in the section 5 2 2 7 2 Examples In this section we present to examples of use of the polychaos class 7 2 1 Product of two random variables In this section we present the polynomial expansion of the product of two random variables We analyse the Scilab script and present the methods which are available to perform the sensi 59 Constructors pc polychaos new file pc polychaos new L srv ny pc polychaos new pc nopt varopt Methods polychaos_setsizetarget pc np polychaos settarget pc output polychaos setinput pc invalue polychaos setdimoutput pc ny polychaos setdegree pc no polychaos getvariance pc polychaos getmean pc Destructor polychaos destroy pc Static methods tokenmatrix polychaos tokens nb polychaos size Figure 7 1 Outline of the methods of the polychaos class tivity analysis T his script is based on the NISP methodology which has been presented in the Introduction chapter We will use
45. new Uniforme vu2 1 vu3 randvar_new Uniforme vu3 2 There are now 3 objects in current use as indicated by the following statements The tokenmatrix is a row matrix containing regular double precision floating point numbers nb randvar size nb 3 tokenmatrix randvar_tokens tokenmatrix 0 1 2 We assume that we have now made our job with the random variables so that it is time to destroy the random variables We call the randvar destroy functions which destroys the variables randvar destroy vui randvar destroy vu2 randvar destroy vu3 We can finally check that there are no random variables left in the memory space 34 gt nb randvar size nb 0 gt tokenmatrix randvar_tokens tokenmatrix Scilab is a wonderful tool to experiment algorithms and make simulations It happens some times that we are managing many variables at the same time and it may happen that at some point we are lost The static methods provides tools to be able to recover from such a situation without closing our Scilab session In the following session we create two random variables gt vui randvar_new Uniforme vul 3 vu2 randvar_new Uniforme vu2 4 Assume now that we have lost the token associated with the variable vu2 We can easily simulate this situation by using the clear which destroys a variable from Scilab s memory space clear
46. om number in the interval 0 1 is generated and then scaled with X a b a X 5 6 Let us now analyse how the uniform random number X 0 1 is computed The uniform random generator is based on the C function rand which returns an integer n in the interval 32 0 RAND MAX The value of the RAND MAX variable is defined in the file stdlib h and is compiler dependent For example with the Visual Studio C 2008 compiler the value is RAND MAX 2 1 32767 5 7 A uniform value X in the range 0 1 is computed from ze N 7 where N RAND MAX and n 0 RAN D M AX 5 8 5 2 Methods In this section we give an overview of the methods which are available in the randvar class 5 2 1 Overview The figure 5 3 presents the methods available in the randvar class The inline help contains the detailed calling sequence for each function and will not be repeated here Constructors rv randvar_new type options Methods value randvar_getvalue rv options randvar_getlog rv Destructor randvar_destroy rv Static methods rvlist randvar_tokens nbrv randvar_size Figure 5 3 Outline of the methods of the randvar class 5 2 2 The Oriented Object system In this section we present the token system which allows to emulate an oriented object program ming with Scilab We also present the naming convention we used to create the names of the f
47. p Compilation of sci polychaos setgroupaddvar cpp Compilation of sci polychaos computeoutput cpp Compilation of sci polychaos setinput cpp Compilation of sci polychaos propagateinput cpp Compilation of sci_polychaos_getanova cpp Compilation of sci polychaos setanova cpp Compilation of sci polychaos getanovaord cpp Compilation of sci polychaos getanovaordco cpp Compilation of sci polychaos realisation cpp Compilation of sci polychaos save cpp Compilation of sci polychaos generatecode cpp Building shared library be patient Generate a cleaner file Generating loader gateway sce Building help Building the master document C tbxnisp help en_US Building the manual file javaHelp C tbxnisp help en_US Please wait building Generating loader sce in this can take a while The following script presents the messages which are generated when the loader of the toolbox 28 is launched The loader script performs the following steps e load the gateway and the NISP library e load the help e load the demo exec C tbxnisp loader sce Start NISP Toolbox Load gateways Load help Load demos It is now necessary to setup your Scilab system so that the toolbox is loaded automatically at startup The way to do this is to configure the Scilab startup configuration file The directory where this file is located is stored in the Scilab variable SCIHOME In the following Scilab session we use Scilab v5 2 0 beta 1 in order to
48. position Y f A zit DI films ti bha alftn gt Tp 2 70 i 1 2 p 1 lt i lt j lt p where fo is a constant and the function of the decomposition satisfy the equalities 1 Fi ye s ECH Ti dEi 0 2 71 0 2 8 Ishigami function TODO 22 Chapter 3 Installation In this section we present the installation process for the toolbox We present the steps which are required to have a running version of the toolbox and presents the several checks which can be performed before using the toolbox 3 1 Introduction There are two possible ways of installing the NISP toolbox in Scilab e use the ATOMS system and get a binary version of the toolbox e build the toolbox from the sources The next two sections present these two ways of using the toolbox Before getting into the installation process let us present some details of the the internal components of the toolbox The following list is an overview of the content of the directories tbanisp demos demonstration scripts tbanisp doc the documentation tbanisp doc usermanual the TEXsources of this manual tbznisp etc startup and shutdow scripts for the toolbox tbrnisp help inline help pages tbanisp macros Scilab macros files sci tbrnisp sci gateway the sources of the gateway tbrnisp src the sources of the NISP library tbinisp tests tests tbxnisp tests nonreg tests tests after some bug has been identified tbxnisp tests unit tests unit t
49. randvar_getsample srvu ans 0 1339746 1 1690525 0 1339746 1 75 0 1339746 2 3309475 1 1 1690525 1 1 75 1 2 3309475 1 8660254 1 1690525 1 8660254 1 75 1 8660254 2 3309475 These two matrices are a 9x2 matrices where each line represents an experiment and each column represents an input random variable The stochastic normalized srvx DOE has been created first then the srvu has been deduced from srvx based on random variable transformations We now use the polychaos_new function and create a new polynomial pc The number of input variables corresponds to the number of variables in the stochastic collection srvx that is 2 and the number of output variables is given as the input argument ny In this particular case the number of experiments to perform is equal to np 9 as returned by the setrandvar_getsize function This parameter is passed to the polynomial pc with the polychaos_setsizetarget function ny 1 pc polychaos new srvx ny np setrandvar getsize srvx polychaos setsizetarget pc np 98 In the next step we perform the simulations prescribed by the DOE We perform this loop in the Scilab language and make a loop over the index k which represents the index of the current experiment while np is the total number of experiments to perform For each loop we get the input from the uncertain collection srvu with the setrandvar getsample function pass it to the Exemple function get back the output which is the
50. ression coefficients of affine models In this section we present the standardized regression coefficients of an affine model Assume that the random variables X are independent with mean E X and finite variances V X for i 1 2 p Let us consider the random variable Y as an affine function of the variables X Y bo y PX 2 4 i 1 2 p where 5 are real parameters for i 1 2 p The expectation of the random variable Y is E Y bo y B E X 2 5 i 1 2 p Since the variables X are independent the variance of the sum of variables is the sum of the variances Hence V Y V B Y V GXI 2 6 i 1 2 p 10 which leads to the equality V Y gt EV 2 7 1 1 2 Hence each term 82V X is the part of the total variance V Y which is caused by the variable Xi The standardized regression coefficient is defined as BV X SRC e 2 8 VY 2 8 for i 1 2 p Hence the sum of the standardized regression coefficients is one SRO SRO A FSRO 1 2 9 2 3 Link with the linear correlation coefficients In this section we present the link between the linear correlation coefficients of an affine model and the standardized regression coefficients Assume that the random variables X are independent with mean E X and finite variances V X for i 1 2 p Let us consider the random variable Y which depends linearily on the variables X by the relationsh
51. ript produces the figures 6 6 and 6 7 Empirical histogram of X1 Empirical Exact Figure 6 6 Monte Carlo Sampling Histogram and exact distribution functions for the first variable 6 2 8 A LHS design In this section we present the creation of a Latin Hypercube Sampling In our example the DOE is based on two random variables the first being Normal with mean 1 0 and standard deviation 0 5 and the second being Uniform in the interval 2 3 We begin by defining two random variables with the randvar new function vul randvar_new Normale 1 0 0 5 vu2 randvar_new Uniforme 2 0 3 0 Then we create a collection of random variables with the setrandvar_new function which creates here an empty collection of random variables Then we add the two random variables to the collection 47 Empirical histogram of X2 Empirical Exact Figure 6 7 Monte Carlo Sampling Histogram and exact distribution functions for the second variable srv setrandvar new setrandvar addrandvar srv vul setrandvar_addrandvar srv vu2 We can now build the DOE so that it is a LHS sampling with 1000
52. s interpreted The figure 6 1 presents the methods available in the setrandvar class A complete description of the input and output arguments of each function is available in the inline help and will not be repeated here More informations about the Oriented Object system used in this toolbox can be found in the section 5 2 2 6 2 Examples In this section we present examples of use of the setrandvar class In the first example we present a Scilab session where we create a Latin Hypercube Sampling In the second part we present various types of DOE which can be generated with this class 6 2 1 A Monte Carlo design with 2 variables In the following example we build a Monte Carlo design of experiments with 2 input random variables The first variable is associated with a Normal distribution function and the second 41 Constructors srv setrandvar new srv setrandvar new n srv setrandvar new file Methods setrandvar_setsample setrandvar_setsample setrandvar_setsample srv k value setrandvar setsample srv value setrandvar save srv file np setrandvar getsize srv sample setrandvar getsample srv k sample setrandvar getsample srv k sample setrandvar getsample srv setrandvar getlog srv nx setrandvar getdimension srv setrandvar freememory srv setrandvar buildsample srv srv2 setrandvar_buildsample srv name np 7 setrandvar buildsample srv name np
53. ses on the behaviour of the mapping in the neighbourhood of a particular point X Dx toward a particular point Y Dy The global sensisitivity analysis models the propagation of uncertainties transforming the whole set Dx into the set Dy In the following we assume that there is only one output variable so that Y R There are two types of analysis that we are going to perform that is uncertainty analysis and sensitivity analysis In uncertainty analysis we assume that fx is the probability density function of the variable X and we are searching for the probability density function fy of the variable Y and by its cumulated density function Fy y P Y lt y This problem is difficult in the general case and this is why we often are looking for an estimate of the expectation of Y as defined by EQ ufr 2 2 and an estimate of its variance Var Y v BO few 2 3 We might also be interested in the computation of the probability over a threshold In sensitivity analysis we focus on the relative importance of the input variable X on the uncertainty of Y This way we can order the input variables so that we can reduce the variability of the most important input variables in order to finally reduce the variability of Y More details on this topic can be found in the papers of Homma and Saltelli 1 or in the work of Sobol 4 The Thesis by Jacques 2 presents an overview of sensitivity analysis 2 2 Standardized reg
54. setrandvar save cpp Sci setrandvar buildsample cpp Sci polychaos new cpp Sci polychaos destroy cpp Sci polychaos tokens cpp Sci polychaos size cpp Sci polychaos setdegree cpp Sci polychaos getdegree cpp Sci polychaos freememory cpp 27 Compilation of sci polychaos getdimoutput cpp Compilation of sci polychaos setdimoutput cpp Compilation of sci polychaos getsizetarget cpp Compilation of sci polychaos setsizetarget cpp Compilation of sci polychaos freememtarget cpp Compilation of sci polychaos settarget cpp Compilation of sci polychaos gettarget cpp Compilation of sci polychaos getdiminput cpp Compilation of sci polychaos getdimexp cpp Compilation of sci polychaos getlog cpp Compilation of sci polychaos computeexp cpp Compilation of sci polychaos getmean cpp Compilation of sci polychaos getvariance cpp Compilation of sci polychaos getcovariance cpp Compilation of sci polychaos getcorrelation cpp Compilation of sci polychaos getindexfirst cpp Compilation of sci polychaos getindextotal cpp Compilation of sci polychaos getmultind cpp Compilation of sci polychaos getgroupind cpp Compilation of sci polychaos setgroupempty cpp Compilation of sci polychaos getgroupinter cpp Compilation of sci polychaos getinvquantile cpp Compilation of sci polychaos buildsample cpp Compilation of sci polychaos getoutput cpp Compilation of sci polychaos getquantile cpp Compilation of sci polychaos getquantwilks cpp Compilation of sci polychaos getsample cp
55. t we create a Latin Hypercube Sampling made of 10 000 points Then get the output of the polynomial on these inputs and plot the histogram of the output polychaos buildsample pc Lhs 10000 0 sample output polychaos getsample pc scf histplot 50 sample_output xtitle Product function Empirical Histogram X P X The previous script produces the figure 7 3 We may explore the following topics e Perform the same analysis where the variable X gt is a normal variable with mean 2 and standard deviation 2 e Check that the development in polynomial chaos on a Hermite Hermite basis does not allow to get exact results See that the convergence can be obtained by increasing the degree e Check that the development on a basis Hermite Legendre allows to get exact results with degree 2 7 2 2 The Ishigami test case In this section we present the Ishigami test case The function Exemple is the model that we consider in this numerical experiment This function takes a vector of size 3 in input and returns a scalar output function y Exemple x a 7 60 Product function Empirical Histogram 0 45 7 0 40 7 0 35 0 30 7 0 25 7 0 20 7 0 15 7 0 10 7 0 05 0 00 y T g T z T 3 T T y T T T
56. ted we use the Scilab function mean to check that the mean is close to 1 which is the expected value of the Normale law when the number of samples is infinite Finally we use the randvar destroy function to destroy our random variable Once done we plot the empirical distribution function of this sample with 50 classes nisp initseed O mu 1 0 sigma 0 5 rv randvar_new Normale mu sigma nbshots 1000 values zeros nbshots for i 1 nbshots values i randvar getvalue rv end mymean mean values mysigma st deviation values myvariance variance values mprintf Meanj is j fi expected 2f in mymean mu mprintf Standard deviationy isy yhfyCexpectedy 4f n mysigma sigma mprintf Variance isy 1 4fu expected 0 4f n myvariance sigma 2 randvar_destroy rv histplot 50 values xtitle Histogram of X X P x The previous script produces the following output Mean is 0 988194 expected 1 000000 Standard deviation is 0 505186 expected 0 500000 Variance is 0 255213 expected 0 250000 The previous script also produces the figure 5 4 5 3 2 Variable transformations In this section we present the transformation of uniform random variables into other types of variables The transformations which are available in the randvar class are presented in figure 5 5 We begin the analysis by a presentation of the theory required to perform transformations Then
57. teexp pc srvx Integration Everything is now ready for the sensitivity analysis Indeed the polychaos_getmean returns the mean while the polychaos_getvariance returns the variance average polychaos getmean pc var polychaos getvariance pc mprintf Meany yo u usf n average mprintf Variancey yu y ust n var mprintf Indice de sensibilite du ler ordre n mprintf Guu Variable X1 4f n polychaos_getindexfirst pc 1 mprintf Guu Variable X2 f n polychaos_getindexfirst pc 2 mprintf Indice de sensibilite Totale n mprintf Guu Variable X1 4f n polychaos_getindextotal pc 1 mprintf uuuuVariable X2 fin polychaos_getindextotal pc 2 99 The previous script produces the following output Mean 1 750000 Variance 1 000000 Indice de sensibilite du 1er ordre Variable X1 0 765625 Variable X2 0 187500 Indice de sensibilite Totale Variable X1 0 812500 Variable X2 0 234375 In order to free the memory required for the computation it is necessary to delete all the objects created so far polychaos_destroy pc randvar destroy vu1 randvar destroy vu2 randvar destroy vx1 randvar destroy vx2 setrandvar_destroy srvu setrandvar_destroy srvx Prior to destroying the objects we can inquire a little more about the density of the output of the chaos polynomial In the following scrip
58. unctions The randvar class provides the following functions e The constructor function randvar new allows to create a new random variable and returns a token rv e The method randvar_getvalue takes the token rv as its first argument In fact all methods takes as their first argument the object on which they apply 33 e The destructor randvar destroy allows to delete the current object from the memory of the library e The static methods randvar tokens and randvar size allows to quiery the current object which are in use More specifically the randvar size function returns the number of current randvar objects and the randvar tokens returns the list of current randvar objects In the following Scilab sessions we present these ideas with practical uses of the toolbox Assume that we start Scilab and that the toolbox is automatically loaded At startup there are no objects so that the randvar size function returns 0 and the randvar tokens function returns an empty matrix nb randvar size nb 0 tokenmatrix randvar tokens tokenmatrix We now create 3 new random variables based on the Uniform distribution function We store the tokens in the variables vu1 vu2 and vu3 These variables are regular Scilab double precision floating point numbers Each value is a token which represents a random variable stored in the toolbox memory space gt vul randvar_new Uniforme vul 0 vu2 randvar
59. we present some of the many the transformations which are provided by the library We now present some additionnal details for the function randvar_getvalue rv rv2 value2 This method allows to transform a random variable sample from one law to another The statement 36 Empirical normal probability distribution function P x Figure 5 4 The histogram of a Normal random variable with 1000 samples Source Target Source Target Normale LogNormale Normale Normale Uniforme Uniforme Exponentielle Exponentielle LogNormale LogNormale LogUniforme LogUniforme Source Target Source Target Uniforme LogUniforme Uniforme Uniforme Normale Normale Exponentielle Exponentielle LogNormale LogNormale LogUniforme LogUniforme Source Target Exponentielle Exponentielle Figure 5 5 Variable transformations available in the randvar class 3T value randvar getvalue rv rv2 value2 returns a random value from the distribution function of the random variable rv by transformation of value2 from the distribution function of random variable rv2 In the following session we transform a uniform random variable sample into a LogUniform variable sample We begin to create a random variable rv from a LogUniform law and parameters a 10 b 20 Then we create a second random variable rv2 from a Uniforme law and parameters a 2 b 3 The main loop is based on the transformation of a sample
60. ws to configure the global behaviour of the toolbox These func tions allows to startup and shutdown the toolbox and initialize the seed of the random number generator They are presented in the figure 4 1 nisp startup Starts up the NISP toolbox nisp shutdown Shuts down the NISP toolbox level nisp verboselevelget Returns the current verbose level nisp verboselevelset level Sets the value of the verbose level nisp initseed seed Sets the seed of the uniform random number generator nisp_destroyall Destroy all current objects nisp_getpath Returns the path to the current module nisp_printall Prints all current objects Figure 4 1 Outline of the configuration methods The user has no need to explicitely call the nisp_startup and nisp shutdown 7 func tions Indeed these functions are called automatically by the etc NISP start and etc NISP quit scripts located in the toolbox directory structure The nisp initseed seed is especially useful when we want to have reproductible re sults It allows to set the seed of the generator at a particular value so that the sequence of uniform pseudo random numbers is deterministic When the toolbox is started up the seed is automatically set to 0 which allows to get the same results from session to session 30 Chapter 5 The randvar class In this section we present the randvar class which allows to define a random variable and to generate random num
61. y for Scilab Scipad Scipad 7 20 sndfile toolbox Read amp write sound files stixbox Statistics toolbox for Scilab 5 2 In the following Scilab session we use the atomsShow function to print the details about the NISP toolbox atomsShow NISP Package NISP Title NISP Summary Non Intrusive Spectral Projection Version 2 1 Depend Category ies Optimization Maintainer s Pierre Marechal lt pierre marechal scilab org gt Michael Baudin lt michael baudin scilab org gt 24 Entity CEA DIGITEO WebSite License LGPL Scilab Version gt 5 2 0 Status Not installed Description This toolbox allows to approximate a given model which is associated with input random variables This toolbox has been created in the context of the OPUS project http opus project fr within the workpackage 2 1 1 Construction de m l ta mod les This project has received funding by Agence Nationale de la recherche http www agence nationale recherche fr See in the help provided in the help en_US directory of the toolbox for more information about its use Use cases are presented in the demos directory In the following Scilab session we use the atomsInstall function to download and install the binary version of the toolbox corresponding to the current operating system gt atomsInstall NISP ans INISP 2 1 allusers D Programs SC3623 1 contrib NISP 2 1 I The allusers option of the atomsInstal
Download Pdf Manuals
Related Search