let me just post my code: t is the time array and n is also an array. For every value of time t, these operations are performed on the entire array n. Then, n is summed to a scalar which represents the system response at time t. I would like to eliminate this for loop if possible. Chris #### code #### b = 4.7 f = [] n = arange(1, N+1, 1) for t in timearray: arg1 = {'S': ((b/t) + (1J*n*pi/t))} exec('from numpy import *', arg1) tempval = eval(transform, arg1)*((-1)**n) rsum = tempval.real.sum() arg2 = {'S': b/t} exec('from numpy import *', arg2) tempval2 = eval(transform, arg2)*0.5 fval = (exp(b) / t) * (tempval2 + rsum) f.append(fval) #### /code ##### On Thu, May 7, 2009 at 1:04 PM, Chris Colbert wrote:

unfortunately, the actual function being processes is not so simple, and involves evaluating user functions input from the prompt as strings. So i have no idea how to do it in Cython.

Let me look into this broadcasting.

Thanks Josef!

On Thu, May 7, 2009 at 12:56 PM, wrote:

...

On Thu, May 7, 2009 at 12:39 PM, Chris Colbert wrote:

...

suppose i have two arrays: n and t, both are 1-D arrays.

for each value in t, I need to use it to perform an element wise scalar operation on every value in n and then sum the results into a single scalar to be stored in the output array.

Is there any way to do this without the for loop like below:

for val in t_array:

out = (n / val).sum() # not the actual function being done, but you get the idea

broad casting should work, e.g.

(n[:,np.newaxis] / val[np.newaxis,:]).sum()

but it constructs the full product array, which is memory intensive for a reduce operation, if the 1d arrays are large.

another candidate for a cython loop if the arrays are large?

Josef _______________________________________________ Numpy-discussion mailing list Numpy-discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion

its part of a larger program for designing PID controllers. This particular function numerical calculates the inverse laplace transform using riemann sums. The exec statements, from what i gather, allow the follow eval statement to be executed in the scope of numpy and its functions. I don't get how it works either, but it doesnt work without it. I've just about got something working using broadcasting and will post it soon. chris On Thu, May 7, 2009 at 1:37 PM, wrote:

On Thu, May 7, 2009 at 1:08 PM, Chris Colbert wrote:

...

let me just post my code:

t is the time array and n is also an array.

For every value of time t, these operations are performed on the entire array n. Then, n is summed to a scalar which represents the system response at time t.

I would like to eliminate this for loop if possible.

Chris

#### code ####

b = 4.7 f = [] n = arange(1, N+1, 1)

for t in timearray: arg1 = {'S': ((b/t) + (1J*n*pi/t))} exec('from numpy import *', arg1) tempval = eval(transform, arg1)*((-1)**n) rsum = tempval.real.sum() arg2 = {'S': b/t} exec('from numpy import *', arg2) tempval2 = eval(transform, arg2)*0.5 fval = (exp(b) / t) * (tempval2 + rsum) f.append(fval)

#### /code #####

I don't understand what the exec statements are doing, I never use it. what is transform? Can you use regular functions instead or is there a special reason for the exec and eval?

In these expressions ((b/t) + (1J*n*pi/t)), (exp(b) / t) broadcasting can be used.

Whats the size of t and n?

Josef _______________________________________________ Numpy-discussion mailing list Numpy-discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion

On Thu, May 7, 2009 at 2:11 PM, Chris Colbert wrote:

alright I got it working. Thanks!

This version is an astonishingly 1900x faster than my original implementation which had two for loops. Both versions are below:

thanks again!

### new fast code ####

    b = 4.7     n = arange(1, N+1, 1.0).reshape(N, -1)     n1 = (-1)**n     prefix = exp(b) / timearray

    arg1 = {'S': b / timearray}     exec('from numpy import *', arg1)     term1 = (0.5) * eval(transform, arg1)

    temp1 = b + (1J * pi * n)     temp2 = temp1 / timearray     arg2 = {'S': temp2}     exec('from numpy import *', arg2)     term2 = (eval(transform, arg2) * n1).sum(axis=0).real

    f = prefix * (term1 + term2)

    return f

