Compare commits

...

7 Commits
main ... dev_yi

Author SHA1 Message Date
cf8b98af4c tmp 2025-04-23 21:55:30 +08:00
3cabeba586 update 2.0 2025-04-23 14:05:58 +08:00
00c6a64d83 tmp 2025-04-21 13:18:40 +08:00
807e472dee tmp 2025-04-04 22:36:13 +08:00
b7508642a7 tmp 2025-02-09 19:49:16 +08:00
8492445240 tmp 2025-02-09 16:12:06 +08:00
963ca494ae tmp 2025-02-03 15:44:09 +08:00
33 changed files with 3372 additions and 140 deletions

View File

@ -1,20 +1,18 @@
cmake_minimum_required(VERSION 3.15.2)
#
project(GCTL_MESH VERSION 1.0)
project(GCTL_MESH VERSION 2.0)
#
include(CMakePackageConfigHelpers)
# ExprTKmacOS 15.4
add_compile_options(-Wno-missing-template-arg-list-after-template-kw)
message(STATUS "Platform: " ${CMAKE_HOST_SYSTEM_NAME})
message(STATUS "Install prefix: " ${CMAKE_INSTALL_PREFIX})
message(STATUS "Processor: " ${CMAKE_HOST_SYSTEM_PROCESSOR})
find_package(GCTL REQUIRED)
find_package(GCTL_GRAPHIC REQUIRED)
message(STATUS "GCTL Version: " ${GCTL_VERSION})
message(STATUS "GCTL_GRAPHIC Version: " ${GCTL_GRAPHIC_VERSION})
#if(${GCTL_VERSION} LESS 1.0)
# message(FATAL_ERROR "GCTL's version must be v1.0 or bigger.")
#endif()
option(GCTL_MESH_EXPRTK "Use the exprtk library." ON)
option(GCTL_MESH_WAVELIB "Use the WaveLib library" ON)

View File

@ -12,18 +12,19 @@ macro(add_example name switch)
#
target_link_libraries(${name} PRIVATE ${GCTL_LIB})
target_link_libraries(${name} PRIVATE ${GCTL_GRAPHIC_LIB})
target_link_libraries(${name} PRIVATE gctl_mesh)
endif()
endmacro()
add_example(mesh_ex1 ON)
add_example(mesh_ex2 ON)
add_example(mesh_ex3 ON)
add_example(mesh_ex4 ON)
add_example(mesh_ex1 OFF)
add_example(mesh_ex2 OFF)
add_example(mesh_ex3 OFF)
add_example(mesh_ex4 OFF)
add_example(mesh_ex5 OFF)
add_example(mesh_ex6 ON)
add_example(mesh_ex7 ON)
add_example(mesh_ex8 ON)
add_example(mesh_ex9 ON)
add_example(mesh_ex10 ON)
add_example(mesh_ex6 OFF)
add_example(mesh_ex7 OFF)
add_example(mesh_ex8 OFF)
add_example(mesh_ex9 OFF)
add_example(mesh_ex10 OFF)
add_example(meshio_ex OFF)
add_example(tri2d_meshio_ex ON)

View File

@ -45,7 +45,10 @@ int main(int argc, char *argv[])
std::cout << "data name: " << data3.name_ << std::endl;
data3.datval_.show();
rgd.diff("data-4", "data-1", "data-3");
rgd.remove_data("data-1");
rgd.show_info();
rgd.save_gmsh_withdata("ex1_out", gctl::OverWrite, gctl::NotPacked);
return 0;
}

109
example/meshio_ex.cpp Normal file
View File

@ -0,0 +1,109 @@
/********************************************************
*
*
*
*
*
*
* Geophysical Computational Tools & Library (GCTL)
*
* Copyright (c) 2023 Yi Zhang (yizhang-geo@zju.edu.cn)
*
* GCTL is distributed under a dual licensing scheme. You can redistribute
* it and/or modify it under the terms of the GNU Lesser General Public
* License as published by the Free Software Foundation, either version 2
* of the License, or (at your option) any later version. You should have
* received a copy of the GNU Lesser General Public License along with this
* program. If not, see <http://www.gnu.org/licenses/>.
*
* If the terms and conditions of the LGPL v.2. would prevent you from using
* the GCTL, please consider the option to obtain a commercial license for a
* fee. These licenses are offered by the GCTL's original author. As a rule,
* licenses are provided "as-is", unlimited in time for a one time fee. Please
* send corresponding requests to: yizhang-geo@zju.edu.cn. Please do not forget
* to include some description of your company and the realm of its activities.
* Also add information on how to contact you by electronic and paper mail.
******************************************************/
#include "../lib/mesh/mesh_io.h"
using namespace gctl;
int main(int argc, char const *argv[]) try
{
mesh_io mshio;
mshio.read_gmsh_v2_ascii("tmp/test");
mshio.info();
array<vertex3dc> nodes;
array<triangle> tris;
mshio.select_elements("Bz (pT)", NodeData);
mshio.export_selected_to(tris, nodes);
//meshio_data &d = mshio.get_data("Bz (pT)", NodeData);
//array<double> dat = d.val;
array<double> dat1 = mshio.get_selected_data("Bx (pT)", NodeData);
array<double> dat2 = mshio.get_selected_data("By (pT)", NodeData);
array<double> dat3 = mshio.get_selected_data("Bz (pT)", NodeData);
mesh_io mshio2;
mshio2.import_from(tris, nodes);
mshio2.add_data("Bx (pT)", dat1, NodeData, OverWrite);
mshio2.add_data("By (pT)", dat2, NodeData, OverWrite);
mshio2.add_data("Bz (pT)", dat3, NodeData, OverWrite);
mshio2.save_gmsh_v2_ascii("tmp");
mshio2.info();
/*
mshio.read_tetgen_ascii("tmp/ex1.1");
mshio.edit_group(Disable, GeometryTag, 5);
mshio.edit_group(GeometryTag, 1, PhysicalTag, 1);
mshio.edit_group(GeometryTag, 2, PhysicalTag, 2);
mshio.edit_group(GeometryTag, 3, PhysicalTag, 3);
mshio.edit_group(GeometryTag, 4, PhysicalTag, 4);
mshio.edit_group(GeometryTag, 1, "Boundary");
mshio.edit_group(GeometryTag, 2, "Body1");
mshio.edit_group(GeometryTag, 3, "Body2");
mshio.edit_group(GeometryTag, 4, "Body3");
mshio.save_gmsh_v2_ascii("tmp/ex1.1");
*/
/*
mshio.read_gmsh_v2_ascii("tmp/ex1.1");
mshio.convert_tags_to_data(GeometryTag);
array<double> body_val(mshio.element_size("Body2"), 2.0);
mshio.add_element_data("BodyValue", "Body2", body_val);
array<double> body_val2(mshio.element_size("Body3"), 1.0);
mshio.add_element_data("BodyValue", "Body3", body_val2);
mshio.save_gmsh_v2_ascii("tmp/ex1.2");
//mshio.save_vtk_legacy_ascii("tmp/ex1.1");
mshio.info();
const array<vertex3dc> &nodes = mshio.get_nodes();
array<tetrahedron> body2_tets;
mshio.export_elements_to(body2_tets, "All");
gmshio gio;
gio.init_file("tmp.msh", Output);
gio.set_packed(NotPacked, Output);
gio.save_mesh(body2_tets, nodes);
*/
/*
mshio.read_gmsh_v2_ascii("tmp/wjb.1");
mshio.edit_group(Disable);
mshio.edit_group(Enable, GeometryTag, 3);
mshio.edit_group(Enable, GeometryTag, 8);
mshio.edit_group(Enable, GeometryTag, 9);
mshio.save_gmsh_v2_ascii("tmp/wjb.2");
*/
return 0;
}
catch(std::exception &e)
{
GCTL_ShowWhatError(e.what(), GCTL_ERROR_ERROR, 0, 0, 0);
}

