wiki:ExampleOldWrapperScilab

Wrapper Scilab (old version)

Prerequisite

  • Scilab
  • OpenTURNS (this tutorial was done with the 0.10.0 version) installed in the /usr folder (Tutorial made with Ubuntu 7.10).
  • Make sure you have automake, build-essential and libtool:
    aptitude install build-essential libtool automake
    

Modification of the source files

  • Copy the folder /usr/share/openturns/WrapperTemplate/ in your home folder for example
  • Change the WrapperTemplate folder to ScilabWrapper

  • Go to the ScilabWrapper folder
  • Save wcode.f as ScilabWrappedCode.f :
    c                                              -*- Fortran -*-
    c
    c   @file  ScilabWrappedCode.f
    c   @brief The code that actually does the computation
    c 
    c   (C) Copyright 2005-2007 EDF-EADS-Phimeca
    c 
    c   This library is free software; you can redistribute it and/or
    c   modify it under the terms of the GNU Lesser General Public
    c   License as published by the Free Software Foundation; either
    c   version 2.1 of the License.
    c 
    c   This library is distributed in the hope that it will be useful
    c   but WITHOUT ANY WARRANTY; without even the implied warranty of
    c   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
    c   Lesser General Public License for more details.
    c 
    c   You should have received a copy of the GNU Lesser General Public
    c   License along with this library; if not, write to the Free Software
    c   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
    c 
    c   @author: $LastChangedBy: dutka $
    c   @date:   $LastChangedDate: 2007-10-04 12:23:01 +0200 (jeu, 04 oct 2007) $
    c   Id:      $Id: wcode.f 554 2007-10-04 10:23:01Z dutka $
    c
    
          SUBROUTINE SCILAB(R,F,G)
          IMPLICIT NONE
    
    c     INPUT PARAMETERS
          REAL*8 R,F
    
    c     OUTPUT PARAMETERS
          REAL*8 G
          REAL*8 RES
    
          INTEGER STATUS
    
          CALL WRITEINPUT(R,F)
          STATUS = SYSTEM('./RunScilab.sh')
          CALL GETRESULT(RES)
    
          G = RES
    
          RETURN
          END
    
          SUBROUTINE WRITEINPUT(I1,I2)
          IMPLICIT NONE
    
    c     INPUT PARAMETERS
          REAL*8  I1,I2
    
          OPEN(10,file='input.txt',status='UNKNOWN')
    
          WRITE(10,20) 'R= ',I1
          WRITE(10,20) 'F= ',I2
    
          CLOSE(10)
    
     20   FORMAT(A,F12.3,A)
    
          RETURN
          END
    
          SUBROUTINE GETRESULT(RES)
          IMPLICIT NONE
    
    c     OUTPUT PARAMETERS
          REAL*8  RES
    
    c     LOCAL PARAMETERS
          INTEGER NB_LIGNS, I
          CHARACTER*256 STRING
          REAL*8 TEMP
    
          NB_LIGNS = 0
          OPEN(10,file='output.txt',status='UNKNOWN')
          DO I=1,NB_LIGNS
             READ(10,'(A)') STRING
          ENDDO
          READ(10,*) STRING, RES
          CLOSE(10)
    
          RETURN
          END
    
  • Save wrapper.c as ScilabWrapper.c :
    /*                                               -*- C -*- */
    /**
     *  @file  ScilabWrapper.c
     *  @brief The wrapper adapts the interface of OpenTURNS and of the wrapped code
     *
     *  (C) Copyright 2005-2007 EDF-EADS-Phimeca
     *
     *  This library is free software; you can redistribute it and/or
     *  modify it under the terms of the GNU Lesser General Public
     *  License as published by the Free Software Foundation; either
     *  version 2.1 of the License.
     *
     *  This library is distributed in the hope that it will be useful
     *  but WITHOUT ANY WARRANTY; without even the implied warranty of
     *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
     *  Lesser General Public License for more details.
     *
     *  You should have received a copy of the GNU Lesser General Public
     *  License along with this library; if not, write to the Free Software
     *  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
     *
     *  @author: $LastChangedBy: dutka $
     *  @date:   $LastChangedDate: 2007-10-05 23:15:55 +0200 (ven, 05 oct 2007) $
     *  Id:      $Id: wrapper.c 556 2007-10-05 21:15:55Z dutka $
     */
    
    #include <stdio.h>
    #include "WrapperInterface.h"
    #include "WrapperCommon.h"
    
    /* WARNING : Please read the following lines
     *
     * In this program, we make the assumption that the end user wishes to
     * call a function (aka NumericalMathFunction) named "wcode" (wcode stands
     * for wrapped code).
     * In order to individualize the wrapper to the user's needs, we encourage
     * you (as the developer of this wrapper) to renamed any occurence of "wcode"
     * (in either case) to the real name of the function.
     * It will also avoid any confusion with other "wcode"s written by other entities
     * or developers.
     */
    
    
    /* If you plan to link this wrapper against a FORTRAN library,
     * remember to change the name WCODE for your actual FORTRAN subroutine name.
     * Otherwise you can ignore these lines.
     *
     * Remember that FORTRAN passes its arguments by reference, not by value as C and C++
     * usually do. So you need to pass the pointer to the arguments rather than its value.
     * This is true for single values (integers, reals, etc.) but not for arrays that are
     * already pointers in the C/C++ environment. Those ones can directly be passed as
     * "values" though they are pointers indeed.
     * Be careful that C and C++ arrays start from 0 and FORTRAN from 1!
     * Be also very careful with the size of the value your plan to pass. Integers in C/C++
     * are not INTEGER*8 in many cases. Float or doubles do not necessarily match REAL*4
     * or REAL*8 in FORTRAN.
     *
     * FORTRAN gives no clue for preventing const values from being altered. So you need to
     * protect them by copying them before calling the FORTRAN subroutine if this import
     * to you.
     *
     * Summary: there are only exceptions to the rule and you need to
     * know exactly what you are doing! You may be desappointed at first, but it will keep
     * you away from segmentation violation and other similar fantasies. ;-)
     */
    #define SCILAB_F77 F77_FUNC(scilab,SCILAB)
    #ifdef __cplusplus
    extern "C"  /* prevent C++ name mangling */
    #endif
    void SCILAB_F77(double * R, double * F, double * G);
    
    
    /* NOTE: all the 'dbg_' function are for debugging needs */
    
    
    /*
     *  This is the declaration of function named 'SCILAB' into the wrapper.
     */
    
    #ifdef __cplusplus
    extern "C" {
    #endif
    
      /*
    *********************************************************************************
    *                                                                               *
    *                             SCILAB function                                    *
    *                                                                               *
    *********************************************************************************
    */
    
    
      /* The wrapper information informs the NumericalMathFunction object that loads the wrapper of the
       * signatures of the wrapper functions. In particular, it hold the size of the input NumericalPoint
       * (inSize_) and of the output NumericalPoint (outSize_).
       * Those information are also used by the gradient and hessian functions to set the correct size
       * of the returned matrix and tensor.
       */
    
      /* The getInfo function is mandatory */
      enum WrapperErrorCode func_getInfo_scilab(void * p_state, struct WrapperInformation * p_info)
      {
        dbg_printEntrance("func_getInfo_scilab");
        dbg_printState("func_getInfo_scilab", p_state);
    
        struct WrapperExchangedData * p_exchangedData = p_state;
    
        p_info->inSize_  = getNumberOfVariables(p_exchangedData, WRAPPER_IN);
        p_info->outSize_ = getNumberOfVariables(p_exchangedData, WRAPPER_OUT);
    
        dbg_printExit("func_getInfo_scilab");
    
        return WRAPPER_OK;
      }
    
    
      /* The state creation/deletion functions allow the wrapper to create or delete a memory location
       * that it will manage itself. It can save in this location any information it needs. The OpenTURNS
       * platform only ensures that the wrapper will receive the state (= the memory location) it works
       * with. If many wrappers are working simultaneously or if the same wrapper is called concurrently,
       * this mechanism will avoid any collision or confusion.
       * The consequence is that NO STATIC DATA should be used in the wrapper OR THE WRAPPER WILL BREAKE
       * one day. You may think that you can't do without static data, but in general this is the footprint
       * of a poor design. But if you persist to use static data, do your work correctly and make use
       * of mutex (for instance) to protect your data against concurrent access. But don't complain about
       * poor computational performance! 
       */
    
      /* The createState function is optional */
      enum WrapperErrorCode func_createState_scilab(void ** p_p_state, const struct WrapperExchangedData * p_exchangedData)
      {
        dbg_printEntrance("func_createState_scilab");
        dbg_printWrapperExchangedData("func_createState_scilab", p_exchangedData);
    
        /* Here we copy the data read in the WML description file into the state */
        copyWrapperExchangedData( (struct WrapperExchangedData **) p_p_state, p_exchangedData);
    
        dbg_printState("func_createState_scilab", *p_p_state);
        dbg_printExit("func_createState_scilab");
    
        return WRAPPER_OK;
      }
    
      /* The deleteState function is optional */
      enum WrapperErrorCode func_deleteState_scilab(void * p_state)
      {
        dbg_printEntrance("func_deleteState_scilab");
        dbg_printState("func_deleteState_scilab", p_state);
    
        /* We free the memory allocated for the state*/
        freeWrapperExchangedData( (struct WrapperExchangedData *) p_state);
    
        dbg_printExit("func_deleteState_scilab");
    
        return WRAPPER_OK;
      }
    
    
    
    
    
    
      /* Any function declared into the wrapper may declare three actual function prefixed with
       * 'init_', 'exec_' and 'finalize_' folowed by the name of the function, here 'wcode'.
       *
       * The 'init_' function is only called once when the NumericalMathFunction object is created.
       * It allows the wrapper to set some internal state, read some external file, prepare the function
       * to run, etc. It takes only one argument, the internal state as created by the 
       *
       * The 'exec_' function is intended to execute what the wrapper is done for: compute an mathematical
       * function or anything else. It takes the internal state pointer as its first argument, the input
       * NumericalPoint pointer as the second and the output NumericalPoint pointer as the third.
       *
       * The 'finalize_' function is only called once when the NumericalMathFunction object is destroyed.
       * It allows the wrapper to flush anything before unloading.
       *
       * Only the 'exec_' function is mandatory.
       */
    
    
      /**
       * Initialization function
       * This function is called once just before the wrapper first called to initialize
       * it, ie create a temparary subdirectory (remember that the wrapper may be called
       * concurrently), store exchanged data in some internal repository, do some
       * pre-computational operation, etc.
       */
      enum WrapperErrorCode func_init_scilab(void * p_state)
      {
        dbg_printEntrance("func_init_scilab");
        dbg_printState("func_init_scilab", p_state);
    
        /* Do something if you need to initialize the wrapper or its environment.
         * If not, you can remove the entire function.
         */
    
        dbg_printExit("func_init_scilab");
    
       return WRAPPER_OK;
      }
    
    
    
    
      /**
       * Execution function
       * This function is called by the platform to do the real work of the wrapper. It may be
       * called concurrently, so be aware of not using shared or global data not protected by
       * a critical section.
       * This function has a mathematical meaning. It operates on one vector (aka point) and
       * returns another vector.
       */
      enum WrapperErrorCode func_exec_scilab(void * p_state, const struct point * inPoint, struct point * outPoint)
      {
        dbg_printEntrance("func_exec_scilab");
        dbg_printState("func_exec_scilab", p_state);
        dbg_printPoint("func_exec_scilab", inPoint);
    
    
    
    
        /* The real computation is here */
    
        /* Call to the FORTRAN subroutine named SCILAB.
         * We have made the assumption that the inputs values won't be touched by the call,
         * so we pass the arguments directly without any protection (copy).
         */
        SCILAB_F77(&(inPoint->data_[0])  /* R */ ,
    	      &(inPoint->data_[1])  /* F */ ,
    	      &(outPoint->data_[0]) /* G */ );
    
    
    
        dbg_printPoint("func_exec_scilab", outPoint);
        dbg_printExit("func_exec_scilab");
    
        return WRAPPER_OK;
      }
    
    
    
    
      /**
       * Finalization function
       * This function is called once just before the wrapper is unloaded. It is the place to flush
       * any output file or free any allocated memory. When this function returns, the wrapper is supposed
       * to have all its work done, so it is not possible to get anymore information from it after that.
       */
      enum WrapperErrorCode func_finalize_scilab(void * p_state)
      {
        dbg_printEntrance("func_finalize_scilab");
        dbg_printState("func_finalize_scilab", p_state);
    
        /* Do something if you need to clean the wrapper or its environment.
         * If not, you can remove the entire function.
         */
    
        dbg_printExit("func_finalize_scilab");
    
        return WRAPPER_OK;
      }
    
    
    
    
    #ifdef __cplusplus
    } /* end extern "C" */
    #endif
    
  • Save wcode.xml.in as ScilabWrapper.xml.in :
    <?xml version="1.0" encoding="ISO-8859-1"?>
    <!-- 
     @file  ScilabWrapper.xml.in
     @brief The description file of the wrapper
    
     (C) Copyright 2005-2007 EDF-EADS-Phimeca
    
     This library is free software; you can redistribute it and/or
     modify it under the terms of the GNU Lesser General Public
     License as published by the Free Software Foundation; either
     version 2.1 of the License.
    
     This library is distributed in the hope that it will be useful
     but WITHOUT ANY WARRANTY; without even the implied warranty of
     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
     Lesser General Public License for more details.
    
     You should have received a copy of the GNU Lesser General Public
     License along with this library; if not, write to the Free Software
     Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
    
     @author: $LastChangedBy: dutka $
     @date:   $LastChangedDate: 2007-10-05 23:15:55 +0200 (ven, 05 oct 2007) $
     Id:      $Id: wcode.xml.in 556 2007-10-05 21:15:55Z dutka $
    -->
    <!DOCTYPE wrapper SYSTEM "@OPENTURNS_WRAPPER_DTD_PATH@/wrapper.dtd">
    
    <wrapper>
      <library>
    
        <!-- The path of the shared object -->
        <path>ScilabWrapper.so</path>
    
    
    
        <!-- This section describes all exchanges data between the wrapper and the platform -->
        <description>
    
          <!-- Those variables are substituted in the files above -->
          <!-- The order of variables is the order of the arguments of the function -->
          <variable-list>
    
            <!-- The definition of a variable -->
            <variable id="R" type="in">
              <comment>Yield strength</comment>
              <unit>MPa</unit>
              <regexp>R=.*</regexp>
              <format>R=%10.5g</format>
            </variable>
    
            <!-- The definition of a variable -->
            <variable id="F" type="in">
              <comment>Axial Force</comment>
              <unit>N</unit>
              <regexp>F=.*</regexp>
              <format>F=%10.5g</format>
            </variable>
    
            <!-- The definition of a variable -->
            <variable id="G" type="out">
              <comment>Limit State</comment>
              <unit>MPa</unit>
              <regexp>G=.*</regexp>
              <format>G=%10.5g</format>
            </variable>
    
          </variable-list>
          
    
    
          <!-- The function that we try to execute through the wrapper -->
          <function provided="yes">scilab</function>
    
          <!-- the gradient is  defined  -->
          <gradient provided="no"></gradient>
    
          <!--  the hessian is  defined  -->
          <hessian provided="no"></hessian>
    
        </description>
    
    
      </library>
    
      <external-code>
        <!-- Those data are external to the platform (input files, etc.)-->
        <data></data>
    
        <wrap-mode type="static-link" state="shared">
          <in-data-transfer mode="arguments" />
          <out-data-transfer mode="arguments" />
        </wrap-mode>
    
        <command># no command</command>
    
      </external-code>
    
    </wrapper>
    
  • Save test_wcode.py as test_ScilabWrapper.py :
    #                                               -*- Python -*-
    #
    #  @file  test_ScilabWrapper.py
    #  @brief A test file for the wrapper code
    #
    #  (C) Copyright 2005-2007 EDF-EADS-Phimeca
    #
    #  This library is free software; you can redistribute it and/or
    #  modify it under the terms of the GNU Lesser General Public
    #  License as published by the Free Software Foundation; either
    #  version 2.1 of the License.
    #
    #  This library is distributed in the hope that it will be useful
    #  but WITHOUT ANY WARRANTY; without even the implied warranty of
    #  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
    #  Lesser General Public License for more details.
    #
    #  You should have received a copy of the GNU Lesser General Public
    #  License along with this library; if not, write to the Free Software
    #  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
    #
    #  @author: $LastChangedBy: dutka $
    #  @date:   $LastChangedDate: 2007-10-05 23:15:55 +0200 (ven, 05 oct 2007) $
    #  Id:      $Id: test_wcode.py 556 2007-10-05 21:15:55Z dutka $
    #
    #  This script imports the OpenTURNS environment in Python, loads
    #  ths wrapper and the wrapped code, and then calls it.
    #
    from openturns import *
    
    f = NumericalMathFunction("ScilabWrapper")
    print "f = ", f
    
    inP = NumericalPoint(f.getInputNumericalPointDimension())
    for i in range(inP.getDimension()):
        inP[i] = 10. * (i + 1)
    print "inP = ", inP
    
    outP = f(inP)
    print "outP = ", outP
    
  • Modify configure.ac :
    #                                               -*- Autoconf -*-
    #
    #  configure.ac
    #
    #  (C) Copyright 2005-2007 EDF-EADS-Phimeca
    #
    #  This library is free software; you can redistribute it and/or
    #  modify it under the terms of the GNU Lesser General Public
    #  License as published by the Free Software Foundation; either
    #  version 2.1 of the License.
    #
    #  This library is distributed in the hope that it will be useful
    #  but WITHOUT ANY WARRANTY; without even the implied warranty of
    #  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
    #  Lesser General Public License for more details.
    #
    #  You should have received a copy of the GNU Lesser General Public
    #  License along with this library; if not, write to the Free Software
    #  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
    #
    #  @author: $LastChangedBy: dutka $
    #  @date:   $LastChangedDate: 2007-10-04 12:23:01 +0200 (jeu, 04 oct 2007) $
    #  Id:      $Id: configure.ac 554 2007-10-04 10:23:01Z dutka $
    #
    #  Process this file with autoconf to produce a configure script.
    #  Set minimal version of Autoconf for configure to work
    #
    AC_INIT(ScilabWrapper.c)
    AC_CONFIG_AUX_DIR(.)
    AM_INIT_AUTOMAKE(ScilabWrapper, 0.0)
    AC_DISABLE_STATIC
    AC_PROG_LIBTOOL
    
    # You need C compiler to build the wrapper
    AC_PROG_CC
    
    # If you don't need FORTRAN, comment the lines below
    AC_PROG_F77
    AC_F77_LIBRARY_LDFLAGS
    AC_F77_WRAPPERS
    
    OT_CHECK_OPENTURNS
    
    AC_OUTPUT([Makefile ScilabWrapper.xml])
    
  • Modify Makefile.am :
    #                                               -*- Makefile -*-
    #
    #  Makefile.am
    #
    #  (C) Copyright 2005-2007 EDF-EADS-Phimeca
    #
    #  This library is free software; you can redistribute it and/or
    #  modify it under the terms of the GNU Lesser General Public
    #  License as published by the Free Software Foundation; either
    #  version 2.1 of the License.
    #
    #  This library is distributed in the hope that it will be useful
    #  but WITHOUT ANY WARRANTY; without even the implied warranty of
    #  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
    #  Lesser General Public License for more details.
    #
    #  You should have received a copy of the GNU Lesser General Public
    #  License along with this library; if not, write to the Free Software
    #  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
    #
    #  @author: $LastChangedBy: dutka $
    #  @date:   $LastChangedDate: 2007-10-05 23:15:55 +0200 (ven, 05 oct 2007) $
    #  Id:      $Id: Makefile.am 556 2007-10-05 21:15:55Z dutka $
    #
    ACLOCAL_AMFLAGS = -I m4
    
    wrapperdir          = $(prefix)/wrappers
    
    wrapper_LTLIBRARIES = ScilabWrapper.la
    ScilabWrapper_la_SOURCES    = ScilabWrapper.c ScilabWrappedCode.f
    ScilabWrapper_la_CPPFLAGS   = $(OPENTURNS_WRAPPER_CPPFLAGS)
    ScilabWrapper_la_LDFLAGS    = -module -no-undefined -version-info 0:0:0
    ScilabWrapper_la_LDFLAGS   += $(OPENTURNS_WRAPPER_LDFLAGS)
    ScilabWrapper_la_LIBADD     = $(OPENTURNS_WRAPPER_LIBS)
    
    # This line is useful only if you compile FORTRAN files
    ScilabWrapper_la_LDFLAGS   += $(FLIBS)
    
    XMLWRAPPERFILE      = ScilabWrapper.xml
    wrapper_DATA        = $(XMLWRAPPERFILE)
    EXTRA_DIST          = $(XMLWRAPPERFILE).in test_ScilabWrapper.py
    

