#include "pch.h" #define XOR_RAND(state, result_v...

创建于:2025年10月18日

使用 GPT-5 Thinking 回答 Chat01

问题

#include "pch.h"

#define XOR_RAND(state, result_var)
do {
int s = state;
s ^= s << 13;
s ^= s >> 17;
s ^= s << 5;
state = s;
result_var = state * 0x1.0p-32f;
} while(0)

#define XOR_RAND_GRSH(state, result_var)
do {
int s = state;
s ^= s << 13;
s ^= s >> 17;
s ^= s << 5;
state = s;
result_var = fmaf(state, 0x1.0p-31f, -1.0f);
} while(0)

#define FABE13_COS(x, result_var)
do {
float abs_val = fabsf(x);
float reduced = fmodf(abs_val, 6.28318530718f);
if(reduced > 3.14159265359f) {
reduced = 6.28318530718f - reduced;
}
if(reduced < 1.57079632679f) {
float val2 = reduced * reduced;
result_var = fmaf(val2, fmaf(val2, 0.0416666667f, -0.5f), 1.0f);
} else {
reduced = 3.14159265359f - reduced;
float val2 = reduced * reduced;
result_var = -fmaf(val2, fmaf(val2, 0.0416666667f, -0.5f), 1.0f);
}
} while(0)

#define FABE13_SINCOS(in, sin_out, cos_out, n)
do {
int i = 0;
int limit = n & ~7;
if(n >= 8) {
static __declspec(align(32)) const __m256 VEC_TWOPI = _mm256_set1_ps(6.28318530718f);
static __declspec(align(32)) const __m256 VEC_PI = _mm256_set1_ps(3.14159265359f);
static __declspec(align(32)) const __m256 VEC_PI_2 = _mm256_set1_ps(1.57079632679f);
static __declspec(align(32)) const __m256 INV_TWOPI = _mm256_set1_ps(0.15915494309189535f);
static __declspec(align(32)) const __m256 BIAS = _mm256_set1_ps(12582912.0f);
static __declspec(align(32)) const __m256 VEC_COS_P3 = _mm256_set1_ps(0.0416666667f);
static __declspec(align(32)) const __m256 VEC_COS_P1 = _mm256_set1_ps(-0.5f);
static __declspec(align(32)) const __m256 VEC_COS_P0 = _mm256_set1_ps(1.0f);
static __declspec(align(32)) const __m256 VEC_SIN_P1 = _mm256_set1_ps(-0.16666666f);
static __declspec(align(32)) const __m256 VEC_SIN_P0 = _mm256_set1_ps(1.0f);
static __declspec(align(32)) const __m256 VEC_ZERO = _mm256_setzero_ps();
while(i < limit) {
__m256 vx = _mm256_load_ps(&in[i]);
__m256 vax = _mm256_andnot_ps(_mm256_set1_ps(-0.0f), vx);
__m256 q = _mm256_fmadd_ps(vax, INV_TWOPI, BIAS);
q = _mm256_sub_ps(q, BIAS);
__m256 r = _mm256_fnmadd_ps(VEC_TWOPI, q, vax);
__m256 r1 = _mm256_min_ps(r, _mm256_sub_ps(VEC_TWOPI, r));
__m256 r2 = _mm256_min_ps(r1, _mm256_sub_ps(VEC_PI, r1));
__m256 t2 = _mm256_mul_ps(r2, r2);
__m256 cosv = _mm256_fmadd_ps(t2, _mm256_fmadd_ps(t2, VEC_COS_P3, VEC_COS_P1), VEC_COS_P0);
__m256 sinv = _mm256_mul_ps(_mm256_fmadd_ps(t2, VEC_SIN_P1, VEC_SIN_P0), r2);
__m256 cflip = _mm256_cmp_ps(r1, VEC_PI_2, _CMP_GT_OQ);
__m256 sflip = _mm256_xor_ps(_mm256_cmp_ps(vx, VEC_ZERO, _CMP_LT_OQ), _mm256_cmp_ps(r, VEC_PI, _CMP_GT_OQ));
cosv = _mm256_blendv_ps(cosv, _mm256_sub_ps(VEC_ZERO, cosv), cflip);
sinv = _mm256_blendv_ps(sinv, _mm256_sub_ps(VEC_ZERO, sinv), sflip);
_mm256_store_ps(&cos_out[i], cosv);
_mm256_store_ps(&sin_out[i], sinv);
i += 8;
}
}
while(i < n) {
float x = in[i];
float ax = fabsf(x);
float q = fmaf(ax, 0.15915494309189535f, 12582912.0f);
q -= 12582912.0f;
float r = fmaf(-6.28318530718f, q, ax);
bool sflip = r > 3.14159265359f;
if(sflip) r = 6.28318530718f - r;
bool cflip = r > 1.57079632679f;
if(cflip) r = 3.14159265359f - r;
float t2 = r * r;
float c = fmaf(t2, fmaf(t2, 0.0416666667f, -0.5f), 1.0f);
float s = fmaf(t2, -0.16666666f, 1.0f) * r;
cos_out[i] = cflip ? -c : c;
sin_out[i] = ((x < 0.0f) ^ sflip) ? -s : s;
++i;
}
} while(0)

static __declspec(noalias) __forceinline float ShekelFunc(const float x, const float seed) noexcept
{
int i = 0;
float current_state = seed, current_res, current_res2, res = 0.0f;
while (i < 10) {
XOR_RAND(current_state, current_res);
const float x_part = fmaf(-current_res, 10.0f, x);
XOR_RAND(current_state, current_res);
XOR_RAND(current_state, current_res2);
float delimiter = fmaf(fmaf(current_res, 20.0f, 5.0f), x_part * x_part, fmaf(current_res2, 0.2f, 1.0f));
delimiter = copysignf(fmaxf(fabsf(delimiter), FLT_MIN), delimiter);
res -= 1.0f / delimiter;
++i;
}
return res;
}

static __declspec(noalias) __forceinline float RastriginFunc(const float x1, const float x2) noexcept
{
const float term1 = fmaf(x1, x1, x2 * x2);
float cos1, cos2;
FABE13_COS(6.28318530717958647692f * x1, cos1);
FABE13_COS(6.28318530717958647692f * x2, cos2);
return (term1 - fmaf(cos1 + cos2, 10.0f, -14.6f)) * fmaf(-term1, 0.25f, 18.42f);
}

