OpenMP support for the code from real-space fc3 to reciprocal-space fc3

This commit is contained in:
Atsushi Togo 2016-03-10 10:57:45 +09:00
parent c6b9d06963
commit 1d0e2a55dd
10 changed files with 107 additions and 59 deletions

View File

@ -36,7 +36,8 @@ from optparse import OptionParser
def get_parser():
parser = OptionParser()
parser.set_defaults(band_indices=None,
parser.set_defaults(average_pp_interaction=False,
band_indices=None,
band_paths=None,
band_points=None,
cell_filename=None,
@ -91,9 +92,13 @@ def get_parser():
mesh_numbers=None,
mesh_divisors=None,
no_kappa_stars=False,
output_filename=None,
phonon_supercell_dimension=None,
pinv_cutoff=1.0e-8,
pp_unit_conversion=None,
primitive_axis=None,
qpoints=None,
quiet=False,
q_direction=None,
read_amplitude=False,
read_collision=None,
@ -101,10 +106,6 @@ def get_parser():
read_fc3=False,
read_gamma=False,
run_with_g=True,
output_filename=None,
pp_unit_conversion=None,
qpoints=None,
quiet=False,
scattering_event_class=None,
show_num_triplets=False,
sigma=None,
@ -116,8 +117,8 @@ def get_parser():
tstep=None,
tsym_type=None,
uplo='L',
average_pp_interaction=False,
verbose=False,
with_g_zero=False,
write_amplitude=False,
write_collision=False,
write_detailed_gamma=False,
@ -385,6 +386,13 @@ def get_parser():
parser.add_option(
"-v", "--verbose", dest="verbose", action="store_true",
help="Detailed run-time information is displayed")
parser.add_option(
"--with_g_zero", dest="with_g_zero",
action="store_true",
help=("Avoid calculating ph-ph interaction for RTA conductivity when "
"integration weight is zero. This is activated when only using "
"tetrahedron method and no smearing method is used. Average "
"ph-ph interaction is not obtained with this option."))
parser.add_option(
"--wgp", "--write_grid_points", dest="write_grid_points",
action="store_true",

View File

@ -46,6 +46,7 @@ class Phono3pySettings(Settings):
self._pp_conversion_factor = None
self._scattering_event_class = None # scattering event class 1 or 2
self._temperatures = None
self._with_g_zero = False
self._write_amplitude = False
self._write_collision = False
self._write_detailed_gamma = False
@ -292,6 +293,12 @@ class Phono3pySettings(Settings):
def get_average_pp_interaction(self):
return self._average_pp_interaction
def set_with_g_zero(self, with_g_zero):
self._with_g_zero = with_g_zero
def get_with_g_zero(self):
return self._with_g_zero
def set_write_amplitude(self, write_amplitude):
self._write_amplitude = write_amplitude
@ -488,6 +495,10 @@ class Phono3pyConfParser(ConfParser):
if self._options.average_pp_interaction:
self._confs['average_pp_interaction'] = '.true.'
if opt.dest == 'with_g_zero':
if self._options.with_g_zero:
self._confs['with_g_zero'] = '.true.'
if opt.dest == 'write_amplitude':
if self._options.write_amplitude:
self._confs['write_amplitude'] = '.true.'
@ -710,6 +721,10 @@ class Phono3pyConfParser(ConfParser):
if confs['average_pp_interaction'] == '.true.':
self.set_parameter('average_pp_interaction', True)
if conf_key == 'with_g_zero':
if confs['with_g_zero'] == '.true.':
self.set_parameter('with_g_zero', True)
if conf_key == 'write_amplitude':
if confs['write_amplitude'] == '.true.':
self.set_parameter('write_amplitude', True)
@ -905,6 +920,10 @@ class Phono3pyConfParser(ConfParser):
self._settings.set_average_pp_interaction(
params['average_pp_interaction'])
# Avoid calculating useless ph-ph calculation for RTA conductivity
if params.has_key('with_g_zero'):
self._settings.set_with_g_zero(params['with_g_zero'])
# Write phonon-phonon interaction amplitudes to hdf5
if params.has_key('write_amplitude'):
self._settings.set_write_amplitude(params['write_amplitude'])

View File

@ -203,7 +203,7 @@ def show_phono3py_settings(settings,
if (settings.get_average_pp_interaction() and
(settings.get_is_bterta() or settings.get_is_lbte())):
print("Use averaged ph-ph interaction")
if log_level > 1:
print("Frequency factor to THz: %s" % frequency_factor_to_THz)
if frequency_step is not None:

View File

@ -428,6 +428,7 @@ class Phono3py:
is_kappa_star=True,
gv_delta_q=None, # for group velocity
run_with_g=True, # integration weights for smearing method, too
with_g_zero=False,
pinv_cutoff=1.0e-8, # for pseudo-inversion of collision matrix
write_gamma=False,
read_gamma=False,
@ -476,6 +477,7 @@ class Phono3py:
is_kappa_star=is_kappa_star,
gv_delta_q=gv_delta_q,
run_with_g=run_with_g,
with_g_zero=with_g_zero,
write_gamma=write_gamma,
read_gamma=read_gamma,
input_filename=input_filename,

View File

@ -24,6 +24,7 @@ def get_thermal_conductivity_RTA(
is_kappa_star=True,
gv_delta_q=1e-4,
run_with_g=True, # integration weights from gaussian smearing function
with_g_zero=False,
write_gamma=False,
read_gamma=False,
input_filename=None,
@ -52,6 +53,7 @@ def get_thermal_conductivity_RTA(
is_kappa_star=is_kappa_star,
gv_delta_q=gv_delta_q,
run_with_g=run_with_g,
with_g_zero=with_g_zero,
log_level=log_level)
if read_gamma:

View File

@ -300,7 +300,8 @@ static void real_to_normal(double *fc3_normal_squared,
shortest_vectors,
multiplicity,
p2s_map,
s2p_map);
s2p_map,
openmp_at_bands);
if (openmp_at_bands) {
#ifdef MEASURE_R2N

View File

@ -59,17 +59,21 @@ void real_to_reciprocal(lapack_complex_double *fc3_reciprocal,
const Darray *shortest_vectors,
const Iarray *multiplicity,
const int *p2s_map,
const int *s2p_map)
const int *s2p_map,
const int openmp_at_bands)
{
int i, j, k, num_patom;
int i, j, k, jk, num_patom;
double pre_phase;
lapack_complex_double pre_phase_factor;
num_patom = multiplicity->dims[1];
for (i = 0; i < num_patom; i++) {
for (j = 0; j < num_patom; j++) {
for (k = 0; k < num_patom; k++) {
if (openmp_at_bands) {
#pragma omp parallel for private(j, k)
for (jk = 0; jk < num_patom * num_patom; jk++) {
j = jk / num_patom;
k = jk % num_patom;
real_to_reciprocal_elements(fc3_reciprocal +
i * 27 * num_patom * num_patom +
j * 27 * num_patom +
@ -83,6 +87,23 @@ void real_to_reciprocal(lapack_complex_double *fc3_reciprocal,
i, j, k);
}
} else {
for (j = 0; j < num_patom; j++) {
for (k = 0; k < num_patom; k++) {
real_to_reciprocal_elements(fc3_reciprocal +
i * 27 * num_patom * num_patom +
j * 27 * num_patom +
k * 27,
q,
fc3,
shortest_vectors,
multiplicity,
p2s_map,
s2p_map,
i, j, k);
}
}
}
pre_phase = 0;

View File

@ -32,6 +32,7 @@
/* ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE */
/* POSSIBILITY OF SUCH DAMAGE. */
#include <stdlib.h>
#include <lapacke.h>
#include <math.h>
#include <phonoc_utils.h>
@ -134,7 +135,8 @@ void reciprocal_to_normal_squared_openmp
const int num_band,
const double cutoff_frequency)
{
int i, j, k, jk, bi, num_atom;
int i, j, k, l, ijk, jk, bi, num_atom, count;
int *n;
#ifdef MEASURE_R2N
double loopTotalCPUTime, loopTotalWallTime;
@ -144,61 +146,52 @@ void reciprocal_to_normal_squared_openmp
num_atom = num_band / 3;
for (i = 0; i < num_band0; i++) {
bi = band_indices[i];
if (freqs0[bi] > cutoff_frequency) {
n = (int*)malloc(sizeof(int) * num_band0 * num_band * num_band);
count = 0;
for (ijk = 0; ijk < num_band0 * num_band * num_band; ijk++) {
if (g_zero[ijk]) {
fc3_normal_squared[ijk] = 0;
} else {
n[count] = ijk;
count++;
}
}
#ifdef MEASURE_R2N
loopStartWallTime = time(NULL);
loopStartCPUTime = clock();
loopStartWallTime = time(NULL);
loopStartCPUTime = clock();
#endif
#pragma omp parallel for private(j, k)
for (jk = 0; jk < num_band * num_band; jk++) {
if (g_zero[i * num_band * num_band + jk]) {
fc3_normal_squared[i * num_band * num_band + jk] = 0;
continue;
}
j = jk / num_band;
k = jk % num_band;
fc3_normal_squared[i * num_band * num_band + jk] =
get_fc3_sum(i, j, k, bi,
freqs0, freqs1, freqs2,
eigvecs0, eigvecs1, eigvecs2,
fc3_reciprocal,
masses,
num_atom,
num_band,
cutoff_frequency);
}
#pragma omp parallel for private(ijk, i, j, k, jk, bi)
for (l = 0; l < count; l++) {
ijk = n[l];
i = ijk / (num_band * num_band);
jk = ijk % (num_band * num_band);
j = jk / num_band;
k = jk % num_band;
bi = band_indices[i];
if (freqs0[bi] > cutoff_frequency) {
fc3_normal_squared[n[l]] = get_fc3_sum(i, j, k, bi,
freqs0, freqs1, freqs2,
eigvecs0, eigvecs1, eigvecs2,
fc3_reciprocal,
masses,
num_atom,
num_band,
cutoff_frequency);
} else {
fc3_normal_squared[n[l]] = 0;
}
}
#ifdef MEASURE_R2N
loopTotalCPUTime = (double)(clock() - loopStartCPUTime) / CLOCKS_PER_SEC;
loopTotalWallTime = difftime(time(NULL), loopStartWallTime);
printf(" Band index %d/%d %1.3fs (%1.3fs CPU)\n",
i + 1, num_band0, loopTotalWallTime, loopTotalCPUTime);
/* #else */
/* printf("*"); */
/* if (i == (num_band0 - 1)) { */
/* printf("\n"); */
/* } */
/* if ((i % 20) == 0 && i != 0) { */
/* printf("\n"); */
/* } */
printf(" %1.3fs (%1.3fs CPU)\n", loopTotalWallTime, loopTotalCPUTime);
#endif
} else {
for (j = 0; j < num_band * num_band; j++) {
fc3_normal_squared[i * num_band * num_band + j] = 0;
}
#ifdef MEASURE_R2N
printf(" Band index %d/%d skipped due to frequency cutoff...\n",
i + 1, num_band0);
#endif
}
}
free(n);
}
static double get_fc3_sum

View File

@ -44,6 +44,7 @@ void real_to_reciprocal(lapack_complex_double *fc3_reciprocal,
const Darray *shortest_vectors,
const Iarray *multiplicity,
const int *p2s_map,
const int *s2p_map);
const int *s2p_map,
const int openmp_at_bands);
#endif

View File

@ -644,6 +644,7 @@ elif run_mode == "conductivity-RTA" or run_mode == "conductivity-LBTE":
is_kappa_star=settings.get_is_kappa_star(),
gv_delta_q=settings.get_group_velocity_delta_q(),
run_with_g=settings.get_run_with_g(),
with_g_zero=settings.get_with_g_zero(),
pinv_cutoff=settings.get_pinv_cutoff(),
write_gamma=settings.get_write_gamma(),
read_gamma=settings.get_read_gamma(),