From 2a4f54cf333c011c7dce75ed233c47be4e582abd Mon Sep 17 00:00:00 2001 From: Yi Zhang Date: Mon, 13 Jan 2025 12:08:35 +0800 Subject: [PATCH] tmp --- lib/mesh/mesh.cpp | 10 ++ lib/mesh/mesh.h | 40 ++++-- lib/mesh/regular_grid.cpp | 120 +++++++++--------- lib/mesh/regular_grid.h | 257 +++++++++++++++++++++++++++++++++++--- 4 files changed, 332 insertions(+), 95 deletions(-) diff --git a/lib/mesh/mesh.cpp b/lib/mesh/mesh.cpp index 5ae0907..3fccb81 100644 --- a/lib/mesh/mesh.cpp +++ b/lib/mesh/mesh.cpp @@ -78,6 +78,16 @@ gctl::meshdata &gctl::base_mesh::get_data(std::string datname) return datalist_[data_index(datname)]; } +const gctl::meshdata &gctl::base_mesh::get_data(std::string datname) const +{ + return datalist_[data_index(datname)]; +} + +gctl::meshdata *gctl::base_mesh::get_data_ptr(std::string datname) +{ + return &datalist_[data_index(datname)]; +} + void gctl::base_mesh::remove_data(std::string datname) { int idx = data_index(datname); diff --git a/lib/mesh/mesh.h b/lib/mesh/mesh.h index 9347c74..175a464 100644 --- a/lib/mesh/mesh.h +++ b/lib/mesh/mesh.h @@ -151,7 +151,7 @@ namespace gctl void set_meshinfo(std::string in_info); /** - * @brief 添加网格数据。 + * @brief 添加网格数据，并返回数据对象的引用。 * * @param in_loctype 网格数据类型（块体或顶点数据） * @param in_valtype 网格数据值类型 @@ -160,13 +160,13 @@ namespace gctl * @param if_output 数据是否可用于输出 * @param nan_val 数据无效值标记 * - * @return 数据对象 + * @return 数据对象的引用（非常量的） */ meshdata &add_data(mesh_data_type_e in_loctype, mesh_data_value_e in_valtype, std::string name, double init_val, bool if_output = true, double nan_val = GCTL_BDL_MAX); /** - * @brief 添加网格数据。 + * @brief 添加网格数据，并返回数据对象的引用。 * * @param in_loctype 网格数据类型（块体或顶点数据） * @param name 数据名称 @@ -174,13 +174,13 @@ namespace gctl * @param if_output 数据是否可用于输出 * @param nan_val 数据无效值标记 * - * @return 数据对象 + * @return 数据对象的引用（非常量的） */ meshdata &add_data(mesh_data_type_e in_loctype, std::string name, const array &init_arr, bool if_output = true, double nan_val = GCTL_BDL_MAX); /** - * @brief 添加网格数据。 + * @brief 添加网格数据，并返回数据对象的引用。 * * @param in_loctype 网格数据类型（块体或顶点数据） * @param name 数据名称 @@ -188,13 +188,13 @@ namespace gctl * @param if_output 数据是否可用于输出 * @param nan_val 数据无效值标记 * - * @return 数据对象 + * @return 数据对象的引用（非常量的） */ meshdata &add_data(mesh_data_type_e in_loctype, std::string name, const array &init_arr, bool if_output = true, double nan_val = GCTL_BDL_MAX); /** - * @brief 添加网格数据。 + * @brief 添加网格数据，并返回数据对象的引用。 * * @param in_loctype 网格数据类型（块体或顶点数据） * @param name 数据名称 @@ -202,19 +202,35 @@ namespace gctl * @param if_output 数据是否可用于输出 * @param nan_val 数据无效值标记 * - * @return 数据对象 + * @return 数据对象的引用（非常量的） */ meshdata &add_data(mesh_data_type_e in_loctype, std::string name, const array &init_arr, bool if_output = true, double nan_val = GCTL_BDL_MAX); /** - * @brief 返回对应名称的数据对象指针。 + * @brief 返回对应名称的数据对象引用。 * * @param datname 数据名 - * @return 数据指针 + * @return 数据对象的引用（非常量的） */ meshdata &get_data(std::string datname); + /** + * @brief 返回对应名称的数据对象常量引用。 + * + * @param datname 数据名 + * @return 数据对象的引用（常量的） + */ + const meshdata &get_data(std::string datname) const; + + /** + * @brief 返回对应名称的数据对象指针。 + * + * @param datname 数据名 + * @return 数据对象的指针 + */ + meshdata *get_data_ptr(std::string datname); + /** * @brief 删除对应名称的数据对象。 * @@ -242,7 +258,7 @@ namespace gctl void save_gmsh_withdata(std::string filename, std::string datname, output_type_e out_mode, index_packed_e packed = Packed); /** - * @brief （基类虚函数）保存网格到Gmsh文件。 + * @brief （基类纯虚函数）保存网格到Gmsh文件。 * * @param filename 文件名 * @param packed 索引是否为紧凑（从0开始） @@ -250,7 +266,7 @@ namespace gctl virtual void save_gmsh(std::string filename, index_packed_e packed = Packed) = 0; /** - * @brief 显示网格的维度信息 + * @brief （基类纯虚函数）显示网格的维度信息 * */ virtual void show_mesh_dimension(std::ostream &os) const = 0; diff --git a/lib/mesh/regular_grid.cpp b/lib/mesh/regular_grid.cpp index 0e1e5d8..c9bbc97 100644 --- a/lib/mesh/regular_grid.cpp +++ b/lib/mesh/regular_grid.cpp @@ -190,7 +190,7 @@ void gctl::regular_grid::clear_regular_grid() int gctl::regular_grid::view(std::string datname) { - meshdata data = get_data(datname); + meshdata &data = get_data(datname); plt.range(rg_xmin, rg_xmin + (rg_xnum - 1)*rg_dx, rg_ymin, rg_ymin + (rg_ynum - 1)*rg_dy); plt.demension(rg_xnum, rg_ynum); @@ -422,36 +422,34 @@ void gctl::regular_grid::save_netcdf_grid(std::string filename, mesh_data_type_e check_initiated(); // check to see if the mesh is not initiated bool if_append = false; - meshdata curr_data; for (size_t i = 0; i < datalist_.size(); i++) { - curr_data = datalist_[i]; - if (curr_data.loctype_ == d_type && d_type == NodeData && - curr_data.valtype_ == Scalar && curr_data.output_ok_) + if (datalist_[i].loctype_ == d_type && d_type == NodeData && + datalist_[i].valtype_ == Scalar && datalist_[i].output_ok_) { if (!if_append) { - gctl::save_netcdf_grid(filename, curr_data.datval_, rg_xnum, rg_ynum, rg_xmin, - rg_dx, rg_ymin, rg_dy, xname, yname, curr_data.name_); + gctl::save_netcdf_grid(filename, datalist_[i].datval_, rg_xnum, rg_ynum, rg_xmin, + rg_dx, rg_ymin, rg_dy, xname, yname, datalist_[i].name_); if_append = true; } else { - gctl::append_netcdf_grid(filename, curr_data.datval_, xname, yname, curr_data.name_); + gctl::append_netcdf_grid(filename, datalist_[i].datval_, xname, yname, datalist_[i].name_); } } - else if (curr_data.loctype_ == d_type && d_type == ElemData && - curr_data.valtype_ == Scalar && curr_data.output_ok_) + else if (datalist_[i].loctype_ == d_type && d_type == ElemData && + datalist_[i].valtype_ == Scalar && datalist_[i].output_ok_) { if (!if_append) { - gctl::save_netcdf_grid(filename, curr_data.datval_, rg_xnum-1, rg_ynum-1, rg_xmin+0.5*rg_dx, - rg_dx, rg_ymin+0.5*rg_dy, rg_dy, xname, yname, curr_data.name_); + gctl::save_netcdf_grid(filename, datalist_[i].datval_, rg_xnum-1, rg_ynum-1, rg_xmin+0.5*rg_dx, + rg_dx, rg_ymin+0.5*rg_dy, rg_dy, xname, yname, datalist_[i].name_); if_append = true; } else { - gctl::append_netcdf_grid(filename, curr_data.datval_, xname, yname, curr_data.name_); + gctl::append_netcdf_grid(filename, datalist_[i].datval_, xname, yname, datalist_[i].name_); } } } @@ -461,7 +459,7 @@ void gctl::regular_grid::save_netcdf_grid(std::string filename, mesh_data_type_e void gctl::regular_grid::save_netcdf_grid(std::string filename, std::string datname, std::string xname, std::string yname) { - meshdata curr_data = get_data(datname); + meshdata &curr_data = get_data(datname); if (!curr_data.output_ok_) { @@ -593,7 +591,7 @@ void gctl::regular_grid::load_surfer_grid(std::string filename, std::string datn void gctl::regular_grid::save_surfer_grid(std::string filename, std::string datname, surfer_file_type_e grid_type) { - meshdata curr_data = get_data(datname); + meshdata &curr_data = get_data(datname); if (!curr_data.output_ok_) { @@ -673,13 +671,10 @@ void gctl::regular_grid::save_gmsh(std::string filename, std::string datname, ou void gctl::regular_grid::save_text(std::string filename, const array &datname) { - mesh_data_type_e d_type; - meshdata curr_data = get_data(datname[0]); - d_type = curr_data.loctype_; - for (size_t i = 1; i < datname.size(); i++) + meshdata &curr_data = get_data(datname[0]); + for (size_t i = 1; i < datalist_.size(); i++) { - curr_data = get_data(datname[i]); - if (d_type != curr_data.loctype_) + if (curr_data.loctype_ != datalist_[i].loctype_) { throw std::runtime_error("[gctl::regular_grid] multiple data location types found."); } @@ -700,10 +695,9 @@ void gctl::regular_grid::save_text(std::string filename, const arrayx + elements[i].ur->x) << " " << 0.5*(elements[i].dl->y + elements[i].ur->y); - for (size_t n = 0; n < datname.size(); n++) + for (size_t n = 0; n < datalist_.size(); n++) { - curr_data = get_data(datname[n]); - ofile << " " << curr_data.datval_[i]; + ofile << " " << datalist_[n].datval_[i]; } ofile << "\n"; } @@ -767,7 +760,7 @@ void gctl::regular_grid::load_data_cloud(const array &in_posi, const a in_val.get(), in_posi.size(), search_xlen, search_ylen, search_deg); // create a new data object - meshdata curr_data = add_data(d_type, Scalar, datname, 0.0, true, GCTL_BDL_MAX); + meshdata &curr_data = add_data(d_type, Scalar, datname, 0.0, true, GCTL_BDL_MAX); for (int i = 0; i < curr_data.datval_.size(); i++) { curr_data.datval_[i] = inte_data->at(i); @@ -785,7 +778,7 @@ void gctl::regular_grid::extract_points(std::string datname, const array parser; if (!parser.compile(expression_str, expression)) throw std::runtime_error("[gctl::regular_grid] Fail to compile the math expression."); - array datas(var.size()); + array datas_ptr(var.size()); mesh_data_value_e val_type1, val_type2; mesh_data_type_e data_type1, data_type2; // Remeber we put the output at the first place. for (size_t i = 1; i < var.size(); i++) { - datas[i] = get_data(data_list[i]); + datas_ptr[i] = get_data_ptr(data_list[i]); } - val_type1 = datas[1].valtype_; - data_type1= datas[1].loctype_; + val_type1 = datas_ptr[1]->valtype_; + data_type1= datas_ptr[1]->loctype_; for (size_t i = 2; i < var.size(); i++) { - val_type2 = datas[i].valtype_; - data_type2= datas[i].loctype_; + val_type2 = datas_ptr[i]->valtype_; + data_type2= datas_ptr[i]->loctype_; if (val_type1 != val_type2 || data_type1 != data_type2) { throw std::runtime_error("[gctl::regular_grid] Can't preform the calculation on selected data."); @@ -1101,22 +1091,24 @@ void gctl::regular_grid::calculator(std::string expression_str, const arrayvaltype_; + data_type2= datas_ptr[0]->loctype_; if (val_type1 != val_type2 || data_type1 != data_type2) { throw std::runtime_error("[gctl::regular_grid] Can't preform the calculation on selected data."); } - for (size_t i = 0; i < datas[0].datval_.size(); i++) + for (size_t i = 0; i < datas_ptr[0]->datval_.size(); i++) { for (size_t v = 0; v < var.size(); v++) { - var[v] = datas[v].datval_[i]; + var[v] = datas_ptr[v]->datval_[i]; } - datas[0].datval_[i] = expression.value(); + datas_ptr[0]->datval_[i] = expression.value(); } return; } @@ -1127,7 +1119,7 @@ void gctl::regular_grid::calculator(std::string expression_str, const arrayJ), 0.0, true, GCTL_BDL_MAX); + meshdata &new_data = add_data(data_type, value_type, datname+"_A"+std::to_string(wt->J), 0.0, true, GCTL_BDL_MAX); idwt2(wt, wavecoeffs, new_data.datval_.get()); int coeff_idx; diff --git a/lib/mesh/regular_grid.h b/lib/mesh/regular_grid.h index ccd3392..b575609 100644 --- a/lib/mesh/regular_grid.h +++ b/lib/mesh/regular_grid.h @@ -48,91 +48,310 @@ namespace gctl /** * 以下为mesh类型的虚函数实现 */ + + /** + * @brief （基类纯虚函数的具体实现）保存网格到Gmsh文件。 + * + * @param filename 文件名 + * @param packed 索引是否为紧凑（从0开始） + */ + void save_gmsh(std::string filename, index_packed_e packed = Packed); - void init(std::string in_name, std::string in_info, int xnum, int ynum, - double xmin, double ymin, double dx, double dy); - + /** + * @brief （基类纯虚函数的具体实现）显示网格的维度信息 + */ void show_mesh_dimension(std::ostream &os) const; + /** + * @brief （基类纯虚函数的具体实现）读入二进制网格文件 + * + * @param filename 文件名 + */ void load_binary(std::string filename); + + /** + * @brief （基类纯虚函数的具体实现）写入二进制网格文件 + * + * @param filename 文件名 + */ void save_binary(std::string filename); /** * 以下为regular_grid的专有函数 */ - regular_grid(); + regular_grid(); ///< 构造函数 + + /** + * @brief 构造函数 + * + * @param in_name 网格名称 + * @param in_info 网格说明 + * @param xnum 网格x方向单元数量 + * @param ynum 网格y方向单元数量 + * @param xmin 网格x方向最小值 + * @param ymin 网格y方向最小值 + * @param dx 网格x方向单元长度 + * @param dy 网格y方向单元长度 + */ regular_grid(std::string in_name, std::string in_info, int xnum, int ynum, double xmin, double ymin, double dx, double dy); - virtual ~regular_grid(); + virtual ~regular_grid(); ///< 析构函数 + + /** + * @brief 初始化规则网格 + * + * @param in_name 网格名称 + * @param in_info 网格说明 + * @param xnum 网格x方向单元数量 + * @param ynum 网格y方向单元数量 + * @param xmin 网格x方向最小值 + * @param ymin 网格y方向最小值 + * @param dx 网格x方向单元长度 + * @param dy 网格y方向单元长度 + */ + void init(std::string in_name, std::string in_info, int xnum, int ynum, + double xmin, double ymin, double dx, double dy); + + /** + * @brief 清除规则网格 + */ void clear_regular_grid(); + + /** + * @brief 使用mathgl显示规则网格 + * + * @param datname 数据名称 + * @return 状态代码 + */ int view(std::string datname); + + /** + * @brief 使用gmt绘制规则网格 + * + * @param datname 数据名称 + */ void plot(std::string datname); - int get_xdim() const; - int get_ydim() const; - double get_xmin() const; - double get_ymin() const; - double get_dx() const; - double get_dy() const; + int get_xdim() const; ///< 获取规则网格的x方向单元数量 + int get_ydim() const; ///< 获取规则网格的y方向单元数量 + double get_xmin() const; ///< 获取规则网格的x方向最小值 + double get_ymin() const; ///< 获取规则网格的y方向最小值 + double get_dx() const; ///< 获取规则网格的x方向单元长度 + double get_dy() const; ///< 获取规则网格的y方向单元长度 #ifdef GCTL_NETCDF + /** + * @brief 从netcdf文件中读取规则网格数据 + * + * @param filename 文件名 + * @param d_type 数据类型 + * @param xname x坐标变量名 + * @param yname y坐标变量名 + * @param zname z坐标变量名 + */ void load_netcdf_grid(std::string filename, mesh_data_type_e d_type, std::string xname = "x", std::string yname = "y", std::string zname = "null"); + + /** + * @brief 保存规则网格数据到netcdf文件中 + * + * @param filename 文件名 + * @param d_type 数据类型 + * @param xname x坐标变量名 + * @param yname y坐标变量名 + */ void save_netcdf_grid(std::string filename, mesh_data_type_e d_type, std::string xname = "x", std::string yname = "y"); + + /** + * @brief 保存规则网格数据到netcdf文件中 + * + * @param filename 文件名 + * @param datname 数据名称 + * @param xname x坐标变量名 + * @param yname y坐标变量名 + */ void save_netcdf_grid(std::string filename, std::string datname, std::string xname = "x", std::string yname = "y"); #endif // GCTL_NETCDF + /** + * @brief 从surfer文件中读取规则网格数据 + * + * @param filename 文件名 + * @param datname 数据名称 + * @param d_type 数据类型 + * @param grid_type 网格类型 + */ void load_surfer_grid(std::string filename, std::string datname, mesh_data_type_e d_type, surfer_file_type_e grid_type = Surfer7Grid); + + /** + * @brief 保存规则网格数据到surfer文件中 + * + * @param filename 文件名 + * @param datname 数据名称 + * @param grid_type 网格类型 + */ void save_surfer_grid(std::string filename, std::string datname, surfer_file_type_e grid_type = Surfer7Grid); - void save_gmsh(std::string filename, index_packed_e packed = Packed); + /** + * @brief 保存规则网格数据到Gmsh文件中 + * + * @param filename 文件名 + * @param out_mode 输出模式 + * @param packed 索引是否为紧凑（从0开始） + */ void save_gmsh(std::string filename, output_type_e out_mode, index_packed_e packed = Packed); + + /** + * @brief 保存规则网格数据到Gmsh文件中 + * + * @param filename 文件名 + * @param datname 数据名称 + * @param out_mode 输出模式 + * @param packed 索引是否为紧凑（从0开始） + */ void save_gmsh(std::string filename, std::string datname, output_type_e out_mode, index_packed_e packed = Packed); + /** + * @brief 保存规则网格数据到文本文件中 + * + * @param filename 文件名 + * @param datname 数据名称 + */ void save_text(std::string filename, const array &datname); + /** + * @brief 从点云数据中提取规则网格数据 + * + * @param in_posi 点云位置数组 + * @param in_val 点云数据数组 + * @param search_xlen 搜索范围x方向长度 + * @param search_ylen 搜索范围y方向长度 + * @param search_deg 搜索范围角度 + * @param datname 新建数据名称 + * @param d_type 新建数据类型 + */ void load_data_cloud(const array &in_posi, const array &in_val, double search_xlen, double search_ylen, double search_deg, std::string datname, mesh_data_type_e d_type = NodeData); + /** + * @brief 从规则网格数据中提取点云数据 + * + * @param datname 数据名称 + * @param in_posi 点云位置数组 + * @param out_val 输出的点云数据数组 + */ void extract_points(std::string datname, const array &in_posi, array &out_val); + + /** + * @brief 从规则网格数据中提取剖面数据 + * + * @param datname 数据名称 + * @param start_p 剖面起点 + * @param end_p 剖面终点 + * @param size_p 剖面点数 + * @param out_posi 输出的剖面位置数组 + * @param out_val 输出的剖面数据数组 + */ void extract_profile(std::string datname, const point2dc &start_p, const point2dc &end_p, int size_p, array &out_posi, array &out_val); + /** + * @brief 计算规则网格数据的梯度数据 + * + * @param datname 数据名称 + * @param gradname 梯度数据名称 + * @param d_type 梯度数据类型 + * @param order 梯度阶数 + */ void gradient(std::string datname, std::string gradname, gradient_type_e d_type, int order = 1); + + /** + * @brief 计算两个规则网格数据的和 + * + * @param newname 新数据名称 + * @param datname 数据名称 + * @param datname2 数据名称2 + */ void sum(std::string newname, std::string datname, std::string datname2); + + /** + * @brief 计算两个规则网格数据的差 + * + * @param newname 新数据名称 + * @param datname 数据名称 + * @param datname2 数据名称2 + */ void diff(std::string newname, std::string datname, std::string datname2); + + /** + * @brief 计算规则网格数据的布尔数据 + * + * @param newname 新数据名称 + * @param datname 数据名称 + * @param maskname 掩码数据名称 + * @param reverse 是否反转掩码 + */ void boolean(std::string newname, std::string datname, std::string maskname, bool reverse = false); #ifdef GCTL_MESH_EXPRTK + + /** + * @brief 计算规则网格数据的函数值，其中f为节点位置（x，y）的函数 + * + * @param expression_str 表达式字符串 + * @param newname 新数据名称 + * @param x_str x变量名称 + * @param y_str y变量名称 + * @param f_str f变量名称 + */ void function(std::string expression_str, std::string newname, std::string x_str = "x", std::string y_str = "y", std::string f_str = "f"); + + /** + * @brief 计算规则网格数据的函数值，计算不同网格数据合成的函数值 + * + * @param expression_str 表达式字符串 + * @param var_list 变量名称列表 + * @param data_list 数据名称列表 + */ void calculator(std::string expression_str, const array &var_list, const array &data_list); + #endif // GCTL_MESH_EXPRTK #ifdef GCTL_MESH_WAVELIB + + /** + * @brief 计算规则网格数据的小波数据 + * + * @param datname 数据名称 + * @param wavename 小波数据名称 + * @param order 小波阶数 + * @param summary 是否输出小波处理摘要 + */ void wavelet(std::string datname, std::string wavename, int order, bool summary = true); + #endif // GCTL_MESH_WAVELIB protected: - int rg_xnum, rg_ynum; - double rg_xmin, rg_ymin; - double rg_dx, rg_dy; + int rg_xnum, rg_ynum; ///< 规则网格的x方向单元数量和y方向单元数量 + double rg_xmin, rg_ymin; ///< 规则网格的x方向最小值和y方向最小值 + double rg_dx, rg_dy; ///< 规则网格的x方向单元长度和y方向单元长度 #ifdef GCTL_GRAPHIC_MATHGL - mathgl_dens plt; + mathgl_dens plt; ///< mathgl绘图对象 #endif // GCTL_GRAPHIC_MATHGL #ifdef GCTL_GRAPHIC_GMT - gmt_JX_single pic; + gmt_JX_single pic; ///< gmt绘图对象 #endif // GCTL_GRAPHIC_GMT - array nodes; - array elements; + array nodes; ///< 规则网格的节点数组 + array elements; ///< 规则网格的单元数组 }; }