Comparison of dislocation density tensor fields derived from discrete dislocation dynamics and crystal plasticity simulations of torsion

The importance of accurate simulation of the plastic deformation of ductile metals to the design of structures and components is well-known. Many techniques exist that address the length scales relevant to deformation pro- cesses, including dislocation dynamics (DD), which models the interaction and evolution of discrete dislocation line segments, and crystal plasticity (CP), which incorporates the crystalline nature and restricted motion of dis- locations into a higher scale continuous field framework. While these two methods are conceptually related, there have been only nominal efforts focused on the system-level material response that use DD-generated information to enhance the fidelity of plasticity models. To ascertain to what degree the predictions of CP are consistent with those of DD, we compare their global and microstructural response in a number of deformation modes. After using nominally homogeneous compression and shear deformation dislocation dynamics simulations to calibrate crystal plasticity flow rule parameters, we compare not only the system-level stress-strain response of prismatic wires in torsion but also the resulting geometrically necessary dislocation density fields. To establish a connection between explicit description of dislocations and the continuum assumed with crystal plasticity simulations, we ascertain the minimum length-scale at which meaningful dislocation density fields appear. Our results show that, for the case of torsion, the two material models can produce comparable spatial dislocation density distributions.


Introduction
The importance of accurate simulation of the plastic deformation of ductile metals in the design of structures and components to performance and failure criteria is well known. Plasticity models describe the influence of elastic and inelastic deformation on stress within a body that undergoes a specific displacement/loading path. These models are conventionally constructed with parameters, such as elastic constants, yield stress and work hardening, fitted to experimentally measured data.
Crystal plasticity (CP) is a particular form of a plasticity model that takes into account some details of the underlying crystal structure. In CP, characteristics of the motion of dislocations (nanometer scale line imperfections in the crystal structure of a material) influence the deformation experienced by the material. Early CP models connected plastic strain rate to Burgers vectors, velocity and line length density of inherent material dislocations [1][2][3]. More recent models cast this response in a finite deformation framework, and decompose the plastic portion of the velocity gradient into separate contributions from various families of dislocations, each associated with a specific slip system [4][5][6].
Given the nanoscale nature of the underlying defects that govern plastic deformation in metals, it is logical to consider whether simulation methods that explicitly model these defects, such as atomistic simulation and discrete dislocation dynamics (DD) see e.g. [7] and [8], can be used to improve the fidelity of CP constitutive relations. The connection between dislocation interaction and plasticity was pioneered by [9], [10], and [11,12]. Since that time, a number of researchers have made connections between dislocation mechanisms and crystal plasticity, see e.g. [13][14][15].
A number of recent efforts have made comparisons of molecular dynamics results with continuum plasticity simulations. [16] compared atomistic simulation of FCC nickel with both CP and macroscopic internal state variable theories for the case of simple shear. This comparison revealed qualitative similarities between the different methods, but also quantitative differences due to the presence of thermal vibrations that occur for atomistic simulation at room temperature. In particular, it was observed that a much narrower stress distribution arose for the finite element analyses than for that extracted from atomistics. This connection to thermal vibration was verified through additional simulations at low temperature (∼10K). Efforts have also been made by Vitek [17,18], Gröger et al. [19][20][21], and Weinberger et al. [22][23][24] to directly use atomistic simulation results to parameterize microscale yield response and CP models in body-centered cubic (BCC) metals. Vitek and coworkers examined the non-Schmidt behavior of the slip systems of group VB and VIB systems using a Finnis-Sinclair inter-atomic potential. 1 In their work, Gröger et al. used bond order potentials to model the behavior of screw dislocations in Mo and W, quantifying a stress/loading-direction dependency to the critical resolved shear stress and characterizing the resulting non-Schmid behavior. Weinberger et al. leveraged this approach to fit parameters of a single crystal yield constitutive model to zero-temperature molecular statics simulation results. They then combined this parameterization with experimental data to develop a temperature and stress dependent model for the evolving plastic strain rate. This approach is particularly well-suited for BCC metals, where lattice friction for individual dislocations represents a significant barrier to glide.
Dislocation dynamics has also been recently employed to make comparisons between continuum plasticity and lower-level simulation methods. [25] built on the insights of Vitek et al. and Gröger et al. in developing a threedimensional discrete DD model to characterize the relationship between dislocation glide behavior and macro-scale plastic slip in single crystal BCC Ta. For FCC metals, where interactions between dislocations and with other obstacles, e.g. precipitates, strongly dictate the speed of dislocation motion and resulting rate of plastic deformation, atomistic simulation is limited in its ability to simulate a significant quantity of dislocations. In this regime, DD can play a useful role in providing relevant characteristics used in the fitting of CP model parameters. Such a connection was exploited by [26], who employed simple boundary value problems to compare the predictions of nonlocal plasticity theories with DD results. [27] later used both atomistic simulation and DD in a hierarchy of modeling methods, with results from atomistics used to quantify individual dislocation mobilities for use in DD simulations, and DD results to quantify parameters related to work hardening within a CP framework. These researchers successfully used this combined approach to predict the mechanical response of an aluminum single crystal deformed under uniaxial compressive loading along the [421] crystal direction. The computed strain-stress response agreed well with experimental data that was not used in model calibration.
In addition to MD and DD, other methods such as phase field crystal methods using diffusive dynamics have also proven useful in simulating plastic behavior with near atomic detail, see, e.g. [28] and [29].
While these efforts have achieved some degree of success, they have been limited in how they use the data from the evolution of explicit dislocations to enhance the fidelity of the plasticity models considered. In general, only the global stress-strain response of a system containing a distribution of dislocations is used to parameterize the hardening behavior; no local information is used to connect the two models. In particular, a direct comparison of the dislocation density tensor field (DDT) resulting from the evolution of the dislocation distribution in the case of DD and from the plastic deformation in the case of CP has not been considered thus far. However, the inclusion of the DDT as an independent field or a parameter in the flow rules and yield criteria of single-crystal plasticity models has been explored, see [30][31][32][33][34][35], and can be used to model length-scale effects [36,37], predict different dislocation microstructures [38], and allows comparison with EBSD and X-ray micro-diffraction data [39][40][41][42][43].
A comparison of the dislocation density tensor field resulting from corresponding DD and CP simulations is the primary contribution of this work. An attendant development is the means of extracting the DDT field from both methods and the associated length-scale at which one may resolve the DD-based densities. After a review of the basic theory and the two simulation methods in Sec. 2 and Sec. 3, we compare corresponding simulations resulting from both methods in Sec. 4. To make the comparison meaningful we make the two models as consistent as possible given the different formulations. First, by using nominally homogeneous compression and shear loading simulations, we calibrate a simple, representative flow rule for the CP model based on the DD response in Sec. 4.1. To eliminate dynamic effects we consider loading rates where both the DD and CP results are relatively insensitive to changes of rate. Another key predicate for consistency is the appropriate scale for the comparison based on the averaging volume used to translate discrete dislocations in DD into a density field. We hypothesize that there exists a range of volumes where a field is resolved and remains practically indistinguishable for small changes in the averaging volume. This notion of an intermediate asymptotic scale, refer to [44], for the dislocation density field is crucial to this work. Then, after parameter calibration and estimation of minimal asymptotic scale, we compare the spatially-varying DDT field, as well as the global response, resulting from DD and CP simulations (Sec. 4.2). Such comparison is performed for torsion of wires with square cross section. The rationale behind this choice is that, on the one hand, torsion induces a deformation state which is inhomogeneous within the cross section, and therefore its plastic component can be associated with geometrically necessary dislocation microstructures. On the other hand, the deformation state is homogeneous along the twist axis, therefore allowing meaningful averages along this direction for all fields extracted from DD and CP. This comparison of inhomogeneous loading in DD is enabled by method of [45] and the MODEL code [46], which treats the boundary effects of a finite system with interacting dislocations. (Other efforts along these lines can be found in [47].) We conclude with a discussion of the results and insights on how aspects of CP models can be improved using DD simulation in the final section. To our knowledge, comparison of DD and CP response at the levels of both overall stress-strain response and microstructural fields has not been investigated before.

Theory
The dislocation density tensor (DDT, α) introduced by [48] is a measure of the geometrically necessary 2 dislocation density in a region of a crystal. The DDT is a fundamental quantity connecting the individual dislocations at the micro-scale to the macroscopic elastic and plastic deformation gradient fields, which are the essential kinematic ingredients in the macroscopic phenomenological theory of plasticity, see, e.g., [49].
The dislocation density tensor α is defined by its relationship with the (total) Burgers vector b for oriented area where κ 0 is the dislocation-free reference configuration and κ is the intermediate (or lattice) configuration where the loop closure b for area A is measured. 3 This definition relies on Stokes theorem and the usual multiplicative decomposition, originally due to [51], of the deformation gradient F, where the deformation gradient F = ∇ X χ is the derivative of the motion x = χ(X, t) with respect to the reference coordinates X, refer to Fig. 1. As denoted in Eq. (1), the DDT α, being the curl of the plastic deformation maps (oriented) elements of area dA in the reference configuration κ 0 to infinitesimal Burgers vectors db in the incompatible lattice configuration. The primary kinematic assumption of crystal plasticity is the rate of plastic deformation is of the form: whereγ a is the slip rate on plane a and the Schmid tensors P a = s a ⊗ n a reside in the intermediate configuration and are associated with the crystallographic planes in which glissile dislocations reside. (Here, s a is the slip direction and n a is the slip plane normal.) This leads to an overall deformation rate given by in the finite deformation context. Fig. 1 shows one such slip plane for an edge dislocation with the slip direction s aligned with the Burgers vector b. In truth, Fig. 1 is a very much simplified picture of plasticity with the one dislocation shown representing a family of parallel dislocation lines as Nye envisaged; however, the essential irreversibility of macro-plasticity is a manifestation of both dislocation motion and entanglement. Under the assumption of linearized kinematics (x ≈ X, small strain and rotation), the total deformation gradient can be written as 2 The geometrically necessary portion of the total dislocation density is the part ensuring compatibility of the overall deformation from reference lattice to current configuration, as opposed to the statistically stored portion. 3 The adjugate F * = det(F)F −T is defined by Nanson's formula da = F * dA. As in [50], the curl of a tensor field F is defined where β e = F e − I and β p = F p − I are the elastic and plastic distortions, respectively. In particular, the process of imposing a relative displacement b (x) to the two sides of a cut in surface A (refer to Fig. 1) results in the plastic distortion tensor, as in [52] If the relative displacement is constant over A , the defect is a Volterra dislocation and it models a crystal dislocation with Burgers vector b . In this case, the dislocation density tensor is concentrated on the closed line bounding A : where Particularly germane to the present study is volume average of the dislocation density tensor (8). The average of Eq. (8) over a region Ω with volume V centered at y can easily be obtained by virtue of the sifting property of the Dirac-δ, that is where the line integral now extends over the portion of the dislocation line inside Ω. Here,α(z) is at the scale of dislocations embedded in lattice and the field α(y) is at a larger length scale where continuous fields can be obtained. This definition classifies dislocation lines into geometrically necessary or statistically stored depending on whether the line is open or closed within the volume Ω. In previous studies [53] we derived a generalization of Eq. (9), which uses bell-shaped smoothing kernels ∆ (x − y) that are not necessarily piecewise constant. This definition is evaluated using a quadrature based on straight segments x ab = x a − x b of length ab = ds dλ with tangents ξ ab = x ab / ab . In a similar fashion the local line density ρ(y) can be estimated: Given that the kernel function has compact support, we can use its radius of influence to estimate the range of length scales for the estimated dislocation density to be smooth and yet resolved the variations in the system. In fact, we will estimate the minimum of the range of intermediate asymptotic scales in order to maximize the resolution of trends in the dislocation density tensor field. In averaging via Eq. (10) to find the minimum number of dislocation lines in a given volume to obtain smooth density field and to compare with CP results, we are also making the conjecture that this is the scale of CP. It bears mentioning that in our previous work [53] we demonstrated that the definition of α, Eq. (9), satisfies the fundamental extensive property of the Burgers vector, namely that the resultant Burgers vector of a Burgers circuit is the sum of Burgers vectors of all individual dislocation lines crossing the enclosed surface. We also showed that Nye's original definition of dislocation density and Eq. (9) become equivalent when certain conditions regarding the spacing between the dislocations and curvature of dislocations are obeyed. Consistency, for example whereρ is the average density, implies the requirement that given (a) the localization ∆(y I ) = 1 V I N I is based on a partition of unity N I (x) | I N I (x) = 1 ∀ x, (b) exact integration of the line integral, and (c) the integration weights equal to kernel volumes w I = V I . For a piece-wise constant basis N I with compact support on a cubical region, exact integration is simply given the two intersection points x 0|1 = λ 0|1 x ab + x b of the segment with the boundary of the cube where λ 0 < λ 1 ∈ [0, 1]. If either end of the segment does not intersect the cube and instead is in the interior of the cube then λ 0|1 = 0|1.