static __declspec(noalias) __forceinline float HillFunc(const float x, const float seed) noexcept
{
int j = 0;
__declspec(align(32)) float angles[14u];
const float start_angle = 6.28318530717958647692f * x;
#pragma loop(ivdep)
while (j < 14) {
angles[j] = start_angle * static_cast<float>(j + 1);
++j;
}
__declspec(align(32)) float sin_vals[14u];
__declspec(align(32)) float cos_vals[14u];
FABE13_SINCOS(angles, sin_vals, cos_vals, 14u);
float current_state = seed, current_res, current_res2;
XOR_RAND(current_state, current_res);
float res = fmaf(current_res, 2.0f, -1.1f);
--j;
while (j >= 0) {
XOR_RAND(current_state, current_res);
XOR_RAND(current_state, current_res2);
res += fmaf(fmaf(current_res, 2.0f, -1.1f), sin_vals[j], fmaf(current_res2, 2.0f, -1.1f) * cos_vals[j]);
--j;
}
return res;
}

static __declspec(noalias) __forceinline float GrishaginFunc(const float x1, const float x2, const float seed) noexcept
{
int j = 0;
__declspec(align(32)) float angles_j[8u];
__declspec(align(32)) float angles_k[8u];
#pragma loop(ivdep)
while (j < 8) {
const float pj_mult = 3.14159265358979323846f * static_cast<float>(j + 1);
angles_j[j] = pj_mult * x1;
angles_k[j] = pj_mult * x2;
++j;
}
__declspec(align(32)) float sin_j[8u], cos_j[8u];
__declspec(align(32)) float sin_k[8u], cos_k[8u];
FABE13_SINCOS(angles_j, sin_j, cos_j, 8u);
FABE13_SINCOS(angles_k, sin_k, cos_k, 8u);
--j;
float part1 = 0.0f;
float part2 = 0.0f;
float current_state = seed, current_res, current_res2;
while (j >= 0) {
size_t k = 0u;
while (k < 8u) {
const float sin_term = sin_j[j] * sin_j[j];
const float cos_term = cos_k[k] * cos_k[k];
XOR_RAND_GRSH(current_state, current_res);
XOR_RAND_GRSH(current_state, current_res2);
part1 = fmaf(current_res, sin_term, fmaf(current_res2, cos_term, part1));
XOR_RAND_GRSH(current_state, current_res);
XOR_RAND_GRSH(current_state, current_res2);
part2 = fmaf(-current_res, cos_term, fmaf(current_res2, sin_term, part2));
++k;
}
--j;
}
return -sqrtf(fmaf(part1, part1, part2 * part2));
}

static __declspec(noalias) __forceinline float Shag(const float _m, const float x1, const float x2, const float y1,
const float y2, const float _N, const float _r) noexcept
{
const float diff = y2 - y1;
const float sign_mult = _mm_cvtss_f32(_mm_castsi128_ps(_mm_set1_epi32(
0x3F800000u | ((reinterpret_cast<const uint32_t>(&diff) & 0x80000000u) ^ 0x80000000u))));
return _N == 1.0f
? fmaf(-(1.0f / _m), diff, x1 + x2) * 0.5f
: _N == 2.0f
? fmaf(sign_mult / (_m * _m), diff * diff * _r, x1 + x2) * 0.5f
: fmaf(sign_mult / powf(_m, _N), powf(diff, _N) * _r, x1 + x2) * 0.5f;
}

__declspec(align(16)) struct Slab final {
char* __restrict base;
char* __restrict current;
char* __restrict end;

text
__declspec(noalias) __forceinline Slab(void* __restrict const memory, const size_t usable_size) noexcept : base(static_cast<char* __restrict>(memory)) , current(base) , end(base + (usable_size & ~static_cast<size_t>(63u))) { }

};

static tbb::enumerable_thread_specific<Slab*> tls( noexcept {
void* __restrict const memory = _aligned_malloc(16777216u, 16u);
Slab* __restrict const slab = static_cast<Slab*>(_aligned_malloc(32u, 16u));
*slab = Slab(memory, 16777216u);

char* __restrict p = slab->base;

#pragma loop(ivdep)
while (p < slab->end) {
*p = 0u;
p += 4096u;
}
return slab;
}());

__declspec(align(64)) struct Interval final {
const float x1;
const float x2;
const float y1;
const float y2;
const float delta_y;
const float ordinate_factor;
const float N_factor;
const float quadratic_term;
const float M;
float R;

text
__declspec(noalias) __forceinline void* operator new(const size_t) noexcept { Slab* __restrict s = tls.local(); char* __restrict result = s->current; s->current += 64u; return result; } __declspec(noalias) __forceinline Interval(const float _x1, const float _x2, const float _y1, const float _y2, const float _N) noexcept : x1(_x1), x2(_x2), y1(_y1), y2(_y2) , delta_y(_y2 - _y1) , ordinate_factor(-(y1 + y2) * 2.0f) , N_factor(_N == 1.0f ? _x2 - _x1 : _N == 2.0f ? sqrtf(_x2 - _x1) : powf(_x2 - _x1, 1.0f / _N)) , quadratic_term((1.0f / N_factor)* delta_y* delta_y) , M((1.0f / N_factor)* fabsf(delta_y)) { } __declspec(noalias) __forceinline void ChangeCharacteristic(const float _m) noexcept { R = fmaf(1.0f / _m, quadratic_term, fmaf(_m, N_factor, ordinate_factor)); }

};

static __declspec(noalias) __forceinline bool ComparePtr(const Interval* __restrict a, const Interval* __restrict b) noexcept
{
return a->R < b->R;
}

const enum List : uint8_t { Top = 0b00u, Down = 0b01u, Left = 0b10u, Right = 0b11u };

__declspec(align(16)) struct Peano2DMap final {
const int levels;
const float a, b, c, d;
const float lenx, leny;
const float inv_lenx;
const uint32_t scale;
const uint8_t start;

text
__declspec(noalias) __forceinline Peano2DMap( const int L, const float _a, const float _b, const float _c, const float _d, const uint8_t startType ) noexcept : levels(L) , a(_a), b(_b), c(_c), d(_d) , lenx(_b - _a) , leny(_d - _c) , inv_lenx(1.0f / (_b - _a)) , scale(static_cast<uint32_t>(1u) << (L << 1)) , start(startType) { }

};

__declspec(align(4)) struct Step final {
const uint8_t next;
const uint8_t dx;
const uint8_t dy;
};

__declspec(align(4)) struct InvStep final {
const uint8_t q;
const uint8_t next;
};

__declspec(align(64)) static const Step g_step_tbl[4][4] = {
{ {Right,0u,0u}, {Top,0u,1u}, {Top,1u,1u}, {Left,1u,0u} },
{ {Left,1u,1u}, {Down,1u,0u}, {Down,0u,0u}, {Right,0u,1u} },
{ {Down,1u,1u}, {Left,0u,1u}, {Left,0u,0u}, {Top,1u,0u} },
{ {Top,0u,0u}, {Right,1u,0u}, {Right,1u,1u}, {Down,0u,1u} }
};

__declspec(align(64)) static const InvStep g_inv_tbl[4][4] = {
{ {0u,Right}, {1u,Top}, {3u,Left}, {2u,Top} },
{ {2u,Down}, {3u,Right}, {1u,Down}, {0u,Left} },
{ {2u,Left}, {1u,Left}, {3u,Top}, {0u,Down} },
{ {0u,Top}, {3u,Down}, {1u,Right}, {2u,Right} }
};

