Add option to rotate material coordinate frame

dune/microstructure/matrix_operations.hh

@@ -105,6 +105,64 @@ namespace MatrixOperations {
 	    return sum;
 	}
 
+
+	/**
+	 * @brief Determines rotation matrix based on an axis and an angle
+	 * 
+	 * @param axis 
+	 * @param angle 
+	 * @return MatrixRT 
+	 */
+	static MatrixRT rotationMatrix(int axis, double angle){
+
+		switch (axis)
+		{
+			case 0:
+			{
+				return {{1.0,          0,            0   },
+						{  0, cos(angle), -1.0*sin(angle)},
+						{  0, sin(angle),      cos(angle)}};
+			}
+			case 1:
+			{
+				return {{     cos(angle),   0, sin(angle)},
+						{              0, 1.0,      0    },
+						{-1.0*sin(angle),   0, cos(angle)}};
+			}
+			case 2:
+			{
+				return {{cos(angle), -1.0*sin(angle),    0},
+						{sin(angle),      cos(angle),    0},
+						{         0,               0,  1.0}};
+			}
+			default:
+				DUNE_THROW(Dune::Exception, " axis not feasible. rotationMatrix is only implemented for 3x3-matrices.");
+		}
+	}
+
+	/**
+	 * @brief 6x6 matrix that transforms the strain tensor. This matrix is used to 
+	 * 		  transform the compliance matrix (given in Voigt notation) into another frame.
+	 * 		  see 'https://en.wikipedia.org/wiki/Orthotropic_material#Condition_for_material_symmetry_2' 
+	 * 		  for details.
+	 * 
+	 * @param axis 
+	 * @param angle 
+	 * @return MatrixRT 
+	 */
+	static Dune::FieldMatrix<double,6,6> rotationMatrixCompliance(int axis, double angle){
+		
+		MatrixRT R = rotationMatrix(axis,angle);
+
+        return {{    R[0][0]*R[0][0],     R[0][1]*R[0][1],     R[0][2]*R[0][2],                 R[0][1]*R[0][2],                 R[0][0]*R[0][2],                 R[0][0]*R[0][1]},
+                {    R[1][0]*R[1][0],     R[1][1]*R[1][1],     R[1][2]*R[1][2],                 R[1][1]*R[1][2],                 R[1][0]*R[1][2],                 R[1][0]*R[1][1]},
+                {    R[2][0]*R[2][0],     R[2][1]*R[2][1],     R[2][2]*R[2][2],                 R[2][1]*R[2][2],                 R[2][0]*R[2][2],                 R[2][0]*R[2][1]},
+                {2.0*R[1][0]*R[2][0], 2.0*R[1][1]*R[2][1], 2.0*R[1][2]*R[2][2], R[1][1]*R[2][2]+R[1][2]*R[2][1], R[1][0]*R[2][2]+R[1][2]*R[2][0], R[1][0]*R[2][1]+R[1][1]*R[2][0]},
+                {2.0*R[0][0]*R[2][0], 2.0*R[0][1]*R[2][1], 2.0*R[0][2]*R[2][2], R[0][1]*R[2][2]+R[0][2]*R[2][1], R[0][0]*R[2][2]+R[0][2]*R[2][0], R[0][0]*R[2][1]+R[0][1]*R[2][0]},
+                {2.0*R[0][0]*R[0][0], 2.0*R[0][1]*R[1][1], 2.0*R[0][2]*R[1][2], R[0][1]*R[1][2]+R[0][2]*R[1][1], R[0][0]*R[1][2]+R[0][2]*R[1][0], R[0][0]*R[1][1]+R[0][1]*R[1][0]}
+               };
+	}
+
 // 	static double linearizedStVenantKirchhoffDensity(double mu, double lambda, MatrixRT E1, MatrixRT E2){    // ?? Check with Robert
 // 		E1= sym(E1);
 // 		E2 = sym(E2);

materials/material_orthotropic.py

@@ -41,7 +41,19 @@ materialParameters_phase1 = [11.2e3,630,1190,700,230,960,0.63 ,0.49,0.37]    # w
 
 #- PHASE 2
 phase2_type="orthotropic"
-materialParameters_phase2 = [10.7e3,430,710,620,23,500, 0.51 ,0.38,0.31]   # Norway spruce parameters (values for compliance matrix) see [Dimwoodie; Timber its nature and behavior p.109]
+# materialParameters_phase2 = [10.7e3,430,710,620,23,500, 0.51 ,0.38,0.31]   # Norway spruce parameters (values for compliance matrix) see [Dimwoodie; Timber its nature and behavior p.109]
+materialParameters_phase2 = [11.2e3,630,1190,700,230,960,0.63 ,0.49,0.37] 
+
+# Pass a set of FrameVectors to transform material properties to this Frame
+# phase2_FrameVector1 = [1, 0 ,0]
+# phase2_FrameVector2 = [0, 5, 0]
+# phase2_FrameVector3 = [0, 0, 1]
+
+phase2_axis = 2
+phase2_angle = np.pi/2.0
+# phase2_angle = 2*np.pi/12
+
+
 
 #- PHASE 3
 phase3_type="isotropic"