feat(nmNum): 增加 Monte-Carlo 压力误差计算

- 使用求解快照中的实测压力,插值计算每组 RMSE。
- 失败、取消或无有效数据时标记误差无效。
- 结果暂存内存,供后续绘图和筛选使用。
feature/Sensitivity-20261012
lvjunjie 1 week ago
parent 15b684a329
commit 04a8759b64

@ -103,12 +103,13 @@ public:
// 运行结果属于分析数据,求解器与界面共用同一结构,避免重复复制整批曲线。
struct RunResult
{
RunResult() : m_nRunIndex(0), m_bSucceeded(false), m_bCancelled(false), m_bVisible(true) {}
RunResult() : m_nRunIndex(0), m_dPressureRmse(-1.0), m_bSucceeded(false), m_bCancelled(false), m_bVisible(true) {}
quint64 m_nRunIndex;
QVector<double> m_vecParameters;
QVector<QVector<double> > m_vecPressure;
QVector<QVector<double> > m_vecLogLog;
QVector<QVector<double> > m_vecSemiLog;
double m_dPressureRmse; // Monte-Carlo 压力均方根误差(MPa),-1 表示尚无有效误差。
bool m_bSucceeded;
bool m_bCancelled;
bool m_bVisible;

@ -211,6 +211,46 @@ bool isFiniteSolverNumber(double value)
#endif
}
double calculateSensitivityPressureRmse(const QVector<QVector<double> >& measured,
const QVector<QVector<double> >& simulated,
const QAtomicInt* cancelRequested)
{
if(measured.size() < 2 || simulated.size() < 2 || simulated[0].isEmpty() ||
measured[0].size() != measured[1].size() ||
simulated[0].size() != simulated[1].size()) return -1.0;
const QVector<double>& times = simulated[0];
const QVector<double>& pressures = simulated[1];
for(int i = 0; i < times.size(); ++i) {
if((i % 256) == 0 && isTaskCancellationRequested(cancelRequested)) return -1.0;
if(!isFiniteSolverNumber(times[i]) || !isFiniteSolverNumber(pressures[i]) ||
(i > 0 && times[i] <= times[i - 1])) return -1.0;
}
// 在实测时间上线性插值原始压力,只统计覆盖范围内的有效点,不外推。
double squaredError = 0.0;
int count = 0;
for(int i = 0; i < measured[0].size(); ++i) {
if((i % 256) == 0 && isTaskCancellationRequested(cancelRequested)) return -1.0;
const double time = measured[0][i];
const double pressure = measured[1][i];
if(!isFiniteSolverNumber(time) || !isFiniteSolverNumber(pressure) ||
time < times.first() || time > times.last()) continue;
const int right = static_cast<int>(
std::lower_bound(times.constBegin(), times.constEnd(), time) - times.constBegin());
double predicted = pressures[right];
if(times[right] != time) {
const double fraction = (time - times[right - 1]) / (times[right] - times[right - 1]);
predicted = pressures[right - 1] + fraction * (pressures[right] - pressures[right - 1]);
}
const double residual = predicted - pressure;
squaredError += residual * residual;
++count;
}
// 普通 Monte-Carlo 的各点权重均为 1;无有效点不能当作完全拟合。
return count > 0 && isFiniteSolverNumber(squaredError)
? std::sqrt(squaredError / count) : -1.0;
}
bool savedPseudoPressure(const nmDataBinaryTools::NM_PEBI_SCENE& scene,
double pressure, double& result)
{
@ -1012,6 +1052,13 @@ bool nmCalculationDllPebiSolverTask::executeSensitivityBatch()
result.m_vecPressure = m_autoFitResultPressure;
result.m_vecLogLog = m_autoFitResultLogLog;
result.m_vecSemiLog = m_autoFitResultSemiLog;
if(m_pSensitivityResult->m_calcType == nmDataSensitive::CALC_MONTE_CARLO) {
// 使用求解时冻结的目标井实测压力,不读取界面后来修改的数据。
const int target = m_pInputSnapshot->m_oSensitivityBaseline.m_nTarget;
result.m_dPressureRmse = calculateSensitivityPressureRmse(
m_pInputSnapshot->m_oGridInput.m_vecWellInputs[target].m_vecHistoryPressure,
result.m_vecPressure, &m_nCancelRequested);
}
} else {
result.m_sFailureReason = m_sFailureReason;
}
@ -1023,6 +1070,7 @@ bool nmCalculationDllPebiSolverTask::executeSensitivityBatch()
result.m_vecPressure.clear();
result.m_vecLogLog.clear();
result.m_vecSemiLog.clear();
result.m_dPressureRmse = -1.0;
} else if(!result.m_bSucceeded) {
++failed;
}

Loading…
Cancel
Save