static Peano2DMap gActiveMap(0, 0.0f, 0.0f, 0.0f, 0.0f, 0b00u);

__declspec(noalias) __forceinline void HitTest2D_analytic(const float x_param, float& out_x1, float& out_x2) noexcept
{
const float a = gActiveMap.a;
const float inv_lenx = gActiveMap.inv_lenx;
const uint32_t scale = gActiveMap.scale;
const uint32_t scale_minus_1 = scale - 1u;
const float lenx = gActiveMap.lenx;
const float leny = gActiveMap.leny;
const float c = gActiveMap.c;
const uint8_t start = gActiveMap.start;
const int levels = gActiveMap.levels;

text
float norm = (x_param - a) * inv_lenx; norm = fminf(fmaxf(norm, 0.0f), 0x1.fffffep-1f); uint32_t idx = static_cast<uint32_t>(norm * static_cast<float>(scale)); idx = idx > scale_minus_1 ? scale_minus_1 : idx; float sx = lenx, sy = leny; float x1 = a, x2 = c; uint8_t type = start; int l = levels - 1;

#pragma loop(ivdep)
while (l >= 0) {
const uint32_t q = (idx >> (l * 2)) & 3u;
const Step s = g_step_tbl[type][q];
type = s.next;
sx *= 0.5f; sy *= 0.5f;
x1 += s.dx ? sx : 0.0f;
x2 += s.dy ? sy : 0.0f;
--l;
}
out_x1 = x1 + sx * 0.5f;
out_x2 = x2 + sy * 0.5f;
}

__declspec(noalias) __forceinline float FindX2D_analytic(const float px, const float py) noexcept
{
const float a = gActiveMap.a;
const float b = gActiveMap.b;
const float c = gActiveMap.c;
const float d = gActiveMap.d;
const float lenx = gActiveMap.lenx;
const float leny = gActiveMap.leny;
const uint32_t scale = gActiveMap.scale;
const uint8_t start = gActiveMap.start;
const int levels = gActiveMap.levels;

text
const float clamped_px = fminf(fmaxf(px, a), b); const float clamped_py = fminf(fmaxf(py, c), d); float sx = lenx, sy = leny; float x0 = a, y0 = c; uint32_t idx = 0u; uint8_t type = start; int l = 0;

#pragma loop(ivdep)
while (l < levels) {
sx *= 0.5f; sy *= 0.5f;
const float mx = x0 + sx;
const float my = y0 + sy;

text
const uint32_t tr = static_cast<uint32_t>((clamped_px > mx) & (clamped_py > my)); const uint32_t tl = static_cast<uint32_t>((clamped_px < mx) & (clamped_py > my)); const uint32_t dl = static_cast<uint32_t>((clamped_px < mx) & (clamped_py < my)); const uint32_t none = static_cast<uint32_t>(1u ^ (tr | tl | dl)); const uint32_t dd = (tr << 1) | tr | tl | (none << 1); const InvStep inv = g_inv_tbl[type][dd]; type = inv.next; idx = (idx << 2) | inv.q; const uint32_t dx = dd >> 1; const uint32_t dy = dd & 1u; x0 += dx ? sx : 0.0f; y0 += dy ? sy : 0.0f; ++l; } const float scale_reciprocal = 1.0f / static_cast<float>(scale); return fmaf(static_cast<float>(idx) * scale_reciprocal, lenx, a);

}

static boost::mpi::environment* __restrict g_env;
static boost::mpi::communicator* __restrict g_world;
static boost::mpi::request recv_req;

extern "C" __declspec(dllexport) __declspec(noalias) __forceinline int AgpInit(const int peanoLevel, const float a, const float b, const float c, const float d) noexcept
{
g_env = new boost::mpi::environment();
g_world = new boost::mpi::communicator();
const int rank = g_world->rank();

text
new(&gActiveMap) Peano2DMap(peanoLevel, a, b, c, d, rank ? static_cast<uint8_t>(Down) : static_cast<uint8_t>(Top)); return rank;

}

__declspec(align(16)) struct CrossMsg final {
float s_x1, s_x2;
float e_x1, e_x2;
float Rtop;
template<typename Archive>
__declspec(noalias) __forceinline void serialize(Archive& __restrict const ar, const unsigned int) noexcept { ar& s_x1& s_x2& e_x1& e_x2& Rtop; }
};

static __declspec(noalias) __forceinline void RecomputeR_ConstM_AVX2(Interval* const* __restrict arr, const size_t n, const float m) noexcept {
const __m256 vm = _mm256_set1_ps(m);

text
__m256 vinvm = _mm256_rcp_ps(vm); vinvm = _mm256_mul_ps(vinvm, _mm256_fnmadd_ps(vm, vinvm, _mm256_set1_ps(2.0f))); size_t i = 0; const int limit = static_cast<int>(n & ~7u); if (n >= 8u) { while (i < static_cast<size_t>(limit)) { __declspec(align(32)) float q[8], nf[8], ord[8]; int k = 0;

#pragma loop(ivdep)
while (k < 8) {
const Interval* const __restrict p = arr[i + k];
q[k] = p->quadratic_term;
nf[k] = p->N_factor;
ord[k] = p->ordinate_factor;
++k;
}
__m256 vq = _mm256_load_ps(q);
__m256 vnf = _mm256_load_ps(nf);
__m256 vod = _mm256_load_ps(ord);

text
__m256 t = _mm256_fmadd_ps(vm, vnf, vod); __m256 res = _mm256_fmadd_ps(vq, vinvm, t); __declspec(align(32)) float out[8]; _mm256_store_ps(out, res); k = 0;

#pragma loop(ivdep)
while (k < 8) {
arr[i + k]->R = out[k];
++k;
}
i += 8u;
}
}
while (i < n) {
arr[i]->ChangeCharacteristic(m);
++i;
}
}