Building the wrapper

  • Run the bootstrap :
    ./bootstrap
    
  • Run :
    ./configure --prefix=/home/%your home folder%/openturns --with-openturns=/usr
    
  • Build the wrapper :
    make
    
  • Install your wrapper in /home/%your home folder%/openturns/wrappers :
    make install
    
  • The wrapper will be seen by OpenTURNS automatically, as it was configured to be placed in a standard folder, where OpenTURNS checks by default the presence of wrappers.
  • To solve a bug of the 0.10.0 version, you also have to copy the file /etc/openturns/openturns.conf to the folder /home/%your home folder%/openturns/etc/.

Running a stochastic analysis using the wrapper

The same physic example as here is used and programmed on a Scilab script, StressedBeam.sci :

// Input file
u=mopen('input.txt','r')
R=mfscanf(u,'%s')
R=mfscanf(u,'%f')
F=mfscanf(u,'%s')
F=mfscanf(u,'%f')
mclose(u);

// Output calculation
G=R-F/(3.14159*100.0) ;

// Output file
u=file('open','output.txt',"unknown")
fprintf(u,'G= %6.9f',G)
file('close',u) 

exit

Some scripts are needed to launch Scilab :

  • RunScilab.sh :
    #!/bin/sh
    scilab -nw -f StressedBeam.sci
    
  • Clean.sh :
    #!/bin/sh
    rm input.txt
    rm output.txt
    

