[ot-users] TR: SpaceFillingC2 speed

Anita Laera anita.laera87 at gmail.com
Thu Jun 15 11:35:39 CEST 2017


This is the script:

from openturns import *


class MomentBeam(OpenTURNSPythonFunction):

    def __init__(self):
        OpenTURNSPythonFunction.__init__(self, 2, 1)
        self.iteration_results = []

    def _exec(self, x):
        P = x[0]
        L = x[1]
        result = P*L
        return [result]

model = MomentBeam()
bending_moment = NumericalMathFunction(model)

force_distribution = Normal(50, 5.0)
L_distribution = Normal(10, 1.5)

collectionMarginals = [force_distribution, L_distribution]
copuleNormal = IndependentCopula(2)

inputDistribution = ComposedDistribution(collectionMarginals, copuleNormal)
inputRandomVector = RandomVector(inputDistribution)
outputVariableOfInterest = RandomVector(bending_moment, inputRandomVector)

threshold = 650
myEvent = Event(outputVariableOfInterest, Greater(), threshold)

optim_solver = AbdoRackwitz()

epsilon = [0.5] * 2
steps = ConstantStep(epsilon)
gradient = CenteredFiniteDifferenceGradient(steps,
bending_moment.getEvaluation())
bending_moment.setGradient(gradient)

meanVector = inputRandomVector.getMean()
myAlgoFORM = FORM(optim_solver, myEvent, meanVector)
myAlgoFORM.run()


I just read the reply from Regis, thank you.
I don't want to change the value in the ResourceMap. I would like to assign
the gradient after creating the function, but if I do as in my script (the
same as in your email) I get these values for the input variables:

input_variables = [*[50.0, 10.0], [50.00001, 10.0], [49.99999, 10.0],
[50.0, 10.00001], [50.0, 9.99999]*, [54.61538461609919,
12.076923076536275], [54.615394616099195, 12.076923076536275],
[54.61537461609919, 12.076923076536275], [54.61538461609919,
12.076933076536275], [54.61538461609919, 12.076913076536275],
[54.65186974033122, 11.893340613204874], [54.651879740331225,
11.893340613204874], [54.65185974033122, 11.893340613204874],
[54.65186974033122, 11.893350613204873], [54.65186974033122,
11.893330613204874], [54.60372712296079, 11.903939968340785], [50.0, 10.0]]

while, if I change the ResourceMap setting the default epsilon to 0.5, then
I get:

input_variables = [[*50.0, 10.0], [50.5, 10.0], [49.5, 10.0], [50.0, 10.5],
[50.0, 9.5]*, [54.61538461538462, 12.076923076923075], [55.11538461538462,
12.076923076923075], [54.11538461538462, 12.076923076923075],
[54.61538461538462, 12.576923076923075], [54.61538461538462,
11.576923076923075], [54.65186973899018, 11.893340613512574],
[55.15186973899018, 11.893340613512574], [54.15186973899018,
11.893340613512574], [54.65186973899018, 12.393340613512574],
[54.65186973899018, 11.393340613512574], [54.60372712265884,
11.903939968405464], [50.0, 10.0]]


2017-06-15 11:18 GMT+02:00 Julien Schueller | Phimeca <schueller at phimeca.com
>:

