diff --git a/Bin/Config/Lang/cn/nmNum_cn.qm b/Bin/Config/Lang/cn/nmNum_cn.qm index 5408661a..b6dbcf60 100644 Binary files a/Bin/Config/Lang/cn/nmNum_cn.qm and b/Bin/Config/Lang/cn/nmNum_cn.qm differ diff --git a/Bin/Config/Lang/cn/nmNum_cn.ts b/Bin/Config/Lang/cn/nmNum_cn.ts index 947891b9..302d0278 100644 --- a/Bin/Config/Lang/cn/nmNum_cn.ts +++ b/Bin/Config/Lang/cn/nmNum_cn.ts @@ -673,8 +673,8 @@ Reason: %1 LM 采样加密：%1 / %2 个目标点；全目标点误差：%3 - LM sampling: layered target points (%1 / %2); acceptance uses all valid target points. - LM 采样方式：目标点分层采样（%1 / %2 点），使用全部有效目标点验收。 + LM sampling: layered target points (%1 / %2); total-stage acceptance uses all valid target points. + LM 采样方式：目标点分层采样（%1 / %2 点），整体阶段使用全部有效目标点验收。 LM sampling: fixed 80 points (original mode). @@ -925,8 +925,8 @@ Reason: %1 === LM 自动拟合 - 局部最优 === - Max iterations reached. Best error: %1, Iterations: %2 - 达到最大迭代次数。最佳误差：%1，迭代次数：%2 + Total-stage budget reached. Best error: %1, Cumulative iterations: %2 + 整体阶段预算已用尽。最佳误差：%1，累计迭代次数：%2 === LM AUTOMATIC FITTING - MAX ITERATIONS === @@ -968,10 +968,6 @@ Reason: %1 === Starting LM Main Loop === === 开始 LM 主循环 === - - LM starting point error: %1; evaluation budget: %2 - LM 起点误差：%1；最大评估次数：%2 - No effective improvement for %1 consecutive steps; rebuilding sensitivity model for confirmation 连续 %1 次无有效改善，正在重建灵敏度模型进行确认 @@ -1065,8 +1061,8 @@ Reason: %1 检测到局部最优 - Maximum iterations reached - 达到最大迭代次数 + Total-stage iteration or evaluation budget reached + 整体阶段迭代或评估预算已用尽 Stopped by user request @@ -1124,10 +1120,6 @@ Reason: %1 Permeability alignment: k=%1, height error=%2, shape error=%3, result=%4 渗透率对齐：k=%1，上下误差=%2，形状误差=%3，结果=%4 - - LM stage 2: optimize pressure and derivative shape; stop after 3 ineffective steps. - LM 阶段二：优先调整压力和压力导数形状，连续 3 步无明显改善后切换。 - LM stage 3: original LM fitting; accept by total error only. LM 阶段三：按原有 LM 拟合，仅依据整体误差接受调整。 @@ -1137,24 +1129,12 @@ Reason: %1 采用灵敏度试算点：整体误差=%1 - Shape stage ended: %1 - 形状阶段结束：%1 - - - 3 consecutive steps without effective shape improvement - 连续 3 步形状没有明显改善 - - - reserve remaining iterations and evaluations for total fitting - 为整体拟合保留剩余迭代和求解预算 - - - no valid shape sensitivity model - 没有有效的形状灵敏度模型 + Sensitivity probe accepted: early relative-slope matching error=%1 + 灵敏度试算点已接受：前期相对斜率匹配误差=%1 - no feasible shape descent step - 没有满足约束的形状下降步 + Shape stage ended: %1 + 形状阶段结束：%1 Permeability alignment ended: %1 @@ -1232,6 +1212,122 @@ Reason: %1 inaccurate model at minimum trust radius 最小信赖半径下仍连续预测失准 + + no remaining shape parameters + 没有其他勾选参数需要调整形状 + + + Early adjustment: adjust wellbore storage, then skin, to match the early pressure-derivative slope difference of the target. + 前期预调整：先调整井储，再调整表皮，使模拟压力与导数的前期斜率差接近目标曲线。 + + + Early adjustment: adjust skin to match the early pressure-derivative slope difference of the target. + 前期预调整：调整表皮，使模拟压力与导数的前期斜率差接近目标曲线。 + + + Wellbore storage adjustment ended: %1; now adjusting skin. + 井储调整结束：%1；开始单独调整表皮。 + + + Early wellbore adjustment ended: %1 + 井储和表皮前期预调整结束：%1 + + + Shape fitting continues with storage and skin fixed; a conditional wellbore recheck follows. + 继续调整形状，暂时固定井储和表皮；随后按条件进行一次井储表皮回检。 + + + 2 consecutive steps without effective early improvement + 连续 2 步前期误差没有明显改善 + + + no feasible early adjustment direction or parameter at bound + 没有明确可行的前期调整方向，或参数已到边界 + + + wellbore storage relative-slope matching + 井储前期相对斜率匹配 + + + skin relative-slope matching + 表皮前期相对斜率匹配 + + + Wellbore storage trial: relative-slope matching error=%1, C=%2 -> %3 + 井储试调：相对斜率匹配误差=%1，井储=%2 → %3 + + + Wellbore recheck: early error=%1, total limit=%2, shape limit=%3. + 井储表皮回检：前期相对斜率匹配误差=%1，整体误差上限=%2，形状误差上限=%3。 + + + early relative-slope error restored within tolerance + 前期相对斜率匹配误差已恢复到容差范围内 + + + new sensitivity for wellbore recheck + 为井储表皮回检重建灵敏度 + + + wellbore recheck completed during sensitivity evaluation + 井储表皮回检在灵敏度评价期间完成 + + + total target reached during shape fitting + 形状调整期间已达到总误差目标 + + + Half-length exploration accepted: L=%1, shape error=%2, total error=%3. + 裂缝半长探索结果已接受：半长=%1，形状误差=%2，整体误差=%3。 + + + fresh sensitivity after half-length exploration + 裂缝半长探索后重建灵敏度 + + + fresh sensitivity at joint shape entry + 进入联合形状调整时重建灵敏度 + + + Adaptive fitting counts: %1 completed shape rounds, %2 total-stage iterations, %3 total evaluations. + 自适应拟合统计：已完成 %1 轮形状搜索，整体阶段迭代 %2 次，累计评估 %3 次。 + + + LM stage 2: adaptive shape search; confirm stagnation after 2 complete rounds without significant improvement. + LM 阶段2：自适应形状搜索；连续两轮完整搜索无显著改善后确认停滞。 + + + LM starting point error: %1; independent total-stage evaluation budget: %2 + LM 起点误差：%1；整体阶段独立评估预算：%2 + + + Shape search round %1: shape error=%2, improvement=%3, required=%4. + 形状搜索第 %1 轮：形状误差=%2，改善量=%3，所需改善量=%4。 + + + Total-stage budget starts now: %1 iterations, %2 evaluations; pre-adjustment is counted separately. + 整体阶段预算开始计数：%1 次迭代、%2 次评估；预调整单独计数。 + + + Unable to build a valid shape sensitivity model. + 无法建立有效的形状灵敏度模型。 + + + confirm stagnation after 2 complete shape rounds + 完成两轮形状搜索后确认停滞 + + + full sensitivity at total-stage entry + 进入整体阶段时建立完整灵敏度 + + + no significant improvement in a complete round after fresh sensitivity confirmation + 重建灵敏度确认后，完整一轮搜索仍无显著改善 + + + no valid early sensitivity model + 无有效的早期灵敏度模型 + nmCalculationSolver @@ -4499,6 +4595,10 @@ Supported types: Vertical, Vertical Fractured, and Horizontal Multi-Fractured We Fitting Curve 拟合曲线 + + Pre-adjustment + 预调整中 + nmWxChangeAnal diff --git a/Include/nmNum/nmCalculation/nmCalculationAutoFitLM.h b/Include/nmNum/nmCalculation/nmCalculationAutoFitLM.h index c700c7ab..579abe33 100644 --- a/Include/nmNum/nmCalculation/nmCalculationAutoFitLM.h +++ b/Include/nmNum/nmCalculation/nmCalculationAutoFitLM.h @@ -48,6 +48,13 @@ struct AutoFitObjectiveBreakdownLM { bool registrationAmbiguous; // 固定 log-time 网格上的双曲线斜率残差，独立于数值残差的采样层级。 QVector shapeResiduals; + // 前期数值残差为 81 点，平行程度残差为 73 个斜率区间，均含第一窗口权重。 + // 平行误差比较模拟与目标各自的压力—导数斜率差，不要求模拟自身斜率差为零。 + QVector earlyValueResiduals; + QVector earlyParallelResiduals; + double earlyValueLoss; + double earlyParallelLoss; + double earlyParallelBias; double pressureVerticalBias; double derivativeVerticalBias; double shapeLoss; @@ -70,6 +77,9 @@ struct AutoFitObjectiveBreakdownLM { , horizontalLoss(std::numeric_limits::quiet_NaN()) , horizontalReliable(false) , registrationAmbiguous(false) + , earlyValueLoss(std::numeric_limits::quiet_NaN()) + , earlyParallelLoss(std::numeric_limits::quiet_NaN()) + , earlyParallelBias(std::numeric_limits::quiet_NaN()) , pressureVerticalBias(0.0) , derivativeVerticalBias(0.0) , shapeLoss(std::numeric_limits::quiet_NaN()) diff --git a/Src/nmNum/nmCalculation/nmCalculationAutoFitLM.cpp b/Src/nmNum/nmCalculation/nmCalculationAutoFitLM.cpp index b932eaa7..7ff0db33 100644 --- a/Src/nmNum/nmCalculation/nmCalculationAutoFitLM.cpp +++ b/Src/nmNum/nmCalculation/nmCalculationAutoFitLM.cpp @@ -316,14 +316,25 @@ static bool trustRegionResidualsValid( for(int i = 0; i < breakdown.shapeResiduals.size(); ++i) { if(!isFiniteNumber(breakdown.shapeResiduals[i])) return false; } + if(breakdown.earlyValueResiduals.size() != 162 || breakdown.earlyParallelResiduals.size() != 146 || + !isFiniteNumber(breakdown.earlyValueLoss) || !isFiniteNumber(breakdown.earlyParallelLoss) || + !isFiniteNumber(breakdown.earlyParallelBias)) return false; + for(int i = 0; i < breakdown.earlyValueResiduals.size(); ++i) { + if(!isFiniteNumber(breakdown.earlyValueResiduals[i])) return false; + } + for(int i = 0; i < breakdown.earlyParallelResiduals.size(); ++i) { + if(!isFiniteNumber(breakdown.earlyParallelResiduals[i])) return false; + } return true; } -// 两类残差在同一次真实评价中获得，联合缓存使阶段切换不必重算灵敏度。 +// 数值、形状和前期平行程度残差共用一次求解，联合缓存供各子阶段复用。 static QVector trustRegionFullResidual(const AutoFitObjectiveBreakdownLM& objective) { QVector residual = objective.residualVector; residual += objective.shapeResiduals; + residual += objective.earlyValueResiduals; + residual += objective.earlyParallelResiduals; return residual; } @@ -337,6 +348,70 @@ static double trustRegionSquaredNorm(const QVector& values) return sum; } +// 第一窗口的形状能量与局部 Fisher 使用相同的中心时间和重叠权重， +// 压力、导数同时参与；不重新插值或调用求解器。 +static double trustRegionEarlyShapeEnergy(const QVector& residual) +{ + const int pointCount = residual.size() / 2; + double energy = 0.0; + for(int row = 0; row < residual.size(); ++row) { + const double coordinate = (row % pointCount + 4.0) / 80.0; + energy += autoFitTimeWindowWeight(coordinate, 0) * residual[row] * residual[row]; + } + return energy; +} + +// 比较模拟与目标各自的压力—导数斜率差；误差为零表示两组曲线的相对走势一致。 +// 单独平移任一曲线不改变此指标；前期数值误差仅作诊断，不参与井储、表皮验收。 +static void populateEarlyWellboreMetrics(AutoFitObjectiveBreakdownLM* objective, + const QVector targetLogs[2], const QVector resultLogs[2], double logTimeSpan) +{ + const int count = 81; + QVector weights(count, 0.0); + double weightSum = 0.0; + for(int i = 0; i < count; ++i) { + weights[i] = autoFitTimeWindowWeight(i / 80.0, 0) * ((i == 0 || i == count - 1) ? 0.5 : 1.0); + weightSum += weights[i]; + } + objective->earlyValueResiduals.fill(0.0, 2 * count); + for(int i = 0; i < count; ++i) { + const double weight = weights[i] / weightSum; + for(int component = 0; component < 2; ++component) { + const int row = component * count + i; + const double scale = qSqrt(0.5 * weight); + objective->earlyValueResiduals[row] = scale * (resultLogs[component][i] - targetLogs[component][i]); + } + } + objective->earlyValueLoss = qSqrt(trustRegionSquaredNorm(objective->earlyValueResiduals)); + + // 复用形状指标的 8 点跨度，窗口权重取斜率区间中心，减少相邻点噪声。 + const int lag = 8; + const int slopeCount = count - lag; + const double logTimeStep = logTimeSpan * lag / (count - 1); + QVector parallelWeights(slopeCount, 0.0); + double parallelWeightSum = 0.0; + for(int i = 0; i < slopeCount; ++i) { + parallelWeights[i] = autoFitTimeWindowWeight((i + lag * 0.5) / (count - 1), 0); + parallelWeightSum += parallelWeights[i]; + } + objective->earlyParallelResiduals.fill(0.0, 2 * slopeCount); + objective->earlyParallelBias = 0.0; + for(int i = 0; i < slopeCount; ++i) { + const double resultSlopeDifference = ((resultLogs[0][i + lag] - resultLogs[0][i]) - + (resultLogs[1][i + lag] - resultLogs[1][i])) / logTimeStep; + const double targetSlopeDifference = ((targetLogs[0][i + lag] - targetLogs[0][i]) - + (targetLogs[1][i + lag] - targetLogs[1][i])) / logTimeStep; + const double relativeSlopeError = resultSlopeDifference - targetSlopeDifference; + const double weight = parallelWeights[i] / parallelWeightSum; + objective->earlyParallelBias += weight * relativeSlopeError; + // 两半各占一半能量以兼容窗口 Fisher；保留相对目标的斜率差残差符号。 + const double residual = qSqrt(0.5 * weight) * relativeSlopeError; + objective->earlyParallelResiduals[i] = residual; + objective->earlyParallelResiduals[slopeCount + i] = residual; + } + objective->earlyParallelLoss = qSqrt(trustRegionSquaredNorm(objective->earlyParallelResiduals)); +} + // 计算同维向量内积；维度不一致表示局部模型无效，返回零让调用方放弃修正。 static double trustRegionDotProduct(const QVector& left, const QVector& right) @@ -421,6 +496,25 @@ static bool solveTrustRegionLinearSystem( return true; } +// 半长候选只由当前值和用户边界生成；投影后去重，不把某个半长写成目标值。 +static QVector trustRegionShapeLengthTrials(double value, double lower, double upper) +{ + QVector trials; + if(!isFiniteNumber(value) || !isFiniteNumber(lower) || !isFiniteNumber(upper) || upper <= lower) + return trials; + const double raw[] = {value * 2.0, value * 0.5, value * 2.0 < upper ? upper : lower}; + for(int i = 0; i < 3; ++i) { + const double trial = qBound(lower, raw[i], upper); + const double tolerance = 1.0e-12 * qMax(1.0, qAbs(trial)); + if(qAbs(trial - value) <= tolerance) continue; + bool duplicate = false; + for(int j = 0; j < trials.size(); ++j) + if(qAbs(trial - trials[j]) <= tolerance) duplicate = true; + if(!duplicate) trials.append(trial); + } + return trials; +} + // Fisher 仅作为当前归一化坐标下的局部信息矩阵，不用于统计置信区间。 struct TrustRegionFisher { @@ -595,43 +689,76 @@ static bool buildTrustRegionFisherStep( (*step)[selected[i]] = solution[i]; } } - double norm = qSqrt(trustRegionSquaredNorm(*step)); - if(!solved || !isFiniteNumber(norm) || norm < minimumStep) { - step->fill(0.0, dimensions); - for(int i = 0; i < count; ++i) { - const int p = selected[i]; - double direction = -global.gradient[p]; - if((coordinates[p] <= minimumStep && direction < 0.0) || - (coordinates[p] >= 1.0 - minimumStep && direction > 0.0)) { - direction = 0.0; - } - (*step)[p] = direction; + auto projectAndPredict = [&]() -> bool { + const double norm = qSqrt(trustRegionSquaredNorm(*step)); + if(!isFiniteNumber(norm) || norm < minimumStep) return false; + if(norm > trustRadius) { + for(int p = 0; p < dimensions; ++p) (*step)[p] *= trustRadius / norm; } - norm = qSqrt(trustRegionSquaredNorm(*step)); - if(norm > minimumStep) { - for(int p = 0; p < dimensions; ++p) { - (*step)[p] *= trustRadius / norm; - } - norm = trustRadius; + for(int p = 0; p < dimensions; ++p) { + (*step)[p] = qBound(0.0, coordinates[p] + (*step)[p], 1.0) - coordinates[p]; } - } - if(norm > trustRadius) { + // 只比较投影后可执行步长的下降，阻尼项不属于真实拟合目标。 + *predictedReduction = -trustRegionDotProduct(global.gradient, *step); for(int p = 0; p < dimensions; ++p) { - (*step)[p] *= trustRadius / norm; + for(int q = 0; q < dimensions; ++q) { + *predictedReduction -= 0.5 * (*step)[p] * global.matrix[p][q] * (*step)[q]; + } } + return isFiniteNumber(*predictedReduction) && *predictedReduction > 1.0e-14 && + qSqrt(trustRegionSquaredNorm(*step)) >= minimumStep; + }; + if(solved && projectAndPredict()) return true; + + // 联合解可能被边界投影破坏；此时尝试可行梯度方向，不直接判为没有下降方向。 + step->fill(0.0, dimensions); + for(int i = 0; i < count; ++i) { + const int p = selected[i]; + double direction = -global.gradient[p]; + if((coordinates[p] <= minimumStep && direction < 0.0) || + (coordinates[p] >= 1.0 - minimumStep && direction > 0.0)) direction = 0.0; + (*step)[p] = direction; } + const double norm = qSqrt(trustRegionSquaredNorm(*step)); + if(!isFiniteNumber(norm) || norm < minimumStep) return false; for(int p = 0; p < dimensions; ++p) { - (*step)[p] = qBound(0.0, coordinates[p] + (*step)[p], 1.0) - coordinates[p]; + (*step)[p] = qBound(0.0, coordinates[p] + (*step)[p] * trustRadius / norm, 1.0) - coordinates[p]; } - // 只比较实际可执行步长的全局预测下降，阻尼项不属于真实拟合目标。 - *predictedReduction = -trustRegionDotProduct(global.gradient, *step); + const double linearReduction = -trustRegionDotProduct(global.gradient, *step); + double curvature = 0.0; for(int p = 0; p < dimensions; ++p) { for(int q = 0; q < dimensions; ++q) { - *predictedReduction -= 0.5 * (*step)[p] * global.matrix[p][q] * (*step)[q]; + curvature += (*step)[p] * global.matrix[p][q] * (*step)[q]; } } - return isFiniteNumber(*predictedReduction) && *predictedReduction > 1.0e-14 && - qSqrt(trustRegionSquaredNorm(*step)) >= minimumStep; + if(!isFiniteNumber(linearReduction) || !isFiniteNumber(curvature) || linearReduction <= 0.0) return false; + // 沿投影梯度最小化局部二次模型，只缩步，不突破信赖域或参数边界。 + const double scale = curvature > 0.0 ? qMin(1.0, linearReduction / curvature) : 1.0; + for(int p = 0; p < dimensions; ++p) (*step)[p] *= scale; + return projectAndPredict(); +} + +// 放大已被真实结果验证可靠的联合方向，仍受最大半径和参数边界限制。 +static bool buildExpandedTrustRegionStep(const TrustRegionFisher& information, + const QVector& coordinates, const QVector& originalStep, + double maximumRadius, QVector* expandedStep, double* prediction) +{ + const double norm = qSqrt(trustRegionSquaredNorm(originalStep)); + if(norm <= 1.0e-12) return false; + const double scale = qMin(2.0, maximumRadius / norm); + if(scale <= 1.01) return false; + expandedStep->resize(originalStep.size()); + double difference = 0.0; + for(int i = 0; i < originalStep.size(); ++i) { + (*expandedStep)[i] = qBound(0.0, coordinates[i] + scale * originalStep[i], 1.0) - coordinates[i]; + difference += qAbs((*expandedStep)[i] - originalStep[i]); + } + if(difference <= 1.0e-10) return false; + *prediction = -trustRegionDotProduct(information.gradient, *expandedStep); + for(int i = 0; i < expandedStep->size(); ++i) + for(int j = 0; j < expandedStep->size(); ++j) + *prediction -= 0.5 * (*expandedStep)[i] * information.matrix[i][j] * (*expandedStep)[j]; + return isFiniteNumber(*prediction) && *prediction > 1.0e-14; } // 每次得到有效真实候选后，使用满足最新割线条件的秩一修正更新完整残差 @@ -1004,7 +1131,8 @@ void nmCalculationAutoFitLM::writeTraceHeader() cols << "sampling_mode" << "sampling_stride" << "sampling_points" << "full_target_points" << "layer_objective" - << "pressure_vertical_bias" << "derivative_vertical_bias"; + << "pressure_vertical_bias" << "derivative_vertical_bias" << "first_window_shape_loss" + << "early_value_loss" << "early_parallel_loss" << "early_parallel_bias"; QTextStream out(&m_traceFile); out << cols.join(",") << "\n"; } @@ -1041,11 +1169,29 @@ void nmCalculationAutoFitLM::writeTraceMetaFile() QTextStream out(&metaFile); out << "{\n"; - out << " \"schema_version\": 10,\n"; - out << " \"strategy\": \"permeability_height_then_shape_then_original_lm\",\n"; + out << " \"schema_version\": 20,\n"; + out << " \"strategy\": \"permeability_height_then_shape_then_joint_lm\",\n"; + out << " \"shape_priority\": \"storage_then_skin_then_shape_without_wellbore_then_optional_wellbore_recheck\",\n"; out << " \"shape_metric\": \"pressure_and_derivative_log_slopes_81_points_lag_8\",\n"; + out << " \"early_parallel_metric\": \"pressure_derivative_log_slope_difference_vs_target_81_points_lag_8\",\n"; + out << " \"early_wellbore_acceptance\": \"relative_slope_matching_decrease_only\",\n"; + out << " \"early_wellbore_total_constraint\": false,\n"; out << " \"shape_stage_total_tolerance\": \"max(0.02, 25% of stage entry total)\",\n"; + out << " \"shape_wellbore_parameters_frozen\": true,\n"; + out << " \"shape_sensitivity_refresh\": \"full Jacobian at joint shape entry and after accepted length exploration\",\n"; + out << " \"shape_length_exploration\": \"each complete shape round: after 3 ordinary candidates or before round completion; bounded 2x, 0.5x and far-bound trials with optional permeability correction\",\n"; + out << " \"shape_round_extra_evaluation_limit\": 8,\n"; + out << " \"shape_expansion_policy\": \"Dfc joint step up to 2x when actual rho exceeds 0.75; retain ordinary step on failure\",\n"; + out << " \"stage2_budget\": \"adaptive, no fixed iteration or evaluation quota; full window/global/length rounds, 2 ineffective rounds then fresh-J confirmation\",\n"; + out << " \"shape_round_improvement_threshold\": \"max(0.0001, 1% of round entry shape loss)\",\n"; + out << " \"iteration_count_scope\": \"max_iterations applies to total-stage LM; trace iteration remains cumulative\",\n"; + out << " \"wellbore_recheck_trigger\": \"early loss exceeds initial early-stage exit by max(0.0001, 10%)\",\n"; + out << " \"wellbore_recheck_acceptance\": \"early decrease with fixed entry total and shape limits, each max(0.0001, 5%)\",\n"; + out << " \"wellbore_recheck_budget\": \"once when early shape degrades; each wellbore parameter stops after 2 ineffective steps, independent of total-stage budget\",\n"; out << " \"total_stage_shape_constraint\": false,\n"; + out << " \"total_parameter_selection\": \"all_valid_free_columns\",\n"; + out << " \"skin_difference_policy\": \"local_scale_in_all_stages\",\n"; + out << " \"difference_failure_policy\": \"halve_before_opposite_direction\",\n"; out << " \"trace_type\": \"finite_difference_lm_trust_region\",\n"; out << " \"run_id\": " << jsonEscape(m_traceRunId) << ",\n"; out << " \"created_at\": " @@ -1163,7 +1309,12 @@ void nmCalculationAutoFitLM::writeTraceRow( << (objectiveBreakdown ? QString::number(objectiveBreakdown->fullPointCount) : QString()) << (objectiveBreakdown ? traceNumber(objectiveBreakdown->layerError) : QString()) << (objectiveBreakdown ? traceNumber(objectiveBreakdown->pressureVerticalBias) : QString()) - << (objectiveBreakdown ? traceNumber(objectiveBreakdown->derivativeVerticalBias) : QString()); + << (objectiveBreakdown ? traceNumber(objectiveBreakdown->derivativeVerticalBias) : QString()) + << (objectiveBreakdown && objectiveBreakdown->valid + ? traceNumber(qSqrt(trustRegionEarlyShapeEnergy(objectiveBreakdown->shapeResiduals))) : QString()) + << (objectiveBreakdown && objectiveBreakdown->valid ? traceNumber(objectiveBreakdown->earlyValueLoss) : QString()) + << (objectiveBreakdown && objectiveBreakdown->valid ? traceNumber(objectiveBreakdown->earlyParallelLoss) : QString()) + << (objectiveBreakdown && objectiveBreakdown->valid ? traceNumber(objectiveBreakdown->earlyParallelBias) : QString()); QTextStream out(&m_traceFile); out << cols.join(",") << "\n"; m_traceFile.flush(); @@ -1628,7 +1779,7 @@ bool nmCalculationAutoFitLM::startAutoFitting() emit logMessageGenerated(tr("=== LM AUTOMATIC FITTING - LOCAL OPTIMUM ===")); } else if(finalReason == LM_MAX_ITERATIONS) { success = true; - message = QString(tr("Max iterations reached. Best error: %1, Iterations: %2")) + message = QString(tr("Total-stage budget reached. Best error: %1, Cumulative iterations: %2")) .arg(m_globalBestFitness, 0, 'e', 4) .arg(m_currentIteration + 1); emit logMessageGenerated(tr("=== LM AUTOMATIC FITTING - MAX ITERATIONS ===")); @@ -1658,7 +1809,6 @@ bool nmCalculationAutoFitLM::startAutoFitting() } emitRunSummary(success, finalReason); - emit progressUpdated(m_maxIterations, m_globalBestFitness); QApplication::processEvents(); msleep(200); QApplication::processEvents(); @@ -1801,11 +1951,11 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() return LM_OPTIMIZATION_FAILED; } - // 真实求解次数比“外层迭代次数”更能反映耗时。预算至少允许完成一次全参数 - // 灵敏度和两次候选评价，同时避免连续重建 Jacobian 导致运行时间失控。 - const int maximumEvaluations = qMax( - m_totalEvaluations + dimensions + 2, - qMax(20, m_maxIterations * 3)); + // 用户迭代设置只约束第三阶段；真实评价额度也在进入该阶段时独立起算。 + // 前期调整按改善情况结束，不能提前消耗掉最终联合 LM 的迭代和求解机会。 + const int totalEvaluationBudget = qMax(dimensions + 2, qMax(20, m_maxIterations * 3)); + int maximumEvaluations = m_totalEvaluations + totalEvaluationBudget; + int totalIterations = 0; // 下列步长均位于归一化内部坐标：0.04 表示参数范围的 4%，信赖半径 // 限制一次联合移动的二范数，相关性门槛用于排除响应近乎共线的参数。 const double sensitivityStep = 0.04; @@ -1938,11 +2088,12 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() } restoreEvaluationState(current); + emit progressUpdated(-1, current.fitness); emit logMessageGenerated(tr("=== Starting LM Main Loop ===")); emit logMessageGenerated( - tr("LM starting point error: %1; evaluation budget: %2") + tr("LM starting point error: %1; independent total-stage evaluation budget: %2") .arg(current.fitness, 0, 'e', 4) - .arg(maximumEvaluations)); + .arg(totalEvaluationBudget)); // 高度预调整只改变用户勾选的渗透率。ln(k) 的初始变化由有符号高度差 // 给出，真实求解后用割线估计修正；拒绝时缩步，不让其他参数补偿高度。 @@ -2010,18 +2161,57 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() if(m_shouldStop) return LM_USER_STOPPED; - // 形状阶段不设达标阈值；连续三步无有效改善或耗用约三分之一预算即切换。 - // 允许整体误差适度回升，但上限固定在高度对齐后，禁止逐步放宽导致漂移。 - bool shapeStage = m_maxIterations >= 3 && current.fitness >= m_targetError; - const double shapeValueLimit = current.fitness + qMax(0.02, 0.25 * current.fitness); - const int shapeEvaluationDeadline = m_totalEvaluations + qMax(0, - (maximumEvaluations - m_totalEvaluations - dimensions - 2) / 3); - const int shapeIterationLimit = qMax(1, m_maxIterations / 3); + // 第二阶段不按迭代比例截断；普通形状按完整搜索轮次判断是否仍有显著改善。 + // 井储、表皮保持冻结，总误差上限仍固定于普通形状阶段入口。 + bool shapeStage = current.fitness >= m_targetError; + double shapeValueLimit = current.fitness + qMax(0.02, 0.25 * current.fitness); + bool jointShapeStarted = false; + double shapeRoundBaseline = current.breakdown.shapeLoss; + int shapeRoundNumber = 0; + int ineffectiveShapeRounds = 0; + bool shapeConfirmationRequested = false; + int shapeCandidateCount = 0; + const int fractureLengthColumn = m_enabledParamIndices.indexOf(6); + const int conductivityColumn = m_enabledParamIndices.indexOf(5); + // 每轮探索有限个候选，但不会因全流程累计次数截断半长与渗透率组合。 + const int maximumShapeExplorationEvaluations = 8; + int shapeExplorationEvaluations = 0; + bool shapeExplorationDone = fractureLengthColumn < 0; + const int storageColumn = m_enabledParamIndices.indexOf(2); + const int skinColumn = m_enabledParamIndices.indexOf(1); + bool hasRemainingShapeParameters = false; + for(int column = 0; column < dimensions; ++column) { + if(m_enabledParamIndices[column] != 1 && m_enabledParamIndices[column] != 2) + hasRemainingShapeParameters = true; + } + const int wellboreParameterCount = (storageColumn >= 0 ? 1 : 0) + (skinColumn >= 0 ? 1 : 0); + bool wellboreRecheck = false; + bool wellboreRecheckDone = false; + double earlyWellboreBaseline = current.breakdown.earlyParallelLoss; + double recheckTotalLimit = 0.0; + double recheckShapeLimit = 0.0; + // 先井储、再表皮；回检沿用相同顺序和前期指标，不回到形状阶段反复循环。 + bool earlyShapeStage = shapeStage && wellboreParameterCount > 0; + int earlyShapeParameter = storageColumn >= 0 ? 2 : 1; + auto parameterAllowedInStage = [&](int column) -> bool { + const int parameterIndex = m_enabledParamIndices[column]; + // 候选和缓存灵敏度探针共用此限制，防止探针绕过形状阶段的冻结。 + if(!shapeStage) return true; + if(earlyShapeStage) return parameterIndex == earlyShapeParameter; + return parameterIndex != 1 && parameterIndex != 2; + }; auto stageError = [&](const TrustRegionEvaluation& point) -> double { + if(earlyShapeStage) return point.breakdown.earlyParallelLoss; return shapeStage ? point.breakdown.shapeLoss : point.fitness; }; auto acceptable = [&](const TrustRegionEvaluation& point, const TrustRegionEvaluation& base) -> bool { if(!point.valid) return false; + if(earlyShapeStage) { + // 初调保持原规则；回检还需保护整轮入口的总误差和形状，不能逐步放宽。 + return point.breakdown.earlyParallelLoss < base.breakdown.earlyParallelLoss && + (!wellboreRecheck || (point.fitness <= recheckTotalLimit && + point.breakdown.shapeLoss <= recheckShapeLimit)); + } return shapeStage ? point.breakdown.shapeLoss < base.breakdown.shapeLoss && point.fitness <= shapeValueLimit // 预调整结束后恢复原 LM：有效候选只按整体误差下降接受。 @@ -2033,8 +2223,13 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() restoreEvaluationState(current); }; emit logMessageGenerated(shapeStage - ? tr("LM stage 2: optimize pressure and derivative shape; stop after 3 ineffective steps.") + ? tr("LM stage 2: adaptive shape search; confirm stagnation after 2 complete rounds without significant improvement.") : tr("LM stage 3: original LM fitting; accept by total error only.")); + if(earlyShapeStage) { + emit logMessageGenerated(earlyShapeParameter == 2 + ? tr("Early adjustment: adjust wellbore storage, then skin, to match the early pressure-derivative slope difference of the target.") + : tr("Early adjustment: adjust skin to match the early pressure-derivative slope difference of the target.")); + } // 有效改善始终相对“上一次有效改善后的误差”累计判断，避免一连串微小 // 下降每次都清零计数；累计达到门槛后才开始新的有效改善基准。 @@ -2085,21 +2280,18 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() }; emit logMessageGenerated(m_layeredSampling - ? tr("LM sampling: layered target points (%1 / %2); acceptance uses all valid target points.") + ? tr("LM sampling: layered target points (%1 / %2); total-stage acceptance uses all valid target points.") .arg(current.breakdown.residualVector.size() / 2).arg(current.breakdown.fullPointCount) : tr("LM sampling: fixed 80 points (original mode).")); - auto enterTotalStage = [&](const QString& reason) { - // 保留当前曲线和完整 J，只重置阶段停滞状态、阻尼及信赖半径。 + auto reuseSensitivityForNextStage = [&]() { + // 更换目标时复用完整 J，清除旧目标下的拒绝和预测失准计数。 rebuildRequested = jacobian.isEmpty() || consecutiveSolverFailures > 0 || acceptedSinceRebuild >= 10; if(jacobian.isEmpty()) rebuildReason = QT_TR_NOOP("no valid sensitivity model"); else if(consecutiveSolverFailures > 0) rebuildReason = QT_TR_NOOP("solver failure before stage switch"); else if(acceptedSinceRebuild >= 10) rebuildReason = QT_TR_NOOP("10 accepted steps since last rebuild"); - // 阶段目标已改变，不把形状目标下的预测失准累计到整体阶段。 consecutivePoorPredictions = 0; modelRebuiltAtMinimumRadius = false; - shapeStage = false; - effectiveImprovementBaseline = current.fitness; consecutiveIneffectiveSteps = 0; consecutiveRejectedSteps = 0; consecutiveSolverFailures = 0; @@ -2108,13 +2300,205 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() globalFallbackAttempted = false; trustRadius = 0.12; damping = 0.01; + }; + auto exploreShapeLengths = [&]() -> bool { + if(shapeExplorationDone || !jointShapeStarted || earlyShapeStage || wellboreRecheck) + return false; + const TrustRegionEvaluation base = current; + TrustRegionEvaluation best; + const QVector trials = trustRegionShapeLengthTrials(base.parameters[fractureLengthColumn], + m_parameterLower[6], m_parameterUpper[6]); + const int explorationLimit = qMin(maximumShapeExplorationEvaluations, shapeExplorationEvaluations + 6); + auto canEvaluate = [&]() -> bool { + return shapeExplorationEvaluations < explorationLimit && processPauseAndStop(); + }; + auto keepBetter = [&](const TrustRegionEvaluation& point) { + if(acceptable(point, base) && (!best.valid || point.breakdown.shapeLoss < best.breakdown.shapeLoss || + (point.breakdown.shapeLoss == best.breakdown.shapeLoss && point.fitness < best.fitness))) best = point; + }; + for(int i = 0; i < trials.size() && canEvaluate(); ++i) { + TrustRegionEvaluation probe; + probe.parameters = base.parameters; + probe.parameters[fractureLengthColumn] = trials[i]; + probe.coordinates = coordinatesFromParameters(probe.parameters); + ++shapeExplorationEvaluations; + probe.valid = evaluateTrustRegionPoint(probe.parameters, &probe.fitness, + &probe.breakdown, &probe.curve, &probe.elapsedMs); + writeTraceRow(m_currentIteration, fractureLengthColumn, "shape_length_probe", probe.parameters, + probe.fitness, probe.valid, probe.elapsedMs, probe.valid ? "probe_valid" : "solver_invalid", + probe.valid ? &probe.breakdown : nullptr); + if(probe.valid) { + keepBetter(probe); + // 半长点即使暂时超出总误差上限，仍允许一次渗透率补偿后验收组合。 + // 沿用高度预调的有符号偏差与初始斜率 -1，远点修正必须真实求解验证。 + if(permeabilityColumn >= 0 && probe.breakdown.verticalReliable && + qAbs(probe.breakdown.verticalCommonBias) > 0.01 && m_parameterLower[0] > 0.0 && + m_parameterUpper[0] > m_parameterLower[0] && canEvaluate()) { + TrustRegionEvaluation corrected; + corrected.parameters = probe.parameters; + const double logRange = qLn(m_parameterUpper[0]) - qLn(m_parameterLower[0]); + const double change = qBound(-0.30 * logRange, probe.breakdown.verticalCommonBias, 0.30 * logRange); + const double k = probe.parameters[permeabilityColumn]; + corrected.parameters[permeabilityColumn] = qBound(m_parameterLower[0], k * qExp(change), m_parameterUpper[0]); + corrected.coordinates = coordinatesFromParameters(corrected.parameters); + if(qAbs(qLn(corrected.parameters[permeabilityColumn] / k)) > 1.0e-5) { + ++shapeExplorationEvaluations; + corrected.valid = evaluateTrustRegionPoint(corrected.parameters, &corrected.fitness, + &corrected.breakdown, &corrected.curve, &corrected.elapsedMs); + writeTraceRow(m_currentIteration, permeabilityColumn, "shape_length_height_probe", corrected.parameters, + corrected.fitness, corrected.valid, corrected.elapsedMs, corrected.valid ? "probe_valid" : "solver_invalid", + corrected.valid ? &corrected.breakdown : nullptr); + if(corrected.valid) keepBetter(corrected); + } + } + } + restoreEvaluationState(base); + } + restoreEvaluationState(base); + // 停止时不把未完成的探索记为完成；正常情况下至多六次调用即可覆盖本轮候选。 + if(m_shouldStop) return false; + shapeExplorationDone = true; + if(!best.valid) return false; + acceptPoint(best); + // 不用跨越不同半长基点的割线修补旧模型；在最终选中的新点重新测完整 J。 + jacobian.clear(); + reuseSensitivityForNextStage(); + rebuildRequested = true; + rebuildReason = QT_TR_NOOP("fresh sensitivity after half-length exploration"); + effectiveImprovementBaseline = current.breakdown.shapeLoss; + writeTraceRow(m_currentIteration, fractureLengthColumn, "shape_exploration_accept", current.parameters, + current.fitness, true, 0, "length_permeability_combination_accepted", ¤t.breakdown); + emit logMessageGenerated(tr("Half-length exploration accepted: L=%1, shape error=%2, total error=%3.") + .arg(current.parameters[fractureLengthColumn], 0, 'g', 6) + .arg(current.breakdown.shapeLoss, 0, 'g', 6).arg(current.fitness, 0, 'g', 6)); + return true; + }; + + auto enterTotalStage = [&](const QString& reason) { + const bool finishedRecheck = wellboreRecheck; + totalIterations = 0; + maximumEvaluations = m_totalEvaluations + totalEvaluationBudget; + shapeStage = false; + earlyShapeStage = false; + wellboreRecheck = false; + effectiveImprovementBaseline = current.fitness; + reuseSensitivityForNextStage(); + // 初调、回检只测井储/表皮；所有出口都在当前点重建完整 J，不能把缺列模型带入整体 LM。 + jacobian.clear(); + rebuildRequested = true; + rebuildReason = QT_TR_NOOP("full sensitivity at total-stage entry"); + emit progressUpdated(0, current.fitness); + emit logMessageGenerated(tr("Total-stage budget starts now: %1 iterations, %2 evaluations; pre-adjustment is counted separately.") + .arg(m_maxIterations).arg(totalEvaluationBudget)); emit logMessageGenerated(tr("Shape stage ended: %1").arg(reason)); emit logMessageGenerated(tr("LM stage 3: original LM fitting; accept by total error only.")); writeTraceRow(m_currentIteration, -1, "stage_switch", current.parameters, current.fitness, - true, 0, "shape_to_total", ¤t.breakdown); + true, 0, finishedRecheck ? "wellbore_recheck_to_total" : "shape_to_total", ¤t.breakdown); + }; + + auto finishShapeStage = [&](const QString& reason) { + if(m_shouldStop) return; + // 只检查一次；用绝对加相对容差区分有意义的退化和接近零时的比例放大。 + if(wellboreRecheckDone || wellboreParameterCount == 0 || !hasRemainingShapeParameters) { + enterTotalStage(reason); + return; + } + wellboreRecheckDone = true; + const double earlyTolerance = qMax(1.0e-4, 0.10 * earlyWellboreBaseline); + if(current.breakdown.earlyParallelLoss <= earlyWellboreBaseline + earlyTolerance) { + writeTraceRow(m_currentIteration, -1, "stage_switch", current.parameters, current.fitness, + true, 0, "wellbore_recheck_skipped_no_degradation", ¤t.breakdown); + enterTotalStage(reason); + return; + } + wellboreRecheck = true; + earlyShapeStage = true; + earlyShapeParameter = storageColumn >= 0 ? 2 : 1; + recheckTotalLimit = current.fitness + qMax(1.0e-4, 0.05 * current.fitness); + recheckShapeLimit = current.breakdown.shapeLoss + qMax(1.0e-4, 0.05 * current.breakdown.shapeLoss); + effectiveImprovementBaseline = current.breakdown.earlyParallelLoss; + // 其他参数已改变，回检入口重测两列，不能沿用初调时的响应方向。 + jacobian.clear(); + reuseSensitivityForNextStage(); + rebuildRequested = true; + rebuildReason = QT_TR_NOOP("new sensitivity for wellbore recheck"); + emit logMessageGenerated(tr("Wellbore recheck: early error=%1, total limit=%2, shape limit=%3.") + .arg(current.breakdown.earlyParallelLoss, 0, 'g', 6) + .arg(recheckTotalLimit, 0, 'g', 6).arg(recheckShapeLimit, 0, 'g', 6)); + writeTraceRow(m_currentIteration, -1, "stage_switch", current.parameters, current.fitness, + true, 0, "shape_to_wellbore_recheck", ¤t.breakdown); + }; + + auto finishEarlyShapeStage = [&](const QString& reason) { + if(earlyShapeParameter == 2 && skinColumn >= 0) { + earlyShapeParameter = 1; + effectiveImprovementBaseline = current.breakdown.earlyParallelLoss; + reuseSensitivityForNextStage(); + emit logMessageGenerated(tr("Wellbore storage adjustment ended: %1; now adjusting skin.").arg(reason)); + writeTraceRow(m_currentIteration, -1, "stage_switch", current.parameters, current.fitness, + true, 0, wellboreRecheck ? "recheck_storage_to_skin" : "storage_to_skin", ¤t.breakdown); + return; + } + if(wellboreRecheck) { + enterTotalStage(reason); + return; + } + earlyShapeStage = false; + earlyWellboreBaseline = current.breakdown.earlyParallelLoss; + emit logMessageGenerated(tr("Early wellbore adjustment ended: %1").arg(reason)); + if(!hasRemainingShapeParameters) { + enterTotalStage(tr("no remaining shape parameters")); + return; + } + effectiveImprovementBaseline = current.breakdown.shapeLoss; + // 井储、表皮可使总误差上升；联合形状阶段从当前工作点建立固定上限。 + shapeValueLimit = current.fitness + qMax(0.02, 0.25 * current.fitness); + reuseSensitivityForNextStage(); + emit logMessageGenerated(tr("Shape fitting continues with storage and skin fixed; a conditional wellbore recheck follows.")); + writeTraceRow(m_currentIteration, -1, "stage_switch", current.parameters, current.fitness, + true, 0, "early_wellbore_to_shape", ¤t.breakdown); + }; + + auto completeShapeSearchRound = [&]() { + if(!shapeStage || earlyShapeStage || !globalFallbackAttempted) return; + // 即使所有局部方向都不可行，也先完成半长组合搜索。大步接受后需在新点重走窗口。 + if(!shapeExplorationDone && exploreShapeLengths()) return; + if(m_shouldStop || !shapeExplorationDone) return; + const double required = qMax(1.0e-4, 0.01 * shapeRoundBaseline); + const double improvement = shapeRoundBaseline - current.breakdown.shapeLoss; + const bool improved = improvement >= required; + ++shapeRoundNumber; + if(improved) { + ineffectiveShapeRounds = 0; + shapeConfirmationRequested = false; + } else ++ineffectiveShapeRounds; + writeTraceRow(m_currentIteration, -1, "shape_round_end", current.parameters, current.fitness, + true, 0, QString(improved ? "improved_round_%1" : "ineffective_round_%1").arg(shapeRoundNumber), ¤t.breakdown); + emit logMessageGenerated(tr("Shape search round %1: shape error=%2, improvement=%3, required=%4.") + .arg(shapeRoundNumber).arg(current.breakdown.shapeLoss, 0, 'g', 6) + .arg(improvement, 0, 'g', 6).arg(required, 0, 'g', 6)); + if(!improved && shapeConfirmationRequested) { + finishShapeStage(tr("no significant improvement in a complete round after fresh sensitivity confirmation")); + return; + } + // 下一轮重新覆盖窗口、全局和半长；单步的小改善不会把当前一轮无限延长。 + shapeRoundBaseline = current.breakdown.shapeLoss; + attemptedWindows.fill(false); + globalFallbackAttempted = false; + shapeCandidateCount = 0; + shapeExplorationEvaluations = 0; + shapeExplorationDone = fractureLengthColumn < 0; + if(!improved && ineffectiveShapeRounds >= 2) { + jacobian.clear(); + reuseSensitivityForNextStage(); + rebuildRequested = true; + rebuildReason = QT_TR_NOOP("confirm stagnation after 2 complete shape rounds"); + shapeConfirmationRequested = true; + } }; auto registerEffectiveImprovement = [&](double fitness) -> bool { + if(shapeStage && !earlyShapeStage) return false; const double requiredImprovement = qMax( shapeStage ? 1.0e-4 : effectiveAbsoluteImprovement, qAbs(effectiveImprovementBaseline) * @@ -2136,9 +2520,15 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() // 无效尝试仍累计，但所有窗口及一次全局回退尝试完毕前，不触发停滞收敛。 // 完整一轮仍无改善时沿用原有重建确认，避免某个难处理窗口提前终止拟合。 auto recordIneffectiveStep = [&]() -> bool { + if(shapeStage && !earlyShapeStage) { + completeShapeSearchRound(); + return false; + } ++consecutiveIneffectiveSteps; - if(shapeStage) { - if(consecutiveIneffectiveSteps >= maximumIneffectiveSteps) enterTotalStage(tr("3 consecutive steps without effective shape improvement")); + if(earlyShapeStage) { + // 一次拒绝只缩步，再给当前参数一次真实尝试；连续无改善才交给下一段。 + if(consecutiveIneffectiveSteps >= 2) + finishEarlyShapeStage(tr("2 consecutive steps without effective early improvement")); return false; } if(!globalFallbackAttempted) { @@ -2172,14 +2562,16 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() .arg(effectiveRelativeImprovement * 100.0, 0, 'f', 2) .arg(maximumIneffectiveSteps)); - if(current.fitness < m_targetError) { + if(!shapeStage && current.fitness < m_targetError) { promoteSampling(true); return samplingRefreshFailed ? LM_OPTIMIZATION_FAILED : LM_TARGET_ACHIEVED; } - // 在同一个真实工作点逐参数做单边差分。首选可用空间更大的方向；只有该方向 - // 求解失败时才补算反方向，因此初次建模通常每个参数只增加一次真实求解。 + // 在同一个真实工作点逐参数做单边差分。失败先在原方向缩步，再反向尝试； + // 正常情况下每列仍只需一次真实求解，重试也计入总评价预算。 auto rebuildSensitivity = [&]() -> bool { + // 第二阶段差分不受第三阶段额度限制；每列仍只有有限的缩步和反向重试。 + const int evaluationLimit = shapeStage ? (std::numeric_limits::max)() : maximumEvaluations; const TrustRegionEvaluation base = current; const QVector baseResidual = trustRegionFullResidual(base.breakdown); const int residualCount = baseResidual.size(); @@ -2202,105 +2594,118 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() for(int column = 0; column < dimensions && - m_totalEvaluations < maximumEvaluations && + m_totalEvaluations < evaluationLimit && processPauseAndStop(); ++column) { - // 单边差分优先选择离边界空间更大的方向；首方向求解无效时才反向 - // 补算，因此正常情况下每个参数只消耗一次真实求解。 + const int parameterIndex = m_enabledParamIndices[column]; + if(earlyShapeStage && parameterIndex != 1 && parameterIndex != 2) continue; + const double lower = m_parameterLower[parameterIndex]; + const double upper = m_parameterUpper[parameterIndex]; + if(upper <= lower) continue; + // 表皮在所有阶段均按当前物理尺度扰动，避免整体阶段退回范围的 4%。 + double localStep = finiteDifferenceStep; + if(parameterIndex == 1) { + localStep = qMin(localStep, 0.02 * qMax(0.1, qAbs(base.parameters[column])) / (upper - lower)); + } else if(shapeStage && useTrustRegionLogScale(parameterIndex, lower, upper)) { + localStep = qMin(localStep, qLn(1.05) / (qLn(upper) - qLn(lower))); + } double positiveRoom = 1.0 - base.coordinates[column]; double negativeRoom = base.coordinates[column]; double preferredSign = positiveRoom >= negativeRoom ? 1.0 : -1.0; bool columnBuilt = false; + double probeStep = localStep; for(int directionAttempt = 0; directionAttempt < 2 && !columnBuilt && - m_totalEvaluations < maximumEvaluations; + m_totalEvaluations < evaluationLimit && processPauseAndStop(); ++directionAttempt) { double direction = directionAttempt == 0 ? preferredSign : -preferredSign; double availableRoom = direction > 0.0 ? positiveRoom : negativeRoom; - const int parameterIndex = m_enabledParamIndices[column]; - const double lower = m_parameterLower[parameterIndex]; - const double upper = m_parameterUpper[parameterIndex]; - double localStep = finiteDifferenceStep; - if(shapeStage && upper > lower && parameterIndex == 1) { - localStep = qMin(localStep, 0.02 * qMax(0.1, qAbs(base.parameters[column])) / (upper - lower)); - } else if(shapeStage && useTrustRegionLogScale(parameterIndex, lower, upper)) { - localStep = qMin(localStep, qLn(1.05) / (qLn(upper) - qLn(lower))); - } - double deltaMagnitude = qMin(localStep, availableRoom); - if(deltaMagnitude < minimumCoordinateStep) { - continue; - } - - TrustRegionEvaluation probe; - probe.coordinates = base.coordinates; - probe.coordinates[column] += direction * deltaMagnitude; - probe.parameters = parametersFromCoordinates(probe.coordinates); - probe.valid = evaluateTrustRegionPoint( - probe.parameters, - &probe.fitness, - &probe.breakdown, - &probe.curve, - &probe.elapsedMs); - - QString decision = probe.valid - ? "sensitivity_valid" - : (directionAttempt == 0 - ? "sensitivity_retry_opposite" - : "sensitivity_invalid"); - writeTraceRow(m_currentIteration, - column, - "trust_region_sensitivity", - probe.parameters, - probe.fitness, - probe.valid, - probe.elapsedMs, - (shapeStage ? "shape_" : "total_") + decision, - probe.valid ? &probe.breakdown : nullptr); - - if(!probe.valid) { - restoreEvaluationState(base); - continue; - } + double deltaMagnitude = qMin(probeStep, availableRoom); + for(int shrinkAttempt = 0; shrinkAttempt < 3 && !columnBuilt && + m_totalEvaluations < evaluationLimit && processPauseAndStop(); ++shrinkAttempt) { + if(deltaMagnitude < minimumCoordinateStep) { + break; + } + const bool canShrink = shrinkAttempt < 2 && deltaMagnitude * 0.5 >= minimumCoordinateStep; + + TrustRegionEvaluation probe; + probe.coordinates = base.coordinates; + probe.coordinates[column] += direction * deltaMagnitude; + probe.parameters = parametersFromCoordinates(probe.coordinates); + probe.valid = evaluateTrustRegionPoint( + probe.parameters, + &probe.fitness, + &probe.breakdown, + &probe.curve, + &probe.elapsedMs); + + QString decision = probe.valid + ? "sensitivity_valid" + : (canShrink ? "sensitivity_retry_smaller" : (directionAttempt == 0 + ? "sensitivity_retry_opposite" + : "sensitivity_invalid")); + writeTraceRow(m_currentIteration, + column, + "trust_region_sensitivity", + probe.parameters, + probe.fitness, + probe.valid, + probe.elapsedMs, + (wellboreRecheck ? "wellbore_recheck_" : (shapeStage ? "shape_" : "total_")) + decision, + probe.valid ? &probe.breakdown : nullptr); + + if(!probe.valid) { + restoreEvaluationState(base); + if(!canShrink) break; + deltaMagnitude *= 0.5; + probeStep = deltaMagnitude; + continue; + } - double delta = probe.coordinates[column] - - base.coordinates[column]; - if(qAbs(delta) < minimumCoordinateStep || - trustRegionFullResidual(probe.breakdown).size() != residualCount) { - restoreEvaluationState(base); - continue; - } + double delta = probe.coordinates[column] - + base.coordinates[column]; + if(qAbs(delta) < minimumCoordinateStep || + trustRegionFullResidual(probe.breakdown).size() != residualCount) { + restoreEvaluationState(base); + break; + } - // 第 column 列是固定残差向量相对内部参数坐标的有限差分： - // J[:,column] = (r_probe-r_base)/delta。 - const QVector probeResidual = trustRegionFullResidual(probe.breakdown); - for(int row = 0; row < residualCount; ++row) { - jacobian[row][column] = (probeResidual[row] - baseResidual[row]) / delta; - } + // 第 column 列是固定残差向量相对内部参数坐标的有限差分： + // J[:,column] = (r_probe-r_base)/delta。 + const QVector probeResidual = trustRegionFullResidual(probe.breakdown); + for(int row = 0; row < residualCount; ++row) { + jacobian[row][column] = (probeResidual[row] - baseResidual[row]) / delta; + } - jacobianColumnValid[column] = true; - columnBuilt = true; + jacobianColumnValid[column] = true; + columnBuilt = true; - if(acceptable(probe, base) && - (!bestProbe.valid || stageError(probe) < stageError(bestProbe))) { - bestProbe = probe; - bestProbeColumn = column; - bestProbeDelta = delta; + // 其他参数仍建立灵敏度供后续复用，但不能绕过当前子阶段的选参限制。 + if(parameterAllowedInStage(column) && acceptable(probe, base) && + (!bestProbe.valid || stageError(probe) < stageError(bestProbe))) { + bestProbe = probe; + bestProbeColumn = column; + bestProbeDelta = delta; + } + restoreEvaluationState(base); } - restoreEvaluationState(base); } } int validColumnCount = 0; + bool hasActiveSensitivity = false; for(int i = 0; i < jacobianColumnValid.size(); ++i) { if(jacobianColumnValid[i]) { ++validColumnCount; + if(parameterAllowedInStage(i)) hasActiveSensitivity = true; } } - if(validColumnCount == 0 || m_shouldStop) { + // 冻结列成功不能掩盖所有可调形状列的求解失败；有效零梯度仍由停滞逻辑处理。 + if(validColumnCount == 0 || (shapeStage && !earlyShapeStage && !hasActiveSensitivity) || m_shouldStop) { restoreEvaluationState(base); return false; } @@ -2323,11 +2728,14 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() current.fitness, true, 0, - shapeStage ? "shape_accepted_cached_probe" : "total_accepted_cached_probe", + wellboreRecheck ? "wellbore_recheck_accepted_cached_probe" : + (earlyShapeStage ? "early_shape_accepted_cached_probe" : + (shapeStage ? "shape_accepted_cached_probe" : "total_accepted_cached_probe")), ¤t.breakdown); emit logMessageGenerated( - tr("Sensitivity probe accepted: total error=%1") - .arg(current.fitness, 0, 'e', 4)); + (earlyShapeStage ? tr("Sensitivity probe accepted: early relative-slope matching error=%1") + : tr("Sensitivity probe accepted: total error=%1")) + .arg(earlyShapeStage ? current.breakdown.earlyParallelLoss : current.fitness, 0, 'e', 4)); } else { restoreEvaluationState(current); } @@ -2350,9 +2758,8 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() int completedIterations = 0; for(int iteration = 0; - iteration < m_maxIterations && - m_totalEvaluations < maximumEvaluations && - !m_shouldStop; + !m_shouldStop && (shapeStage || + (totalIterations < m_maxIterations && m_totalEvaluations < maximumEvaluations)); ++iteration) { m_currentIteration = iteration; completedIterations = iteration + 1; @@ -2360,19 +2767,31 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() if(!processPauseAndStop()) { break; } - if(shapeStage && (iteration >= shapeIterationLimit || - m_totalEvaluations + (rebuildRequested ? dimensions : 0) >= shapeEvaluationDeadline)) { - enterTotalStage(tr("reserve remaining iterations and evaluations for total fitting")); + if(wellboreRecheck && current.breakdown.earlyParallelLoss <= + earlyWellboreBaseline + qMax(1.0e-4, 0.10 * earlyWellboreBaseline)) + enterTotalStage(tr("early relative-slope error restored within tolerance")); + if(shapeStage && !earlyShapeStage && !jointShapeStarted) { + jointShapeStarted = true; + shapeRoundBaseline = current.breakdown.shapeLoss; + // 初调已改变井储/表皮，普通形状必须在新的工作点测完整响应。 + jacobian.clear(); + reuseSensitivityForNextStage(); + rebuildRequested = true; + rebuildReason = QT_TR_NOOP("fresh sensitivity at joint shape entry"); } + if(shapeStage && !earlyShapeStage && shapeCandidateCount >= 3) + exploreShapeLengths(); + if(m_shouldStop) break; + if(!shapeStage && (totalIterations >= m_maxIterations || m_totalEvaluations >= maximumEvaluations)) break; if(!shapeStage && m_layeredSampling && m_samplingStride > 1) { const int remainingEvaluations = maximumEvaluations - m_totalEvaluations; // 已准备重建时，提前计入本轮差分及候选的开销，避免紧接着因预算再加密。 const bool reserveFinalBudget = remainingEvaluations <= 2 * (dimensions + 1) || (rebuildRequested && remainingEvaluations <= 3 * (dimensions + 1)); const int layerDeadline = qMax(1, m_maxIterations * (m_samplingStride == 4 ? 1 : 2) / 3); - const bool mergeRefinement = rebuildRequested && iteration < layerDeadline && - (iteration + 1 >= layerDeadline || reserveFinalBudget); - if(reserveFinalBudget || iteration >= layerDeadline || mergeRefinement) { + const bool mergeRefinement = rebuildRequested && totalIterations < layerDeadline && + (totalIterations + 1 >= layerDeadline || reserveFinalBudget); + if(reserveFinalBudget || totalIterations >= layerDeadline || mergeRefinement) { const bool promoted = promoteSampling(reserveFinalBudget); if(samplingRefreshFailed) break; if(promoted && mergeRefinement) { @@ -2384,14 +2803,19 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() emit logMessageGenerated(tr("Rebuilding sensitivity model: %1").arg(tr(rebuildReason))); writeTraceRow(m_currentIteration, -1, "sensitivity_rebuild", current.parameters, current.fitness, true, 0, rebuildReason, ¤t.breakdown); - const bool confirmingStagnation = - stagnationConfirmationRequested; if(!rebuildSensitivity()) { if(shapeStage && !m_shouldStop) { - enterTotalStage(tr("no valid shape sensitivity model")); - rebuildRequested = true; - rebuildReason = QT_TR_NOOP("no valid sensitivity model"); - continue; + if(earlyShapeStage) { + finishEarlyShapeStage(tr("no valid early sensitivity model")); + rebuildRequested = true; + rebuildReason = QT_TR_NOOP("no valid sensitivity model"); + continue; + } + // 真实差分全失败不能冒充停滞收敛；独立探索若找到有效新点则重试。 + if(exploreShapeLengths()) continue; + m_lastError = tr("Unable to build a valid shape sensitivity model."); + stopReason = LM_OPTIMIZATION_FAILED; + break; } if(promoteSampling(false)) continue; stopReason = m_shouldStop @@ -2399,31 +2823,40 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() : LM_LOCAL_OPTIMUM; break; } - if(current.fitness < m_targetError) { + if(shapeStage && !earlyShapeStage && current.fitness < m_targetError) { + finishShapeStage(tr("total target reached during shape fitting")); + --iteration; + continue; + } + if(!shapeStage && current.fitness < m_targetError) { promoteSampling(true); stopReason = LM_TARGET_ACHIEVED; break; } - if(m_totalEvaluations >= maximumEvaluations) { + if(!shapeStage && m_totalEvaluations >= maximumEvaluations) { stopReason = LM_MAX_ITERATIONS; break; } - const bool rebuildEffective = - registerEffectiveImprovement(stageError(current)); - if(confirmingStagnation && !rebuildEffective) { - if(promoteSampling(false)) continue; - emit logMessageGenerated( - tr("Sensitivity rebuild produced no effective improvement; " - "local convergence detected")); - stopReason = LM_LOCAL_OPTIMUM; - break; + if(wellboreRecheck && current.breakdown.earlyParallelLoss <= + earlyWellboreBaseline + qMax(1.0e-4, 0.10 * earlyWellboreBaseline)) { + enterTotalStage(tr("wellbore recheck completed during sensitivity evaluation")); + --iteration; + continue; } + // 单列探针没有改善不代表联合步无效；重建后继续真实候选评价，再确认停滞。 + registerEffectiveImprovement(stageError(current)); } + // 第三阶段独立计数；方向不可行也消耗一次局部尝试，不能无限缩步循环。 + if(!shapeStage) ++totalIterations; // 每轮从最新 J 和当前残差重算窗口 Fisher，包含有限差分与割线更新的变化。 - const QVector objectiveResidual = shapeStage - ? current.breakdown.shapeResiduals : current.breakdown.residualVector; - const int rowOffset = shapeStage ? current.breakdown.residualVector.size() : 0; + const QVector objectiveResidual = earlyShapeStage + ? current.breakdown.earlyParallelResiduals + : (shapeStage ? current.breakdown.shapeResiduals : current.breakdown.residualVector); + const int rowOffset = earlyShapeStage + ? current.breakdown.residualVector.size() + current.breakdown.shapeResiduals.size() + + current.breakdown.earlyValueResiduals.size() + : (shapeStage ? current.breakdown.residualVector.size() : 0); const QVector > objectiveJacobian = jacobian.mid(rowOffset, objectiveResidual.size()); QVector objectiveCoordinates = current.breakdown.sampleCoordinates; if(shapeStage) { @@ -2437,15 +2870,37 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() const QVector information = buildTrustRegionFisher( objectiveJacobian, objectiveResidual, jacobianColumnValid, objectiveCoordinates); const TrustRegionFisher& global = information[kAutoFitTimeWindowCount]; + // 前期残差已含第一窗口权重，使用全局和，不能再次乘局部窗口权重。 + const TrustRegionFisher& stepInformation = global; + QVector stageColumnValid = jacobianColumnValid; + for(int column = 0; column < dimensions; ++column) { + if(!parameterAllowedInStage(column)) stageColumnValid[column] = false; + } QVector selectedColumns; QVector coordinateStep; double predictedReduction = 0.0; int selectedWindow = -1; - // 最大能量窗口优先。某窗口无可行全局下降步时立即换下一个，不调用 DLL。 - // 全部窗口都处理不动后，再用全局 Fisher 作一次补充选参。 - while(selectedColumns.isEmpty()) { - selectedWindow = nextTrustRegionWindow( + if(!shapeStage) { + // 整体阶段联合求解全部有效自由列，不按窗口得分、三参数上限或共线阈值删列。 + // 边界列也参与耦合求解，是否向内移动由联合解及其边界投影决定。 + for(int column = 0; column < dimensions; ++column) { + // 仅跳过完全没有总残差响应的列，保留弱敏感列及相互相关的列。 + if(stageColumnValid[column] && global.matrix[column][column] > 0.0) + selectedColumns.append(column); + } + // 当前候选已使用全局信息，停滞判断无需再逐一尝试窗口子集。 + attemptedWindows.fill(true); + globalFallbackAttempted = true; + if(selectedColumns.isEmpty() || !buildTrustRegionFisherStep( + global, selectedColumns, current.coordinates, damping, + trustRadius, minimumCoordinateStep, &coordinateStep, &predictedReduction)) { + selectedColumns.clear(); + } + } + // 预调整阶段保留原窗口筛选及子集比较，整体阶段直接采用上面的联合步。 + while(shapeStage && selectedColumns.isEmpty()) { + selectedWindow = earlyShapeStage ? 0 : nextTrustRegionWindow( objectiveWindows, attemptedWindows); if(selectedWindow >= 0) { attemptedWindows[selectedWindow] = true; @@ -2454,9 +2909,17 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() } const TrustRegionFisher& local = selectedWindow >= 0 ? information[selectedWindow] : global; - const QVector proposed = selectTrustRegionFisherColumns( - local, global, jacobianColumnValid, - current.coordinates, minimumCoordinateStep); + QVector proposed; + if(earlyShapeStage) { + // 当前只开放井储或表皮；其余阶段仍使用原 Fisher 组合筛选。 + for(int column = 0; column < dimensions; ++column) { + if(stageColumnValid[column]) proposed.append(column); + } + } else { + proposed = selectTrustRegionFisherColumns( + local, global, stageColumnValid, + current.coordinates, minimumCoordinateStep); + } // 最多三个推荐参数，比较其全部非空子集（最多七组），避免首参数必选。 // 这里只做小矩阵运算，真正的候选评价每轮仍至多一次。 for(int mask = 1; mask < (1 << proposed.size()); ++mask) { @@ -2468,15 +2931,12 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() } QVector step; double reduction = 0.0; - if(!buildTrustRegionFisherStep( - global, columns, current.coordinates, - damping, trustRadius, minimumCoordinateStep, - &step, &reduction)) { - continue; - } - // 仅形状预调整限制整体偏离；第三阶段直接使用原 LM 步长， - // 不预测形状上限，也不因形状变化缩短候选步长。 - if(shapeStage) { + // 井储、表皮均由相对斜率残差的实际灵敏度确定方向，不预设参数增减。 + const bool stepBuilt = buildTrustRegionFisherStep(stepInformation, columns, current.coordinates, + damping, trustRadius, minimumCoordinateStep, &step, &reduction); + if(!stepBuilt) continue; + // 仅联合形状阶段预测总误差上限；井储、表皮步长不受总误差裁剪。 + if(shapeStage && !earlyShapeStage) { auto predictedGuard = [&](double scale) -> double { double energy = 0.0; for(int row = 0; row < current.breakdown.residualVector.size(); ++row) { @@ -2499,10 +2959,10 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() const double stepScale = low; for(int column = 0; column < dimensions; ++column) step[column] *= stepScale; if(qSqrt(trustRegionSquaredNorm(step)) < minimumCoordinateStep) continue; - reduction = -trustRegionDotProduct(global.gradient, step); + reduction = -trustRegionDotProduct(stepInformation.gradient, step); for(int a = 0; a < dimensions; ++a) for(int b = 0; b < dimensions; ++b) - reduction -= 0.5 * step[a] * global.matrix[a][b] * step[b]; + reduction -= 0.5 * step[a] * stepInformation.matrix[a][b] * step[b]; if(!isFiniteNumber(reduction) || reduction <= 1.0e-14) continue; } } @@ -2515,14 +2975,24 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() predictedReduction = reduction; } } - if(selectedWindow < 0) { + if(earlyShapeStage || selectedWindow < 0) { break; } } - // 只有窗口候选与全局回退均无方向，才收缩半径并进入原有重建/收敛处理。 + // 当前阶段没有可行下降步时，收缩半径并进入原有重建/收敛处理。 if(selectedColumns.isEmpty()) { - if(shapeStage) { enterTotalStage(tr("no feasible shape descent step")); continue; } + if(earlyShapeStage) { + finishEarlyShapeStage(tr("no feasible early adjustment direction or parameter at bound")); + // 尚未求解候选，原迭代留给其余参数，不额外占用迭代或求解预算。 + --iteration; + continue; + } + if(shapeStage) { + completeShapeSearchRound(); + --iteration; + continue; + } if(trustRadius <= minimumTrustRadius * 1.01 && modelRebuiltAtMinimumRadius) { if(promoteSampling(false)) continue; @@ -2541,7 +3011,7 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() continue; } - const double stepNorm = qSqrt(trustRegionSquaredNorm(coordinateStep)); + double stepNorm = qSqrt(trustRegionSquaredNorm(coordinateStep)); QVector candidateCoordinates = current.coordinates; for(int column = 0; column < dimensions; ++column) { candidateCoordinates[column] += coordinateStep[column]; @@ -2551,19 +3021,27 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() for(int i = 0; i < selectedColumns.size(); ++i) { selectedParameterIndices << QString::number(m_enabledParamIndices[selectedColumns[i]]); } - const QString selectionName = QString(shapeStage ? "shape_" : "total_") + (selectedWindow >= 0 + QString selectionName = (wellboreRecheck ? QString("recheck_") : QString()) + QString(earlyShapeStage + ? (earlyShapeParameter == 2 ? "storage_parallel_" : "skin_parallel_") : (shapeStage ? "shape_" : "total_")) + (selectedWindow >= 0 ? QString("window_%1").arg(selectedWindow + 1) : QString("global")) + "_params_" + selectedParameterIndices.join("_"); TrustRegionEvaluation candidate; candidate.coordinates = candidateCoordinates; candidate.parameters = parametersFromCoordinates(candidate.coordinates); + if(earlyShapeStage && earlyShapeParameter == 2) { + emit logMessageGenerated(tr("Wellbore storage trial: relative-slope matching error=%1, C=%2 -> %3") + .arg(current.breakdown.earlyParallelLoss, 0, 'g', 5) + .arg(current.parameters[storageColumn], 0, 'g', 6) + .arg(candidate.parameters[storageColumn], 0, 'g', 6)); + } candidate.valid = evaluateTrustRegionPoint( candidate.parameters, &candidate.fitness, &candidate.breakdown, &candidate.curve, &candidate.elapsedMs); + if(shapeStage && !earlyShapeStage) ++shapeCandidateCount; if(!candidate.valid) { // 求解失败的候选不能改变 current。先完整恢复上一个已接受参数和 @@ -2599,6 +3077,49 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() continue; } + // 普通联合步真实表现可靠时，给含导流能力的方向一次较大步试算。 + // 先保存原步；扩展失败或不如原步时仍保留原候选，本轮额外调用共用八次上限。 + if(shapeStage && !earlyShapeStage && shapeExplorationDone && conductivityColumn >= 0 && + selectedColumns.contains(conductivityColumn) && acceptable(candidate, current) && + shapeExplorationEvaluations < maximumShapeExplorationEvaluations && processPauseAndStop()) { + const double ordinaryEnergy = 0.5 * trustRegionSquaredNorm(objectiveResidual); + const double ordinaryReduction = ordinaryEnergy - 0.5 * trustRegionSquaredNorm(candidate.breakdown.shapeResiduals); + const double predictionFloor = qMax(1.0e-14, ordinaryEnergy * 1.0e-8); + QVector expandedStep; + double expandedPrediction = 0.0; + if(predictedReduction > predictionFloor && ordinaryReduction / predictedReduction > 0.75 && + buildExpandedTrustRegionStep(stepInformation, current.coordinates, coordinateStep, + maximumTrustRadius, &expandedStep, &expandedPrediction)) { + writeTraceRow(m_currentIteration, -1, "shape_step_base", candidate.parameters, + candidate.fitness, true, candidate.elapsedMs, "eligible_for_expansion", &candidate.breakdown); + TrustRegionEvaluation expanded; + expanded.coordinates = current.coordinates; + for(int i = 0; i < dimensions; ++i) expanded.coordinates[i] += expandedStep[i]; + expanded.parameters = parametersFromCoordinates(expanded.coordinates); + ++shapeExplorationEvaluations; + expanded.valid = evaluateTrustRegionPoint(expanded.parameters, &expanded.fitness, + &expanded.breakdown, &expanded.curve, &expanded.elapsedMs); + const bool useExpanded = acceptable(expanded, current) && expanded.breakdown.shapeLoss < candidate.breakdown.shapeLoss; + writeTraceRow(m_currentIteration, -1, "shape_step_expanded", expanded.parameters, + expanded.fitness, expanded.valid, expanded.elapsedMs, useExpanded ? "selected_expanded" : "retain_ordinary", + expanded.valid ? &expanded.breakdown : nullptr); + if(useExpanded) { + candidate = expanded; + coordinateStep = expandedStep; + predictedReduction = expandedPrediction; + stepNorm = qSqrt(trustRegionSquaredNorm(coordinateStep)); + selectionName += "_expanded"; + } + // 两次真实调用已分别记录，下方候选行只记录最终决策，避免重复统计耗时。 + candidate.elapsedMs = -1; + restoreEvaluationState(current); + } + } + if(m_shouldStop) { + restoreEvaluationState(current); + stopReason = LM_USER_STOPPED; + break; + } consecutiveSolverFailures = 0; // 有效候选即使最终被拒绝，也提供了一条真实割线，可用于修正下一轮 // 局部模型；是否成为新工作点由当前阶段的目标和约束共同决定。 @@ -2611,9 +3132,12 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() // reductionRatio 衡量局部线性模型的可信度：接近 1 表示预测准确； // 值较小表示虽然可能下降，但模型低估了非线性，需要收紧下一步。 const double objectiveEnergy = 0.5 * trustRegionSquaredNorm(objectiveResidual); - double actualReduction = objectiveEnergy - 0.5 * ( - trustRegionSquaredNorm(shapeStage ? candidate.breakdown.shapeResiduals : candidate.breakdown.residualVector)); - double reductionRatio = actualReduction / predictedReduction; + const QVector candidateResidual = earlyShapeStage + ? candidate.breakdown.earlyParallelResiduals + : (shapeStage ? candidate.breakdown.shapeResiduals : candidate.breakdown.residualVector); + const double candidateEnergy = 0.5 * trustRegionSquaredNorm(candidateResidual); + double actualReduction = objectiveEnergy - candidateEnergy; + double reductionRatio = predictedReduction > 1.0e-14 ? actualReduction / predictedReduction : 0.0; // 使用同一采样层、同一阶段目标比较预测和实际改善。接近收敛时的微小 // 预测量交给停滞逻辑处理，避免比例数值波动反复触发昂贵的全参数重建。 const double predictionFloor = qMax(1.0e-14, objectiveEnergy * 1.0e-8); @@ -2631,7 +3155,8 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() } else { fullDataRejections = 0; } - QString componentName = shapeStage ? "shape" : "total"; + QString componentName = earlyShapeStage + ? (earlyShapeParameter == 2 ? "storage_parallel" : "skin_parallel") : (shapeStage ? "shape" : "total"); if(accepted) { // 当前阶段接受候选后同步发布参数、曲线和诊断。 @@ -2669,14 +3194,6 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() restoreEvaluationState(current); } - // 候选只要更优就继续作为 current 保存；是否足以解除停滞，则统一 - // 相对上一次有效改善基准判断。拒绝和微小改善都会累计无效次数。 - const bool effectiveImprovement = - registerEffectiveImprovement(stageError(current)); - if(!effectiveImprovement && recordIneffectiveStep()) { - stopReason = LM_LOCAL_OPTIMUM; - } - writeTraceRow(m_currentIteration, -1, "trust_region_candidate", @@ -2694,6 +3211,10 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() componentDisplayName = tr("vertical deviation"); } else if(componentName == "horizontal") { componentDisplayName = tr("horizontal deviation"); + } else if(componentName == "storage_parallel") { + componentDisplayName = tr("wellbore storage relative-slope matching"); + } else if(componentName == "skin_parallel") { + componentDisplayName = tr("skin relative-slope matching"); } else if(componentName == "shape") { componentDisplayName = tr("shape deviation"); } else if(componentName == "total") { @@ -2705,9 +3226,13 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() .arg(iteration + 1) .arg(componentDisplayName) .arg(selectedColumns.size()) - .arg(candidate.fitness, 0, 'e', 4) + .arg(earlyShapeStage ? candidate.breakdown.earlyParallelLoss : candidate.fitness, 0, 'e', 4) .arg(accepted ? tr("accepted") : tr("rejected"))); - emit progressUpdated(iteration + 1, m_globalBestFitness); + emit progressUpdated(shapeStage ? -1 : totalIterations, m_globalBestFitness); + + // 先记录当前试调结果，再发布阶段切换，避免日志显示为未试调就结束。 + const bool effectiveImprovement = registerEffectiveImprovement(stageError(current)); + if(!effectiveImprovement && recordIneffectiveStep()) stopReason = LM_LOCAL_OPTIMUM; if(stopReason == LM_LOCAL_OPTIMUM || fullDataRejections >= 2) { if(promoteSampling(false)) { @@ -2717,12 +3242,15 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() if(stopReason == LM_LOCAL_OPTIMUM || samplingRefreshFailed) break; } - if(current.fitness < m_targetError) { + if(shapeStage && !earlyShapeStage && current.fitness < m_targetError) { + finishShapeStage(tr("total target reached during shape fitting")); + } + if(!shapeStage && current.fitness < m_targetError) { promoteSampling(true); stopReason = LM_TARGET_ACHIEVED; break; } - if(selectedWindow < 0 && trustRadius <= minimumTrustRadius * 1.01 && + if(!shapeStage && selectedWindow < 0 && trustRadius <= minimumTrustRadius * 1.01 && consecutiveRejectedSteps >= 2 && consecutivePoorPredictions >= 2) { if(modelRebuiltAtMinimumRadius) { if(promoteSampling(false)) continue; @@ -2738,6 +3266,8 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting() m_currentIteration = completedIterations - 1; } restoreEvaluationState(current); + emit logMessageGenerated(tr("Adaptive fitting counts: %1 completed shape rounds, %2 total-stage iterations, %3 total evaluations.") + .arg(shapeRoundNumber).arg(totalIterations).arg(m_totalEvaluations)); if(m_shouldStop) { return LM_USER_STOPPED; @@ -3581,6 +4111,7 @@ double nmCalculationAutoFitLM::calculateLogLogCurveError( const int lag = 8; const double span = qLn(overlapMaxX) - qLn(overlapMinX); QVector residuals[2]; + QVector targetLogs[2], resultLogs[2]; double biases[2] = {0.0, 0.0}; for(int component = 0; component < 2; ++component) { for(int i = 0; i < count; ++i) { @@ -3593,6 +4124,8 @@ double nmCalculationAutoFitLM::calculateLogLogCurveError( time, &resultValue)) return false; const double residual = resultValue - targetValue; residuals[component].append(residual); + targetLogs[component].append(targetValue); + resultLogs[component].append(resultValue); // 对数时间梯形权重等价于互补窗口加权求和，避免密集段主导高度。 biases[component] += residual * ((i == 0 || i == count - 1) ? 0.5 : 1.0) / (count - 1); @@ -3614,6 +4147,7 @@ double nmCalculationAutoFitLM::calculateLogLogCurveError( } } objective->shapeLoss = qSqrt(trustRegionSquaredNorm(objective->shapeResiduals)); + populateEarlyWellboreMetrics(objective, targetLogs, resultLogs, span); return isFiniteNumber(objective->shapeLoss); }; @@ -4625,7 +5159,7 @@ QString nmCalculationAutoFitLM::getStopReasonDescription(StopReasonLM reason) co return tr("Local optimum detected"); case LM_MAX_ITERATIONS: - return tr("Maximum iterations reached"); + return tr("Total-stage iteration or evaluation budget reached"); case LM_USER_STOPPED: return tr("Stopped by user request"); diff --git a/Src/nmNum/nmSubWxs/nmWxAutomaticFittingStart.cpp b/Src/nmNum/nmSubWxs/nmWxAutomaticFittingStart.cpp index 4a508097..9de73cd7 100644 --- a/Src/nmNum/nmSubWxs/nmWxAutomaticFittingStart.cpp +++ b/Src/nmNum/nmSubWxs/nmWxAutomaticFittingStart.cpp @@ -689,16 +689,19 @@ void nmWxAutomaticfittingStart::setPseudoPressureMode(bool enabled) void nmWxAutomaticfittingStart::onFittingProgress(int iteration, double fitness) { - // 确保iteration在合理范围内 - int displayIteration = qMax(1, qMin(iteration, m_maxIterations)); - - // 更新进度条 - progressBar->setValue(displayIteration); - double progress = (double)displayIteration / m_maxIterations * 100; - progressBar->setFormat(QString("%1%").arg(progress, 0, 'f', 1)); + // LM 预调整没有固定总步数，单独显示阶段；整体阶段从局部迭代 0 开始计进度。 + if (m_autoFitterLM && iteration < 0) { + progressBar->setValue(0); + progressBar->setFormat(tr("Pre-adjustment")); + currentIterationValue->setText(tr("Pre-adjustment")); + } else { + int displayIteration = qMax(m_autoFitterLM ? 0 : 1, qMin(iteration, m_maxIterations)); + progressBar->setValue(displayIteration); + double progress = (double)displayIteration / m_maxIterations * 100; + progressBar->setFormat(QString("%1%").arg(progress, 0, 'f', 1)); + currentIterationValue->setText(QString::number(displayIteration)); + } - // 更新参数显示 - currentIterationValue->setText(QString::number(displayIteration)); currentComfortValue->setText(formatScientific(fitness)); // 更新最佳适应度