TITLE: How to find finite trigonometric products
QUESTION [10 upvotes]: I wonder how to prove ?
$$\prod_{k=1}^{n}\left(1+2\cos\frac{2\pi 3^k}{3^n+1} \right)=1$$
give me a tip

REPLY [3 votes]: Let $S_n = \sum_{k=0}^n 3^k = \frac{3^{n+1}-1}{2}$.  Then
$$3^{n}- S_{n-1} = 3^{n} - \frac{3^{n}-1}{2} = \frac{3^{n}+1}{2} = S_{n-1}+1.
$$
Now by induction we have the following product identity for $n \geq 0$:
$$
\begin{eqnarray}
\prod_{k=0}^{n}\left(z^{3^k}+1+z^{-3^k}\right) &=& \left(z^{3^{n}}+1+z^{-3^{n}}\right)\prod_{k=0}^{n-1}\left(z^{3^k}+1+z^{-3^k}\right) \\
&=& \left(z^{3^{n}}+1+z^{-3^{n}}\right) \left(\sum_{k=-S_{n-1}}^{S_{n-1}} z^k\right) \\
&=&\sum_{k=S_{n-1}+1}^{S_n}z^k + \sum_{k=-S_{n-1}}^{S_{n-1}}z^k+\sum_{k=-S_n}^{-S_{n-1}-1} z^k \\
&=& \sum_{k=-S_n}^{S_n} z^k
\end{eqnarray}
$$
Now take $z = \exp\left(\frac{\pi \, i}{3^n + 1}\right)$ and use that $z^{3^n+1}=-1$ to get
$$\begin{eqnarray}
\prod_{k=0}^n\left(1 + 2 \cos \left(\frac{2 \pi \,3^k}{3^n+1}\right)\right) &=& \sum_{k=-S_n}^{S_n}z^{2k} = \frac{z^{2S_n+1}-z^{-2S_n-1}}{z-z^{-1}} = \frac{z^{3^{n+1}}-z^{-3^{n+1}}}{z-z^{-1}} \\
&=& \frac{z^{3(3^n+1)-3} - z^{-3(3^n+1)+3}}{z-z^{-1}} = \frac{z^3-z^{-3}}{z-z^{-1}} = z^2 + 1 + z^{-2} \\
&=& 1 + 2\cos\left(\frac{2\pi}{3^n+1}\right)
\end{eqnarray}$$
and your identity follows.


