Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Z matrix to Cartesian converter uses average atomic weights #10

Open
AlexB67 opened this issue May 22, 2020 · 0 comments
Open

Z matrix to Cartesian converter uses average atomic weights #10

AlexB67 opened this issue May 22, 2020 · 0 comments

Comments

@AlexB67
Copy link

AlexB67 commented May 22, 2020

Hello
Thanks for the great youtube courses and utilities. I watched a many :)

I think the weights should be that of the most abundant isotope or relevant isotope in question, not the average weights. This is how other packages do it for ab-initio.

I use this reference:
https://physics.nist.gov/cgi-bin/Compositions/stand_alone.pl?ele=&ascii=ascii

I noticed while testing my C++ Z matrix code v psi4 which agree, but yours code gives very slightly different centre of mass coordinates when I add the following code (to your utility) to get the centre of mass adjusted Cartesian coordinates. Ultimately, this is caused by the weights you use.

Sorry. I don't know python to write code in normally, but this works okay.

print xyz coordinates to screen

def print_coords(self, comment):
    angstrom_to_bohr = 1.889725989
    total_mass = 0.0
    
    for i in range(self.n_atoms):
        total_mass += at_masses[self.atoms[i].attype]

    x_com = 0.0
    y_com = 0.0
    z_com = 0.0

    for i in range(self.n_atoms):
        x_com += self.atoms[i].coords[0] * at_masses[self.atoms[i].attype]
        y_com += self.atoms[i].coords[1] * at_masses[self.atoms[i].attype]
        z_com += self.atoms[i].coords[2] * at_masses[self.atoms[i].attype]
    
    x_com /= total_mass
    y_com /= total_mass
    z_com /= total_mass

    for i in range(self.n_atoms):
        self.atoms[i].coords[0] -= x_com
        self.atoms[i].coords[1] -= y_com
        self.atoms[i].coords[2] -= z_com
        
    print('\n')
    print('%s\n' % (comment), end='')
    for i in range(self.n_atoms):
        if self.atoms[i].attype == 'X':
                continue
        
        print('%-2s' % (self.atoms[i].attype), end='')
        
        for j in range(3):
            print(' %14.10f' % (self.atoms[i].coords[j] * angstrom_to_bohr), end='')
        
        print('\n', end='')

I think a bit higher precision in the printed output and definitions of weights would also be useful due to getting rounding errors you get from the cartesians you produce when used in calculations when comparing with others.

o' course weight are not used anyway in the utility to convert, so it is minor. The coordinates without COM adjustment are correct up to that point to within rounding errors.

Best regards
...

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

1 participant