The anomalous behavior of coefficient of normal restitution in the oblique impact
The coefficient of normal restitution in an oblique impact is %numerically and theoretically studied. Using a two-dimensional lattice model for an elastic disk and an elastic wall, we demonstrate %the impact of the elastic disk on the elastic wall %at various incident angles. %In the condition that the translational energy of the disk %is fixed, we measure the coefficient of normal restitution for %each incident angle. %We find that the coefficient of normal restitution can exceed unity and has a peak against the incident angle %between an elastic disk and an elastic wall in our simulation. %In addition, we find the coefficient of friction also has a peak %against the incident angle. We also explain this behavior based upon a phenomenological theory.
main.tex
% Group addresses by affiliation; use superscriptaddress for long
% author lists, or if there are many overlapping affiliations.
% For Phys. Rev. appearance, change preprint to twocolumn.
% Choose pra, prb, prc, prd, pre, prl, prstab, or rmp for journal
%  Add 'draft' option to mark overfull boxes with black boxes
%  Add 'showpacs' option to make PACS codes appear
%  Add 'showkeys' option to make keywords appear
%\documentclass[aps,prl,preprint,superscriptaddress,showpacs]{revtex4}
\documentclass[aps,prl,twocolumn,superscriptaddress]{revtex4}
\usepackage{graphicx}
%\usepackage{amsmath}
% You should use BibTeX and apsrev.bst for references
% Choosing a journal automatically selects the correct APS
% BibTeX style file (bst file), so only uncomment the line
% below if necessary.
\bibliographystyle{apsrev}

\begin{document}

% Use the \preprint command to place your local institutional report
% number in the upper righthand corner of the title page in preprint mode.
% Multiple \preprint commands are allowed.
% Use the 'preprintnumbers' class option to override journal defaults
% to display numbers if necessary
%\preprint{}

%Title of paper
\title{The anomalous behavior of coefficient of normal restitution in 
the oblique impact}

% repeat the \author .. \affiliation  etc. as needed
% \email, \thanks, \homepage, \altaffiliation all apply to the current
% author. Explanatory text should go in the []'s, actual e-mail
% address or url should go in the {}'s for \email and \homepage.
% Please use the appropriate macro foreach each type of information

