Skip to content

Fix: issues caused by workRange #101

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.

Sign up for GitHub

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

Merged
HPDell merged 5 commits into from
Feb 24, 2025
Merged

Fix: issues caused by workRange #101

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
2ba66f8
fix: min distance remove 0
HPDell Feb 23, 2025
ede8c62
edit: make mWorkRange as pair
HPDell Feb 23, 2025
6839162
edit: move gwr fit process to fitInitial
HPDell Feb 23, 2025
a806a80
edit: handle min and max distance
HPDell Feb 23, 2025
d745453
fix: init mWorkRange in GWRRobust
HPDell Feb 24, 2025
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
2 changes: 1 addition & 1 deletion include/gwmodelpp/GWRBasic.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -763,7 +763,7 @@ class GWRBasic : public GWRBase, public IBandwidthSelectable, public IVarialbeSe
int mWorkerId = 0;
int mWorkerNum = 1;
arma::uword mWorkRangeSize = 0;
std::optional<std::pair<arma::uword, arma::uword>> mWorkRange;
std::pair<arma::uword, arma::uword> mWorkRange = std::make_pair(arma::uword(0), arma::uword(0));

arma::mat mBetasSE; //!< \~english Standard errors of coefficient estimates. \~chinese 回归系数估计值的标准差。
arma::vec mSHat; //!< \~english A vector of \f$tr(S)\f$ and \f$tr(SS^T)\f$. \~chinese 由 \f$tr(S)\f$ 和 \f$tr(SS^T)\f$ 组成的向量。
Expand Down
4 changes: 3 additions & 1 deletion include/gwmodelpp/GWRMultiscale.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -645,6 +645,8 @@ class GWRMultiscale : public SpatialMultiscaleAlgorithm,
return mSpatialWeights[i].weight<BandwidthWeight>();
}

arma::mat fitInitial();

