0

I wrote a function to create a domain to collect the intenisty of light which comes into a probe. When I use a single value for the starting point is fine, but when I set an array of poitns (because I need to describe a surface and the integrate the intensity onto this surface) it raises the ValueError: setting an array element with a sequence.

I tried to verify the shape of the different arrays and it (3,) and (3,), so in theory there is no inconsistence... Here the function

def probe_intensity(cone_FOV, n_circles):
    probeI = 0
    # probe location information
    #probe_angle = 5.0/4.0*np.pi
    #tvar=np.linspace(4.09,4.186,30)
    ray_cone_angle = cone_FOV/(2.0*n_circles)
   # T0 = np.array([0, 0, 1])
    r = 6.2
    x = 5.5
    probe_angle = np.linspace(1.2334*np.pi-1,1.2334*np.pi+1,30)
    probe_direction_angle = probe_angle - np.pi
    #R1=np.ones(len(probe_angle))
    # probe direction
    for k in probe_angle:

# getting the estimated cone for a single ray
#ray_cone_angle = cone_FOV/(2.0*n_circles)

# probe position information
#r = 6.2
#x = 5.5

        R0 = np.array([5.5,
               r*np.sin(probe_angle)+0.0001,
               r*np.cos(probe_angle)+0.0001])*0.001

        gamma = probe_direction_angle

        ROT1 = np.matrix([[1, 0, 0],
                  [0, np.cos(gamma), np.sin(gamma)],
                  [0, -np.sin(gamma), np.cos(gamma)]])


# vector T0 before any rotation
        T0 = np.array([0, 0, 1])

# vector T1 - probe axis vector
        T1 = np.array(ROT1*np.matrix(T0).T).flatten()

    # rotation for the cone_half_angle
        for ray_carrier_angle in [(i + 0.5)*ray_cone_angle for i in range(n_circles)]:

            alpha = ray_carrier_angle

            ROT2 = np.matrix([[np.cos(alpha), 0, np.sin(alpha)],
                          [0, 1, 0],
                          [-np.sin(alpha), 0, np.cos(alpha)]])

        # rotating to make a full circle of the cone
            for beta in np.linspace(0, 2*np.pi, int(np.round(n_circles*np.pi))):
                ROT3 = np.matrix([[np.cos(beta), np.sin(beta), 0],
                              [-np.sin(beta), np.cos(beta), 0],
                              [0, 0, 1]])

            # second rotaiton for the cone after the probe direction rotation
                T3 = np.array(ROT1*ROT3*ROT2*np.matrix(T0).T).flatten()
                print R0
                print R0.shape
                print T3
                print T3.shape
                T3=list(T3)
                R0=list(R0)


                ray_solution = ray_trace_solve_ivp(R0, T3, 0.0005)
                t, I = field.integrate_trace(ray_solution, 0.0005)
            #s,e,alpha,I = field.getspectra

        # Add the intensity of the ray to the intensity gathered at the probe. Dot product takes care of the projection
            probeI += I[-1]*np.dot(T0, T3)
            print probeI

return probeI

the other function which is called is

def ray_trace_solve_ivp(R0, T0, optical=False, dt=np.inf, atol=1e-6, rtol=1e-2):
    y0 = np.r_[R0, T0]

    if optical:
        differential_equation = eikonalODE1system_optical
    else:
        differential_equation = eikonalODE1system_physical

    sol = solve_ivp(differential_equation, [0, 0.5], y0,
                    events=limit_functions,
                    dense_output=True,
                    max_step=dt,
                    atol=atol,
                    rtol=rtol)
    return sol

I get the following error message:

  File "/home/tont_fe/anaconda2/lib/python2.7/site-packages/spyder_kernels/customize/spydercustomize.py", line 827, in runfile
    execfile(filename, namespace)

  File "/home/tont_fe/anaconda2/lib/python2.7/site-packages/spyder_kernels/customize/spydercustomize.py", line 102, in execfile
    builtins.execfile(filename, *where)

  File "/home/tont_fe/data/cfd/OH_analysis/OHstar_ray_tracing/oh_ray_trace_candidate_cone_federica.py", line 1109, in 
    print("probe intensity=", probe_intensity(np.pi/20.0, 10))

  File "/home/tont_fe/data/cfd/OH_analysis/OHstar_ray_tracing/oh_ray_trace_candidate_cone_federica.py", line 878, in probe_intensity
    ray_solution = ray_trace_solve_ivp(R0, T3, 0.0005)

  File "/home/tont_fe/data/cfd/OH_analysis/OHstar_ray_tracing/oh_ray_trace_candidate_cone_federica.py", line 690, in ray_trace_solve_ivp
    rtol=rtol)

  File "/home/tont_fe/anaconda2/lib/python2.7/site-packages/scipy/integrate/_ivp/ivp.py", line 456, in solve_ivp
    solver = method(fun, t0, y0, tf, vectorized=vectorized, **options)

  File "/home/tont_fe/anaconda2/lib/python2.7/site-packages/scipy/integrate/_ivp/rk.py", line 96, in __init__
    support_complex=True)

  File "/home/tont_fe/anaconda2/lib/python2.7/site-packages/scipy/integrate/_ivp/base.py", line 120, in __init__
    self._fun, self.y = check_arguments(fun, y0, support_complex)

  File "/home/tont_fe/anaconda2/lib/python2.7/site-packages/scipy/integrate/_ivp/base.py", line 15, in check_arguments
    y0 = y0.astype(dtype, copy=False)