View File

@ -0,0 +1,55 @@
/********************************************************
*
*
*
*
*
*
* Geophysical Computational Tools & Library (GCTL)
*
* Copyright (c) 2023 Yi Zhang (yizhang-geo@zju.edu.cn)
*
* GCTL is distributed under a dual licensing scheme. You can redistribute
* it and/or modify it under the terms of the GNU Lesser General Public
* License as published by the Free Software Foundation, either version 2
* of the License, or (at your option) any later version. You should have
* received a copy of the GNU Lesser General Public License along with this
* program. If not, see <http://www.gnu.org/licenses/>.
*
* If the terms and conditions of the LGPL v.2. would prevent you from using
* the GCTL, please consider the option to obtain a commercial license for a
* fee. These licenses are offered by the GCTL's original author. As a rule,
* licenses are provided "as-is", unlimited in time for a one time fee. Please
* send corresponding requests to: yizhang-geo@zju.edu.cn. Please do not forget
* to include some description of your company and the realm of its activities.
* Also add information on how to contact you by electronic and paper mail.
******************************************************/
#include "../lib/mesh/tri2d_mesh.h"
using namespace gctl;
int main(int argc, char const *argv[]) try
{
mesh_io mshio;
mshio.read_gmsh_v2_ascii("tmp/tri2d");
mshio.info();
array<vertex2dc> nodes;
array<triangle2d> tris;
mshio.select_elements(_3NodeTriangle);
mshio.export_selected_to(tris, nodes);
_1d_array data = mshio.get_selected_data("Model slowness (s/km)", ElemData);
triangle2d_mesh tri2dmesh;
tri2dmesh.init("Test", "This is a test.", nodes, tris);
tri2dmesh.add_data(ElemData, "Test Data", data);
tri2dmesh.show_info();
tri2dmesh.save_gmsh_withdata("tri2d_out", OverWrite);
return 0;
}
catch(std::exception &e)
{
GCTL_ShowWhatError(e.what(), GCTL_ERROR_ERROR, 0, 0, 0);
}

View File

@ -27,9 +27,7 @@ set_target_properties(gctl_mesh_static PROPERTIES CXX_STANDARD 17 CXX_STANDARD_R
#
target_link_libraries(gctl_mesh PUBLIC ${GCTL_LIB})
target_link_libraries(gctl_mesh PUBLIC ${GCTL_GRAPHIC_LIB})
target_link_libraries(gctl_mesh_static ${GCTL_LIB})
target_link_libraries(gctl_mesh_static ${GCTL_GRAPHIC_LIB})
if(GCTL_MESH_WAVELIB)
target_link_libraries(gctl_mesh PUBLIC ${WaveLib_LIB})

View File

@ -29,6 +29,8 @@
#define _GCTL_MESH_H
#include "mesh/gctl_mesh_config.h"
#include "mesh/meshdata.h"
#include "mesh/meshio.h"
#include "mesh/mesh.h"
#include "mesh/regular_grid.h"
#include "mesh/regular_mesh_2d.h"
@ -36,6 +38,7 @@
#include "mesh/linear_mesh_2d.h"
#include "mesh/linear_mesh_3d.h"
#include "mesh/tri2d_mesh.h"
#include "mesh/tri_mesh.h"
#include "mesh/tet_mesh.h"
#include "mesh/regular_mesh_sph_3d.h"

View File

@ -276,9 +276,8 @@ const gctl::array<double> *gctl::linear_mesh_2d::get_ysizes() const
void gctl::linear_mesh_2d::save_gmsh(std::string filename, index_packed_e packed)
{
std::ofstream outfile;
gctl::open_outfile(outfile, filename, ".msh");
gctl::save2gmsh(outfile, elements, nodes, packed);
outfile.close();
meshio_.init_file(filename, Output);
meshio_.set_packed(packed, Output);
meshio_.save_mesh(elements, nodes);
return;
}

View File

@ -77,6 +77,6 @@ namespace gctl
array<vertex2dc> nodes;
array<rectangle2d> elements;
};
}
};
#endif //_GCTL_LINEAR_MESH_2D_H

View File

@ -366,9 +366,8 @@ const gctl::array<double> *gctl::linear_mesh_3d::get_zsizes() const
void gctl::linear_mesh_3d::save_gmsh(std::string filename, index_packed_e packed)
{
std::ofstream outfile;
gctl::open_outfile(outfile, filename, ".msh");
gctl::save2gmsh(outfile, elements, nodes, packed);
outfile.close();
meshio_.init_file(filename, Output);
meshio_.set_packed(packed, Output);
meshio_.save_mesh(elements, nodes);
return;
}

View File

@ -84,6 +84,6 @@ namespace gctl
array<vertex3dc> nodes;
array<block> elements;
};
}
};
#endif //_GCTL_LINEAR_MESH_3D_H

View File

@ -35,6 +35,15 @@ gctl::base_mesh::base_mesh()
meshinfo_ = "Undefined";
node_num_ = ele_num_ = 0;
initialized_ = false;
// 注意这里我们一定要预先为datalist_分配空间因为
// datalist_为vector类型所以元素增加时可能造成
// 原有的迭代器失效造成访问失败。一个简单的例子是做
// 网格加法时我们通过get_data获取两个网格数据并将
// 计算结果通过add_data保存这时候就可能出现vector
// 扩容造成get_data获取的数据无法访问访问失败。
// 默认可以保存100个网格数据应该是够用了。同时
// 在add_data中添加相应的错误提示。
datalist_.reserve(100);
}
gctl::base_mesh::~base_mesh()
@ -51,6 +60,7 @@ void gctl::base_mesh::clear()
node_num_ = ele_num_ = 0;
initialized_ = false;
destroy_vector(datalist_);
datalist_.reserve(100);
return;
}
@ -182,6 +192,8 @@ gctl::meshdata &gctl::base_mesh::add_data(mesh_data_type_e in_loctype, mesh_data
std::string name, double init_val, bool if_output, double nan_val)
{
check_initiated();
if (datalist_.size() == 100) throw std::runtime_error("[gctl::base_mesh] Maximal data number reached.");
meshdata new_data(in_loctype, in_valtype, 0, name, if_output, nan_val);
if (in_loctype == NodeData && in_valtype == Scalar) new_data.datval_.resize(node_num_, init_val);
@ -200,6 +212,8 @@ gctl::meshdata &gctl::base_mesh::add_data(mesh_data_type_e in_loctype, std::stri
const array<double> &init_arr, bool if_output, double nan_val)
{
check_initiated();
if (datalist_.size() == 100) throw std::runtime_error("[gctl::base_mesh] Maximal data number reached.");
meshdata new_data(in_loctype, Scalar, 0, name, if_output, nan_val);
if (in_loctype == NodeData)
@ -222,6 +236,8 @@ gctl::meshdata &gctl::base_mesh::add_data(mesh_data_type_e in_loctype, std::stri
const array<point3dc> &init_arr, bool if_output, double nan_val)
{
check_initiated();
if (datalist_.size() == 100) throw std::runtime_error("[gctl::base_mesh] Maximal data number reached.");
meshdata new_data(in_loctype, Vector, 0, name, if_output, nan_val);
if (in_loctype == NodeData)
@ -258,6 +274,8 @@ gctl::meshdata &gctl::base_mesh::add_data(mesh_data_type_e in_loctype, std::stri
const array<tensor> &init_arr, bool if_output, double nan_val)
{
check_initiated();
if (datalist_.size() == 100) throw std::runtime_error("[gctl::base_mesh] Maximal data number reached.");
meshdata new_data(in_loctype, Tensor, 0, name, if_output, nan_val);
if (in_loctype == NodeData)

View File

@ -29,8 +29,7 @@
#define _GCTL_BASE_MESH_H
#include "meshdata.h"
#include "gctl/io.h"
#include "gctl/algorithm.h"
#include "meshio.h"
namespace gctl
{
@ -377,6 +376,6 @@ namespace gctl
*/
void save_datablock(std::ofstream &outfile);
};
}
};
#endif //_GCTL_BASE_MESH_H

View File

@ -63,6 +63,8 @@ void gctl::meshdata::create(mesh_data_type_e in_loctype,
mesh_data_value_e in_valtype, size_t size,
std::string name, bool if_output, double nan_val)
{
if (name == "") throw std::runtime_error("[gctl::meshdata::create] Invalid data name.");
name_ = name;
loctype_ = in_loctype;
valtype_ = in_valtype;

View File

@ -29,9 +29,8 @@
#define _GCTL_MESHDATA_H
#include "gctl_mesh_config.h"
#include "gctl/core.h"
#include "gctl/geometry/point3c.h"
#include "gctl/geometry/tensor.h"
#include "gctl/poly/point3c.h"
#include "gctl/poly/tensor.h"
namespace gctl
{
@ -102,7 +101,7 @@ namespace gctl
bool if_output = true, double nan_val = GCTL_BDL_MAX);
/**
* @brief
* @brief
*/
array<point3dc> export_vector() const;
@ -139,6 +138,6 @@ namespace gctl
*/
void save_binary(std::ofstream &outfile);
};
}
};
#endif // _GCTL_MESHDATA_H

