Re: algopy question

[Sebastian Walter]

Hello Andrey, 

I'm not sure I fully understand your question.

Could you provide a full example that can be run in Python,

or a mathematical description?

best regards,

Sebastian

On Wed, May 15, 2013 at 10:52 AM, Yachmenev, Andrey <a.yac...@ucl.ac.uk> wrote:

Dear Sebastian,

I am very new to python, I switched to it from fortran since I've been told that it is much easier to treat problems related with automatic differentiation. Among few python modules for automatic differentiation I found AlgoPy the easiest to employ. With simple test problems I did it works fine, but I wonder if you have a wrapper to simplify the computation of mixed partial derivatives.

Like, if I provide the following test matrix function, then wrapper function returns a list of requested partial derivatives (with respect to r[:]) of each matrix element.

def ddt(x, x0, terms, coefs):
   """ Computes elements of D.D^T where D is expanded in Taylor series terms of x-x0 """
    nx = len(x)
    nterms = len(terms)
    dimen = len(coefs[0][:,0])
    mat = numpy.zeros((dimen,dimen), 'float64')
    for iterm in range(nterms):
        prod = 1.0
        for ix in range(nx):
            prod = prod * (x[ix] - x0[ix])**terms[iterm][ix]
        mat[:,:] = mat[:,:] + coefs[iterm][:,:] * prod
    return numpy.dot(mat, numpy.transpose(mat))

Thank you in advance.

With Best Regards,
Andrey


---

Andrey Yachmenev, Dr.
Department of Physics & Astronomy
University College London
Gower Street, London WC1E 6BT, UK
Tel (work): +(44) 020 7679 0172
Email: a.yac...@ucl.ac.uk

[Sebastian Walter]

Hi Andrey,

yes, that's certainly possible.

Since Algopy only computes in univariate Taylor arithmetic, one has to convert several univariate Taylor polynomials to one multivariate Taylor polynomial.

Basically, the idea is

1) INPUT: multivariate polynomial A(r1,..,rn)

2) CONVERT: to set of univariate polynomials A1(r), A2(r), ... 

3) COMPUTE: in univariate Taylor arithmetic

4) CONVERT: set of univariate polynomials f(A1(r)), f(A2(r)), ...  back to multivariate polynomial

The functionality is implemented in the submodule `algopy.exact_interpolation`.

The mathematical details are described on page 315 of the book "Evaluating Derivatives" by Andreas Griewank,

Chapter 13, Subsection: Multivariate Tensors via Univariate Tensors.

A more detailed and easier to understand description can be found in the original paper "Evaluating higher derivative tensors by forward propagation of univariate Taylor series"

by  Andreas Griewank, Jean Utke, Andrea Walther.

Example (taken from `test_exact_interpolation.py`):

It assumes the special case when the input polynomial is

A(r1,...,rn) = A_{0} + A_{1,0,..} r1 + A_{0,1,0...} r2 + ... + A_{0,...0,1} rn .

    def test_interpolation(self):

        def f(x):

            return x[0] + x[1] + 3.*x[0]*x[1] + 7.*x[1]*x[1] + 17.*x[0]*x[0]*x[0]

        

        N = 2

        D = 5

        deg_list = [0,1,2,3,4]

        coeff_list = []

        for n,deg in enumerate(deg_list):

            Gamma, rays = generate_Gamma_and_rays(N,deg)   # 1) INPUT

            x = UTPM(numpy.zeros((D,) + rays.shape)) # 2) CONVERT

            x.data[1,:,:] = rays

            y = f(x) # 3) COMPUTE

            coeff_list.append(numpy.dot(Gamma, y.data[deg])) # 4) CONVERT

            

        assert_array_almost_equal([0], coeff_list[0])

        assert_array_almost_equal([1,1], coeff_list[1])

        assert_array_almost_equal([0,3,7], coeff_list[2])

        assert_array_almost_equal([17,0,0,0], coeff_list[3])    

If your input polynomial is not of the special form as above you'll have to adapt the conversion step.

regards,

Sebastian

On Sat, May 18, 2013 at 1:19 PM, Yachmenev, Andrey <a.yac...@ucl.ac.uk> wrote:

Yes, exactly, but involving cross-terms as well like A_{2,2,...}*r1^2*r2^2. These cross-term partial derivatives are the main problem for me in terms of algopy.


From: Sebastian Walter [sebastia...@gmail.com]

Sent: 18 May 2013 11:16
To: Yachmenev, Andrey
Subject: Re: algopy question

so your A is of the form

A = A_0 + A_{1,0,...,0} r1 + A_{0,1,0,...,0} r2 + ... + A_{0,...,1} rn + A_{2,0,...,0) r1^2 + ... + A_{0,...,2} rn^2 + terms in r1^3 ... rn^3 etc. ?

On Fri, May 17, 2013 at 4:45 PM, Yachmenev, Andrey <a.yac...@ucl.ac.uk> wrote:

Dear Sebastian,

thank you for helping me.
I would like to compute high-order partial derivatives over r1, r2, r3, ... rn variables of symmetric matrix M
M(r1,r2,r3,...rn) = A(r1,r2,r3,...rn) * A(r1,r2,r3,...rn)^T,
where A(r1,r2,r3,...rn) is expanded in Taylor series in r1, r2, r3, ... rn (around some point r1_0, r2_0, r3_0, ... rn_0).
Basically I know Taylor expansion of symmetric matrix A and would like to convert it into expansion of matrix M=A*A^T.

With Best Regards,
Andrey

---

Andrey Yachmenev, Dr.
Department of Physics & Astronomy
University College London
Gower Street, London WC1E 6BT, UK
Tel (work): +(44) 020 7679 0172
Email: a.yac...@ucl.ac.uk

---

From: Sebastian Walter [sebastia...@gmail.com]
Sent: 17 May 2013 15:22
To: Yachmenev, Andrey; alg...@googlegroups.com
Subject: Re: algopy question