Method
We employed the discrete dislocation dynamics (DD) method developed by [45], as implemented in the computer program MODEL [46]. A brief explanation of the method is given here; a full account of the method can be found in [45]. The crystal plasticity model was implemented in the finite deformation code Albany [54]. Given the aim of this work to compare and correlate DD and CP, we have maximized the commonalities of the two modeling techniques. For both models, we use: (a) isotropic elasticity with the slip planes of a cubic lattice, 4 (b) acceleration is ignored in the balance of linear momentum, however there are dynamic/viscous effects in each model, (c) the momentum balance is solved on finite element meshes with corresponding boundary conditions. In this mode, both the response of CP, governed by its flow rule, and that of DD, governed by its mobility equation, are relatively insensitive to loading rate.
The formula (10) allows to construct a continuum dislocation density field from DD simulations, which can be directly compared to the corresponding quantity defined in Sec. 3.2 in the CP framework.

Dislocation dynamics
In general, DD approaches treat dislocation lines explicitly, with lines decomposed into discrete edges (segments) bounded by vertices. Given the formulas in Sec. 2, the dislocation density field α(x) can be computed for a network of dislocations; however, their evolution from a given initial configuration must be computed from a closure equation yielding the dislocation velocity w as a function of the local stress state. Such an equation of motion for the discrete dislocation configuration can be obtained from thermodynamic considerations [45]. For FCC metals and by neglecting inertia terms, the equation of motion can be written as: Eq. (15) expresses the balance of the total configurational force per unit line of dislocation, and is comprised of: a viscous drag contribution, 5 −Bw, and a mechanical contribution, (σb) × ξ, the Peach-Koehler force. In contrast to CP, the assumption of linearized kinematics adopted in DD allows to write the (small-strain) stress field σ appearing in the Peach-Koehler force as the sum of two terms: (a) a dislocation-dislocation interaction term which can be computed accurately using the Peach-Koehler stress equation [58] and (b) a correction term which accounts for the mechanical boundary conditions imposed to the simulation domain.
The DD formulation used in this work solves Eq. (15) in weak/variational form for the dislocation velocity field using the finite element method, refer to [45]. The numerical implementation is based on discretization of dislocation lines into segments connecting pairs of nodes. Segment shape functions and nodal degrees of freedom are then used for the unknown dislocation velocity field w, so that Eq. (15), upon assembly over the dislocation network, is transformed into a discrete system of equations for the nodal velocities. Finally, the nodal velocities are used to update the nodal positions and evolve the dislocation configuration in time.

