added license information to all files

updated documentation
added example file
added graphene structure
This commit is contained in:
pmla 2019-02-15 18:49:23 -05:00
parent a8cee136fe
commit 35d2de298f
35 changed files with 2856 additions and 1524 deletions

View File

@ -38,20 +38,20 @@ ico (icosahedral) = 4
sc (simple cubic) = 5
dcub (diamond cubic) = 6
dhex (diamond hexagonal) = 7
other = 8 :ul
graphene = 8 :ul
The value of the PTM structure will be 0 for atoms not in the specified
The value of the PTM structure will be 0 for unknown types and -1 for atoms not in the specified
compute group. The choice of structures to search for can be specified using the "structures"
argument, which is a hyphen-separated list of structure keywords.
Two convenient pre-set options are provided:
default: fcc-hcp-bcc-ico
all: fcc-hcp-bcc-ico-sc-dcub-dhex :ul
all: fcc-hcp-bcc-ico-sc-dcub-dhex-graphene :ul
The 'default' setting detects the same structures as the Common Neighbor Analysis method.
The 'all' setting searches for all structure types. A small performance penalty is
incurred for the diamond structures, so it is not recommended to use this option if
it is known that the simulation does not contain diamond structures.
The 'all' setting searches for all structure types. A performance penalty is
incurred for the diamond and graphene structures, so it is not recommended to use this option if
it is known that the simulation does not contain these structures.
PTM identifies structures using two steps. First, a graph isomorphism test is used

View File

@ -0,0 +1,47 @@
# LAMMPS Input File for Grain Boundaries
# Mark Tschopp, Dec2009
# Modified by PM Larsen to demonstrate the use of PTM in LAMMPS
# This file will generate a single Sigma5(310) STGB and run PTM
# ---------- Initialize Simulation ---------------------
clear
units metal
dimension 3
boundary p p p
atom_style atomic
# ---------- Create Sigma5(310) Structure ---------
lattice fcc 4.05
region whole block 0.000000 12.807225 -64.0361225 64.0361225 0.000000 4.050000 units box
create_box 2 whole
region upper block INF INF 0.000000 64.0361225 INF INF units box
lattice fcc 4.05 orient x 0 3 1 orient y 0 -1 3 orient z 1 0 0
create_atoms 1 region upper
region lower block INF INF -64.0361225 0.000000 INF INF units box
lattice fcc 4.05 orient x 0 3 -1 orient y 0 1 3 orient z 1 0 0
create_atoms 2 region lower
group upper type 1
group lower type 2
mass 1 1.0
mass 2 1.0
# ---------- Define Interatomic Potential ---------------------
pair_style lj/cut 2.5
pair_coeff * * 1 1
pair_coeff 1 1 1 1.1 2.8
# ---------- Displace atoms and delete overlapping atoms ---------------------
displace_atoms upper move 0 0 0 units lattice
delete_atoms overlap 0.35 lower upper
# ---------- Define PTM settings (default structures, RMSD threshold of 0.1) -------------------
compute ptm all ptm/atom default 0.1
# ---------- Dump data into Data file -------------
reset_timestep 0
dump 1 all cfg 10000 dump.ptm_example_*.cfg mass type xs ys zs c_ptm[1] c_ptm[2] c_ptm[3] c_ptm[4] c_ptm[5] c_ptm[6] c_ptm[7]
dump_modify 1 element Al Al
run 1
print "All done"

View File

@ -36,10 +36,9 @@ under
#include "ptm_functions.h"
#define MAX_NEIGHBORS 30
#define NUM_COLUMNS 7
#define UNKNOWN 0
#define OTHER 8
#define PTM_LAMMPS_UNKNOWN -1
#define PTM_LAMMPS_OTHER 0
using namespace LAMMPS_NS;
@ -69,7 +68,9 @@ ComputePTMAtom::ComputePTMAtom(LAMMPS *lmp, int narg, char **arg)
char *ptr = structures;
const char *strings[] = {"fcc", "hcp", "bcc", "ico", "sc",
"dcub", "dhex", "all", "default"};
"dcub", "dhex", "graphene", "all", "default"};
int num_strings = sizeof(strings) / sizeof(const char*);
int32_t flags[] = {
PTM_CHECK_FCC,
PTM_CHECK_HCP,
@ -78,6 +79,7 @@ ComputePTMAtom::ComputePTMAtom(LAMMPS *lmp, int narg, char **arg)
PTM_CHECK_SC,
PTM_CHECK_DCUB,
PTM_CHECK_DHEX,
PTM_CHECK_GRAPHENE,
PTM_CHECK_ALL,
PTM_CHECK_FCC | PTM_CHECK_HCP | PTM_CHECK_BCC | PTM_CHECK_ICO};
@ -85,7 +87,7 @@ ComputePTMAtom::ComputePTMAtom(LAMMPS *lmp, int narg, char **arg)
while (*ptr != '\0') {
bool found = false;
for (int i = 0; i < 9; i++) {
for (int i = 0; i < num_strings; i++) {
int len = strlen(strings[i]);
if (strncmp(ptr, strings[i], len) == 0) {
input_flags |= flags[i];
@ -156,47 +158,83 @@ void ComputePTMAtom::init_list(int id, NeighList *ptr) { list = ptr; }
/* ---------------------------------------------------------------------- */
typedef struct
{
double **x;
int *numneigh;
int **firstneigh;
int *ilist;
int nlocal;
} ptmnbrdata_t;
typedef struct {
int index;
double d;
} ptmnbr_t;
static bool sorthelper_compare(ptmnbr_t const &a, ptmnbr_t const &b) {
return a.d < b.d;
}
static int get_neighbors(double *pos, int jnum, int *jlist, double **x,
double (*nbr)[3]) {
static int get_neighbours(void* vdata, size_t central_index, size_t atom_index, int num, size_t* nbr_indices, int32_t* numbers, double (*nbr_pos)[3])
{
ptmnbrdata_t* data = (ptmnbrdata_t*)vdata;
double **x = data->x;
double *pos = x[atom_index];
int *jlist = NULL;
int jnum = 0;
if (atom_index < data->nlocal) {
jlist = data->firstneigh[atom_index];
jnum = data->numneigh[atom_index];
}
else {
jlist = data->firstneigh[central_index];
jnum = data->numneigh[central_index];
}
ptmnbr_t *nbr_order = new ptmnbr_t[jnum];
int count = 0;
for (int jj = 0; jj < jnum; jj++) {
int j = jlist[jj];
j &= NEIGHMASK;
if (j == atom_index)
continue;
double dx = pos[0] - x[j][0];
double dy = pos[1] - x[j][1];
double dz = pos[2] - x[j][2];
double rsq = dx * dx + dy * dy + dz * dz;
nbr_order[jj].index = j;
nbr_order[jj].d = rsq;
nbr_order[count].index = j;
nbr_order[count].d = rsq;
count++;
}
std::sort(nbr_order, nbr_order + jnum, &sorthelper_compare);
int num_nbrs = std::min(MAX_NEIGHBORS, jnum);
std::sort(nbr_order, nbr_order + count, &sorthelper_compare);
int num_nbrs = std::min(num - 1, count);
nbr[0][0] = nbr[0][1] = nbr[0][2] = 0;
nbr_pos[0][0] = nbr_pos[0][1] = nbr_pos[0][2] = 0;
nbr_indices[0] = atom_index;
numbers[0] = 0;
for (int jj = 0; jj < num_nbrs; jj++) {
int j = nbr_order[jj].index;
nbr[jj + 1][0] = x[j][0] - pos[0];
nbr[jj + 1][1] = x[j][1] - pos[1];
nbr[jj + 1][2] = x[j][2] - pos[2];
nbr_pos[jj + 1][0] = x[j][0] - pos[0];
nbr_pos[jj + 1][1] = x[j][1] - pos[1];
nbr_pos[jj + 1][2] = x[j][2] - pos[2];
nbr_indices[jj + 1] = j;
numbers[jj + 1] = 0;
}
delete[] nbr_order;
return num_nbrs;
return num_nbrs + 1;
}
void ComputePTMAtom::compute_peratom() {
@ -229,50 +267,29 @@ void ComputePTMAtom::compute_peratom() {
double **x = atom->x;
int *mask = atom->mask;
int nlocal = atom->nlocal;
ptmnbrdata_t nbrlist = {x, numneigh, firstneigh, ilist, atom->nlocal};
for (int ii = 0; ii < inum; ii++) {
int i = ilist[ii];
output[i][0] = UNKNOWN;
output[i][0] = PTM_LAMMPS_UNKNOWN;
if (!(mask[i] & groupbit))
continue;
double *pos = x[i];
int *jlist = firstneigh[i];
int jnum = numneigh[i];
if (jnum <= 0)
continue;
// get neighbours ordered by increasing distance
double nbr[MAX_NEIGHBORS + 1][3];
int num_nbrs = get_neighbors(pos, jnum, jlist, x, nbr);
// check that we have enough neighbours for the desired structure types
int32_t flags = 0;
if (num_nbrs >= PTM_NUM_NBRS_SC && (input_flags & PTM_CHECK_SC))
flags |= PTM_CHECK_SC;
if (num_nbrs >= PTM_NUM_NBRS_FCC && (input_flags & PTM_CHECK_FCC))
flags |= PTM_CHECK_FCC;
if (num_nbrs >= PTM_NUM_NBRS_HCP && (input_flags & PTM_CHECK_HCP))
flags |= PTM_CHECK_HCP;
if (num_nbrs >= PTM_NUM_NBRS_ICO && (input_flags & PTM_CHECK_ICO))
flags |= PTM_CHECK_ICO;
if (num_nbrs >= PTM_NUM_NBRS_BCC && (input_flags & PTM_CHECK_BCC))
flags |= PTM_CHECK_BCC;
if (num_nbrs >= PTM_NUM_NBRS_DCUB && (input_flags & PTM_CHECK_DCUB))
flags |= PTM_CHECK_DCUB;
if (num_nbrs >= PTM_NUM_NBRS_DHEX && (input_flags & PTM_CHECK_DHEX))
flags |= PTM_CHECK_DHEX;
// now run PTM
int8_t mapping[MAX_NEIGHBORS + 1];
int32_t type, alloy_type;
double scale, rmsd, interatomic_distance, lattice_constant;
double q[4], F[9], F_res[3], U[9], P[9];
ptm_index(local_handle, flags, num_nbrs + 1, nbr, NULL, true, &type,
&alloy_type, &scale, &rmsd, q, F, F_res, U, P, mapping,
&interatomic_distance, &lattice_constant);
double scale, rmsd, interatomic_distance;
double q[4];
bool standard_orientations = false;
ptm_index(local_handle, i, get_neighbours, (void*)&nbrlist,
input_flags, standard_orientations,
&type, &alloy_type, &scale, &rmsd, q,
NULL, NULL, NULL, NULL, &interatomic_distance, NULL, NULL);
if (rmsd > rmsd_threshold) {
type = PTM_MATCH_NONE;
@ -280,8 +297,10 @@ void ComputePTMAtom::compute_peratom() {
// printf("%d type=%d rmsd=%f\n", i, type, rmsd);
if (type == PTM_MATCH_NONE)
type = OTHER;
if (type == PTM_MATCH_NONE) {
type = PTM_LAMMPS_OTHER;
rmsd = INFINITY;
}
output[i][0] = type;
output[i][1] = rmsd;

View File

@ -1,3 +1,12 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#include <algorithm>
#include "ptm_constants.h"
#include "ptm_initialize_data.h"
@ -6,8 +15,8 @@ namespace ptm {
#define NUM_ALLOY_TYPES 3
static uint32_t typedata[NUM_ALLOY_TYPES][3] = {
{PTM_MATCH_FCC, PTM_ALLOY_L10, 0x000001fe},
{PTM_MATCH_FCC, PTM_ALLOY_L12_CU, 0x0000001e},
{PTM_MATCH_FCC, PTM_ALLOY_L10, 0x00000db6},
{PTM_MATCH_FCC, PTM_ALLOY_L12_CU, 0x00000492},
{PTM_MATCH_FCC, PTM_ALLOY_L12_AU, 0x00001ffe},
};
@ -78,6 +87,10 @@ static int32_t canonical_alloy_representation(const refdata_t* ref, int8_t* mapp
int32_t find_alloy_type(const refdata_t* ref, int8_t* mapping, int32_t* numbers)
{
for (int i=0;i<ref->num_nbrs+1;i++)
if (numbers[i] == -1)
return PTM_ALLOY_NONE;
if (test_pure(ref->num_nbrs, numbers))
return PTM_ALLOY_PURE;
@ -97,6 +110,11 @@ int32_t find_alloy_type(const refdata_t* ref, int8_t* mapping, int32_t* numbers)
if (test_shell_structure(ref, mapping, numbers, 4))
return PTM_ALLOY_SIC;
if (ref->type == PTM_MATCH_GRAPHENE)
if (test_shell_structure(ref, mapping, numbers, 3))
return PTM_ALLOY_BN;
return PTM_ALLOY_NONE;
}

View File

@ -1,3 +1,12 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#ifndef PTM_ALLOY_TYPES_H
#define PTM_ALLOY_TYPES_H

View File

@ -0,0 +1,127 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#ifndef PTM_ALT_TEMPLATES_H
#define PTM_ALT_TEMPLATES_H
#include <cmath>
const double ptm_template_hcp_alt1[PTM_NUM_POINTS_HCP][3] = {
{ 0, 0, 0 },
{ 1, 0, 0 },
{ -0.5, -sqrt(3)/2, 0 },
{ -0.5, -sqrt(3)/6, -sqrt(6)/3 },
{ 0, sqrt(3)/3, -sqrt(6)/3 },
{ 0.5, -sqrt(3)/6, -sqrt(6)/3 },
{ -1, 0, 0 },
{ -0.5, sqrt(3)/2, 0 },
{ 0.5, sqrt(3)/2, 0 },
{ 0.5, -sqrt(3)/2, 0 },
{ 0.5, -sqrt(3)/6, sqrt(6)/3 },
{ 0, sqrt(3)/3, sqrt(6)/3 },
{ -0.5, -sqrt(3)/6, sqrt(6)/3 },
};
const double ptm_template_dcub_alt1[PTM_NUM_POINTS_DCUB][3] = {
{ 0, 0, 0 },
{ 4/(sqrt(3)+6*sqrt(2)), -4/(sqrt(3)+6*sqrt(2)), 4/(sqrt(3)+6*sqrt(2)) },
{ 4/(sqrt(3)+6*sqrt(2)), 4/(sqrt(3)+6*sqrt(2)), -4/(sqrt(3)+6*sqrt(2)) },
{ -4/(sqrt(3)+6*sqrt(2)), -4/(sqrt(3)+6*sqrt(2)), -4/(sqrt(3)+6*sqrt(2)) },
{ -4/(sqrt(3)+6*sqrt(2)), 4/(sqrt(3)+6*sqrt(2)), 4/(sqrt(3)+6*sqrt(2)) },
{ 8/(sqrt(3)+6*sqrt(2)), 0, 8/(sqrt(3)+6*sqrt(2)) },
{ 0, -8/(sqrt(3)+6*sqrt(2)), 8/(sqrt(3)+6*sqrt(2)) },
{ 8/(sqrt(3)+6*sqrt(2)), -8/(sqrt(3)+6*sqrt(2)), 0 },
{ 0, 8/(sqrt(3)+6*sqrt(2)), -8/(sqrt(3)+6*sqrt(2)) },
{ 8/(sqrt(3)+6*sqrt(2)), 0, -8/(sqrt(3)+6*sqrt(2)) },
{ 8/(sqrt(3)+6*sqrt(2)), 8/(sqrt(3)+6*sqrt(2)), 0 },
{ -8/(sqrt(3)+6*sqrt(2)), 0, -8/(sqrt(3)+6*sqrt(2)) },
{ 0, -8/(sqrt(3)+6*sqrt(2)), -8/(sqrt(3)+6*sqrt(2)) },
{ -8/(sqrt(3)+6*sqrt(2)), -8/(sqrt(3)+6*sqrt(2)), 0 },
{ -8/(sqrt(3)+6*sqrt(2)), 8/(sqrt(3)+6*sqrt(2)), 0 },
{ -8/(sqrt(3)+6*sqrt(2)), 0, 8/(sqrt(3)+6*sqrt(2)) },
{ 0, 8/(sqrt(3)+6*sqrt(2)), 8/(sqrt(3)+6*sqrt(2)) },
};
const double ptm_template_dhex_alt1[PTM_NUM_POINTS_DHEX][3] = {
{ 0, 0, 0 },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), -4*sqrt(3)/(3*sqrt(3)+18*sqrt(2)) },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), -4*sqrt(3)/(3*sqrt(3)+18*sqrt(2)) },
{ 0, 8*sqrt(6)/(3*sqrt(3)+18*sqrt(2)), -4*sqrt(3)/(3*sqrt(3)+18*sqrt(2)) },
{ 0, 0, 4*sqrt(3)/(sqrt(3)+6*sqrt(2)) },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(sqrt(3)+6*sqrt(2)), 0 },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), -16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ -8*sqrt(2)/(sqrt(3)+6*sqrt(2)), 0, 0 },
{ 8*sqrt(2)/(sqrt(3)+6*sqrt(2)), 0, 0 },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), -16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(sqrt(3)+6*sqrt(2)), 0 },
{ 0, 8*sqrt(6)/(3*sqrt(3)+18*sqrt(2)), -16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(sqrt(3)+6*sqrt(2)), 0 },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(sqrt(3)+6*sqrt(2)), 0 },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), 16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ 0, 8*sqrt(6)/(3*sqrt(3)+18*sqrt(2)), 16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), 16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
};
const double ptm_template_dhex_alt2[PTM_NUM_POINTS_DHEX][3] = {
{ 0, 0, 0 },
{ 0, -8*sqrt(6)/(3*sqrt(3)+18*sqrt(2)), 4*sqrt(3)/(3*sqrt(3)+18*sqrt(2)) },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), 4*sqrt(3)/(3*sqrt(3)+18*sqrt(2)) },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), 4*sqrt(3)/(3*sqrt(3)+18*sqrt(2)) },
{ 0, 0, -4*sqrt(3)/(sqrt(3)+6*sqrt(2)) },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(sqrt(3)+6*sqrt(2)), 0 },
{ 0, -8*sqrt(6)/(3*sqrt(3)+18*sqrt(2)), 16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(sqrt(3)+6*sqrt(2)), 0 },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(sqrt(3)+6*sqrt(2)), 0 },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), 16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ -8*sqrt(2)/(sqrt(3)+6*sqrt(2)), 0, 0 },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), 16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ 8*sqrt(2)/(sqrt(3)+6*sqrt(2)), 0, 0 },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(sqrt(3)+6*sqrt(2)), 0 },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), -16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), -16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ 0, -8*sqrt(6)/(3*sqrt(3)+18*sqrt(2)), -16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
};
const double ptm_template_dhex_alt3[PTM_NUM_POINTS_DHEX][3] = {
{ 0, 0, 0 },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), 4*sqrt(3)/(3*sqrt(3)+18*sqrt(2)) },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), 4*sqrt(3)/(3*sqrt(3)+18*sqrt(2)) },
{ 0, 8*sqrt(6)/(3*sqrt(3)+18*sqrt(2)), 4*sqrt(3)/(3*sqrt(3)+18*sqrt(2)) },
{ 0, 0, -4*sqrt(3)/(sqrt(3)+6*sqrt(2)) },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(sqrt(3)+6*sqrt(2)), 0 },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), 16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ 8*sqrt(2)/(sqrt(3)+6*sqrt(2)), 0, 0 },
{ -8*sqrt(2)/(sqrt(3)+6*sqrt(2)), 0, 0 },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), 16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(sqrt(3)+6*sqrt(2)), 0 },
{ 0, 8*sqrt(6)/(3*sqrt(3)+18*sqrt(2)), 16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(sqrt(3)+6*sqrt(2)), 0 },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(sqrt(3)+6*sqrt(2)), 0 },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), -16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ 0, 8*sqrt(6)/(3*sqrt(3)+18*sqrt(2)), -16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), -16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
};
const double ptm_template_graphene_alt1[PTM_NUM_POINTS_GRAPHENE][3] = {
{ 0, 0, 0 },
{ 3*sqrt(3)/22-9./11, -3./22+3*sqrt(3)/11, 0 },
{ 9./11-3*sqrt(3)/22, -3./22+3*sqrt(3)/11, 0 },
{ 0, -6*sqrt(3)/11+3./11, 0 },
{ -18./11+3*sqrt(3)/11, 0, 0 },
{ 3*sqrt(3)/22-9./11, -9./22+9*sqrt(3)/11, 0 },
{ 9./11-3*sqrt(3)/22, -9./22+9*sqrt(3)/11, 0 },
{ -3*sqrt(3)/11+18./11, 0, 0 },
{ 9./11-3*sqrt(3)/22, -9*sqrt(3)/11+9./22, 0 },
{ 3*sqrt(3)/22-9./11, -9*sqrt(3)/11+9./22, 0 },
};
#endif

View File

@ -1,4 +1,13 @@
#include <cstring>
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#include <string.h>
#include <climits>
#include <algorithm>
#include "ptm_graph_tools.h"

View File

@ -1,3 +1,12 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#ifndef PTM_CANONICAL_COLOURED_H
#define PTM_CANONICAL_COLOURED_H

View File