If you don't do code generation and have control over transform, then, I think, it would be more readable to replace the exec and eval by a function call to transform. I haven't found a case yet where eval is necessary, except for code generation as in sympy. Josef

##### old slow code ######     b = 4.7     f = []

    for t in timearray:         rsum = 0.0         for n in range(1, N+1):             arg1 = {'S': ((b/t) + (1J*n*pi/t))}             exec('from numpy import *', arg1)             tempval = eval(transform, arg1)*((-1)**n)             rsum = rsum + tempval.real         arg2 = {'S': b/t}         exec('from numpy import *', arg2)         tempval2 = eval(transform, arg2)*0.5         fval = (exp(b) / t) * (tempval2 + rsum)         f.append(fval)

   return f

On Thu, May 7, 2009 at 1:41 PM, Chris Colbert wrote:

...

its part of a larger program for designing PID controllers. This particular function numerical calculates the inverse laplace transform using riemann sums.

The exec statements, from what i gather, allow the follow eval statement to be executed in the scope of numpy and its functions. I don't get how it works either, but it doesnt work without it.

I've just about got something working using broadcasting and will post it soon.

chris

On Thu, May 7, 2009 at 1:37 PM, wrote:

...

On Thu, May 7, 2009 at 1:08 PM, Chris Colbert wrote:

...

let me just post my code:

t is the time array and n is also an array.

For every value of time t, these operations are performed on the entire array n. Then, n is summed to a scalar which represents the system response at time t.

I would like to eliminate this for loop if possible.

Chris

#### code ####

b = 4.7 f = [] n = arange(1, N+1, 1)

for t in timearray:         arg1 = {'S': ((b/t) + (1J*n*pi/t))}         exec('from numpy import *', arg1)         tempval = eval(transform, arg1)*((-1)**n)         rsum = tempval.real.sum()         arg2 = {'S': b/t}         exec('from numpy import *', arg2)         tempval2 = eval(transform, arg2)*0.5         fval = (exp(b) / t) * (tempval2 + rsum)         f.append(fval)

#### /code #####

I don't understand what the exec statements are doing, I never use it. what is transform? Can you use regular functions instead or is there a special reason for the exec and eval?

In these expressions ((b/t) + (1J*n*pi/t)),  (exp(b) / t) broadcasting can be used.

Whats the size of t and n?

Josef _______________________________________________ Numpy-discussion mailing list Numpy-discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion

_______________________________________________ Numpy-discussion mailing list Numpy-discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion

On Thu, May 7, 2009 at 3:10 PM, Chris Colbert wrote:

the user of the program inputs the transform in a text field. So I have no way of know the function apriori.

that doesn't mean I still couldn't throw the exec and eval commands into another function just to clean things up.

Chris

No, I think this is ok then, it is similar to what sympy.lambdify does. Now you call exec only twice, which might have been the main slowdown in the loop version. In this case, users need to write vectorized transform functions that handle the full n x t array. Josef

On Thu, May 7, 2009 at 2:45 PM, wrote:

...

On Thu, May 7, 2009 at 2:11 PM, Chris Colbert wrote:

...

alright I got it working. Thanks!

This version is an astonishingly 1900x faster than my original implementation which had two for loops. Both versions are below:

thanks again!

### new fast code ####

    b = 4.7     n = arange(1, N+1, 1.0).reshape(N, -1)     n1 = (-1)**n     prefix = exp(b) / timearray

    arg1 = {'S': b / timearray}     exec('from numpy import *', arg1)     term1 = (0.5) * eval(transform, arg1)

    temp1 = b + (1J * pi * n)     temp2 = temp1 / timearray     arg2 = {'S': temp2}     exec('from numpy import *', arg2)     term2 = (eval(transform, arg2) * n1).sum(axis=0).real

    f = prefix * (term1 + term2)

    return f

