YXS
/
nmWTAI-Platform
8
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

refactor(nmNum): 调整直井拟合流程并清理窗口与 Fisher 遗留逻辑

Browse Source 
- 普通直井跳过形状调整，直接进入总误差 LM
- 移除时间窗口及无用诊断，保留 LM 整体矩阵计算
feature/AutoFit-Optimize-20260914
lvjunjie 2 weeks ago
parent 1b5d535df6
commit b501816f6b
2 changed files with 46 additions and 308 deletions
Show Stats Download Patch File Download Diff File

28
Include/nmNum/nmCalculation/nmCalculationAutoFitLM.h
Unescape Escape View File

@ -12,29 +12,14 @@
#include "nmCalculation_global.h"
// 固定对数时间窗口的误差诊断。时间边界为不含重叠区的基础边界；
// rmsError 是局部加权均方根，energy 是对整体残差能量的贡献。
struct AutoFitTimeWindowLM {
double timeMin;
double timeMax;
double weightSum;
double rmsError;
double energy;
AutoFitTimeWindowLM()
: timeMin(0.0), timeMax(0.0), weightSum(0.0), rmsError(0.0), energy(0.0)
{}
};
// 双对数曲线误差分解。预调整使用纵向偏差，形状阶段使用双曲线斜率残差，
// 整体阶段只使用 total 接受候选；时间窗口和残差用于 Fisher 选参。
// 整体阶段只使用 total 接受候选；两个 LM 目标均联合调整有效参数。
struct AutoFitObjectiveBreakdownLM {
bool valid;
double total;
double pressureLoss;
double derivativeLoss;
QVector<double> residualVector;
QVector<AutoFitTimeWindowLM> timeWindows;
double verticalCommonBias;
double verticalLoss;
bool verticalReliable;
@ -44,14 +29,6 @@ struct AutoFitObjectiveBreakdownLM {
bool registrationAmbiguous;
// 与整体误差共用 log-time 采样网格，斜率跨度约为完整对数时域的 10%。
QVector<double> shapeResiduals;
// 前期数值和形状残差沿用同一网格与斜率跨度，均含第一窗口权重。
// earlyParallelResiduals/Loss 沿用历史字段名，现为压力、导数分别匹配目标的形状误差。
// earlyParallelBias 仍记录相对斜率偏差，仅供诊断。
QVector<double> earlyValueResiduals;
QVector<double> earlyParallelResiduals;
double earlyValueLoss;
double earlyParallelLoss;
double earlyParallelBias;
double pressureVerticalBias;
double derivativeVerticalBias;
double shapeLoss;
@ -71,9 +48,6 @@ struct AutoFitObjectiveBreakdownLM {
, horizontalLoss(std::numeric_limits<double>::quiet_NaN())
, horizontalReliable(false)
, registrationAmbiguous(false)
, earlyValueLoss(std::numeric_limits<double>::quiet_NaN())
, earlyParallelLoss(std::numeric_limits<double>::quiet_NaN())
, earlyParallelBias(std::numeric_limits<double>::quiet_NaN())
, pressureVerticalBias(0.0)
, derivativeVerticalBias(0.0)
, shapeLoss(std::numeric_limits<double>::quiet_NaN())

320
Src/nmNum/nmCalculation/nmCalculationAutoFitLM.cpp
Unescape Escape View File

@ -27,8 +27,6 @@
#endif
static const bool kAutoFitDiagnosticTraceEnabled = true;
static const int kAutoFitTimeWindowCount = 4;
static const double kAutoFitTimeWindowOverlapRatio = 0.20;
static const int kAutoFitIntervalsPerDecade = 20;
static int autoFitShapeLag(int pointCount)
@ -37,57 +35,6 @@ static int autoFitShapeLag(int pointCount)
return qMax(1, qRound((pointCount - 1) * 0.10));
}
static double autoFitTimeWindowWeight(double coordinate, int windowIndex)
{
// 重叠区以基础边界为中心，总宽度为窗口宽度的 20%；相邻权重互补。
const double width = 1.0 / kAutoFitTimeWindowCount;
const double overlap = width * kAutoFitTimeWindowOverlapRatio;
const double halfOverlap = overlap * 0.5;
const double pi = 3.14159265358979323846;
double weight = 1.0;
if(windowIndex > 0) {
double u = qBound(0.0,
(coordinate - windowIndex * width + halfOverlap) / overlap,
1.0);
weight *= 0.5 * (1.0 - qCos(pi * u));
}
if(windowIndex + 1 < kAutoFitTimeWindowCount) {
double u = qBound(0.0,
(coordinate - (windowIndex + 1) * width + halfOverlap) / overlap,
1.0);
weight *= 0.5 * (1.0 + qCos(pi * u));
}
return weight;
}
static QVector<AutoFitTimeWindowLM> calculateAutoFitTimeWindows(
const QVector<double>& residualVector, double timeMin, double timeMax)
{
// 残差前后两半分别为压力和导数，已包含各占一半及时间采样权重的归一化。
const int pointCount = residualVector.size() / 2;
const double logMin = qLn(timeMin);
const double logSpan = qLn(timeMax) - logMin;
QVector<AutoFitTimeWindowLM> windows(kAutoFitTimeWindowCount);
for(int k = 0; k < windows.size(); ++k) {
AutoFitTimeWindowLM& window = windows[k];
window.timeMin = k == 0 ? timeMin
: qExp(logMin + logSpan * k / kAutoFitTimeWindowCount);
window.timeMax = k + 1 == windows.size() ? timeMax
: qExp(logMin + logSpan * (k + 1) / kAutoFitTimeWindowCount);
for(int i = 0; i < pointCount; ++i) {
const double weight = autoFitTimeWindowWeight(
static_cast<double>(i) / (pointCount - 1), k);
window.weightSum += weight;
window.energy += weight *
(residualVector[i] * residualVector[i] +
residualVector[pointCount + i] * residualVector[pointCount + i]);
}
window.rmsError = window.weightSum > 0.0
? qSqrt(window.energy * pointCount / window.weightSum) : 0.0;
}
return windows;
}
static inline bool isFiniteNumber(double value)
{
#ifdef Q_OS_WIN
@ -279,26 +226,14 @@ static bool trustRegionResidualsValid(
for(int i = 0; i < breakdown.shapeResiduals.size(); ++i) {
if(!isFiniteNumber(breakdown.shapeResiduals[i])) return false;
}
if(breakdown.earlyValueResiduals.size() != 2 * pointCount ||
breakdown.earlyParallelResiduals.size() != 2 * slopeCount ||
!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;
}
// 数值、整体形状和前期形状残差共用一次求解，联合缓存供各子阶段复用。
// 数值和整体形状残差共用一次求解，联合缓存供两个 LM 目标复用。
static QVector<double> trustRegionFullResidual(const AutoFitObjectiveBreakdownLM& objective)
{
QVector<double> residual = objective.residualVector;
residual += objective.shapeResiduals;
residual += objective.earlyValueResiduals;
residual += objective.earlyParallelResiduals;
return residual;
}
@ -312,108 +247,6 @@ static double trustRegionSquaredNorm(const QVector<double>& values)
return sum;
}
// 数值残差已含 sqrt(0.5 * 窗口权重)，相减并换底后得到 log10 间距的加权 RMS。
// 不减去平均间距，保留“两条曲线整体离得太远”的信息；共同纵移会自然抵消。
static double trustRegionEarlyGapLoss(const AutoFitObjectiveBreakdownLM& objective)
{
const int count = objective.earlyValueResiduals.size() / 2;
double energy = 0.0;
for(int i = 0; i < count; ++i) {
const double difference = objective.earlyValueResiduals[i] - objective.earlyValueResiduals[count + i];
energy += 2.0 * difference * difference;
}
return qSqrt(energy) / qLn(10.0);
}
// 前期带符号的平均间距差：正值表示模拟比目标分得更开，负值表示更靠近。
// 复用现有加权数值残差，恢复第一窗口的归一化权重，不用 RMS 推断偏差符号。
static double trustRegionEarlyGapBias(const AutoFitObjectiveBreakdownLM& objective)
{
const int count = objective.earlyValueResiduals.size() / 2;
QVector<double> weights(count, 0.0);
double weightSum = 0.0;
for(int i = 0; i < count; ++i) {
weights[i] = autoFitTimeWindowWeight(static_cast<double>(i) / (count - 1), 0) *
((i == 0 || i == count - 1) ? 0.5 : 1.0);
weightSum += weights[i];
}
double bias = 0.0;
for(int i = 0; i < count; ++i) {
bias += qSqrt(2.0 * weights[i] / weightSum) *
(objective.earlyValueResiduals[i] - objective.earlyValueResiduals[count + i]);
}
return bias / qLn(10.0);
}
// 第一窗口的形状能量与局部 Fisher 使用相同的中心时间和重叠权重，
// 压力、导数同时参与；不重新插值或调用求解器。
static double trustRegionEarlyShapeEnergy(const AutoFitObjectiveBreakdownLM& objective)
{
const QVector<double>& residual = objective.shapeResiduals;
const int pointCount = objective.residualVector.size() / 2;
const int lag = autoFitShapeLag(pointCount);
const int slopeCount = pointCount - lag;
double energy = 0.0;
for(int row = 0; row < residual.size(); ++row) {
const double coordinate = (row % slopeCount + lag * 0.5) / (pointCount - 1);
energy += autoFitTimeWindowWeight(coordinate, 0) * residual[row] * residual[row];
}
return energy;
}
// 前期复用整体形状的双曲线斜率残差，只按第一窗口重新加权并归一化。
// 单独平移任一曲线不改变此指标；前期数值误差仅作诊断，不参与井储、表皮验收。
static void populateEarlyWellboreMetrics(AutoFitObjectiveBreakdownLM* objective,
const QVector<double> targetLogs[2], const QVector<double> resultLogs[2], double logTimeSpan)
{
const int count = targetLogs[0].size();
QVector<double> weights(count, 0.0);
double weightSum = 0.0;
for(int i = 0; i < count; ++i) {
weights[i] = autoFitTimeWindowWeight(static_cast<double>(i) / (count - 1), 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));
// 复用整体形状的斜率跨度，窗口权重取区间中心，随采样数量同步变化。
const int lag = autoFitShapeLag(count);
const int slopeCount = count - lag;
const double logTimeStep = logTimeSpan * lag / (count - 1);
QVector<double> 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;
// 整体残差已含 sqrt(0.5 / slopeCount)，换成第一窗口权重后仍保持双曲线等权。
// 压力、导数分别匹配目标；相对斜率偏差只保留为诊断，不参与前期验收。
const double scale = qSqrt(slopeCount * weight);
objective->earlyParallelResiduals[i] = scale * objective->shapeResiduals[i];
objective->earlyParallelResiduals[slopeCount + i] = scale * objective->shapeResiduals[slopeCount + i];
}
objective->earlyParallelLoss = qSqrt(trustRegionSquaredNorm(objective->earlyParallelResiduals));
}
// 计算同维向量内积；维度不一致表示局部模型无效，返回零让调用方放弃修正。
static double trustRegionDotProduct(const QVector<double>& left,
const QVector<double>& right)
@ -498,65 +331,45 @@ static bool solveTrustRegionLinearSystem(
return true;
}
// Fisher 仅作为当前归一化坐标下的局部信息矩阵，不用于统计置信区间。
struct TrustRegionFisher
// LM 正规方程所需的 J^T*J 与 J^T*r，使用当前归一化参数坐标。
struct TrustRegionNormalEquations
{
QVector<QVector<double> > matrix;
QVector<double> gradient;
explicit TrustRegionFisher(int dimensions = 0)
explicit TrustRegionNormalEquations(int dimensions = 0)
: matrix(dimensions, QVector<double>(dimensions, 0.0))
, gradient(dimensions, 0.0)
{}
};
static QVector<TrustRegionFisher> buildTrustRegionFisher(
static TrustRegionNormalEquations buildTrustRegionNormalEquations(
const QVector<QVector<double> >& jacobian,
const QVector<double>& residual,
const QVector<bool>& columnValid,
const QVector<double>& coordinates = QVector<double>())
const QVector<bool>& columnValid)
{
// 残差和 J 已包含压力/导数及时间采样权重，只再乘一次窗口权重。
// 前半段是压力，后半段是导数；同一时间点的两行使用相同权重。
// 残差和 J 已包含压力、导数及采样权重，直接累加整体矩阵，避免重复计权。
const int dimensions = columnValid.size();
const int pointCount = residual.size() / 2;
QVector<TrustRegionFisher> information(
kAutoFitTimeWindowCount + 1, TrustRegionFisher(dimensions));
for(int k = 0; k < kAutoFitTimeWindowCount; ++k) {
TrustRegionFisher& local = information[k];
TrustRegionNormalEquations equations(dimensions);
for(int row = 0; row < residual.size(); ++row) {
const double weight = autoFitTimeWindowWeight(
coordinates.isEmpty() ? static_cast<double>(row % pointCount) / (pointCount - 1)
: coordinates[row % pointCount], k);
for(int p = 0; p < dimensions; ++p) {
if(!columnValid[p]) {
continue;
}
local.gradient[p] += weight * jacobian[row][p] * residual[row];
if(!columnValid[p]) continue;
equations.gradient[p] += jacobian[row][p] * residual[row];
for(int q = 0; q <= p; ++q) {
if(columnValid[q]) {
local.matrix[p][q] +=
weight * jacobian[row][p] * jacobian[row][q];
}
if(columnValid[q])
equations.matrix[p][q] += jacobian[row][p] * jacobian[row][q];
}
}
}
// 互补窗口之和恢复全局 J^T*J 和 J^T*r，不能对各窗口再单独平均。
TrustRegionFisher& global = information[kAutoFitTimeWindowCount];
for(int p = 0; p < dimensions; ++p) {
global.gradient[p] += local.gradient[p];
for(int q = 0; q <= p; ++q) {
local.matrix[q][p] = local.matrix[p][q];
global.matrix[p][q] += local.matrix[p][q];
global.matrix[q][p] = global.matrix[p][q];
}
}
for(int q = 0; q < p; ++q)
equations.matrix[q][p] = equations.matrix[p][q];
}
return information;
return equations;
}
static bool buildTrustRegionFisherStep(
const TrustRegionFisher& global,
static bool buildTrustRegionStep(
const TrustRegionNormalEquations& global,
const QVector<int>& selected,
const QVector<double>& coordinates,
double damping, double trustRadius, double minimumStep,
@ -997,15 +810,8 @@ void nmCalculationAutoFitLM::writeTraceHeader()
<< "late_trend_reliable"
<< "registration_ambiguous";
for(int k = 0; k < kAutoFitTimeWindowCount; ++k) {
const QString prefix = QString("window_%1_").arg(k + 1);
cols << prefix + "time_min" << prefix + "time_max"
<< prefix + "weight_sum" << prefix + "rms_error" << prefix + "energy";
}
cols << "sampling_mode" << "sampling_points"
<< "pressure_vertical_bias" << "derivative_vertical_bias" << "first_window_shape_loss"
<< "early_value_loss" << "early_parallel_loss" << "early_parallel_bias" << "early_gap_loss" << "early_gap_bias";
<< "pressure_vertical_bias" << "derivative_vertical_bias";
QTextStream out(&m_traceFile);
out << cols.join(",") << "\n";
}
@ -1042,8 +848,8 @@ void nmCalculationAutoFitLM::writeTraceMetaFile()
QTextStream out(&metaFile);
out << "{\n";
out << " \"schema_version\": 54,\n";
out << " \"strategy\": \"permeability_height_then_joint_lm_shape_then_total\",\n";
out << " \"schema_version\": 55,\n";
out << " \"strategy\": \"permeability_height_then_joint_lm_total_for_vertical_well_else_shape_then_total\",\n";
out << " \"height_acceptance\": \"reliable_vertical_loss_decrease; no_shape_or_total_loss_constraint\",\n";
out << " \"shape_metric\": \"pressure_and_derivative_log_slopes_shared_sampling_grid_lag_round_10_percent_min_1\",\n";
out << " \"shape_parameter_selection\": \"all_valid_free_columns_joint_LM; no_single_parameter_sweep_or_wellbore_recheck\",\n";
@ -1076,10 +882,8 @@ void nmCalculationAutoFitLM::writeTraceMetaFile()
out << " \"max_iterations\": " << m_maxIterations << ",\n";
out << " \"target_error\": " << jsonNumber(m_targetError) << "\n";
out << " },\n";
// 区间尚未冻结时写 null；首次有效评价后重写元数据，保存实际窗口基准。
out << " \"time_windows\": {\n";
out << " \"count\": " << kAutoFitTimeWindowCount << ",\n";
out << " \"overlap_ratio\": " << jsonNumber(kAutoFitTimeWindowOverlapRatio) << ",\n";
// 区间尚未冻结时写 null；首次有效评价后记录固定的误差比较范围。
out << " \"comparison_time_range\": {\n";
out << " \"comparison_time_min\": "
<< (m_comparisonTimeMin > 0.0 ? jsonNumber(m_comparisonTimeMin) : QString("null")) << ",\n";
out << " \"comparison_time_max\": "
@ -1155,32 +959,10 @@ void nmCalculationAutoFitLM::writeTraceRow(
}
}
// 无效候选也补齐窗口列，保证轨迹每行结构一致。
for(int k = 0; k < kAutoFitTimeWindowCount; ++k) {
if(objectiveBreakdown && objectiveBreakdown->valid &&
k < objectiveBreakdown->timeWindows.size()) {
const AutoFitTimeWindowLM& window = objectiveBreakdown->timeWindows[k];
cols << traceNumber(window.timeMin) << traceNumber(window.timeMax)
<< traceNumber(window.weightSum) << traceNumber(window.rmsError)
<< traceNumber(window.energy);
} else {
for(int i = 0; i < 5; ++i) {
cols << QString();
}
}
}
cols << "log_uniform_per_decade"
<< (objectiveBreakdown ? QString::number(objectiveBreakdown->residualVector.size() / 2) : QString())
<< (objectiveBreakdown ? traceNumber(objectiveBreakdown->pressureVerticalBias) : QString())
<< (objectiveBreakdown ? traceNumber(objectiveBreakdown->derivativeVerticalBias) : QString())
<< (objectiveBreakdown && objectiveBreakdown->valid
? traceNumber(qSqrt(trustRegionEarlyShapeEnergy(*objectiveBreakdown))) : QString())
<< (objectiveBreakdown && objectiveBreakdown->valid ? traceNumber(objectiveBreakdown->earlyValueLoss) : QString())
<< (objectiveBreakdown && objectiveBreakdown->valid ? traceNumber(objectiveBreakdown->earlyParallelLoss) : QString())
<< (objectiveBreakdown && objectiveBreakdown->valid ? traceNumber(objectiveBreakdown->earlyParallelBias) : QString())
<< (objectiveBreakdown && objectiveBreakdown->valid ? traceNumber(trustRegionEarlyGapLoss(*objectiveBreakdown)) : QString())
<< (objectiveBreakdown && objectiveBreakdown->valid ? traceNumber(trustRegionEarlyGapBias(*objectiveBreakdown)) : QString());
<< (objectiveBreakdown ? traceNumber(objectiveBreakdown->derivativeVerticalBias) : QString());
QTextStream out(&m_traceFile);
out << cols.join(",") << "\n";
m_traceFile.flush();
@ -1824,6 +1606,10 @@ bool nmCalculationAutoFitLM::evaluateTrustRegionPoint(
StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting()
{
// 本次拟合开始时按目标井的实际井型确定流程，普通直井跳过形状调整。
nmDataAnalyzeManager* dataManager = nmDataAnalyzeManager::getCurrentInstance();
nmDataWellBase* targetWell = dataManager ? dataManager->findWellByName(m_targetWellName) : nullptr;
const bool skipShapeStage = targetWell && targetWell->getWellType() == NM_WELL_MODEL::Vertical_Well;
const int dimensions = getEnabledParameterCount();
if(dimensions <= 0 || m_enabledParamIndices.size() != dimensions) {
m_lastError = tr("No valid parameters are available for trust-region fitting");
@ -1863,11 +1649,10 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting()
const char* rebuildReason = QT_TR_NOOP("initial sensitivity model");
bool modelRebuiltAtMinimumRadius = false;
bool stagnationConfirmationRequested = false;
bool globalFallbackAttempted = false;
StopReasonLM stopReason = LM_MAX_ITERATIONS;
// jacobian 的行对应固定采样的残差，列对应用户勾选的参数。
// Fisher 按当前目标取数值或形状残差行，两段均使用所有勾选参数。
// LM 按当前目标取数值或形状残差行，两段均使用所有勾选参数。
QVector<QVector<double> > jacobian;
QVector<bool> jacobianColumnValid(dimensions, false);
@ -2058,8 +1843,8 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting()
if(m_shouldStop) return LM_USER_STOPPED;
if(completedIterations >= m_maxIterations) return LM_MAX_ITERATIONS;
// 原逐参数第二阶段已移除；同一个联合 LM 先匹配整体形状，再匹配总误差。
bool shapeStage = true;
// 直井直接匹配总误差，其他井型仍先匹配整体形状，再匹配总误差。
bool shapeStage = !skipShapeStage;
auto stageError = [&](const TrustRegionEvaluation& point) -> double {
return shapeStage ? point.breakdown.shapeLoss : point.fitness;
};
@ -2086,7 +1871,6 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting()
}
// 有效改善后更新累计基准，并清除停滞计数。
globalFallbackAttempted = false;
effectiveImprovementBaseline = fitness;
consecutiveIneffectiveSteps = 0;
stagnationConfirmationRequested = false;
@ -2096,10 +1880,6 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting()
// 两个目标均沿用原 LM 的有效改善门槛及重建后停滞确认。
auto recordIneffectiveStep = [&]() -> bool {
++consecutiveIneffectiveSteps;
if(!globalFallbackAttempted) {
return false;
}
globalFallbackAttempted = false;
if(consecutiveIneffectiveSteps < maximumIneffectiveSteps) {
return false;
}
@ -2308,7 +2088,8 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting()
return true;
};
for(int phase = 0; phase < 2 && !m_shouldStop && completedIterations < m_maxIterations; ++phase) {
for(int phase = skipShapeStage ? 1 : 0;
phase < 2 && !m_shouldStop && completedIterations < m_maxIterations; ++phase) {
shapeStage = phase == 0;
phaseIterations = 0;
maximumEvaluations = m_totalEvaluations + phaseEvaluationBudget;
@ -2323,7 +2104,6 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting()
consecutivePoorPredictions = 0;
modelRebuiltAtMinimumRadius = false;
stagnationConfirmationRequested = false;
globalFallbackAttempted = false;
jacobian.clear();
rebuildRequested = true;
rebuildReason = shapeStage ? QT_TR_NOOP("full sensitivity at shape LM entry")
@ -2336,7 +2116,8 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting()
DEBUG_OUT(tr("LM phase budget: %1 remaining shared iterations, %2 evaluations.")
.arg(m_maxIterations - completedIterations).arg(phaseEvaluationBudget));
writeTraceRow(completedIterations, -1, "stage_switch", current.parameters, current.fitness,
true, 0, shapeStage ? "height_to_shape_lm" : "shape_to_total", &current.breakdown);
true, 0, shapeStage ? "height_to_shape_lm" : (skipShapeStage ? "height_to_total" : "shape_to_total"),
&current.breakdown);
while(!m_shouldStop && completedIterations < m_maxIterations && m_totalEvaluations < maximumEvaluations) {
m_currentIteration = completedIterations;
@ -2374,18 +2155,8 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting()
? current.breakdown.shapeResiduals : current.breakdown.residualVector;
const int rowOffset = shapeStage ? current.breakdown.residualVector.size() : 0;
const QVector<QVector<double> > objectiveJacobian = jacobian.mid(rowOffset, objectiveResidual.size());
QVector<double> objectiveCoordinates;
if(shapeStage) {
// Fisher 窗口与形状残差使用同一组斜率区间中心。
const int pointCount = current.breakdown.residualVector.size() / 2;
const int lag = autoFitShapeLag(pointCount);
for(int i = 0; i < pointCount - lag; ++i)
objectiveCoordinates.append((i + lag * 0.5) / (pointCount - 1));
}
const QVector<TrustRegionFisher> information = buildTrustRegionFisher(
objectiveJacobian, objectiveResidual, jacobianColumnValid, objectiveCoordinates);
// 两段都使用当前目标的全局信息，不按时间窗口另行选参。
const TrustRegionFisher& global = information[kAutoFitTimeWindowCount];
const TrustRegionNormalEquations global = buildTrustRegionNormalEquations(
objectiveJacobian, objectiveResidual, jacobianColumnValid);
QVector<int> selectedColumns;
QVector<double> coordinateStep;
double predictedReduction = 0.0;
@ -2395,8 +2166,7 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting()
if(jacobianColumnValid[column] && global.matrix[column][column] > 0.0)
selectedColumns.append(column);
}
globalFallbackAttempted = true;
if(selectedColumns.isEmpty() || !buildTrustRegionFisherStep(
if(selectedColumns.isEmpty() || !buildTrustRegionStep(
global, selectedColumns, current.coordinates, damping,
trustRadius, minimumCoordinateStep, &coordinateStep, &predictedReduction)) {
selectedColumns.clear();
@ -2426,7 +2196,7 @@ StopReasonLM nmCalculationAutoFitLM::runTrustRegionFitting()
for(int column = 0; column < dimensions; ++column) {
candidateCoordinates[column] += coordinateStep[column];
}
// 窗口和参数索引仅写入已有跟踪文件，不向拟合日志窗口增加分段信息。
// 参数索引保留在跟踪文件中，便于检查每次联合调整的参数。
QStringList selectedParameterIndices;
for(int i = 0; i < selectedColumns.size(); ++i) {
selectedParameterIndices << QString::number(m_enabledParamIndices[selectedColumns[i]]);
@ -3300,8 +3070,8 @@ double nmCalculationAutoFitLM::calculateLogLogCurveError(
const QVector<QVector<double> >& target,
const QVector<QVector<double> >& result)
{
// 首次有效评价后固定公共时间范围。窗口与上下、左右、形状只负责诊断，
// 主目标仍由完整压力和导数残差计算，避免窗口重叠造成重复计权。
// 首次有效评价后固定公共时间范围，保证各候选在相同时间点比较。
// 总误差由完整压力和导数残差计算，形状阶段单独使用斜率残差。
// 整个计算过程均位于 log(time)-log(value) 坐标。
const double invalidLoss = 1.0e10;
const double valueFloor = 1.0e-12;
@ -3427,7 +3197,7 @@ double nmCalculationAutoFitLM::calculateLogLogCurveError(
}
// 仅首次有效评价使用交集建立基准；失败试算不能冻结区间，后续候选
// 必须覆盖完整基准，不允许靠丢失首尾点缩小误差或改变窗口位置。
// 必须覆盖完整基准，不允许靠丢失首尾点缩小误差。
const bool comparisonRangeFixed = m_comparisonTimeMin > 0.0;
double overlapMinX = comparisonRangeFixed
? m_comparisonTimeMin : qMax(targetMinX, resultMinX);
@ -3453,7 +3223,6 @@ double nmCalculationAutoFitLM::calculateLogLogCurveError(
const int lag = autoFitShapeLag(count);
const double span = qLn(overlapMaxX) - qLn(overlapMinX);
QVector<double> residuals[2];
QVector<double> 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) {
@ -3466,9 +3235,7 @@ 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);
}
@ -3489,7 +3256,6 @@ double nmCalculationAutoFitLM::calculateLogLogCurveError(
}
}
objective->shapeLoss = qSqrt(trustRegionSquaredNorm(objective->shapeResiduals));
populateEarlyWellboreMetrics(objective, targetLogs, resultLogs, span);
return isFiniteNumber(objective->shapeLoss);
};
@ -4119,8 +3885,6 @@ double nmCalculationAutoFitLM::calculateLogLogCurveError(
isFiniteNumber(breakdown.total) &&
breakdown.total >= 0.0;
if(breakdown.valid && breakdown.total < 1.0e9) {
breakdown.timeWindows = calculateAutoFitTimeWindows(
breakdown.residualVector, overlapMinX, overlapMaxX);
if(!comparisonRangeFixed) {
// 整个目标及诊断均有效后才提交基准，中点回退可重新建立区间。
m_comparisonTimeMin = overlapMinX;

Write Preview
Loading…
Cancel
Save