Chapter 20, Problem 16P

Soft tissue follows an exponential deformation behavior in uniaxial tension while it is in the physiologic or normal range of elongation. This can be expressed as

$\sigma =\frac{E_0}{a}\left(e^{a\epsilon }-1\right)$

where $\sigma =$ stress, $\epsilon =$ strain, and $E_0$ anda are material constants that are determined experimentally. To evaluate the two material constants, the above equation is differentiated with respect to $\epsilon$. Using the above equation establishes the fundamental relationship for soft tissue

$\frac{d\sigma }{d\epsilon }=E_0+a\sigma$

To evaluate $E_0$ and a, stress-strain data are plotted as $d\sigma /d\epsilon$ versus $\sigma$ and the intercept and slope of this plot are the two material constants, respectively.

In the following table is stress-strain data for heart chordate tendineae (small tendons used to hold heart valves closed during contraction of the heart muscle; these data are from loading the tissue, while different curves are produced on unloading).

$\sigma ,1{\text{0}}^{\text{3}}{\text{N/m}}^{\text{2}}$	87.8	96.6	176	263	351	571	834	1229	1624	2107	2678	3380	4258
$\epsilon ,{\text{ 10}}^{-\text{3}}\text{m/m}$	153	204	255	306	357	408	459	510	561	612	663	714	765

Calculate the derivative $d\sigma /d\epsilon$ using finite differences. Plot these data and eliminate the data points near the zero points that appear not to follow the straight-line relationship. The error in these data comes from the inability of the instrumentation to read the small values in this region. Perform a regression analysis of the remaining data points to determine the values of $E_0$ and a.

Plot the stress versus strain data points along with the analytic curve expressed by the first equation. This will indicate how well the analytic curve matches these data.

Many times this does not work well because the value of $E_0$ is difficult to evaluate using this technique. To solve this problem $E_0$ is not used. A data point is selected $\left(\bar{\sigma },\bar{\epsilon }\right)$ that is in the middle of the regression analysis range. These values are substituted into the first equation and a value for $E_0/a$ is determined and substituted into the first equation, which becomes

$\sigma =\left(\frac{\bar{\sigma }}{e^{a\bar{\epsilon }}-1}\right)\left(e^{a\epsilon -1}\right)$

Using this approach, experimental data that are well defined will produce a good match of the data points and the analytic curve. Use this new relationship and again plot the stress versus strain data points and this new analytic curve.

To determine

To calculate: The values of $E_0$ and a using a regression analysis on the derivatives based on finite differences of the provided data points and then obtain a function of the form $\sigma =\frac{E_0}{a}\left(e^{a\epsilon }-1\right)$ that fits the data points correctly.

Answer to Problem 16P

Solution: The function that fits the data perfectly is $\sigma =804.5591945812\left(e^{0.0035499643\epsilon }-1\right)$.

Explanation of Solution

Given Information:

The provided values of the stress $\sigma$ in ${\text{10}}^3\text{N}/{\text{m}}^2$ based on the strain $\epsilon$ in ${10}^{-3}\text{m}/\text{m}$ are:

$\sigma$, ${\text{10}}^3\text{N}/{\text{m}}^2$	87.8	96.6	176	263	351	571
$\epsilon$, ${10}^{-3}\text{m}/\text{m}$	153	204	255	306	357	408

And,

$\sigma$, ${\text{10}}^3\text{N}/{\text{m}}^2$	834	1229	1624	2107	2678	3380	4258
$\epsilon$, ${10}^{-3}\text{m}/\text{m}$	459	510	561	612	663	714	765

Formula Used:

The formula for the derivative of a function at a point c using forward difference is:

