mcodev31 / libmsym

molecular point group symmetry lib
MIT License
73 stars 33 forks source link

Symmetrization of a molecule to the desired symmetry with THRESHOLD. #5

Closed ghost closed 7 years ago

ghost commented 7 years ago

Dear libmsym developers.

I have a slightly asymmetric molecule in XYZ-file

16
Coordinates from ORCA-job cubane E -308.476504210567
H       2.066906     -1.290195      0.022380
C       1.078819     -0.807735     -0.052882
H       2.067165      1.289619      0.017269
C       1.078771      0.807357     -0.054339
H       0.000865     -1.501089     -2.071487
C       0.000492     -1.331808     -0.985550
H      -0.001572      1.502336     -2.071168
C      -0.000676      1.331920     -0.985471
H      -0.000515     -1.457760      1.968468
C      -0.000379     -0.778815      1.103109
H       0.001330      1.457400      1.968659
C       0.000580      0.778691      1.103119
H      -2.067067     -1.289825      0.017976
C      -1.078730     -0.807503     -0.054085
H      -2.067142      1.289825      0.020021
C      -1.078846      0.807581     -0.053185

And I want to symmetrize it to C2V point group as I can do it with NWCHEM http://www.nwchem-sw.org/index.php/Release66:Geometry#SYMMETRY_--_Symmetry_Group_Input

I try to run example function from msym_example.c with different thresholds parameter, especially I with different DEFAULT_PERMUTATION_THRESHOLD, but I receive on

Error Invalid equivalence set: Cannot find permutation for C2 when determining equivalence set from point group C2v

What should I do to get symmetrized molecule.

Best, Vladimir.

mcodev31 commented 7 years ago

Hi Vladimir

That typically happens when the permutation threshold is too low. Changing all thresholds to 1.0e-2 gives C2v:

16

H 2.067069977 -1.289866032 0.018250820 C 1.078791503 -0.807544001 -0.054783425 H 2.067069977 1.289866032 0.018250820 C 1.078791503 0.807544001 -0.054783425 H -0.000000000 -1.501712408 -2.072488160 C -0.000000000 -1.331863956 -0.986671172 H 0.000000000 1.501712408 -2.072488160 C 0.000000000 1.331863956 -0.986671172 H -0.000000000 -1.457580052 1.967402839 C -0.000000000 -0.778753029 1.101953328 H 0.000000000 1.457580052 1.967402839 C 0.000000000 0.778753029 1.101953328 H -2.067069977 -1.289866032 0.018250820 C -1.078791503 -0.807544001 -0.054783425 H -2.067069977 1.289866032 0.018250820 C -1.078791503 0.807544001 -0.054783425

ghost commented 7 years ago

Hi, Marcus. I have found in source code default threshold values. https://github.com/mcodev31/libmsym/blob/master/src/context.h#L19 Does msym_example.c use it?

mcodev31 commented 7 years ago

Yes msym_example uses a NULL threshold, which means default. Note that this is just example code, libmsym is meant as a library for other projects to integrate (e.g. Avogadro and Luscus have versions with symmetry from libmsym)

If you want to use your own threshold do this:

#define MY_ZERO_THRESHOLD 1.0e-2
#define MY_GEOMETRY_THRESHOLD 1.0e-2
#define MY_ANGLE_THRESHOLD 1.0e-2
#define MY_EQUIVALENCE_THRESHOLD 1.0e-2
#define MY_EIGFACT_THRESHOLD 1.0e-2
#define MY_PERMUTATION_THRESHOLD 1.0e-2
#define MY_ORTHOGONALIZATION_THRESHOLD 1.0e-2

const msym_thresholds_t thresholds = {
    .zero = MY_ZERO_THRESHOLD,
    .geometry = MY_GEOMETRY_THRESHOLD,
    .angle = MY_ANGLE_THRESHOLD,
    .equivalence = MY_EQUIVALENCE_THRESHOLD,
    .eigfact = MY_EIGFACT_THRESHOLD,
    .permutation = MY_PERMUTATION_THRESHOLD,
    .orthogonalization = MY_ORTHOGONALIZATION_THRESHOLD
};

int main(int argc, const char * argv[]) {
    int ret = 1;
    if(argc == 2){
        ret = example(argv[1],&thresholds);
        fflush(stdout);
    } else {
        printf("usage msym_example <xyz-file>");
    }
    return ret;
}
ghost commented 7 years ago

Oh, thanks a lot. I have no questions on it.