# Symbolic matrix differentiation with Sympy !

After a few days spent computing strange Jacobians for my Ph.D. thesis, I figured out that my computer could actually do most of the computations for me : all I needed was an automatic matrix differentiator, i.e. an algorithm that would tell me how a function or matrices $F(X,Y,...)$ varies when the matrices $(X,Y...)$ are changed into $(X + \partial X, Y + \partial Y,...)$ where the $(\partial X,\partial Y,...)$ are small. The very useful Matrix Cookbook gives for instance

$\partial (X^T) = (\partial X)^T$

$\partial (XY) = (\partial X)Y + X (\partial Y)$

$\partial(X^{-1}) = - X^{-1}(\partial X) X^{-1}$

… And so on. We can see that it is not completely unlike classic differentiation, but the non-commutativity slightly complicates things, and the computed formulae can be monstruous : imagine differentiating by hand the hat matrix $H = X \left(X^\top \Sigma^{-1} X\right)^{-1} X^\top \Sigma^{-1}$, where $\Sigma$ is a constant ðŸ˜¯

Maybe I didn’t look properly, but I didn’t find any program to automatize the rules described above, so I went on building something around Sympy, a young yet sophisticated symbolic mathematics module for the Python programmation language, whose authors claim that it is easily extendable. I put that claim to the test and indeed I found it pretty easy to work my way through.

In what follows I explain what I did, in tutorial-like fashion, with the hope that this can help understanding how one can tweak Sympy to meet one’s own needs. But I am no Sympy expert, and I take any tips and comments in the comments section below ! If you are in a hurry, you’ll find the complete source code of the solution at the bottom of this blog.

First, a matrix is a symbol, which is lucky because sympy already has a Symbol class. You only have to specify that these symbols are not commutative:


def matrices(names):
''' Call with  A,B,C = matrix('A B C') '''
return symbols(names,commutative=False)



#### Matrix-specific transformations

All the regular operations like addition, substraction, multiplication… are already implemented in the Symbol class, so we can focus on matrix-specific transformations, like the operator $\partial$ and the inverse function $X \mapsto X^{-1}$ :

d = Function("d",commutative=False)
inv = Function("inv",commutative=False)


That was tough ! Now for the transpose function $X \mapsto X^T$ we want the general rules

$(X+Y...)^T = X^T + Y^T ...$
$(XY...)^T = ...Y^TX^T$

to automatically applyÂ  when transposing an expression. Here is the code I wrote for that. To understand it you needÂ  to know that when I call some Sympy Function F over a bunch of arguments X,Y,Z, the result is an object of class “F” with a list of arguments (.args)Â  equal to [X,Y,Z]. Moreover successive additions or multiplications are flattened into one single object: for instance
$X^2 + 2*X*Y + 3*Z + 5$
will be an object of class “Add” with arguments
$[X^2,2*X*Y,3*Z,5]$,
while
$X^2*Y*X*Z$
will be of class “Mul” with arguments
$[X^2,Y,X,Z]$.

class t(Function):
''' The transposition, with special rules
t(A+B) = t(A) + t(B) and t(AB) = t(B)t(A) '''
is_commutative = False
def __new__(cls,arg):
return Add(*[t(A) for A in arg.args])
elif arg.is_Mul:
L = len(arg.args)
return Mul(*[t(arg.args[L-i-1]) for i in range(L)])
else:
return Function.__new__(cls,arg)


That’s enough transformations for now, let us explain to sympy how to differentiate a matrix expression !

#### Matrix differentiation

We write all the just enough differentiation rules of the Matrix Cookbook into one dictionnary:

MATRIX_DIFF_RULES = {
# e =expression, s = a list of symbols respsect to which
# we want to differentiate

Symbol : lambda e,s : d(e) if (e in s) else 0,
Mul :  lambda e,s : Mul(matDiff(e.args[0],s),Mul(*e.args[1:]))
+  Mul(e.args[0],matDiff(Mul(*e.args[1:]),s)) ,
t :   lambda e,s : t( matDiff(e.args[0],s) ),
inv : lambda e,s : - e * matDiff(e.args[0],s) * e
}


and apply them recursively to the expression we want to treat :

def matDiff(expr,symbols):
if expr.__class__ in MATRIX_DIFF_RULES.keys():
return  MATRIX_DIFF_RULES[expr.__class__](expr,symbols)
else:
return 0


