Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Consolidate some code #66

Merged
merged 1 commit into from
Oct 31, 2024
Merged
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
16 changes: 2 additions & 14 deletions src/CALPHADFreeEnergyFunctionsBinary.cc
Original file line number Diff line number Diff line change
Expand Up @@ -251,27 +251,15 @@ bool CALPHADFreeEnergyFunctionsBinary::computeCeqT(
//=======================================================================

void CALPHADFreeEnergyFunctionsBinary::computePhasesFreeEnergies(
const double temperature, const double* const hphi, const double conc,
const double temperature, const double* const phi, const double conc,
double& fl, double& fa)
{
// std::cout<<"CALPHADFreeEnergyFunctionsBinary::computePhasesFreeEnergies()"<<endl;

double c[2] = { conc, conc };

// evaluate temperature dependent parameters
CalphadDataType fA[2];
CalphadDataType fB[2];
int ret = computePhaseConcentrations(temperature, &conc, phi, c);

CalphadDataType Lmix_L[4];
CalphadDataType Lmix_A[4];
computeTdependentParameters(temperature, Lmix_L, Lmix_A, fA, fB);

double RT = GASCONSTANT_R_JPKPMOL * temperature;

CALPHADConcSolverBinary solver;
solver.setup(conc, hphi[0], 1. - hphi[0], RT, Lmix_L, Lmix_A, fA, fB);
int ret = solver.ComputeConcentration(
c, newton_tol_, newton_maxits_, newton_alpha_);
if (ret < 0)
{
#if 0
Expand Down
34 changes: 1 addition & 33 deletions src/CALPHADFreeEnergyFunctionsBinary2Ph1Sl.cc
Original file line number Diff line number Diff line change
Expand Up @@ -316,39 +316,7 @@ void CALPHADFreeEnergyFunctionsBinary2Ph1Sl::computePhasesFreeEnergies(

double c[2] = { conc, conc };

// evaluate temperature dependent parameters
CalphadDataType fA[2];
CalphadDataType fB[2];

CalphadDataType Lmix_L[4];
CalphadDataType Lmix_A[4];
computeTdependentParameters(temperature, Lmix_L, Lmix_A, fA, fB);

double RTinv = 1.0 / (GASCONSTANT_R_JPKPMOL * temperature);

// Get the sublattice stoichiometry
int p[2];
int q[2];
p[0] = sublattice_stoichiometry_phaseL_[0];
p[1] = sublattice_stoichiometry_phaseA_[0];
q[0] = sublattice_stoichiometry_phaseL_[1];
q[1] = sublattice_stoichiometry_phaseA_[1];

CALPHADConcSolverBinary2Ph1Sl solver;
solver.setup(conc, hphi[0], RTinv, Lmix_L, Lmix_A, fA, fB, p, q);
int ret = solver.ComputeConcentration(
c, newton_tol_, newton_maxits_, newton_alpha_);
if (ret < 0)
{
#if 0
std::cerr << "ERROR in "
"CALPHADFreeEnergyFunctionsBinary2Ph1Sl::"
"computePhasesFreeEnergies()"
" ---"
<< "conc=" << conc << ", hphi=" << hphi << std::endl;
abort();
#endif
}
int ret = computePhaseConcentrations(temperature, &conc, hphi, c);

// assert(c[0] >= 0.);
fl = computeFreeEnergy(temperature, &c[0], PhaseIndex::phaseL, false);
Expand Down
90 changes: 1 addition & 89 deletions src/CALPHADFreeEnergyFunctionsBinary3Ph2Sl.cc
Original file line number Diff line number Diff line change
Expand Up @@ -299,95 +299,7 @@ void CALPHADFreeEnergyFunctionsBinary3Ph2Sl::computePhasesFreeEnergies(

double c[3] = { conc, conc, conc };

// evaluate temperature dependent parameters
CalphadDataType fA[3];
CalphadDataType fB[3];

CalphadDataType Lmix_L[4];
CalphadDataType Lmix_A[4];
CalphadDataType Lmix_B[4];

computeTdependentParameters(temperature, Lmix_L, Lmix_A, Lmix_B, fA, fB);

double RTinv = 1.0 / (GASCONSTANT_R_JPKPMOL * temperature);

// Get the sublattice stoichiometry
int p[3];
int q[3];
p[0] = sublattice_stoichiometry_phaseL_[0];
p[1] = sublattice_stoichiometry_phaseA_[0];
p[2] = sublattice_stoichiometry_phaseB_[0];
q[0] = sublattice_stoichiometry_phaseL_[1];
q[1] = sublattice_stoichiometry_phaseA_[1];
q[2] = sublattice_stoichiometry_phaseB_[1];

CALPHADConcSolverBinary3Ph2Sl solver;
solver.setup(conc, hphi[0], hphi[1], hphi[2], RTinv, Lmix_L, Lmix_A, Lmix_B,
fA, fB, p, q);

// Loop that changes the initial conditions if the solver doesn't converge
int ret = -1;
int reset_index = 0;

#ifndef HAVE_OPENMP_OFFLOAD
std::random_device rd;
std::mt19937 prng(rd());
std::uniform_real_distribution<double> dist_generator(0.0, 1.0);
#endif

while (ret == -1)
{
ret = solver.ComputeConcentration(
c, newton_tol_, newton_maxits_, newton_alpha_);

#ifndef HAVE_OPENMP_OFFLOAD
if (ret == -1)
{
if (reset_index >= newton_max_resets_)
{
std::cout << "WARNING: Maximum number of restarts ("
<< newton_max_resets_
<< ") reached in "
"Newton solve without convergence."
<< std::endl;
break;
}
else
{
// Randomly pick sublattice occupations
double y0 = dist_generator(prng);
double y1 = dist_generator(prng);
double y2 = dist_generator(prng);

c[0] = (y0 + p[0]) / (p[0] + q[0]);
c[1] = (y1 + p[1]) / (p[1] + q[1]);
c[2] = (y0 + p[2]) / (p[2] + q[2]);

// std::cout << "New initial conditions: " << c[0] << " " <<
// c[1]
// << " " << c[2] << std::endl;
}
}
#else
break;
#endif

reset_index++;
}

if (ret < 0)
{
#ifndef HAVE_OPENMP_OFFLOAD
std::cerr << "ERROR in "
"CALPHADFreeEnergyFunctionsBinary3Ph2Sl::"
"computePhasesFreeEnergies()"
" ---"
<< conc << ", hphi0=" << hphi[0] << ", hphi1=" << hphi[1]
<< ", hphi2=" << hphi[2] << "x=" << c[0] << ", " << c[1]
<< ", " << c[2] << std::endl;
// abort();
#endif
}
int ret = computePhaseConcentrations(temperature, &conc, hphi, c);

// assert(c[0] >= 0.);
fl = computeFreeEnergy(temperature, &c[0], PhaseIndex::phaseL, false);
Expand Down
19 changes: 1 addition & 18 deletions src/CALPHADFreeEnergyFunctionsBinaryThreePhase.cc
Original file line number Diff line number Diff line change
Expand Up @@ -236,25 +236,8 @@ void CALPHADFreeEnergyFunctionsBinaryThreePhase::computePhasesFreeEnergies(

double c[3] = { conc, conc, conc };

// evaluate temperature dependent parameters
CalphadDataType fA[3];
CalphadDataType fB[3];

CalphadDataType Lmix_L[4];
CalphadDataType Lmix_A[4];
CalphadDataType Lmix_B[4];
int ret = computePhaseConcentrations(temperature, &conc, hphi, c);

computeTdependentParameters(temperature, Lmix_L, Lmix_A, Lmix_B, fA, fB);
assert(fA[0] == fA[0]);
assert(Lmix_L[0] == Lmix_L[0]);

double RT = GASCONSTANT_R_JPKPMOL * temperature;

CALPHADConcSolverBinaryThreePhase solver;
solver.setup(
conc, hphi[0], hphi[1], hphi[2], RT, Lmix_L, Lmix_A, Lmix_B, fA, fB);
int ret = solver.ComputeConcentration(
c, newton_tol_, newton_maxits_, newton_alpha_);
if (ret < 0)
{
#if 0
Expand Down
40 changes: 2 additions & 38 deletions src/CALPHADFreeEnergyFunctionsTernary.cc
Original file line number Diff line number Diff line change
Expand Up @@ -513,45 +513,9 @@ void CALPHADFreeEnergyFunctionsTernary::computePhasesFreeEnergies(

double cauxilliary[4] = { conc0, conc1, conc0, conc1 };

CalphadDataType L_AB_L[4];
CalphadDataType L_AC_L[4];
CalphadDataType L_BC_L[4];

CalphadDataType L_ABC_L[3];

CalphadDataType L_AB_S[4];
CalphadDataType L_AC_S[4];
CalphadDataType L_BC_S[4];

CalphadDataType L_ABC_S[3];
double conc[2] = { conc0, conc1 };
computePhaseConcentrations(temperature, conc, hphi, cauxilliary);

CalphadDataType fA[2];
CalphadDataType fB[2];
CalphadDataType fC[2];

computeTdependentParameters(temperature, L_AB_L, L_AC_L, L_BC_L, L_ABC_L,
L_AB_S, L_AC_S, L_BC_S, L_ABC_S, fA, fB, fC);

double RTinv = 1.0 / (GASCONSTANT_R_JPKPMOL * temperature);
CALPHADConcSolverTernary solver;
solver.setup(conc0, conc1, hphi[0], RTinv, L_AB_L, L_AC_L, L_BC_L, L_AB_S,
L_AC_S, L_BC_S, L_ABC_L, L_ABC_S, fA, fB, fC);
int ret
= solver.ComputeConcentration(cauxilliary, newton_tol_, newton_maxits_);
#ifndef HAVE_OPENMP_OFFLOAD
if (ret < 0)
{
std::cerr << "ERROR in "
"CALPHADFreeEnergyFunctionsTernary::"
"computePhasesFreeEnergies() "
"---"
<< "conc0=" << conc0 << ", conc1=" << conc1
<< ", hphi=" << hphi[0] << std::endl;
abort();
}

assert(conc0 >= 0.);
#endif
double concl[2] = { cauxilliary[0], cauxilliary[1] };
fl = computeFreeEnergy(temperature, &concl[0], PhaseIndex::phaseL, false);

Expand Down
Loading