2024
lib/mesh/meshio.cpp Normal file

File diff suppressed because it is too large Load Diff

734
lib/mesh/meshio.h Normal file
View File

@ -0,0 +1,734 @@
/********************************************************
*
*
*
*
*
*
* Geophysical Computational Tools & Library (GCTL)
*
* Copyright (c) 2023 Yi Zhang (yizhang-geo@zju.edu.cn)
*
* GCTL is distributed under a dual licensing scheme. You can redistribute
* it and/or modify it under the terms of the GNU Lesser General Public
* License as published by the Free Software Foundation, either version 2
* of the License, or (at your option) any later version. You should have
* received a copy of the GNU Lesser General Public License along with this
* program. If not, see <http://www.gnu.org/licenses/>.
*
* If the terms and conditions of the LGPL v.2. would prevent you from using
* the GCTL, please consider the option to obtain a commercial license for a
* fee. These licenses are offered by the GCTL's original author. As a rule,
* licenses are provided "as-is", unlimited in time for a one time fee. Please
* send corresponding requests to: yizhang-geo@zju.edu.cn. Please do not forget
* to include some description of your company and the realm of its activities.
* Also add information on how to contact you by electronic and paper mail.
******************************************************/
#ifndef _GCTL_MESH_IO_H
#define _GCTL_MESH_IO_H
#include <unordered_set>
#include <map>
#include "gctl/math/refellipsoid.h"
#include "gctl/io/triangle_io.h"
#include "gctl/io/tetgen_io.h"
#include "gctl/io/gmsh_io.h"
namespace gctl
{
/**
* @brief
*/
#define DEFAULT_INVALID_TAG -9999
/**
* @brief
*/
enum element_type_enum
{
NotSet,
_2NodeLine,
_3NodeTriangle,
_4NodeQuadrangle,
_4NodeTetrahedron,
_8NodeHexahedron,
_6NodePrism,
_5NodePyramid,
_3NodeSecondOrderLine,
_6NdoeSecondOrderLine,
_9NodeSecondOrderQuadrangle,
_10NodeSecondOrderTetrahedron,
_27NodeSecondOrderHexahedron,
_18NodeSecondOrderPrism,
_14NodeSecondOrderPyramid,
_1NodePoint,
_8NodeSecondOrderQuadrangle,
_20NdoeSecondOrderHexahedron,
_15NodeSecondOrderPrism,
_13NodeSecondOrderPyramid,
_9NodeThirdOrderIncompleteTriangle,
_10NdoeThirdOrderTriangle,
_12NodeFourthOrderIncompleteTriangle,
_15NodeFourthOrderTriangle,
_15NodeFifthOrderCompleteTriangle,
_21NodeFifthOrderCompleteTriangle,
_4NodeThirdOrderEdge,
_5NodeFourthOrderEdge,
_6NodeFifthOrderEdge,
_20NodeThirdOrderTetrahedron,
_35NodeFourthOrderTetrahedron,
_56NodeFifithOrderTetrahedron,
_64NodeThirdOrderHexahedron,
_125NodeFourthOrderHexahedron,
};
/**
* @brief
*
*/
enum element_tag_enum
{
PhysicalTag, // 元素的物理分组标签
GeometryTag, // 元素的几何分组标签
PartitionTag, // 元素的剖分分组标签
NodeTag, // 顶点的标签(仅用于输出顶点标签数据)
};
/**
* @brief
*
*/
struct meshio_element
{
bool enabled; // 单元体是否有效
int id; // 单元体编号
element_type_enum type; // 单元体类型
array<vertex3dc*> vert_ptrs; // 顶点指针数组
array<meshio_element*> neigh_ptrs; // 相邻单元体指针数组
meshio_element();
};
/**
* @brief
*
*/
struct meshio_data
{
bool enabled; // 数据体是否有效
mesh_data_type_e d_type; // 数据类型
array<std::string> str_tag; // 字符串类型的标签(默认为一个,即为数据的名称)
array<double> real_tag; // 实数类型的标签默认为一个等于0.0
array<int> int_tag; // 整数类型的标签(最少三个,最后一个为数据的长度)
array<void*> tar_ptrs; // 数据连接的对象指针数组 具体类型为vertex3dc*或meshio_element*
array<double> val; // 数据值(目前仅支持标量数据,后续再添加对矢量数据的支持)
meshio_data();
/**
* @brief
*
*/
void clear();
/**
* @brief
*
*/
bool pass_check();
};
/**
* @brief
*
*/
struct meshio_element_group
{
bool enabled; // 组是否有效
element_type_enum type; // 组内单元体类型
int phys_group; // 物理分组标签
int geom_group; // 几何分组标签
int part_group; // 剖分分组标签
std::string name; // 组名
std::vector<meshio_element*> elem_ptrs; // 组内单元体指针数组
meshio_element_group();
/**
* @brief
*
*/
void enable_elements();
/**
* @brief
*
*/
void disable_elements();
};
/**
* @brief
*
*/
class mesh_io
{
public:
mesh_io();
virtual ~mesh_io();
/**
* @brief
*
*/
void reset();
/**
* @brief
*
* @param ss clog
*/
void info(std::ostream &ss = std::clog);
/**
* @brief
*
* @param swt 使Enable或Disable
* @param e_type NotSet
*/
void edit_group(switch_type_e swt, element_type_enum e_type = NotSet);
/**
* @brief
*
* @param swt 使Enable或Disable
* @param grp_name
*/
void edit_group(switch_type_e swt, std::string grp_name);
/**
* @brief
*
* @param swt 使Enable或Disable
* @param tag_type PhysicalTagGeometryTag或者PartitionTag
* @param tag
*/
void edit_group(switch_type_e swt, element_tag_enum tag_type, int tag);
/**
* @brief
*
* @param anchor_type PhysicalTagGeometryTag或者PartitionTag
* @param anchor_group
* @param new_name
*/
void edit_group(element_tag_enum anchor_type, int anchor_group, std::string new_name);
/**
* @brief
*
* @param anchor_type PhysicalTagGeometryTag或者PartitionTag
* @param anchor_group
* @param tar_type PhysicalTagGeometryTag或者PartitionTag
* @param tar_group
*/
void edit_group(element_tag_enum anchor_type, int anchor_group, element_tag_enum tar_type, int tar_group);
/**
* @brief
*
* @param tag_type PhysicalTagGeometryTag或者PartitionTag
* @param group_name
* @return
*/
int group_tag(element_tag_enum tag_type, std::string group_name);
/**
* @brief
*
* @return
*/
const array<vertex3dc> &get_all_nodes();
/**
* @brief
*
* @return
*/
const array<meshio_element> &get_all_elems();
/**
* @brief
*
* @param e_type NotSet
*/
void select_elements(element_type_enum e_type = NotSet);
/**
* @brief
*
* @param tag_type
* @param tag
*/
void select_elements(element_tag_enum tag_type, int tag);
/**
* @brief
*
* @param phys_name
*/
void select_elements(std::string phys_name);
/**
* @brief
*
* @param datname
* @param dtype
*/
void select_elements(std::string dat_name, mesh_data_type_e dtype);
/**
* @brief 使select_elements函数选择
*
* @return
*/
size_t selected_node_size();
/**
* @brief 使select_elements函数选择
*
* @return
*/
size_t selected_element_size();
/**
* @brief -9999
*
* @return
*/
const array<int> &selected_node_tags();
/**
* @brief -9999
*
* @return
*/
const array<int> &selected_element_tags(element_tag_enum tag_type);
/**
* @brief gmsh格式分组表
*
* @param g_groups gmsh格式表
*/
void get_gmsh_physical_groups(std::vector<gmsh_physical_group> &g_groups);
/**
* @brief
*
* @param tag_type
*/
void convert_tags_to_data(element_tag_enum tag_type);
/**
* @brief name的数据
*
* @param name
* @param type
*
* @return -1
*/
int if_saved_data(std::string name, mesh_data_type_e type);
/**
* @brief
*
* @param name
* @param type
* @return
*/
meshio_data &get_data(std::string name, mesh_data_type_e type);
/**
* @brief
*
* @param name
* @param type
* @return
*/
meshio_data *get_data_ptr(std::string name, mesh_data_type_e type);
/**
* @brief
*
* @param name
* @param type
*
* @return
*/
array<double> get_selected_data(std::string name, mesh_data_type_e type);
/**
* @brief
*
* @param data
* @param name
* @param dtype
* @param op OverWrite或Append
*/
void add_data(std::string name, const array<double> &data,
mesh_data_type_e dtype, output_type_e op = OverWrite);
/**
* @brief
*
* @param swt 使Enable或Disable
* @param dataname null
*/
void edit_data(switch_type_e swt, std::string dataname = "null");
/**
* @brief triangle软件输出的网格剖分文件
*
* @param filename .node和.ele文件必须在同一路径下.neigh文件不是必须的
* @param is_packed 0
*/
void read_triangle_ascii(std::string filename, index_packed_e is_packed = Packed);
/**
* @brief tetgen软件输出的网格剖分文件
*
* @param filename .node和.ele文件必须在同一路径下.neigh文件不是必须的
* @param is_packed 0
*/
void read_tetgen_ascii(std::string filename, index_packed_e is_packed = Packed);
/**
* @brief Gmsh软件输出的网格剖分文件v2.2ASCII文件
*
* @param filename
* @param is_packed 0
*/
void read_gmsh_v2_ascii(std::string filename, index_packed_e is_packed = NotPacked);
/**
* @brief Gmsh软件格式的网格剖分文件v2.2ASCII文件
*
* @param filename
* @param is_packed 0
*/
void save_gmsh_v2_ascii(std::string filename, index_packed_e is_packed = NotPacked);
/**
* @brief Paraview软件格式的网格剖分文件
*
* @param filename
*/
void save_vtk_legacy_ascii(std::string filename);
/**
* @brief
*
* @param filename
* @param dataname
* @param out_coor Cartesian或Spherical
* @param refr out_coor为Cartesian时无效
* @param refR out_coor为Cartesian时无效
*/
void save_data_to_xyz(std::string filename, std::string dataname = "null", coordinate_system_e out_coor = Cartesian, double refr = GCTL_Earth_Radius, double refR = GCTL_Earth_Radius);
/**
* @brief
*
* @tparam T
* @param elems
* @param nodes
*/
template <typename T>
void export_selected_to(array<T> &elems, array<vertex3dc> &nodes);
/**
* @brief
*
* @tparam T
* @param elems
* @param nodes
* @param x_id x索引0
* @param y_id y索引1
*
* @note xyz坐标中的对应项分别赋值给二维坐标的xy值{0,1}xoy平面
* {1,0}yox平面{0,2}xoz平面{1,2}yoz平面
*/
template <typename T>
void export_selected_to(array<T> &elems, array<vertex2dc> &nodes, int x_id = 0, int y_id = 1);
/**
* @brief
*
* @tparam T
* @param elems
* @param nodes
*/
template <typename T>
void import_from(const array<T> &elems, const array<vertex3dc> &nodes);
/**
* @brief
*
* @tparam T
* @param elems
* @param nodes
* @param x_id x坐标到三维坐标系统的索引0
* @param y_id y坐标到三维坐标系统的索引1
*
* @note (x,y)(x,y,z){0,1}(x,y)(x,y,0)
* {1,0}(x,y)(y,x,0){0,2}(x,y)(x,0,z){1,2}(x,y)(0,y,z)
*/
template <typename T>
void import_from(const array<T> &elems, const array<vertex2dc> &nodes, int x_id = 0, int y_id = 1);
protected:
bool initialized_; // 类型是否已经初始化完成
// 有效的顶点、单元体和单元体组的数量
size_t valid_node_size_, valid_elem_size_, valid_group_size_;
array<vertex3dc> nodes_; // 网格顶点 当顶点索引为无效值时将不会被输出
array<meshio_element> elems_; // 网格元素
std::vector<meshio_data> datas_; // 网格数据
std::vector<meshio_element_group> groups_; // 网格单元体组
array<int> nodes_tag_; // 顶点标签
array<int> selected_node_tag_; // 被选择的顶点的标签
array<int> selected_elem_tag_; // 被选择的单元体的标签
std::vector<vertex3dc*> selected_nodes_; // 被选择的顶点
std::vector<meshio_element*> selected_elems_; // 被选择的单元体
std::vector<meshio_element_group*> selected_groups_; // 被选择的单元体组
element_type_enum elem_gmshcode2type_[94]; // gmsh的单元体类型数组 数组索引为gmsh的单元体类型码值
element_type_enum elem_vtkcode2type_[14]; // vtk的单元体类型数组 数组索引为vtk的单元体类型码值
std::map<element_type_enum, int> elem_type2gmshcode_; // 单元体类型到gmsh类型码值的映射
std::map<element_type_enum, int> elem_type2vtkcode_; // 单元体类型到vtk类型码值的映射
std::map<element_type_enum, int> elem_type2size_; // 单元体类型到单元体顶点数量的映射
std::map<element_type_enum, std::string> elem_type2name_; // 单元体类型到单元体名称的映射
std::string elem_name(element_type_enum e_type); // 获取单元体名称字符串
int elem_gmsh_code(element_type_enum e_type); // 获取单元体gmsh类型码值
int elem_vtk_code(element_type_enum e_type); // 获取单元体vtk类型码值
int elem_size(element_type_enum e_type); // 获取单元体顶点数量
element_type_enum elem_gmsh_type(int code); // 获取对应gmsh类型码的单元体类型
element_type_enum elem_vtk_type(int code); // 获取对应vtk类型码的单元体类型
void update_indexing(); // 更新索引(对网格的元素进行操作后需调用)
void sort_groups(); // 对单元体组进行梳理(对网格的元素进行操作后需调用)
element_type_enum match_type(const std::type_info &tinfo); // 根据type_info返回单元体类型
};
template <typename T>
void gctl::mesh_io::export_selected_to(array<T> &elems, array<vertex3dc> &nodes)
{
const std::type_info &tinfo = typeid(T);
element_type_enum oe_type = match_type(tinfo);
size_t n = selected_nodes_.size();
size_t s = selected_elems_.size();
std::map<int, int> node_map;
nodes.resize(n);
for (size_t i = 0; i < n; i++)
{
node_map[selected_nodes_[i]->id] = i; // 原网格顶点索引到新网格顶点索引的映射
nodes[i].id = i;
nodes[i].x = selected_nodes_[i]->x;
nodes[i].y = selected_nodes_[i]->y;
nodes[i].z = selected_nodes_[i]->z;
}
int vnum = elem_size(oe_type);
elems.resize(s);
for (size_t i = 0; i < s; i++)
{
if (selected_elems_[i]->type != oe_type)
throw gctl::invalid_argument("[gctl::mesh_io::export_to] The selected element is not " + elem_name(oe_type) + ".");
elems[i].id = i;
for (size_t v = 0; v < vnum; v++)
elems[i].vert[v] = nodes.get(node_map[selected_elems_[i]->vert_ptrs[v]->id]);
}
node_map.clear();
return;
}
template <typename T>
void gctl::mesh_io::export_selected_to(array<T> &elems, array<vertex2dc> &nodes, int x_id, int y_id)
{
const std::type_info &tinfo = typeid(T);
element_type_enum oe_type = match_type(tinfo);
size_t n = selected_nodes_.size();
size_t s = selected_elems_.size();
std::map<int, int> node_map;
nodes.resize(n);
double xyz_ref[3];
for (size_t i = 0; i < n; i++)
{
node_map[selected_nodes_[i]->id] = i; // 原网格顶点索引到新网格顶点索引的映射
nodes[i].id = i;
xyz_ref[0] = selected_nodes_[i]->x;
xyz_ref[1] = selected_nodes_[i]->y;
xyz_ref[2] = selected_nodes_[i]->z;
nodes[i].x = xyz_ref[x_id];
nodes[i].y = xyz_ref[y_id];
}
int vnum = elem_size(oe_type);
elems.resize(s);
for (size_t i = 0; i < s; i++)
{
if (selected_elems_[i]->type != oe_type)
throw gctl::invalid_argument("[gctl::mesh_io::export_to] The selected element is not " + elem_name(oe_type) + ".");
elems[i].id = i;
for (size_t v = 0; v < vnum; v++)
elems[i].vert[v] = nodes.get(node_map[selected_elems_[i]->vert_ptrs[v]->id]);
}
node_map.clear();
return;
}
template <typename T>
void gctl::mesh_io::import_from(const array<T> &elems, const array<vertex3dc> &nodes)
{
reset(); // 重置网格数据
const std::type_info &tinfo = typeid(T);
element_type_enum oe_type = match_type(tinfo);
int vnum = elem_size(oe_type);
valid_node_size_ = nodes.size();
valid_elem_size_ = elems.size();
nodes_.resize(valid_node_size_);
selected_nodes_.resize(valid_node_size_);
for (size_t i = 0; i < valid_node_size_; i++)
{
nodes_[i].id = i;
nodes_[i].x = nodes[i].x;
nodes_[i].y = nodes[i].y;
nodes_[i].z = nodes[i].z;
selected_nodes_[i] = nodes_.get(i);
}
elems_.resize(valid_elem_size_);
selected_elems_.resize(valid_elem_size_);
for (size_t i = 0; i < valid_elem_size_; i++)
{
elems_[i].id = i;
elems_[i].type = oe_type;
elems_[i].vert_ptrs.resize(vnum);
for (size_t v = 0; v < vnum; v++)
{
elems_[i].vert_ptrs[v] = nodes_.get(elems[i].vert[v]->id);
}
selected_elems_[i] = elems_.get(i);
}
// 所有单元体都属于同一个组
valid_group_size_ = 1;
groups_.resize(1);
groups_[0].enabled = true;
groups_[0].type = oe_type;
groups_[0].phys_group = 0; // 默认组别为0
groups_[0].geom_group = 0; // 默认组别为0
groups_[0].part_group = 3;
for (size_t i = 0; i < elems_.size(); i++)
{
groups_[0].elem_ptrs.push_back(elems_.get(i));
}
groups_[0].enable_elements();
initialized_ = true;
selected_groups_.resize(1);
selected_groups_[0] = &groups_[0];
return;
}
template <typename T>
void gctl::mesh_io::import_from(const array<T> &elems, const array<vertex2dc> &nodes, int x_id, int y_id)
{
reset(); // 重置网格数据
const std::type_info &tinfo = typeid(T);
element_type_enum oe_type = match_type(tinfo);
int vnum = elem_size(oe_type);
valid_node_size_ = nodes.size();
valid_elem_size_ = elems.size();
int xyz_ref[3];
nodes_.resize(valid_node_size_);
selected_nodes_.resize(valid_node_size_);
for (size_t i = 0; i < valid_node_size_; i++)
{
nodes_[i].id = i;
xyz_ref[0] = 0.0;
xyz_ref[1] = 0.0;
xyz_ref[2] = 0.0;
xyz_ref[x_id] = nodes[i].x;
xyz_ref[y_id] = nodes[i].y;
nodes_[i].x = xyz_ref[0];
nodes_[i].y = xyz_ref[1];
nodes_[i].z = xyz_ref[2];
selected_nodes_[i] = nodes_.get(i);
}
elems_.resize(valid_elem_size_);
selected_elems_.resize(valid_elem_size_);
for (size_t i = 0; i < valid_elem_size_; i++)
{
elems_[i].id = i;
elems_[i].type = oe_type;
elems_[i].vert_ptrs.resize(vnum);
for (size_t v = 0; v < vnum; v++)
{
elems_[i].vert_ptrs[v] = nodes_.get(elems[i].vert[v]->id);
}
selected_elems_[i] = elems_.get(i);
}
// 所有单元体都属于同一个组
valid_group_size_ = 1;
groups_.resize(1);
groups_[0].enabled = true;
groups_[0].type = oe_type;
groups_[0].phys_group = 0; // 默认组别为0
groups_[0].geom_group = 0; // 默认组别为0
groups_[0].part_group = 3;
for (size_t i = 0; i < elems_.size(); i++)
{
groups_[0].elem_ptrs.push_back(elems_.get(i));
}
groups_[0].enable_elements();
initialized_ = true;
selected_groups_.resize(1);
selected_groups_[0] = &groups_[0];
return;
}
};
#endif // _GCTL_MESH_IO_H

