Skip to content

Commit

Permalink
threshold scaling and fix possible null pointer deref
Browse files Browse the repository at this point in the history
threshold scaling based on found vectors
fix possible null pointer dereference in context
minor fixes
  • Loading branch information
mcodev31 committed Nov 9, 2015
1 parent 4adbebf commit 29d7dba
Show file tree
Hide file tree
Showing 8 changed files with 75 additions and 59 deletions.
22 changes: 12 additions & 10 deletions examples/example.c
Original file line number Diff line number Diff line change
Expand Up @@ -70,13 +70,11 @@ int example(const char* in_file, msym_thresholds_t *thresholds){
int length = read_xyz(in_file, &elements);
if(length <= 0) return -1;

//orbitalsl = sizeof(orbitals)/sizeof(char*);
//bfsl = orbitalsl*length;
//bfs = calloc(bfsl, sizeof(msym_basis_function_t));
double (*psalcs)[bfsl] = NULL; //calloc(bfsl, sizeof(*salcs)); // SALCs in matrix form, and input for symmetrization
double *pcmem = NULL; // calloc(bfsl, sizeof(double)); // Some temporary memory
int *pspecies = NULL; // calloc(bfsl, sizeof(*species));
msym_partner_function_t *ppf = NULL; //calloc(bfsl, sizeof(*pf));

double (*psalcs)[bfsl] = NULL; // SALCs in matrix form, and input for symmetrization
double *pcmem = NULL; // Some temporary memory
int *pspecies = NULL;
msym_partner_function_t *ppf = NULL;

/* Create a context */
msym_context ctx = msymCreateContext();
Expand Down Expand Up @@ -306,9 +304,13 @@ int example(const char* in_file, msym_thresholds_t *thresholds){
do{

int sel_ss = 0, sel_salc = 0;
for(int i = 0; i < msrsl;i++) printf("\t [%d] %s (%d SALCs with %d partner functions each)\n",i, mct->s[msrs[i].s].name, msrs[i].salcl, mct->s[msrs[i].s].d);

do {printf("\nChoose subspace [0-%d]:",msrsl-1);} while(scanf(" %d", &sel_ss) <= 0 || sel_ss < 0 || sel_ss >= msrsl);
for(int i = 0; i < msrsl;i++){
if(msrs[i].salcl > 0)
printf("\t [%d] %s (%d SALCs with %d partner functions each)\n",i, mct->s[msrs[i].s].name, msrs[i].salcl, mct->s[msrs[i].s].d);
else
printf("\t [-] %s (no SALCs of this symmetry species)\n",mct->s[msrs[i].s].name);
}
do {printf("\nChoose subspace [0-%d]:",msrsl-1);} while(scanf(" %d", &sel_ss) <= 0 || sel_ss < 0 || sel_ss >= msrsl || msrs[sel_ss].salcl <= 0);

int salcl = msrs[sel_ss].salcl;

Expand Down
53 changes: 27 additions & 26 deletions src/context.c
Original file line number Diff line number Diff line change
Expand Up @@ -117,7 +117,7 @@ msym_error_t msymSetThresholds(msym_context ctx, const msym_thresholds_t *thresh

msym_error_t msymGetThresholds(msym_context ctx, const msym_thresholds_t **thresholds){
msym_error_t ret = MSYM_SUCCESS;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->thresholds == NULL) ret = MSYM_INVALID_THRESHOLD;
*thresholds = ctx->thresholds;
err:
Expand All @@ -126,7 +126,7 @@ msym_error_t msymGetThresholds(msym_context ctx, const msym_thresholds_t **thres

msym_error_t msymSetElements(msym_context ctx, int length, msym_element_t *elements){
msym_error_t ret = MSYM_SUCCESS;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
/* Allow manual setting of point group before elements */
if(NULL != ctx->es) ctxDestroyPointGroup(ctx);
return ctxSetElements(ctx, length, elements);
Expand All @@ -137,7 +137,7 @@ msym_error_t msymSetElements(msym_context ctx, int length, msym_element_t *eleme
msym_error_t msymGetElements(msym_context ctx, int *length, msym_element_t **elements){
msym_error_t ret = MSYM_SUCCESS;
msym_element_t *relements = NULL;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->elements == NULL || ctx->ext.elements == NULL) {ret = MSYM_INVALID_ELEMENTS;goto err;}

*elements = ctx->ext.elements;
Expand All @@ -164,7 +164,7 @@ msym_error_t msymGetEquivalenceSets(msym_context ctx, int *length, const msym_eq

msym_error_t msymGetBasisFunctions(msym_context ctx, int *length, msym_basis_function_t **basis){
msym_error_t ret = MSYM_SUCCESS;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->basis == NULL) {
msymSetErrorDetails("Found no basis functions");
ret = MSYM_INVALID_BASIS_FUNCTIONS;
Expand All @@ -179,7 +179,7 @@ msym_error_t msymGetBasisFunctions(msym_context ctx, int *length, msym_basis_fun

msym_error_t msymSetBasisFunctions(msym_context ctx, int length, msym_basis_function_t *basis){
msym_error_t ret = MSYM_SUCCESS;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->elements == NULL) {ret = MSYM_INVALID_ELEMENTS;goto err;}
ctxDestroyBasisFunctions(ctx);
ctx->basis = malloc(sizeof(msym_basis_function_t[length]));
Expand Down Expand Up @@ -229,7 +229,7 @@ msym_error_t msymSetBasisFunctions(msym_context ctx, int length, msym_basis_func

msym_error_t msymGetPointGroupName(msym_context ctx, int l, char buf[l]){
msym_error_t ret = MSYM_SUCCESS;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->pg == NULL) {ret = MSYM_INVALID_POINT_GROUP;goto err;}
snprintf(buf, l, "%s",ctx->pg->name);
err:
Expand All @@ -238,7 +238,7 @@ msym_error_t msymGetPointGroupName(msym_context ctx, int l, char buf[l]){

msym_error_t msymGetPointGroupType(msym_context ctx, msym_point_group_type_t *t, int *n){
msym_error_t ret = MSYM_SUCCESS;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->pg == NULL) {ret = MSYM_INVALID_POINT_GROUP;goto err;}

*t = ctx->pg->type;
Expand All @@ -254,7 +254,7 @@ msym_error_t msymGetSubgroups(msym_context ctx, int *sgl, const msym_subgroup_t
msym_subgroup_t *gsg = NULL;
int gsgl = 0;

if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->pg == NULL) {ret = MSYM_INVALID_POINT_GROUP;goto err;}
if(ctx->pg->perm == NULL && !(isLinearPointGroup(ctx->pg) && !isLinearSubgroup(ctx->pg))) {
ret = MSYM_INVALID_PERMUTATION;
Expand Down Expand Up @@ -337,7 +337,7 @@ msym_error_t msymGetCharacterTable(msym_context ctx, const msym_character_table_

msym_error_t msymGetCenterOfMass(msym_context ctx, double v[3]){
msym_error_t ret = MSYM_SUCCESS;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->elements == NULL) {ret = MSYM_INVALID_ELEMENTS;goto err;}
vcopy(ctx->cm, v);
err:
Expand All @@ -355,7 +355,7 @@ msym_error_t msymSetCenterOfMass(msym_context ctx, double cm[3]){

msym_error_t msymGetGeometry(msym_context ctx, msym_geometry_t *geometry){
msym_error_t ret = MSYM_SUCCESS;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->elements == NULL) {ret = MSYM_INVALID_ELEMENTS;goto err;}
if(ctx->geometry == MSYM_GEOMETRY_UNKNOWN) {ret = MSYM_INVALID_GEOMETRY;goto err;}
*geometry = ctx->geometry;
Expand All @@ -365,15 +365,15 @@ msym_error_t msymGetGeometry(msym_context ctx, msym_geometry_t *geometry){

msym_error_t msymGetPrincipalMoments(msym_context ctx, double eigval[3]){
msym_error_t ret = MSYM_SUCCESS;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->elements == NULL) {ret = MSYM_INVALID_ELEMENTS;goto err;}
vcopy(ctx->eigval, eigval);
err:
return ret;
}
msym_error_t msymGetPrincipalAxes(msym_context ctx, double eigvec[3][3]){
msym_error_t ret = MSYM_SUCCESS;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->elements == NULL) {ret = MSYM_INVALID_ELEMENTS;goto err;}
mcopy(ctx->eigvec, eigvec);
err:
Expand All @@ -382,9 +382,9 @@ msym_error_t msymGetPrincipalAxes(msym_context ctx, double eigvec[3][3]){

msym_error_t msymGetRadius(msym_context ctx, double *radius){
msym_error_t ret = MSYM_SUCCESS;
double r = 0.0;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->elements == NULL) {ret = MSYM_INVALID_ELEMENTS;goto err;}
double r = 0.0;
for(int i = 0;i < ctx->elementsl;i++){
double abs = vabs(ctx->elements[i].v);
r = r > abs ? r : abs;
Expand All @@ -398,7 +398,7 @@ msym_error_t msymGetRadius(msym_context ctx, double *radius){
msym_error_t msymGetSymmetryOperations(msym_context ctx, int *sopsl, const msym_symmetry_operation_t **sops){
msym_error_t ret = MSYM_SUCCESS;
msym_symmetry_operation_t *rsops = NULL;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->pg == NULL || ctx->pg->sops == NULL) {ret = MSYM_INVALID_POINT_GROUP;goto err;}

*sops = ctx->pg->sops;
Expand All @@ -413,12 +413,12 @@ msym_error_t msymGetSymmetryOperations(msym_context ctx, int *sopsl, const msym_

msym_error_t msymReleaseContext(msym_context ctx){
msym_error_t ret = MSYM_SUCCESS;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
free(ctx->thresholds);
ctxDestroyElements(ctx);
ctxDestroyPointGroup(ctx);
free(ctx);
err:
//err:
return ret;
}

Expand All @@ -433,7 +433,7 @@ msym_error_t msymReleaseContext(msym_context ctx){
msym_error_t ctxSetElements(msym_context ctx, int length, msym_element_t elements[length]){
msym_error_t ret = MSYM_SUCCESS;
msym_thresholds_t *thresholds = NULL;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}

ctxDestroyElements(ctx);

Expand Down Expand Up @@ -466,6 +466,7 @@ msym_error_t ctxSetElements(msym_context ctx, int length, msym_element_t element

return ret;
err:

free(ctx->elements);
free(ctx->pelements);
free(ctx->ext.elements);
Expand All @@ -480,8 +481,8 @@ msym_error_t ctxSetElements(msym_context ctx, int length, msym_element_t element

msym_error_t ctxGetThresholds(msym_context ctx, msym_thresholds_t **thresholds){
msym_error_t ret = MSYM_SUCCESS;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx->thresholds == NULL) ret = MSYM_INVALID_THRESHOLD;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->thresholds == NULL) {ret = MSYM_INVALID_THRESHOLD; goto err;}
msym_thresholds_t *t = ctx->thresholds;
if(t->angle < 1.0 && !signbit(t->angle) &&
t->equivalence < 1.0 && !signbit(t->equivalence) &&
Expand All @@ -501,7 +502,7 @@ msym_error_t ctxGetThresholds(msym_context ctx, msym_thresholds_t **thresholds){

msym_error_t ctxGetElements(msym_context ctx, int *l, msym_element_t **elements){
msym_error_t ret = MSYM_SUCCESS;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->elements == NULL) {ret = MSYM_INVALID_ELEMENTS; goto err;}
*elements = (msym_element_t *) ctx->elements;
*l = ctx->elementsl;
Expand All @@ -511,7 +512,7 @@ msym_error_t ctxGetElements(msym_context ctx, int *l, msym_element_t **elements)

msym_error_t ctxGetExternalElements(msym_context ctx, int *l, msym_element_t **elements){
msym_error_t ret = MSYM_SUCCESS;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->ext.elements == NULL) {ret = MSYM_INVALID_ELEMENTS; goto err;}
*elements = (msym_element_t *) ctx->ext.elements;
*l = ctx->elementsl;
Expand All @@ -537,7 +538,7 @@ msym_error_t ctxUpdateExternalElementCoordinates(msym_context ctx){

msym_error_t ctxGetInternalElement(msym_context ctx, msym_element_t *ext, msym_element_t **element){
msym_error_t ret = MSYM_SUCCESS;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->ext.elements == NULL) {ret = MSYM_INVALID_ELEMENTS;goto err;}
if(ext < ctx->ext.elements || ext >= ctx->ext.elements+ctx->elementsl){
msymSetErrorDetails("Element pointer (%p) outside memory block (%p -> %p)", ext, ctx->ext.elements, ctx->ext.elements + ctx->elementsl);
Expand All @@ -551,7 +552,7 @@ msym_error_t ctxGetInternalElement(msym_context ctx, msym_element_t *ext, msym_e

msym_error_t ctxGetSubgroups(msym_context ctx, int *sgl, msym_subgroup_t **sg){
msym_error_t ret = MSYM_SUCCESS;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->sg == NULL) {ret = MSYM_INVALID_SUBGROUPS;goto err;}
*sg = ctx->sg;
*sgl = ctx->sgl;
Expand All @@ -561,7 +562,7 @@ msym_error_t ctxGetSubgroups(msym_context ctx, int *sgl, msym_subgroup_t **sg){

msym_error_t ctxSetSubgroups(msym_context ctx, int sgl, msym_subgroup_t *sg){
msym_error_t ret = MSYM_SUCCESS;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
ctxDestroySubgroups(ctx);
ctx->sg = sg;
ctx->sgl = sgl;
Expand All @@ -572,7 +573,7 @@ msym_error_t ctxSetSubgroups(msym_context ctx, int sgl, msym_subgroup_t *sg){

msym_error_t ctxGetElementPtrs(msym_context ctx, int *l, msym_element_t ***pelements){
msym_error_t ret = MSYM_SUCCESS;
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT;goto err;}
if(ctx == NULL) {ret = MSYM_INVALID_CONTEXT; return ret;}
if(ctx->pelements == NULL) {ret = MSYM_INVALID_ELEMENTS; goto err;}
*pelements = (msym_element_t **) ctx->pelements;
*l = ctx->elementsl;
Expand Down
2 changes: 1 addition & 1 deletion src/context.h
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,7 @@
#define DEFAULT_EQUIVALENCE_THRESHOLD 5.0e-4
#define DEFAULT_EIGFACT_THRESHOLD 1.0e-3
#define DEFAULT_PERMUTATION_THRESHOLD 5.0e-3
#define DEFAULT_ORTHOGONALIZATION_THRESHOLD 0.01
#define DEFAULT_ORTHOGONALIZATION_THRESHOLD 1.0e-2


msym_error_t ctxGetThresholds(msym_context ctx, msym_thresholds_t **thresholds);
Expand Down
10 changes: 9 additions & 1 deletion src/equivalence_set.c
Original file line number Diff line number Diff line change
Expand Up @@ -48,6 +48,13 @@ msym_error_t generateEquivalenceSet(msym_point_group_t *pg, int length, msym_ele
msym_equivalence_set_t *ges = calloc(length,sizeof(msym_equivalence_set_t));
int gel = 0;
int gesl = 0;

if(pg->order <= 0){
ret = MSYM_INVALID_POINT_GROUP;
msymSetErrorDetails("Point group of zero order when determining equivalence set");
goto err;
}

for(int i = 0;i < length;i++){
double ev[3];
vsub(elements[i].v, cm, ev);
Expand Down Expand Up @@ -83,7 +90,7 @@ msym_error_t generateEquivalenceSet(msym_point_group_t *pg, int length, msym_ele
}
}

if(pg->order % aes->length != 0){
if(!aes->length || (pg->order % aes->length != 0)){
msymSetErrorDetails("Equivalence set length (%d) not a divisor of point group order (%d)",pg->order);
ret = MSYM_INVALID_EQUIVALENCE_SET;
goto err;
Expand Down Expand Up @@ -185,6 +192,7 @@ msym_error_t findPointGroupEquivalenceSets(msym_point_group_t *pg, int length, m
free(pelements);
return ret;
err:
free(ges);
free(pelements);
return ret;

Expand Down
3 changes: 3 additions & 0 deletions src/linalg.c
Original file line number Diff line number Diff line change
Expand Up @@ -535,6 +535,7 @@ int mgs(int l, const double m[l][l], double o[l][l], int n, double t){
int mgs2(int l, int lm, const double m[l][l], double o[l][l], int n, double t){

int nm = n + lm + MGS_ADD;
double ts = l/(1.0 + l);
for(int i = 0; i < l && n < nm;i++){
if(vlabs(l, m[i]) < t){
continue;
Expand All @@ -550,8 +551,10 @@ int mgs2(int l, int lm, const double m[l][l], double o[l][l], int n, double t){
for(int k = 0;k < l;k++) o[n][k] -= c*o[j][k];
}
if(vlabs(l, o[n]) >= t){
t *= ts;
vlnorm(l, o[n]);
n++;

}
}
}
Expand Down
2 changes: 1 addition & 1 deletion src/point_group.c
Original file line number Diff line number Diff line change
Expand Up @@ -793,8 +793,8 @@ msym_error_t transformAxes(msym_point_group_type_t type, int n, msym_symmetry_op
case (MSYM_POINT_GROUP_TYPE_Ih) : {
msym_symmetry_operation_t *dprimary = NULL;
msym_symmetry_operation_t *sop;
double z = -2.0;
for(sop = sops; sop < (sops + sopsl); sop++){
double z = -2.0;
if(sop->type == PROPER_ROTATION && sop->order == 2){
double v[3];
vnorm2(sop->v,v);
Expand Down
Loading

0 comments on commit 29d7dba

Please sign in to comment.