static __declspec(noalias) __forceinline void RecomputeR_AffineM_AVX2(Interval* const* __restrict arr, const size_t n, const float GLOBAL_FACTOR, const float alpha) noexcept {
const __m256 vGF = _mm256_set1_ps(GLOBAL_FACTOR);
const __m256 va = _mm256_set1_ps(alpha);

text
size_t i = 0; int limit = static_cast<int>(n & ~7u); if (n >= 8u) { while (i < static_cast<size_t>(limit)) { __declspec(align(32)) float len[8], Mv[8], q[8], nf[8], ord[8]; int k = 0;

#pragma loop(ivdep)
while (k < 8) {
const Interval* p = arr[i + k];
len[k] = p->x2 - p->x1;
Mv[k] = p->M;
q[k] = p->quadratic_term;
nf[k] = p->N_factor;
ord[k] = p->ordinate_factor;
++k;
}
__m256 vlen = _mm256_load_ps(len);
__m256 vM = _mm256_load_ps(Mv);
__m256 vq = _mm256_load_ps(q);
__m256 vnf = _mm256_load_ps(nf);
__m256 vod = _mm256_load_ps(ord);

text
__m256 vm = _mm256_fmadd_ps(vGF, vlen, _mm256_mul_ps(va, vM)); __m256 vinvm = _mm256_rcp_ps(vm); vinvm = _mm256_mul_ps(vinvm, _mm256_fnmadd_ps(vm, vinvm, _mm256_set1_ps(2.0f))); __m256 t = _mm256_fmadd_ps(vm, vnf, vod); __m256 res = _mm256_fmadd_ps(vq, vinvm, t); __declspec(align(32)) float out[8]; _mm256_store_ps(out, res); k = 0;

#pragma loop(ivdep)
while (k < 8) {
arr[i + k]->R = out[k];
++k;
}
i += 8u;
}
}
while (i < n) {
const Interval* const __restrict p = arr[i];
const float mi = fmaf(GLOBAL_FACTOR, p->x2 - p->x1, p->M * alpha);
arr[i]->R = fmaf(1.0f / mi, p->quadratic_term, fmaf(mi, p->N_factor, p->ordinate_factor));
++i;
}
}

extern "C" __declspec(dllexport) __declspec(noalias)
void AGP_1D(const float global_iterations,
const float a, const float b, const float r,
const bool mode, const float epsilon, const float seed,
float** __restrict out_data, size_t* __restrict out_len) noexcept
{
int schetchick = 0;
const float initial_length = b - a;
float dmax = initial_length;
const float threshold_03 = 0.3f * initial_length;
const float inv_threshold_03 = 1.0f / threshold_03;
const float start_val = ShekelFunc(a, seed);
float best_f = ShekelFunc(b, seed);
float x_Rmax_1 = a;
float x_Rmax_2 = b;
float y_Rmax_1 = start_val;
float y_Rmax_2 = best_f;

text
std::vector<float, boost::alignment::aligned_allocator<float, 16u>> Extr; std::vector<Interval* __restrict, boost::alignment::aligned_allocator<Interval* __restrict, 64u>> R; Extr.clear(); Extr.reserve(static_cast<size_t>(global_iterations) << 2u); R.clear(); R.reserve(static_cast<size_t>(global_iterations) << 1u); R.emplace_back(new Interval(a, b, start_val, best_f, 1.0f)); float Mmax = R.front()->M; float m = r * Mmax; while (true) { float new_point = Shag(m, x_Rmax_1, x_Rmax_2, y_Rmax_1, y_Rmax_2, 1.0f, r); float new_value = ShekelFunc(new_point, seed); if (new_value < best_f) { best_f = new_value; Extr.emplace_back(best_f); Extr.emplace_back(new_point); } std::pop_heap(R.begin(), R.end(), ComparePtr); Interval* __restrict promejutochny_otrezok = R.back(); float new_x1 = promejutochny_otrezok->x1; float new_x2 = promejutochny_otrezok->x2; float len2 = new_x2 - new_point; float len1 = new_point - new_x1; float interval_len = len1 < len2 ? len1 : len2; if (interval_len < epsilon || ++schetchick == static_cast<int>(global_iterations)) { Extr.emplace_back(static_cast<float>(schetchick)); Extr.emplace_back(interval_len); *out_len = Extr.size(); *out_data = reinterpret_cast<float* __restrict>(CoTaskMemAlloc(sizeof(float) * (*out_len))); memcpy(*out_data, Extr.data(), sizeof(float) * (*out_len)); return; } Interval* __restrict curr = new Interval(new_x1, new_point, promejutochny_otrezok->y1, new_value, 1.0f); Interval* __restrict curr1 = new Interval(new_point, new_x2, new_value, promejutochny_otrezok->y2, 1.0f); float currM = curr->M > curr1->M ? curr->M : curr1->M; size_t r_size = R.size(); if (mode) { if (len2 + len1 == dmax) { dmax = len2 > len1 ? len2 : len1; size_t i = 0u;

#pragma loop(ivdep)
while (i < r_size) {
float len_item = R[i]->x2 - R[i]->x1;
if (len_item > dmax) dmax = len_item;
++i;
}
}

text
if (threshold_03 > dmax && schetchick % 3 == 0 || 10.0f * dmax < initial_length) { if (currM > Mmax) { Mmax = currM; m = r * Mmax; } float progress = fmaf(-inv_threshold_03, dmax, 1.0f); float alpha = fmaf(progress, progress, 1.0f); float betta = 2.0f - alpha; float MULTIPLIER = (1.0f / dmax) * Mmax; float global_coeff = fmaf(MULTIPLIER, r, -MULTIPLIER); float GLOBAL_FACTOR = betta * global_coeff; curr->ChangeCharacteristic(fmaf(GLOBAL_FACTOR, len1, curr->M * alpha)); curr1->ChangeCharacteristic(fmaf(GLOBAL_FACTOR, len2, curr1->M * alpha)); if (r_size < 64u) { size_t i = 0u;

#pragma loop(ivdep)
while (i < r_size) {
R[i]->ChangeCharacteristic(fmaf(GLOBAL_FACTOR, R[i]->x2 - R[i]->x1, R[i]->M * alpha));
++i;
}
}
else {
RecomputeR_AffineM_AVX2(R.data(), r_size, GLOBAL_FACTOR, alpha);
}
std::make_heap(R.begin(), R.end(), ComparePtr);
}
else {
if (currM > Mmax) {
if (currM - Mmax < Mmax * 0.15f) {
Mmax = currM;
m = r * Mmax;
curr->ChangeCharacteristic(m);
curr1->ChangeCharacteristic(m);
}
else {
Mmax = currM;
m = r * Mmax;
curr->ChangeCharacteristic(m);
curr1->ChangeCharacteristic(m);
if (r_size < 64u) {
size_t i = 0u;
#pragma loop(ivdep)
while (i < r_size) {
R[i]->ChangeCharacteristic(m);
++i;
}
}
else {
RecomputeR_ConstM_AVX2(R.data(), r_size, m);
}
std::make_heap(R.begin(), R.end(), ComparePtr);
}
}
else {
curr->ChangeCharacteristic(m);
curr1->ChangeCharacteristic(m);
}
}
}
else {
if (currM > Mmax) {
if (currM - Mmax < Mmax * 0.15f) {
Mmax = currM;
m = r * Mmax;
curr->ChangeCharacteristic(m);
curr1->ChangeCharacteristic(m);
}
else {
Mmax = currM;
m = r * Mmax;
curr->ChangeCharacteristic(m);
curr1->ChangeCharacteristic(m);
if (r_size < 64u) {
size_t i = 0u;
#pragma loop(ivdep)
while (i < r_size) {
R[i]->ChangeCharacteristic(m);
++i;
}
}
else {
RecomputeR_ConstM_AVX2(R.data(), r_size, m);
}
std::make_heap(R.begin(), R.end(), ComparePtr);
}
}
else {
curr->ChangeCharacteristic(m);
curr1->ChangeCharacteristic(m);
}
}

text
R.back() = curr; std::push_heap(R.begin(), R.end(), ComparePtr); R.emplace_back(curr1); std::push_heap(R.begin(), R.end(), ComparePtr); Interval* __restrict top_ptr = R.front(); x_Rmax_1 = top_ptr->x1; x_Rmax_2 = top_ptr->x2; y_Rmax_1 = top_ptr->y1; y_Rmax_2 = top_ptr->y2; }

}

