diff --git a/src/test_rotation.m b/src/test_rotation.m index 590b04e..1bc728a 100644 --- a/src/test_rotation.m +++ b/src/test_rotation.m @@ -1,3 +1,6 @@ +%% +clear variables + %% density = 3.291 @@ -28,89 +31,18 @@ tensor_olivine_ref = stiffnessTensor(Cij_Ol,... 'unit', 'GPa',... cs_olivine,... - 'density', density); - -%% Totation 90 degrees axe 3: x3 or z - -C11 = 182.1; -C22 = 280.2; -C33 = 207.6; -C44 = 68.8; -C55 = 56.8; -C66 = 68.5; -C12 = 71.9; -C13 = 70.1; -C23 = 67.2; - -% Elastic stiffness tensor (in GPa) -Cij_Ol_90x3 =... - [[C11 C12 C13 0.0 0.0 0.0];... - [ C12 C22 C23 0.0 0.0 0.0];... - [ C13 C23 C33 0.0 0.0 0.0];... - [ 0.0 0.0 0.0 C44 0.0 0.0];... - [ 0.0 0.0 0.0 0.0 C55 0.0];... - [ 0.0 0.0 0.0 0.0 0.0 C66]]; - -% generate stiffness tensor -tensor_olivine_90x3 = stiffnessTensor(Cij_Ol_90x3,... - 'unit', 'GPa',... - cs_olivine,... - 'density', density); - -%% Totation 90 degrees axe 1: x1 or x? - -C11 = 280.2; -C22 = 207.6; -C33 = 182.1; -C44 = 68.8; -C55 = 56.8; -C66 = 68.5; -C12 = 67.2; -C13 = 71.9; -C23 = 70.1; - -% Elastic stiffness tensor (in GPa) -Cij_Ol_90x1 =... - [[C11 C12 C13 0.0 0.0 0.0];... - [ C12 C22 C23 0.0 0.0 0.0];... - [ C13 C23 C33 0.0 0.0 0.0];... - [ 0.0 0.0 0.0 C44 0.0 0.0];... - [ 0.0 0.0 0.0 0.0 C55 0.0];... - [ 0.0 0.0 0.0 0.0 0.0 C66]]; - -% generate stiffness tensor -tensor_olivine_90x1 = stiffnessTensor(Cij_Ol_90x1,... - 'unit', 'GPa',... - cs_olivine,... - 'density', density); - - -%% Totation 90 degrees axe 2: x2 or y? + 'density', density) -C11 = 207.6; -C22 = 182.1; -C33 = 280.2; -C44 = 68.5; -C55 = 68.8; -C66 = 56.8; -C12 = 70.1; -C13 = 67.2; -C23 = 71.9; +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +%% +rot3 = rotation.byAxisAngle(zvector, 90*degree); +tensor_olivine_rot90zvector = rotate(tensor_olivine_ref, rot3) -% Elastic stiffness tensor (in GPa) -Cij_Ol_90x2 =... - [[C11 C12 C13 0.0 0.0 0.0];... - [ C12 C22 C23 0.0 0.0 0.0];... - [ C13 C23 C33 0.0 0.0 0.0];... - [ 0.0 0.0 0.0 C44 0.0 0.0];... - [ 0.0 0.0 0.0 0.0 C55 0.0];... - [ 0.0 0.0 0.0 0.0 0.0 C66]]; +rot1 = rotation.byAxisAngle(xvector, 90*degree); +tensor_olivine_rot90xvector = rotate(tensor_olivine_ref, rot1) -% generate stiffness tensor -tensor_olivine_90x2 = stiffnessTensor(Cij_Ol_90x2,... - 'unit', 'GPa',... - cs_olivine,... - 'density', density); +rot2 = rotation.byAxisAngle(yvector, 90*degree); +tensor_olivine_rot90yvector = rotate(tensor_olivine_ref, rot2) %% % estimate vp velocities @@ -121,17 +53,17 @@ [minVp, minVpPos] = min(vp); % Rotation 90 x3 -[vp_x3, ~, ~, ~, ~, ~] = velocity(tensor_olivine_90x3); +[vp_x3, ~, ~, ~, ~, ~] = velocity(tensor_olivine_rot90zvector); [maxVp_x3, maxVpPos_x3] = max(vp_x3); [minVp_x3, minVpPos_x3] = min(vp_x3); % Rotation 90 x1 -[vp_x1, ~, ~, ~, ~, ~] = velocity(tensor_olivine_90x1); +[vp_x1, ~, ~, ~, ~, ~] = velocity(tensor_olivine_rot90xvector); [maxVp_x1, maxVpPos_x1] = max(vp_x1); [minVp_x1, minVpPos_x1] = min(vp_x1); % Rotation 90 x2 -[vp_x2, ~, ~, ~, ~, ~] = velocity(tensor_olivine_90x2); +[vp_x2, ~, ~, ~, ~, ~] = velocity(tensor_olivine_rot90yvector); [maxVp_x2, maxVpPos_x2] = max(vp_x2); [minVp_x2, minVpPos_x2] = min(vp_x2); @@ -198,191 +130,16 @@ colormap(flipud(brewermap(256, 'Spectral'))) -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -%% Totation 45 degrees axe 3: x3 or z - -C11 = 220.03; -C22 = 220.03; -C33 = 207.6; -C44 = 62.8; -C55 = 62.8; -C66 = 79.63; -C12 = 83.03; -C13 = 68.65; -C23 = 68.65; -C16 = -24.53; -C26 = -24.53; -C36 = 1.45; -C45 = -6.0; - -% Elastic stiffness tensor (in GPa) -Cij_Ol_45x3 =... - [[C11 C12 C13 0.0 0.0 C16];... - [ C12 C22 C23 0.0 0.0 C26];... - [ C13 C23 C33 0.0 0.0 C36];... - [ 0.0 0.0 0.0 C44 C45 0.0];... - [ 0.0 0.0 0.0 C45 C55 0.0];... - [ C16 C26 C36 0.0 0.0 C66]]; - -% generate stiffness tensor -tensor_olivine_45x3 = stiffnessTensor(Cij_Ol_45x3,... - 'unit', 'GPa',... - cs_olivine,... - 'density', density); - -%% Totation 45 degrees axe 1: x1 or x? - -C11 = 280.2; -C22 = 200.98; -C33 = 200.98; -C44 = 62.8; -C55 = 62.8; -C66 = 62.38; -C12 = 69.55; -C13 = 69.55; -C23 = 63.98; -C16 = -2.35; -C26 = 6.38; -C36 = 6.38; -C45 = -6.0; - -% Elastic stiffness tensor (in GPa) -Cij_Ol_45x1 =... - [[C11 C12 C13 0.0 0.0 C16];... - [ C12 C22 C23 0.0 0.0 C26];... - [ C13 C23 C33 0.0 0.0 C36];... - [ 0.0 0.0 0.0 C44 C45 0.0];... - [ 0.0 0.0 0.0 C45 C55 0.0];... - [ C16 C26 C36 0.0 0.0 C66]]; - -% generate stiffness tensor -tensor_olivine_45x1 = stiffnessTensor(Cij_Ol_45x1,... - 'unit', 'GPa',... - cs_olivine,... - 'density', density); - - -%% Totation 45 degrees axe 2: x2 or y? - -C11 = 224.35; -C22 = 182.1; -C33 = 224.35; -C44 = 62.65; -C55 = 88.35; -C66 = 62.65; -C12 = 71.0; -C13 = 86.75; -C23 = 71.0; -C15 = -18.15; -C25 = -0.9; -C35 = -18.15; -C46 = 5.85; - -% Elastic stiffness tensor (in GPa) -Cij_Ol_45x2 =... - [[C11 C12 C13 0.0 C15 0.0];... - [ C12 C22 C23 0.0 C25 0.0];... - [ C13 C23 C33 0.0 C35 0.0];... - [ 0.0 0.0 0.0 C44 0.0 C46];... - [ C15 C25 C35 0.0 C55 0.0];... - [ 0.0 0.0 0.0 C46 0.0 C66]]; - -% generate stiffness tensor -tensor_olivine_45x2 = stiffnessTensor(Cij_Ol_45x2,... - 'unit', 'GPa',... - cs_olivine,... - 'density', density); - -%% -% estimate vp velocities - -% Rotation 45 x3 -[vp_x3, ~, ~, ~, ~, ~] = velocity(tensor_olivine_45x3); -[maxVp_x3, maxVpPos_x3] = max(vp_x3); -[minVp_x3, minVpPos_x3] = min(vp_x3); - -% Rotation 45 x1 -[vp_x1, ~, ~, ~, ~, ~] = velocity(tensor_olivine_45x1); -[maxVp_x1, maxVpPos_x1] = max(vp_x1); -[minVp_x1, minVpPos_x1] = min(vp_x1); - -% Rotation 45 x2 -[vp_x2, ~, ~, ~, ~, ~] = velocity(tensor_olivine_45x2); -[maxVp_x2, maxVpPos_x2] = max(vp_x2); -[minVp_x2, minVpPos_x2] = min(vp_x2); - - -%% figure Vp comparative -blackMarker = {'Marker','s','MarkerSize',14,'antipodal',... - 'MarkerEdgeColor','white','MarkerFaceColor','black','doNotDraw'}; -whiteMarker = {'Marker','o','MarkerSize',14,'antipodal',... - 'MarkerEdgeColor','black','MarkerFaceColor','white','doNotDraw'}; - -mtexFig = newMtexFigure('figSize', 'huge', 'layout', [1 4]); - -plot(vp, 'contourf', 'complete', 'doNotDraw', 'upper') -mtexTitle('Olivine reference', 'doNotDraw') -hold on -plot(maxVpPos(1), blackMarker{:}) -hold on -plot(minVpPos(1), whiteMarker{:}) -hold off -% percentage anisotropy -AVp = 200*(maxVp - minVp)./(maxVp + minVp); -xlabel({['Anisotropy = ', num2str(AVp,'%6.1f')]; ['min ', num2str(minVp,'%6.1f'), ' - ', num2str(maxVp,'%6.1f'), ' max']}) - - -nextAxis(1,2) -plot(vp_x3, 'contourf', 'complete', 'doNotDraw', 'upper') -mtexTitle('Olivine 45 degrees z(x3)', 'doNotDraw') -hold on -plot(maxVpPos_x3(1), blackMarker{:}) -hold on -plot(minVpPos_x3(1), whiteMarker{:}) -hold off -% percentage anisotropy -AVp_x3 = 200*(maxVp_x3 - minVp_x3)./(maxVp_x3 + minVp_x3); -xlabel({['Anisotropy = ', num2str(AVp_x3,'%6.1f')]; ['min ', num2str(minVp_x3,'%6.1f'), ' - ', num2str(maxVp_x3,'%6.1f'), ' max']}) - - -nextAxis(1,3) -plot(vp_x1, 'contourf', 'complete', 'doNotDraw', 'upper') -mtexTitle('Olivine 45 degrees x(x1)', 'doNotDraw') -hold on -plot(maxVpPos_x1(1), blackMarker{:}) -hold on -plot(minVpPos_x1(1), whiteMarker{:}) -hold off -% percentage anisotropy -AVp_x1 = 200*(maxVp_x1 - minVp_x1)./(maxVp_x1 + minVp_x1); -xlabel({['Anisotropy = ', num2str(AVp_x1,'%6.1f')]; ['min ', num2str(minVp_x1,'%6.1f'), ' - ', num2str(maxVp_x1,'%6.1f'), ' max']}) - -nextAxis(1,4) -plot(vp_x2, 'contourf', 'complete', 'doNotDraw', 'upper') -mtexTitle('Olivine 45 degrees y(x2)', 'doNotDraw') -hold on -plot(maxVpPos_x2(1), blackMarker{:}) -hold on -plot(minVpPos_x2(1), whiteMarker{:}) -hold off -% percentage anisotropy -AVp_x2 = 200*(maxVp_x2 - minVp_x2)./(maxVp_x2 + minVp_x2); -xlabel({['Anisotropy = ', num2str(AVp_x2,'%6.1f')]; ['min ', num2str(minVp_x2,'%6.1f'), ' - ', num2str(maxVp_x2,'%6.1f'), ' max']}) - -setColorRange('equal') -mtexColorbar('figSize', 'huge', 'title', 'km/s', 'FontSize', 24) -colormap(flipud(brewermap(256, 'Spectral'))) - %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% rot3 = rotation.byAxisAngle(zvector, 45*degree); -tensor_olivine_rot45zvector = rotate(tensor_olivine_ref, rot3); +tensor_olivine_rot45zvector = rotate(tensor_olivine_ref, rot3) rot1 = rotation.byAxisAngle(xvector, 45*degree); -tensor_olivine_rot45xvector = rotate(tensor_olivine_ref, rot1); +tensor_olivine_rot45xvector = rotate(tensor_olivine_ref, rot1) rot2 = rotation.byAxisAngle(yvector, 45*degree); -tensor_olivine_rot45yvector = rotate(tensor_olivine_ref, rot2); +tensor_olivine_rot45yvector = rotate(tensor_olivine_ref, rot2) diff --git a/src/test_rotation2.ipynb b/src/test_rotation2.ipynb index 884d697..0e23efd 100644 --- a/src/test_rotation2.ipynb +++ b/src/test_rotation2.ipynb @@ -4,7 +4,9 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "Rotating a rank-4 tensor involves applying the rotation to each of the tensor’s indices. For a stiffness tensor $C_{abcd}$ and a rotation matrix $R$, the transformation rule for rotating is given by:\n", + "# General procedure for rotating the stiffness tensor\n", + "\n", + "Rotating a rank-4 tensor involves applying the rotation to each of the tensor’s indices. For a stiffness tensor $C_{ijkl}$ and a rotation matrix $R$, the transformation rule for rotating is given by (we use $C_{abcd}$ for the original tensor):\n", "\n", "$$\n", "C'_{ijkl} = R_{ia} R_{jb} R_{kc} R_{ld} C_{abcd}\n", @@ -23,7 +25,7 @@ "source": [ "import numpy as np\n", "from scipy.spatial.transform import Rotation as R\n", - "from seismic_tools import _rearrange_tensor, _tensor_in_voigt\n", + "from tensor_tools import _rearrange_tensor, _tensor_in_voigt\n", "\n", "np.set_printoptions(suppress=True)\n", "\n", @@ -32,6 +34,12 @@ " stiffness_tensor: np.ndarray, rotation_matrix: np.ndarray\n", ") -> np.ndarray:\n", " \"\"\"Rotate a 3x3x3x3 symmetric stiffness tensor using numpy.einsum.\n", + " The operation is as follows:\n", + "\n", + " C'ijkl = Ria x Rjb x Rkc x Rld x Cabcd\n", + "\n", + " where Cabcd and C'ijkl are the original and the rotated tensor,\n", + " respectively, and R the rotation matrix\n", "\n", " Parameters\n", " ----------\n", @@ -152,13 +160,23 @@ "output_type": "stream", "text": [ "Rotated Stiffness Tensor:\n", - " [[220.02 83.02 68.65 0. 0. 24.52]\n", - " [ 83.02 220.02 68.65 0. 0. 24.53]\n", - " [ 68.65 68.65 207.6 0. 0. -1.45]\n", - " [ 0. 0. 0. 62.8 6. 0. ]\n", - " [ 0. 0. 0. 6. 62.8 0. ]\n", - " [ 24.52 24.53 -1.45 0. 0. 79.62]]\n" + "\n" ] + }, + { + "data": { + "text/plain": [ + "array([[220.02, 83.02, 68.65, 0. , 0. , 24.52],\n", + " [ 83.02, 220.02, 68.65, 0. , 0. , 24.53],\n", + " [ 68.65, 68.65, 207.6 , 0. , 0. , -1.45],\n", + " [ 0. , 0. , 0. , 62.8 , 6. , 0. ],\n", + " [ 0. , 0. , 0. , 6. , 62.8 , 0. ],\n", + " [ 24.52, 24.53, -1.45, 0. , 0. , 79.62]])" + ] + }, + "execution_count": 5, + "metadata": {}, + "output_type": "execute_result" } ], "source": [ @@ -172,7 +190,24 @@ "\n", "Cij_rot = _tensor_in_voigt(rotated_tensor)\n", "\n", - "print(\"Rotated Stiffness Tensor:\\n\", np.around(Cij_rot, 2))" + "print(\"Rotated Stiffness Tensor:\\n\")\n", + "np.around(Cij_rot, 2)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "**MTEX result**:\n", + "```\n", + "tensor in Voigt matrix representation:\n", + " 220.02 83.02 68.65 0 0 24.53\n", + " 83.02 220.03 68.65 0 0 24.53\n", + " 68.65 68.65 207.6 0 0 -1.45\n", + " 0 0 0 62.8 6 0\n", + " 0 0 0 6 62.8 0\n", + " 24.52 24.53 -1.45 0 0 79.62\n", + "```" ] }, { @@ -184,14 +219,23 @@ "name": "stdout", "output_type": "stream", "text": [ - "Rotated Stiffness Tensor:\n", - " [[280.2 69.55 69.55 2.35 0. 0. ]\n", - " [ 69.55 189.28 75.68 -6.38 0. 0. ]\n", - " [ 69.55 75.68 189.27 -6.37 0. 0. ]\n", - " [ 2.35 -6.38 -6.37 62.38 0. 0. ]\n", - " [ 0. 0. 0. 0. 68.65 -0.15]\n", - " [ 0. 0. 0. 0. -0.15 68.65]]\n" + "Rotated Stiffness Tensor:\n" ] + }, + { + "data": { + "text/plain": [ + "array([[280.2 , 69.55, 69.55, 2.35, 0. , 0. ],\n", + " [ 69.55, 189.28, 75.68, -6.38, 0. , 0. ],\n", + " [ 69.55, 75.68, 189.27, -6.37, 0. , 0. ],\n", + " [ 2.35, -6.38, -6.37, 62.38, 0. , 0. ],\n", + " [ 0. , 0. , 0. , 0. , 68.65, -0.15],\n", + " [ 0. , 0. , 0. , 0. , -0.15, 68.65]])" + ] + }, + "execution_count": 6, + "metadata": {}, + "output_type": "execute_result" } ], "source": [ @@ -203,7 +247,24 @@ "\n", "Cij_rot = _tensor_in_voigt(rotated_tensor)\n", "\n", - "print(\"Rotated Stiffness Tensor:\\n\", np.around(Cij_rot, 2))" + "print(\"Rotated Stiffness Tensor:\")\n", + "np.around(Cij_rot, 2)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "**MTEX result**:\n", + "```\n", + "tensor in Voigt matrix representation:\n", + " 280.2 69.55 69.55 2.35 0 0\n", + " 69.55 189.27 75.67 -6.38 0 0\n", + " 69.55 75.67 189.27 -6.38 0 0\n", + " 2.35 -6.38 -6.38 62.38 0 0\n", + " 0 0 0 0 68.65 -0.15\n", + " 0 0 0 0 -0.15 68.65\n", + "```" ] }, { @@ -215,14 +276,23 @@ "name": "stdout", "output_type": "stream", "text": [ - "Rotated Stiffness Tensor:\n", - " [[224.35 71. 86.75 0. -18.15 0. ]\n", - " [ 71. 182.1 71. 0. -0.9 0. ]\n", - " [ 86.75 71. 224.35 0. -18.15 0. ]\n", - " [ 0. 0. 0. 62.65 0. -5.85]\n", - " [-18.15 -0.9 -18.15 0. 88.35 0. ]\n", - " [ 0. 0. 0. -5.85 0. 62.65]]\n" + "Rotated Stiffness Tensor:\n" ] + }, + { + "data": { + "text/plain": [ + "array([[224.35, 71. , 86.75, 0. , -18.15, 0. ],\n", + " [ 71. , 182.1 , 71. , 0. , -0.9 , 0. ],\n", + " [ 86.75, 71. , 224.35, 0. , -18.15, 0. ],\n", + " [ 0. , 0. , 0. , 62.65, 0. , -5.85],\n", + " [-18.15, -0.9 , -18.15, 0. , 88.35, 0. ],\n", + " [ 0. , 0. , 0. , -5.85, 0. , 62.65]])" + ] + }, + "execution_count": 7, + "metadata": {}, + "output_type": "execute_result" } ], "source": [ @@ -234,14 +304,216 @@ "\n", "Cij_rot = _tensor_in_voigt(rotated_tensor)\n", "\n", - "print(\"Rotated Stiffness Tensor:\\n\", np.around(Cij_rot, 2))" + "print(\"Rotated Stiffness Tensor:\")\n", + "np.around(Cij_rot, 2)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "**MTEX result**:\n", + "```\n", + "tensor in Voigt matrix representation:\n", + " 224.35 71 86.75 0 -18.15 0\n", + " 71 182.1 71 0 -0.9 0\n", + " 86.75 71 224.35 0 -18.15 0\n", + " 0 0 0 62.65 0 -5.85\n", + " -18.15 -0.9 -18.15 0 88.35 0\n", + " 0 0 0 -5.85 0 62.65\n", + "```" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "MTEX pole figure plot\n", + "\n", + "![90 rot](https://raw.githubusercontent.com/marcoalopez/PyRockWave/main/src/img/mtex_ol_rotated.png)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "___" ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 8, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Rotated Stiffness Tensor:\n" + ] + }, + { + "data": { + "text/plain": [ + "array([[182.1, 71.9, 70.1, 0. , 0. , 0. ],\n", + " [ 71.9, 280.2, 67.2, 0. , 0. , 0. ],\n", + " [ 70.1, 67.2, 207.6, 0. , 0. , 0. ],\n", + " [ 0. , 0. , 0. , 68.8, 0. , 0. ],\n", + " [ 0. , 0. , 0. , 0. , 56.8, 0. ],\n", + " [ 0. , 0. , 0. , 0. , 0. , 68.5]])" + ] + }, + "execution_count": 8, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "# Define a rotation (for example, a 90-degree rotation around the z-axis)\n", + "rotation = R.from_euler('z', 90, degrees=True)\n", + "rotated_tensor = rotate_stiffness_tensor(Cijkl, rotation.as_matrix())\n", + "Cij_rot = _tensor_in_voigt(rotated_tensor)\n", + "\n", + "print(\"Rotated Stiffness Tensor:\") \n", + "np.around(Cij_rot, 2)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "**MTEX result**:\n", + "```\n", + " tensor in Voigt matrix representation:\n", + " 182.1 71.9 70.1 0 0 0\n", + " 71.9 280.2 67.2 0 0 0\n", + " 70.1 67.2 207.6 0 0 0\n", + " 0 0 0 68.8 0 0\n", + " 0 0 0 0 56.8 0\n", + " 0 0 0 0 0 68.5\n", + "```" + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Rotated Stiffness Tensor:\n" + ] + }, + { + "data": { + "text/plain": [ + "array([[280.2, 67.2, 71.9, 0. , 0. , 0. ],\n", + " [ 67.2, 207.6, 70.1, 0. , 0. , 0. ],\n", + " [ 71.9, 70.1, 182.1, 0. , 0. , 0. ],\n", + " [ 0. , 0. , 0. , 56.8, 0. , 0. ],\n", + " [ 0. , 0. , 0. , 0. , 68.5, 0. ],\n", + " [ 0. , 0. , 0. , 0. , 0. , 68.8]])" + ] + }, + "execution_count": 9, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "# Define a rotation (for example, a 90-degree rotation around the x-axis)\n", + "rotation = R.from_euler('x', 90, degrees=True)\n", + "rotated_tensor = rotate_stiffness_tensor(Cijkl, rotation.as_matrix())\n", + "Cij_rot = _tensor_in_voigt(rotated_tensor)\n", + "\n", + "print(\"Rotated Stiffness Tensor:\") \n", + "np.around(Cij_rot, 2)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "**MTEX result**:\n", + "```\n", + "tensor in Voigt matrix representation:\n", + " 280.2 67.2 71.9 0 0 0\n", + " 67.2 207.6 70.1 0 0 0\n", + " 71.9 70.1 182.1 0 0 0\n", + " 0 0 0 56.8 0 0\n", + " 0 0 0 0 68.5 0\n", + " 0 0 0 0 0 68.8\n", + "```" + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Rotated Stiffness Tensor:\n" + ] + }, + { + "data": { + "text/plain": [ + "array([[207.6, 70.1, 67.2, 0. , 0. , 0. ],\n", + " [ 70.1, 182.1, 71.9, 0. , 0. , 0. ],\n", + " [ 67.2, 71.9, 280.2, 0. , 0. , 0. ],\n", + " [ 0. , 0. , 0. , 68.5, 0. , 0. ],\n", + " [ 0. , 0. , 0. , 0. , 68.8, 0. ],\n", + " [ 0. , 0. , 0. , 0. , 0. , 56.8]])" + ] + }, + "execution_count": 10, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "# Define a rotation (for example, a 90-degree rotation around the y-axis)\n", + "rotation = R.from_euler('y', 90, degrees=True)\n", + "rotated_tensor = rotate_stiffness_tensor(Cijkl, rotation.as_matrix())\n", + "Cij_rot = _tensor_in_voigt(rotated_tensor)\n", + "\n", + "print(\"Rotated Stiffness Tensor:\") \n", + "np.around(Cij_rot, 2)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "**MTEX result**:\n", + "```\n", + "tensor in Voigt matrix representation:\n", + " 207.6 70.1 67.2 0 0 0\n", + " 70.1 182.1 71.9 0 0 0\n", + " 67.2 71.9 280.2 0 0 0\n", + " 0 0 0 68.5 0 0\n", + " 0 0 0 0 68.8 0\n", + " 0 0 0 0 0 56.8\n", + "```" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "MTEX pole figure plot\n", + "\n", + "![90 rot](https://raw.githubusercontent.com/marcoalopez/PyRockWave/main/src/img/mtex_ol.png)" + ] + }, + { + "cell_type": "markdown", "metadata": {}, - "outputs": [], "source": [] } ],