Modified to speed up.

git-svn-id: https://svn.code.sf.net/p/dct/code/trunk@598 4c865b51-4357-4376-afb4-474e03ccb993

zUtil_Maths/gtMathsRotationTensor.m

@@ -24,17 +24,31 @@ function Srot = gtMathsRotationTensor(om, rotcomp)
 if (length(om) == 1)
     Srot = rotcomp.const + rotcomp.cos*cosd(om) + rotcomp.sin*sind(om);
 else
-    n = length(om);
     
-    omcol(1,1,:) = om;
-    cosOmcol = cosd(omcol);
-    sinOmcol = sind(omcol);
+    % Expand cos(om)
+    cosom  = reshape(cosd(om),1,[]);
+    cosomr = cosom(ones(9,1),:);
+    cosomr = reshape(cosomr,3,[]);
+
+    % Expand sin(om)
+    sinom  = reshape(sind(om),1,[]);
+    sinomr = sinom(ones(9,1),:);
+    sinomr = reshape(sinomr,3,[]);
+    
+    % Expand matrices
+    ind1 = [1; 2; 3];
+    ind2 = ind1(:,ones(1,length(om)));
+    
+    constr = rotcomp.const(ind1, ind2);
+    cosr   = rotcomp.cos(ind1, ind2);
+    sinr   = rotcomp.sin(ind1, ind2);
+    
+    % Multiply element-wise
+    Srot = constr + cosr.*cosomr + sinr.*sinomr;
+    
+    % Reshape final output
+    Srot = reshape(Srot,3,3,[]);
     
-    Srot = rotcomp.const(:, :, ones(n,1)) ...
-           + rotcomp.cos(:, :, ones(n,1)) ...
-           .* cosOmcol(ones(3, 1), ones(3,1), :) ...
-           + rotcomp.sin(:, :, ones(n,1)) ...
-           .* sinOmcol(ones(3, 1), ones(3,1), :);
 end
     
 end

zUtil_Strain2/gtStrainPlaneNormals.m

-function [pldef, drel] = gtStrainPlaneNormals(pl, G)
+function [pldef, drel] = gtStrainPlaneNormals(pl, defT)
 % GTSTRAINPLANENORMALS Given a deformation tensor, calculates new plane 
 % normals with no approximations.
 % 
-%   [pldef, drel] = gtStrainPlaneNormals(pl, G)
+%   [pldef, drel] = gtStrainPlaneNormals(pl, deft)
 %
 % Returns the new orientations of plane normals 'pl' after a deformation 
-% described by the deformation tensor 'G' (9 distinct components).
+% described by the deformation tensor 'defT' (9 distinct components).
 % Also returns the elongations 'drel' due to the deformation in the 
 % directions of the plane normals. In a crystal this equals to the ratios 
 % of the d-spacings after and before the deformation.
 % The results are exact, no approximation is used.
-% The 'pl' and 'G' coordinates must be given in the same reference frame.
+% The 'pl' and 'defT' coordinates must be given in the same reference frame.
 %
 % INPUT
-%   pl - normalised coordinates of plane normals (3xn)
-%   G  - deformation tensor (3x3)
+%   pl    - normalised coordinates of plane normals (3xn)
+%   defT  - deformation tensor (3x3xn)
 %
 % OUTPUT
 %   pldef - plane normals in deformed state
@@ -61,9 +61,22 @@ e2 = [pl(2,:).*e1(3,:) - pl(3,:).*e1(2,:);
 % sum(ch.*pl, 1) < 0
 
 
-% e-s in deformed state:
-e1def = G*e1 + e1;
-e2def = G*e2 + e2;
+% e1 and e2 vectors in deformed state:
+if (size(defT,3) == 1)
+    e1def = defT*e1 + e1;
+    e2def = defT*e2 + e2;
+else
+    % Expand vectors for multiplication element-wise
+    e1r   = reshape(e1,1,[]);
+    defTe = reshape(defT,3,[],1).*[e1r; e1r; e1r];
+    defTe = defTe(:,1:3:end) + defTe(:,2:3:end) + defTe(:,3:3:end);
+    e1def = defTe + e1;
+
+    e2r   = reshape(e2,1,[]);
+    defTe = reshape(defT,3,[],1).*[e2r; e2r; e2r];
+    defTe = defTe(:,1:3:end) + defTe(:,2:3:end) + defTe(:,3:3:end);
+    e2def = defTe + e2;
+end
 
 
 % new deformed pl
@@ -84,14 +97,23 @@ pldef     = pldef./normpldef([1 1 1],:);
 %sum(pldef.*pl, 1) >= 0
 
 
+% New length along original plane normal
+% (exact calculation, non-linear in the defT components)
+
+if (size(defT,3) == 1)
+    fvec = defT*pl + pl;
+else
+    plr   = reshape(pl,1,[]);
+    defTe = reshape(defT,3,[],1).*[plr; plr; plr];
+    defTe = defTe(:,1:3:end) + defTe(:,2:3:end) + defTe(:,3:3:end);
+    
+    fvec = defTe + pl;   
+end
 
-% new relative d-spacing (exact calculation, non-linear in the G components)
-fvec = pl + G*pl;   % change along original plane normal
 
 % new normalised d-spacing
 drel = sum(pldef.*fvec, 1);
 
 
-
 end % of function