@ -1,30 +1,43 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#ifndef PTM_CONSTANTS_H
#define PTM_CONSTANTS_H
#include <cmath>
//------------------------------------
// definitions
//------------------------------------
#define PTM_NO_ERROR 0
#define PTM_CHECK_FCC (1 << 0)
#define PTM_CHECK_HCP (1 << 1)
#define PTM_CHECK_BCC (1 << 2)
#define PTM_CHECK_ICO (1 << 3)
#define PTM_CHECK_SC (1 << 4)
#define PTM_CHECK_DCUB (1 << 5)
#define PTM_CHECK_DHEX (1 << 6)
#define PTM_CHECK_NONDIAMOND (PTM_CHECK_SC | PTM_CHECK_FCC | PTM_CHECK_HCP | PTM_CHECK_ICO | PTM_CHECK_BCC)
#define PTM_CHECK_ALL (PTM_CHECK_SC | PTM_CHECK_FCC | PTM_CHECK_HCP | PTM_CHECK_ICO | PTM_CHECK_BCC | PTM_CHECK_DCUB | PTM_CHECK_DHEX)
#define PTM_CHECK_FCC (1 << 0)
#define PTM_CHECK_HCP (1 << 1)
#define PTM_CHECK_BCC (1 << 2)
#define PTM_CHECK_ICO (1 << 3)
#define PTM_CHECK_SC (1 << 4)
#define PTM_CHECK_DCUB (1 << 5)
#define PTM_CHECK_DHEX (1 << 6)
#define PTM_CHECK_GRAPHENE (1 << 7)
#define PTM_CHECK_DEFAULT (PTM_CHECK_FCC | PTM_CHECK_HCP | PTM_CHECK_ICO | PTM_CHECK_BCC)
#define PTM_CHECK_ALL (PTM_CHECK_SC | PTM_CHECK_FCC | PTM_CHECK_HCP | PTM_CHECK_ICO | PTM_CHECK_BCC | PTM_CHECK_DCUB | PTM_CHECK_DHEX | PTM_CHECK_GRAPHENE)
#define PTM_MATCH_NONE 0
#define PTM_MATCH_FCC 1
#define PTM_MATCH_HCP 2
#define PTM_MATCH_BCC 3
#define PTM_MATCH_ICO 4
#define PTM_MATCH_SC 5
#define PTM_MATCH_DCUB 6
#define PTM_MATCH_DHEX 7
#define PTM_MATCH_NONE 0
#define PTM_MATCH_FCC 1
#define PTM_MATCH_HCP 2
#define PTM_MATCH_BCC 3
#define PTM_MATCH_ICO 4
#define PTM_MATCH_SC 5
#define PTM_MATCH_DCUB 6
#define PTM_MATCH_DHEX 7
#define PTM_MATCH_GRAPHENE 8
#define PTM_ALLOY_NONE 0
#define PTM_ALLOY_PURE 1
@ -33,9 +46,10 @@
#define PTM_ALLOY_L12_AU 4
#define PTM_ALLOY_B2 5
#define PTM_ALLOY_SIC 6
#define PTM_ALLOY_BN 7
#define PTM_MAX_INPUT_POINTS 35
#define PTM_MAX_INPUT_POINTS 19
#define PTM_MAX_NBRS 16
#define PTM_MAX_POINTS (PTM_MAX_NBRS + 1)
#define PTM_MAX_FACETS 28 //2 * PTM_MAX_NBRS - 4
@ -52,16 +66,18 @@
#define PTM_NUM_NBRS_SC 6
#define PTM_NUM_NBRS_DCUB 16
#define PTM_NUM_NBRS_DHEX 16
#define PTM_NUM_NBRS_GRAPHENE 9
#define PTM_NUM_POINTS_FCC (PTM_NUM_NBRS_FCC + 1)
#define PTM_NUM_POINTS_HCP (PTM_NUM_NBRS_HCP + 1)
#define PTM_NUM_POINTS_BCC (PTM_NUM_NBRS_BCC + 1)
#define PTM_NUM_POINTS_ICO (PTM_NUM_NBRS_ICO + 1)
#define PTM_NUM_POINTS_SC (PTM_NUM_NBRS_SC + 1)
#define PTM_NUM_POINTS_DCUB (PTM_NUM_NBRS_DCUB + 1)
#define PTM_NUM_POINTS_DHEX (PTM_NUM_NBRS_DHEX + 1)
#define PTM_NUM_POINTS_DCUB (PTM_NUM_NBRS_DCUB + 1)
#define PTM_NUM_POINTS_DHEX (PTM_NUM_NBRS_DHEX + 1)
#define PTM_NUM_POINTS_GRAPHENE (PTM_NUM_NBRS_GRAPHENE + 1)
const int ptm_num_nbrs[8] = {0, PTM_NUM_NBRS_FCC, PTM_NUM_NBRS_HCP, PTM_NUM_NBRS_BCC, PTM_NUM_NBRS_ICO, PTM_NUM_NBRS_SC, PTM_NUM_NBRS_DCUB, PTM_NUM_NBRS_DHEX};
const int ptm_num_nbrs[9] = {0, PTM_NUM_NBRS_FCC, PTM_NUM_NBRS_HCP, PTM_NUM_NBRS_BCC, PTM_NUM_NBRS_ICO, PTM_NUM_NBRS_SC, PTM_NUM_NBRS_DCUB, PTM_NUM_NBRS_DHEX, PTM_NUM_NBRS_GRAPHENE};
//------------------------------------
// template structures
@ -69,106 +85,134 @@ const int ptm_num_nbrs[8] = {0, PTM_NUM_NBRS_FCC, PTM_NUM_NBRS_HCP, PTM_NUM_NBRS
//these point sets have barycentre {0, 0, 0} and are scaled such that the mean neighbour distance is 1
const double ptm_template_fcc[PTM_NUM_POINTS_FCC][3] = { { 0. , 0. , 0. },
{ 0. , 0.707106781187, 0.707106781187 },
{ 0. , -0.707106781187, -0.707106781187 },
{ 0. , 0.707106781187, -0.707106781187 },
{ 0. , -0.707106781187, 0.707106781187 },
{ 0.707106781187, 0. , 0.707106781187 },
{ -0.707106781187, 0. , -0.707106781187 },
{ 0.707106781187, 0. , -0.707106781187 },
{ -0.707106781187, 0. , 0.707106781187 },
{ 0.707106781187, 0.707106781187, 0. },
{ -0.707106781187, -0.707106781187, 0. },
{ 0.707106781187, -0.707106781187, 0. },
{ -0.707106781187, 0.707106781187, 0. } };
const double ptm_template_fcc[PTM_NUM_POINTS_FCC][3] = {
{ 0, 0, 0 },
{ sqrt(2)/2, sqrt(2)/2, 0 },
{ 0, sqrt(2)/2, sqrt(2)/2 },
{ sqrt(2)/2, 0, sqrt(2)/2 },
{ -sqrt(2)/2, -sqrt(2)/2, 0 },
{ 0, -sqrt(2)/2, -sqrt(2)/2 },
{ -sqrt(2)/2, 0, -sqrt(2)/2 },
{ -sqrt(2)/2, sqrt(2)/2, 0 },
{ 0, -sqrt(2)/2, sqrt(2)/2 },
{ -sqrt(2)/2, 0, sqrt(2)/2 },
{ sqrt(2)/2, -sqrt(2)/2, 0 },
{ 0, sqrt(2)/2, -sqrt(2)/2 },
{ sqrt(2)/2, 0, -sqrt(2)/2 },
};
const double ptm_template_hcp[PTM_NUM_POINTS_HCP][3] = { { 0. , 0. , 0. },
{ 0.707106781186, 0. , 0.707106781186 },
{ -0.235702260395, -0.942809041583, -0.235702260395 },
{ 0.707106781186, 0.707106781186, 0. },
{ -0.235702260395, -0.235702260395, -0.942809041583 },
{ 0. , 0.707106781186, 0.707106781186 },
{ -0.942809041583, -0.235702260395, -0.235702260395 },
{ -0.707106781186, 0.707106781186, 0. },
{ 0. , 0.707106781186, -0.707106781186 },
{ 0.707106781186, 0. , -0.707106781186 },
{ 0.707106781186, -0.707106781186, 0. },
{ -0.707106781186, 0. , 0.707106781186 },
{ 0. , -0.707106781186, 0.707106781186 } };
const double ptm_template_hcp[PTM_NUM_POINTS_HCP][3] = {
{ 0, 0, 0 },
{ 0.5, -sqrt(3)/2, 0 },
{ -1, 0, 0 },
{ -0.5, sqrt(3)/6, -sqrt(6)/3 },
{ 0.5, sqrt(3)/6, -sqrt(6)/3 },
{ 0, -sqrt(3)/3, -sqrt(6)/3 },
{ -0.5, sqrt(3)/2, 0 },
{ 0.5, sqrt(3)/2, 0 },
{ 1, 0, 0 },
{ -0.5, -sqrt(3)/2, 0 },
{ 0, -sqrt(3)/3, sqrt(6)/3 },
{ 0.5, sqrt(3)/6, sqrt(6)/3 },
{ -0.5, sqrt(3)/6, sqrt(6)/3 },
};
const double ptm_template_bcc[PTM_NUM_POINTS_BCC][3] = { { 0. , 0. , 0. },
{ -0.541451884327, -0.541451884327, -0.541451884327 },
{ 0.541451884327, 0.541451884327, 0.541451884327 },
{ 0.541451884327, -0.541451884327, -0.541451884327 },
{ -0.541451884327, 0.541451884327, 0.541451884327 },
{ -0.541451884327, 0.541451884327, -0.541451884327 },
{ 0.541451884327, -0.541451884327, 0.541451884327 },
{ -0.541451884327, -0.541451884327, 0.541451884327 },
{ 0.541451884327, 0.541451884327, -0.541451884327 },
{ 0. , 0. , -1.082903768655 },
{ 0. , 0. , 1.082903768655 },
{ 0. , -1.082903768655, 0. },
{ 0. , 1.082903768655, 0. },
{ -1.082903768655, 0. , 0. },
{ 1.082903768655, 0. , 0. } };
const double ptm_template_bcc[PTM_NUM_POINTS_BCC][3] = {
{ 0, 0, 0 },
{ 7*sqrt(3)/3-7./2, 7*sqrt(3)/3-7./2, 7*sqrt(3)/3-7./2 },
{ 7./2-7*sqrt(3)/3, 7*sqrt(3)/3-7./2, 7*sqrt(3)/3-7./2 },
{ 7*sqrt(3)/3-7./2, 7*sqrt(3)/3-7./2, 7./2-7*sqrt(3)/3 },
{ 7./2-7*sqrt(3)/3, 7./2-7*sqrt(3)/3, 7*sqrt(3)/3-7./2 },
{ 7*sqrt(3)/3-7./2, 7./2-7*sqrt(3)/3, 7*sqrt(3)/3-7./2 },
{ 7./2-7*sqrt(3)/3, 7*sqrt(3)/3-7./2, 7./2-7*sqrt(3)/3 },
{ 7./2-7*sqrt(3)/3, 7./2-7*sqrt(3)/3, 7./2-7*sqrt(3)/3 },
{ 7*sqrt(3)/3-7./2, 7./2-7*sqrt(3)/3, 7./2-7*sqrt(3)/3 },
{ 14*sqrt(3)/3-7, 0, 0 },
{ 7-14*sqrt(3)/3, 0, 0 },
{ 0, 14*sqrt(3)/3-7, 0 },
{ 0, 7-14*sqrt(3)/3, 0 },
{ 0, 0, 14*sqrt(3)/3-7 },
{ 0, 0, 7-14*sqrt(3)/3 },
};
const double ptm_template_ico[PTM_NUM_POINTS_ICO][3] = { { 0. , 0. , 0. },
{ 0. , 0.525731112119, 0.850650808352 },
{ 0. , -0.525731112119, -0.850650808352 },
{ 0. , 0.525731112119, -0.850650808352 },
{ 0. , -0.525731112119, 0.850650808352 },
{ -0.525731112119, -0.850650808352, 0. },
{ 0.525731112119, 0.850650808352, 0. },
{ 0.525731112119, -0.850650808352, 0. },
{ -0.525731112119, 0.850650808352, 0. },
{ -0.850650808352, 0. , -0.525731112119 },
{ 0.850650808352, 0. , 0.525731112119 },
{ 0.850650808352, 0. , -0.525731112119 },
{ -0.850650808352, 0. , 0.525731112119 } };
const double ptm_template_ico[PTM_NUM_POINTS_ICO][3] = {
{ 0, 0, 0 },
{ 0, 0, 1 },
{ 0, 0, -1 },
{ -sqrt((5-sqrt(5))/10), (5+sqrt(5))/10, -sqrt(5)/5 },
{ sqrt((5-sqrt(5))/10), -(5+sqrt(5))/10, sqrt(5)/5 },
{ 0, -2*sqrt(5)/5, -sqrt(5)/5 },
{ 0, 2*sqrt(5)/5, sqrt(5)/5 },
{ sqrt((5+sqrt(5))/10), -(5-sqrt(5))/10, -sqrt(5)/5 },
{ -sqrt((5+sqrt(5))/10), (5-sqrt(5))/10, sqrt(5)/5 },
{ -sqrt((5+sqrt(5))/10), -(5-sqrt(5))/10, -sqrt(5)/5 },
{ sqrt((5+sqrt(5))/10), (5-sqrt(5))/10, sqrt(5)/5 },
{ sqrt((5-sqrt(5))/10), (5+sqrt(5))/10, -sqrt(5)/5 },
{ -sqrt((5-sqrt(5))/10), -(5+sqrt(5))/10, sqrt(5)/5 },
};
const double ptm_template_sc[PTM_NUM_POINTS_SC][3] = { { 0. , 0. , 0. },
{ 0. , 0. , -1. },
{ 0. , 0. , 1. },
{ 0. , -1. , 0. },
{ 0. , 1. , 0. },
{ -1. , 0. , 0. },
{ 1. , 0. , 0. } };
const double ptm_template_sc[PTM_NUM_POINTS_SC][3] = {
{ 0, 0, 0 },
{ 0, 0, -1 },
{ 0, 0, 1 },
{ 0, -1, 0 },
{ 0, 1, 0 },
{ -1, 0, 0 },
{ 1, 0, 0 },
};
const double ptm_template_dcub[PTM_NUM_POINTS_DCUB][3] = { { 0. , 0. , 0. },
{ -0.391491627053, 0.391491627053, 0.391491627053 },
{ -0.391491627053, -0.391491627053, -0.391491627053 },
{ 0.391491627053, -0.391491627053, 0.391491627053 },
{ 0.391491627053, 0.391491627053, -0.391491627053 },
{ -0.782983254107, 0. , 0.782983254107 },
{ -0.782983254107, 0.782983254107, 0. },
{ 0. , 0.782983254107, 0.782983254107 },
{ -0.782983254107, -0.782983254107, 0. },
{ -0.782983254107, 0. , -0.782983254107 },
{ 0. , -0.782983254107, -0.782983254107 },
{ 0. , -0.782983254107, 0.782983254107 },
{ 0.782983254107, -0.782983254107, 0. },
{ 0.782983254107, 0. , 0.782983254107 },
{ 0. , 0.782983254107, -0.782983254107 },
{ 0.782983254107, 0. , -0.782983254107 },
{ 0.782983254107, 0.782983254107, 0. } };
const double ptm_template_dcub[PTM_NUM_POINTS_DCUB][3] = {
{ 0, 0, 0 },
{ 4/(sqrt(3)+6*sqrt(2)), 4/(sqrt(3)+6*sqrt(2)), 4/(sqrt(3)+6*sqrt(2)) },
{ 4/(sqrt(3)+6*sqrt(2)), -4/(sqrt(3)+6*sqrt(2)), -4/(sqrt(3)+6*sqrt(2)) },
{ -4/(sqrt(3)+6*sqrt(2)), -4/(sqrt(3)+6*sqrt(2)), 4/(sqrt(3)+6*sqrt(2)) },
{ -4/(sqrt(3)+6*sqrt(2)), 4/(sqrt(3)+6*sqrt(2)), -4/(sqrt(3)+6*sqrt(2)) },
{ 8/(sqrt(3)+6*sqrt(2)), 8/(sqrt(3)+6*sqrt(2)), 0 },
{ 0, 8/(sqrt(3)+6*sqrt(2)), 8/(sqrt(3)+6*sqrt(2)) },
{ 8/(sqrt(3)+6*sqrt(2)), 0, 8/(sqrt(3)+6*sqrt(2)) },
{ 0, -8/(sqrt(3)+6*sqrt(2)), -8/(sqrt(3)+6*sqrt(2)) },
{ 8/(sqrt(3)+6*sqrt(2)), -8/(sqrt(3)+6*sqrt(2)), 0 },
{ 8/(sqrt(3)+6*sqrt(2)), 0, -8/(sqrt(3)+6*sqrt(2)) },
{ -8/(sqrt(3)+6*sqrt(2)), -8/(sqrt(3)+6*sqrt(2)), 0 },
{ 0, -8/(sqrt(3)+6*sqrt(2)), 8/(sqrt(3)+6*sqrt(2)) },
{ -8/(sqrt(3)+6*sqrt(2)), 0, 8/(sqrt(3)+6*sqrt(2)) },
{ -8/(sqrt(3)+6*sqrt(2)), 0, -8/(sqrt(3)+6*sqrt(2)) },
{ -8/(sqrt(3)+6*sqrt(2)), 8/(sqrt(3)+6*sqrt(2)), 0 },
{ 0, 8/(sqrt(3)+6*sqrt(2)), -8/(sqrt(3)+6*sqrt(2)) },
};
const double ptm_template_dhex[PTM_NUM_POINTS_DHEX][3] = {
{ 0, 0, 0 },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), -4*sqrt(3)/(3*sqrt(3)+18*sqrt(2)) },
{ 0, -8*sqrt(6)/(3*sqrt(3)+18*sqrt(2)), -4*sqrt(3)/(3*sqrt(3)+18*sqrt(2)) },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), -4*sqrt(3)/(3*sqrt(3)+18*sqrt(2)) },
{ 0, 0, 4*sqrt(3)/(sqrt(3)+6*sqrt(2)) },
{ -8*sqrt(2)/(sqrt(3)+6*sqrt(2)), 0, 0 },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), -16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(sqrt(3)+6*sqrt(2)), 0 },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(sqrt(3)+6*sqrt(2)), 0 },
{ 0, -8*sqrt(6)/(3*sqrt(3)+18*sqrt(2)), -16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), -4*sqrt(6)/(sqrt(3)+6*sqrt(2)), 0 },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), -16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(sqrt(3)+6*sqrt(2)), 0 },
{ 8*sqrt(2)/(sqrt(3)+6*sqrt(2)), 0, 0 },
{ 0, -8*sqrt(6)/(3*sqrt(3)+18*sqrt(2)), 16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ 4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), 16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
{ -4*sqrt(2)/(sqrt(3)+6*sqrt(2)), 4*sqrt(6)/(3*(sqrt(3)+6*sqrt(2))), 16*sqrt(3)/(3*(sqrt(3)+6*sqrt(2))) },
};
const double ptm_template_graphene[PTM_NUM_POINTS_GRAPHENE][3] = {
{ 0, 0, 0 },
{ 0, -3./11+6*sqrt(3)/11, 0 },
{ -3*sqrt(3)/22+9./11, -3*sqrt(3)/11+3./22, 0 },
{ -9./11+3*sqrt(3)/22, -3*sqrt(3)/11+3./22, 0 },
{ -9./11+3*sqrt(3)/22, -9./22+9*sqrt(3)/11, 0 },
{ -3*sqrt(3)/22+9./11, -9./22+9*sqrt(3)/11, 0 },
{ -3*sqrt(3)/11+18./11, 0, 0 },
{ -3*sqrt(3)/22+9./11, -9*sqrt(3)/11+9./22, 0 },
{ -9./11+3*sqrt(3)/22, -9*sqrt(3)/11+9./22, 0 },
{ -18./11+3*sqrt(3)/11, 0, 0 },
};
const double ptm_template_dhex[PTM_NUM_POINTS_DHEX][3] = { { 0. , 0. , 0. },
{ -0.391491627053, -0.391491627053, -0.391491627053 },
{ 0.391491627053, -0.391491627053, 0.391491627053 },
{ -0.391491627053, 0.391491627053, 0.391491627053 },
{ 0.391491627053, 0.391491627053, -0.391491627053 },
{ -0.260994418036, -1.043977672142, -0.260994418036 },
{ -1.043977672142, -0.260994418036, -0.260994418036 },
{ -0.260994418036, -0.260994418036, -1.043977672142 },
{ 0.782983254107, 0. , 0.782983254107 },
{ 0.782983254107, -0.782983254107, 0. },
{ 0. , -0.782983254107, 0.782983254107 },
{ 0. , 0.782983254107, 0.782983254107 },
{ -0.782983254107, 0.782983254107, 0. },
{ -0.782983254107, 0. , 0.782983254107 },
{ 0.782983254107, 0.782983254107, 0. },
{ 0. , 0.782983254107, -0.782983254107 },
{ 0.782983254107, 0. , -0.782983254107 } };
#endif

View File

@ -1,6 +1,15 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#include <cmath>
#include <cfloat>
#include <cstring>
#include <string.h>
#include <cassert>
#include <algorithm>
#include "ptm_convex_hull_incremental.h"

View File

@ -1,3 +1,12 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#ifndef PTM_CONVEX_HULL_INCREMENTAL_H
#define PTM_CONVEX_HULL_INCREMENTAL_H

View File

