feat(nmCalculation): 建立 PEBI 静态属性计算基础

- 新增静态属性值对象、属性目录和独立 Kriging 计算服务
- 统一手工求解与自动拟合的 k/phi/h 输入捕获和单位转换
- 校验原始单元共同索引域并完善失败事务清理
- 统一建网、Kriging 和模型求解的授权路径快照
feature/nmNum-postprocessing-view-20260910
lh 1 month ago
parent bb11b64542
commit 36f17ae321

@ -0,0 +1,53 @@
#ifndef NMPEBISTATICPROPERTYCALCULATOR_H
#define NMPEBISTATICPROPERTYCALCULATOR_H
#include "nmCalculation_global.h"
#include "nmPebiStaticPropertyValues.h"
#include "nmPropertyInterpolationData.h"
#include <QAtomicInt>
#include <QPointF>
#include <QString>
#include <QVector>
enum nmPebiStaticPropertyCalculationStatus
{
NM_PebiStaticPropertyCalculation_Succeeded = 0,
NM_PebiStaticPropertyCalculation_Cancelled,
NM_PebiStaticPropertyCalculation_Failed
};
/** @brief 一次静态属性计算独占的不可变值输入。 */
struct NMCALCULATION_EXPORT nmPebiStaticPropertyCalculationInput
{
nmPebiStaticPropertyCalculationInput();
bool m_bCaptured; ///< 全部字段已经在 DataManager 所属线程完成冻结。
double m_dReferencePermeabilityMilliDarcy;
double m_dReferencePorosity;
double m_dReferenceThicknessMeters;
QVector<nmPropertyInterpolationDataSet> m_vecDataSets;
QString m_sLicensePath;
};
/**
* @brief 计算按原始 PEBI 求解单元顺序排列的 k/phi/h 静态属性。
*
* 服务只解释旧插值键 k、phi、h,不读取 DataManager、DLL 网格结构、VTK
* 或界面对象。失败和取消都不会向调用方提交部分数组。
*/
class NMCALCULATION_EXPORT nmPebiStaticPropertyCalculator
{
public:
static nmPebiStaticPropertyCalculationStatus calculate(
const QVector<QPointF>& vecOriginalCellPropertyLocations,
const nmPebiStaticPropertyCalculationInput& oInput,
nmPebiStaticPropertyValues& oOutput,
QString* pErrorMessage = NULL,
const QAtomicInt* pCancelRequested = NULL);
private:
nmPebiStaticPropertyCalculator();
};
#endif // NMPEBISTATICPROPERTYCALCULATOR_H

