forked from deepmodeling/abacus-develop
-
Notifications
You must be signed in to change notification settings - Fork 0
Expand file tree
/
Copy pathread_stru.cpp
More file actions
92 lines (83 loc) · 2.36 KB
/
read_stru.cpp
File metadata and controls
92 lines (83 loc) · 2.36 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
#include "read_stru.h"
#include "source_base/timer.h"
#include "source_base/vector3.h"
#include "source_base/mathzone.h"
bool unitcell::check_tau(const Atom* atoms,
const int& ntype,
const double& lat0)
{
ModuleBase::TITLE("UnitCell","check_tau");
ModuleBase::timer::start("UnitCell","check_tau");
ModuleBase::Vector3<double> diff = 0.0;
double norm = 0.0;
double tolerence_bohr = 1.0e-3;
for(int T1=0; T1< ntype; T1++)
{
for(int I1=0; I1< atoms[T1].na; I1++)
{
double shortest_norm = 10000.0; // a large number
for(int T2=0; T2<ntype; T2++)
{
for(int I2=0; I2<atoms[T2].na; I2++)
{
if(T1==T2 && I1==I2)
{
shortest_norm = 0.0;
}
else
{
diff = atoms[T1].tau[I1] - atoms[T2].tau[I2];
norm = diff.norm() * lat0;
if( shortest_norm > norm )
{
shortest_norm = norm;
}
if( norm < tolerence_bohr ) // unit is Bohr
{
GlobalV::ofs_warning << " two atoms are too close!" << std::endl;
GlobalV::ofs_warning << " type:" << atoms[T1].label << " atom " << I1 + 1 << std::endl;
GlobalV::ofs_warning << " type:" << atoms[T2].label << " atom " << I2 + 1 << std::endl;
GlobalV::ofs_warning << " distance = " << norm << " Bohr" << std::endl;
ModuleBase::timer::end("UnitCell","check_tau");
return false;
}
}
}
}
}
}
ModuleBase::timer::end("UnitCell","check_tau");
return true;
}
void unitcell::check_dtau(Atom* atoms,
const int& ntype,
const double& lat0,
ModuleBase::Matrix3& latvec)
{
for(int it=0; it<ntype; it++)
{
Atom* atom1 = &atoms[it];
for(int ia=0; ia<atoms[it].na; ia++)
{
// mohan add 2011-04-07
// fmod(x,1.0) set the result between the [0,1.0),
// while the x may be the negtivate value,thus we add 10000.
atom1->taud[ia].x=fmod(atom1->taud[ia].x + 10000,1.0);
atom1->taud[ia].y=fmod(atom1->taud[ia].y + 10000,1.0);
atom1->taud[ia].z=fmod(atom1->taud[ia].z + 10000,1.0);
double cx2=0.0;
double cy2=0.0;
double cz2=0.0;
ModuleBase::Mathzone::Direct_to_Cartesian(
atom1->taud[ia].x, atom1->taud[ia].y, atom1->taud[ia].z,
latvec.e11, latvec.e12, latvec.e13,
latvec.e21, latvec.e22, latvec.e23,
latvec.e31, latvec.e32, latvec.e33,
cx2, cy2, cz2);
atom1->tau[ia].x = cx2;
atom1->tau[ia].y = cy2;
atom1->tau[ia].z = cz2;
}
}
return;
}