Crystal plasticity
In our CP model we use Kirchhoff-St. Venant elasticity S = CE e with elastic Lagrange strain E e = 1 2 F T e F e − I . This elastic constitutive relation reduces to σ = Cε in the limit of small strain and rotation. The elastic modulus C with cubic symmetry is where e ij is the dyad e i ⊗ e j , which we reduce to isotropy by setting C 44 = G, C 11 = 2G 1−ν 1−2ν , C 12 = 2G 2ν 1−2ν for correspondence with the DD simulations.
For the plastic behavior, the kinematic assumption in Eq. (4) can be seen as a plastic flow rule of the form found in [59,Ex. 2.6]. A particularly simple hardening rule for the rate of slipγ a for slip system a iṡ subject to a τ crit a threshold and (isotropic) hardening modulus H a . Here τ a = P a · σ is the shear stress resolved on the slip plane with normal n a and the Cauchy stress is simply σ = 1 det F FS(E e )F T . The implicit update algorithm largely follows [59].
Given a solution to momentum balance, we recover a dislocation density field by first using a local, element-wise L 2 projection to move F p data at the integration points to the element nodes. Then using the finite element basis 5 In general, the drag coefficient B is a function of velocity itself, i.e. the relation between drag force and dislocation velocity is non-linear; however, for the simulations presented here we assume that B is a constant material coefficient. Also of note is the fact that the drag coefficient is in general a second-order tensor exhibiting anisotropy depending on the crystal structure, and that it may also depend on the local dislocation character (i.e. edge versus screw), and on the type of dislocation motion (i.e. glide versus climb). As described in [45], the method used here utilizes Lagrange multipliers to constrain climb processes, which are often only active at very high temperatures. Approaches that do incorporate climb include the works by [55], [56] and [57].
functions N I to interpolate F p on each element, which allows differentiation with respect the reference coordinates X, we obtain: where I indexes nodes, is the permutation tensor, e and E are basis vectors in the spatial and reference frames (respectively), and we have dropped the subscript p for clarity. Chain rule provides the spatial gradient N I,B in terms of the partial derivatives of the basis N I with respect to the local coordinates ζ a where ∂ ζ X is the Jacobian of the local ζ-to-referential X coordinate map. Lastly, we average the local node values to obtain a unique value for every node and hence a continuous plastic deformation curl field. Alternately a global L 2 projection could be used.