Don't forget to make the scripts executables :

chmod +x Scilab.sh
chmod +x Clean.sh

The analysis is done using the TUI :

from openturns import *
from math import * 

stochasticDimension = 2
RandomGenerator().SetSeed(0)

# Definition of the function to use : Scilab Coupling
LimitState = NumericalMathFunction("ScilabWrapper")
dim = LimitState.getInputNumericalPointDimension()

# Mean
mean = NumericalPoint(dim, 0.0)
mean[0] = 300.0
mean[1] = 75000.0

# Standard deviation
sigma = NumericalPoint(dim, 0.0)
sigma[0] = 30.0
sigma[1] = 5000.0

# Parameters for the lognormal distribution:
MuLog = NumericalPoint(dim, 0)
SigmaLog = NumericalPoint(dim, 0)

MuLog[0] = log(mean[0]*mean[0]/(sqrt(mean[0]*mean[0]+sigma[0]*sigma[0])))
SigmaLog[0] = sqrt(log(1.0+(sigma[0]*sigma[0]/(mean[0]*mean[0]))))
print "MuLog=" , MuLog
print "SigmaLog=", SigmaLog

component = Description(1)
BorneInf = 0.0

# Initialization of the distribution collection:
aCollection = DistributionCollection()

# Create a first marginal : distribution 1D
marginal = LogNormal(MuLog[0], SigmaLog[0], BorneInf)
marginal.setName("Yield strength")
component[0] = "R"
marginal.setDescription(component)
# Graphical output of the PDF
pdfgraph=marginal.drawPDF()
pdfgraph.draw("/tmp", "pdf_R", 640, 480)
# Fill the first marginal of aCollection
aCollection.add( Distribution(marginal, "Yield strength") )