extern "C" __declspec(dllexport) __declspec(noalias)
void AGP_2D(
const float N, const float global_iterations,
const float a, const float b, const float c, const float d, const float r,
const bool mode, const float epsilon, const float seed,
float** __restrict out_data, size_t* __restrict out_len) noexcept
{
const int rank = g_world->rank();
const int partner = rank ^ 1;
int dummy;
if (rank) recv_req = g_world->irecv(1, 0, dummy);

text
const float inv_divider = ldexpf(1.0f, -((gActiveMap.levels << 1) + 1)); const float x_addition = (b - a) * inv_divider; const float y_addition = (d - c) * inv_divider; const float true_start = a + x_addition; const float true_end = b - x_addition; float x_Rmax_1 = true_start; float x_Rmax_2 = true_end; const float initial_length = true_end - true_start; float dmax = initial_length; const float threshold_03 = 0.3f * initial_length; const float inv_threshold_03 = 1.0f / threshold_03; const float start_val = rank ? RastriginFunc(true_end, d - y_addition) : RastriginFunc(true_start, c + y_addition); float y_Rmax_1 = start_val; float best_f = rank ? RastriginFunc(true_start, d - y_addition) : RastriginFunc(true_end, c + y_addition); float y_Rmax_2 = best_f; std::vector<float, boost::alignment::aligned_allocator<float, 16u>> Extr; std::vector<Interval* __restrict, boost::alignment::aligned_allocator<Interval* __restrict, 64u>> R; Extr.clear(); Extr.reserve(static_cast<size_t>(global_iterations) << 2u); R.clear(); R.reserve(static_cast<size_t>(global_iterations) << 1u); R.emplace_back(new Interval(true_start, true_end, start_val, best_f, 2.0f)); float Mmax = R.front()->M; float m = r * Mmax; int schetchick = 0; float interval_len = 0.0f; while (true) { if (g_world->iprobe(1, 0)) { Extr.emplace_back(static_cast<float>(schetchick)); Extr.emplace_back(interval_len); *out_len = Extr.size(); *out_data = reinterpret_cast<float* __restrict>(CoTaskMemAlloc(sizeof(float) * (*out_len))); memcpy(*out_data, Extr.data(), sizeof(float) * (*out_len)); recv_req.wait(); recv_req = g_world->irecv(1, 0, dummy); return; } float cooling = ldexpf(1.0f, -(1.0f / 138.63f) * ++schetchick); int T = static_cast<int>(fmaf(20.0f, cooling, 10.0f)); float k = fmaf(0.2f, cooling, 0.7f); Interval* __restrict top_ptr = R.front(); if (schetchick % T == 0) { float s_x1, s_x2, e_x1, e_x2; HitTest2D_analytic(top_ptr->x1, s_x1, s_x2); HitTest2D_analytic(top_ptr->x2, e_x1, e_x2); CrossMsg outbound = CrossMsg{ s_x1, s_x2, e_x1, e_x2, top_ptr->R }; CrossMsg inbound; g_world->sendrecv(partner, 0, outbound, partner, 0, inbound); const float sx = FindX2D_analytic(inbound.s_x1, inbound.s_x2); const float ex = FindX2D_analytic(inbound.e_x1, inbound.e_x2); Interval* __restrict injected = new Interval(sx, ex, RastriginFunc(inbound.s_x1, inbound.s_x2), RastriginFunc(inbound.e_x1, inbound.e_x2), 2.0f); injected->R = inbound.Rtop * k; R.emplace_back(injected); std::push_heap(R.begin(), R.end(), ComparePtr); } float new_point = Shag(m, x_Rmax_1, x_Rmax_2, y_Rmax_1, y_Rmax_2, 2.0f, r); float new_x1_val, new_x2_val; HitTest2D_analytic(new_point, new_x1_val, new_x2_val); float new_value = RastriginFunc(new_x1_val, new_x2_val); if (new_value < best_f) { best_f = new_value; Extr.emplace_back(best_f); Extr.emplace_back(new_x1_val); Extr.emplace_back(new_x2_val); } std::pop_heap(R.begin(), R.end(), ComparePtr); Interval* __restrict promejutochny_otrezok = R.back(); float segment_x1 = promejutochny_otrezok->x1; float segment_x2 = promejutochny_otrezok->x2; float len2 = segment_x2 - new_point; float len1 = new_point - segment_x1; interval_len = len1 < len2 ? len1 : len2; if (interval_len < epsilon || schetchick == static_cast<int>(global_iterations)) { g_world->send(0, 0, 0); return; } Interval* __restrict curr = new Interval(segment_x1, new_point, promejutochny_otrezok->y1, new_value, 2.0f); Interval* __restrict curr1 = new Interval(new_point, segment_x2, new_value, promejutochny_otrezok->y2, 2.0f); float currM = curr->M > curr1->M ? curr->M : curr1->M; size_t r_size = R.size(); if (mode) { if (len2 + len1 == dmax) { dmax = len2 > len1 ? len2 : len1; size_t i = 0u;

#pragma loop(ivdep)
while (i < r_size) {
float len_item = R[i]->x2 - R[i]->x1;
if (len_item > dmax) dmax = len_item;
++i;
}
}

text
if (threshold_03 > dmax && schetchick % 3 == 0 || 10.0f * dmax < initial_length) { if (currM > Mmax) { Mmax = currM; m = r * Mmax; } float progress = fmaf(-inv_threshold_03, dmax, 1.0f); float alpha = fmaf(progress, progress, 1.0f); float betta = 2.0f - alpha; float MULTIPLIER = (1.0f / dmax) * Mmax; float global_coeff = fmaf(MULTIPLIER, r, -MULTIPLIER); float GLOBAL_FACTOR = betta * global_coeff; curr->ChangeCharacteristic(fmaf(GLOBAL_FACTOR, len1, curr->M * alpha)); curr1->ChangeCharacteristic(fmaf(GLOBAL_FACTOR, len2, curr1->M * alpha)); if (r_size < 64u) { size_t i = 0u;

#pragma loop(ivdep)
while (i < r_size) {
R[i]->ChangeCharacteristic(fmaf(GLOBAL_FACTOR, R[i]->x2 - R[i]->x1, R[i]->M * alpha));
++i;
}
}
else {
RecomputeR_AffineM_AVX2(R.data(), r_size, GLOBAL_FACTOR, alpha);
}
std::make_heap(R.begin(), R.end(), ComparePtr);
}
else {
if (currM > Mmax) {
if (currM - Mmax < Mmax * 0.15f) {
Mmax = currM;
m = r * Mmax;
curr->ChangeCharacteristic(m);
curr1->ChangeCharacteristic(m);
}
else {
Mmax = currM;
m = r * Mmax;
curr->ChangeCharacteristic(m);
curr1->ChangeCharacteristic(m);
if (r_size < 64u) {
size_t i = 0u;
#pragma loop(ivdep)
while (i < r_size) {
R[i]->ChangeCharacteristic(m);
++i;
}
}
else {
RecomputeR_ConstM_AVX2(R.data(), r_size, m);
}
std::make_heap(R.begin(), R.end(), ComparePtr);
}
}
else {
curr->ChangeCharacteristic(m);
curr1->ChangeCharacteristic(m);
}
}
}
else {
if (currM > Mmax) {
if (currM - Mmax < Mmax * 0.15f) {
Mmax = currM;
m = r * Mmax;
curr->ChangeCharacteristic(m);
curr1->ChangeCharacteristic(m);
}
else {
Mmax = currM;
m = r * Mmax;
curr->ChangeCharacteristic(m);
curr1->ChangeCharacteristic(m);
if (r_size < 64u) {
size_t i = 0u;
#pragma loop(ivdep)
while (i < r_size) {
R[i]->ChangeCharacteristic(m);
++i;
}
}
else {
RecomputeR_ConstM_AVX2(R.data(), r_size, m);
}
std::make_heap(R.begin(), R.end(), ComparePtr);
}
}
else {
curr->ChangeCharacteristic(m);
curr1->ChangeCharacteristic(m);
}
}

text
R.back() = curr; std::push_heap(R.begin(), R.end(), ComparePtr); R.emplace_back(curr1); std::push_heap(R.begin(), R.end(), ComparePtr); top_ptr = R.front(); x_Rmax_1 = top_ptr->x1; x_Rmax_2 = top_ptr->x2; y_Rmax_1 = top_ptr->y1; y_Rmax_2 = top_ptr->y2; }

}

