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

Calculating energy #38

Closed
anneb101 opened this issue Aug 23, 2016 · 9 comments
Closed

Calculating energy #38

anneb101 opened this issue Aug 23, 2016 · 9 comments

Comments

@anneb101
Copy link

Hello,

I am trying to understand how exactly the energy of a group of atoms is calculated given an input of the atomic structure and the potential. From what I understand, after the potential is defined (pot = Potential(...)), this is executed as the pot.calc(a,'energy=energy') method (or something similar).

Would you be able to please point me to which exact file has this subroutine with the actual lines that do this calculation at the lowest level?

Thanks,
Anne

@albapa
Copy link
Member

albapa commented Aug 23, 2016

Hi Anne,

I think the syntax should be

pot.calc(a,energy=True)

It's difficult to say which subroutine is executed, because it depends on
the potential that you initialised. For example, for SW that will invoke a
chain like this:

potential_calc() in Potential.f95
Potential_Simple_Calc() in Potential_Simple.f95
IP_calc() in IP.f95
IPModel_SW_Calc() in IPModel_SW.f95

and that's the lowest level, apart from some auxiliary
routines IPModel_SW_Calc() might call.

I hope this helps.

Albert

On 23 August 2016 at 19:12, astronomydomine2 notifications@github.com
wrote:

Hello,

I am trying to understand how exactly the energy of a group of atoms is
calculated given an input of the atomic structure and the potential. From
what I understand, after the potential is defined (pot = Potential(...)),
this is executed as the pot.calc(a,'energy=energy') method (or something
similar).

Would you be able to please point me to which exact file has this
subroutine with the actual lines that do this calculation at the lowest
level?

Thanks,
Anne


You are receiving this because you are subscribed to this thread.
Reply to this email directly, view it on GitHub
#38, or mute the thread
https://github.com/notifications/unsubscribe-auth/AEKWQuGwhgdaDXC0BU7pwc4QPNts1lUhks5qizgVgaJpZM4JrOTa
.

Dr Albert Bartok-Partay
Leverhulme Early Career Fellow

Department of Engineering, University of Cambridge, Cambridge, CB21PZ
Magdalene College, Cambridge, CB30AG

T: +44 1223 748232
W: http://www.bartokpartay.com
E: ab686@cam.ac.uk
PGP public key: http://pgp.mit.edu (ID:FC646FB1)

@anneb101
Copy link
Author

Thanks a lot; is there any publication or documentation on how you got the exact formula to compute energy based on potential, for the different potentials?

Thanks
Anne

@albapa
Copy link
Member

albapa commented Aug 24, 2016

Most are based on the appropriate publications - sometimes this is
documented in the code, but very often it is evident given the name of the
potential. Which ones are you interested in?

Albert

On 24 August 2016 at 14:44, astronomydomine2 notifications@github.com
wrote:

Thanks a lot; is there any publication or documentation on how you got the
exact formula to compute energy based on potential, for the different
potentials?

Thanks
Anne


You are receiving this because you commented.
Reply to this email directly, view it on GitHub
#38 (comment), or mute
the thread
https://github.com/notifications/unsubscribe-auth/AEKWQkuTBlS5IaXPk9BTxbS7Uv2Me6dxks5qjEqzgaJpZM4JrOTa
.

@anneb101
Copy link
Author

Hi Albert,
I am interested in potentials for quartz and borosilicate glass.
From what I understand quartz uses: p = Potential('IP TS', param_filename='TS_params.xml')
Do you know what would be most appropriate for borosilicate glass?

I looked up the Yukawa potential/screened Coulomb potential but could not determine what sort of algorithm computes the energy. My goal is to get Ebulk and Eslab to calculate gamma (surface energy density)

Thanks
Anne

@jameskermode
Copy link
Member

The TS potential implemented in QUIP uses a short-range Coulomb interaction with the Yukawa method. The exact formulation as well as parameters for SiO2 (same as those in the .xml file you mention) are reported in this article:

http://scitation.aip.org/content/aip/journal/jcp/133/9/10.1063/1.3475565

It might be possible to use it for borosilicates too but would require reparameterisation.

@anneb101
Copy link
Author

Thank you James. I looked through the equations in the paper; is it the Eqq+Eqd+Edd sum for every atom pair, as listed in equations 4-6? I am a bit unclear still on how energy computed.

I am looking at surface.py. It looks like after the input of surface atoms, bulk atoms, and potential, the surface.energy is computed using pot.calc(surface...) and bulk.energy is calculated using pot.calc(bulk...). I am interested in understanding how exactly these quantities are computed before being plugged into bulk_energy_per_sio2 and the final returned surface energy value.

def surface_energy(pot, bulk, surface, dir=2):
if not hasattr(bulk, 'energy'):
bulk.calc_connect()
pot.calc(bulk, energy=True, force=True)
bulk_energy_per_sio2 = bulk.energy/(bulk.z == 14).sum()
surface.calc_connect()
pot.calc(surface, energy=True, force=True)
area = surface.cell_volume()/surface.lattice[dir,dir]
return (surface.energy - (surface.z == 14).sum()*bulk_energy_per_sio2)/(2.0*area)

@jameskermode
Copy link
Member

Those equations give the electrostatic contribution to the potential energy. There is also a short-range Morse-Stretch term. You can find all the details in the Calc interface in src/Potetials/IPModel_TS.f95:

https://github.com/libAtoms/QUIP/blob/public/src/Potentials/IPModel_TS.f95#L417

@gabor1
Copy link
Contributor

gabor1 commented Mar 23, 2017

Can I close this now?

@albapa
Copy link
Member

albapa commented Mar 23, 2017 via email

@gabor1 gabor1 closed this as completed Mar 23, 2017
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

4 participants