% \affiliation command applies to all authors since the last
% \affiliation command. The \affiliation command should follow the
% other information
% \affiliation can be followed by \email, \homepage, \thanks as well.
\author{Hiroto Kuninaka}
\email[E-mail: ]{kuninaka@kuchem.kyoto-u.ac.jp}
%\homepage[URL: ]{http://ace.phys.h.kyoto-u.ac.jp/~kuninaka/}
\thanks{Present address: Department of Chemistry, Kyoto University, 
Sakyo-ku, Kyoto, Japan, 606-8502}
%\altaffiliation{}
\affiliation{Graduate School of Human and Environmental Studies, 
Kyoto University, Sakyo-ku, Kyoto, Japan, 606-8501}
\author{Hisao Hayakawa}
\affiliation{Department of Physics, Yoshida-south campus, 
Kyoto University, Sakyo-ku, Kyoto, Japan, 606-8501
 }
%Collaboration name if desired (requires use of superscriptaddress
%option in \documentclass). \noaffiliation is required (may also be
%used with the \author command).
%\collaboration can be followed by \email, \homepage, \thanks as well.
%\collaboration{Hisao Hayakawa}
%\noaffiliation

\date{\today}

\begin{abstract}
The coefficient of normal restitution in an oblique impact is
%numerically and 
theoretically studied. 
Using a two-dimensional lattice model for an elastic disk 
and an elastic wall, 
we demonstrate
%the impact of the elastic disk on the elastic wall 
%at various incident angles.
%In the condition that the translational energy of the disk
%is fixed, we measure the coefficient of normal restitution for
%each incident angle.
%We find 
that the coefficient of normal restitution can exceed unity 
and has a peak against the incident angle 
%between an elastic disk and an elastic wall 
in our simulation. 
%In addition, we find the coefficient of friction also has a peak 
%against the incident angle. 
We also explain this behavior based upon 
a phenomenological theory.
\end{abstract}

% insert suggested PACS numbers in braces on next line
%\pacs{PACS number(s): 45.10.-b, 45.50.Tn, 45.70.-n, 62.20.-x, 81.40.Pq}
\pacs{45.50.-j, 45.50.Tn, 45.70.-n, 62.20.-x, 81.40.Pq}
%PACS number(s): 45.10.-b, 45.50.Tn, 45.70.-n, 62.20.-x, 81.40.Pq
% insert suggested keywords - APS authors don't need to do this
%\keywords{}

\maketitle

The coefficient of normal restitution $e$ 
is introduced to determine the normal component 
of the post-collisional velocity in the collision 
of two materials. 
The coefficient $e$ is defined by 
\begin{equation}\label{COR}
{\bf v}(\tau_{c}) \cdot {\bf n}=-e {\bf v}(0) \cdot {\bf n},
\end{equation} 
where ${\bf v}(\tau)$ is the relative velocity of 
the centers of mass of two colliding materials 
at time $\tau$ measured from the initial contact, 
$\tau_{c}$ is the duration of a collision, 
and ${\bf n}$ is the unit vector normal to the contact plane. 
Though some text books of elementary physics state 
that $e$ is a material constant, 
many experiments and simulations show that $e$ decreases 
with increasing impact velocity\cite{stronge}. 
The dependence of $e$ on the low impact velocity is theoretically 
treated by the quasi-static theory
\cite{kuwabara,morgado,brilliantov96}. 
We also recognize that $e$ can be less than unity 
for the normal impact without the introduction of any explicit 
dissipation, because the macroscopic inelasticities 
originate in the transfer of the energy 
from the translational mode to the internal modes 
such as the vibrations\cite{morgado,gerl,ces}. 
 
While $e$ has been believed to be less than unity 
in most situations, 
it is recently reported that $e$ can exceed unity 
in oblique impacts\cite{smith,calsamiglia,louge}. 
In particular, Louge and Adams\cite{louge} observed 
oblique impacts of a hard aluminum oxide sphere 
on a thick elastoplastic polycarbonate plate 
%with the incident angle $\gamma$. 
and found that $e$ grows monotonically with the magnitude of 
% (Fig.\ref{fig1}).
the tangent of the incident angle $\gamma$. 
In their experiment, 
Young's modulus of the plate is $100$ times smaller 
than that of the aluminum oxide sphere. 
They also suggested that $e$ can exceed unity 
for the most oblique impacts. 
 
In this letter, 
we demonstrate 
that our two-dimensional simulation of the oblique impact 
based on Hamilton's equation 
has yielded an increasing $e$ with $\tan\gamma$ and 
$e$ exceeds unity at the critical incident angle. 
Finally, we explain our results 
%with the aid of 
by our phenomenological theory. 
%\cite{galin}.

%%%%%%%%%%%% FIGURES %%%%%%%%%%%%
\begin{figure}[htbp]
\begin{center}
\includegraphics[width=0.5\textwidth]{fig1.eps}
\end{center}
\caption{The elastic disk and wall consisted of random lattice.}
\label{fig1}
\end{figure}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Let us introduce our numerical model\cite{kuninaka_jpsj2003}. 
Our model consists of an elastic disk and an elastic wall 
(Fig. \ref{fig1}). 
The width and the height of the wall are $8 R$ and $2 R$, 
respectively, where $R$ is the mean radius of the undeformed disk. 
%We adopt the fixed boundary condition 
%for the both side ends and the bottom of the wall. 
The both side ends and the bottom of the wall are fixed. 
We place $800$ mass points at random in a disk with the radius $R$ 
and $4000$ mass points at random in a wall 
for the disk and the wall, respectively.  
%the position 
%$(x_j, y_j)(1 \le j \le 2000)$, 
%where $x_j$ and $y_j$ are uniform random numbers 
%in the range of $-4R \le x_j \le 4R$ and $-2R \le y_j \le 0$ 
%respectively. 
%After that, 
We connect each mass point with its neighbor mass points 
by the Delaunay triangulation algorithm\cite{delaunay}, 
and undeformed nonlinear springs are placed on all the connections.
%all the connections are replaced 
%by undeformed nonlinear springs.
%We place an undeformed nonlinear spring for each connected line. 

%The spring interaction between connected mass points 
%is given by 
%\begin{equation}\label{interaction}
% V^{(i)}(x)=
%\frac{1}{2} k^{(i)}_{a} x^{2}+\frac{1}{4} k^{(i)}_{b} x^{4},
%\hspace{3mm}i=d\text{(disk)},w\text{(wall)},
%\end{equation} 

Each mass point $i$ on the lower half boundary of the disk 
feels the force, 
${\bf F}(l_{s}^{(i)})=aV_0\exp(-a l_{s}^{(i)}){\bf n}_{s}^{(i)}$, 
where $l_{s}^{(i)}$ is the distance 
between $i$-th surface mass point of the disk and 
the nearest surface spring of the wall, $a=300/R$, $V_0=amc^2R/2$, 
$m$ is the mass of each mass point, $c$ is the one-dimensional 
velocity of sound, 
and ${\bf n}_{s}^{(i)}$ is the unit vector normal to the connection 
between two surface mass points of the wall
\cite{kuninaka_jpsj2003}. 
%Then, two particles of both ends of the surface spring are applied by 
%the reaction forces which are decided by the balance of the torques 
%in the opposite direction. 
We should note that the strong repulsion ${\bf F}(l_{s}^{(i)})$ 
is introduced to inhibit the penetration of the disk 
to the surface of the wall\cite{gerl}. 
Thus, the dynamical equation of motion for each mass point $i$ 
of the lower half boundary of the disk is described by
\begin{equation}\label{keq}
m \frac{d^{2}{\bf r}_{i}}{d \tau^{2}}=
\sum^{N_{i}}_{j=1} \left\{-k_{a} {\bf x}_{ij} - k_{b} {\bf x}_{ij}^{3}
\right\}+\Theta(l_{th}-l_{s}^{(i)}) a V_0 \exp(-a l_{s}^{(i)})
{\bf n}_{s}^{(i)},
\end{equation}
where ${\bf r}_{i}$ is the position of $i$-th mass point, 
$N_{i}$ is the number of mass points connected to $i$-th mass point, 
${\bf x}_{ij}$ is the relative deformation of the spring 
between $i$-th and $j$-th connected mass points, 
$k_{a}$ and $k_{b}=k_{a}\times10^{-3}/R^{2}$ are the spring constants. 
Here we introduce the step function $\Theta(x)$, 
i.e. $\Theta(x)=1$ for $x \ge 0$ and $\Theta(x)=0$ for $x < 0$, 
and the threshold length $l_{th}$ which is the average of 
the natural lengths of the springs of the disk. 
%$\theta(x)=1 \hspace{2mm} \text{for} \hspace{2mm} x \geq 0
%\hspace{2mm} \text{and } 
%\hspace{2mm} 0\hspace{2mm}\text{for} \hspace{2mm} x<0$, 
%and $l_{th}$ is the threshold length which is assumed to 
%be the average of the natural lengths of the springs of the disk.   
For internal mass points, the last term of 
the right hand side of eq.(\ref{keq}) is omitted. 
In most of our simulations, 
we adopt $k_{a} = k^{(d)}_{a} =1.0 \times m c^2/R^2$ 
for the disk and 
$k_{a} = k^{(w)}_{a} =1.0 \times 10^{-2} m c^2/R^2$ for the wall. 
%where $m$ and $c$ are the mass of each mass point and 
%the one-dimensional velocity of sound, respectively. 
We do not introduce any dissipative mechanism in this model. 
Thus, during a collision, 
a part of initial translational energy of the disk is 
distributed into the vibrational energy of the disk and the wall. 
It should be noted that the macroscopic dissipation 
can be interpreted as the irreversible transfer of the energy 
from the translational motion to the internal vibration. 
When we introduce explicit dissipations in the model 
and add the gravity to the disk, we have confirmed 
that the compression of the disk can be described 
by two-dimensional Hertzian contact theory
\cite{kuninaka_doctor}.  

In this model, the roughness of the surfaces is important 
to make the disk rotate after a collision\cite{kuninaka_jpsj2003}. 
We modify the positions of the surface mass points of the flat 
wall and the smooth disk by using normal random numbers 
whose average and standard deviation are 
$0$ and $3 \times 10^{-2}R$, respectively. 

Poisson's ratio $\nu$ and Young's modulus $E$ of this model 
can be evaluated by adding the external force to 
stretch the rectangle of random lattice numerically. 
%For the relaxation of the vibrational motion of the band, 
%we add the viscous force term in the spring force 
%which is proportional to the relative velocity 
%of two connected mass points.
We obtain Poisson's ratio $\nu=(7.50 \pm 0.11) \times 10^{-2}$ 
and Young's modulus $E=(9.54 \pm 0.231)\times 10^{3}mc^{2}/R^{2}$, 
respectively\cite{kuninaka_jpsj2003}. 
%by pulling the band mechanically with $10$ magnitudes 
%of pulling force from $2.0 \times 10^{2} mc^{2}/R$ 
%to $3.0 \times 10^{2} mc^{2}/R$ 
%and averaging those $10$ results. 

%In our simulation, we define the incident angle $\gamma$ 
%as the angle between the normal vector of the wall and the initial 
%velocity vector of the disk(see Fig. \ref{fig1}). 
We solve the dynamical equation of motion (\ref{keq}) 
for each mass point with the initial speed $|{\bf v}(0)|=0.1 c$ 
and the incident angle $\gamma$,  
and determine $e$ for each $\gamma$ according to eq.(\ref{COR}). 
All the results in this letter are obtained by averaging 
the results of 100 disks with different configurations of 
mass points. 
We use the fourth order symplectic integrator 
with the time step $\Delta \tau = 10^{-3}R/c$. 

%%%%%%%%%%%% FIGURES %%%%%%%%%%%%
\begin{figure}[htbp]
\begin{center}
\includegraphics[width=0.4\textwidth]{fig2.eps}
\end{center}
\caption{Numerical and theoretical results of the relation between
 $\Psi_1$ and $e$.}
\label{fig2}
\end{figure}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Figure \ref{fig2} is the relation between $e$ and 
$\Psi_1 \equiv -({\bf v}(0)\cdot{\bf t})/({\bf v}(0)\cdot{\bf n})
=\tan\gamma$, 
where ${\bf t}$ is the unit vector vertical to ${\bf n}$. 
%tangential to the contact plane. 
The cross point is the mean value and the error bar is  
the standard deviation of 100 samples for each $\gamma$. 
This result shows that $e$ increases with increasing $\Psi_1$ 
to exceed unity, and has a peak around $\Psi_1 = 6.0$. 
The behavior of $e$ having the peak is contrast to 
that in the experiment by Louge and Adams\cite{louge}. 

Here, let us explain 
%clarify the mechanism of 
our results. 
Louge and Adams\cite{louge} suggest 
that their results can be explained by the rotation of 
the normal unit vector ${\bf n}$ arising from 
the local deformation of the wall's surface. 
Thus, we aim to determine the angle of rotation of the unit vector 
$\alpha$ at each $\gamma$ from the theory of elasticity. 

%%%%%%%%%%%% FIGURES %%%%%%%%%%%%
\begin{figure}[htbp]
\begin{center}
\includegraphics[width=0.4\textwidth]{fig3.eps}
\end{center}
\caption{The schematic figure of a hard disk sliding on a soft wall. 
$x$ coordinates of both ends of the contact area AB are $x=x_a$ and 
$x=x_b$.}
\label{fig3}
\end{figure}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Figure \ref{fig3} is the schematic figure of a hard disk moving 
from the left to the right on a wall,  
%$x$ coordinates of both ends of the contact area AB are $x=x_a$ and 
%$x=x_b$. 
where $l\equiv|x_b-x_a|$. 
Assuming that $l$ is small compared to $R$, 
the perimeter of the contact is approximated by 
\begin{equation}\label{contact}
 f(x)=(x-x_c)^2/2R-y_c,
\end{equation}
where $(x_c,y_c)$ is the lowest position of the disk. 
To calculate $\tan\alpha \equiv (f(x_b)-f(x_a))/l$, 
we need to know the ratio of $|x_c-x_a|$ to $l$.  
%because the left and the right sides of the contact area 
%are not equal due to the frictional force. 
From the theory of elasticity\cite{galin}, 
this ratio can be estimated as 
\begin{equation}\label{ratio1}
 \frac{x_c-x_a}{l} = 1-\theta \hspace{5mm} \text{with} \hspace{5mm}
 \theta = \frac{1}{\pi}\arctan \frac{1-2\nu}{\mu (2-2\nu)},
\end{equation}
where $\mu$ is the coefficient of the friction. 
We evaluate $\mu$ 
by $\mu \equiv |{\bf J} \cdot {\bf t}|/|{\bf J} \cdot {\bf n}|$ 
with ${\bf J}=M({\bf v}(\tau)-{\bf v}(0))$, where $M$ is the mass 
of the disk. 
The cross points in Fig.\ref{fig4} represent $\mu$ calculated 
from our simulation, where $\mu$ has a peak 
around $\Psi_1=3.0$. 
%Note that Eq.(\ref{ratio1}) equals to $1/2$ when $\mu=0$.
From eqs.(\ref{contact}) and (\ref{ratio1}),  
$\tan\alpha$ is given by 
\begin{equation}\label{tanalpha2}
\tan \alpha = \frac{2\theta-1}{2-2\theta} \frac{|x_c-x_a|}{R},
\end{equation}
where we evaluate $|x_c-x_a|=0.55R$ 
from the maximum value of the compression of the disk. 
% through its definition 
%$\mu = |J_t|/|J_n|$, where $J_{t}$ and $J_{n}$ are respectively 
%tangential and normal components of impulse ${\bf J}$. 
From eqs.(\ref{ratio1}) and (\ref{tanalpha2}), 
we obtain the relation between $\Psi_1$ and $\tan\alpha$.

Next, let us calculate $e$ 
from the relation between $\Psi_1$ and $\tan\alpha$. 
We introduce the rotated unit vectors, 
${\bf n}_{\alpha}$ and ${\bf t}_{\alpha}$, as 
${\bf n}_{\alpha} = \cos\alpha {\bf n}-\sin\alpha {\bf t}$ and 
${\bf t}_{\alpha} = \sin\alpha {\bf n}+\cos\alpha {\bf t}$, 
respectively. 
%The impulse ${\bf J}$ is represented as
%\begin{equation}
%{\bf J} = J_n {\bf n}+J_t {\bf t}
%= J_n^{\alpha} {\bf n^{\alpha}}+J_t^{\alpha} {\bf t^{\alpha}}.
%\end{equation}
%Thus, $J_n$ and $J_t$ are expressed as 
%${\bf n}^{\alpha}=\cos\alpha {\bf n}-\sin\alpha {\bf t}$ and  
%${\bf t}^{\alpha}=\sin\alpha {\bf n}+\cos\alpha {\bf t}$, 
%because ${\bf n}=(0,1)$ and  ${\bf t}=(1,0)$, 
%we obtain
%\begin{eqnarray}
%$J_n = J^{\alpha}_{n}\cos\alpha+J^{\alpha}_{t}\sin\alpha$ and 
%$J_t = -J^{\alpha}_{n}\sin\alpha+J^{\alpha}_{t}\cos\alpha$.
%\end{eqnarray}
%By expressing the velocities of the contact point of the disk 
%before and after collision as ${\bf u}$ and ${\bf u}^{'}$ respectively , 
%\begin{eqnarray}
%${\bf u} = u_n {\bf n}+u_t {\bf t}
%= u_n^{\alpha} {\bf n}^{\alpha}+u_t^{\alpha}{\bf t}^{\alpha}$ and
%${\bf u}^{'} = u_n^{'} {\bf n}+u_t^{'} {\bf t}
%= u_n^{\alpha'} {\bf n}^{\alpha}+u_t^{\alpha'}{\bf t}^{\alpha}$.
%\end{eqnarray}
%Note that $u_n=u_n^{\alpha}\cos\alpha+u_t^{\alpha}\sin\alpha$ and 
%$u_t=-u_n^{\alpha}\sin\alpha+u_t^{\alpha}\cos\alpha$.
%From these expressions, 
%$u_n$ and $u_n^{'}$ are represented as 
%$u_n=u_n^{\alpha}\cos\alpha+u_t^{\alpha}\sin\alpha$ and 
%$u_n^{'}=u_n^{\alpha'}\cos\alpha+u_t^{\alpha'}\sin\alpha$.
%Here, 
By introducing 
$e_{\alpha} \equiv -({\bf v}(\tau_{c}) \cdot {\bf n}_{\alpha})/({\bf v}(0) \cdot
{\bf n}_{\alpha})$,  
% as COR for the collision of an inclined surface. 
we can express $e$ in terms of $e_{\alpha}$ as 
\begin{equation}\label{CORmod2}
e=\frac{e_{\alpha}+\Psi_2^{\alpha} \tan\alpha}{1-\Psi_1^{\alpha} \tan\alpha},
\end{equation}
where 
$\Psi_1^{\alpha}=-({\bf v}(0) \cdot {\bf t}_{\alpha})/({\bf v}(0) \cdot 
{\bf n}_{\alpha})$ 
and $\Psi_2^{\alpha}=-({\bf v}(\tau_{c}) \cdot {\bf t}_{\alpha})/({\bf v}(0) \cdot 
{\bf n}_{\alpha})$. 
%$e^{\alpha}$ is interpreted a constant to $\Psi_1$ less than $1$. 
%$\Psi_1^{\alpha}$ and $\Psi_2^{\alpha}$ in Eq.(\ref{CORmod2}) can 
%be expressed by $\mu$, $\Psi_1$, and $\tan\alpha$. 
%By noting that 
%$u_t^{\alpha}=u_n\sin\alpha+u_t\cos\alpha$ and
%$u_n^{\alpha}=u_n\cos\alpha-u_t\sin\alpha$, 
$\Psi_1^{\alpha}$ also can be rewritten as 
\begin{equation}\label{psi1}
\Psi_{1}^{\alpha}=(\Psi_1-\tan\alpha)/(1+\Psi_1 \tan\alpha).
\end{equation}
%from the definitions of ${\bf n}_{\alpha}$ and ${\bf t}_{\alpha}$. 
On the other hand, in the oblique impact of slipping disks, 
$\Psi_2^{\alpha}$ is given by  
\begin{equation}\label{from_wal}
\Psi_2^{\alpha}=\Psi_1^{\alpha}-3(1+e_{\alpha})\mu_{\alpha}
\end{equation}
in the two-dimensional situation\cite{walton_rheol}.  
%$\mu^{\alpha}$ in Eq.(\ref{from_wal}) is a function of $\mu$ 
%and $\tan\alpha$.
In eq.(\ref{from_wal}), $\mu_{\alpha}$, defined by 
$\mu_{\alpha}=|{\bf J} \cdot {\bf t}_{\alpha}|/|{\bf J} \cdot {\bf
n}_{\alpha}|$, is given by 
\begin{equation}\label{mualpha}
\mu_{\alpha}=\frac{\mu+\tan\alpha}{1-\mu\tan\alpha}.
\end{equation}
%Here, $\mu_{\alpha}$ is the coefficient of sliding friction 
%without deformation of the wall. 
%From eqs.(\ref{CORmod2})-(\ref{mualpha}), we can determine 
%$e$ as a function of $\alpha$, $\mu$, $e_{\alpha}$, and $\Psi_{1}$.  
To draw the solid line in Fig.\ref{fig2}, at first, 
we calculate $\mu$ and $\tan\alpha$ for each $\Psi_1$. 
By choosing a fitting parameter $e_{\alpha}=0.95$ and 
substituting eqs.(\ref{psi1})-(\ref{mualpha}) into eq.(\ref{CORmod2}) 
we obtain $e$ as a function of $\alpha$, $\mu$, and $\Psi_{1}$. 
%With the aid of  $\Psi_1^{\alpha}$ and $\Psi_2^{\alpha}$ 
%in eqs.(\ref{psi1}) and (\ref{from_wal}), we obtain $e$ 
%by substituting them into eq.(\ref{CORmod2}) for each $\Psi_1$, 
%where $e_{\alpha}=0.95$ is a fitting parameter. 
%The solid line of Fig.\ref{fig2} is Eq.(\ref{CORmod2}) with 
%$e_{\alpha}=0.95$ which is the value in the normal impact. 
%and is assumed to be constant in all range of $\Psi_1$. 
All points are interpolated with the cubic spline interpolation method 
to draw the solid curve. 
Such the theoretical description of $e$ 
is qualitatively consistent with our numerical result, 
though the theoretical value is a little smaller 
than the observed value. 
%As can be seen, COR $e$ depends on the relation between 
%$\mu$ and $\Psi_1$.  

%%%%%%%%%%%% FIGURES %%%%%%%%%%%%
\begin{figure}[htbp]
\begin{center}
\includegraphics[width=0.4\textwidth]{fig4.eps}
\end{center}
\caption{Numerical and theoretical results of the relation between
 $\Psi_1$ and $\mu$.}
\label{fig4}
\end{figure}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


Now, let us consider how $\mu$ depends on $\Psi_1$ 
based on a phenomenological argument. 
For simplicity, we replace the roughness on the surface 
of the wall by a periodic array of asperities.
When the disk hits one asperity,  
 a fraction of the energy is absorbed in the wall. We assume that 
the tangential velocity $v_t^{(i-1)}$ changes to
$v_t^{(i)}=(1-\eta)v_t^{(i-1)}$ when the disk hits the $i$-th asperity.
Thus, if the disk contacts $N$ asperities during the impact, the
tangential speed at the release point 
becomes $v_t(\tau_{c})=v_{t}^{(N)}=(1-\eta)^Nv_t(0)$.
Here we estimate the number of contacted asperities during the
collision as $N=\rho l_{sl}$, where $\rho$ is the number of 
the asperities 
in a unit length on the surface and $l_{sl}$
is the length of sliding  which can be evaluated as 
$l_{sl}=v_t(0)\tau_{c}$ with $\tau_{c}=\pi (R/c)
\sqrt{\ln\left(4 c/v_n(0)\right)}$\cite{ces}.
Thus, the tangential impulse $J_t\equiv {\bf J}\cdot {\bf t}$ 
is approximated by
\begin{equation}
J_t=M\{(1-\eta)^{\rho l_{sl}}-1\}v_t(0)\simeq -\eta M \rho l_{sl} 
v_t(0),
\end{equation}
where we assume small $\eta$. 
Now, we should answer the question how to determine $\eta$.
When the impact velocity is large enough, $\eta$ becomes smaller 
because the asperities are broken down when the disk hits them.
Therefore we may assume the form 
$\eta=\eta_0/(1+\beta v^{2}_{t}(0)/v^{2}_{n}(0))
=\eta_0/(1+\beta \Psi^{2}_{1})$ with 
the introduction of the dimensionless parameters $\eta_0$ and $\beta$.
Taking into account 
$J_n \equiv {\bf J}\cdot{\bf n}=-M(1+e)v_n(0)$ and the
definition of $\mu=J_t/J_n$ we obtain
\begin{equation}
\mu=\frac{\pi\eta_{0} \rho R}{1+e}
\frac{\Psi_{1}}{1+\beta \Psi^{2}_{1}}
\sqrt{\frac{\Psi^{2}_{1}}{1+\Psi^{2}_{1}}
\ln\left(40\sqrt{1+\Psi^{2}_{1}}\right)} \frac{|{\bf v}(0)|}{c}.
\label{toy}
\end{equation}
Here we use $\cos\gamma=\sqrt{1/(1+\Psi^{2}_{1})}$,
$\sin\gamma=\sqrt{\Psi^{2}_{1}/(1+\Psi^{2}_{1})}$. 
%In this expression, 
%$e$ is determined by the simulation data. 
We adopt the numerical results for the value of $e$ 
at each $\Psi_1$ (Fig.\ref{fig2}) and $|{\bf v}(0)|/c = 0.1$. 
%$\beta$ and $\eta_{0} \rho$ are fitting parameters. 
%Assuming that the disk rotates scarcely after collision,  
%we calculate Eq.(\ref{toy}) by the contact point velocity ${\bf u}$ 
%instead of ${\bf v}$. 
The solid curve in Fig.\ref{fig4} is eq.(\ref{toy}) 
with fitting parameters $\beta=0.21$ and $\eta_0 \rho R=0.18$, 
which reproduces our numerical result. 
We do not claim that our simple argument explains 
the experimental result because of two fitting parameters. 
However, we would like to emphasize that our picture 
captures the essence of the behavior of $\mu$ in the oblique impact.

Let us discuss our result. 
First, we emphasize that the novel phenomena of $e$ exceeding unity 
% observed by Louge and 
%Adams 
are obtained from the local deformation of the soft wall 
for the oblique impacts of a hard disk. 
When we simulate the impact between a disk and a hard wall,  
%refer to the origin of the relation between $e$ and $\Psi_1$. 
%As described previously, for $e$ to exceed $1$, 
%it is essential for the wall surface to deform locally, which 
%produces the rotation of the normal unit vector to the wall surface. 
%If $\alpha=0$ in Eq.(\ref{CORmod2}), 
%$e$ is equal to $e^{\alpha}$ which is a constant smaller than $1$. 
%We have also carried out the simulation when 
for $k^{(w)}_a = 10 \times k^{(d)}_a$, 
$e$ fluctuates around a constant to exceed unity abruptly 
around $\Psi_1=4.5$\cite{kuninaka_doctor}. 
This tendency resembles the experimental results 
by Calsamiglia {\it et al.}\cite{calsamiglia}. 
%This result can be attributed to the insufficient local deformation 
%originated in the hardness of the wall. 
Thus, for smooth increase of $e$ to exceed unity, 
the wall should be softer than the disk. 
%In fact, the Young's modulus of the wall is $100$ times softer 
%than that of the sphere in the experiment by Louge and Adams\cite{louge}.
In addition, 
it is important to fix the initial kinetic energy of the disk. 
So far, we have confirmed that $e$ does not exceed unity when $\Psi_1$ is 
controlled by changing $v_t$ with fixed
$v_n$\cite{kuninaka_jpsj2003}. 
%the initial kinetic energy of the disk becomes large 
%for the large $\Psi_1$ so that the local deformation of the wall is 
%easy to collapse in the grazing collision. 

The second, the initial velocity of the disk and the local deformation 
of the wall are much larger than those in the experimental ones 
in ref.\cite{louge}. 
% that the local dissipation in springs and 
%gravity have not affected our numerical results. 
%Our simple estimation of $\alpha$ have made the solid curve 
%in Fig.\ref{fig2} less than simulation data globally.
%And it became clear that the solid curve can fit our data more 
%correctly when $|x_c-x_a|$ in Eq.(\ref{tanalpha2}) takes 
%larger value. 
%We may need to estimate the size of the contact area more precisely.
They cause the most significant difference 
between our result and their result\cite{louge}.  
Because of the high speed impact in our case, 
there is a peak of $e$ for small $\gamma$. 
In fact, our simulation with 
$|{\bf v}(0)|=0.01c$ shows the shift of the peak for larger $\gamma$. 
Therefore, we expect that our model reproduces the result of 
ref.\cite{louge} for the low impact speed. 
In addition,  we have carried out simulations 
with a disk of $400$ mass points and a wall of $2000$ mass points 
to check the effect of the system size. 
Although there is a slight difference between the results, 
the data are also well reproducible by our phenomenological theory.  
%from the relation between $\mu$ and $\Psi_1$.  

The third, the local deformation of the wall also affects 
the relation between $\mu$ and $\Psi_1$. 
In early studies, 
it has been shown that $\mu$ depends on the impact velocity
\cite{gorham,louge}. 
%They also have determined the relation between $\mu$ and $\Psi_1$ 
%and shown that the relation can be fitted by an exponential function.
In our simulation, $\mu$ has a peak around $\Psi_1=3.0$. 
This behavior is interpreted as 
that the asperities are flattened for large $v_t$. 
%This results in the peak of $\mu$ to $\Psi_1$. 
%We derived the relation between $\mu$ and $\Psi_1$ 
%in eq.(\ref{toy}). 
%By controlling $\beta$ and $\eta_0 \rho$, 
%$\mu$ can increase with the increasing $\gamma$ 
%as in ref.\cite{louge}. 
Equation (\ref{toy}) indicates that $\mu$ can increase 
with increasing $\gamma$ as in ref.\cite{louge} 
if we choose a suitable set of $\beta$ and $\eta_{0} \rho$. 
The difference between the results of their experiment and our
simulation may be explained by the choice of these parameters.

The fourth, we adopt the static theory of elasticity to explain our
numerical results in this letter. 
However, it is important to solve the time-dependent equation 
of the deformation of the wall's surface 
to analyze the dynamics of impact phenomena. 
The dynamical analysis will be our future task.

In final, we indicate that the friction coefficient is derived from
our Hamiltonian model. This friction comes from the irreversible 
energy transfer of the macroscopic translational motion 
to the internal vibration. The irreversibility is indeed related 
to the second law of thermodynamics 
if the number of the internal degrees of freedom is infinite. However, 
the irreversibility in our system which includes only $10^3$ mass points
in a disk is incomplete as we can see that
the Hertzian contact theory cannot be recovered 
without introduction of the explicit dissipation. 
We believe that the complete treatment of inelastic collision
of macroscopic materials without introduction of dissipations will
be a fundamental subject of nonequilibrium statistical mechanics.

In summary, we have carried out the two-dimensional simulation 
of the oblique impact of an elastic disk on an elastic wall. 
We have found that  $e$ can exceed unity in the oblique impact, 
which is attributed to the local deformation of the wall. 
We have estimated the magnitude of the local defomation $\alpha$ 
based on the static theory of elasticity and derived the relation 
between $e$ and $\Psi_1$ 
by taking into account the rotation of the normal unit vector of 
the wall's surface. The relation between $\mu$ and $\Psi_1$ is also 
related to the local deformation and is explained 
by the simple argument. 
%assuming the decrease of the tangential impulse proportional to 
%the initial tangential velocity of the disk.

\begin{acknowledgments}
We would like to thank Y. Tanaka and S. Nagahiro for their 
valuable comments. 
We would also like to thank M. Y. Louge for his detailed 
introduction of their experiment to us and some valuable comments. 
Parts of numerical computation in this work were carried out 
at Yukawa Institute Computer Facility. 
This study is partially supported by the Grant-in-Aid of 
Ministry of Education, Culture, Sports, Science and Technology (MEXT), 
Japan (Grant No. 15540393) and the Grant-in-Aid for the 21st century COE 
"Center for Diversity and Universality in Physics" from MEXT, Japan. 
\end{acknowledgments}

%\bibliographystyle{apsrev}
\bibliography{impact}
%\newpage
%%%%%%%%%%%% FIGURES %%%%%%%%%%%%
%\begin{figure}[htbp]
%\begin{center}
%\includegraphics[width=1.\textwidth]{fig1.eps}
%\end{center}
%\caption{The elastic disk and wall consisted of random lattice.}
%\label{fig1}
%\end{figure}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\newpage
%%%%%%%%%%%% FIGURES %%%%%%%%%%%%
%\begin{figure}[htbp]
%\begin{center}
%\includegraphics[width=1.\textwidth]{fig2.eps}
%\end{center}
%\caption{Numerical and theoretical results of the relation between
% $\Psi_1$ and $e$.}
%\label{fig2}
%\end{figure}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\newpage
%%%%%%%%%%%% FIGURES %%%%%%%%%%%%
%\begin{figure}[htbp]
%\begin{center}
%\includegraphics[width=1.\textwidth]{fig3.eps}
%\end{center}
%\caption{The schematic figure of a hard disk sliding on a soft wall. 
%$x$ coordinates of both ends of the contact area AB are $x=x_a$ and 
%$x=x_b$.}
%\label{fig3}
%\end{figure}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\newpage
%%%%%%%%%%%% FIGURES %%%%%%%%%%%%
%\begin{figure}[htbp]
%\begin{center}
%\includegraphics[width=1.\textwidth]{fig4.eps}
%\end{center}
%\caption{Numerical and theoretical results of the relation between
% $\Psi_1$ and $\mu$.}
%\label{fig4}
%\end{figure}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\end{document}


main.bbl
\begin{thebibliography}{15}
\expandafter\ifx\csname natexlab\endcsname\relax\def\natexlab#1{#1}\fi
\expandafter\ifx\csname bibnamefont\endcsname\relax
  \def\bibnamefont#1{#1}\fi
\expandafter\ifx\csname bibfnamefont\endcsname\relax
  \def\bibfnamefont#1{#1}\fi
\expandafter\ifx\csname citenamefont\endcsname\relax
  \def\citenamefont#1{#1}\fi
\expandafter\ifx\csname url\endcsname\relax
  \def\url#1{\texttt{#1}}\fi
\expandafter\ifx\csname urlprefix\endcsname\relax\def\urlprefix{URL }\fi
\providecommand{\bibinfo}[2]{#2}
\providecommand{\eprint}[2][]{\url{#2}}

\bibitem[{\citenamefont{Stronge}(2000)}]{stronge}
\bibinfo{author}{\bibfnamefont{W.~J.} \bibnamefont{Stronge}},
  \emph{\bibinfo{title}{Impact Mechanics}} (\bibinfo{publisher}{Cambridge
  University Press}, \bibinfo{address}{Cambridge}, \bibinfo{year}{2000}).

\bibitem[{\citenamefont{Kuwabara and Kono}(1987)}]{kuwabara}
\bibinfo{author}{\bibfnamefont{G.}~\bibnamefont{Kuwabara}} \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{K.}~\bibnamefont{Kono}},
  \bibinfo{journal}{Jpn. J. Appl. Phys.} \textbf{\bibinfo{volume}{26}},
  \bibinfo{pages}{1230} (\bibinfo{year}{1987}).

\bibitem[{\citenamefont{Morgado and Oppenheim}(1997)}]{morgado}
\bibinfo{author}{\bibfnamefont{W.~A.~M.} \bibnamefont{Morgado}}
  \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{I.}~\bibnamefont{Oppenheim}},
  \bibinfo{journal}{Phys.\ Rev.\ E} \textbf{\bibinfo{volume}{55}},
  \bibinfo{pages}{1940} (\bibinfo{year}{1997}).

\bibitem[{\citenamefont{Brilliantov et~al.}(1996)\citenamefont{Brilliantov,
  Spahn, Hertzsch, and P{\"o}schel}}]{brilliantov96}
\bibinfo{author}{\bibfnamefont{N.~V.} \bibnamefont{Brilliantov}},
  \bibinfo{author}{\bibfnamefont{F.}~\bibnamefont{Spahn}},
  \bibinfo{author}{\bibfnamefont{J.-M.} \bibnamefont{Hertzsch}},
  \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{T.}~\bibnamefont{P{\"o}schel}},
  \bibinfo{journal}{Phys.\ Rev.\ E} \textbf{\bibinfo{volume}{53}},
  \bibinfo{pages}{5382} (\bibinfo{year}{1996}).

\bibitem[{\citenamefont{Gerl and Zippelius}(1999)}]{gerl}
\bibinfo{author}{\bibfnamefont{F.}~\bibnamefont{Gerl}} \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{A.}~\bibnamefont{Zippelius}},
  \bibinfo{journal}{Phys.\ Rev.\ E} \textbf{\bibinfo{volume}{59}},
  \bibinfo{pages}{2361} (\bibinfo{year}{1999}).

\bibitem[{\citenamefont{Hayakawa and Kuninaka}(2002)}]{ces}
\bibinfo{author}{\bibfnamefont{H.}~\bibnamefont{Hayakawa}} \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{H.}~\bibnamefont{Kuninaka}},
  \bibinfo{journal}{Chem. Eng. Sci.} \textbf{\bibinfo{volume}{57}},
  \bibinfo{pages}{239} (\bibinfo{year}{2002}).

\bibitem[{\citenamefont{Smith and Liu}(1992)}]{smith}
\bibinfo{author}{\bibfnamefont{C.~E.} \bibnamefont{Smith}} \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{P.-P.} \bibnamefont{Liu}},
  \bibinfo{journal}{J. of Appl. Mech.} \textbf{\bibinfo{volume}{59}},
  \bibinfo{pages}{963} (\bibinfo{year}{1992}).

\bibitem[{\citenamefont{Calsamiglia et~al.}(1997)\citenamefont{Calsamiglia,
  Kennedy, Chatterjee, Ruina, and Jenkins}}]{calsamiglia}
\bibinfo{author}{\bibfnamefont{J.}~\bibnamefont{Calsamiglia}},
  \bibinfo{author}{\bibfnamefont{S.~W.} \bibnamefont{Kennedy}},
  \bibinfo{author}{\bibfnamefont{A.}~\bibnamefont{Chatterjee}},
  \bibinfo{author}{\bibfnamefont{A.}~\bibnamefont{Ruina}}, \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{J.~T.} \bibnamefont{Jenkins}},
  \bibinfo{journal}{J. of Appl. Mech.} \textbf{\bibinfo{volume}{66}},
  \bibinfo{pages}{146} (\bibinfo{year}{1997}).

\bibitem[{\citenamefont{Louge and Adams}(2002)}]{louge}
\bibinfo{author}{\bibfnamefont{M.~Y.} \bibnamefont{Louge}} \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{M.~E.} \bibnamefont{Adams}},
  \bibinfo{journal}{Phys.\ Rev.\ E} \textbf{\bibinfo{volume}{65}},
  \bibinfo{pages}{021303} (\bibinfo{year}{2002}).

\bibitem[{\citenamefont{Kuninaka and Hayakawa}(2003)}]{kuninaka_jpsj2003}
\bibinfo{author}{\bibfnamefont{H.}~\bibnamefont{Kuninaka}} \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{H.}~\bibnamefont{Hayakawa}},
  \bibinfo{journal}{J. Phys. Soc. Jpn.} \textbf{\bibinfo{volume}{72}},
  \bibinfo{pages}{1655} (\bibinfo{year}{2003}).

\bibitem[{\citenamefont{Sugihara}(2001)}]{delaunay}
\bibinfo{author}{\bibfnamefont{K.}~\bibnamefont{Sugihara}},
  \emph{\bibinfo{title}{Data Structure and Algorithms(in Japanese)}}
  (\bibinfo{publisher}{Kyoritsu}, \bibinfo{address}{Japan},
  \bibinfo{year}{2001}).

\bibitem[{\citenamefont{Kuninaka}(2004)}]{kuninaka_doctor}
\bibinfo{author}{\bibfnamefont{H.}~\bibnamefont{Kuninaka}}, Ph.D. thesis,
  \bibinfo{school}{Kyoto University} (\bibinfo{year}{2004}); 
\bibinfo{author}{\bibfnamefont{H.}~\bibnamefont{Hayakawa}} \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{H.}~\bibnamefont{Kuninaka}},
  \bibinfo{journal}{Phase Trans.} \textbf{\bibinfo{volume}{77}},
  \bibinfo{pages}{889} (\bibinfo{year}{2004}).

\bibitem[{\citenamefont{Galin}(1953)}]{galin}
\bibinfo{author}{\bibfnamefont{L.~A.} \bibnamefont{Galin}},
  \emph{\bibinfo{title}{Contact problems in the theory of elasticity}}
  (\bibinfo{publisher}{Gostekhizdat}, \bibinfo{address}{Moskow},
  \bibinfo{year}{1953}).

\bibitem[{\citenamefont{Walton and Braun}(1986)}]{walton_rheol}
\bibinfo{author}{\bibfnamefont{O.~R.} \bibnamefont{Walton}} \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{R.~L.} \bibnamefont{Braun}},
  \bibinfo{journal}{J. Rheol.} \textbf{\bibinfo{volume}{30}},
  \bibinfo{pages}{949} (\bibinfo{year}{1986}).

\bibitem[{\citenamefont{Gorham and Kharaz}(2000)}]{gorham}
\bibinfo{author}{\bibfnamefont{D.~A.} \bibnamefont{Gorham}} \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{A.~H.} \bibnamefont{Kharaz}},
  \bibinfo{journal}{Powder Technology} \textbf{\bibinfo{volume}{112}},
  \bibinfo{pages}{193} (\bibinfo{year}{2000}).

\end{thebibliography}