Simple as pie, and we are done ! Let’s play around with our new toy:

X,S = matrices("X S")
H= X*inv(t(X)*inv(S)*X)*t(X)*inv(S)
print  mdiff(H,X)

>>> X*(inv(t(X)*inv(S)*X)*t(d(X))*inv(S) - inv(t(X)*inv(S)*X)*(t(X)*inv(S)*d(X) + t(d(X))*inv(S)*X)*inv(t(X)*inv(S)*X)*t(X)*inv(S)) + d(X)*inv(t(X)*inv(S)*X)*t(X)*inv(S)


It works ! But we must concede that it is barely readable… so let us put some style in these expressions !

## Cosmetics

You don’t need latex to have nice-looking formulae. Here is what we want to say to sympy:

• Don’t write inv(X), write XÂ¯Â¹
• Don’t write t(X), write X’
• Don’t write d(X), write âˆ‚X
• In a general way, don’t put parenthesis if the transformation applies to one symbol only.

The programmers of Sympy have made it all easy to customize the default printing method. You just write a printing method with complements for the functions of our own. First we need to include this line at the top of our source file to explain that we are going to use strange characters:

# -*- coding: utf-8 -*-


Then we write our own printer class that herits from the default one but has a few more methods :

class matStrPrinter(StrPrinter):
''' Nice printing for console mode : XÂ¯Â¹, X', âˆ‚X '''

def _print_inv(self, expr):
if expr.args[0].is_Symbol:
return  self._print(expr.args[0]) +'Â¯Â¹'
else:
return '(' +  self._print(expr.args[0]) + ')Â¯Â¹'

def _print_t(self, expr):
return  self._print(expr.args[0]) +"'"

def _print_d(self, expr):
if expr.args[0].is_Symbol:
return 'âˆ‚'+  self._print(expr.args[0])
else:
return 'âˆ‚('+  self._print(expr.args[0]) +')'


A this point the people who wrote the Sympy documentation want us to directly replace the default printer with this command:

Basic.__str__ = lambda self: matStrPrinter().doprint(self)


It indeed seems to be the only way to ensure that all the expressions in the recursion will print well with the standard “print” method, but it is a little to definitive for me (what if I have several customized printers to use in the same script ?). Moreover I don’t like the ‘*’ between the matrices and I didn’t find any smart way to put them away. For all these reasons I wrote my own matrix printing command:

def matPrint(m):
mem = Basic.__str__
Basic.__str__ = lambda self: matStrPrinter().doprint(self)
print str(m).replace('*','')
Basic.__str__ = mem


Now let’s try again:

matPrint(  matDiff(H,X) )

