From b4e5f07f7340637403abe012df7427a8d0a0fc0e Mon Sep 17 00:00:00 2001 From: Thomas White Date: Fri, 22 Jul 2011 17:57:22 +0200 Subject: Add left coset decomposition --- src/symmetry.c | 220 +++++++++++++++++++++++++++++++++++++++++++++++-- tests/symmetry_check.c | 24 ++++++ 2 files changed, 238 insertions(+), 6 deletions(-) diff --git a/src/symmetry.c b/src/symmetry.c index efd96b5f..e9c7d642 100644 --- a/src/symmetry.c +++ b/src/symmetry.c @@ -198,6 +198,36 @@ static void add_symop(SymOpList *ops, } +/* Add a operation to a SymOpList */ +static void add_copied_op(SymOpList *ops, struct sym_op *copyme) +{ + int n; + signed int *h, *k, *l; + + if ( ops->n_ops == ops->max_ops ) { + ops->max_ops += 16; + alloc_ops(ops); + } + + n = ops->n_ops; + + h = malloc(3*sizeof(signed int)); + k = malloc(3*sizeof(signed int)); + l = malloc(3*sizeof(signed int)); + + memcpy(h, copyme->h, 3*sizeof(signed int)); + memcpy(k, copyme->k, 3*sizeof(signed int)); + memcpy(l, copyme->l, 3*sizeof(signed int)); + + ops->ops[n].h = h; + ops->ops[n].k = k; + ops->ops[n].l = l; + ops->ops[n].order = copyme->order; + + ops->n_ops++; +} + + /** * num_equivs: * @ops: A %SymOpList @@ -974,6 +1004,33 @@ static int is_inversion(const struct sym_op *op) } +static int is_identity(const struct sym_op *op) +{ + if ( (op->h[0]!=1) || (op->h[1]!=0) || (op->h[2]!=0) ) return 0; + if ( (op->k[0]!=0) || (op->k[1]!=1) || (op->k[2]!=0) ) return 0; + if ( (op->l[0]!=0) || (op->l[1]!=0) || (op->l[2]!=1) ) return 0; + return 1; +} + + +static signed int determinant(const struct sym_op *op) +{ + signed int det = 0; + + det += op->h[0] * (op->k[1]*op->l[2] - op->k[2]*op->l[1]); + det -= op->h[1] * (op->k[0]*op->l[2] - op->k[2]*op->l[0]); + det += op->h[2] * (op->k[0]*op->l[1] - op->k[1]*op->l[0]); + + return det; +} + + +/** + * is_centrosymmetric: + * @s: A %SymOpList + * + * Returns: non-zero if @s contains an inversion operation + */ int is_centrosymmetric(const SymOpList *s) { int i, n; @@ -987,27 +1044,178 @@ int is_centrosymmetric(const SymOpList *s) } +static int ops_equal(const struct sym_op *op, + signed int *h, signed int *k, signed int *l) +{ + if ( (op->h[0]!=h[0]) || (op->h[1]!=h[1]) || (op->h[2]!=h[2]) ) + return 0; + if ( (op->k[0]!=k[0]) || (op->k[1]!=k[1]) || (op->k[2]!=k[2]) ) + return 0; + if ( (op->l[0]!=l[0]) || (op->l[1]!=l[1]) || (op->l[2]!=l[2]) ) + return 0; + return 1; +} + + +/* Return true if a*b = ans */ +static int check_mult(const struct sym_op *ans, + const struct sym_op *a, const struct sym_op *b) +{ + signed int *ans_h, *ans_k, *ans_l; + int val; + + ans_h = malloc(3*sizeof(signed int)); + ans_k = malloc(3*sizeof(signed int)); + ans_l = malloc(3*sizeof(signed int)); + + combine_ops(a->h, a->k, a->l, b->h, b->k, b->l, ans_h, ans_k, ans_l); + val = ops_equal(ans, ans_h, ans_k, ans_l); + + free(ans_h); + free(ans_k); + free(ans_l); + + return val; +} + + /** - * get_twins: - * - * Calculate twinning laws. - * + * get_ambiguities: + * @source: The "source" symmetry, a %SymOpList + * @target: The "target" symmetry, a %SymOpList + + * Calculates twinning laws. Returns a %SymOpList containing the twinning + * operators, which are the symmetry operations which can be added to @target + * to generate @source. Only rotations are allowable - no mirrors nor + * inversions. * To count the number of possibilities, use num_ops() on the result. + * + * Returns: A %SymOpList containing the twinning operators */ -SymOpList *get_twins(const SymOpList *source, const SymOpList *target) +SymOpList *get_ambiguities(const SymOpList *source, const SymOpList *target) { int n_src, n_tgt; int i; - SymOpList *twins = new_symoplist(); + SymOpList *twins; + SymOpList *src_reordered; + SymOpMask *used; + char *name; n_src = num_ops(source); n_tgt = num_ops(target); + src_reordered = new_symoplist(); + used = new_symopmask(source); + + /* Find identity */ for ( i=0; imask[i] == 0 ) continue; + if ( is_identity(&source->ops[i]) ) { + add_copied_op(src_reordered, &source->ops[i]); + used->mask[i] = 0; + } + } + /* Find binary options (order=2) of first kind (determinant positive) */ + for ( i=0; imask[i] == 0 ) continue; + if ( (source->ops[i].order == 2) + && (determinant(&source->ops[i]) > 0) ) { + add_copied_op(src_reordered, &source->ops[i]); + used->mask[i] = 0; + } + } + + /* Find other operations of first kind (determinant positive) */ + for ( i=0; imask[i] == 0 ) continue; + if ( determinant(&source->ops[i]) > 0 ) { + add_copied_op(src_reordered, &source->ops[i]); + used->mask[i] = 0; + } + } + + /* Find inversion */ + for ( i=0; imask[i] == 0 ) continue; + if ( is_inversion(&source->ops[i]) ) { + add_copied_op(src_reordered, &source->ops[i]); + used->mask[i] = 0; + } + } + /* Find binary options of second kind (determinant negative) */ + for ( i=0; imask[i] == 0 ) continue; + if ( (source->ops[i].order == 2) + && (determinant(&source->ops[i]) < 0) ) { + add_copied_op(src_reordered, &source->ops[i]); + used->mask[i] = 0; + } } + /* Find other operations of second kind (determinant negative) */ + for ( i=0; imask[i] == 0 ) continue; + if ( determinant(&source->ops[i]) < 0 ) { + add_copied_op(src_reordered, &source->ops[i]); + used->mask[i] = 0; + } + } + + int n_left_over = 0; + for ( i=0; imask[i] == 0 ) continue; + n_left_over++; + } + if ( n_left_over != 0 ) { + ERROR("%i operations left over after rearranging for" + " left coset decomposition.\n", n_left_over); + } + + if ( num_ops(src_reordered) != num_ops(source) ) { + ERROR("%i ops went to %i after rearranging.\n", + num_ops(src_reordered), num_ops(source)); + } + + free_symopmask(used); + used = new_symopmask(src_reordered); + + for ( i=0; imask[i] == 0 ) continue; + + for ( j=1; jops[k], + &src_reordered->ops[i], + &target->ops[j]) ) + { + used->mask[k] = 0; + } + } + + } + + } + + twins = new_symoplist(); + for ( i=0; imask[i] == 0 ) continue; + add_copied_op(twins, &src_reordered->ops[i]); + } + + free_symopmask(used); + free_symoplist(src_reordered); + + name = malloc(64); + snprintf(name, 63, "%s -> %s", symmetry_name(source), + symmetry_name(target)); + twins->name = name; + return twins; } diff --git a/tests/symmetry_check.c b/tests/symmetry_check.c index 5a2c91d5..6934a3d9 100644 --- a/tests/symmetry_check.c +++ b/tests/symmetry_check.c @@ -65,6 +65,28 @@ static void check_pg_props(const char *pg, int answer, int centro, int *fail) } +static void check_subgroup(const char *ssource, const char *starget, + int *fail) +{ + SymOpList *source; + SymOpList *target; + SymOpList *twins; + + source = get_pointgroup(ssource); + target = get_pointgroup(starget); + if ( (source == NULL) || (target == NULL) ) { + *fail = 1; + return; + } + + twins = get_ambiguities(source, target); + describe_symmetry(twins); + + free_symoplist(target); + free_symoplist(source); +} + + int main(int argc, char *argv[]) { int fail = 0; @@ -127,5 +149,7 @@ int main(int argc, char *argv[]) check_pg_props( "m-3m", 48, 1, &fail); STATUS("\n"); + check_subgroup("2/m", "m", &fail); + return fail; } -- cgit v1.2.3