Skip to content

Commit

Permalink
make d3 optional again
Browse files Browse the repository at this point in the history
Signed-off-by: Conrad Hübler <Conrad.Huebler@gmx.net>
  • Loading branch information
conradhuebler committed May 28, 2024
1 parent f33e2b9 commit dabc733
Show file tree
Hide file tree
Showing 6 changed files with 53 additions and 4 deletions.
15 changes: 14 additions & 1 deletion src/core/dftd3interface.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -47,23 +47,29 @@ DFTD3Interface::DFTD3Interface(const json& controller)
m_damping = parameter["d_damping"];
m_functional = parameter["d_func"];
CreateParameter();
#ifdef USE_D3
m_error = dftd3_new_error();
#endif
}

DFTD3Interface::DFTD3Interface()
{
#ifdef USE_D3
m_error = dftd3_new_error();
#endif
}

DFTD3Interface::~DFTD3Interface()
{
#ifdef USE_D3
dftd3_delete_model(&m_disp);
dftd3_delete_structure(&m_mol);
dftd3_delete_error(&m_error);
dftd3_delete_param(&m_param);

delete m_coord;
delete m_attyp;
#endif
}

void DFTD3Interface::PrintParameter() const
Expand All @@ -73,6 +79,7 @@ void DFTD3Interface::PrintParameter() const

void DFTD3Interface::CreateParameter()
{
#ifdef USE_D3
if (m_d3_a1 > 1e-8 || m_d3_a2 > 1e-8 || m_d3_s6 > 1e-8 || m_d3_s8 > 1e-8 || m_d3_s9 > 1e-8) {
if (m_damping.compare("bj")) {
m_param = dftd3_new_rational_damping(m_error, m_d3_s6, m_d3_s8, m_d3_s9, m_d3_a1, m_d3_a2, m_d3_alp);
Expand Down Expand Up @@ -101,6 +108,7 @@ void DFTD3Interface::CreateParameter()
}
delete[] cstr;
}
#endif
}

void DFTD3Interface::UpdateParameters(const json& controller)
Expand Down Expand Up @@ -148,9 +156,10 @@ bool DFTD3Interface::InitialiseMolecule(const std::vector<int>& atomtypes)
m_coord[3 * i + 2] = (3 * i + 2) / au;
m_attyp[i] = atomtypes[i];
}

#ifdef USE_D3
m_mol = dftd3_new_structure(m_error, atomtypes.size(), m_attyp, m_coord, NULL, NULL);
m_disp = dftd3_new_d3_model(m_error, m_mol);
#endif

return true;
}
Expand All @@ -167,15 +176,19 @@ double DFTD3Interface::DFTD3Calculation(double* grad)
double energy = 0;
double sigma[9];

#ifdef USE_D3
dftd3_update_structure(m_error, m_mol, m_coord, NULL);
dftd3_get_dispersion(m_error, m_mol, m_disp, m_param, &energy, grad, sigma);
#endif

return energy;
}

void DFTD3Interface::clear()
{
#ifdef USE_D3
dftd3_delete_model(&m_disp);
dftd3_delete_structure(&m_mol);
dftd3_delete_error(&m_error);
#endif
}
6 changes: 5 additions & 1 deletion src/core/dftd3interface.h
Original file line number Diff line number Diff line change
Expand Up @@ -20,8 +20,9 @@
#pragma once

#include "src/tools/general.h"

#ifdef USE_D3
#include "s-dftd3.h"
#endif

#include "src/core/molecule.h"

Expand Down Expand Up @@ -85,8 +86,11 @@ class DFTD3Interface {

double m_bet = 8;
bool m_atm = false;

#ifdef USE_D3
dftd3_error m_error;
dftd3_structure m_mol;
dftd3_model m_disp;
dftd3_param m_param;
#endif
};
10 changes: 9 additions & 1 deletion src/core/forcefieldthread.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -423,16 +423,21 @@ D3Thread::D3Thread(int thread, int threads)
setAutoDelete(false);
m_final_factor = 1 / 2625.15 * 4.19;
m_d = 1e-7;
#ifdef USE_D3
m_d3 = new DFTD3Interface();
#endif
}

D3Thread::~D3Thread()
{
#ifdef USE_D3
delete m_d3;
#endif
}

