Problem with vectorization

[poncev]

I want to integrate a very complicated function but I have not succeed yet. Here I tried to define a function that exemplifies my difficulties, I guess the problem has to do with vectorization. It'd be extremely hard to redefine my original function to handle vectorizations

import numpy as np
import quadpy

def f(x):
    v = np.array([np.cos(x[1]), np.cos(x[0])])
    w = np.inner(v, np.array([1, 0]))
    return w

domain = quadpy.c2.rectangle_points([0, 1], [0, 1])
scheme = quadpy.c2.get_good_scheme(10)
val = scheme.integrate(f, domain)
print(val)

I get the error: ValueError: shapes (2,28) and (2,) not aligned: 28 (dim 1) != 2 (dim 0)

Any suggestions? Thank you.

[poncev]

I have zeroed in on the first point where the integrate method fails. During the evaluation of the function I want to integrate, I compute a matrix. Let us say t is one of the arguments, a real number, but scheme.integrate passes an array of points of integration, let us say t = [t_0, ..., t_7]. Then I define matrix similar to

R = np.array([[1, 0, 0],
              [0, np.cos(t), np.sin(t)],
              [0, -np.sin(t), np.cos(t)]])

Now, since t is an array, not a real value as I intended, the final matrix is misshaped, we get

array([1, 0, 0],
      [0, array([cos(t_0), ..., cos(t_7)]), array([sin(t_0), ..., sin(t_7)])], 
      ...)

Since the function to integrate has many more terms like this, I wonder whether there is a method to fix this kind of problems efficiently.

Thank you!

[nschloe]

I get the error: ValueError: shapes (2,28) and (2,) not aligned: 28 (dim 1) != 2 (dim 0)

Well, the error is in your own code. The line

    w = np.inner(v, np.array([1, 0]))

fails.

I don't understand your second post. You're wondering why quadpy gives you multiple input points?

[poncev]

Thank you for replying. I also posted the question in stackoverflow; hopefully, it is more clear there. I'm not asking why quadpy gives multiple input points.

The function in my first post is f(x, y) = cos(y). The function works as expected if I pass [x_0, x_1], but I wanted to ilustrate the problems that may arise in the library. I admit that the example is artificial. I know that that line fails in quadpy, but it does not fail when you use the function as originally intended, with two real numbers.

In my second post I wrote a more realistic example. When I defined the function, I was thinking in a matrix with real valued entries (I didn't know I would use quadpy). Quadpy inputs a vector, so the matrix has mixed vector and real valued entries, and when I try to multiply this matrix with other marices or vectors, numpy does not know how to deal with the situation and throws errors. When I input a real value everything goes well.

I want to learn how I should define my functions so that they are compatible with quadpy. I have been able to integrate basic functions with quadpy, but when my functions entail matrices and operations on vectors, everything crashes For example, in my code, should I replace constant*vector by np.tensordot(constant, vector)? The problem is that quadpy replaces the constant by a vector, and the vector by an array of arrays.

I suppose that if I learn how to write vectorized versions of my functions, then the code will be more efficient.

[nschloe]

Perhaps this helps:

https://github.com/sigma-py/quadpy/wiki/Dimensionality-of-input-and-output-arrays

[poncev]

Thank you! I will check next week.

[nschloe]

Closing for lack of feedback. Feel free to reopen at any time!