@ -1,3 +1,12 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#include "ptm_deformation_gradient.h"
namespace ptm {

View File

@ -1,3 +1,12 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#ifndef PTM_DEFORMATION_GRADIENT_H
#define PTM_DEFORMATION_GRADIENT_H
@ -8,139 +17,247 @@ namespace ptm {
void calculate_deformation_gradient(int num_points, const double (*ideal_points)[3], int8_t* mapping, double (*normalized)[3], const double (*penrose)[3], double* F, double* res);
//sc
#define k_sc 0.5
const double penrose_sc[PTM_NUM_POINTS_SC][3] = {
{0, 0, 0},
{0, 0, -k_sc},
{0, 0, k_sc},
{0, -k_sc, 0},
{0, k_sc, 0},
{-k_sc, 0, 0},
{k_sc, 0, 0},
};
const double penrose_sc[PTM_NUM_POINTS_SC][3] = {
{ 0, 0, 0 },
{ 0, 0, -0.5 },
{ 0, 0, 0.5 },
{ 0, -0.5, 0 },
{ 0, 0.5, 0 },
{ -0.5, 0, 0 },
{ 0.5, 0, 0 },
};
//fcc
#define k_fcc 0.17677669529663678216
const double penrose_fcc[PTM_NUM_POINTS_FCC][3] = {
{0, 0, 0},
{0, k_fcc, k_fcc},
{0, -k_fcc, -k_fcc},
{0, k_fcc, -k_fcc},
{0, -k_fcc, k_fcc},
{k_fcc, 0, k_fcc},
{-k_fcc, 0, -k_fcc},
{k_fcc, 0, -k_fcc},
{-k_fcc, 0, k_fcc},
{k_fcc, k_fcc, -0},
{-k_fcc, -k_fcc, 0},
{k_fcc, -k_fcc, 0},
{-k_fcc, k_fcc, -0},
};
{ 0, 0, 0 },
{ sqrt(2)/8, sqrt(2)/8, 0 },
{ 0, sqrt(2)/8, sqrt(2)/8 },
{ sqrt(2)/8, 0, sqrt(2)/8 },
{ -sqrt(2)/8, -sqrt(2)/8, 0 },
{ 0, -sqrt(2)/8, -sqrt(2)/8 },
{ -sqrt(2)/8, 0, -sqrt(2)/8 },
{ -sqrt(2)/8, sqrt(2)/8, 0 },
{ 0, -sqrt(2)/8, sqrt(2)/8 },
{ -sqrt(2)/8, 0, sqrt(2)/8 },
{ sqrt(2)/8, -sqrt(2)/8, 0 },
{ 0, sqrt(2)/8, -sqrt(2)/8 },
{ sqrt(2)/8, 0, -sqrt(2)/8 },
};
//hcp
#define k_hcp 0.17677669529663678216
const double penrose_hcp[PTM_NUM_POINTS_HCP][3] = {
{0, 0, 0},
{k_hcp, 0, k_hcp},
{-k_hcp/3, -4*k_hcp/3, -k_hcp/3},
{k_hcp, k_hcp, 0},
{-k_hcp/3, -k_hcp/3, -4*k_hcp/3},
{0, k_hcp, k_hcp},
{-4*k_hcp/3, -k_hcp/3, -k_hcp/3},
{-k_hcp, k_hcp, -0},
{0, k_hcp, -k_hcp},
{k_hcp, 0, -k_hcp},
{k_hcp, -k_hcp, 0},
{-k_hcp, 0, k_hcp},
{0, -k_hcp, k_hcp},
};
{ 0, 0, 0 },
{ 1./8, -sqrt(3)/8, 0 },
{ -1./4, 0, 0 },
{ -1./8, sqrt(3)/24, -sqrt(6)/12 },
{ 1./8, sqrt(3)/24, -sqrt(6)/12 },
{ 0, -sqrt(3)/12, -sqrt(6)/12 },
{ -1./8, sqrt(3)/8, 0 },
{ 1./8, sqrt(3)/8, 0 },
{ 1./4, 0, 0 },
{ -1./8, -sqrt(3)/8, 0 },
{ 0, -sqrt(3)/12, sqrt(6)/12 },
{ 1./8, sqrt(3)/24, sqrt(6)/12 },
{ -1./8, sqrt(3)/24, sqrt(6)/12 },
};
const double penrose_hcp_alt1[PTM_NUM_POINTS_HCP][3] = {
{ 0, 0, 0 },
{ 1./4, 0, 0 },
{ -1./8, -sqrt(3)/8, 0 },
{ -1./8, -sqrt(3)/24, -sqrt(6)/12 },
{ 0, sqrt(3)/12, -sqrt(6)/12 },
{ 1./8, -sqrt(3)/24, -sqrt(6)/12 },
{ -1./4, 0, 0 },
{ -1./8, sqrt(3)/8, 0 },
{ 1./8, sqrt(3)/8, 0 },
{ 1./8, -sqrt(3)/8, 0 },
{ 1./8, -sqrt(3)/24, sqrt(6)/12 },
{ 0, sqrt(3)/12, sqrt(6)/12 },
{ -1./8, -sqrt(3)/24, sqrt(6)/12 },
};
//ico
#define k_ico 0.13143277802974323576
#define phi 1.61803398874989490253
//((1.0 + sqrt(5)) / 2)
const double penrose_ico[PTM_NUM_POINTS_ICO][3] = {
{0, 0, 0},
{0, k_ico, phi*k_ico},
{0, -k_ico, -phi*k_ico},
{0, k_ico, -phi*k_ico},
{0, -k_ico, phi*k_ico},
{-k_ico, -phi*k_ico, -0},
{k_ico, phi*k_ico, 0},
{k_ico, -phi*k_ico, 0},
{-k_ico, phi*k_ico, -0},
{-phi*k_ico, 0, -k_ico},
{phi*k_ico, 0, k_ico},
{phi*k_ico, 0, -k_ico},
{-phi*k_ico, 0, k_ico},
};
{ 0, 0, 0 },
{ 0, 0, 0.25 },
{ 0, 0, -0.25 },
{ -sqrt(-10*sqrt(5) + 50)/40, sqrt(5)/40 + 1./8, -sqrt(5)/20 },
{ sqrt(-10*sqrt(5) + 50)/40, -1./8 - sqrt(5)/40, sqrt(5)/20 },
{ 0, -sqrt(5)/10, -sqrt(5)/20 },
{ 0, sqrt(5)/10, sqrt(5)/20 },
{ sqrt(10*sqrt(5) + 50)/40, -1./8 + sqrt(5)/40, -sqrt(5)/20 },
{ -sqrt(10*sqrt(5) + 50)/40, -sqrt(5)/40 + 1./8, sqrt(5)/20 },
{ -sqrt(10*sqrt(5) + 50)/40, -1./8 + sqrt(5)/40, -sqrt(5)/20 },
{ sqrt(10*sqrt(5) + 50)/40, -sqrt(5)/40 + 1./8, sqrt(5)/20 },
{ sqrt(-10*sqrt(5) + 50)/40, sqrt(5)/40 + 1./8, -sqrt(5)/20 },
{ -sqrt(-10*sqrt(5) + 50)/40, -1./8 - sqrt(5)/40, sqrt(5)/20 },
};
//bcc
#define k_bcc 0.11543038598460284017
const double penrose_bcc[PTM_NUM_POINTS_BCC][3] = {
{0, 0, 0},
{-k_bcc, -k_bcc, -k_bcc},
{k_bcc, k_bcc, k_bcc},
{k_bcc, -k_bcc, -k_bcc},
{-k_bcc, k_bcc, k_bcc},
{-k_bcc, k_bcc, -k_bcc},
{k_bcc, -k_bcc, k_bcc},
{-k_bcc, -k_bcc, k_bcc},
{k_bcc, k_bcc, -k_bcc},
{0, 0, -2*k_bcc},
{0, 0, 2*k_bcc},
{0, -2*k_bcc, 0},
{0, 2*k_bcc, 0},
{-2*k_bcc, 0, 0},
{2*k_bcc, 0, -0},
};
{ 0, 0, 0 },
{ 3./56 + sqrt(3)/28, 3./56 + sqrt(3)/28, 3./56 + sqrt(3)/28 },
{ -sqrt(3)/28 - 3./56, 3./56 + sqrt(3)/28, 3./56 + sqrt(3)/28 },
{ 3./56 + sqrt(3)/28, 3./56 + sqrt(3)/28, -sqrt(3)/28 - 3./56 },
{ -sqrt(3)/28 - 3./56, -sqrt(3)/28 - 3./56, 3./56 + sqrt(3)/28 },
{ 3./56 + sqrt(3)/28, -sqrt(3)/28 - 3./56, 3./56 + sqrt(3)/28 },
{ -sqrt(3)/28 - 3./56, 3./56 + sqrt(3)/28, -sqrt(3)/28 - 3./56 },
{ -sqrt(3)/28 - 3./56, -sqrt(3)/28 - 3./56, -sqrt(3)/28 - 3./56 },
{ 3./56 + sqrt(3)/28, -sqrt(3)/28 - 3./56, -sqrt(3)/28 - 3./56 },
{ 3./28 + sqrt(3)/14, 0, 0 },
{ -sqrt(3)/14 - 3./28, 0, 0 },
{ 0, 3./28 + sqrt(3)/14, 0 },
{ 0, -sqrt(3)/14 - 3./28, 0 },
{ 0, 0, 3./28 + sqrt(3)/14 },
{ 0, 0, -sqrt(3)/14 - 3./28 },
};
//dcub
#define kdcub 0.07095369570691034689
const double penrose_dcub[PTM_NUM_POINTS_DCUB][3] = {
{ 0, 0, 0 },
{ -kdcub, kdcub, kdcub },
{ -kdcub, -kdcub, -kdcub },
{ kdcub, -kdcub, kdcub },
{ kdcub, kdcub, -kdcub },
{ -2 * kdcub, 0, 2 * kdcub },
{ -2 * kdcub, 2 * kdcub, 0 },
{ 0, 2 * kdcub, 2 * kdcub },
{ -2 * kdcub, -2 * kdcub, 0 },
{ -2 * kdcub, 0, -2 * kdcub },
{ 0, -2 * kdcub, -2 * kdcub },
{ 0, -2 * kdcub, 2 * kdcub },
{ 2 * kdcub, -2 * kdcub, 0 },
{ 2 * kdcub, 0, 2 * kdcub },
{ 0, 2 * kdcub, -2 * kdcub },
{ 2 * kdcub, 0, -2 * kdcub },
{ 2 * kdcub, 2 * kdcub, 0 },
};
{ 0, 0, 0 },
{ 23./(48*(-sqrt(3) + 6*sqrt(2))), 23./(48*(-sqrt(3) + 6*sqrt(2))), 23./(48*(-sqrt(3) + 6*sqrt(2))) },
{ 23./(48*(-sqrt(3) + 6*sqrt(2))), -23./(-48*sqrt(3) + 288*sqrt(2)), -23./(-48*sqrt(3) + 288*sqrt(2)) },
{ -23./(-48*sqrt(3) + 288*sqrt(2)), -23./(-48*sqrt(3) + 288*sqrt(2)), 23./(48*(-sqrt(3) + 6*sqrt(2))) },
{ -23./(-48*sqrt(3) + 288*sqrt(2)), 23./(48*(-sqrt(3) + 6*sqrt(2))), -23./(-48*sqrt(3) + 288*sqrt(2)) },
{ 23./(24*(-sqrt(3) + 6*sqrt(2))), 23./(24*(-sqrt(3) + 6*sqrt(2))), 0 },
{ 0, 23./(24*(-sqrt(3) + 6*sqrt(2))), 23./(24*(-sqrt(3) + 6*sqrt(2))) },
{ 23./(24*(-sqrt(3) + 6*sqrt(2))), 0, 23./(24*(-sqrt(3) + 6*sqrt(2))) },
{ 0, -23./(-24*sqrt(3) + 144*sqrt(2)), -23./(-24*sqrt(3) + 144*sqrt(2)) },
{ 23./(24*(-sqrt(3) + 6*sqrt(2))), -23./(-24*sqrt(3) + 144*sqrt(2)), 0 },
{ 23./(24*(-sqrt(3) + 6*sqrt(2))), 0, -23./(-24*sqrt(3) + 144*sqrt(2)) },
{ -23./(-24*sqrt(3) + 144*sqrt(2)), -23./(-24*sqrt(3) + 144*sqrt(2)), 0 },
{ 0, -23./(-24*sqrt(3) + 144*sqrt(2)), 23./(24*(-sqrt(3) + 6*sqrt(2))) },
{ -23./(-24*sqrt(3) + 144*sqrt(2)), 0, 23./(24*(-sqrt(3) + 6*sqrt(2))) },
{ -23./(-24*sqrt(3) + 144*sqrt(2)), 0, -23./(-24*sqrt(3) + 144*sqrt(2)) },
{ -23./(-24*sqrt(3) + 144*sqrt(2)), 23./(24*(-sqrt(3) + 6*sqrt(2))), 0 },
{ 0, 23./(24*(-sqrt(3) + 6*sqrt(2))), -23./(-24*sqrt(3) + 144*sqrt(2)) },
};
const double penrose_dcub_alt1[PTM_NUM_POINTS_DCUB][3] = {
{ 0, 0, 0 },
{ 23./(48*(-sqrt(3) + 6*sqrt(2))), -23./(-48*sqrt(3) + 288*sqrt(2)), 23./(48*(-sqrt(3) + 6*sqrt(2))) },
{ 23./(48*(-sqrt(3) + 6*sqrt(2))), 23./(48*(-sqrt(3) + 6*sqrt(2))), -23./(-48*sqrt(3) + 288*sqrt(2)) },
{ -23./(-48*sqrt(3) + 288*sqrt(2)), -23./(-48*sqrt(3) + 288*sqrt(2)), -23./(-48*sqrt(3) + 288*sqrt(2)) },
{ -23./(-48*sqrt(3) + 288*sqrt(2)), 23./(48*(-sqrt(3) + 6*sqrt(2))), 23./(48*(-sqrt(3) + 6*sqrt(2))) },
{ 23./(24*(-sqrt(3) + 6*sqrt(2))), 0, 23./(24*(-sqrt(3) + 6*sqrt(2))) },
{ 0, -23./(-24*sqrt(3) + 144*sqrt(2)), 23./(24*(-sqrt(3) + 6*sqrt(2))) },
{ 23./(24*(-sqrt(3) + 6*sqrt(2))), -23./(-24*sqrt(3) + 144*sqrt(2)), 0 },
{ 0, 23./(24*(-sqrt(3) + 6*sqrt(2))), -23./(-24*sqrt(3) + 144*sqrt(2)) },
{ 23./(24*(-sqrt(3) + 6*sqrt(2))), 0, -23./(-24*sqrt(3) + 144*sqrt(2)) },
{ 23./(24*(-sqrt(3) + 6*sqrt(2))), 23./(24*(-sqrt(3) + 6*sqrt(2))), 0 },
{ -23./(-24*sqrt(3) + 144*sqrt(2)), 0, -23./(-24*sqrt(3) + 144*sqrt(2)) },
{ 0, -23./(-24*sqrt(3) + 144*sqrt(2)), -23./(-24*sqrt(3) + 144*sqrt(2)) },
{ -23./(-24*sqrt(3) + 144*sqrt(2)), -23./(-24*sqrt(3) + 144*sqrt(2)), 0 },
{ -23./(-24*sqrt(3) + 144*sqrt(2)), 23./(24*(-sqrt(3) + 6*sqrt(2))), 0 },
{ -23./(-24*sqrt(3) + 144*sqrt(2)), 0, 23./(24*(-sqrt(3) + 6*sqrt(2))) },
{ 0, 23./(24*(-sqrt(3) + 6*sqrt(2))), 23./(24*(-sqrt(3) + 6*sqrt(2))) },
};
#define kdhex 0.04730246380471011397
const double penrose_dhex[PTM_NUM_POINTS_DHEX][3] = {
{ 0, 0, 0 },
{ -kdcub, -kdcub, -kdcub },
{ kdcub, -kdcub, kdcub },
{ -kdcub, kdcub, kdcub },
{ kdcub, kdcub, -kdcub },
{ -kdhex, -4 * kdhex, -kdhex },
{ -4 * kdhex, -kdhex, -kdhex },
{ -kdhex, -kdhex, -4 * kdhex },
{ 2 * kdcub, 0, 2 * kdcub },
{ 2 * kdcub, -2 * kdcub, 0 },
{ 0, -2 * kdcub, 2 * kdcub },
{ 0, 2 * kdcub, 2 * kdcub },
{ -2 * kdcub, 2 * kdcub, 0 },
{ -2 * kdcub, 0, 2 * kdcub },
{ 2 * kdcub, 2 * kdcub, 0 },
{ 0, 2 * kdcub, -2 * kdcub },
{ 2 * kdcub, 0, -2 * kdcub },
};
{ 0, 0, 0 },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), 23*sqrt(6)/(144*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(3)/(-144*sqrt(3) + 864*sqrt(2)) },
{ 0, -23*sqrt(6)/(-72*sqrt(3) + 432*sqrt(2)), -23*sqrt(3)/(-144*sqrt(3) + 864*sqrt(2)) },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(6)/(144*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(3)/(-144*sqrt(3) + 864*sqrt(2)) },
{ 0, 0, 23*sqrt(3)/(48*(-sqrt(3) + 6*sqrt(2))) },
{ -23*sqrt(2)/(-24*sqrt(3) + 144*sqrt(2)), 0, 0 },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), 23*sqrt(6)/(144*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(3)/(-36*sqrt(3) + 216*sqrt(2)) },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), 23*sqrt(6)/(48*(-sqrt(3) + 6*sqrt(2))), 0 },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(6)/(-48*sqrt(3) + 288*sqrt(2)), 0 },
{ 0, -23*sqrt(6)/(-72*sqrt(3) + 432*sqrt(2)), -23*sqrt(3)/(-36*sqrt(3) + 216*sqrt(2)) },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), -23*sqrt(6)/(-48*sqrt(3) + 288*sqrt(2)), 0 },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(6)/(144*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(3)/(-36*sqrt(3) + 216*sqrt(2)) },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(6)/(48*(-sqrt(3) + 6*sqrt(2))), 0 },
{ 23*sqrt(2)/(24*(-sqrt(3) + 6*sqrt(2))), 0, 0 },
{ 0, -23*sqrt(6)/(-72*sqrt(3) + 432*sqrt(2)), 23*sqrt(3)/(36*(-sqrt(3) + 6*sqrt(2))) },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(6)/(144*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(3)/(36*(-sqrt(3) + 6*sqrt(2))) },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), 23*sqrt(6)/(144*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(3)/(36*(-sqrt(3) + 6*sqrt(2))) },
};
const double penrose_dhex_alt1[PTM_NUM_POINTS_DHEX][3] = {
{ 0, 0, 0 },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), -23*sqrt(6)/(-144*sqrt(3) + 864*sqrt(2)), -23*sqrt(3)/(-144*sqrt(3) + 864*sqrt(2)) },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(6)/(-144*sqrt(3) + 864*sqrt(2)), -23*sqrt(3)/(-144*sqrt(3) + 864*sqrt(2)) },
{ 0, 23*sqrt(6)/(72*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(3)/(-144*sqrt(3) + 864*sqrt(2)) },
{ 0, 0, 23*sqrt(3)/(48*(-sqrt(3) + 6*sqrt(2))) },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), -23*sqrt(6)/(-48*sqrt(3) + 288*sqrt(2)), 0 },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), -23*sqrt(6)/(-144*sqrt(3) + 864*sqrt(2)), -23*sqrt(3)/(-36*sqrt(3) + 216*sqrt(2)) },
{ -23*sqrt(2)/(-24*sqrt(3) + 144*sqrt(2)), 0, 0 },
{ 23*sqrt(2)/(24*(-sqrt(3) + 6*sqrt(2))), 0, 0 },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(6)/(-144*sqrt(3) + 864*sqrt(2)), -23*sqrt(3)/(-36*sqrt(3) + 216*sqrt(2)) },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(6)/(-48*sqrt(3) + 288*sqrt(2)), 0 },
{ 0, 23*sqrt(6)/(72*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(3)/(-36*sqrt(3) + 216*sqrt(2)) },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), 23*sqrt(6)/(48*(-sqrt(3) + 6*sqrt(2))), 0 },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(6)/(48*(-sqrt(3) + 6*sqrt(2))), 0 },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(6)/(-144*sqrt(3) + 864*sqrt(2)), 23*sqrt(3)/(36*(-sqrt(3) + 6*sqrt(2))) },
{ 0, 23*sqrt(6)/(72*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(3)/(36*(-sqrt(3) + 6*sqrt(2))) },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), -23*sqrt(6)/(-144*sqrt(3) + 864*sqrt(2)), 23*sqrt(3)/(36*(-sqrt(3) + 6*sqrt(2))) },
};
const double penrose_dhex_alt2[PTM_NUM_POINTS_DHEX][3] = {
{ 0, 0, 0 },
{ 0, -23*sqrt(6)/(-72*sqrt(3) + 432*sqrt(2)), 23*sqrt(3)/(144*(-sqrt(3) + 6*sqrt(2))) },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), 23*sqrt(6)/(144*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(3)/(144*(-sqrt(3) + 6*sqrt(2))) },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(6)/(144*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(3)/(144*(-sqrt(3) + 6*sqrt(2))) },
{ 0, 0, -23*sqrt(3)/(-48*sqrt(3) + 288*sqrt(2)) },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), -23*sqrt(6)/(-48*sqrt(3) + 288*sqrt(2)), 0 },
{ 0, -23*sqrt(6)/(-72*sqrt(3) + 432*sqrt(2)), 23*sqrt(3)/(36*(-sqrt(3) + 6*sqrt(2))) },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(6)/(-48*sqrt(3) + 288*sqrt(2)), 0 },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), 23*sqrt(6)/(48*(-sqrt(3) + 6*sqrt(2))), 0 },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), 23*sqrt(6)/(144*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(3)/(36*(-sqrt(3) + 6*sqrt(2))) },
{ -23*sqrt(2)/(-24*sqrt(3) + 144*sqrt(2)), 0, 0 },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(6)/(144*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(3)/(36*(-sqrt(3) + 6*sqrt(2))) },
{ 23*sqrt(2)/(24*(-sqrt(3) + 6*sqrt(2))), 0, 0 },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(6)/(48*(-sqrt(3) + 6*sqrt(2))), 0 },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), 23*sqrt(6)/(144*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(3)/(-36*sqrt(3) + 216*sqrt(2)) },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(6)/(144*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(3)/(-36*sqrt(3) + 216*sqrt(2)) },
{ 0, -23*sqrt(6)/(-72*sqrt(3) + 432*sqrt(2)), -23*sqrt(3)/(-36*sqrt(3) + 216*sqrt(2)) },
};
const double penrose_dhex_alt3[PTM_NUM_POINTS_DHEX][3] = {
{ 0, 0, 0 },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(6)/(-144*sqrt(3) + 864*sqrt(2)), 23*sqrt(3)/(144*(-sqrt(3) + 6*sqrt(2))) },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), -23*sqrt(6)/(-144*sqrt(3) + 864*sqrt(2)), 23*sqrt(3)/(144*(-sqrt(3) + 6*sqrt(2))) },
{ 0, 23*sqrt(6)/(72*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(3)/(144*(-sqrt(3) + 6*sqrt(2))) },
{ 0, 0, -23*sqrt(3)/(-48*sqrt(3) + 288*sqrt(2)) },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(6)/(-48*sqrt(3) + 288*sqrt(2)), 0 },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(6)/(-144*sqrt(3) + 864*sqrt(2)), 23*sqrt(3)/(36*(-sqrt(3) + 6*sqrt(2))) },
{ 23*sqrt(2)/(24*(-sqrt(3) + 6*sqrt(2))), 0, 0 },
{ -23*sqrt(2)/(-24*sqrt(3) + 144*sqrt(2)), 0, 0 },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), -23*sqrt(6)/(-144*sqrt(3) + 864*sqrt(2)), 23*sqrt(3)/(36*(-sqrt(3) + 6*sqrt(2))) },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), -23*sqrt(6)/(-48*sqrt(3) + 288*sqrt(2)), 0 },
{ 0, 23*sqrt(6)/(72*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(3)/(36*(-sqrt(3) + 6*sqrt(2))) },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), 23*sqrt(6)/(48*(-sqrt(3) + 6*sqrt(2))), 0 },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), 23*sqrt(6)/(48*(-sqrt(3) + 6*sqrt(2))), 0 },
{ -23*sqrt(2)/(-48*sqrt(3) + 288*sqrt(2)), -23*sqrt(6)/(-144*sqrt(3) + 864*sqrt(2)), -23*sqrt(3)/(-36*sqrt(3) + 216*sqrt(2)) },
{ 0, 23*sqrt(6)/(72*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(3)/(-36*sqrt(3) + 216*sqrt(2)) },
{ 23*sqrt(2)/(48*(-sqrt(3) + 6*sqrt(2))), -23*sqrt(6)/(-144*sqrt(3) + 864*sqrt(2)), -23*sqrt(3)/(-36*sqrt(3) + 216*sqrt(2)) },
};
const double penrose_graphene[PTM_NUM_POINTS_GRAPHENE][3] = {
{ 0, 0, 0 },
{ 0, 2./63 + 4*sqrt(3)/63, 0 },
{ sqrt(3)/63 + 2./21, -2*sqrt(3)/63 - 1./63, 0 },
{ -2./21 - sqrt(3)/63, -2*sqrt(3)/63 - 1./63, 0 },
{ -2./21 - sqrt(3)/63, 1./21 + 2*sqrt(3)/21, 0 },
{ sqrt(3)/63 + 2./21, 1./21 + 2*sqrt(3)/21, 0 },
{ 2*sqrt(3)/63 + 4./21, 0, 0 },
{ sqrt(3)/63 + 2./21, -2*sqrt(3)/21 - 1./21, 0 },
{ -2./21 - sqrt(3)/63, -2*sqrt(3)/21 - 1./21, 0 },
{ -4./21 - 2*sqrt(3)/63, 0, 0 },
};
const double penrose_graphene_alt1[PTM_NUM_POINTS_GRAPHENE][3] = {
{ 0, 0, 0 },
{ -2./21 - sqrt(3)/63, 1./63 + 2*sqrt(3)/63, 0 },
{ sqrt(3)/63 + 2./21, 1./63 + 2*sqrt(3)/63, 0 },
{ 0, -4*sqrt(3)/63 - 2./63, 0 },
{ -4./21 - 2*sqrt(3)/63, 0, 0 },
{ -2./21 - sqrt(3)/63, 1./21 + 2*sqrt(3)/21, 0 },
{ sqrt(3)/63 + 2./21, 1./21 + 2*sqrt(3)/21, 0 },
{ 2*sqrt(3)/63 + 4./21, 0, 0 },
{ sqrt(3)/63 + 2./21, -2*sqrt(3)/21 - 1./21, 0 },
{ -2./21 - sqrt(3)/63, -2*sqrt(3)/21 - 1./21, 0 },
};
}
#endif

View File

