Skip to content

Commit bba46c1

Browse files
authored
BUG: Fix Monoclinic Symmertry Operator. Validate all other symmetry operators (#23)
* Validate all LaueOps Classes. * New SymmetryOperator values were generated starting from Quaternions that are defined in EMsoft. That code is in the gen_sym_code.cpp class. Output was hand copied to the LaueOps classes and then optimized from there. * Added API to return the Rotation Point Group value as a string Signed-off-by: Michael Jackson <mike.jackson@bluequartz.net> --------- Signed-off-by: Michael Jackson <mike.jackson@bluequartz.net>
1 parent d7db8a3 commit bba46c1

29 files changed

+1100
-432
lines changed

CMakeLists.txt

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -15,7 +15,7 @@ set(CMAKE_CXX_EXTENSIONS OFF)
1515
set(CMAKE_EXPORT_COMPILE_COMMANDS ON)
1616

1717
# set project's name
18-
project(EbsdLibProj VERSION 1.0.30)
18+
project(EbsdLibProj VERSION 1.0.31)
1919

2020
# ---------- Setup output Directories -------------------------
2121
set(CMAKE_LIBRARY_OUTPUT_DIRECTORY

Docs/Index.md

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,8 @@
1+
# Various Bits of Documentation for EbsdLib
2+
3+
EbsdLib is primarily used in the [DREAM3D](https://www.dream3d.io) family of applications and libraries.
4+
5+
## Rotation Point Groups
6+
7+
The PDF is courtesy of Dr. Anthony Rollett from Carnegie Mellon University. The original
8+
URL is [http://pajarito.materials.cmu.edu/lectures/L3-OD_symmetry-21Jan16-slide_50-operators.pdf](http://pajarito.materials.cmu.edu/lectures/L3-OD_symmetry-21Jan16-slide_50-operators.pdf)
Binary file not shown.

Source/Apps/SourceList.cmake

Lines changed: 4 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -6,7 +6,10 @@ add_executable(make_ipf ${EbsdLibProj_SOURCE_DIR}/Source/Apps/make_ipf.cpp)
66
target_link_libraries(make_ipf PUBLIC EbsdLib)
77
target_include_directories(make_ipf PUBLIC ${EbsdLibProj_SOURCE_DIR}/Source)
88

9-
109
add_executable(convert_orientations ${EbsdLibProj_SOURCE_DIR}/Source/Apps/ConvertOrientations.cpp)
1110
target_link_libraries(convert_orientations PUBLIC EbsdLib)
1211
target_include_directories(convert_orientations PUBLIC ${EbsdLibProj_SOURCE_DIR}/Source)
12+
13+
add_executable(gen_sym_code ${EbsdLibProj_SOURCE_DIR}/Source/Apps/gen_sym_code.cpp)
14+
target_link_libraries(gen_sym_code PUBLIC EbsdLib)
15+
target_include_directories(gen_sym_code PUBLIC ${EbsdLibProj_SOURCE_DIR}/Source)

Source/Apps/gen_sym_code.cpp

Lines changed: 240 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,240 @@
1+
2+
#include "EbsdLib/Core/EbsdLibConstants.h"
3+
4+
#include "EbsdLib/Core/OrientationTransformation.hpp"
5+
#include "EbsdLib/LaueOps/LaueOps.h"
6+
#include "EbsdLib/Math/EbsdLibMath.h"
7+
#include "EbsdLib/Math/Matrix3X1.hpp"
8+
#include "EbsdLib/Math/Matrix3X3.hpp"
9+
10+
#include <cmath>
11+
#include <fstream>
12+
#include <iomanip>
13+
#include <iostream>
14+
#include <limits>
15+
#include <map>
16+
#include <string>
17+
#include <vector>
18+
19+
const double sq22 = 0.7071067811865475244; // sqrt(2)/2
20+
static double sq32 = std::sqrt(3.0) / 2.0;
21+
22+
const double half = 0.50; // 1/2
23+
24+
std::vector<QuatD> SYM_Qsymop = {
25+
QuatD(0.0, 0.0, 0.0, 0.0), // 0: JUNK.
26+
QuatD(1.0, 0.0, 0.0, 0.0), // 1: identity operator
27+
QuatD(0.0, 1.0, 0.0, 0.0), // 2: 180@[100]
28+
QuatD(0.0, 0.0, 1.0, 0.0), // 3: 180@[010]
29+
QuatD(0.0, 0.0, 0.0, 1.0), // 4: 180@[001]
30+
QuatD(sq22, sq22, 0.0, 0.0), // 5: 90@[100]
31+
QuatD(sq22, 0.0, sq22, 0.0), // 6: 90@[010]
32+
QuatD(sq22, 0.0, 0.0, sq22), // 7: 90@[001]
33+
QuatD(sq22, -sq22, 0.0, 0.0), // 8: 270@[100]
34+
QuatD(sq22, 0.0, -sq22, 0.0), // 9: 270@[010]
35+
QuatD(sq22, 0.0, 0.0, -sq22), // 10: 270@[001]
36+
QuatD(0.0, sq22, sq22, 0.0), // 11: 180@[110]
37+
QuatD(0.0, -sq22, sq22, 0.0), // 12: 180@[-110]
38+
QuatD(0.0, 0.0, sq22, sq22), // 13: 180@[011]
39+
QuatD(0.0, 0.0, -sq22, sq22), // 14: 180@[0-11]
40+
QuatD(0.0, sq22, 0.0, sq22), // 15: 180@[101]
41+
QuatD(0.0, -sq22, 0.0, sq22), // 16: 180@[-101]
42+
QuatD(half, half, half, half), // 17: 120@[111]
43+
QuatD(half, -half, -half, -half), // 18: 120@[-1-1-1]
44+
QuatD(half, half, -half, half), // 19: 120@[1-11]
45+
QuatD(half, -half, half, -half), // 20: 120@[-11-1]
46+
QuatD(half, -half, half, half), // 21: 120@[-111]
47+
QuatD(half, half, -half, -half), // 22: 120@[1-1-1]
48+
QuatD(half, -half, -half, half), // 23: 120@[-1-11]
49+
QuatD(half, half, half, -half), // 24: 120@[11-1]
50+
QuatD(sq32, 0.0, 0.0, half), // 25: 60@[001] (hexagonal/trigonal operators start here)
51+
QuatD(half, 0.0, 0.0, sq32), // 26: 120@[001]
52+
QuatD(0.0, 0.0, 0.0, 1.0), // 27: 180@[001] (duplicate from above, but useful to keep it here)
53+
QuatD(-half, 0.0, 0.0, sq32), // 28: 240@[001]
54+
QuatD(-sq32, 0.0, 0.0, half), // 29: 300@[001]
55+
QuatD(0.0, 1.0, 0.0, 0.0), // 30: 180@[100]
56+
QuatD(0.0, sq32, half, 0.0), // 31: 180@[xxx]
57+
QuatD(0.0, half, sq32, 0.0), // 32: 180@[xxx]
58+
QuatD(0.0, 0.0, 1.0, 0.0), // 33: 180@[010]
59+
QuatD(0.0, -half, sq32, 0.0), // 34: 180@[xxx]
60+
QuatD(0.0, -sq32, half, 0.0), // 35: 180@[xxx]
61+
};
62+
63+
std::vector<QuatD> InitRotationPointGroup(int rotPointGroup)
64+
{
65+
std::vector<QuatD> sym_ops;
66+
67+
// identity operator is part of all point groups
68+
// sym_ops.at = 0.D0 // initialize all entries to zero
69+
sym_ops.push_back(SYM_Qsymop[1]);
70+
size_t i = 0;
71+
// select statement for each individual rotational point group (see typedefs.f90 for SYM_Qsymop definitions)
72+
switch(rotPointGroup)
73+
{
74+
case(1): // 1 (no additional symmetry elements)
75+
sym_ops.resize(1);
76+
// sym_ops.at(2-1) = -sym_ops.at(1)
77+
break;
78+
case(2): // 2 (we'll assume that the two-fold axis lies along the e_y-axis)
79+
sym_ops.resize(2);
80+
sym_ops.at(2 - 1) = SYM_Qsymop[3];
81+
break;
82+
case(222): // 222
83+
sym_ops.resize(4);
84+
for(i = 2; i <= 4; i++)
85+
sym_ops.at(i - 1) = SYM_Qsymop[i];
86+
87+
break;
88+
case(4): // 4
89+
sym_ops.resize(4);
90+
sym_ops.at(2 - 1) = SYM_Qsymop[4];
91+
sym_ops.at(3 - 1) = SYM_Qsymop[7];
92+
sym_ops.at(4 - 1) = SYM_Qsymop[10];
93+
break;
94+
case(422): // 422
95+
sym_ops.resize(8);
96+
sym_ops.at(2 - 1) = SYM_Qsymop[4];
97+
sym_ops.at(3 - 1) = SYM_Qsymop[7];
98+
sym_ops.at(4 - 1) = SYM_Qsymop[10];
99+
sym_ops.at(5 - 1) = SYM_Qsymop[2];
100+
sym_ops.at(6 - 1) = SYM_Qsymop[3];
101+
sym_ops.at(7 - 1) = SYM_Qsymop[11];
102+
sym_ops.at(8 - 1) = SYM_Qsymop[12];
103+
break;
104+
case(3): // 3
105+
sym_ops.resize(3);
106+
sym_ops.at(2 - 1) = SYM_Qsymop[26];
107+
sym_ops.at(3 - 1) = SYM_Qsymop[28];
108+
// call FatalError('InitDictionaryIndexing','this symmetry has not yet been implemented (pg 3)')
109+
break;
110+
case(32): // 32 (needs special handling)
111+
sym_ops.resize(6);
112+
sym_ops.at(2 - 1) = SYM_Qsymop[26];
113+
sym_ops.at(3 - 1) = SYM_Qsymop[28];
114+
sym_ops.at(4 - 1) = SYM_Qsymop[30];
115+
sym_ops.at(5 - 1) = SYM_Qsymop[32];
116+
sym_ops.at(6 - 1) = SYM_Qsymop[34];
117+
// call FatalError('InitDictionaryIndexing','this symmetry has not yet been implemented (pg 32)')
118+
break;
119+
case(6): // 6
120+
sym_ops.resize(6);
121+
for(i = 25; i <= 29; i++) // i=25,29
122+
sym_ops.at(i - 23 - 1) = SYM_Qsymop[i];
123+
124+
break;
125+
case(622): // 622
126+
sym_ops.resize(12);
127+
for(i = 25; i <= 35; i++) // do i=25,35
128+
sym_ops.at(i - 23 - 1) = SYM_Qsymop[i];
129+
130+
break;
131+
case(23): // 23
132+
sym_ops.resize(12);
133+
for(i = 2; i <= 4; i++) // do i=2,4
134+
sym_ops.at(i - 1) = SYM_Qsymop[i];
135+
136+
for(i = 17; i <= 24; i++) // do i=17,24
137+
sym_ops.at(4 + (i - 16) - 1) = SYM_Qsymop[i];
138+
139+
break;
140+
case(432): // 432
141+
sym_ops.resize(24);
142+
for(i = 2; i <= 24; i++) // do i=2,24
143+
sym_ops.at(i - 1) = SYM_Qsymop[i];
144+
145+
break;
146+
147+
default: // this should never happen ...
148+
throw std::runtime_error("Rotation Point Group is not recognized");
149+
}
150+
151+
return sym_ops;
152+
}
153+
154+
void Print3x3(const OrientationD& m)
155+
{
156+
std::cout << std::fixed;
157+
std::cout << std::setprecision(16);
158+
size_t i = 0;
159+
std::cout << " {{" << m[i++] << ", " << m[i++] << ", " << m[i++] << "},\n";
160+
std::cout << " {" << m[i++] << ", " << m[i++] << ", " << m[i++] << "},\n";
161+
std::cout << " {" << m[i++] << ", " << m[i++] << ", " << m[i++] << "}},\n";
162+
}
163+
164+
void PrintQuat(const QuatD& q, bool sv)
165+
{
166+
std::cout << std::fixed;
167+
std::cout << std::setprecision(16);
168+
size_t i = 0;
169+
if(sv)
170+
{
171+
std::cout << " QuatD(" << q[1] << ", " << q[2] << ", " << q[3] << ", " << q[0] << ")";
172+
}
173+
{
174+
std::cout << " QuatD(" << q[i++] << ", " << q[i++] << ", " << q[i++] << ", " << q[i++] << ")";
175+
}
176+
}
177+
178+
void PrintRod(const OrientationD& q)
179+
{
180+
std::cout << std::fixed;
181+
std::cout << std::setprecision(16);
182+
183+
std::cout << " {" << (std::isinf(q[0]) ? 10.0E12 : q[0]) << ", " << (std::isinf(q[1]) ? 10.0E12 : q[1]) << ", " << (std::isinf(q[2]) ? 10.0E12 : q[2]) << ", "
184+
<< (std::isinf(q[3]) ? 10.0E12 : q[3]) << "}";
185+
}
186+
187+
// -----------------------------------------------------------------------------
188+
int main(int argc, char* argv[])
189+
{
190+
std::cout << "Starting SymOP GenCode" << std::endl;
191+
192+
std::vector<int> rotPointGroup = {1, 2, 222, 4, 422, 3, 32, 6, 622, 23, 432};
193+
std::vector<std::string> names = {"TriclinicOps", "MonoclinicOps", "OrthorhombicOps", "TetragonalLowOps", "TetragonalOps", "TrigonalLowOps",
194+
"TrigonalOps", "HexagonalLowOps", "HexagonalOps", "CubicLowOps", "CubicOps"};
195+
196+
for(size_t i = 0; i < rotPointGroup.size(); i++)
197+
{
198+
std::cout << "###########################################################################\n";
199+
std::cout << names[i] << "\n// Rotation Point Group: " << rotPointGroup[i] << std::endl;
200+
auto symOPs = InitRotationPointGroup(rotPointGroup[i]);
201+
202+
for(auto& quat : symOPs)
203+
{
204+
quat = QuatD(quat[1], quat[2], quat[3], quat[0]);
205+
}
206+
207+
std::cout << std::setprecision(8);
208+
std::cout << "// clang-format off\n";
209+
std::cout << "static const std::vector<QuatD> QuatSym ={\n";
210+
for(const auto& symOp : symOPs)
211+
{
212+
PrintQuat(symOp, false);
213+
std::cout << ",\n";
214+
}
215+
std::cout << "};\n\n";
216+
217+
std::cout << "static const std::vector<OrientationD> RodSym = {\n";
218+
for(const auto& symOp : symOPs)
219+
{
220+
auto ro = OrientationTransformation::qu2ro<QuatD, OrientationD>(symOp);
221+
PrintRod(ro);
222+
std::cout << ",\n";
223+
}
224+
std::cout << "};\n\n";
225+
226+
std::cout << "static const double MatSym[k_SymOpsCount][3][3] = {\n";
227+
for(const auto& symOp : symOPs)
228+
{
229+
auto om = OrientationTransformation::qu2om<QuatD, OrientationD>(symOp);
230+
Print3x3(om);
231+
std::cout << " \n";
232+
}
233+
std::cout << "};\n";
234+
std::cout << "// clang-format on\n\n";
235+
}
236+
237+
return -1;
238+
239+
return 0;
240+
}

0 commit comments

Comments
 (0)