View File

@ -652,10 +652,9 @@ void gctl::regular_grid::save_surfer_grid(std::string filename, std::string datn
void gctl::regular_grid::save_gmsh(std::string filename, index_packed_e packed)
{
std::ofstream outfile;
gctl::open_outfile(outfile, filename, ".msh");
gctl::save2gmsh(outfile, elements, nodes, packed);
outfile.close();
meshio_.init_file(filename, Output);
meshio_.set_packed(packed, Output);
meshio_.save_mesh(elements, nodes);
return;
}
@ -951,7 +950,6 @@ void gctl::regular_grid::diff(std::string newname, std::string datname, std::str
}
meshdata &new_data = add_data(data_type1, value_type1, newname, 0.0, true, GCTL_BDL_MAX);
for (size_t i = 0; i < new_data.datval_.size(); i++)
{
new_data.datval_[i] = data_1.datval_[i] - data_2.datval_[i];

View File

@ -28,7 +28,10 @@
#ifndef _GCTL_REGULAR_GRID_H
#define _GCTL_REGULAR_GRID_H
#include "gctl/graphic.h"
#include "gctl/math/interpolate.h"
#include "gctl/graphic/gmt.h"
#include "gctl/io/surfer_io.h"
#include "gctl/io/netcdf_io.h"
#include "mesh.h"
#ifdef GCTL_MESH_EXPRTK
@ -328,6 +331,6 @@ namespace gctl
array<vertex2dc> nodes; ///< 规则网格的节点数组
array<rectangle2d> elements; ///< 规则网格的单元数组
};
}
};
#endif //_GCTL_REGULAR_GRID_H