/**
* \~english
* @brief The backfitting algorithm.
Expand Down Expand Up @@ -817,7 +819,7 @@ class GWRMultiscale : public SpatialMultiscaleAlgorithm,
#ifdef ENABLE_MPI
arma::uword mWorkRangeSize = 0;
#endif // ENABLE_MPI
std::optional<std::pair<arma::uword, arma::uword>> mWorkRange;
std::pair<arma::uword, arma::uword> mWorkRange = std::make_pair(arma::uword(0), arma::uword(0));

public:
static int treeChildCount; //!< \~english \~chinese
Expand Down
63 changes: 27 additions & 36 deletions src/gwmodelpp/GWRBasic.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -63,6 +63,8 @@ mat GWRBasic::fit()
{
GWM_LOG_STAGE("Initializing");
uword nDp = mCoords.n_rows, nVars = mX.n_cols;
mWorkRange = make_pair(uword(0), nDp);

#ifdef ENABLE_MPI
// Sync x, y, coords with other processes
if (mParallelType & ParallelType::MPI)
Expand Down Expand Up @@ -250,12 +252,11 @@ mat GWRBasic::fitCoreSerial(const mat &x, const vec &y, const SpatialWeight &sw)
mBetasSE = mat(nVar, nDp, fill::zeros);
mSHat = vec(2, fill::zeros);
mQDiag = vec(nDp, fill::zeros);
std::pair<uword, uword> workRange = mWorkRange.value_or(make_pair(0, nDp));
uword rangeSize = workRange.second - workRange.first;
uword rangeSize = mWorkRange.second - mWorkRange.first;
mS = mat(isStoreS() ? rangeSize : 1, nDp, fill::zeros);
mC = cube(nVar, nDp, isStoreC() ? rangeSize : 1, fill::zeros);
// cout << mWorkerId << " process work range: [" << workRange.first << "," << workRange.second << "]\n";
for (uword i = workRange.first; i < workRange.second; i++)
// cout << mWorkerId << " process work range: [" << mWorkRange.first << "," << mWorkRange.second << "]\n";
for (uword i = mWorkRange.first; i < mWorkRange.second; i++)
{
GWM_LOG_STOP_BREAK(mStatus);
vec w = mSpatialWeight.weightVector(i);
Expand All @@ -274,8 +275,8 @@ mat GWRBasic::fitCoreSerial(const mat &x, const vec &y, const SpatialWeight &sw)
vec p = - si.t();
p(i) += 1.0;
mQDiag += p % p;
mS.row(isStoreS() ? (i - workRange.first) : 0) = si;
mC.slice(isStoreC() ? (i - workRange.first) : 0) = ci;
mS.row(isStoreS() ? (i - mWorkRange.first) : 0) = si;
mC.slice(isStoreC() ? (i - mWorkRange.first) : 0) = ci;
}
catch (const exception& e)
{
Expand All @@ -292,8 +293,7 @@ mat GWRBasic::fitCoreCVSerial(const mat& x, const vec& y, const SpatialWeight& s
{
uword nDp = mCoords.n_rows, nVar = x.n_cols;
mat betas(nVar, nDp, fill::zeros);
std::pair<uword, uword> workRange = mWorkRange.value_or(make_pair(0, nDp));
for (uword i = workRange.first; i < workRange.second; i++)
for (uword i = mWorkRange.first; i < mWorkRange.second; i++)
{
GWM_LOG_STOP_BREAK(mStatus);
vec w = sw.weightVector(i);
Expand All @@ -320,8 +320,7 @@ mat GWRBasic::fitCoreSHatSerial(const mat& x, const vec& y, const SpatialWeight&
uword nDp = mCoords.n_rows, nVar = x.n_cols;
mat betas(nVar, nDp, fill::zeros);
shat = vec(2, arma::fill::zeros);
std::pair<uword, uword> workRange = mWorkRange.value_or(make_pair(0, nDp));
for (uword i = workRange.first; i < workRange.second; i++)
for (uword i = mWorkRange.first; i < mWorkRange.second; i++)
{
GWM_LOG_STOP_BREAK(mStatus);
vec w = sw.weightVector(i);
Expand Down Expand Up @@ -460,16 +459,15 @@ mat GWRBasic::fitCoreOmp(const mat& x, const vec& y, const SpatialWeight& sw)
uword nDp = mCoords.n_rows, nVar = x.n_cols;
mat betas(nVar, nDp, fill::zeros);
mBetasSE = mat(nVar, nDp, fill::zeros);
std::pair<uword, uword> workRange = mWorkRange.value_or(make_pair(0, nDp));
uword rangeSize = workRange.second - workRange.first;
uword rangeSize = mWorkRange.second - mWorkRange.first;
mS = mat(isStoreS() ? rangeSize : 1, nDp, fill::zeros);
mC = cube(nVar, nDp, isStoreC() ? rangeSize : 1, fill::zeros);
mat shat_all(2, mOmpThreadNum, fill::zeros);
mat qDiag_all(nDp, mOmpThreadNum, fill::zeros);
bool success = true;
std::exception except;
#pragma omp parallel for num_threads(mOmpThreadNum)
for (uword i = workRange.first; i < workRange.second; i++)
for (uword i = mWorkRange.first; i < mWorkRange.second; i++)
{
GWM_LOG_STOP_CONTINUE(mStatus);
if (success)
Expand All @@ -491,8 +489,8 @@ mat GWRBasic::fitCoreOmp(const mat& x, const vec& y, const SpatialWeight& sw)
vec p = - si.t();
p(i) += 1.0;
qDiag_all.col(thread) += p % p;
mS.row(isStoreS() ? (i - workRange.first) : 0) = si;
mC.slice(isStoreC() ? (i - workRange.first) : 0) = ci;
mS.row(isStoreS() ? (i - mWorkRange.first) : 0) = si;
mC.slice(isStoreC() ? (i - mWorkRange.first) : 0) = ci;
}
catch (const exception& e)
{
Expand All @@ -518,9 +516,8 @@ arma::mat GWRBasic::fitCoreCVOmp(const arma::mat& x, const arma::vec& y, const S
uword nDp = mCoords.n_rows, nVar = x.n_cols;
mat betas(nVar, nDp, fill::zeros);
bool flag = true;
std::pair<uword, uword> workRange = mWorkRange.value_or(make_pair(0, nDp));
#pragma omp parallel for num_threads(mOmpThreadNum)
for (uword i = workRange.first; i < workRange.second; i++)
for (uword i = mWorkRange.first; i < mWorkRange.second; i++)
{
GWM_LOG_STOP_CONTINUE(mStatus);
if (flag)
Expand Down Expand Up @@ -551,9 +548,8 @@ arma::mat GWRBasic::fitCoreSHatOmp(const arma::mat& x, const arma::vec& y, const
mat betas(nVar, nDp, fill::zeros);
mat shat_all(2, mOmpThreadNum, fill::zeros);
int flag = true;
std::pair<uword, uword> workRange = mWorkRange.value_or(make_pair(0, nDp));
#pragma omp parallel for num_threads(mOmpThreadNum)
for (uword i = workRange.first; i < workRange.second; i++)
for (uword i = mWorkRange.first; i < mWorkRange.second; i++)
{
GWM_LOG_STOP_CONTINUE(mStatus);
if (flag)
Expand Down Expand Up @@ -597,8 +593,7 @@ mat GWRBasic::fitCoreCuda(const mat& x, const vec& y, const SpatialWeight& sw)
mQDiag = vec(nDp, arma::fill::zeros);
mS = mat(isStoreS() ? nDp : 1, nDp, arma::fill::zeros);
mC = cube(nVar, nDp, isStoreC() ? nDp : 1, arma::fill::zeros);
std::pair<uword, uword> workRange = mWorkRange.value_or(make_pair(0, nDp));
uword rangeSize = workRange.second - workRange.first;
uword rangeSize = mWorkRange.second - mWorkRange.first;
size_t groups = rangeSize / mGroupLength + (rangeSize % mGroupLength == 0 ? 0 : 1);
cumat u_xt(xt), u_y(y);
cumat u_betas(nVar, nDp), u_betasSE(nVar, nDp);
Expand All @@ -616,7 +611,7 @@ mat GWRBasic::fitCoreCuda(const mat& x, const vec& y, const SpatialWeight& sw)
for (size_t i = 0; i < groups; i++)
{
GWM_LOG_STOP_BREAK(mStatus);
size_t begin = workRange.first + i * mGroupLength, length = (begin + mGroupLength > workRange.second) ? (workRange.second - begin) : mGroupLength;
size_t begin = mWorkRange.first + i * mGroupLength, length = (begin + mGroupLength > mWorkRange.second) ? (mWorkRange.second - begin) : mGroupLength;
for (size_t j = 0, e = begin + j; j < length; j++, e++)
{
checkCudaErrors(sw.weightVector(e, u_dists.dmem(), u_weights.dmem()));
Expand Down Expand Up @@ -676,8 +671,7 @@ arma::mat GWRBasic::predictCuda(const mat& locations, const mat& x, const vec& y
{
uword nRp = locations.n_rows, nDp = mCoords.n_rows, nVar = x.n_cols;
mat betas(nVar, nRp, fill::zeros);
std::pair<uword, uword> workRange = mWorkRange.value_or(make_pair(0, nDp));
uword rangeSize = workRange.second - workRange.first;
uword rangeSize = mWorkRange.second - mWorkRange.first;
size_t groups = rangeSize / mGroupLength + (rangeSize % mGroupLength == 0 ? 0 : 1);
cumat u_xt(x.t()), u_y(y);
custride u_xtwx(nVar, nVar, mGroupLength), u_xtwy(nVar, 1, mGroupLength), u_xtwxI(nVar, nVar, mGroupLength);
Expand All @@ -690,7 +684,7 @@ arma::mat GWRBasic::predictCuda(const mat& locations, const mat& x, const vec& y
for (size_t i = 0; i < groups; i++)
{
GWM_LOG_STOP_BREAK(mStatus);
size_t begin = workRange.first + i * mGroupLength, length = (begin + mGroupLength > workRange.second) ? (workRange.second - begin) : mGroupLength;
size_t begin = mWorkRange.first + i * mGroupLength, length = (begin + mGroupLength > mWorkRange.second) ? (mWorkRange.second - begin) : mGroupLength;
for (size_t j = 0, e = begin + j; j < length; j++, e++)
{
checkCudaErrors(mSpatialWeight.weightVector(e, u_dists.dmem(), u_weights.dmem()));
Expand Down Expand Up @@ -731,13 +725,12 @@ arma::mat GWRBasic::fitCoreCVCuda(const arma::mat& x, const arma::vec& y, const
p_info = new int[mGroupLength];
checkCudaErrors(cudaMalloc(&d_info, sizeof(int) * mGroupLength));
bool success = true;
std::pair<uword, uword> workRange = mWorkRange.value_or(make_pair(0, nDp));
uword rangeSize = workRange.second - workRange.first;
uword rangeSize = mWorkRange.second - mWorkRange.first;
size_t groups = rangeSize / mGroupLength + (rangeSize % mGroupLength == 0 ? 0 : 1);
for (size_t i = 0; i < groups && success; i++)
{
GWM_LOG_STOP_BREAK(mStatus);
size_t begin = workRange.first + i * mGroupLength, length = (begin + mGroupLength > workRange.second) ? (workRange.second - begin) : mGroupLength;
size_t begin = mWorkRange.first + i * mGroupLength, length = (begin + mGroupLength > mWorkRange.second) ? (mWorkRange.second - begin) : mGroupLength;
for (size_t j = 0, e = begin + j; j < length; j++, e++)
{
checkCudaErrors(sw.weightVector(e, u_dists.dmem(), u_weights.dmem()));
Expand Down Expand Up @@ -782,13 +775,12 @@ arma::mat GWRBasic::fitCoreSHatCuda(const arma::mat& x, const arma::vec& y, cons
p_info = new int[mGroupLength];
checkCudaErrors(cudaMalloc(&d_info, sizeof(int) * mGroupLength));
bool success = true;
std::pair<uword, uword> workRange = mWorkRange.value_or(make_pair(0, nDp));
uword rangeSize = workRange.second - workRange.first;
uword rangeSize = mWorkRange.second - mWorkRange.first;
size_t groups = rangeSize / mGroupLength + (rangeSize % mGroupLength == 0 ? 0 : 1);
for (size_t i = 0; i < groups && success; i++)
{
GWM_LOG_STOP_BREAK(mStatus);
size_t begin = workRange.first + i * mGroupLength, length = (begin + mGroupLength > workRange.second) ? (workRange.second - begin) : mGroupLength;
size_t begin = mWorkRange.first + i * mGroupLength, length = (begin + mGroupLength > mWorkRange.second) ? (mWorkRange.second - begin) : mGroupLength;
for (size_t j = 0, e = begin + j; j < length; j++, e++)
{
checkCudaErrors(sw.weightVector(e, u_dists.dmem(), u_weights.dmem()));
Expand Down Expand Up @@ -991,11 +983,10 @@ arma::mat gwm::GWRBasic::fitMpi()
// }
GWM_MPI_MASTER_END
GWM_MPI_WORKER_BEGIN
std::pair<uword, uword> workRange = mWorkRange.value_or(make_pair(0, nDp));
// printf("%d process work range: [%lld, %lld]\n", mWorkerId, workRange.first, workRange.second);
// cout << mWorkerId << " process work range: [" << workRange.first << "," << workRange.second << "]\n";
mat betas = mBetas.cols(workRange.first, workRange.second - 1);
mat betasSE = mBetasSE.cols(workRange.first, workRange.second - 1);
// printf("%d process work range: [%lld, %lld]\n", mWorkerId, mWorkRange.first, mWorkRange.second);
// cout << mWorkerId << " process work range: [" << mWorkRange.first << "," << mWorkRange.second << "]\n";
mat betas = mBetas.cols(mWorkRange.first, mWorkRange.second - 1);
mat betasSE = mBetasSE.cols(mWorkRange.first, mWorkRange.second - 1);
MPI_Send(betas.memptr(), betas.n_elem, MPI_DOUBLE, 0, int(GWRBasicFitMpiTags::Betas), MPI_COMM_WORLD);
MPI_Send(betasSE.memptr(), betasSE.n_elem, MPI_DOUBLE, 0, int(GWRBasicFitMpiTags::BetasSE), MPI_COMM_WORLD);
MPI_Send(mSHat.memptr(), mSHat.n_elem, MPI_DOUBLE, 0, int(GWRBasicFitMpiTags::SHat), MPI_COMM_WORLD);
Expand Down
Loading