10
  • Fix the indents. Show us (and yourself) exactly where the error occurs. Tell us about the variables at that point (type, shape,dtype, etc). In other words all the useful stuff needed to debug such an error. Commented Nov 2, 2019 at 19:24
  • I think that the error occurs when calling the ray_solve_ivp function, where at the beginning I pack the two arrays in a single array, so when the np.r_(R0,T0) is called. y0 = y0.astype(dtype, copy=False). Could maybe help specifying a dtype in y0? and if yes, which one? Commented Nov 2, 2019 at 20:12
  • The error message should show exactly which line has the problem! Commented Nov 2, 2019 at 20:15
  • I put the error messagee in the post! Commented Nov 2, 2019 at 20:23
  • Have you checked your inputs to the ivp function agaist the documented requirements? type, shape dtype etc? Commented Nov 2, 2019 at 21:17

1 Answer 1

0

So the relevant parts of your code are:

probe_angle = np.linspace(1.2334*np.pi-1,1.2334*np.pi+1,30)
R0 = np.array([5.5,
               r*np.sin(probe_angle)+0.0001,
               r*np.cos(probe_angle)+0.0001])*0.001
T0 = np.array([0, 0, 1])
y0 = np.r_[R0, T0]

Or simplified a bit:

In [96]: R0 = np.array([1, np.linspace(0,1,3), np.linspace(0,3,3)])             
In [97]: R0                                                                     
Out[97]: array([1, array([0. , 0.5, 1. ]), array([0. , 1.5, 3. ])], dtype=object)
In [98]: T0 = np.array([0,0,1])                                                 
In [99]: np.r_[R0,T0]                                                           
Out[99]: 
array([1, array([0. , 0.5, 1. ]), array([0. , 1.5, 3. ]), 0, 0, 1],
      dtype=object)
In [100]: _.astype(float)                                                       
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-100-b12cc607d111> in <module>
----> 1 _.astype(float)

ValueError: setting an array element with a sequence.

Your R0, and then y0 constructor, makes an object dtype array - an array that contains arrays.

But solve_ivp requires a (n,) float array, where n is the number of values that differential_equation produces.

The immediate problem is in creating a numeric array from that y0. Since I don't know what you need or want, I'm not going to propose an fix.

Also I don't know anything about the differential_equation. I don't know how many value it returns. Just one, three? Or does it vary with size of the y0 input?

Step back and work on a simpler problem. Make sure you understand the requirements of solve_ivp, both the function and the y0.

And when you come back with new questions, please make it easier to help!

Sign up to request clarification or add additional context in comments.

3 Comments

So the differential equations gives as output an array of 6 elements. I think I found a way to solve the error anyway. I just take each element of the array R0 (made of multiple arrays ) and I extract an array from this values. Essentially I separate the different arrays and put them in a list. What I want is essentially to have a cycle which produces a results for each pair of arrays R0 and T3. I would like to post here the solution I adopted, so you can tell me if it can be a useful workaround.
for k in probe_angle: R0 = np.array([x, rnp.sin(probe_angle)+0.0001, rnp.cos(probe_angle)+0.0001])*0.001 for i in range(0,len(probe_angle)): R_i=np.array([x,R0[1][i],R0[2][i]]) Rlist.append(R_i)
T3 = np.array(ROT1*ROT3*ROT2*np.matrix(T0).T).flatten() print T3 for idx in range(0,len(probe_angle)): T3_idx=np.array([x,T3[1][idx],T3[2][idx]]) T3list.append(T3_idx) for elem in Rlist: for item in T3list: ray_solution = ray_trace_solve_ivp(elem, item, 0.0005) t, I = field.integrate_trace(ray_solution, 0.0005)

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.