View File

@ -213,9 +213,8 @@ double gctl::regular_mesh_2d::get_ysize() const
void gctl::regular_mesh_2d::save_gmsh(std::string filename, index_packed_e packed)
{
std::ofstream outfile;
gctl::open_outfile(outfile, filename, ".msh");
gctl::save2gmsh(outfile, elements, nodes, packed);
outfile.close();
meshio_.init_file(filename, Output);
meshio_.set_packed(packed, Output);
meshio_.save_mesh(elements, nodes);
return;
}

View File

@ -74,6 +74,6 @@ namespace gctl
array<vertex2dc> nodes;
array<rectangle2d> elements;
};
}
};
#endif //_GCTL_REGULAR_MESH_2D_H

View File

@ -260,9 +260,8 @@ double gctl::regular_mesh_3d::get_zsize() const
void gctl::regular_mesh_3d::save_gmsh(std::string filename, index_packed_e packed)
{
std::ofstream outfile;
gctl::open_outfile(outfile, filename, ".msh");
gctl::save2gmsh(outfile, elements, nodes, packed);
outfile.close();
meshio_.init_file(filename, Output);
meshio_.set_packed(packed, Output);
meshio_.save_mesh(elements, nodes);
return;
}

View File

@ -77,6 +77,6 @@ namespace gctl
array<vertex3dc> nodes;
array<block> elements;
};
}
};
#endif //_GCTL_REGULAR_MESH_3D_H