extern "C" __declspec(dllexport) __declspec(noalias) __forceinline void AgpWaitStartAndRun() noexcept
{
int dummy;
recv_req = g_world->irecv(0, 0, dummy);
float* __restrict buf;
size_t len;
while (true) {
if (recv_req.test()) {
AGP_2D(2.0f, 10000.0f, -2.2f, 1.8f, -2.2f, 1.8f, 2.5f, false, 0.0001f, GetTickCount(), &buf, &len);
recv_req = g_world->irecv(0, 0, dummy);
}
}
}

extern "C" __declspec(dllexport) __declspec(noalias) __forceinline void AgpStartWorkers() noexcept
{
g_world->isend(1, 0, 0);
}

extern "C" __declspec(dllexport) __declspec(noalias) __forceinline void AGP_Free(float* __restrict p) noexcept
{
CoTaskMemFree(p);
} мой одномерный/двумерный АГП, проблема двумерного АГП: алгоритм завершается сразу же после первого же обмена между процессами, после примерно 30 итераций, теоретически могло быть что сразу же после первого же обмена были найдены очень перспективные области, но когда менял параметры обмена и установил чтобы обмен был через каждые 10 итераций то тогда алгоритм завешался примерно на 10-11 итерациях - чего быть не должно, если поставить обмен через каждые 50 итераций - то алгоритм зависнет и вызовет ошибку: C:\Users\maxim>mpiexec -n 2 D:\TEST_AGP\x64\Release\TEST_AGP.exe

job aborted:
[ranks] message

[0] terminated

[1] process exited without calling finalize

---- error analysis -----

[1] on DESKTOP-C1VT6JF
D:\TEST_AGP\x64\Release\TEST_AGP.exe ended prematurely and may have crashed. exit code 0xc0000409

---- error analysis ----- - это dll вызываемая в управляемом winforms коде: MyForm.cpp #include "MyForm.h"
#include <float.h>

using namespace System;
using namespace System::Windows::Forms;

typedef int(__cdecl* PInit)(int, float, float, float, float);
typedef void(__cdecl* PStartWorkers)();

[STAThread]
int main() {
HMODULE h = LoadLibraryW(L"TEST_FUNC.dll");
auto AgpInit = (PInit)GetProcAddress(h, "AgpInit");
auto AgpWaitStartAndRun = (PStartWorkers)GetProcAddress(h, "AgpWaitStartAndRun");

text
const int rank = AgpInit(12, -2.2f, 1.8f, -2.2f, 1.8f); if (rank == 0) { Application::EnableVisualStyles(); Application::SetCompatibleTextRenderingDefault(false); Application::Run(gcnew TESTAGP::MyForm(h)); } else { AgpWaitStartAndRun(); } return 0;

} MyForm.h #pragma once

#define WIN32_LEAN_AND_MEAN
#include <Windows.h>

#include <algorithm>

#include <vector>

#include <utility>