If you don't do code generation and have control over transform, then, I think, it would be more readable to replace the exec and eval by a function call to transform.

I haven't found a case yet where eval is necessary, except for code generation as in sympy.

Josef

...

##### old slow code ######     b = 4.7     f = []

    for t in timearray:         rsum = 0.0         for n in range(1, N+1):             arg1 = {'S': ((b/t) + (1J*n*pi/t))}             exec('from numpy import *', arg1)             tempval = eval(transform, arg1)*((-1)**n)             rsum = rsum + tempval.real         arg2 = {'S': b/t}         exec('from numpy import *', arg2)         tempval2 = eval(transform, arg2)*0.5         fval = (exp(b) / t) * (tempval2 + rsum)         f.append(fval)

   return f

On Thu, May 7, 2009 at 1:41 PM, Chris Colbert wrote:

...

its part of a larger program for designing PID controllers. This particular function numerical calculates the inverse laplace transform using riemann sums.

The exec statements, from what i gather, allow the follow eval statement to be executed in the scope of numpy and its functions. I don't get how it works either, but it doesnt work without it.

I've just about got something working using broadcasting and will post it soon.

chris

On Thu, May 7, 2009 at 1:37 PM, wrote:

...

On Thu, May 7, 2009 at 1:08 PM, Chris Colbert wrote:

...

let me just post my code:

t is the time array and n is also an array.

For every value of time t, these operations are performed on the entire array n. Then, n is summed to a scalar which represents the system response at time t.

I would like to eliminate this for loop if possible.

Chris

#### code ####

b = 4.7 f = [] n = arange(1, N+1, 1)

for t in timearray:         arg1 = {'S': ((b/t) + (1J*n*pi/t))}         exec('from numpy import *', arg1)         tempval = eval(transform, arg1)*((-1)**n)         rsum = tempval.real.sum()         arg2 = {'S': b/t}         exec('from numpy import *', arg2)         tempval2 = eval(transform, arg2)*0.5         fval = (exp(b) / t) * (tempval2 + rsum)         f.append(fval)

#### /code #####

I don't understand what the exec statements are doing, I never use it. what is transform? Can you use regular functions instead or is there a special reason for the exec and eval?

In these expressions ((b/t) + (1J*n*pi/t)),  (exp(b) / t) broadcasting can be used.

Whats the size of t and n?

Josef _______________________________________________ Numpy-discussion mailing list Numpy-discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion

_______________________________________________ Numpy-discussion mailing list Numpy-discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion

_______________________________________________ Numpy-discussion mailing list Numpy-discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion

_______________________________________________ Numpy-discussion mailing list Numpy-discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion

that's essentially what the eval statement does. On Thu, May 7, 2009 at 4:22 PM, wrote:

On Thu, May 7, 2009 at 3:39 PM, wrote:

...

On Thu, May 7, 2009 at 3:10 PM, Chris Colbert wrote:

...

the user of the program inputs the transform in a text field. So I have no way of know the function apriori.

that doesn't mean I still couldn't throw the exec and eval commands into another function just to clean things up.

Chris

No, I think this is ok then, it is similar to what sympy.lambdify does. Now you call exec only twice, which might have been the main slowdown in the loop version.

In this case, users need to write vectorized transform functions that handle the full n x t array.

Josef

this would be an alternative, which might also work better in a loop, requires numpy.* in local scope

...

...

...

transform = 'sqrt(x)' exec('def fun(x): return ' + transform) from numpy import * fun(5) 2.2360679774997898 fun(np.arange(5)) array([ 0. , 1. , 1.41421356, 1.73205081, 2. ])

fun

Josef _______________________________________________ Numpy-discussion mailing list Numpy-discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion