mcocdawc / chemcoord

A python module for manipulating cartesian and internal coordinates.
GNU Lesser General Public License v3.0
72 stars 19 forks source link

changing ONLY the position of the origin #57

Closed oferfi closed 2 years ago

oferfi commented 5 years ago

I'm working now with some linear organic chains. it would be very usefull to tell CC to start from .e.g the terminal carbon atom, and then move automatically along the carbon chain. I thoung of doing so by chaning the "Definition of the construction table" in a particular way: I want only to define the origin and let CC cunstruct the rest.

is it possible?

attached file: ht-cis.txt

Code Sample

#!/usr/bin/env python
import chemcoord as cc
import time
import pandas as pd
import sympy
import argparse
import numpy as np
parser = argparse.ArgumentParser()
parser.add_argument('-i', metavar = "input xyz file", dest='xyz', help='XYZ',type=str)
input = parser.parse_args()

xyz_coord = cc.Cartesian.read_xyz(input.xyz, start_index=1)
zcoord = xyz_coord.get_zmat()
new_lst=np.array([],dtype=int)

for i in range(1,len(xyz_coord)+1):
 new_lst=np.append(new_lst,i)
print(zcoord.change_numbering(new_lst))
print(zcoord.loc[:,'bond'])

Problem description

please note at the output details that the C2 and C3 refer to the same carbon which (in the case of HEXATRIEN) means that the origin is somewhere in the middle of the chain.

Expected Output

my goal isto get output somewhat like this: 1 C origin 1.001936 e_z 42.471043 e_x -0.0 2 C 1 1.448768 e_z 147.275520 e_x -0.0 3 C 2 RRR ... 4 C 3 RRR ... 5 C 4 RRR ... 6 C 5 RRR ... 7 H 1 RRR ... 8 H 1 RRR ... the rest of H atoms

Output of cc.show_versions()

atom b bond a angle d dihedral 1 C origin 1.001936 e_z 42.471043 e_x -0.0 2 C 1 1.448768 e_z 147.275520 e_x -0.0 3 C 1 1.353050 2 126.182617 e_x -0.0 4 H 1 1.084522 2 115.933171 3 -180.0 5 C 2 1.342678 1 122.754077 3 -180.0 6 H 2 1.082337 1 118.529921 3 -0.0 7 C 3 1.448768 1 126.182617 2 -0.0 8 H 3 1.084522 1 117.884211 2 -180.0 9 H 5 1.082067 2 120.862352 1 -0.0 10 H 5 1.079835 2 121.431689 1 -180.0 11 C 7 1.342678 3 122.754077 1 -180.0 12 H 7 1.082337 3 118.529921 1 -0.0 13 H 11 1.079835 7 121.431689 3 -180.0 14 H 11 1.082067 7 120.862352 3 -0.0
mcocdawc commented 5 years ago

If I understand you correctly and you are doing this for a problem regarding internal coordinates (e.g. vibrations), I would forget about the upper right triangle data with ['origin', 'e_x', 'e_y', 'e_z'] etc. This is just implemented to control the translational and rotational degrees of freedom of a molecule.

The following is the canonical way of supplying a partially predefined construction table.

molecule = cc.Cartesian.read_xyz('./ht-cis.xyz', start_index=1)

zmolecule = xyz_coord.get_zmat()

c_table_start = pd.DataFrame([['origin', 'e_z', 'e_x'],
                              [6, 'e_z', 'e_x'],
                              [4, 6, 'e_x']],
                              columns=['b', 'a', 'd'],
                              index=[6, 4, 2])

molecule.get_zmat(molecule.get_construction_table(fragment_list=[(molecule, frag_c_table)]))

PS: Arrays vs lists: I would not use arrays for dynamically appending elements. python lists are far better tailored for this task.

new_lst = []
for i in ...:
    new_lst.append(...)

click vs argparse: I would use click for constructing CLI tools from python.

#!/usr/bin/env python
import click
import chemcoord as cc
import pandas as pd

@click.command()
@click.option('--file', type=str, help="input xyz file")
def read_convert(file):
    molecule = cc.Cartesian.read_xyz(file, start_index=1)
    c_table_start = pd.DataFrame([['origin', 'e_z', 'e_x'],
                                  [6, 'e_z', 'e_x'],
                                  [4, 6, 'e_x']],
                                 columns=['b', 'a', 'd'],
                                 index=[6, 4, 2])
    c_table = molecule.get_construction_table(
        fragment_list=[(molecule, c_table_start)])
    print(molecule.get_zmat(c_table))

read_convert()

PPS: In the documentation it says:

  1. A list of tuples is given. Each tuple contains the fragment with its corresponding construction table in the form of:

[(frag1, c_table1), (frag2, c_table2)...] If the construction table of a fragment is not complete, the rest of each fragment’s construction table is calculated automatically.