Results
For our comparison studies we chose face-centered cubic single crystals of Cu. Its 12 slip systems are defined by the Schmid tensors {s a ⊗ n a } = { 110 ⊗ {111}}. As mentioned, we employ an isotropic elasticity model with effective Poisson ratio ν = 0.34 and shear modulus G = 48 GPa. Given that MODEL is a dimensionless code, we adopt the same normalization in our results. Lengths are normalized to the material's Burgers vector b = |b| = 2.556Å, stress-like quantities are normalized by G, the drag coefficient is set at B = 10 −4 Pa-s, and time is normalized by τ = B G = 2.08 fs. After preliminary studies, over a wide range initial (fully gissile) dislocation densities and loading rates, we chose an initial dislocation density of ρ 0 = 10 13 m −2 = 6.533 × 10 −7 b −2 , consisting of prismatic dislocation loops with 110 -type Burgers vectors lying in {111} planes, and loading rate of 10 −11 dimensionless inverse-time units (corresponding to 4.8 × 10 3 sec −1 ) to approximate the rate-independent loading regime. From the initial dislocation density we employ a length-scale = 1 √ ρ0 = 1.29 × 10 3 b to non-dimensionalize α fields. Given that initially the dislocation network consists entirely of loops, the total dislocation density tensor starts at zero. 6 Note, Eq. (15) introduces a time-scale to DD and Eq. (17) introduces time-scale to CP, yet neither has inertia so scaling time can lead to self-similar solutions dependent only on loading and not the rate of loading. For DD, this limit is achieved when the dislocation interactions and collisions are sufficiently resolved. For CP, time-scale and loading rate invariance is apparent only as m → ∞ and also for small strains with a linear flow rule, as in [60] and [61]. We examine the response of DD and CP for strains up to 1%, well below the limit of validity of the linear elastic model in DD.
In Sec. 4.1 we use homogeneous compression and shear loading results from DD to calibrate the CP flow rule and to estimate a minimum asymptotic scale. Then, in Sec. 4.2 we compare the dislocation density fields and global response resulting from inhomogeneous loading due to torsion. For all comparisons an ensemble of 10 replicas of the DD systems are used to alleviate sensitivity to initial dislocation arrangement.