> Hello Anita,
>
>
> Are you also using openturns from conda on osx ?
>
> Could you show us the script about your default epsilon ? That should work
> even without loading the xml defaults.
>
>
> j
> ------------------------------
> *De :* Anita Laera <anita.laera87 at gmail.com>
> *Envoyé :* jeudi 15 juin 2017 10:56:15
> *À :* Julien Schueller | Phimeca
> *Cc :* roy; haddad at imacs.polytechnique.fr; users
> *Objet :* Re: [ot-users] SpaceFillingC2 speed
>
> I have the same message every time I start a calculation.
>
> Also, when I specify a certain value for the epsilon of the centered
> gradient I want to use for AbdoRackwitz() in FORM, it keeps on using the
> default epsilon 1e-5 (as in the ResourceMap).
>
>
> 2017-06-15 10:44 GMT+02:00 Julien Schueller | Phimeca <
> schueller at phimeca.com>:
>
>> Hi @roy
>>
>>
>> The message: "WRN - The configuration file has not been found, using
>> default parameters. "
>>
>> Is due to an error of the xml configuration loading code specific to osx.
>>
>> I tried to debug it once, the openturns.conf file was really in the
>> correct location though.
>>
>> Sofiane, do you have this message when you compile on osx box ?. I wonder
>> if it's related to conda.
>>
>>
>> j
>>
>>
>>
>>
>>
>>
>>
>>
>> ------------------------------
>> *De :* users-bounces at openturns.org <users-bounces at openturns.org> de la
>> part de roy <roy at cerfacs.fr>
>> *Envoyé :* jeudi 15 juin 2017 10:15:39
>> *À :* D. Barbier
>> *Cc :* users
>> *Objet :* Re: [ot-users] SpaceFillingC2 speed
>>
>> Hello Denis,
>>
>> Indeed now OT is faster using ot.Sample(sample).
>>
>> Regarding numba, it has to be pure python and not numpy for it to work
>> efficiently.
>>
>> import numpy as np
>> import timeit
>> import openturns as ot
>> from numba import jit, njit
>>
>>
>> def discrepancy(sample):
>>     n_sample = len(sample)
>>     dim = sample.shape[1]
>>
>>     abs_ = abs(sample - 0.5)
>>     disc1 = np.sum(np.prod(1 + 0.5 * abs_ - 0.5 * abs_ ** 2, axis=1))
>>
>>     prod_arr = 1
>>     for i in range(dim):
>>         s0 = sample[:, i]
>>         prod_arr *= (1 +
>>                      0.5 * abs(s0[:, None] - 0.5) + 0.5 * abs(s0 - 0.5) -
>>                      0.5 * abs(s0[:, None] - s0))
>>     disc2 = prod_arr.sum()
>>
>>     c2 = (13 / 12) ** dim - 2 / n_sample * disc1 + 1 / (n_sample ** 2) *
>> disc2
>>     return np.sqrt(c2)
>>
>> @jit
>> def discrepancy_numba(sample):
>>     n_sample = len(sample)
>>     dim = sample.shape[1]
>>
>>     abs_ = abs(sample - 0.5)
>>     disc1 = np.sum(np.prod(1 + 0.5 * abs_ - 0.5 * abs_ ** 2, axis=1))
>>
>>     prod_arr = 1
>>     for i in range(dim):
>>         s0 = sample[:, i]
>>         prod_arr *= (1 +
>>                      0.5 * abs(s0[:, None] - 0.5) + 0.5 * abs(s0 - 0.5) -
>>                      0.5 * abs(s0[:, None] - s0))
>>     disc2 = prod_arr.sum()
>>
>>     c2 = (13 / 12) ** dim - 2 / n_sample * disc1 + 1 / (n_sample ** 2) *
>> disc2
>>     return np.sqrt(c2)
>>
>> @njit
>> def discrepancy_faster_numba(sample):
>>     disc1 = 0
>>     n_sample = len(sample)
>>     dim = sample.shape[1]
>>
>>     for i in range(n_sample):
>>         prod = 1
>>         for item in sample[i]:
>>             sub = abs(item - 0.5)
>>             prod *= 1 + 0.5 * sub - 0.5 * sub ** 2
>>         disc1 += prod
>>
>>     disc2 = 0
>>     for i in range(n_sample):
>>         for j in range(n_sample):
>>             prod = 1
>>             for k in range(dim):
>>                 a = 0.5 * abs(sample[i,k] - 0.5)
>>                 b = 0.5 * abs(sample[j,k] - 0.5)
>>                 c = 0.5 * abs(sample[i,k] - sample[j,k])
>>                 prod *= 1 + a + b - c
>>             disc2 += prod
>>
>>     c2 = (13 / 12) ** dim - 2 / n_sample * disc1 + 1 / (n_sample ** 2) *
>> disc2
>>     return np.sqrt(c2)
>>
>>
>> sample = np.random.random_sample((500, 2))
>> ot_sample = ot.Sample(sample)
>> print(discrepancy(sample))
>> print(discrepancy_numba(sample))
>> print(discrepancy_faster_numba(sample))
>> print(ot.SpaceFillingC2().evaluate(sample))
>>
>> print('Function time: ', timeit.repeat('discrepancy(sample)',
>> number=500, repeat=4, setup="from __main__ import discrepancy, sample"))
>> print('numba time: ', timeit.repeat('discrepancy_numba(sample)',
>> number=500, repeat=4, setup="from __main__ import discrepancy_numba,
>> sample"))
>> print('Fast numba time: ', timeit.repeat('discrepancy_faster_numba(sample)',
>> number=500, repeat=4, setup="from __main__ import discrepancy_faster_numba,
>> sample"))
>> print('OT time: ', timeit.repeat('ot.SpaceFillingC2().evaluate(ot_sample)',
>> number=500, repeat=4, setup="from __main__ import ot_sample, ot"))
>>
>>
>>
>> [34m [1mWRN - The configuration file has not been found, using default
>> parameters. [0m      #### IF YOU HAPPEN TO KNOW HOW TO REMOVE THIS BY THE
>> WAY
>> 0.0181493670249
>> 0.0181493670249
>> 0.018149367024149737
>> 0.018149367024149737
>> Function time:  [4.525451728957705, 4.541200206964277,
>> 4.4143504980020225, 4.56408092204947]
>> numba time:  [4.3976798499934375, 4.876463262015022, 5.385470865992829,
>> 5.138608552981168]
>> Fast numba time:  [0.6634743280010298, 0.6538278009975329,
>> 0.7077985780197196, 0.6579875709721819]
>> OT time:  [0.7988348260405473, 0.7220299079781398, 0.7797102630138397,
>> 0.7526425909600221]
>> [Finished in 53.8s]
>>
>>
>> So using numba is here again faster. Even if I use a large sample (1000)
>> numba is slightly faster.
>>
>>
>> Sincerely,
>>
>> *Pamphile ROY*
>> Chercheur doctorant en Quantification d’Incertitudes
>> CERFACS - Toulouse (31) - France
>> +33 (0) 5 61 19 31 57 <+33%205%2061%2019%2031%2057>
>> +33 (0) 7 86 43 24 22 <+33%207%2086%2043%2024%2022>
>>
>>
>>
>> Le 14 juin 2017 à 23:20, D. Barbier <bouzim at gmail.com> a écrit :
>>
>> Hello Pamphile,
>>
>> The problem is that your sample case is small, so the conversion from
>> a numpy array into an OT Sample has a significant cost.
>> If you rerun your benchmark on
>>  otsample = ot.Sample(sample)
>> (or directly generate a random sample with OT), you will see that our
>> version is much faster.
>>
>> BTW I was intrigued by your results with numba, but could not achieve
>> the same speedup, my gain is almost negligible.  Can you please show
>> your test case with numba?  Did you use a GPU?
>>
>> Regards,
>> Denis
>>
>> 2017-06-14 10:32 GMT+02:00 roy <roy at cerfacs.fr>:
>>
>> Hi,
>>
>> Thanks for the feedback, indeed that could explain the behaviours.
>>
>>
>> Pamphile ROY
>> Chercheur doctorant en Quantification d’Incertitudes
>> CERFACS - Toulouse (31) - France
>> +33 (0) 5 61 19 31 57 <+33%205%2061%2019%2031%2057>
>> +33 (0) 7 86 43 24 22 <+33%207%2086%2043%2024%2022>
>>
>>
>>
>> Le 14 juin 2017 à 10:15, HADDAD Sofiane <sofiane_haddad at yahoo.fr> a
>> écrit :
>>
>> Hi,
>>
>> It also depends on sample size
>>
>> With sample's size=1000, I get this :
>>
>> 0.00975831343631
>> 0.009758313432154839
>> Function time:  [19.408187157008797, 21.296883990988135,
>> 19.92589810100617]
>> OT time:  [4.125010760006262, 4.1429947539872956, 4.138353090995224]
>>
>> For small samples, maybe we spend more time for the generation of small
>> objects than the evaluation itself
>>
>> Regards,
>> Sofiane
>>
>>
>> Le Mercredi 14 juin 2017 0h22, D. Barbier <bouzim at gmail.com> a écrit :
>>
>>
>> On 2017-06-13 12:01 GMT+02:00 roy wrote:
>>
>> Hi everyone,
>>
>> I was playing with Centered discrepancy and wrote my function before I saw
>> the class SpaceFillingC2.
>> There is no issue except that I get 2x speedup with my python version.
>> There
>> might be room for improvement as I can even get a 10x on my version using
>> numba.
>>
>> [...]
>>
>> Hello Pamphile,
>>
>> I will have a look, thanks a lot for your feedback.
>> Regards,
>>
>> Denis
>>
>> _______________________________________________
>> OpenTURNS users mailing list
>> users at openturns.org
>> http://openturns.org/mailman/listinfo/users
>> <http://smtp.datacenter.ibo.fr:32224/?dmVyPTEuMDAxJiY2MWQ5YjA4OTI1NmIxYWY0Nj01OTQyNEJCMV8zMzk1OF8xMjc2MV8xJiYyMGNjMjM5Njg2ZjM1Nzc9MTIzMyYmdXJsPWh0dHAlM0ElMkYlMkZvcGVudHVybnMlMkVvcmclMkZtYWlsbWFuJTJGbGlzdGluZm8lMkZ1c2Vycw==>
>>
>>
>>
>>
>>
>> _______________________________________________
>> OpenTURNS users mailing list
>> users at openturns.org
>> http://openturns.org/mailman/listinfo/users
>> <http://smtp.datacenter.ibo.fr:32224/?dmVyPTEuMDAxJiY2MWQ5YjA4OTI1NmIxYWY0Nj01OTQyNEJCMV8zMzk1OF8xMjc2MV8xJiYyMGNjMjM5Njg2ZjM1Nzc9MTIzMyYmdXJsPWh0dHAlM0ElMkYlMkZvcGVudHVybnMlMkVvcmclMkZtYWlsbWFuJTJGbGlzdGluZm8lMkZ1c2Vycw==>
>>
>>
>
> _______________________________________________
> OpenTURNS users mailing list
> users at openturns.org
> http://openturns.org/mailman/listinfo/users
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://openturns.org/pipermail/users/attachments/20170615/0c6cc827/attachment.html>


More information about the users mailing list