Skip to content

Commit 8f71e5e

Browse files
slipchenkoslipchenko
authored
Merge pull request #14 from libefp2/qchem
synchronized with libefp in Q-Chem 6.1
2 parents 79d6cd8 + 136a9d1 commit 8f71e5e

8 files changed

Lines changed: 625 additions & 241 deletions

File tree

‎.DS_Store‎

2 KB
Binary file not shown.

‎Doxyfile‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -58,7 +58,7 @@ PROJECT_LOGO =
5858
# entered, it will be relative to the location where doxygen was started. If
5959
# left blank the current directory will be used.
6060

61-
OUTPUT_DIRECTORY =
61+
OUTPUT_DIRECTORY = ../WEB
6262

6363
# If the CREATE_SUBDIRS tag is set to YES then doxygen will create 4096 sub-
6464
# directories (in 2 levels) under the output directory of each output format and
@@ -790,7 +790,7 @@ WARN_LOGFILE =
790790
# spaces. See also FILE_PATTERNS and EXTENSION_MAPPING
791791
# Note: If this tag is empty the current directory is searched.
792792

793-
INPUT = src/efp.h
793+
INPUT = src efpmd/src
794794

795795
# This tag can be used to specify the character encoding of the source files
796796
# that doxygen parses. Internally doxygen uses the UTF-8 encoding. Doxygen uses

‎src/efp.c‎