Calibration with homogeneous loading simulations
With the given elastic constants and symmetries, we use compression and (simple) shear loading DD and CP simulations to calibrate the hardening parameters of the flow rule specified in Eq. (17). For these nominally homogeneous loading cases we employed a cube of dimensions 2000b × 2000b × 2000b oriented with the 100 crystal directions. In compression, the top face of the cube (the face with outward normal in the +e 3 direction) is displaced in e 3 and the bottom face x 3 (face with normal in the −e 3 direction) is fixed; the +x 1 and +x 2 lateral faces are constrained in their normal directions to promote homogeneous deformation and allow for expansion. In shear, the top face is displaced in x 1 while fixed in x 2 and x 3 ; the bottom face is fixed in x 1 , x 2 and x 3 and all lateral faces are fixed in x 2 and x 3 . Fig. 2 shows the stress-strain response for compression. The initial slope of the DD response, 2.54G, compares well with the expected Young's modulus of 2(1 + ν)G= 2.68G. Likewise, Fig. 3 shows the stress-strain response for simple shear loading. In this case, the initial slope is 0.92G near the expected value 1.00G. These results imply that the initial distribution of dislocations minimally affects the material's elastic response. After approximately 0.001 strain in compression and 0.004 strain in shear the total line density starts to rise rapidly which corresponds to a distinct decrease in the slope of the stress-strain response. Using direct fits to the DD response we determined that the slope of this plastic regime was 0.12G in the compression case and 0.42G in the shear case and that the onset of plastic response was at 0.00125 and 0.0031 strain, respectively. These values were used as initial guesses in calibrating the CP model. Fig. 2 and 3 also show that the calibrated CP stress-strain response compares well to that of DD using the hardening parameters: C a =1 (inverse time units), τ crit a = 0.11 GPa = 0.0022 G, and H a = 30.5 GPa = 0.635 G, for all slip systems a ∈ [1,12]. We employed the exponent m=4 which was the highest value that offered numerical stability. (The issue of instability upon approaching the rate independent limit m → ∞ is well-known, as discussed by [60], [61] and [62].) Lastly, by using a kernel estimator, Eq. (10), positioned near the center of a 4000b × 4000b × 4000b system we were able to estimate the asymptotic region for the dislocation density tensor field using the compression case. Fig.  4 shows that three relevant components of α converge to near constant, non-trivial values for kernel radii above ≈ 500 b, which is also well below the system size used for this study.

