NISP Toolbox Manual - Forge
Contents
1. Figure 2 6 The histogram of a Log Uniform random variable with 1000 samples 20 Chapter 3 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 3 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 is interpreted The figure 3 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 2 2 2 3 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 wi2. mprintf Varianceuuuucu 4fAn var mprintf Firstiorder sensitivityuindex n mprintf u uuuyVariable X1 fAn polychaos_getindexfirst pc 1 mprintf u uuyVariable X2 fin polychaos_getindexfirst pc 2 mprintf uuuuVariablewX3u u4fAn polychaos_getindexfirst pc 3 mprintf Total sensitivity index n mprintf Guu Variable X1 Y 4f n polychaos_getindextotal pc 1 mprintf Guu Variable X2 4f n polychaos_getindextotal pc 2 mprintf uuuuVariablewX3u u4fAn polychaos_getindextotal pc 3 The previous script produces the following output 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 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 ofivariables n mprintf uuuuGroupeuXiuwetuX2u 4fAn polychaos_getgroupind pc The previous script produces the following output Fraction of the vari
3. 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 92 Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation Compilation of of of of of of of of of of of nisp_msg cpp nisp_conf cpp nisp_ort cpp nisp_pc cpp nisp_polyrule cpp 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 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 Compilati
4. end for i 1 n 1 plot 0 1 cut i cut i end The previous script produces the figure 3 13 32 LHS Max Min Design 1 O 0 9 0 8 0 7 0 6 0 5 x2 0 4 0 3 0 2 0 1 0 0 0 0 0 1 0 2 0 3 0 4 0 5 0 6 0 7 0 8 0 9 1 0 x1 Figure 3 12 Latin Hypercube Sampling Computed with setrandvar_buildsample and the LhsMaxMin option 0 97 0 8 0 77 0 6 0 57 0 4 O 0 3 7 0 2 0 17 0 0 T T T T T 0 0 0 1 0 2 0 3 0 4 0 5 0 6 0 7 0 8 0 9 1 0 Figure 3 13 Latin Hypercube Sampling Computed with nisp_buildlhs 33 3 2 5 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 3 14 presents this sampling and the figures 3 15 and 3 16 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 Monte Carlo Sampling Variable 2 Uniforme 1 0 1 0 0 2 0 0 0 2 0 4 0 6 0 8 1 0 Variable 1 Uniforme 1 0 1 0 Figure
5. Name f x E X V X Normale VE exp 1 e u oi T z b a if x a b b a b a Uniforme 0 ie ad 5 1 8 BE Xexp Ax ifa gt 0 1 1 Exponentielle d 0 in lt 0 i TW 1 In z p 7 LogNormale ca var XP 2 de exp u ei exp o 1 exp 2u 0 0 ifa lt 0 1 17 if x exp a exp b 2_ 2 Wi 4 b a 1 P EXP exp b exp a 1 exp b exp a 2 LogUniforme 0 Fu 5 ha E x Figure 2 1 Distributions functions of the randvar class The expected value is denoted by E X and the variance is denoted by V X 11 Name Parameter 1 a Parameter 2 b Conditions Normale a U 01 c gt 0 Uniforme a 0 b 1 a lt b Exponentielle a 1 LogNormale u 0 o 1 a U LogUniforme a 0 b 1 a lt b Figure 2 2 Default parameters for distributions functions 2 1 2 The Log Normal distribution A log normal distribution is a probability distribution of a random variable whose logarithm is normally distributed If X is a random variable with a normal distribution then Y exp X has a log normal distribution The LogNormale law is defined by the expected value y and the standard deviation o of the underlying Normal random variable In other words the random variable log X has mean y and variance 0 The expected value and the variance of the Log Normal law are given by A ES exp 2 zc 2 1 explo
6. legend Total First jorder pos 1 xtitle Ishigami Sensitivity indices Fonction Ishigami Indice de sensibilit EST totale CU premier ordre Figure 4 5 Ishigami function Sensitivity indices 48 Chapter 5 Thanks Many thanks to Allan Cornet 49 Bibliography 1 Michael Baudin Programming in scilab 2010 2 Michael Baudin and Jean Marc Martinez Polyn mes de chaos sous scilab via la librairie nisp In 2 mes Journ es de Statistique 24 au 28 mai 2010 Marseille France 2010 3 Didier Pelat Bases et m thodes pour le traitement des donn es Bruits et Signaux Master M2 Recherche Astronomie astrophysique 2006 90 Appendix A Installation A 1 Architecture of the toolbox 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 tbrnisp doc the documentation tbinisp doc usermanual the ETEXsources of this manual tbrnisp etc startup and shutdow scripts for the toolbox tbrnisp help inline help pages tbrnisp macros Scilab macros files sci tbrnisp sci_gateway the sources of the gateway tbanisp src the sources of the NISP library tbanisp tests tests tbanisp tests nonreg_tests tests after some bug has been identified tbanisp tests unit_tests unit tests A 2 Installing the toolbox from the sources In this section we present th
7. 1 exp 2u 0 2 2 p e 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 1 Wa 2 3 o ln 1 nx 2 4 In particular the expected value u of with the Normal random variable satisfies the equation u In E X a 2 5 Caution These are the parameters for NISP v2 5 Earlier versions of the toolbox used other parameters See the changelog txt file in your toolbox to see how to update your scripts 2 1 3 The Log Uniform distribution A log uniform distribution is a probability distribution of a random variable whose logarithm has a uniform distribution If X is a random variable with a uniform distribution then Y exp X has a log uniform distribution The LogUniforme law is defined by the minimum a and the maximum b of the underlying Uniform random variable In other words the random variable log X has minimum a and maximum b 12 If X is a LogUniforme random variable with minimum A and maximum B then a log A 2 6 log B 2 7 Caution These are the parameters for NISP v2 5 Earlier versions of the toolbox used other parameters See the changelog txt file in your toolbox to see how to update your scripts 2 1 4 Uniform random number generation In this section we present the generation of uni
8. 3 gt tokenmatrix randvar_tokens tokenmatrix O 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 14 gt randvar_destroy vul gt randvar_destroy vu2 gt randvar_destroy vu3 We can finally check that there are no random variables left in the memory space gt nb randvar_size nb O gt tokenmatrix randvar_tokens tokenmatrix L 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 vul randvar_new Uniforme vul 3 gt 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 function which destroys a variable from Scilab s memory space gt clear vu2 gt randvar_getvalue vu2 l 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 th
9. 6 X Y 1 X Y 2 vu Figure 1 1 The NISP methodology 13 The 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 e Manage various types of random variables uniform normal exponential log normal Generate random numbers from a given random variable e 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 sta
10. 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 Scilr oC 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 1 5 Configuration functions In this section we present functions which allow to configure the NISP toolbox The nisp_ functions allows 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 1 5 Public functions 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 NISP random number generator nisp_destroyall Destroy all current NISP objects nisp_getpath Returns the path to the current module nisp_printall Prints all current NISP objects Private functions nisp_startup Starts up the NISP toolbox nisp_shutdown Shuts down the NISP toolbox Figure 1 5 Outline of the configuration methods The nisp_initseed se
11. groupe de variables aleatoires srv setrandvar_new setrandvar_addrandvar srv rvi setrandvar_addrandvar srv rv2 np 1000 setrandvar_buildsample srv MonteCarlo np 77 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 26 x linspace 1 3 1000 p dnorm x 1 0 5 plot x p r 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 Empirical histogram of X2 X P X legend Empirical Exact Clean up setrandvar_destroy srv randvar_destroy rv1 randvar_destroy rv2 The previous script produces the figures 3 6 and 3 7 Empirical histogram of X1 Empirical Exact Figure 3 6 Monte Carlo Sampling Histogram and exact distribution functions for the first variable 3 2 3 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 mea
12. o o o o o o o cog pee eee o o D e D D D e o o o e oor besse o o ec o o o 004 besse o o o ccc poo o o o o eo pee eee ec o D e D D e e o e e ood be ee o ee ce e o o e o o eed Peee o o o o poo o o o o o o eed eee ooo eee e o 9 o o ee ip ee S S S S E E S o o e UN o 0 bd d o o e S S S S 5 0 4 0 5 0 6 0 7 Variable 1 Uniforme 0 0 1 0 0 3 0 2 Figure 3 18 Quadrature sampling Two uniform variables in the interval 1 1 36 Sobol Sampling 1 0 F e as a ION PS ds D I Sek s AO Aa D Je A E AE 084 DEA 48 ea U AS a RCA PO ge DN w Wi Be H D 067 AC E e e si Ba ef Sy Tos 7 ce Re 04 4 a Sc B S Se ef sg Ee 1 P wf E Eat ES sen ee s Y e e 8 a AER Z 00 N E H Ze o 02e x oe ef SE KSM 5 ie e tete e gt 047 Ra Vi EN Wi FE 5 T ee e ot fe DNA 2 061 7 A ne te qe DN SE e K N Ea 083 9 H KEE Lo si e ao e D Ee Ta T R a T 1 0 AAA a E 1 0 0 8 0 6 0 4 0 2 0 0 Variable 1 Uniforme 1 0 1 0 Figure 3 19 Sobol sampling Two uniform variables in the interval 1 1 37 Chapter 4 The polychaos class 4 1 Introduction The polychaos class allows to manage a polynomial chaos expansio
13. 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 5 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 While the previous script is perfectly correct it can be very slow when the number of ex periments is large This is because the interpreter has to perform a large number of loops with matrices of small size In general this produces much slower script and should be avoided More details on this topic are presented in 1 44 4 2 3 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 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 r m 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_n
14. 1000 samples 17 are provided by the library Any outcome from a random variable x can be transformed into the outcome of another variable y This can be done by computing u Fx x where Fx is the cumulated distribution function of X This transforms the random variable X into the random variable U which has a uniform distribution in 0 1 Then we invert the CDF by computing y Fy u where Fy is the cumulated distribution function of Y We now present some examples of the function randvar_getvalue rv rv2 value2 which performs this transformation The statement 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 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 0 a 10 0 b 20 0 rv randvar_new LogUniforme a b rv2 randvar_new Uniforme 2 3 nbshots 1000 valuesLou ze
15. 3 14 Monte Carlo Sampling Two uniform variables in the interval 1 1 With the setrandvar_buildsample function we can change the type of sampling by changing the second argument This way we can create the Petras Quadrature and Sobol sampling presented in figures 3 17 3 18 and 3 19 34 Variable 1 Uniforme 1 0 1 0 0 8 0 7 0 6 7 0 57 0 4 0 3 0 2 0 1 0 0 Figure 3 15 Latin Hypercube Sampling First uniform variable in 1 1 Variable 2 Uniforme 1 0 1 0 0 8 0 7 7 0 6 7 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 3 16 Latin Hypercube Sampling Second uniform variable in 1 1 39 Petras Sampling 0 5 kee o o o pe 0 4 kee o o o o e Coon ZH ejgeue Variable 1 Uniforme 0 0 1 0 Figure 3 17 Petras sampling Two uniform variables in the interval 1 1 Quadrature Sampling 0 L 0 0 BUMOJUN ZH oigeue o S A hooo o o
16. NISP Toolbox Manual Michael Baudin EDF Jean Marc Martinez CEA Version 0 4 February 2013 Abstract This document is an introduction to the NISP module v2 5 We present the installation process of the module in binary from ATOMS or from the sources Then we present the configuration functions and the randvar setrandvar and polychaos classes Several examples are provided for each class which provides an overview of the use of NISP in practical situations Contents 1 Introduction Et ARS UG POE cacon RARA Ae Oe AO O TS The ODP DDI 52214 Sea YR RAS RAR ER AE RD AAA AA lay The NIGP maane MM 1 4 Installing the toolbox from ATOMS Parc MM 2 The randvar class 241 The distart on MEMBRES 2 ek ds we de pk ra A A 2 1 2 The Log Normal distribution 2 2 44 REF ers dau RR 2 1 3 The Log Uniform distribution 4 4 4 se RRRRRRRFER GR eS 2 1 4 Uniform random number generation Oo a Lai Ni ni nie in e AF a ia SS BS e NEI ON YN Ori a Na Dee ads ie 2 242 The Cee Ob Ject systemi 2 544 FR A E EENHEETEN CAN T ROA WEE fog a e c secor e De RR UE MES MES 3 The setrandvar class dl TOME cse eia O R caidas dl IATA e e e Y RW RAEE A e YI HY A 3 2 1 A Monte Carlo design with 2 variables 3 2 2 A Monte Carlo design with 2 variables 3 2 3 A LHS design cierra A 3 2 4 A note on the LHS samplings k SF FRAU
17. UR AA 3 2 5 Other types of DOES RR RR ee RE FR RR ER ba 4 The polychaos class EL GE Y RE AAA Pe A Boe we LS PED ea pa a af iu Be md EA Wb ACB dn MAB AR E A21 Praduet ol two random variables 2 24 45 44 5454 4245845454 42 2 A cC ST ee SRE ROR eS A R d 423 The leen test case o coc Las die dos de AE OER dE dot de 5 Thanks Jarre 10 11 11 12 12 13 13 13 13 16 16 16 21 21 21 21 25 27 31 34 38 38 38 38 43 45 49 Bibliography A Installation A 1 Architecture of the toolbox A 2 Installing the toolbox from the sources Copyright 2012 2013 Michael Baudin Copyright 2008 2011 INRIA Michael Baudin Copyright 2008 2011 CEA Jean Marc Martinez This file must be used under the terms of the Creative Commons Attribution ShareAlike 3 0 Unported License http creativecommons org licenses by sa 3 0 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 Licen
18. _ 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 4 1 Outline of the methods of the polychaos class tivity analysis This script is based on the NISP methodology which has been presented in the Introduction chapter We will use 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 vx1 randvar_new Normale vx2 randvar_new Uniforme srvx setrandvar_new setrandvar_addrandvar srvx vxi 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 deviati
19. ance 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 1 0 0 0 313953 O 1 O 0 442325 1 1 0 1 55229e 009 001 8 08643e 031 46 O 1 0 243721 1 1 7 26213e 031 1 1 1 6007e 007 EB Oor 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 the outputs We plot it with the histplot Scilab graphic function which produces the figure 4 4 polychaos_buildsample pc Lhs 10000 0 sample_output polychaos_getsample pc scf histplot 50 sample_output xtitle Ishigami Histogram Fonction Ishigami Histogramme normalis 0 12 0 10 0 08 0 06 0 04 0 02 0 00 1 7 f IE 10 5 0 5 10 15 20 Figure 4 4 Ishigami function Histogram of the output of the chaos polynomial on a LHS design with 10 000 experiments We can plot a bar graph of the sensitivity indices as presented in figure 4 5 for i 1 nx indexfirst i polychaos_getindexfirst pc i indextotal i polychaos_getindextotal pc i end scf 47 bar indextotal 0 2 blue bar indexfirst 0 15 yellow
20. andvar_buildsample srvx Guadrature 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 O 0 1127017 O 0 5 O 0 8872983 1 7320508 0 1127017 1 7320508 0 5 1 7320508 0 8872983 gt setrandvar_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
21. argument the object on which they apply 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 gt nb randvar_size nb O gt tokenmatrix randvar_tokens tokenmatrix L 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 vui randvar_new Uniforme vul O gt vu2 randvar_new Uniforme vu2 1 gt vu3 randvar_new Uniforme vu3 2 There are now 3 objects in current use as indicated by the following statements The token matrix is a row matrix containing regular double precision floating point numbers gt nb randvar_size nb
22. ce LGPL as all components of the OPUS project This module was presented in the 42 mes Journ es de Statistique du 24 au 28 mai 2010 2 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 4 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
23. computed from the coefficients The figure 1 1 presents the NISP methodology The process requires that the user 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 Jr Ku Depending on the type of random variable the polynomials Vy 4 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 reguired the DOE define the collection of normalized random variables JC ze 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 Y i i where Y f X 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 gt Parameter Solver gt Projection
24. d 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 30 Variable 2 Uniforme 2 0 3 0 Figure 3 10 Latin Hypercube Sampling Uniform random variable Our numerical simulation is now finished but we must destroy the objects so that the memory managed by the toolbox is deleted randvar_destroy vul randvar_destroy vu2 setrandvar_destroy srv 3 2 4 A note on the LHS samplings We emphasize that the LHS sampling which is provided by the setrandvar_buildsample function is so that the points are centered within their cells In the following script we create a LHS sampling with 10 points srv setrandvar_new 2 n
25. e 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 A gt randvar_tokens ans L 15 2 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 2 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 computed 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 plo
26. e 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 51 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 tbznisp load sce script to load all commands and setup the documentation 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 tbznisp rundemos sce script to run all demonstration scripts and get a 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 gt exec C tbxnisp builder sce Building sources
27. ed is useful when we want to have reproductible results It allows to set the seed of the generator at a particular value so that the seguence 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 The user has no need to explicitely call the nisp_startup and nisp_shutdown functions Indeed these functions are called automatically by the etc NISP start and etc NISP quit scripts located in the toolbox directory structure See the section A 1 for details on this topic 10 Chapter 2 The randvar class In this section we present the randvar class which allows to define a random variable and to generate random numbers from a given distribution function 2 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 2 1 1 Overview The table 2 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 2 2 which also presents the conditions which must be satisfied by the parameters
28. es 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 gt 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 polychaosl passed ref created 003 004 tbxnisp randvari passed ref created 004 004 tbxnisp setrandvarl passed ref created Summary tests 4 100 passed 0 o failed 0 0 skipped 0 0 length 3 84 sec Tests ending the 2009 11 18 at 12 47 48 96
29. ew 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 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_getsize srvu polychaos_setsizetarget pc np 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 45 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 uAf n average
30. form random numbers Since v2 5 the library uses the Mersenne Twister pseudo random number generator This generator is from M Matsumoto and T Nishimura Mersenne Twister A 623 dimensionally equidistributed uniform pseudorandom number generator ACM Trans on Modeling and Com puter Simulation Vol 8 No 1 January pp 3 30 1998 2 2 Methods In this section we give an overview of the methods which are available in the randvar class 2 2 1 Overview The figure 2 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 L options Methods value randvar_getvalue rv options randvar_getlog rv Destructor randvar_destroy rv Static methods rvlist randvar_tokens nbrv randvar_size Figure 2 3 Outline of the methods of the randvar class 2 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 functions The randvar class provides the following functions 13 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
31. haos_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 54 Building the manual file javaHelp in C tbxnisp help en_US Please wait building this can take a while Generating loader sce The following script presents the messages which are generated when the loader of the toolbox 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 gt 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 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 do
32. ion is not provided by Scilab but is provided by the Distfun module Indeed the distfun_normpdf function of the Distfun module computes the Normal probability 23 Empirical histogram of X1 3 a o o 3 a Figure 3 2 Monte Carlo Sampling Normal random variable Empirical histogram of X2 0 0 2 0 2 1 2 2 2 3 24 2 5 2 6 2 7 2 8 2 9 3 0 Figure 3 3 Monte Carlo Sampling Uniform random variable 24 distribution function In order to install this module we can run the atomsInstall function as in the following script atomsInstall distfun 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 3 4 y Empirical Exact ee lu 0 05 4 3 n ot o 3 a Figure 3 4 Monte Carlo Sampling Histogram and exact distribution functions for the fir
33. m 0 45 0 40 0 35 0 30 0 25 7 S a 0 20 7 0 15 7 0 10 7 0 05 7 0 00 T T T m T T T T T T y T T ee T T 1 3 2 1 0 1 2 3 4 5 6 7 Figure 4 3 Product function Histogram of the output on a LHS design with 10000 experiments We may explore the following topics e Perform the same analysis where the variable X 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 4 2 2 A note on performance In this section we emphasize vectorization which can be used to improve the performance of a script when we compute the output of a function on a given sampling 43 In order to use vectorization the core feature that we used in the Exemple is the use of the elementwise multiplication denoted by In the Exemple function below reproduced here for simplicity the input x is a np by 2 matrix of doubles where np is the number of experiments and y is a np by 1 matrix of doubles function y Exemple x y 1 x 1 x 2 endfunc
34. mtarget 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_getinvguantile cpp Compilation of sci_polychaos_buildsample cpp Compilation of sci_polychaos_getoutput cpp Compilation of sci_polychaos_getguantile cpp Compilation of sci_polychaos_getguantwilks cpp Compilation of sci_polychaos_getsample cpp 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_polyc
35. 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 Variableuflw uNormale 1 0 0 5 my_handle children y_label text Variable 2 Uniforme 2 0 3 0 28 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 3 8 Latin Hypercube Sampling The first variable is Normal the second variable is Uniform 29 The following script allows to plot the histogram of the two variables which are presented in figures 3 9 and 3 10 Plot Var Hi my_handle scf clf my_handle reset histplot 50 sampling 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 1 0 0 5 0 0 0 5 1 0 1 5 2 0 2 5 3 0 Figure 3 9 Latin Hypercube Sampling Normal random variable We can use the mean an
36. n 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 and Scilab It can be used as a meta model The figure 4 1 presents the most commonly used methods available in the polychaos class More methods are presented in figure 4 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 2 2 2 4 2 Examples In this section we present to examples of use of the polychaos class 4 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 38 Constructors pc polychaos_ new file pc polychaos_ new srv ny pc polychaos
37. n 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 27 Empirical histogram of X2 Empirical Exact Figure 3 7 Monte Carlo Sampling Histogram and exact distribution functions for the second variable 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 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 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 3 8 my_handle scf clf
38. n appendix in section A 2 The ATOMS component is the Scilab tool which allows to search download install and load toolboxes ATOMS comes with Scilab v5 2 The 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 64 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 reguires to run the atomsInstall function then to restart Scilab 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 atomsInstall function can be used to install the toolbox for all the users of this computer Then we restart Scilab and the NISP toolbox is automaticall loaded T 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
39. on 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 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_setrandvar_save cpp sci_setrandvar_buildsample cpp sci_polychaos_new cpp sci_polychaos_destroy cpp 53 Compilation of sci_polychaos_tokens cpp Compilation of sci_polychaos_size cpp Compilation of sci_polychaos_setdegree cpp Compilation of sci_polychaos_getdegree cpp Compilation of sci_polychaos_freememory cpp 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_freeme
40. on 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 39 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_getmultind pc polychaos_getlog pc polychaos_getinvguantile pc threshold polychaos_getindextotal pc polychaos_getindexfirst pc ny polychaos_getdimoutput pc nx polychaos_getdiminput pc p polychaos_getdimexp pc no polychaos_getdegree pc polychaos_getcovariance pc polychaos_getcorrelation pc polychaos_getanova pc polychaos_generatecode pc filename funname polychaos_computeoutput pc polychaos_computeexp L pc srv method polychaos_computeexp pc pc2 invalue varopt polychaos_buildsample pc type np order Figure 4 2 More methods from the polychaos class 40 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 setr
41. p 10 setrandvar_buildsample srv Lhs np sampling setrandvar_getsample srv scf plot sampling 1 sampling 2 bo xtitle LHS Design X1 X2 Add the cuts cut linspace 0 1 np 1 for 1 np 1 plot cut i cut i 0 1 end for sj zi npo 1 plot 0 1 cut i cut i end setrandvar_destroy srv 31 The previous script produces the figure 3 11 LHS Design O 0 0 0 1 0 2 0 3 0 4 0 5 0 6 0 7 0 8 0 9 1 0 x1 Figure 3 11 Latin Hypercube Sampling Computed with setrandvar_buildsample and the Lhs option The LhsMaxMin sampling provided by the setrandvar_buildsample function tries to max imize the minimum distance between the points in the sampling The ntry parameter is the number of random points generated before the best is accepted in the sampling np 10 ntry 100 setrandvar_buildsample srv LhsMaxMin np ntry sampling setrandvar_getsample srv The previous script produces the figure 3 12 On the other hand the nisp_buildlhs function produces a more classical LHS sampling where the points are randomly picked within their cells n 5 s 2s sampling nisp_buildlhs s n scf plot sampling 1 sampling 2 bo Add the cuts cut linspace 0 1 n 1 for i 1 n 1 plot cut i cut i 0 1
42. ros 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 histogram 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 2 5 and 2 6 18 Empirical histogram Uniform variable Ee amp 08 0 6 7 0 47 0 0 2 0 2 1 2 2 2 3 2 4 2 5 2 6 2 7 2 8 2 9 3 0 Figure 2 5 The histogram of a Uniform random variable with 1000 samples 19 Empirical histogram Log Uniform variable 0 18 0 12 7
43. st 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 3 5 3 2 2 A Monte Carlo design with 2 variables In this section we create a Monte Carlo design with 2 variables 25 Empirical histogram of X2 Figure 3 5 Monte Carlo Sampling Histogram and exact distribution functions for the second variable 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 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 nisp_initseed O rvi randvar_new Normale 1 0 0 5 rv2 randvar_new Exponentielle 5 Definition d un
44. t of random 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 setrand var_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 i np 5000 setrandvar_buildsample srvu MonteCarlo np sampling setrandvar_getsample srvu Check sampling of random variable 1 mean sampling 1 Expectation 1 Check sampling of random variable 2 mean sampling 2 Expectation 2 5 ii scf histplot 50 sampling 1 xtitle Empirical histogram of X1 scf histplot 50 sampling 2 xtitle Empirical histogram of X2 Clean up setrandvar_destroy srvu randvar_destroy rvul randvar_destroy rvu2 The previous script produces the following output gt mean sampling 1 Expectation 1 ans 1 0064346 gt mean sampling 2 Expectation 2 5 ans 2 5030984 The prevous script also produces the figures 3 2 and 3 3 We may now want to add the exact distribution to these histograms and compare The Normal distribution funct
45. t 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 Mean f expected y f 1n mymean mu mprint Std dev u f expected ff in mysigma sigma mprintf Variancey Af Cexpectedy A4f 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 0 988194 expected 1 000000 Std dev 0 505186 expected 0 500000 Variance 0 255213 expected 0 250000 The previous script also produces the figure 2 4 2 3 2 Variable transformations In this section we present the transformation of a random variable with given distribution into a variable with another distribution Then we present some of the many the transformations which 16 Histogram of X 1 0 0 5 0 0 0 5 1 0 1 5 2 0 2 5 3 0 Figure 2 4 The histogram of a Normal random variable with
46. tements in the Scilab console help randvar The previous statements opens the Help Browser 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 See in the section A 1 for details on the internal structure of the toolbox File Help Browser gt l S NISP Toolbox 1 NISP Toolbox A Table of Contents i yd Ew nisp Functions to configure the NISP library LS eted overview An overview of the NISP toolbox polychaos A class to manage a Polynomial Chaos expansion randvar x randvar A class to manage a random variable ak A setrandvar A class to manage a set of random variables Figure 1 2 The NISP inline help 1 4 Installing the toolbox from ATOMS 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 In this section we present the method to install NISP from ATOMS The installation of the toolbox from the sources in presented i
47. th this class 3 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 21 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 setrandvar_buildsample srv name np 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 7 Figure 3 1 Outline of the methods of the setrandvar class 22 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 results can be reproduced The setrandvar_new function is used to create a new se
48. the polychaos_setsizetarget function ny 1 pc polychaos_ new srvx ny np setrandvar_getsize srvx polychaos_setsizetarget pc np 41 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 then transferred to the polynomial pc by the polychaos_settarget 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_computeexp 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 you Uhsf n average mprintf Varianceyyuu sf n var mprintf Indice de sensibilite du ler ordre n mprintf GuuuVariable X1 4f n pol
49. tion The elementwise multiplication allows to multiply the two first columns of x and sets the result into the output y in one single statement Since Scilab uses optimized numerical libraries this allows to get the best performance in most situations In the previous section we have shown that we can compute the output of the Exemple function in one single call to the function outputdata Exemple inputdata This call allows to produce all the outputs as fast as possible and is the recommended method The reason is that the previous script lets Scilab perform computations with large matrices In fact there is another slower method to perform the same computation We 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 then transferred to the polynomial pc by the polychaos_settarget function This is slow for k 1 np inputdata setrandvar_getsample srvu k outputdata Exemple inputdata mprintf Experiment d inputy f 4f uoutputu u f n 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
50. ychaos_getindexfirst pc 1 mprintf Guu Variable X2 4f n polychaos_getindexfirst pc 2 mprintf Indice de sensibilite Totale n mprintf Guu uVariable X1 4f n polychaos_getindextotal pc 1 mprintf Guu Variable X2 f n polychaos_getindextotal pc 2 The previous script produces the following output Mean 1 750000 Variance 1 000000 Indice de sensibilite du ler 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 vul randvar_destroy vu2 randvar_destroy vx1 randvar_destroy vx2 setrandvar_destroy srvu setrandvar_destroy srvx 42 Prior to destroying the objects we can inquire a little more about the density of the output of the chaos polynomial In the following script 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 4 3 Product function Empirical Histogra
Download Pdf Manuals
Related Search