Skip to content
Closed
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
15 changes: 13 additions & 2 deletions source/source_cell/read_atoms_helper.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -391,7 +391,7 @@ bool parse_atom_properties(std::ifstream& ifpos,
else if ( tmpid == "mag" || tmpid == "magmom")
{
set_element_mag_zero = true;
double tmpamg=0;
double tmpamg=0.0;
ifpos >> tmpamg;
tmp=ifpos.get();
while (tmp==' ')
Expand All @@ -409,12 +409,17 @@ bool parse_atom_properties(std::ifstream& ifpos,
+pow(atom.m_loc_[ia].z,2));
input_vec_mag=true;


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

// print mag
std::cout << " read mag for atom " << ia+1
<< ": " << atom.mag[ia] << std::endl;
}
else if ( tmpid == "angle1")
{
Expand All @@ -432,7 +437,7 @@ bool parse_atom_properties(std::ifstream& ifpos,
}
else if ( tmpid == "lambda")
{
double tmplam=0;
double tmplam=0.0;
ifpos >> tmplam;
tmp=ifpos.get();
while (tmp==' ')
Expand All @@ -444,6 +449,12 @@ bool parse_atom_properties(std::ifstream& ifpos,
ifpos.putback(tmp);
ifpos >> atom.lambda[ia].y>>atom.lambda[ia].z;
atom.lambda[ia].x=tmplam;

// print lambda
std::cout << " read lambda for atom " << ia+1
<< ": " << atom.lambda[ia].x
<< " " << atom.lambda[ia].y
<< " " << atom.lambda[ia].z << std::endl;
}
else
{
Expand Down
13 changes: 9 additions & 4 deletions source/source_lcao/module_deltaspin/basic_funcs.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -142,8 +142,13 @@ void print_2d(const std::string info, const std::vector<ModuleBase::Vector3<doub
for (const auto &row : array)
{
iat += 1;
if (nspin == 2) { ofs << FmtCore::format("ATOM %6d %20.10f\n", iat, row.z*unit_convert);
} else if (nspin == 4) { ofs << FmtCore::format("ATOM %6d %20.10f %20.10f %20.10f\n", iat, row.x*unit_convert, row.y*unit_convert, row.z*unit_convert);
}
if (nspin == 2)
{
ofs << FmtCore::format(" atom %6d %20.10f\n", iat, row.z*unit_convert);
}
else if (nspin == 4)
{
ofs << FmtCore::format(" atom %6d %20.10f %20.10f %20.10f\n", iat, row.x*unit_convert, row.y*unit_convert, row.z*unit_convert);
}
}
}
}
28 changes: 18 additions & 10 deletions source/source_lcao/module_deltaspin/lambda_loop.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -100,24 +100,32 @@ void spinconstrain::SpinConstrain<std::complex<double>>::run_lambda_loop(int out

template <>
void spinconstrain::SpinConstrain<std::complex<double>>::run_lambda_loop(
int outer_step,
int outer_step,
bool rerun)
{
// init controlling parameters
int nat = this->get_nat();
int ntype = this->get_ntype();
const int nat = this->get_nat();
const int ntype = this->get_ntype();

std::vector<ModuleBase::Vector3<double>> initial_lambda(nat,0.0);
std::vector<ModuleBase::Vector3<double>> delta_lambda(nat,0.0);

// set nu, dnu and dnu_last_step
std::vector<ModuleBase::Vector3<double>> dnu(nat, 0.0), dnu_last_step(nat, 0.0);

// two controlling temp variables
std::vector<ModuleBase::Vector3<double>> temp_1(nat, 0.0);
std::vector<ModuleBase::Vector3<double>> spin(nat, 0.0), delta_spin(nat, 0.0);
std::vector<ModuleBase::Vector3<double>> search(nat, 0.0), search_old(nat, 0.0);
std::vector<ModuleBase::Vector3<double>> new_spin(nat, 0.0), spin_plus(nat, 0.0);

double alpha_opt, alpha_plus;
double beta = 0.0, g = 0.0, mean_error = 0.0, mean_error_old = 0.0, rms_error = 0.0;
double alpha_opt = 0.0;
double alpha_plus = 0.0;
double beta = 0.0;
double g = 0.0;
double mean_error = 0.0;
double mean_error_old = 0.0;
double rms_error = 0.0;

double alpha_trial = this->alpha_trial_;

Expand All @@ -143,9 +151,9 @@ void spinconstrain::SpinConstrain<std::complex<double>>::run_lambda_loop(
this->cal_mw_from_lambda(i_step);
spin = this->Mi_;
where_fill_scalar_else_2d(this->constrain_, 0, zero, this->lambda_, initial_lambda);
print_2d("initial lambda (eV/uB): ", initial_lambda, this->nspin_, ModuleBase::Ry_to_eV);
print_2d("initial spin (uB): ", spin, this->nspin_);
print_2d("target spin (uB): ", this->target_mag_, this->nspin_);
print_2d(" initial lambda (eV/uB): ", initial_lambda, this->nspin_, ModuleBase::Ry_to_eV);
print_2d(" initial mag (uB): ", spin, this->nspin_);
print_2d(" target mag (uB): ", this->target_mag_, this->nspin_);
i_step++;
}
else
Expand Down Expand Up @@ -221,10 +229,10 @@ void spinconstrain::SpinConstrain<std::complex<double>>::run_lambda_loop(
}
mean_error = sum_2d(temp_1) / nat;
rms_error = std::sqrt(mean_error);
std::cout<<"Current RMS: "<<rms_error<<std::endl;
std::cout<<" current RMS "<<rms_error<<std::endl;
if(rms_error > this->current_sc_thr_ * 10 && rerun == true && this->higher_mag_prec == true)
{
std::cout<<"Error: RMS error is too large, rerun the loop"<<std::endl;
std::cout<<" Error: RMS error is too large, rerun the loop"<<std::endl;
this->run_lambda_loop(outer_step, false);
}
}
Expand Down
33 changes: 17 additions & 16 deletions source/source_lcao/module_deltaspin/lambda_loop_helper.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -4,31 +4,32 @@
template <>
void spinconstrain::SpinConstrain<std::complex<double>>::print_termination()
{
print_2d("after-optimization spin (uB): (print in the inner loop): ", this->Mi_, this->nspin_);
print_2d("after-optimization lambda (eV/uB): (print in the inner loop): ", this->lambda_, this->nspin_, ModuleBase::Ry_to_eV);
std::cout << "Inner optimization for lambda ends." << std::endl;
std::cout << "===============================================================================" << std::endl;
print_2d(" mag (uB): ", this->Mi_, this->nspin_);
print_2d(" lambda (eV/uB): ", this->lambda_, this->nspin_, ModuleBase::Ry_to_eV);
std::cout << " Inner optimization for lambda ends." << std::endl;
std::cout << " ----------------------------------------------------------------\n" << std::endl;
}

template <>
bool spinconstrain::SpinConstrain<std::complex<double>>::check_rms_stop(int outer_step,
int i_step,
double rms_error,
double duration,
double total_duration)
int i_step,
double rms_error,
double duration,
double total_duration)
{
std::cout << "Step (Outer -- Inner) = " << outer_step << " -- " << std::left << std::setw(5) << i_step + 1
<< " RMS = " << rms_error << " TIME(s) = " << std::setw(11) << duration << std::endl;
// std::cout << "Step (Outer -- Inner) = " << outer_step << " -- " << std::left << std::setw(5)
std::cout << " step "<< i_step + 1
<< " RMS " << rms_error << " TIME(s) " << std::setw(11) << duration << std::endl;
if (rms_error < this->current_sc_thr_ || i_step == this->nsc_ - 1)
{
if (rms_error < this->current_sc_thr_)
{
std::cout << "Meet convergence criterion ( < " << this->current_sc_thr_ << " ), exit.";
std::cout << " Meet convergence criterion ( < " << this->current_sc_thr_ << " ), exit.";
std::cout << " Total TIME(s) = " << total_duration << std::endl;
}
else if (i_step == this->nsc_ - 1)
{
std::cout << "Reach maximum number of steps ( " << this->nsc_ << " ), exit.";
std::cout << " Reach maximum number of steps ( " << this->nsc_ << " ), exit.";
std::cout << " Total TIME(s) = " << total_duration << std::endl;
}
this->print_termination();
Expand All @@ -41,9 +42,9 @@ bool spinconstrain::SpinConstrain<std::complex<double>>::check_rms_stop(int oute
template <>
void spinconstrain::SpinConstrain<std::complex<double>>::print_header()
{
std::cout << "===============================================================================" << std::endl;
std::cout << "Inner optimization for lambda begins ..." << std::endl;
std::cout << "Covergence criterion for the iteration: " << this->sc_thr_ << std::endl;
std::cout << " ----------------------------------------------------------------" << std::endl;
std::cout << " Inner optimization for lambda begins ..." << std::endl;
std::cout << " Covergence threshold " << this->sc_thr_ << std::endl;
}

/// check restriction
Expand Down Expand Up @@ -180,4 +181,4 @@ bool spinconstrain::SpinConstrain<std::complex<double>>::check_gradient_decay(
}
}
return false;
}
}
31 changes: 30 additions & 1 deletion source/source_lcao/module_deltaspin/spin_constrain.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -197,6 +197,8 @@ void SpinConstrain<TK>::set_sc_lambda()
template <typename TK>
void SpinConstrain<TK>::set_target_mag()
{
ModuleBase::TITLE("SpinConstrain", "set_target_mag_1");

this->check_atomCounts();
int nat = this->get_nat();
this->target_mag_.resize(nat, 0.0);
Expand All @@ -208,6 +210,9 @@ void SpinConstrain<TK>::set_target_mag()
int index = element_data.index;
int iat = this->get_iat(itype, index);
ModuleBase::Vector3<double> mag(0.0, 0.0, 0.0);

std::cout << " Here, mag_type is " << element_data.mag_type << std::endl;

if (element_data.mag_type == 0)
{
mag.x = element_data.target_mag[0];
Expand All @@ -222,11 +227,17 @@ void SpinConstrain<TK>::set_target_mag()
mag.y = element_data.target_mag_val * std::sin(radian_angle1) * std::sin(radian_angle2);
mag.z = element_data.target_mag_val * std::cos(radian_angle1);
if (std::abs(mag.x) < 1e-14)
{
mag.x = 0.0;
}
if (std::abs(mag.y) < 1e-14)
{
mag.y = 0.0;
}
if (std::abs(mag.z) < 1e-14)
{
mag.z = 0.0;
}
}
this->target_mag_[iat] = mag;
}
Expand Down Expand Up @@ -286,6 +297,10 @@ void SpinConstrain<TK>::set_sc_lambda(const ModuleBase::Vector3<double>* lambda_
template <typename TK>
void SpinConstrain<TK>::set_target_mag(const ModuleBase::Vector3<double>* target_mag_in, int nat_in)
{
ModuleBase::TITLE("SpinConstrain", "set_target_mag_2");

std::cout << " set_target_mag_2 " << std::endl;

this->check_atomCounts();
int nat = this->get_nat();
if (nat_in != nat)
Expand All @@ -302,8 +317,22 @@ void SpinConstrain<TK>::set_target_mag(const ModuleBase::Vector3<double>* target
template <typename TK>
void SpinConstrain<TK>::set_target_mag(const std::vector<ModuleBase::Vector3<double>>& target_mag_in)
{
int nat = this->get_nat();
ModuleBase::TITLE("SpinConstrain", "set_target_mag_3");

std::cout << " set_target_mag_3 " << std::endl;

const int nat = this->get_nat();
assert(target_mag_in.size() == nat);

std::cout << " target_mag_in size: " << nat << std::endl;
for (int iat = 0; iat < nat; iat++)
{
std::cout << " atom=" << iat+1
<< ": " << target_mag_in[iat].x
<< " " << target_mag_in[iat].y
<< " " << target_mag_in[iat].z << std::endl;
}

if (this->nspin_ == 2)
{
this->target_mag_.resize(nat, 0.0);
Expand Down
Loading