int D3Thread::execute()
{
#ifdef USE_D3
for (int i = 0; i < m_atom_types.size(); ++i) {
m_d3->UpdateAtom(i, m_geometry(i, 0), m_geometry(i, 1), m_geometry(i, 2));
}
Expand All @@ -447,7 +452,10 @@ int D3Thread::execute()
}
} else
m_vdw_energy = m_d3->DFTD3Calculation(0);

#else
std::cerr << "D3 is not included, sorry for that" << std::endl;
exit(1);
#endif
return 0;
}

Expand Down
6 changes: 6 additions & 0 deletions src/core/forcefieldthread.h
Original file line number Diff line number Diff line change
Expand Up @@ -168,18 +168,24 @@ class D3Thread : public ForceFieldThread {
~D3Thread();
void setParamater(const json& parameter)
{
#ifdef USE_D3
m_d3->UpdateParametersD3(parameter);
#endif
}

void Initialise(const std::vector<int>& atom_types)
{
m_atom_types = atom_types;
#ifdef USE_D3
m_d3->InitialiseMolecule(m_atom_types);
#endif
}
virtual int execute() override;

private:
#ifdef USE_D3
DFTD3Interface* m_d3;
#endif
std::vector<int> m_atom_types;
};

Expand Down
17 changes: 16 additions & 1 deletion src/core/qmdff.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -665,8 +665,10 @@ QMDFF::QMDFF(const json& controller)
d3settings["d_a1"] = a1;
d3settings["d_a2"] = a2;
d3settings["d_alp"] = 1;

#ifdef USE_D3
m_d3 = new DFTD3Interface(d3settings);
#endif

m_h4correction.allocate(m_atom_types.size());

readFF(parameter);
Expand All @@ -685,7 +687,10 @@ QMDFF::QMDFF(const json& controller)

QMDFF::~QMDFF()
{
#ifdef USE_D3
delete m_d3;
#endif

delete m_threadpool;
for (int i = 0; i < m_stored_threads.size(); ++i)
delete m_stored_threads[i];
Expand Down Expand Up @@ -805,7 +810,9 @@ void QMDFF::setMolecule(const std::vector<int>& atom_types, const Matrix& geomet
m_atom_types = atom_types;
m_geometry = geometry;
m_gradient = Eigen::MatrixXd::Zero(m_atom_types.size(), 3);
#ifdef USE_D3
m_d3->InitialiseMolecule(m_atom_types);
#endif
}

void QMDFF::writeParameterFile(const std::string& file) const
Expand Down Expand Up @@ -1182,7 +1189,9 @@ double QMDFF::Calculate(bool grd, bool verbose)
m_h4correction.GradientHH()[i].x = 0;
m_h4correction.GradientHH()[i].y = 0;
m_h4correction.GradientHH()[i].z = 0;
#ifdef USE_D3
m_d3->UpdateAtom(i, m_geometry(i, 0), m_geometry(i, 1), m_geometry(i, 2));
#endif
}

double energy = 0.0;
Expand Down Expand Up @@ -1230,7 +1239,10 @@ double QMDFF::Calculate(bool grd, bool verbose)

if (grd) {
double grad[3 * m_atom_types.size()];
#ifdef USE_D3
d3_energy = m_d3->DFTD3Calculation(grad);
#endif

for (int i = 0; i < m_atom_types.size(); ++i) {
double val = grad[3 * i + 0] * au;
if (!std::isnan(val) && std::abs(val) < 1e10)
Expand All @@ -1243,7 +1255,10 @@ double QMDFF::Calculate(bool grd, bool verbose)
m_gradient(i, 2) += val;
}
} else
#ifdef USE_D3
d3_energy = m_d3->DFTD3Calculation(0);
#endif

energy += m_final_factor * m_h4_scaling * energy_h4 + m_final_factor * m_hh_scaling * energy_hh + d3_energy;
}
energy = bond_energy + angle_energy + dihedral_energy + inversion_energy; // + vdw_energy;
Expand Down
3 changes: 3 additions & 0 deletions src/core/qmdff.h
Original file line number Diff line number Diff line change
Expand Up @@ -344,5 +344,8 @@ class QMDFF {
std::vector<QMDFFThread*> m_stored_threads;
CxxThreadPool* m_threadpool;
hbonds4::H4Correction m_h4correction;

#ifdef USE_D3
DFTD3Interface* m_d3;
#endif
};

0 comments on commit dabc733

Please sign in to comment.