@ -1,3 +1,12 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#ifndef PTM_FUNCTIONS_H
#define PTM_FUNCTIONS_H
@ -15,8 +24,11 @@ extern "C" {
#endif
int ptm_index( ptm_local_handle_t local_handle, int32_t flags, int num_points, double (*atomic_positions)[3], int32_t* atomic_numbers, bool topological_ordering, //inputs
int32_t* p_type, int32_t* p_alloy_type, double* p_scale, double* p_rmsd, double* q, double* F, double* F_res, double* U, double* P, int8_t* mapping, double* p_interatomic_distance, double* p_lattice_constant); //outputs
int ptm_index( ptm_local_handle_t local_handle,
size_t atom_index, int (get_neighbours)(void* vdata, size_t central_index, size_t atom_index, int num, size_t* nbr_indices, int32_t* numbers, double (*nbr_pos)[3]), void* nbrlist,
int32_t flags, bool output_conventional_orientation, //inputs
int32_t* p_type, int32_t* p_alloy_type, double* p_scale, double* p_rmsd, double* q, double* F, double* F_res, double* U, double* P, double* p_interatomic_distance, double* p_lattice_constant,
size_t* output_indices); //outputs
#ifdef __cplusplus

View File

@ -1,3 +1,12 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#ifndef PTM_FUNDAMENTAL_MAPPINGS_H
#define PTM_FUNDAMENTAL_MAPPINGS_H
@ -11,172 +20,309 @@ namespace ptm {
#define NUM_DCUB_MAPPINGS 12
#define NUM_DHEX_MAPPINGS 3
#define NUM_GRAPHENE_MAPPINGS 6
#define NUM_CONVENTIONAL_GRAPHENE_MAPPINGS 12
#define NUM_CONVENTIONAL_HEX_MAPPINGS 12
#define NUM_CONVENTIONAL_DCUB_MAPPINGS 24
#define NUM_CONVENTIONAL_DHEX_MAPPINGS 12
//-----------------------------------------------------------------------------
// sc
//-----------------------------------------------------------------------------
const int8_t mapping_sc[NUM_CUBIC_MAPPINGS][PTM_MAX_POINTS] = {
{0, 1, 2, 3, 4, 5, 6},
{0, 2, 1, 4, 3, 5, 6},
{0, 2, 1, 3, 4, 6, 5},
{0, 1, 2, 4, 3, 6, 5},
{0, 3, 4, 5, 6, 1, 2},
{0, 5, 6, 2, 1, 4, 3},
{0, 6, 5, 1, 2, 4, 3},
{0, 4, 3, 5, 6, 2, 1},
{0, 5, 6, 1, 2, 3, 4},
{0, 4, 3, 6, 5, 1, 2},
{0, 3, 4, 6, 5, 2, 1},
{0, 6, 5, 2, 1, 3, 4},
{0, 3, 4, 2, 1, 5, 6},
{0, 6, 5, 3, 4, 1, 2},
{0, 1, 2, 5, 6, 4, 3},
{0, 4, 3, 1, 2, 5, 6},
{0, 5, 6, 3, 4, 2, 1},
{0, 1, 2, 6, 5, 3, 4},
{0, 2, 1, 5, 6, 3, 4},
{0, 5, 6, 4, 3, 1, 2},
{0, 3, 4, 1, 2, 6, 5},
{0, 2, 1, 6, 5, 4, 3},
{0, 6, 5, 4, 3, 2, 1},
{0, 4, 3, 2, 1, 6, 5} };
{ 0, 1, 2, 3, 4, 5, 6},
{ 0, 3, 4, 2, 1, 5, 6},
{ 0, 6, 5, 3, 4, 1, 2},
{ 0, 1, 2, 5, 6, 4, 3},
{ 0, 1, 2, 6, 5, 3, 4},
{ 0, 5, 6, 3, 4, 2, 1},
{ 0, 4, 3, 1, 2, 5, 6},
{ 0, 3, 4, 5, 6, 1, 2},
{ 0, 6, 5, 2, 1, 3, 4},
{ 0, 5, 6, 2, 1, 4, 3},
{ 0, 3, 4, 6, 5, 2, 1},
{ 0, 6, 5, 1, 2, 4, 3},
{ 0, 4, 3, 6, 5, 1, 2},
{ 0, 4, 3, 5, 6, 2, 1},
{ 0, 5, 6, 1, 2, 3, 4},
{ 0, 2, 1, 4, 3, 5, 6},
{ 0, 2, 1, 5, 6, 3, 4},
{ 0, 5, 6, 4, 3, 1, 2},
{ 0, 6, 5, 4, 3, 2, 1},
{ 0, 2, 1, 6, 5, 4, 3},
{ 0, 2, 1, 3, 4, 6, 5},
{ 0, 3, 4, 1, 2, 6, 5},
{ 0, 4, 3, 2, 1, 6, 5},
{ 0, 1, 2, 4, 3, 6, 5},
};
//-----------------------------------------------------------------------------
// fcc
//-----------------------------------------------------------------------------
const int8_t mapping_fcc[NUM_CUBIC_MAPPINGS][PTM_MAX_POINTS] = {
{0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12},
{0, 2, 1, 4, 3, 7, 8, 5, 6, 11, 12, 9, 10},
{0, 3, 4, 1, 2, 6, 5, 8, 7, 12, 11, 10, 9},
{0, 4, 3, 2, 1, 8, 7, 6, 5, 10, 9, 12, 11},
{0, 9, 10, 11, 12, 1, 2, 4, 3, 5, 6, 8, 7},
{0, 7, 8, 6, 5, 11, 12, 10, 9, 2, 1, 4, 3},
{0, 8, 7, 5, 6, 10, 9, 11, 12, 4, 3, 2, 1},
{0, 11, 12, 9, 10, 2, 1, 3, 4, 7, 8, 6, 5},
{0, 5, 6, 8, 7, 9, 10, 12, 11, 1, 2, 3, 4},
{0, 10, 9, 12, 11, 4, 3, 1, 2, 8, 7, 5, 6},
{0, 12, 11, 10, 9, 3, 4, 2, 1, 6, 5, 7, 8},
{0, 6, 5, 7, 8, 12, 11, 9, 10, 3, 4, 1, 2},
{0, 3, 4, 2, 1, 9, 10, 11, 12, 7, 8, 5, 6},
{0, 12, 11, 9, 10, 8, 7, 5, 6, 1, 2, 4, 3},
{0, 5, 6, 7, 8, 4, 3, 2, 1, 11, 12, 10, 9},
{0, 4, 3, 1, 2, 11, 12, 9, 10, 5, 6, 7, 8},
{0, 9, 10, 12, 11, 7, 8, 6, 5, 3, 4, 2, 1},
{0, 8, 7, 6, 5, 1, 2, 3, 4, 12, 11, 9, 10},
{0, 7, 8, 5, 6, 3, 4, 1, 2, 9, 10, 12, 11},
{0, 11, 12, 10, 9, 5, 6, 8, 7, 4, 3, 1, 2},
{0, 1, 2, 4, 3, 12, 11, 10, 9, 8, 7, 6, 5},
{0, 6, 5, 8, 7, 2, 1, 4, 3, 10, 9, 11, 12},
{0, 10, 9, 11, 12, 6, 5, 7, 8, 2, 1, 3, 4},
{0, 2, 1, 3, 4, 10, 9, 12, 11, 6, 5, 8, 7} };
{ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12},
{ 0, 12, 11, 1, 9, 8, 4, 6, 2, 7, 3, 5, 10},
{ 0, 2, 7, 9, 5, 10, 12, 11, 4, 6, 8, 1, 3},
{ 0, 10, 3, 8, 7, 6, 11, 1, 9, 2, 4, 12, 5},
{ 0, 7, 9, 2, 10, 12, 5, 4, 3, 8, 1, 6, 11},
{ 0, 11, 1, 12, 8, 4, 9, 2, 10, 3, 5, 7, 6},
{ 0, 3, 8, 10, 6, 11, 7, 9, 5, 4, 12, 2, 1},
{ 0, 3, 1, 2, 6, 4, 5, 12, 7, 11, 9, 10, 8},
{ 0, 11, 6, 7, 8, 3, 10, 5, 9, 4, 2, 12, 1},
{ 0, 5, 12, 10, 2, 9, 7, 11, 3, 1, 8, 6, 4},
{ 0, 6, 7, 11, 3, 10, 8, 9, 1, 2, 12, 4, 5},
{ 0, 8, 9, 4, 11, 12, 1, 2, 6, 7, 5, 3, 10},
{ 0, 9, 4, 8, 12, 1, 11, 6, 10, 5, 3, 7, 2},
{ 0, 12, 10, 5, 9, 7, 2, 3, 4, 8, 6, 1, 11},
{ 0, 2, 3, 1, 5, 6, 4, 8, 12, 10, 11, 9, 7},
{ 0, 10, 5, 12, 7, 2, 9, 4, 11, 6, 1, 8, 3},
{ 0, 1, 12, 11, 4, 9, 8, 10, 6, 5, 7, 3, 2},
{ 0, 8, 10, 3, 11, 7, 6, 5, 1, 12, 2, 4, 9},
{ 0, 5, 4, 6, 2, 1, 3, 8, 7, 9, 11, 10, 12},
{ 0, 4, 6, 5, 1, 3, 2, 7, 12, 11, 10, 9, 8},
{ 0, 7, 11, 6, 10, 8, 3, 1, 5, 12, 4, 2, 9},
{ 0, 9, 2, 7, 12, 5, 10, 3, 11, 1, 6, 8, 4},
{ 0, 6, 5, 4, 3, 2, 1, 12, 8, 10, 9, 11, 7},
{ 0, 4, 8, 9, 1, 11, 12, 10, 2, 3, 7, 5, 6},
};
//-----------------------------------------------------------------------------
// bcc
//-----------------------------------------------------------------------------
const int8_t mapping_bcc[NUM_CUBIC_MAPPINGS][PTM_MAX_POINTS] = {
{0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14},
{0, 4, 3, 2, 1, 7, 8, 5, 6, 10, 9, 12, 11, 13, 14},
{0, 6, 5, 7, 8, 2, 1, 3, 4, 10, 9, 11, 12, 14, 13},
{0, 8, 7, 5, 6, 3, 4, 2, 1, 9, 10, 12, 11, 14, 13},
{0, 1, 2, 7, 8, 3, 4, 5, 6, 11, 12, 13, 14, 9, 10},
{0, 4, 3, 7, 8, 5, 6, 2, 1, 13, 14, 10, 9, 12, 11},
{0, 8, 7, 3, 4, 2, 1, 5, 6, 14, 13, 9, 10, 12, 11},
{0, 4, 3, 5, 6, 2, 1, 7, 8, 12, 11, 13, 14, 10, 9},
{0, 1, 2, 5, 6, 7, 8, 3, 4, 13, 14, 9, 10, 11, 12},
{0, 8, 7, 2, 1, 5, 6, 3, 4, 12, 11, 14, 13, 9, 10},
{0, 6, 5, 3, 4, 7, 8, 2, 1, 11, 12, 14, 13, 10, 9},
{0, 6, 5, 2, 1, 3, 4, 7, 8, 14, 13, 10, 9, 11, 12},
{0, 7, 8, 6, 5, 1, 2, 4, 3, 11, 12, 10, 9, 13, 14},
{0, 3, 4, 6, 5, 8, 7, 1, 2, 14, 13, 11, 12, 9, 10},
{0, 5, 6, 1, 2, 8, 7, 4, 3, 9, 10, 13, 14, 12, 11},
{0, 5, 6, 8, 7, 4, 3, 1, 2, 12, 11, 9, 10, 13, 14},
{0, 7, 8, 1, 2, 4, 3, 6, 5, 13, 14, 11, 12, 10, 9},
{0, 3, 4, 8, 7, 1, 2, 6, 5, 9, 10, 14, 13, 11, 12},
{0, 7, 8, 4, 3, 6, 5, 1, 2, 10, 9, 13, 14, 11, 12},
{0, 5, 6, 4, 3, 1, 2, 8, 7, 13, 14, 12, 11, 9, 10},
{0, 3, 4, 1, 2, 6, 5, 8, 7, 11, 12, 9, 10, 14, 13},
{0, 2, 1, 6, 5, 4, 3, 8, 7, 10, 9, 14, 13, 12, 11},
{0, 2, 1, 8, 7, 6, 5, 4, 3, 14, 13, 12, 11, 10, 9},
{0, 2, 1, 4, 3, 8, 7, 6, 5, 12, 11, 10, 9, 14, 13} };
{ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14},
{ 0, 3, 6, 8, 2, 1, 7, 4, 5, 9, 10, 14, 13, 11, 12},
{ 0, 2, 6, 1, 7, 4, 3, 8, 5, 13, 14, 11, 12, 10, 9},
{ 0, 5, 1, 8, 2, 4, 3, 6, 7, 12, 11, 9, 10, 13, 14},
{ 0, 2, 4, 6, 5, 1, 7, 8, 3, 11, 12, 10, 9, 13, 14},
{ 0, 3, 1, 6, 5, 8, 2, 4, 7, 14, 13, 11, 12, 9, 10},
{ 0, 5, 4, 1, 7, 8, 2, 6, 3, 9, 10, 13, 14, 12, 11},
{ 0, 1, 3, 5, 6, 2, 8, 7, 4, 13, 14, 9, 10, 11, 12},
{ 0, 6, 7, 3, 4, 2, 8, 5, 1, 11, 12, 14, 13, 10, 9},
{ 0, 8, 3, 7, 1, 5, 6, 2, 4, 12, 11, 14, 13, 9, 10},
{ 0, 6, 2, 7, 1, 3, 4, 5, 8, 14, 13, 10, 9, 11, 12},
{ 0, 4, 2, 5, 6, 7, 1, 3, 8, 12, 11, 13, 14, 10, 9},
{ 0, 4, 7, 2, 8, 5, 6, 3, 1, 13, 14, 10, 9, 12, 11},
{ 0, 8, 5, 3, 4, 7, 1, 2, 6, 14, 13, 9, 10, 12, 11},
{ 0, 1, 5, 2, 8, 3, 4, 7, 6, 11, 12, 13, 14, 9, 10},
{ 0, 8, 7, 5, 6, 3, 4, 2, 1, 9, 10, 12, 11, 14, 13},
{ 0, 3, 8, 1, 7, 6, 5, 4, 2, 11, 12, 9, 10, 14, 13},
{ 0, 5, 8, 4, 3, 1, 7, 6, 2, 13, 14, 12, 11, 9, 10},
{ 0, 7, 4, 8, 2, 6, 5, 1, 3, 14, 13, 12, 11, 10, 9},
{ 0, 7, 6, 4, 3, 8, 2, 1, 5, 12, 11, 10, 9, 14, 13},
{ 0, 6, 3, 2, 8, 7, 1, 5, 4, 10, 9, 11, 12, 14, 13},
{ 0, 2, 1, 4, 3, 6, 5, 8, 7, 10, 9, 13, 14, 11, 12},
{ 0, 7, 8, 6, 5, 4, 3, 1, 2, 10, 9, 14, 13, 12, 11},
{ 0, 4, 5, 7, 1, 2, 8, 3, 6, 10, 9, 12, 11, 13, 14},
};
//-----------------------------------------------------------------------------
// ico
//-----------------------------------------------------------------------------
const int8_t mapping_ico[NUM_ICO_MAPPINGS][PTM_MAX_POINTS] = {
{0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12},
{0, 10, 9, 8, 7, 5, 6, 2, 1, 12, 11, 3, 4},
{0, 1, 2, 9, 10, 7, 8, 11, 12, 5, 6, 3, 4},
{0, 4, 3, 8, 7, 2, 1, 11, 12, 9, 10, 6, 5},
{0, 6, 5, 9, 10, 4, 3, 7, 8, 12, 11, 2, 1},
{0, 12, 11, 3, 4, 7, 8, 10, 9, 2, 1, 6, 5},
{0, 4, 3, 6, 5, 9, 10, 2, 1, 8, 7, 11, 12},
{0, 8, 7, 2, 1, 4, 3, 10, 9, 5, 6, 11, 12},
{0, 10, 9, 3, 4, 12, 11, 5, 6, 8, 7, 2, 1},
{0, 12, 11, 6, 5, 2, 1, 7, 8, 3, 4, 10, 9},
{0, 1, 2, 11, 12, 9, 10, 5, 6, 3, 4, 7, 8},
{0, 8, 7, 11, 12, 5, 6, 4, 3, 2, 1, 10, 9},
{0, 6, 5, 2, 1, 12, 11, 4, 3, 9, 10, 7, 8},
{0, 3, 4, 5, 6, 1, 2, 10, 9, 12, 11, 7, 8},
{0, 3, 4, 7, 8, 12, 11, 1, 2, 5, 6, 10, 9},
{0, 6, 5, 7, 8, 9, 10, 12, 11, 2, 1, 4, 3},
{0, 9, 10, 11, 12, 4, 3, 1, 2, 7, 8, 6, 5},
{0, 11, 12, 9, 10, 1, 2, 4, 3, 8, 7, 5, 6},
{0, 8, 7, 5, 6, 10, 9, 11, 12, 4, 3, 2, 1},
{0, 10, 9, 2, 1, 8, 7, 12, 11, 3, 4, 5, 6},
{0, 12, 11, 2, 1, 10, 9, 6, 5, 7, 8, 3, 4},
{0, 9, 10, 6, 5, 7, 8, 4, 3, 11, 12, 1, 2},
{0, 8, 7, 10, 9, 2, 1, 5, 6, 11, 12, 4, 3},
{0, 6, 5, 12, 11, 7, 8, 2, 1, 4, 3, 9, 10},
{0, 11, 12, 8, 7, 4, 3, 5, 6, 1, 2, 9, 10},
{0, 4, 3, 11, 12, 8, 7, 9, 10, 6, 5, 2, 1},
{0, 4, 3, 9, 10, 11, 12, 6, 5, 2, 1, 8, 7},
{0, 12, 11, 10, 9, 3, 4, 2, 1, 6, 5, 7, 8},
{0, 5, 6, 8, 7, 11, 12, 10, 9, 3, 4, 1, 2},
{0, 7, 8, 6, 5, 12, 11, 9, 10, 1, 2, 3, 4},
{0, 10, 9, 12, 11, 2, 1, 3, 4, 5, 6, 8, 7},
{0, 7, 8, 1, 2, 9, 10, 3, 4, 12, 11, 6, 5},
{0, 5, 6, 1, 2, 3, 4, 11, 12, 8, 7, 10, 9},
{0, 7, 8, 12, 11, 3, 4, 6, 5, 9, 10, 1, 2},
{0, 1, 2, 5, 6, 11, 12, 3, 4, 7, 8, 9, 10},
{0, 11, 12, 1, 2, 5, 6, 9, 10, 4, 3, 8, 7},
{0, 5, 6, 3, 4, 10, 9, 1, 2, 11, 12, 8, 7},
{0, 5, 6, 10, 9, 8, 7, 3, 4, 1, 2, 11, 12},
{0, 3, 4, 12, 11, 10, 9, 7, 8, 1, 2, 5, 6},
{0, 9, 10, 7, 8, 1, 2, 6, 5, 4, 3, 11, 12},
{0, 9, 10, 1, 2, 11, 12, 7, 8, 6, 5, 4, 3},
{0, 7, 8, 3, 4, 1, 2, 12, 11, 6, 5, 9, 10},
{0, 11, 12, 5, 6, 8, 7, 1, 2, 9, 10, 4, 3},
{0, 1, 2, 7, 8, 3, 4, 9, 10, 11, 12, 5, 6},
{0, 3, 4, 10, 9, 5, 6, 12, 11, 7, 8, 1, 2},
{0, 2, 1, 4, 3, 8, 7, 6, 5, 12, 11, 10, 9},
{0, 2, 1, 12, 11, 6, 5, 10, 9, 8, 7, 4, 3},
{0, 9, 10, 4, 3, 6, 5, 11, 12, 1, 2, 7, 8},
{0, 11, 12, 4, 3, 9, 10, 8, 7, 5, 6, 1, 2},
{0, 2, 1, 10, 9, 12, 11, 8, 7, 4, 3, 6, 5},
{0, 5, 6, 11, 12, 1, 2, 8, 7, 10, 9, 3, 4},
{0, 10, 9, 5, 6, 3, 4, 8, 7, 2, 1, 12, 11},
{0, 12, 11, 7, 8, 6, 5, 3, 4, 10, 9, 2, 1},
{0, 7, 8, 9, 10, 6, 5, 1, 2, 3, 4, 12, 11},
{0, 2, 1, 8, 7, 10, 9, 4, 3, 6, 5, 12, 11},
{0, 8, 7, 4, 3, 11, 12, 2, 1, 10, 9, 5, 6},
{0, 6, 5, 4, 3, 2, 1, 9, 10, 7, 8, 12, 11},
{0, 2, 1, 6, 5, 4, 3, 12, 11, 10, 9, 8, 7},
{0, 3, 4, 1, 2, 7, 8, 5, 6, 10, 9, 12, 11},
{0, 4, 3, 2, 1, 6, 5, 8, 7, 11, 12, 9, 10} };
{ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12},
{ 0, 6, 5, 2, 1, 12, 11, 4, 3, 9, 10, 7, 8},
{ 0, 6, 5, 9, 10, 4, 3, 7, 8, 12, 11, 2, 1},
{ 0, 8, 7, 2, 1, 4, 3, 10, 9, 5, 6, 11, 12},
{ 0, 10, 9, 3, 4, 12, 11, 5, 6, 8, 7, 2, 1},
{ 0, 8, 7, 11, 12, 5, 6, 4, 3, 2, 1, 10, 9},
{ 0, 1, 2, 11, 12, 9, 10, 5, 6, 3, 4, 7, 8},
{ 0, 1, 2, 9, 10, 7, 8, 11, 12, 5, 6, 3, 4},
{ 0, 10, 9, 8, 7, 5, 6, 2, 1, 12, 11, 3, 4},
{ 0, 12, 11, 3, 4, 7, 8, 10, 9, 2, 1, 6, 5},
{ 0, 4, 3, 6, 5, 9, 10, 2, 1, 8, 7, 11, 12},
{ 0, 12, 11, 6, 5, 2, 1, 7, 8, 3, 4, 10, 9},
{ 0, 4, 3, 8, 7, 2, 1, 11, 12, 9, 10, 6, 5},
{ 0, 3, 4, 5, 6, 1, 2, 10, 9, 12, 11, 7, 8},
{ 0, 11, 12, 9, 10, 1, 2, 4, 3, 8, 7, 5, 6},
{ 0, 3, 4, 7, 8, 12, 11, 1, 2, 5, 6, 10, 9},
{ 0, 8, 7, 5, 6, 10, 9, 11, 12, 4, 3, 2, 1},
{ 0, 10, 9, 2, 1, 8, 7, 12, 11, 3, 4, 5, 6},
{ 0, 11, 12, 8, 7, 4, 3, 5, 6, 1, 2, 9, 10},
{ 0, 6, 5, 7, 8, 9, 10, 12, 11, 2, 1, 4, 3},
{ 0, 6, 5, 12, 11, 7, 8, 2, 1, 4, 3, 9, 10},
{ 0, 9, 10, 11, 12, 4, 3, 1, 2, 7, 8, 6, 5},
{ 0, 12, 11, 2, 1, 10, 9, 6, 5, 7, 8, 3, 4},
{ 0, 4, 3, 11, 12, 8, 7, 9, 10, 6, 5, 2, 1},
{ 0, 7, 8, 6, 5, 12, 11, 9, 10, 1, 2, 3, 4},
{ 0, 8, 7, 10, 9, 2, 1, 5, 6, 11, 12, 4, 3},
{ 0, 10, 9, 12, 11, 2, 1, 3, 4, 5, 6, 8, 7},
{ 0, 9, 10, 6, 5, 7, 8, 4, 3, 11, 12, 1, 2},
{ 0, 4, 3, 9, 10, 11, 12, 6, 5, 2, 1, 8, 7},
{ 0, 12, 11, 10, 9, 3, 4, 2, 1, 6, 5, 7, 8},
{ 0, 7, 8, 1, 2, 9, 10, 3, 4, 12, 11, 6, 5},
{ 0, 5, 6, 8, 7, 11, 12, 10, 9, 3, 4, 1, 2},
{ 0, 5, 6, 1, 2, 3, 4, 11, 12, 8, 7, 10, 9},
{ 0, 11, 12, 5, 6, 8, 7, 1, 2, 9, 10, 4, 3},
{ 0, 3, 4, 12, 11, 10, 9, 7, 8, 1, 2, 5, 6},
{ 0, 9, 10, 7, 8, 1, 2, 6, 5, 4, 3, 11, 12},
{ 0, 7, 8, 3, 4, 1, 2, 12, 11, 6, 5, 9, 10},
{ 0, 3, 4, 10, 9, 5, 6, 12, 11, 7, 8, 1, 2},
{ 0, 1, 2, 7, 8, 3, 4, 9, 10, 11, 12, 5, 6},
{ 0, 1, 2, 5, 6, 11, 12, 3, 4, 7, 8, 9, 10},
{ 0, 11, 12, 1, 2, 5, 6, 9, 10, 4, 3, 8, 7},
{ 0, 5, 6, 3, 4, 10, 9, 1, 2, 11, 12, 8, 7},
{ 0, 5, 6, 10, 9, 8, 7, 3, 4, 1, 2, 11, 12},
{ 0, 9, 10, 1, 2, 11, 12, 7, 8, 6, 5, 4, 3},
{ 0, 7, 8, 12, 11, 3, 4, 6, 5, 9, 10, 1, 2},
{ 0, 2, 1, 12, 11, 6, 5, 10, 9, 8, 7, 4, 3},
{ 0, 2, 1, 4, 3, 8, 7, 6, 5, 12, 11, 10, 9},
{ 0, 9, 10, 4, 3, 6, 5, 11, 12, 1, 2, 7, 8},
{ 0, 7, 8, 9, 10, 6, 5, 1, 2, 3, 4, 12, 11},
{ 0, 2, 1, 8, 7, 10, 9, 4, 3, 6, 5, 12, 11},
{ 0, 11, 12, 4, 3, 9, 10, 8, 7, 5, 6, 1, 2},
{ 0, 10, 9, 5, 6, 3, 4, 8, 7, 2, 1, 12, 11},
{ 0, 8, 7, 4, 3, 11, 12, 2, 1, 10, 9, 5, 6},
{ 0, 3, 4, 1, 2, 7, 8, 5, 6, 10, 9, 12, 11},
{ 0, 2, 1, 10, 9, 12, 11, 8, 7, 4, 3, 6, 5},
{ 0, 12, 11, 7, 8, 6, 5, 3, 4, 10, 9, 2, 1},
{ 0, 4, 3, 2, 1, 6, 5, 8, 7, 11, 12, 9, 10},
{ 0, 2, 1, 6, 5, 4, 3, 12, 11, 10, 9, 8, 7},
{ 0, 5, 6, 11, 12, 1, 2, 8, 7, 10, 9, 3, 4},
{ 0, 6, 5, 4, 3, 2, 1, 9, 10, 7, 8, 12, 11},
};
//-----------------------------------------------------------------------------
// hcp
//-----------------------------------------------------------------------------
const int8_t mapping_hcp[NUM_HEX_MAPPINGS][PTM_MAX_POINTS] = {
{0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12},
{0, 5, 6, 1, 2, 3, 4, 9, 10, 12, 11, 8, 7},
{0, 3, 4, 5, 6, 1, 2, 12, 11, 7, 8, 10, 9},
{0, 4, 3, 2, 1, 6, 5, 11, 12, 10, 9, 7, 8},
{0, 2, 1, 6, 5, 4, 3, 8, 7, 11, 12, 9, 10},
{0, 6, 5, 4, 3, 2, 1, 10, 9, 8, 7, 12, 11} };
{ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12},
{ 0, 2, 7, 4, 5, 3, 8, 1, 9, 6, 12, 10, 11},
{ 0, 7, 1, 5, 3, 4, 9, 2, 6, 8, 11, 12, 10},
{ 0, 6, 9, 10, 11, 12, 1, 8, 7, 2, 3, 4, 5},
{ 0, 8, 6, 12, 10, 11, 2, 9, 1, 7, 4, 5, 3},
{ 0, 9, 8, 11, 12, 10, 7, 6, 2, 1, 5, 3, 4},
};
const int8_t template_indices_hcp[NUM_CONVENTIONAL_HEX_MAPPINGS] = {0, 1, 1, 0, 0, 1, 0, 0, 1, 1, 0, 1};
const int8_t mapping_hcp_conventional[NUM_CONVENTIONAL_HEX_MAPPINGS][PTM_MAX_POINTS] = {
{ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12},
{ 0, 2, 7, 4, 5, 3, 8, 1, 9, 6, 12, 10, 11},
{ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12},
{ 0, 2, 7, 4, 5, 3, 8, 1, 9, 6, 12, 10, 11},
{ 0, 7, 1, 5, 3, 4, 9, 2, 6, 8, 11, 12, 10},
{ 0, 8, 6, 12, 10, 11, 2, 9, 1, 7, 4, 5, 3},
{ 0, 6, 9, 10, 11, 12, 1, 8, 7, 2, 3, 4, 5},
{ 0, 8, 6, 12, 10, 11, 2, 9, 1, 7, 4, 5, 3},
{ 0, 6, 9, 10, 11, 12, 1, 8, 7, 2, 3, 4, 5},
{ 0, 9, 8, 11, 12, 10, 7, 6, 2, 1, 5, 3, 4},
{ 0, 9, 8, 11, 12, 10, 7, 6, 2, 1, 5, 3, 4},
{ 0, 7, 1, 5, 3, 4, 9, 2, 6, 8, 11, 12, 10},
};
//-----------------------------------------------------------------------------
// dcub
//-----------------------------------------------------------------------------
const int8_t mapping_dcub[NUM_DCUB_MAPPINGS][PTM_MAX_POINTS] = {
{0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16},
{0, 2, 1, 4, 3, 9, 8, 10, 6, 5, 7, 14, 16, 15, 11, 13, 12},
{0, 4, 3, 2, 1, 15, 16, 14, 12, 13, 11, 10, 8, 9, 7, 5, 6},
{0, 3, 4, 1, 2, 13, 12, 11, 16, 15, 14, 7, 6, 5, 10, 9, 8},
{0, 4, 2, 1, 3, 14, 15, 16, 9, 10, 8, 6, 5, 7, 12, 11, 13},
{0, 4, 1, 3, 2, 16, 14, 15, 7, 6, 5, 13, 11, 12, 9, 8, 10},
{0, 1, 4, 2, 3, 6, 7, 5, 14, 16, 15, 9, 10, 8, 13, 12, 11},
{0, 3, 1, 2, 4, 11, 13, 12, 5, 7, 6, 8, 9, 10, 16, 14, 15},
{0, 3, 2, 4, 1, 12, 11, 13, 10, 8, 9, 15, 14, 16, 5, 6, 7},
{0, 2, 4, 3, 1, 10, 9, 8, 15, 14, 16, 12, 13, 11, 6, 7, 5},
{0, 1, 3, 4, 2, 7, 5, 6, 13, 11, 12, 16, 15, 14, 8, 10, 9},
{0, 2, 3, 1, 4, 8, 10, 9, 11, 12, 13, 5, 7, 6, 15, 16, 14} };
{ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16},
{ 0, 1, 3, 4, 2, 7, 5, 6, 11, 13, 12, 14, 15, 16, 8, 10, 9},
{ 0, 4, 1, 3, 2, 16, 14, 15, 7, 6, 5, 12, 13, 11, 9, 8, 10},
{ 0, 2, 3, 1, 4, 8, 10, 9, 13, 12, 11, 6, 7, 5, 15, 16, 14},
{ 0, 4, 2, 1, 3, 14, 15, 16, 9, 10, 8, 7, 5, 6, 12, 13, 11},
{ 0, 3, 2, 4, 1, 12, 13, 11, 10, 8, 9, 16, 14, 15, 5, 6, 7},
{ 0, 3, 1, 2, 4, 13, 11, 12, 5, 7, 6, 10, 9, 8, 16, 14, 15},
{ 0, 2, 4, 3, 1, 10, 9, 8, 15, 14, 16, 13, 11, 12, 6, 7, 5},
{ 0, 1, 4, 2, 3, 6, 7, 5, 14, 16, 15, 8, 10, 9, 11, 12, 13},
{ 0, 2, 1, 4, 3, 9, 8, 10, 6, 5, 7, 15, 16, 14, 13, 11, 12},
{ 0, 4, 3, 2, 1, 15, 16, 14, 12, 11, 13, 9, 8, 10, 7, 5, 6},
{ 0, 3, 4, 1, 2, 11, 12, 13, 16, 15, 14, 5, 6, 7, 10, 9, 8},
};
const int8_t template_indices_dcub[NUM_CONVENTIONAL_DCUB_MAPPINGS] = {0, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 0, 1, 1, 0};
const int8_t mapping_dcub_conventional[NUM_CONVENTIONAL_DCUB_MAPPINGS][PTM_MAX_POINTS] = {
{ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16},
{ 0, 2, 1, 4, 3, 9, 8, 10, 6, 5, 7, 15, 16, 14, 13, 11, 12},
{ 0, 4, 1, 3, 2, 16, 14, 15, 7, 6, 5, 12, 13, 11, 9, 8, 10},
{ 0, 1, 3, 4, 2, 7, 5, 6, 11, 13, 12, 14, 15, 16, 8, 10, 9},
{ 0, 4, 2, 1, 3, 14, 15, 16, 9, 10, 8, 7, 5, 6, 12, 13, 11},
{ 0, 2, 3, 1, 4, 8, 10, 9, 13, 12, 11, 6, 7, 5, 15, 16, 14},
{ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16},
{ 0, 1, 3, 4, 2, 7, 5, 6, 11, 13, 12, 14, 15, 16, 8, 10, 9},
{ 0, 4, 1, 3, 2, 16, 14, 15, 7, 6, 5, 12, 13, 11, 9, 8, 10},
{ 0, 2, 3, 1, 4, 8, 10, 9, 13, 12, 11, 6, 7, 5, 15, 16, 14},
{ 0, 4, 2, 1, 3, 14, 15, 16, 9, 10, 8, 7, 5, 6, 12, 13, 11},
{ 0, 3, 2, 4, 1, 12, 13, 11, 10, 8, 9, 16, 14, 15, 5, 6, 7},
{ 0, 3, 1, 2, 4, 13, 11, 12, 5, 7, 6, 10, 9, 8, 16, 14, 15},
{ 0, 2, 4, 3, 1, 10, 9, 8, 15, 14, 16, 13, 11, 12, 6, 7, 5},
{ 0, 1, 4, 2, 3, 6, 7, 5, 14, 16, 15, 8, 10, 9, 11, 12, 13},
{ 0, 2, 1, 4, 3, 9, 8, 10, 6, 5, 7, 15, 16, 14, 13, 11, 12},
{ 0, 2, 4, 3, 1, 10, 9, 8, 15, 14, 16, 13, 11, 12, 6, 7, 5},
{ 0, 1, 4, 2, 3, 6, 7, 5, 14, 16, 15, 8, 10, 9, 11, 12, 13},
{ 0, 3, 2, 4, 1, 12, 13, 11, 10, 8, 9, 16, 14, 15, 5, 6, 7},
{ 0, 3, 1, 2, 4, 13, 11, 12, 5, 7, 6, 10, 9, 8, 16, 14, 15},
{ 0, 4, 3, 2, 1, 15, 16, 14, 12, 11, 13, 9, 8, 10, 7, 5, 6},
{ 0, 4, 3, 2, 1, 15, 16, 14, 12, 11, 13, 9, 8, 10, 7, 5, 6},
{ 0, 3, 4, 1, 2, 11, 12, 13, 16, 15, 14, 5, 6, 7, 10, 9, 8},
{ 0, 3, 4, 1, 2, 11, 12, 13, 16, 15, 14, 5, 6, 7, 10, 9, 8},
};
//-----------------------------------------------------------------------------
// dhex
//-----------------------------------------------------------------------------
const int8_t mapping_dhex[NUM_DHEX_MAPPINGS][PTM_MAX_POINTS] = {
{0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16},
{0, 1, 3, 4, 2, 6, 7, 5, 11, 13, 12, 14, 16, 15, 8, 9, 10},
{0, 1, 4, 2, 3, 7, 5, 6, 14, 15, 16, 8, 10, 9, 11, 13, 12} };
{ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16},
{ 0, 3, 1, 2, 4, 12, 11, 13, 5, 6, 7, 9, 8, 10, 16, 14, 15},
{ 0, 2, 3, 1, 4, 8, 9, 10, 12, 11, 13, 6, 5, 7, 15, 16, 14},
};
const int8_t template_indices_dhex[NUM_CONVENTIONAL_DHEX_MAPPINGS] = {0, 1, 1, 0, 0, 3, 2, 2, 3, 3, 2, 1};
const int8_t mapping_dhex_conventional[NUM_CONVENTIONAL_DHEX_MAPPINGS][PTM_MAX_POINTS] = {
{ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16},
{ 0, 3, 1, 2, 4, 12, 11, 13, 5, 6, 7, 9, 8, 10, 16, 14, 15},
{ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16},
{ 0, 3, 1, 2, 4, 12, 11, 13, 5, 6, 7, 9, 8, 10, 16, 14, 15},
{ 0, 2, 3, 1, 4, 8, 9, 10, 12, 11, 13, 6, 5, 7, 15, 16, 14},
{ 0, 2, 3, 1, 4, 8, 9, 10, 12, 11, 13, 6, 5, 7, 15, 16, 14},
{ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16},
{ 0, 2, 3, 1, 4, 8, 9, 10, 12, 11, 13, 6, 5, 7, 15, 16, 14},
{ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16},
{ 0, 3, 1, 2, 4, 12, 11, 13, 5, 6, 7, 9, 8, 10, 16, 14, 15},
{ 0, 3, 1, 2, 4, 12, 11, 13, 5, 6, 7, 9, 8, 10, 16, 14, 15},
{ 0, 2, 3, 1, 4, 8, 9, 10, 12, 11, 13, 6, 5, 7, 15, 16, 14},
};
//-----------------------------------------------------------------------------
// graphene
//-----------------------------------------------------------------------------
const int8_t mapping_graphene[NUM_GRAPHENE_MAPPINGS][PTM_MAX_POINTS] = {
{ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9},
{ 0, 2, 3, 1, 6, 7, 8, 9, 4, 5},
{ 0, 3, 1, 2, 8, 9, 4, 5, 6, 7},
{ 0, 2, 1, 3, 7, 6, 5, 4, 9, 8},
{ 0, 3, 2, 1, 9, 8, 7, 6, 5, 4},
{ 0, 1, 3, 2, 5, 4, 9, 8, 7, 6},
};
const int8_t template_indices_graphene[NUM_CONVENTIONAL_GRAPHENE_MAPPINGS] = {0, 1, 1, 0, 0, 1, 0, 0, 1, 1, 0, 1};
const int8_t mapping_graphene_conventional[NUM_CONVENTIONAL_GRAPHENE_MAPPINGS][PTM_MAX_POINTS] = {
{ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9},
{ 0, 2, 3, 1, 6, 7, 8, 9, 4, 5},
{ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9},
{ 0, 2, 3, 1, 6, 7, 8, 9, 4, 5},
{ 0, 3, 1, 2, 8, 9, 4, 5, 6, 7},
{ 0, 3, 2, 1, 9, 8, 7, 6, 5, 4},
{ 0, 2, 1, 3, 7, 6, 5, 4, 9, 8},
{ 0, 3, 2, 1, 9, 8, 7, 6, 5, 4},
{ 0, 2, 1, 3, 7, 6, 5, 4, 9, 8},
{ 0, 1, 3, 2, 5, 4, 9, 8, 7, 6},
{ 0, 1, 3, 2, 5, 4, 9, 8, 7, 6},
{ 0, 3, 1, 2, 8, 9, 4, 5, 6, 7},
};
}
#endif

File diff suppressed because it is too large Load Diff

View File

@ -1,3 +1,12 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#ifndef PTM_GRAPH_DATA_H
#define PTM_GRAPH_DATA_H

View File

@ -1,4 +1,13 @@
#include <cstring>
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#include <string.h>
#include <algorithm>
#include "ptm_graph_tools.h"
#include "ptm_constants.h"

View File

@ -1,3 +1,12 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#ifndef PTM_GRAPH_TOOLS_H
#define PTM_GRAPH_TOOLS_H

View File

@ -1,125 +1,180 @@
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cmath>
#include <cfloat>
#include <cassert>
#include <algorithm>
#include "ptm_convex_hull_incremental.h"
#include "ptm_graph_data.h"
#include "ptm_deformation_gradient.h"
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#include "ptm_alloy_types.h"
#include "ptm_constants.h"
#include "ptm_convex_hull_incremental.h"
#include "ptm_deformation_gradient.h"
#include "ptm_functions.h"
#include "ptm_graph_data.h"
#include "ptm_initialize_data.h"
#include "ptm_neighbour_ordering.h"
#include "ptm_normalize_vertices.h"
#include "ptm_quat.h"
#include "ptm_polar.h"
#include "ptm_initialize_data.h"
#include "ptm_quat.h"
#include "ptm_structure_matcher.h"
#include "ptm_functions.h"
#include "ptm_constants.h"
#include <algorithm>
#include <cassert>
#include <cfloat>
#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <string.h>
static double calculate_interatomic_distance(int type, double scale) {
assert(type >= 1 && type <= 8);
//todo: verify that c == norm(template[1])
static double calculate_interatomic_distance(int type, double scale)
{
assert(type >= 1 && type <= 7);
double c[8] = {0, 1, 1, (7. - 3.5 * sqrt(3)), 1, 1, sqrt(3) * 4. / (6 * sqrt(2) + sqrt(3)), sqrt(3) * 4. / (6 * sqrt(2) + sqrt(3))};
// these values should be equal to norm(template[1])
double c[9] = {0,
1,
1,
(7. - 3.5 * sqrt(3)),
1,
1,
sqrt(3) * 4. / (6 * sqrt(2) + sqrt(3)),
sqrt(3) * 4. / (6 * sqrt(2) + sqrt(3)),
-3. / 11 + 6 * sqrt(3) / 11};
return c[type] / scale;
}
static double calculate_lattice_constant(int type, double interatomic_distance)
{
assert(type >= 1 && type <= 7);
double c[8] = {0, 2 / sqrt(2), 2 / sqrt(2), 2. / sqrt(3), 2 / sqrt(2), 1, 4 / sqrt(3), 4 / sqrt(3)};
static double calculate_lattice_constant(int type,
double interatomic_distance) {
assert(type >= 1 && type <= 8);
double c[9] = {0, 2 / sqrt(2), 2 / sqrt(2), 2. / sqrt(3), 2 / sqrt(2),
1, 4 / sqrt(3), 4 / sqrt(3), sqrt(3)};
return c[type] * interatomic_distance;
}
static int rotate_into_fundamental_zone(int type, double* q)
{
if (type == PTM_MATCH_SC) return ptm::rotate_quaternion_into_cubic_fundamental_zone(q);
if (type == PTM_MATCH_FCC) return ptm::rotate_quaternion_into_cubic_fundamental_zone(q);
if (type == PTM_MATCH_BCC) return ptm::rotate_quaternion_into_cubic_fundamental_zone(q);
if (type == PTM_MATCH_ICO) return ptm::rotate_quaternion_into_icosahedral_fundamental_zone(q);
if (type == PTM_MATCH_HCP) return ptm::rotate_quaternion_into_hcp_fundamental_zone(q);
if (type == PTM_MATCH_DCUB) return ptm::rotate_quaternion_into_diamond_cubic_fundamental_zone(q);
if (type == PTM_MATCH_DHEX) return ptm::rotate_quaternion_into_diamond_hexagonal_fundamental_zone(q);
static int rotate_into_fundamental_zone(int type,
bool output_conventional_orientation,
double *q) {
if (type == PTM_MATCH_SC)
return ptm::rotate_quaternion_into_cubic_fundamental_zone(q);
if (type == PTM_MATCH_FCC)
return ptm::rotate_quaternion_into_cubic_fundamental_zone(q);
if (type == PTM_MATCH_BCC)
return ptm::rotate_quaternion_into_cubic_fundamental_zone(q);
if (type == PTM_MATCH_ICO)
return ptm::rotate_quaternion_into_icosahedral_fundamental_zone(q);
if (type == PTM_MATCH_HCP || type == PTM_MATCH_GRAPHENE) {
if (output_conventional_orientation) {
return ptm::rotate_quaternion_into_hcp_conventional_fundamental_zone(q);
} else {
return ptm::rotate_quaternion_into_hcp_fundamental_zone(q);
}
}
if (type == PTM_MATCH_DCUB) {
if (output_conventional_orientation) {
return ptm::rotate_quaternion_into_cubic_fundamental_zone(q);
} else {
return ptm::rotate_quaternion_into_diamond_cubic_fundamental_zone(q);
}
}
if (type == PTM_MATCH_DHEX) {
if (output_conventional_orientation) {
return ptm::rotate_quaternion_into_hcp_conventional_fundamental_zone(q);
} else {
return ptm::rotate_quaternion_into_diamond_hexagonal_fundamental_zone(q);
}
}
return -1;
}
static void order_points(ptm_local_handle_t local_handle, int num_points, double (*unpermuted_points)[3], int32_t* unpermuted_numbers, bool topological_ordering,
int8_t* ordering, double (*points)[3], int32_t* numbers)
{
if (topological_ordering)
{
double normalized_points[PTM_MAX_INPUT_POINTS][3];
ptm::normalize_vertices(num_points, unpermuted_points, normalized_points);
int ret = ptm::calculate_neighbour_ordering((void*)local_handle, num_points, (const double (*)[3])normalized_points, ordering);
if (ret != 0)
topological_ordering = false;
}
if (!topological_ordering)
for (int i=0;i<num_points;i++)
ordering[i] = i;
for (int i=0;i<num_points;i++)
{
memcpy(points[i], &unpermuted_points[ordering[i]], 3 * sizeof(double));
if (unpermuted_numbers != NULL)
numbers[i] = unpermuted_numbers[ordering[i]];
}
}
static void output_data(ptm::result_t* res, int num_points, int32_t* unpermuted_numbers, double (*points)[3], int32_t* numbers, int8_t* ordering,
int32_t* p_type, int32_t* p_alloy_type, double* p_scale, double* p_rmsd, double* q, double* F, double* F_res,
double* U, double* P, int8_t* mapping, double* p_interatomic_distance, double* p_lattice_constant)
{
*p_type = PTM_MATCH_NONE;
if (p_alloy_type != NULL)
*p_alloy_type = PTM_ALLOY_NONE;
if (mapping != NULL)
memset(mapping, -1, num_points * sizeof(int8_t));
const ptm::refdata_t* ref = res->ref_struct;
static void output_data(ptm::result_t *res, double (*points)[3],
int32_t *numbers, size_t *ordering,
bool output_conventional_orientation, int32_t *p_type,
int32_t *p_alloy_type, double *p_scale, double *p_rmsd,
double *q, double *F, double *F_res, double *U,
double *P, double *p_interatomic_distance,
double *p_lattice_constant, size_t *output_indices) {
const ptm::refdata_t *ref = res->ref_struct;
if (ref == NULL)
return;
*p_type = ref->type;
if (p_alloy_type != NULL && unpermuted_numbers != NULL)
if (p_alloy_type != NULL)
*p_alloy_type = ptm::find_alloy_type(ref, res->mapping, numbers);
int bi = rotate_into_fundamental_zone(ref->type, res->q);
int8_t temp[PTM_MAX_POINTS];
for (int i=0;i<ref->num_nbrs+1;i++)
temp[ref->mapping[bi][i]] = res->mapping[i];
memcpy(res->mapping, temp, (ref->num_nbrs+1) * sizeof(int8_t));
if (F != NULL && F_res != NULL)
memset(temp, -1, PTM_MAX_POINTS * sizeof(int8_t));
int bi = rotate_into_fundamental_zone(ref->type,
output_conventional_orientation, res->q);
//todo: return if bi == -1
if (bi != -1)
{
if (output_conventional_orientation & (ref->type == PTM_MATCH_HCP || ref->type == PTM_MATCH_GRAPHENE || ref->type == PTM_MATCH_DCUB || ref->type == PTM_MATCH_DHEX))
{
for (int i = 0; i < ref->num_nbrs + 1; i++)
temp[ref->mapping_conventional[bi][i]] = res->mapping[i];
}
else
{
for (int i = 0; i < ref->num_nbrs + 1; i++)
temp[ref->mapping[bi][i]] = res->mapping[i];
}
}
memcpy(res->mapping, temp, (ref->num_nbrs + 1) * sizeof(int8_t));
if (F != NULL && F_res != NULL) {
double scaled_points[PTM_MAX_INPUT_POINTS][3];
ptm::subtract_barycentre(ref->num_nbrs + 1, points, scaled_points);
for (int i = 0;i<ref->num_nbrs + 1;i++)
{
for (int i = 0; i < ref->num_nbrs + 1; i++) {
scaled_points[i][0] *= res->scale;
scaled_points[i][1] *= res->scale;
scaled_points[i][2] *= res->scale;
}
ptm::calculate_deformation_gradient(ref->num_nbrs + 1, ref->points, res->mapping, scaled_points, ref->penrose, F, F_res);
const double (*ref_template)[3] = ref->points;
const double (*ref_penrose)[3] = ref->penrose;
if (output_conventional_orientation & (ref->type == PTM_MATCH_HCP || ref->type == PTM_MATCH_GRAPHENE || ref->type == PTM_MATCH_DCUB || ref->type == PTM_MATCH_DHEX))
{
if (ref->template_indices[bi] == 1)
{
ref_template = ref->points_alt1;
ref_penrose = ref->penrose_alt1;
}
else if (ref->template_indices[bi] == 2)
{
ref_template = ref->points_alt2;
ref_penrose = ref->penrose_alt2;
}
else if (ref->template_indices[bi] == 3)
{
ref_template = ref->points_alt3;
ref_penrose = ref->penrose_alt3;
}
}
ptm::calculate_deformation_gradient(ref->num_nbrs + 1, ref_template,
res->mapping, scaled_points, ref_penrose,
F, F_res);
if (ref->type == PTM_MATCH_GRAPHENE) // hack for pseudo-2d structures
F[8] = 1;
if (P != NULL && U != NULL)
ptm::polar_decomposition_3x3(F, false, U, P);
}
if (mapping != NULL)
for (int i=0;i<ref->num_nbrs + 1;i++)
mapping[i] = ordering[res->mapping[i]];
if (output_indices != NULL)
for (int i = 0; i < ref->num_nbrs + 1; i++)
output_indices[i] = ordering[res->mapping[i]];
double interatomic_distance = calculate_interatomic_distance(ref->type, res->scale);
double lattice_constant = calculate_lattice_constant(ref->type, interatomic_distance);
double interatomic_distance =
calculate_interatomic_distance(ref->type, res->scale);
double lattice_constant =
calculate_lattice_constant(ref->type, interatomic_distance);
if (p_interatomic_distance != NULL)
*p_interatomic_distance = interatomic_distance;
@ -132,67 +187,73 @@ static void output_data(ptm::result_t* res, int num_points, int32_t* unpermuted_
memcpy(q, res->q, 4 * sizeof(double));
}
extern bool ptm_initialized;
int ptm_index( ptm_local_handle_t local_handle, int32_t flags,
int num_points, double (*unpermuted_points)[3], int32_t* unpermuted_numbers, bool topological_ordering,
int32_t* p_type, int32_t* p_alloy_type, double* p_scale, double* p_rmsd, double* q, double* F, double* F_res,
double* U, double* P, int8_t* mapping, double* p_interatomic_distance, double* p_lattice_constant)
{
int ptm_index(ptm_local_handle_t local_handle, size_t atom_index,
int(get_neighbours)(void *vdata, size_t central_index, size_t atom_index, int num,
size_t *nbr_indices, int32_t *numbers,
double (*nbr_pos)[3]),
void *nbrlist, int32_t flags,
bool output_conventional_orientation, int32_t *p_type,
int32_t *p_alloy_type, double *p_scale, double *p_rmsd, double *q,
double *F, double *F_res, double *U, double *P,
double *p_interatomic_distance, double *p_lattice_constant,
size_t *output_indices) {
assert(ptm_initialized);
assert(num_points <= PTM_MAX_INPUT_POINTS);
if (flags & PTM_CHECK_SC)
assert(num_points >= PTM_NUM_POINTS_SC);
if (flags & PTM_CHECK_BCC)
assert(num_points >= PTM_NUM_POINTS_BCC);
if (flags & (PTM_CHECK_FCC | PTM_CHECK_HCP | PTM_CHECK_ICO))
assert(num_points >= PTM_NUM_POINTS_FCC);
if (flags & (PTM_CHECK_DCUB | PTM_CHECK_DHEX))
assert(num_points >= PTM_NUM_POINTS_DCUB);
int ret = 0;
ptm::result_t res;
res.ref_struct = NULL;
res.rmsd = INFINITY;
int8_t ordering[PTM_MAX_INPUT_POINTS];
double points[PTM_MAX_POINTS][3];
int32_t numbers[PTM_MAX_POINTS];
size_t ordering[PTM_MAX_INPUT_POINTS];
int32_t numbers[PTM_MAX_INPUT_POINTS];
double points[PTM_MAX_INPUT_POINTS][3];
int8_t dordering[PTM_MAX_INPUT_POINTS];
double dpoints[PTM_MAX_POINTS][3];
int32_t dnumbers[PTM_MAX_POINTS];
size_t dordering[PTM_MAX_INPUT_POINTS];
int32_t dnumbers[PTM_MAX_INPUT_POINTS];
double dpoints[PTM_MAX_INPUT_POINTS][3];
size_t gordering[PTM_MAX_INPUT_POINTS];
int32_t gnumbers[PTM_MAX_INPUT_POINTS];
double gpoints[PTM_MAX_INPUT_POINTS][3];
ptm::convexhull_t ch;
double ch_points[PTM_MAX_INPUT_POINTS][3];
int num_lpoints = 0;
if (flags & (PTM_CHECK_SC | PTM_CHECK_FCC | PTM_CHECK_HCP | PTM_CHECK_ICO | PTM_CHECK_BCC))
{
int num_lpoints = std::min(std::min(PTM_MAX_POINTS, 20), num_points);
order_points(local_handle, num_lpoints, unpermuted_points, unpermuted_numbers, topological_ordering, ordering, points, numbers);
ptm::normalize_vertices(num_lpoints, points, ch_points);
ch.ok = false;
if (flags & PTM_CHECK_SC)
ret = match_general(&ptm::structure_sc, ch_points, points, &ch, &res);
if (flags & (PTM_CHECK_SC | PTM_CHECK_FCC | PTM_CHECK_HCP | PTM_CHECK_ICO |
PTM_CHECK_BCC)) {
int min_points = PTM_NUM_POINTS_SC;
if (flags & (PTM_CHECK_FCC | PTM_CHECK_HCP | PTM_CHECK_ICO))
ret = match_fcc_hcp_ico(ch_points, points, flags, &ch, &res);
min_points = PTM_NUM_POINTS_FCC;
if (flags & PTM_CHECK_BCC)
ret = match_general(&ptm::structure_bcc, ch_points, points, &ch, &res);
min_points = PTM_NUM_POINTS_BCC;
num_lpoints = ptm::calculate_neighbour_ordering(
local_handle, atom_index, min_points, get_neighbours, nbrlist, false,
ordering, points, numbers);
if (num_lpoints >= min_points) {
ptm::normalize_vertices(num_lpoints, points, ch_points);
ch.ok = false;
if (flags & PTM_CHECK_SC)
ret = match_general(&ptm::structure_sc, ch_points, points, &ch, &res);
if (flags & (PTM_CHECK_FCC | PTM_CHECK_HCP | PTM_CHECK_ICO))
ret = match_fcc_hcp_ico(ch_points, points, flags, &ch, &res);
if (flags & PTM_CHECK_BCC)
ret = match_general(&ptm::structure_bcc, ch_points, points, &ch, &res);
}
}
if (flags & (PTM_CHECK_DCUB | PTM_CHECK_DHEX))
{
ret = ptm::calculate_diamond_neighbour_ordering(num_points, unpermuted_points, unpermuted_numbers, dordering, dpoints, dnumbers);
if (ret == 0)
{
if (flags & (PTM_CHECK_DCUB | PTM_CHECK_DHEX)) {
const int num_inner = 4, num_outer = 3, max_snbrs = 12;
ret = ptm::calculate_two_shell_neighbour_ordering(
(void *)local_handle, atom_index, get_neighbours, nbrlist, num_inner,
num_outer, max_snbrs, false, dordering, dpoints, dnumbers);
if (ret == 0) {
ptm::normalize_vertices(PTM_NUM_NBRS_DCUB + 1, dpoints, ch_points);
ch.ok = false;
@ -200,19 +261,42 @@ int ptm_index( ptm_local_handle_t local_handle, int32_t flags,
}
}
if (res.ref_struct != NULL && (res.ref_struct->type == PTM_MATCH_DCUB || res.ref_struct->type == PTM_MATCH_DHEX))
{
output_data( &res, num_points, unpermuted_numbers, dpoints, dnumbers, dordering,
p_type, p_alloy_type, p_scale, p_rmsd, q, F, F_res,
U, P, mapping, p_interatomic_distance, p_lattice_constant);
if (flags & PTM_CHECK_GRAPHENE) {
const int num_inner = 3, num_outer = 2, max_snbrs = 12;
ret = ptm::calculate_two_shell_neighbour_ordering(
(void *)local_handle, atom_index, get_neighbours, nbrlist, num_inner,
num_outer, max_snbrs, true, gordering, gpoints, gnumbers);
if (ret == 0) {
ret = match_graphene(gpoints, &res);
}
}
else
{
output_data( &res, num_points, unpermuted_numbers, points, numbers, ordering,
p_type, p_alloy_type, p_scale, p_rmsd, q, F, F_res,
U, P, mapping, p_interatomic_distance, p_lattice_constant);
*p_type = PTM_MATCH_NONE;
if (p_alloy_type != NULL)
*p_alloy_type = PTM_ALLOY_NONE;
if (output_indices != NULL)
memset(output_indices, -1, PTM_MAX_INPUT_POINTS * sizeof(size_t));
if (res.ref_struct == NULL)
return PTM_NO_ERROR;
if (res.ref_struct->type == PTM_MATCH_DCUB ||
res.ref_struct->type == PTM_MATCH_DHEX) {
output_data(&res, dpoints, dnumbers, dordering,
output_conventional_orientation, p_type, p_alloy_type, p_scale,
p_rmsd, q, F, F_res, U, P, p_interatomic_distance,
p_lattice_constant, output_indices);
} else if (res.ref_struct->type == PTM_MATCH_GRAPHENE) {
output_data(&res, gpoints, gnumbers, gordering,
output_conventional_orientation, p_type, p_alloy_type, p_scale,
p_rmsd, q, F, F_res, U, P, p_interatomic_distance,
p_lattice_constant, output_indices);
} else {
output_data(&res, points, numbers, ordering, output_conventional_orientation,
p_type, p_alloy_type, p_scale, p_rmsd, q, F, F_res, U, P,
p_interatomic_distance, p_lattice_constant, output_indices);
}
return PTM_NO_ERROR;
}

View File

@ -1,6 +1,15 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <string.h>
#include <cmath>
#include <cfloat>
#include <cassert>

View File

@ -1,3 +1,12 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#ifndef PTM_INITIALIZE_DATA_H
#define PTM_INITIALIZE_DATA_H
@ -9,6 +18,7 @@
#include "ptm_neighbour_ordering.h"
#include "ptm_canonical_coloured.h"
#include "ptm_convex_hull_incremental.h"
#include "ptm_alt_templates.h"
namespace ptm {
@ -20,23 +30,181 @@ typedef struct
int num_facets;
int max_degree;
int num_graphs;
int num_mappings;
graph_t* graphs;
const double (*points)[3];
const double (*points_alt1)[3];
const double (*points_alt2)[3];
const double (*points_alt3)[3];
const double (*penrose)[3];
const double (*penrose_alt1)[3];
const double (*penrose_alt2)[3];
const double (*penrose_alt3)[3];
int num_mappings;
const int8_t (*mapping)[PTM_MAX_POINTS];
const int8_t (*mapping_conventional)[PTM_MAX_POINTS];
const int8_t *template_indices;
} refdata_t;
//refdata_t structure_sc = { .type = PTM_MATCH_SC, .num_nbrs = 6, .num_facets = 8, .max_degree = 4, .num_graphs = NUM_SC_GRAPHS, .graphs = graphs_sc, .points = ptm_template_sc, .penrose = penrose_sc , .mapping = mapping_sc };
const refdata_t structure_sc = { PTM_MATCH_SC, 6, 8, 4, NUM_SC_GRAPHS, NUM_CUBIC_MAPPINGS, graphs_sc, ptm_template_sc, penrose_sc, mapping_sc };
const refdata_t structure_fcc = { PTM_MATCH_FCC, 12, 20, 6, NUM_FCC_GRAPHS, NUM_CUBIC_MAPPINGS, graphs_fcc, ptm_template_fcc, penrose_fcc, mapping_fcc };
const refdata_t structure_hcp = { PTM_MATCH_HCP, 12, 20, 6, NUM_HCP_GRAPHS, NUM_HEX_MAPPINGS, graphs_hcp, ptm_template_hcp, penrose_hcp, mapping_hcp };
const refdata_t structure_ico = { PTM_MATCH_ICO, 12, 20, 6, NUM_ICO_GRAPHS, NUM_ICO_MAPPINGS, graphs_ico, ptm_template_ico, penrose_ico, mapping_ico };
const refdata_t structure_bcc = { PTM_MATCH_BCC, 14, 24, 8, NUM_BCC_GRAPHS, NUM_CUBIC_MAPPINGS, graphs_bcc, ptm_template_bcc, penrose_bcc, mapping_bcc };
const refdata_t structure_dcub = { PTM_MATCH_DCUB, 16, 28, 8, NUM_DCUB_GRAPHS, NUM_DCUB_MAPPINGS, graphs_dcub, ptm_template_dcub, penrose_dcub, mapping_dcub };
const refdata_t structure_dhex = { PTM_MATCH_DHEX, 16, 28, 8, NUM_DHEX_GRAPHS, NUM_DHEX_MAPPINGS, graphs_dhex, ptm_template_dhex, penrose_dhex, mapping_dhex };
const refdata_t structure_sc = { PTM_MATCH_SC, //.type
6, //.num_nbrs
8, //.num_facets
4, //.max_degree
NUM_SC_GRAPHS, //.num_graphs
graphs_sc, //.graphs
ptm_template_sc, //.points
NULL, //.points_alt1
NULL, //.points_alt2
NULL, //.points_alt3
penrose_sc, //.penrose
NULL, //.penrose_alt1
NULL, //.penrose_alt2
NULL, //.penrose_alt3
NUM_CUBIC_MAPPINGS, //.num_mappings
mapping_sc, //.mapping
NULL, //.mapping_conventional
NULL, //.template_indices
};
const refdata_t structure_fcc = { PTM_MATCH_FCC, //.type
12, //.num_nbrs
20, //.num_facets
6, //.max_degree
NUM_FCC_GRAPHS, //.num_graphs
graphs_fcc, //.graphs
ptm_template_fcc, //.points
NULL, //.points_alt1
NULL, //.points_alt2
NULL, //.points_alt3
penrose_fcc, //.penrose
NULL, //.penrose_alt1
NULL, //.penrose_alt2
NULL, //.penrose_alt3
NUM_CUBIC_MAPPINGS, //.num_mappings
mapping_fcc, //.mapping
NULL, //.mapping_conventional
NULL, //.template_indices
};
const refdata_t structure_hcp = { PTM_MATCH_HCP, //.type
12, //.num_nbrs
20, //.num_facets
6, //.max_degree
NUM_HCP_GRAPHS, //.num_graphs
graphs_hcp, //.graphs
ptm_template_hcp, //.points
ptm_template_hcp_alt1, //.points_alt1
NULL, //.points_alt2
NULL, //.points_alt3
penrose_hcp, //.penrose
penrose_hcp_alt1, //.penrose_alt1
NULL, //.penrose_alt2
NULL, //.penrose_alt3
NUM_HEX_MAPPINGS, //.num_mappings
mapping_hcp, //.mapping
mapping_hcp_conventional, //.mapping_conventional
template_indices_hcp, //.template_indices
};
const refdata_t structure_ico = { PTM_MATCH_ICO, //.type
12, //.num_nbrs
20, //.num_facets
6, //.max_degree
NUM_ICO_GRAPHS, //.num_graphs
graphs_ico, //.graphs
ptm_template_ico, //.points
NULL, //.points_alt1
NULL, //.points_alt2
NULL, //.points_alt3
penrose_ico, //.penrose
NULL, //.penrose_alt1
NULL, //.penrose_alt2
NULL, //.penrose_alt3
NUM_ICO_MAPPINGS, //.num_mappings
mapping_ico, //.mapping
NULL, //.mapping_conventional
NULL, //.template_indices
};
const refdata_t structure_bcc = { PTM_MATCH_BCC, //.type
14, //.num_nbrs
24, //.num_facets
8, //.max_degree
NUM_BCC_GRAPHS, //.num_graphs
graphs_bcc, //.graphs
ptm_template_bcc, //.points
NULL, //.points_alt1
NULL, //.points_alt2
NULL, //.points_alt3
penrose_bcc, //.penrose
NULL, //.penrose_alt1
NULL, //.penrose_alt2
NULL, //.penrose_alt3
NUM_CUBIC_MAPPINGS, //.num_mappings
mapping_bcc, //.mapping
NULL, //.mapping_conventional
NULL, //.template_indices
};
const refdata_t structure_dcub = { PTM_MATCH_DCUB, //.type
16, //.num_nbrs
28, //.num_facets
8, //.max_degree
NUM_DCUB_GRAPHS, //.num_graphs
graphs_dcub, //.graphs
ptm_template_dcub, //.points
ptm_template_dcub_alt1, //.points_alt1
NULL, //.points_alt2
NULL, //.points_alt3
penrose_dcub, //.penrose
penrose_dcub_alt1, //.penrose_alt1
NULL, //.penrose_alt2
NULL, //.penrose_alt3
NUM_DCUB_MAPPINGS, //.num_mappings
mapping_dcub, //.mapping
mapping_dcub_conventional, //.mapping_conventional
template_indices_dcub, //.template_indices
};
const refdata_t structure_dhex = { PTM_MATCH_DHEX, //.type
16, //.num_nbrs
28, //.num_facets
8, //.max_degree
NUM_DHEX_GRAPHS, //.num_graphs
graphs_dhex, //.graphs
ptm_template_dhex, //.points
ptm_template_dhex_alt1, //.points_alt1
ptm_template_dhex_alt2, //.points_alt2
ptm_template_dhex_alt3, //.points_alt3
penrose_dhex, //.penrose
penrose_dhex_alt1, //.penrose_alt1
penrose_dhex_alt2, //.penrose_alt2
penrose_dhex_alt3, //.penrose_alt3
NUM_DHEX_MAPPINGS, //.num_mappings
mapping_dhex, //.mapping
mapping_dhex_conventional, //.mapping_conventional
template_indices_dhex, //.template_indices
};
const refdata_t structure_graphene = { PTM_MATCH_GRAPHENE, //.type
9, //.num_nbrs
-1, //.num_facets
-1, //.max_degree
-1, //.num_graphs
NULL, //.graphs
ptm_template_graphene, //.points
ptm_template_graphene_alt1, //.points_alt1
NULL, //.points_alt2
NULL, //.points_alt3
penrose_graphene, //.penrose
penrose_graphene_alt1, //.penrose_alt1
NULL, //.penrose_alt2
NULL, //.penrose_alt3
-1, //.num_mappings
mapping_graphene, //.mapping
mapping_graphene_conventional, //.mapping_conventional
template_indices_graphene, //.template_indices
};
}
#ifdef __cplusplus

View File

@ -1,11 +1,24 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
//todo: normalize vertices
#include <cstdlib>
#include <cmath>
#include <cstring>
#include <cassert>
#include <algorithm>
#include <set>
#include "ptm_constants.h"
#include "ptm_voronoi_cell.h"
#include "ptm_neighbour_ordering.h"
#include "ptm_normalize_vertices.h"
namespace ptm {
@ -15,8 +28,20 @@ typedef struct
double area;
double dist;
int index;
int inner;
int32_t number;
double offset[3];
} sorthelper_t;
typedef struct
{
size_t index;
int32_t number;
double area;
double offset[3];
} solidnbr_t;
static bool sorthelper_compare(sorthelper_t const& a, sorthelper_t const& b)
{
if (a.area > b.area)
@ -31,10 +56,37 @@ static bool sorthelper_compare(sorthelper_t const& a, sorthelper_t const& b)
return false;
}
//todo: change voronoi code to return errors rather than exiting
static int calculate_voronoi_face_areas(int num_points, const double (*_points)[3], double* normsq, double max_norm, ptm_voro::voronoicell_neighbor* v, std::vector<int>& nbr_indices, std::vector<double>& face_areas)
static double dot_product(double* a, double* b)
{
const double k = 1000 * max_norm; //todo: reduce this constant
return a[0] * b[0] + a[1] * b[1] + a[2] * b[2];
}
static void cross_product(double* a, double* b, double* c)
{
c[0] = a[1] * b[2] - a[2] * b[1];
c[1] = a[2] * b[0] - a[0] * b[2];
c[2] = a[0] * b[1] - a[1] * b[0];
}
static double calculate_solid_angle(double* R1, double* R2, double* R3) //norms of R1-R3 must be 1
{
double R2R3[3];
cross_product(R2, R3, R2R3);
double numerator = dot_product(R1, R2R3);
double r1r2 = dot_product(R1, R2);
double r2r3 = dot_product(R2, R3);
double r3r1 = dot_product(R3, R1);
double denominator = 1 + r1r2 + r3r1 + r2r3;
return fabs(2 * atan2(numerator, denominator));
}
//todo: change voronoi code to return errors rather than exiting
static int calculate_voronoi_face_areas(int num_points, const double (*_points)[3], double* normsq, double max_norm, ptm_voro::voronoicell_neighbor* v, bool calc_solid_angles,
std::vector<int>& nbr_indices, std::vector<double>& face_areas)
{
const double k = 10 * max_norm;
v->init(-k,k,-k,k,-k,k);
for (int i=1;i<num_points;i++)
@ -46,11 +98,84 @@ static int calculate_voronoi_face_areas(int num_points, const double (*_points)[
}
v->neighbors(nbr_indices);
v->face_areas(face_areas);
return 0;
//v->face_areas(face_areas);
if (!calc_solid_angles)
{
v->face_areas(face_areas);
return 0;
}
else
{
std::vector<int> face_vertices;
std::vector<double> vertices;
v->face_vertices(face_vertices);
v->vertices(0, 0, 0, vertices);
size_t num_vertices = vertices.size() / 3;
for (size_t i=0;i<num_vertices;i++)
{
double x = vertices[i * 3 + 0];
double y = vertices[i * 3 + 1];
double z = vertices[i * 3 + 2];
double s = sqrt(x*x + y*y + z*z);
vertices[i * 3 + 0] /= s;
vertices[i * 3 + 1] /= s;
vertices[i * 3 + 2] /= s;
}
int num_faces = v->number_of_faces();
#ifdef DEBUG
printf("number of voronoi faces: %d\n", num_faces);
#endif
//std::vector<double> solids(face_areas.size()+1);
size_t c = 0;
for (int current_face=0;current_face<num_faces;current_face++)
{
int num = face_vertices[c++];
int point_index = nbr_indices[current_face];
if (point_index > 0)
{
double solid_angle = 0;
int u = face_vertices[c];
int v = face_vertices[c+1];
for (int i=2;i<num;i++)
{
int w = face_vertices[c+i];
double omega = calculate_solid_angle(&vertices[u*3], &vertices[v*3], &vertices[w*3]);
solid_angle += omega;
v = w;
}
face_areas[current_face] = solid_angle;
//solids[current_face] = solid_angle;
//face_areas[point_index] = solid_angle;
}
c += num;
}
#ifdef DEBUG
printf("\n");
for (int i=0;i<solids.size();i++)
{
printf("%d\t%f\t%f\n", i, solids[i], face_areas[i]);
}
#endif
assert(c == face_vertices.size());
return 0;
}
}
int calculate_neighbour_ordering(void* _voronoi_handle, int num_points, const double (*_points)[3], int8_t* ordering)
static int _calculate_neighbour_ordering(void* _voronoi_handle, int num_points, double (*_points)[3], bool calc_solid_angles, sorthelper_t* data)
{
assert(num_points <= PTM_MAX_INPUT_POINTS);
@ -59,7 +184,8 @@ int calculate_neighbour_ordering(void* _voronoi_handle, int num_points, const do
double max_norm = 0;
double points[PTM_MAX_INPUT_POINTS][3];
double normsq[PTM_MAX_INPUT_POINTS];
for (int i = 0;i<num_points;i++)
for (int i=0;i<num_points;i++)
{
double x = _points[i][0] - _points[0][0];
double y = _points[i][1] - _points[0][1];
@ -70,16 +196,13 @@ int calculate_neighbour_ordering(void* _voronoi_handle, int num_points, const do
normsq[i] = x*x + y*y + z*z;
max_norm = std::max(max_norm, normsq[i]);
#ifdef DEBUG
printf("point %d: %f\t%f\t%f\t%f\n", i, x, y, z, x*x + y*y + z*z);
#endif
}
max_norm = sqrt(max_norm);
std::vector<int> nbr_indices(num_points + 6);
std::vector<double> face_areas(num_points + 6);
int ret = calculate_voronoi_face_areas(num_points, points, normsq, max_norm, voronoi_handle, nbr_indices, face_areas);
int ret = calculate_voronoi_face_areas(num_points, points, normsq, max_norm, voronoi_handle, calc_solid_angles, nbr_indices, face_areas);
if (ret != 0)
return ret;
@ -93,24 +216,61 @@ int calculate_neighbour_ordering(void* _voronoi_handle, int num_points, const do
areas[index] = face_areas[i];
}
sorthelper_t data[PTM_MAX_INPUT_POINTS];
for (int i=0;i<num_points;i++)
{
assert(areas[i] == areas[i]);
data[i].area = areas[i];
data[i].dist = normsq[i];
data[i].index = i;
memcpy(data[i].offset, points[i], 3 * sizeof(double));
}
std::sort(data, data + num_points, &sorthelper_compare);
return ret;
}
#ifdef DEBUG
for (int i=0;i<num_points;i++)
printf("%d %f\n", data[i].index, data[i].area);
#endif
int calculate_neighbour_ordering( void* _voronoi_handle, size_t atom_index, int min_points, int (get_neighbours)(void* vdata, size_t central_index, size_t atom_index, int num, size_t* nbr_indices, int32_t* numbers,
double (*nbr_pos)[3]), void* nbrlist, bool calc_solid_angles,
size_t* indices, double (*points)[3], int32_t* numbers)
{
size_t nbr_indices[PTM_MAX_INPUT_POINTS];
int32_t nbr_numbers[PTM_MAX_INPUT_POINTS];
double nbr_pos[PTM_MAX_INPUT_POINTS][3];
int num_points = get_neighbours(nbrlist, atom_index, atom_index, PTM_MAX_INPUT_POINTS, nbr_indices, nbr_numbers, nbr_pos);
if (num_points < min_points)
return -1;
sorthelper_t data[PTM_MAX_INPUT_POINTS];
int ret = _calculate_neighbour_ordering(_voronoi_handle, num_points, nbr_pos, calc_solid_angles, data);
if (ret != 0)
return ret;
for (int i=0;i<num_points;i++)
ordering[i] = data[i].index;
{
int index = data[i].index;
indices[i] = nbr_indices[index];
numbers[i] = nbr_numbers[index];
memcpy(&points[i], nbr_pos[index], 3 * sizeof(double));
}
return num_points;
}
static int find_diamond_neighbours(void* _voronoi_handle, int num_points, double (*_points)[3], size_t* nbr_indices, int32_t* nbr_numbers, int num_solid_nbrs, bool calc_solid_angles, solidnbr_t* nbrlist)
{
sorthelper_t data[PTM_MAX_INPUT_POINTS];
int ret = _calculate_neighbour_ordering(_voronoi_handle, num_points, _points, calc_solid_angles, data);
if (ret != 0)
return ret;
int n = std::min(num_solid_nbrs, num_points - 1);
for (int i=0;i<n;i++)
{
nbrlist[i].index = nbr_indices[data[i+1].index];
nbrlist[i].area = data[i+1].area;
nbrlist[i].number = nbr_numbers[data[i+1].index];
memcpy(nbrlist[i].offset, data[i+1].offset, 3 * sizeof(double));
}
return ret;
}
@ -127,78 +287,112 @@ void voronoi_uninitialize_local(void* _ptr)
delete ptr;
}
#define MAX_INNER 4
#define MAX_SNBRS 12
typedef struct
int calculate_two_shell_neighbour_ordering( void* _voronoi_handle, size_t atom_index, int (get_neighbours)(void* vdata, size_t central_index, size_t atom_index, int num, size_t* nbr_indices, int32_t* numbers, double (*nbr_pos)[3]), void* nbrlist,
int num_inner, int num_outer, int max_snbrs, bool calc_solid_angles,
size_t* nbr_indices, double (*points)[3], int32_t* numbers)
{
double dist;
int p;
int index;
} diamond_t;
assert(num_inner <= MAX_INNER);
static bool diamond_compare(diamond_t const& a, diamond_t const& b)
{
return a.dist < b.dist;
}
size_t central_nbr_indices[MAX_SNBRS + 1];
int32_t central_nbr_numbers[MAX_SNBRS + 1];
double central_nbr_pos[MAX_SNBRS + 1][3];
int num_points = get_neighbours(nbrlist, atom_index, atom_index, max_snbrs + 1, central_nbr_indices, central_nbr_numbers, central_nbr_pos);
if (num_points < num_inner + 1)
return -1;
int calculate_diamond_neighbour_ordering( int num_points, double (*unpermuted_points)[3], int32_t* unpermuted_numbers,
int8_t* ordering, double (*points)[3], int32_t* numbers)
{
assert(num_points <= PTM_MAX_INPUT_POINTS);
solidnbr_t central_solid[MAX_SNBRS];
int ret = find_diamond_neighbours(_voronoi_handle, num_points, central_nbr_pos, central_nbr_indices, central_nbr_numbers, max_snbrs, calc_solid_angles, central_solid);
if (ret != 0)
return ret;
diamond_t data[4 * (PTM_MAX_INPUT_POINTS - 5)];
int index = 0;
for (int i=5;i<num_points;i++)
sorthelper_t data[MAX_INNER * 6];
nbr_indices[0] = atom_index;
numbers[0] = central_nbr_numbers[0];
memset(&points[0], 0, 3 * sizeof(double));
std::set<size_t> claimed;
claimed.insert(atom_index);
for (int i=0;i<num_inner;i++)
{
for (int j=1;j<5;j++)
nbr_indices[i+1] = central_solid[i].index;
numbers[i+1] = central_solid[i].number;
memcpy(&points[i+1], central_solid[i].offset, 3 * sizeof(double));
claimed.insert(central_solid[i].index);
}
int index = 0;
for (int i=0;i<num_inner;i++)
{
size_t inner_index = central_solid[i].index;
//----------------------------------------------
size_t inner_nbr_indices[MAX_SNBRS + 1];
int32_t inner_nbr_numbers[MAX_SNBRS + 1];
double inner_nbr_pos[MAX_SNBRS + 1][3];
int num_points = get_neighbours(nbrlist, atom_index, inner_index, max_snbrs + 1, inner_nbr_indices, inner_nbr_numbers, inner_nbr_pos);
if (num_points < num_inner + 1)
return -1;
solidnbr_t inner_solid[MAX_SNBRS];
ret = find_diamond_neighbours(_voronoi_handle, num_points, inner_nbr_pos, inner_nbr_indices, inner_nbr_numbers, max_snbrs, calc_solid_angles, inner_solid);
if (ret != 0)
return ret;
//----------------------------------------------
int n = std::min(6, num_points);
for (int j=0;j<n;j++)
{
double dx = unpermuted_points[i][0] - unpermuted_points[j][0];
double dy = unpermuted_points[i][1] - unpermuted_points[j][1];
double dz = unpermuted_points[i][2] - unpermuted_points[j][2];
bool already_claimed = claimed.find(inner_solid[j].index) != claimed.end();
if (already_claimed)
continue;
double d = dx*dx + dy*dy + dz*dz;
data[index].inner = i;
data[index].index = inner_solid[j].index;
data[index].area = inner_solid[j].area;
data[index].dist = 0;
data[index].number = inner_solid[j].number;
memcpy(data[index].offset, inner_solid[j].offset, 3 * sizeof(double));
for (int k=0;k<3;k++)
data[index].offset[k] += central_solid[i].offset[k];
data[index].p = j - 1;
data[index].index = i;
data[index].dist = d;
index++;
}
}
int n = index;
std::sort(data, data + n, &diamond_compare);
for (index=0;index<5;index++)
ordering[index] = index;
std::sort(data, data + n, &sorthelper_compare);
int num_found = 0;
bool hit[PTM_MAX_INPUT_POINTS] = {0};
int counts[4] = {0};
int counts[MAX_INNER] = {0};
for (int i=0;i<n;i++)
{
int p = data[i].p;
int q = data[i].index;
if (hit[q] || counts[p] >= 3)
int inner = data[i].inner;
int nbr_atom_index = data[i].index;
bool already_claimed = claimed.find(nbr_atom_index) != claimed.end();
if (counts[inner] >= num_outer || already_claimed)
continue;
ordering[1 + 4 + 3 * p + counts[p]] = q;
counts[p]++;
index++;
nbr_indices[1 + num_inner + num_outer * inner + counts[inner]] = nbr_atom_index;
numbers[1 + num_inner + num_outer * inner + counts[inner]] = data[i].number;
memcpy(points[1 + num_inner + num_outer * inner + counts[inner]], &data[i].offset, 3 * sizeof(double));
claimed.insert(nbr_atom_index);
counts[inner]++;
num_found++;
if (num_found >= 12)
if (num_found >= num_inner * num_outer)
break;
}
if (num_found != 12)
if (num_found != num_inner * num_outer)
return -1;
for (int i=0;i<PTM_NUM_NBRS_DCUB+1;i++)
{
memcpy(points[i], &unpermuted_points[ordering[i]], 3 * sizeof(double));
if (unpermuted_numbers != NULL)
numbers[i] = unpermuted_numbers[ordering[i]];
}
return 0;
}

View File

@ -1,14 +1,26 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#ifndef PTM_NEIGHBOUR_ORDERING_H
#define PTM_NEIGHBOUR_ORDERING_H
#include <inttypes.h>
#include <cstddef>
namespace ptm {
int calculate_neighbour_ordering(void* voronoi_handle, int num_points, const double (*_points)[3], int8_t* ordering);
int calculate_neighbour_ordering( void* _voronoi_handle, size_t atom_index, int min_points, int (get_neighbours)(void* vdata, size_t central_index, size_t atom_index, int num, size_t* nbr_indices, int32_t* numbers,
double (*nbr_pos)[3]), void* nbrlist, bool calc_solid_angles,
size_t* nbr_indices, double (*points)[3], int32_t* numbers);
int calculate_diamond_neighbour_ordering( int num_points, double (*unpermuted_points)[3], int32_t* unpermuted_numbers,
int8_t* ordering, double (*points)[3], int32_t* numbers);
int calculate_two_shell_neighbour_ordering( void* _voronoi_handle, size_t atom_index, int (get_neighbours)(void* vdata, size_t central_index, size_t atom_index, int num, size_t* nbr_indices, int32_t* numbers, double (*nbr_pos)[3]), void* nbrlist,
int num_inner, int num_outer, int max_snbrs, bool calc_solid_angles,
size_t* nbr_indices, double (*points)[3], int32_t* numbers);
void* voronoi_initialize_local();
void voronoi_uninitialize_local(void* ptr);

View File

@ -1,3 +1,12 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#include <cmath>
namespace ptm {

View File

@ -1,3 +1,12 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#ifndef PTM_NORMALIZE_VERTICES_H
#define PTM_NORMALIZE_VERTICES_H

View File

@ -1,5 +1,5 @@
/*******************************************************************************
* -/_|:|_|_\-
* -/_|:|_|_\-
*
* This code is a modification of D.L. Theobald's QCP rotation code.
* It has been adapted to calculate the polar decomposition of a 3x3 matrix
@ -14,7 +14,7 @@
* USA
*
* dtheobald@brandeis.edu
*
*
* Pu Liu
* Johnson & Johnson Pharmaceutical Research and Development, L.L.C.
* 665 Stockton Drive
@ -22,7 +22,7 @@
* USA
*
* pliu24@its.jnj.com
*
*
*
* If you use this QCP rotation calculation method in a publication, please
* reference:
@ -33,7 +33,7 @@
* Acta Crystallographica A 61(4):478-480.
*
* Pu Liu, Dmitris K. Agrafiotis, and Douglas L. Theobald (2009)
* "Fast determination of the optimal rotational matrix for macromolecular
* "Fast determination of the optimal rotational matrix for macromolecular
* superpositions."
* Journal of Computational Chemistry 31(7):1561-1563.
*
@ -63,7 +63,7 @@
* DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
* THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
* (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
* OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
* OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*
* Source: started anew.
*
@ -82,14 +82,13 @@
*
* 2016/05/29 QCP method adapted for polar decomposition of a 3x3 matrix,
* for use in Polyhedral Template Matching.
*
*
******************************************************************************/
#include <cmath>
#include <algorithm>
#include <cstring>
#include <string.h>
#include "ptm_quat.h"
#include "ptm_polar.h"
namespace ptm {
@ -327,7 +326,7 @@ void InnerProduct(double *A, int num, const double (*coords1)[3], double (*coord
A[6] += z1 * x2;
A[7] += z1 * y2;
A[8] += z1 * z2;
A[8] += z1 * z2;
}
}

View File

@ -1,9 +1,17 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#ifndef PTM_POLAR_H
#define PTM_POLAR_H
#include <stdint.h>
#include <stdbool.h>
#include <inttypes.h>
namespace ptm {

View File

@ -1,138 +1,159 @@
#include <cstring>
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#include <string.h>
#include <cmath>
#include <cfloat>
#include <algorithm>
namespace ptm {
#define SIGN(x) (x >= 0 ? 1 : -1)
#define MIN(X, Y) (((X) < (Y)) ? (X) : (Y))
#define MAX(X, Y) (((X) > (Y)) ? (X) : (Y))
#define SQRT_2 1.4142135623730951454746218587388284504414
#define HALF_SQRT_2 0.7071067811865474617150084668537601828575
const double generator_cubic[24][4] = {
{ 1, 0, 0, 0 },
{ sqrt(2)/2, sqrt(2)/2, 0, 0 },
{ sqrt(2)/2, 0, sqrt(2)/2, 0 },
{ sqrt(2)/2, 0, 0, sqrt(2)/2 },
{ sqrt(2)/2, 0, 0, -sqrt(2)/2 },
{ sqrt(2)/2, 0, -sqrt(2)/2, 0 },
{ sqrt(2)/2, -sqrt(2)/2, -0, -0 },
{ 0.5, 0.5, 0.5, 0.5 },
{ 0.5, 0.5, 0.5, -0.5 },
{ 0.5, 0.5, -0.5, 0.5 },
{ 0.5, 0.5, -0.5, -0.5 },
{ 0.5, -0.5, 0.5, 0.5 },
{ 0.5, -0.5, 0.5, -0.5 },
{ 0.5, -0.5, -0.5, 0.5 },
{ 0.5, -0.5, -0.5, -0.5 },
{ 0, 1, 0, 0 },
{ 0, sqrt(2)/2, sqrt(2)/2, 0 },
{ 0, sqrt(2)/2, 0, sqrt(2)/2 },
{ 0, sqrt(2)/2, 0, -sqrt(2)/2 },
{ 0, sqrt(2)/2, -sqrt(2)/2, 0 },
{ 0, 0, 1, 0 },
{ 0, 0, sqrt(2)/2, sqrt(2)/2 },
{ 0, 0, sqrt(2)/2, -sqrt(2)/2 },
{ 0, 0, 0, 1 },
};
#define PHI 1.6180339887498949025257388711906969547272
#define HALF_PHI 0.8090169943749474512628694355953484773636
const double generator_diamond_cubic[12][4] = {
{ 1, 0, 0, 0 },
{ 0.5, 0.5, 0.5, 0.5 },
{ 0.5, 0.5, 0.5, -0.5 },
{ 0.5, 0.5, -0.5, 0.5 },
{ 0.5, 0.5, -0.5, -0.5 },
{ 0.5, -0.5, 0.5, 0.5 },
{ 0.5, -0.5, 0.5, -0.5 },
{ 0.5, -0.5, -0.5, 0.5 },
{ 0.5, -0.5, -0.5, -0.5 },
{ 0, 1, 0, 0 },
{ 0, 0, 1, 0 },
{ 0, 0, 0, 1 },
};
#define INV_PHI 0.6180339887498947915034364086750429123640
#define HALF_INV_PHI 0.3090169943749473957517182043375214561820
#define SQRT_5_ 2.23606797749978969640917366873127623544061835961152572427089
#define SQRT_2_3 0.8164965809277260344600790631375275552273
#define SQRT_1_6 0.4082482904638630172300395315687637776136
const double generator_hcp[6][4] = {
{ 1, 0, 0, 0 },
{ 0.5, 0, 0, sqrt(3)/2 },
{ 0.5, 0, 0, -sqrt(3)/2 },
{ 0, sqrt(3)/2, 0.5, 0 },
{ 0, sqrt(3)/2, -0.5, 0 },
{ 0, 0, 1, 0 },
};
double generator_cubic[24][4] = { {1, 0, 0, 0 },
{0, 1, 0, 0 },
{0, 0, 1, 0 },
{0, 0, 0, 1 },
{0.5, 0.5, 0.5, 0.5 },
{0.5, 0.5, -0.5, 0.5 },
{0.5, -0.5, 0.5, 0.5 },
{0.5, -0.5, -0.5, 0.5 },
{-0.5, 0.5, 0.5, 0.5 },
{-0.5, 0.5, -0.5, 0.5 },
{-0.5, -0.5, 0.5, 0.5 },
{-0.5, -0.5, -0.5, 0.5 },
{HALF_SQRT_2, HALF_SQRT_2, 0, 0 },
{HALF_SQRT_2, 0, HALF_SQRT_2, 0 },
{HALF_SQRT_2, 0, 0, HALF_SQRT_2 },
{-HALF_SQRT_2, HALF_SQRT_2, 0, 0 },
{-HALF_SQRT_2, 0, HALF_SQRT_2, 0 },
{-HALF_SQRT_2, 0, 0, HALF_SQRT_2 },
{0, HALF_SQRT_2, HALF_SQRT_2, 0 },
{0, HALF_SQRT_2, 0, HALF_SQRT_2 },
{0, 0, HALF_SQRT_2, HALF_SQRT_2 },
{0, -HALF_SQRT_2, HALF_SQRT_2, 0 },
{0, -HALF_SQRT_2, 0, HALF_SQRT_2 },
{0, 0, -HALF_SQRT_2, HALF_SQRT_2 } };
const double generator_hcp_conventional[12][4] = {
{ 1, 0, 0, 0 },
{ sqrt(3)/2, 0, 0, 0.5 },
{ sqrt(3)/2, 0, 0, -0.5 },
{ 0.5, 0, 0, sqrt(3)/2 },
{ 0.5, 0, 0, -sqrt(3)/2 },
{ 0, 1, 0, 0 },
{ 0, sqrt(3)/2, 0.5, 0 },
{ 0, sqrt(3)/2, -0.5, 0 },
{ 0, 0.5, sqrt(3)/2, 0 },
{ 0, 0.5, -sqrt(3)/2, 0 },
{ 0, 0, 1, 0 },
{ 0, 0, 0, 1 },
};
double generator_diamond_cubic[12][4] = { {1, 0, 0, 0 },
{0, 1, 0, 0 },
{0, 0, 1, 0 },
{0, 0, 0, 1 },
{0.5, 0.5, 0.5, 0.5 },
{0.5, 0.5, -0.5, 0.5 },
{0.5, -0.5, 0.5, 0.5 },
{0.5, -0.5, -0.5, 0.5 },
{-0.5, 0.5, 0.5, 0.5 },
{-0.5, 0.5, -0.5, 0.5 },
{-0.5, -0.5, 0.5, 0.5 },
{-0.5, -0.5, -0.5, 0.5 } };
const double generator_diamond_hexagonal[3][4] = {
{ 1, 0, 0, 0 },
{ 0.5, 0, 0, sqrt(3)/2 },
{ 0.5, 0, 0, -sqrt(3)/2 },
};
double generator_hcp[6][4] = { {1, 0, 0, 0},
{0.5, 0.5, 0.5, 0.5},
{0.5, -0.5, -0.5, -0.5},
{0, SQRT_2_3, -SQRT_1_6, -SQRT_1_6},
{0, SQRT_1_6, -SQRT_2_3, SQRT_1_6},
{0, SQRT_1_6, SQRT_1_6, -SQRT_2_3} };
const double generator_icosahedral[60][4] = {
{ 1, 0, 0, 0 },
{ (1+sqrt(5))/4, 0.5, sqrt(25-10*sqrt(5))/10, sqrt(50-10*sqrt(5))/20 },
{ (1+sqrt(5))/4, 0.5, -sqrt(25-10*sqrt(5))/10, -sqrt(50-10*sqrt(5))/20 },
{ (1+sqrt(5))/4, 1/(1+sqrt(5)), sqrt(10*sqrt(5)+50)/20, -sqrt(50-10*sqrt(5))/20 },
{ (1+sqrt(5))/4, 1/(1+sqrt(5)), -sqrt(10*sqrt(5)+50)/20, sqrt(50-10*sqrt(5))/20 },
{ (1+sqrt(5))/4, 0, sqrt(50-10*sqrt(5))/10, sqrt(50-10*sqrt(5))/20 },
{ (1+sqrt(5))/4, 0, 0, sqrt(5./8-sqrt(5)/8) },
{ (1+sqrt(5))/4, 0, 0, -sqrt(5./8-sqrt(5)/8) },
{ (1+sqrt(5))/4, 0, -sqrt(50-10*sqrt(5))/10, -sqrt(50-10*sqrt(5))/20 },
{ (1+sqrt(5))/4, -1/(1+sqrt(5)), sqrt(10*sqrt(5)+50)/20, -sqrt(50-10*sqrt(5))/20 },
{ (1+sqrt(5))/4, -1/(1+sqrt(5)), -sqrt(10*sqrt(5)+50)/20, sqrt(50-10*sqrt(5))/20 },
{ (1+sqrt(5))/4, -0.5, sqrt(25-10*sqrt(5))/10, sqrt(50-10*sqrt(5))/20 },
{ (1+sqrt(5))/4, -0.5, -sqrt(25-10*sqrt(5))/10, -sqrt(50-10*sqrt(5))/20 },
{ 0.5, (1+sqrt(5))/4, sqrt(50-10*sqrt(5))/20, -sqrt(25-10*sqrt(5))/10 },
{ 0.5, (1+sqrt(5))/4, -sqrt(50-10*sqrt(5))/20, sqrt(25-10*sqrt(5))/10 },
{ 0.5, 0.5, sqrt((5+2*sqrt(5))/20), sqrt(25-10*sqrt(5))/10 },
{ 0.5, 0.5, sqrt(25-10*sqrt(5))/10, -sqrt((5+2*sqrt(5))/20) },
{ 0.5, 0.5, -sqrt(25-10*sqrt(5))/10, sqrt((5+2*sqrt(5))/20) },
{ 0.5, 0.5, -sqrt((5+2*sqrt(5))/20), -sqrt(25-10*sqrt(5))/10 },
{ 0.5, 1/(1+sqrt(5)), sqrt(10*sqrt(5)+50)/20, sqrt((5+2*sqrt(5))/20) },
{ 0.5, 1/(1+sqrt(5)), -sqrt(10*sqrt(5)+50)/20, -sqrt((5+2*sqrt(5))/20) },
{ 0.5, 0, sqrt((5+sqrt(5))/10), -sqrt(25-10*sqrt(5))/10 },
{ 0.5, 0, sqrt(50-10*sqrt(5))/10, -sqrt((5+2*sqrt(5))/20) },
{ 0.5, 0, -sqrt(50-10*sqrt(5))/10, sqrt((5+2*sqrt(5))/20) },
{ 0.5, 0, -sqrt((5+sqrt(5))/10), sqrt(25-10*sqrt(5))/10 },
{ 0.5, -1/(1+sqrt(5)), sqrt(10*sqrt(5)+50)/20, sqrt((5+2*sqrt(5))/20) },
{ 0.5, -1/(1+sqrt(5)), -sqrt(10*sqrt(5)+50)/20, -sqrt((5+2*sqrt(5))/20) },
{ 0.5, -0.5, sqrt((5+2*sqrt(5))/20), sqrt(25-10*sqrt(5))/10 },
{ 0.5, -0.5, sqrt(25-10*sqrt(5))/10, -sqrt((5+2*sqrt(5))/20) },
{ 0.5, -0.5, -sqrt(25-10*sqrt(5))/10, sqrt((5+2*sqrt(5))/20) },
{ 0.5, -0.5, -sqrt((5+2*sqrt(5))/20), -sqrt(25-10*sqrt(5))/10 },
{ 0.5, -(1+sqrt(5))/4, sqrt(50-10*sqrt(5))/20, -sqrt(25-10*sqrt(5))/10 },
{ 0.5, -(1+sqrt(5))/4, -sqrt(50-10*sqrt(5))/20, sqrt(25-10*sqrt(5))/10 },
{ 1/(1+sqrt(5)), (1+sqrt(5))/4, sqrt(50-10*sqrt(5))/20, sqrt(10*sqrt(5)+50)/20 },
{ 1/(1+sqrt(5)), (1+sqrt(5))/4, -sqrt(50-10*sqrt(5))/20, -sqrt(10*sqrt(5)+50)/20 },
{ 1/(1+sqrt(5)), 0.5, sqrt((5+2*sqrt(5))/20), -sqrt(10*sqrt(5)+50)/20 },
{ 1/(1+sqrt(5)), 0.5, -sqrt((5+2*sqrt(5))/20), sqrt(10*sqrt(5)+50)/20 },
{ 1/(1+sqrt(5)), 0, sqrt((5+sqrt(5))/10), sqrt(10*sqrt(5)+50)/20 },
{ 1/(1+sqrt(5)), 0, 0, sqrt(1-1/(2*sqrt(5)+6)) },
{ 1/(1+sqrt(5)), 0, 0, -sqrt(1-1/(2*sqrt(5)+6)) },
{ 1/(1+sqrt(5)), 0, -sqrt((5+sqrt(5))/10), -sqrt(10*sqrt(5)+50)/20 },
{ 1/(1+sqrt(5)), -0.5, sqrt((5+2*sqrt(5))/20), -sqrt(10*sqrt(5)+50)/20 },
{ 1/(1+sqrt(5)), -0.5, -sqrt((5+2*sqrt(5))/20), sqrt(10*sqrt(5)+50)/20 },
{ 1/(1+sqrt(5)), -(1+sqrt(5))/4, sqrt(50-10*sqrt(5))/20, sqrt(10*sqrt(5)+50)/20 },
{ 1/(1+sqrt(5)), -(1+sqrt(5))/4, -sqrt(50-10*sqrt(5))/20, -sqrt(10*sqrt(5)+50)/20 },
{ 0, 1, 0, 0 },
{ 0, (1+sqrt(5))/4, sqrt(5./8-sqrt(5)/8), 0 },
{ 0, (1+sqrt(5))/4, sqrt(50-10*sqrt(5))/20, -sqrt(50-10*sqrt(5))/10 },
{ 0, (1+sqrt(5))/4, -sqrt(50-10*sqrt(5))/20, sqrt(50-10*sqrt(5))/10 },
{ 0, (1+sqrt(5))/4, -sqrt(5./8-sqrt(5)/8), 0 },
{ 0, 0.5, sqrt((5+2*sqrt(5))/20), sqrt(50-10*sqrt(5))/10 },
{ 0, 0.5, sqrt(25-10*sqrt(5))/10, sqrt((5+sqrt(5))/10) },
{ 0, 0.5, -sqrt(25-10*sqrt(5))/10, -sqrt((5+sqrt(5))/10) },
{ 0, 0.5, -sqrt((5+2*sqrt(5))/20), -sqrt(50-10*sqrt(5))/10 },
{ 0, 1/(1+sqrt(5)), sqrt(1-1/(2*sqrt(5)+6)), 0 },
{ 0, 1/(1+sqrt(5)), sqrt(10*sqrt(5)+50)/20, -sqrt((5+sqrt(5))/10) },
{ 0, 1/(1+sqrt(5)), -sqrt(10*sqrt(5)+50)/20, sqrt((5+sqrt(5))/10) },
{ 0, 1/(1+sqrt(5)), -sqrt(1-1/(2*sqrt(5)+6)), 0 },
{ 0, 0, sqrt((5+sqrt(5))/10), -sqrt(50-10*sqrt(5))/10 },
{ 0, 0, sqrt(50-10*sqrt(5))/10, sqrt((5+sqrt(5))/10) },
};
double generator_diamond_hexagonal[3][4] = { {1, 0, 0, 0},
{0.5, 0.5, 0.5, 0.5},
{0.5, -0.5, -0.5, -0.5} };
double generator_icosahedral[60][4] = { {1, 0, 0, 0},
{HALF_PHI, -HALF_INV_PHI, -0.5, 0},
{HALF_PHI, 0, -HALF_INV_PHI, -0.5},
{HALF_PHI, -0.5, 0, -HALF_INV_PHI},
{HALF_PHI, HALF_INV_PHI, -0.5, 0},
{HALF_PHI, 0, HALF_INV_PHI, -0.5},
{HALF_PHI, -0.5, 0, HALF_INV_PHI},
{HALF_PHI, 0.5, 0, -HALF_INV_PHI},
{HALF_PHI, 0, -HALF_INV_PHI, 0.5},
{HALF_PHI, -HALF_INV_PHI, 0.5, 0},
{HALF_PHI, 0, HALF_INV_PHI, 0.5},
{HALF_PHI, HALF_INV_PHI, 0.5, 0},
{HALF_PHI, 0.5, 0, HALF_INV_PHI},
{0.5, HALF_PHI, -HALF_INV_PHI, 0},
{0.5, HALF_PHI, HALF_INV_PHI, 0},
{0.5, 0.5, 0.5, 0.5},
{0.5, 0.5, 0.5, -0.5},
{0.5, 0.5, -0.5, 0.5},
{0.5, 0.5, -0.5, -0.5},
{0.5, HALF_INV_PHI, 0, HALF_PHI},
{0.5, HALF_INV_PHI, 0, -HALF_PHI},
{0.5, 0, HALF_PHI, -HALF_INV_PHI},
{0.5, 0, HALF_PHI, HALF_INV_PHI},
{0.5, 0, -HALF_PHI, -HALF_INV_PHI},
{0.5, 0, -HALF_PHI, HALF_INV_PHI},
{0.5, -HALF_INV_PHI, 0, HALF_PHI},
{0.5, -HALF_INV_PHI, 0, -HALF_PHI},
{0.5, -0.5, 0.5, 0.5},
{0.5, -0.5, 0.5, -0.5},
{0.5, -0.5, -0.5, 0.5},
{0.5, -0.5, -0.5, -0.5},
{0.5, -HALF_PHI, -HALF_INV_PHI, 0},
{0.5, -HALF_PHI, HALF_INV_PHI, 0},
{HALF_INV_PHI, -HALF_PHI, 0, -0.5},
{HALF_INV_PHI, 0, -0.5, -HALF_PHI},
{HALF_INV_PHI, -0.5, -HALF_PHI, 0},
{HALF_INV_PHI, 0, 0.5, -HALF_PHI},
{HALF_INV_PHI, -HALF_PHI, 0, 0.5},
{HALF_INV_PHI, 0.5, -HALF_PHI, 0},
{HALF_INV_PHI, HALF_PHI, 0, -0.5},
{HALF_INV_PHI, -0.5, HALF_PHI, 0},
{HALF_INV_PHI, 0, -0.5, HALF_PHI},
{HALF_INV_PHI, HALF_PHI, 0, 0.5},
{HALF_INV_PHI, 0, 0.5, HALF_PHI},
{HALF_INV_PHI, 0.5, HALF_PHI, 0},
{0, 1, 0, 0},
{0, HALF_PHI, -0.5, HALF_INV_PHI},
{0, HALF_PHI, -0.5, -HALF_INV_PHI},
{0, HALF_PHI, 0.5, HALF_INV_PHI},
{0, HALF_PHI, 0.5, -HALF_INV_PHI},
{0, 0.5, HALF_INV_PHI, -HALF_PHI},
{0, 0.5, HALF_INV_PHI, HALF_PHI},
{0, 0.5, -HALF_INV_PHI, -HALF_PHI},
{0, 0.5, -HALF_INV_PHI, HALF_PHI},
{0, HALF_INV_PHI, -HALF_PHI, 0.5},
{0, HALF_INV_PHI, -HALF_PHI, -0.5},
{0, HALF_INV_PHI, HALF_PHI, 0.5},
{0, HALF_INV_PHI, HALF_PHI, -0.5},
{0, 0, 1, 0},
{0, 0, 0, 1} };
static void quat_rot(double* r, double* a, double* b)
{
@ -142,13 +163,13 @@ static void quat_rot(double* r, double* a, double* b)
b[3] = (r[0] * a[3] + r[1] * a[2] - r[2] * a[1] + r[3] * a[0]);
}
static int rotate_quaternion_into_fundamental_zone(int num_generators, double (*generator)[4], double* q)
static int rotate_quaternion_into_fundamental_zone(int num_generators, const double (*generator)[4], double* q)
{
double max = 0.0;
int i = 0, bi = -1;
for (i=0;i<num_generators;i++)
{
double* g = generator[i];
const double* g = generator[i];
double t = fabs(q[0] * g[0] - q[1] * g[1] - q[2] * g[2] - q[3] * g[3]);
if (t > max)
{
@ -158,7 +179,7 @@ static int rotate_quaternion_into_fundamental_zone(int num_generators, double (*
}
double f[4];
quat_rot(q, generator[bi], f);
quat_rot(q, (double*)generator[bi], f);
memcpy(q, &f, 4 * sizeof(double));
if (q[0] < 0)
{
@ -191,6 +212,11 @@ int rotate_quaternion_into_hcp_fundamental_zone(double* q)
return rotate_quaternion_into_fundamental_zone(6, generator_hcp, q);
}
int rotate_quaternion_into_hcp_conventional_fundamental_zone(double* q)
{
return rotate_quaternion_into_fundamental_zone(12, generator_hcp_conventional, q);
}
int rotate_quaternion_into_diamond_hexagonal_fundamental_zone(double* q)
{
return rotate_quaternion_into_fundamental_zone(3, generator_diamond_hexagonal, q);
@ -219,6 +245,26 @@ void normalize_quaternion(double* q)
q[3] /= size;
}
void quaternion_to_rotation_matrix(double* q, double* u)
{
double a = q[0];
double b = q[1];
double c = q[2];
double d = q[3];
u[0] = a*a + b*b - c*c - d*d;
u[1] = 2*b*c - 2*a*d;
u[2] = 2*b*d + 2*a*c;
u[3] = 2*b*c + 2*a*d;
u[4] = a*a - b*b + c*c - d*d;
u[5] = 2*c*d - 2*a*b;
u[6] = 2*b*d - 2*a*c;
u[7] = 2*c*d + 2*a*b;
u[8] = a*a - b*b - c*c + d*d;
}
void rotation_matrix_to_quaternion(double* u, double* q)
{
double r11 = u[0];
@ -236,14 +282,14 @@ void rotation_matrix_to_quaternion(double* u, double* q)
q[2] = (1.0 - r11 + r22 - r33) / 4.0;
q[3] = (1.0 - r11 - r22 + r33) / 4.0;
q[0] = sqrt(MAX(0, q[0]));
q[1] = sqrt(MAX(0, q[1]));
q[2] = sqrt(MAX(0, q[2]));
q[3] = sqrt(MAX(0, q[3]));
q[0] = sqrt(std::max(0., q[0]));
q[1] = sqrt(std::max(0., q[1]));
q[2] = sqrt(std::max(0., q[2]));
q[3] = sqrt(std::max(0., q[3]));
double m0 = MAX(q[0], q[1]);
double m1 = MAX(q[2], q[3]);
double max = MAX(m0, m1);
double m0 = std::max(q[0], q[1]);
double m1 = std::max(q[2], q[3]);
double max = std::max(m0, m1);
int i = 0;
for (i=0;i<4;i++)
@ -278,30 +324,10 @@ void rotation_matrix_to_quaternion(double* u, double* q)
normalize_quaternion(q);
}
void quaternion_to_rotation_matrix(double* q, double* u)
{
double a = q[0];
double b = q[1];
double c = q[2];
double d = q[3];
u[0] = a*a + b*b - c*c - d*d;
u[1] = 2*b*c - 2*a*d;
u[2] = 2*b*d + 2*a*c;
u[3] = 2*b*c + 2*a*d;
u[4] = a*a - b*b + c*c - d*d;
u[5] = 2*c*d - 2*a*b;
u[6] = 2*b*d - 2*a*c;
u[7] = 2*c*d + 2*a*b;
u[8] = a*a - b*b - c*c + d*d;
}
double quat_quick_misorientation(double* q1, double* q2)
{
double t = quat_dot(q1, q2);
t = MIN(1, MAX(-1, t));
t = std::min(1., std::max(-1., t));
return 2 * t * t - 1;
}
@ -310,91 +336,5 @@ double quat_misorientation(double* q1, double* q2)
return acos(quat_quick_misorientation(q1, q2));
}
double quat_quick_disorientation_cubic(double* q0, double* q1)
{
double qrot[4];
double qinv[4] = {q0[0], -q0[1], -q0[2], -q0[3]};
quat_rot(qinv, q1, qrot);
rotate_quaternion_into_cubic_fundamental_zone(qrot);
double t = qrot[0];
t = MIN(1, MAX(-1, t));
return 2 * t * t - 1;
}
double quat_disorientation_cubic(double* q0, double* q1)
{
return acos(quat_quick_disorientation_cubic(q0, q1));
}
double quat_quick_disorientation_diamond_cubic(double* q0, double* q1)
{
double qrot[4];
double qinv[4] = {q0[0], -q0[1], -q0[2], -q0[3]};
quat_rot(qinv, q1, qrot);
rotate_quaternion_into_diamond_cubic_fundamental_zone(qrot);
double t = qrot[0];
t = MIN(1, MAX(-1, t));
return 2 * t * t - 1;
}
double quat_disorientation_diamond_cubic(double* q0, double* q1)
{
return acos(quat_quick_disorientation_diamond_cubic(q0, q1));
}
double quat_quick_disorientation_hcp(double* q0, double* q1)
{
double qrot[4];
double qinv[4] = {q0[0], -q0[1], -q0[2], -q0[3]};
quat_rot(qinv, q1, qrot);
rotate_quaternion_into_hcp_fundamental_zone(qrot);
double t = qrot[0];
t = MIN(1, MAX(-1, t));
return 2 * t * t - 1;
}
double quat_disorientation_hcp(double* q0, double* q1)
{
return acos(quat_quick_disorientation_hcp(q0, q1));
}
double quat_quick_disorientation_diamond_hexagonal(double* q0, double* q1)
{
double qrot[4];
double qinv[4] = {q0[0], -q0[1], -q0[2], -q0[3]};
quat_rot(qinv, q1, qrot);
rotate_quaternion_into_diamond_hexagonal_fundamental_zone(qrot);
double t = qrot[0];
t = MIN(1, MAX(-1, t));
return 2 * t * t - 1;
}
double quat_disorientation_diamond_hexagonal(double* q0, double* q1)
{
return acos(quat_quick_disorientation_diamond_hexagonal(q0, q1));
}
double quat_quick_disorientation_icosahedral(double* q0, double* q1)
{
double qrot[4];
double qinv[4] = {q0[0], -q0[1], -q0[2], -q0[3]};
quat_rot(qinv, q1, qrot);
rotate_quaternion_into_icosahedral_fundamental_zone(qrot);
double t = qrot[0];
t = MIN(1, MAX(-1, t));
return 2 * t * t - 1;
}
double quat_disorientation_icosahedral(double* q0, double* q1)
{
return acos(quat_quick_disorientation_icosahedral(q0, q1));
}
}

View File

@ -1,3 +1,12 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#ifndef PTM_QUAT_H
#define PTM_QUAT_H
@ -7,26 +16,15 @@ int rotate_quaternion_into_cubic_fundamental_zone(double* q);
int rotate_quaternion_into_diamond_cubic_fundamental_zone(double* q);
int rotate_quaternion_into_icosahedral_fundamental_zone(double* q);
int rotate_quaternion_into_hcp_fundamental_zone(double* q);
int rotate_quaternion_into_hcp_conventional_fundamental_zone(double* q);
int rotate_quaternion_into_diamond_hexagonal_fundamental_zone(double* q);
void normalize_quaternion(double* q);
void quaternion_to_rotation_matrix(double* q, double* U);
void rotation_matrix_to_quaternion(double* u, double* q);
double quat_dot(double* a, double* b);
double quat_quick_misorientation(double* q1, double* q2);
double quat_misorientation(double* q1, double* q2);
double quat_quick_disorientation_cubic(double* q0, double* q1);
double quat_disorientation_cubic(double* q0, double* q1);
double quat_quick_disorientation_diamond_cubic(double* q0, double* q1);
double quat_disorientation_diamond_cubic(double* q0, double* q1);
double quat_quick_disorientation_hcp(double* q0, double* q1);
double quat_disorientation_hcp(double* q0, double* q1);
double quat_quick_disorientation_diamond_hexagonal(double* q0, double* q1);
double quat_disorientation_diamond_hexagonal(double* q0, double* q1);
double quat_quick_disorientation_icosahedral(double* q0, double* q1);
double quat_disorientation_icosahedral(double* q0, double* q1);
}
#endif

View File

@ -1,6 +1,15 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <string.h>
#include <cmath>
#include <cfloat>
#include <cassert>
@ -294,5 +303,74 @@ int match_dcub_dhex(double (*ch_points)[3], double (*points)[3], int32_t flags,
return PTM_NO_ERROR;
}
static void check_graphs_graphene( const refdata_t* s,
int num_points,
const double (*ideal_points)[3],
double (*normalized)[3],
int8_t* mapping,
result_t* res)
{
double G1 = 0, G2 = 0;
for (int i=0;i<num_points;i++)
{
double x1 = ideal_points[i][0];
double y1 = ideal_points[i][1];
double z1 = ideal_points[i][2];
double x2 = normalized[i][0];
double y2 = normalized[i][1];
double z2 = normalized[i][2];
G1 += x1 * x1 + y1 * y1 + z1 * z1;
G2 += x2 * x2 + y2 * y2 + z2 * z2;
}
double E0 = (G1 + G2) / 2;
double q[4], scale = 0;
double rmsd = calc_rmsd(num_points, ideal_points, normalized, mapping, G1, G2, E0, q, &scale);
if (rmsd < res->rmsd)
{
res->rmsd = rmsd;
res->scale = scale;
res->ref_struct = s;
memcpy(res->q, q, 4 * sizeof(double));
memcpy(res->mapping, mapping, sizeof(int8_t) * num_points);
}
}
int match_graphene(double (*points)[3], result_t* res)
{
int num_nbrs = structure_graphene.num_nbrs;
int num_points = num_nbrs + 1;
const double (*ideal_points)[3] = structure_graphene.points;
double normalized[PTM_MAX_POINTS][3];
subtract_barycentre(num_points, points, normalized);
int8_t mapping[PTM_MAX_POINTS];
for (int i=0;i<num_points;i++)
mapping[i] = i;
for (int i=0;i<2;i++)
{
std::swap(mapping[4], mapping[5]);
for (int j=0;j<2;j++)
{
std::swap(mapping[6], mapping[7]);
for (int k=0;k<2;k++)
{
std::swap(mapping[8], mapping[9]);
check_graphs_graphene(&structure_graphene, num_points, ideal_points, normalized, mapping, res);
}
}
}
return PTM_NO_ERROR;
}
}

View File

@ -1,3 +1,12 @@
/*Copyright (c) 2016 PM Larsen
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#ifndef PTM_STRUCTURE_MATCHER_H
#define PTM_STRUCTURE_MATCHER_H
@ -19,6 +28,7 @@ typedef struct
int match_general(const refdata_t* s, double (*ch_points)[3], double (*points)[3], convexhull_t* ch, result_t* res);
int match_fcc_hcp_ico(double (*ch_points)[3], double (*points)[3], int32_t flags, convexhull_t* ch, result_t* res);
int match_dcub_dhex(double (*ch_points)[3], double (*points)[3], int32_t flags, convexhull_t* ch, result_t* res);
int match_graphene(double (*points)[3], result_t* res);
}

View File

@ -1,3 +1,45 @@
/*
Voro++ Copyright (c) 2008, The Regents of the University of California, through
Lawrence Berkeley National Laboratory (subject to receipt of any required
approvals from the U.S. Dept. of Energy). All rights reserved.
Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are met:
(1) Redistributions of source code must retain the above copyright notice, this
list of conditions and the following disclaimer.
(2) Redistributions in binary form must reproduce the above copyright notice,
this list of conditions and the following disclaimer in the documentation
and/or other materials provided with the distribution.
(3) Neither the name of the University of California, Lawrence Berkeley
National Laboratory, U.S. Dept. of Energy nor the names of its contributors may
be used to endorse or promote products derived from this software without
specific prior written permission.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR
ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
(INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
You are under no obligation whatsoever to provide any bug fixes, patches, or
upgrades to the features, functionality or performance of the source code
("Enhancements") to anyone; however, if you choose to make your Enhancements
available either publicly, or directly to Lawrence Berkeley National
Laboratory, without imposing a separate written license agreement for such
Enhancements, then you hereby grant the following license: a non-exclusive,
royalty-free perpetual license to install, use, modify, prepare derivative
works, incorporate into other computer software, distribute, and sublicense
such enhancements or derivative works thereof, in binary and source code form.
*/
// Voro++, a 3D cell-based Voronoi library
//
// Author : Chris H. Rycroft (LBL / UC Berkeley)
@ -267,6 +309,7 @@ inline bool voronoicell_base::search_for_outside_edge(vc_class &vc,int &up) {
* \param[in,out] stackp2 a pointer to the end of the stack entries. */
template<class vc_class>
inline void voronoicell_base::add_to_stack(vc_class &vc,int lp,int *&stackp2) {
(void)vc;
for(int *k(ds2);k<stackp2;k++) if(*k==lp) return;
if(stackp2==stacke2) add_memory_ds2(stackp2);
*(stackp2++)=lp;
@ -1360,6 +1403,71 @@ void voronoicell_neighbor::neighbors(std::vector<int> &v) {
reset_edges();
}
/** Returns the number of faces of a computed Voronoi cell.
* \return The number of faces. */
int voronoicell_base::number_of_faces() {
int i,j,k,l,m,s=0;
for(i=1;i<p;i++) for(j=0;j<nu[i];j++) {
k=ed[i][j];
if(k>=0) {
s++;
ed[i][j]=-1-k;
l=cycle_up(ed[i][nu[i]+j],k);
do {
m=ed[k][l];
ed[k][l]=-1-m;
l=cycle_up(ed[k][nu[k]+l],m);
k=m;
} while (k!=i);
}
}
reset_edges();
return s;
}
/** Returns a vector of the vertex vectors in the global coordinate system.
* \param[out] v the vector to store the results in.
* \param[in] (x,y,z) the position vector of the particle in the global
* coordinate system. */
void voronoicell_base::vertices(double x,double y,double z,std::vector<double> &v) {
v.resize(3*p);
double *ptsp=pts;
for(int i=0;i<3*p;i+=3) {
v[i]=x+*(ptsp++)*0.5;
v[i+1]=y+*(ptsp++)*0.5;
v[i+2]=z+*(ptsp++)*0.5;
}
}
/** For each face, this routine outputs a bracketed sequence of numbers
* containing a list of all the vertices that make up that face.
* \param[out] v the vector to store the results in. */
void voronoicell_base::face_vertices(std::vector<int> &v) {
int i,j,k,l,m,vp(0),vn;
v.clear();
for(i=1;i<p;i++) for(j=0;j<nu[i];j++) {
k=ed[i][j];
if(k>=0) {
v.push_back(0);
v.push_back(i);
ed[i][j]=-1-k;
l=cycle_up(ed[i][nu[i]+j],k);
do {
v.push_back(k);
m=ed[k][l];
ed[k][l]=-1-m;
l=cycle_up(ed[k][nu[k]+l],m);
k=m;
} while (k!=i);
vn=v.size();
v[vp]=vn-vp-1;
vp=vn;
}
}
reset_edges();
}
// Explicit instantiation
template bool voronoicell_base::nplane(voronoicell_neighbor&,double,double,double,double,int);
template void voronoicell_base::check_memory_for_copy(voronoicell_neighbor&,voronoicell_base*);

View File

@ -1,3 +1,46 @@
/*
Voro++ Copyright (c) 2008, The Regents of the University of California, through
Lawrence Berkeley National Laboratory (subject to receipt of any required
approvals from the U.S. Dept. of Energy). All rights reserved.
Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are met:
(1) Redistributions of source code must retain the above copyright notice, this
list of conditions and the following disclaimer.
(2) Redistributions in binary form must reproduce the above copyright notice,
this list of conditions and the following disclaimer in the documentation
and/or other materials provided with the distribution.
(3) Neither the name of the University of California, Lawrence Berkeley
National Laboratory, U.S. Dept. of Energy nor the names of its contributors may
be used to endorse or promote products derived from this software without
specific prior written permission.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR
ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
(INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
You are under no obligation whatsoever to provide any bug fixes, patches, or
upgrades to the features, functionality or performance of the source code
("Enhancements") to anyone; however, if you choose to make your Enhancements
available either publicly, or directly to Lawrence Berkeley National
Laboratory, without imposing a separate written license agreement for such
Enhancements, then you hereby grant the following license: a non-exclusive,
royalty-free perpetual license to install, use, modify, prepare derivative
works, incorporate into other computer software, distribute, and sublicense
such enhancements or derivative works thereof, in binary and source code form.
*/
// Voro++, a 3D cell-based Voronoi library
//
// Author : Chris H. Rycroft (LBL / UC Berkeley)

View File

@ -1,3 +1,46 @@
/*
Voro++ Copyright (c) 2008, The Regents of the University of California, through
Lawrence Berkeley National Laboratory (subject to receipt of any required
approvals from the U.S. Dept. of Energy). All rights reserved.
Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are met:
(1) Redistributions of source code must retain the above copyright notice, this
list of conditions and the following disclaimer.
(2) Redistributions in binary form must reproduce the above copyright notice,
this list of conditions and the following disclaimer in the documentation
and/or other materials provided with the distribution.
(3) Neither the name of the University of California, Lawrence Berkeley
National Laboratory, U.S. Dept. of Energy nor the names of its contributors may
be used to endorse or promote products derived from this software without
specific prior written permission.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR
ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
(INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
You are under no obligation whatsoever to provide any bug fixes, patches, or
upgrades to the features, functionality or performance of the source code
("Enhancements") to anyone; however, if you choose to make your Enhancements
available either publicly, or directly to Lawrence Berkeley National
Laboratory, without imposing a separate written license agreement for such
Enhancements, then you hereby grant the following license: a non-exclusive,
royalty-free perpetual license to install, use, modify, prepare derivative
works, incorporate into other computer software, distribute, and sublicense
such enhancements or derivative works thereof, in binary and source code form.
*/
// Voro++, a 3D cell-based Voronoi library
//
// Author : Chris H. Rycroft (LBL / UC Berkeley)