diff --git a/SO3Fun/@SO3FunHarmonic/quadrature.m b/SO3Fun/@SO3FunHarmonic/quadrature.m
index 3d4b23d0a..b682f631e 100644
--- a/SO3Fun/@SO3FunHarmonic/quadrature.m
+++ b/SO3Fun/@SO3FunHarmonic/quadrature.m
@@ -60,8 +60,12 @@
 % Curtis quadrature grid in fundamental region. 
 % Therefore adjust the bandwidth to crystal and specimen symmetry.
 bw = adjustBandwidth(N,SRight,SLeft);
-SO3G = quadratureSO3Grid(bw,'ClenshawCurtis',SRight,SLeft,'ABG');
-  
+if check_option(varargin,'GaussLegendre')
+  SO3G = quadratureSO3Grid(bw,'GaussLegendre',SRight,SLeft,'ABG');
+else
+  SO3G = quadratureSO3Grid(bw,'ClenshawCurtis',SRight,SLeft,'ABG');
+end
+
 % Only evaluate unique orientations
 values = f.eval(SO3G);
   
diff --git a/TensorAnalysis/@strainTensor/calcTaylor.m b/TensorAnalysis/@strainTensor/calcTaylor.m
index bb6ae466a..c376da030 100644
--- a/TensorAnalysis/@strainTensor/calcTaylor.m
+++ b/TensorAnalysis/@strainTensor/calcTaylor.m
@@ -46,7 +46,8 @@
   bw = get_option(varargin,'bandwidth',32);
   numOut = nargout;
   F = SO3FunHandle(@(rot) calcTaylorFun(rot,eps,sS,numOut,varargin{:}),sS.CS,eps.CS);
-  SO3F = SO3FunHarmonic(F,'bandwidth',bw);
+  % Use Gauss-Legendre quadrature, since the evaluation process is very expansive
+  SO3F = SO3FunHarmonic(F,'bandwidth',bw,'GaussLegendre');
   M = SO3F(1);
   if nargout>1
     b = [];