Comparison of inhomogeneous loading simulations
For comparing the predictions of the DD and CP models we simulated the torsion of a L 1 ×L 2 ×L 3 = 2000b×2000b× 4000b rectangular wire. The torsion was effected by a twisting displacement of the top relative to the bottom faces, with the bottom face fully fixed and the lateral faces traction-free. Both 100 and 110 crystallographic orientations of the lateral faces of the wire are explored. Since the elastic response is isotropic, only slip plane orientation relative to the loading distinguishes these two cases. The general solution for the stress and displacement fields in an isotropic elastic wire due to torsion are: and in terms of a warping function w = w(x 1 , x 2 ). Given the warping function characteristic of our prismatic wire with square cross-section, the torsional modulus K relating the torsional moment M e 3 = x × σ dA to the twist per axial distance θ/L 3 is where L 1 = L 2 is the side length of the square cross-section. The torque-twist response for both models are shown in Fig. 5. The correspondence between the calibrated CP model and the DD response for the 100 orientation of the wire is similar to that shown in Fig. 3 for shear. For the 110 oriented wire, however, the DD response shows a delayed onset of significant plasticity, as evidenced both by the torque-twist response as well as the line density evolution, whereas the CP torque-twist for this orientation is nearly the same as for 100 . Fig. 6 shows the evolution of the dislocation network with twist for the 100 case. The dislocations are primarily screw dislocations on different glide planes. They reach their shown locations after a series of cross slip events from a few initial common sources. Clearly, between twist θ = 0.02 radians and 0.04 radians, the initial random network multiplied through Frank-Read mechanisms (no nucleation model is employed) and aligned with the diagonals of the square cross-section. This regime also corresponds to the apparent onset of plasticity shown in Fig. 5. Fig. 7a shows the dislocation network and resulting line density at distinctly plastic state, twist angle θ = 0.06 radians for both the 100 and 110 cases. In the 100 case, dislocation lines spanning the diagonals of the wire's square cross-section pile up; whereas, in the 110 case shown, the majority of the lines are vertical and span the opposing face mid-lines while a minority span one pair of the diagonal corners as in the 100 case. Clearly from the line networks themselves and the resulting line density, the dislocation lines in the 110 case are more dispersed spatially and less numerous in total than in the 100 case. (Note these networks are still dense relative to the nominally homogeneous cases where there is no tendency for dislocations to accumulate in any particular region of the body.) In DD, dislocation lines are generated in high stress regions near the corners of the wire in torsion and migrate to the low/zero stress basins shown in green in Fig. 7b, given the Peach-Kohler forces in Eq. (15). For the two cases, the orientations of the slip planes relative to the wire are different, but elastic stress field is same, up until yield, due to the isotropic moduli assumption. The resolved shear stress (RSS) is plotted in Fig. 7b for both cases. In the 100 case the diagonal basin in RSS is dominant and in the 110 case the vertical basin is dominant but broader than in the 100 .
Since our system size and total line density are at the threshold of being able to be resolved given the estimated minimum asymptotic scale, we average the DD data through the x 3 -axis since slices of the solution are nominally equivalent along the twist axis. Fig. 8 and Fig. 9 show the two dimensional projection of the in-plane components of the dislocation density tensor for the 100 and 110 orientations of the wire, respectively. In general, the α fields generated with DD and CP are only comparable up to a divergence-free/solenoidal field i.e.  scalar function does not contribute to α (which follows directly from a Helmholtz decomposition of F p ). Also, the solutions have nominal 90 degree rotation and inversion symmetry. The fields resulting from both DD and CP have comparable magnitudes. And yet, the α fields generated by DD have significant components in the center of wire, unlike CP. In the 100 orientation case, the magnitudes of α components are similar; but in the 110 case the weak pattern correspondence breaks down due to more dispersed nature of dislocation line density. As mentioned, in DD lines are generated at high stress regions at the lateral edges of the prismatic wire and migrate to the RSS basins either diagonal or crossing the midpoints of the wire's cross-section. In contrast, in a CP simulation, a background of mobile, geometrically necessary dislocations are assumed to be present everywhere and plastic slip occurs in regions of high stress. In CP there is no mechanism for explicit generation of dislocations nor compatibility of plastic slip induced by the motion of the implicit dislocation network. Hence, the dislocation density field produced by CP only responds to the regions of high RSS, and there is no driving force for migration/flow of dislocation density to low RSS basins.  As a last comparison, Fig. 10 shows the sensitivity of the dislocation density field for the 100 case predicted by CP to the flow rule (17) exponent m. Clearly the significant dislocation density approaches the center of the crosssection and the "×" pattern becomes more prominent and sharper with increased m showing better correspondence with the DD results. An explanation for this better correspondence is discussed in the concluding section. Due to numerical issues with the form the simple flow rule (17) it was not possible to increase m to further improve the correspondence with the DD results.