@ -13,6 +13,7 @@
#include "nmDataNumericalAnalysisCase.h"
#include "nmNumericalResultPersistence.h"
#include "nmDataGaugeRecord.h"
#include "nmPropertyInterpolationData.h"
#include "nmGridDefine.h"
@ -59,58 +60,6 @@ struct BackgroundImageInfo {
QImage objImage; // 具体的QImage对象
};
/* 插值测点数据 */
struct NM_DATA_EXPORT nmPropertyInterpolationPointData {
nmPropertyInterpolationPointData()
: x(0.0), y(0.0), value(0.0)
{
}
nmPropertyInterpolationPointData(double xValue, double yValue, double propertyValue)
: x(xValue), y(yValue), value(propertyValue)
{
}
double x;
double y;
double value;
};
/* 一组属性插值数据及其Kriging参数 */
struct NM_DATA_EXPORT nmPropertyInterpolationDataSet {
nmPropertyInterpolationDataSet()
: property("k"),
xDisplayUnit("m"),
yDisplayUnit("m"),
valueDisplayUnit("mD"),
rangeDisplayUnit("m"),
useForCalculation(false),
showPoints(true),
showLabels(true),
nugget(0.01),
sill(100.0),
range(1000.0),
model(0)
{
}
QString name;
QString property;
// 仅记录界面显示单位;坐标、Range及属性值始终按基准单位保存.
QString xDisplayUnit;
QString yDisplayUnit;
QString valueDisplayUnit;
QString rangeDisplayUnit;
bool useForCalculation;
QVector<nmPropertyInterpolationPointData> points;
bool showPoints;
bool showLabels;
double nugget;
double sill;
double range;
int model;
};
enum NM_SOLVER_MODEL_TYPE {
SMT_Oil_ConstPvt = 1,
SMT_Oil_VariablePvt = 2,

@ -0,0 +1,70 @@
#ifndef NMPEBIRESULTPROPERTYCATALOG_H
#define NMPEBIRESULTPROPERTYCATALOG_H
#include "nmData_global.h"
#include <QList>
#include <QString>
/** @brief 结果属性的数据时态。 */
enum NM_PEBI_RESULT_PROPERTY_TIME_TYPE
{
NM_PebiResultProperty_Static = 0,
NM_PebiResultProperty_Dynamic
};
/** @brief 结果属性的生产来源,不与界面文字或 VTK 类型绑定。 */
enum NM_PEBI_RESULT_PROPERTY_SOURCE
{
NM_PebiResultPropertySource_StaticCalculation = 0,
NM_PebiResultPropertySource_SolverOutput
};
/** @brief 结果属性与几何元素的逻辑关联。 */
enum NM_PEBI_RESULT_PROPERTY_ASSOCIATION
{
NM_PebiResultPropertyAssociation_CellData = 0
};
/** @brief 结果属性建议采用的数值范围策略。 */
enum NM_PEBI_RESULT_PROPERTY_RANGE_POLICY
{
NM_PebiResultPropertyRange_CurrentArray = 0,
NM_PebiResultPropertyRange_AllFrames
};
/** @brief 一项稳定结果属性的非界面契约。 */
struct NM_DATA_EXPORT nmPebiResultPropertyDescriptor
{
nmPebiResultPropertyDescriptor();
QString m_sPropertyId; ///< 不参与翻译的稳定业务 ID。
QString m_sDisplayNameSourceText; ///< 供后续界面翻译的英文源文本。
QString m_sBaseUnit; ///< 数据层基准单位;无量纲时为空。
QString m_sDefaultDisplayUnit; ///< 默认显示单位;无量纲时为空。
NM_PEBI_RESULT_PROPERTY_TIME_TYPE m_eTimeType;
NM_PEBI_RESULT_PROPERTY_SOURCE m_eSource;
NM_PEBI_RESULT_PROPERTY_ASSOCIATION m_eAssociation;
NM_PEBI_RESULT_PROPERTY_RANGE_POLICY m_eRangePolicy;
};
/**
* @brief 提供第一批 PEBI 结果属性的稳定只读目录。
*
* 目录只描述数据契约,不创建 LUT、VTK 数组或界面控件。旧工程中的
* 插值键 k/phi/h 由计算服务独立解释,不能使用本目录放宽其合法范围。
*/
class NM_DATA_EXPORT nmPebiResultPropertyCatalog
{
public:
static QList<nmPebiResultPropertyDescriptor> descriptors();
static bool findDescriptor(
const QString& sPropertyId,
nmPebiResultPropertyDescriptor& oDescriptor);
static bool isKnownProperty(const QString& sPropertyId);
private:
nmPebiResultPropertyCatalog();
};
#endif // NMPEBIRESULTPROPERTYCATALOG_H

@ -0,0 +1,30 @@
#ifndef NMPEBISTATICPROPERTYVALUES_H
#define NMPEBISTATICPROPERTYVALUES_H
#include "nmData_global.h"
#include <QVector>
/**
* @brief 按原始 PEBI 求解单元索引保存的三项静态属性。
*
* 该类型不包含 VTK、求解器或界面对象。三个数组共享同一个原始单元
* 索引域,其中渗透率使用数据层单位 mD,孔隙度无量纲,厚度使用 m。
*/
struct NM_DATA_EXPORT nmPebiStaticPropertyValues
{
/** @brief 清空三项属性,使对象回到未计算状态。 */
void clear();
/**
* @brief 检查数组长度和全部数值是否可用于指定原始单元集合。
* @param nOriginalCellCount 原始 PEBI 求解单元数量,必须大于零。
*/
bool isValidForCellCount(int nOriginalCellCount) const;
QVector<double> m_vecPermeabilityMilliDarcy; ///< 渗透率,单位 mD。
QVector<double> m_vecPorosity; ///< 孔隙度,无量纲。
QVector<double> m_vecThicknessMeters; ///< 厚度,单位 m。
};
#endif // NMPEBISTATICPROPERTYVALUES_H

@ -0,0 +1,65 @@
#ifndef NMPROPERTYINTERPOLATIONDATA_H
#define NMPROPERTYINTERPOLATIONDATA_H
#include "nmData_global.h"
#include <QString>
#include <QVector>
/** @brief 一条按数据层基准单位保存的属性插值测点。 */
struct NM_DATA_EXPORT nmPropertyInterpolationPointData
{
nmPropertyInterpolationPointData()
: x(0.0), y(0.0), value(0.0)
{
}
nmPropertyInterpolationPointData(double xValue,
double yValue,
double propertyValue)
: x(xValue), y(yValue), value(propertyValue)
{
}
double x; ///< X 坐标,使用数据层长度基准单位。
double y; ///< Y 坐标,使用数据层长度基准单位。
double value; ///< 属性值,使用对应属性的数据层基准单位。
};
/** @brief 一组属性插值测点、Kriging 参数和界面显示偏好。 */
struct NM_DATA_EXPORT nmPropertyInterpolationDataSet
{
nmPropertyInterpolationDataSet()
: property("k"),
xDisplayUnit("m"),
yDisplayUnit("m"),
valueDisplayUnit("mD"),
rangeDisplayUnit("m"),
useForCalculation(false),
showPoints(true),
showLabels(true),
nugget(0.01),
sill(100.0),
range(1000.0),
model(0)
{
}
QString name;
QString property;
// 以下显示单位只保存界面偏好;坐标、Range 和属性值始终保存基准值。
QString xDisplayUnit;
QString yDisplayUnit;
QString valueDisplayUnit;
QString rangeDisplayUnit;
bool useForCalculation;
QVector<nmPropertyInterpolationPointData> points;
bool showPoints;
bool showLabels;
double nugget;
double sill;
double range;
int model;
};
#endif // NMPROPERTYINTERPOLATIONDATA_H

@ -15,6 +15,7 @@
#include "nmDataFault.h"
#include "nmDataTimeStepSetting.h"
#include "nmPebiResultSnapshotBuilder.h"
#include "nmPebiStaticPropertyCalculator.h"
#include "nmCalculationPebiGrid.h"
#include "nmCalculationUtils.h"
@ -35,6 +36,7 @@
#include <stdexcept>
#include <cmath>
#include <float.h>
#include <limits.h>
#include <string>
#include <vtkSmartPointer.h>
@ -61,7 +63,8 @@ struct nmPebiSolverInputSnapshot
nmPebiGridInputSnapshot m_oGridInput; ///< 网格、场景、井顺序和授权路径值快照。
nmPebiGridResult m_oGridResult; ///< 已有网格副本或后台新生成的局部网格。
QVector<nmPropertyInterpolationDataSet> m_vecPropertyDataSets; ///< 属性插值值副本。
/** @brief 在 DataManager 线程冻结的静态属性参考值、插值数据和授权路径。 */
nmPebiStaticPropertyCalculationInput m_oStaticPropertyInput;
vtkSmartPointer<vtkUnstructuredGrid> m_pBaseGrid; ///< 完整成果使用的基础 VTK 网格。
vtkSmartPointer<vtkUnstructuredGrid> m_pSourceBaseGrid; ///< 主线程捕获引用,后台只读并深拷贝。
int m_nSolverType; ///< PEBI 求解器实现类型。
@ -258,10 +261,21 @@ bool captureResultWellMetadata(
return true;
}
bool captureResultParameters(
/**
* @brief 从同一个储层对象冻结静态属性输入及可选的正式结果储层参数。
*
* 自动拟合也必须捕获静态属性输入,但不需要正式结果元数据。只有全部字段都复制
* 成功后才设置 m_bCaptured,防止默认构造的零值被误认为有效业务输入。
*/
bool captureStaticPropertyInput(
nmDataAnalyzeManager* pDataManager,
nmPebiSolverInputSnapshot& oSnapshot)
nmPebiSolverInputSnapshot& oSnapshot,
bool bCaptureResultReservoir,
nmDataReservoir** ppReservoir)
{
if(ppReservoir != NULL) {
*ppReservoir = NULL;
}
if(pDataManager == NULL ||
QThread::currentThread() != pDataManager->thread()) {
return false;
@ -272,26 +286,61 @@ bool captureResultParameters(
return false;
}
nmPebiResultReservoirParameters& oReservoir =
oSnapshot.m_oResultReservoirParameters;
oReservoir.m_dInitialPressure = pReservoir->getInitialPressure()
.getValue().toDouble();
oReservoir.m_dPermeability = pReservoir->getPermeability()
.getValue().toDouble();
oReservoir.m_dThickness = pReservoir->getThickness()
.getValue().toDouble();
oReservoir.m_dPorosity = pReservoir->getPorosity()
.getValue().toDouble();
oReservoir.m_dTotalCompressibility = pReservoir->getCt()
.getValue().toDouble();
oReservoir.m_dRockCompressibility = pReservoir->getCf()
.getValue().toDouble();
oReservoir.m_dOilSaturation = pReservoir->getSoi()
.getValue().toDouble();
oReservoir.m_dGasSaturation = pReservoir->getSgi()
.getValue().toDouble();
oReservoir.m_dWaterSaturation = pReservoir->getSwi()
.getValue().toDouble();
nmPebiStaticPropertyCalculationInput oStaticInput;
// 数据层保持显示和持久化基准单位;渗透率只在组装求解器 Base.k 时转成 D。
oStaticInput.m_dReferencePermeabilityMilliDarcy =
pReservoir->getPermeability().getValue().toDouble();
oStaticInput.m_dReferencePorosity =
pReservoir->getPorosity().getValue().toDouble();
oStaticInput.m_dReferenceThicknessMeters =
pReservoir->getThickness().getValue().toDouble();
oStaticInput.m_vecDataSets =
pDataManager->getPropertyInterpolationDataSets();
// 授权路径随网格输入最先冻结;Kriging、建网和模型求解必须复用同一份值,
// 不能在分阶段捕获期间从 DataManager 再读出另一个时刻的路径。
oStaticInput.m_sLicensePath =
oSnapshot.m_oGridInput.m_sLicensePath;
oStaticInput.m_bCaptured = true;
if(bCaptureResultReservoir) {
nmPebiResultReservoirParameters oReservoir;
oReservoir.m_dInitialPressure = pReservoir->getInitialPressure()
.getValue().toDouble();
oReservoir.m_dPermeability =
oStaticInput.m_dReferencePermeabilityMilliDarcy;
oReservoir.m_dThickness =
oStaticInput.m_dReferenceThicknessMeters;
oReservoir.m_dPorosity = oStaticInput.m_dReferencePorosity;
oReservoir.m_dTotalCompressibility = pReservoir->getCt()
.getValue().toDouble();
oReservoir.m_dRockCompressibility = pReservoir->getCf()
.getValue().toDouble();
oReservoir.m_dOilSaturation = pReservoir->getSoi()
.getValue().toDouble();
oReservoir.m_dGasSaturation = pReservoir->getSgi()
.getValue().toDouble();
oReservoir.m_dWaterSaturation = pReservoir->getSwi()
.getValue().toDouble();
oSnapshot.m_oResultReservoirParameters = oReservoir;
}
oSnapshot.m_oStaticPropertyInput = oStaticInput;
if(ppReservoir != NULL) {
*ppReservoir = pReservoir;
}
return true;
}
/** @brief 捕获完整结果专用的 PVT 和求解设置,不重复查询储层对象。 */
bool captureResultParameters(
nmDataAnalyzeManager* pDataManager,
nmDataReservoir* pReservoir,
nmPebiSolverInputSnapshot& oSnapshot)
{
if(pDataManager == NULL || pReservoir == NULL ||
QThread::currentThread() != pDataManager->thread()) {
return false;
}
nmPebiResultPvtParameters& oPvt =
oSnapshot.m_oResultPvtParameters;
@ -359,122 +408,54 @@ bool captureResultParameters(
return true;
}
// 将启用的数据组插值到全部网格单元中心,并覆盖对应的求解器属性数组.
bool applyPropertyInterpolation(
HX_NWTM_MODEL_INPUT& modelInput,
const QVector<nmPropertyInterpolationDataSet>& dataSets,
const QString& licensePath,
QString& errorMessage,
/**
* @brief 校验 DLL 三个原始单元数组共享的结构索引域并提取属性计算位置。
*
* 数量一致只能证明三个数组可按共同下标访问。下标 i 是否表示同一个求解单元
* 由 DLL 接口契约保证,应用层不重新排序,也不把 Trinodexy[i] 宣称为几何质心。
*/
bool extractStaticPropertyLocations(
const HX_NWTM_GRID_OUTPUT1& oGridOutput1,
const HX_NWTM_GRID_OUTPUT2& oGridOutput2,
QVector<QPointF>& vecLocations,
QString& sErrorMessage,
const QAtomicInt* pCancelRequested)
{
bool hasEnabledDataSet = false;
for(int i = 0; i < dataSets.size(); ++i) {
if(dataSets[i].useForCalculation) {
hasEnabledDataSet = true;
break;
}
}
if(!hasEnabledDataSet) {
return true;
vecLocations.clear();
sErrorMessage.clear();
const size_t nLocationCount = oGridOutput2.Trinodexy.size();
const size_t nIsPlotCount = oGridOutput1.PEBI_cell.isplot.size();
const size_t nPolygonIndexCount = oGridOutput1.PEBI_cell.pindex.size();
if(nLocationCount == 0 ||
nLocationCount != nIsPlotCount ||
nLocationCount != nPolygonIndexCount ||
nLocationCount > static_cast<size_t>(INT_MAX)) {
sErrorMessage = QLatin1String(
"PEBI property index arrays are empty or inconsistent.");
return false;
}
QVector<QPointF> targetPoints;
targetPoints.reserve(static_cast<int>(modelInput.GRID.Trinodexy.size()));
for(size_t cellIndex = 0;
cellIndex < modelInput.GRID.Trinodexy.size(); ++cellIndex) {
if((cellIndex % 256) == 0 &&
vecLocations.reserve(static_cast<int>(nLocationCount));
for(size_t nIndex = 0; nIndex < nLocationCount; ++nIndex) {
if((nIndex % 256) == 0 &&
isTaskCancellationRequested(pCancelRequested)) {
vecLocations.clear();
return false;
}
const dVec1& cellPosition = modelInput.GRID.Trinodexy[cellIndex];
if(cellPosition.size() < 2) {
errorMessage = QString("Grid cell %1 has no valid center coordinate.")
.arg(static_cast<qulonglong>(cellIndex));
return false;
}
targetPoints.append(QPointF(cellPosition[0], cellPosition[1]));
}
bool kApplied = false;
bool phiApplied = false;
bool hApplied = false;
for(int dataSetIndex = 0; dataSetIndex < dataSets.size(); ++dataSetIndex) {
if(isTaskCancellationRequested(pCancelRequested)) {
// Trinodexy 提供属性计算位置;本阶段只读取 X/Y,不涉及 TRI_cell。
const dVec1& vecPosition = oGridOutput2.Trinodexy[nIndex];
if(vecPosition.size() < 2 ||
!isFiniteSolverNumber(vecPosition[0]) ||
!isFiniteSolverNumber(vecPosition[1])) {
sErrorMessage = QString::fromLatin1(
"PEBI property location %1 is invalid.")
.arg(static_cast<qulonglong>(nIndex));
vecLocations.clear();
return false;
}
const nmPropertyInterpolationDataSet& dataSet = dataSets[dataSetIndex];
if(!dataSet.useForCalculation) {
continue;
}
dVec1* solverValues = NULL;
bool* propertyApplied = NULL;
if(dataSet.property == "k") {
solverValues = &modelInput.Base.k;
propertyApplied = &kApplied;
} else if(dataSet.property == "phi") {
solverValues = &modelInput.Base.phi;
propertyApplied = &phiApplied;
} else if(dataSet.property == "h") {
solverValues = &modelInput.Base.h;
propertyApplied = &hApplied;
} else {
errorMessage = QString("Dataset '%1' has an unknown property.")
.arg(dataSet.name);
return false;
}
if(*propertyApplied) {
errorMessage = QString(
"More than one dataset is enabled for property %1.")
.arg(dataSet.property);
return false;
}
QVector<QPointF> measurementPoints;
QVector<double> measurementValues;
measurementPoints.reserve(dataSet.points.size());
measurementValues.reserve(dataSet.points.size());
for(int pointIndex = 0; pointIndex < dataSet.points.size(); ++pointIndex) {
const nmPropertyInterpolationPointData& point = dataSet.points[pointIndex];
measurementPoints.append(QPointF(point.x, point.y));
measurementValues.append(point.value);
}
QVector<double> interpolationValues;
QString calculationError;
if(!nmCalculationUtils::calculateKriging(
targetPoints,
measurementPoints,
measurementValues,
dataSet.nugget,
dataSet.sill,
dataSet.range,
dataSet.model,
licensePath,
interpolationValues,
&calculationError,
pCancelRequested)) {
errorMessage = QString("Dataset '%1': %2")
.arg(dataSet.name)
.arg(calculationError);
return false;
}
solverValues->resize(interpolationValues.size());
for(int valueIndex = 0; valueIndex < interpolationValues.size(); ++valueIndex) {
if((valueIndex % 256) == 0 &&
isTaskCancellationRequested(pCancelRequested)) {
return false;
}
// 属性插值数据层以mD保存;渗透率在写入PEBI数组时转换为D。
(*solverValues)[valueIndex] = dataSet.property == "k"
? nmCalculationUtils::milliDarcyToDarcy(
interpolationValues[valueIndex])
: interpolationValues[valueIndex];
}
*propertyApplied = true;
vecLocations.append(QPointF(vecPosition[0], vecPosition[1]));
}
return true;
@ -485,7 +466,9 @@ bool applyPropertyInterpolation(
*/
bool buildModelInputFromSnapshot(
const nmPebiSolverInputSnapshot& oSnapshot,
const QVector<QPointF>& vecPropertyLocations,
HX_NWTM_MODEL_INPUT& oModelInput,
nmPebiStaticPropertyValues& oStaticProperties,
QString& sErrorMessage,
const QAtomicInt* pCancelRequested)
{
@ -535,7 +518,7 @@ bool buildModelInputFromSnapshot(
oModelInput.PVT.Sw = oScene.PVT.Sw;
oModelInput.PVT.Krw = oScene.PVT.Krw;
// 第三步:先用储层参考值填满全部网格单元,再按捕获的属性数据组覆盖。
// 第三步:恢复与静态属性无关的储层和时间步参数。
oModelInput.Base.Pi = oScene.Base.Pi;
oModelInput.Base.Cti = oScene.Base.Cti;
oModelInput.Base.Cf = oScene.Base.Cf;
@ -546,17 +529,45 @@ bool buildModelInputFromSnapshot(
oModelInput.Base.dt_Min = oScene.Base.dt_Min;
oModelInput.Base.dt_Max = oScene.Base.dt_Max;
const size_t nCellCount = oModelInput.GRID.Trinodexy.size();
oModelInput.Base.k = dVec1(nCellCount, oScene.Base.k_ref);
oModelInput.Base.phi = dVec1(nCellCount, oScene.Base.phi_ref);
oModelInput.Base.h = dVec1(nCellCount, oScene.Base.h_ref);
return applyPropertyInterpolation(
oModelInput,
oSnapshot.m_vecPropertyDataSets,
oSnapshot.m_oGridInput.m_sLicensePath,
sErrorMessage,
// 第四步:Calculator 只生成原始单元顺序的业务值,不接触 DLL 模型结构。
const nmPebiStaticPropertyCalculationStatus eStatus =
nmPebiStaticPropertyCalculator::calculate(
vecPropertyLocations,
oSnapshot.m_oStaticPropertyInput,
oStaticProperties,
&sErrorMessage,
pCancelRequested);
if(eStatus != NM_PebiStaticPropertyCalculation_Succeeded) {
return false;
}
const int nCellCount = vecPropertyLocations.size();
dVec1 vecPermeability(static_cast<size_t>(nCellCount));
dVec1 vecPorosity(static_cast<size_t>(nCellCount));
dVec1 vecThickness(static_cast<size_t>(nCellCount));
for(int nIndex = 0; nIndex < nCellCount; ++nIndex) {
if((nIndex % 256) == 0 &&
isTaskCancellationRequested(pCancelRequested)) {
return false;
}
// 静态属性层统一保存 mD;求解器 Base.k 的唯一写入边界转换为 D。
vecPermeability[static_cast<size_t>(nIndex)] =
nmCalculationUtils::milliDarcyToDarcy(
oStaticProperties.m_vecPermeabilityMilliDarcy[nIndex]);
vecPorosity[static_cast<size_t>(nIndex)] =
oStaticProperties.m_vecPorosity[nIndex];
vecThickness[static_cast<size_t>(nIndex)] =
oStaticProperties.m_vecThicknessMeters[nIndex];
}
if(isTaskCancellationRequested(pCancelRequested)) {
return false;
}
// 三个临时数组全部完成后再一起交换,失败和取消不会留下半套求解输入。
oModelInput.Base.k.swap(vecPermeability);
oModelInput.Base.phi.swap(vecPorosity);
oModelInput.Base.h.swap(vecThickness);
return true;
}
bool isReasonableLogLogValue(double value)
@ -688,6 +699,11 @@ void nmCalculationDllPebiSolverTask::run()
discardPendingFullResult();
bSucceeded = false;
}
if(!bSucceeded) {
// execute() 的所有普通失败也必须释放压力帧、网格别名和候选快照。
// 该清理函数允许在取消或异常分支再次调用。
discardPendingFullResult();
}
if(isCancelRequested()) {
// 停止只丢弃任务局部结果,主线程不会进入 commitResult()。
m_nWasCancelled.fetchAndStoreOrdered(1);
@ -899,15 +915,20 @@ bool nmCalculationDllPebiSolverTask::captureManualInputStep(
if(pState->m_ePhase ==
nmPebiManualCaptureState::CapturePhase_SolverSettings) {
reportStage(tr("Capturing solver settings..."), 0, 0);
m_pInputSnapshot->m_vecPropertyDataSets =
m_pDataManager->getPropertyInterpolationDataSets();
m_pInputSnapshot->m_nSolverType =
m_pDataManager->getPebiSolverType();
m_pInputSnapshot->m_nOmpThreads =
m_pDataManager->getPebiOmpThreads();
m_pInputSnapshot->m_nIluReuseSteps =
m_pDataManager->getPebiIluReuseSteps();
if(!captureResultParameters(m_pDataManager,
nmDataReservoir* pCapturedReservoir = NULL;
// 手工求解无条件冻结静态输入;正式结果参数复用同一次储层对象读取。
if(!captureStaticPropertyInput(m_pDataManager,
*m_pInputSnapshot,
true,
&pCapturedReservoir) ||
!captureResultParameters(m_pDataManager,
pCapturedReservoir,
*m_pInputSnapshot)) {
return fail(tr("Failed to capture solver settings."));
}
@ -1059,18 +1080,25 @@ bool nmCalculationDllPebiSolverTask::captureInputSnapshot(
}
}
// 第三步:复制属性插值及 DLL 求解配置。这些设置只影响模型求解结果,
// 不应在工作线程中再次从 DataManager 查询。
m_pInputSnapshot->m_vecPropertyDataSets =
m_pDataManager->getPropertyInterpolationDataSets();
// 第三步:冻结 DLL 求解配置和静态属性输入。静态输入必须在自动拟合分支
// 判断之外捕获,使每个粒子任务都读取候选参数写入后的最新储层参考值。
m_pInputSnapshot->m_nSolverType =
m_pDataManager->getPebiSolverType();
m_pInputSnapshot->m_nOmpThreads =
m_pDataManager->getPebiOmpThreads();
m_pInputSnapshot->m_nIluReuseSteps =
m_pDataManager->getPebiIluReuseSteps();
if(!m_pInputSnapshot->m_bAutoFitTargetOnly &&
!captureResultParameters(m_pDataManager, *m_pInputSnapshot)) {
nmDataReservoir* pCapturedReservoir = NULL;
const bool bCaptureFullResult =
!m_pInputSnapshot->m_bAutoFitTargetOnly;
if(!captureStaticPropertyInput(m_pDataManager,
*m_pInputSnapshot,
bCaptureFullResult,
&pCapturedReservoir) ||
(bCaptureFullResult &&
!captureResultParameters(m_pDataManager,
pCapturedReservoir,
*m_pInputSnapshot))) {
return fail(tr("Failed to capture solver settings."));
}
@ -1293,24 +1321,50 @@ bool nmCalculationDllPebiSolverTask::execPebiMode()
: fail(tr("The generated grid is invalid."));
}
// 第三步:由场景值快照重建完整模型输入。属性插值可以耗时,但它只读取任务
// 自己的点集副本,因此放在工作线程不会阻塞或竞争界面数据。
// 第三步:先验证 DLL 原始单元数组的共同结构索引域。自动拟合不创建 VTK,
// 因此这项约束必须位于求解任务内,不能只依赖 VTK 转换路径。
QVector<QPointF> vecPropertyLocations;
QString sPreparationError;
if(!extractStaticPropertyLocations(oGridResult.m_oGridOutput1,
oGridResult.m_oGridOutput2,
vecPropertyLocations,
sPreparationError,
pCancelRequested)) {
if(isCancelRequested()) {
return false;
}
const QString sLogMessage =
QString::fromLatin1("Invalid PEBI property index domain: %1")
.arg(sPreparationError);
qWarning() << sLogMessage;
zxLogInstance::getInstance()->writeLogF(sLogMessage);
return fail(tr("The generated grid is invalid."));
}
if(isCancelRequested()) {
return false;
}
// 第四步:由任务值快照重建完整模型输入。静态属性计算只读取任务自己的
// 位置和测点副本,因此工作线程不会访问或竞争 DataManager 对象。
// 直接初始化可避免先构造临时模型输入、再复制完整 GRID 的额外内存峰值。
HX_NWTM_MODEL_INPUT oModelInput(oGridResult.m_oGridOutput2);
QString sInterpolationError;
nmPebiStaticPropertyValues oStaticProperties;
reportStage(tr("Preparing model input..."), 0, 0);
if(!buildModelInputFromSnapshot(*m_pInputSnapshot,
vecPropertyLocations,
oModelInput,
sInterpolationError,
oStaticProperties,
sPreparationError,
pCancelRequested)) {
if(isCancelRequested()) {
return false;
}
const QString sLogMessage =
QString("Property interpolation failed: %1")
.arg(sInterpolationError);
QString::fromLatin1("Static property preparation failed: %1")
.arg(sPreparationError);
qWarning() << sLogMessage;
zxLogInstance::getInstance()->writeLogF(sLogMessage);
// 沿用既有已翻译提示;详细原因已经写入日志,不扩展本阶段 UI 文案。
return fail(tr("Property interpolation failed."));
}
@ -1318,7 +1372,8 @@ bool nmCalculationDllPebiSolverTask::execPebiMode()
return false;
}
// 第四步:建网、模型求解和 Kriging 共用 DLL 全局状态,整个 DLL 调用必须串行。
// 第五步:模型求解与前面的建网、Kriging 共用 DLL 全局状态。Kriging 已在
// calculateKriging() 内独立持锁,这里只串行保护模型求解入口。
reportStage(tr("Waiting for solver..."), 0, 0);
nmSolverDllMutexLocker oDllLocker;
if(!oDllLocker.lock(nmCalculationUtils::getHxNwtmDllMutex(),
@ -1910,11 +1965,13 @@ void nmCalculationDllPebiSolverTask::discardPendingFullResult()
m_pPendingResultSnapshot.clear();
m_vecPendingWellResults.clear();
m_mapPendingTimeSteps.clear();
m_autoFitResultPressure.clear();
m_autoFitResultLogLog.clear();
m_autoFitResultSemiLog.clear();
m_bPendingFullResultReady = false;
// 候选可能与输入基础网格共享 VTK 对象,丢弃或发布后必须清除任务侧别名。
if(m_pInputSnapshot != NULL &&
!m_pInputSnapshot->m_bAutoFitTargetOnly) {
if(m_pInputSnapshot != NULL) {
m_pInputSnapshot->m_pBaseGrid = NULL;
m_pInputSnapshot->m_pSourceBaseGrid = NULL;
m_pInputSnapshot->m_oGridResult.m_pUnstructuredGrid = NULL;

@ -0,0 +1,260 @@
#include "nmPebiStaticPropertyCalculator.h"
#include "nmCalculationUtils.h"
#include <QtGlobal>
#include <float.h>
#include <math.h>
#if !defined(_MSC_VER)
#include <cmath>
#endif
namespace {
enum StaticPropertyKind
{
StaticProperty_Permeability = 0,
StaticProperty_Porosity,
StaticProperty_Thickness,
StaticProperty_Invalid
};
bool isCancellationRequested(const QAtomicInt* pCancelRequested)
{
return pCancelRequested != NULL &&
static_cast<int>(*pCancelRequested) != 0;
}
bool isFiniteValue(double dValue)
{
#if defined(_MSC_VER)
return _finite(dValue) != 0;
#else
return std::isfinite(dValue);
#endif
}
StaticPropertyKind propertyKindForLegacyKey(const QString& sPropertyKey)
{
// 旧工程持久化键必须保持封闭;正式属性目录中的 pressure 不能进入此路径。
if(sPropertyKey == QLatin1String("k")) {
return StaticProperty_Permeability;
}
if(sPropertyKey == QLatin1String("phi")) {
return StaticProperty_Porosity;
}
if(sPropertyKey == QLatin1String("h")) {
return StaticProperty_Thickness;
}
return StaticProperty_Invalid;
}
nmPebiStaticPropertyCalculationStatus cancelled(QString* pErrorMessage)
{
// 主动停止不是插值错误,交由任务状态单独报告。
if(pErrorMessage != NULL) {
pErrorMessage->clear();
}
return NM_PebiStaticPropertyCalculation_Cancelled;
}
nmPebiStaticPropertyCalculationStatus failed(
QString* pErrorMessage,
const QString& sErrorMessage)
{
if(pErrorMessage != NULL) {
*pErrorMessage = sErrorMessage;
}
return NM_PebiStaticPropertyCalculation_Failed;
}
QVector<double>* valuesForKind(
nmPebiStaticPropertyValues& oValues,
StaticPropertyKind eKind)
{
switch(eKind) {
case StaticProperty_Permeability:
return &oValues.m_vecPermeabilityMilliDarcy;
case StaticProperty_Porosity:
return &oValues.m_vecPorosity;
case StaticProperty_Thickness:
return &oValues.m_vecThicknessMeters;
default:
return NULL;
}
}
}
nmPebiStaticPropertyCalculationInput::
nmPebiStaticPropertyCalculationInput()
: m_bCaptured(false),
m_dReferencePermeabilityMilliDarcy(0.0),
m_dReferencePorosity(0.0),
m_dReferenceThicknessMeters(0.0)
{
}
nmPebiStaticPropertyCalculationStatus
nmPebiStaticPropertyCalculator::calculate(
const QVector<QPointF>& vecOriginalCellPropertyLocations,
const nmPebiStaticPropertyCalculationInput& oInput,
nmPebiStaticPropertyValues& oOutput,
QString* pErrorMessage,
const QAtomicInt* pCancelRequested)
{
oOutput.clear();
if(pErrorMessage != NULL) {
pErrorMessage->clear();
}
if(isCancellationRequested(pCancelRequested)) {
return cancelled(pErrorMessage);
}
if(!oInput.m_bCaptured) {
return failed(pErrorMessage,
QLatin1String("Static property input was not captured."));
}
if(vecOriginalCellPropertyLocations.isEmpty()) {
return failed(pErrorMessage,
QLatin1String("No PEBI property locations are available."));
}
if(!isFiniteValue(oInput.m_dReferencePermeabilityMilliDarcy) ||
!isFiniteValue(oInput.m_dReferencePorosity) ||
!isFiniteValue(oInput.m_dReferenceThicknessMeters)) {
return failed(pErrorMessage,
QLatin1String("A reference static property value is invalid."));
}
for(int nIndex = 0;
nIndex < vecOriginalCellPropertyLocations.size(); ++nIndex) {
if((nIndex % 256) == 0 &&
isCancellationRequested(pCancelRequested)) {
return cancelled(pErrorMessage);
}
const QPointF& oLocation =
vecOriginalCellPropertyLocations[nIndex];
if(!isFiniteValue(oLocation.x()) || !isFiniteValue(oLocation.y())) {
return failed(pErrorMessage,
QString::fromLatin1(
"PEBI property location %1 contains an invalid coordinate.")
.arg(nIndex));
}
}
// 先完整检查启用键及唯一性,避免完成一次昂贵 Kriging 后才发现重复数据组。
bool aPropertyEnabled[3] = { false, false, false };
for(int nDataSetIndex = 0;
nDataSetIndex < oInput.m_vecDataSets.size(); ++nDataSetIndex) {
if(isCancellationRequested(pCancelRequested)) {
return cancelled(pErrorMessage);
}
const nmPropertyInterpolationDataSet& oDataSet =
oInput.m_vecDataSets[nDataSetIndex];
if(!oDataSet.useForCalculation) {
continue;
}
const StaticPropertyKind eKind =
propertyKindForLegacyKey(oDataSet.property);
if(eKind == StaticProperty_Invalid) {
return failed(pErrorMessage,
QString::fromLatin1(
"Dataset '%1' has an unknown static property key.")
.arg(oDataSet.name));
}
if(aPropertyEnabled[static_cast<int>(eKind)]) {
return failed(pErrorMessage,
QString::fromLatin1(
"More than one dataset is enabled for property %1.")
.arg(oDataSet.property));
}
aPropertyEnabled[static_cast<int>(eKind)] = true;
}
const int nCellCount = vecOriginalCellPropertyLocations.size();
nmPebiStaticPropertyValues oCandidate;
oCandidate.m_vecPermeabilityMilliDarcy.fill(
oInput.m_dReferencePermeabilityMilliDarcy, nCellCount);
oCandidate.m_vecPorosity.fill(
oInput.m_dReferencePorosity, nCellCount);
oCandidate.m_vecThicknessMeters.fill(
oInput.m_dReferenceThicknessMeters, nCellCount);
for(int nDataSetIndex = 0;
nDataSetIndex < oInput.m_vecDataSets.size(); ++nDataSetIndex) {
if(isCancellationRequested(pCancelRequested)) {
return cancelled(pErrorMessage);
}
const nmPropertyInterpolationDataSet& oDataSet =
oInput.m_vecDataSets[nDataSetIndex];
if(!oDataSet.useForCalculation) {
continue;
}
QVector<QPointF> vecMeasurementPoints;
QVector<double> vecMeasurementValues;
vecMeasurementPoints.reserve(oDataSet.points.size());
vecMeasurementValues.reserve(oDataSet.points.size());
for(int nPointIndex = 0;
nPointIndex < oDataSet.points.size(); ++nPointIndex) {
const nmPropertyInterpolationPointData& oPoint =
oDataSet.points[nPointIndex];
vecMeasurementPoints.append(QPointF(oPoint.x, oPoint.y));
vecMeasurementValues.append(oPoint.value);
}
QVector<double> vecInterpolationValues;
QString sCalculationError;
if(!nmCalculationUtils::calculateKriging(
vecOriginalCellPropertyLocations,
vecMeasurementPoints,
vecMeasurementValues,
oDataSet.nugget,
oDataSet.sill,
oDataSet.range,
oDataSet.model,
oInput.m_sLicensePath,
vecInterpolationValues,
&sCalculationError,
pCancelRequested)) {
// calculateKriging() 只有 bool 返回值,取消必须读取原子状态判定。
if(isCancellationRequested(pCancelRequested)) {
return cancelled(pErrorMessage);
}
return failed(pErrorMessage,
QString::fromLatin1("Dataset '%1': %2")
.arg(oDataSet.name)
.arg(sCalculationError));
}
const StaticPropertyKind eKind =
propertyKindForLegacyKey(oDataSet.property);
QVector<double>* pTargetValues = valuesForKind(oCandidate, eKind);
if(pTargetValues == NULL ||
vecInterpolationValues.size() != nCellCount) {
return failed(pErrorMessage,
QString::fromLatin1(
"Dataset '%1' returned an invalid value count.")
.arg(oDataSet.name));
}
*pTargetValues = vecInterpolationValues;
}
if(!oCandidate.isValidForCellCount(nCellCount)) {
return failed(pErrorMessage,
QLatin1String("The calculated static properties are invalid."));
}
if(isCancellationRequested(pCancelRequested)) {
return cancelled(pErrorMessage);
}
// 三个 qSwap 均为容器所有权交换;到达此处后不再执行可能失败的操作。
qSwap(oOutput.m_vecPermeabilityMilliDarcy,
oCandidate.m_vecPermeabilityMilliDarcy);
qSwap(oOutput.m_vecPorosity, oCandidate.m_vecPorosity);
qSwap(oOutput.m_vecThicknessMeters, oCandidate.m_vecThicknessMeters);
return NM_PebiStaticPropertyCalculation_Succeeded;
}

@ -0,0 +1,88 @@
#include "nmPebiResultPropertyCatalog.h"
namespace {
nmPebiResultPropertyDescriptor makeDescriptor(
const char* pPropertyId,
const char* pDisplayNameSourceText,
const char* pBaseUnit,
NM_PEBI_RESULT_PROPERTY_TIME_TYPE eTimeType,
NM_PEBI_RESULT_PROPERTY_SOURCE eSource,
NM_PEBI_RESULT_PROPERTY_RANGE_POLICY eRangePolicy)
{
nmPebiResultPropertyDescriptor oDescriptor;
oDescriptor.m_sPropertyId = QString::fromLatin1(pPropertyId);
oDescriptor.m_sDisplayNameSourceText =
QString::fromLatin1(pDisplayNameSourceText);
oDescriptor.m_sBaseUnit = QString::fromLatin1(pBaseUnit);
oDescriptor.m_sDefaultDisplayUnit = oDescriptor.m_sBaseUnit;
oDescriptor.m_eTimeType = eTimeType;
oDescriptor.m_eSource = eSource;
oDescriptor.m_eAssociation =
NM_PebiResultPropertyAssociation_CellData;
oDescriptor.m_eRangePolicy = eRangePolicy;
return oDescriptor;
}
}
nmPebiResultPropertyDescriptor::nmPebiResultPropertyDescriptor()
: m_eTimeType(NM_PebiResultProperty_Static),
m_eSource(NM_PebiResultPropertySource_StaticCalculation),
m_eAssociation(NM_PebiResultPropertyAssociation_CellData),
m_eRangePolicy(NM_PebiResultPropertyRange_CurrentArray)
{
}
QList<nmPebiResultPropertyDescriptor>
nmPebiResultPropertyCatalog::descriptors()
{
// 每次返回值副本,避免 VS2010 下函数静态容器的并发初始化问题。
QList<nmPebiResultPropertyDescriptor> listDescriptors;
listDescriptors.append(makeDescriptor(
"permeability", "Permeability", "mD",
NM_PebiResultProperty_Static,
NM_PebiResultPropertySource_StaticCalculation,
NM_PebiResultPropertyRange_CurrentArray));
listDescriptors.append(makeDescriptor(
"porosity", "Porosity", "",
NM_PebiResultProperty_Static,
NM_PebiResultPropertySource_StaticCalculation,
NM_PebiResultPropertyRange_CurrentArray));
listDescriptors.append(makeDescriptor(
"thickness", "Thickness", "m",
NM_PebiResultProperty_Static,
NM_PebiResultPropertySource_StaticCalculation,
NM_PebiResultPropertyRange_CurrentArray));
listDescriptors.append(makeDescriptor(
"pressure", "Pressure", "MPa",
NM_PebiResultProperty_Dynamic,
NM_PebiResultPropertySource_SolverOutput,
NM_PebiResultPropertyRange_AllFrames));
return listDescriptors;
}
bool nmPebiResultPropertyCatalog::findDescriptor(
const QString& sPropertyId,
nmPebiResultPropertyDescriptor& oDescriptor)
{
const QList<nmPebiResultPropertyDescriptor> listDescriptors =
descriptors();
for(int nIndex = 0; nIndex < listDescriptors.size(); ++nIndex) {
if(listDescriptors[nIndex].m_sPropertyId == sPropertyId) {
oDescriptor = listDescriptors[nIndex];
return true;
}
}
// 查询失败时清除调用方旧值,避免把上一次命中的描述误当成当前属性。
oDescriptor = nmPebiResultPropertyDescriptor();
return false;
}
bool nmPebiResultPropertyCatalog::isKnownProperty(
const QString& sPropertyId)
{
nmPebiResultPropertyDescriptor oDescriptor;
return findDescriptor(sPropertyId, oDescriptor);
}

@ -0,0 +1,53 @@
#include "nmPebiStaticPropertyValues.h"
#include <float.h>
#include <math.h>
#if !defined(_MSC_VER)
#include <cmath>
#endif
namespace {
bool isFinitePropertyValue(double dValue)
{
#if defined(_MSC_VER)
return _finite(dValue) != 0;
#else
return std::isfinite(dValue);
#endif
}
bool containsOnlyFiniteValues(const QVector<double>& vecValues)
{
for(int nIndex = 0; nIndex < vecValues.size(); ++nIndex) {
if(!isFinitePropertyValue(vecValues[nIndex])) {
return false;
}
}
return true;
}
}
void nmPebiStaticPropertyValues::clear()
{
m_vecPermeabilityMilliDarcy.clear();
m_vecPorosity.clear();
m_vecThicknessMeters.clear();
}
bool nmPebiStaticPropertyValues::isValidForCellCount(
int nOriginalCellCount) const
{
if(nOriginalCellCount <= 0 ||
m_vecPermeabilityMilliDarcy.size() != nOriginalCellCount ||
m_vecPorosity.size() != nOriginalCellCount ||
m_vecThicknessMeters.size() != nOriginalCellCount) {
return false;
}
return containsOnlyFiniteValues(m_vecPermeabilityMilliDarcy) &&
containsOnlyFiniteValues(m_vecPorosity) &&
containsOnlyFiniteValues(m_vecThicknessMeters);
}
Loading…
Cancel
Save