View File

@ -271,16 +271,16 @@ void gctl::regular_mesh_sph_3d::save_gmsh(std::string filename, index_packed_e p
gctl::array<vertex3dc> tmp_nodes(nodes.size());
for (int i = 0; i < nodes.size(); i++)
{
tmp_c = nodes[i].s2c();
tmp_c = s2c(nodes[i]);
tmp_nodes[i].id = i;
tmp_nodes[i].x = tmp_c.x;
tmp_nodes[i].y = tmp_c.y;
tmp_nodes[i].z = tmp_c.z;
}
std::ofstream outfile;
gctl::open_outfile(outfile, filename, ".msh");
gctl::save2gmsh(outfile, elements, tmp_nodes, packed);
outfile.close();
meshio_.init_file(filename, Output);
meshio_.set_packed(packed, Output);
meshio_.save_mesh(elements, tmp_nodes);
tmp_nodes.clear();
return;
}

View File

@ -76,6 +76,6 @@ namespace gctl
array<vertex3ds> nodes;
array<tesseroid> elements;
};
}
};
#endif // _GCTL_REGULAR_MESH_SPH_3D_H

View File

@ -197,9 +197,8 @@ void gctl::tetrahedron_mesh::load_gmsh(std::string filename, index_packed_e pack
void gctl::tetrahedron_mesh::save_gmsh(std::string filename, index_packed_e packed)
{
std::ofstream outfile;
gctl::open_outfile(outfile, filename, ".msh");
gctl::save2gmsh(outfile, elements, nodes, packed);
outfile.close();
meshio_.init_file(filename, Output);
meshio_.set_packed(packed, Output);
meshio_.save_mesh(elements, nodes);
return;
}

View File

@ -69,6 +69,6 @@ namespace gctl
array<vertex3dc> nodes;
array<tetrahedron> elements;
};
}
};
#endif //_GCTL_TET_MESH_H

View File

