#include "nmCalculationAutoFitGA.h" //#include "nmCalculationExeSolverTask.h" #include "nmCalculationDllPebiSolverTask.h" #include "nmDataAnalyzeManager.h" #include "nmDataWellBase.h" #include "nmDataVerticalWell.h" #include "nmDataVerticalFracturedWell.h" #include "nmDataHorizontalFracturedWell.h" #include "nmDataReservoir.h" #include "nmDataAutomaticFitting.h" #include #include #include #include #include #include #include #ifdef Q_OS_WIN #include #include #define DEBUG_OUT(msg) OutputDebugStringA(QString("[AutoFit] %1\n").arg(msg).toLocal8Bit().data()) #endif // 常量定义 const double nmCalculationAutoFitGA::MIN_FITNESS_IMPROVEMENT = 1e-8; const double nmCalculationAutoFitGA::MUTATION_STRENGTH = 0.1; const int nmCalculationAutoFitGA::CONVERGENCE_CHECK_INTERVAL = 10; const int nmCalculationAutoFitGA::MAX_STAGNATION_GENERATIONS = 20; // 无穷大和NaN检查 static inline bool isFiniteNumber(double value) { #ifdef Q_OS_WIN return _finite(value) != 0 && !_isnan(value); #else return std::isfinite(value); #endif } // sleep函数 static inline void msleep(int ms) { #ifdef Q_OS_WIN Sleep(ms); #endif } // 构造函数 nmCalculationAutoFitGA::nmCalculationAutoFitGA(QObject* parent) : QObject(parent) , m_isRunning(false) , m_shouldStop(false) , m_isPaused(false) , m_currentGeneration(0) , m_bestFitness(1e10) , m_worstFitness(-1e10) , m_averageFitness(1e10) , m_previousBestFitness(1e10) , m_populationSize(40) , m_maxGenerations(100) , m_targetError(0.001) , m_crossoverRate(0.8) , m_mutationRate(0.1) , m_elitismRate(0.1) , m_tournamentSize(3) , m_useUniformCrossover(true) , m_totalEvaluations(0) , m_successfulEvaluations(0) , m_evaluationInProgress(0) , m_consecutiveFailures(0) , m_progressTimer(0) , m_userInitialFitness(1e10) , m_consecutiveFailedGenerations(0) , m_maxConsecutiveFailures(3) , m_hasValidUserSolution(false) , m_improvementThreshold(0.05) , m_diversityThreshold(0.05) , m_convergenceVarianceThreshold(1e-8) , m_trueConvergenceWindow(15) , m_localOptimumWindow(8) , m_nearTargetFactor(2.0) , m_farTargetFactor(10.0) , m_targetWellName("") { DEBUG_OUT(QString("GA Constructor: this=0x%1").arg((quintptr)this, 0, 16)); // 初始化随机数种子 qsrand(QTime::currentTime().msec()); // 初始化临时目录用于DLL求解器 initializeTemporaryDirectory(); // 初始化最优个体 m_bestIndividual.fitness = 1e10; m_bestIndividual.isEvaluated = false; // 创建进度更新定时器 m_progressTimer = new QTimer(this); connect(m_progressTimer, SIGNAL(timeout()), this, SLOT(updateProgress())); DEBUG_OUT("AutoFit GA calculator initialized (data-driven mode)"); DEBUG_OUT("GA Constructor completed"); } // 析构函数 nmCalculationAutoFitGA::~nmCalculationAutoFitGA() { DEBUG_OUT(QString("GA Destructor: this=0x%1").arg((quintptr)this, 0, 16)); // 首先停止算法 if(m_isRunning) { m_shouldStop = true; // 立即设置停止标志 // 等待当前操作完成,增加超时时间 int waitCount = 0; while(m_isRunning && waitCount < 100) { QApplication::processEvents(QEventLoop::ExcludeUserInputEvents, 50); msleep(50); waitCount++; } // 如果仍在运行则强制停止 if(m_isRunning) { DEBUG_OUT("Force stopping GA - timeout reached"); m_isRunning = false; } } // 清理临时目录 cleanupTemporaryDirectory(); // 断开所有信号连接,防止回调已销毁的对象 disconnect(this, nullptr, nullptr, nullptr); DEBUG_OUT("GA Destructor completed"); } // ==================== 目录处理方法 ==================== void nmCalculationAutoFitGA::initializeTemporaryDirectory() { QString timestamp = QDateTime::currentDateTime().toString("yyyyMMdd_hhmmss_zzz"); QString processId = QString::number(QCoreApplication::applicationPid()); m_tempDirectory = QApplication::applicationDirPath() + "/autofit_temp_" + processId + "_" + timestamp; // 确保目录不存在 int counter = 0; QString originalPath = m_tempDirectory; while(QDir(m_tempDirectory).exists() && counter < 100) { m_tempDirectory = originalPath + "_" + QString::number(counter); counter++; } if(QDir().mkpath(m_tempDirectory)) { DEBUG_OUT(QString("Initialized temp directory: %1").arg(m_tempDirectory)); } else { DEBUG_OUT(QString("Warning: Failed to create temp directory: %1").arg(m_tempDirectory)); m_tempDirectory = QApplication::applicationDirPath(); } } void nmCalculationAutoFitGA::cleanupTemporaryDirectory() { if(QDir(m_tempDirectory).exists()) { if(removeDirectoryRecursively(m_tempDirectory)) { DEBUG_OUT("Temp directory cleaned up successfully"); } else { DEBUG_OUT("Warning: Failed to clean up temp directory completely"); } } } bool nmCalculationAutoFitGA::removeDirectoryRecursively(const QString& path) { QDir dir(path); if(!dir.exists()) { return true; } // 递归删除子目录和文件 QFileInfoList entries = dir.entryInfoList(QDir::NoDotAndDotDot | QDir::AllEntries | QDir::Hidden); bool allRemoved = true; for(int i = 0; i < entries.size(); ++i) { const QFileInfo& entry = entries[i]; if(entry.isDir()) { if(!removeDirectoryRecursively(entry.absoluteFilePath())) { allRemoved = false; } } else { QFile file(entry.absoluteFilePath()); // 处理只读文件 if(!file.permissions().testFlag(QFile::WriteUser)) { file.setPermissions(file.permissions() | QFile::WriteUser); } if(!file.remove()) { DEBUG_OUT(QString("Failed to remove file: %1").arg(entry.absoluteFilePath())); allRemoved = false; } } } // 删除目录本身 if(allRemoved) { return dir.rmdir(path); } return false; } // ==================== 公共接口方法 ==================== void nmCalculationAutoFitGA::setTargetLogLogData(const QVector>& targetData) { m_targetLogLogData = targetData; DEBUG_OUT(QString("Target LogLog data set: %1 arrays").arg(targetData.size())); if(targetData.size() >= 3) { DEBUG_OUT(QString("LogLog data points: X=%1, Y1=%2, Y2=%3") .arg(targetData[0].size()) .arg(targetData[1].size()) .arg(targetData[2].size())); } } bool nmCalculationAutoFitGA::startAutoFitting() { if(m_isRunning) { m_lastError = "GA fitting is already running"; return false; } DEBUG_OUT("=== GA AUTO FITTING START ==="); // 发送初始化日志 emit logMessageGenerated(tr("=== GA Automatic Fitting Started ===")); emit logMessageGenerated(tr("Algorithm: Genetic Algorithm")); try { // 从数据管理器加载所有配置 if(!loadAllConfigFromDataManager()) { emit logMessageGenerated(tr("ERROR: Failed to load configuration from data manager")); return false; } int enabledParams = getEnabledParameterCount(); emit logMessageGenerated(tr("Enabled parameters count: %1").arg(enabledParams)); if(enabledParams == 0) { m_lastError = "No parameters enabled for optimization"; emit logMessageGenerated(tr("ERROR: No parameters enabled for optimization")); return false; } if(m_targetLogLogData.isEmpty() || m_targetLogLogData.size() < 3) { m_lastError = "Target LogLog data is empty or insufficient"; emit logMessageGenerated(tr("ERROR: Target LogLog data is empty or insufficient")); return false; } // 检查数据一致性 if(m_targetLogLogData[0].size() != m_targetLogLogData[1].size() || m_targetLogLogData[0].size() != m_targetLogLogData[2].size()) { m_lastError = "Target LogLog data arrays have inconsistent sizes"; emit logMessageGenerated(tr("ERROR: Target LogLog data arrays have inconsistent sizes")); return false; } emit logMessageGenerated(tr("Target data validation passed (%1 data points)").arg(m_targetLogLogData[0].size())); // 使用保存的初始值进行精英保护 QVector savedInitialValues = m_initialValues; // 重置状态 resetOptimizer(); m_isRunning = true; m_shouldStop = false; m_isPaused = false; m_currentGeneration = 0; m_consecutiveFailures = 0; m_consecutiveFailedGenerations = 0; // 精英保护:评估用户初始解 if(!savedInitialValues.isEmpty()) { m_userInitialSolution = savedInitialValues; emit logMessageGenerated(tr("=== Evaluating Initial Solution (Elite Protection) ===")); // 输出初始参数值 QString paramStr = tr("Initial parameters: "); for(int i = 0; i < m_userInitialSolution.size(); ++i) { paramStr += QString("[%1]=%2 ").arg(i).arg(m_userInitialSolution[i], 0, 'f', 6); } emit logMessageGenerated(paramStr); try { emit logMessageGenerated(tr("Starting initial solution evaluation...")); DEBUG_OUT(QString("Before evaluateFitness: m_userInitialSolution size = %1").arg(m_userInitialSolution.size())); if(m_userInitialSolution.isEmpty()) { emit logMessageGenerated(tr("ERROR: m_userInitialSolution is empty!")); return false; } for(int i = 0; i < m_userInitialSolution.size(); ++i) { emit logMessageGenerated(tr("Initial param[%1] = %2").arg(i).arg(m_userInitialSolution[i], 0, 'f', 6)); } m_userInitialFitness = evaluateGenes(m_userInitialSolution); emit logMessageGenerated(tr("evaluateGenes returned: %1").arg(m_userInitialFitness, 0, 'e', 10)); if(m_userInitialFitness < 1e9) { emit logMessageGenerated(tr("Taking SUCCESS branch (fitness < 1e9)")); m_hasValidUserSolution = true; m_bestFitness = m_userInitialFitness; m_bestIndividual.genes = m_userInitialSolution; m_bestIndividual.fitness = m_userInitialFitness; m_bestIndividual.isEvaluated = true; emit logMessageGenerated(tr("Initial solution evaluation successful")); emit logMessageGenerated(tr("Initial fitness (error): %1").arg(m_userInitialFitness, 0, 'e', 4)); emit logMessageGenerated(tr("Elite protection activated - initial solution will be preserved if no significant improvement found")); } else { emit logMessageGenerated(tr("Taking FAILURE branch (fitness >= 1e9)")); m_hasValidUserSolution = false; emit logMessageGenerated(tr("Initial solution evaluation failed - starting with random initialization")); } } catch(...) { m_hasValidUserSolution = false; emit logMessageGenerated(tr("Exception during initial solution evaluation")); } // 恢复初始值供种群初始化使用 m_initialValues = savedInitialValues; } // 初始化种群 initializePopulation(); emit logMessageGenerated(tr("Population initialized: %1 individuals, %2 dimensions").arg(m_populationSize).arg(getEnabledParameterCount())); // GA主循环 emit logMessageGenerated(tr("=== Starting GA Main Loop ===")); for(m_currentGeneration = 0; m_currentGeneration < m_maxGenerations && !m_shouldStop; ++m_currentGeneration) { // 每次迭代都输出标题,或者只在重要迭代输出详细信息 bool shouldOutputDetail = (m_currentGeneration % qMax(1, m_maxGenerations / 10) == 0) || (m_currentGeneration < 5) || (m_currentGeneration >= m_maxGenerations - 2); // 每次迭代都输出标题 emit logMessageGenerated(tr("--- Generation %1/%2 ---").arg(m_currentGeneration + 1).arg(m_maxGenerations)); // 只在特定迭代输出详细统计信息 if(shouldOutputDetail) { emit logMessageGenerated(tr("Current best error: %1").arg(m_bestFitness, 0, 'e', 4)); emit logMessageGenerated(tr("Total evaluations: %1 (successful: %2, failures: %3)") .arg(m_totalEvaluations).arg(m_successfulEvaluations).arg(m_totalEvaluations - m_successfulEvaluations)); } // 检查暂停状态 while(m_isPaused && !m_shouldStop) { QApplication::processEvents(); //msleep(100); } if(m_shouldStop) { emit logMessageGenerated(tr("Optimization stopped by user request")); break; } // 1. 评估种群 evaluatePopulation(); if(m_shouldStop) break; // 2. 更新统计信息 updatePopulationStatistics(); // 3. 发射进度信号 emit progressUpdated(m_currentGeneration, m_bestFitness); if(shouldOutputDetail) { emit logMessageGenerated(tr("Generation %1 completed: best = %2, avg = %3, worst = %4") .arg(m_currentGeneration + 1).arg(m_bestFitness, 0, 'e', 4) .arg(m_averageFitness, 0, 'e', 4).arg(m_worstFitness, 0, 'e', 4)); } // 记录收敛历史 m_convergenceHistory.append(m_bestFitness); // 更新收敛指标 updateConvergenceMetrics(); // 智能收敛判断 StopReasonGA stopReason = analyzeOptimizationStatus(); if(stopReason == GA_TARGET_ACHIEVED) { emit logMessageGenerated(tr("=== TARGET ACHIEVED ===")); emit logMessageGenerated(tr("Target error achieved! Current error: %1 < Target: %2") .arg(m_bestFitness, 0, 'e', 4).arg(m_targetError, 0, 'e', 4)); emit logMessageGenerated(tr("Optimization completed successfully after %1 generations") .arg(m_currentGeneration + 1)); break; } else if(stopReason == GA_TRUE_CONVERGENCE) { emit logMessageGenerated(tr("=== TRUE CONVERGENCE DETECTED ===")); emit logMessageGenerated(tr("Algorithm has converged to a stable solution")); emit logMessageGenerated(tr("Final error: %1 after %2 generations") .arg(m_bestFitness, 0, 'e', 4).arg(m_currentGeneration + 1)); emit logMessageGenerated(tr("Solution quality: %1 (1.0 = target achieved)") .arg(m_targetError / qMax(1e-10, m_bestFitness), 0, 'f', 3)); break; } else if(stopReason == GA_LOCAL_OPTIMUM) { emit logMessageGenerated(tr("=== LOCAL OPTIMUM DETECTED ===")); emit logMessageGenerated(tr("Algorithm appears to be trapped in local optimum")); emit logMessageGenerated(tr("Current error: %1 after %2 generations") .arg(m_bestFitness, 0, 'e', 4).arg(m_currentGeneration + 1)); emit logMessageGenerated(tr("Suggestion: Try restarting with different parameters or larger search space")); break; } else if(stopReason == GA_CONSECUTIVE_FAILURES) { emit logMessageGenerated(tr("=== CONSECUTIVE FAILURES ===")); emit logMessageGenerated(tr("Too many consecutive failed generations (%1/%2)") .arg(m_consecutiveFailedGenerations).arg(m_maxConsecutiveFailures)); break; } else if(stopReason == GA_CONTINUE_OPTIMIZATION) { // 继续优化,每10次迭代输出一次状态 if(m_currentGeneration > 0 && m_currentGeneration % 10 == 0) { double diversity = calculatePopulationDiversity(); emit logMessageGenerated(tr(" Optimization status: diversity=%1") .arg(diversity, 0, 'f', 4)); } } // 4. 创建新一代种群 if(m_currentGeneration < m_maxGenerations - 1) { QVector newPopulation; newPopulation.reserve(m_populationSize); // 精英保留 applyElitism(newPopulation); // 生成新个体直到填满种群 while(newPopulation.size() < m_populationSize && !m_shouldStop) { // 选择父代 int parent1Index = tournamentSelection(); int parent2Index = tournamentSelection(); // 确保父代不同 while(parent1Index == parent2Index && m_population.size() > 1) { parent2Index = tournamentSelection(); } GAIndividual offspring1, offspring2; // 交叉 if(random01() < m_crossoverRate) { crossover(m_population[parent1Index], m_population[parent2Index], offspring1, offspring2); } else { offspring1 = m_population[parent1Index]; offspring2 = m_population[parent2Index]; } // 变异 if(random01() < m_mutationRate) { mutate(offspring1); } if(random01() < m_mutationRate) { mutate(offspring2); } // 边界约束 clampToLimits(offspring1.genes); clampToLimits(offspring2.genes); // 添加到新种群 if(newPopulation.size() < m_populationSize) { newPopulation.append(offspring1); } if(newPopulation.size() < m_populationSize) { newPopulation.append(offspring2); } } // 替换种群 m_population = newPopulation; // 自适应参数调整 adaptiveParameterUpdate(m_currentGeneration); } // 输出迭代结束标记 if(shouldOutputDetail) { emit logMessageGenerated(tr("Generation %1 completed - Current best: %2") .arg(m_currentGeneration + 1).arg(m_bestFitness, 0, 'e', 4)); } // 强制处理事件,保持界面响应 QApplication::processEvents(); // 在世代间稍作停顿,减少系统负载 if(m_currentGeneration % 3 == 2) { //msleep(200); } } // 最终结果验证和保护 validateAndProtectFinalResult(); } catch(const std::exception& e) { m_lastError = QString("Critical exception in GA main loop: %1").arg(e.what()); emit logMessageGenerated(tr("CRITICAL ERROR: %1").arg(e.what())); cleanupTemporaryDirectory(); m_isRunning = false; if(m_progressTimer) m_progressTimer->stop(); emit fittingFinished(false, m_lastError); return false; } catch(...) { m_lastError = "Unknown critical exception in GA main loop"; emit logMessageGenerated(tr("CRITICAL ERROR: Unknown exception in GA main loop")); cleanupTemporaryDirectory(); m_isRunning = false; if(m_progressTimer) m_progressTimer->stop(); emit fittingFinished(false, m_lastError); return false; } m_isRunning = false; if(m_progressTimer) m_progressTimer->stop(); // 应用最终参数 if(!m_bestIndividual.genes.isEmpty()) { try { emit logMessageGenerated(tr("Applying optimized parameters to model...")); applyParametersToDataManager(m_bestIndividual.genes); saveOptimizationResult(); // 输出最终优化结果 emit logMessageGenerated(tr("=== Optimization Results ===")); emit logMessageGenerated(tr("Final error: %1").arg(m_bestFitness, 0, 'e', 4)); emit logMessageGenerated(tr("Total generations: %1").arg(m_currentGeneration + 1)); emit logMessageGenerated(tr("Total evaluations: %1 (successful: %2)") .arg(m_totalEvaluations).arg(m_successfulEvaluations)); // 输出最优参数值 QString finalParams = tr("Optimized parameters: "); for(int i = 0; i < m_bestIndividual.genes.size(); ++i) { finalParams += QString("[%1]=%2 ").arg(i).arg(m_bestIndividual.genes[i], 0, 'f', 6); } emit logMessageGenerated(finalParams); emit logMessageGenerated(tr("Parameters applied successfully to data manager")); } catch(const std::exception& e) { emit logMessageGenerated(tr("ERROR: Failed to apply final parameters: %1").arg(e.what())); m_lastError = QString("Failed to apply final parameters: %1").arg(e.what()); } catch(...) { emit logMessageGenerated(tr("ERROR: Unknown error applying final parameters")); m_lastError = "Failed to apply final parameters due to unknown error"; } } // 判断系统确定最终结果 bool success; QString message; StopReasonGA finalReason = analyzeOptimizationStatus(); if(finalReason == GA_TARGET_ACHIEVED) { success = true; //message = QString("GA optimization completed successfully. Target achieved. Best error: %1, Generations: %2") // .arg(m_bestFitness, 0, 'e', 4).arg(m_currentGeneration + 1); emit logMessageGenerated(tr("=== GA OPTIMIZATION SUCCESSFUL ===")); } else if(finalReason == GA_TRUE_CONVERGENCE) { success = true; //message = QString("GA optimization converged to stable solution. Best error: %1, Generations: %2") // .arg(m_bestFitness, 0, 'e', 4).arg(m_currentGeneration + 1); emit logMessageGenerated(tr("=== GA OPTIMIZATION CONVERGED ===")); } else if(finalReason == GA_LOCAL_OPTIMUM) { success = true; //message = QString("GA optimization trapped in local optimum. Best error: %1, Generations: %2") // .arg(m_bestFitness, 0, 'e', 4).arg(m_currentGeneration + 1); emit logMessageGenerated(tr("=== GA OPTIMIZATION - LOCAL OPTIMUM ===")); } else if(finalReason == GA_MAX_ITERATIONS) { success = true; //message = QString("GA optimization completed. Max generations reached. Best error: %1, Generations: %2") // .arg(m_bestFitness, 0, 'e', 4).arg(m_currentGeneration + 1); emit logMessageGenerated(tr("=== GA OPTIMIZATION - MAX GENERATIONS ===")); } else if(finalReason == GA_USER_STOPPED) { success = true; //message = QString("GA optimization stopped by user. Best error: %1, Generations: %2") // .arg(m_bestFitness, 0, 'e', 4).arg(m_currentGeneration + 1); emit logMessageGenerated(tr("=== GA OPTIMIZATION STOPPED BY USER ===")); } else if(finalReason == GA_CONSECUTIVE_FAILURES) { success = false; //message = QString("GA optimization failed due to consecutive failures. Best error: %1, Generations: %2") // .arg(m_bestFitness, 0, 'e', 4).arg(m_currentGeneration + 1); emit logMessageGenerated(tr("=== GA OPTIMIZATION FAILED ===")); } else { success = false; //message = QString("GA optimization ended unexpectedly. Best error: %1, Generations: %2") // .arg(m_bestFitness, 0, 'e', 4).arg(m_currentGeneration + 1); emit logMessageGenerated(tr("=== GA OPTIMIZATION - UNKNOWN END ===")); } emit logMessageGenerated(tr("Result: %1").arg(success ? "SUCCESS" : "FAILED")); emit fittingFinished(success, message); cleanupTemporaryDirectory(); return success; } void nmCalculationAutoFitGA::stopFitting() { if(m_isRunning) { DEBUG_OUT("Stop request received, setting stop flag..."); // 添加停止日志 emit logMessageGenerated(tr("=== User Stop Request Received ===")); emit logMessageGenerated(tr("Gracefully stopping GA optimization...")); m_shouldStop = true; if(m_progressTimer) { m_progressTimer->stop(); } // 等待当前评估完成,缩短超时时间 int waitCount = 0; while(m_evaluationInProgress > 0 && waitCount < 30) { // 减少等待时间 QApplication::processEvents(QEventLoop::ExcludeUserInputEvents, 50); msleep(50); waitCount++; } // 超时时强制重置 if(m_evaluationInProgress > 0) { DEBUG_OUT("Force resetting evaluation counter"); emit logMessageGenerated(tr("Force stopping current evaluation...")); m_evaluationInProgress = 0; } // 确保运行标志被清除 m_isRunning = false; emit logMessageGenerated(tr("GA optimization stop request processed")); DEBUG_OUT("Stop request processed"); } else { DEBUG_OUT("Stop request received but GA is not running"); emit logMessageGenerated(tr("Stop request received but optimization is not running")); } cleanupTemporaryDirectory(); } bool nmCalculationAutoFitGA::isRunning() const { return m_isRunning; } int nmCalculationAutoFitGA::getCurrentGeneration() const { return m_currentGeneration; } QVector nmCalculationAutoFitGA::getBestSolution() const { return m_bestIndividual.genes; } double nmCalculationAutoFitGA::getBestFitness() const { return m_bestFitness; } QString nmCalculationAutoFitGA::getLastError() const { return m_lastError; } void nmCalculationAutoFitGA::resetOptimizer() { m_population.clear(); m_eliteIndividuals.clear(); m_bestIndividual = GAIndividual(); m_bestFitness = 1e10; m_worstFitness = -1e10; m_averageFitness = 1e10; m_previousBestFitness = 1e10; m_currentGeneration = 0; m_totalEvaluations = 0; m_successfulEvaluations = 0; m_convergenceHistory.clear(); m_lastError.clear(); m_initialValues.clear(); m_userInitialSolution.clear(); m_userInitialFitness = 1e10; m_hasValidUserSolution = false; m_diversityHistory.clear(); DEBUG_OUT("GA optimizer reset"); } void nmCalculationAutoFitGA::setGATargetWellName(const QString& wellName) { m_targetWellName = wellName; } void nmCalculationAutoFitGA::updateProgress() { // 这个槽函数在定时器触发时被调用,可以用来更新界面或执行周期性任务 if(m_isRunning) { QApplication::processEvents(); } } // ==================== 数据加载方法 ==================== bool nmCalculationAutoFitGA::loadAllConfigFromDataManager() { try { loadOptimizationConfig(); loadParameterBounds(); extractUserInitialValues(); return true; } catch(...) { m_lastError = "Failed to load configuration from data manager"; return false; } } void nmCalculationAutoFitGA::loadOptimizationConfig() { nmDataAnalyzeManager* dataManager = nmDataAnalyzeManager::getCurrentInstance(); nmDataAutomaticFitting fittingData = dataManager->getAutomaticFittingDataCopy(); // 加载基本配置 m_maxGenerations = fittingData.getIterationCount().getValue().toInt(); m_targetError = fittingData.getErrorTolerance().getValue().toDouble(); // 设置GA默认参数 m_populationSize = 40; m_crossoverRate = 0.8; m_mutationRate = 0.1; m_elitismRate = 0.15; m_tournamentSize = 3; m_useUniformCrossover = true; DEBUG_OUT(QString("Loaded GA optimization config: generations=%1, error=%2, population=%3") .arg(m_maxGenerations).arg(m_targetError).arg(m_populationSize)); } void nmCalculationAutoFitGA::loadParameterBounds() { nmDataAnalyzeManager* dataManager = nmDataAnalyzeManager::getCurrentInstance(); nmDataAutomaticFitting fittingData = dataManager->getAutomaticFittingDataCopy(); // 获取参数选择状态 m_parameterSelected.resize(11); m_parameterSelected[0] = fittingData.getPermeabilitySelected(); m_parameterSelected[1] = fittingData.getSkinSelected(); m_parameterSelected[2] = fittingData.getWellboreStorageSelected(); m_parameterSelected[3] = fittingData.getPorositySelected(); m_parameterSelected[4] = fittingData.getInitialPressureSelected(); m_parameterSelected[5] = fittingData.getThicknessSelected(); m_parameterSelected[6] = fittingData.getCtSelected(); m_parameterSelected[7] = fittingData.getCfSelected(); m_parameterSelected[8] = fittingData.getSoiSelected(); m_parameterSelected[9] = fittingData.getSwiSelected(); m_parameterSelected[10] = fittingData.getSgiSelected(); // 获取参数边界 m_parameterLower.resize(11); m_parameterUpper.resize(11); m_parameterLower[0] = fittingData.getPermeabilityMin().getValue().toDouble(); m_parameterUpper[0] = fittingData.getPermeabilityMax().getValue().toDouble(); m_parameterLower[1] = fittingData.getSkinMin().getValue().toDouble(); m_parameterUpper[1] = fittingData.getSkinMax().getValue().toDouble(); m_parameterLower[2] = fittingData.getWellboreStorageMin().getValue().toDouble(); m_parameterUpper[2] = fittingData.getWellboreStorageMax().getValue().toDouble(); m_parameterLower[3] = fittingData.getPorosityMin().getValue().toDouble(); m_parameterUpper[3] = fittingData.getPorosityMax().getValue().toDouble(); m_parameterLower[4] = fittingData.getInitialPressureMin().getValue().toDouble(); m_parameterUpper[4] = fittingData.getInitialPressureMax().getValue().toDouble(); m_parameterLower[5] = fittingData.getThicknessMin().getValue().toDouble(); m_parameterUpper[5] = fittingData.getThicknessMax().getValue().toDouble(); m_parameterLower[6] = fittingData.getCtMin().getValue().toDouble(); m_parameterUpper[6] = fittingData.getCtMax().getValue().toDouble(); m_parameterLower[7] = fittingData.getCfMin().getValue().toDouble(); m_parameterUpper[7] = fittingData.getCfMax().getValue().toDouble(); m_parameterLower[8] = fittingData.getSoiMin().getValue().toDouble(); m_parameterUpper[8] = fittingData.getSoiMax().getValue().toDouble(); m_parameterLower[9] = fittingData.getSwiMin().getValue().toDouble(); m_parameterUpper[9] = fittingData.getSwiMax().getValue().toDouble(); m_parameterLower[10] = fittingData.getSgiMin().getValue().toDouble(); m_parameterUpper[10] = fittingData.getSgiMax().getValue().toDouble(); // 更新启用参数索引 m_enabledParamIndices.clear(); for(int i = 0; i < m_parameterSelected.size(); ++i) { if(m_parameterSelected[i]) { m_enabledParamIndices.append(i); } } DEBUG_OUT(QString("Loaded parameter bounds: %1 enabled parameters") .arg(m_enabledParamIndices.size())); } void nmCalculationAutoFitGA::extractUserInitialValues() { nmDataAnalyzeManager* dataManager = nmDataAnalyzeManager::getCurrentInstance(); nmDataReservoir reservoirData = dataManager->getReservoirDataCopy(); //QVector wells = dataManager->getWellDataList(); nmDataWellBase* pTargetWell = dataManager->findWellByName(m_targetWellName); m_initialValues.clear(); // 按照启用参数的顺序提取初始值 for(int i = 0; i < m_enabledParamIndices.size(); ++i) { int paramIndex = m_enabledParamIndices[i]; double initialValue = 0.0; switch(paramIndex) { case 0: // 渗透率 initialValue = reservoirData.getPermeability().getValue().toDouble(); break; case 1: // 表皮系数 if(pTargetWell) { initialValue = pTargetWell->getPerforation(0)->getSkin().getValue().toDouble(); } break; case 2: // 井筒储集系数 if(pTargetWell) { initialValue = pTargetWell->getWellboreStorage().getValue().toDouble(); } break; case 3: // 孔隙度 initialValue = reservoirData.getPorosity().getValue().toDouble(); break; case 4: // 初始压力 initialValue = reservoirData.getInitialPressure().getValue().toDouble(); break; case 5: // 储层厚度 initialValue = reservoirData.getThickness().getValue().toDouble(); break; case 6: // 综合压缩系数 initialValue = reservoirData.getCt().getValue().toDouble(); break; case 7: // 岩石压缩系数 initialValue = reservoirData.getCf().getValue().toDouble(); break; case 8: // 初始含油饱和度 initialValue = reservoirData.getSoi().getValue().toDouble(); break; case 9: // 初始含水饱和度 initialValue = reservoirData.getSwi().getValue().toDouble(); break; case 10: // 初始含气饱和度 initialValue = reservoirData.getSgi().getValue().toDouble(); break; } m_initialValues.append(initialValue); } DEBUG_OUT(QString("Extracted %1 user initial values").arg(m_initialValues.size())); for(int i = 0; i < m_initialValues.size(); ++i) { DEBUG_OUT(QString(" Initial[%1] = %2").arg(i).arg(m_initialValues[i], 0, 'e', 3)); } // 验证初始值 if(!validateInitialValues()) { DEBUG_OUT("Warning: Some initial values are outside parameter bounds"); } } // ==================== 遗传算法核心方法 ==================== void nmCalculationAutoFitGA::initializePopulation() { int dimensions = getEnabledParameterCount(); if(dimensions == 0) return; m_population.clear(); m_population.resize(m_populationSize); bool hasValidInitials = !m_initialValues.isEmpty() && m_initialValues.size() >= dimensions; int guidedCount = hasValidInitials ? qMax(2, m_populationSize / 2) : qMax(1, m_populationSize / 3); DEBUG_OUT(QString("Enhanced population initialization: %1 individuals, %2 guided, %3 random") .arg(m_populationSize).arg(guidedCount).arg(m_populationSize - guidedCount)); for(int i = 0; i < m_populationSize; ++i) { GAIndividual& individual = m_population[i]; individual.genes.resize(dimensions); individual.fitness = 1e10; individual.isEvaluated = false; // 基因初始化 for(int j = 0; j < dimensions; ++j) { int paramIndex = m_enabledParamIndices[j]; double range = m_parameterUpper[paramIndex] - m_parameterLower[paramIndex]; if(i == 0 && hasValidInitials) { // 第一个个体:使用用户初始值 individual.genes[j] = m_initialValues[j]; } else if(i < guidedCount && hasValidInitials) { // 引导搜索策略 double searchRadius; if(i <= guidedCount / 3) { searchRadius = range * 0.03; } else if(i <= guidedCount * 2 / 3) { searchRadius = range * 0.08; } else { searchRadius = range * 0.15; } double offset = (random01() - 0.5) * searchRadius; individual.genes[j] = m_initialValues[j] + offset; } else { // 随机初始化 individual.genes[j] = m_parameterLower[paramIndex] + random01() * range; } } // 边界约束 clampToLimits(individual.genes); } } double nmCalculationAutoFitGA::evaluateGenes(const QVector& genes) { const QString funcName = QString("evaluateGenes[Gen%1]").arg(m_currentGeneration); static int callCount = 0; callCount++; try { DEBUG_OUT(QString("%1: Call #%2 - Starting evaluation with %3 genes") .arg(funcName).arg(callCount).arg(genes.size())); // 打印参数值用于对比 QString paramStr = "Parameters: "; for(int i = 0; i < genes.size(); ++i) { paramStr += QString("[%1]=%2 ").arg(i).arg(genes[i], 0, 'f', 6); } DEBUG_OUT(QString("%1: %2").arg(funcName).arg(paramStr)); // 1. 参数有效性检查 if(!validateParameters(genes)) { DEBUG_OUT(QString("%1: Call #%2 - Invalid parameters").arg(funcName).arg(callCount)); return 1e10; } // 2. 检查数据管理器状态 nmDataAnalyzeManager* dataManager = nmDataAnalyzeManager::getCurrentInstance(); if(!dataManager) { DEBUG_OUT(QString("%1: Call #%2 - DataManager is null").arg(funcName).arg(callCount)); return 1e10; } // 3. 应用参数到数据管理器 try { DEBUG_OUT(QString("%1: Call #%2 - Applying parameters to DataManager").arg(funcName).arg(callCount)); applyParametersToDataManager(genes); DEBUG_OUT(QString("%1: Call #%2 - Parameters applied successfully").arg(funcName).arg(callCount)); } catch(const std::exception& e) { DEBUG_OUT(QString("%1: Call #%2 - Failed to apply parameters: %3").arg(funcName).arg(callCount).arg(e.what())); return 1e10; } catch(...) { DEBUG_OUT(QString("%1: Call #%2 - Unknown error applying parameters").arg(funcName).arg(callCount)); return 1e10; } // 4. 运行求解器 QVector > solverResult; const int maxRetries = 2; bool solverSuccess = false; for(int retry = 0; retry <= maxRetries; ++retry) { if(m_shouldStop) return 1e10; try { DEBUG_OUT(QString("%1: Call #%2 - Solver attempt %3/%4") .arg(funcName).arg(callCount).arg(retry + 1).arg(maxRetries + 1)); // 在求解器调用前添加短暂延迟,确保状态稳定 if(retry > 0) { DEBUG_OUT(QString("%1: Call #%2 - Retry delay before solver attempt") .arg(funcName).arg(callCount)); msleep(1000); // 增加延迟时间 } solverResult = runSolver(); if(!solverResult.isEmpty() && validateSolverResult(solverResult)) { DEBUG_OUT(QString("%1: Call #%2 - Solver successful on attempt %3, result size: %4") .arg(funcName).arg(callCount).arg(retry + 1).arg(solverResult[0].size())); solverSuccess = true; break; } else { DEBUG_OUT(QString("%1: Call #%2 - Solver failed on attempt %3 - empty or invalid result") .arg(funcName).arg(callCount).arg(retry + 1)); if(retry < maxRetries) { DEBUG_OUT(QString("%1: Call #%2 - Will retry solver").arg(funcName).arg(callCount)); } } } catch(const std::exception& e) { DEBUG_OUT(QString("%1: Call #%2 - Solver exception on attempt %3: %4") .arg(funcName).arg(callCount).arg(retry + 1).arg(e.what())); } catch(...) { DEBUG_OUT(QString("%1: Call #%2 - Unknown solver exception on attempt %3") .arg(funcName).arg(callCount).arg(retry + 1)); } } if(!solverSuccess) { DEBUG_OUT(QString("%1: Call #%2 - All solver attempts failed").arg(funcName).arg(callCount)); return 1e10; } // 5. 获取双对数结果数据 QVector > resultLogLogData; try { nmDataWellBase* pTargetWell = dataManager->findWellByName(m_targetWellName); if(pTargetWell) { resultLogLogData = pTargetWell->getResultLogLog(); if(!validateLogLogData(resultLogLogData)) { DEBUG_OUT(QString("%1: Call #%2 - Invalid result LogLog data").arg(funcName).arg(callCount)); return 1e10; } DEBUG_OUT(QString("%1: Call #%2 - LogLog result data obtained, size: %3") .arg(funcName).arg(callCount).arg(resultLogLogData[0].size())); } else { DEBUG_OUT(QString("%1: Call #%2 - No wells found after solver").arg(funcName).arg(callCount)); return 1e10; } } catch(const std::exception& e) { DEBUG_OUT(QString("%1: Call #%2 - Error getting LogLog result: %3").arg(funcName).arg(callCount).arg(e.what())); return 1e10; } catch(...) { DEBUG_OUT(QString("%1: Call #%2 - Unknown error getting LogLog result").arg(funcName).arg(callCount)); return 1e10; } // 6. 双对数曲线对齐和误差计算 double error; try { error = calculateLogLogCurveError(m_targetLogLogData, resultLogLogData); if(!isFiniteNumber(error) || error < 0) { DEBUG_OUT(QString("%1: Call #%2 - Invalid error value: %3").arg(funcName).arg(callCount).arg(error)); return 1e10; } DEBUG_OUT(QString("%1: Call #%2 - Evaluation successful, error = %3") .arg(funcName).arg(callCount).arg(error, 0, 'e', 6)); } catch(const std::exception& e) { DEBUG_OUT(QString("%1: Call #%2 - LogLog error calculation failed: %3").arg(funcName).arg(callCount).arg(e.what())); return 1e10; } catch(...) { DEBUG_OUT(QString("%1: Call #%2 - Unknown error in LogLog error calculation").arg(funcName).arg(callCount)); return 1e10; } return error; } catch(const std::exception& e) { DEBUG_OUT(QString("%1: Call #%2 - Top-level exception: %3").arg(funcName).arg(callCount).arg(e.what())); return 1e10; } catch(...) { DEBUG_OUT(QString("%1: Call #%2 - Unknown top-level exception").arg(funcName).arg(callCount)); return 1e10; } } double nmCalculationAutoFitGA::evaluateIndividual(GAIndividual& individual) { // 如果已经评估过,直接返回 if(individual.isEvaluated) { return individual.fitness; } // 调用新的评估函数 individual.fitness = evaluateGenes(individual.genes); individual.isEvaluated = true; return individual.fitness; } void nmCalculationAutoFitGA::evaluatePopulation() { int successfulEvaluations = 0; int totalEvaluations = 0; int currentGenerationFailed = 0; for(int i = 0; i < m_population.size() && !m_shouldStop; ++i) { GAIndividual& individual = m_population[i]; if(!individual.isEvaluated) { try { double previousFitness = individual.fitness; individual.fitness = evaluateGenes(individual.genes); individual.isEvaluated = true; totalEvaluations++; m_totalEvaluations++; if(individual.fitness < 1e9) { successfulEvaluations++; m_successfulEvaluations++; // 更新全局最优 if(individual.fitness < m_bestFitness) { m_bestFitness = individual.fitness; m_bestIndividual = individual; emit logMessageGenerated(tr(" Individual %1 improved: %2 -> %3") .arg(i + 1).arg(previousFitness, 0, 'e', 3).arg(individual.fitness, 0, 'e', 3)); } } else { currentGenerationFailed++; emit logMessageGenerated(tr(" Individual %1: evaluation failed").arg(i + 1)); } // 每评估2个个体处理一次事件 if(i % 2 == 0) { QApplication::processEvents(); } } catch(const std::exception& e) { emit logMessageGenerated(tr(" Individual %1: Exception: %2").arg(i + 1).arg(e.what())); individual.fitness = 1e10; individual.isEvaluated = true; totalEvaluations++; currentGenerationFailed++; m_totalEvaluations++; } catch(...) { emit logMessageGenerated(tr(" Individual %1: Unknown exception").arg(i + 1)); individual.fitness = 1e10; individual.isEvaluated = true; totalEvaluations++; currentGenerationFailed++; m_totalEvaluations++; } } } if(m_shouldStop) return; double currentSuccessRate = totalEvaluations > 0 ? (double)successfulEvaluations / totalEvaluations : 0.0; // 更新连续失败代数计数 if(successfulEvaluations == 0) { m_consecutiveFailedGenerations++; emit logMessageGenerated(tr("WARNING: No successful evaluations in generation %1 (consecutive failures: %2)") .arg(m_currentGeneration + 1).arg(m_consecutiveFailedGenerations)); } else { m_consecutiveFailedGenerations = 0; // 重置连续失败计数 } // 输出当前代统计 emit logMessageGenerated(tr("Current generation stats: %1 successful, %2 failed out of %3 individuals (success rate: %4%)") .arg(successfulEvaluations).arg(currentGenerationFailed) .arg(totalEvaluations).arg(currentSuccessRate * 100, 0, 'f', 1)); // 检查是否需要停止优化 if(m_consecutiveFailedGenerations >= m_maxConsecutiveFailures) { m_lastError = QString("Too many consecutive failed generations (%1)").arg(m_consecutiveFailedGenerations); emit logMessageGenerated(tr("ERROR: Too many consecutive failed generations (%1/%2) - stopping optimization") .arg(m_consecutiveFailedGenerations).arg(m_maxConsecutiveFailures)); return; } // 警告低成功率但不立即停止 if(currentSuccessRate < 0.5 && m_currentGeneration > 3) { emit logMessageGenerated(tr("WARNING: Low success rate (%1%) in generation %2, but continuing optimization") .arg(currentSuccessRate * 100, 0, 'f', 1).arg(m_currentGeneration + 1)); } } int nmCalculationAutoFitGA::tournamentSelection() { int bestIndex = qrand() % m_population.size(); double bestFitness = m_population[bestIndex].fitness; for(int i = 1; i < m_tournamentSize; ++i) { int candidateIndex = qrand() % m_population.size(); double candidateFitness = m_population[candidateIndex].fitness; if(candidateFitness < bestFitness) { bestIndex = candidateIndex; bestFitness = candidateFitness; } } return bestIndex; } int nmCalculationAutoFitGA::rouletteWheelSelection() { // 计算适应度总和(使用倒数,因为我们要最小化) double totalFitness = 0.0; double maxFitness = -1e10; // 找到最大适应度值 for(int i = 0; i < m_population.size(); ++i) { if(m_population[i].fitness > maxFitness) { maxFitness = m_population[i].fitness; } } // 计算转换后的适应度总和 for(int i = 0; i < m_population.size(); ++i) { double transformedFitness = maxFitness - m_population[i].fitness + 1e-6; totalFitness += transformedFitness; } // 轮盘赌选择 double randomValue = random01() * totalFitness; double cumulativeFitness = 0.0; for(int i = 0; i < m_population.size(); ++i) { double transformedFitness = maxFitness - m_population[i].fitness + 1e-6; cumulativeFitness += transformedFitness; if(cumulativeFitness >= randomValue) { return i; } } return m_population.size() - 1; // 备用选择 } void nmCalculationAutoFitGA::crossover(const GAIndividual& parent1, const GAIndividual& parent2, GAIndividual& offspring1, GAIndividual& offspring2) { if(m_useUniformCrossover) { uniformCrossover(parent1, parent2, offspring1, offspring2); } else { singlePointCrossover(parent1, parent2, offspring1, offspring2); } } void nmCalculationAutoFitGA::singlePointCrossover(const GAIndividual& parent1, const GAIndividual& parent2, GAIndividual& offspring1, GAIndividual& offspring2) { int dimensions = parent1.genes.size(); if(dimensions == 0) return; // 初始化子代 offspring1.genes.resize(dimensions); offspring2.genes.resize(dimensions); offspring1.fitness = 1e10; offspring2.fitness = 1e10; offspring1.isEvaluated = false; offspring2.isEvaluated = false; // 选择交叉点 int crossoverPoint = qrand() % dimensions; // 执行交叉 for(int i = 0; i < dimensions; ++i) { if(i < crossoverPoint) { offspring1.genes[i] = parent1.genes[i]; offspring2.genes[i] = parent2.genes[i]; } else { offspring1.genes[i] = parent2.genes[i]; offspring2.genes[i] = parent1.genes[i]; } } } void nmCalculationAutoFitGA::uniformCrossover(const GAIndividual& parent1, const GAIndividual& parent2, GAIndividual& offspring1, GAIndividual& offspring2) { int dimensions = parent1.genes.size(); if(dimensions == 0) return; // 初始化子代 offspring1.genes.resize(dimensions); offspring2.genes.resize(dimensions); offspring1.fitness = 1e10; offspring2.fitness = 1e10; offspring1.isEvaluated = false; offspring2.isEvaluated = false; // 均匀交叉 for(int i = 0; i < dimensions; ++i) { if(random01() < 0.5) { offspring1.genes[i] = parent1.genes[i]; offspring2.genes[i] = parent2.genes[i]; } else { offspring1.genes[i] = parent2.genes[i]; offspring2.genes[i] = parent1.genes[i]; } } } void nmCalculationAutoFitGA::mutate(GAIndividual& individual) { // 使用高斯变异作为主要方式 gaussianMutation(individual); } void nmCalculationAutoFitGA::gaussianMutation(GAIndividual& individual) { for(int i = 0; i < individual.genes.size(); ++i) { if(random01() < m_mutationRate) { int paramIndex = m_enabledParamIndices[i]; double range = m_parameterUpper[paramIndex] - m_parameterLower[paramIndex]; double sigma = range * MUTATION_STRENGTH; // 高斯变异 double mutation = gaussianRandom(0.0, sigma); individual.genes[i] += mutation; // 边界处理 individual.genes[i] = qMax(m_parameterLower[paramIndex], qMin(m_parameterUpper[paramIndex], individual.genes[i])); } } // 标记为未评估 individual.isEvaluated = false; } void nmCalculationAutoFitGA::polynomialMutation(GAIndividual& individual) { const double eta = 20.0; // 分布指数 for(int i = 0; i < individual.genes.size(); ++i) { if(random01() < m_mutationRate) { int paramIndex = m_enabledParamIndices[i]; double lower = m_parameterLower[paramIndex]; double upper = m_parameterUpper[paramIndex]; double y = individual.genes[i]; double delta1 = (y - lower) / (upper - lower); double delta2 = (upper - y) / (upper - lower); double rnd = random01(); double mut_pow = 1.0 / (eta + 1.0); double deltaq; if(rnd <= 0.5) { double xy = 1.0 - delta1; double val = 2.0 * rnd + (1.0 - 2.0 * rnd) * pow(xy, eta + 1.0); deltaq = pow(val, mut_pow) - 1.0; } else { double xy = 1.0 - delta2; double val = 2.0 * (1.0 - rnd) + 2.0 * (rnd - 0.5) * pow(xy, eta + 1.0); deltaq = 1.0 - pow(val, mut_pow); } y += deltaq * (upper - lower); individual.genes[i] = qMax(lower, qMin(upper, y)); } } // 标记为未评估 individual.isEvaluated = false; } void nmCalculationAutoFitGA::applyElitism(QVector& newPopulation) { int eliteCount = static_cast(m_populationSize * m_elitismRate); if(eliteCount == 0) return; // 对种群按适应度排序 QVector sortedPopulation = m_population; // 冒泡排序(适应度从小到大) for(int i = 0; i < sortedPopulation.size() - 1; ++i) { for(int j = 0; j < sortedPopulation.size() - 1 - i; ++j) { if(sortedPopulation[j].fitness > sortedPopulation[j + 1].fitness) { GAIndividual temp = sortedPopulation[j]; sortedPopulation[j] = sortedPopulation[j + 1]; sortedPopulation[j + 1] = temp; } } } // 复制精英个体到新种群 for(int i = 0; i < eliteCount && i < sortedPopulation.size(); ++i) { newPopulation.append(sortedPopulation[i]); } DEBUG_OUT(QString("Applied elitism: %1 elite individuals preserved").arg(eliteCount)); } void nmCalculationAutoFitGA::updatePopulationStatistics() { if(m_population.isEmpty()) return; m_previousBestFitness = m_bestFitness; double sum = 0.0; double minFitness = 1e10; double maxFitness = -1e10; int validCount = 0; for(int i = 0; i < m_population.size(); ++i) { const GAIndividual& individual = m_population[i]; if(individual.isEvaluated && individual.fitness < 1e9) { sum += individual.fitness; validCount++; if(individual.fitness < minFitness) { minFitness = individual.fitness; } if(individual.fitness > maxFitness) { maxFitness = individual.fitness; } // 更新全局最优 if(individual.fitness < m_bestFitness) { // 只有显著改进时才更新 double improvement = (m_bestFitness - individual.fitness); double relativeImprovement = improvement / qMax(1e-10, qAbs(m_bestFitness)); if(relativeImprovement > m_improvementThreshold) { m_bestFitness = individual.fitness; m_bestIndividual = individual; DEBUG_OUT(QString("Global best updated with %1% improvement: %2") .arg(relativeImprovement * 100, 0, 'f', 3) .arg(m_bestFitness, 0, 'e', 4)); } } } } if(validCount > 0) { m_averageFitness = sum / validCount; m_worstFitness = maxFitness; } else { m_averageFitness = 1e10; m_worstFitness = 1e10; } // 在没有找到更好解时才检查精英保护 if(m_hasValidUserSolution && m_userInitialFitness < m_bestFitness) { DEBUG_OUT("Elite protection: No significant improvement found, checking initial solution"); // 检查初始解是否仍然是最优的 double improvement = m_bestFitness - m_userInitialFitness; double relativeImprovement = improvement / qMax(1e-10, qAbs(m_bestFitness)); if(relativeImprovement > m_improvementThreshold * 0.5) { // 使用更宽松的阈值 DEBUG_OUT("Elite protection: Restoring user initial solution"); m_bestFitness = m_userInitialFitness; m_bestIndividual.genes = m_userInitialSolution; m_bestIndividual.fitness = m_userInitialFitness; m_bestIndividual.isEvaluated = true; } } } bool nmCalculationAutoFitGA::checkConvergence() { if(m_convergenceHistory.size() < CONVERGENCE_CHECK_INTERVAL) { return false; } // 检查最近几次世代的改进 double recentBest = m_convergenceHistory.last(); double oldBest = m_convergenceHistory[m_convergenceHistory.size() - CONVERGENCE_CHECK_INTERVAL]; double improvement = oldBest - recentBest; // 如果连续多代没有显著改进,认为已收敛 if(improvement < MIN_FITNESS_IMPROVEMENT) { // 检查是否连续停滞 int stagnationCount = 0; for(int i = m_convergenceHistory.size() - 1; i >= qMax(0, m_convergenceHistory.size() - MAX_STAGNATION_GENERATIONS); --i) { if(i > 0) { double diff = m_convergenceHistory[i - 1] - m_convergenceHistory[i]; if(diff < MIN_FITNESS_IMPROVEMENT) { stagnationCount++; } else { break; } } } return stagnationCount >= MAX_STAGNATION_GENERATIONS; } return false; } void nmCalculationAutoFitGA::adaptiveParameterUpdate(int generation) { // 自适应调整变异率 double progress = static_cast(generation) / m_maxGenerations; // 早期探索,后期开发 m_mutationRate = 0.2 * (1.0 - progress) + 0.05 * progress; // 自适应调整交叉率 if(generation > 0) { double improvement = m_previousBestFitness - m_bestFitness; if(improvement < MIN_FITNESS_IMPROVEMENT) { // 如果改进很小,增加探索性 m_mutationRate = qMin(0.3, m_mutationRate * 1.1); m_crossoverRate = qMax(0.6, m_crossoverRate * 0.95); } else { // 如果有明显改进,增加开发性 m_mutationRate = qMax(0.05, m_mutationRate * 0.9); m_crossoverRate = qMin(0.9, m_crossoverRate * 1.05); } } } void nmCalculationAutoFitGA::validateAndProtectFinalResult() { if(!m_hasValidUserSolution) { emit logMessageGenerated(tr("No initial solution for elite protection")); return; } emit logMessageGenerated(tr("=== Final Result Validation (Elite Protection) ===")); // 使用已有的评估结果 double finalFitness = m_bestFitness; double initialFitness = m_userInitialFitness; emit logMessageGenerated(tr("Comparing results: Initial=%1, Final=%2") .arg(initialFitness, 0, 'e', 4).arg(finalFitness, 0, 'e', 4)); // 计算改进程度 double improvement = initialFitness - finalFitness; double relativeImprovement = improvement / qMax(1e-10, qAbs(initialFitness)); emit logMessageGenerated(tr("Improvement: %1 (%2%)") .arg(improvement, 0, 'e', 4).arg(relativeImprovement * 100, 0, 'f', 2)); if(relativeImprovement < m_improvementThreshold) { emit logMessageGenerated(tr("Elite protection triggered: insufficient improvement")); emit logMessageGenerated(tr("Threshold: %1%, Actual: %2%") .arg(m_improvementThreshold * 100, 0, 'f', 2) .arg(relativeImprovement * 100, 0, 'f', 4)); emit logMessageGenerated(tr("Restoring initial solution as final result")); m_bestFitness = initialFitness; m_bestIndividual.genes = m_userInitialSolution; m_bestIndividual.fitness = initialFitness; m_bestIndividual.isEvaluated = true; emit logMessageGenerated(tr("Initial solution restored successfully")); } else { emit logMessageGenerated(tr("Final result validated - significant improvement achieved")); } } StopReasonGA nmCalculationAutoFitGA::analyzeOptimizationStatus() { // 1. 检查用户停止 if(m_shouldStop) { return GA_USER_STOPPED; } // 2. 检查连续失败 if(m_consecutiveFailedGenerations >= m_maxConsecutiveFailures) { return GA_CONSECUTIVE_FAILURES; } // 3. 检查是否达到目标精度 if(m_bestFitness < m_targetError) { return GA_TARGET_ACHIEVED; } // 4. 检查是否达到最大迭代数 if(m_currentGeneration >= m_maxGenerations - 1) { return GA_MAX_ITERATIONS; } // 5. 需要足够的历史数据才能判断收敛 if(m_convergenceHistory.size() < m_localOptimumWindow) { return GA_CONTINUE_OPTIMIZATION; } // 6. 检查真正的收敛 if(m_convergenceHistory.size() >= m_trueConvergenceWindow && checkTrueConvergence()) { return GA_TRUE_CONVERGENCE; } // 7. 检查局部最优陷阱 if(checkLocalOptimumTrap()) { return GA_LOCAL_OPTIMUM; } return GA_CONTINUE_OPTIMIZATION; } bool nmCalculationAutoFitGA::checkTrueConvergence() const { if(m_convergenceHistory.size() < m_trueConvergenceWindow) { return false; } // 1. 检查解质量 - 如果已经接近目标,小改进可能是真收敛 bool nearTarget = (m_bestFitness < m_targetError * m_nearTargetFactor); // 2. 检查适应度稳定性 - 长期小幅波动 double recentVariance = calculateFitnessVariance(10); double recentMean = 0.0; int windowSize = qMin(10, m_convergenceHistory.size()); // 计算最近窗口的均值 for(int i = m_convergenceHistory.size() - windowSize; i < m_convergenceHistory.size(); ++i) { recentMean += m_convergenceHistory[i]; } recentMean /= windowSize; double relativeVariance = recentVariance / qMax(1e-10, recentMean * recentMean); bool stableError = (relativeVariance < m_convergenceVarianceThreshold); // 3. 检查种群多样性 - 应该收敛到同一区域 double currentDiversity = calculatePopulationDiversity(); bool lowDiversity = (currentDiversity < m_diversityThreshold); // 4. 检查长期改进趋势 double longTermImprovement = calculateLongTermImprovement(m_trueConvergenceWindow); bool minimalLongTermImprovement = (longTermImprovement < 1e-4); // 0.01% // 真收敛的判断条件 bool isConverged = nearTarget || (stableError && lowDiversity && minimalLongTermImprovement); if(isConverged) { DEBUG_OUT("=== TRUE CONVERGENCE ANALYSIS ==="); DEBUG_OUT(QString("nearTarget=%1 (fitness=%2, target*factor=%3)") .arg(nearTarget).arg(m_bestFitness, 0, 'e', 4) .arg(m_targetError * m_nearTargetFactor, 0, 'e', 4)); DEBUG_OUT(QString("stableError=%1 (relativeVariance=%2)") .arg(stableError).arg(relativeVariance, 0, 'e', 6)); DEBUG_OUT(QString("lowDiversity=%1 (diversity=%2, threshold=%3)") .arg(lowDiversity).arg(currentDiversity, 0, 'f', 6).arg(m_diversityThreshold)); DEBUG_OUT(QString("longTermImprovement=%1%") .arg(longTermImprovement * 100, 0, 'f', 4)); } return isConverged; } bool nmCalculationAutoFitGA::checkLocalOptimumTrap() const { if(m_convergenceHistory.size() < m_localOptimumWindow) { return false; } // 1. 检查解质量 - 如果距离目标还很远,停滞就可能是局部最优 bool farFromTarget = (m_bestFitness > m_targetError * m_farTargetFactor); // 2. 检查短期改进 - 近期改进非常小 double shortTermImprovement = calculateLongTermImprovement(m_localOptimumWindow); bool poorShortTermImprovement = (shortTermImprovement < 1e-5); // 0.001% // 3. 检查种群多样性 - 可能过早聚集或无效分散 double currentDiversity = calculatePopulationDiversity(); bool problematicDiversity = (currentDiversity < m_diversityThreshold * 0.1) || (currentDiversity > m_diversityThreshold * 5.0); // 4. 检查适应度方差 - 可能卡在平坦区域 double fitnessVariance = calculateFitnessVariance(m_localOptimumWindow); bool flatFitnessLandscape = (fitnessVariance < m_convergenceVarianceThreshold * 0.1); // 局部最优的判断条件 bool isLocalOptimum = farFromTarget && poorShortTermImprovement && (problematicDiversity || flatFitnessLandscape); if(isLocalOptimum) { DEBUG_OUT("=== LOCAL OPTIMUM ANALYSIS ==="); DEBUG_OUT(QString("farFromTarget=%1 (fitness=%2, target*factor=%3)") .arg(farFromTarget).arg(m_bestFitness, 0, 'e', 4) .arg(m_targetError * m_farTargetFactor, 0, 'e', 4)); DEBUG_OUT(QString("poorShortTermImprovement=%1 (improvement=%2%)") .arg(poorShortTermImprovement).arg(shortTermImprovement * 100, 0, 'f', 4)); DEBUG_OUT(QString("problematicDiversity=%1 (diversity=%2)") .arg(problematicDiversity).arg(currentDiversity, 0, 'f', 6)); DEBUG_OUT(QString("fitnessVariance=%1") .arg(fitnessVariance, 0, 'e', 6)); } return isLocalOptimum; } double nmCalculationAutoFitGA::calculatePopulationDiversity() const { if(m_population.size() < 2) return 0.0; int dimensions = getEnabledParameterCount(); if(dimensions == 0) return 0.0; double totalDiversity = 0.0; for(int dim = 0; dim < dimensions; ++dim) { // 计算该维度上所有个体的均值 double mean = 0.0; for(int i = 0; i < m_population.size(); ++i) { mean += m_population[i].genes[dim]; } mean /= m_population.size(); // 计算该维度上的方差 double variance = 0.0; for(int i = 0; i < m_population.size(); ++i) { double diff = m_population[i].genes[dim] - mean; variance += diff * diff; } variance /= m_population.size(); // 归一化到参数范围 int paramIndex = m_enabledParamIndices[dim]; double range = m_parameterUpper[paramIndex] - m_parameterLower[paramIndex]; double normalizedStd = sqrt(variance) / qMax(1e-10, range); totalDiversity += normalizedStd; } return totalDiversity / dimensions; } double nmCalculationAutoFitGA::calculateFitnessVariance(int windowSize) const { if(m_convergenceHistory.size() < windowSize) { return 1e10; // 数据不足,返回大值 } // 计算最近windowSize次迭代的方差 double mean = 0.0; int startIdx = m_convergenceHistory.size() - windowSize; for(int i = startIdx; i < m_convergenceHistory.size(); ++i) { mean += m_convergenceHistory[i]; } mean /= windowSize; double variance = 0.0; for(int i = startIdx; i < m_convergenceHistory.size(); ++i) { double diff = m_convergenceHistory[i] - mean; variance += diff * diff; } variance /= windowSize; return variance; } double nmCalculationAutoFitGA::calculateLongTermImprovement(int windowSize) const { if(m_convergenceHistory.size() < windowSize) { return 1.0; // 数据不足,假设有改进 } double oldFitness = m_convergenceHistory[m_convergenceHistory.size() - windowSize]; double improvement = (oldFitness - m_bestFitness) / qMax(1e-10, qAbs(oldFitness)); return improvement; } void nmCalculationAutoFitGA::updateConvergenceMetrics() { // 更新多样性历史 m_diversityHistory.append(calculatePopulationDiversity()); // 保持历史长度合理(最多保留50个数据点) const int maxHistorySize = 50; while(m_diversityHistory.size() > maxHistorySize) { m_diversityHistory.remove(0); } } // ==================== 参数应用方法 ==================== bool nmCalculationAutoFitGA::validateParameters(const QVector& parameters) const { if(parameters.size() != getEnabledParameterCount()) { return false; } for(int i = 0; i < parameters.size(); ++i) { if(!isFiniteNumber(parameters[i])) { return false; } // 检查参数范围 if(i < m_enabledParamIndices.size()) { int paramIndex = m_enabledParamIndices[i]; if(paramIndex >= 0 && paramIndex < m_parameterLower.size()) { if(parameters[i] < m_parameterLower[paramIndex] || parameters[i] > m_parameterUpper[paramIndex]) { return false; } } } } // 直接拦截会导致求解器数值崩溃的参数值 for(int i = 0; i < parameters.size() && i < m_enabledParamIndices.size(); ++i) { int paramIndex = m_enabledParamIndices[i]; double value = parameters[i]; switch(paramIndex) { case 0: // 渗透率:必须大于零 if(value <= 1e-8) { DEBUG_OUT(QString("Rejecting near-zero permeability: %1").arg(value)); return false; } break; case 2: // 井筒储集系数:必须大于零 if(value <= 1e-10) { DEBUG_OUT(QString("Rejecting near-zero wellbore storage: %1").arg(value)); return false; } break; case 3: // 孔隙度:必须在合理范围 if(value <= 1e-6 || value >= 0.99) { DEBUG_OUT(QString("Rejecting unrealistic porosity: %1").arg(value)); return false; } break; case 6: // 综合压缩系数:必须大于零 if(value <= 1e-8) { DEBUG_OUT(QString("Rejecting near-zero total compressibility: %1").arg(value)); return false; } break; } } return true; } bool nmCalculationAutoFitGA::validateLogLogData(const QVector>& logLogData) const { // 检查基本结构 if(logLogData.size() < 3) { DEBUG_OUT("LogLog data has less than 3 arrays"); return false; } // 检查数组大小一致性 int size = logLogData[0].size(); if(size == 0) { DEBUG_OUT("Empty LogLog data"); return false; } if(logLogData[1].size() != size || logLogData[2].size() != size) { DEBUG_OUT(QString("LogLog data size mismatch: X=%1, Y1=%2, Y2=%3") .arg(logLogData[0].size()) .arg(logLogData[1].size()) .arg(logLogData[2].size())); return false; } // 检查最小数据点数 if(size < 5) { DEBUG_OUT(QString("Too few LogLog data points: %1").arg(size)); return false; } // 数据有效性检查 for(int i = 0; i < size; ++i) { if(!isFiniteNumber(logLogData[0][i]) || !isFiniteNumber(logLogData[1][i]) || !isFiniteNumber(logLogData[2][i])) { DEBUG_OUT(QString("Invalid LogLog data at index %1").arg(i)); return false; } } return true; } bool nmCalculationAutoFitGA::validateInitialValues() const { if(m_initialValues.size() != m_enabledParamIndices.size()) { DEBUG_OUT("Initial values count mismatch with enabled parameters"); return false; } bool allValid = true; for(int i = 0; i < m_initialValues.size(); ++i) { int paramIndex = m_enabledParamIndices[i]; double value = m_initialValues[i]; if(!isFiniteNumber(value)) { DEBUG_OUT(QString("Initial value[%1] is not finite: %2").arg(i).arg(value)); allValid = false; continue; } if(paramIndex < m_parameterLower.size() && paramIndex < m_parameterUpper.size()) { double minVal = m_parameterLower[paramIndex]; double maxVal = m_parameterUpper[paramIndex]; if(value < minVal || value > maxVal) { DEBUG_OUT(QString("Initial value[%1] = %2 is outside bounds [%3, %4]") .arg(i).arg(value).arg(minVal).arg(maxVal)); allValid = false; } } } return allValid; } void nmCalculationAutoFitGA::applyParametersToDataManager(const QVector& parameters) { if(parameters.size() != getEnabledParameterCount()) { return; } updateReservoirParameters(parameters); updateWellParameters(parameters); } void nmCalculationAutoFitGA::updateReservoirParameters(const QVector& parameters) { nmDataAnalyzeManager* dataManager = nmDataAnalyzeManager::getCurrentInstance(); nmDataReservoir reservoirData = dataManager->getReservoirDataCopy(); int paramIndex = 0; for(int i = 0; i < m_parameterSelected.size() && i < 11; ++i) { if(m_parameterSelected[i] && paramIndex < parameters.size()) { double value = parameters[paramIndex]; switch(i) { case 0: // 渗透率 reservoirData.getPermeability().setValue(value); break; case 3: // 孔隙度 reservoirData.getPorosity().setValue(value); break; case 4: // 初始压力 reservoirData.getInitialPressure().setValue(value); break; case 5: // 储层厚度 reservoirData.getThickness().setValue(value); break; case 6: // 综合压缩系数 reservoirData.getCt().setValue(value); break; case 7: // 岩石压缩系数 reservoirData.getCf().setValue(value); break; case 8: // 初始含油饱和度 reservoirData.getSoi().setValue(value); break; case 9: // 初始含水饱和度 reservoirData.getSwi().setValue(value); break; case 10: // 初始含气饱和度 reservoirData.getSgi().setValue(value); break; } paramIndex++; } } // 更新数据管理器 dataManager->updateReservoirData(reservoirData); } void nmCalculationAutoFitGA::updateWellParameters(const QVector& parameters) { nmDataAnalyzeManager* dataManager = nmDataAnalyzeManager::getCurrentInstance(); // 获取井列表,更新第一口井的参数 //QVector wells = dataManager->getWellDataList(); nmDataWellBase* pWell = dataManager->findWellByName(m_targetWellName); if(!pWell) return; //nmDataWellBase* pWell = wells[0]; // 使用第一口井 int paramIndex = 0; for(int i = 0; i < m_parameterSelected.size() && i < 9; ++i) { if(m_parameterSelected[i] && paramIndex < parameters.size()) { double value = parameters[paramIndex]; switch(i) { case 1: { // 表皮系数 nmDataAttribute skinAttr = pWell->getPerforation(0)->getSkin(); skinAttr.setValue(value); pWell->setRateDependentSkin(skinAttr); } break; case 2: { // 井筒储集系数 nmDataAttribute wellboreAttr = pWell->getWellboreStorage(); wellboreAttr.setValue(value); pWell->setWellboreStorage(wellboreAttr); } break; } paramIndex++; } } // 根据井类型更新到数据管理器 updateWellToDataManager(pWell); } void nmCalculationAutoFitGA::updateWellToDataManager(nmDataWellBase* pWell) { if(!pWell) return; nmDataAnalyzeManager* dataManager = nmDataAnalyzeManager::getCurrentInstance(); NM_WELL_MODEL wellType = pWell->getWellType(); switch(wellType) { case NM_WELL_MODEL::Vertical_Well: { nmDataVerticalWell* pVerticalWell = dynamic_cast(pWell); if(pVerticalWell) { QVector wells; wells.append(*pVerticalWell); dataManager->updateVerticalWells(wells); } break; } case NM_WELL_MODEL::Vertical_Fractured_Well: { nmDataVerticalFracturedWell* pVFracturedWell = dynamic_cast(pWell); if(pVFracturedWell) { QVector wells; wells.append(*pVFracturedWell); dataManager->updateVerticalFracturedWells(wells); } break; } case NM_WELL_MODEL::Horizontal_Fractured_Well: { nmDataHorizontalFracturedWell* pHFracturedWell = dynamic_cast(pWell); if(pHFracturedWell) { QVector wells; wells.append(*pHFracturedWell); dataManager->updateHorizontalFracturedWells(wells); } break; } default: break; } } void nmCalculationAutoFitGA::clampToLimits(QVector& parameters) const { for(int i = 0; i < parameters.size() && i < m_enabledParamIndices.size(); ++i) { int paramIndex = m_enabledParamIndices[i]; if(paramIndex >= 0 && paramIndex < m_parameterLower.size()) { parameters[i] = qMax(m_parameterLower[paramIndex], qMin(m_parameterUpper[paramIndex], parameters[i])); } } } // ==================== 求解器相关方法 ==================== QVector> nmCalculationAutoFitGA::runSolver() { //return runSolverExe(); return runSolverDll(); } bool nmCalculationAutoFitGA::validateSolverResult(const QVector>& result) const { // 基本检查 if(result.size() < 2) { DEBUG_OUT("Solver result has less than 2 arrays"); return false; } if(result[0].size() != result[1].size()) { DEBUG_OUT(QString("Size mismatch: X=%1, Y=%2").arg(result[0].size()).arg(result[1].size())); return false; } if(result[0].size() == 0) { DEBUG_OUT("Empty solver result"); return false; } // 检查最小数据点数 if(result[0].size() < 10) { DEBUG_OUT(QString("Too few data points: %1").arg(result[0].size())); return false; } // 数据有效性检查 for(int i = 0; i < result[0].size(); ++i) { if(!isFiniteNumber(result[0][i]) || !isFiniteNumber(result[1][i])) { DEBUG_OUT(QString("Invalid data at index %1: X=%2, Y=%3") .arg(i).arg(result[0][i]).arg(result[1][i])); return false; } } return true; } QVector> nmCalculationAutoFitGA::runSolverDll() { DEBUG_OUT("SOLVER DLL START"); if(m_evaluationInProgress > 0) { DEBUG_OUT("DLL Solver already running, skipping"); return QVector>(); } ++m_evaluationInProgress; QVector> result; nmCalculationDllPebiSolverTask* dllTask = nullptr; try { DEBUG_OUT("Creating DLL solver task"); dllTask = new nmCalculationDllPebiSolverTask(m_tempDirectory); if(m_shouldStop) { DEBUG_OUT("Should stop - cleaning up and returning empty result"); delete dllTask; --m_evaluationInProgress; return result; } DEBUG_OUT("Starting DLL solver execution..."); // 异步执行 dllTask->start(); // 等待完成 int waitTime = 0; const int maxWait = 30000; // 30秒超时 const int checkInterval = 500; while(waitTime < maxWait) { QApplication::processEvents(QEventLoop::ExcludeUserInputEvents, 100); if(!dllTask->isRunning()) { DEBUG_OUT("DLL solver task completed"); break; } if(m_shouldStop) { DEBUG_OUT("DLL solver task terminated by user"); dllTask->terminate(); break; } msleep(checkInterval); waitTime += checkInterval; } // 超时处理 if(dllTask->isRunning()) { DEBUG_OUT("DLL solver task timeout, terminating..."); dllTask->terminate(); dllTask->wait(2000); delete dllTask; dllTask = nullptr; --m_evaluationInProgress; m_consecutiveFailures++; return result; } // 验证结果数据是否已更新 nmDataAnalyzeManager* dataManager = nmDataAnalyzeManager::getCurrentInstance(); //QVector wells = dataManager->getWellDataList(); nmDataWellBase* pTargetWell = dataManager->findWellByName(m_targetWellName); if(!pTargetWell) { DEBUG_OUT("No wells found in data manager after DLL execution"); delete dllTask; --m_evaluationInProgress; return result; } // 验证结果数据 QVector> pressureResult = pTargetWell->getResultPressure(); QVector> logLogResult = pTargetWell->getResultLogLog(); DEBUG_OUT(QString("DLL result verification - Pressure arrays: %1, LogLog arrays: %2") .arg(pressureResult.size()).arg(logLogResult.size())); if(pressureResult.size() >= 2) { DEBUG_OUT(QString("Pressure result - Time points: %1, Pressure points: %2") .arg(pressureResult[0].size()).arg(pressureResult[1].size())); if(pressureResult[0].size() > 0) { DEBUG_OUT(QString("Sample pressure data - Time[0]: %1, Time[last]: %2, P[0]: %3, P[last]: %4") .arg(pressureResult[0][0]) .arg(pressureResult[0][pressureResult[0].size()-1]) .arg(pressureResult[1][0]) .arg(pressureResult[1][pressureResult[1].size()-1])); } } // 数据有效性检查 if(pressureResult.size() >= 2 && pressureResult[0].size() > 0 && pressureResult[1].size() > 0) { result = pressureResult; DEBUG_OUT(QString("Got DLL solver result: %1 points").arg(result[0].size())); m_consecutiveFailures = 0; // 检查结果是否与之前不同 static QVector lastPressureResult; bool isDifferentFromLast = false; if(lastPressureResult.isEmpty() || lastPressureResult.size() != pressureResult[1].size()) { isDifferentFromLast = true; } else { for(int i = 0; i < qMin(5, pressureResult[1].size()); ++i) { if(qAbs(lastPressureResult[i] - pressureResult[1][i]) > 1e-12) { isDifferentFromLast = true; break; } } } if(isDifferentFromLast) { DEBUG_OUT("RESULT VERIFICATION: Got NEW result data from DLL"); lastPressureResult = pressureResult[1]; } else { DEBUG_OUT("!!! WARNING: Result data appears to be identical to previous run !!!"); } } else { DEBUG_OUT("DLL solver result is empty or invalid"); DEBUG_OUT(QString("Pressure result size: %1, Array sizes: %2, %3") .arg(pressureResult.size()) .arg(pressureResult.size() > 0 ? pressureResult[0].size() : 0) .arg(pressureResult.size() > 1 ? pressureResult[1].size() : 0)); m_consecutiveFailures++; } } catch(const std::bad_alloc& e) { DEBUG_OUT(QString("Memory allocation failed in DLL solver: %1").arg(e.what())); m_consecutiveFailures++; } catch(const std::exception& e) { DEBUG_OUT(QString("Exception in DLL solver: %1").arg(e.what())); m_consecutiveFailures++; } catch(...) { DEBUG_OUT("Unknown exception in DLL solver"); m_consecutiveFailures++; } // 清理DLL任务 if(dllTask) { if(dllTask->isRunning()) { dllTask->terminate(); dllTask->wait(3000); } DEBUG_OUT("Cleaning up DLL solver task..."); delete dllTask; dllTask = nullptr; } QApplication::processEvents(QEventLoop::AllEvents, 100); --m_evaluationInProgress; DEBUG_OUT(QString("SOLVER DLL END - ResultPoints: %1") .arg(result.isEmpty() ? 0 : result[0].size())); return result; } //QVector> nmCalculationAutoFitGA::runSolverExe() //{ // DEBUG_OUT("SOLVER EXE START"); // // if(m_evaluationInProgress > 0) { // DEBUG_OUT("EXE Solver already running, skipping"); // return QVector>(); // } // // ++m_evaluationInProgress; // QVector> result; // nmCalculationExeSolverTask* exeTask = nullptr; // // try { // DEBUG_OUT("Creating EXE solver task"); // exeTask = new nmCalculationExeSolverTask(QString()); // // if(m_shouldStop) { // DEBUG_OUT("Should stop - cleaning up and returning empty result"); // delete exeTask; // --m_evaluationInProgress; // return result; // } // // DEBUG_OUT("Starting EXE solver execution..."); // // // 启动 // bool executeSuccess = exeTask->execute(); // // // 检查执行状态 // if(!executeSuccess) { // DEBUG_OUT(QString("EXE solver execution failed: %1").arg(exeTask->getLastError())); // DEBUG_OUT(QString("EXE solver exit code: %1").arg(exeTask->getExitCode())); // // // 清理并返回空结果 // delete exeTask; // exeTask = nullptr; // --m_evaluationInProgress; // m_consecutiveFailures++; // return result; // } // // DEBUG_OUT("EXE solver execution completed successfully"); // // // 验证结果数据是否已更新 // nmDataAnalyzeManager* dataManager = nmDataAnalyzeManager::getCurrentInstance(); // QVector wells = dataManager->getWellDataList(); // // if(wells.isEmpty()) { // DEBUG_OUT("No wells found in data manager after EXE execution"); // delete exeTask; // --m_evaluationInProgress; // return result; // } // // // 验证结果数据的时间戳或唯一性 // QVector> pressureResult = wells[0]->getResultPressure(); // QVector> logLogResult = wells[0]->getResultLogLog(); // // DEBUG_OUT(QString("Raw result verification - Pressure arrays: %1, LogLog arrays: %2") // .arg(pressureResult.size()).arg(logLogResult.size())); // // if(pressureResult.size() >= 2) { // DEBUG_OUT(QString("Pressure result - Time points: %1, Pressure points: %2") // .arg(pressureResult[0].size()).arg(pressureResult[1].size())); // // if(pressureResult[0].size() > 0) { // DEBUG_OUT(QString("Sample pressure data - Time[0]: %1, Time[last]: %2, P[0]: %3, P[last]: %4") // .arg(pressureResult[0][0]) // .arg(pressureResult[0][pressureResult[0].size() - 1]) // .arg(pressureResult[1][0]) // .arg(pressureResult[1][pressureResult[1].size() - 1])); // } // } // // if(logLogResult.size() >= 3) { // DEBUG_OUT(QString("LogLog result - X points: %1, Y1 points: %2, Y2 points: %3") // .arg(logLogResult[0].size()).arg(logLogResult[1].size()).arg(logLogResult[2].size())); // // if(logLogResult[0].size() > 0) { // DEBUG_OUT(QString("Sample LogLog data - X[0]: %1, X[last]: %2, Y1[0]: %3, Y2[0]: %4") // .arg(logLogResult[0][0]) // .arg(logLogResult[0][logLogResult[0].size() - 1]) // .arg(logLogResult[1][0]) // .arg(logLogResult[2][0])); // } // } // // // 数据有效性检查 // if(pressureResult.size() >= 2 && pressureResult[0].size() > 0 && pressureResult[1].size() > 0) { // result = pressureResult; // DEBUG_OUT(QString("Got EXE solver result: %1 points").arg(result[0].size())); // m_consecutiveFailures = 0; // // // 检查结果是否与之前不同 // static QVector lastPressureResult; // bool isDifferentFromLast = false; // // if(lastPressureResult.isEmpty() || lastPressureResult.size() != pressureResult[1].size()) { // isDifferentFromLast = true; // } else { // for(int i = 0; i < qMin(5, pressureResult[1].size()); ++i) { // if(qAbs(lastPressureResult[i] - pressureResult[1][i]) > 1e-12) { // isDifferentFromLast = true; // break; // } // } // } // // if(isDifferentFromLast) { // DEBUG_OUT("RESULT VERIFICATION: Got NEW result data from EXE"); // lastPressureResult = pressureResult[1]; // } else { // DEBUG_OUT("!!! WARNING: Result data appears to be identical to previous run !!!"); // } // // } else { // DEBUG_OUT("EXE solver result is empty or invalid"); // DEBUG_OUT(QString("Pressure result size: %1, Array sizes: %2, %3") // .arg(pressureResult.size()) // .arg(pressureResult.size() > 0 ? pressureResult[0].size() : 0) // .arg(pressureResult.size() > 1 ? pressureResult[1].size() : 0)); // m_consecutiveFailures++; // } // // } catch(const std::exception& e) { // DEBUG_OUT(QString("Exception in EXE solver: %1").arg(e.what())); // m_consecutiveFailures++; // } catch(...) { // DEBUG_OUT("Unknown exception in EXE solver"); // m_consecutiveFailures++; // } // // // 清理EXE任务 // if(exeTask) { // DEBUG_OUT("Cleaning up EXE solver task..."); // delete exeTask; // exeTask = nullptr; // } // // QApplication::processEvents(QEventLoop::AllEvents, 100); // --m_evaluationInProgress; // // DEBUG_OUT(QString("SOLVER EXE END - ResultPoints: %1") // .arg(result.isEmpty() ? 0 : result[0].size())); // // return result; //} // ==================== 数据处理方法 ==================== QVector nmCalculationAutoFitGA::interpolateData( const QVector& source, const QVector& targetX) const { QVector result; if(source.isEmpty() || targetX.isEmpty()) { DEBUG_OUT("Warning: Empty data in interpolation"); return result; } // 数据清理和验证合并 QVector validSource; const double MAX_REASONABLE_VALUE = 1e12; for(int i = 0; i < source.size(); ++i) { const QPointF& point = source[i]; // 检查数值有效性 if(!isFiniteNumber(point.x()) || !isFiniteNumber(point.y())) { DEBUG_OUT(QString("Skipping invalid data point at index %1: X=%2, Y=%3") .arg(i).arg(point.x()).arg(point.y())); continue; } // 检查极大值 - 跳过而不是失败 if(qAbs(point.y()) > MAX_REASONABLE_VALUE) { DEBUG_OUT(QString("Skipping extremely large Y value at index %1: %2") .arg(i).arg(point.y())); continue; } // 只保留有效的数据点 validSource.append(point); } // 检查清理后的数据是否足够 if(validSource.size() < 3) { DEBUG_OUT(QString("Insufficient valid data points after cleaning: %1") .arg(validSource.size())); return result; // 返回空结果,但不算失败 } DEBUG_OUT(QString("Data cleaning: %1 -> %2 valid points") .arg(source.size()).arg(validSource.size())); // 对清理后的数据进行排序 QVector sortedSource = validSource; for(int i = 0; i < sortedSource.size() - 1; ++i) { for(int j = 0; j < sortedSource.size() - 1 - i; ++j) { if(sortedSource[j].x() > sortedSource[j + 1].x()) { QPointF temp = sortedSource[j]; sortedSource[j] = sortedSource[j + 1]; sortedSource[j + 1] = temp; } } } double sourceMinX = sortedSource.first().x(); double sourceMaxX = sortedSource.last().x(); double targetMinX = targetX[0]; double targetMaxX = targetX[0]; for(int i = 1; i < targetX.size(); ++i) { if(targetX[i] < targetMinX) targetMinX = targetX[i]; if(targetX[i] > targetMaxX) targetMaxX = targetX[i]; } DEBUG_OUT(QString("Source X range: [%1, %2], Target X range: [%3, %4]") .arg(sourceMinX).arg(sourceMaxX).arg(targetMinX).arg(targetMaxX)); // 安全插值算法 for(int i = 0; i < targetX.size(); ++i) { double x = targetX[i]; if(!isFiniteNumber(x)) continue; double y = 0.0; // 插值逻辑(保持原有逻辑,但使用sortedSource) if(x <= sourceMinX) { if(sortedSource.size() >= 2) { double dx = sortedSource[1].x() - sortedSource[0].x(); if(qAbs(dx) > 1e-10) { double slope = (sortedSource[1].y() - sortedSource[0].y()) / dx; slope = qMax(-1e6, qMin(1e6, slope)); y = sortedSource[0].y() + slope * (x - sortedSource[0].x()); } else { y = sortedSource[0].y(); } } else { y = sortedSource[0].y(); } } else if(x >= sourceMaxX) { if(sortedSource.size() >= 2) { int lastIdx = sortedSource.size() - 1; double dx = sortedSource[lastIdx].x() - sortedSource[lastIdx - 1].x(); if(qAbs(dx) > 1e-10) { double slope = (sortedSource[lastIdx].y() - sortedSource[lastIdx - 1].y()) / dx; slope = qMax(-1e6, qMin(1e6, slope)); y = sortedSource[lastIdx].y() + slope * (x - sortedSource[lastIdx].x()); } else { y = sortedSource[lastIdx].y(); } } else { y = sortedSource.last().y(); } } else { // 内插 bool found = false; for(int j = 0; j < sortedSource.size() - 1; ++j) { if(x >= sortedSource[j].x() && x <= sortedSource[j + 1].x()) { double dx = sortedSource[j + 1].x() - sortedSource[j].x(); if(qAbs(dx) > 1e-10) { double ratio = (x - sortedSource[j].x()) / dx; y = sortedSource[j].y() + ratio * (sortedSource[j + 1].y() - sortedSource[j].y()); } else { y = sortedSource[j].y(); } found = true; break; } } if(!found) { // 使用最近点 double minDist = 1e10; for(int k = 0; k < sortedSource.size(); ++k) { double dist = qAbs(sortedSource[k].x() - x); if(dist < minDist) { minDist = dist; y = sortedSource[k].y(); } } } } // 最终数值检查 if(!isFiniteNumber(y)) { y = sortedSource.size() > 0 ? sortedSource[0].y() : 1.0; } y = qMax(-1e12, qMin(1e12, y)); result.append(QPointF(x, y)); } if(result.isEmpty()) { DEBUG_OUT("LogLog interpolation failed"); return result; } DEBUG_OUT(QString("Interpolation completed: %1 -> %2 points") .arg(validSource.size()).arg(result.size())); return result; } double nmCalculationAutoFitGA::calculateLogLogCurveError( const QVector>& target, const QVector>& result) const { // 验证数据 if(!validateLogLogData(target) || !validateLogLogData(result)) { return 1e10; } try { // 数据对齐:找到X值的重叠区域 double targetMinX = target[0][0]; double targetMaxX = target[0][0]; for(int i = 1; i < target[0].size(); ++i) { if(target[0][i] < targetMinX) targetMinX = target[0][i]; if(target[0][i] > targetMaxX) targetMaxX = target[0][i]; } double resultMinX = result[0][0]; double resultMaxX = result[0][0]; for(int i = 1; i < result[0].size(); ++i) { if(result[0][i] < resultMinX) resultMinX = result[0][i]; if(result[0][i] > resultMaxX) resultMaxX = result[0][i]; } double overlapMinX = qMax(targetMinX, resultMinX); double overlapMaxX = qMin(targetMaxX, resultMaxX); if(overlapMinX >= overlapMaxX) { DEBUG_OUT("No overlap between target and result LogLog curves"); return 1e10; } // 生成公共X网格进行插值 QVector commonX; int numPoints = 50; if(overlapMinX > 0 && overlapMaxX > 0) { // 对数空间均匀分布 double logMin = qLn(overlapMinX); double logMax = qLn(overlapMaxX); for(int i = 0; i < numPoints; ++i) { double logX = logMin + i * (logMax - logMin) / (numPoints - 1); double x = qExp(logX); // 数值保护 if(!isFiniteNumber(x) || x <= 0) { continue; } commonX.append(x); } DEBUG_OUT("Using log-uniform grid for better early-time coverage"); } if(commonX.isEmpty()) { DEBUG_OUT("Failed to generate common X grid"); return 1e10; } // 插值目标曲线 QVector targetCurve1, targetCurve2; for(int i = 0; i < target[0].size(); ++i) { // 检查数据有效性 if(isFiniteNumber(target[0][i]) && isFiniteNumber(target[1][i]) && isFiniteNumber(target[2][i])) { targetCurve1.append(QPointF(target[0][i], target[1][i])); targetCurve2.append(QPointF(target[0][i], target[2][i])); } } if(targetCurve1.isEmpty() || targetCurve2.isEmpty()) { DEBUG_OUT("Target curves are empty after filtering"); return 1e10; } QVector alignedTarget1 = interpolateData(targetCurve1, commonX); QVector alignedTarget2 = interpolateData(targetCurve2, commonX); // 插值结果曲线 QVector resultCurve1, resultCurve2; for(int i = 0; i < result[0].size(); ++i) { // 检查数据有效性 if(isFiniteNumber(result[0][i]) && isFiniteNumber(result[1][i]) && isFiniteNumber(result[2][i])) { resultCurve1.append(QPointF(result[0][i], result[1][i])); resultCurve2.append(QPointF(result[0][i], result[2][i])); } } if(resultCurve1.isEmpty() || resultCurve2.isEmpty()) { DEBUG_OUT("Result curves are empty after filtering"); return 1e10; } QVector alignedResult1 = interpolateData(resultCurve1, commonX); QVector alignedResult2 = interpolateData(resultCurve2, commonX); // 检查插值结果 if(alignedTarget1.isEmpty() || alignedTarget2.isEmpty() || alignedResult1.isEmpty() || alignedResult2.isEmpty()) { DEBUG_OUT("LogLog interpolation failed"); return 1e10; } if(alignedTarget1.size() != alignedResult1.size() || alignedTarget2.size() != alignedResult2.size()) { DEBUG_OUT("LogLog interpolation size mismatch"); return 1e10; } // 计算两条曲线的误差 double error1 = calculateCurveError(alignedTarget1, alignedResult1); double error2 = calculateCurveError(alignedTarget2, alignedResult2); // 检查个别误差是否有效 if(!isFiniteNumber(error1) || error1 > 1e9) { DEBUG_OUT(QString("Curve1 error is invalid: %1").arg(error1)); error1 = 1e10; } if(!isFiniteNumber(error2) || error2 > 1e9) { DEBUG_OUT(QString("Curve2 error is invalid: %1").arg(error2)); error2 = 1e10; } // 组合误差 - 添加保护 double combinedError; if(error1 > 1e9 && error2 > 1e9) { combinedError = 1e10; } else if(error1 > 1e9) { combinedError = error2; } else if(error2 > 1e9) { combinedError = error1; } else { combinedError = 0.5 * error1 + 0.5 * error2; } DEBUG_OUT(QString("LogLog errors: Curve1=%1, Curve2=%2, Combined=%3") .arg(error1, 0, 'e', 4).arg(error2, 0, 'e', 4).arg(combinedError, 0, 'e', 4)); return qMin(1e9, combinedError); } catch(const std::exception& e) { DEBUG_OUT(QString("Exception in LogLog error calculation: %1").arg(e.what())); return 1e10; } catch(...) { DEBUG_OUT("Unknown exception in LogLog error calculation"); return 1e10; } } double nmCalculationAutoFitGA::calculateCurveError( const QVector& curve1, const QVector& curve2) const { if(curve1.size() != curve2.size() || curve1.isEmpty()) { return 1e10; } // 预检查:确保没有无穷大值 for(int i = 0; i < curve1.size(); ++i) { if(!isFiniteNumber(curve1[i].y()) || !isFiniteNumber(curve2[i].y())) { DEBUG_OUT(QString("Infinite value detected at index %1: Y1=%2, Y2=%3") .arg(i).arg(curve1[i].y()).arg(curve2[i].y())); return 1e10; } if(qAbs(curve1[i].y()) > 1e12 || qAbs(curve2[i].y()) > 1e12) { DEBUG_OUT(QString("Extremely large value detected at index %1") .arg(i)); return 1e10; } } double totalError = 0.0; double totalWeight = 0.0; int validPoints = 0; for(int i = 0; i < curve1.size(); ++i) { double y1 = curve1[i].y(); double y2 = curve2[i].y(); // 跳过异常值 if(!isFiniteNumber(y1) || !isFiniteNumber(y2)) { continue; } // 自适应权重:根据Y值大小调整,添加上限 double weightFactor = qMin(100.0, qAbs(y1) * 0.01); double weight = 1.0 / (1.0 + weightFactor); // 相对误差和绝对误差的组合 double yMax = qMax(qAbs(y1), qAbs(y2)); yMax = qMax(1e-12, yMax); // 防止除零 double relativeError = qAbs(y1 - y2) / yMax; double absoluteError = qAbs(y1 - y2); // 限制误差值 relativeError = qMin(1e6, relativeError); absoluteError = qMin(1e6, absoluteError); // 误差组合:相对误差为主,绝对误差为辅 double pointError = 0.7 * relativeError + 0.3 * absoluteError; if(isFiniteNumber(pointError) && pointError < 1e10) { totalError += weight * pointError * pointError; totalWeight += weight; validPoints++; } } if(totalWeight > 0 && validPoints > 0) { double result = sqrt(totalError / totalWeight); // 最终检查 if(!isFiniteNumber(result)) { DEBUG_OUT("Final error calculation produced infinite result"); return 1e10; } return qMin(1e9, result); // 限制最大误差值 } else { DEBUG_OUT(QString("No valid points for error calculation: validPoints=%1") .arg(validPoints)); return 1e10; } } // ==================== 工具方法 ==================== double nmCalculationAutoFitGA::random01() const { return static_cast(qrand()) / RAND_MAX; } double nmCalculationAutoFitGA::gaussianRandom(double mean, double stddev) const { static bool hasSpare = false; static double spare; if(hasSpare) { hasSpare = false; return spare * stddev + mean; } hasSpare = true; double u = qMax(random01(), 1e-12); double v = random01(); double mag = stddev * sqrt(-2.0 * log(u)); spare = mag * cos(2.0 * 3.14159265359 * v); return mag * sin(2.0 * 3.14159265359 * v) + mean; } int nmCalculationAutoFitGA::getEnabledParameterCount() const { int count = 0; for(int i = 0; i < m_parameterSelected.size(); ++i) { if(m_parameterSelected[i]) count++; } return count; } void nmCalculationAutoFitGA::saveOptimizationResult() { DEBUG_OUT(QString("GA optimization result: fitness=%1, evaluations=%2/%3") .arg(m_bestFitness, 0, 'e', 4) .arg(m_successfulEvaluations) .arg(m_totalEvaluations)); }