Skip to content

Commit

Permalink
Merge branch 'develop' into warning
Browse files Browse the repository at this point in the history
  • Loading branch information
hongriTianqi committed Jan 18, 2024
2 parents da261ad + 93c9ce3 commit 8be30ee
Show file tree
Hide file tree
Showing 15 changed files with 445 additions and 374 deletions.
16 changes: 13 additions & 3 deletions .github/workflows/test.yml
Original file line number Diff line number Diff line change
Expand Up @@ -2,27 +2,37 @@ name: Integration Test and Unit Test

on:
pull_request:

concurrency:
group: ${{ github.workflow }}-${{ github.ref }}
cancel-in-progress: true

jobs:
test:
name: Test
runs-on: self-hosted
if: github.repository_owner == 'deepmodeling'
container: ghcr.io/deepmodeling/abacus-gnu
container:
image: ghcr.io/deepmodeling/abacus-gnu
volumes:
- /tmp/ccache:/github/home/.ccache
steps:
- name: Checkout
uses: actions/checkout@v4
with:
submodules: recursive

- name: Install Ccache
run: |
sudo apt-get update
sudo apt-get install -y ccache
- name: Build
run: |
cmake -B build -DBUILD_TESTING=ON -DENABLE_DEEPKS=ON -DENABLE_LIBXC=ON -DENABLE_LIBRI=ON -DENABLE_PAW=ON
cmake --build build -j8
cmake --install build
- name: Test
env:
GTEST_COLOR: 'yes'
Expand Down
49 changes: 34 additions & 15 deletions source/module_base/math_sphbes.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -808,7 +808,7 @@ void Sphbes::dsphbesj(const int n,
}
}

void Sphbes::sphbes_zeros(const int l, const int n, double* const zeros)
void Sphbes::sphbes_zeros(const int l, const int n, double* const zeros, const bool return_all)
{
assert( n > 0 );
assert( l >= 0 );
Expand All @@ -818,10 +818,22 @@ void Sphbes::sphbes_zeros(const int l, const int n, double* const zeros)
// This property enables us to use bracketing method recursively
// to find all zeros of j_l from the zeros of j_0.

// if l is odd , j_0 --> j_1 --> j_3 --> j_5 --> ...
// if l is even, j_0 --> j_2 --> j_4 --> j_6 --> ...

int nz = n + (l+1)/2; // number of effective zeros in buffer
// If return_all is true, zeros of j_0, j_1, ..., j_l will all be returned
// such that zeros[l*n+i] is the i-th zero of j_l. As such, it is required
// that the array "zeros" has a size of (l+1)*n.
//
// If return_all is false, only the zeros of j_l will be returned
// and "zeros" is merely required to have a size of n.
// Note that in this case the bracketing method can be applied with a stride
// of 2 instead of 1:
// j_0 --> j_1 --> j_3 --> j_5 --> ... --> j_l (odd l)
// j_0 --> j_2 --> j_4 --> j_6 --> ... --> j_l (even l)

// Every recursion step reduces the number of zeros by 1.
// If return_all is true, one needs to start with n+l zeros of j_0
// to ensure n zeros of j_l; otherwise with a stride of 2 one only
// needs to start with n+(l+1)/2 zeros of j_0
int nz = n + ( return_all ? l : (l+1)/2 );
double* buffer = new double[nz];

// zeros of j_0 = sin(x)/x is just n*pi
Expand All @@ -831,27 +843,34 @@ void Sphbes::sphbes_zeros(const int l, const int n, double* const zeros)
buffer[i] = (i+1) * PI;
}

int ll = 1;
int ll; // active l
auto jl = [&ll] (double x) { return sphbesj(ll, x); };

