468 lines
8.8 KiB
C
468 lines
8.8 KiB
C
/*
|
|
Space Group Info's (c) 1994-96 Ralf W. Grosse-Kunstleve
|
|
*/
|
|
|
|
#include <stdio.h>
|
|
#include <stdlib.h>
|
|
|
|
|
|
#ifdef APP_INCLUDE
|
|
#include APP_INCLUDE
|
|
#endif
|
|
|
|
#ifndef AppMalloc
|
|
#define AppMalloc(ptr, n) (ptr) = malloc((n) * sizeof (*(ptr)))
|
|
#endif
|
|
#ifndef AppFree
|
|
#define AppFree(ptr, n) free(ptr)
|
|
#endif
|
|
|
|
|
|
#include "sginfo.h"
|
|
|
|
|
|
/* Non elegant way to get s.i. vectors and moduli:
|
|
1. Build field with legal reference points marked (TestField)
|
|
2. Go through list of possible s.i. vects and mods:
|
|
Verify with TestField
|
|
*/
|
|
|
|
|
|
void MarkLegalOrigins(const T_SgInfo * SgInfo, int *TestField)
|
|
{
|
|
int O[3], V[3], lx, ly, lz, mx, my, mz, i;
|
|
int IsFine, iList, iLoopInv, nLoopInv;
|
|
int BufMx[9];
|
|
const T_RTMx *lsmx;
|
|
int nTrV, iTrV;
|
|
const int *TrV;
|
|
|
|
|
|
nLoopInv = Sg_nLoopInv(SgInfo);
|
|
|
|
nTrV = SgInfo->LatticeInfo->nTrVector;
|
|
|
|
switch (SgInfo->LatticeInfo->Code) {
|
|
default:
|
|
case 'P':
|
|
lx = ly = lz = 12;
|
|
break;
|
|
case 'A':
|
|
lx = ly = 12;
|
|
lz = 6;
|
|
break;
|
|
case 'B':
|
|
ly = lz = 12;
|
|
lx = 6;
|
|
break;
|
|
case 'C':
|
|
lz = lx = 12;
|
|
ly = 6;
|
|
break;
|
|
case 'I':
|
|
lx = ly = 12;
|
|
lz = 6;
|
|
break;
|
|
case 'R':
|
|
lx = ly = 12;
|
|
lz = 4;
|
|
break;
|
|
case 'S':
|
|
lz = lx = 12;
|
|
ly = 4;
|
|
break;
|
|
case 'T':
|
|
ly = lz = 12;
|
|
lx = 4;
|
|
break;
|
|
case 'F':
|
|
lx = 12;
|
|
ly = lz = 6;
|
|
break;
|
|
}
|
|
|
|
for (O[0] = 0; O[0] < 12; O[0]++)
|
|
for (O[1] = 0; O[1] < 12; O[1]++)
|
|
for (O[2] = 0; O[2] < 12; O[2]++) {
|
|
IsFine = 1;
|
|
|
|
for (iList = 0; IsFine && iList < SgInfo->nList; iList++) {
|
|
lsmx = &SgInfo->ListSeitzMx[iList];
|
|
|
|
for (iLoopInv = 0; IsFine && iLoopInv < nLoopInv; iLoopInv++) {
|
|
if (iLoopInv == 0)
|
|
for (i = 0; i < 9; i++) {
|
|
if (i % 4)
|
|
BufMx[i] = lsmx->s.R[i];
|
|
else
|
|
BufMx[i] = lsmx->s.R[i] - 1;
|
|
} else
|
|
for (i = 0; i < 9; i++) {
|
|
if (i % 4)
|
|
BufMx[i] = -lsmx->s.R[i];
|
|
else
|
|
BufMx[i] = -lsmx->s.R[i] - 1;
|
|
}
|
|
|
|
RotMx_t_Vector(V, BufMx, O, 12);
|
|
|
|
TrV = SgInfo->LatticeInfo->TrVector;
|
|
|
|
for (iTrV = 0; iTrV < nTrV; iTrV++) {
|
|
mx = (V[0] * (STBF / 12) + *TrV++) % STBF;
|
|
my = (V[1] * (STBF / 12) + *TrV++) % STBF;
|
|
mz = (V[2] * (STBF / 12) + *TrV++) % STBF;
|
|
|
|
if (mx == 0 && my == 0 && mz == 0)
|
|
break;
|
|
}
|
|
|
|
if (iTrV == nTrV)
|
|
IsFine = 0;
|
|
}
|
|
}
|
|
|
|
if (!(O[0] < lx && O[1] < ly && O[2] < lz))
|
|
IsFine = -IsFine;
|
|
|
|
*TestField++ = IsFine;
|
|
|
|
#if DEBUG_MarkLegalOrigins
|
|
if (IsFine == 1)
|
|
putc(' ', stdout);
|
|
else if (IsFine == -1)
|
|
putc('#', stdout);
|
|
if (IsFine != 0)
|
|
fprintf(stdout, " %2d %2d %2d\n", O[0], O[1], O[2]);
|
|
#endif
|
|
}
|
|
}
|
|
|
|
|
|
#define IsArbitraryShift(iShift) \
|
|
( (iShift) == 1 || (iShift) == 5 \
|
|
|| (iShift) == 7 || (iShift) == 11)
|
|
|
|
|
|
int Verify_si(int h, int k, int l, const int *TestField)
|
|
{
|
|
int O[3], TH;
|
|
|
|
|
|
for (O[0] = 0; O[0] < 12; O[0]++)
|
|
for (O[1] = 0; O[1] < 12; O[1]++)
|
|
for (O[2] = 0; O[2] < 12; O[2]++) {
|
|
if (*TestField++) {
|
|
TH = h * O[0] + k * O[1] + l * O[2];
|
|
TH %= 12;
|
|
if (TH)
|
|
return 0;
|
|
|
|
if (IsArbitraryShift(O[0]))
|
|
TH += h;
|
|
if (IsArbitraryShift(O[1]))
|
|
TH += k;
|
|
if (IsArbitraryShift(O[2]))
|
|
TH += l;
|
|
if (TH)
|
|
return 0;
|
|
}
|
|
}
|
|
|
|
return 1;
|
|
}
|
|
|
|
|
|
int Is_si(const T_SgInfo * SgInfo, int h, int k, int l)
|
|
{
|
|
int i_si_v, u;
|
|
const int *si_v, *si_m;
|
|
|
|
|
|
si_v = SgInfo->si_Vector;
|
|
si_m = SgInfo->si_Modulus;
|
|
|
|
for (i_si_v = 0; i_si_v < SgInfo->n_si_Vector; i_si_v++) {
|
|
u = *si_v++ * h;
|
|
u += *si_v++ * k;
|
|
u += *si_v++ * l;
|
|
|
|
if (*si_m) {
|
|
if (u % (*si_m))
|
|
return 0;
|
|
} else {
|
|
if (u)
|
|
return 0;
|
|
}
|
|
|
|
si_m++;
|
|
}
|
|
|
|
return 1;
|
|
}
|
|
|
|
|
|
int Set_si(T_SgInfo * SgInfo)
|
|
{
|
|
static const int TabTrial_si[] = {
|
|
0,
|
|
|
|
1, 0, 2, -1, 4, /* I -4 */
|
|
1, 2, -1, 0, 4,
|
|
1, -1, 0, 2, 4,
|
|
|
|
1, 2, 4, 3, 6, /* P 3 2 */
|
|
1, 4, 3, 2, 6,
|
|
1, 3, 2, 4, 6,
|
|
|
|
1, 1, 1, 1, 4,
|
|
1, 1, 1, 1, 2,
|
|
1, 1, 1, 1, 0,
|
|
|
|
1, 0, 0, 1, 2,
|
|
1, 0, 1, 0, 2,
|
|
1, 1, 0, 0, 2,
|
|
|
|
1, 0, 0, 1, 0,
|
|
1, 0, 1, 0, 0,
|
|
1, 1, 0, 0, 0,
|
|
|
|
2, 1, -1, 0, 3,
|
|
0, 0, 1, 0,
|
|
2, -1, 0, 1, 3,
|
|
0, 1, 0, 0,
|
|
2, 0, 1, -1, 3,
|
|
1, 0, 0, 0,
|
|
|
|
2, 0, 1, 1, 4, /* F 2x */
|
|
1, 0, 0, 0,
|
|
2, 1, 0, 1, 4, /* F 2y */
|
|
0, 1, 0, 0,
|
|
2, 1, 1, 0, 4, /* F 2z */
|
|
0, 0, 1, 0,
|
|
|
|
2, 1, 0, 0, 2,
|
|
0, 0, 1, 2,
|
|
2, 0, 1, 0, 2,
|
|
0, 0, 1, 2,
|
|
2, 1, 0, 0, 2,
|
|
0, 1, 0, 2,
|
|
|
|
2, 1, 1, 0, 2,
|
|
0, 0, 1, 2,
|
|
2, 1, 0, 1, 2,
|
|
0, 1, 0, 2,
|
|
2, 0, 1, 1, 2,
|
|
1, 0, 0, 2,
|
|
|
|
2, 1, 0, 0, 2,
|
|
0, 0, 1, 0,
|
|
2, 0, 1, 0, 2,
|
|
0, 0, 1, 0,
|
|
2, 1, 0, 0, 2,
|
|
0, 1, 0, 0,
|
|
|
|
2, 1, 0, 0, 0,
|
|
0, 0, 1, 2,
|
|
2, 0, 1, 0, 0,
|
|
0, 0, 1, 2,
|
|
2, 1, 0, 0, 0,
|
|
0, 1, 0, 2,
|
|
|
|
2, 1, 1, 0, 2,
|
|
0, 0, 1, 0,
|
|
2, 1, 0, 1, 2,
|
|
0, 1, 0, 0,
|
|
2, 0, 1, 1, 2,
|
|
1, 0, 0, 0,
|
|
|
|
2, 1, 0, 0, 0,
|
|
0, 0, 1, 0,
|
|
2, 0, 1, 0, 0,
|
|
0, 0, 1, 0,
|
|
2, 1, 0, 0, 0,
|
|
0, 1, 0, 0,
|
|
|
|
3, 1, 0, 0, 2,
|
|
0, 1, 0, 2,
|
|
0, 0, 1, 2,
|
|
|
|
3, 1, 0, 0, 0,
|
|
0, 1, 0, 2,
|
|
0, 0, 1, 2,
|
|
|
|
3, 1, 0, 0, 2,
|
|
0, 1, 0, 0,
|
|
0, 0, 1, 2,
|
|
|
|
3, 1, 0, 0, 2,
|
|
0, 1, 0, 2,
|
|
0, 0, 1, 0,
|
|
|
|
3, 1, 0, 0, 2,
|
|
0, 1, 0, 0,
|
|
0, 0, 1, 0,
|
|
|
|
3, 1, 0, 0, 0,
|
|
0, 1, 0, 2,
|
|
0, 0, 1, 0,
|
|
|
|
3, 1, 0, 0, 0,
|
|
0, 1, 0, 0,
|
|
0, 0, 1, 2,
|
|
|
|
3, 1, 0, 0, 0,
|
|
0, 1, 0, 0,
|
|
0, 0, 1, 0,
|
|
|
|
3, -1, 0, 0, 2, /* -A 1 */
|
|
0, -1, 1, 4,
|
|
0, 1, 1, 4,
|
|
|
|
3, -1, 0, 1, 4, /* -B 1 */
|
|
0, -1, 0, 2,
|
|
1, 0, 1, 4,
|
|
|
|
3, 1, 1, 0, 4, /* -C 1 */
|
|
1, -1, 0, 4,
|
|
0, 0, -1, 2,
|
|
|
|
3, -1, 1, 1, 4, /* -I 1 */
|
|
1, -1, 1, 4,
|
|
1, 1, -1, 4,
|
|
|
|
3, 0, 1, 1, 4, /* -F 1 */
|
|
1, 0, 1, 4,
|
|
1, 1, 0, 4,
|
|
|
|
3, -1, 0, 0, 0, /* A 2x */
|
|
0, -1, 1, 4,
|
|
0, 1, 1, 4,
|
|
|
|
3, -1, 0, 1, 4, /* B 2y */
|
|
0, -1, 0, 0,
|
|
1, 0, 1, 4,
|
|
|
|
3, 1, 1, 0, 4, /* C 2z */
|
|
1, -1, 0, 4,
|
|
0, 0, -1, 0,
|
|
|
|
-1
|
|
};
|
|
|
|
int h, k, l, iList;
|
|
int Maxh, Maxk, Maxl;
|
|
int Minh, Mink, Minl;
|
|
int nTestField, *TestField;
|
|
int nProperty, *Property, *pp;
|
|
int IsFine, would_be, is;
|
|
int i_si, *si_v;
|
|
const int *trial_si;
|
|
|
|
|
|
SgInfo->n_si_Vector = -1;
|
|
|
|
nTestField = 12 * 12 * 12;
|
|
AppMalloc(TestField, nTestField);
|
|
if (TestField == NULL) {
|
|
SetSgError("Not enough core");
|
|
return -1;
|
|
}
|
|
|
|
MarkLegalOrigins(SgInfo, TestField);
|
|
|
|
Maxh = Maxk = Maxl = 7;
|
|
SetListMin_hkl(SgInfo, Maxk, Maxl, &Minh, &Mink, &Minl);
|
|
|
|
nProperty = (Maxh - Minh + 1)
|
|
* (Maxk - Mink + 1)
|
|
* (Maxl - Minl + 1);
|
|
AppMalloc(Property, nProperty);
|
|
if (Property == NULL) {
|
|
SetSgError("Not enough core");
|
|
AppFree(TestField, nTestField);
|
|
return -1;
|
|
}
|
|
|
|
pp = Property;
|
|
for (h = Minh; h <= Maxh; h++)
|
|
for (k = Mink; k <= Maxk; k++)
|
|
for (l = Minl; l <= Maxl; l++) {
|
|
iList = IsSysAbsent_hkl(SgInfo, h, k, l, NULL);
|
|
if (SgError != NULL) {
|
|
AppFree(Property, nProperty);
|
|
AppFree(TestField, nTestField);
|
|
return -1;
|
|
}
|
|
|
|
if (iList == 0)
|
|
*pp++ = Verify_si(h, k, l, TestField);
|
|
else
|
|
*pp++ = -1;
|
|
}
|
|
|
|
trial_si = TabTrial_si;
|
|
while (*trial_si >= 0) {
|
|
SgInfo->n_si_Vector = *trial_si++;
|
|
si_v = SgInfo->si_Vector;
|
|
for (i_si = 0; i_si < SgInfo->n_si_Vector; i_si++) {
|
|
*si_v++ = *trial_si++;
|
|
*si_v++ = *trial_si++;
|
|
*si_v++ = *trial_si++;
|
|
SgInfo->si_Modulus[i_si] = *trial_si++;
|
|
}
|
|
|
|
IsFine = 1;
|
|
|
|
pp = Property;
|
|
for (h = Minh; IsFine && h <= Maxh; h++)
|
|
for (k = Mink; IsFine && k <= Maxk; k++)
|
|
for (l = Minl; IsFine && l <= Maxl; l++) {
|
|
is = *pp++;
|
|
|
|
if (is >= 0) {
|
|
would_be = Is_si(SgInfo, h, k, l);
|
|
if (is != would_be)
|
|
IsFine = 0;
|
|
}
|
|
}
|
|
|
|
if (IsFine) {
|
|
AppFree(Property, nProperty);
|
|
AppFree(TestField, nTestField);
|
|
return 0;
|
|
}
|
|
}
|
|
|
|
SgInfo->n_si_Vector = -1;
|
|
SetSgError("Internal Error: Can't determine s.i. vectors and moduli");
|
|
|
|
AppFree(Property, nProperty);
|
|
AppFree(TestField, nTestField);
|
|
|
|
return -1;
|
|
}
|
|
|
|
|
|
void Set_uvw(const T_SgInfo * SgInfo, int h, int k, int l, int *uvw)
|
|
{
|
|
int i_si_v, u;
|
|
const int *si_v, *si_m;
|
|
|
|
|
|
si_v = SgInfo->si_Vector;
|
|
si_m = SgInfo->si_Modulus;
|
|
|
|
for (i_si_v = 0; i_si_v < SgInfo->n_si_Vector; i_si_v++) {
|
|
u = *si_v++ * h;
|
|
u += *si_v++ * k;
|
|
u += *si_v++ * l;
|
|
|
|
if (*si_m)
|
|
u %= (*si_m);
|
|
si_m++;
|
|
|
|
uvw[i_si_v] = u;
|
|
}
|
|
}
|