当前位置： 首页 > news >正文

C++实现MATLAB矩阵计算程序

news 2025/8/10 5:29:04

使用GCC C++代码，整合了MATIO库、Eigen矩阵计算、多线程处理及日志系统，支持解析M代码中的矩阵运算并输出多种格式结果：

#include <iostream>
#include <fstream>
#include <vector>
#include <map>
#include <regex>
#include <thread>
#include <future>
#include <queue>
#include <mutex>
#include <condition_variable>
#include <filesystem>
#include <sstream>
#include <cctype>
#include <stdexcept>
#include <chrono>
#include <<iomanip>
#include <backward.hpp>
#include <matio.h>
#include <Eigen/Dense>
#include <spdlog/spdlog.h>
#include <spdlog/sinks/stdout_color_sinks.h>
#include <spdlog/sinks/basic_file_sink.h>namespace fs = std::filesystem;
namespace backward { backward::SignalHandling sh; }using Eigen::MatrixXd;
using Eigen::VectorXd;
using EigenMap = std::map<std::string, MatrixXd>;// 全局日志器初始化
std::shared_ptr<spdlog::logger> init_logger() {auto console_sink = std::make_shared<spdlog::sinks::stdout_color_sink_mt>();console_sink->set_pattern("%^[%Y-%m-%d %H:%M:%S.%e] [%l]%$ %v");auto file_sink = std::make_shared<spdlog::sinks::basic_file_sink_mt>("matrix_calc.log", true);file_sink->set_pattern("[%Y-%m-%d %H:%M:%S.%e] [%l] %v");std::vector<spdlog::sink_ptr> sinks{console_sink, file_sink};auto logger = std::make_shared<spdlog::logger>("mat_exec", sinks.begin(), sinks.end());logger->set_level(spdlog::level::debug);spdlog::register_logger(logger);return logger;
}auto logger = init_logger();// 异常处理宏定义
#define THROW_EXCEPTION(msg) \do { \std::stringstream ss; \ss << msg << " [文件: " << __FILE__ << " 行: " << __LINE__ << "]"; \logger->error(ss.str()); \throw std::runtime_error(ss.str()); \} while(0)// 矩阵运算工具类
namespace MatrixOps {MatrixXd add(const MatrixXd& a, const MatrixXd& b) {if (a.rows() != b.rows() || a.cols() != b.cols())THROW_EXCEPTION("矩阵加法维度不匹配");return a + b;}MatrixXd multiply(const MatrixXd& a, const MatrixXd& b) {if (a.cols() != b.rows())THROW_EXCEPTION("矩阵乘法维度不匹配");return a * b;}MatrixXd inv(const MatrixXd& m) {if (m.rows() != m.cols()) THROW_EXCEPTION("求逆要求方阵");return m.inverse();}MatrixXd pinv(const MatrixXd& m, double tol = 1e-6) {Eigen::JacobiSVD<MatrixXd> svd(m, Eigen::ComputeThinU | Eigen::ComputeThinV);VectorXd sv = svd.singularValues();MatrixXd S_inv = MatrixXd::Zero(m.cols(), m.rows());for (int i = 0; i < sv.size(); ++i) {if (sv(i) > tol) S_inv(i, i) = 1.0 / sv(i);}return svd.matrixV() * S_inv * svd.matrixU().transpose();}MatrixXd eig(const MatrixXd& m) {if (m.rows() != m.cols()) THROW_EXCEPTION("特征值要求方阵");Eigen::EigenSolver<MatrixXd> solver(m);return solver.eigenvalues().real();}MatrixXd svd(const MatrixXd& m) {Eigen::JacobiSVD<MatrixXd> svd(m, Eigen::ComputeThinU | Eigen::ComputeThinV);return svd.singularValues();}MatrixXd chol(const MatrixXd& m) {if (m.rows() != m.cols()) THROW_EXCEPTION("Cholesky分解要求方阵");Eigen::LLT<MatrixXd> llt(m);if (llt.info() == Eigen::NumericalIssue)THROW_EXCEPTION("矩阵非正定，无法Cholesky分解");return llt.matrixL();}
}// MAT文件处理类
class MatFileHandler {
public:static EigenMap load(const std::string& path) {if (!fs::exists(path)) THROW_EXCEPTION("文件不存在: " + path);mat_t* mat = Mat_Open(path.c_str(), MAT_ACC_RDONLY);if (!mat) THROW_EXCEPTION("无法打开MAT文件: " + path);EigenMap variables;matvar_t* var;while ((var = Mat_VarReadNext(mat)) != nullptr) {if (var->data_type == MAT_T_DOUBLE && var->rank == 2) {MatrixXd matrix(var->dims[0], var->dims[1]);memcpy(matrix.data(), var->data, var->nbytes);variables[var->name] = matrix;logger->debug("加载变量: {} [{}x{}]", var->name, var->dims[0], var->dims[1]);}Mat_VarFree(var);}Mat_Close(mat);return variables;}static void save(const std::string& path, const EigenMap& variables) {mat_t* mat = Mat_CreateVer(path.c_str(), nullptr, MAT_FT_MAT5);if (!mat) THROW_EXCEPTION("无法创建MAT文件: " + path);for (const auto& [name, matrix] : variables) {size_t dims[2] = {static_cast<size_t>(matrix.rows()),static_cast<size_t>(matrix.cols())};matvar_t* var = Mat_VarCreate(name.c_str(), MAT_C_DOUBLE, MAT_T_DOUBLE, 2, dims,(void*)matrix.data(), MAT_F_DONT_COPY_DATA);if (!var || Mat_VarWrite(mat, var, MAT_COMPRESSION_NONE) != 0) {Mat_VarFree(var);Mat_Close(mat);THROW_EXCEPTION("写入变量失败: " + name);}Mat_VarFree(var);logger->debug("保存变量: {} [{}x{}]", name, dims[0], dims[1]);}Mat_Close(mat);}
};// 其他格式输出类
class FormatExporter {
public:static void save_csv(const std::string& path, const EigenMap& variables) {std::ofstream file(path);if (!file.is_open()) THROW_EXCEPTION("无法创建CSV文件: " + path);for (const auto& [name, matrix] : variables) {file << "=== " << name << " (" << matrix.rows() << "x" << matrix.cols() << ") ===" << std::endl;for (int i = 0; i < matrix.rows(); ++i) {for (int j = 0; j < matrix.cols(); ++j) {file << matrix(i, j);if (j < matrix.cols() - 1) file << ",";}file << std::endl;}file << std::endl;}}// 可扩展实现Excel和Parquet格式static void save_excel(const std::string& path, const EigenMap& variables) {logger->warn("Excel格式暂未实现: {}", path);}static void save_parquet(const std::string& path, const EigenMap& variables) {logger->warn("Parquet格式暂未实现: {}", path);}
};// M代码解析器
class MCodeParser {
public:explicit MCodeParser(const std::string& script_path) {parse_script(script_path);}void execute(EigenMap& variables) const {for (size_t i = 0; i < instructions.size(); ++i) {try {const auto& [lhs, rhs] = instructions[i];MatrixXd result = evaluate_expression(rhs, variables);variables[lhs] = result;logger->info("执行成功: {} = [{}x{}矩阵] (行 {})", lhs, result.rows(), result.cols(), i + 1);} catch (const std::exception& e) {logger->error("执行失败 (行 {}): {}\n{}", i + 1, instructions[i].first + " = " + instructions[i].second, e.what());}}}private:std::vector<std::pair<std::string, std::string>> instructions;void parse_script(const std::string& path) {std::ifstream file(path);if (!file.is_open()) THROW_EXCEPTION("无法打开M文件: " + path);std::string line;int line_num = 0;while (std::getline(file, line)) {line_num++;// 移除注释和空白size_t comment_pos = line.find('%');if (comment_pos != std::string::npos) {line = line.substr(0, comment_pos);}line = trim(line);if (line.empty()) continue;// 解析赋值语句size_t eq_pos = line.find('=');if (eq_pos == std::string::npos) {logger->warn("忽略无效语句 (行 {}): {}", line_num, line);continue;}std::string lhs = trim(line.substr(0, eq_pos));std::string rhs = trim(line.substr(eq_pos + 1));if (lhs.empty() || rhs.empty()) {logger->warn("忽略不完整赋值 (行 {}): {}", line_num, line);continue;}instructions.emplace_back(lhs, rhs);}logger->info("解析完成: {} 条有效指令", instructions.size());}static std::string trim(const std::string& s) {auto start = std::find_if_not(s.begin(), s.end(), ::isspace);auto end = std::find_if_not(s.rbegin(), s.rend(), ::isspace).base();return (start < end) ? std::string(start, end) : "";}MatrixXd evaluate_expression(const std::string& expr, const EigenMap& vars) const {// 处理函数调用static std::regex func_re(R"((\w+)\s*\(\s*([^)]*)\s*\))");std::smatch match;if (std::regex_match(expr, match, func_re)) {std::string func = match[1];std::vector<std::string> args = split_arguments(match[2]);return call_function(func, args, vars);}// 处理二元运算static std::vector<std::pair<std::regex, std::function<MatrixXd(const MatrixXd&, const MatrixXd&)>>> ops = {{std::regex(R"(^\s*(\w+)\s*\+\s*(\w+)\s*$)"), MatrixOps::add},{std::regex(R"(^\s*(\w+)\s*\*\s*(\w+)\s*$)"), MatrixOps::multiply}};for (const auto& [op_re, op_func] : ops) {if (std::regex_match(expr, match, op_re)) {std::string a_name = match[1];std::string b_name = match[2];if (!vars.count(a_name)) THROW_EXCEPTION("未定义变量: " + a_name);if (!vars.count(b_name)) THROW_EXCEPTION("未定义变量: " + b_name);return op_func(vars.at(a_name), vars.at(b_name));}}// 处理变量引用if (vars.count(expr)) {return vars.at(expr);}THROW_EXCEPTION("无法解析表达式: " + expr);}static std::vector<std::string> split_arguments(const std::string& args_str) {std::vector<std::string> args;std::istringstream iss(args_str);std::string arg;while (std::getline(iss, arg, ',')) {arg = trim(arg);if (!arg.empty()) args.push_back(arg);}return args;}MatrixXd call_function(const std::string& name, const std::vector<std::string>& args, const EigenMap& vars) const {if (args.empty()) THROW_EXCEPTION("函数调用缺少参数: " + name);// 检查所有参数是否存在for (const auto& arg : args) {if (!vars.count(arg)) THROW_EXCEPTION("函数参数未定义: " + arg);}// 调用对应矩阵函数if (name == "inv" && args.size() == 1) return MatrixOps::inv(vars.at(args[0]));if (name == "pinv" && args.size() == 1) return MatrixOps::pinv(vars.at(args[0]));if (name == "eig" && args.size() == 1) return MatrixOps::eig(vars.at(args[0]));if (name == "svd" && args.size() == 1) return MatrixOps::svd(vars.at(args[0]));if (name == "chol" && args.size() == 1) return MatrixOps::chol(vars.at(args[0]));THROW_EXCEPTION("不支持的函数或参数数量错误: " + name);}
};// 线程池处理器
class ThreadPool {
public:explicit ThreadPool(size_t threads = std::thread::hardware_concurrency()) : stop(false) {for (size_t i = 0; i < threads; ++i) {workers.emplace_back([this] {while (true) {std::function<void()> task;{std::unique_lock<std::mutex> lock(queue_mutex);condition.wait(lock, [this] { return stop || !tasks.empty(); });if (stop && tasks.empty()) return;task = std::move(tasks.front());tasks.pop();}task();}});}logger->info("线程池初始化完成，线程数: {}", workers.size());}template<class F>void enqueue(F&& f) {{std::unique_lock<std::mutex> lock(queue_mutex);if (stop) throw std::runtime_error("线程池已停止，无法添加任务");tasks.emplace(std::forward<F>(f));}condition.notify_one();}~ThreadPool() {{std::unique_lock<std::mutex> lock(queue_mutex);stop = true;}condition.notify_all();for (std::thread& worker : workers) {if (worker.joinable()) worker.join();}}private:std::vector<std::thread> workers;std::queue<std::function<void()>> tasks;std::mutex queue_mutex;std::condition_variable condition;bool stop;
};// 文件处理函数
void process_file(const std::string& mat_path, const std::string& m_script,const std::vector<std::string>& output_formats) {try {logger->info("开始处理: {}", mat_path);// 加载MAT文件变量EigenMap variables = MatFileHandler::load(mat_path);logger->info("加载变量数量: {}", variables.size());// 解析并执行M代码MCodeParser parser(m_script);parser.execute(variables);// 生成输出文件名（移除原扩展名）std::string base_name = fs::path(mat_path).stem().string();// 保存结果到指定格式for (const auto& format : output_formats) {std::string output_path = base_name + "." + format;if (format == "mat") {MatFileHandler::save(output_path, variables);} else if (format == "csv") {FormatExporter::save_csv(output_path, variables);} else if (format == "xlsx") {FormatExporter::save_excel(output_path, variables);} else if (format == "parquet") {FormatExporter::save_parquet(output_path, variables);} else {logger->warn("不支持的输出格式: {}", format);continue;}logger->info("已保存到: {}", output_path);}logger->info("处理完成: {}", mat_path);
} catch (const std::exception& e) {
logger->error("处理文件失败: {}\n{}", mat_path, e.what());
}
}// 目录处理函数
void process_directory(const std::string& dir_path,
const std::string& m_script,
const std::vectorstd::string& output_formats,
ThreadPool& pool) {
if (!fs::is_directory(dir_path)) {
logger->error("不是有效目录: {}", dir_path);
return;
}for (const auto& entry : fs::directory_iterator(dir_path)) {
if (entry.is_regular_file() && entry.path().extension() == ".mat") {
pool.enqueue([=] {
process_file(entry.path().string(), m_script, output_formats);
});
}
}
logger->info("目录任务已全部加入队列: {}", dir_path);
}// 命令行参数解析
struct CmdArgs {
std::string input_path;
std::string script_path;
std::vectorstd::string output_formats;
size_t threads = std::thread::hardware_concurrency();
};CmdArgs parse_args(int argc, char* argv[]) {
if (argc < 4) {
std::cerr << "用法: " << argv[0]
<< " <输入路径(.mat或目录)> <M脚本路径> <输出格式(用逗号分隔)> [线程数]\n"
<< "示例: " << argv[0]
<< " data.mat calc.m mat,csv 4\n"
<< "支持格式: mat, csv, xlsx(待实现), parquet(待实现)\n";
exit(1);
}CmdArgs args;
args.input_path = argv[1];
args.script_path = argv[2];// 解析输出格式
std::istringstream iss(argv[3]);
std::string format;
while (std::getline(iss, format, ',')) {
if (!format.empty()) args.output_formats.push_back(format);
}// 解析线程数
if (argc >= 5) {
try {
args.threads = std::stoul(argv[4]);
if (args.threads == 0) args.threads = 1;
} catch (...) {
logger->warn("无效线程数，使用默认值: {}", args.threads);
}
}// 验证输入
if (!fs::exists(args.input_path)) {
THROW_EXCEPTION("输入路径不存在: " + args.input_path);
}
if (!fs::exists(args.script_path) || fs::path(args.script_path).extension() != ".m") {
THROW_EXCEPTION("无效M脚本文件: " + args.script_path);
}
if (args.output_formats.empty()) {
THROW_EXCEPTION("未指定输出格式");
}return args;
}int main(int argc, char* argv[]) {
try {
logger->info("程序启动");// 解析命令行参数
CmdArgs args = parse_args(argc, argv);
logger->info("输入路径: {}", args.input_path);
logger->info("M脚本: {}", args.script_path);
logger->info("输出格式: {}", fmt::join(args.output_formats, ","));
logger->info("线程数: {}", args.threads);// 初始化线程池
ThreadPool pool(args.threads);// 处理输入（文件或目录）
if (fs::is_regular_file(args.input_path) &&
fs::path(args.input_path).extension() == ".mat") {
// 处理单个MAT文件
pool.enqueue([&] {
process_file(args.input_path, args.script_path, args.output_formats);
});
} else if (fs::is_directory(args.input_path)) {
// 处理目录下所有MAT文件
process_directory(args.input_path, args.script_path, args.output_formats, pool);
} else {
THROW_EXCEPTION("输入必须是.mat文件或目录: " + args.input_path);
}// 等待所有任务完成（通过线程池析构实现）
logger->info("等待所有任务完成...");
} catch (const std::exception& e) {
logger->critical("程序异常终止: {}", e.what());
return 1;
}logger->info("程序正常结束");
return 0;
}

代码说明

  1. 核心功能实现：

    • 使用MATIO库完成.mat文件的读写操作
    • 基于Eigen库实现矩阵的基本运算（加、乘、逆、伪逆等）
    • 支持解析M代码中的简单赋值语句和函数调用
    • 通过线程池实现多文件并行处理

  2. 日志与异常处理：

    • 集成spdlog实现带时间戳的日志记录（控制台+文件）
    • 详细的错误信息包含文件名、行号和调用堆栈
    • 异常处理覆盖文件操作、矩阵运算和参数解析等场景

  3. 扩展性设计：

    • 采用模块化结构，便于添加新的矩阵运算函数
    • 预留Excel和Parquet格式输出接口
    • 命令行参数支持灵活配置输入输出和线程数

  4. 使用方法：

    • 编译时需链接matio、eigen、spdlog和线程库
    • 运行时通过命令行指定输入（文件/目录）、M脚本、输出格式和线程数
http://www.lryc.cn/news/614857.html

相关文章：

  • 【Redis7.x】docker配置主从+sentinel监控遇到的问题与解决
  • Debian 系统更新命令
  • PDF 转 HTML API 数据接口
  • 免费PDF编辑软件 pdf24-creator 及其安装包
  • 力扣-74.搜索二维矩阵
  • MyBatis联合查询 - 注解篇
  • 【洛谷题单】--分支结构(三)
  • JAVA基础-使用BIO / NIO实现聊天室功能
  • 一文详解 C++ 继承体系
  • AI_RAG
  • 本地连接跳板机
  • 10. 怎么实现深拷贝？
  • ABP VNext + Apache Kafka Exactly-Once 语义：金融级消息一致性实战
  • VSCode添加Python、Java注释技巧、模板
  • 笔试——Day33
  • java web项目入门了解
  • 微信原生小程序 Timeline 组件实现
  • 在Word和WPS文字中快速拆分、合并表格
  • JavaWeb03——javascript基础语法
  • C++-AVL树
  • 微软将于 10 月停止混合 Exchange 中的共享 EWS 访问
  • SOLi-LABS Page-3 (Stacked injections) --39-53关
  • 使用 Vuepress + GitHub Pages 搭建项目文档（2）- 使用 GitHub Actions 工作流自动部署
  • 如何解决 Vue 项目启动时出现的 “No such module: http_parser” 错误问题
  • 2G内存的服务器用宝塔安装php的fileinfo拓展时总是卡死无法安装成功的解决办法
  • 企业级web应用服务器TOMCAT入门详解
  • kettle插件-kettle MinIO插件，轻松解决文件上传到MinIO服务器
  • 解决本地连接服务器ollama的错误
  • 大语言模型提示工程与应用：大语言模型对抗性提示安全防御指南
  • LLVM编译器入门