TITLE: Why do mathematicians understand structure preservation in a such a way?
QUESTION [1 upvotes]: This question is about homomorphism and it's definition, which still feels kind of weird and unreliable to me.
To be more obvious, let's take graph $G = (Vertices(G), Edges(G))$. Then homomorphism of $G$ to somewhat graph $G'$ consists of two independent functions: $vrt : Vertices(G) \mapsto Vertices(G')$, $edg : Edges(G) \mapsto Edges(G')$. As such, in this particular case, homomorphism itself becomes a complicated object: it can be decomposed into some simpler primitives and thus one has to define certain rule how to apply those primitives to $G$ in order to get $G'$. Formally it looks somewhat like:
$$(\forall v_G \in Vertices(G) : vrt(v_G) \in Vertices(G')) \land (\forall e_G \in Edges(G) : edg(e_G) \in Edges(G')) \tag{1}$$
Observation #1: I can't see how it follows from the homomorphism's definition. There is no $(\circ)$ given explicitly; what kind of operation is it? Which are the elements of that operation? Can someone show me why $(vrt, edg)$ definitions follow from that?
Observation #2: I see homomorphism cares not to degrade existing structure, at the same time saying literally nothing about upgrading it? Since homomorphism is complicated object on it's own right, it can preserve existing structure and add something, which does not exist in the original one. And (1) requirement from the above would hold indeed. In order to ensure no upgrading took place, one should extend (1) to:
$$\text{vertices of the G did not degrade}
\land
\text{edges of the G did not degrade}
\land
\text{all the vertices of the G' came from the G}
\land
\text{all the edges of the G' came from the the G} \tag{2}$$
where:

vertices of the G did not degrade = $(\forall v_G \in Vertices(G) :
   vrt(v_G) \in Vertices(G'))$;
edges of the G did not degrade = $(\forall e_G \in Edges(G) : edg(e_G) \in Edges(G'))$
all the vertices of the G' came from the G = $\not \exists v_{G'} \in vrt(G') : \forall v_G \in vrt(G) : vrt(v_G) \not = v_{G'}$
all the edges of the G' came from the G = $\not \exists e_{G'} \in edg(G') : \forall e_G \in edg(G) : edg(e_G) \not = e_{G'}$

As a conclusion. I don't have any teacher to learn math, I do it only with help of the books published online and this platform. And for me it feels like homomorphism is both: 1) extraordinary important, 2) poorly defined idea. At the first sight it seems rigorous, but the deeper one dig, the sloppier it becomes. It feels like mathematicians just don't bother to redefine loosed things: whenever it comes down to somewhat specific structures, like graphs or monoid, or whatever else, homomorphism emerges not as direct consequence of the $F(g \circ f) = F(g) \bullet F(f)$, but as something intuitive, something like "we all understand that to preserve such-and-such particular structure, one has to fulfill following requirements". And yet I haven't seen "no-upgrade" requirement (observation #2), which seems quite valuable.
So, please, help me get deeper into math. What's wrong with my understanding?

REPLY [3 votes]: In general, perhaps the best advice is to be skeptical but open minded.  Observe what the standard notion of morphism does for you in each mathematical context, and observe the analogies between them.  In fact, the analogies are so strong that eventually people invented whole subjects to deal with them, category theory, universal algebra.  But these general abstract theories are definitely not the place to start learning mathematics.
I think if you just take the examples of theories that seem to concern you at the moment, you can identify two rather different sorts: There are algebraic structures like groups, rings, modules, etc.  and there are relational structures such as graphs and partially ordered sets.  Perhaps it would help you to make this distinction.
For the algebraic examples, a morphism is always a map, not necessarily injective or surjective, that preserves the algebraic operation or operations.  For example, for rings there are two operations $\cdot$ and $+$, and  the requirement is $\phi(x\cdot y) = \phi(x)\cdot \phi(y)$ and  $\phi(x + y) = \phi(x) + \phi(y)$.  (Frequently, one works with rings with multiplicative identity, and requires in addition that $\phi$ take the identity of one ring to the identity of the other.  You already see from this that there are choices to be made; there are  different theories depending on these choices.)
For graphs, there are also many slightly different flavors: graphs with only single edges allowed, graphs with multiple edges allowed, directed graphs, and more.  Each of these would require a slightly different notion of structure preserving map.  The simplest flavor would be graphs with only single edges between vertices (and possibly loops).  This is the flavor treated in the wikipedia article https://en.wikipedia.org/wiki/Graph_homomorphism.  In this case the graph consists of the vertex set and a relation $e(x,y)$ which is symmetric but neither transitive nor reflexive.  Namely $e(x, y)$ holds if there is an edge connecting vertices $x$ and $y$.  Now a graph homomorphism is a map on vertices which preserves the relation:  $e(x, y)$ implies $e(\phi(x), \phi(y))$.
As soon as one allows more complicated setups, say graphs with multiple edges, then one has to invent a more complicated way both to describe what a graph is and what a morphism is.  But the idea should still be that edges get mapped to edges.  You could just start with this intuitive idea yourself and work out how to formalize this.
Now if you really want to make an effort, you could bring the algebraic examples and the relational examples under one roof, e.g by viewing an operation as a relation,  $R(x, y, z)$ holds if $x\cdot y = z$.  Now the homomorphism property would translate to $R(x, y, z)$ implies $R(\phi(x), \phi(y), \phi(z))$.  For me, this does absolutely nothing to help me understand groups or rings.


TITLE: About eigenvalues, eigenvectors, and their relation with the root square of a matrix.
QUESTION [0 upvotes]: I was trying to obtain the square root of a matrix through the eigenvalues and eigenvectors, and there is something that doesn't add up in some of the demonstrations that I observed after getting stuck.
So, being $Q$ the eigenvector column matrix of $A$ and  $\Lambda$ a diagonal matrix with the eigenvalues of $A$, we can put $A$ as:
$$ A = Q \Lambda Q^{-1} $$
From this is easily deducible that:
$$ A^n = (Q \Lambda Q^{-1})(Q \Lambda Q^{-1})...(Q \Lambda Q^{-1}) = Q \Lambda^n Q^{-1} $$
But the demonstration of this only seems valid for all the natural numbers and I cannot see how it can extend directly for the integer or rational numbers.

REPLY [0 votes]: Not every matrix has a square root. For example, consider
$$
            A = \left[\begin{array}{cc}0 & 1 \\ 0 & 0\end{array}\right]
$$
This has no $2\times 2$ square root $B$. If it did, then $B^4=A^2=0$ would hold, which would mean that the minimal polynomial for $B$ would either be $\lambda$, or $\lambda^2$. However, $B^2=A\ne 0$. So $B=0$ would have to hold, which is a contradiction.


TITLE: How to approach on this - finding minimum distance of point on the ellipse from the centre of it.
QUESTION [0 upvotes]: Question
The minimum distance of any point on the ellipse $$x^2+3y^2+4xy=4$$ from its centre is ______.
Attempt
Converted the given expression into $$(x+2y)^2-y^2=4$$. But, this becomes equation of hyperbola of the form $$\frac{x^2}{a^2}-\frac{y^2}{b^2}=1$$.
Then how it is ellipse?
Any hints or suggestion?

REPLY [1 votes]: The matrix
$$
\begin{pmatrix} 1 & 2 \\ 2 & 3 \end{pmatrix}
$$
defining the hyperbola $(x,y)B(x,y)^T=4$ has eigenvalues $2\pm\sqrt 5$. The eigenvector associated to the positive eigenvalue is $(\sqrt5-1,2)$. Therefore the vertices are obtained by solving
$$
\left\{\begin{array}{l}
x^2+4xy+3y^2=4 \\
(x,y) = t(\sqrt5-1,2) \\
\end{array}\right.
$$
The solutions are
$$
x_0 = \sqrt{\frac{14}{\sqrt{5
   }}-6}
\qquad
y_0 = \sqrt{\frac{6}{\sqrt{
   5}}-2}
$$
and the opposite point $-(x_0,y_0)$.
Therefore the distance from the center (which is the origin) is
$$
|(x_0,y_0)| = \sqrt{x_0^2+y_0^2} = 2 \sqrt{\sqrt5-2} .
$$


TITLE: Possible Word Combinations
QUESTION [0 upvotes]: How many possible word combinations can exist for a alphabet of 4 vowels and 8 consonants, 12 letters altogether where a word can start with a vowel or consonant, but every word must end in a vowel?  Moreover, every vowel must be followed by a consonant (except at the end). Each consonant can be followed by another consonant, or a vowel. Double consonants must be followed by a vowel.  Max letters per word are six.  How many possible words can exist?

REPLY [1 votes]: First write down strings like VCCVCV that represent allowable words where V is a vowel and C is a consonant.  Without the restrictions there would be $2^6=64$ of them.  Ending in a vowel cuts that in half, so there are less than $32$ of them.  For each one use the multiplication principle.  Each V gives $4$ possibilities, each C gives $8$, so VCCVCV gives $4^38^3=2^{15}=32768$ words.  You just need to find how many patters there are with $3$ vowels and $3$ consonants, then mulitply by $32768$ and add to the results with other numbers of vowels/consonants.

REPLY [0 votes]: Find, by chasing cases, the allowed words of length $6$ over the alphabet $\{V,C\}$. Sort them into types according to the number of vowels occurring. Then for each type determine in how many ways you can distribute the available $4$ vowels and $8$ consonants on a $\{V,C\}$-word of this type.
A word can begin with $VC$, $CVC$, or $CCVC$, and has to end with $CV$. This leaves the four types
$$VCCVCV,\quad VCVCCV,\quad CVCVCV,\quad CCVCCV\ ,$$
three of them of type $3V\ 3C$, and one of them of type $2V\ 4C$. It follows that there are
$$3\cdot4^3\cdot8^3+ 4^2\cdot8^4=163\,840$$
admissible words.


TITLE: What is your favorite group?
QUESTION [6 upvotes]: I would like to know about your favorite group(s). Since groups do appear everywhere in mathematics and there are plenty of them, which ones have drawn your attention the most or surprised you? Please not just name the group, but also provide some facts about it why you find this one particularly interesting.
I'll start by mentioning Grigorchuks group.
Because is was the first group i encountered, wich is finitely generated but not finitely presented. Also it was the first group discovered with intermediate growth. It has a lot of "strange" properties like:

it's infinite but residually finite
it's amenable but not elementary amenable
every proper quotient group is finite
every maximal subgroup has finite index

Also Grigorchuks group acts as a key-counterexample in infinite group theory. My professor once told me: "If you have a conjecture about infinite groups, try it one Grigorchuks group. If it holds, it might be worth trying to prove it."
edit: flagged for community wiki.

REPLY [1 votes]: What about the smallest non-trivial group $\mathbb{Z}/2\mathbb{Z}$ ? See the discussion Fantastic properties of $\mathbb{Z}/2\mathbb{Z}$ for many convincing arguments.


TITLE: What do Bethuel et al. mean by "diagonal argument" in this context?
QUESTION [1 upvotes]: In their paper (with full text here) Bethuel et al. use an argument on page 636 that (simplified) reads:
Assume you have a sequence of functions $(v^n)_{n \in \mathbb{N}}$ with domain $\mathbb{R}^N$, that are uniformly bounded in $C^k$, i.e. for every $k \in \mathbb{N}$ there is a constant $K(k)$, s.t. \begin{equation} \Vert v^n \Vert_{C^k(\mathbb{R}^N)} \leq K(k).\end{equation}
They then use Ascoli's theorem to assert that for every compact ball $B(0,j)$ there is a smooth map $\mathfrak{v}^j$ s.t. for $n \to \infty$ \begin{equation} v^n \to \mathfrak{v}^j ~~~ \mbox{ in } C^k(B(0,j)), \forall k \in \mathbb{N}.\end{equation}
This part I do understand. The next step is to show that there is a function $\mathfrak{v}$, s.t. $v^n \to \mathfrak{v}$ in $C^k(K)$ for any $K \subset \subset \mathbb{R}^N$. Their (full) argument is:

"To conclude we let $j \to \infty$, and we invoke a diagonal argument, so that in particular $\mathfrak{v} = \mathfrak{v}^j$ does not depend on the ball $B(0,j)$."

What do they mean by "diagonal argument"? Can anyone provide me with any source containing an example for such an argument?

REPLY [1 votes]: The argument the closest to being "diagonal" that I can think of is that $v^n(x)\to\mathfrak v^j(x)$ for every $j$ and every $x$ in $B(0,j)$. By the uniqueness of the pointwise limit, for every $i$ and $j$, the functions $\mathfrak v^i:B(0,i)\to\mathbb R$ and $\mathfrak v^j:B(0,j)\to\mathbb R$ coincide on $B(0,\min(i,j))$.
One can define a function $\mathfrak v:\mathbb R^N\to\mathbb R$ unambiguously by $\mathfrak v(x)=\mathfrak v^j(x)$ for every $x$ and every $j\gt|x|$. Then $v^n\to\mathfrak v$ in $C^k(B(0,j))$, for every $j$, hence $v^n\to\mathfrak v$ in $C^k(K)$, for every bounded subset $K\subset\mathbb R^N$.


TITLE: During $200$ test demands zero failures were discovered. What is the posterior probability that the software meets the business requirement?
QUESTION [0 upvotes]: A business software system is required to achieve at most $1$ failure in $1000$ demands (this is equivalent to a probability on demand: $pfd=10^{-3}$). The business manager wants to to evaluate the system before delivery using a combination of expert judgement about the process applied in its development and failure data collected during testing.
From an analysis of the development process applied, the business manager believes that there is a $70$% chance that the system will meet the requirement and a $30$% chance it is an order of magnitude worse than that ($pfd=10^{-2}$). If the $pfd$ can take two values only, this is expressed as the prior:
$P(pfd=10^{-3})=0.7$
$P(pfd=10^{-2})=0.3$
The testing process is assumed to compromise of a sequence of independent test demands, each of which can result in failure or success.
The likelihood of observing $f$ failures in $d$ tests is defined by the binomial distribution conditioned on the number of demands, $d$, and the $pfd$. Thus:
$P(f|d,pfd)=\frac{d!}{f!(d-f)!}(pfd)^{f}(1-pfd)^{d-f}$
The posterior distribution for $pfd$, using the Bayes theorem, is:
$P(pfd|d,f)=\frac{P(f|pfd,d)P(pfd)}{P(f)}$
$a)$ State the formula needed to calculate the marginal probability of $P(f).$
$b)$ During $200$ test demands zero failures were discovered. What is the posterior probability that the software meets the business requirement?
$My$ $Working:$
Let $F$ be random variable representing number of failures in $d$ number of tests. then:
$Pr(success)=0.7$ ; $Pr(failure)=0.3$
$a)$ Honestly I don't even know how to begin with. There is only one random variable involved as far as I understand, but the question is asking for marginal distribution which exists when there are at least two random variables involved.
$b)$ Since there are zero failures, this means $pfd=0.7$ and From the given binomial distribution we obtain: $P(0|200,0.7)=\frac{200!}{0!(200-0)!}0.7^{0}(1-0.7)^{200-0} = 0.3^{200}$; Plugging this value in posterior distribution gives:
$P(0.7|200,0)=\frac{P(0|0.7,200)P(pfd)}{P(0)}$
I think I am jumping here and there. I am not sure whether my calculations are right. Any help and guidance will be appreciated

REPLY [1 votes]: From the binomial distribution, you have $P(f=0 \mid pfd, n=200) = (1-pfd)^{200}$
and a prior of $P(pfd=10^{-3})=0.7$ and  $P(pfd=10^{-2})=0.3$
Since $P(f) = \sum\limits_k P(f \mid pfd_k, n=200) P(pfd_k)$
you get a posterior distribution of $$P(pfd=10^{-3} \mid f=0, n=200)=\dfrac{0.7 \times (1-10^{-3})^{200}}{0.7 \times (1-10^{-3})^{200}+0.3 \times (1-10^{-2})^{200}}$$ which is about $0.934$


text:
\section{Experiments}
\label{sec:xp}

The code and data are available online\footnote{\url{https://gitlab.com/nnadisic/ssnmf}}.

All experiments have been performed on a personal computer with an i5
processor, with a clock frequency of 2.30GHz.
All algorithms are single-threaded.
All are implemented in Matlab, except the sparse NNLS solver \arbo{}, which is implemented in C++ with a Matlab MEX interface.

As far as we know, no algorithm other than \brassens{} can tackle SSNMF with more than one interior point (see Section~\ref{sec:relwork}) hence comparisons with existing works are unfortunately limited.
For example, separable NMF algorithms can only identify the exterior vertices; see Section~\ref{sec:intro}.
However, we will compare \brassens{} to SNPA on a real multispectral image in Section~\ref{sec:xp-real}, to show the advantages of the SSNMF model over separable NMF.
In Section~\ref{sec:xp-synth}, we illustrate the correctness and efficiency of \brassens{} on synthetic data sets.

\subsection{Synthetic Data Sets}\label{sec:xp-synth}

In this section, we illustrate the behaviour of \brassens{} in different experimental setups.
The generation of a synthetic data set is done as follows: for a given number of dimensions $m$, number of vertices $r$, number of data points $n$, and data sparsity $k$, we generate matrices $W \in \mathbb{R}_+^{m \times r}$ such that the last $r-m$ columns of $W$ are linear combinations of the first $m$ columns, and $H \in \mathbb{R}_+^{r \times n}$ such that $H = [ I_r , H' ]$ and $\|H(:,j)\|_0 \leq k$ for all $j$.
We use the uniform distribution in the interval [0,1] to generate random numbers (columns of $W$ and columns of $H$), and then normalize the columns of $W$ and $H$ to have unit $\ell_1$ norm.
We then compute $M = W H$.
This way, the matrix $M$ is $k$-sparse $r$-separable, with $r$ vertices, of which $r-m$ are interior vertices (in fact, the first $m$ columns of $W$ are linearly independent with probability one as they are generated randomly).
We then run \brassens{} on $M$, with the parameter $k$, and no noise-tolerance.
For a given setup, we perform 30 rounds of generation and solving, and we measure the median of the running time and the median of the number of candidates extracted by kSSNPA.
This number of candidates corresponds to $| \mathcal{J}' |$ in \cref{alg:ssnmf}, that is, the number of interior points selected by kSSNPA as potential interior vertices.
Note that a larger number of candidates only results in an increased computation time, and does not change the output of the algorithm which is guaranteed to extract all vertices (\cref{th:recovery}).

\cref{fig:xpnvar} shows the behaviour of \brassens{} when $n$ varies, with fixed $m=3$, $k=2$, and $r=5$.
To the best of our knowledge, this case is not handled by any other algorithm in the literature.
Both the number of candidates and the run time grow slower than linear.

The irregularities in the plot are due to the high variance between runs.
Indeed, if vertices are generated in a way that some segments between vertices are very close to each other, \brassens{} typically selects more candidates before identifying all columns of $W$.
\begin{figure}[ht!]
    \centering
    \input{img/xp-nvar}
    \caption{Results for \brassens{} on synthetic data sets for different values of $n$, with fixed $m=3$, $k=2$, and $r=5$ with 3 exterior and 2 interior vertices. The values showed are the medians over 30 experiments.}\label{fig:xpnvar}
\end{figure}

In \cref{tab:synth2} we compare the performance of \brassens{} for several sets of parameters.
The number of candidates grows relatively slowly as the dimensions $(m,n)$ of the problem increase, showing the efficiency of the screening performed by kSSNPA.
However, the run time grows rather fast when the dimensions $(m,n)$ grow.
This is because, not only the number of NNLS subproblems to solve increase, but also their size.
\setlength{\tabcolsep}{6pt}
\begin{table}[ht]
\centering
\caption{Results for \brassens{} on synthetic data sets (median over 30 experiments).}
\label{tab:synth2}
\begin{tabular}{c|c|c|c|c|c}
m & n  & r  & k & Number of candidates & Run time in seconds \\ \hline
3 & 25 & 5  & 2 & 5.5                  & 0.26                \\
4 & 30 & 6  & 3 & 8.5                  & 3.30                \\
5 & 35 & 7  & 4 & 9.5                  & 38.71               \\
6 & 40 & 8  & 5 & 13                   & 395.88              \\

\end{tabular}
\end{table}
In all cases, as guaranteed by Theorem~\ref{th:recovery}, \brassens{} was able to correctly identify the columns of $W$. Again, as far as we know, no existing algorithms in the literature can perform this task.

To summarize, our synthetic experiments show the efficiency of the screening done by kSSNPA, and the capacity of \brassens{} to handle medium-scale data sets.

\subsection{Blind Multispectral Unmixing}\label{sec:xp-real}

A multispectral image is an image composed of various wavelength ranges, called \emph{spectral bands}, where every pixel is described by its \emph{spectral signature}.
This signature is a vector representing the amount of energy measured for this pixel in every considered spectral band.
Multispectral images usually have a small number of bands (between 3 and 15).
These bands can be included or not in the spectrum of visible light.
In the NMF model, if the columns of $M$ are the $n$ pixels of the image, then its rows represent the $m$ spectral bands.

The unmixing of a multispectral image consists in identifying the different materials present in that image.
When the spectral signatures of the materials present in the image are unknown, it is referred to as \emph{blind unmixing}.
The use of NMF for blind unmixing of multispectral images relies on the linear mixing model, that is, the assumption that the spectral signature of a pixel is the linear combination of the spectral signatures of the material present in this pixel.
This corresponds exactly to the NMF model, which is therefore able to identify both the materials present in the image ($W$) and the proportions/abundances of materials present in every pixel ($H$); see~\cite{bioucas2012hyperspectral,ma2013signal} for more details.

Let us apply \brassens{} to the unmixing of the well-known Urban satellite image~\cite{fyzhu2014hyperspectraldata}, composed of
$309 \times 309$ pixels. The original cleaned image has 162 bands, but we only keep 3 bands, namely the bands 2, 80, and 133 -- these were obtained by selecting different bands with SPA applied on $M^T$ -- to obtain a data set of size $3 \times \num{94249}$. The question is: can we still recover materials by using only 3 bands? (The reason for this choice is that this data set is well known and the ground truth is available, which is not the case of most multispectral images with only 3 bands.)
We first normalize all columns of $M$ so that they sum to one.
Then, we run \brassens{} with a sparsity constraint $k = 2$ (this means that we assume that a pixel can be composed of at most 2 materials, which is reasonable for this relatively high resolution image) and a noise-tolerance threshold of 4\%; this means that we stop SNPA and kSSNPA when
$\|M - M(:,\mathcal{J}) H\|_F \leq 0.04 \|M\|_F$.
\brassens{} extracts 5 columns of the input matrix.
For comparison, we run SNPA with $r = 5$.
Note that this setup corresponds to underdetermined blind unmixing, because $m = 3 < r = 5$.
It would not be possible to tackle this problem using standard NMF algorithms (that would return a trivial solution such as $M = I_3 M$).
It can be solved with SNPA, but SNPA cannot identify interior vertices.

SNPA extracts the 5 vertices in \num{3.8} seconds.
\brassens{} extracts 5 vertices, including one interior vertex, in \num{33} seconds.
The resulting abundance maps are showed in \cref{fig:mapbrassens}.
They correspond to the reshaped rows of $H$, hence they show which pixel contains which extracted material (they are more easily interpretable than the spectral signatures contained in the columns of $W$).
The materials they contain are given in \cref{tab:unmixing}, using the ground truth from~\cite{zhu2017hyperspectral}.
We see that \brassens{} produces a better solution, as the materials present in the image are better separated: the first three abundance maps of \brassens{} are sparser and correspond to well-defined materials. The last two abundances maps of SNPA and of \brassens{} are similar but extracted in a different order.
The running time of \brassens{} is reasonable, although ten times higher than SNPA.

\begin{figure}[ht!]
    \centering
    \includegraphics[width=\textwidth]{img/mapsnpa.png}
    \vspace{0.1cm}\\
    \includegraphics[width=\textwidth]{img/mapbrassens.png}
    \caption{Abundances maps of materials (that is, reshaped rows of $H$) extracted by SNPA (top) and by \brassens{} (bottom) in the Urban image with only 3 spectral bands.}\label{fig:mapbrassens}
\end{figure}

\begin{table}[ht!]
\centering
\caption{Interpretation of the unimixing results from Figure~\ref{fig:mapbrassens}.}\label{tab:unmixing}
\begin{tabular}{c|c|c}
Image & Materials extracted by SNPA & Materials extracted by \brassens{} \\ \hline
1     & Grass + trees + roof tops   & Grass + trees                   \\
2     & Roof tops 1                 & Roof tops 1                     \\
3     & Dirt + road + roof tops     & Road                            \\
4     & Dirt + grass                & Roof tops 1 and 2 + road        \\
5     & Roof tops 1 + dirt + road   & Dirt + grass
\end{tabular}
\end{table}

config: arxiv

file: 2006.07553/arXiv SSNMF/60xp.tex


text:
\begin{document}

\title[Nonholonomic Constraints: a New Viewpoint]{Nonholonomic
Constraints: a New Viewpoint}

\author[J. Grabowski]{J. Grabowski}
\address{Janusz Grabowski: Polish Academy of Sciences, Institute of Mathematics,
 \'Sniadeckich 8, P.O. Box 21, 00-956 Warszawa, Poland
} \email{jagrab@impan.gov.pl}

\author[M. de Le\'on]{M. de Le\'on}
\address{Manuel de Le\'on:
Instituto de Ciencias Matem\'aticas (CSIC-UAM-UC3M-UCM),  Serrano 123, 28006
Madrid, Spain} \email{mdeleon@imaff.cfmac.csic.es}

\author[J.\ C.\ Marrero]{J.\ C.\ Marrero}
\address{Juan Carlos Marrero:
Departamento de Matem\'atica Fundamental y Unidad Asociada
ULL-CSIC Geo\-metr\'{\i}a Diferencial y Mec\'anica Geom\'etrica,
Facultad de Ma\-te\-m\'a\-ti\-cas, Universidad de la Laguna, La
Laguna, Tenerife, Canary Islands, Spain} \email{jcmarrer@ull.es}

\author[D.\ Mart\'{\i}n de Diego]{D. Mart\'{\i}n de Diego}
\address{David Mart\'{\i}n de Diego:
Instituto de Ciencias Matem\'aticas (CSIC-UAM-UC3M-UCM),  Serrano 123, 28006
Madrid, Spain} \email{d.martin@imaff.cfmac.csic.es}

\keywords{General algebroids, Lie algebroids, double vector
bundle, Nonholonomic mechanics, Lagrange-d'Alembert's equations,
nonholonomic bracket.}

\subjclass[2000]{70F25, 53D17, 70G45, 17B66, 70H03, 70H45}

\date{}

\thanks{This work has been partially supported by the Polish
Ministry of Science and Higher Education
under the grant No. N201 005 31/0115, MEC (Spain)
Grants MTM 2006-03322, MTM 2007-62478, project Ingenio
Mathematica (i-MATH) No. CSD 2006-00032 (Consolider-Ingenio 2010)
and S-0505/ESP/0158 of the CAM}

\begin{abstract}
The purpose of this paper is to show that, at least for
Lagrangians of mechanical type, nonholonomic Euler-Lagrange
equations for a nonholonomic linear constraint $D$ may be viewed
as non-constrained Euler-Lagrange equations but on a new
(generally not Lie) algebroid structure on $D$. The proposed novel
formalism allows us to treat in a unified way a variety of
situations in nonholonomic mechanics and gives rise to a version
of Neoether Theorem producing actual first integrals in case of
symmetries.
\end{abstract}

\maketitle

\tableofcontents
\section{Introduction}
There are many approaches to geometric mechanics in the
literature. We will work with a natural generalization of the
framework for studying mechanical systems proposed by
W.~M.~Tulczyjew \cite{Tul1,Tul2} (see also \cite{TU} and
references therein). In the simplest form, the phase dynamics of
the system is understand as the lagrangian submanifold $\zG$ of
the symplectic manifold $(\sT\sT^*M,\dt\omega_M)$ equipped with
tangent lift $\dt\omega_M$ of the canonical symplectic form
$\omega_M$ of $\sT^*M$. Here, $M$ represents the configuration
manifold of the system and $\zG$ is obtained from the lagrangian
submanifold $\xd L(M)\subset \sT^*\sT M$ induced by the Lagrangian
$L:\sT M\ra\R$ via the canonical isomorphism $\ze_M:\sT^*\sT M\ra
\sT\sT^*M$. In other words, the phase dynamics, as well as the
Euler-Lagrange equations, are obtained in a simple way by means of
the {\it Tulczyjew differential} $\zL_L=\ze_M\circ\xd
L:M\ra\sT\sT^* M$. It is important to observe that both
$\sT\sT^*M$ and $\sT^*\sT M$ are double vector bundles over
$\sT^*M$ and $\sT M$ (see \cite{GU2} and references therein). The
resulting submanifold $\zG$ of $\sT\sT^*M$ is a particular case of
modelling dynamical systems as implicit differential equations
defined by differential inclusions (see \cite{MMT1,MMT2}).

This framework admits an immediate generalization for more general
morphisms $\ze:\sT^* E\ra\sT E^*$ of canonical double vector
bundles associated with a vector bundle $E$ (see \cite{GGU,GG}),
inducing the Tulczyjew differentials $\zL_L=\ze\circ\xd L:M\ra\sT
E^*$. This generalization includes as a particular case a theory
of mechanical systems based on Lie algebroids, as proposed by
A.~Weinstein \cite{We} and developed by many authors, however in a
different geometrical setting (see, for instance,
\cite{CoLeMaMaMa,CoLeMaMa,LMM,Li,mart,Medina}). The motivation for
study systems on Lie algebroids is that they often appear
naturally as results of some reduction procedures. This is is a
situation similar to the one known in the theory of Hamiltonian
systems: reductions may lead from a symplectic to a Poisson
structure.

An additional challenge and one of the most fascinating topics in
geometric mechanics is the study of constraints in this context.
Of course, a general problem of putting constraints for the system
in a variational setting involves constraints for velocities as
well as constraints for virtual displacements, as was noticed
already in \cite{Tu2}. In some cases, however, one assumes that
the constraints can be determined from a constraint subset $D$ of
$\sT M$ (or, of $E$ in the algebroid case) by certain
well-described procedures. The best known approaches of this type
refer to the so called {\it vakonomic} and {\it nonholonomic}
constraints. In the simplest situation, for $D$ being a linear
nonholonomic constraint, i.e. just a vector subbundle of $\sT M$
(or, of $E$ in the algebroid context), this procedure describes
the nonholonomic Euler-Lagrange equations by means of the {\it
d'Alembert principle}, having analogs also in the algebroid case
\cite{CoLeMaMa,LMMdD,M4,GG}. We should stress that our
nonholonomic constraints are linear in the broader sense, i.e.
they are subbundles over submanifolds of the original base
manifold. The nonholonomic Euler-Lagrange equations are commonly
viewed as being not variational equations. In \cite{GG} it has
been pointed out that it is not exactly the case, if we extend
slightly our understanding of Variational Calculus.

In this paper we continue studying nonholonomic constraints on
algebroids and showing that, at least for Lagrangians of
mechanical type, the nonholonomic Euler-Lagrange equations are
just non-constrained Euler-Lagrange equations but for a special
algebroid structure on the constraint subbundle (see also
\cite{LMMdD}). This shows that mechanical systems based on general
(not necessary Lie) algebroids appear naturally in the presence of
nonholonomic constraints and gives a powerful geometrical tool
when dealing with constrained systems. In particular, we get a
version of Noether Theorem with true first integrals for
nonholonomic systems. We do not get all possible algebroids on $D$
applying our procedure. In particular, if the original structure
was a Lie algebroid, then the new algebroid bracket is
automatically skew- symmetric, so we deal with a {\it quasi-Lie
algebroid}. One can the associate with the sequence of procedures,
like reduction by symmetries and passing to a nonholonomic
constraint, the sequence of the corresponding novel structures
serving as appropriate geometrical tools in describing the
systems:
$$ \sT M\quad \xrightarrow[]{\text{reduction by symmetries}}\quad \text{Lie algebroid} \quad
\xrightarrow[]{\text{nonholonomic constraint}}\quad\text{quasi-Lie
algebroid}\,.
$$
All this is of course closely related to the discovery of the role
of the nonholonomic quasi-Poisson brackets
\cite{Mr,SM,IbLeMaMa,CaLeMa}, this time not in the Hamilton but in
the Lagrange picture.

The paper is organized as follows. In the next section, we recall
after \cite{GGU} the basic ideas of developing mechanics on a
general algebroid $E$, in particular, the Euler-Lagrange
equations. Then, we construct in this setting an analog of the
Tulczyjew differential for linear nonholonomic constraints,
together with the corresponding nonholonomic Euler-Lagrange
equations. In Section 4 we study reductions of a general algebroid
to an algebroid on a nonholonomic constraint $D$ (satisfying a
natural admissibility condition) along a given projection. We
discuss also the problem, what algebroids can be obtained in this
way, if we start with Lie algebroids.

Section 5 is the most important part of our paper. For Lagrangian
functions $L$ of mechanical type (`kinetic energy - potential'),
we show that the non-constrained Euler-Lagrange equations on $D$,
derived for the reduced Lagrangian $l=L_{\mid D}$ and for the
reduced algebroid structure on $D$ along the orthogonal projection
associated with the kinetic energy, coincide with the nonholonomic
Euler-Lagrange equations. Moreover, for this nonholonomic case, we
can apply therefore the generalization of the Noether Theorem
proved in \cite{GGU} to obtain actual first integrals. Passing to
the nonholonomic constraints does not requires therefore any
change in our unified algebroid approach to mechanics:
nonholonomic Euler-Lagrange equations are included in our
framework.

We end up with two well-known examples of nonholonomic
constraints, the Chaplygin sleigh and the snakebord, to show how
simply the corresponding equations of motions can be derived by
means of our method.
\section{Geometric Mechanics on general algebroids}

Let $M$ be a smooth manifold and $(x^i)$, $i=1, \ldots, n$, be a
coordinate system in $M$. Denote by $\tau_M\colon \sT M\rightarrow
M$ the tangent vector bundle, with induced coordinates $(x^i,
\dot{x}^j)$, and by $\pi_M\colon \sT^\ast M \rightarrow M$ the
contangent bundle, with induced coordinates $(x^i, p_j)$.

Let $\tau: E\rightarrow M$ be a vector bundle and $\pi: E^\ast\rightarrow M$ its dual
bundle. Taking a local basis $\{e_a\}$ of sections of $E$, then we have the
corresponding local coordinates $(x^i, y^a)$ on $E$, where $y^a(e)$ is the
$a$th-coordinate of $e\in E$ in the given basis. We denote by $(x^i, \zx_a)$ the
corresponding coordinates of the dual bundle $\pi: E^\ast\rightarrow M$. One can also
say that $\zx_a$ is the fiber-wise linear local function $\zi(e_a)$ on $E^\ast$
corresponding to the local section $e_a$ of $E$. We have also adapted local
coordinates:
    \beas
     (x^i, y^a,{\dot x}^j, {\dot y}^b ) \quad \text{in} \ \sT E\,,& \qquad
    (x^i, \zx_a, {\dot x}^j, {\dot \zx}_b) \quad \text{in} \ \sT E^\* \,,\\
    (x^i, y^a, p_j, \zp_b)\quad \text{in}\ \sT^\*E\,,&\qquad
    (x^i, \zx_a, p_j, \zf^b)\quad \text{in}\ \sT^\* E^\* \,.
    \eeas

It is well known (cf. \cite{KU}) that the cotangent bundles $\sT^\*E$ and $\sT^\*E^\*$
are examples of double vector bundles:
$$\xymatrix{
\sT^\ast E^\ast\ar[rr]^{\sT^\ast\zp} \ar[d]_{\zp_{E^\ast}} && E\ar[d]^{\zt} \\
E^\ast\ar[rr]^{\zp} && M } \qquad {,}\qquad \xymatrix{
\sT^\ast E\ar[rr]^{\sT^\ast\zt} \ar[d]_{\zp_{E}} && E^\ast\ar[d]^{\zp} \\
E\ar[rr]^{\zt} && M }.
$$
The concept of a double vector bundle is due to J.~Pradines
\cite{Pr1,Pr2}, see also \cite{Ma,KU}. In particular, all arrows
correspond to vector bundle structures and all pairs of vertical
and horizontal arrows are vector bundle morphisms. The double
vector bundles have been recently characterized \cite{GR} in a
simple way as two vector bundle structures whose Euler vector
fields commute. The above double vector bundles are canonically
isomorphic with the isomorphism
    \be\label{iso}\cR_\zt \colon \sT^\*E \rightarrow \sT^\* E^\*
                    \ee
  being simultaneously an anti-symplectomorphism
  (cf. \cite{KU,GU2}). In adapted local coordinates, $\cR_\zt$ is given by
    $$\cR_\zt(x^i, y^a, p_j, \zp_b) = (x^i, \zp_a, -p_j,y^b).
                              $$
This means that we can identify $\zp_b$ and $\zx_b$, as well as $y^b$ with $\zf^b$, so
we will use local coordinates $(x,y,p,\zx)$ on $\sT^\ast E$ and local coordinates
$(x,\zx,p,y)$ on $\sT^\ast E^\ast$ in full agreement with the isomorphism (\ref{iso}).

For describing the dynamics of a Lagrangian or Hamiltonian systems
it is necessary to give an additional ingredient, typically, a Lie
algebroid structure or, more generally, the algebroid structure.
In several papers, many authors have studied Lie algebroid
structures as a unified geometric framework, general enough to
account for different mechanical systems (defined on tangent
bundles, on Lie algebras, quotients by Lie groups actions, etc.).

\medskip It is well known that Lie algebroid structures on a vector bundle $E$
correspond to linear Poisson tensors on $E^\*$. A 2-contravariant tensor $\zP$ on
$E^\*$ is called {\it linear} if the corresponding mapping $\widetilde{\zP} \colon
\sT^\* E^\* \rightarrow \sT E^\*$ induced by contraction is a morphism of double vector
bundles. This is the same as to say that the corresponding bracket of functions is
closed on (fiber-wise) linear functions. The commutative diagram
$$\xymatrix{
\sT^\ast E^\ast\ar[r]^{\widetilde\Pi}  & \sT E^\ast \\
\sT^\ast E\ar[u]_{\cR_\tau}\ar[ur]^{\ze} & },
$$
composed with (\ref{iso}), describes a one-to-one correspondence between linear
2-contravariant tensors $\zP$ on $E^\*$ and homomorphisms of double vector bundles (cf.
\cite{KU, GU2}) covering the identity on $E^\*$:

\be\xymatrix{
 & \sT^\ast E \ar[rrr]^{\varepsilon} \ar[dr]^{\pi_E}
 \ar[ddl]_{\sT^\ast\tau}
 & & & \sT E^\ast\ar[dr]^{\sT\pi}\ar[ddl]_/-20pt/{\tau_{E^\ast}}
 & \\
 & & E\ar[rrr]^/-20pt/{\zr}\ar[ddl]_/-20pt/{\tau}
 & & & \sT M \ar[ddl]_{\tau_M}\\
 E^\ast\ar[rrr]^/-20pt/{id}\ar[dr]^{\pi}
 & & & E^\ast\ar[dr]^{\pi} & &  \\
 & M\ar[rrr]^{id}& & & M &
}\label{F1.3} \ee In local coordinates, every  $\ze$ as above is
of the form
\be\label{F1.4} \ze(x^i,y^a,p_j,\zp_b) = (x^i, \zp_a,
\sum_c\zr^j_c(x)y^c, \sum_{a,c}{\mathcal C}^c_{ab}(x) y^a\zp_c +
\sum_i\zs^i_b(x) p_i) \ee
and it corresponds to the linear tensor
\be\label{tensor} \zP_\ze =\sum_{a,b,c}{\mathcal C}^c_{ab}(x)\zx_c
\partial _{\zx_a}\otimes
\partial _{\zx_b} + \sum_{i,b}\zr^i_b(x) \partial _{\zx_b} \otimes
\partial _{x^i} - \sum_{a,j}\zs^j_a(x)\partial _{x^j} \otimes
\partial _{\zx_a}\,.
\ee
In \cite{GU2} by {\it algebroids} we meant the morphisms
(\ref{F1.3}) of double vector bundles covering the identity on
$E^\*$, while {\it Lie algebroids} were those algebroids for which
the tensor $\zP_\ze$ is a Poisson tensor. We can consider the {\it
adjoint tensor} $\zP_\ze^+$, i.e. the 2-contravariant tensor
obtained from $\zP_\ze$ by transposition:
$$ \zP_\ze^{+} =\sum_{a,b,c}{\mathcal C}^c_{ba}(x)\zx_c \partial _{\zx_a}\otimes
\partial _{\zx_b} + \sum_{i,b}\zr^i_b(x) \partial _{x^i}\otimes\partial _{\zx_b}
 - \sum_{a,j}\zs^j_a(x) \partial _{\zx_a}\otimes\partial _{x^j}\,,
$$
and the corresponding {\it adjoint algebroid} structure $\ze^+$.
Algebroids $\ze$ corresponding to skew-symmetric $\zP_\ze$
(anti-symmetric brackets $[\cdot,\cdot]_\ze$), i.e. such that
$\ze^{+}=-\ze$, we will call {\it quasi-Lie algebroids}. The
relation to the canonical definition of Lie algebroid is given by
the following theorem (cf. { \cite{GU3, GU2}}).

\begin{theorem}
An algebroid structure $(E,\ze)$ can be equivalently defined as a bilinear bracket $[\cdot ,\cdot]_\ze $ on
sections of $\zt\colon E\rightarrow M$, together with vector bundle morphisms $\zr,\, \zs \colon E\rightarrow
\sT M$ (left and right anchors), such that
$$ [fX,gY]_\ze = f\zr(X)(g)Y -g\zs(Y)(f) X
+fg [X,Y]_\ze
$$
         for $f,g \in \cC^\infty (M)$, $X,Y\in\Sec{E}$.
The bracket and anchors are related to the 2-contravariant tensor $\zP_\ze$ by the
formulae
\beas
        \zi([X,Y]_\ze)&= \{\zi(X), \zi(Y)\}_{\zP_\ze},  \\
        \zp^\*(\zr(X)(f))       &= \{\zi(X), \zp^\*f\}_{\zP_\ze}, \\
        \zp^\*(\zs(X)(f))       &= \{\zp^\* f, \zi(X)\}_{\zP_\ze}.
                                                   \eeas
        The algebroid $(E,\ze)$ is a quasi-Lie algebroid if and only if the tensor $\zP_\ze$ is
        skew-symmetric and it is a Lie algebroid if and only if the tensor
$\zP_\ze$ is a Poisson tensor.
\end{theorem}

The canonical example of a mapping $\ze$ in the case of $E=\sT M$
is given by $\ze = \ze_M = \za^{-1}_M$ -- the inverse to the
Tulczyjew isomorphism $\za_M:\sT\sT^\*M\ra\sT^\*\sT M$. In
general, the algebroid structure map $\ze$ is not an isomorphism
and, consequently, its dual $\zk = \ze^{\ast}$ with respect to the
right projection is a relation and not a mapping. Some authors
(see \cite{LMM} for the case of Lie algebroids) have introduced
the concept of prolongation of a Lie algebroid to maintain some of
the original properties of the Tulzcyjew mapping  (isomorphism and
symplecticity), but, for the purposes of this paper, it is not
necessary to use this formalism.

The double vector bundle morphism (\ref{F1.3}) can be extended to the following
algebroid analogue of the so called Tulczyjew triple
\be\xymatrix@C-5pt{
 &\sT^\ast E^\ast \ar[rrr]^{\widetilde{\zP}}
\ar[ddl]_{\pi_{E^\ast}} \ar[dr]^{\sT^\ast\pi}
 &  &  & \sT E^\ast \ar[ddl]_/-25pt/{\zt_{E^\ast}} \ar[dr]^{\sT\pi}
 &  &  & \sT^\ast E \ar[ddl]_/-25pt/{\sT^\ast\zt} \ar[dr]^{\pi_E}
\ar[lll]_{\varepsilon}
 & \\
 & & E \ar[rrr]^/-20pt/{\zr}\ar[ddl]_/-20pt/{\zt}
 & & & \sT M\ar[ddl]_/-20pt/{\zt_M}
 & & & E\ar[lll]_/+20pt/{\zr}\ar[ddl]_{\zt}
 \\
 E^\ast\ar[rrr]^/-20pt/{id} \ar[dr]^{\pi}
 & & & E^\ast\ar[dr]^{\pi}
 & & & E^\ast\ar[dr]^{\pi}\ar[lll]_/-20pt/{id}
 & & \\
 & M\ar[rrr]^{id}
 & & & M & & & M\ar[lll]_{id} &
}\label{F1.3b}\ee The left-hand side is Hamiltonian, the right-hand side is Lagrangian,
and the dynamics lives in the middle.

We introduce now the dynamics through a  Lagrangian $L:E\ra{\R}$
which defines two smooth maps: the {\em Legendre mapping}: {
$\lambda_L:E\rightarrow E^\ast$,
$\lambda_L=\tau_{E^\ast}\circ\ze\circ\xd L$,} which is covered by
the {\em Tulczyjew differential} { $\Lambda_{L}: E\rightarrow \sT
E^\ast$, $\Lambda_{L}=\ze\circ\xd L$}:
$$\xymatrix{
\sT^\ast E\ar[rr]^{\ze}  && \sT E^\ast\ar[d]^{\zt_{E^\ast}} \\
E\ar[rr]^{\zl_L}\ar[u]^{\xd L}\ar[rru]^{\Lambda_{L}} && E^\ast }.
$$
The lagrangian function defines {\em the phase dynamics} {
${\zG_L}=\Lambda_{L}(E)\subset \sT E^\ast$} which can be
understood as an implicit differential equation on $E^\*$,
solutions of which are `phase trajectories' of the system
$\zb:\R\ra E^\ast$ and satisfy $\bt{\zb}(t)\in {\zG_L}$, where
$\bt{\zb}$ is the tangent prolongation of $\zb$,
$\bt\zb(t)=(\zb(t),\dot{\zb}(t))$. An analog of the {\it
Euler-Lagrange equations} for curves $\gamma: \R\rightarrow E$ is
in turn:
$$(E_L)\qquad \bt(\zl_L\circ\gamma)=\Lambda_L\circ\gamma.$$
In local coordinates, ${\zG_L}$ has the parametrization by
$(x^i,y^a)$ via $\zL_L$ in the form (cf. (\ref{F1.4}))
\be\zL_L(x^i,y^a)= (x^a,\frac{\partial L}{\partial
y^a}(x,y), \sum_c\zr^j_c(x)y^c, \sum_{a,c}{\mathcal C}^c_{ab}(x) y^a\frac{\partial L}{\partial y^c}(x,y) +\sum_i
\zs^i_b(x)\frac{\partial L}{\partial x^i}(x,y)) \label{F1.4a}\ee and the equation $(E_L)$, for
$\zg(t)=(x^i(t),y^a(t))$, reads \be (E_L):\qquad\frac{\xd x^i}{\xd t}=\sum_c\zr^i_c(x)y^c, \qquad
\frac{\xd}{\xd t}\left(\frac{\partial L}{\partial y^b}\right)= \sum_{a,c}{\mathcal C}^c_{ab}(x) y^a\frac{\partial
L}{\partial y^c}(x,y) + \sum_i\zs^i_b(x)\frac{\partial L}{\partial x^i}(x,y)\,.\label{EL2}
\ee
As one can see from (\ref{EL2}), the solutions are automatically
{\it admissible curves} in $E$, i.e. the velocity
$\bt(\zt\circ\zg)(t)$ is $\zr(\zg(t))$.

With this framework it is possible to write in a unified point of
view many equations of different mechanical systems that usually,
in the literature, appears as different ones (Classical
Euler-Lagrange equations, Lagrange-Poincar\'e equations after
reduction by the action of a Lie group, Euler-Poincar\'e
equations, etc.). (See \cite{CoLeMaMaMa,LMM} for applications of
the theory in the case when $(E,\ze)$ is  a Lie algebroid).

\section{Nonholonomic mechanics}\label{koi}
We can start from a general Lie algebroid $(E,\ze)$ keeping in
mind the standard case $E=\sT M$. A (linear) nonholonomic
Lagrangian system is determined by a Lagrangian function $L: E\ra
\R$ and a vector subbundle $D$, $\hbox{rank}\, D=n-r$, of the
bundle $E$. We will accept subbundles over a submanifold, so let
us denote $D_M=\zt(D)$. By $i_D: D\hookrightarrow E$ let us denote
the inclusion and by $i^*_D: E^*_{\mid D_M}\rightarrow D^*$ the
dual map.

Because the solutions of the dynamics in $E$ should be admissible
curves, we need an {\it admissibility condition} ensuring that
there are admissible curves through every point of $D$. The
natural condition we take is $\zr(D)\subset \sT D_M$. By {\it
strong admissibility condition} we will mean that $D$ satisfies
the integrabilty condition with respect to both algebroid
structures: $\ze$ and $\ze^+$, i.e.
\be\label{ic}\zr(D)\subset \sT
D_M\quad\text{and}\quad\zs(D)\subset \sT D_M\,.
\ee
The Lagrangian function $L: E\rightarrow \R$ and the
vector subbundle $D$ define also the smooth map -- the {\it
constrained Tulczyjew differential} \be\label{nhc}\zL_L^D:
{D}\rightarrow \sT D^\ast,\qquad \zL_L^D=\sT
i_D^\ast\circ\ze\circ\xd L\ee covering the {\it constrained
Legendre map} \be\label{nh0}\zl_L^D:D\ra D^*,\qquad
\zl_L^D=i_{D}^\ast \circ \zl_L\,.
\ee
The diagram picture is the following
\be\label{nhc1} \xymatrix{ \sT^\ast E_{\mid D}\ar[rr]^{\ze}  && \sT (E^\ast_{\mid
D_M})\ar[rr]^{\sT{i_D^\ast}}&&
\sT D^\ast\ar[d]^{\zt_{D^\ast}}\\
D\ar[rr]^{\zl_L}\ar[u]^{\xd
L}\ar@{.>}[rrrru]^{\zL_L^D}\ar[rr]^{\zl_{L}}\ar@/d3ex/@{.>}[rrrr]_{\zl_L^D} &&
E^\ast_{\mid D_M}\ar[rr]^{i_D^\ast} && D^\ast}
\ee

\begin{definition}{\rm The bundle $D^\*$ we call
{\it the phase space} of the nonholonomic system and the implicit
differential equation ${\zG^D_L}=\zL_L^D(D)$, being a subset of
$\sT D^*$, we interpret as the {\it phase dynamics} of the
nonholonomic system. A curve $\gamma: I\rightarrow D$ is a
solution of the {\it nonholonomic Euler-Lagrange equation} if and
only if $\zL_L^D\circ \gamma$ is an admissible curve in $\sT
D^\*$, i.e. the nonholonomic Euler-Lagrange equation reads
\be\label{nEL}(nE_L):\qquad\qquad\bt(\zl^D_L\circ\zg)=\zL_L^D\circ \gamma\,.\ee}
\end{definition}

\bigskip To find the explicit form of the nonholonomic Euler-Lagrange equation (\ref{nEL}),
consider local coordinates $(x^I)=(x^i,x^\zi)$ on a open set ${U}$
of $M$ such that $D_M$ is
 determined by the constraint $x^\zi=0$.
A local basis $\{e_a\}_{a=1,\dots,n-r}$ of sections of $D$ we can
extend to local sections of $E$ and complete them to a local basis
of sections $\{e_a, e_{\alpha}\}$ of the vector bundle $E$. Then,
in coordinates $(x^I,y^A)=(x^i,x^\zi,y^a,y^\alpha)$ adapted to
this bases, the local equations defining the constrained subbundle
$D$ as a vector subbundle of $E$ over $D_M$ are $x^\zi=0$,
$y^\alpha=0$. Note that admissibility of the constraint $D$ means
that $\zr_a^\zi(x^i,0)=0$.

Taking local coordinates $(x^i, y^a)$ on $D$, we may write then
$i_D: D\hookrightarrow E$ as $i_D(x^i, y^a)=(x^i, 0, y^a, 0)$ and
$i_D^*(x^i,0,\zx_a,\zx_\za)=(x^i,\zx_a)$, so
$$\sT i_D^*(x^i,0,\zx_A,\dot{x}^j,0,\dot{\zx}_B)=(x^i,\zx_a,\dot{x}^j,\dot{\zx}_b)\,.$$
For the adapted local coordinates
$(x^i,x^\zi,y^a,y^\za,p_j,p_\zg,\zx_b,\zx_\zb)$ in $\sT^*E$, the
map $\ze$ reduced to $(\sT^*E)_{\mid D}$ takes values in
$\sT(E^*_{\mid D_M})$ (admissibility) and reads
$$\ze(x^i,0,y^a,0,p_J,\zx_B)=(x^i,0,\zx_A,\zr_d^j(x^i,0)y^d,0,{\mathcal
C}^C_{dB}(x^i,0)y^d\zx_C+\zs^J_B(x^i,0)p_J)\,.
$$
Therefore
\be\label{yy}\sT i_D^*\circ\ze(x^i,0,y^a,0,p_J,\zx_B)=(x^i,\zx_a,\zr_d^j(x^i,0)y^d,{\mathcal
C}^C_{db}(x^i,0)y^d\zx_C+\zs^J_b(x^i,0)p_J)
\ee
and
\begin{eqnarray*}
&&\zL_L^D(x^i, y^a)\\
&&=Ti_D^*(\ze(x^i, 0,y^a, 0, \frac{\partial L}{\partial
x^J}(x^i,0, y^a, 0), \frac{\partial L}{\partial
y^B}(x^i,0, y^a, 0))\\
&&=(x^i, \frac{\partial L}{\partial y^b}(x^i, 0,y^a, 0),
\rho^j_d(x^i,0) y^d, {{\mathcal C}}^C_{db}(x^i,0) y^d
\frac{\partial L}{\partial y^C}(x^i, 0,y^a,
0)+\sigma^J_b(x^i,0)\frac{\partial L}{\partial x^J}(x^i,0, y^a,
0))
\end{eqnarray*}
Therefore, locally,  the nonholonomic Euler-Lagrange equations read:
\begin{eqnarray}\label{nELe1}&\qquad
x^\zi=0,\quad y^\za=0,\quad \frac{d x^j}{dt}=\rho^j_d(x^i,0)y^d,\\
 &\frac{d}{dt}\left(\frac{\partial L}{\partial
y^b}(x^i, 0,y^a, 0)\right)= {{\mathcal C}}^C_{db}(x^i,0) y^d
\frac{\partial L}{\partial y^C}(x^i,0, y^a,
0)+\sigma^J_b(x^i,0)\frac{\partial L}{\partial x^J}(x^i, 0,y^a,
0))\,.\label{nELe2}
\end{eqnarray}
In the case of a Lie algebroid $\sigma^i_A= \rho^i_A$, and if the
subbundle $D$ is over the whole base manifold $M$, $D_M=M$, the
previous equations are precisely the nonholonomic equations
obtained in \cite{CoLeMaMa} (see Equations 3.8).

\section{The nonholonomic reduction of algebroids}

It has been recognized a long time ago \cite{Mr,SM} that
nonholonomic constraints may lead to certain {\it nonholonomic
brackets} that do not satisfy the Jacobi identity. We will show
that linear nonholonomic systems of mechanical type on general
algebroids are again systems on general algebroids. But even if we
start with a Lie algebroid, the new algebroid is, in general, no
longer a Lie algebroid but certain quasi-Lie algebroid associated
with a linear bi-vector field. This observation, made already in
\cite{GG} (see also \cite{LMMdD}), puts new light to the role of
{\it quasi-Poisson brackets}, i.e. the brackets represented by
arbitrary bi-vector fields and not satisfying, in general, the
Jacobi identity (see \cite{GG,LMMdD}). This shows, on the other
hand, that developing Mechanics on general algebroids makes sense,
as the reduction to a nonholonomic constraint will move us, in
general, from the Lie algebroid picture into a more general one.

Let $\ze$ be an algebroid structure on a vector bundle $E$ over
$M$ associated with the tensor $\zP_\ze$ with the local form
(\ref{tensor}). For a linear subbundle $D$ in $E$ over a
submanifold $D_M\subset M$, satisfying the strong admissibility
condition (\ref{ic}), consider a decomposition
\be\label{dec}E_{\mid D_M}=D\oplus_{D_M} D^\perp
\ee
and the associated projection $P:{E}_{\mid D_M}\ra{D}$. We can
construct the morphism ${\ze_P}: \sT^*D\rightarrow \sT D^*$ of
double vector bundles as follows.

To have the corresponding expressions in local coordinates,
consider local coordinates $(x^I)=(x^i,x^\zi)$ on a open set ${U}$
of $M$ such that $D_M$ is determined by the constraint $x^\zi=0$.
Local bases $\{e_a\}_{a=1,\dots,n-r}$ and
$\{e_\za\}_{\za=n-r+1,\dots,n}$ of sections of $D$ and $D^\perp$,
respectively, can be extended to a basis $e_A$ of local sections
of $E$. Then, we get the coordinates
$(x^I,y^A)=(x^i,x^\zi,y^a,y^\alpha)$ adapted to these bases and
the adapted coordinates
$(x^I,y^A,p_J,\zx_B)=(x^i,x^\zi,y^a,y^\alpha,p_j,p_\zg,\zx_b,\zx_\zb)$
in $\sT^\ast E$. We get automatically coordinates $(x^i,y^a)$ in
$D$ and $(x^i,y^a,y^\za)$ in $E_{\mid D_M}$, and the adapted
coordinates $(x^i,y^a,p_j,\zx_a)$ and $(x^i,y^A,p_j,\zx_B)$ in
$\sT^\ast D$ and $\sT^\ast(E_{\mid D_M})$, respectively.

Let us consider the phase lift $\sT^*P:\sT^*D\rel\sT^*(E_{\mid
D_M})$ which, as often happens which phase lifts, is not a map but
only a relation. In our local coordinates,
$$\sT^*P(x^i,y^a,p_j,\zx_b)=(x^i,y^A,p_j,\zx_b,0)\,.$$
We have also the embedding $i_{E_{\mid D_M}}:E_{\mid D_M}\ra E$
whose phase lift, restricted to $\sT^*E_{\mid D}$,
$$\sT^*(i_{E_{{\mid D}_{M}}}):\sT^*E_{\mid D}
\rel\sT^*(E_{\mid D_M})
$$
is the relation
$$\sT^*(i_{E_{{\mid
D}_{M}}})(x^i,0,y^a,0,p_J,\zx_B)=(x^i,y^A,p_j,\zx_B)\,.
$$
The composition of relations
$$\sT^*(i_{E_{{\mid D}_{M}}})^{-1}\circ\sT^*P:\sT^*D\rel\sT^*E_{\mid D}
$$
has the local form
$$\left(\sT^*(i_{E_{{\mid D}_{M}}})^{-1}\circ\sT^*P\right)(x^i,y^a,p_j,\zx_b)=
(x^i,0,y^a,0,p_J,\zx_b,0)\,.
$$
But, as $\zs_b^\zi(x^i,0)=0$ (strong admissibility condition), the
value of $\sT i_D^*\circ\ze$ in (\ref{yy}) does not depend on
$p_\zg$ and the composition
\be\label{nhp} {\ze_P}= \sT i^*_D\circ
\ze\circ \sT^*(i_{E_{{\mid D}_{M}}})^{-1}\circ\sT^*P: \sT^*D\ra
\sT D^* \ee
is a well-defined map which, in local coordinates,
reads
\be\label{yyy}{\ze_P}(x^i,y^a,p_j,\zx_b)=(x^i,\zx_a,\zr_d^j(x^i,0)y^d,{\mathcal
C}^c_{db}(x^i,0)y^d\zx_c+\zs^k_b(x^i,0)p_k)\,. \ee This is of
course an algebroid structure on the bundle $D$. This algebroid
structure can be described in a more straightforward way as
follows.

The decomposition (\ref{dec}) gives the dual decomposition
\be\label{dec1}E_{\mid D_M}^*=D^*\oplus_{D_M} (D^\perp)^* \ee and
the corresponding projection $P^*:E_{\mid D_M}^*\ra D^*$. This
projection does not give us a canonical projection of vectors and,
more generally, contravariant tensors on $E$ at points of $E_{\mid
D_M}^*$, unless they are tangent to $E_{\mid D_M}^*$. But the
tensor $\zP_\ze$ with the local form (\ref{tensor}) is a sum of
tensor product with at least one part in the product being
vertical. Moreover, due to strong admissibility condition for $D$,
there is a unique decomposition
$\zP_\ze=\zP_{\ze_D}+\zP_{\ze_{D^\perp}}$, where $\zP_{\ze_D}$ is
a linear tensor tangent to $D^*$ and the vertical parts of
$\zP_{\ze_{D^\perp}}$ are tangent to $(D^\perp)^*$, so they will
be killed by the projection, independently how the other part of
the product is. This gives a well-defined projection
$\sT(P^*)(\zP_\ze)=\zP_{\ze_P}$, when we restrict to the points of
$D^*$.

To put it differently, we can take any smooth projection
${\mathcal P}$ from a neighbourhood of $E_{\mid D_M}^*$ onto
$E_{\mid D_M}^*$. Then,
\be\label{gg}\sT(P^*)(\zP_\ze):=\sT(P^*\circ{\mathcal P})(\zP_\ze)
\ee
is a linear 2-contravariant tensor on $D$ which does not depend on
the choice of $\mathcal P$. This tensor defines on $D^*$ the {\it
nonholonomic bracket} $\{\cdot,\cdot\}_{\ze_P}$ associated with
the projection $P$.

On the level of the algebroid bracket $[\cdot,\cdot]_\ze$ this
procedure is the following. The strong admissibility condition for
$D$ implies that the bracket $[X,Y]_\ze$ of sections of $D$ is a
section of $E_{\mid D_M}$. Projecting this section to $D$ along
$P$ gives us a bracket
\be\label{nbra}[X,Y]_{\ze_P}=P[X,Y]_\ze
\ee
on sections of $D$ -- the {\it nonholonomic restriction of
$[\cdot,\cdot]$ along $P$}. This is an algebroid bracket with the
original anchors.

Of course, if the constraint bundle $D$ is a subalgebroid of $E$,
then it is a {\it holonomic} constraint and the projection $P$
plays no role: the tensor $\zP_\ze$ is tangent to $D$ and we just
take the restriction. It is however clear that, in the true
nonholonomic case, the nonholonomic bracket need not satisfy the
Jacobi identity, even when the original bracket does. In local
coordinates,
$$ \zP_{\ze_{P}} =\sum_{a,b,c}{\mathcal C}^c_{ab}(x^i,0)\zx_c \partial _{\zx_a}\otimes \partial _{\zx_b} +
\sum_{j,b}\zr^j_b(x^i,0) \partial _{\zx_b} \otimes \partial _{x^j}
- \sum_{a,j}\zs^j_a(x^i,0)\partial _{x^j} \otimes \partial
_{\zx_a}.
$$
and the corresponding nonholonomic bracket reads:
\begin{eqnarray*}
\{\zx_a, \zx_b\}_{\ze_{P}}&=&{\mathcal C}^c_{ab}(x^i,0)\zx_c, \\
\{\zx_b, x^j\}_{\ze_{P}}&=&\zr^j_b(x^i,0), \\
\{x^j, \zx_a\}_{\ze_{P}}&=&-\zs^j_a(x^i,0),\\
\{x^k, x^j\}_{\ze_{P}}&=& 0.
\end{eqnarray*}
In the case of a Lie algebroid, we have ${\mathcal
C}^c_{ab}=-{\mathcal C}^c_{ba}$, $\sigma^j_A=\rho^j_A$, and the
above bracket corresponds to the one introduced by \cite{SM} (see
also \cite{CoLeMaMa} and references therein) for the subbundles
over the total base $M$.

A natural question arises here: what Lie algebroid structures we
can obtain as nonholonomic restrictions of Lie algebroid brackets?
Of course, as the brackets must be skew-symmetric automatically,
the algebroids must be necessarily quasi-Lie algebroids. We will
call them {\it nonholonomic quasi-Lie algebroids}. In other words:
which linear bi-vector fields on a vector bundle $D \to M$ can be
obtained by projections of Poisson tensors from a bigger vector
bundle $E \to M$. For nonholonomic quasi-Lie algebras the answer
is simple.
\begin{theorem} Any quasi-Lie algebroid $D$ with the trivial anchor
is a nonholonomic quasi-Lie algebroid.
\end{theorem}
\begin{proof} As the anchor is trivial, we can clearly reduce to finite-dimensional
quasi-Lie algebras. Take a basis $e_i$ of $D$ and the
corresponding linear coordinates $y^i$, $i=1,\dots,n$. The
algebroid bracket in $D$ is then determined by the structure
constants $c_{ij}^k$,
$$[e_i,e_j]=c_{ij}^ke_k\,,$$
satisfying $c_{ij}^k=-c_{ji}^k$. Let us consider a new algebra $E$
with the basis $e_i,f_j$, $i=1,\dots,n$ and the bracket
$[\cdot,\cdot]_E$ for which $f_j$ are central elements and
$$[e_i,e_j]_E=c_{ij}^kf_k\,.$$
Since the algebra is 2-step nilpotent, it satisfies the Jacobi
identity and it is a Lie algebra. We can them embed $D$ in $E$ by
putting $\zf(e_i)=(e_i+f_i)$ and take the complementary subspace
$D^\perp$ as spanned by elements $e_i-f_i$. The projection $P:E\ra
D$ is therefore given by
$$P(e_i)=e_i+f_i,\quad P(f_i)=e_i+f_i\,.$$
The nonholonomic bracket on the embedded submanifold $D$ is
therefore
$$[e_i+f_i,e_j+f_j]=P(c_{ij}^kf_k)=c_{ij}^k(e_k+f_k)\,,$$
thus the original bracket in $D$.
\end{proof}

\medskip
For a general quasi-Lie algebroid the situation is much more
complicated and we do not know a full characterization of
nonholonomic quasi-Lie algebroids. Note however that not all
quasi-Lie algebroids are nonholonomic reductions of Lie
algebroids.
\begin{example}{\rm Take a quasi-Lie algebroid structure $\ze$ on
$D=\sT M$ with the anchor $\zr=\id_{\sT M}$. The corresponding
tensor on $\sT^* M$ has therefore the local form
$$\zP_\ze=\frac{1}{2}c_{ij}^k(x)p_k\partial_{p_i}\wedge\partial_{p_j}
+\partial_{p_i}\wedge\partial_{x^i}\,.
$$
Suppose that the corresponding bracket is a nonholonomic reduction
of a Lie algebroid bracket $[\cdot,\cdot]_E$ of a bigger vector
bundle $E$ along a projection $P$,
$$[X,Y]_\ze=P([X,Y]_E)\,.
$$
Since the anchor $\zr$ of $D$ is the restriction of the anchor
$\zr_E$ to $D$, the Lie algebroid $E$ satisfies $\zr_E(E_{\mid
M})=\sT M$, thus $E_{\mid M}=D\oplus_M K$, where $K$ is the kernel
of $\zr_E$ over $M$. For sections $X,Y$ of $D$ we have therefore
$\zr_E([X,Y]_E)=\zr([X,Y])$, thus, as the anchor map is for Lie
algebroids a homomorphism of brackets,
$$\zr([X,Y])=\zr_E([X,Y]_E)=[\zr_E(X),\zr_E(Y)]_{vf}
=[\zr(X),\zr(Y)]_{vf}\,,
$$
where $[\cdot,\cdot]_{vf}$ is the bracket of vector fields. We get
that the anchor $\zr$ maps the algebroid bracket into the bracket
of vector fields which is not the case for generally non-zero
structure functions $c_{ij}^k(x)$. }
\end{example}

\section{Nonholonomic systems of mechanical type}
\noindent Let us consider now on our general algebroid $(E,\ze)$ a
Lagrangian $L: E\rightarrow \R$ of mechanical type, i.e. \be\label
{ml} L(e)=\frac{1}{2}{G} (e, e)-V(\tau(e))\,, \ee where ${\mathcal
G}:E\times_M E\to \R$ is a bundle metric on $E$ and $V:
M\rightarrow \R$ is the potential function. Let us consider also a
a vector subbundle $D$ of $E$ over $D_M$ satisfying the strong
admissibility condition. Having the metric ${\mathcal G}$ in $E$,
we have the natural decomposition $E_{\mid
D_M}=D\oplus_{D_M}D^\perp$ with $D^\perp$ being the orthogonal
complement of $D$ with respect to the metric ${\mathcal G}$,
accompanied with the corresponding (this time -- orthogonal)
projection $P:E_{\mid D_M}\ra D$. The fundamental observation in
this case is the following.
\begin{theorem}\label{Theorem1}The nonholonomic Tulczyjew
differential ${\zL^D_L}$ coincides with the Tulczyjew differential
$\zL_l$ associated with the restricted Lagrangian $l=L_{\mid D}:
D\rightarrow \R$ and the nonholonomic algebroid structure $\ze_P$
on $D$. In other words, the nonholonomic Euler-Lagrange equations
associated with $L$ and the nonholonomic constraint $D$ in $E$
coincide with the non-constrained Euler-Lagrange equations for the
Lagrangian function $l$ on the algebroid $(D,{\ze_P})$.
\end{theorem}
\begin{proof}
Let us consider local coordinates $(x^I)=(x^i,x^\zi)$ on a open
set ${U}$ of $M$ such that $D_M$ is determined by the constraint
$x^\zi=0$. Take a local basis of orthonormal sections $\{e_a,
e_{\alpha}\}$ of $E$ adapted to the orthogonal decomposition
$E_{\mid D_M}=D\oplus_{\mid D_M} D^{\perp}$, i.e, $\hbox{span}
\{e_a\}=D$ and $\hbox{span} \{e_\alpha\}=D^\perp$. We get the
induced coordinates $(x^I,y^A)$ in $E$ in which the Lagrangian
takes the form
$$L(x,y)=\frac{1}{2}\sum _A(y^A)^2-V(x^i,x^\zi)\,.$$
Note that the `mass' is 1 in these coordinates. On $E_{\mid D_M}$
we have the induced coordinates $(x^i, y^a, y^{\alpha})$, so $D$
is locally determined by the constraints $x^\zi=0,y^{\alpha}=0$
and the restricted Lagrangian is:
\[
l(x^i, y^a)=\frac{1}{2}\sum_{a=1}^{r} (y^a)^2-V(x^i,0).
\]
In our coordinates, the projector ${P}$ reads
\[
{P}(x^i, y^a, y^{\alpha})=(x^i, y^a)
\]
and we deduce that $${\ze_P}(x^i,
y^a,p_j,\zx_b)=(x^i,\zx^a,\zr^j_d(x^i,0)y^d,{\mathcal
C}^c_{db}(x^i,0)y^d\zx_c+\zs^k_b(x^i,0)p_k)\,,
$$
so that
$$\ze_{P}\circ \xd l(x^i, y^a)
=\left(x^i, y^a, \rho^j_d(x^i,0)y^d, \sum_{c}{{\mathcal
C}}^c_{db}(x^i,0)y^d y^c-\sigma^k_b\frac{\partial V}{\partial
x^k}(x^i,0)\right)\,.
$$
We get therefore the Euler-Lagrange equations for the Lagrangian
$l$ on $D$ of the form
\begin{eqnarray*}
  \frac{d x^j}{dt}&=&\rho^j_a(x^i,0)y^a,\\
 \frac{d}{dt}\left(\frac{\partial l}{\partial
y^b}(x^i,y^a)\right)&=& {{\mathcal C}}^c_{db}(x^i,0) y^d
y^c+\sigma^j_b(x)\frac{\partial l}{\partial x^j}(x^i,y^a),
\end{eqnarray*}
which are exactly the nonholonomic equations (\ref{nELe1}) and
(\ref{nELe2}), if we take to account that
$$\frac{\partial L}{\partial y^C}(x^i,0, y^a, 0)=\zd^C_ay^a$$
and
$$\sigma^\zg_b(x^i,0)=0\,.$$
\end{proof}
The above theorem makes it clear that, at least for Lagrangians of
mechanical type, the nonholonomic Euler-Lagrange equations are
just Euler-Lagrange equations, but on the nonholonomic restriction
of the algebroid. In any case, passing to a liner nonholonomic
constraint does not move us out of the Mechanics on algebroids.
Our theory is therefore complete with respect to passing to the
nonholonomic case, what was the main problem in understanding the
nonholonomic constraint in the Lie algebroid (thus canonical)
case.

\medskip
Having interpreted the nonholonomic Euler-Lagrange equations as
just Euler-Lagrange equations, but on a reduced algebroid, we can
use the generalization of Noether Theorem formulated for general
algebroids in \cite{GGU} to obtain a nonholonomic Noether Theorem,
at least for Lagrangians of mechanical type.

The generalized Noether Theorem for a general algebroid $(E,\ze)$
is based on the concept of the complete lift $\dt^\ze(K)$ of
tensor fields $K$ being sections of the tensor products
$E^{\otimes k}$ to the corresponding contravariant tensor fields
on $E$. This is a natural generalization of the standard concept
of the tangent lift $\dt$ which lifts contravariant tensors on a
manifold $M$ to the corresponding tensor fields on $\sT M$, (cf.
\cite{GU1,YI}).

For a vector bundle $\zt:E\ra M$, let $\otimes ^k(\zt)$ be the
space of sections of the tensor-product bundle $E^{\otimes k}$
over $M$. With any tensor field $K\in \otimes ^k(\zt)$ we can
associate the linear function $\zi(K)$ on the dual bundle
$(E^{\otimes k})^*=(E^*)^{\otimes k}$ and the vertical lift
$\textsf{v}_\zt(K) \in \otimes ^k(\ezt_E)$  In local coordinates,
$$ \zi(f^{a_1\dotsi a_k}(x) e_{a_1}\otimes \cdots \otimes
e_{a_k}) = f^{a_1\dotsi a_k}(x)\zx_{a_1\dotsi a_k}
$$
and
$$ \textsf{v}_\zt(f^{a_1\dotsi a_k}(x) e_{a_1}\otimes \cdots \otimes
e_{a_k}) = f^{a_1\dotsi a_k}(x) \partial _{y^{a_1}} \otimes \cdots
\otimes \partial _{y^{a_k}}\,.
$$
A particular case of the vertical
lift is the lift $\vt(K)$ of a contravariant tensor field $K$ on
$M$ into a contravariant tensor field on $\sT M$. It is well known
(see \cite{YI,GU1}) that in the case of $E=\sT M$ we have also the
tangent lift $\dt \colon \otimes (\zt_M) \rightarrow \otimes
(\zt_{\ss T M})$ which is a $\vt$-derivation. What has been done
for the tangent bundle, can be repeated in the case of an
arbitrary algebroid $(E,\ze)$. Note first that we can extend $\ze$
naturally to mappings (cf. \cite{GU3, GU2})
$$ \ze^{\otimes r} \colon \otimes ^r_E \sT^\* E \longrightarrow
\sT \otimes ^r_M E^\*,\quad r\ge 0.
$$

\begin{theorem}{\em \cite{GU3, GU2}}
Let $(E,\ze)$ be an algebroid. For $K\in \otimes ^k(\zt)$, $k\ge
0$, the equality
\be \zi(\dt^\ze(K)) = \dt (\zi(K))\circ \ze^{\otimes k} \ee
defines the  tensor field $\dt^\ze(K) \in \otimes ^k(\zt_E)$ which
is linear and the mapping
$$ \dt^\ze \colon \oplus_k\otimes^k (\zt) \longrightarrow \oplus_k\otimes^k (\zt_E) $$
is a $\vt$-derivation of degree $0$. In local coordinates, the
lifts of functions on $M$ and sections of $E$ read:
\bea
\dt^\ze(f(x))& =&  y^a \zr^i_a(x) \frac{\partial f}{\partial
x^i}(x),\\
\dt^\ze(f^a(x)e_a) &=& f^a(x)\zs^i_a(x)\partial _{x^i} + \left(
y^a \zr^i_a(x) \frac{\partial f^c}{\partial x^i}(x) + {\mathcal
C}^c_{ab}(x)y^af^b(x) \right) \partial _{y^c}\,.
\eea
Conversely, if $D\colon \oplus_k\otimes^k (\zt) \longrightarrow
\oplus_k\otimes^k (\zt_E)$ is a $\vt$-derivation of degree $0$
such that $D(K)$ is linear for each $K \in \otimes ^1(\zt)$, then
there is an algebroid structure $\ze$ on $\zt\colon E\rightarrow
M$ such that $D = \dt^\ze$. This algebroid structure is a Lie
algebroid if and only if
$$\dt^\ze([X,Y]_\ze) = [\dt^\ze(X), \dt^\ze(Y)]_{vf}
$$
for all $X,Y\in\otimes ^1(\zt)$.
\end{theorem}
Suppose now that we are dealing with a Lagrangian $L$ of
mechanical type and a nonholonomic constraint subbundle $D$ of
$(E,\ze)$ satisfying the strong admissibility assumption, so that
we have also the nonoholomic reduced algebroid $(D,\ze_P)$. We
will say that a pair $(X,f)$ consisting of a section $X$ of $D$
and a function $f$ on $D_M$ is a {\it symmetry of the nonholonomic
problem} associated with the Lagrangian $L$, if
\be\label{Ne}
{\dt^{\ze_P}(X)}(l)=\dt^{\ze_{P}}(f)\,,
\ee
where $l=L_{\mid D}$ is the restriction of the Lagrangian to the
constraint.
\begin{theorem}{[Nonholonomic Noether Theorem]}
The following are equivalent:
\begin{itemize}
\item[(a)] the pair $(X,f)$ is a symmetry of the nonholonomic
problem associated with the Lagrangian $l$;

\item[(b)] the function $\zi_{X}-f\circ\zt_{D^*}$ on $D^*$ is a
constant of the motion of the nonholonomic phase dynamics on
$D^*$, i.e. $\dt(\zi_{X}-f\circ\zt_{D^*})$ vanishes on
${\zG^D_L}$;

\item[(c)] the function $$(\zi_{X}-f\circ\zt_{D^*})\circ\zl_L^D
$$
is a constant of the motion for the nonholonomic Euler-Lagrange
equation.
\end{itemize}
\end{theorem}
\begin{proof}
The proof is completely analogous to that for  Theorem 4 in
\cite{GGU}.
\end{proof}
For a section $X$ of $D$,  not verifying necessarily condition (\ref{Ne}), we obtain
\be\label{nme}
{\dt^{\ze_P}(X)}(l)\circ\zg(t)=\frac{\xd}{\xd t}\left(
\zi_{X}\circ\zl_L^D\circ\zg(t)\right)
\ee
for any solution $\zg:\R\ra D$ of the nonholonomic Euler-Lagrange
equation. This last equation can be interpreted as a general
version of the {\it nonholonomic momentum equation} studied for
several authors (see \cite{CoLeMaMa,BKMM} and references therein)
for nonholonomic systems with symmetry.

\section{Examples}

\begin{example}{The  Chaplygin sleigh.}

{\rm As an example of nonholonomic system on a Lie algebra, we
study   the  Chaplygin sleigh which describes a rigid body sliding
on a plane. The body is supported in three points, two of which
slides freely without friction while the third point is a knife
edge.  This imposes the constraint of no motion orthogonal to this
edge (see \cite{Cha,NeFu}).

The configuration space before reduction is the Lie group
$G=SE(2)$ of the Euclidean motions of the 2-dimensional plane
$\R^2$. We will need in the sequel to fix some notation about the
Lie algebra $\mathfrak{se}(2)$. First of all its elements are
matrices of the form
\[
\hat{\xi}=\left(
\begin{array}{ccc}
0&\xi_3&\xi_1\\
-\xi_3&0&\xi_2\\
0&0&0
\end{array}
\right)
\]
and a basis of the Lie algebra $\mathfrak{se}(2)\cong \R ^3$
is given by
\[
E_1=\left(
\begin{array}{ccc}
0&0&1\\
0&0&0\\
0&0&0
\end{array}
\right) ,\qquad E_2=\left(
\begin{array}{ccc}
0&0&0\\
0&0&1\\
0&0&0
\end{array}
\right) ,\quad
E_3=\left(
\begin{array}{ccc}
0&-1&0\\
1&0&0\\
0&0&0
\end{array}
\right).
\]
We have that
\[
{}[E_3,E_1]=E_2,\quad [E_2, E_3]=E_1,\quad [E_1, E_2]=0.
\]
An element $\xi \in \mathfrak{se}(2)$ is of the form
\[
\xi =v_1\, E_1+v_2\, E_2+\omega \, E_3.
\]

The Chaplygin system is described by the kinetic Lagrangian function
\[
\begin{array}{rrcl}
L:&\mathfrak{se}(2)&\longrightarrow&\R\\
  &(v_1, v_2, \omega)&\longmapsto&\frac{1}{2}\left[ (J+m(a^2+b^2))\omega^2 + mv_1^2+m v_2^2-2bm\omega v_1-2am\omega v_2\right]
  \end{array}
\]
where   $m$ and $J$ denotes the mass and moment of inertia of the
sleigh relative to the contact point and $(a, b)$ represents the
position of the center of mass with respect to the body frame
determined placing the origin at the contact point and the first
coordinate axis in the direction of the knife axis. Additionally,
the   system is subjected to the nonholonomic constraint
determined by the linear subspace of $\mathfrak{se}(2)$:
\[
D=\{(v_1, v_2, \omega)\in se(2)\; |\; v_2=0\}\, .
\]
Instead of $\{E_1, E_2, E_3\}$ we take the basis of $\mathfrak{se}(2)$:
\[
\{e_1=E_3, e_2=E_1, e_3= -ma E_3-mab E_1+(J+ma^2) E_2\}
\]
which is  a basis adapted to the decomposition $D\oplus D^\perp$; $D=\hbox{span }\{ e_1, e_2\}$ and $D^\perp=\hbox{span }\{ e_3\}$.

In the induced coordinates $(y^1, y^2)$ on $D$ the restricted lagrangian is
\[
l(y^1, y^2)=\frac{1}{2}\left[ (J+m(a^2+b^2))(y^1)^2 +
m(y^2)^2-2bmy^1y^2\right]\; ,
\]
and moreover,
\[
{}[e_1, e_2]_{\ze_P}=\frac{ma}{J+ma^2} e_1+\frac{mab}{J+ma^2}e_2\;
,\] Therefore, ${\mathcal C}^1_{12}=\frac{ma}{J+ma^2}$ and
${\mathcal C}^2_{12}=\frac{mab}{J+ma^2}$.

Then,
 $$\ze_P\circ \xd l(y^1, y^2))
=\left( (J+m(a^2+b^2))y^1 -bm y^2, my^2-bmy^1, -ma y^1y^2,
ma(y^1)^2\right).
$$
In consequence, the equations of motion are:
\begin{eqnarray*}
(J+m(a^2+b^2))\dot{y}^1 -bm \dot{y}^2&=& -ma y^1y^2\\
m\dot{y}^2-bm\dot{y}^1&=&ma(y^1)^2
\end{eqnarray*}
or,
\begin{eqnarray*}
\dot{y}^1 &=& \frac{ma}{J+ma^2}y^1\left(by^1-y^2\right)\\
\dot{y}^2&=&
\frac{ma}{J+ma^2}y^1\left((J+m(a^2+b^2))y^1-by^2\right).
\end{eqnarray*}

}

\end{example}

\begin{example}{The snakeboard}
{\rm

 As a mechanical system the snakeboard has as configuration
space $Q=SE(2)\times T^2$ with coordinates $(x, y, \theta,
\psi, \phi)$ (see ~\cite{BuZe,LMMdD,MaKo}).

\medskip
\begin{center}
\includegraphics{snakeboard.pdf}
\end{center}

\medskip
The  nonholonomic dynamics is described by
\begin{itemize}
\item The \emph{Lagrangian}
\begin{eqnarray*}
 L(q, \dot{q})= \frac{1}{2}m (\dot{x}^2+\dot{y}^2) + \frac{1}{2}
 (J+2J_1)\dot{\theta}^2
 +\frac{1}{2}J_0(\dot{\theta}+\dot{\psi})^2
 +J_1\dot{\phi}^2,
 \end{eqnarray*}
where $m$ is the total mass of the board, $J>0$ is the moment of
inertia of the board, $J_0>0$ is the moment of inertia of the
rotor of the snakeboard mounted on the body's center of mass and
$J_1>0$  is the moment of inertia of each wheel axles. The
distance between the center of the board and the wheels is denoted
by $r$. For simplicity (see  \cite{MaKo}), we assume that $J+J_0+2J_1=mr^2$.

\item The \emph{nonholonomic constraints} induced by the non sliding condition in the
sideways direction of the wheels:
\begin{align*}
-\dot{x}\sin (\theta+\phi) +\dot{y}\cos (\theta+\phi)
-r\dot{\theta}\cos\phi&=0\\
-\dot{x}\sin (\theta-\phi) +\dot{y}\cos (\theta-\phi)
+r\dot{\theta}\cos\phi&=0.
\end{align*}
\end{itemize}

Observe that the Lagrangian is induced by the riemannian  metric ${\mathcal G}$ on $Q$,
\[
{\mathcal G}=mdx^2+mdy^2+mr^2d\theta^2+J_0 d\theta\otimes\psi+J_0 d\psi\otimes
d\theta + J_0d\psi^2+2J_1 d\phi^2.
\]

 The constraint  subbundle $\tau_D: D\longmapsto Q$ is
 \[
{D}=\text{span}\left\{ e_1=\frac{\partial}{\partial
\psi}, e_2=\frac{\partial}{\partial \phi},
e_3=a\frac{\partial}{\partial x}+b\frac{\partial}{\partial
y}+c\frac{\partial}{\partial \theta}\right\}.
\]
where
\begin{eqnarray*}
a&=&-r(\cos\phi\cos(\theta-\phi)+\cos\phi\cos(\theta+\phi))=-2r\cos^2\phi\cos\theta\\
b&=&-r(\cos\phi\sin(\theta-\phi)+\cos\phi\sin(\theta+\phi))=-2r\cos^2\phi\sin\theta\\
c&=&\sin(2\phi).
\end{eqnarray*}

The orthogonal complement of $D$ is spanned by
\[
{D}^\perp=\text{span}\left\{ e_4=-b\frac{\partial}{\partial x}+a\frac{\partial}{\partial y},
e_5=(cJ-cmr^2)\frac{\partial}{\partial x}
+am\frac{\partial}{\partial \theta}-am\frac{\partial}{\partial \psi}\right\}.
\]
In the induced coordinates $(x, y, \theta,
\psi, \phi, y^1, y^2, y^3)$ on $D$  the restricted lagrangian is
\[
l((x, y, \theta,
\psi, \phi, y^1, y^2, y^3)=2mr^2\cos^2\phi (y^3)^2+J_0 c y^1y^3+\frac{1}{2}J_0 (y^1)^2+J_1 (y^2)^2.
\]
where now the nonholonomic constraints are rewritten as: $y^4=0$ and $y^5=0$.
After some straightforward computations we deduce that
\begin{eqnarray*}
{}[e_1, e_2]_{\ze_P}&=&0,
\\
{}[e_1, e_3]_{\ze_P}&=&0,\\
{}[e_2, e_3]_{\ze_P}&=&
\frac{2mr^2\cos^2\phi}{mr^2-J_0\sin^2\phi} e_1
-\frac{(mr^2+\cos 2\phi)\tan\phi}{mr^2-J_0\sin^2\phi} e_3.
\end{eqnarray*}
Thus, the unique non vanishing structure functions are:
\[
{\mathcal C}_{23}^1=-{\mathcal C}_{32}^1=\frac{2mr^2\cos^2\phi}{mr^2-J_0\sin^2\phi}, \qquad
{\mathcal C}_{23}^3=-{\mathcal C}_{32}^3=-\frac{(mr^2+\cos 2\phi)\tan\phi}{mr^2-J_0\sin^2\phi}  \; .
\]

Therefore, the equation of motion of the snakeboard are
\begin{eqnarray*}
&\dot{x}=ay^3,\qquad \dot{y}=by^3, \qquad \dot{\theta}=cy^3,\qquad \dot{\psi}=y^1, \qquad \dot{\phi}=y^2,&\\
&\frac{d}{dt}(y^1+cy^3)=0,\qquad \frac{d}{dt}(y^2)=0,&\\
&\frac{d}{dt}(4mr^2\cos^2\phi y^3+J_0 cy^1)= 2J_0\cos(2\phi)y^1y^2-2mr^2\sin(2\phi)y^2y^3&\; .
\end{eqnarray*}

}
\end{example}