# Create a second marginal : distribution 1D
marginal = Normal(mean[1], sigma[1])
marginal.setName("Traction_load")
component[0] = "F"
marginal.setDescription(component)
# Graphical output of the PDF
pdfgraph=marginal.drawPDF()
pdfgraph.draw("/tmp", "pdf_F", 640, 480)
# Fill the second marginal of aCollection
aCollection.add( Distribution(marginal, "Traction_load") )

# Create a copula : IndependentCopula (no correlation)
aCopula = IndependentCopula(aCollection.getSize())
aCopula.setName("Independent copula")

# Instanciate one distribution object
myDistribution = ComposedDistribution(aCollection, Copula(aCopula))
myDistribution.setName("myDist")

# We create a 'usual' RandomVector from the Distribution
vect = RandomVector(Distribution(myDistribution))

# We create a composite random vector
G = RandomVector(LimitState, vect)

# We create an Event from this RandomVector
myEvent = Event(G, ComparisonOperator(Less()), 0.0, "Event 1")

cv = 0.10
NbSim = 100

myAlgo = MonteCarlo(myEvent)
myAlgo.setMaximumOuterSampling(NbSim)
myAlgo.setBlockSize(1)
myAlgo.setMaximumCoefficientOfVariation(cv)

myAlgo.run()

result = myAlgo.getResult()
probability = result.getProbabilityEstimate()
print "MonteCarlo result=" , result.str()
print "Number of executed iterations =" , result.getOuterSampling()
print "Pf = " , probability
print "CV =" , sqrt(result.getVarianceEstimate())
Last modified 8 years ago Last modified on 07/23/09 14:25:46