Module 1: MATLAB Project - mathematical modeling of plasticity

Project Brief

In this project you will write MATLAB code defining several functions from the plastic behavior of materials, and use curve-fitting approaches to extract parameters from experimental measurements.

Successful completion will demonstrate competence in MATLAB data manipulation, analysis, and visualization. These skills are invaluable in analyzing the output of scientific computing software.


You should submit your scripts by creating a subdirectory called
and copying your code into that directory by 11:59pm on 4 April 2022. Late submissions will not be accepted; let me know in advance if you will have difficulty with completion..

I will give you feedback on the expectations listed below and for the overall script. Plots possessing unlabeled axes (variable and units!), illegibly small text, or inappropriately scaled axes will be penalized.

Specific Expectations

  1. Multiaxial loading. Two of the equations to convert a uniaxial yield stress \(\sigma_\text{Y}\) into the yield “surface” for a general stress state are the Tresca criterion, \[\max{\left|\sigma_i - \sigma_j\right|} \ge \sigma_\text{Y},\] and the von Mises criterion, \[\left(\sigma_1-\sigma_2\right)^2 + \left(\sigma_2-\sigma_3\right)^2 + \left(\sigma_3-\sigma_1\right)^2 \ge 2\sigma_\text{Y}^2\] where yield occurs when the criterion is true. Make two three dimensional plots showing the yield surface (where the inequality is a strict equality) in \(\sigma_1\), \(\sigma_2\), \(\sigma_3\): one plot for the Tresca criterion (use \(\sigma_\text{Y}\)=100 MPa) and one plot for the von Mises criterion (use \(\sigma_\text{Y}\) = 100 MPa). Comment on the differences between the two surfaces.

  2. Solute/dislocation interaction. The hydrostatic pressure of an edge dislocation is given by the equation \(p = -(Kb/4\pi r)\sin\theta\) for bulk modulus \(K\), Burgers vector \(b\), distance from slip plane \(r\) and angle from slip plane \(\theta\). If a solute is introduced at some location, that changes the volume by \(\Delta V = a_0^3 \epsilon\), the interaction energy is given by \(p\Delta V\); this interaction depends on the solute position, because the pressure is position dependent. Consider a solute of magnesium in aluminum (\(G\) = 29 GPa, \(\nu\) = 0.33, \(a_0\) = 0.405 nm, and for Mg in Al, \(\epsilon\) = 0.2). If the solute is placed at \(\sqrt{3} a_0\) above the slip system, (1) make a plot of the interaction energy as the dislocation glides in its slip plane, (2) plot the interaction force in the glide plane as a function of position, and (3) determine the maximum interaction force. Comment on the size of the interaction energies and forces.

  3. Work-hardening model. Kocks ( suggested a physically-motivated model for the evolution of dislocation density with two contributions: dislocation production due to jogs, that is proportional to the inverse distance between dislocations divided by a mean-free path \(\Lambda\) of dislocation motion, and dislocation annihilation that assumes a constant density of “recovery sites,” so that the rate of annihilation is proportional to the dislocation density, \[\frac{d\rho}{d\varepsilon_\text{T}} = \frac{\rho^{1/2}}{\Lambda} - \frac{2\rho}{\varepsilon_\text{r}}\] for true strain \(\varepsilon_\text{T}\), and where \(\varepsilon_\text{r}\) (a unitless strain) is an empirical parameter governing the recovery of dislocations. Combined with the equation for work-hardening, \[\sigma_\text{T}(\rho) = \sigma_0 + \bar m\alpha Gb\left(\sqrt{\rho} - \sqrt{\rho_0}\right)\] with \(\bar m = 3.06\) (FCC-averaged Schmid factor) and \(\alpha = 0.2\) for the dislocation hardening, we can determine a true-stress / true-strain relationship for work-hardening that includes two empirical parameters.

    1. First, consider nickel (\(G\) = 125 GPa, \(b\) = 0.25 nm) which has a yield stress of 15 MPa at an initial (well-annealed) dislocation density of \(10^{10} \text{m}^{-2}\). Using \(\Lambda\) = 30 nm and \(\varepsilon_\text{r}\) = 0.7, use MATLAB to numerically integrate your dislocation density vs. strain from 0 to a true strain of 1, and plot it.
    2. Next, with your values of dislocation density, make a plot of true stress vs. true strain over the same range.
    3. Convert your expressions into a parameterized model of stress and strain, and determine the interaction parameters corresponding to the data in /class/mse404pla.