Lines changed: 132 additions & 76 deletions
Original file line numberDiff line numberDiff line change
@@ -502,6 +502,9 @@ compute_two_body_range(struct efp *efp, size_t frag_from, size_t frag_to,
502502

503503
(void)data;
504504

505+
// ligand is a fragment (not QM)
506+
bool if_pairwise = efp->opts.enable_pairwise && efp->opts.ligand > -1;
507+
505508
#ifdef _OPENMP
506509
#pragma omp parallel for schedule(dynamic) reduction(+:e_elec,e_disp,e_xr,e_cp)
507510
#endif
@@ -533,7 +536,7 @@ compute_two_body_range(struct efp *efp, size_t frag_from, size_t frag_to,
533536
e_cp += ecp;
534537

535538
/* */
536-
if (efp->opts.enable_pairwise) {
539+
if (if_pairwise) {
537540
if (i == efp->ligand_index) {
538541
efp->pair_energies[fr_j].exchange_repulsion = exr;
539542
efp->pair_energies[fr_j].charge_penetration = ecp;
@@ -547,11 +550,10 @@ compute_two_body_range(struct efp *efp, size_t frag_from, size_t frag_to,
547550
}
548551
if (do_elec(&efp->opts) && efp->frags[i].n_multipole_pts > 0 &&
549552
efp->frags[fr_j].n_multipole_pts > 0) {
550-
e_elec_tmp = efp_frag_frag_elec(efp,
551-
i, fr_j);
553+
e_elec_tmp = efp_frag_frag_elec(efp, i, fr_j);
552554
e_elec += e_elec_tmp;
553555
/* */
554-
if (efp->opts.enable_pairwise) {
556+
if (if_pairwise) {
555557
if (i == efp->ligand_index)
556558
efp->pair_energies[fr_j].electrostatic = e_elec_tmp;
557559
if (fr_j == efp->ligand_index)
@@ -564,7 +566,7 @@ compute_two_body_range(struct efp *efp, size_t frag_from, size_t frag_to,
564566
i, fr_j, s, ds);
565567
e_disp += e_disp_tmp;
566568
/* */
567-
if (efp->opts.enable_pairwise) {
569+
if (if_pairwise) {
568570
if (i == efp->opts.ligand)
569571
efp->pair_energies[fr_j].dispersion = e_disp_tmp;
570572
if (fr_j == efp->opts.ligand)
@@ -582,6 +584,12 @@ compute_two_body_range(struct efp *efp, size_t frag_from, size_t frag_to,
582584
efp->energy.dispersion += e_disp;
583585
efp->energy.exchange_repulsion += e_xr;
584586
efp->energy.charge_penetration += e_cp;
587+
588+
if (efp->opts.print > 0) {
589+
printf(" In compute_two_body_range() \n");
590+
print_ene(&efp->energy);
591+
print_energies(efp);
592+
}
585593
}
586594

587595
EFP_EXPORT enum efp_result
@@ -1208,18 +1216,6 @@ efp_prepare(struct efp *efp)
12081216
efp->n_polarizable_pts += efp->frags[i].n_polarizable_pts;
12091217
}
12101218

1211-
/*
1212-
efp->n_fragment_field_pts = 0;
1213-
if (efp->opts.enable_pairwise && efp->opts.ligand != -1) {
1214-
size_t ligand_idx = efp->opts.ligand;
1215-
for (size_t i = 0; i < efp->n_frag; i++) {
1216-
efp->frags[i].fragment_field_offset = efp->n_fragment_field_pts;
1217-
efp->n_fragment_field_pts += efp->frags[ligand_idx].n_polarizable_pts;
1218-
}
1219-
}
1220-
efp->fragment_field = (vec_t *)calloc(efp->n_fragment_field_pts, sizeof(vec_t));
1221-
*/
1222-
12231219
efp->grad = (six_t *)calloc(efp->n_frag, sizeof(six_t));
12241220
efp->skiplist = (char *)calloc(efp->n_frag * efp->n_frag, 1);
12251221
efp->pair_energies = (struct efp_energy *)calloc(efp->n_frag, sizeof(struct efp_energy));
@@ -1421,38 +1417,28 @@ efp_get_frag_multipole_count(struct efp *efp, size_t frag_idx, size_t *n_mult)
14211417

14221418
/*
14231419
EFP_EXPORT enum efp_result
1424-
efp_get_frag_mult_rank(struct efp *efp, size_t frag_idx, size_t mult_idx, size_t *rank)
1420+
efp_get_frag_rank(struct efp *efp, size_t frag_idx, size_t *rank)
14251421
{
14261422
assert(efp);
14271423
assert(rank);
14281424
assert(frag_idx < efp->n_frag);
1429-
assert(mult_idx < efp->frags[frag_idx].n_multipole_pts);
14301425
14311426
struct frag *frag = efp->frags + frag_idx;
1432-
struct multipole_pt *pt = frag->multipole_pts + mult_idx;
1433-
1434-
*rank = -1;
1435-
1436-
for (size_t t = 0; t < 10; t++)
1437-
if (pt->octupole[t] != 0) {
1438-
*rank = 3;
1439-
return EFP_RESULT_SUCCESS;
1440-
}
14411427
1442-
for (size_t t = 0; t < 6; t++)
1443-
if (pt->quadrupole[t] != 0) {
1444-
*rank = 2;
1445-
return EFP_RESULT_SUCCESS;
1446-
}
1447-
1448-
if (pt->dipole.x != 0 || pt->dipole.y != 0 || pt->dipole.z != 0) {
1449-
*rank = 1;
1450-
return EFP_RESULT_SUCCESS;
1451-
}
1452-
1453-
if (pt->monopole != 0) {
1454-
*rank = 0;
1455-
return EFP_RESULT_SUCCESS;
1428+
*rank = 0;
1429+
for (size_t i=0; i<frag->n_multipole_pts; i++) {
1430+
struct multipole_pt *pt = frag->multipole_pts + i;
1431+
size_t rank_tmp = 0;
1432+
if (pt->if_dip)
1433+
rank_tmp = 1;
1434+
if (pt->if_quad)
1435+
rank_tmp = 2;
1436+
if (pt->if_oct)
1437+
rank_tmp = 3;
1438+
if (rank_tmp > *rank)
1439+
*rank = rank_tmp;
1440+
if (*rank == 3)
1441+
break;
14561442
}
14571443
14581444
return EFP_RESULT_SUCCESS;
@@ -1678,7 +1664,6 @@ efp_get_induced_dipole_conj_values(struct efp *efp, double *dip)
16781664

16791665
for (size_t i = 0; i < efp->n_frag; i++) {
16801666
struct frag *frag = efp->frags + i;
1681-
16821667
for (size_t j = 0; j < frag->n_polarizable_pts; j++) {
16831668
struct polarizable_pt *pt = frag->polarizable_pts + j;
16841669

@@ -1836,7 +1821,6 @@ efp_shutdown(struct efp *efp)
18361821
free(efp->ai_orbital_energies);
18371822
free(efp->ai_dipole_integrals);
18381823
free(efp->skiplist);
1839-
// free(efp->fragment_field);
18401824
free(efp->pair_energies);
18411825
free(efp->symmlist);
18421826
free(efp);
@@ -2007,27 +1991,6 @@ efp_create(void)
20071991
return efp;
20081992
}
20091993

2010-
EFP_EXPORT enum efp_result
2011-
efp_set_electron_density_field_fn(struct efp *efp,
2012-
efp_electron_density_field_fn fn)
2013-
{
2014-
assert(efp);
2015-
2016-
efp->get_electron_density_field = fn;
2017-
2018-
return EFP_RESULT_SUCCESS;
2019-
}
2020-
2021-
EFP_EXPORT enum efp_result
2022-
efp_set_electron_density_field_user_data(struct efp *efp, void *user_data)
2023-
{
2024-
assert(efp);
2025-
2026-
efp->get_electron_density_field_user_data = user_data;
2027-
2028-
return EFP_RESULT_SUCCESS;
2029-
}
2030-
20311994
EFP_EXPORT enum efp_result
20321995
efp_get_frag_count(struct efp *efp, size_t *n_frag)
20331996
{
@@ -2182,6 +2145,65 @@ efp_get_frag_mult_pt(struct efp *efp, size_t frag_idx, size_t pt_idx,
21822145
return EFP_RESULT_SUCCESS;
21832146
}
21842147

2148+
EFP_EXPORT enum efp_result
2149+
efp_get_frag_pol_pt(struct efp *efp, size_t frag_idx, size_t pt_idx,
2150+
struct efp_pol_pt *pol_pt)
2151+
{
2152+
assert(efp);
2153+
assert(pol_pt);
2154+
assert(frag_idx < efp->n_frag);
2155+
assert(pt_idx < efp->frags[frag_idx].n_polarizable_pts);
2156+
2157+
struct frag *frag;
2158+
frag = efp->frags + frag_idx;
2159+
struct polarizable_pt *pt;
2160+
pt = frag->polarizable_pts + pt_idx;
2161+
2162+
pol_pt->x = pt->x;
2163+
pol_pt->y = pt->y;
2164+
pol_pt->z = pt->z;
2165+
pol_pt->indip[0] = pt->indip.x;
2166+
pol_pt->indip[1] = pt->indip.y;
2167+
pol_pt->indip[2] = pt->indip.z;
2168+
pol_pt->indipconj[0] = pt->indipconj.x;
2169+
pol_pt->indipconj[1] = pt->indipconj.y;
2170+
pol_pt->indipconj[2] = pt->indipconj.z;
2171+
pol_pt->indip_gs[0] = pt->indip_gs.x;
2172+
pol_pt->indip_gs[1] = pt->indip_gs.y;
2173+
pol_pt->indip_gs[2] = pt->indip_gs.z;
2174+
pol_pt->indipconj_gs[0] = pt->indipconj_gs.x;
2175+
pol_pt->indipconj_gs[1] = pt->indipconj_gs.y;
2176+
pol_pt->indipconj_gs[2] = pt->indipconj_gs.z;
2177+
pol_pt->ai_field[0] = pt->elec_field_wf.x;
2178+
pol_pt->ai_field[1] = pt->elec_field_wf.y;
2179+
pol_pt->ai_field[2] = pt->elec_field_wf.z;
2180+
2181+
return EFP_RESULT_SUCCESS;
2182+
}
2183+
2184+
EFP_EXPORT enum efp_result
2185+
save_ai_field_pol_pt(struct efp *efp, struct efp_pol_pt *pol_pt, size_t frag_idx, size_t pt_idx) {
2186+
assert(efp);
2187+
assert(pol_pt);
2188+
assert(frag_idx < efp->n_frag);
2189+
assert(pt_idx < efp->frags[frag_idx].n_polarizable_pts);
2190+
2191+
struct frag *frag;
2192+
frag = efp->frags + frag_idx;
2193+
struct polarizable_pt *pt;
2194+
pt = frag->polarizable_pts + pt_idx;
2195+
2196+
pt->elec_field_wf.x = pol_pt->ai_field[0];
2197+
pt->elec_field_wf.y = pol_pt->ai_field[1];
2198+
pt->elec_field_wf.z = pol_pt->ai_field[2];
2199+
2200+
if (efp->opts.print > 1)
2201+
print_pol_pt(efp,frag_idx,pt_idx);
2202+
2203+
return EFP_RESULT_SUCCESS;
2204+
}
2205+
2206+
21852207
EFP_EXPORT void
21862208
efp_torque_to_derivative(const double *euler, const double *torque,
21872209
double *deriv)
@@ -2437,6 +2459,8 @@ print_pol_pt(struct efp *efp, size_t frag_index, size_t pol_index) {
24372459
printf(" conjugated induced dipole %lf %lf %lf\n", pt->indipconj.x, pt->indipconj.y, pt->indipconj.z);
24382460
printf(" old induced dipole %lf %lf %lf\n", pt->indip_old.x, pt->indip_old.y, pt->indip_old.z);
24392461
printf(" old conjugated induced dipole %lf %lf %lf\n", pt->indipconj_old.x, pt->indipconj_old.y, pt->indipconj_old.z);
2462+
printf(" gr state induced dipole %lf %lf %lf\n", pt->indip_gs.x, pt->indip_gs.y, pt->indip_gs.z);
2463+
printf(" gr state conjug induced dipole %lf %lf %lf\n", pt->indipconj_gs.x, pt->indipconj_gs.y, pt->indipconj_gs.z);
24402464

24412465
printf("\n");
24422466
}
@@ -2470,17 +2494,6 @@ print_frag_info(struct efp *efp, size_t frag_index) {
24702494
}
24712495

24722496
print_ligand(efp, frag_index);
2473-
2474-
for (int i=0; i < fr->n_multipole_pts; i++) {
2475-
struct efp_mult_pt *pt;
2476-
pt = (struct efp_mult_pt *)malloc(sizeof(struct efp_mult_pt));
2477-
2478-
efp_get_frag_mult_pt(efp, frag_index, i, pt);
2479-
print_efp_mult_pt(pt);
2480-
2481-
free(pt);
2482-
}
2483-
24842497
printf("\n");
24852498
}
24862499

@@ -2498,3 +2511,46 @@ print_efp_mult_pt(struct efp_mult_pt *pt) {
24982511
printf(" screen0 = %lf, if_screen0 = %s\n", pt->screen0, pt->if_screen ? "true" : "false");
24992512
printf("\n");
25002513
}
2514+
2515+
void
2516+
print_efp_pol_pt(struct efp_pol_pt *pt) {
2517+
printf("\n Polarizability point \n");
2518+
printf(" Coordinates %lf %lf %lf\n", pt->x, pt->y, pt->z);
2519+
printf(" induced dipole %lf %lf %lf\n", pt->indip[0], pt->indip[1], pt->indip[2]);
2520+
printf(" conjug induced dipole %lf %lf %lf\n", pt->indipconj[0], pt->indipconj[1], pt->indipconj[2]);
2521+
printf(" gr state induced dipole %lf %lf %lf\n", pt->indip_gs[0], pt->indip_gs[1], pt->indip_gs[2]);
2522+
printf(" gr state conjug induced dipole %lf %lf %lf\n", pt->indipconj_gs[0], pt->indipconj_gs[1], pt->indipconj_gs[2]);
2523+
printf(" AI field %lf %lf %lf\n", pt->ai_field[0], pt->ai_field[1], pt->ai_field[2]);
2524+
printf("\n");
2525+
}
2526+
2527+
void print_ene(struct efp_energy *energy) {
2528+
2529+
printf(" EFP ENERGY COMPONENTS \n");
2530+
2531+
printf(" ELECTROSTATIC ENERGY %lf \n", energy->electrostatic);
2532+
printf(" AI ELECTROSTATIC ENERGY %lf \n", energy->ai_electrostatic);
2533+
printf(" CHARGE PENETRATION ENERGY %lf \n", energy->charge_penetration);
2534+
printf(" POINT CHARGES ENERGY %lf \n", energy->electrostatic_point_charges);
2535+
printf(" POLARIZATION ENERGY %lf \n", energy->polarization);
2536+
printf(" EXC STATE POLARIZATION ENERGY %lf \n", energy->exs_polarization);
2537+
printf(" AI POLARIZATION ENERGY %lf \n", energy->ai_polarization);
2538+
printf(" DISPERSION ENERGY %lf \n", energy->dispersion);
2539+
printf(" AI DISPERSION ENERGY %lf \n", energy->ai_dispersion);
2540+
printf(" EXCHANGE-REPULSION ENERGY %lf \n", energy->exchange_repulsion);
2541+
double ene_sum = energy->electrostatic + energy->charge_penetration +
2542+
energy->electrostatic_point_charges + energy->polarization +
2543+
energy->dispersion + energy->exchange_repulsion;
2544+
printf(" SUM ENERGY %lf \n", ene_sum);
2545+
printf(" TOTAL ENERGY %lf \n", energy->total);
2546+
printf("\n");
2547+
}
2548+
2549+
void print_energies(struct efp *efp) {
2550+
printf(" --- PAIRWISE ENERGIES --- \n");
2551+
for (size_t i=0; i<efp->n_frag; i++) {
2552+
printf(" PAIR ENERGY on FRAGMENT %d %s \n", i, efp->frags[i].name);
2553+
print_ene(&efp->pair_energies[i]);
2554+
}
2555+
printf("\n");
2556+
}

0 commit comments

Comments
 (0)