@ -25,8 +25,8 @@
* Also add information on how to contact you by electronic and paper mail.
******************************************************/
#ifndef _GCTL_TRI_MESH_H
#define _GCTL_TRI_MESH_H
#ifndef _GCTL_TRI2D_MESH_H
#define _GCTL_TRI2D_MESH_H
#include "mesh.h"
@ -69,7 +69,7 @@ namespace gctl
void load_gmsh_groups();
/**
* @brief
* @brief
*
* @param datname
* @param phynames
@ -82,7 +82,9 @@ namespace gctl
array<triangle2d> elems_;
array<gmsh_physical_group> groups_;
_2i_vector elems_tag_;
};
}
#endif //_GCTL_TRI_MESH_H
mesh_io fio_;
};
};
#endif //_GCTL_TRI2D_MESH_H

254
lib/mesh/tri_mesh.cpp Normal file
View File

@ -0,0 +1,254 @@
/********************************************************
*
*
*
*
*
*
* Geophysical Computational Tools & Library (GCTL)
*
* Copyright (c) 2023 Yi Zhang (yizhang-geo@zju.edu.cn)
*
* GCTL is distributed under a dual licensing scheme. You can redistribute
* it and/or modify it under the terms of the GNU Lesser General Public
* License as published by the Free Software Foundation, either version 2
* of the License, or (at your option) any later version. You should have
* received a copy of the GNU Lesser General Public License along with this
* program. If not, see <http://www.gnu.org/licenses/>.
*
* If the terms and conditions of the LGPL v.2. would prevent you from using
* the GCTL, please consider the option to obtain a commercial license for a
* fee. These licenses are offered by the GCTL's original author. As a rule,
* licenses are provided "as-is", unlimited in time for a one time fee. Please
* send corresponding requests to: yizhang-geo@zju.edu.cn. Please do not forget
* to include some description of your company and the realm of its activities.
* Also add information on how to contact you by electronic and paper mail.
******************************************************/
#include "tri_mesh.h"
void gctl::triangle_mesh::init(std::string in_name, std::string in_info, const array<vertex3dc> &in_nodes_,
const array<triangle> &in_triangles)
{
check_initiated(true); // 检查是否已经初始化
base_mesh::init(TRI_TET_MESH, MESH_3D, in_name, in_info);
node_num_ = in_nodes_.size();
ele_num_ = in_triangles.size();
nodes_.resize(node_num_);
for (int i = 0; i < node_num_; i++)
{
nodes_[i] = in_nodes_[i];
}
elems_.resize(ele_num_);
for (int i = 0; i < ele_num_; i++)
{
elems_[i].id = i;
for (int j = 0; j < 3; j++)
{
elems_[i].vert[j] = nodes_.get(in_triangles[i].vert[j]->id);
}
}
initialized_ = true;
return;
}
void gctl::triangle_mesh::show_mesh_dimension(std::ostream &os) const
{
os << "Node num: " << node_num_ << std::endl;
os << "Elem num: " << ele_num_ << std::endl;
return;
}
void gctl::triangle_mesh::load_binary(std::string filename)
{
check_initiated(true); // 检查是否已经初始化
std::ifstream infile;
gctl::open_infile(infile, filename, ".2m", std::ios::in|std::ios::binary);
// 读入网格头信息
load_headinfo(infile, TRI_TET_MESH, MESH_3D);
// 读入网格信息
infile.read((char*)&node_num_, sizeof(int));
infile.read((char*)&ele_num_, sizeof(int));
nodes_.resize(node_num_);
elems_.resize(ele_num_);
for (int i = 0; i < node_num_; i++)
{
infile.read((char*)nodes_.get(i), sizeof(gctl::vertex3dc));
}
int in_index;
for (int i = 0; i < ele_num_; i++)
{
elems_[i].id = i;
for (int j = 0; j < 3; j++)
{
infile.read((char*)&in_index, sizeof(int));
elems_[i].vert[j] = nodes_.get(in_index);
}
}
// 读入模型数据单元
load_datablock(infile);
infile.close();
initialized_ = true;
return;
}
void gctl::triangle_mesh::save_binary(std::string filename)
{
check_initiated(); // 检查是否已经初始化
std::ofstream outfile;
gctl::open_outfile(outfile, filename, ".2m", std::ios::out|std::ios::binary);
// 首先输出网格的头信息
save_headinfo(outfile);
// 输出网格信息
outfile.write((char*)&node_num_, sizeof(int));
outfile.write((char*)&ele_num_, sizeof(int));
for (int i = 0; i < node_num_; i++)
{
outfile.write((char*)nodes_.get(i), sizeof(gctl::vertex3dc));
}
int in_index;
for (int i = 0; i < ele_num_; i++)
{
for (int j = 0; j < 3; j++)
{
in_index = elems_[i].vert[j]->id;
outfile.write((char*)&in_index, sizeof(int));
}
}
// 输出的模型数据单元
save_datablock(outfile);
outfile.close();
return;
}
gctl::triangle_mesh::triangle_mesh() : base_mesh::base_mesh(){}
gctl::triangle_mesh::triangle_mesh(std::string in_name, std::string in_info, const array<vertex3dc> &in_nodes_,
const array<triangle> &in_triangles)
{
init(in_name, in_info, in_nodes_, in_triangles);
}
gctl::triangle_mesh::~triangle_mesh(){}
const gctl::array<gctl::vertex3dc> &gctl::triangle_mesh::get_nodes() const
{
check_initiated(); // 检查是否已经初始化
return nodes_;
}
const gctl::array<gctl::triangle> &gctl::triangle_mesh::get_elements() const
{
check_initiated(); // 检查是否已经初始化
return elems_;
}
void gctl::triangle_mesh::load_gmsh(std::string filename, index_packed_e packed)
{
meshio_.init_file(filename, Input);
meshio_.set_packed(packed, Input);
meshio_.read_mesh(elems_, nodes_, &elems_tag_);
// 设置名称与信息等
node_num_ = nodes_.size();
ele_num_ = elems_.size();
base_mesh::init(TRI_TET_MESH, MESH_3D, filename, "Imported from a .msh file.");
initialized_ = true;
return;
}
void gctl::triangle_mesh::load_gmsh_groups()
{
check_initiated(); // 检查是否已经初始化
meshio_.read_physical_groups(groups_);
return;
}
void gctl::triangle_mesh::save_gmsh(std::string filename, index_packed_e packed)
{
meshio_.init_file(filename, Output);
meshio_.set_packed(packed, Output);
meshio_.save_mesh(elems_, nodes_);
return;
}
void gctl::triangle_mesh::groups2data(std::string datname, _1s_vector phynames, _1d_vector phyvals)
{
check_initiated(); // 检查是否已经初始化
if (phynames.size() != phyvals.size() || phynames.empty())
{
throw std::invalid_argument("The size of phynames and phyvals must be the same and non-empty.");
}
if (groups_.empty() || elems_tag_.empty())
{
throw std::invalid_argument("No physical groups found.");
}
if (elems_tag_.empty()) throw std::invalid_argument("No physical tags found.");
for (size_t i = 0; i < ele_num_; i++)
{
if (elems_tag_[i].empty()) throw std::invalid_argument("No physical tags found.");
}
array<int> tag_idx(phynames.size());
for (size_t i = 0; i < phynames.size(); i++)
{
tag_idx[i] = meshio_.physical_name2tag(groups_, phynames[i]);
}
if (saved(datname))
{
meshdata &dat = get_data(datname);
if (dat.loctype_ != ElemData) throw std::invalid_argument("The data type must be ElemData.");
for (size_t i = 0; i < ele_num_; i++)
{
for (size_t e = 0; e < phynames.size(); e++)
{
if (elems_tag_[i][0] == tag_idx[e])
{
dat.datval_[i] = phyvals[e];
break;
}
}
}
}
else
{
meshdata &newdat = add_data(ElemData, Scalar, datname, 0.0);
for (size_t i = 0; i < ele_num_; i++)
{
for (size_t e = 0; e < phynames.size(); e++)
{
if (elems_tag_[i][0] == tag_idx[e])
{
newdat.datval_[i] = phyvals[e];
break;
}
}
}
}
return;
}

85
lib/mesh/tri_mesh.h Normal file
View File

@ -0,0 +1,85 @@
/********************************************************
*
*
*
*
*
*
* Geophysical Computational Tools & Library (GCTL)
*
* Copyright (c) 2023 Yi Zhang (yizhang-geo@zju.edu.cn)
*
* GCTL is distributed under a dual licensing scheme. You can redistribute
* it and/or modify it under the terms of the GNU Lesser General Public
* License as published by the Free Software Foundation, either version 2
* of the License, or (at your option) any later version. You should have
* received a copy of the GNU Lesser General Public License along with this
* program. If not, see <http://www.gnu.org/licenses/>.
*
* If the terms and conditions of the LGPL v.2. would prevent you from using
* the GCTL, please consider the option to obtain a commercial license for a
* fee. These licenses are offered by the GCTL's original author. As a rule,
* licenses are provided "as-is", unlimited in time for a one time fee. Please
* send corresponding requests to: yizhang-geo@zju.edu.cn. Please do not forget
* to include some description of your company and the realm of its activities.
* Also add information on how to contact you by electronic and paper mail.
******************************************************/
#ifndef _GCTL_TRI_MESH_H
#define _GCTL_TRI_MESH_H
#include "mesh.h"
namespace gctl
{
class triangle_mesh : public base_mesh
{
public:
/**
* mesh类型的虚函数实现
*/
void load_gmsh(std::string filename, index_packed_e packed = Packed);
void save_gmsh(std::string filename, index_packed_e packed = Packed);
void load_binary(std::string filename);
void save_binary(std::string filename);
void show_mesh_dimension(std::ostream &os) const;
/**
* triangle_mesh的专有函数
*/
triangle_mesh();
triangle_mesh(std::string in_name, std::string in_info, const array<vertex3dc> &in_nodes,
const array<triangle> &in_triangles);
virtual ~triangle_mesh();
void init(std::string in_name, std::string in_info, const array<vertex3dc> &in_nodes,
const array<triangle> &in_triangles);
const array<vertex3dc> &get_nodes() const;
const array<triangle> &get_elements() const;
void load_gmsh_groups();
/**
* @brief
*
* @param datname
* @param phynames
* @param phyvals
*/
void groups2data(std::string datname, _1s_vector phynames, _1d_vector phyvals);
protected:
array<vertex3dc> nodes_;
array<triangle> elems_;
array<gmsh_physical_group> groups_;
_2i_vector elems_tag_;
};
};
#endif //_GCTL_TRI2D_MESH_H

View File