namespace TESTAGP {
using namespace System;
using namespace System::IO;

text
ref class MyForm : public System::Windows::Forms::Form { public: MyForm(HMODULE hLib) : hLib(hLib) { // Передаем дескриптор из main InitializeComponent(); // Загружаем функции из DLL f = (agp_c)GetProcAddress(hLib, "AGP_2D"); pStart = (start_workers)GetProcAddress(hLib, "AgpStartWorkers"); pFree = (free_agp)GetProcAddress(hLib, "AGP_Free"); pInit = (PInit)GetProcAddress(hLib, "AgpInit"); } protected: ~MyForm() { delete components; } private: HINSTANCE hLib = nullptr; typedef void(__cdecl* agp_c)(float, float, float, float, float, float, float, bool, float, float, float**, size_t*); typedef void(__cdecl* start_workers)(); typedef void(__cdecl* free_agp)(float*); typedef int(__cdecl* PInit)(int, float, float, float, float); agp_c f = nullptr; start_workers pStart = nullptr; free_agp pFree = nullptr; PInit pInit = nullptr; System::Windows::Forms::Button^ button1; System::Windows::Forms::TextBox^ textBox2; System::Windows::Forms::DataVisualization::Charting::Chart^ chart2; System::Windows::Forms::Label^ label2; System::Windows::Forms::Label^ label3; System::Windows::Forms::TextBox^ textBox1; System::Windows::Forms::TextBox^ textBox3; System::Windows::Forms::TextBox^ textBox4; System::Windows::Forms::TextBox^ textBox5; System::Windows::Forms::TextBox^ textBox6; System::Windows::Forms::Label^ label6; System::Windows::Forms::Label^ label7; System::Windows::Forms::Label^ label8; System::Windows::Forms::Label^ label9; System::Windows::Forms::Label^ label10; System::Windows::Forms::Label^ label1; System::Windows::Forms::TextBox^ textBox7; System::Windows::Forms::TextBox^ textBox8; System::ComponentModel::Container^ components; System::Void button1_Click(System::Object^ sender, System::EventArgs^ e) { chart2->Series[0]->Points->Clear(); chart2->Series[1]->Points->Clear(); chart2->Series[2]->Points->Clear(); chart2->Series[3]->Points->Clear(); static LARGE_INTEGER start, end; QueryPerformanceCounter(&start); float* buf = nullptr; size_t len = 0u; pStart(); // разбудили rank != 0 f(2.0f, 10000.0f, -2.2f, 1.8f, -2.2f, 1.8f, 2.5f, false, 0.0001f, GetTickCount(), &buf, &len); if (buf != nullptr) { std::vector<float> Extr_2D(buf, buf + len); pFree(buf); // Для 2D ветки: извлекаем schetchick и interval_len textBox7->Text = Convert::ToString(Extr_2D.back()); // schetchick Extr_2D.pop_back(); textBox6->Text = Convert::ToString(Extr_2D.back()); // interval_len Extr_2D.pop_back(); // Извлекаем последнюю тройку (x2, x1, f) float x2 = Extr_2D.back(); Extr_2D.pop_back(); float x1 = Extr_2D.back(); Extr_2D.pop_back(); float f_val = Extr_2D.back(); Extr_2D.pop_back(); // Выводим координаты и значение функции textBox4->Text = Convert::ToString(x1); textBox3->Text = Convert::ToString(x2); textBox2->Text = Convert::ToString(f_val); // Новое поле для значения функции // Отображаем точку на графике chart2->Series[2]->Points->AddXY(x1, x2); // Обрабатываем остальные точки (в порядке: f, x1, x2) while (!Extr_2D.empty()) { float x2_point = Extr_2D.back(); Extr_2D.pop_back(); float x1_point = Extr_2D.back(); Extr_2D.pop_back(); float f_point = Extr_2D.back(); Extr_2D.pop_back(); // Добавляем точку на график (только координаты) chart2->Series[0]->Points->AddXY(x1_point, x2_point); } } QueryPerformanceCounter(&end); LARGE_INTEGER freq; QueryPerformanceFrequency(&freq); const int multiplier = ((1ULL << 32) / freq.QuadPart) * 1'000'000ULL; textBox5->Text = Convert::ToString( ((end.QuadPart - start.QuadPart) * multiplier) >> 32 ) + " microseconds"; //textBox2->Text = Convert::ToString(Extr_2D.back()); //textBox1->Text = Convert::ToString(Extr_2D_LNA.back()); } Странно также то что и при завершении через 10 итераций и при завершении через 30 итераций минимум находится как -74.38 тогда как истинный минимум примерно -99, сама схема алгоритма не считая обмена между процессами точно корректна - проверял без MPI - минимум находится, я не стал делать схема с асинхронным обменом так как предполагаю что процессы работают примерно с одинаковой скоростью, ещё возможно сомнительный момент это пересчёт всех характеристик только когда новое M не просто лучше старого, но и лучше на определённый процент, думай как можно дольше чтобы установить в чём проблема и предложи мне конкретные шаги или как обнаружить проблему или конкретные исправления в коде, которые от неё избавятся

思考

Identifying bugs in code

So, there’s a discrepancy. The code uses ierre to receive from rank 1 even when rank is 0, which seems wrong. The partner exchange between ranks should be the focus, not just using iprobe() for a hardcoded rank (1). This creates a bug by not checking for the correct partner when probing for messages.

Examining termination logic

Okay, at the end of the loop, the code sends to rank 0, but rank 0 sending to itself is strange since iprobe only checks for messages from rank 1. This is problematic since iprobe won’t trigger properly. The issue arises because rank 0 only gets the irecv for rank 1, but not for itself.

Analyzing termination issue

Okay, in the termination check, both ranks use iprobe, but only rank 1 has posted a pending irecv. This makes things inconsistent, especially since iprobe on rank 0 could still check for messages from rank 1, even though rank 0 never posts an irecv. Also, it's unclear who sends the 'stop' signal.

Checking message communication flow

In the function MyForm::button1_Click, pStart() calls AgpStartWorkers() which sends a message to rank 1. Meanwhile, AgpWaitStartAndRun() listens for a message from rank 0. However, inside AGP_2D, there's an issue where the receive request (irecv) is set incorrectly on rank 1—it's expecting from source 1 instead of partner 0.

Analyzing early termination issue

The termination trigger seems to stem from iprobe(1,0), which only becomes true if a message from rank 1 is waiting—though the sendrecv matching process should have cleared any corresponding message. The issue arises from incorrect irecv expectations, where it checks source 1 instead of the partner rank. This mismatch causes an early detection of a message.

Identifying potential bugs and termination issue

The termination issue appears to be due to an improperly handled message, where iprobe(1, 0) on rank 0 detects a message that may not be intended, potentially from the first 'start' message. There's also a race condition involving recv_req and iprobe potentially using uninitialized data on rank 0, leading to undefined behavior or a crash. A fix involves replacing hardcoded values with the correct partner rank for message communication.

Identifying MPI-related bugs and improvements

The algorithm prematurely exits due to iprobe mistakenly matching exchange messages as termination signals, triggered by using the same tag. This causes issues where recv_req is uninitialized. Fixes include using distinct tags for different message types, ensuring correct initialization of recv_req, and using local recv_req for AGP_2D. Additionally, avoid global tag reuse and ensure MPI_Finalize is properly called.

Identifying bugs and fixes in AGP_2D

The issue with AGP_2D arises from misused global variables and improperly initialized recv_req, causing premature termination. Changes to improve the code include replacing hardcoded tags, differentiating message types with unique tags, adding persistent irecv for correct message matching, and refactoring the busy-wait into blocking recv.

Addressing termination bugs and memory alignment

The proposed changes focus on fixing the termination bug by ensuring proper handling of recv_req and symmetrical operations between ranks. Key updates include correcting message passing logic, preventing blocking on rank 0, and ensuring proper memory alignment when using SIMD operations. Additionally, adding a finalization function can solve the exit issue with MPI.

Thought for 2m 31s

回答

Коротко: у тебя пересекаются MPI-теги и перепутан источник в iprobe/irecv. В результате процесс видит «чужое» сообщение обмена (CrossMsg) как сигнал остановки и вываливается как раз «после первого обмена». Плюс есть гонка/UB: на ранге 0 ты вызываешь recv_req.wait() по неинициализированному request, что и даёт зависоны/0xc0000409. Ниже — минимальный набор правок, после которых поведение станет корректным.


Что именно ломает работу

  1. Один и тот же тег 0 используется для всего подряд:
  • старт воркера (AgpStartWorkers),
  • регулярный обмен sendrecv в AGP_2D,
  • «сигнал остановки».

На каждой итерации цикла ты делаешь в начале:

cpp
if (g_world->iprobe(1, 0)) { /* ... */ return; }

А чуть позже — обычный обмен данными:

cpp
g_world->sendrecv(partner, 0, outbound, partner, 0, inbound);

Если партнёр ушёл вперёд на пол-итерации, его «обменное» сообщение (тоже c тегом 0) уже висит в очереди, и твоё iprobe(1,0) срабатывает — алгоритм «думает», что пришла команда завершаться, и выходит. Именно поэтому при обмене каждые ~30 итераций выходишь на ~30-й, при «каждые 10» — на 10-й/11-й и т.п.

  1. Неправильный источник в iprobe/irecv:
    В начале AGP_2D:
cpp
if (rank) recv_req = g_world->irecv(1, 0, dummy); // ← source=1 захардкожен

и проверка:

cpp
if (g_world->iprobe(1, 0)) { ... recv_req.wait(); ... }
  • Для ранга 0 recv_req вообще не инициализируется, но в ветке iprobe ты всё равно вызываешь recv_req.wait() → UB/зависон/краш.
  • Для ранга 1 ожидание сообщения от самого себя (source=1) — очевидная логическая ошибка.
  1. Уведомление о завершении посылается всегда на ранг 0:
cpp
g_world->send(0, 0, 0);

Если завершился ранг 0, он шлёт сообщение самому себе; партнёр останавливается только «случайно», когда поймает чужое сообщение с тем же тегом. Это ещё один источник гонок.

  1. В «локальном» завершении ты не заполняешь out_data/out_len
    — на ранге 0 WinForms ждёт буфер, но AGP_2D может вернуться вообще без данных.

Минимально инвазивные исправления

Введи разные теги и всегда адресуйся к партнёру:

cpp
enum : int { TAG_START = 100, TAG_STOP = 101, TAG_XCHG = 102 };

Старт/пробуждение воркера

cpp
// AgpWaitStartAndRun() recv_req = g_world->irecv(0, TAG_START, dummy); // ждём только START ... if (recv_req.test()) { AGP_2D(...); recv_req = g_world->irecv(0, TAG_START, dummy); }
cpp
// AgpStartWorkers() g_world->isend(1, TAG_START, 0);

Обмен в 2D-алгоритме

В AGP_2D убери вообще глобальный recv_req и всю логику с irecv(1,0,...):

cpp
const int rank = g_world->rank(); const int partner = rank ^ 1; int dummy; // В НАЧАЛЕ каждой итерации — проверяем ТОЛЬКО TAG_STOP и ТОЛЬКО от партнёра: boost::mpi::status st; if (g_world->iprobe(partner, TAG_STOP, st)) { g_world->recv(partner, TAG_STOP, dummy); // поглощаем стоп // Отдаём результаты UI ПЕРЕД выходом Extr.emplace_back(static_cast<float>(schetchick)); Extr.emplace_back(interval_len); *out_len = Extr.size(); *out_data = (float*)CoTaskMemAlloc(sizeof(float) * (*out_len)); std::memcpy(*out_data, Extr.data(), sizeof(float) * (*out_len)); return; }

Обмен треками делай через свой тег:

cpp
g_world->sendrecv(partner, TAG_XCHG, outbound, partner, TAG_XCHG, inbound);

Локальное завершение — всегда формируем вывод и уведомляем партнёра:

cpp
if (interval_len < epsilon || schetchick == static_cast<int>(global_iterations)) { Extr.emplace_back(static_cast<float>(schetchick)); Extr.emplace_back(interval_len); *out_len = Extr.size(); *out_data = (float*)CoTaskMemAlloc(sizeof(float) * (*out_len)); std::memcpy(*out_data, Extr.data(), sizeof(float) * (*out_len)); g_world->isend(partner, TAG_STOP, 0); // не блокируемся return; }

Обрати внимание: здесь нет recv_req.wait() — он был источником зависаний/крашей на ранге 0. Мы либо «поглощаем» стоп немедленно через recv(...) после iprobe, либо сами рассылаем стоп партнёру и выходим.

(Опционально, но очень желательно) Завершение MPI

Ты создаёшь boost::mpi::environment и communicator, но нигде их не уничтожаешь → «process exited without calling finalize». Добавь:

cpp
extern "C" __declspec(dllexport) void AgpFinalize() noexcept { delete g_world; g_world = nullptr; delete g_env; g_env = nullptr; }

И вызови это один раз при закрытии приложения (из WinForms, перед FreeLibrary).


Почему после правок исчезнут «ранние выходы» и зависания

  • iprobe больше не видит сообщения обмена (они на TAG_XCHG) и реагирует только на явные «стопы» (TAG_STOP) от правильного источника.
  • Никаких wait() по неинициализированному request.
  • Процесс, который действительно решил завершиться, сообщает это партнёру, а не рангу 0 «вообще».
  • В любой ветке завершения ранга 0 ты возвращаешь буфер в UI, так что парсер WinForms получает валидные данные (и не падает).

Небольшие дополнительные замечания

  • В AGP_2D ты кладёшь schetchick/interval_len в буфер только в ветке с iprobe. Из-за этого раньше UI иногда получал «странные» пары значений. С правкой выше — порядок заполнения данных становится предсказуемым.
  • 0xC0000409 — это stack buffer overrun (защитный куки-чек), но в твоём случае он с высокой вероятностью возникал из-за UB вокруг recv_req.wait() и «сам себе послал, сам не принял» при send(0,0,0).
  • Твой расчёт T (период обмена) действительно выдаёт ~30 на старте (T ≈ 20*0.995 + 10), поэтому «после первого обмена» ≈ «на ~30-й итерации».

Быстрый чек-лист для локализации, если что-то останется:

  1. Включи логирование с рангом/тегом вокруг всех send/recv/sendrecv/iprobe (буквально printf("[r%d] iprobe src=%d tag=%d hit\n",rank,src,tag)), чтобы увидеть порядок.
  2. Временно замени iprobe+recv на блокирующий recv под условием «мы заранее точно знаем, что партнёр пошлёт стоп» — так проверишь, что обменные сообщения больше не задевают стоп-канал.
  3. Запусти с -mca pml ob1 -mca btl self,sm,tcp (OpenMPI) или аналог в MS-MPI и посмотрись, не «залипает» ли блокирующий send где-то ещё (после наших правок не должен).
  4. На время отладки дай уникальные теги из разных диапазонов (100/101/102 как выше) — это наглядно.

Если хочешь, могу прислать патч-дифф по конкретным файлам — но принципиально нужны именно эти четыре изменения: теги, источники, убрать recv_req.wait() из 2D и посылать STOP партнёру.

分享这个问答