if (l % 2 == 1)
int stride;
std::function<void()> copy_if_needed;
int offset = 0; // keeps track of the position in zeros for next copy (used when return_all == true)
if (return_all)
{
for (int i = 0; i < nz-1; i++)
{
buffer[i] = illinois(jl, buffer[i], buffer[i+1], 1e-15, 50);
}
--nz;
copy_if_needed = [&](){ std::copy(buffer, buffer + n, zeros + offset); offset += n; };
stride = 1;
ll = 1;
}
else
{
copy_if_needed = [](){};
stride = 2;
ll = 2 - l % 2;
}

for (ll = 2 + l%2; ll <= l; ll += 2, --nz)
for (; ll <= l; ll += stride, --nz)
{
copy_if_needed();
for (int i = 0; i < nz-1; i++)
{
buffer[i] = illinois(jl, buffer[i], buffer[i+1], 1e-15, 50);
}
}

std::copy(buffer, buffer + n, zeros);
std::copy(buffer, buffer + n, zeros + offset);
delete[] buffer;
}

Expand Down
13 changes: 9 additions & 4 deletions source/module_base/math_sphbes.h
Original file line number Diff line number Diff line change
Expand Up @@ -126,13 +126,18 @@ class Sphbes
* This function computes the first n positive zeros of the l-th order
* spherical Bessel function of the first kind.
*
* @param[in] l order of the spherical Bessel function
* @param[in] n number of zeros to be computed
* @param[out] zeros on exit, contains the first n positive zeros in ascending order
* @param[in] l (maximum) order of the spherical Bessel function
* @param[in] n number of zeros to be computed (for each j_l if return_all is true)
* @param[out] zeros on exit, contains the positive zeros.
* @param[in] return_all if true, return all zeros from j_0 to j_l such that zeros[l*n+i]
* is the i-th zero of j_l. If false, return only the first n zeros of j_l.
*
* @note The size of array "zeros" must be at least (l+1)*n if return_all is true, and n otherwise.
*/
static void sphbes_zeros(const int l,
const int n,
double* const zeros
double* const zeros,
bool return_all = false
);

private:
Expand Down
16 changes: 14 additions & 2 deletions source/module_base/test/math_sphbes_test.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -352,15 +352,27 @@ TEST_F(Sphbes, Zeros)

int lmax = 20;
int nzeros = 500;
double* zeros = new double[nzeros];
double* zeros = new double[nzeros*(lmax+1)];
for (int l = 0; l <= lmax; ++l)
{
ModuleBase::Sphbes::sphbes_zeros(l, nzeros, zeros);
ModuleBase::Sphbes::sphbes_zeros(l, nzeros, zeros, false);
for (int i = 0; i < nzeros; ++i)
{
EXPECT_LT(std::abs(ModuleBase::Sphbes::sphbesj(l, zeros[i])), 1e-14);
}
}


ModuleBase::Sphbes::sphbes_zeros(lmax, nzeros, zeros, true);
for (int l = 0; l <= lmax; ++l)
{
for (int i = 0; i < nzeros; ++i)
{
EXPECT_LT(std::abs(ModuleBase::Sphbes::sphbesj(l, zeros[l*nzeros+i])), 1e-14);
}
}

delete[] zeros;
}

