netgen/libsrc/geom2d/python_geom2d.cpp

277 lines
10 KiB
C++
Raw Normal View History

2014-09-10 17:16:09 +06:00
#ifdef NG_PYTHON
2014-10-08 21:48:48 +06:00
#include <../general/ngpython.hpp>
2014-09-10 17:16:09 +06:00
#include <meshing.hpp>
#include <geometry2d.hpp>
using namespace netgen;
2016-08-26 17:33:57 +05:00
2015-08-08 15:58:41 +05:00
namespace netgen
{
extern std::shared_ptr<NetgenGeometry> ng_geometry;
}
2016-11-04 16:14:52 +05:00
DLL_HEADER void ExportGeom2d(py::module &m)
2014-09-10 17:16:09 +06:00
{
2016-11-04 16:14:52 +05:00
py::class_<SplineGeometry2d, shared_ptr<SplineGeometry2d>>
(m, "SplineGeometry",
"a 2d boundary representation geometry model by lines and splines")
2016-11-04 16:14:52 +05:00
.def(py::init<>())
.def("__init__",
[](SplineGeometry2d *instance, const string & filename)
2015-08-08 15:58:41 +05:00
{
cout << "load geometry";
ifstream ist(filename);
2016-11-04 16:14:52 +05:00
new (instance) SplineGeometry2d();
instance->Load (filename.c_str());
ng_geometry = shared_ptr<SplineGeometry2d>(instance, NOOP_Deleter);
})
2015-08-08 15:58:41 +05:00
2014-09-10 17:16:09 +06:00
.def("Load",&SplineGeometry2d::Load)
.def("AppendPoint", FunctionPointer
2016-07-08 20:40:59 +05:00
([](SplineGeometry2d &self, double px, double py, double maxh, bool hpref)
{
Point<2> p;
p(0) = px;
p(1) = py;
GeomPoint<2> gp(p);
gp.hmax = maxh;
2016-07-08 20:40:59 +05:00
gp.hpref = hpref;
self.geompoints.Append(gp);
return self.geompoints.Size()-1;
}),
2016-11-04 16:14:52 +05:00
py::arg("x"), py::arg("y"), py::arg("maxh") = 1e99, py::arg("hpref")=false)
.def("Append", FunctionPointer([](SplineGeometry2d &self, py::list segment, int leftdomain, int rightdomain,
py::object bc, py::object copy, double maxh, bool hpref)
2014-09-10 17:16:09 +06:00
{
2016-11-04 16:14:52 +05:00
py::extract<std::string> segtype(segment[0]);
SplineSegExt * seg;
if (segtype().compare("line") == 0)
{
2016-11-04 16:14:52 +05:00
py::extract<int> point_index1(segment[1]);
py::extract<int> point_index2(segment[2]);
//point_index1.check()
LineSeg<2> * l = new LineSeg<2>(self.GetPoint(point_index1()), self.GetPoint(point_index2()));
seg = new SplineSegExt(*l);
}
else if (segtype().compare("spline3") == 0)
{
2016-11-04 16:14:52 +05:00
py::extract<int> point_index1(segment[1]);
py::extract<int> point_index2(segment[2]);
py::extract<int> point_index3(segment[3]);
SplineSeg3<2> * seg3 = new SplineSeg3<2>(self.GetPoint(point_index1()), self.GetPoint(point_index2()), self.GetPoint(point_index3()));
seg = new SplineSegExt(*seg3);
}
else
{
cout << "Appended segment is not a line or a spline3" << endl;
}
seg->leftdom = leftdomain;
seg->rightdom = rightdomain;
2015-11-01 16:07:26 +05:00
seg->hmax = maxh;
2016-07-08 20:40:59 +05:00
seg->hpref_left = hpref;
seg->hpref_right = hpref;
seg->reffak = 1;
seg->copyfrom = -1;
2016-11-04 16:14:52 +05:00
if (py::extract<int>(copy).check())
seg->copyfrom = py::extract<int>(copy)()+1;
2016-10-20 15:28:51 +05:00
2016-11-04 16:14:52 +05:00
if (py::extract<int>(bc).check())
seg->bc = py::extract<int>(bc)();
else if (py::extract<string>(bc).check())
{
2016-11-04 16:14:52 +05:00
string bcname = py::extract<string>(bc)();
int bcnum = self.GetBCNumber(bcname);
if (bcnum == 0)
bcnum = self.AddBCName(bcname);
seg->bc = bcnum;
}
else
seg->bc = self.GetNSplines()+1;
self.AppendSegment(seg);
2016-10-20 15:28:51 +05:00
return self.GetNSplines()-1;
2016-11-04 18:55:15 +05:00
}), py::arg("point_indices"), py::arg("leftdomain") = 1, py::arg("rightdomain") = py::int_(0),
2016-11-04 16:14:52 +05:00
py::arg("bc")=NGDummyArgument(), py::arg("copy")=NGDummyArgument(), py::arg("maxh")=1e99, py::arg("hpref")=false
)
2014-09-16 22:39:54 +06:00
2016-11-04 16:14:52 +05:00
.def("AppendSegment", FunctionPointer([](SplineGeometry2d &self, py::list point_indices, int leftdomain, int rightdomain)
{
2016-11-04 16:14:52 +05:00
int npts = py::len(point_indices);
2014-09-16 22:39:54 +06:00
SplineSegExt * seg;
2016-11-04 16:14:52 +05:00
//int a = py::extract<int>(point_indices[0]);
2014-09-16 22:39:54 +06:00
if (npts == 2)
{
2016-11-04 16:14:52 +05:00
LineSeg<2> * l = new LineSeg<2>(self.GetPoint(py::extract<int>(point_indices[0])()), self.GetPoint(py::extract<int>(point_indices[1])()));
2014-09-16 22:39:54 +06:00
seg = new SplineSegExt(*l);
}
else if (npts == 3)
{
2016-11-04 16:14:52 +05:00
SplineSeg3<2> * seg3 = new SplineSeg3<2>(self.GetPoint(py::extract<int>(point_indices[0])()), self.GetPoint(py::extract<int>(point_indices[1])()), self.GetPoint(py::extract<int>(point_indices[2])()));
2014-09-16 22:39:54 +06:00
seg = new SplineSegExt(*seg3);
}
seg->leftdom = leftdomain;
seg->rightdom = rightdomain;
2014-09-10 17:16:09 +06:00
seg->hmax = 1e99;
seg->reffak = 1;
seg->copyfrom = -1;
self.AppendSegment(seg);
2016-11-04 18:55:15 +05:00
}), py::arg("point_indices"), py::arg("leftdomain") = 1, py::arg("rightdomain") = py::int_(0))
2014-09-16 22:39:54 +06:00
//.def("AppendSegment", FunctionPointer([](SplineGeometry2d &self, int point_index1, int point_index2)//, int leftdomain, int rightdomain)
// {
// LineSeg<2> * l = new LineSeg<2>(self.GetPoint(point_index1), self.GetPoint(point_index2));
// SplineSegExt * seg = new SplineSegExt(*l);
// seg->leftdom = 1;// leftdomain;
// seg->rightdom = 0;// rightdomain;
// seg->hmax = 1e99;
// seg->reffak = 1;
// seg->copyfrom = -1;
// self.AppendSegment(seg);
2016-11-04 16:14:52 +05:00
// }))//, (py::arg("self"), py::arg("point_index1"), py::arg("point_index2"), py::arg("leftdomain") = 1, py::arg("rightdomain") = 0) )
2014-09-16 22:39:54 +06:00
//.def("AppendSegment", FunctionPointer([](SplineGeometry2d &self, int point_index1, int point_index2, int point_index3)//, int leftdomain, int rightdomain)
// {
// SplineSeg3<2> * seg3 = new SplineSeg3<2>(self.GetPoint(point_index1), self.GetPoint(point_index2), self.GetPoint(point_index3));
// SplineSegExt * seg = new SplineSegExt(*seg3);
// seg->leftdom = 1;// leftdomain;
// seg->rightdom = 0;// rightdomain;
// seg->hmax = 1e99;
// seg->reffak = 1;
// seg->copyfrom = -1;
// self.AppendSegment(seg);
2016-11-04 16:14:52 +05:00
// }))//, (py::arg("self"), py::arg("point_index1"), py::arg("point_index2"), py::arg("point_index3"), py::arg("leftdomain") = 1, py::arg("rightdomain") = 0 ) )
2015-09-12 18:02:56 +05:00
.def("SetMaterial", &SplineGeometry2d::SetMaterial)
.def("SetDomainMaxH", &SplineGeometry2d::SetDomainMaxh)
2015-09-12 18:02:56 +05:00
2014-09-10 17:16:09 +06:00
.def("PlotData", FunctionPointer([](SplineGeometry2d &self)
{
Box<2> box(self.GetBoundingBox());
double xdist = box.PMax()(0) - box.PMin()(0);
double ydist = box.PMax()(1) - box.PMin()(1);
2016-11-04 16:14:52 +05:00
py::tuple xlim = py::make_tuple(box.PMin()(0) - 0.1*xdist, box.PMax()(0) + 0.1*xdist);
py::tuple ylim = py::make_tuple(box.PMin()(1) - 0.1*ydist, box.PMax()(1) + 0.1*ydist);
2014-09-10 17:16:09 +06:00
2016-11-04 16:14:52 +05:00
py::list xpoints, ypoints;
2014-09-10 17:16:09 +06:00
for (int i = 0; i < self.splines.Size(); i++)
{
2016-11-04 16:14:52 +05:00
py::list xp, yp;
2014-09-10 17:16:09 +06:00
if (self.splines[i]->GetType().compare("line")==0)
{
2014-09-16 22:39:54 +06:00
GeomPoint<2> p1 = self.splines[i]->StartPI();
GeomPoint<2> p2 = self.splines[i]->EndPI();
2016-11-04 16:14:52 +05:00
xp.append(py::cast(p1(0)));
xp.append(py::cast(p2(0)));
yp.append(py::cast(p1(1)));
yp.append(py::cast(p2(1)));
2014-09-10 17:16:09 +06:00
}
else if (self.splines[i]->GetType().compare("spline3")==0)
{
double len = self.splines[i]->Length();
int n = floor(len/(0.05*min(xdist,ydist)));
2014-09-16 22:39:54 +06:00
for (int j = 0; j <= n; j++)
2014-09-10 17:16:09 +06:00
{
GeomPoint<2> point = self.splines[i]->GetPoint(j*1./n);
2016-11-04 16:14:52 +05:00
xp.append(py::cast(point(0)));
yp.append(py::cast(point(1)));
2014-09-10 17:16:09 +06:00
}
}
else
{
cout << "spline is neither line nor spline3" << endl;
}
xpoints.append(xp);
ypoints.append(yp);
2014-09-10 17:16:09 +06:00
}
2016-11-04 16:14:52 +05:00
return py::tuple(py::make_tuple(xlim, ylim, xpoints, ypoints));
2014-09-10 17:16:09 +06:00
}))
.def("PointData", FunctionPointer([](SplineGeometry2d &self)
{
2016-11-04 16:14:52 +05:00
py::list xpoints, ypoints, pointindex;
2014-09-10 17:16:09 +06:00
for (int i = 0; i < self.geompoints.Size(); i++)
{
2016-11-04 16:14:52 +05:00
pointindex.append(py::cast(i));
xpoints.append(py::cast(self.geompoints[i][0]));
ypoints.append(py::cast(self.geompoints[i][1]));
2014-09-10 17:16:09 +06:00
}
2016-11-04 16:14:52 +05:00
return py::tuple(py::make_tuple(xpoints, ypoints, pointindex));
2014-09-10 17:16:09 +06:00
}))
.def("SegmentData", FunctionPointer([](SplineGeometry2d &self)
{
2016-11-04 16:14:52 +05:00
py::list leftpoints, rightpoints, leftdom, rightdom;
2014-09-10 17:16:09 +06:00
for (int i = 0; i < self.splines.Size(); i++)
{
GeomPoint<2> point = self.splines[i]->GetPoint(0.5);
Vec<2> normal = self.GetSpline(i).GetTangent(0.5);
double temp = normal(0);
normal(0) = normal(1);
normal(1) = -temp;
2016-11-04 16:14:52 +05:00
leftdom.append(py::cast(self.GetSpline(i).leftdom));
rightdom.append(py::cast(self.GetSpline(i).rightdom));
2014-09-10 17:16:09 +06:00
2016-11-04 16:14:52 +05:00
rightpoints.append(py::make_tuple(point(0), point(1), normal(0)<0, normal(1)<0));
leftpoints.append(py::make_tuple(point(0), point(1), normal(0)<0, normal(1)<0));
2014-09-10 17:16:09 +06:00
}
2016-11-04 16:14:52 +05:00
return py::tuple(py::make_tuple(leftpoints, rightpoints, leftdom, rightdom));
2014-09-10 17:16:09 +06:00
}))
.def("Print", FunctionPointer([](SplineGeometry2d &self)
{
for (int i = 0; i < self.geompoints.Size(); i++)
{
cout << i << " : " << self.geompoints[i][0] << " , " << self.geompoints[i][1] << endl;
}
//Box<2> box(self.GetBoundingBox());
//cout << box.PMin() << endl;
//cout << box.PMax() << endl;
cout << self.splines.Size() << endl;
for (int i = 0; i < self.splines.Size(); i++)
{
cout << self.splines[i]->GetType() << endl;
//cout << i << " : " << self.splines[i]->GetPoint(0.1) << " , " << self.splines[i]->GetPoint(0.5) << endl;
}
}))
2015-08-08 15:58:41 +05:00
.def("GenerateMesh", FunctionPointer([](shared_ptr<SplineGeometry2d> self, MeshingParameters & mparam)
2014-09-10 17:16:09 +06:00
{
2015-08-08 15:58:41 +05:00
shared_ptr<Mesh> mesh = make_shared<Mesh> ();
2015-08-08 22:10:48 +05:00
mesh->SetGeometry(self);
SetGlobalMesh (mesh);
2015-08-08 15:58:41 +05:00
ng_geometry = self;
self->GenerateMesh(mesh, mparam);
2014-09-26 16:10:16 +06:00
return mesh;
2014-09-10 17:16:09 +06:00
}))
;
}
2016-11-04 16:14:52 +05:00
PYBIND11_PLUGIN(libgeom2d) {
py::module m("geom2d", "pybind geom2d");
ExportGeom2d(m);
return m.ptr();
2014-09-26 15:34:25 +06:00
}
2014-09-10 17:16:09 +06:00
#endif