@ -443,7 +443,7 @@ void save_text(const std::vector<std::string> &cmd_units)
void data_cloud(const std::vector<std::string> &cmd_units)
{
// data-cloud \<new-data-name\> node|cell \<data-could-file\> \<search-x\> \<search-y\> \<search-deg\> [\<text_descriptor\>]
// data-cloud \<new-data-name\> node|cell \<data-could-file\> \<search-x\> \<search-y\> \<search-deg\>
if (cmd_units.size() < 7) throw std::runtime_error("data-cloud: insufficient parameters.");
gctl::array<std::string> copy_str(7, "null");
@ -456,34 +456,16 @@ void data_cloud(const std::vector<std::string> &cmd_units)
if (copy_str[1] == "node") d_type = NodeData;
else if (copy_str[1] == "cell") d_type = ElemData;
else throw std::runtime_error("open-surfer: invalid grid data type.");
text_descriptor desc;
desc.file_name_ = copy_str[2];
desc.file_ext_ = ".txt";
desc.col_str_ = "0,1,2";
desc.delimiter_ = ' ';
desc.att_sym_ = '#';
desc.head_num_ = 0;
if (copy_str[6] != "null")
{
parse_string_to_value(copy_str[6], '/', true, desc.col_str_,
desc.delimiter_, desc.att_sym_, desc.head_num_);
}
std::vector<double> posix_vec, posiy_vec, data_vec;
read_text2vectors(desc, posix_vec, posiy_vec, data_vec);
geodsv_io text_in;
text_in.load_text(copy_str[2]);
array<point2dc> posi_arr(posix_vec.size());
array<double> posi_val;
array<point2dc> posi_arr;
_1d_array data_vec;
text_in.get_column_point2dc(posi_arr, 0, 1);
text_in.get_column(data_vec, 2);
posi_val.input(data_vec);
for (size_t i = 0; i < posi_arr.size(); i++)
{
posi_arr[i].x = posix_vec[i];
posi_arr[i].y = posiy_vec[i];
}
rg.load_data_cloud(posi_arr, posi_val, atof(copy_str[3].c_str()), atof(copy_str[4].c_str()), atof(copy_str[5].c_str()), copy_str[0], d_type);
rg.load_data_cloud(posi_arr, data_vec, atof(copy_str[3].c_str()), atof(copy_str[4].c_str()), atof(copy_str[5].c_str()), copy_str[0], d_type);
return;
}
@ -530,31 +512,19 @@ void wavelet(const std::vector<std::string> &cmd_units)
void sum_data(const std::vector<std::string> &cmd_units)
{
// sum-data \<new-data-name\> \<data-name1\> \<data-name2\>
// sum-data \<new-data\> \<data1\> \<data2\>
if (cmd_units.size() < 4) throw std::runtime_error("sum: insufficient parameters.");
gctl::array<std::string> copy_str(3, "null");
for (size_t i = 1; i <= GCTL_MIN(cmd_units.size() - 1, 3); i++)
{
copy_str[i - 1] = cmd_units[i];
}
rg.sum(copy_str[0], copy_str[1], copy_str[2]);
rg.sum(cmd_units[1], cmd_units[2], cmd_units[3]);
return;
}
void diff_data(const std::vector<std::string> &cmd_units)
{
// sum-data \<new-data-name\> \<data-name1\> \<data-name2\>
// diff \<new-data\> \<data1\> \<data2\>
if (cmd_units.size() < 4) throw std::runtime_error("diff: insufficient parameters.");
gctl::array<std::string> copy_str(3, "null");
for (size_t i = 1; i <= GCTL_MIN(cmd_units.size() - 1, 3); i++)
{
copy_str[i - 1] = cmd_units[i];
}
rg.diff(copy_str[0], copy_str[1], copy_str[2]);
rg.diff(cmd_units[1], cmd_units[2], cmd_units[3]);
return;
}
@ -706,7 +676,7 @@ void get_stats(const std::vector<std::string> &cmd_units)
void get_profile(const std::vector<std::string> &cmd_units)
{
// profile \<data-name\> \<out-file\> {\<x\>,\<y\> \<x\>,\<y\> \<num\>}|{\<xy-file\> [\<file-format\>]}
// profile \<data-name\> \<out-file\> {\<x\>,\<y\> \<x\>,\<y\> \<num\>}|{\<xy-file\>}
if (cmd_units.size() < 4) throw std::runtime_error("profile: insufficient parameters.");
gctl::array<std::string> copy_str(5, "null");
@ -726,29 +696,9 @@ void get_profile(const std::vector<std::string> &cmd_units)
}
else // File exist
{
text_descriptor desc;
desc.file_name_ = copy_str[2];
desc.file_ext_ = ".txt";
desc.att_sym_ = '#';
desc.col_str_ = "0,1";
desc.delimiter_ = ' ';
desc.head_num_ = 0;
if (copy_str[3] != "null")
{
gctl::parse_string_to_value(copy_str[3], '/', true, desc.col_str_, desc.delimiter_,
desc.att_sym_, desc.head_num_);
}
std::vector<double> xs, ys;
gctl::read_text2vectors(desc, xs, ys);
xys.resize(xs.size());
for (size_t i = 0; i < xys.size(); i++)
{
xys[i].x = xs[i];
xys[i].y = ys[i];
}
geodsv_io fio;
fio.load_text(copy_str[2]);
fio.get_column_point2dc(xys, 0, 1);
rg.extract_points(copy_str[0], xys, p_data);
}

View File

@ -28,6 +28,8 @@
#ifndef GCTL_MESH_GRIDMANAGER_H
#define GCTL_MESH_GRIDMANAGER_H
#include "gctl/io/term_io.h"
#include "gctl/io/dsv_io.h"
#include "../../lib/mesh.h"
using namespace gctl;

View File

@ -22,8 +22,8 @@ Mathematic expresstions could be used to specify the grid ranges.
4. `save gmsh <file>` Save the grid using the Gmsh legacy format. This command will write all data that are allowed for outputing.
5. `save text <file>` Save the grid data to text file.
#### cloud \<new-data-name\> node|cell \<data-could-file\> \<search-x\> \<search-y\> \<search-deg\> [\<text_descriptor\>]
Import a randomly distributed data could (namely a xyz data). Note that a grid must exist before the use this command. 'node' or 'cell' indicates whether the input could data should be node or cell registed. \<data-could-file\> is the input data file. \<search-x\>,\<search-y\> and \<search-deg\> defines an oval of which the command is used to interpolate the could data to grid points. Additional options for reading the \<data-could-file\> could be given by the [\<text_descriptor\>] argument which has a format of \<order_str\>/\<delimiter\>/\<annotate\>/\<head_record\>. For example 0,1,2/,/#/0 means that the intput data could is stored int the first three columns separated by commas. The file has no head records and any line starts with '#' is an annotation.
#### cloud \<new-data-name\> node|cell \<data-could-file\> \<search-x\> \<search-y\> \<search-deg\>
Import a randomly distributed data could (namely a xyz data). Note that a grid must exist before the use this command. 'node' or 'cell' indicates whether the input could data should be node or cell registed. \<data-could-file\> is the input data file. \<search-x\>,\<search-y\> and \<search-deg\> defines an oval of which the command is used to interpolate the could data to grid points.
#### gradient \<data-name\> \<new-data-name\> dx|dy \<order\>
Calculate directional gradient(s) of the selected data. The new gradient data will have the same
@ -32,10 +32,10 @@ data type of the used data. 'dx' and 'dy' give the calculation directions. highe
#### wavelet \<data-name\> \<wavelet-name\> \<order\> [\<show-summary\>]
Using a wavelet \<wavelet-name\> to decompose the selecte data. Avaiable wavelets include 'bd2'... .The order of decomposition could be set using \<order\>. And set \<show-summary\> to yes to see the summary.
#### sum \<new-data-name\> \<data-name1\> \<data-name2\>
#### sum \<new-data\> \<data1\> \<data2\>
Calculate the sum of the two data. The arguments' format should be clear enough.
#### diff \<new-data-name\> \<data-name1\> \<data-name2\>
#### diff \<new-data\> \<data1\> \<data2\>
Calculate the difference of the two data. The arguments' format should be clear enough.
#### bool \<new-data-name\> \<data-name1> \<bool-data\> [reverse]