TEST_F(Sphbes, ZerosOld)
Expand Down
189 changes: 95 additions & 94 deletions source/module_cell/read_atoms.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -535,100 +535,101 @@ bool UnitCell::read_atom_positions(std::ifstream &ifpos, std::ofstream &ofs_runn
ModuleBase::GlobalFunc::ZEROS(atoms[it].mag,na);
for (int ia = 0;ia < na; ia++)
{
// modify the reading of frozen ions and velocities -- Yuanbo Li 2021/8/20
ifpos >> v.x >> v.y >> v.z;
mv.x = true ;
mv.y = true ;
mv.z = true ;
atoms[it].vel[ia].set(0,0,0);
atoms[it].mag[ia]=magnet.start_magnetization[it];//if this line is used, default startmag_type would be 2
atoms[it].angle1[ia]=0;
atoms[it].angle2[ia]=0;
atoms[it].m_loc_[ia].set(0,0,0);

std::string tmpid;
tmpid = ifpos.get();

if( (int)tmpid[0] < 0 )
{
std::cout << "read_atom_positions, mismatch in atom number for atom type: " << atoms[it].label << std::endl;
exit(1);
}

bool input_vec_mag=false;
bool input_angle_mag=false;
while ( (tmpid != "\n") && (ifpos.eof()==false) && (tmpid !="#") )
{
tmpid = ifpos.get() ;
// old method of reading frozen ions
char tmp = (char)tmpid[0];
if ( tmp >= 48 && tmp <= 57 )
{
mv.x = std::stoi(tmpid);
ifpos >> mv.y >> mv.z ;
}
// new method of reading frozen ions and velocities
if ( tmp >= 'a' && tmp <='z')
{
ifpos.putback(tmp);
ifpos >> tmpid;
}
if ( tmpid == "m" )
{
ifpos >> mv.x >> mv.y >> mv.z ;
}
else if ( tmpid == "v" ||tmpid == "vel" || tmpid == "velocity" )
{
ifpos >> atoms[it].vel[ia].x >> atoms[it].vel[ia].y >> atoms[it].vel[ia].z;
}
else if ( tmpid == "mag" || tmpid == "magmom")
{
set_element_mag_zero = true;
double tmpamg=0;
ifpos >> tmpamg;
tmp=ifpos.get();
while (tmp==' ')
{
tmp=ifpos.get();
}

if((tmp >= 48 && tmp <= 57) or tmp=='-')
{
ifpos.putback(tmp);
ifpos >> atoms[it].m_loc_[ia].y>>atoms[it].m_loc_[ia].z;
atoms[it].m_loc_[ia].x=tmpamg;
atoms[it].mag[ia]=sqrt(pow(atoms[it].m_loc_[ia].x,2)+pow(atoms[it].m_loc_[ia].y,2)+pow(atoms[it].m_loc_[ia].z,2));
input_vec_mag=true;

}
else
{
ifpos.putback(tmp);
atoms[it].mag[ia]=tmpamg;
}

// atoms[it].mag[ia];
}
else if ( tmpid == "angle1")
{
ifpos >> atoms[it].angle1[ia];
atoms[it].angle1[ia]=atoms[it].angle1[ia]/180 *ModuleBase::PI;
input_angle_mag=true;
set_element_mag_zero = true;
}
else if ( tmpid == "angle2")
{
ifpos >> atoms[it].angle2[ia];
atoms[it].angle2[ia]=atoms[it].angle2[ia]/180 *ModuleBase::PI;
input_angle_mag=true;
set_element_mag_zero = true;
}

}
while ( (tmpid != "\n") && (ifpos.eof()==false) )
{
tmpid = ifpos.get();
}
// modify the reading of frozen ions and velocities -- Yuanbo Li 2021/8/20
ifpos >> v.x >> v.y >> v.z;
mv.x = true ;
mv.y = true ;
mv.z = true ;
atoms[it].vel[ia].set(0,0,0);
atoms[it].mag[ia]=magnet.start_magnetization[it];//if this line is used, default startmag_type would be 2
atoms[it].angle1[ia]=0;
atoms[it].angle2[ia]=0;
atoms[it].m_loc_[ia].set(0,0,0);

std::string tmpid;
tmpid = ifpos.get();

if( (int)tmpid[0] < 0 )
{
std::cout << "read_atom_positions, mismatch in atom number for atom type: " << atoms[it].label << std::endl;
exit(1);
}

bool input_vec_mag=false;
bool input_angle_mag=false;
// read if catch goodbit before "\n" and "#"
while ( (tmpid != "\n") && (ifpos.good()) && (tmpid !="#") )
{
tmpid = ifpos.get() ;
// old method of reading frozen ions
char tmp = (char)tmpid[0];
if ( tmp >= 48 && tmp <= 57 )
{
mv.x = std::stoi(tmpid);
ifpos >> mv.y >> mv.z ;
}
// new method of reading frozen ions and velocities
if ( tmp >= 'a' && tmp <='z')
{
ifpos.putback(tmp);
ifpos >> tmpid;
}
if ( tmpid == "m" )
{
ifpos >> mv.x >> mv.y >> mv.z ;
}
else if ( tmpid == "v" ||tmpid == "vel" || tmpid == "velocity" )
{
ifpos >> atoms[it].vel[ia].x >> atoms[it].vel[ia].y >> atoms[it].vel[ia].z;
}
else if ( tmpid == "mag" || tmpid == "magmom")
{
set_element_mag_zero = true;
double tmpamg=0;
ifpos >> tmpamg;
tmp=ifpos.get();
while (tmp==' ')
{
tmp=ifpos.get();
}

if((tmp >= 48 && tmp <= 57) or tmp=='-')
{
ifpos.putback(tmp);
ifpos >> atoms[it].m_loc_[ia].y>>atoms[it].m_loc_[ia].z;
atoms[it].m_loc_[ia].x=tmpamg;
atoms[it].mag[ia]=sqrt(pow(atoms[it].m_loc_[ia].x,2)+pow(atoms[it].m_loc_[ia].y,2)+pow(atoms[it].m_loc_[ia].z,2));
input_vec_mag=true;

}
else
{
ifpos.putback(tmp);
atoms[it].mag[ia]=tmpamg;
}

// atoms[it].mag[ia];
}
else if ( tmpid == "angle1")
{
ifpos >> atoms[it].angle1[ia];
atoms[it].angle1[ia]=atoms[it].angle1[ia]/180 *ModuleBase::PI;
input_angle_mag=true;
set_element_mag_zero = true;
}
else if ( tmpid == "angle2")
{
ifpos >> atoms[it].angle2[ia];
atoms[it].angle2[ia]=atoms[it].angle2[ia]/180 *ModuleBase::PI;
input_angle_mag=true;
set_element_mag_zero = true;
}
}
// move to next line
while ( (tmpid != "\n") && (ifpos.good()) )
{
tmpid = ifpos.get();
}
std::string mags;
//cout<<"mag"<<atoms[it].mag[ia]<<"angle1"<<atoms[it].angle1[ia]<<"angle2"<<atoms[it].angle2[ia]<<'\n';

Expand Down
2 changes: 1 addition & 1 deletion source/module_hamilt_lcao/module_deltaspin/cal_mw.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -33,7 +33,7 @@ ModuleBase::matrix SpinConstrain<std::complex<double>, psi::DEVICE_CPU>::cal_MW_
const char N_char = 'N';
const int one_int = 1;
const std::complex<double> one_float = {1.0, 0.0}, zero_float = {0.0, 0.0};
pzgemm_(&T_char,
pzgemm_(&N_char,
&T_char,
&nw,
&nw,
Expand Down
4 changes: 2 additions & 2 deletions source/module_io/mulliken_charge.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -44,7 +44,7 @@ ModuleBase::matrix ModuleIO::cal_mulliken(const std::vector<std::vector<double>>
const char N_char = 'N';
const int one_int = 1;
const double one_float = 1.0, zero_float = 0.0;
pdgemm_(&T_char,
pdgemm_(&N_char,
&T_char,
&GlobalV::NLOCAL,
&GlobalV::NLOCAL,
Expand Down Expand Up @@ -156,7 +156,7 @@ ModuleBase::matrix ModuleIO::cal_mulliken(const std::vector<std::vector<std::com
const char N_char = 'N';
const int one_int = 1;
const std::complex<double> one_float = {1.0, 0.0}, zero_float = {0.0, 0.0};
pzgemm_(&T_char,
pzgemm_(&N_char,
&T_char,
&GlobalV::NLOCAL,
&GlobalV::NLOCAL,
Expand Down
Loading

0 comments on commit 8be30ee

Please sign in to comment.