On Mon, Nov 15, 2010 at 9:39 AM, MJ <mlj.o...@gmail.com> wrote:
> Hi!!
>
> It's Maria again... :D
>
> First of all I have to tell you that the module it's working perfectly
> and I get my matrices and the corresponding derivatives without a
> problem.
I'm glad to hear that.
>
> But...of course...as I told you I'm using this program to model
> optical systems, and therefore one of the main variables (let's say a
> certain x[1], for example, in algopy) I have to use is the index of
> refraction of the different materials light passes through, which is a
> complex number. To test if algopy would work with them I set this
> small example and I got an error.
complex numbers are not (yet) officially supported.
Since there are no tests that include complex numbers,
it is possible that there are bugs in algopy that either may result in an error,
or worse, may result in an incorrect computation.
In conclusion: You can give it a try, but you'll have to be extra careful.
If you have small test problems that I could use to verify the
correctness of the computation,
I'd be glad to incorporate them into algopy.
>
>>>> import numpy
>>>> from numpy import matrix
>>>> import cmath
>>>> from cmath import sin, cos, sinh, cosh
>>>> import algopy
>>>> from algopy import UTPM, zeros
>>>>
>>>> def f(x):
> ... M = zeros((2,2), dtype=x)
> ... M[0,0] = M[1,1] = 1
> ... M[0,1] = sin(x[0]) + cos(x[1])
> ... M[1,0] = (sin(x[0]) + cos(x[1]))/2
> ... return M
> ...
>>>> x = UTPM(zeros((2,2,2)))
>>>> x.data[0,:] = [2+2j, 0+1j]
> Traceback (most recent call last):
> File "<stdin>", line 1, in <module>
> ValueError: setting an array element with a sequence.
>
> Any ideas on how can I introduce this type of numbers in the process??
I had to make some small modifications of algopy to get your example to work.
They will only be available in the next release 0.3.0.
If you want to give it a try now, you can get the bleeding edge
version from https://github.com/b45ch1/algopy ,
for instance on Linux with `$git clone git://github.com/b45ch1/algopy.git`.
You need to initialize the UTPM instance with a complex zeros array.
Below a code that should work
-------------------- start code ----------------------------
import numpy
import algopy
from algopy import UTPM, zeros
def f(x):
M = zeros((2,2), dtype=x)
M[0,0] = M[1,1] = 1
M[0,1] = numpy.sin(x[0]) + numpy.cos(x[1])
M[1,0] = (numpy.sin(x[0]) + numpy.cos(x[1]))/2
return M
x = UTPM(zeros((2,2,2),dtype=complex))
x.data[0,:] = [2+2j, 0+1j]
y = f(x)
-------------------- start code ----------------------------
In case you want the Jacobian of the function f you can use
----------------------- start code ------------------------
import numpy
import algopy
from algopy import UTPM, zeros
def f(x):
M = zeros((2,2), dtype=x)
M[0,0] = M[1,1] = 1
M[0,1] = numpy.sin(x[0]) + numpy.cos(x[1])
M[1,0] = (numpy.sin(x[0]) + numpy.cos(x[1]))/2
return M
x = UTPM.init_jacobian([2+2j, 0+1j])
y = f(x)
print UTPM.extract_jacobian(y)
----------------------- end code ------------------------
and get as output
python test_experimental.py[[[ 0.00000000+0.j 0.00000000+0.j ]
[-1.56562584-3.29789484j 0.00000000-1.17520119j]]
[[-0.78281292-1.64894742j 0.00000000-0.5876006j ]
[ 0.00000000+0.j 0.00000000+0.j ]]]
which is hopefully correct...
regards,
Sebastian
>
> Thanks in advance!!!!
>
> Maria.
your function maps from R^2 --> R^{2 \times 2} and therefore the Jacobian is a
R^{2 \times 2 \times 2} matrix.
You can also reformulate your function such that it is a mapping
R^2 --> R^{4}
----------- start code -----------
import numpy
import algopy
from algopy import UTPM, zeros
def f(x):
M = zeros((2,2), dtype=x)
M[0,0] = M[1,1] = 1
M[0,1] = numpy.sin(x[0]) + numpy.cos(x[1])
M[1,0] = (numpy.sin(x[0]) + numpy.cos(x[1]))/2
return M.reshape(4)
x = UTPM.init_jacobian([2+2j, 0+1j])
y = f(x)
print UTPM.extract_jacobian(y)
------------ end code -------------------
then you get the R^{4 \times 2} Jacobian matrix as output: