---
abstract: 'A three-dimensional continuum dislocation theory for single crystals containing curved dislocations is proposed. A set of governing equations and boundary conditions is derived for the true placement, plastic slips, and loop functions in equilibrium that minimize the free energy of crystal among all admissible functions, provided the resistance to the dislocation motion is negligible. For the non-vanishing resistance to dislocation motion the governing equations are derived from the variational equation that includes the dissipation function. A simplified theory for small strains is also provided. An asymptotic solution is found for the two-dimensional problem of a single crystal beam deforming in single slip and simple shear.'
address: |
    Lehrstuhl für Mechanik - Materialtheorie, Ruhr-Universität Bochum,\
    D-44780 Bochum, Germany
author:
- 'K. C. Le[^1]'
title: 'Three-dimensional continuum dislocation theory'
---

dislocations (A) ,crystal plasticity (B) ,finite strain (B) ,variational calculus (C) .

Introduction {#sec:introduction}
============

In view of a huge number of dislocations appearing in plastically deformed crystals (which typically lies in the range $10^8\div 10^{15}$ dislocations per square meter) the necessity of developing a physically meaningful continuum dislocation theory (CDT) to describe the evolution of dislocation network and predict the formation of microstructure in terms of mechanical and thermal loading conditions becomes clear to all researchers in crystal plasticity. One of the main guiding principles in seeking such a continuum dislocation theory has first been proposed by @Hansen1986 in form of the so-called LEDS-hypothesis: the true dislocation structure in the final state of deformation minimizes the energy of crystal among all admissible dislocation configurations. In view of numerous experimental evidences supporting this hypothesis (see, i.e., [@Hughes1997; @Kuhlmann1989; @Kuhlmann2001; @Laird1986]), its use in constructing the continuum dislocation theory seems to be quite reasonable and appealing. For the practical realization one needs to i) specify the whole set of unknown functions and state variables of the continuum dislocation theory, and ii) lay down the free energy of crystals as their functional to be minimized. Such program has been implemented by @Berdichevsky06a in the linear, and by @Le2014nonlinear in the nonlinear setting of CDT for networks of dislocations, whose lines are straight and remain so during the whole deformation process (see also [@Ortiz99; @Ortiz00]). The developed CDT has been successfully applied to various two-dimensional problems of dislocation pileups, bending, torsion, as well as formation of dislocation patterns in single crystals (see [@Berdichevsky-Le07; @Kaluza2011torsion; @Kochmann08b; @Kochmann08a; @Kochmann09; @Koster2015; @Le08a; @Le08b; @Le2009plane; @Le2012polygonization; @Le2013on]). Let us mention here the similar approaches suggested in [@Acharya2000; @Acharya2001; @Engels2012; @Gurtin2002; @Gurtin2007; @Mayeur2014; @Oztop2013] which do not use the LEDS-hypothesis explicitly but employ instead the extended principle of virtual work for the gradient plasticity. However, as experiences and experiments show, dislocation lines are in general loops that, as a rule, can change their directions and curvatures depending on the loading condition and crystal’s geometry. Therefore the extension of CDT to networks of dislocations whose lines are curves in the slip planes is inevitable. To the best of author’s knowledge, such three-dimensional continuum dislocation theory based on the LEDS-hypothesis for curved dislocations has not been developed until now. It became also clear to him that the latter’s absence was due to the missing scalar dislocation densities for the network of curved dislocations.

The first attempt at constructing a continuum theory that can predict in principle not only the dislocation densities but also the direction and curvature of the dislocation lines has been made by @Hochrainer2007 in form of the so-called continuum dislocation dynamics. Their theory starts with the definition of the dislocation density that contains also the information about the orientation and curvature of the dislocation lines. Then the set of kinematic equations is derived for the dislocation density and curvature that requires the knowledge about the dislocation velocity. The relation between the dislocation density and the macroscopic plastic slip rate via the dislocation velocity is postulated in form of Orowan’s equation. The couple system of crystal plasticity and continuum dislocation dynamics becomes closed by the constitutive equation of a flow rule type (see [@Hochrainer2014; @Sandfeld2011; @Sandfeld2015; @Wulfinghoff2015]). In addition to the heavy computational cost of such theory, the relation to thermodynamics of crystal plasticity and to the LEDS-hypothesis is completely lost: the equilibrium solution found in this theory may not minimize the energy of crystal among all admissible dislocation configurations. Let us mention also a continuum approach proposed recently by @Zhu2013 [@Zhu2014] in which the three-dimensional dislocation structure is characterized by two families of disregistry functions that may take only integer values. The dislocation density can then be expressed in their terms. The coupled system of equations is derived from the underlying discrete dislocation dynamics for the displacement and disregistry functions. This approach is subject to the same critics as that proposed in [@Hochrainer2007].

The aim of this paper is to extend the nonlinear continuum dislocation theory (CDT) developed recently by @Le2014nonlinear to the case of crystals containing curved dislocations. Provided the dislocation network is regular in the sense that nearby dislocations have nearly the same direction and orientation, we introduce a loop function whose level curves coincide with the dislocation lines. Taking an infinitesimal area perpendicular to the dislocation line at some point of the crystal, we express the densities of edge and screw dislocations at that point through the resultant Burgers vectors of dislocations whose lines cross this area at right angle. Such scalar densities contain not only the information about the number of dislocations, but also the information about the orientation and curvature of the dislocation lines. In case of dislocation motion we introduce the vector of normal velocity of dislocation line through the time derivative of the loop function. Following @Kroener92 and [@Berdichevsky06] we require that the free energy density of crystal depends only on the elastic strain tensor and on the above scalar densities of dislocations. Then we formulate a new variational principle of CDT according to which the placement, the plastic slip, and the loop function in the final state of equilibrium minimize the free energy functional among all admissible functions. We derive from this variational principle a new set of equilibrium equations, boundary conditions, and constitutive equations for these unknown functions. In case the resistance to dislocation motion is significant, the variational principle must be replaced by the variational equation that takes the dissipation into account. The constructed theory is generalized for single crystals having a finite number of active slip systems. We provide also the simplifications of the theory for small strains. As compared to the continuum dislocation dynamics proposed in [@Hochrainer2007; @Zhu2014] our theory is advantageous not only in the computational cost due to its simplicity, but also in its full consistency with the LEDS-hypothesis. In the problem of single crystal beam having only one active slip system and deforming in simple shear, the energy minimization problem reduces to the two-dimensional variational problem. We solve this problem analytically for the circular cross section and asymptotically for the rectangular cross section. We will show that this solution reduces to that found in [@Berdichevsky-Le07] for the crystals with thin and long cross-section.

The paper is organized as follows. After this short introduction we present in Section 2 the three-dimensional kinematics for single crystals deforming in single slip. Section 3 formulates the variational principles of the three-dimensional CDT and derives its governing equations. Section 4 extends this nonlinear theory to the case of single crystals with $n$ active slip system. Section 5 studies the three-dimensional small strain CDT. Section 6 is devoted to the analytical and asymptotic solutions of the two-dimensional energy minimization problem of a single crystal beam deforming in simple shear. Finally, Section 7 concludes the paper.

3-D kinematics for single crystals deforming in single slip {#sec:ncdt}
===========================================================

Nonlinear CDT starts from the basic kinematic resolution of the deformation gradient $\mathbf{F}=\partial \mathbf{y}/\partial \mathbf{x}$ into elastic and plastic parts [@Bilby57] $$\label{eq:2.1}
    \mathbf{F}=\mathbf{F}^e\cdot \mathbf{F}^p.$$ We attribute an active role to the plastic deformation: $\mathbf{F}^p$ is the deformation [*creating*]{} dislocations (either inside or at the boundary of the volume element) or [*changing*]{} their positions in the crystal without distorting the lattice parallelism (see Fig. \[fig:Nhslip\]). On the contrary, the elastic deformation $\mathbf{F}^e$ deforms the crystal lattice having [*frozen*]{} dislocations [@Le2014nonlinear]. Note that the lattice vectors remain unchanged when the plastic deformation is applied, while they change together with the shape vectors by the elastic deformation.

![Multiplicative decomposition[]{data-label="fig:Nhslip"}](Fp){height="6.5cm"}

We consider first a single crystal deforming in single slip. In this case let us denote the right-handed triad of unit lattice vectors of the active slip system by $\mathbf{s}$, $\mathbf{p}$, and $\mathbf{m}$, where $\mathbf{s}$ points to the slip direction, $\mathbf{p}$ lies in the slip plane and is perpendicular to $\mathbf{s}$, and $\mathbf{m}$ is normal to the slip plane. Without restricting generality we may choose the rectangular cartesian coordinate system $(x_1,x_2,x_3)$ in the reference configuration such that its basis vectors coincide with these lattice vectors (see Fig. \[fig:loop\]) $$\mathbf{e}_1=\mathbf{s},\quad \mathbf{e}_2=\mathbf{p}, \quad \mathbf{e}_3=\mathbf{m}.$$ The plastic deformation is then given by $$\label{eq:2.0}
\mathbf{F}^p=\mathbf{I}+\beta (\mathbf{x})\mathbf{s}\otimes \mathbf{m}=\mathbf{I}+\beta (\mathbf{x})\mathbf{e}_1\otimes \mathbf{e}_3,$$ with $\beta $ being the plastic slip. We assume that all dislocations causing this plastic deformation lie completely in the slip planes and the dislocation network is regular in the sense that nearby dislocations have nearly the same direction and orientation. This enables one to introduce a scalar function $l(x_1,x_2,x_3)$ (called a loop function) such that its level curves $$\label{eq:2.a}
l(x_1,x_2,c_3)=c,$$ with $c_3$ and $c$ being constants, coincide with the dislocation lines. Thus, in this three-dimensional kinematics we admit, according to equation , only the conservative motion of dislocations and exclude from consideration the dislocation climb which is an important mechanism of temperature-dependent creep. We denote by ${\boldsymbol{\nu}}$ and ${\boldsymbol{\tau}}$ the plane unit vectors normal and tangential to the dislocation line. From equation follow $${\boldsymbol{\nu}}=\frac{1}{\sqrt{l_{,1}^2+l_{,2}^2}}(l_{,1}\mathbf{e}_1+l_{,2}\mathbf{e}_2), \quad {\boldsymbol{\tau}}=\frac{1}{\sqrt{l_{,1}^2+l_{,2}^2}}(-l_{,2}\mathbf{e}_1+l_{,1}\mathbf{e}_2),$$ where the comma before an index denotes the partial derivative with respect to the corresponding coordinate. Note that ${\boldsymbol{\nu}}$, ${\boldsymbol{\tau}}$, $\mathbf{m}$ form a right-handed basis vectors of the three-dimensional space (see Fig. \[fig:loop\]).

![A dislocation loop in the chosen coordinate system[]{data-label="fig:loop"}](loop){height="6.5cm"}

@Ortiz99 introduced the resultant Burgers vector of excess dislocations, whose lines cross the area $\mathcal{A}$ in the reference configuration, in the following way $$\label{eq:2.2}
\mathbf{b}_r=\oint_{\mathcal{C}} \mathbf{F}^p\cdot d\mathbf{x} ,$$ where $\mathcal{C}$ is the close contour surrounding $\mathcal{A}$. @Le2014nonlinear have shown that, in the continuum limit, when the atomic distance goes to zero at the fixed sizes of the representative volume element and the fixed density of dislocations per area of unit cell, integral gives the total closure failure induced by $\mathbf{F}^p$ which must be equal to the resultant Burgers vector. It is natural to assume $\mathbf{F}^p$ continuously differentiable in this continuum limit, so, applying Stoke’s theorem we get from $$\mathbf{b}_r=- \int _{\mathcal{A}}(\mathbf{F}^p\times \nabla ) \cdot \mathbf{n}da,$$ where $\times $ denotes the vector product, $da$ the surface element, and $\mathbf{n}$ the unit vector normal to $\mathcal{A}$. This legitimates the introduction of the dislocation density tensor $$\mathbf{T}=-\mathbf{F}^p\times \nabla .$$ For the plastic deformation taken from $$\mathbf{T}=-\mathbf{F}^p\times \nabla =\mathbf{s}\otimes (\nabla \beta \times \mathbf{m}).$$ If we choose now an infinitesimal area $da$ with the unit normal vector ${\boldsymbol{\tau}}$, then the resultant Burgers vector of all excess dislocations, whose dislocation lines cross this area at right angle is given by $$\mathbf{b}_r=\mathbf{T}\cdot {\boldsymbol{\tau}}\, da=-\mathbf{s}\, (\nabla \beta \cdot {\boldsymbol{\nu}})da=-\mathbf{s}\, \partial _\nu \beta \, da.$$ This resultant Burgers vector can be decomposed into the sum of two vectors $$\mathbf{b}_r=\mathbf{b}_{r\perp}+\mathbf{b}_{r\parallel}=-({\boldsymbol{\nu}}s_\nu +{\boldsymbol{\tau}}s_\tau )\partial _\nu \beta\, da,$$ where $s_\nu =\mathbf{s}\cdot {\boldsymbol{\nu}}=\nu _1$ and $s_\tau =\mathbf{s}\cdot {\boldsymbol{\tau}}=\tau _1$ are the projections of the slip vector onto the normal and tangential direction to the dislocation line, respectively. This allows us to define two scalar densities (or the numbers of excess dislocations per unit area) of edge and screw dislocations $$\label{eq:2.3a}
\begin{split}
\rho_\perp =\frac{|\mathbf{b}_{r\perp}|}{b}=\frac{1}{b}|s_\nu \partial _\nu \beta| =\frac{1}{b}\left| \frac{l_{,1}(\beta _{,1}l_{,1}+\beta _{,2}l_{,2})}{l_{,1}^2+l_{,2}^2}\right| , 
\\
\rho_\parallel =\frac{|\mathbf{b}_{r\parallel}|}{b}=\frac{1}{b}|s_\tau \partial _\nu \beta|=\frac{1}{b}\left| \frac{l_{,2}(\beta _{,1}l_{,1}+\beta _{,2}l_{,2})}{l_{,1}^2+l_{,2}^2}\right| ,
\end{split}$$ with $b$ the magnitude of Burgers vector. We see that the three-dimensional dislocation densities $\rho_\perp $ and $\rho_\parallel $ depend on both the gradient of the plastic slip and the gradient of the loop function $l(\mathbf{x})$ through the vectors ${\boldsymbol{\nu}}$ and ${\boldsymbol{\tau}}$.

Consider now the case of motion of dislocation loops in the slip plane. In this case we allow the loop function to depend explicitly on time $t$ such that equation $$\label{eq:2.4}
l(x_1,x_2,c_3,t)=c$$ with fixed constants $c_3$ and $c$ describes one and the same dislocation line during its motion in the slip plane. Letting $\zeta $ be the variable along the dislocation line, we may represent the level curve defined by in the parametric form $$x_1=x_1(\zeta ,t),\quad x_2=x_2(\zeta ,t).$$ Fixing $\zeta $ and taking the differential of we obtain $$l_{,1}dx_1+l_{,2}dx_2+l_{,t}dt=0$$ that yields $$l_{,1}\frac{dx_1}{dt}+l_{,2}\frac{dx_2}{dt}=-l_{,t}.$$ Since $\mathbf{v}=\frac{dx_1}{dt}\mathbf{e}_1+\frac{dx_2}{dt}\mathbf{e}_2$ is the velocity of the fixed point on the dislocation line with coordinate $\zeta $, we define the normal velocity of the dislocation line as follows $$\label{eq:2.5}
v_\nu =\mathbf{v}\cdot {\boldsymbol{\nu}}=\frac{1}{\sqrt{l_{,1}^2+l_{,2}^2}}(l_{,1}\frac{dx_1}{dt}+l_{,2}\frac{dx_2}{dt})=-\frac{\dot{l}}{\sqrt{l_{,1}^2+l_{,2}^2}},$$ with $\dot{l}=l_{,t}$. This kinematic quantity will be used in the case of non-vanishing resistance to dislocation motion.

Governing equations for single crystals deforming in single slip
================================================================

According to @Kroener92, the elastic deformation $\mathbf{F}^e$ and the dislocation densities $\rho_\perp $ and $\rho_\parallel $ characterize the current state of the crystal, so these quantities are the state variables of the continuum dislocation theory. The reason why the plastic deformation $\mathbf{F}^p$ cannot be qualified for the state variable is that it depends on the cut surfaces and consequently on the whole history of creating dislocations. Likewise, the gradient of plastic strain tensor $\mathbf{C}^p$ cannot be used as the state variable by the same reason. In contrary, the dislocation densities depend only on the characteristics of dislocations in the current state (Burgers vector and positions of dislocation lines) and not on how they are created, so $\rho_\perp $ and $\rho_\parallel $, in addition to $\mathbf{F}^e$, are the proper state variables. Thus, if we consider isothermal processes of deformation, then the free energy per unit volume of crystal (assumed as macroscopically homogeneous) must be a function of $\mathbf{F}^e$, $\rho_\perp $, and $\rho_\parallel $ $$\psi = \psi (\mathbf{F}^e,\rho_\perp ,\rho_\parallel ).$$ Now, if we superimpose an elastic rotation $\mathbf{R}$ onto the actual deformation of the body, then the total and elastic deformation change according to $$\mathbf{F}^{\ast}=\mathbf{R}\cdot \mathbf{F},\quad \mathbf{F}^{e\ast}=\mathbf{R}\cdot \mathbf{F}^e.$$ At the same time, the dislocation densities $\rho_\perp $ and $\rho_\parallel $ remain unchanged. As such superimposed elastic rotation does not change the elastic strain and the dislocation densities, we expect that the energy remains unchanged. The standard argument (see, e.g., [@Gurtin1981]) leads then to $$\psi = \psi (\mathbf{C}^e,\rho_\perp ,\rho_\parallel ),$$ where $\mathbf{C}^e$ is the elastic strain defined by $$\mathbf{C}^e=\mathbf{F}^{eT}\cdot \mathbf{F}^e.$$

Let the undeformed single crystal occupy some region $\mathcal{V}$ of the three-dimensional space. The boundary of this region, $\partial \mathcal{V}$, is assumed to be the closure of union of two non-intersecting surfaces, $\partial _k$ and $\partial _s$. Let the placement be a given smooth function of coordinates $$\label{eq:2.6}
\mathbf{y}(\mathbf{x})=\mathbf{x}+\mathbf{u}_0(\mathbf{x}) \quad \text{at $\partial _k$},$$ where $\mathbf{u}_0(\mathbf{x})=\mathbf{y}(\mathbf{x})-\mathbf{x}$ is the given displacement vector. Such condition does not admit dislocations to reach this part $\partial _k$ of the boundary, so we set $$\label{eq:2.6a}
\beta (\mathbf{x})=0, \quad l(\mathbf{x})=0 \quad \text{at $\partial _k$}.$$ At the remaining part $\partial _s$ the “dead” load (traction) $\mathbf{t}$ is specified. Note that, in case the whole boundary is free, we do not have any kinematic constraint at $\partial \mathcal{V}$. However, as the energy density is invariant with respect to the shift of the loop function $l(\mathbf{x})$ on an arbitrary constant which does not change the dislocation densities, we can impose on this scalar function the following constraint $$\int_{\mathcal{V}}l(\mathbf{x})\, dx=0,$$ where $dx=dx_1dx_2dx_3$ denotes the volume element. If no body force acts on this crystal, then its energy functional is defined as $$\label{eq:2.7}
I[\mathbf{y}(\mathbf{x}),\beta (\mathbf{x}),l(\mathbf{x})]=\int_{\mathcal{V}}w(\mathbf{F},\beta ,\nabla \beta ,\nabla l)\, dx-\int_{\partial _s} \mathbf{t} \cdot \mathbf{y}\, da,$$ where $$\label{eq:2.8}
w(\mathbf{F},\beta ,\nabla \beta ,\nabla l)=\psi (\mathbf{C}^e,\rho_\perp ,\rho_\parallel ),$$ Provided the resistance to the dislocation motion is negligibly small and no surfaces of discontinuity occur inside crystals, then the following variational principle is valid for single crystals with one active slip system: the true placement vector $\check{\mathbf{y}}(\mathbf{x})$, the true plastic slip $\check{\beta }(\mathbf{x})$, and the true loop function $\check{l}(\mathbf{x})$ in the [*final*]{} equilibrium state of deformation minimize energy functional among all continuously differentiable fields $\mathbf{y}(\mathbf{x})$, $\beta (\mathbf{x})$, and $l(\mathbf{x})$ satisfying constraints and .

Let us derive the equilibrium equations from this variational principle. We compute the first variation of functional $$\delta I= \int_{\mathcal{V}}
\left( \mathbf{P}\mathbf{:}\delta \mathbf{y}\nabla +\frac{\partial w}{\partial \beta }\delta \beta +\frac{\partial w}{\partial \nabla \beta }\cdot \nabla \delta \beta +\frac{\partial w}{\nabla l}\cdot \nabla \delta l \right) dx
- \int_{\partial _s} \mathbf{t} \cdot \delta \mathbf{y}\, da ,$$ where $\mathbf{P}=\partial w /\partial \mathbf{F}$. Integrating the first, third, and fourth term by parts with the help of Gauss’ theorem and taking the conditions and into account, we obtain $$\begin{gathered}
\label{eq:2.9}
\delta I = \int_{\mathcal{V}}[ -\delta \mathbf{y}\cdot (\mathbf{P}\cdot \nabla ) +(w_\beta -\nabla \cdot w_{\nabla \beta })\delta \beta - (\nabla \cdot w_{\nabla l})\delta l ] \, dx
\\
+ \int_{\partial _s} [(\mathbf{P}\cdot \mathbf{n}-\mathbf{t} )\cdot \delta \mathbf{y} +w_{\nabla \beta }\cdot \mathbf{n}\, \delta \beta]\, da + \int_{\partial \mathcal{V}}w_{\nabla l }\cdot \mathbf{n}\, \delta l \, da=0.\end{gathered}$$ Equation implies that the minimizer must satisfy in $\mathcal{V}$ the equilibrium equations $$\label{eq:2.10}
\mathbf{P}\cdot \nabla =0,
\quad
-w_\beta +\nabla \cdot w_{\nabla \beta }=0, \quad \nabla \cdot w_{\nabla l} =0,$$ subjected to the kinematic boundary conditions and at $\partial _k$, and the following natural boundary conditions $$\label{eq:2.11}
\mathbf{P}\cdot \mathbf{n}=\mathbf{t},
\quad
w_{\nabla \beta }\cdot \mathbf{n}=0 , \quad w_{\nabla l }\cdot \mathbf{n}=0 \quad \text{at $\partial _s$}.$$ We call $\mathbf{P}$ the first Piola-Kirchhoff stress tensor, $\tau _r=-w_\beta $ the resolved shear stress (or Schmid stress), and $\varsigma =-\nabla \cdot w_{\nabla \beta }$ the back stress. The first equation of is nothing else but the equilibrium of macro-forces acting on the crystal, the second equation represents the equilibrium of micro-forces acting on dislocations, while the last one expresses the equilibrium condition for the curved dislocation lines.

The constitutive equations for $\mathbf{P}=w_{\mathbf{F}}$, $-w_\beta $, $w_{\nabla \beta }$, and $w_{\nabla l}$ can easily be obtained from the free energy density . First, we express $\mathbf{F}^e$ in terms of $\mathbf{F}$ and $\beta $ with the use of and $$\mathbf{F}^e=\mathbf{F}\cdot \mathbf{F}^{p-1}=\mathbf{F}\cdot (\mathbf{I}-\beta \mathbf{s}\otimes \mathbf{m}).$$ Now, the standard differentiation using the chain rule and the above relation yields the first Piola-Kirchhoff stress tensor $$\label{eq:2.12}
\mathbf{P}=w_{\mathbf{F}}=2\mathbf{F}^e\cdot \psi _{\mathbf{C}^e} \cdot \mathbf{F}^{p-T}.$$ For the resolved shear stress (Schmid stress) we get $$\label{eq:2.13}
\tau _r=-w_\beta =2\mathbf{s}\cdot \mathbf{F}^{p-T}\cdot \mathbf{C}^e \cdot \psi _{\mathbf{C}^e}\cdot \mathbf{m}.$$ Likewise, from follows $$\label{eq:2.14}
w_{\nabla \beta }=\frac{1}{b}[\psi _{\rho_\perp }{\text{sign}}(s_\nu \partial _\nu \beta )s_\nu +\psi _{\rho_\parallel }{\text{sign}}(s_\tau \partial _\nu \beta )s_\tau ] {\boldsymbol{\nu}}.$$ Thus, the vector $w_{\nabla \beta }$ is two-dimensional. Finally, we compute $w_{\nabla l }$ directly in components using formulas . Since $\rho _\perp $ and $\rho _\parallel $ do not depend on $l_{,3}$, so $w_{l_{,3}}=0$, and the vector $w_{\nabla l}$ is also two-dimensional. For its first two components we have $$\label{eq:2.15}
\begin{split}
w_{l_{,1}}=\frac{1}{b}\left[ \psi _{\rho _\perp }{\text{sign}}(s_\nu \partial _\nu \beta )\left( -\frac{2l_{,1}^2(\beta _{,1}l_{,1}+\beta _{,2}l_{,2})}{(l_{,1}^2+l_{,2}^2)^2}+\frac{2\beta _{,1}l_{,1}+\beta _{,2}l_{,2}}{l_{,1}^2+l_{,2}^2} \right) \right.
\\
+\left. \psi _{\rho _\parallel }{\text{sign}}(s_\tau \partial _\nu \beta )\left( -\frac{2l_{,1}l_{,2}(\beta _{,1}l_{,1}+\beta _{,2}l_{,2})}{(l_{,1}^2+l_{,2}^2)^2}+\frac{\beta _{,1}l_{,2}}{l_{,1}^2+l_{,2}^2} \right) \right] ,
\\
w_{l_{,2}}=\frac{1}{b}\left[ \psi _{\rho _\perp }{\text{sign}}(s_\nu \partial _\nu \beta )\left( -\frac{2l_{,1}l_{,2}(\beta _{,1}l_{,1}+\beta _{,2}l_{,2})}{(l_{,1}^2+l_{,2}^2)^2}+\frac{\beta _{,2}l_{,1}}{l_{,1}^2+l_{,2}^2} \right) \right.
\\
+\left. \psi _{\rho _\parallel }{\text{sign}}(s_\tau \partial _\nu \beta )\left( -\frac{2l_{,2}^2(\beta _{,1}l_{,1}+\beta _{,2}l_{,2})}{(l_{,1}^2+l_{,2}^2)^2}+\frac{\beta _{,1}l_{,1}+2\beta _{,2}l_{,2}}{l_{,1}^2+l_{,2}^2} \right) \right] .
\end{split}$$ Substituting the constitutive equations - into - we get the completely new system of equations and boundary conditions which, together with and , enable one to determine $\check{\mathbf{y}}(\mathbf{x})$, $\check{\beta }(\mathbf{x})$, and $\check{l}(\mathbf{x})$. Note that equations $_1$ and $_2$ are coupled via the first Piola-Kirchhoff stress tensor and the Schmid stress containing both $\mathbf{F}$ and $\beta $, while equations $_2$ and $_3$ are coupled because both contain the gradients of $\beta $ and $l$. All equations are strongly nonlinear partial differential equations.

The above theory has been developed for the case of negligibly small resistance to dislocation motion and plastic slip. In real crystals there is however always the resistance to the dislocation motion and plastic slip causing the energy dissipation that changes the above variational principle as well as the equilibrium conditions. We assume that the dissipation function depends on the plastic slip rate $\dot{\beta }$ and on the normal velocity of the dislocation loop $v_\nu $ given by (or, equivalently, on $\dot{l}$). Thus, $$D=D(\dot{\beta },\dot{l}),$$ When the dissipation is taken into account, the above formulated variational principle must be modified. Following [@Sedov65; @Berdichevsky1967] we require that the true placement $\check{\mathbf{y}}(\mathbf{x},t)$, the true plastic slips $\check{\beta }(\mathbf{x},t)$, and the true loop function $\check{l}(\mathbf{x},t)$ obey the variational equation $$\label{eq:2.15b}
\delta I+\int_{\mathcal{V}} (\frac{\partial D}{\partial \dot{\beta }}\delta \beta +\frac{\partial D}{\partial \dot{l}}\delta l)\, dx=0$$ for all variations of admissible fields $\mathbf{y}(\mathbf{x},t)$, $\beta (\mathbf{x},t)$, and $l(\mathbf{x},t)$ satisfying the constraints and . Together with the above formula for $\delta I$ and the arbitrariness of $\delta \mathbf{y}$, $\delta \beta $, and $\delta l$ in $\mathcal{V}$ as well as at $\partial _s$, equation yields $$\label{eq:2.15c}
\mathbf{P}\cdot \nabla =0,
\quad
-w_\beta +\nabla \cdot w_{\nabla \beta }=\frac{\partial D}{\partial \dot{\beta }}, \quad \nabla \cdot w_{\nabla l}=\frac{\partial D}{\partial \dot{l}},$$ which are subjected to the kinematic boundary conditions and , and the natural boundary conditions . The constitutive equations remain exactly the same as -. For the rate-independent theory the dissipation function can be assumed in a simple form $$D=K_1|\dot{\beta }|+K_2|\dot{l}|,$$ with $K_1$ and $K_2$ being positive constants. We call $K_1$ the critical resolved shear stress and $K_2$ the Peierls threshold. In this case equations $_{2,3}$ become $$-w_\beta +\nabla \cdot w_{\nabla \beta }=K_1\, {\text{sign}}\dot{\beta }, \quad \nabla \cdot w_{\nabla l}=K_2\, {\text{sign}}\dot{l}$$ for non-vanishing $\dot{\beta}$ and $\dot{l}$. These are the yield conditions for $\beta $ and $l$: $\dot{\beta}$ and $\dot{l}$ are non-zero if and only if $$|-w_\beta +\nabla \cdot w_{\nabla \beta }|=K_1, \quad |\nabla \cdot w_{\nabla l}|=K_2.$$ On the contrary, if the expressions on the left-hand sides are less than $K_1$ and $K_2$, the plastic slip cannot evolve and the dislocation lines cannot move: $\dot{\beta}=0$ and $\dot{l}=0$. Thus, they are frozen in the crystal.

Extension to multiple slip
==========================

The extension to the case of single crystals having $n$ active slip systems can be done straightforwardly under the assumption[^2] $$\label{eq:2.16}
\mathbf{F}^p=\mathbf{I}+\sum_{\mathfrak{a}=1}^n \beta ^\mathfrak{a}(\mathbf{x}) \mathbf{s}^\mathfrak{a}\otimes \mathbf{m}^\mathfrak{a} ,$$ with $\beta ^\mathfrak{a}$ being the plastic slip, where the pair of constant and mutually orthogonal unit vectors $\mathbf{s}^\mathfrak{a}$ and $\mathbf{m}^\mathfrak{a}$ is used to denote the slip direction and the normal to the slip planes of the corresponding $\mathfrak{a}$-th slip system, respectively. Here and later, the Gothic upper index $\mathfrak{a}$ running from 1 to $n$ numerates the slip systems, so one could clearly distinguish $\beta ^{\mathfrak{a}}$ from the power function. We denote by $\mathbf{p}^{\mathfrak{a}}$ the unit vector lying in the slip plane such that $\mathbf{s}^{\mathfrak{a}}$, $\mathbf{p}^{\mathfrak{a}}$, and $\mathbf{m}^{\mathfrak{a}}$ form a right-handed basis vectors. For each slip system we can introduce the coordinates associated with these basis vectors $$\label{eq:2.17}
\xi ^{\mathfrak{a}}_1=\mathbf{s}^{\mathfrak{a}}\cdot \mathbf{x}, \quad \xi ^{\mathfrak{a}}_2=\mathbf{p}^{\mathfrak{a}}\cdot \mathbf{x},\quad \xi ^{\mathfrak{a}}_3=\mathbf{m}^{\mathfrak{a}}\cdot \mathbf{x}.$$ Equations can be regarded as the one-to-one linear transformation relating ${\boldsymbol{\xi}}^{\mathfrak{a}}$ and $\mathbf{x}$ according to $${\boldsymbol{\xi}}^{\mathfrak{a}}=\mathbf{M}^{\mathfrak{a}}\cdot \mathbf{x}, \quad \mathbf{x}=\mathbf{M}^{{\mathfrak{a}}-1}{\boldsymbol{\xi}}^{\mathfrak{a}},$$ where $\mathbf{M}^{\mathfrak{a}}$ is the $3\times 3$ matrix whose rows are basis vectors $\mathbf{s}^{\mathfrak{a}}$, $\mathbf{p}^{\mathfrak{a}}$, and $\mathbf{m}^{\mathfrak{a}}$. Thus, any function of $\mathbf{x}$ can be expressed as function of ${\boldsymbol{\xi}}^{\mathfrak{a}}$ and [*vice versa*]{}. For the plastic slip $\beta ^\mathfrak{a}$ caused by dislocations of the slip system ${\mathfrak{a}}$ we assume that their lines lie completely in the slip planes parallel to the $(\xi^{\mathfrak{a}}_1,\xi^{\mathfrak{a}}_2)$-plane. To describe the latter we introduce the loop function $l^{\mathfrak{a}}(\xi^{\mathfrak{a}}_1,\xi^{\mathfrak{a}}_2, \xi^{\mathfrak{a}}_3)$ such that its level curves $$\label{eq:2.18}
l^{\mathfrak{a}}(\xi ^{\mathfrak{a}}_1,\xi ^{\mathfrak{a}}_2,c _3)=c,$$ where $c_3$ and $c$ are constants, coincide with the dislocation lines. We denote by ${\boldsymbol{\nu}}^{\mathfrak{a}}$ and ${\boldsymbol{\tau}}^{\mathfrak{a}}$ the plane unit vectors normal and tangential to the dislocation line. From equation follow $${\boldsymbol{\nu}}^{\mathfrak{a}}=\frac{1}{\sqrt{(l^{\mathfrak{a}}_{;1})^2+(l^{\mathfrak{a}}_{;2})^2}}(l^{\mathfrak{a}}_{;1}\mathbf{s}^{\mathfrak{a}}+l^{\mathfrak{a}}_{;2}\mathbf{p}^{\mathfrak{a}}), \quad {\boldsymbol{\tau}}^{\mathfrak{a}}=\frac{1}{\sqrt{(l^{\mathfrak{a}}_{;1})^2+(l^{\mathfrak{a}}_{;2})^2}}(-l^{\mathfrak{a}}_{;2}\mathbf{s}^{\mathfrak{a}}+l^{\mathfrak{a}}_{;1}\mathbf{p}^{\mathfrak{a}}),$$ where the semicolon in indices denotes the partial derivatives of the loop function with respect to $\xi ^{\mathfrak{a}}_1$, $\xi ^{\mathfrak{a}}_2$, so these vectors lie in the slip planes parallel to the $(\xi _1,\xi _2)$-plane as expected.

For the plastic deformation the dislocation density tensor becomes $$\label{eq:2.20}
\mathbf{T}=-\mathbf{F}^p\times \nabla =\sum_{\mathfrak{a}=1}^n \mathbf{s}^{\mathfrak{a}}\otimes (\nabla \beta ^{\mathfrak{a}}\times \mathbf{m}^{\mathfrak{a}}).$$ To characterize the geometrically necessary dislocations belonging to one slip system we consider one term $\mathbf{T}^{\mathfrak{a}}=\mathbf{s}^{\mathfrak{a}}\otimes (\nabla \beta ^{\mathfrak{a}}\times \mathbf{m}^{\mathfrak{a}})$ in the sum . Let us choose an infinitesimal area $da$ with the unit normal vector ${\boldsymbol{\tau}}^{\mathfrak{a}}$ and compute the resultant Burgers vector of all excess dislocations of the system ${\mathfrak{a}}$, whose dislocation lines cross this area at right angle $$\mathbf{b}^{\mathfrak{a}}_r=\mathbf{T}^{\mathfrak{a}}\cdot {\boldsymbol{\tau}}^{\mathfrak{a}}\, da=-\mathbf{s}^{\mathfrak{a}}\, (\nabla \beta ^{\mathfrak{a}}\cdot {\boldsymbol{\nu}}^{\mathfrak{a}})da=-\mathbf{s}^{\mathfrak{a}}\, \partial _{\nu ^{\mathfrak{a}}}\beta ^{\mathfrak{a}}\, da.$$ This resultant Burgers vector can be decomposed into the sum of two vectors $$\mathbf{b}^{\mathfrak{a}}_r=\mathbf{b}^{\mathfrak{a}}_{r\perp}+\mathbf{b}^{\mathfrak{a}}_{r\parallel}=-({\boldsymbol{\nu}}^{\mathfrak{a}}s^{\mathfrak{a}}_\nu +{\boldsymbol{\tau}}^{\mathfrak{a}}s^{\mathfrak{a}}_\tau )\partial _{\nu ^{\mathfrak{a}}}\beta^{\mathfrak{a}}\, da,$$ where $s^{\mathfrak{a}}_\nu =\mathbf{s}^{\mathfrak{a}}\cdot {\boldsymbol{\nu}}^{\mathfrak{a}}$ and $s^{\mathfrak{a}}_\tau =\mathbf{s}^{\mathfrak{a}}\cdot {\boldsymbol{\tau}}^{\mathfrak{a}}$ are the projections of the slip vector onto the normal and tangential direction to the dislocation line, respectively. This allows us to define two scalar densities of edge and screw dislocations of the corresponding slip system $$\label{eq:2.21}
\begin{split}
\rho^{\mathfrak{a}}_\perp =\frac{|\mathbf{b}^{\mathfrak{a}}_{r\perp}|}{b} =\frac{1}{b}|s^{\mathfrak{a}}_\nu \partial _{\nu ^{\mathfrak{a}}}\beta^{\mathfrak{a}}| =\frac{1}{b}\left| \frac{l^{\mathfrak{a}}_{;1}(\beta ^{\mathfrak{a}}_{;1}l^{\mathfrak{a}}_{;1}+\beta ^{\mathfrak{a}}_{;2}l^{\mathfrak{a}}_{;2})}{(l^{\mathfrak{a}}_{;1})^2+(l^{\mathfrak{a}}_{;2})^2}\right| , 
\\
\rho^{\mathfrak{a}}_\parallel =\frac{|\mathbf{b}^{\mathfrak{a}}_{r\parallel }|}{b} =\frac{1}{b}|s^{\mathfrak{a}}_\tau \partial _{\nu ^{\mathfrak{a}}}\beta^{\mathfrak{a}}|=\frac{1}{b}\left| \frac{l^{\mathfrak{a}}_{;2}(\beta ^{\mathfrak{a}}_{;1}l^{\mathfrak{a}}_{;1}+\beta ^{\mathfrak{a}}_{;2}l^{\mathfrak{a}}_{;2})}{(l^{\mathfrak{a}}_{;1})^2+(l^{\mathfrak{a}}_{;2})^2}\right| ,
\end{split}$$ We see that the dislocation densities $\rho^{\mathfrak{a}}_\perp $ and $\rho^{\mathfrak{a}}_\parallel $ depend only on the partial derivatives $\beta ^{\mathfrak{a}}_{;\alpha }$ and $l^{\mathfrak{a}}_{;\alpha }$, $\alpha =1,2$. For the moving dislocations we allow the loop functions to depend on time $t$ such that the level curves $$l^{\mathfrak{a}}(\xi ^{\mathfrak{a}}_1,\xi ^{\mathfrak{a}}_2,c _3,t)=c,$$ with $c_3$ and $c$ being constants, coincide with the dislocation lines during their motion. Similar to the single slip we introduce the normal velocities of dislocation lines as follows $$v^{\mathfrak{a}}_\nu =\mathbf{v}^{\mathfrak{a}}\cdot {\boldsymbol{\nu}}^{\mathfrak{a}}=-\frac{\dot{l}^{\mathfrak{a}}}{\sqrt{(l^{\mathfrak{a}}_{;1})^2+(l^{\mathfrak{a}}_{;2})^2}},$$ with $\dot{l}=l_{,t}$. These kinematic quantities will be used in the model with dissipation.

From the above discussion of kinematics we see that a single crystal with $n$ active slip systems is a generalized continuum with $3+2n$ degrees of freedom at each point: $\mathbf{y}(\mathbf{x})$, $\beta ^{\mathfrak{a}}(\mathbf{x})$, and $l^{\mathfrak{a}}({\boldsymbol{\xi}}^{\mathfrak{a}}(\mathbf{x}))$, ${\mathfrak{a}}=1,\ldots ,n$. We require that the free energy per unit volume of crystal (assumed as macroscopically homogeneous) must be a function of $\mathbf{C}^e=\mathbf{F}^{eT}\cdot \mathbf{F}^e$ (where $\mathbf{F}^e=\mathbf{F}\cdot \mathbf{F}^{p-1}$), $\rho^{\mathfrak{a}}_\perp $, and $\rho^{\mathfrak{a}}_\parallel $ $$\psi = \psi (\mathbf{C}^e,\rho^{\mathfrak{a}}_\perp ,\rho^{\mathfrak{a}}_\parallel ).$$ Under the same loading condition as for the crystal with single slip we write down the energy functional $$\label{eq:2.22}
I[\mathbf{y}(\mathbf{x}),\beta ^{\mathfrak{a}}(\mathbf{x}),l^{\mathfrak{a}}({\boldsymbol{\xi}}^{\mathfrak{a}}(\mathbf{x}))]=\int_{\mathcal{V}}w(\mathbf{F},\beta ^{\mathfrak{a}},\nabla \beta ^{\mathfrak{a}},\nabla l^{\mathfrak{a}})\, dx-\int_{\partial _s} \mathbf{t} \cdot \mathbf{y}\, da,$$ where $$w(\mathbf{F},\beta ^{\mathfrak{a}},\nabla \beta ^{\mathfrak{a}},\nabla l^{\mathfrak{a}})=\psi (\mathbf{C}^e,\rho^{\mathfrak{a}}_\perp ,\rho^{\mathfrak{a}}_\parallel ).$$ Provided the resistance to the dislocation motion is negligibly small, we formulate the following variational principle for single crystals with $n$ active slip systems: the true placement vector $\check{\mathbf{y}}(\mathbf{x})$, the true plastic slips $\check{\beta }^{\mathfrak{a}}(\mathbf{x})$, and the true loop functions $\check{l}^{\mathfrak{a}}({\boldsymbol{\xi}}^{\mathfrak{a}}(\mathbf{x}))$ in the [*final*]{} equilibrium state of deformation minimize energy functional among all continuously differentiable fields $\mathbf{y}(\mathbf{x})$, $\beta ^{\mathfrak{a}}(\mathbf{x})$, and $l^{\mathfrak{a}}({\boldsymbol{\xi}}^{\mathfrak{a}}(\mathbf{x}))$ satisfying the constraints $$\label{eq:2.23}
\mathbf{y}(\mathbf{x})=\mathbf{x}+\mathbf{u}_0(\mathbf{x}),\quad \beta ^{\mathfrak{a}}(\mathbf{x})=0, \quad l^{\mathfrak{a}}(\mathbf{x})=0 \quad \text{at $\partial _k$},$$

Applying the same calculus of variation and taking into account the arbitrariness of the variations of $\mathbf{y}(\mathbf{x})$, $\beta (\mathbf{x})$, and $l({\boldsymbol{\xi}}^{\mathfrak{a}}(\mathbf{x}))$ in $\mathcal{V}$ as well as at $\partial _s$, one can show that the minimizer must satisfy in $\mathcal{V}$ the equilibrium equations $$\label{eq:2.25}
\mathbf{P}\cdot \nabla =0,
\quad
-w_\beta ^{\mathfrak{a}}+\nabla \cdot w_{\nabla \beta ^{\mathfrak{a}}}=0, \quad \nabla \cdot w_{\nabla l^{\mathfrak{a}}} =0,$$ subjected to the kinematic boundary conditions at $\partial _k$ and the following natural boundary conditions $$\label{eq:2.26}
\mathbf{P}\cdot \mathbf{n}=\mathbf{t},
\quad
w_{\nabla \beta ^{\mathfrak{a}}}\cdot \mathbf{n}=0, \quad w_{\nabla l^{\mathfrak{a}}}\cdot \mathbf{n}=0 \quad \text{at $\partial _s$}.$$ The constitutive equations for $\mathbf{P}=w_{\mathbf{F}}$, $-w_{\beta ^{\mathfrak{a}}}$, $w_{\nabla \beta ^{\mathfrak{a}}}$, and $w_{\nabla l^{\mathfrak{a}}}$ can easily be obtained from the above free energy density by standard differentiation. For the first Piola-Kirchhoff stress tensor and the Schmid stresses we have $$\label{eq:2.27}
\mathbf{P}=w_{\mathbf{F}}=2\mathbf{F}^e\cdot \psi _{\mathbf{C}^e} \cdot \mathbf{F}^{p-T}.$$ $$\label{eq:2.28}
\tau ^{\mathfrak{a}}_r=-w_{\beta ^{\mathfrak{a}}}=2\mathbf{s}^{\mathfrak{a}}\cdot \mathbf{F}^{p-T}\cdot \mathbf{C}^e \cdot \psi _{\mathbf{C}^e}\cdot \mathbf{m}^{\mathfrak{a}}.$$ Likewise, from follows $$\label{eq:2.29}
w_{\nabla \beta ^{\mathfrak{a}}}=\frac{1}{b}[\psi _{\rho^{\mathfrak{a}}_\perp }{\text{sign}}(s^{\mathfrak{a}}_\nu \partial _{\nu ^{\mathfrak{a}}} \beta ^{\mathfrak{a}})s^{\mathfrak{a}}_\nu +\psi _{\rho^{\mathfrak{a}}_\parallel }{\text{sign}}(s^{\mathfrak{a}}_\tau \partial _{\nu ^{\mathfrak{a}}} \beta ^{\mathfrak{a}})s^{\mathfrak{a}}_\tau ] {\boldsymbol{\nu}}^{\mathfrak{a}}.$$ Thus, the vectors $w_{\nabla \beta ^{\mathfrak{a}}}$ are two-dimensional. Finally, for $w_{\nabla l^{\mathfrak{a}}}$ we have $$\label{eq:2.30}
w_{\nabla l^{\mathfrak{a}}}=\left( \psi _{\rho _\perp ^{\mathfrak{a}}}\frac{\partial \rho _\perp ^{\mathfrak{a}}}{\partial l^{\mathfrak{a}}_{;1}}+\psi _{\rho _\parallel ^{\mathfrak{a}}}\frac{\partial \rho _\parallel ^{\mathfrak{a}}}{\partial l^{\mathfrak{a}}_{;1}}\right)\mathbf{s}^{\mathfrak{a}}+\left( \psi _{\rho _\perp ^{\mathfrak{a}}}\frac{\partial \rho _\perp ^{\mathfrak{a}}}{\partial l^{\mathfrak{a}}_{;2}}+\psi _{\rho _\parallel ^{\mathfrak{a}}}\frac{\partial \rho _\parallel ^{\mathfrak{a}}}{\partial l^{\mathfrak{a}}_{;2}}\right)\mathbf{p}^{\mathfrak{a}}.$$ Differentiating formulas for the dislocation densities $\rho^{\mathfrak{a}}_\perp $ and $\rho^{\mathfrak{a}}_\parallel $ with respect to $l^{\mathfrak{a}}_{;1}$ and $l^{\mathfrak{a}}_{;2}$, we get $$\begin{aligned}
\frac{\partial \rho _\perp ^{\mathfrak{a}}}{\partial l^{\mathfrak{a}}_{;1}}&=\frac{1}{b} {\text{sign}}(s^{\mathfrak{a}}_\nu \partial _{\nu ^{\mathfrak{a}}}\beta ^{\mathfrak{a}})\left( -\frac{2(l^{\mathfrak{a}}_{;1})^2(\beta ^{\mathfrak{a}}_{;1}l^{\mathfrak{a}}_{;1}+\beta ^{\mathfrak{a}}_{;2}l^{\mathfrak{a}}_{;2})}{((l^{\mathfrak{a}}_{;1})^2+(l^{\mathfrak{a}}_{;2})^2)^2}+\frac{2\beta ^{\mathfrak{a}}_{;1}l^{\mathfrak{a}}_{;1}+\beta ^{\mathfrak{a}}_{;2}l^{\mathfrak{a}}_{;2}}{(l^{\mathfrak{a}}_{;1})^2+(l^{\mathfrak{a}}_{;2})^2} \right) ,
\\
\frac{\partial \rho _\parallel ^{\mathfrak{a}}}{\partial l^{\mathfrak{a}}_{;1}}&=\frac{1}{b}{\text{sign}}(s^{\mathfrak{a}}_\tau \partial _{\nu ^{\mathfrak{a}}}\beta ^{\mathfrak{a}})\left( -\frac{2l^{\mathfrak{a}}_{;1}l^{\mathfrak{a}}_{;2}(\beta ^{\mathfrak{a}}_{;1}l^{\mathfrak{a}}_{;1}+\beta ^{\mathfrak{a}}_{;2}l^{\mathfrak{a}}_{;2})}{((l^{\mathfrak{a}}_{;1})^2+(l^{\mathfrak{a}}_{;2})^2)^2}+\frac{\beta ^{\mathfrak{a}}_{;1}l^{\mathfrak{a}}_{;2}}{(l^{\mathfrak{a}}_{;1})^2+(l^{\mathfrak{a}}_{;2})^2} \right) ,
\\
\frac{\partial \rho _\perp ^{\mathfrak{a}}}{\partial l^{\mathfrak{a}}_{;2}}&=\frac{1}{b}{\text{sign}}(s^{\mathfrak{a}}_\nu \partial _{\nu ^{\mathfrak{a}}}\beta ^{\mathfrak{a}})\left( -\frac{2l^{\mathfrak{a}}_{;1}l^{\mathfrak{a}}_{;2}(\beta ^{\mathfrak{a}}_{;1}l^{\mathfrak{a}}_{;1}+\beta ^{\mathfrak{a}}_{;2}l^{\mathfrak{a}}_{;2})}{((l^{\mathfrak{a}}_{;1})^2+(l^{\mathfrak{a}}_{;2})^2)^2}+\frac{\beta ^{\mathfrak{a}}_{;2}l^{\mathfrak{a}}_{;1}}{(l^{\mathfrak{a}}_{;1})^2+(l^{\mathfrak{a}}_{;2})^2} \right) ,
\\
\frac{\partial \rho _\parallel ^{\mathfrak{a}}}{\partial l^{\mathfrak{a}}_{;1}}&=\frac{1}{b} {\text{sign}}(s^{\mathfrak{a}}_\tau \partial _{\nu ^{\mathfrak{a}}}\beta ^{\mathfrak{a}})\left( -\frac{2(l^{\mathfrak{a}}_{;2})^2(\beta ^{\mathfrak{a}}_{;1}l^{\mathfrak{a}}_{;1}+\beta ^{\mathfrak{a}}_{;2}l^{\mathfrak{a}}_{;2})}{((l^{\mathfrak{a}}_{;1})^2+(l^{\mathfrak{a}}_{;2})^2)^2}+\frac{\beta ^{\mathfrak{a}}_{;1}l^{\mathfrak{a}}_{;1}+2\beta ^{\mathfrak{a}}_{;2}l^{\mathfrak{a}}_{;2}}{(l^{\mathfrak{a}}_{;1})^2+(l^{\mathfrak{a}}_{;2})^2} \right) .\end{aligned}$$ Substituting the constitutive equations - into - we get the completely new system of equations and boundary conditions which, together with , enables one to determine $\check{\mathbf{y}}(\mathbf{x})$, $\check{\beta }^{\mathfrak{a}}(\mathbf{x})$, and $\check{l}^{\mathfrak{a}}(\mathbf{x})$.

For the case of non-zero resistance to dislocation motion leading to the energy dissipation we take the dissipation function in the form $$D=D(\dot{\beta }^{\mathfrak{a}},\dot{l}^{\mathfrak{a}}),$$ We require that the true placement $\check{\mathbf{y}}(\mathbf{x},t)$, the true plastic slips $\check{\beta }^{\mathfrak{a}}(\mathbf{x},t)$, and the true loop functions $\check{l}^{\mathfrak{a}}(\mathbf{x},t)$ obey the variational equation $$\label{eq:2.32}
\delta I+\int_{\mathcal{V}} \sum_{{\mathfrak{a}}=1}^n (\frac{\partial D}{\partial \dot{\beta }^{\mathfrak{a}}}\delta \beta ^{\mathfrak{a}}+\frac{\partial D}{\partial \dot{l}^{\mathfrak{a}}}\delta l^{\mathfrak{a}})\, dx=0$$ for all variations of admissible fields $\mathbf{y}(\mathbf{x},t)$, $\beta ^{\mathfrak{a}}(\mathbf{x},t)$, and $l^{\mathfrak{a}}(\mathbf{x},t)$ satisfying the constraints . It is then easy to show by exactly the same arguments like those used at the end of the previous Section that equation yields $$\mathbf{P}\cdot \nabla =0,
\quad
-w_{\beta ^{\mathfrak{a}}}+\nabla \cdot w_{\nabla \beta ^{\mathfrak{a}}}=\frac{\partial D}{\partial \dot{\beta }^{\mathfrak{a}}}, \quad \nabla \cdot w_{\nabla l^{\mathfrak{a}}} =\frac{\partial D}{\partial \dot{l}^{\mathfrak{a}}},$$ which are subjected to the boundary conditions and . The constitutive equations remain exactly the same as -.

Small strain theory
===================

Let us simplify the above theory for small strains. In this case, instead of the placement $\mathbf{y}(\mathbf{x})$ we regards the displacement $\mathbf{u}(\mathbf{x})$ that is related to the former by $$\mathbf{u}(\mathbf{x})=\mathbf{y}(\mathbf{x})-\mathbf{x}$$ as the unknown function. Thus, the total compatible deformation is $$\mathbf{F}=\frac{\partial \mathbf{y}}{\partial \mathbf{x}} =\mathbf{I}+\mathbf{u}\nabla .$$ We assume that the displacement gradient $\mathbf{u}\nabla $ (called distortion) is small compared with $\mathbf{I}$. Concerning the plastic deformation given by we also assume that the plastic slips $\beta ^{\mathfrak{a}}$ are much smaller than 1. Using the the multiplicative resolution to express $\mathbf{F}^e$ through $\mathbf{F}$ and $\mathbf{F}^{p-1}$ and neglecting the small nonlinear terms in it, we get $$\mathbf{F}^e=\mathbf{F}\cdot \mathbf{F}^{p-1}=\mathbf{I}+{\boldsymbol{\beta}}^e$$ where $${\boldsymbol{\beta}}^e =\mathbf{u}\nabla -\sum_{{\mathfrak{a}}=1}^n \beta ^{\mathfrak{a}}\mathbf{s}^{\mathfrak{a}}\otimes \mathbf{m}^{\mathfrak{a}}.$$ The last equation can be interpreted as the additive resolution of the total distortion into the plastic and elastic parts. The total compatible strain tensor field can be obtained from the displacement field according to $${\boldsymbol{\varepsilon}}=\frac{1}{2}(\mathbf{u}\nabla +\nabla \mathbf{u}).$$ The incompatible plastic strain tensor field is the symmetric part of the plastic distortion field $${\boldsymbol{\varepsilon}}^p =\frac{1}{2}({\boldsymbol{\beta}}+{\boldsymbol{\beta}}^T)=\frac{1}{2}\sum_{{\mathfrak{a}}=1}^n \beta ^{\mathfrak{a}}(\mathbf{s}^{\mathfrak{a}}\otimes \mathbf{m}^{\mathfrak{a}}+\mathbf{m}^{\mathfrak{a}}\otimes \mathbf{s}^{\mathfrak{a}}).$$ Accordingly, the elastic strain tensor field is equal to $${\boldsymbol{\varepsilon}}^e={\boldsymbol{\varepsilon}}-{\boldsymbol{\varepsilon}}^p.$$

The dislocation densities remain exactly the same as in the finite strain theory. They are given by the formulas . Concerning the free energy density we will assume that it depends on the elastic strain ${\boldsymbol{\varepsilon}}^e$ and on the dislocation densities $\rho^{\mathfrak{a}}_\perp $ and $\rho^{\mathfrak{a}}_\parallel $ $$\psi = \psi ({\boldsymbol{\varepsilon}}^e,\rho^{\mathfrak{a}}_\perp ,\rho^{\mathfrak{a}}_\parallel ).$$ The energy of crystal containing dislocations reads $$\label{eq:3.5}
I[\mathbf{u}(\mathbf{x}),\beta (\mathbf{x}),l(\mathbf{x})]=\int_{\mathcal{V}}w(\mathbf{u}\nabla ,\beta ^{\mathfrak{a}},\nabla \beta ^{\mathfrak{a}},\nabla l^{\mathfrak{a}})\, dx-\int_{\partial _s} \mathbf{t}\cdot \mathbf{u}\, da,$$ where $w(\mathbf{u}\nabla ,\beta ^{\mathfrak{a}},\nabla \beta ^{\mathfrak{a}},\nabla l^{\mathfrak{a}})=\psi ({\boldsymbol{\varepsilon}}^e,\rho_\perp ,\rho_\parallel )$. Provided the resistance to the dislocation motion can be neglected, then the following variational principle is valid for single crystals: the true displacement field $\check{\mathbf{u}}(\mathbf{x})$, the true plastic slips $\check{\beta }^{\mathfrak{a}}(\mathbf{x})$, and the true loop functions $\check{l}^{\mathfrak{a}}(\mathbf{x})$ in the [*final*]{} state of deformation in equilibrium minimize energy functional among all admissible fields satisfying the constraints $$\label{eq:3.6}
\mathbf{u}(\mathbf{x})=\mathbf{u}_0(\mathbf{x}),\quad \beta ^{\mathfrak{a}}(\mathbf{x})=0,\quad l^{\mathfrak{a}}(\mathbf{x})=0  \quad \text{at $\partial _k$}.$$

The standard calculus of variation similar to the previous case leads to the equilibrium equations $$\label{eq:3.7}
{\boldsymbol{\sigma}}\cdot \nabla =0, \quad w_\beta ^{\mathfrak{a}}-\nabla \cdot w_{\nabla \beta ^{\mathfrak{a}}}=0, \quad \nabla \cdot w_{\nabla l^{\mathfrak{a}}} =0,$$ subjected to the kinematic boundary conditions at $\partial _k$ and the following natural boundary conditions $$\label{eq:3.8}
{\boldsymbol{\sigma}}\cdot \mathbf{n}=\mathbf{t},\quad w_{\nabla \beta ^{\mathfrak{a}}}\cdot \mathbf{n}=0, \quad w_{\nabla l^{\mathfrak{a}}}\cdot \mathbf{n}=0 \quad \text{at $\partial _s$}.$$ The constitutive equation for the Cauchy stress tensor becomes $$\label{eq:3.9}
{\boldsymbol{\sigma}}=\frac{\partial \psi }{\partial {\boldsymbol{\varepsilon}}^e}.$$ For the Schmid stress we obtain $$\label{eq:3.10}
\tau ^{\mathfrak{a}}_r=-w_{\beta ^{\mathfrak{a}}}=\mathbf{s}^{\mathfrak{a}}\cdot {\boldsymbol{\sigma}}\cdot \mathbf{m}^{\mathfrak{a}}.$$ The constitutive equations for $w_{\nabla \beta ^{\mathfrak{a}}}$ and $w_{\nabla l^{\mathfrak{a}}}$ remain unchanged as compared with and . The new set of governing equations and boundary conditions are obtained by substituting , , , and into the equilibrium equations and boundary conditions . Note that even for small strain the system of governing equations remain as a whole nonlinear.

It is a simple matter to modify the theory for the case of non-zero resistance to dislocation motion and plastic slip leading to the energy dissipation.

Simple shear deformation of a single crystal beam {#sec:pldeformation}
=================================================

![Simple shear deformation of a beam of rectangular cross section[]{data-label="fig:plane"}](planestrain){width="8cm"}

Let us consider now the simple shear deformation of a single crystal beam having only one active slip system. The crystal occupies in its initial configuration a long cylinder of an arbitrary cross section such that $(x_1,x_2)\in \mathcal{A}$ and $0\le x_3\le L$ (see Fig. \[fig:plane\] for the beam of rectangular cross section). As before, the slip system is chosen such that the vectors $\mathbf{s}$, $\mathbf{p}$, and $\mathbf{m}$ coincide with $\mathbf{e}_1$, $\mathbf{e}_2$, and $\mathbf{e}_3$, respectively. We realize the simple shear deformation by placing this crystal beam in a “hard” device with the prescribed displacements at the boundary of the crystal such that $$\label{eq:5.1}
y_1=x_1+\gamma x_3,\quad y_2=x_2,\quad y_3=x_3.$$ We assume that the length of the crystal $L$ is large enough compared with the sizes of the cross section to guarantee the uniform simple shear deformation state for which equation is valid everywhere in the beam. If the overall shear $\gamma $ is sufficiently small, then it is natural to expect that the crystal deforms elastically and the plastic slip as well as the loop function vanish. If this parameter exceeds some critical threshold, then dislocation loops may appear (see one dislocation loop in Fig. \[fig:plane\]). Due to the almost translational invariance in $x_3$-direction we may assume that $\beta (\mathbf{x})$ and $l(\mathbf{x})$ depend only on $x_1$ and $x_2$ (except perhaps the neighborhoods of $x_3=0$ and $x_3=L$).

For simplicity let us consider the small strain theory. Then the only non-zero components of the total strain tensor, under the condition that is valid everywhere, are $$\varepsilon _{13}=\varepsilon _{31}=\frac{1}{2}\gamma .$$ Since the plastic distortion tensor ${\boldsymbol{\beta}}$ has only one non-zero component $\beta _{13}=\beta (x_1,x_2)$, the non-zero components of the elastic strain tensor read $$\label{eq:5.2}
\varepsilon ^e_{13}=\varepsilon ^e_{31}=\frac{1}{2}(\gamma -\beta (x_1,x_2)).$$ With the loop function being $l(x_1,x_2)$ the dislocation densities are given by $$\label{eq:5.3}
\begin{split}
\rho_\perp  =\frac{1}{b}\left| \frac{l_{,1}(\beta _{,1}l_{,1}+\beta _{,2}l_{,2})}{l_{,1}^2+l_{,2}^2}\right| , 
\\
\rho_\parallel =\frac{1}{b}\left| \frac{l_{,2}(\beta _{,1}l_{,1}+\beta _{,2}l_{,2})}{l_{,1}^2+l_{,2}^2}\right| .
\end{split}$$

For the small strain theory we propose the free energy per unit volume of the undeformed crystal in the most simple form $$\label{eq:5.4}
\psi ({\boldsymbol{\varepsilon}}^e,\rho _\perp, \rho _\parallel ) = \frac{1}{2}\lambda ({\text{tr}}{\boldsymbol{\varepsilon}}^e)^2 + \mu {\text{tr}}({\boldsymbol{\varepsilon}}^e \cdot {\boldsymbol{\varepsilon}}^e) + \frac{1}{2}\mu k_1 \frac{\rho_\perp ^2}{\rho _s^2}+ \frac{1}{2}\mu k_2 \frac{\rho_\parallel ^2}{\rho _s^2} .$$ Here ${\boldsymbol{\varepsilon}}^e$ is the elastic strain tensor, $\lambda $ and $\mu $ are the Lamé constants, $k_1$ and $k_2$ are material constants, while $\rho _s$ can be interpreted as the saturated dislocation density. The first two terms in represent the free energy of the crystal due to the macroscopic elastic strain, where we assume that the crystal is elastically isotropic. The last two terms in correspond to the energy of the dislocation network for moderate dislocation densities [@Gurtin2002; @Gurtin2007]. Note that, for the small or extremely large dislocation densities close to the saturated value, the logarithmic energy proposed by @Berdichevsky06 [@Berdichevsky06a] turns out to be more appropriate. Substituting formulas and into we obtain $$\begin{gathered}
w(\beta ,\nabla \beta ,\nabla l)=\frac{1}{2}\mu (\gamma -\beta )^2+\frac{1}{2}\frac{\mu }{b^2\rho _s^2}\left[ k_1\frac{l_{,1}^2(\beta _{,1}l_{,1}+\beta _{,2}l_{,2})^2}{(l_{,1}^2+l_{,2}^2)^2} \right.
\\
\left. +k_2 \frac{l_{,2}^2(\beta _{,1}l_{,1}+\beta _{,2}l_{,2})^2}{(l_{,1}^2+l_{,2}^2)^2} \right] .\end{gathered}$$ Since in this case the side boundary does not allow dislocations to reach it, we can pose on both functions $\beta (x_1,x_2)$ and $l(x_1,x_2)$ the Dirichlet boundary conditions $$\label{eq:5.6}
\beta (x_1,x_2)=0,\quad l(x_1,x_2)=0 \quad \text{for $(x_1,x_2)\in \partial \mathcal{A}$},$$ The variational problem reduces to minimizing the two-dimensional functional $$\label{eq:5.7}
I[\beta (x_1,x_2),l(x_1,x_2)]=L\int_{\mathcal{A}} w(\beta ,\nabla \beta ,\nabla l)dx_1dx_2$$ among all admissible functions $\beta (x_1,x_2)$ and $l(x_1,x_2)$ satisfying the boundary conditions .

It is convenient to simplify the functional and minimize it in the dimensionless form. Introducing the dimensionless variables and quantity $$\bar{x}_1=b\rho _sx_1, \quad \bar{x}_2=b\rho _sx_2,\quad (\bar{x}_1,\bar{x}_2)\in \bar{\mathcal{A}},\quad \bar{I}=\frac{Ib^2\rho _s^2}{\mu L},$$ we can write functional in the form $$\label{eq:5.8}
I=\frac{1}{2} \int_{\mathcal{A}}\left[ (\gamma -\beta )^2+ k_1 \frac{l_{,1}^2(\beta _{,1}l_{,1}+\beta _{,2}l_{,2})^2}{(l_{,1}^2+l_{,2}^2)^2} + k_2 \frac{l_{,2}^2(\beta _{,1}l_{,1}+\beta _{,2}l_{,2})^2}{(l_{,1}^2+l_{,2}^2)^2} \right] dx_1dx_2,$$ where the bar over the quantities are dropped for short. The problem is to minimize functional among all admissible functions $\beta (x_1,x_2)$ and $l(x_1,x_2)$ satisfying Dirichlet boundary conditions . Since the integrand of is positive definite, the existence of the minimizer in this variational problem is guaranteed.

In one special case the problem degenerates and admits an analytical solution. Indeed, if we choose $k_1=k_2=k$ and take $\mathcal{A}$ to be a circular cross section whose boundary is given in the polar coordinates by $r=R$, then due to the symmetry we may assume that both $\beta $ and $l$ are functions of $r=\sqrt{x_1^2+x_2^2}$ only. It is now a simple matter to show that $$\frac{l_{,1}^2(\beta _{,1}l_{,1}+\beta _{,2}l_{,2})^2}{(l_{,1}^2+l_{,2}^2)^2} + \frac{l_{,2}^2(\beta _{,1}l_{,1}+\beta _{,2}l_{,2})^2}{(l_{,1}^2+l_{,2}^2)^2}=\beta _{,r}^2.$$ Thus, functional (normalized by $2\pi $) does not depend on $l(r)$ and takes the form $$I=\int_0^R [ \frac{1}{2} (\gamma -\beta )^2+\frac{1}{2} k \beta _{,r}^2 ] r dr,$$ which leads to the Euler equation $$(\gamma -\beta )+\frac{k}{r}(\beta _{,r}r)_{,r}=0.$$ This is nothing else but the inhomogeneous modified Bessel equation that yields the following solution (regular at $r=0$) $$\beta (r)=\gamma +C I_0(r/\sqrt{k}),$$ with $I_0(x)$ being the modified Bessel function of the first kind. The coefficient $C$ must be found from the boundary condition $\beta (R)=0$ giving $$C=-\gamma /I_0(R/\sqrt{k}).$$ Thus, $$\label{eq:5.9}
\beta (r)=\gamma [1-I_0(r/\sqrt{k})/I_0(R/\sqrt{k})].$$ The plot of solution for different values of $\gamma $ is shown in Fig. \[fig:sbeta\], where we choose the radius such that $\bar{R}=b\rho _sR=1$ and the parameter $k=10^{-4}$. As $\gamma $ increases, the plastic slip increases too. It is seen from this Figure that the plastic slip is nearly constant in the middle of the cross section and changes strongly only in the thin layer in form of ring near the boundary. Since $l$ is a function of $r$, dislocations in form of circles pile up against the boundary of the cross section, leaving the middle of the cross section almost dislocation-free. Since the dislocation loops are circles, they have the purely edge character at $\varphi =0$ and $\varphi =\pi $ and the purely screw character at $\varphi =\pm \pi /2$. For all other angles the dislocation loops have the mixed character. Note that if $k_1\ne k_2$, the strictly axi-symmetric solution is no longer valid because the contributions of the edge and screw components to the energy of dislocation network are not equal. Note also that, if the logarithmic energy is used instead of the quadratic energy, $\beta $ becomes non-zero only if $\gamma >\gamma _c$, so there is a threshold stress for the dislocation nucleation (see [@Berdichevsky-Le07]). Besides, the existence of the dislocation-free zone in the middle of the cross-section can be established.

![The plastic slip $\beta (r)$: i) $\gamma =0.001$, ii) $\gamma =0.005$, iii) $\gamma =0.01$.[]{data-label="fig:sbeta"}](Sbeta){width="8cm"}

For an arbitrary cross section and for $k_1\ne k_2$ the problem does not admit exact analytical solution. However, based on the character of solution that changes strongly only in the normal direction to the boundary observed in the previous case, we may use the asymptotic method to find the solution in the thin boundary layer. Take for example the rectangular cross-section $(x_1,x_2)\in (0,W)\times (0,H)$, with $W$ and $H$ being its width and height. In this case let us assume that there are two boundary layers near the left and right boundaries $x_1=0$ and $x_1=W$, and two other boundary layers near the bottom and top boundaries $x_2=0$ and $x_2=H$. Near the boundaries parallel to the $x_2$ axis the derivative with respect to $x_2$ can be neglect as compared to the derivative with respect to $x_1$, while the opposite is true near the boundaries parallel to $x_1$-axis. In the middle of the cross section the plastic slip remains constant. Functional reduces to the sum of four integrals, and all of them do not depend on $l$. Consider for instance the integral along the left boundary layer $$I_1[\beta (x_1)]=\int_0^\lambda [ \frac{1}{2} (\gamma -\beta )^2+\frac{1}{2}k_1\beta _{,1}^2]dx_1,$$ with $\lambda $ being still an unknown length. We omit here the integration over $x_2$ (because in this boundary layer it plays just the role of a parameter) and try to find the minimum among $\beta $. The standard variational calculus leads to the Euler equation $$\label{eq:5.11}
\gamma -\beta +k_1\beta _{,11}=0$$ which must be subjected to the boundary conditions $$\label{eq:5.12}
\beta (0)=0.$$ The solution of and that does not grow exponentially as $x_1\to \infty $ reads $$\beta (x_1)=\gamma (1-e^{-x_1/\sqrt{k_1}}).$$ Due to the mirror symmetry of the problem, the plastic slip in the boundary layer near $x_1=W$ must be $$\beta (x)=\gamma (1-e^{-(W-x_1)/\sqrt{k_1}}).$$ Similarly, the solution in the boundary layers parallel to the $x_1$ axis equals $$\label{eq:5.14}
\beta (x_2)=\begin{cases}
   \gamma (1-e^{-x_2/\sqrt{k_2}})   & \text{near $x_2=0$}, \\
    \gamma (1-e^{-(H-x_2)/\sqrt{k_2}})   & \text{near $x_2=H$}.
\end{cases}$$ Thus, the width of the boundary layers parallel to the $x_1$ axis must be of the order $\sqrt{k_2}$, while that of the boundary layers parallel to the $x_2$-axis must be of the order $\sqrt{k_1}$. Since $l$ depends only on the normal coordinate to the boundary, the dislocation lines must be parallel to the boundary of the cross section except at four corners where the asymptotic solution becomes no longer valid (the numerical solution should lead to a smoothing of the corners of dislocation loops). We see that near the vertical boundaries the dislocation loops have the edge character, while near the horizontal boundaries they have the screw character. This agrees well with the widths of the boundary layers determined by the corresponding contributions of the edge and screw components to the energy of the dislocation network.

For very thin rectangular cross section with $H\ll W$ we may neglect the influence of the edges near $x_1= 0$ and $x_1=W$ by considering the dislocation network in the central part of the beam. In this case we will have only screw dislocations which pile up against two obstacle at $x_2=0$ and $x_2=H$. The solution reduces to that found in [@Berdichevsky-Le07].

Conclusion
==========

In this paper we have developed the nonlinear CDT for crystals containing curved dislocations based on the LEDS-hypothesis. The completely new set of equilibrium equations, boundary conditions and constitutive equations have been derived from the principle of minimum free energy. As the outcome, we have obtained the system of strongly nonlinear partial differential equations for the placement, the plastic slip, and the loop function. In the case of non-vanishing resistance to dislocation motion we have derived the governing equations from the variational equation that takes the dissipation into account. We have extended the theory to the case of multiple slip and simplified it for small strains. The application of the theory has been illustrated on the problem of single crystal beam having one slip system and deforming in simple shear. Under the simplified assumption $k_1=k_2$ the analytical solution of this problem has been found for the circular cross section. For arbitrary cross sections the problem has been solved by the asymptotic method. We have shown that the asymptotic solution found for the rectangular cross section reduces to the well-known solution in [@Berdichevsky-Le07] if it is thin and long.

[*Acknowledgments*]{}

The financial support by the German Science Foundation (DFG) through the research projects LE 1216/4-2 and GP01-G within the Collaborative Research Center 692 (SFB692) is gratefully acknowledged.

[42]{} Acharya, A., 2001. A model of crystal plasticity based on the theory of continuously distributed dislocations. J. Mech. Phys. Solids 49, 761-784. Acharya, A., Bassani, J.L., 2000. Lattice incompatibility and a gradient theory of crystal plasticity. J. Mech. Phys. Solids 48, 1565-1595. Berdichevsky, V.L., 2006a. Continuum theory of dislocations revisited. Continuum Mech. Therm. 18, 195-222. Berdichevsky, V.L., 2006b. On thermodynamics of crystal plasticity. Scripta Materialia 54, 711-716. Berdichevsky, V.L., Le, K.C., 2007. Dislocation nucleation and work hardening in anti-planed constrained shear. Continuum Mech. Therm. 18, 455-467. Berdichevsky, V.L., Sedov, L.I., 1967. Dynamic theory of continuously distributed dislocations. Its relation to plasticity theory. Appl. Math. Mech. (PMM) 31, 989-1006. Bilby, B.A., Gardner, L.R.T., Stroh, A.N., 1957. Continuous distributions of dislocations and the theory of plasticity, in: Extrait des actes du $IX^e$ congrès international de mécanique appliquée, pp. 35-44. Engels, P., Ma, A., Hartmaier, A., 2012. Continuum simulation of the evolution of dislocation densities during nanoindentation. Int. J. Plasticity 38, 159-169. Gurtin, M.E., 1981. An Introduction to Continuum Mechanics. Academic Press, New York. Gurtin, M. E., 2002. A gradient theory of single-crystal viscoplasticity that accounts for geometrically necessary dislocations. J. Mech. Phys. Solids 50, 5-32. Gurtin, M.E., Anand, L., Lele, S.P., 2007. Gradient single-crystal plasticity with free energy dependent on dislocation densities. J. Mech. Phys. Solids 55, 1853-1878. Hansen, N., Kuhlmann-Wilsdorf, D., 1986. Low energy dislocation structures due to unidirectional deformation at low temperatures. Mater. Sci. Eng. 81, 141-161. Hochrainer, T., Sandfeld, S., Zaiser, M., Gumbsch, P., 2014. Continuum dislocation dynamics: towards a physical theory of crystal plasticity. J. Mech. Phys. Solids 63, 167-178. Hochrainer, T., Zaiser, M., Gumbsch, P., 2007. A three-dimensional continuum theory of dislocation systems: kinematics and mean-field formulation. Phil. Mag. 87, 1261-1282. Hughes, D.A., Hansen, N., 1997. High angle boundaries formed by grain subdivision mechanisms. Acta Mater. 45, 3871-3886. Kaluza, M., Le, K.C., 2011. On torsion of a single crystal rod. Int. J. Plasticity 27, 460-469. Kochmann, D.M., Le, K.C., 2008a. Dislocation pile-ups in bicrystals within continuum dislocation theory. Int. J. Plasticity 24, 2125-2147. Kochmann, D.M., Le, K.C., 2008b. Plastic deformation of bicrystals within continuum dislocation theory. Math. Mech. Solids 14, 540-563. Kochmann, D.M., Le, K.C., 2009. A continuum model for initiation and evolution of deformation twinning. J. Mech. Phys. Solids 57, 987-1002. Koster, M., Le, K.C., Nguyen, B.D., 2015. Formation of grain boundaries in ductile single crystals at finite plastic deformations. Int. J. Plasticity 69, 134-151. Kr[ö]{}ner, E., 1992. Mikrostrukturmechanik. GAMM-Mitteilungen 15, 104-119. Kuhlmann-Wilsdorf, D., 1989. Theory of plastic deformation: - properties of low energy dislocation structures. Mater. Sci. Eng. A 113, 1-41. Kuhlmann-Wilsdorf, D., 2001. Q: Dislocations structures - how far from equilibrium? A: Very close indeed. Mater. Sci. Eng. A 315, 211-216. Kuhlmann-Wilsdorf, D., Hansen, N., 1991. Geometrically necessary, incidental and subgrain boundaries. Scripta Metall. Mater. 25, 1557-1562. Laird, C., Charsley, P., Mughrabi, H., 1986. Low energy dislocation structures produced by cyclic deformation. Mater. Sci. Eng. 81, 433-450. Le, K.C., 2010. Introduction to Micromechanics. Nova Science, New York. Le, K.C., G[ü]{}nther, C., 2014. Nonlinear continuum dislocation theory revisited. Int. J. Plasticity 53, 164-178. Le, K.C., Nguyen, B.D., 2012. Polygonization: Theory and comparison with experiments. Int. J. Eng. Sci. 59, 211-218. Le, K.C., Nguyen, B.D., 2013. On bending of single crystal beam with continuously distributed dislocations. Int. J. Plasticity 48, 152-167. Le, K.C., Sembiring, P., 2008a. Analytical solution of plane constrained shear problem for single crystals within continuum dislocation theory. Arch. Appl. Mech. 78, 587–597. Le, K.C., Sembiring, P., 2008b. Plane-constrained shear of a single crystal strip with two active slip-systems. J. Mech. Phys. Solids 56, 2541-2554. Le, K.C., Sembiring, P., 2009. Plane constrained uniaxial extension of a single crystal strip. Int. J. Plasticity 25, 1950-1969. Mayeur, J.R., McDowell, D.L., 2014. A comparison of Gurtin type and micropolar theories of generalized single crystal plasticity. Int. J. Plasticity 57, 29-51. ztop, M.S., Niordson, C.F., Kysar, J.W., 2013. Length-scale effect due to periodic variation of geometrically necessary dislocation densities. Int. J. Plasticity 41, 189-201. Ortiz, M., Repetto, E.A., 1999. Nonconvex energy minimization and dislocation structures in ductile single crystals. J. Mech. Phys. Solids 47, 397-462. Ortiz, M., Repetto, E.A., Stainier, L., 2000. A theory of subgrain dislocation structures. J. Mech. Phys. Solids 48, 2077-2114. Sandfeld, S., Hochrainer, T., Zaiser, M., Gumbsch, P., 2011. Continuum modeling of dislocation plasticity: Theory, numerical implementation, and validation by discrete dislocation simulations. Journal of Materials Research 26, 623-632. Sandfeld, S., Thawinan, E., Wieners, C., 2015. A link between microstructure evolution and macroscopic response in elasto-plasticity: formulation and numerical approximation of the higher-dimensional continuum dislocation dynamics theory. Int. J. Plasticity 72, 1-20. Sedov, L.I., 1965. Mathematical methods of constructing models of continuum media. Usp. Matem. Nauk. 20, 123-182. Wulfinghoff, S., B[ö]{}hlke, T., 2015. Gradient crystal plasticity including dislocation-based work-hardening and dislocation transport. Int. J. Plasticity 69, 152-169. Zhu, X., Dai, S., Xiang, Y., 2013. Numerical simulation of dynamics of dislocation arrays and long-range stress fields of nonplanar dislocation arrays. Int. J. Plasticity 43, 85-100. Zhu, X., Xiang, Y., 2014. Continuum framework for dislocation structure, energy and dynamics of dislocation arrays and low angle grain boundaries. J. Mech. Phys. Solids 69, 175-194.

[^1]: phone: +49 234 32-26033, email: chau.le@rub.de. The paper is dedicated to the 70th birthday of my teacher V. Berdichevsky.

[^2]: In conventional crystal plasticity the kinematic equation for $\mathbf{F}^p$ is usually formulated in rate form that does not always reduces to [@Ortiz99; @Ortiz00].