$f\text{'}\left(c\right)=\frac{f\left(c+h\right)-f\left(c\right)}{h}$

The formula for the derivative of a function at a point c using backward difference is:

$f\text{'}\left(c\right)=\frac{f\left(c\right)-f\left(c-h\right)}{h}$

The formula for the derivative of a function at a point c using central difference is:

$f\text{'}\left(c\right)=\frac{f\left(c+h\right)-f\left(c-h\right)}{2h}$

Calculation:

Consider the formula for the stress:

$\sigma =\frac{E_0}{a}\left(e^{a\epsilon }-1\right)$

Differentiate both sides with respect to $\epsilon$ to obtain:

$\frac{d\sigma }{d\epsilon }=E_0+a\sigma$

This is a linear equation.

First obtain the value of the derivative of the stress $\sigma$ at the first data point using the forward difference. Note that the data is equally spaced with a step-size of 51.

$\begin{array}{c}\sigma \text{'}\left(153\right)=\frac{\epsilon \left(204\right)-\epsilon \left(153\right)}{51}\\ =\frac{96.6-87.8}{51}\\ =0.1725\end{array}$

Now the derivative at the second data point can be obtained using the central difference. Note that the data is equally spaced with a step-size of 51.

$\begin{array}{c}\sigma \text{'}\left(204\right)=\frac{\epsilon \left(255\right)-\epsilon \left(153\right)}{51}\\ =\frac{176-87.8}{2\left(51\right)}\\ =0.8647\end{array}$

Now do the same for all the data points except the last one:

$\begin{array}{c}\sigma \text{'}\left(255\right)=\frac{263-96.6}{2\left(51\right)}\\ =1.6314\\ \sigma \text{'}\left(306\right)=\frac{351-176}{2\left(51\right)}\\ =1.7157\end{array}$

And,

$\begin{array}{c}\sigma \text{'}\left(357\right)=\frac{571-263}{2\left(51\right)}\\ =3.0196\\ \sigma \text{'}\left(408\right)=\frac{834-351}{2\left(51\right)}\\ =4.7353\end{array}$

And,

$\begin{array}{c}\sigma \text{'}\left(459\right)=\frac{1229-571}{2\left(51\right)}\\ =6.4510\\ \sigma \text{'}\left(510\right)=\frac{1624-834}{2\left(51\right)}\\ =7.7451\end{array}$

And,

$\begin{array}{c}\sigma \text{'}\left(561\right)=\frac{2107-1229}{2\left(51\right)}\\ =8.6078\\ \sigma \text{'}\left(612\right)=\frac{2678-1624}{2\left(51\right)}\\ =10.3333\end{array}$

And,

$\begin{array}{c}\sigma \text{'}\left(663\right)=\frac{3380-2107}{2\left(51\right)}\\ =12.4804\\ \sigma \text{'}\left(714\right)=\frac{4258-2678}{2\left(51\right)}\\ =15.4902\end{array}$

Now the derivative at the last data point can be obtained using the backward difference. Note that the data is equally spaced with a step-size of 51.

$\begin{array}{c}\sigma \text{'}\left(153\right)=\frac{\epsilon \left(765\right)-\epsilon \left(714\right)}{51}\\ =\frac{4258-3380}{51}\\ =17.2157\end{array}$

Hence, the derivatives have been obtained as:

$\begin{array}{ccc}\sigma & \epsilon & \sigma \text{'}\left(\epsilon \right)\\ 87.8& 153& 0.1725\\ 96.6& 204& 0.8647\\ 176& 255& 1.6314\\ 263& 306& 1.7157\\ 351& 357& 3.0196\\ 571& 408& 4.7353\\ 834& 459& 6.4510\\ 1229& 510& 7.7451\\ 1624& 561& 8.6078\\ 2107& 612& 10.3333\\ 2678& 663& 12.4804\\ 3380& 714& 15.4902\\ 4258& 765& 17.2157\end{array}$

Now perform a regression analysis on the derivative with respect to $\sigma$ and obtain a linear equation for the derivative in terms of $\sigma$.

This can be done using MS-Excel.

Step 1: Enter the data points in a blank workbook.

Step 2: Select the data points, go to insert and the select scatter plot.

This gives:

It can be seen that the first four points do not follow a linear relationship. So disregard the first four points and obtain a line of best fit for the remaining points.

This can be done using MATLAB.

Code:

clc;

clear all;

% Enter the data points

Stress=[87.896.6176263351571834122916242107267833804258];

Strain=153:51:765;

% Enter the derivative values calculated

Diff=[0.17250.86471.63141.71573.01964.73536.5417.74518.607810.333312.480415.490217.2157];

% Select the new data excluding the first four data points

Newdata=[Stress(5:end);Diff(5:end)];

% Perform a linear regression on these points

polyfit(Newdata(1, :),Newdata(2, :),1);

Output:

This implies that the linear model for the derivative is:

$\frac{d\sigma }{d\epsilon }=2.8561564180+0.0035499643\sigma$

This gives the value of a as 0.0035499643 and the value of $E_0$ as 2.8561564180.

Now, substitute this in the formula for stress to obtain:

$\begin{array}{c}\sigma =\frac{2.8561564180}{0.0035499643}\left(e^{0.0035499643\epsilon }-1\right)\\ =804.5591945812\left(e^{0.0035499643\epsilon }-1\right)\end{array}$

To check whether this model is a good fit for the data, plot this model along with the provided data points. This can be done using MATLAB.

Code:

clc;

clear all;

% Enter the data points

Stress=[87.896.6176263351571834122916242107267833804258];

Strain=153:51:765;

% Enter the derivative values calculated

Diff=[0.17250.86471.63141.71573.01964.73536.5417.74518.607810.333312.480415.490217.2157];

% Select the new data excluding the first four data points

Newdata=[Stress(5:end);Diff(5:end)];

% Perform a linear regression on these points

polyfit(Newdata(1, :),Newdata(2, :),1);

% Enter the formula for Stress obtained using the values of E0 and a

eqn=804.5591945812*(exp(0.0035499643*Strain)-1);

% Now plot the points along with the data

scatter(Strain,Stress);hold on

plot(Strain,eqn)

Output:

Clearly, the obtained function is not a good fit for the data.

Rewrite the expression for the original functionto obtain:

$\frac{E_0}{a}=\frac{\sigma }{\left(e^{a\epsilon }-1\right)}$

Use this to obtain:

$\sigma =\frac{\sigma \text{'}}{e^{a\epsilon \text{'}}-1}\left(e^{a\epsilon }-1\right)$

Now consider a data point in the middle of the data $\left(\sigma \text{'},\epsilon \text{'}\right)=\left(2107,612\right)$.

Use this to get the formula for the strain as:

$\begin{array}{c}\sigma =\frac{2107}{e^{0.0035499643\left(612\right)}-1}\left(e^{0.0035499643\epsilon }-1\right)\\ =\frac{2107}{e^{0.0035499643\left(612\right)}-1}\left(e^{0.0035499643\epsilon }-1\right)\\ =270.79152781\left(e^{0.0035499643\epsilon }-1\right)\end{array}$

To check whether this model is a good fit for the data, plot this model along with the provided data points. This can be done using MATLAB.

Code:

clc;

clear all;

% Enter the data points

Stress=[87.896.6176263351571834122916242107267833804258];

Strain=153:51:765;

% Enter the derivative values calculated

Diff=[0.17250.86471.63141.71573.01964.73536.5417.74518.607810.333312.480415.490217.2157];

% Select the new data excluding the first four data points

Newdata=[Stress(5:end);Diff(5:end)];

% Perform a linear regression on these points

polyfit(Newdata(1, :),Newdata(2, :),1);

% Enter the formula for Stress obtained using the values of E0 and a

eqn=804.5591945812*(exp(0.0035499643*Strain)-1);

eqn_new=270.79152781*(exp(0.0035499643*Strain)-1);

% Now plot the points along with the data

scatter(Strain,Stress);hold on

plot(Strain,eqn_new)

Output:

It is clear that the function obtained using the data point from the middle is a much better fit for the provided data points.

Hence, the desired function that fits the data points is $\sigma =804.5591945812\left(e^{0.0035499643\epsilon }-1\right)$.