>>> X((X'SÂ¯Â¹X)Â¯Â¹âˆ‚X'SÂ¯Â¹ - (X'SÂ¯Â¹X)Â¯Â¹(X'SÂ¯Â¹âˆ‚X + âˆ‚X'SÂ¯Â¹X)(X'SÂ¯Â¹X)Â¯Â¹X'SÂ¯Â¹) + âˆ‚X(X'SÂ¯Â¹X)Â¯Â¹X'SÂ¯Â¹


Way better ! Notice, however that the first X is a factor or two terms of a sum, one of these being in turn a factor involving a sum. We can break these imbricated expressions into a sum of products (which is simpler to analyze) by calling twice Sympy’s ‘expand’ function :

matPrint(  expand ( expand ( matDiff(H,X) ) ) )

>>> X(X'SÂ¯Â¹X)Â¯Â¹âˆ‚X'SÂ¯Â¹ + âˆ‚X(X'SÂ¯Â¹X)Â¯Â¹X'SÂ¯Â¹ - X(X'SÂ¯Â¹X)Â¯Â¹X'SÂ¯Â¹âˆ‚X(X'SÂ¯Â¹X)Â¯Â¹X'SÂ¯Â¹ - X(X'SÂ¯Â¹X)Â¯Â¹âˆ‚X'SÂ¯Â¹X(X'SÂ¯Â¹X)Â¯Â¹X'SÂ¯Â¹


That’s nice enough for the console ! Now what if I want to report this fundamental result to the scientific community ? Let’s see how to generate the LaTeX code to embed the formula in a document.

## Automatic LaTeX code generation

One of the nice features of Sympy is its ability to generate LaTeX code (and even display compiled formulae). For matrix computations we would like Sympy to follow these rules:

• $X^{-1}$ writes X^{-1}
• $X^{T}$ writes X^{T}
• $\partial X$ writes \partial X
• In a general way, don’t put parenthesis if the transformation applies to one symbol only.

Like in the previous section we create a new ‘printer’ class with these features, and then write a function inpired by Sympy’s latex method (which is simply called latex() ):

class matLatPrinter(LatexPrinter):
''' Printing instructions for latex : X^{-1},  X^T, \partial X '''

def _print_inv(self, expr):
if expr.args[0].is_Symbol:
return self._print(expr.args[0]) +'^{-1}'
else:
return '(' + self._print(expr.args[0]) + ')^{-1}'
def _print_t(self, expr):
return  self._print(expr.args[0]) +'^T'

def _print_d(self, expr):
if expr.args[0].is_Symbol:
return '\partial '+ self._print(expr.args[0])
else:
return '\partial ('+ self._print(expr.args[0]) +')'

def matLatex(expr, profile=None, **kargs):
if profile is not None:
profile.update(kargs)
else:
profile = kargs
return matLatPrinter(profile).doprint(expr)



We now try it on the derivative of H:

print matLatex( matDiff(H,X) )

>>> $X \left((X^T S^{-1} X)^{-1} \partial X^T S^{-1} - (X^T S^{-1} X)^{-1} \left(X^T S^{-1} \partial X + \partial X^T S^{-1} X\right) (X^T S^{-1} X)^{-1} X^T S^{-1}\right) + \partial X (X^T S^{-1} X)^{-1} X^T S^{-1}$


which once compiled yields
$X \left((X^T S^{-1} X)^{-1} \partial X^T S^{-1} - (X^T S^{-1} X)^{-1} \left(X^T S^{-1} \partial X + \partial X^T S^{-1} X\right) (X^T S^{-1} X)^{-1} X^T S^{-1}\right) + \partial X (X^T S^{-1} X)^{-1} X^T S^{-1}$

Yeah !

## Source Code

I hope this helped. If you just came here to get a matrix differentiator, here is the code with the example. Have fun ðŸ™‚


# Declaration

# -*- coding: utf-8 -*-

#----------------------------------------------------------------------
#
# FUNCTIONS FOR THE AUTOMATIC DIFFERENTIATION  OF MATRICES WITH SYMPY
#
#----------------------------------------------------------------------

from sympy import *
from sympy.printing.str import StrPrinter
from sympy.printing.latex import LatexPrinter

#####  M  E  T  H  O  D  S

def matrices(names):
''' Call with  A,B,C = matrix('A B C') '''
return symbols(names,commutative=False)

# Transformations

d = Function("d",commutative=False)
inv = Function("inv",commutative=False)

class t(Function):
''' The transposition, with special rules
t(A+B) = t(A) + t(B) and t(AB) = t(B)t(A) '''
is_commutative = False
def __new__(cls,arg):
return Add(*[t(A) for A in arg.args])
elif arg.is_Mul:
L = len(arg.args)
return Mul(*[t(arg.args[L-i-1]) for i in range(L)])
else:
return Function.__new__(cls,arg)

# Differentiation

MATRIX_DIFF_RULES = {

# e =expression, s = a list of symbols respsect to which
# we want to differentiate

Symbol : lambda e,s : d(e) if (e in s) else 0,
Mul :  lambda e,s : Mul(matDiff(e.args[0],s),Mul(*e.args[1:]))
+  Mul(e.args[0],matDiff(Mul(*e.args[1:]),s)) ,
t :   lambda e,s : t( matDiff(e.args[0],s) ),
inv : lambda e,s : - e * matDiff(e.args[0],s) * e
}

def matDiff(expr,symbols):
if expr.__class__ in MATRIX_DIFF_RULES.keys():
return  MATRIX_DIFF_RULES[expr.__class__](expr,symbols)
else:
return 0

#####  C  O  S  M  E  T  I  C  S

# Console mode

class matStrPrinter(StrPrinter):
''' Nice printing for console mode : XÂ¯Â¹, X', âˆ‚X '''

def _print_inv(self, expr):
if expr.args[0].is_Symbol:
return  self._print(expr.args[0]) +'Â¯Â¹'
else:
return '(' +  self._print(expr.args[0]) + ')Â¯Â¹'

def _print_t(self, expr):
return  self._print(expr.args[0]) +"'"

def _print_d(self, expr):
if expr.args[0].is_Symbol:
return 'âˆ‚'+  self._print(expr.args[0])
else:
return 'âˆ‚('+  self._print(expr.args[0]) +')'

def matPrint(m):
mem = Basic.__str__
Basic.__str__ = lambda self: matStrPrinter().doprint(self)
print str(m).replace('*','')
Basic.__str__ = mem

# Latex mode

class matLatPrinter(LatexPrinter):
''' Printing instructions for latex : X^{-1},  X^T, \partial X '''

def _print_inv(self, expr):
if expr.args[0].is_Symbol:
return self._print(expr.args[0]) +'^{-1}'
else:
return '(' + self._print(expr.args[0]) + ')^{-1}'
def _print_t(self, expr):
return  self._print(expr.args[0]) +'^T'

def _print_d(self, expr):
if expr.args[0].is_Symbol:
return '\partial '+ self._print(expr.args[0])
else:
return '\partial ('+ self._print(expr.args[0]) +')'

def matLatex(expr, profile=None, **kargs):
if profile is not None:
profile.update(kargs)
else:
profile = kargs
return matLatPrinter(profile).doprint(expr)

#####    T  E  S  T  S

X,S = matrices("X S")
H= X*inv(t(X)*inv(S)*X)*t(X)*inv(S)

matPrint(  expand( expand( matDiff(H,X) ) ) )

print matLatex( matDiff(H,X) )


# Incommensuration with Platonions

A set of incommensurate (or incommensurable) numbers is a collection of reals such that it is not possible to obtain one of them by multiplying the others by rational coefficients and summing the results. A more formal definition states that the familyÂ $(I_i)_{i \leq N}$ is incommensurateÂ  if for any set of integer coefficients $(\lambda_i )_{i \leq N}$ such that

$\sum_{i \leq N}\lambda_i\,I_i\,=0$

it follows that the $\lambda_i$ are all null ($N$ can be infinite). Such numbers are quite useful : for instance if $(\omega_1,\omega_2,...\omega_N)$ are incommensurate, then the vector $V(t)=\big(cos(\omega_1 t),cos(\omega_2 t),...cos( \omega_N t) \big)$ can take any value in $[-1;1]^N$ as long as we set the variable $t$ appropriately. I recently found this trick used in Fourier amplitude sensitivity testing where it enables to collapse multiple integrals into single-dimensional ones.

Now that we have an idea why we could need them for, how do we find incommensurate numbers ? I have been surprised of how few one can find on this subject on the internet (links are welcome in the comments section !). It is clear that irrationals are good candidates, and among irrationals, square-roots of integers seem likely to provide incommensurate families. But not all of them : for instance $\sqrt{2}$ and $\sqrt{8}$ are not incommensurable ( since $2\sqrt{2}=\sqrt{8}$) ! In a very pleasant and erudite text, Kevin Brown focuses on the numbers of the form $\sqrt{N}$ where the integer $N$ is not a multiple of the square of an integer bigger than 1:

Â $P = \{\sqrt{n} \in \mathbb{N}, n\geq 1 \} - \{\sqrt{k*a^2} | (k,a) \in \mathbb{N}^2, k\geq 1,a \geq2 \}$

The first elements of this set are

$\sqrt{1},\sqrt{2},\sqrt{3},\sqrt{5},\sqrt{7},\sqrt{10},\sqrt{11},\sqrt{13},\sqrt{14},\sqrt{15}, \\ \sqrt{17},\sqrt{19},\sqrt{21},\sqrt{22},\sqrt{23},\sqrt{26},\sqrt{29},...$

Any finite subset of $P$, Brown shows, is incommensurate. The reals that you can produce with elements ofÂ $P$ are named PlatonionsÂ  (because of Platon and the quarternions… yep), but for the sake of simplicity I will call Platonions the elements of $P$ (which, in Brown’s vocable, would be called a Platonic field basis…).Â  These numbers provide a very cheap way of generating incommensurate sets : just take an array of the first N integers, discards the integers of the form $k*a^2$ (necessarily we will have $2\leq a\leq \sqrt{N}$ and $1\leq k \leq N/a^2$ ) and return the square roots of the remaining values. Here is the Python code:

def platonions_up_to(N):
""" Returns numbers of the form sqrt(I) where I is not a
multiple of a square. Such numbers are incommensurable."""
squareMultiples = [k*(a**2) for a in range(2,int(math.sqrt(N))+2)
for k in range(1,1+int((N+1)/(a**2)))]
return [math.sqrt(i) for i in range(1,N+1) if not i in squareMultiples]


Now what if we want exactly $K$ incommensurate numbers ? How big should $N$ be to be sure to find at least $K$ Platonions inferior to $N$ ? The good news is that Platonions are plentyfull. More precisely, the fraction of integers smaller than $N$ that are Platonions is pretty much constant:

import matplotlib.pyplot as plt
Nlist = range(50,5051,500) # N = 50,550,1050,1550... 5050
proportionslist = [1.0*len(platonions_up_to(N))/N for N in Nlist]
plt.plot(Nlist,proportionslist);plt.show()


Fraction of integers smaller than N (x-axis) whose square-root is a Platonion

We see that the fraction integers smaller than N whose square-root is a Platonion is always bigger than 0.60. To be sure to find at least $K$ Platonions we thus need to consider the $K/0.60$ first integers only. In Python:

def findKincommensurates(K):
return platonions_up_to(int(K/0.60))[:K]


That’s as far as I will go today ! You can now reread the title and assess the progress you made in a topic that you had no idea existed a few minutes ago.

# Gipsy-like piano music (with sheet music)

A few days ago I recorded this :

I am aware that the piano has a most un-gipsic sound, but I never had the courage to learn guitar, so I am trying to bring the wonderful repertory of gipsy jazz to my favorite instrument. As I have been asked for the sheet music, I wrote a transcription :

Une Routine Manouche (sheet music)

For the fun I placed this piece of music under a Creative Common licence, which basically lets you free to do whatever you want with it, even make money ( but don’t be too confident about that ! ), without needing my permission.Â  All I ask is that you cite a certain “A. de la Marmotte” as the author, so that I can google it in ten years and see what has become of my little routine.

In case you are wondering, I edited the music with the free software Lilypond, which is somehow an equivalent of LaTeX for sheet music.Â  If you want to modify the score, here is the source file I wrote (copy it in a file like “myFile.ly” and run with Lilypond):

 % OPEN WITH LILYPOND

\version "2.12.3"

%{
*
* U N EÂ Â  R O U T I N EÂ Â  M A N O U C H E
*
* Version 0 - 31/03/2012
*
*
* The full text of the licence can be found here:
*
* In a nutshell, YOU ARE FREE to share (copy, distribute, transmit),
* to remix or adapt the work, and to make commercial use of the work.
*
* BUT I ask you to cite "A. de la Marmotte" as the original author
* if you do any of the things above, just so that I can see what this
* little piece is becoming.
*
* Ideas for improvements are most welcome !
*
%}

title = "Une routine manouche"
subtitle = "by A. de Lamarmotte"
subsubtitle = "After the chorus from V. Monti's famous Czardas"
arranger="Grenoble, France, March 2012"
composer = "This is a free work ( licence CC BY 3.0)"
}

% In what follows, one line represents one bar.

global = {
\key d \minor
\time 4/4
}

right = \new Voice \with {
}{
\global
d''4 \grace{b''16[ c''' cis''']} d'''8-5 cis''' c''' b''4-2 bes''-4
a'' gis''8 a'' bes'' a'' f'' d''4.-.-1 a'4.-3 r4
r2. <g-1 cis'-4>8 <a-2 d'-5>4
f'8-3 e' g' f'-1 a' g' bes'
a'-1 d''-3 cis'' e'' d''-1 f'' e'' g''
f''4 g'' f''8 g''16 f''16 e''8 bes'4. r2 r8 bes''-3
\times 2/3 {a''-2 bes'' a''} g'' bes''Â  \times 2/3 {a'' bes'' a''} g'' bes''
\times 2/3 {a'' bes'' a''} gis'' a'' bes'' a''-1 c''' bes''
\times 2/3 {a'' bes''-3 a''} gis''-2 a'' f'''-. a'' e'''-. a''
es'''-. a'' d'''-. a'' cis'''-. a'' c'''-. a''
b''-4 bes''-3 f'' d'' b'-4 bes' f' d'
e'-2 b' gis' e' d' gis' e' d'-1
cis'-2 e'-1 g' bes' e'-1 g' bes' cis''
g'-1 bes' cis'' e'' bes'-1 cis'' e'' g''
d''4-1 d''' \acciaccatura{gis''16} a''2-3
\times 2/3 {g''8 a'' g''} f''-1 g''-2 \acciaccatura{gis''16-3} a''8-5 f''8-4 d''4-2
<<
{ r4 <c'' d'' f''> <b' d'' f''>2Â  <c'' d'' f''>8 <c'' d'' f''>4
<b' d'' f''>4.} \\
{ a'1~a'2. }
>>Â  r8 \times 2/3 {a''16-5 f'' e''}
d''8-1 cis''-3 c'' b' bes'-2 a' d''-4 cis''
c'' b' bes'-2 a' d'' cis'' c'' b'-1
<g'-2 bes'-4>4.-. <g' bes'>4 g'8 bes'4
<g' a'>2. r8 bes''-4

\times 2/3 {g''8-2 bes'' g''} e'' bes'' \times 2/3 {g'' bes'' g''} e'' bes''
\times 2/3 {g'' bes'' g''} e'' g''-4 \times 2/3 {e''-2 g'' e''} cis'' bes''-4
a''-^-3 f'' d'' g''-^-4 e''-2 cis''-1 f''-^-4 d''-2
b'-1 e''-^ cis'' bes' d''-^ a' f' cis''-3
c'' b'-1 bes'-4 a' gis' g' c''-4 b'
bes' a'-1 gis'-2 g' c''-5 b' bes' a'
<e' a' d''>8-^ <e' a' d''>4.-^ <e' a' d''>4-^Â Â  <e' a' d''>8-^ <e' a' d''>4.-^

r4 e''-^-3 \times 2/3 {dis''8-^ e'' dis''}
d''4.-^-1 c''8-2 d''-1 e''-2 g''-3 bes''-4
d'''4.-5 \times 2/3 {bes''16 g'' e''} d''4 r4
\ottava #1
r4 \times 2/3 {a'''8-3 c''''-5 a'''} f'''-1 a'''-5 d'''-2 f'''-4
\ottava #0
c'''-1 d'''-5 a'' c''' f''-1 a''-5 d''-1 f''-2
gis''-4 a'' g'' f'' e'' d'' cis''-4 bes'
a' g' f' e' d' cis'-4 bes-2 a-1
e'-4( cis'-2 d'-1) g'-4( e'-2 f'-1) bes' gis'
a' e'' cis'' d'' g'' e'' f'' g''-4
f''4-3 g'' f''8 g''16 f'' e''4
<bes' bes''>4-. r8 <bes' e'' bes''>4. r4
a'8 d'' e'' f'' e''4 d''
<b' a''>8 f''4-. <b' d''>4 f''8 d'' b'-1
bes'4-2 f'' d''8 bes' a'4-1
<< {r4 <d''-4 f''-5> r8 <a' cis'' e''>4 <a' d''>8~<a' d''>4} \\ {gis'2-2 g'4.-1 f'8-1~f'4}>>
r2.
r2 e'4-^ \times 2/3 {dis'8-^ e' dis'}
d'4. d'''8 \times 2/3 { bes'' d''' bes'' } g'' bes''
d''2 e'4-^ \times 2/3 {dis'8-^ e' dis'}
d'4 c'8 d' f' a' d'' f''
d'''2 a'4-^ \times 2/3 {gis'8-^ a' gis'}
r8 <bes'' f'''>8 <cis''' g''> r8 <bes' f''> <cis'' g'> r8 gis'
<< {r4 <cis'' f''> r8 <a' cis'' e''>4 r8 <a' d'' f''>4.-^ <a' d'' f''>4.-^ <a' d'' f''>4-^}
\\ {a'2-^ g'4.-^ f'8-^~f'1}>>
r4 d''4 e'' fis''
\times 2/3 {g''8 a'' g''} fis''8 a''8 \times 2/3 {g''8 a'' g''} fis''8 a''8
g'' a'' ais'' b'' c''' d''' b'' d'''
\times 2/3 {c'''8 d''' c'''} b'' d''' \times 2/3 {c'''8 d''' c'''} b'' d'''
\times 4/5 {c'''8-3 d''' c''' a'' f''} \times 4/5 {c''8-3 d'' c'' a' f'}
d'-2 f'-1 g'-3 aes'-4 g' f' d'-1 c'-2
d'-1 f' g' aes' g' f' d'-1 c'-2
<< {r4 <f' a' d''>4.<f' a' d''>4.}\\{d'2.-1 d'4}>>
<f' a' d''>4. <f' a' d''>4. <f' a' d''>4
r8 <d''' f''>8 <gis'' d''> r8 <d'' f'>8 <gis' d'> r8 <d'' f'>
r8 <d''' e''>8 <g'' d''> r8 <d'' e'>8 <g' d'> r4
\times 8/12 {d'16_\markup{ \small \bold Freely } [( f' a' b'Â Â  d'' f'' a'' b''
\ottava #1
d''' f''' a''' b'''])} d''''2
\ottava #0
r1 \bar "|."

}

left = \new Voice \with {
}{
\global
d4( <f a b>-.) a,( <f a b>-.)
e <g bes d'> cis <e g a>
d-^ <f a b> e-^ <g bes d'>
f-^ <a b d'> a, <cis e>
d <f a b> a, <f a b>
d <f a b> c <f a b>
bes, <g bes d'> g,Â  <g bes d'>
bes, <g bes d'> g,Â  <g bes d'>
e <g bes d'> bes,Â  <g bes d'>
<< {r4 <g bes cis'>2Â  <g bes cis'>4} \\ {ees1-^} >>
d4-^ <f a b> a, <f a b>
d <f a b> a, <f a b>
<< {r4 <f bes>2 <f bes>4} \\ {d1} >>
e4-^ <gis b> d-^ <gis b>
cis-^ <e g bes> e-^ <g bes cis'>
g-^ <bes cis' e'> bes-^ <cis' e' g'>
\clef treble
d'-^Â  <f' a' > cis'-^Â  <f' a' >
c'-^Â  <f' a' > a-^Â  <cis' e' g'>
\clef bass
d <f a b> a, <f a b>
d <f a b> a, <f a b>
d-^ <f a b> e-^ < bes d'>
f-^ <bes d'> fis <bes d'>
g <bes d' e'> e <bes d' e'>
g <a d' e'> e <a d' e'>
g <bes cis' e'> e <bes cis' e'>
cis <bes cis' e'>Â  e <bes cis' e'>
d <f a b> a, <f a b>
d <f a b> a, <f a b>
e <g bes c'> cis <g bes c'>
e <g bes c'> cis <g bes c'>
d <f a b> d <f a b>
d <f a b> <d, d>-^ <dis, dis>-^
<e, e> <g bes c'> c <g bes c'>
d <g bes c'> e <g bes c'>
f < a c' d'> c < a c' d'>
f < a c' d'> c < a c' d'>
cis <e g a> a, <e g a>
b, <e g a> cis, <e g a>
d-^ <f a> cis-^ <f a>
c-^ <f a> b,-^ <f a>
bes,-^ <g bes d' e'> g, <g bes d' e'>
bes,Â Â  <g bes d' e'> g, <g bes d' e'>
d <f a b> a, <f a b>
d <f a b> a, <f a b>
bes, <f bes d'> f, <f bes d'>
d <e gis b> cis <e g a>
d <f a b> a, <f a b>
d <f a b> a, <f a b>
e <g bes c'> c <g bes c'>
e <g bes c'> c <g bes c'>
f <a> c <a c'>
f <a c' d'> c <a c' d'>
\clef treble
<g bes cis' g'>4.-^\arpeggio <bes' cis''>4.-. <bes cis'>4-.
\clef bass
a4 <cis' f'> g <a cis' e'>
d <f a b> a, <f a b>
d <f a b> a, <f a b>
e <g bes c'> c <g bes c'>
e <g bes c'> c <g bes c'>
f <a c'> c <a c'>
f <a c' d'> c <a c' d'>
<bes, gis>4 <gis bes d'> <bes, gis> <gis bes d'>
<a, g>4 <g a c'> <a, g> <g a c'>
d <f a b> a, <f a b>
d <f a b> a, <f a b>
<bes, gis>4. <bes gis'>4. <bes, gis>4
<a, g>4. <a g'>4. <a, g>4
<d f a b>2.\arpeggio

\ottava #-1
\acciaccatura {cis,8}Â  d,4 d,,1
\bar "|."

}
#(set-global-staff-size 18)

\score {
\new PianoStaff \with {
instrumentName = "Piano"
} <<
\new Staff = "right" \with {
midiInstrument = "acoustic grand"
} \right
\new Staff = "left" \with {
midiInstrument = "acoustic grand"
} { \clef bass \left }
>>
\layout { }
\midi {
\context {
\Score
tempoWholesPerMinute = #(ly:make-moment 100 4)
}
}
}