Discussion
In this work, we have established that a simple CP model calibrated to DD data can produce similar microstructural fields as well as overall response. This similarity is predicated on resolving the dislocation density field at a scale above the intermediate asymptotic length scale intrinsic to the material with its native dislocation population. For the material we examined, Cu with line density ρ = 10 13 m −2 , we established that the minimum asymptotic scale is approximately 500 b. Physically, the asymptotic scale is determined by how densely dislocation lines pack after being generated at high RSS regions and migrating to low RSS regions. In particular, we estimated the scale at which continuum dislocation density fields should appear with nominally homogeneous simulations where there is little difference between high and low RSS regions. This estimate proved robust in constructing dislocation density fields in the inhomogeneous loading cases we considered. For the finite systems we examined, it appears that the dislocation pile-ups on particular slip planes reach a characteristic spacing due to the fact that forces between dislocations scale as b 2 , while forces due to the applied stress scale as b (refer to Chap. 21 in [63]). 7 Hence, we conjecture that the asymptotic scale for the dislocation density is strongly dependent on material and less dependent on the particular loading case. By examining Eq. (10) and Eq. (11), the connection between the dislocation line density ρ and dislocation density tensor α can be seen. Clearly ρ indicates the general magnitude of α, i.e. there are no α contributions where ρ is zero; but, α adds vectorial information from the line direction as well as from the Burgers vector, see Eq. (10) In contrast, with a CP model the dislocation density α results from the curl of the plastic deformation gradient through the CP flow rule, Eq. (17). So, unlike the explicit representation of dislocation lines in DD, CP, through its flow rule, conflates generation and migration and only creates mobile slip in high stress/deformation gradient regions. Moreover it simulates the hardening due to entanglement of dislocation lines with the flow rule.
A loose correspondence can be made between the CP flow rule and the dynamics of the dislocation lines that governs DD evolution. The Peach-Kohler forces driving the dislocation dynamics are analogous to the resolved shear stresses driving plastic evolution in the flow rule so that the dislocation line velocity is proportional to the Peach-Kohler force, and the plastic slip rate is proportional to the local resolved shear stress. From our results, is it clear that, although the global response of corresponding DD and CP models can be made the same, the microstructural dislocation density tensor field is much more sensitive to details of the calibrated flow rule. Although algorithmic stability issues prevented us from exploring the response of flow rules with high exponents, we were able to show that increased exponent tends to improve the qualitative similarity of the in-plane α fields for a prismatic wire in torsion where there is a strong tendency for the dislocation network to align with geometric features. It stands to reason that this improvement will saturate before the m → ∞ limit is reached, where plastic slip, and by analogy dislocation velocity, is insensitive to RSS magnitude. In fact, it is hard to imagine how a flow rule could be constructed to perfectly represent the α patterns seen with DD unless some compatibility of the plastic slip representing dislocation motion is enforced. This modeling direction is currently being developed within the framework of continuum dislocation theory, as in [64], [65] and [66], where dislocation densities are independent field quantities evolving according to transport equations, and plastic slip is computed from dislocation fluxes.
In general, the continuum field, α given by Eqs. (3) and (10), captures the (geometrically necessary) dislocation density and is a quantitative measure useful in understanding the microstructures that result from large strains in the plastic regime. This tensor can be used in the yield function of macroscopic phenomenological theories of plasticity to reproduce realistic dislocation microstructures. Such models form alternatives to the crystal-plasticity formulations based on active slip systems. For example, [38] include the DDT in the yield function, and develop a rate-independent macroscopic theory based on material symmetry. As the authors demonstrate, inclusion of DDT in the yield function introduces an additional length-scale into the macroscopic theory. Using this length-scale as a parameter, [38] have shown the formation of dislocation microstructures with spatially localized DDT similar to our work. For similar structures, the value of the length-scale in the theory are also similar to the range of intermediate asymptotic length scales that we estimated are needed to obtain a field description of the DDT from the DD simulations. This suggests that performing DD simulations to evaluate the DDT tensor may serve as a practical tool for calibrating and validating the phenomenological theories of plasticity. Also the correspondence between DD and CP could perhaps be improved by using a flow rule with explicit dependence on the DDT computed from the plastic deformation gradient. An analysis along these lines is left for future work.