branch: devel-tetsuo-xml-binary commit a79688f46f6ce867615dc824bdf18b1e0154400b Author: Tetsuo Koyama <tkoyama...@gmail.com> AuthorDate: Sun May 31 12:02:00 2020 +0000
Add the scripting interface for vtu export --- interface/src/gf_mesh_fem_get.cc | 33 ++++++++++++++++++ interface/src/gf_mesh_get.cc | 25 ++++++++++++++ interface/src/gf_slice_get.cc | 72 ++++++++++++++++++++++++++++++++++++++++ src/getfem/getfem_export.h | 4 +-- src/getfem_export.cc | 16 ++++----- tests/helmholtz.cc | 3 +- 6 files changed, 141 insertions(+), 12 deletions(-) diff --git a/interface/src/gf_mesh_fem_get.cc b/interface/src/gf_mesh_fem_get.cc index 03abb0c..8a179b1 100644 --- a/interface/src/gf_mesh_fem_get.cc +++ b/interface/src/gf_mesh_fem_get.cc @@ -763,6 +763,39 @@ void gf_mesh_fem_get(getfemint::mexargs_in& m_in, } ); + /*@GET ('export to vtu',@str filename, ... ['ascii'], U, 'name'...) + Export a @tmf and some fields to a vtu file. + + The FEM and geometric transformations will be mapped to order 1 + or 2 isoparametric Pk (or Qk) FEMs (as VTK(XML) does not handle higher + order elements). If you need to represent high-order FEMs or + high-order geometric transformations, you should consider + SLICE:GET('export to vtu').@*/ + sub_command + ("export to vtu", 0, -1, 0, 0, + std::string fname = in.pop().to_string(); + bool ascii = false; + while (in.remaining() && in.front().is_string()) { + std::string cmd2 = in.pop().to_string(); + if (cmd_strmatch(cmd2, "ascii")) + ascii = true; + else THROW_BADARG("expecting 'ascii', got " << cmd2); + } + getfem::vtu_export exp(fname, ascii); + exp.exporting(*mf); + exp.write_mesh(); + int count = 1; + while (in.remaining()) { + const getfem::mesh_fem *mf2 = mf; + if (in.remaining() >= 2 && is_meshfem_object(in.front())) + mf2 = to_meshfem_object(in.pop()); + darray U = in.pop().to_darray(); + in.last_popped().check_trailing_dimension(int(mf2->nb_dof())); + exp.write_point_data(*mf2, U, get_vtk_dataset_name(in, count)); + count+=1; + } + ); + /*@GET ('export to dx',@str filename, ...['as', @str mesh_name][,'edges']['serie',@str serie_name][,'ascii'][,'append'], U, 'name'...) Export a @tmf and some fields to an OpenDX file. diff --git a/interface/src/gf_mesh_get.cc b/interface/src/gf_mesh_get.cc index 6c9687f..9fcf011 100644 --- a/interface/src/gf_mesh_get.cc +++ b/interface/src/gf_mesh_get.cc @@ -1247,6 +1247,31 @@ void gf_mesh_get(getfemint::mexargs_in& m_in, ); + /*@GET ('export to vtu', @str filename, ... [,'ascii'][,'quality']) + Exports a mesh to a VTK(XML) file . + + If 'quality' is specified, an estimation of the quality of each + convex will be written to the file. + + See also MESH_FEM:GET('export to vtu'), SLICE:GET('export to vtu').@*/ + sub_command + ("export to vtu", 1, 3, 0, 1, + bool write_q = false; bool ascii = false; + std::string fname = in.pop().to_string(); + while (in.remaining() && in.front().is_string()) { + std::string cmd2 = in.pop().to_string(); + if (cmd_strmatch(cmd2, "ascii")) + ascii = true; + else if (cmd_strmatch(cmd2, "quality")) + write_q = true; + else THROW_BADARG("expecting 'ascii' or 'quality', got " << cmd2); + } + getfem::vtu_export exp(fname, ascii); + exp.exporting(*pmesh); + exp.write_mesh(); if (write_q) exp.write_mesh_quality(*pmesh); + ); + + /*@GET ('export to dx', @str filename, ... [,'ascii'][,'append'][,'as',@str name,[,'serie',@str serie_name]][,'edges']) Exports a mesh to an OpenDX file. diff --git a/interface/src/gf_slice_get.cc b/interface/src/gf_slice_get.cc index 985072c..9a47f09 100644 --- a/interface/src/gf_slice_get.cc +++ b/interface/src/gf_slice_get.cc @@ -438,6 +438,78 @@ void gf_slice_get(getfemint::mexargs_in& m_in, ); + /*@GET ('export to vtu', @str filename, ...) + Export a slice to VTK(XML). + + Following the `filename`, you may use any of the following options: + + - if 'ascii' is not used, the file will contain binary data + (non portable, but fast). + - if 'edges' is used, the edges of the original mesh will be + written instead of the slice content. + + More than one dataset may be written, just list them. Each dataset + consists of either: + + - a field interpolated on the slice (scalar, vector or tensor), + followed by an optional name. + - a mesh_fem and a field, followed by an optional name. + + Examples: + + - SLICE:GET('export to vtu', 'test.vtu', Usl, 'first_dataset', mf, + U2, 'second_dataset') + - SLICE:GET('export to vtu', 'test.vtu', 'ascii', mf,U2) + - SLICE:GET('export to vtu', 'test.vtu', 'edges', 'ascii', Uslice)@*/ + sub_command + ("export to vtu",1, -1, 0, 0, + std::string fname = in.pop().to_string(); + bool ascii = false; + bool edges = false; + while (in.remaining() && in.front().is_string()) { + std::string cmd2 = in.pop().to_string(); + if (cmd_strmatch(cmd2, "ascii")) + ascii = true; + else if (cmd_strmatch(cmd2, "edges")) + edges = true; + else THROW_BADARG("expecting 'ascii' or 'edges', got " << cmd2); + } + getfem::vtu_export exp(fname, ascii); + getfem::stored_mesh_slice sl_edges; + const getfem::stored_mesh_slice *vtu_slice = sl; + getfem::mesh m_edges; + if (edges) { + vtu_slice = &sl_edges; + dal::bit_vector slice_edges; + getfem::mesh_slicer slicer(sl->linked_mesh()); + getfem::slicer_build_edges_mesh action(m_edges,slice_edges); + slicer.push_back_action(action); slicer.exec(*sl); + sl_edges.build(m_edges, getfem::slicer_none()); + } + exp.exporting(*vtu_slice); + exp.write_mesh(); + int count = 1; + if (in.remaining()) { + do { + if (in.remaining() >= 2 && is_meshfem_object(in.front())) { + const getfem::mesh_fem &mf = *to_meshfem_object(in.pop()); + + darray U = in.pop().to_darray(); + in.last_popped().check_trailing_dimension(int(mf.nb_dof())); + + exp.write_point_data(mf,U,get_vtk_dataset_name(in, count)); + } else if (in.remaining()) { + darray slU = in.pop().to_darray(); + in.last_popped().check_trailing_dimension(int(vtu_slice->nb_points())); + + exp.write_sliced_point_data(slU,get_vtk_dataset_name(in, count)); + } else THROW_BADARG("don't know what to do with this argument") + count+=1; + } while (in.remaining()); + } + ); + + /*@GET ('export to pov', @str filename) Export a the triangles of the slice to POV-RAY.@*/ sub_command diff --git a/src/getfem/getfem_export.h b/src/getfem/getfem_export.h index b4436cf..d5888ff 100644 --- a/src/getfem/getfem_export.h +++ b/src/getfem/getfem_export.h @@ -286,8 +286,8 @@ namespace getfem { class vtu_export : public vtk_export { public: - vtu_export(const std::string& fname) : vtk_export(fname, true, false) {} - vtu_export(std::ostream &os_) : vtk_export(os_, true, false) {} + vtu_export(const std::string& fname, bool ascii_) : vtk_export(fname, ascii_, false) {} + vtu_export(std::ostream &os_, bool ascii_) : vtk_export(os_, ascii_, false) {} }; /** @brief A (quite large) class for exportation of data to IBM OpenDX. diff --git a/src/getfem_export.cc b/src/getfem_export.cc index acfca88..63abc0c 100644 --- a/src/getfem_export.cc +++ b/src/getfem_export.cc @@ -189,7 +189,6 @@ namespace getfem } void vtk_export::init() { - GMM_ASSERT1(vtk || ascii, "vtu support only ascii format."); strcpy(header, "Exported by GetFEM"); psl = 0; dim_ = dim_type(-1); static int test_endian = 0x01234567; @@ -432,12 +431,12 @@ namespace getfem } for (size_type ic=0; ic < psl->nb_convex(); ++ic) { const getfem::mesh_slicer::cs_simplexes_ct& s = psl->simplexes(ic); - for (size_type i=0; i < s.size(); ++i) { + for (const auto &val : s) { if (vtk) { - write_val(int(s[i].dim()+1)); + write_val(int(val.dim()+1)); } - for (size_type j=0; j < s[i].dim()+1; ++j) - write_val(int(s[i].inodes[j] + nodes_cnt)); + for (size_type j=0; j < val.dim()+1; ++j) + write_val(int(val.inodes[j] + nodes_cnt)); write_separ(); } nodes_cnt += psl->nodes(ic).size(); @@ -452,14 +451,13 @@ namespace getfem int cnt = 0; for (size_type ic=0; ic < psl->nb_convex(); ++ic) { const getfem::mesh_slicer::cs_simplexes_ct& s = psl->simplexes(ic); - for (size_type i=0; i < s.size(); ++i) { + for (const auto &val : s) if (vtk) { - write_val(int(vtk_simplex_code[s[i].dim()])); + write_val(int(vtk_simplex_code[val.dim()])); } else { - cnt += int(s[i].dim()+1); + cnt += int(val.dim()+1); write_val(cnt); } - } write_separ(); splx_cnt -= s.size(); } diff --git a/tests/helmholtz.cc b/tests/helmholtz.cc index 6b4a8d1..6cf28ac 100644 --- a/tests/helmholtz.cc +++ b/tests/helmholtz.cc @@ -235,7 +235,8 @@ int main(int argc, char *argv[]) { getfem::vtk_export vtk_exp(p.datafilename + ".vtk", p.PARAM.int_value("VTK_EXPORT")==1, true); cout << "export to " << p.datafilename + ".vtu" << "..\n"; - getfem::vtu_export vtu_exp(p.datafilename + ".vtu"); + getfem::vtu_export vtu_exp(p.datafilename + ".vtu", + p.PARAM.int_value("VTK_EXPORT")==0); getfem::stored_mesh_slice sl(p.mesh, p.mesh.nb_convex() < 2000 ? 8 : 6); vtk_exp.exporting(sl); vtk_exp.write_point_data(p.mf_u, gmm::real_part(U), "helmholtz_rfield");