# 3D plot of a anisotropy of young's modulus

Sorry, I know a lot of similar questions have been asked before, but I just can't find one that exactly describes my problem.

How can i draw a 3D plot of anisotrophy of young's modulus from following equation

As he drawn in Mathematica http://demonstrations.wolfram.com/AnisotropicElasticity/

Looking forward to yours help

• Use SphericalPlot3D, it is the easiest way. – Mauricio Fernández Apr 13 '17 at 17:08
• You can download the .nb file from the demonstration and see for yourself ;) [use the link that says Download Author Code »(preview »)] – Sumit Apr 13 '17 at 18:36

I will define Young's modulus viewed in tensile direction $\boldsymbol{n}$ of a material with compliance $\mathbb{S}$, being the inverse of the stiffness $\mathbb{C}$, i.e., $\mathbb{C} = \mathbb{S}^{-1}$ (in the sense of the inverse of second-order symmetric tensors), as done in this publication

http://www.uni-magdeburg.de/ifme/zeitschrift_tm/2001_Heft2/Boehlke_Brueggeman.pdf

$$Y(\boldsymbol{n},\mathbb{S}) = (\mathbb{S} \cdot \boldsymbol{n}^{\otimes 4})^{-1} = (S_{ijkl} n_i n_j n_k n_l)^{-1}$$

In the following code I create the compliance of a cubic material and plot Young's modulus with SphericalPlot3D

(*tensile direction*)
n = {Cos[phi]*Sin[theta], Sin[phi]*Sin[theta], Cos[theta]};
(*Basic tensors*)
I2 = IdentityMatrix@3;
IdI = TensorProduct[I2, I2];
I4 = TensorTranspose[IdI, {1, 3, 2, 4}];
IS = (I4 + TensorTranspose[I4, {1, 2, 4, 3}])/2;
(*Projectors of cubic materials*)
Pc1 = 1/3*IdI;
Dc = ConstantArray[0, {3, 3, 3, 3}];
Dc[[1, 1, 1, 1]] = 1;
Dc[[2, 2, 2, 2]] = 1;
Dc[[3, 3, 3, 3]] = 1;
Pc2 = Dc - Pc1;
Pc3 = IS - Pc1 - Pc2;
(*Spectral represention of cubic compliance*)
comp = 1/l1*Pc1 + 1/l2*Pc2 + 1/l3*Pc3;
l1 = 10;
l2 = 3;
l3 = 8;
(*Young's modulues*)
Y = (Total[comp*TensorProduct[n, n, n, n], Infinity])^(-1) // Simplify;
(*Plot*)
SphericalPlot3D[Y, {theta, 0, Pi}, {phi, 0, 2 Pi}]


• Thanx for your help but still i have some problem could you please help me little more in my case C11 = 354.06, C12 = 100.20, C44 = 170 and S11 = 0.003227305, S12 = -0.000711874, S44 = 0.005882353 (all units in GPa) by using these values plot should be spherical but it not – user48104 Apr 14 '17 at 4:25
• @Mauricio Lobos Fernández, can you also show the % Poisson ratio % 1. Maximum / 2. Minimum / 3. Average and % Shear modulus % 1. Maximum / 2. Minimum / 3. Average, both are also important, thanking you! – ABCDEMMM Jul 20 '18 at 1:45
• That depends on the definition you use and usually the number of variables are too high to plot. – Mauricio Fernández Jul 21 '18 at 8:10
 e[n_, s11_, b1_] := 1/(s11 - b1 Total[Times @@@ Subsets[n^2, {2}]])
SphericalPlot3D[
e[FromSphericalCoordinates[{1, p, t}], 2, 3],
{p, 0, Pi}, {t, -Pi, Pi}]


• Thanx for your help @george2079 after pasting the code in Mathematica i just get empty frame could you be help me little more whats going wrong? – user48104 Apr 14 '17 at 4:34
• if you have an older version that doesn't have FromSphericalCoordinates replace that with {Cos[t] Sin[p], Sin[p] Sin[t], Cos[p]} – george2079 Apr 14 '17 at 14:51
• After replacing the FromSphericalCoordinates into {Cos[t] Sin[p], Sin[p] Sin[t], Cos[p]} still getting the empty frame i dont know why?? – user48104 Apr 14 '17 at 16:41
• @george2079, can you also show the % Poisson ratio % 1. Maximum / 2. Minimum / 3. Average and % Shear modulus % 1. Maximum / 2. Minimum / 3. Average, both are also important, thanking you! – ABCDEMMM Jul 20 '18 at 1:46