Ideally, we would implement these rules symbolically, so that they work
even for something like MatrixSymbol('x', n, n). See
https://code.google.com/p/sympy/issues/detail?id=2759. Getting this working
would be a really nice feature for SymPy.Aaron Meurer On Thu, Oct 31, 2013 at 4:09 AM, Saullo Castro <[email protected]>wrote: > I have the functional: > > > <https://lh5.googleusercontent.com/-YyVcmvfWlB8/UnIq_SuBdjI/AAAAAAAABYE/NkiWxG35LS8/s1600/fig1.png> > > Where A is a function of "v". The non-linear system of equations necessary > to find "v" is obtained doing: > > > <https://lh5.googleusercontent.com/-SElccyNi1xo/UnIrNKzuKxI/AAAAAAAABYM/xd6tHk0SJT0/s1600/fig2.png> > > According to the differentiation rules, this system could be written in > matrix form as: > > > <https://lh6.googleusercontent.com/-Qm0cHaMnEAA/UnIrYLIfs5I/AAAAAAAABYU/ocDPAfGA4ic/s1600/fig3.png> > > So I've implemented in Sympy, with the help of NumPy "ndarray" the > differentiation about a vector. In this case the shape of (\partial A / > \partial v is (n X n X n), where "n" is the length of vector "v". > > The function is simple: > def vdiff(x, vector): > x = np.array(x) > shape = x.shape > ans = [np.array([e.diff(vi) for e in x.ravel()]) for vi in vector] > ans = [a.reshape(shape) for a in ans] > return np.array(ans).swapaxes(0, 1) > > > > And I've tested with two different cases: > > def test_00(): > from sympy.abc import a, b, c, x > A11 = x*a*c**2 > A12 = x**2*a*b*c > A13 = x**2*a**3*b**5 > A21 = x**3*a**2*b*c > A22 = x**4*a*b**2*c**5 > A23 = 5*x**4*a*b**2*c > A31 = x**4*a*b**2*c**4 > A32 = 4*x**4*a*b**2*c**2 > A33 = 4*x**4*a**5*b**2*c > A = np.array([[A11, A12, A13], > [A21, A22, A23], > [A31, A32, A33]]) > v = np.array([a, b, c]) > F = (v.T.dot(A)).dot(v) > Av = vdiff(A, v) > p1 = v.dot(A) > p2 = A.dot(v) > p3 = v.dot(Av.dot(v)) > new = p1 + p2 + p3 > ref = np.array([F.diff(a), F.diff(b), F.diff(c)]) > print sympy.simplify(Matrix(ref-new)) > > def test_01(): > from sympy.abc import a, b, c, x > sympy.var('c1, c2') > A11 = x*a*c**2 > A12 = x**2*a*b*c > A13 = x**2*a**3*b**5 > A21 = x**3*a**2*b*c > A22 = x**4*a*b**2*c**5 > A23 = 5*x**4*a*b**2*c > A = np.array([[A11, A12, A13], > [A21, A22, A23]]) > v = np.array([a, b, c]) > cc = np.array([c1, c2]) > F = cc.dot(A.dot(v)) > ref = np.array([F.diff(a), F.diff(b), F.diff(c)]) > Av = vdiff(A, v) > p1 = cc.dot(A) > p2 = cc.dot(Av.dot(v)) > new = p1 + p2 > print sympy.simplify(Matrix(ref-new)) > > > I'd like to share here and get some feedback in case someone has a more > straightforward way to implement this... > > > Thank you! > > Saullo > > > > > -- > You received this message because you are subscribed to the Google Groups > "sympy" group. > To unsubscribe from this group and stop receiving emails from it, send an > email to [email protected]. > To post to this group, send email to [email protected]. > Visit this group at http://groups.google.com/group/sympy. > For more options, visit https://groups.google.com/groups/opt_out. > -- You received this message because you are subscribed to the Google Groups "sympy" group. To unsubscribe from this group and stop receiving emails from it, send an email to [email protected]. To post to this group, send email to [email protected]. Visit this group at http://groups.google.com/group/sympy. For more options, visit https://groups.google.com/groups/opt_out.
