mirror of https://github.com/phonopy/phonopy.git
Merge branch 'rc'
This commit is contained in:
commit
984914ca4d
|
@ -173,100 +173,6 @@ void py_get_dynamical_matrices_with_dd_openmp_over_qpoints(
|
|||
lambda, use_Wang_NAC);
|
||||
}
|
||||
|
||||
void py_get_dynamical_matrix(nb::ndarray<> py_dynamical_matrix,
|
||||
nb::ndarray<> py_force_constants,
|
||||
nb::ndarray<> py_q, nb::ndarray<> py_svecs,
|
||||
nb::ndarray<> py_multi, nb::ndarray<> py_masses,
|
||||
nb::ndarray<> py_s2p_map, nb::ndarray<> py_p2s_map,
|
||||
long use_openmp) {
|
||||
double(*dm)[2];
|
||||
double *fc;
|
||||
double *q;
|
||||
double(*svecs)[3];
|
||||
double *m;
|
||||
long(*multi)[2];
|
||||
long *s2p_map;
|
||||
long *p2s_map;
|
||||
long num_patom;
|
||||
long num_satom;
|
||||
|
||||
dm = (double(*)[2])py_dynamical_matrix.data();
|
||||
fc = (double *)py_force_constants.data();
|
||||
q = (double *)py_q.data();
|
||||
svecs = (double(*)[3])py_svecs.data();
|
||||
m = (double *)py_masses.data();
|
||||
multi = (long(*)[2])py_multi.data();
|
||||
s2p_map = (long *)py_s2p_map.data();
|
||||
p2s_map = (long *)py_p2s_map.data();
|
||||
num_patom = py_p2s_map.shape(0);
|
||||
num_satom = py_s2p_map.shape(0);
|
||||
|
||||
if (py_q.ndim() == 1) {
|
||||
phpy_get_dynamical_matrix_at_q(dm, num_patom, num_satom, fc, q, svecs,
|
||||
multi, m, s2p_map, p2s_map, NULL,
|
||||
use_openmp);
|
||||
} else {
|
||||
phpy_get_dynamical_matrices_openmp_over_qpoints(
|
||||
dm, num_patom, num_satom, fc, (double(*)[3])q, py_q.shape(0), svecs,
|
||||
multi, m, s2p_map, p2s_map, NULL);
|
||||
}
|
||||
}
|
||||
|
||||
void py_get_nac_dynamical_matrix(
|
||||
nb::ndarray<> py_dynamical_matrix, nb::ndarray<> py_force_constants,
|
||||
nb::ndarray<> py_q, nb::ndarray<> py_svecs, nb::ndarray<> py_multi,
|
||||
nb::ndarray<> py_masses, nb::ndarray<> py_s2p_map, nb::ndarray<> py_p2s_map,
|
||||
nb::ndarray<> py_q_cart, nb::ndarray<> py_born, double factor,
|
||||
long use_openmp) {
|
||||
double(*dm)[2];
|
||||
double *fc;
|
||||
double *q_cart;
|
||||
double *q;
|
||||
double(*svecs)[3];
|
||||
double *m;
|
||||
double(*born)[3][3];
|
||||
long(*multi)[2];
|
||||
long *s2p_map;
|
||||
long *p2s_map;
|
||||
long num_patom;
|
||||
long num_satom;
|
||||
|
||||
long n;
|
||||
double(*charge_sum)[3][3];
|
||||
|
||||
dm = (double(*)[2])py_dynamical_matrix.data();
|
||||
fc = (double *)py_force_constants.data();
|
||||
q_cart = (double *)py_q_cart.data();
|
||||
q = (double *)py_q.data();
|
||||
svecs = (double(*)[3])py_svecs.data();
|
||||
m = (double *)py_masses.data();
|
||||
born = (double(*)[3][3])py_born.data();
|
||||
multi = (long(*)[2])py_multi.data();
|
||||
s2p_map = (long *)py_s2p_map.data();
|
||||
p2s_map = (long *)py_p2s_map.data();
|
||||
num_patom = py_p2s_map.shape(0);
|
||||
num_satom = py_s2p_map.shape(0);
|
||||
|
||||
charge_sum =
|
||||
(double(*)[3][3])malloc(sizeof(double[3][3]) * num_patom * num_patom);
|
||||
n = num_satom / num_patom;
|
||||
|
||||
phpy_get_charge_sum(charge_sum, num_patom, factor / n, q_cart, born);
|
||||
|
||||
if (py_q.ndim() == 1) {
|
||||
phpy_get_dynamical_matrix_at_q(dm, num_patom, num_satom, fc, q, svecs,
|
||||
multi, m, s2p_map, p2s_map, charge_sum,
|
||||
use_openmp);
|
||||
} else {
|
||||
phpy_get_dynamical_matrices_openmp_over_qpoints(
|
||||
dm, num_patom, num_satom, fc, (double(*)[3])q, py_q.shape(0), svecs,
|
||||
multi, m, s2p_map, p2s_map, charge_sum);
|
||||
}
|
||||
|
||||
free(charge_sum);
|
||||
charge_sum = NULL;
|
||||
}
|
||||
|
||||
void py_get_recip_dipole_dipole(
|
||||
nb::ndarray<> py_dd, nb::ndarray<> py_dd_q0, nb::ndarray<> py_G_list,
|
||||
nb::ndarray<> py_q_cart, nb::ndarray<> py_q_direction,
|
||||
|
@ -668,8 +574,6 @@ NB_MODULE(_phonopy, m) {
|
|||
m.def("transpose_compact_fc", &py_transpose_compact_fc);
|
||||
m.def("dynamical_matrices_with_dd_openmp_over_qpoints",
|
||||
&py_get_dynamical_matrices_with_dd_openmp_over_qpoints);
|
||||
m.def("dynamical_matrix", &py_get_dynamical_matrix);
|
||||
m.def("nac_dynamical_matrix", &py_get_nac_dynamical_matrix);
|
||||
m.def("recip_dipole_dipole", &py_get_recip_dipole_dipole);
|
||||
m.def("recip_dipole_dipole_q0", &py_get_recip_dipole_dipole_q0);
|
||||
m.def("derivative_dynmat", &py_get_derivative_dynmat);
|
||||
|
|
25
c/dynmat.c
25
c/dynmat.c
|
@ -193,8 +193,6 @@ static void add_dynmat_dd_at_q(
|
|||
double q_cart[3];
|
||||
double mm;
|
||||
|
||||
q_dir_cart = NULL;
|
||||
|
||||
dd = (double(*)[2])malloc(sizeof(double[2]) * num_patom * num_patom * 9);
|
||||
get_q_cart(q_cart, q, reciprocal_lattice);
|
||||
dym_get_recip_dipole_dipole(dd, dd_q0, G_list, num_G_points, num_patom,
|
||||
|
@ -218,29 +216,6 @@ static void add_dynmat_dd_at_q(
|
|||
dd = NULL;
|
||||
}
|
||||
|
||||
long dym_get_dynamical_matrices_openmp_over_qpoints(
|
||||
double (*dynamical_matrices)[2], // [q-points, num_band, num_band,
|
||||
// (real, imag)]
|
||||
const long num_patom, const long num_satom, const double *fc,
|
||||
const double (*qpoints)[3], const long n_qpoints, const double (*svecs)[3],
|
||||
const long (*multi)[2], const double *mass, const long *s2p_map,
|
||||
const long *p2s_map, const double (*charge_sum)[3][3]) {
|
||||
long i, adrs_shift;
|
||||
|
||||
adrs_shift = num_patom * num_patom * 9;
|
||||
|
||||
#ifdef _OPENMP
|
||||
#pragma omp parallel for
|
||||
#endif
|
||||
for (i = 0; i < n_qpoints; i++) {
|
||||
dym_get_dynamical_matrix_at_q(
|
||||
dynamical_matrices + adrs_shift * i, num_patom, num_satom, fc,
|
||||
qpoints[i], svecs, multi, mass, s2p_map, p2s_map, charge_sum, 0);
|
||||
}
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
||||
/// @brief charge_sum is NULL if G-L NAC or no-NAC.
|
||||
long dym_get_dynamical_matrix_at_q(double (*dynamical_matrix)[2],
|
||||
const long num_patom, const long num_satom,
|
||||
|
|
|
@ -45,11 +45,6 @@ long dym_dynamical_matrices_with_dd_openmp_over_qpoints(
|
|||
const double *q_direction, const double nac_factor,
|
||||
const double (*dd_q0)[2], const double (*G_list)[3],
|
||||
const long num_G_points, const double lambda, const long use_Wang_NAC);
|
||||
long dym_get_dynamical_matrices_openmp_over_qpoints(
|
||||
double (*dynamical_matrices)[2], const long num_patom, const long num_satom,
|
||||
const double *fc, const double (*qpoints)[3], const long n_qpoints,
|
||||
const double (*svecs)[3], const long (*multi)[2], const double *mass,
|
||||
const long *s2p_map, const long *p2s_map, const double (*charge_sum)[3][3]);
|
||||
long dym_get_dynamical_matrix_at_q(double (*dynamical_matrix)[2],
|
||||
const long num_patom, const long num_satom,
|
||||
const double *fc, const double q[3],
|
||||
|
|
24
c/phonopy.c
24
c/phonopy.c
|
@ -97,30 +97,6 @@ long phpy_dynamical_matrices_with_dd_openmp_over_qpoints(
|
|||
num_G_points, lambda, use_Wang_NAC);
|
||||
}
|
||||
|
||||
long phpy_get_dynamical_matrices_openmp_over_qpoints(
|
||||
double (*dynamical_matrices)[2], const long num_patom, const long num_satom,
|
||||
const double *fc, const double (*qpoints)[3], const long n_qpoints,
|
||||
const double (*svecs)[3], const long (*multi)[2], const double *mass,
|
||||
const long *s2p_map, const long *p2s_map,
|
||||
const double (*charge_sum)[3][3]) {
|
||||
return dym_get_dynamical_matrices_openmp_over_qpoints(
|
||||
dynamical_matrices, num_patom, num_satom, fc, qpoints, n_qpoints, svecs,
|
||||
multi, mass, s2p_map, p2s_map, charge_sum);
|
||||
}
|
||||
|
||||
long phpy_get_dynamical_matrix_at_q(double (*dynamical_matrix)[2],
|
||||
const long num_patom, const long num_satom,
|
||||
const double *fc, const double q[3],
|
||||
const double (*svecs)[3],
|
||||
const long (*multi)[2], const double *mass,
|
||||
const long *s2p_map, const long *p2s_map,
|
||||
const double (*charge_sum)[3][3],
|
||||
const long use_openmp) {
|
||||
return dym_get_dynamical_matrix_at_q(dynamical_matrix, num_patom, num_satom,
|
||||
fc, q, svecs, multi, mass, s2p_map,
|
||||
p2s_map, charge_sum, use_openmp);
|
||||
}
|
||||
|
||||
void phpy_get_charge_sum(
|
||||
double (*charge_sum)[3][3], const long num_patom,
|
||||
const double factor, /* 4pi/V*unit-conv and denominator */
|
||||
|
|
13
c/phonopy.h
13
c/phonopy.h
|
@ -58,19 +58,6 @@ long phpy_dynamical_matrices_with_dd_openmp_over_qpoints(
|
|||
const double *q_direction, const double nac_factor,
|
||||
const double (*dd_q0)[2], const double (*G_list)[3],
|
||||
const long num_G_points, const double lambda, const long use_Wang_NAC);
|
||||
long phpy_get_dynamical_matrices_openmp_over_qpoints(
|
||||
double (*dynamical_matrices)[2], const long num_patom, const long num_satom,
|
||||
const double *fc, const double (*qpoints)[3], const long n_qpoints,
|
||||
const double (*svecs)[3], const long (*multi)[2], const double *mass,
|
||||
const long *s2p_map, const long *p2s_map, const double (*charge_sum)[3][3]);
|
||||
long phpy_get_dynamical_matrix_at_q(double (*dynamical_matrix)[2],
|
||||
const long num_patom, const long num_satom,
|
||||
const double *fc, const double q[3],
|
||||
const double (*svecs)[3],
|
||||
const long (*multi)[2], const double *mass,
|
||||
const long *s2p_map, const long *p2s_map,
|
||||
const double (*charge_sum)[3][3],
|
||||
const long use_openmp);
|
||||
void phpy_get_charge_sum(
|
||||
double (*charge_sum)[3][3], const long num_patom,
|
||||
const double factor, /* 4pi/V*unit-conv and denominator */
|
||||
|
|
42
c/rgrid.c
42
c/rgrid.c
|
@ -38,19 +38,13 @@
|
|||
#include <stddef.h>
|
||||
#include <stdio.h>
|
||||
|
||||
static void get_all_grid_addresses(long grid_address[][3], const long mesh[3]);
|
||||
static long get_double_grid_index(const long address_double[3],
|
||||
const long mesh[3]);
|
||||
static long get_grid_index_single_mesh(const long address[3],
|
||||
const long mesh[3]);
|
||||
static void reduce_grid_address(long address[3], const long mesh[3]);
|
||||
static void reduce_double_grid_address(long address[3], const long mesh[3]);
|
||||
static long mat_modulo_l(const long a, const long b);
|
||||
|
||||
void rgd_get_all_grid_addresses(long grid_address[][3], const long mesh[3]) {
|
||||
get_all_grid_addresses(grid_address, mesh);
|
||||
}
|
||||
|
||||
long rgd_get_double_grid_index(const long address_double[3],
|
||||
const long mesh[3]) {
|
||||
return get_double_grid_index(address_double, mesh);
|
||||
|
@ -66,30 +60,6 @@ void rgd_get_double_grid_address(long address_double[3], const long address[3],
|
|||
reduce_double_grid_address(address_double, mesh);
|
||||
}
|
||||
|
||||
static void get_all_grid_addresses(long grid_address[][3], const long mesh[3]) {
|
||||
long i, j, k;
|
||||
long grid_index;
|
||||
long address[3];
|
||||
|
||||
for (i = 0; i < mesh[0]; i++) {
|
||||
address[0] = i;
|
||||
for (j = 0; j < mesh[1]; j++) {
|
||||
address[1] = j;
|
||||
for (k = 0; k < mesh[2]; k++) {
|
||||
address[2] = k;
|
||||
grid_index = get_grid_index_single_mesh(address, mesh);
|
||||
|
||||
assert(mesh[0] * mesh[1] * mesh[2] > grid_index);
|
||||
|
||||
grid_address[grid_index][0] = address[0];
|
||||
grid_address[grid_index][1] = address[1];
|
||||
grid_address[grid_index][2] = address[2];
|
||||
reduce_grid_address(grid_address[grid_index], mesh);
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
static long get_double_grid_index(const long address_double[3],
|
||||
const long mesh[3]) {
|
||||
long i;
|
||||
|
@ -118,18 +88,6 @@ static long get_grid_index_single_mesh(const long address[3],
|
|||
#endif
|
||||
}
|
||||
|
||||
static void reduce_grid_address(long address[3], const long mesh[3]) {
|
||||
long i;
|
||||
|
||||
for (i = 0; i < 3; i++) {
|
||||
#ifndef GRID_BOUNDARY_AS_NEGATIVE
|
||||
address[i] -= mesh[i] * (address[i] > mesh[i] / 2);
|
||||
#else
|
||||
address[i] -= mesh[i] * (address[i] > (mesh[i] - 1) / 2);
|
||||
#endif
|
||||
}
|
||||
}
|
||||
|
||||
static void reduce_double_grid_address(long address[3], const long mesh[3]) {
|
||||
long i;
|
||||
|
||||
|
|
|
@ -70,7 +70,6 @@
|
|||
/* without GRID_BOUNDARY_AS_NEGATIVE, e.g., [-2, -1, 0, 1, 2, 3]. */
|
||||
/* with GRID_BOUNDARY_AS_NEGATIVE, e.g., [-3, -2, -1, 0, 1, 2]. */
|
||||
|
||||
void rgd_get_all_grid_addresses(long grid_address[][3], const long mesh[3]);
|
||||
long rgd_get_double_grid_index(const long address_double[3],
|
||||
const long mesh[3]);
|
||||
void rgd_get_double_grid_address(long address_double[3], const long address[3],
|
||||
|
|
|
@ -2,7 +2,11 @@
|
|||
|
||||
# Change Log
|
||||
|
||||
## Jul-13-2024: Version 2.26.3
|
||||
## Jul-14-2024: Version 2.26.5
|
||||
|
||||
- Refactoring of dynamical matrix code.
|
||||
|
||||
## Jul-13-2024: Version 2.26.4
|
||||
|
||||
- Collection of minor updates to follow updates of external libraries.
|
||||
|
||||
|
|
|
@ -55,7 +55,7 @@ copyright = "2009, Atsushi Togo"
|
|||
# The short X.Y version.
|
||||
version = "2.26"
|
||||
# The full version, including alpha/beta/rc tags.
|
||||
release = "2.26.4"
|
||||
release = "2.26.5"
|
||||
|
||||
# The language for content autogenerated by Sphinx. Refer to documentation
|
||||
# for a list of supported languages.
|
||||
|
|
|
@ -175,6 +175,16 @@ class Phonopy:
|
|||
self._is_symmetry = is_symmetry
|
||||
self._calculator = calculator
|
||||
self._store_dense_svecs = store_dense_svecs
|
||||
if int(self.version.split(".")[0]) > 2 and not store_dense_svecs:
|
||||
warnings.warn(
|
||||
(
|
||||
"store_dense_svecs=False is not supported in Phonopy v3"
|
||||
"and later versions."
|
||||
),
|
||||
DeprecationWarning,
|
||||
stacklevel=2,
|
||||
)
|
||||
self._store_dense_svecs = True
|
||||
self._use_SNF_supercell = use_SNF_supercell
|
||||
self._log_level = log_level
|
||||
|
||||
|
|
|
@ -286,42 +286,7 @@ class DynamicalMatrix:
|
|||
self._force_constants = np.array(fc, dtype="double", order="C")
|
||||
|
||||
def _run_c_dynamical_matrix(self, q):
|
||||
import phonopy._phonopy as phonoc
|
||||
|
||||
fc = self._force_constants
|
||||
mass = self._pcell.masses
|
||||
size_prim = len(mass)
|
||||
dm_shape = size_prim * 3, size_prim * 3
|
||||
dm = np.zeros(dm_shape + (2,), dtype="double")
|
||||
|
||||
if fc.shape[0] == fc.shape[1]: # full-fc
|
||||
s2p_map = self._s2p_map
|
||||
p2s_map = self._p2s_map
|
||||
else: # compact-fc
|
||||
s2p_map = self._s2pp_map
|
||||
p2s_map = np.arange(len(self._p2s_map), dtype="int_")
|
||||
|
||||
phonoc.dynamical_matrix(
|
||||
dm,
|
||||
fc,
|
||||
np.array(q, dtype="double"),
|
||||
self._svecs,
|
||||
self._multi,
|
||||
mass,
|
||||
s2p_map,
|
||||
p2s_map,
|
||||
self._use_openmp * 1,
|
||||
)
|
||||
|
||||
# Data of dm array are stored in memory by the C order of
|
||||
# (size_prim * 3, size_prim * 3, 2), where the last 2 means
|
||||
# real and imaginary parts. This code assumes this memory
|
||||
# order is that expected by numpy. Otherwise, numpy complex array
|
||||
# should be created as follows:
|
||||
# dm_double = dm.view(dtype='double').reshape(size_prim * 3,
|
||||
# size_prim * 3, 2)
|
||||
# dm = dm_double[:, :, 0] + 1j * dm_double[:, :, 1]
|
||||
self._dynamical_matrix = dm.view(dtype=self._dtype_complex).reshape(dm_shape)
|
||||
self._dynamical_matrix = run_dynamical_matrix_solver_c(self, q, is_nac=False)
|
||||
|
||||
def _run_py_dynamical_matrix(self, q):
|
||||
"""Python implementation of building dynamical matrix.
|
||||
|
@ -582,7 +547,7 @@ class DynamicalMatrixGL(DynamicalMatrixNAC):
|
|||
shape=(primitive atoms, supercell atoms, 3, 3) for compact FC.
|
||||
dtype='double'
|
||||
with_full_terms : bool, optional
|
||||
When False, only reciprocal terms are considered for NAC. This is the
|
||||
When False, only reciprocal terms are considered for NAC. False is the
|
||||
default and the reasonable choice.
|
||||
decimals : int, optional, default=None
|
||||
Number of decimals. Use like dm.round(decimals).
|
||||
|
@ -753,18 +718,23 @@ class DynamicalMatrixGL(DynamicalMatrixNAC):
|
|||
self.show_nac_message()
|
||||
|
||||
def _compute_dynamical_matrix(self, q_red, q_direction):
|
||||
if self._Gonze_force_constants is None:
|
||||
self.make_Gonze_nac_dataset()
|
||||
if self._with_full_terms:
|
||||
if self._Gonze_force_constants is None:
|
||||
self.make_Gonze_nac_dataset()
|
||||
|
||||
if self._log_level > 2:
|
||||
print("%d %s" % (self._Gonze_count + 1, q_red))
|
||||
self._Gonze_count += 1
|
||||
fc = self._force_constants
|
||||
self._force_constants = self._Gonze_force_constants
|
||||
self._run(q_red)
|
||||
self._force_constants = fc
|
||||
dm_dd = self._get_Gonze_dipole_dipole(q_red, q_direction)
|
||||
self._dynamical_matrix += dm_dd
|
||||
if self._log_level > 2:
|
||||
print("%d %s" % (self._Gonze_count + 1, q_red))
|
||||
self._Gonze_count += 1
|
||||
fc = self._force_constants
|
||||
self._force_constants = self._Gonze_force_constants
|
||||
self._run(q_red)
|
||||
self._force_constants = fc
|
||||
dm_dd = self._get_Gonze_dipole_dipole(q_red, q_direction)
|
||||
self._dynamical_matrix += dm_dd
|
||||
else:
|
||||
self._dynamical_matrix = run_dynamical_matrix_solver_c(
|
||||
self, q_red, q_direction
|
||||
)
|
||||
|
||||
def _get_Gonze_dipole_dipole(self, q_red, q_direction):
|
||||
num_atom = len(self._pcell)
|
||||
|
@ -1069,19 +1039,22 @@ class DynamicalMatrixWang(DynamicalMatrixNAC):
|
|||
|
||||
def _compute_dynamical_matrix(self, q_red, q_direction):
|
||||
# Wang method (J. Phys.: Condens. Matter 22 (2010) 202201)
|
||||
if q_direction is None:
|
||||
q_cart = np.dot(q_red, self._rec_lat.T)
|
||||
else:
|
||||
q_cart = np.dot(q_direction, self._rec_lat.T)
|
||||
|
||||
constant = self._get_constant_factor(
|
||||
q_cart, self._dielectric, self._pcell.volume, self._unit_conversion
|
||||
)
|
||||
try:
|
||||
import phonopy._phonopy as phonoc # noqa F401
|
||||
|
||||
self._run_c_Wang_dynamical_matrix(q_red, q_cart, constant)
|
||||
self._dynamical_matrix = run_dynamical_matrix_solver_c(
|
||||
self, q_red, q_direction
|
||||
)
|
||||
# self._run_c_Wang_dynamical_matrix(q_red, q_cart, constant)
|
||||
except ImportError:
|
||||
if q_direction is None:
|
||||
q_cart = np.dot(q_red, self._rec_lat.T)
|
||||
else:
|
||||
q_cart = np.dot(q_direction, self._rec_lat.T)
|
||||
|
||||
constant = self._get_constant_factor(
|
||||
q_cart, self._dielectric, self._pcell.volume, self._unit_conversion
|
||||
)
|
||||
num_atom = len(self._pcell)
|
||||
fc_backup = self._force_constants.copy()
|
||||
nac_q = self._get_charge_sum(num_atom, q_cart, self._born) * constant
|
||||
|
@ -1089,47 +1062,6 @@ class DynamicalMatrixWang(DynamicalMatrixNAC):
|
|||
self._run(q_red)
|
||||
self._force_constants = fc_backup
|
||||
|
||||
def _run_c_Wang_dynamical_matrix(self, q_red, q_cart, factor):
|
||||
import phonopy._phonopy as phonoc
|
||||
|
||||
fc = self._force_constants
|
||||
mass = self._pcell.masses
|
||||
size_prim = len(mass)
|
||||
dm = np.zeros((size_prim * 3, size_prim * 3), dtype=self._dtype_complex)
|
||||
|
||||
if fc.shape[0] == fc.shape[1]: # full fc
|
||||
phonoc.nac_dynamical_matrix(
|
||||
dm.view(dtype="double"),
|
||||
fc,
|
||||
np.array(q_red, dtype="double"),
|
||||
self._svecs,
|
||||
self._multi,
|
||||
mass,
|
||||
self._s2p_map,
|
||||
self._p2s_map,
|
||||
np.array(q_cart, dtype="double"),
|
||||
self._born,
|
||||
factor,
|
||||
self._use_openmp * 1,
|
||||
)
|
||||
else:
|
||||
phonoc.nac_dynamical_matrix(
|
||||
dm.view(dtype="double"),
|
||||
fc,
|
||||
np.array(q_red, dtype="double"),
|
||||
self._svecs,
|
||||
self._multi,
|
||||
mass,
|
||||
self._s2pp_map,
|
||||
np.arange(len(self._p2s_map), dtype="int_"),
|
||||
np.array(q_cart, dtype="double"),
|
||||
self._born,
|
||||
factor,
|
||||
self._use_openmp * 1,
|
||||
)
|
||||
|
||||
self._dynamical_matrix = dm
|
||||
|
||||
def _get_charge_sum(self, num_atom, q, born):
|
||||
nac_q = np.zeros((num_atom, num_atom, 3, 3), dtype="double", order="C")
|
||||
A = np.dot(q, born)
|
||||
|
@ -1213,21 +1145,52 @@ def get_dynamical_matrix(
|
|||
def run_dynamical_matrix_solver_c(
|
||||
dm: Union[DynamicalMatrix, DynamicalMatrixWang, DynamicalMatrixGL],
|
||||
qpoints,
|
||||
nac_q_direction: Optional[np.ndarray] = None, # in reduced coordinates
|
||||
nac_q_direction: Optional[np.ndarray] = None,
|
||||
is_nac: Optional[bool] = None,
|
||||
):
|
||||
"""Bulid and solve dynamical matrices on grid in C-API.
|
||||
|
||||
If dynamical matrices at many qpoints are calculated, it is recommended not
|
||||
to use this function one qpoint by one qpoint to avoid overhead in the
|
||||
preparation steps.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
dm : DynamicalMatrix
|
||||
DynamicalMatrix instance.
|
||||
qpoints : array_like,
|
||||
q-points in crystallographic coordinates.
|
||||
shape=(n_qpoints, 3), dtype='double', order='C'
|
||||
q-points in crystallographic coordinates. shape=(n_qpoints, 3),
|
||||
dtype='double', order='C'
|
||||
nac_q_direction : array_like, optional
|
||||
See Interaction.nac_q_direction. Default is None.
|
||||
Direction of q from Gamma point given in reduced coordinates. This is
|
||||
only activated when q-point->[0,0,0] case. For example, this is used for
|
||||
q->[0,0,0] where approaching direction is known, e.g., band structure
|
||||
calculation. Default is None.
|
||||
is_nac : bool, optional
|
||||
True if NAC is considered. Default is None. If None, it is determined
|
||||
from dm.is_nac().
|
||||
|
||||
"""
|
||||
import phonopy._phonopy as phonoc
|
||||
|
||||
if (
|
||||
isinstance(qpoints, np.ndarray)
|
||||
and qpoints.dtype == np.dtype("double")
|
||||
and qpoints.flags.aligned
|
||||
and qpoints.flags.owndata
|
||||
and qpoints.flags.c_contiguous
|
||||
):
|
||||
_qpoints = qpoints
|
||||
else:
|
||||
_qpoints = np.array(qpoints, dtype="double", order="C")
|
||||
qpoints_ndim = _qpoints.ndim
|
||||
_qpoints = _qpoints.reshape(-1, 3)
|
||||
|
||||
if is_nac is None:
|
||||
_is_nac = dm.is_nac()
|
||||
else:
|
||||
_is_nac = is_nac
|
||||
|
||||
(
|
||||
svecs,
|
||||
multi,
|
||||
|
@ -1240,7 +1203,7 @@ def run_dynamical_matrix_solver_c(
|
|||
) = _extract_params(dm)
|
||||
|
||||
use_Wang_NAC = False
|
||||
if isinstance(dm, DynamicalMatrixGL):
|
||||
if _is_nac and dm.nac_method == "gonze":
|
||||
gonze_nac_dataset = dm.Gonze_nac_dataset
|
||||
if gonze_nac_dataset[0] is None:
|
||||
dm.make_Gonze_nac_dataset()
|
||||
|
@ -1258,8 +1221,7 @@ def run_dynamical_matrix_solver_c(
|
|||
G_list = np.zeros(3) # dummy value
|
||||
Lambda = 0.0 # dummy value
|
||||
fc = dm.force_constants
|
||||
if isinstance(dm, DynamicalMatrixWang):
|
||||
use_Wang_NAC = True
|
||||
use_Wang_NAC = _is_nac and dm.nac_method == "wang"
|
||||
|
||||
if nac_q_direction is None:
|
||||
is_nac_q_zero = True
|
||||
|
@ -1268,23 +1230,13 @@ def run_dynamical_matrix_solver_c(
|
|||
is_nac_q_zero = False
|
||||
_nac_q_direction = np.array(nac_q_direction, dtype="double")
|
||||
|
||||
if born is None:
|
||||
_born = np.zeros(9) # dummy variable
|
||||
else:
|
||||
_born = born
|
||||
|
||||
if dielectric is None:
|
||||
_dielectric = np.zeros(9) # dummy variable
|
||||
else:
|
||||
_dielectric = dielectric
|
||||
|
||||
p2s, s2p = _get_fc_elements_mapping(dm, fc)
|
||||
|
||||
dtype_complex = "c%d" % (np.dtype("double").itemsize * 2)
|
||||
dynmat = np.zeros((len(qpoints), len(p2s) * 3, len(p2s) * 3), dtype=dtype_complex)
|
||||
phonoc.dynamical_matrices_with_dd_openmp_over_qpoints(
|
||||
dynmat.view(dtype="double"),
|
||||
np.array(qpoints, dtype="double", order="C"),
|
||||
_qpoints,
|
||||
fc,
|
||||
svecs,
|
||||
multi,
|
||||
|
@ -1293,19 +1245,22 @@ def run_dynamical_matrix_solver_c(
|
|||
s2p,
|
||||
p2s,
|
||||
_nac_q_direction,
|
||||
_born,
|
||||
_dielectric,
|
||||
born,
|
||||
dielectric,
|
||||
rec_lattice,
|
||||
float(nac_factor),
|
||||
nac_factor,
|
||||
dd_q0,
|
||||
G_list,
|
||||
Lambda,
|
||||
dm.is_nac() * 1,
|
||||
_is_nac * 1,
|
||||
is_nac_q_zero * 1,
|
||||
use_Wang_NAC * 1, # use_Wang_NAC
|
||||
)
|
||||
|
||||
return dynmat
|
||||
if qpoints_ndim == 1:
|
||||
return dynmat[0]
|
||||
else:
|
||||
return dynmat
|
||||
|
||||
|
||||
def _extract_params(dm: Union[DynamicalMatrix, DynamicalMatrixNAC]):
|
||||
|
@ -1316,17 +1271,17 @@ def _extract_params(dm: Union[DynamicalMatrix, DynamicalMatrixNAC]):
|
|||
else:
|
||||
_svecs, _multi = sparse_to_dense_svecs(svecs, multi)
|
||||
|
||||
masses = np.array(dm.primitive.masses, dtype="double")
|
||||
masses = dm.primitive.masses
|
||||
rec_lattice = np.array(np.linalg.inv(dm.primitive.cell), dtype="double", order="C")
|
||||
positions = np.array(dm.primitive.positions, dtype="double", order="C")
|
||||
if isinstance(dm, DynamicalMatrixNAC):
|
||||
positions = dm.primitive.positions
|
||||
if dm.is_nac():
|
||||
born = dm.born
|
||||
nac_factor = dm.nac_factor
|
||||
nac_factor = float(dm.nac_factor)
|
||||
dielectric = dm.dielectric_constant
|
||||
else:
|
||||
born = None
|
||||
nac_factor = 0
|
||||
dielectric = None
|
||||
born = np.zeros(9) # dummy variable
|
||||
nac_factor = 0.0 # dummy variable
|
||||
dielectric = np.zeros(9) # dummy variable
|
||||
|
||||
return (
|
||||
_svecs,
|
||||
|
@ -1340,19 +1295,15 @@ def _extract_params(dm: Union[DynamicalMatrix, DynamicalMatrixNAC]):
|
|||
)
|
||||
|
||||
|
||||
def _get_fc_elements_mapping(dm, fc):
|
||||
def _get_fc_elements_mapping(dm: DynamicalMatrix, fc: np.ndarray):
|
||||
p2s_map = dm.primitive.p2s_map
|
||||
s2p_map = dm.primitive.s2p_map
|
||||
if fc.shape[0] == fc.shape[1]: # full fc
|
||||
fc_p2s = p2s_map
|
||||
fc_s2p = s2p_map
|
||||
return np.array(p2s_map, dtype="int_"), np.array(s2p_map, dtype="int_")
|
||||
else: # compact fc
|
||||
primitive = dm.primitive
|
||||
p2p_map = primitive.p2p_map
|
||||
s2pp_map = np.array(
|
||||
[p2p_map[s2p_map[i]] for i in range(len(s2p_map))], dtype="intc"
|
||||
[p2p_map[s2p_map[i]] for i in range(len(s2p_map))], dtype="int_"
|
||||
)
|
||||
fc_p2s = np.arange(len(p2s_map), dtype="intc")
|
||||
fc_s2p = s2pp_map
|
||||
|
||||
return np.array(fc_p2s, dtype="int_"), np.array(fc_s2p, dtype="int_")
|
||||
return np.arange(len(p2s_map), dtype="int_"), s2pp_map
|
||||
|
|
|
@ -211,7 +211,9 @@ class PhonopyAtoms:
|
|||
@property
|
||||
def positions(self):
|
||||
"""Setter and getter of positions in Cartesian coordinates."""
|
||||
return np.dot(self._scaled_positions, self._cell)
|
||||
return np.array(
|
||||
np.dot(self._scaled_positions, self._cell), dtype="double", order="C"
|
||||
)
|
||||
|
||||
@positions.setter
|
||||
def positions(self, positions):
|
||||
|
|
|
@ -34,4 +34,4 @@
|
|||
# ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
|
||||
# POSSIBILITY OF SUCH DAMAGE.
|
||||
|
||||
__version__ = "2.26.4"
|
||||
__version__ = "2.26.5"
|
||||
|
|
Loading…
Reference in New Issue