Skip to content

Commit

Permalink
pull trunk develop
Browse files Browse the repository at this point in the history
  • Loading branch information
VicTrqt committed Oct 17, 2024
2 parents fdf4780 + 2dc6bd2 commit d324222
Show file tree
Hide file tree
Showing 13 changed files with 2,554 additions and 51 deletions.
5 changes: 4 additions & 1 deletion .github/workflows/test.yml
Original file line number Diff line number Diff line change
Expand Up @@ -25,6 +25,9 @@ jobs:
config:
- os: ubuntu-latest
python: '3.11'
config:
- os: ubuntu-latest
python: '3.12'

runs-on: ${{ matrix.config.os }}
env:
Expand All @@ -51,7 +54,7 @@ jobs:
cp abipy/data/managers/gh_manager.yml $HOME/.abinit/abipy/manager.yml
cp abipy/data/managers/gh_scheduler.yml $HOME/.abinit/abipy/scheduler.yml
# FIXME TEMPORARY HACK
pip install git+https://github.com/gmatteo/pymatgen.git@master -U
#pip install git+https://github.com/gmatteo/pymatgen.git@master -U
- name: pytest
run: |
Expand Down
30 changes: 17 additions & 13 deletions abipy/abio/factories.py
Original file line number Diff line number Diff line change
Expand Up @@ -1381,16 +1381,17 @@ def scf_for_phonons(structure, pseudos, kppa=None, ecut=None, pawecutdg=None, nb

def ddkpert_from_gsinput(gs_input, ddk_pert, nband=None, use_symmetries=False, ddk_tol=None, manager=None) -> AbinitInput:
"""
Returns an |AbinitInput| to perform a DDK calculations for a specific perturbation and based on a ground state |AbinitInput|.
Returns an |AbinitInput| to perform a DDK calculations for a specific perturbation based on a ground state |AbinitInput|.
Args:
gs_input: an |AbinitInput| representing a ground state calculation, likely the SCF performed to get the WFK.
gs_input: an |AbinitInput| representing a ground state calculation, likely the SCF performed to get the WFK.
ddk_pert: dict with the Abinit variables defining the perturbation
Example: {'idir': 1, 'ipert': 4, 'qpt': [0.0, 0.0, 0.0]},
use_symmetries: boolean that determines if the irreducible components of the perturbation are used.
Default to False. (TODO: Should be implemented)
ddk_tol: a dictionary with a single key defining the type of tolerance used for the DDK calculations and its value. Default: {"tolvrs": 1.0e-22}.
manager: |TaskManager| of the task. If None, the manager is initialized from the config file.
ddk_tol: a dictionary with a single key defining the type of tolerance used for the DDK calculations and its value.
Default: {"tolvrs": 1.0e-22}.
manager: |TaskManager| of the task. If None, the manager is initialized from the config file.
"""
gs_input = gs_input.deepcopy()
gs_input.pop_irdvars()
Expand All @@ -1409,19 +1410,20 @@ def ddkpert_from_gsinput(gs_input, ddk_pert, nband=None, use_symmetries=False, d

return ddk_inp


def ddepert_from_gsinput(gs_input, dde_pert, use_symmetries=True, dde_tol=None, manager=None) -> AbinitInput:
"""
Returns an |AbinitInput| to perform a DDE calculations for a specific perturbation and based on a ground state |AbinitInput|.
Returns an |AbinitInput| to perform a DDE calculations for a specific perturbation based on a ground state |AbinitInput|.
Args:
gs_input: an |AbinitInput| representing a ground state calculation, likely the SCF performed to get the WFK.
gs_input: an |AbinitInput| representing a ground state calculation, likely the SCF performed to get the WFK.
dde_pert: dict with the Abinit variables defining the perturbation
Example: {'idir': 1, 'ipert': 4, 'qpt': [0.0, 0.0, 0.0]},
Example: {'idir': 1, 'ipert': 4, 'qpt': [0.0, 0.0, 0.0]},
use_symmetries: boolean that determines if the irreducible components of the perturbation are used.
Default to True. Should be set to False for nonlinear coefficients calculation.
dde_tol: a dictionary with a single key defining the type of tolerance used for the DDE calculations and
dde_tol: a dictionary with a single key defining the type of tolerance used for the DDE calculations and
its value. Default: {"tolvrs": 1.0e-22}.
manager: |TaskManager| of the task. If None, the manager is initialized from the config file.
manager: |TaskManager| of the task. If None, the manager is initialized from the config file.
"""
gs_input = gs_input.deepcopy()
gs_input.pop_irdvars()
Expand All @@ -1434,15 +1436,16 @@ def ddepert_from_gsinput(gs_input, dde_pert, use_symmetries=True, dde_tol=None,

return dde_inp


def dtepert_from_gsinput(gs_input, dte_pert, manager=None) -> AbinitInput:
"""
Returns an |AbinitInput| to perform a DTE calculations for a specific perturbation and based on a ground state |AbinitInput|.
Args:
gs_input: an |AbinitInput| representing a ground state calculation, likely the SCF performed to get the WFK.
gs_input: an |AbinitInput| representing a ground state calculation, likely the SCF performed to get the WFK.
dte_pert: dict with the Abinit variables defining the perturbation
Example: {'idir': 1, 'ipert': 4, 'qpt': [0.0, 0.0, 0.0]},
manager: |TaskManager| of the task. If None, the manager is initialized from the config file.
Example: {'idir': 1, 'ipert': 4, 'qpt': [0.0, 0.0, 0.0]},
manager: |TaskManager| of the task. If None, the manager is initialized from the config file.
"""
gs_input = gs_input.deepcopy()
gs_input.pop_irdvars()
Expand All @@ -1452,6 +1455,7 @@ def dtepert_from_gsinput(gs_input, dte_pert, manager=None) -> AbinitInput:

return dte_inp


def dte_from_gsinput(gs_input, use_phonons=True, ph_tol=None, ddk_tol=None, dde_tol=None,
skip_dte_permutations=False, manager=None) -> MultiDataset:
"""
Expand Down Expand Up @@ -1509,7 +1513,7 @@ def dte_from_gsinput(gs_input, use_phonons=True, ph_tol=None, ddk_tol=None, dde_
gs_input.set_vars(nband=nband)
gs_input.pop('nbdbuf', None)
multi_dte = gs_input.make_dte_inputs(phonon_pert=use_phonons, skip_permutations=skip_dte_permutations,
manager=manager)
manager=manager)
multi_dte.add_tags(atags.DTE)
multi.extend(multi_dte)

Expand Down
43 changes: 37 additions & 6 deletions abipy/abio/tests/test_factories.py
Original file line number Diff line number Diff line change
@@ -1,3 +1,4 @@
import json
import abipy.data as abidata
import abipy.abilab as abilab

Expand All @@ -7,8 +8,9 @@
from abipy.abio.factories import *
from abipy.abio.factories import (BandsFromGsFactory, IoncellRelaxFromGsFactory, HybridOneShotFromGsFactory,
ScfForPhononsFactory, PhononsFromGsFactory, PiezoElasticFactory, PiezoElasticFromGsFactory, ShiftMode)
from abipy.abio.factories import _find_nscf_nband_from_gsinput
import json
from abipy.abio.factories import _find_nscf_nband_from_gsinput, minimal_scf_input
from abipy.abio.input_tags import DDK, DDE, PH_Q_PERT, STRAIN, DTE, PH_Q_PERT


write_inputs_to_json = False

Expand Down Expand Up @@ -304,7 +306,6 @@ def test_phonons_from_gsinput(self):
with_bec=False, ph_tol=None, ddk_tol=None, dde_tol=None)
self.abivalidate_multi(multi)

from abipy.abio.input_tags import DDK, DDE, PH_Q_PERT
inp_ddk = multi.filter_by_tags(DDK)[0]
inp_dde = multi.filter_by_tags(DDE)[0]
inp_ph_q_pert_1 = multi.filter_by_tags(PH_Q_PERT)[0]
Expand Down Expand Up @@ -499,7 +500,7 @@ def test_dfpt_from_gsinput(self):
do_dte=True, ph_tol=None, ddk_tol=None, dde_tol=None)
self.abivalidate_multi(multi)

from abipy.abio.input_tags import DDK, DDE, PH_Q_PERT, STRAIN, DTE

inp_ddk = multi.filter_by_tags(DDK)[0]
inp_dde = multi.filter_by_tags(DDE)[0]
inp_ph_q_pert_1 = multi.filter_by_tags(PH_Q_PERT)[0]
Expand Down Expand Up @@ -530,10 +531,40 @@ def test_dfpt_from_gsinput(self):
self.assert_input_equality('dfpt_from_gsinput_dte.json', inp_dte)

def test_minimal_scf_input(self):
from abipy.abio.factories import minimal_scf_input
inp = minimal_scf_input(self.si_structure, self.si_pseudo)

self.abivalidate_input(inp)

self.assertEqual(inp["nband"], 1)
self.assertEqual(inp["nstep"], 0)

def test_ddkpert_from_gsinput(self):
gs_inp = gs_input(self.si_structure, self.si_pseudo, kppa=None, ecut=2, spin_mode="unpolarized")
gs_inp["nband"] = 4
gs_inp["autoparal"] = 1
gs_inp["npfft"] = 10

ddk_pert = {'idir': 1, 'ipert': 3, 'qpt': [0.0, 0.0, 0.0]}
ddk_input = ddkpert_from_gsinput(gs_inp, ddk_pert)
assert ddk_input["tolwfr"] == 1.0e-22
assert "autoparal" not in ddk_input
assert "npfft" not in ddk_input
self.abivalidate_input(ddk_input)

dde_pert = {'idir': 1, 'ipert': 4, 'qpt': [0.0, 0.0, 0.0]}
dde_input = ddepert_from_gsinput(gs_inp, dde_pert)
assert "autoparal" not in dde_input
assert "npfft" not in dde_input
assert dde_input["tolvrs"] == 1.0e-22
self.abivalidate_input(dde_input)

#dte_pert = {'i1dir': 1, 'ipert': 4, 'qpt': [0.0, 0.0, 0.0]}
#dte_input = dtepert_from_gsinput(gs_inp, dte_pert)
#assert "autoparal" not in dte_input
#assert "npfft" not in dte_input
#assert dte_input["tolvrs"] == 1.0e-22
#self.abivalidate_input(dte_input)

#def test_dte_from_gsinput(self):
# gs_inp = gs_input(self.si_structure, self.si_pseudo, kppa=None, ecut=2, spin_mode="unpolarized")
# multi = dte_from_gsinput(gs_inp, use_phonons=True)
# self.abivalidate_input(multi)
1 change: 1 addition & 0 deletions abipy/core/abinit_units.py
Original file line number Diff line number Diff line change
Expand Up @@ -119,6 +119,7 @@ def wlabel_from_units(units: str, unicode=False) -> str:
'cm-1': r'Frequency (cm$^{-1}$)',
'cm^-1': r'Frequency (cm$^{-1}$)',
'thz': r'Frequency (Thz)',
'hbar': r'Angular momentum ($\hbar$)',
}
try:
s = d[units.lower().strip()]
Expand Down
193 changes: 193 additions & 0 deletions abipy/dfpt/deformation_utils.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,193 @@
# deformation_utils.py

import numpy as np
from pymatgen.core import Structure, Lattice, Element
from abipy.core.symmetries import AbinitSpaceGroup
from pymatgen.symmetry.analyzer import SpacegroupAnalyzer
from abipy.abio.inputs import AbinitInput
import re


def generate_deformations_volumic(structure, eps_V=0.02, scales=None):
if scales is None:
scales = [-1, 0, 1, 2, 3]
rprim = structure.lattice.matrix
structures_new = {}

for i in scales:
rprim2 = np.copy(rprim)
rprim2[:, :] = rprim[:, :] * (1.00 + eps_V * i)**(1/3.)

structure2 = structure.copy()
structure2.lattice = Lattice(rprim2)
#structure2.scale_lattice(structure2.volume*(1.00 + eps_V * i))
namei = int(round(1000 * (1.00 + eps_V * i)))
formatted_namei = f"{namei:04d}"
structures_new[formatted_namei] = structure2


return structures_new

def generate_deformations(structure , eps=0.005):
spgrp = AbinitSpaceGroup.from_structure(structure )
print (spgrp)
spgrp_number=spgrp.spgid
rprim= structure.lattice.matrix

rprim2 = np.copy(rprim)
rprim_new = {}
structures_new = {}

if 1 <= spgrp_number <= 2:
disp=[[1,1,1,1,1,1], [0,1,1,1,1,1], [2,1,1,1,1,1], [1,0,1,1,1,1], [1,2,1,1,1,1], [1,1,0,1,1,1],
[1,1,2,1,1,1], [1,1,1,0,1,1], [1,1,1,2,1,1], [1,1,1,1,0,1], [1,1,1,1,2,1], [1,1,1,1,1,0],
[1,1,1,1,1,2], [0,0,1,1,1,1], [1,0,0,1,1,1], [1,1,0,0,1,1], [1,1,1,0,0,1], [1,1,1,1,0,0],
[0,1,0,1,1,1], [0,1,1,0,1,1], [0,1,1,1,0,1], [0,1,1,1,1,0], [1,0,1,0,1,1], [1,0,1,1,0,1],
[1,0,1,1,1,0], [1,1,0,1,0,1], [1,1,0,1,1,0], [1,1,1,0,1,0] , [0 ,0,0,0,0,0]]
if abs(rprim[1, 0]) > 1e-9 or abs(rprim[2, 0]) > 1e-9 or abs(rprim[2, 1]) > 1e-9:
print("Warning: The lattice is oriented such that xz =xy =yz =0 .")
rprim0 = np.copy(rprim)
a=rprim[0, :]
b=rprim[1, :]
c=rprim[2, :]
norm_a = np.linalg.norm(a)
norm_b = np.linalg.norm(b)
norm_c = np.linalg.norm(c)

# Compute angles between vectors
cos_ab = np.dot(a, b) / (norm_a * norm_b)
cos_ac = np.dot(a, c) / (norm_a * norm_c)
cos_bc = np.dot(b, c) / (norm_b * norm_c)

rprim0[0,0] = 1.0
rprim0[0,1] = 0.0
rprim0[0,2] = 0.0
rprim0[1,0] = cos_ab
rprim0[1,1] = np.sqrt(1-cos_ab**2)
rprim0[1,2] = 0.0
rprim0[2,0] = cos_ac
rprim0[2,1] = (cos_bc-rprim0[1,0]*rprim0[2,0])/rprim0[1,1]
rprim0[2,2] = np.sqrt(1.0-rprim0[2,0]**2-rprim0[2,1]**2)
rprim0[0,:] = rprim0[0,:]*norm_a
rprim0[1,:] = rprim0[1,:]*norm_b
rprim0[2,:] = rprim0[2,:]*norm_c
print("Old rprim:")
print(rprim)
print("New rprim:")
print(rprim0)

for pair in disp:
i,j,k,l,m,n = pair
rprim2[ :,0] = rprim0[ :,0] * (1.00 + eps * i) + rprim0[ :,1] * (eps * l) +rprim0[ :,2] * (eps * m)
rprim2[ :,1] = rprim0[ :,1] * (1.00 + eps * j) + rprim0[ :,2] * (eps * n)
rprim2[ :,2] = rprim0[ :,2] * (1.00 + eps * k)

namei = int(round(1000 * (1.00 + eps * i)))
namej = int(round(1000 * (1.00 + eps * j)))
namek = int(round(1000 * (1.00 + eps * k)))
namel = int(round(1000 * (1.00 + eps * l)))
namem = int(round(1000 * (1.00 + eps * m)))
namen = int(round(1000 * (1.00 + eps * n)))
formatted_namei = f"{namei:04d}_{namej:04d}_{namek:04d}_{namel:04d}_{namem:04d}_{namen:04d}"

structure2=structure.copy()
structure2.lattice=Lattice(rprim2)
structures_new[formatted_namei] = structure2

return structures_new
elif 3 <= spgrp_number <= 15:
disp=[[1,1,1,1], [0,1,1,1], [2,1,1,1], [1,0,1,1], [1,2,1,1], [1,1,0,1], [1,1,2,1], [1,1,1,0],
[1,1,1,2], [0,0,1,1], [1,0,0,1], [1,1,0,0], [0,1,0,1], [1,0,1,0], [0,1,1,0]]
if abs(rprim[1, 0]) > 1e-9 or abs(rprim[0, 1]) > 1e-9 or abs(rprim[2, 1]) > 1e-9 or abs(rprim[1, 2]) > 1e-9:
print("Error: Monoclinic structure with yx=xy=0 and yz=zy=0 lattice required.")
elif abs(rprim[0, 2]) > 1e-9 :
print("Warning: The lattice is oriented such that xz = 0.")
rprim0 = np.copy(rprim)
a=rprim[0, :]
b=rprim[1, :]
c=rprim[2, :]
norm_a = np.linalg.norm(a)
norm_b = np.linalg.norm(b)
norm_c = np.linalg.norm(c)

# Compute angles between vectors
cos_ab = np.dot(a, b) / (norm_a * norm_b)
cos_ac = np.dot(a, c) / (norm_a * norm_c)
cos_bc = np.dot(b, c) / (norm_b * norm_c)

rprim0[0,0] = norm_a
rprim0[0,2] = 0.0
rprim0[1,1] = norm_b
rprim0[2,0] = norm_c*cos_ac
rprim0[2,2] = norm_c*np.sqrt(1-cos_ac**2)
print("Old rprim:")
print(rprim)
print("New rprim:")
print(rprim0)

for pair in disp:
i,j,k,l = pair
rprim2[ :,0] = rprim0[ :,0] * (1.00 + eps * i) +rprim0[ :,2] * (eps * l)
rprim2[ :,1] = rprim0[ :,1] * (1.00 + eps * j)
rprim2[ :,2] = rprim0[ :,2] * (1.00 + eps * k)

namei = int(round(1000 * (1.00 + eps * i)))
namej = int(round(1000 * (1.00 + eps * j)))
namek = int(round(1000 * (1.00 + eps * k)))
namel = int(round(1000 * (1.00 + eps * l)))
formatted_namei = f"{namei:04d}_{namej:04d}_{namek:04d}_{namel:04d}"

structure2=structure.copy()
structure2.lattice=Lattice(rprim2)
structures_new[formatted_namei] = structure2

return structures_new
elif 16 <= spgrp_number <= 74:
disp=[[0,0,1],[0,1,0],[1,0,0],[1,1,1],[0,1,1],[2,1,1],[1,0,1],[1,2,1],[1,1,0],[1,1,2]]
for pair in disp:
i,j,k = pair
rprim2[ :,0] = rprim[ :,0] * (1.00 + eps * i)
rprim2[ :,1] = rprim[ :,1] * (1.00 + eps * j)
rprim2[ :,2] = rprim[ :,2] * (1.00 + eps * k)

namei = int(round(1000 * (1.00 + eps * i)))
namej = int(round(1000 * (1.00 + eps * j)))
namek = int(round(1000 * (1.00 + eps * k)))
formatted_namei = f"{namei:04d}_{namej:04d}_{namek:04d}"

structure2=structure.copy()
structure2.lattice=Lattice(rprim2)
structures_new[formatted_namei] = structure2

return structures_new
elif 75 <= spgrp_number <= 194:
disp=[[0,0],[1,1],[0,1],[2,1],[1,0],[1,2]]
for pair in disp:
i, k = pair
rprim2[ :,0] = rprim[ :,0] * (1.00 + eps * i)
rprim2[ :,1] = rprim[ :,1] * (1.00 + eps * i)
rprim2[ :,2] = rprim[ :,2] * (1.00 + eps * k)

namei = int(round(1000 * (1.00 + eps * i)))
namek = int(round(1000 * (1.00 + eps * k)))
formatted_namei = f"{namei:04d}_{namek:04d}"
rprim_new[formatted_namei] = rprim2

structure2=structure.copy()
structure2.lattice=Lattice(rprim2)
structures_new[formatted_namei] = structure2

return structures_new
elif 195 <= spgrp_number <= 230:
for i in range(3):
rprim2[ :,0] = rprim[ :,0] * (1.00 + eps * i)
rprim2[ :,1] = rprim[ :,1] * (1.00 + eps * i)
rprim2[ :,2] = rprim[ :,2] * (1.00 + eps * i)
namei = int(round(1000 * (1.00 + eps * i)))
formatted_namei = f"{namei:04d}"

structure2=structure.copy()
structure2.lattice=Lattice(rprim2)
structures_new[formatted_namei] = structure2
return structures_new

Loading

0 comments on commit d324222

Please sign in to comment.