diff --git a/config.h.in b/config.h.in index 1b506b31..a71eeea9 100644 --- a/config.h.in +++ b/config.h.in @@ -67,9 +67,6 @@ /* Define to the one symbol short name of this package. */ #undef PACKAGE_TARNAME -/* Define to the home page for this package. */ -#undef PACKAGE_URL - /* Define to the version of this package. */ #undef PACKAGE_VERSION diff --git a/libsrc/csg/Makefile.am b/libsrc/csg/Makefile.am index a09eb358..6502caee 100644 --- a/libsrc/csg/Makefile.am +++ b/libsrc/csg/Makefile.am @@ -8,7 +8,7 @@ revolution.hpp spline3d.hpp vscsg.hpp AM_CPPFLAGS = -I$(top_srcdir)/libsrc/include $(TCL_INCLUDES) METASOURCES = AUTO -noinst_LTLIBRARIES = libcsg.la libcsgvis.la +lib_LTLIBRARIES = libcsg.la libcsgvis.la libcsg_la_SOURCES = algprim.cpp brick.cpp \ @@ -18,6 +18,8 @@ manifold.cpp meshsurf.cpp polyhedra.cpp revolution.cpp singularref.cpp \ solid.cpp specpoin.cpp spline3d.cpp surface.cpp triapprox.cpp - libcsgvis_la_SOURCES = vscsg.cpp csgpkg.cpp + +libcsgvis_la_LIBADD = libcsg.la $(top_builddir)/libsrc/geom2d/libgeom2d.la + diff --git a/libsrc/csg/csgparser.cpp b/libsrc/csg/csgparser.cpp index b694cdc9..c2728bcd 100644 --- a/libsrc/csg/csgparser.cpp +++ b/libsrc/csg/csgparser.cpp @@ -7,9 +7,6 @@ namespace netgen { - //using namespace netgen; - - static kwstruct defkw[] = { { TOK_RECO, "algebraic3d" }, @@ -695,6 +692,65 @@ namespace netgen } + template + void LoadSpline (SplineGeometry & spline, CSGScanner & scan) + { + double hd; + Point x; + int nump, numseg; + + //scan.ReadNext(); + scan >> nump >> ';'; + + hd = 1; + spline.geompoints.SetSize(nump); + for(int i = 0; i> x(0) >> ',' >> x(1) >> ';'; + else if(D==3) + scan >> x(0) >> ',' >> x(1) >> ',' >> x(2) >> ';'; + + spline.geompoints[i] = GeomPoint(x,hd); + } + + scan >> numseg;// >> ';'; + + spline.splines.SetSize(numseg); + + int pnums,pnum1,pnum2,pnum3; + + + for(int i = 0; i> ';' >> pnums >> ','; + if (pnums == 2) + { + scan >> pnum1 >> ',' >> pnum2;// >> ';'; + spline.splines[i] = new LineSeg(spline.geompoints[pnum1-1], + spline.geompoints[pnum2-1]); + } + else if (pnums == 3) + { + scan >> pnum1 >> ',' >> pnum2 >> ',' + >> pnum3;// >> ';'; + spline.splines[i] = new SplineSeg3(spline.geompoints[pnum1-1], + spline.geompoints[pnum2-1], + spline.geompoints[pnum3-1]); + } + else if (pnums == 4) + { + scan >> pnum1 >> ',' >> pnum2 >> ',' + >> pnum3;// >> ';'; + spline.splines[i] = new CircleSeg(spline.geompoints[pnum1-1], + spline.geompoints[pnum2-1], + spline.geompoints[pnum3-1]); + } + } + } + + + void ParseFlags (CSGScanner & scan, Flags & flags) { @@ -1118,7 +1174,8 @@ namespace netgen ParseChar (scan, '('); SplineGeometry<2> * newspline = new SplineGeometry<2>; - newspline->CSGLoad(scan); + // newspline->CSGLoad(scan); + LoadSpline (*newspline, scan); ParseChar (scan, ')'); ParseChar (scan, ';'); @@ -1143,7 +1200,8 @@ namespace netgen ParseChar (scan, '('); SplineGeometry<3> * newspline = new SplineGeometry<3>; - newspline->CSGLoad(scan); + // newspline->CSGLoad(scan); + LoadSpline (*newspline, scan); ParseChar (scan, ')'); ParseChar (scan, ';'); diff --git a/libsrc/geom2d/Makefile.am b/libsrc/geom2d/Makefile.am index 761cc52a..db027238 100644 --- a/libsrc/geom2d/Makefile.am +++ b/libsrc/geom2d/Makefile.am @@ -10,4 +10,6 @@ libgeom2d_la_SOURCES = genmesh2d.cpp geom2dmesh.cpp spline.cpp \ libgeom2dvis_la_SOURCES = geom2dpkg.cpp vsgeom2d.cpp +libgeom2dvis_la_LIBADD = libgeom2d.la + diff --git a/libsrc/geom2d/genmesh2d.cpp b/libsrc/geom2d/genmesh2d.cpp index 5897a986..61550840 100644 --- a/libsrc/geom2d/genmesh2d.cpp +++ b/libsrc/geom2d/genmesh2d.cpp @@ -1,7 +1,7 @@ #include -#include +#include #include -#include "meshing.hpp" + namespace netgen { diff --git a/libsrc/geom2d/geom2dmesh.cpp b/libsrc/geom2d/geom2dmesh.cpp index e1d00206..acd84d00 100644 --- a/libsrc/geom2d/geom2dmesh.cpp +++ b/libsrc/geom2d/geom2dmesh.cpp @@ -1,9 +1,8 @@ #include - -#include -#include #include +#include + namespace netgen { diff --git a/libsrc/geom2d/geom2dpkg.cpp b/libsrc/geom2d/geom2dpkg.cpp index ab3746f6..6ceb7ccc 100644 --- a/libsrc/geom2d/geom2dpkg.cpp +++ b/libsrc/geom2d/geom2dpkg.cpp @@ -1,12 +1,13 @@ #include #include #include -#include #include -#include +#include +#include +#include #include "vsgeom2d.hpp" diff --git a/libsrc/geom2d/spline.cpp b/libsrc/geom2d/spline.cpp index 99758c8e..524c3107 100644 --- a/libsrc/geom2d/spline.cpp +++ b/libsrc/geom2d/spline.cpp @@ -5,9 +5,9 @@ Spline curve for Mesh generator */ #include -#include #include #include +#include namespace netgen { diff --git a/libsrc/geom2d/spline.hpp b/libsrc/geom2d/spline.hpp index 156d8d89..44ae3703 100644 --- a/libsrc/geom2d/spline.hpp +++ b/libsrc/geom2d/spline.hpp @@ -11,877 +11,876 @@ namespace netgen { -void CalcPartition (double l, double h, double h1, double h2, - double hcurve, double elto0, Array & points); + void CalcPartition (double l, double h, double h1, double h2, + double hcurve, double elto0, Array & points); -/* - Spline curves for 2D mesh generation -*/ + /* + Spline curves for 2D mesh generation + */ -/// Geometry point -template < int D > -class GeomPoint : public Point -{ -public: - /// refinement factor at point - double refatpoint; - /// max mesh-size at point - double hmax; - /// hp-refinement - bool hpref; - - /// - GeomPoint () { ; } - - /// - GeomPoint (const Point & ap, double aref = 1, bool ahpref=false) - : Point(ap), refatpoint(aref), hpref(ahpref) { ; } -}; - - - - -/// base class for 2d - segment -template < int D > -class SplineSeg -{ -public: - /// left domain - int leftdom; - /// right domain - int rightdom; - /// refinement at line - double reffak; - /// maximal h; - double hmax; - /// boundary condition number - int bc; - /// copy spline mesh from other spline (-1.. do not copy) - int copyfrom; - /// perfrom anisotropic refinement (hp-refinement) to edge - bool hpref_left; - /// perfrom anisotropic refinement (hp-refinement) to edge - bool hpref_right; - /// - int layer; - - - SplineSeg () + /// Geometry point + template < int D > + class GeomPoint : public Point { - layer = 1; + public: + /// refinement factor at point + double refatpoint; + /// max mesh-size at point + double hmax; + /// hp-refinement + bool hpref; + + /// + GeomPoint () { ; } + + /// + GeomPoint (const Point & ap, double aref = 1, bool ahpref=false) + : Point(ap), refatpoint(aref), hpref(ahpref) { ; } + }; + + + + + /// base class for 2d - segment + template < int D > + class SplineSeg + { + public: + /// left domain + int leftdom; + /// right domain + int rightdom; + /// refinement at line + double reffak; + /// maximal h; + double hmax; + /// boundary condition number + int bc; + /// copy spline mesh from other spline (-1.. do not copy) + int copyfrom; + /// perfrom anisotropic refinement (hp-refinement) to edge + bool hpref_left; + /// perfrom anisotropic refinement (hp-refinement) to edge + bool hpref_right; + /// + int layer; + + + SplineSeg () + { + layer = 1; + } + + /// calculates length of curve + virtual double Length () const; + /// returns point at curve, 0 <= t <= 1 + virtual Point GetPoint (double t) const = 0; + /// returns a (not necessarily unit-length) tangent vector for 0 <= t <= 1 + virtual Vec GetTangent (const double t) const + { cerr << "GetTangent not implemented for spline base-class" << endl; Vec dummy; return dummy;} + virtual void GetDerivatives (const double t, + Point & point, + Vec & first, + Vec & second) const {;} + /// partitionizes curve + void Partition (double h, double elto0, + Mesh & mesh, Point3dTree & searchtree, int segnr) const; + /// returns initial point on curve + virtual const GeomPoint & StartPI () const = 0; + /// returns terminal point on curve + virtual const GeomPoint & EndPI () const = 0; + /** writes curve description for fepp: + for implicitly given quadratic curves, the 6 coefficients of + the polynomial + $$ a x^2 + b y^2 + c x y + d x + e y + f = 0 $$ + are written to ost */ + void PrintCoeff (ostream & ost) const; + + virtual void GetCoeff (Vector & coeffs) const = 0; + + virtual void GetPoints (int n, Array > & points); + + /** calculates (2D) lineintersections: + for lines $$ a x + b y + c = 0 $$ the interecting points are calculated + and stored in points */ + virtual void LineIntersections (const double a, const double b, const double c, + Array < Point > & points, const double eps) const + {points.SetSize(0);} + + virtual double MaxCurvature(void) const = 0; + + virtual string GetType(void) const {return "splinebase";} + + virtual void Project (const Point point, Point & point_on_curve, double & t) const + { cerr << "Project not implemented for spline base-class" << endl;} + + virtual void GetRawData (Array & data) const + { cerr << "GetRawData not implemented for spline base-class" << endl;} + + }; + + + /// Straight line form p1 to p2 + template< int D > + class LineSeg : public SplineSeg + { + /// + GeomPoint p1, p2; + public: + /// + LineSeg (const GeomPoint & ap1, const GeomPoint & ap2); + /// + virtual double Length () const; + /// + inline virtual Point GetPoint (double t) const; + /// + virtual Vec GetTangent (const double t) const; + + + virtual void GetDerivatives (const double t, + Point & point, + Vec & first, + Vec & second) const; + /// + virtual const GeomPoint & StartPI () const { return p1; }; + /// + virtual const GeomPoint & EndPI () const { return p2; } + /// + virtual void GetCoeff (Vector & coeffs) const; + + virtual string GetType(void) const {return "line";} + + virtual void LineIntersections (const double a, const double b, const double c, + Array < Point > & points, const double eps) const; + + virtual double MaxCurvature(void) const {return 0;} + + virtual void Project (const Point point, Point & point_on_curve, double & t) const; + + virtual void GetRawData (Array & data) const; + }; + + + /// curve given by a rational, quadratic spline (including ellipses) + template< int D > + class SplineSeg3 : public SplineSeg + { + /// + GeomPoint p1, p2, p3; + + mutable double proj_latest_t; + public: + /// + SplineSeg3 (const GeomPoint & ap1, + const GeomPoint & ap2, + const GeomPoint & ap3); + /// + inline virtual Point GetPoint (double t) const; + /// + virtual Vec GetTangent (const double t) const; + + + virtual void GetDerivatives (const double t, + Point & point, + Vec & first, + Vec & second) const; + /// + virtual const GeomPoint & StartPI () const { return p1; }; + /// + virtual const GeomPoint & EndPI () const { return p3; } + /// + virtual void GetCoeff (Vector & coeffs) const; + + virtual string GetType(void) const {return "spline3";} + + const GeomPoint & TangentPoint (void) const { return p2; } + + virtual void LineIntersections (const double a, const double b, const double c, + Array < Point > & points, const double eps) const; + + virtual double MaxCurvature(void) const; + + virtual void Project (const Point point, Point & point_on_curve, double & t) const; + + virtual void GetRawData (Array & data) const; + }; + + + // Gundolf Haase 8/26/97 + /// A circle + template < int D > + class CircleSeg : public SplineSeg + { + /// + private: + GeomPoint p1, p2, p3; + //const GeomPoint &p1, &p2, &p3; + Point pm; + double radius, w1,w3; + public: + /// + CircleSeg (const GeomPoint & ap1, + const GeomPoint & ap2, + const GeomPoint & ap3); + /// + virtual Point GetPoint (double t) const; + /// + virtual const GeomPoint & StartPI () const { return p1; } + /// + virtual const GeomPoint & EndPI () const { return p3; } + /// + virtual void GetCoeff (Vector & coeffs) const; + /// + double Radius() const { return radius; } + /// + double StartAngle() const { return w1; } + /// + double EndAngle() const { return w3; } + /// + const Point & MidPoint(void) const {return pm; } + + virtual string GetType(void) const {return "circle";} + + virtual void LineIntersections (const double a, const double b, const double c, + Array < Point > & points, const double eps) const; + + virtual double MaxCurvature(void) const {return 1./radius;} + }; + + + + + + + /// + template + class DiscretePointsSeg : public SplineSeg + { + Array > pts; + GeomPoint p1n, p2n; + public: + /// + DiscretePointsSeg (const Array > & apts); + /// + virtual ~DiscretePointsSeg (); + /// + virtual Point GetPoint (double t) const; + /// + virtual const GeomPoint & StartPI () const { return p1n; }; + /// + virtual const GeomPoint & EndPI () const { return p2n; } + /// + virtual void GetCoeff (Vector & coeffs) const {;} + + virtual double MaxCurvature(void) const {return 1;} + }; + + + + + + + + // calculates length of spline-curve + template + double SplineSeg :: Length () const + { + int n = 100; + double dt = 1.0 / n; + + Point pold = GetPoint (0); + + double l = 0; + for (int i = 1; i <= n; i++) + { + Point p = GetPoint (i * dt); + l += Dist (p, pold); + pold = p; + } + + return l; } - /// calculates length of curve - virtual double Length () const; - /// returns point at curve, 0 <= t <= 1 - virtual Point GetPoint (double t) const = 0; - /// returns a (not necessarily unit-length) tangent vector for 0 <= t <= 1 - virtual Vec GetTangent (const double t) const - { cerr << "GetTangent not implemented for spline base-class" << endl; Vec dummy; return dummy;} - virtual void GetDerivatives (const double t, - Point & point, - Vec & first, - Vec & second) const {;} - /// partitionizes curve - void Partition (double h, double elto0, - Mesh & mesh, Point3dTree & searchtree, int segnr) const; - /// returns initial point on curve - virtual const GeomPoint & StartPI () const = 0; - /// returns terminal point on curve - virtual const GeomPoint & EndPI () const = 0; - /** writes curve description for fepp: - for implicitly given quadratic curves, the 6 coefficients of - the polynomial - $$ a x^2 + b y^2 + c x y + d x + e y + f = 0 $$ - are written to ost */ - void PrintCoeff (ostream & ost) const; - - virtual void GetCoeff (Vector & coeffs) const = 0; - - virtual void GetPoints (int n, Array > & points); - - /** calculates (2D) lineintersections: - for lines $$ a x + b y + c = 0 $$ the interecting points are calculated - and stored in points */ - virtual void LineIntersections (const double a, const double b, const double c, - Array < Point > & points, const double eps) const - {points.SetSize(0);} - - virtual double MaxCurvature(void) const = 0; - - virtual string GetType(void) const {return "splinebase";} - - virtual void Project (const Point point, Point & point_on_curve, double & t) const - { cerr << "Project not implemented for spline base-class" << endl;} - - virtual void GetRawData (Array & data) const - { cerr << "GetRawData not implemented for spline base-class" << endl;} - -}; -/// Straight line form p1 to p2 -template< int D > -class LineSeg : public SplineSeg -{ - /// - GeomPoint p1, p2; -public: - /// - LineSeg (const GeomPoint & ap1, const GeomPoint & ap2); - /// - virtual double Length () const; - /// - inline virtual Point GetPoint (double t) const; - /// - virtual Vec GetTangent (const double t) const; + // partitionizes spline curve + template + void SplineSeg :: Partition (double h, double elto0, + Mesh & mesh, Point3dTree & searchtree, int segnr) const + { + int i, j; + double l; // , r1, r2, ra; + double lold, dt, frac; + int n = 100; + Point p, pold, mark, oldmark; + Array curvepoints; + double edgelength, edgelengthold; + l = Length(); - - virtual void GetDerivatives (const double t, - Point & point, - Vec & first, - Vec & second) const; - /// - virtual const GeomPoint & StartPI () const { return p1; }; - /// - virtual const GeomPoint & EndPI () const { return p2; } - /// - virtual void GetCoeff (Vector & coeffs) const; - - virtual string GetType(void) const {return "line";} - - virtual void LineIntersections (const double a, const double b, const double c, - Array < Point > & points, const double eps) const; - - virtual double MaxCurvature(void) const {return 0;} - - virtual void Project (const Point point, Point & point_on_curve, double & t) const; - - virtual void GetRawData (Array & data) const; -}; + double h1 = min (StartPI().hmax, h/StartPI().refatpoint); + double h2 = min (EndPI().hmax, h/EndPI().refatpoint); + double hcurve = min (hmax, h/reffak); -/// curve given by a rational, quadratic spline (including ellipses) -template< int D > -class SplineSeg3 : public SplineSeg -{ - /// - GeomPoint p1, p2, p3; + CalcPartition (l, h, h1, h2, hcurve, elto0, curvepoints); + // cout << "curvepoints = " << curvepoints << endl; - mutable double proj_latest_t; -public: - /// - SplineSeg3 (const GeomPoint & ap1, - const GeomPoint & ap2, - const GeomPoint & ap3); - /// - inline virtual Point GetPoint (double t) const; - /// - virtual Vec GetTangent (const double t) const; + dt = 1.0 / n; - - virtual void GetDerivatives (const double t, - Point & point, - Vec & first, - Vec & second) const; - /// - virtual const GeomPoint & StartPI () const { return p1; }; - /// - virtual const GeomPoint & EndPI () const { return p3; } - /// - virtual void GetCoeff (Vector & coeffs) const; + l = 0; + j = 1; - virtual string GetType(void) const {return "spline3";} - - const GeomPoint & TangentPoint (void) const { return p2; } - - virtual void LineIntersections (const double a, const double b, const double c, - Array < Point > & points, const double eps) const; - - virtual double MaxCurvature(void) const; - - virtual void Project (const Point point, Point & point_on_curve, double & t) const; - - virtual void GetRawData (Array & data) const; -}; - - -// Gundolf Haase 8/26/97 -/// A circle -template < int D > -class CircleSeg : public SplineSeg -{ - /// -private: - GeomPoint p1, p2, p3; - //const GeomPoint &p1, &p2, &p3; - Point pm; - double radius, w1,w3; -public: - /// - CircleSeg (const GeomPoint & ap1, - const GeomPoint & ap2, - const GeomPoint & ap3); - /// - virtual Point GetPoint (double t) const; - /// - virtual const GeomPoint & StartPI () const { return p1; } - /// - virtual const GeomPoint & EndPI () const { return p3; } - /// - virtual void GetCoeff (Vector & coeffs) const; - /// - double Radius() const { return radius; } - /// - double StartAngle() const { return w1; } - /// - double EndAngle() const { return w3; } - /// - const Point & MidPoint(void) const {return pm; } - - virtual string GetType(void) const {return "circle";} - - virtual void LineIntersections (const double a, const double b, const double c, - Array < Point > & points, const double eps) const; - - virtual double MaxCurvature(void) const {return 1./radius;} -}; - - - - - - -/// -template -class DiscretePointsSeg : public SplineSeg -{ - Array > pts; - GeomPoint p1n, p2n; -public: - /// - DiscretePointsSeg (const Array > & apts); - /// - virtual ~DiscretePointsSeg (); - /// - virtual Point GetPoint (double t) const; - /// - virtual const GeomPoint & StartPI () const { return p1n; }; - /// - virtual const GeomPoint & EndPI () const { return p2n; } - /// - virtual void GetCoeff (Vector & coeffs) const {;} - - virtual double MaxCurvature(void) const {return 1;} -}; - - - - - - - -// calculates length of spline-curve -template -double SplineSeg :: Length () const -{ - Point p, pold; - - int i, n = 100; - double dt = 1.0 / n; - - pold = GetPoint (0); - - double l = 0; - for (i = 1; i <= n; i++) - { - p = GetPoint (i * dt); - l += Dist (p, pold); - pold = p; - } - return l; -} - - - -// partitionizes spline curve -template -void SplineSeg :: Partition (double h, double elto0, - Mesh & mesh, Point3dTree & searchtree, int segnr) const -{ - int i, j; - double l; // , r1, r2, ra; - double lold, dt, frac; - int n = 100; - Point p, pold, mark, oldmark; - Array curvepoints; - double edgelength, edgelengthold; - l = Length(); - - double h1 = min (StartPI().hmax, h/StartPI().refatpoint); - double h2 = min (EndPI().hmax, h/EndPI().refatpoint); - double hcurve = min (hmax, h/reffak); - - - CalcPartition (l, h, h1, h2, hcurve, elto0, curvepoints); - // cout << "curvepoints = " << curvepoints << endl; - - dt = 1.0 / n; - - l = 0; - j = 1; - - pold = GetPoint (0); - lold = 0; - oldmark = pold; - edgelengthold = 0; - Array locsearch; - - for (i = 1; i <= n; i++) - { - p = GetPoint (i*dt); - l = lold + Dist (p, pold); - while (j < curvepoints.Size() && (l >= curvepoints[j] || i == n)) - { - frac = (curvepoints[j]-lold) / (l-lold); - edgelength = i*dt + (frac-1)*dt; - // mark = pold + frac * (p-pold); - mark = GetPoint (edgelength); - - // cout << "mark = " << mark << " =?= " << GetPoint (edgelength) << endl; + pold = GetPoint (0); + lold = 0; + oldmark = pold; + edgelengthold = 0; + Array locsearch; + for (i = 1; i <= n; i++) + { + p = GetPoint (i*dt); + l = lold + Dist (p, pold); + while (j < curvepoints.Size() && (l >= curvepoints[j] || i == n)) { - PointIndex pi1 = -1, pi2 = -1; + frac = (curvepoints[j]-lold) / (l-lold); + edgelength = i*dt + (frac-1)*dt; + // mark = pold + frac * (p-pold); + mark = GetPoint (edgelength); - Point3d mark3(mark(0), mark(1), 0); - Point3d oldmark3(oldmark(0), oldmark(1), 0); + // cout << "mark = " << mark << " =?= " << GetPoint (edgelength) << endl; - Vec<3> v (1e-4*h, 1e-4*h, 1e-4*h); - searchtree.GetIntersecting (oldmark3 - v, oldmark3 + v, locsearch); + { + PointIndex pi1 = -1, pi2 = -1; + + Point3d mark3(mark(0), mark(1), 0); + Point3d oldmark3(oldmark(0), oldmark(1), 0); - for (int k = 0; k < locsearch.Size(); k++) - if ( mesh[PointIndex(locsearch[k])].GetLayer() == layer) - pi1 = locsearch[k]; - // if (locsearch.Size()) pi1 = locsearch[0]; + Vec<3> v (1e-4*h, 1e-4*h, 1e-4*h); + searchtree.GetIntersecting (oldmark3 - v, oldmark3 + v, locsearch); + + for (int k = 0; k < locsearch.Size(); k++) + if ( mesh[PointIndex(locsearch[k])].GetLayer() == layer) + pi1 = locsearch[k]; + // if (locsearch.Size()) pi1 = locsearch[0]; - searchtree.GetIntersecting (mark3 - v, mark3 + v, locsearch); - for (int k = 0; k < locsearch.Size(); k++) - if ( mesh[PointIndex(locsearch[k])].GetLayer() == layer) - pi2 = locsearch[k]; - // if (locsearch.Size()) pi2 = locsearch[0]; + searchtree.GetIntersecting (mark3 - v, mark3 + v, locsearch); + for (int k = 0; k < locsearch.Size(); k++) + if ( mesh[PointIndex(locsearch[k])].GetLayer() == layer) + pi2 = locsearch[k]; + // if (locsearch.Size()) pi2 = locsearch[0]; - /* - for (PointIndex pk = PointIndex::BASE; - pk < mesh.GetNP()+PointIndex::BASE; pk++) - { - if (Dist (mesh[pk], oldmark3) < 1e-4 * h) pi1 = pk; - if (Dist (mesh[pk], mark3) < 1e-4 * h) pi2 = pk; - } - */ + /* + for (PointIndex pk = PointIndex::BASE; + pk < mesh.GetNP()+PointIndex::BASE; pk++) + { + if (Dist (mesh[pk], oldmark3) < 1e-4 * h) pi1 = pk; + if (Dist (mesh[pk], mark3) < 1e-4 * h) pi2 = pk; + } + */ - // cout << "pi1 = " << pi1 << endl; - // cout << "pi2 = " << pi2 << endl; + // cout << "pi1 = " << pi1 << endl; + // cout << "pi2 = " << pi2 << endl; - if (pi1 == -1) - { - pi1 = mesh.AddPoint(oldmark3, layer); - searchtree.Insert (oldmark3, pi1); - } - if (pi2 == -1) - { - pi2 = mesh.AddPoint(mark3, layer); - searchtree.Insert (mark3, pi2); - } + if (pi1 == -1) + { + pi1 = mesh.AddPoint(oldmark3, layer); + searchtree.Insert (oldmark3, pi1); + } + if (pi2 == -1) + { + pi2 = mesh.AddPoint(mark3, layer); + searchtree.Insert (mark3, pi2); + } - Segment seg; - seg.edgenr = segnr; - seg.si = bc; // segnr; - seg[0] = pi1; - seg[1] = pi2; - seg.domin = leftdom; - seg.domout = rightdom; - seg.epgeominfo[0].edgenr = segnr; - seg.epgeominfo[0].dist = edgelengthold; - seg.epgeominfo[1].edgenr = segnr; - seg.epgeominfo[1].dist = edgelength; - seg.singedge_left = hpref_left; - seg.singedge_right = hpref_right; - mesh.AddSegment (seg); - } + Segment seg; + seg.edgenr = segnr; + seg.si = bc; // segnr; + seg[0] = pi1; + seg[1] = pi2; + seg.domin = leftdom; + seg.domout = rightdom; + seg.epgeominfo[0].edgenr = segnr; + seg.epgeominfo[0].dist = edgelengthold; + seg.epgeominfo[1].edgenr = segnr; + seg.epgeominfo[1].dist = edgelength; + seg.singedge_left = hpref_left; + seg.singedge_right = hpref_right; + mesh.AddSegment (seg); + } - oldmark = mark; - edgelengthold = edgelength; - j++; - } + oldmark = mark; + edgelengthold = edgelength; + j++; + } - pold = p; - lold = l; + pold = p; + lold = l; + } + } + + + template + void SplineSeg :: GetPoints (int n, Array > & points) + { + points.SetSize (n); + if (n >= 2) + for (int i = 0; i < n; i++) + points[i] = GetPoint(double(i) / (n-1)); + } + + + template + void SplineSeg :: PrintCoeff (ostream & ost) const + { + Vector u(6); + + GetCoeff(u); + + for ( int i=0; i<6; i++) + ost << u[i] << " "; + ost << endl; + } + + + + /* + Implementation of line-segment from p1 to p2 + */ + + + template + LineSeg :: LineSeg (const GeomPoint & ap1, + const GeomPoint & ap2) + : p1(ap1), p2(ap2) + { + ; + } + + + template + inline Point LineSeg :: GetPoint (double t) const + { + return p1 + t * (p2 - p1); + } + + template + Vec LineSeg :: GetTangent (const double t) const + { + return p2-p1; + } + + template + void LineSeg :: GetDerivatives (const double t, + Point & point, + Vec & first, + Vec & second) const + { + first = p2 - p1; + point = p1 + t * first; + second = 0; + } + + + template + double LineSeg :: Length () const + { + return Dist (p1, p2); + } + + + template + void LineSeg :: GetCoeff (Vector & coeffs) const + { + coeffs.SetSize(6); + + double dx = p2(0) - p1(0); + double dy = p2(1) - p1(1); + + coeffs[0] = coeffs[1] = coeffs[2] = 0; + coeffs[3] = -dy; + coeffs[4] = dx; + coeffs[5] = -dx * p1(1) + dy * p1(0); + } + + + + template + void LineSeg :: LineIntersections (const double a, const double b, const double c, + Array < Point > & points, const double eps) const + { + points.SetSize(0); + + double denom = -a*p2(0)+a*p1(0)-b*p2(1)+b*p1(1); + if(fabs(denom) < 1e-20) + return; + + double t = (a*p1(0)+b*p1(1)+c)/denom; + if((t > -eps) && (t < 1.+eps)) + points.Append(GetPoint(t)); + } + + + + template + void LineSeg :: Project (const Point point, Point & point_on_curve, double & t) const + { + Vec v = p2-p1; + double l = v.Length(); + v *= 1./l; + t = (point-p1)*v; + + if(t<0) t = 0; + if(t>l) t = l; + + point_on_curve = p1+t*v; + + t *= 1./l; + } + + + template + void LineSeg :: GetRawData (Array & data) const + { + data.Append(2); + for(int i=0; i + SplineSeg3 :: SplineSeg3 (const GeomPoint & ap1, + const GeomPoint & ap2, + const GeomPoint & ap3) + : p1(ap1), p2(ap2), p3(ap3) + { + proj_latest_t = 0.5; + } + + template + inline Point SplineSeg3 :: GetPoint (double t) const + { + double x, y, w; + double b1, b2, b3; + + b1 = (1-t)*(1-t); + b2 = sqrt(2.0) * t * (1-t); + b3 = t * t; + + x = p1(0) * b1 + p2(0) * b2 + p3(0) * b3; + y = p1(1) * b1 + p2(1) * b2 + p3(1) * b3; + w = b1 + b2 + b3; + + if(D==3) + { + double z = p1(2) * b1 + p2(2) * b2 + p3(2) * b3; + return Point (x/w, y/w, z/w); + } + else + return Point (x/w, y/w); + } + + + + + template + Vec SplineSeg3 :: GetTangent (const double t) const + { + const double b1 = (1.-t)*((sqrt(2.)-2.)*t-sqrt(2.)); + const double b2 = sqrt(2.)*(1.-2.*t); + const double b3 = t*((sqrt(2.)-2)*t+2.); + + + Vec retval; + for(int i=0; i + void SplineSeg3 :: GetCoeff (Vector & u) const + { + DenseMatrix a(6, 6); + DenseMatrix ata(6, 6); + Vector f(6); + + u.SetSize(6); + + // ata.SetSymmetric(1); + + double t = 0; + for (int i = 0; i < 5; i++, t += 0.25) + { + Point p = GetPoint (t); + a(i, 0) = p(0) * p(0); + a(i, 1) = p(1) * p(1); + a(i, 2) = p(0) * p(1); + a(i, 3) = p(0); + a(i, 4) = p(1); + a(i, 5) = 1; + } + a(5, 0) = 1; + + CalcAtA (a, ata); + + u = 0; + u(5) = 1; + a.MultTrans (u, f); + ata.Solve (f, u); + + // the sign + Point p0 = GetPoint(0); + Vec ht = GetTangent(0); + Vec<2> tang(ht(0), ht(1)); + + double gradx = 2.*u(0)*p0(0) + u(2)*p0(1) + u(3); + double grady = 2.*u(1)*p0(1) + u(2)*p0(0) + u(4); + Vec<2> gradn (grady, -gradx); + + if (tang * gradn < 0) u *= -1; + } + + /* + template + double SplineSeg3 :: MaxCurvature(void) const + { + Vec v1 = p1-p2; + Vec v2 = p3-p2; + double l1 = v1.Length(); + double l2 = v2.Length(); + (*testout) << "v1 " << v1 << " v2 " << v2 << endl; + + double cosalpha = v1*v2/(l1*l2); + + (*testout) << "cosalpha " << cosalpha << endl; + + return sqrt(cosalpha + 1.)/(min2(l1,l2)*(1.-cosalpha)); } -} + */ + + + template + void SplineSeg3 :: LineIntersections (const double a, const double b, const double c, + Array < Point > & points, const double eps) const + { + points.SetSize(0); + + double t; + + const double c1 = a*p1(0) - sqrt(2.)*a*p2(0) + a*p3(0) + + b*p1(1) - sqrt(2.)*b*p2(1) + b*p3(1) + + (2.-sqrt(2.))*c; + const double c2 = -2.*a*p1(0) + sqrt(2.)*a*p2(0) -2.*b*p1(1) + sqrt(2.)*b*p2(1) + (sqrt(2.)-2.)*c; + const double c3 = a*p1(0) + b*p1(1) + c; + + if(fabs(c1) < 1e-20) + { + if(fabs(c2) < 1e-20) + return; + + t = -c3/c2; + if((t > -eps) && (t < 1.+eps)) + points.Append(GetPoint(t)); + return; + } + + const double discr = c2*c2-4.*c1*c3; + + if(discr < 0) + return; + + if(fabs(discr/(c1*c1)) < 1e-14) + { + t = -0.5*c2/c1; + if((t > -eps) && (t < 1.+eps)) + points.Append(GetPoint(t)); + return; + } + + t = (-c2 + sqrt(discr))/(2.*c1); + if((t > -eps) && (t < 1.+eps)) + points.Append(GetPoint(t)); + + t = (-c2 - sqrt(discr))/(2.*c1); + if((t > -eps) && (t < 1.+eps)) + points.Append(GetPoint(t)); + } -template -void SplineSeg :: GetPoints (int n, Array > & points) -{ - points.SetSize (n); - if (n >= 2) - for (int i = 0; i < n; i++) - points[i] = GetPoint(double(i) / (n-1)); -} + template < int D > + void SplineSeg3 :: GetRawData (Array & data) const + { + data.Append(3); + for(int i=0; i -void SplineSeg :: PrintCoeff (ostream & ost) const -{ - Vector u(6); + //######################################################################## + // circlesegment - GetCoeff(u); - - for ( int i=0; i<6; i++) - ost << u[i] << " "; - ost << endl; -} - - - -/* - Implementation of line-segment from p1 to p2 -*/ - - -template -LineSeg :: LineSeg (const GeomPoint & ap1, - const GeomPoint & ap2) - : p1(ap1), p2(ap2) -{ - ; -} - - -template -inline Point LineSeg :: GetPoint (double t) const -{ - return p1 + t * (p2 - p1); -} - -template -Vec LineSeg :: GetTangent (const double t) const -{ - return p2-p1; -} - -template -void LineSeg :: GetDerivatives (const double t, - Point & point, - Vec & first, - Vec & second) const -{ - first = p2 - p1; - point = p1 + t * first; - second = 0; -} - - -template -double LineSeg :: Length () const -{ - return Dist (p1, p2); -} - - -template -void LineSeg :: GetCoeff (Vector & coeffs) const -{ - coeffs.SetSize(6); - - double dx = p2(0) - p1(0); - double dy = p2(1) - p1(1); - - coeffs[0] = coeffs[1] = coeffs[2] = 0; - coeffs[3] = -dy; - coeffs[4] = dx; - coeffs[5] = -dx * p1(1) + dy * p1(0); -} - - - -template -void LineSeg :: LineIntersections (const double a, const double b, const double c, - Array < Point > & points, const double eps) const -{ - points.SetSize(0); - - double denom = -a*p2(0)+a*p1(0)-b*p2(1)+b*p1(1); - if(fabs(denom) < 1e-20) - return; - - double t = (a*p1(0)+b*p1(1)+c)/denom; - if((t > -eps) && (t < 1.+eps)) - points.Append(GetPoint(t)); -} - - - -template -void LineSeg :: Project (const Point point, Point & point_on_curve, double & t) const -{ - Vec v = p2-p1; - double l = v.Length(); - v *= 1./l; - t = (point-p1)*v; - - if(t<0) t = 0; - if(t>l) t = l; - - point_on_curve = p1+t*v; - - t *= 1./l; -} - - -template -void LineSeg :: GetRawData (Array & data) const -{ - data.Append(2); - for(int i=0; i -SplineSeg3 :: SplineSeg3 (const GeomPoint & ap1, + template + CircleSeg :: CircleSeg (const GeomPoint & ap1, const GeomPoint & ap2, const GeomPoint & ap3) - : p1(ap1), p2(ap2), p3(ap3) -{ - proj_latest_t = 0.5; -} - -template -inline Point SplineSeg3 :: GetPoint (double t) const -{ - double x, y, w; - double b1, b2, b3; - - b1 = (1-t)*(1-t); - b2 = sqrt(2.0) * t * (1-t); - b3 = t * t; - - x = p1(0) * b1 + p2(0) * b2 + p3(0) * b3; - y = p1(1) * b1 + p2(1) * b2 + p3(1) * b3; - w = b1 + b2 + b3; - - if(D==3) - { - double z = p1(2) * b1 + p2(2) * b2 + p3(2) * b3; - return Point (x/w, y/w, z/w); - } - else - return Point (x/w, y/w); -} - - - - -template -Vec SplineSeg3 :: GetTangent (const double t) const -{ - const double b1 = (1.-t)*((sqrt(2.)-2.)*t-sqrt(2.)); - const double b2 = sqrt(2.)*(1.-2.*t); - const double b3 = t*((sqrt(2.)-2)*t+2.); - - - Vec retval; - for(int i=0; i -void SplineSeg3 :: GetCoeff (Vector & u) const -{ - DenseMatrix a(6, 6); - DenseMatrix ata(6, 6); - Vector f(6); - - u.SetSize(6); - - // ata.SetSymmetric(1); - - double t = 0; - for (int i = 0; i < 5; i++, t += 0.25) - { - Point p = GetPoint (t); - a(i, 0) = p(0) * p(0); - a(i, 1) = p(1) * p(1); - a(i, 2) = p(0) * p(1); - a(i, 3) = p(0); - a(i, 4) = p(1); - a(i, 5) = 1; - } - a(5, 0) = 1; - - CalcAtA (a, ata); - - u = 0; - u(5) = 1; - a.MultTrans (u, f); - ata.Solve (f, u); - - // the sign - Point p0 = GetPoint(0); - Vec ht = GetTangent(0); - Vec<2> tang(ht(0), ht(1)); - - double gradx = 2.*u(0)*p0(0) + u(2)*p0(1) + u(3); - double grady = 2.*u(1)*p0(1) + u(2)*p0(0) + u(4); - Vec<2> gradn (grady, -gradx); + : p1(ap1), p2(ap2), p3(ap3) + { + Vec v1,v2; - if (tang * gradn < 0) u *= -1; -} - -/* -template -double SplineSeg3 :: MaxCurvature(void) const -{ - Vec v1 = p1-p2; - Vec v2 = p3-p2; - double l1 = v1.Length(); - double l2 = v2.Length(); - (*testout) << "v1 " << v1 << " v2 " << v2 << endl; - - double cosalpha = v1*v2/(l1*l2); - - (*testout) << "cosalpha " << cosalpha << endl; - - return sqrt(cosalpha + 1.)/(min2(l1,l2)*(1.-cosalpha)); -} -*/ + v1 = p1 - p2; + v2 = p3 - p2; + Point p1t(p1+v1); + Point p2t(p3+v2); -template -void SplineSeg3 :: LineIntersections (const double a, const double b, const double c, - Array < Point > & points, const double eps) const -{ - points.SetSize(0); - - double t; - - const double c1 = a*p1(0) - sqrt(2.)*a*p2(0) + a*p3(0) - + b*p1(1) - sqrt(2.)*b*p2(1) + b*p3(1) - + (2.-sqrt(2.))*c; - const double c2 = -2.*a*p1(0) + sqrt(2.)*a*p2(0) -2.*b*p1(1) + sqrt(2.)*b*p2(1) + (sqrt(2.)-2.)*c; - const double c3 = a*p1(0) + b*p1(1) + c; - - if(fabs(c1) < 1e-20) - { - if(fabs(c2) < 1e-20) - return; - - t = -c3/c2; - if((t > -eps) && (t < 1.+eps)) - points.Append(GetPoint(t)); - return; - } - - const double discr = c2*c2-4.*c1*c3; - - if(discr < 0) - return; - - if(fabs(discr/(c1*c1)) < 1e-14) - { - t = -0.5*c2/c1; - if((t > -eps) && (t < 1.+eps)) - points.Append(GetPoint(t)); - return; - } - - t = (-c2 + sqrt(discr))/(2.*c1); - if((t > -eps) && (t < 1.+eps)) - points.Append(GetPoint(t)); - - t = (-c2 - sqrt(discr))/(2.*c1); - if((t > -eps) && (t < 1.+eps)) - points.Append(GetPoint(t)); -} - - -template < int D > -void SplineSeg3 :: GetRawData (Array & data) const -{ - data.Append(3); - for(int i=0; i -CircleSeg :: CircleSeg (const GeomPoint & ap1, - const GeomPoint & ap2, - const GeomPoint & ap3) - : p1(ap1), p2(ap2), p3(ap3) -{ - Vec v1,v2; - - v1 = p1 - p2; - v2 = p3 - p2; - - Point p1t(p1+v1); - Point p2t(p3+v2); - - // works only in 2D!!!!!!!!! + // works only in 2D!!!!!!!!! - Line2d g1t,g2t; + Line2d g1t,g2t; - g1t.P1() = Point<2>(p1(0),p1(1)); - g1t.P2() = Point<2>(p1t(0),p1t(1)); - g2t.P1() = Point<2>(p3(0),p3(1)); - g2t.P2() = Point<2>(p2t(0),p2t(1)); + g1t.P1() = Point<2>(p1(0),p1(1)); + g1t.P2() = Point<2>(p1t(0),p1t(1)); + g2t.P1() = Point<2>(p3(0),p3(1)); + g2t.P2() = Point<2>(p2t(0),p2t(1)); - Point<2> mp = CrossPoint (g1t,g2t); + Point<2> mp = CrossPoint (g1t,g2t); - pm(0) = mp(0); pm(1) = mp(1); - radius = Dist(pm,StartPI()); - Vec2d auxv; - auxv.X() = p1(0)-pm(0); auxv.Y() = p1(1)-pm(1); - w1 = Angle(auxv); - auxv.X() = p3(0)-pm(0); auxv.Y() = p3(1)-pm(1); - w3 = Angle(auxv); - if ( fabs(w3-w1) > M_PI ) - { - if ( w3>M_PI ) w3 -= 2*M_PI; - if ( w1>M_PI ) w1 -= 2*M_PI; - } -} + pm(0) = mp(0); pm(1) = mp(1); + radius = Dist(pm,StartPI()); + Vec2d auxv; + auxv.X() = p1(0)-pm(0); auxv.Y() = p1(1)-pm(1); + w1 = Angle(auxv); + auxv.X() = p3(0)-pm(0); auxv.Y() = p3(1)-pm(1); + w3 = Angle(auxv); + if ( fabs(w3-w1) > M_PI ) + { + if ( w3>M_PI ) w3 -= 2*M_PI; + if ( w1>M_PI ) w1 -= 2*M_PI; + } + } -template -Point CircleSeg :: GetPoint (double t) const -{ - if (t >= 1.0) { return p3; } + template + Point CircleSeg :: GetPoint (double t) const + { + if (t >= 1.0) { return p3; } - double phi = StartAngle() + t*(EndAngle()-StartAngle()); - Vec tmp(cos(phi),sin(phi)); + double phi = StartAngle() + t*(EndAngle()-StartAngle()); + Vec tmp(cos(phi),sin(phi)); - return pm + Radius()*tmp; -} + return pm + Radius()*tmp; + } -template -void CircleSeg :: GetCoeff (Vector & coeff) const -{ - coeff[0] = coeff[1] = 1.0; - coeff[2] = 0.0; - coeff[3] = -2.0 * pm[0]; - coeff[4] = -2.0 * pm[1]; - coeff[5] = sqr(pm[0]) + sqr(pm[1]) - sqr(Radius()); -} + template + void CircleSeg :: GetCoeff (Vector & coeff) const + { + coeff[0] = coeff[1] = 1.0; + coeff[2] = 0.0; + coeff[3] = -2.0 * pm[0]; + coeff[4] = -2.0 * pm[1]; + coeff[5] = sqr(pm[0]) + sqr(pm[1]) - sqr(Radius()); + } -template -void CircleSeg :: LineIntersections (const double a, const double b, const double c, - Array < Point > & points, const double eps) const -{ - points.SetSize(0); + template + void CircleSeg :: LineIntersections (const double a, const double b, const double c, + Array < Point > & points, const double eps) const + { + points.SetSize(0); - double px=0,py=0; + double px=0,py=0; - if(fabs(b) > 1e-20) - py = -c/b; - else - px = -c/a; + if(fabs(b) > 1e-20) + py = -c/b; + else + px = -c/a; - const double c1 = a*a + b*b; - const double c2 = 2. * ( a*(py-pm(1)) - b*(px-pm(0))); - const double c3 = pow(px-pm(0),2) + pow(py-pm(1),2) - pow(Radius(),2); + const double c1 = a*a + b*b; + const double c2 = 2. * ( a*(py-pm(1)) - b*(px-pm(0))); + const double c3 = pow(px-pm(0),2) + pow(py-pm(1),2) - pow(Radius(),2); - const double discr = c2*c2 - 4*c1*c3; + const double discr = c2*c2 - 4*c1*c3; - if(discr < 0) - return; + if(discr < 0) + return; - Array t; + Array t; - if(fabs(discr) < 1e-20) - t.Append(-0.5*c2/c1); - else - { - t.Append((-c2+sqrt(discr))/(2.*c1)); - t.Append((-c2-sqrt(discr))/(2.*c1)); - } + if(fabs(discr) < 1e-20) + t.Append(-0.5*c2/c1); + else + { + t.Append((-c2+sqrt(discr))/(2.*c1)); + t.Append((-c2-sqrt(discr))/(2.*c1)); + } - for(int i=0; i p (px-t[i]*b,py+t[i]*a); + for(int i=0; i p (px-t[i]*b,py+t[i]*a); - double angle = atan2(p(1),p(0))+M_PI; + double angle = atan2(p(1),p(0))+M_PI; - if(angle > StartAngle()-eps && angle < EndAngle()+eps) - points.Append(p); - } -} + if(angle > StartAngle()-eps && angle < EndAngle()+eps) + points.Append(p); + } + } -template -DiscretePointsSeg :: DiscretePointsSeg (const Array > & apts) - : pts (apts) -{ - for(int i=0; i + DiscretePointsSeg :: DiscretePointsSeg (const Array > & apts) + : pts (apts) + { + for(int i=0; i -DiscretePointsSeg :: ~DiscretePointsSeg () -{ ; } + template + DiscretePointsSeg :: ~DiscretePointsSeg () + { ; } -template -Point DiscretePointsSeg :: GetPoint (double t) const -{ - double t1 = t * (pts.Size()-1); - int segnr = int(t1); - if (segnr < 0) segnr = 0; - if (segnr >= pts.Size()) segnr = pts.Size()-1; + template + Point DiscretePointsSeg :: GetPoint (double t) const + { + double t1 = t * (pts.Size()-1); + int segnr = int(t1); + if (segnr < 0) segnr = 0; + if (segnr >= pts.Size()) segnr = pts.Size()-1; - double rest = t1 - segnr; + double rest = t1 - segnr; - return pts[segnr] + rest*Vec(pts[segnr+1]-pts[segnr]); -} + return pts[segnr] + rest*Vec(pts[segnr+1]-pts[segnr]); + } -typedef GeomPoint<2> GeomPoint2d; -typedef SplineSeg<2> SplineSegment; -typedef LineSeg<2> LineSegment; -typedef SplineSeg3<2> SplineSegment3; -typedef CircleSeg<2> CircleSegment; -typedef DiscretePointsSeg<2> DiscretePointsSegment; + typedef GeomPoint<2> GeomPoint2d; + typedef SplineSeg<2> SplineSegment; + typedef LineSeg<2> LineSegment; + typedef SplineSeg3<2> SplineSegment3; + typedef CircleSeg<2> CircleSegment; + typedef DiscretePointsSeg<2> DiscretePointsSegment; } diff --git a/libsrc/geom2d/splinegeometry.cpp b/libsrc/geom2d/splinegeometry.cpp index f0ea24ee..bb9c3104 100644 --- a/libsrc/geom2d/splinegeometry.cpp +++ b/libsrc/geom2d/splinegeometry.cpp @@ -6,1279 +6,1277 @@ #include -#include +#include #include -#include "meshing.hpp" namespace netgen { -template -void SplineGeometry :: LoadDataV2 ( ifstream & infile ) -{ - // new parser by Astrid Sinwel + template + void SplineGeometry :: LoadDataV2 ( ifstream & infile ) + { + // new parser by Astrid Sinwel - PrintMessage (1, "Load 2D Geometry V2"); - int nump, leftdom, rightdom; - Point x; - int hi1, hi2, hi3; - double hd; - char buf[50], ch; - int pointnr; + PrintMessage (1, "Load 2D Geometry V2"); + int nump, leftdom, rightdom; + Point x; + int hi1, hi2, hi3; + double hd; + char buf[50], ch; + int pointnr; - string keyword; + string keyword; - Array < GeomPoint > infilepoints (0); - Array pointnrs (0); - nump = 0; - int numdomains = 0; + Array < GeomPoint > infilepoints (0); + Array pointnrs (0); + nump = 0; + int numdomains = 0; - TestComment ( infile ); - // refinement factor - infile >> elto0; - TestComment ( infile ); + TestComment ( infile ); + // refinement factor + infile >> elto0; + TestComment ( infile ); - // test if next ch is a letter, i.e. new keyword starts - bool ischar = false; + // test if next ch is a letter, i.e. new keyword starts + bool ischar = false; - while ( infile.good() ) - { - infile >> keyword; + while ( infile.good() ) + { + infile >> keyword; - ischar = false; + ischar = false; - if ( keyword == "points" ) - { - PrintMessage (3, "load points"); - infile.get(ch); - infile.putback(ch); - - // test if ch is a letter - if ( int(ch) >= 65 && int(ch) <=90 ) - ischar = true; - if ( int(ch) >= 97 && int(ch) <= 122 ) - ischar = true; - - while ( ! ischar ) - { - TestComment ( infile ); - infile >> pointnr; - // pointnrs 1-based - if ( pointnr > nump ) nump = pointnr; - pointnrs.Append(pointnr); - - for(int j=0; j> x(j); - // hd is now optional, default 1 - // infile >> hd; - hd = 1; - - Flags flags; - - - // get flags, - ch = 'a'; - // infile >> ch; - do - { - infile.get (ch); - // if another int-value, set refinement flag to this value - // (corresponding to old files) - if ( int (ch) >= 48 && int(ch) <= 57 ) - { - infile.putback(ch); - infile >> hd; - infile.get(ch); - } - } - while (isspace(ch) && ch != '\n'); - while (ch == '-') - { - char flag[100]; - flag[0]='-'; - infile >> (flag+1); - flags.SetCommandLineFlag (flag); - ch = 'a'; - do { - infile.get (ch); - } while (isspace(ch) && ch != '\n'); - } - if (infile.good()) - infile.putback (ch); - - if ( hd == 1 ) - hd = flags.GetNumFlag ( "ref", 1.0); - // geompoints.Append (GeomPoint(x, hd)); - - infilepoints.Append ( GeomPoint(x, hd) ); - infilepoints.Last().hpref = flags.GetDefineFlag ("hpref"); - infilepoints.Last().hmax = flags.GetNumFlag ("maxh", 1e99); - - TestComment(infile); - infile.get(ch); - infile.putback(ch); - - // test if letter - if ( int(ch) >= 65 && int(ch) <=90 ) - ischar = true; - if ( int(ch) >= 97 && int(ch) <= 122 ) - ischar = true; - } - - // infile.putback (ch); - - geompoints.SetSize(nump); - for ( int i = 0; i < nump; i++ ) - { - geompoints[pointnrs[i] - 1] = infilepoints[i]; - geompoints[pointnrs[i] - 1].hpref = infilepoints[i].hpref; - } - TestComment(infile); - } - - else if ( keyword == "segments" ) - { - PrintMessage (3, "load segments"); - - bcnames.SetSize(0); - infile.get(ch); - infile.putback(ch); - int i = 0; - - // test if ch is a letter - if ( int(ch) >= 65 && int(ch) <=90 ) - ischar = true; - if ( int(ch) >= 97 && int(ch) <= 122 ) - ischar = true; - - while ( !ischar ) //ch != 'p' && ch != 'm' ) - { - i++; - TestComment ( infile ); - - SplineSeg * spline = 0; - TestComment ( infile ); - - infile >> leftdom >> rightdom; - - if ( leftdom > numdomains ) numdomains = leftdom; - if ( rightdom > numdomains ) numdomains = rightdom; - - - infile >> buf; - // type of spline segement - if (strcmp (buf, "2") == 0) - { // a line - infile >> hi1 >> hi2; - spline = new LineSeg(geompoints[hi1-1], - geompoints[hi2-1]); - } - else if (strcmp (buf, "3") == 0) - { // a rational spline - infile >> hi1 >> hi2 >> hi3; - spline = new SplineSeg3 (geompoints[hi1-1], - geompoints[hi2-1], - geompoints[hi3-1]); - } - else if (strcmp (buf, "4") == 0) - { // an arc - infile >> hi1 >> hi2 >> hi3; - spline = new CircleSeg (geompoints[hi1-1], - geompoints[hi2-1], - geompoints[hi3-1]); - break; - } - else if (strcmp (buf, "discretepoints") == 0) - { - int npts; - infile >> npts; - Array< Point > pts(npts); - for (int j = 0; j < npts; j++) - for(int k=0; k> pts[j](k); - - spline = new DiscretePointsSeg (pts); - } - - // infile >> spline->reffak; - spline -> leftdom = leftdom; - spline -> rightdom = rightdom; - splines.Append (spline); - - - // hd is now optional, default 1 - // infile >> hd; - hd = 1; - infile >> ch; - - // get refinement parameter, if it is there - //infile.get (ch); - // if another int-value, set refinement flag to this value - // (corresponding to old files) - - if ( int (ch) >= 48 && int(ch) <= 57 ) - { - infile.putback(ch); - infile >> hd; - infile >> ch ; - } - - // get flags, - Flags flags; - while (ch == '-') - { - char flag[100]; - flag[0]='-'; - infile >> (flag+1); - flags.SetCommandLineFlag (flag); - ch = 'a'; - infile >> ch; - } - - if (infile.good()) - infile.putback (ch); - - splines.Last()->bc = int (flags.GetNumFlag ("bc", i+1)); - splines.Last()->hpref_left = int (flags.GetDefineFlag ("hpref")) || - int (flags.GetDefineFlag ("hprefleft")); - splines.Last()->hpref_right = int (flags.GetDefineFlag ("hpref")) || - int (flags.GetDefineFlag ("hprefright")); - splines.Last()->copyfrom = int (flags.GetNumFlag ("copy", -1)); - splines.Last()->reffak = flags.GetNumFlag ("ref", 1 ); - splines.Last()->hmax = flags.GetNumFlag ("maxh", 1e99 ); - if ( hd != 1 ) splines.Last()->reffak = hd; - - if ( flags.StringFlagDefined("bcname") ) - { - int mybc = splines.Last()->bc-1; - for ( int ii = bcnames.Size(); ii <= mybc; ii++ ) - bcnames.Append ( new string ("default")); - if ( bcnames[mybc] ) delete bcnames[mybc]; - bcnames[mybc] = new string (flags.GetStringFlag("bcname","") ); - } - - TestComment(infile); - infile.get(ch); - infile.putback(ch); - - // test if ch is a letter - if ( int(ch) >= 65 && int(ch) <=90 ) - ischar = true; - if ( int(ch) >= 97 && int(ch) <= 122 ) - ischar = true; - - } - - infile.get(ch); - infile.putback(ch); - - - } - else if ( keyword == "materials" ) - { - TestComment ( infile ); - int domainnr; - char material[100]; - - if ( !infile.good() ) - return; - - materials.SetSize(numdomains) ; - maxh.SetSize ( numdomains ) ; - for ( int i = 0; i < numdomains; i++) - maxh[i] = 1000; - quadmeshing.SetSize ( numdomains ); - quadmeshing = false; - tensormeshing.SetSize ( numdomains ); - tensormeshing = false; - layer.SetSize ( numdomains ); - layer = 1; - - - TestComment ( infile ); - - for ( int i=0; i> domainnr; - infile >> material; - - strcpy (materials[domainnr-1], material); - - Flags flags; - ch = 'a'; - infile >> ch; - while (ch == '-') - { - char flag[100]; - flag[0]='-'; - infile >> (flag+1); - flags.SetCommandLineFlag (flag); - ch = 'a'; - infile >> ch; - } - - if (infile.good()) - infile.putback (ch); - - maxh[domainnr-1] = flags.GetNumFlag ( "maxh", 1000); - if (flags.GetDefineFlag("quad")) quadmeshing[domainnr-1] = true; - if (flags.GetDefineFlag("tensor")) tensormeshing[domainnr-1] = true; - layer[domainnr-1] = int(flags.GetNumFlag ("layer", 1)); - } - } - } - return; -} - - - - - - - -// check if comments in a .in2d file... -// template -// void SplineGeometry :: TestComment ( ifstream & infile ) -// { -// bool comment = true; -// char ch; -// infile.get(ch); -// infile.putback(ch); -// int ii = 0; -// while ( comment == true && ii < 100) -// { -// infile.get(ch); -// if ( ch == '#' ) -// while ( ch != '\n') -// { -// infile.get(ch); -// comment = false; -// } -// else if ( ch == '\n' ) -// { -// comment = true; -// ii ++; -// } -// else -// { -// infile.putback(ch); -// comment = false; -// } -// -// infile.get(ch) ; -// if ( ch == '\n' || ch == '#' ) -// { -// comment = true; -// } -// infile.putback(ch); -// if ( !comment ) break; -// } -// cerr << "** comment done" << endl; -// cerr << " * last char was " << ch << endl; -// return; -// -// } - -// herbert: fixed TestComment -template -void SplineGeometry :: TestComment ( ifstream & infile ) -{ - bool comment = true; - char ch; - while ( comment == true && !infile.eof() ) { - infile.get(ch); - if ( ch == '#' ) { // skip comments - while ( ch != '\n' && !infile.eof() ) { - infile.get(ch); - } - } - else if ( ch == '\n' ) { // skip empty lines - ; - } - else if ( isspace(ch) ) { // skip whitespaces - ; - } - else { // end of comment - infile.putback(ch); - comment = false; - } - } - return; -} - - - -template -SplineGeometry :: ~SplineGeometry() -{ - for(int i=0; i -int SplineGeometry :: Load (const Array & raw_data, const int startpos) -{ - int pos = startpos; - if(raw_data[pos] != D) - throw NgException("wrong dimension of spline raw_data"); - - pos++; - - elto0 = raw_data[pos]; pos++; - - splines.SetSize(int(raw_data[pos])); - pos++; - - Array< Point > pts(3); - - for(int i=0; i= 65 && int(ch) <=90 ) + ischar = true; + if ( int(ch) >= 97 && int(ch) <= 122 ) + ischar = true; + + while ( ! ischar ) + { + TestComment ( infile ); + infile >> pointnr; + // pointnrs 1-based + if ( pointnr > nump ) nump = pointnr; + pointnrs.Append(pointnr); + + for(int j=0; j> x(j); + // hd is now optional, default 1 + // infile >> hd; + hd = 1; + + Flags flags; + + + // get flags, + ch = 'a'; + // infile >> ch; + do + { + infile.get (ch); + // if another int-value, set refinement flag to this value + // (corresponding to old files) + if ( int (ch) >= 48 && int(ch) <= 57 ) + { + infile.putback(ch); + infile >> hd; + infile.get(ch); + } + } + while (isspace(ch) && ch != '\n'); + while (ch == '-') + { + char flag[100]; + flag[0]='-'; + infile >> (flag+1); + flags.SetCommandLineFlag (flag); + ch = 'a'; + do { + infile.get (ch); + } while (isspace(ch) && ch != '\n'); + } + if (infile.good()) + infile.putback (ch); + + if ( hd == 1 ) + hd = flags.GetNumFlag ( "ref", 1.0); + // geompoints.Append (GeomPoint(x, hd)); + + infilepoints.Append ( GeomPoint(x, hd) ); + infilepoints.Last().hpref = flags.GetDefineFlag ("hpref"); + infilepoints.Last().hmax = flags.GetNumFlag ("maxh", 1e99); + + TestComment(infile); + infile.get(ch); + infile.putback(ch); + + // test if letter + if ( int(ch) >= 65 && int(ch) <=90 ) + ischar = true; + if ( int(ch) >= 97 && int(ch) <= 122 ) + ischar = true; + } + + // infile.putback (ch); + + geompoints.SetSize(nump); + for ( int i = 0; i < nump; i++ ) + { + geompoints[pointnrs[i] - 1] = infilepoints[i]; + geompoints[pointnrs[i] - 1].hpref = infilepoints[i].hpref; + } + TestComment(infile); } - if (type == 2) - { - splines[i] = new LineSeg(GeomPoint(pts[0],1), - GeomPoint(pts[1],1)); - //(*testout) << "appending line segment " - // << pts[0] << " -- " << pts[1] << endl; - } - else if (type == 3) - { - splines[i] = new SplineSeg3(GeomPoint(pts[0],1), - GeomPoint(pts[1],1), - GeomPoint(pts[2],1)); - //(*testout) << "appending spline segment " - // << pts[0] << " -- " << pts[1] << " -- " << pts[2] << endl; - - } - else - throw NgException("something wrong with spline raw data"); - - } - return pos; -} - -template -void SplineGeometry :: GetRawData (Array & raw_data) const -{ - raw_data.Append(D); - raw_data.Append(elto0); - - - raw_data.Append(splines.Size()); - for(int i=0; iGetRawData(raw_data); - - -} - -template -void SplineGeometry :: CSGLoad (CSGScanner & scan) -{ - double hd; - Point x; - int nump, numseg; - - //scan.ReadNext(); - scan >> nump >> ';'; - - hd = 1; - geompoints.SetSize(nump); - for(int i = 0; i> x(0) >> ',' >> x(1) >> ';'; - else if(D==3) - scan >> x(0) >> ',' >> x(1) >> ',' >> x(2) >> ';'; - - geompoints[i] = GeomPoint(x,hd); - } - - scan >> numseg;// >> ';'; - - splines.SetSize(numseg); - - int pnums,pnum1,pnum2,pnum3; - - - for(int i = 0; i> ';' >> pnums >> ','; - if (pnums == 2) - { - scan >> pnum1 >> ',' >> pnum2;// >> ';'; - splines[i] = new LineSeg(geompoints[pnum1-1], - geompoints[pnum2-1]); - } - else if (pnums == 3) - { - scan >> pnum1 >> ',' >> pnum2 >> ',' - >> pnum3;// >> ';'; - splines[i] = new SplineSeg3(geompoints[pnum1-1], - geompoints[pnum2-1], - geompoints[pnum3-1]); - } - else if (pnums == 4) - { - scan >> pnum1 >> ',' >> pnum2 >> ',' - >> pnum3;// >> ';'; - splines[i] = new CircleSeg(geompoints[pnum1-1], - geompoints[pnum2-1], - geompoints[pnum3-1]); - - } - - } -} - - - - -template -void SplineGeometry :: Load (const char * filename) -{ - - ifstream infile; - Point x; - char buf[50]; - - - infile.open (filename); - - if ( ! infile.good() ) - throw NgException(string ("Input file '") + - string (filename) + - string ("' not available!")); - - TestComment ( infile ); - - infile >> buf; // file recognition - - tensormeshing.SetSize(0); - quadmeshing.SetSize(0); - - TestComment ( infile ); - if ( strcmp (buf, "splinecurves2dnew") == 0 ) - { - LoadDataNew ( infile ); - } - else if ( strcmp (buf, "splinecurves2dv2") == 0 ) - { - LoadDataV2 ( infile ); - } - else - { - LoadData(infile ); - } - infile.close(); -} - - -template -void SplineGeometry :: LoadDataNew ( ifstream & infile ) -{ - - int nump, numseg, leftdom, rightdom; - Point x; - int hi1, hi2, hi3; - double hd; - char buf[50], ch; - int pointnr; - - - TestComment ( infile ); - infile >> elto0; - TestComment ( infile ); - - infile >> nump; - geompoints.SetSize(nump); - - for (int i = 0; i < nump; i++) - { - TestComment ( infile ); - infile >> pointnr; - if ( pointnr > nump ) - { - throw NgException(string ("Point number greater than total number of points") ); - } - for(int j=0; j> x(j); - - - // hd is now optional, default 1 - // infile >> hd; - hd = 1; - - Flags flags; - - - // get flags, - ch = 'a'; - // infile >> ch; - do - { - - infile.get (ch); - // if another int-value, set refinement flag to this value - // (corresponding to old files) - if ( int (ch) >= 48 && int(ch) <= 57 ) - { - infile.putback(ch); - infile >> hd; - infile.get(ch); - } - } - while (isspace(ch) && ch != '\n'); - while (ch == '-') - { - char flag[100]; - flag[0]='-'; - infile >> (flag+1); - flags.SetCommandLineFlag (flag); - ch = 'a'; - do { - infile.get (ch); - } while (isspace(ch) && ch != '\n'); - } - - if (infile.good()) - infile.putback (ch); - - if ( hd == 1 ) - hd = flags.GetNumFlag ( "ref", 1.0); - // geompoints.Append (GeomPoint(x, hd)); - geompoints[pointnr-1] = GeomPoint(x, hd); - geompoints[pointnr-1].hpref = flags.GetDefineFlag ("hpref"); - } - - TestComment ( infile ); - - infile >> numseg; - bcnames.SetSize(numseg); - for ( int i = 0; i < numseg; i++ ) - bcnames[i] = 0;//new"default"; - - SplineSeg * spline = 0; - for (int i = 0; i < numseg; i++) - { - TestComment ( infile ); - - infile >> leftdom >> rightdom; - - // cout << "add spline " << i << ", left = " << leftdom << endl; - - infile >> buf; - // type of spline segement - if (strcmp (buf, "2") == 0) - { // a line - infile >> hi1 >> hi2; - spline = new LineSeg (geompoints[hi1-1], - geompoints[hi2-1]); - } - else if (strcmp (buf, "3") == 0) - { // a rational spline - infile >> hi1 >> hi2 >> hi3; - spline = new SplineSeg3 (geompoints[hi1-1], - geompoints[hi2-1], - geompoints[hi3-1]); - } - else if (strcmp (buf, "4") == 0) - { // an arc - infile >> hi1 >> hi2 >> hi3; - spline = new CircleSeg (geompoints[hi1-1], - geompoints[hi2-1], - geompoints[hi3-1]); -// break; - } - else if (strcmp (buf, "discretepoints") == 0) - { - int npts; - infile >> npts; - Array< Point > pts(npts); - for (int j = 0; j < npts; j++) - for(int k=0; k> pts[j](k); - - spline = new DiscretePointsSeg (pts); - } - - // infile >> spline->reffak; - spline -> leftdom = leftdom; - spline -> rightdom = rightdom; - splines.Append (spline); - - // hd is now optional, default 1 - // infile >> hd; - hd = 1; - infile >> ch; - - // get refinement parameter, if it is there - // infile.get (ch); - // if another int-value, set refinement flag to this value - // (corresponding to old files) - if ( int (ch) >= 48 && int(ch) <= 57 ) - { - infile.putback(ch); - infile >> hd; - infile >> ch ; - } - - Flags flags; - while (ch == '-') - { - char flag[100]; - flag[0]='-'; - infile >> (flag+1); - flags.SetCommandLineFlag (flag); - ch = 'a'; - infile >> ch; - } - - if (infile.good()) - infile.putback (ch); - - splines.Last()->bc = int (flags.GetNumFlag ("bc", i+1)); - splines.Last()->hpref_left = int (flags.GetDefineFlag ("hpref")) || - int (flags.GetDefineFlag ("hprefleft")); - splines.Last()->hpref_right = int (flags.GetDefineFlag ("hpref")) || - int (flags.GetDefineFlag ("hprefright")); - splines.Last()->copyfrom = int (flags.GetNumFlag ("copy", -1)); - splines.Last()->reffak = flags.GetNumFlag ("ref", 1 ); - splines.Last()->hmax = flags.GetNumFlag ("maxh", 1e99 ); - - if ( flags.StringFlagDefined("bcname") ) - { - int mybc = splines.Last()->bc-1; - if ( bcnames[mybc] ) delete bcnames[mybc]; - bcnames[mybc] = new string (flags.GetStringFlag("bcname","") ); - } - - if ( hd != 1 ) - splines.Last()->reffak = hd; - } - if ( !infile.good() ) - return; - TestComment ( infile ); - int numdomains; - int domainnr; - char material[100]; - - if ( !infile.good() ) - return; - - infile >> numdomains; - materials.SetSize(numdomains) ; - maxh.SetSize ( numdomains ) ; - for ( int i = 0; i < numdomains; i++) - maxh[i] = 1000; - - TestComment ( infile ); - - for ( int i=0; i> domainnr; - infile >> material; - strcpy(materials[domainnr-1], material); - - Flags flags; - ch = 'a'; - infile >> ch; - while (ch == '-') - { - char flag[100]; - flag[0]='-'; - infile >> (flag+1); - flags.SetCommandLineFlag (flag); - ch = 'a'; - infile >> ch; - } - - if (infile.good()) - infile.putback (ch); - - maxh[domainnr-1] = flags.GetNumFlag ( "maxh", 1000); - } - return; -} - - - -template -void SplineGeometry :: LoadData ( ifstream & infile ) -{ - - int nump, numseg, leftdom, rightdom; - Point x; - int hi1, hi2, hi3; - double hd; - char buf[50], ch; - - materials.SetSize(0); - maxh.SetSize(0); - infile >> elto0; - - TestComment ( infile ); - - infile >> nump; - for (int i = 0; i < nump; i++) - { - TestComment ( infile ); - for(int j=0; j> x(j); - infile >> hd; - - Flags flags; - - ch = 'a'; - // infile >> ch; - do { - infile.get (ch); - } while (isspace(ch) && ch != '\n'); - while (ch == '-') - { - char flag[100]; - flag[0]='-'; - infile >> (flag+1); - flags.SetCommandLineFlag (flag); - ch = 'a'; - do { - infile.get (ch); - } while (isspace(ch) && ch != '\n'); - } - - if (infile.good()) - infile.putback (ch); - - geompoints.Append (GeomPoint(x, hd)); - geompoints.Last().hpref = flags.GetDefineFlag ("hpref"); - geompoints.Last().hmax = 1e99; - } - - PrintMessage (3, nump, " points loaded"); - TestComment ( infile ); - - infile >> numseg; - bcnames.SetSize(numseg); - for ( int i = 0; i < numseg; i++ ) - bcnames[i] = 0; // "default"; - - SplineSeg * spline = 0; - - PrintMessage (3, numseg, " segments loaded"); - for (int i = 0; i < numseg; i++) - { - TestComment ( infile ); - - infile >> leftdom >> rightdom; - - // cout << "add spline " << i << ", left = " << leftdom << ", right = " << rightdom << endl; - - infile >> buf; - // type of spline segement - if (strcmp (buf, "2") == 0) - { // a line - infile >> hi1 >> hi2; - spline = new LineSeg(geompoints[hi1-1], - geompoints[hi2-1]); - } - else if (strcmp (buf, "3") == 0) - { // a rational spline - infile >> hi1 >> hi2 >> hi3; - spline = new SplineSeg3 (geompoints[hi1-1], - geompoints[hi2-1], - geompoints[hi3-1]); - } - else if (strcmp (buf, "4") == 0) - { // an arc - infile >> hi1 >> hi2 >> hi3; - spline = new CircleSeg (geompoints[hi1-1], - geompoints[hi2-1], - geompoints[hi3-1]); -// break; - } - else if (strcmp (buf, "discretepoints") == 0) - { - int npts; - infile >> npts; - Array< Point > pts(npts); - for (int j = 0; j < npts; j++) - for(int k=0; k> pts[j](k); - - spline = new DiscretePointsSeg (pts); - } - - infile >> spline->reffak; - spline -> leftdom = leftdom; - spline -> rightdom = rightdom; - spline -> hmax = 1e99; - splines.Append (spline); - - - Flags flags; - ch = 'a'; - infile >> ch; - while (ch == '-') - { - char flag[100]; - flag[0]='-'; - infile >> (flag+1); - flags.SetCommandLineFlag (flag); - ch = 'a'; - infile >> ch; - } - - if (infile.good()) - infile.putback (ch); - - splines.Last()->bc = int (flags.GetNumFlag ("bc", i+1)); - splines.Last()->hpref_left = int (flags.GetDefineFlag ("hpref")) || - int (flags.GetDefineFlag ("hprefleft")); - splines.Last()->hpref_right = int (flags.GetDefineFlag ("hpref")) || - int (flags.GetDefineFlag ("hprefright")); - splines.Last()->copyfrom = int (flags.GetNumFlag ("copy", -1)); - if ( flags.StringFlagDefined("bcname") ) - { - int mybc = splines.Last()->bc-1; - if ( bcnames[mybc] ) delete bcnames[mybc]; - bcnames[mybc] = new string (flags.GetStringFlag("bcname","") ); - } - } -} - - - -void SplineGeometry2d :: PartitionBoundary (double h, Mesh & mesh2d) -{ - enum { D = 2 }; - Box bbox; - GetBoundingBox (bbox); - double dist = Dist (bbox.PMin(), bbox.PMax()); - Point<3> pmin; - Point<3> pmax; - - pmin(2) = -dist; pmax(2) = dist; - for(int j=0;jleftdom : splines[i]->rightdom; - if (dom != 0) splines[i] -> layer = GetDomainLayer (dom); + else if ( keyword == "segments" ) + { + PrintMessage (3, "load segments"); + + bcnames.SetSize(0); + infile.get(ch); + infile.putback(ch); + int i = 0; + + // test if ch is a letter + if ( int(ch) >= 65 && int(ch) <=90 ) + ischar = true; + if ( int(ch) >= 97 && int(ch) <= 122 ) + ischar = true; + + while ( !ischar ) //ch != 'p' && ch != 'm' ) + { + i++; + TestComment ( infile ); + + SplineSeg * spline = 0; + TestComment ( infile ); + + infile >> leftdom >> rightdom; + + if ( leftdom > numdomains ) numdomains = leftdom; + if ( rightdom > numdomains ) numdomains = rightdom; + + + infile >> buf; + // type of spline segement + if (strcmp (buf, "2") == 0) + { // a line + infile >> hi1 >> hi2; + spline = new LineSeg(geompoints[hi1-1], + geompoints[hi2-1]); + } + else if (strcmp (buf, "3") == 0) + { // a rational spline + infile >> hi1 >> hi2 >> hi3; + spline = new SplineSeg3 (geompoints[hi1-1], + geompoints[hi2-1], + geompoints[hi3-1]); + } + else if (strcmp (buf, "4") == 0) + { // an arc + infile >> hi1 >> hi2 >> hi3; + spline = new CircleSeg (geompoints[hi1-1], + geompoints[hi2-1], + geompoints[hi3-1]); + break; + } + else if (strcmp (buf, "discretepoints") == 0) + { + int npts; + infile >> npts; + Array< Point > pts(npts); + for (int j = 0; j < npts; j++) + for(int k=0; k> pts[j](k); + + spline = new DiscretePointsSeg (pts); + } + + // infile >> spline->reffak; + spline -> leftdom = leftdom; + spline -> rightdom = rightdom; + splines.Append (spline); + + + // hd is now optional, default 1 + // infile >> hd; + hd = 1; + infile >> ch; + + // get refinement parameter, if it is there + //infile.get (ch); + // if another int-value, set refinement flag to this value + // (corresponding to old files) + + if ( int (ch) >= 48 && int(ch) <= 57 ) + { + infile.putback(ch); + infile >> hd; + infile >> ch ; + } + + // get flags, + Flags flags; + while (ch == '-') + { + char flag[100]; + flag[0]='-'; + infile >> (flag+1); + flags.SetCommandLineFlag (flag); + ch = 'a'; + infile >> ch; + } + + if (infile.good()) + infile.putback (ch); + + splines.Last()->bc = int (flags.GetNumFlag ("bc", i+1)); + splines.Last()->hpref_left = int (flags.GetDefineFlag ("hpref")) || + int (flags.GetDefineFlag ("hprefleft")); + splines.Last()->hpref_right = int (flags.GetDefineFlag ("hpref")) || + int (flags.GetDefineFlag ("hprefright")); + splines.Last()->copyfrom = int (flags.GetNumFlag ("copy", -1)); + splines.Last()->reffak = flags.GetNumFlag ("ref", 1 ); + splines.Last()->hmax = flags.GetNumFlag ("maxh", 1e99 ); + if ( hd != 1 ) splines.Last()->reffak = hd; + + if ( flags.StringFlagDefined("bcname") ) + { + int mybc = splines.Last()->bc-1; + for ( int ii = bcnames.Size(); ii <= mybc; ii++ ) + bcnames.Append ( new string ("default")); + if ( bcnames[mybc] ) delete bcnames[mybc]; + bcnames[mybc] = new string (flags.GetStringFlag("bcname","") ); + } + + TestComment(infile); + infile.get(ch); + infile.putback(ch); + + // test if ch is a letter + if ( int(ch) >= 65 && int(ch) <=90 ) + ischar = true; + if ( int(ch) >= 97 && int(ch) <= 122 ) + ischar = true; + + } + + infile.get(ch); + infile.putback(ch); + + + } + else if ( keyword == "materials" ) + { + TestComment ( infile ); + int domainnr; + char material[100]; + + if ( !infile.good() ) + return; + + materials.SetSize(numdomains) ; + maxh.SetSize ( numdomains ) ; + for ( int i = 0; i < numdomains; i++) + maxh[i] = 1000; + quadmeshing.SetSize ( numdomains ); + quadmeshing = false; + tensormeshing.SetSize ( numdomains ); + tensormeshing = false; + layer.SetSize ( numdomains ); + layer = 1; + + + TestComment ( infile ); + + for ( int i=0; i> domainnr; + infile >> material; + + strcpy (materials[domainnr-1], material); + + Flags flags; + ch = 'a'; + infile >> ch; + while (ch == '-') + { + char flag[100]; + flag[0]='-'; + infile >> (flag+1); + flags.SetCommandLineFlag (flag); + ch = 'a'; + infile >> ch; + } + + if (infile.good()) + infile.putback (ch); + + maxh[domainnr-1] = flags.GetNumFlag ( "maxh", 1000); + if (flags.GetDefineFlag("quad")) quadmeshing[domainnr-1] = true; + if (flags.GetDefineFlag("tensor")) tensormeshing[domainnr-1] = true; + layer[domainnr-1] = int(flags.GetNumFlag ("layer", 1)); + } + } } + return; + } - for (int i = 0; i < splines.Size(); i++) - if (splines[i]->copyfrom == -1) + + + + + + + // check if comments in a .in2d file... + // template + // void SplineGeometry :: TestComment ( ifstream & infile ) + // { + // bool comment = true; + // char ch; + // infile.get(ch); + // infile.putback(ch); + // int ii = 0; + // while ( comment == true && ii < 100) + // { + // infile.get(ch); + // if ( ch == '#' ) + // while ( ch != '\n') + // { + // infile.get(ch); + // comment = false; + // } + // else if ( ch == '\n' ) + // { + // comment = true; + // ii ++; + // } + // else + // { + // infile.putback(ch); + // comment = false; + // } + // + // infile.get(ch) ; + // if ( ch == '\n' || ch == '#' ) + // { + // comment = true; + // } + // infile.putback(ch); + // if ( !comment ) break; + // } + // cerr << "** comment done" << endl; + // cerr << " * last char was " << ch << endl; + // return; + // + // } + + // herbert: fixed TestComment + template + void SplineGeometry :: TestComment ( ifstream & infile ) + { + bool comment = true; + char ch; + while ( comment == true && !infile.eof() ) { + infile.get(ch); + if ( ch == '#' ) { // skip comments + while ( ch != '\n' && !infile.eof() ) { + infile.get(ch); + } + } + else if ( ch == '\n' ) { // skip empty lines + ; + } + else if ( isspace(ch) ) { // skip whitespaces + ; + } + else { // end of comment + infile.putback(ch); + comment = false; + } + } + return; + } + + + + template + SplineGeometry :: ~SplineGeometry() + { + for(int i=0; i + int SplineGeometry :: Load (const Array & raw_data, const int startpos) + { + int pos = startpos; + if(raw_data[pos] != D) + throw NgException("wrong dimension of spline raw_data"); + + pos++; + + elto0 = raw_data[pos]; pos++; + + splines.SetSize(int(raw_data[pos])); + pos++; + + Array< Point > pts(3); + + for(int i=0; ileftdom ), GetDomainMaxh ( splines[i]->rightdom ) ); - double maximum = max2 ( GetDomainMaxh ( splines[i]->leftdom ), GetDomainMaxh ( splines[i]->rightdom ) ); - minimum = min2 ( minimum, h ); - maximum = min2 ( maximum, h); - if ( minimum > 0 ) - splines[i]->Partition(minimum, elto0, mesh2d, searchtree, i+1); - else if ( maximum > 0 ) - splines[i]->Partition(maximum, elto0, mesh2d, searchtree, i+1); - else - splines[i]->Partition(h, elto0, mesh2d, searchtree, i+1); + int type = int(raw_data[pos]); + pos++; + + for(int j=0; j(GeomPoint(pts[0],1), + GeomPoint(pts[1],1)); + //(*testout) << "appending line segment " + // << pts[0] << " -- " << pts[1] << endl; + } + else if (type == 3) + { + splines[i] = new SplineSeg3(GeomPoint(pts[0],1), + GeomPoint(pts[1],1), + GeomPoint(pts[2],1)); + //(*testout) << "appending spline segment " + // << pts[0] << " -- " << pts[1] << " -- " << pts[2] << endl; + + } + else + throw NgException("something wrong with spline raw data"); + + } + return pos; + } + + template + void SplineGeometry :: GetRawData (Array & raw_data) const + { + raw_data.Append(D); + raw_data.Append(elto0); + + + raw_data.Append(splines.Size()); + for(int i=0; iGetRawData(raw_data); + + + } + + /* + template + void SplineGeometry :: CSGLoad (CSGScanner & scan) + { + double hd; + Point x; + int nump, numseg; + + //scan.ReadNext(); + scan >> nump >> ';'; + + hd = 1; + geompoints.SetSize(nump); + for(int i = 0; i> x(0) >> ',' >> x(1) >> ';'; + else if(D==3) + scan >> x(0) >> ',' >> x(1) >> ',' >> x(2) >> ';'; + + geompoints[i] = GeomPoint(x,hd); + } + + scan >> numseg;// >> ';'; + + splines.SetSize(numseg); + + int pnums,pnum1,pnum2,pnum3; + + + for(int i = 0; i> ';' >> pnums >> ','; + if (pnums == 2) + { + scan >> pnum1 >> ',' >> pnum2;// >> ';'; + splines[i] = new LineSeg(geompoints[pnum1-1], + geompoints[pnum2-1]); + } + else if (pnums == 3) + { + scan >> pnum1 >> ',' >> pnum2 >> ',' + >> pnum3;// >> ';'; + splines[i] = new SplineSeg3(geompoints[pnum1-1], + geompoints[pnum2-1], + geompoints[pnum3-1]); + } + else if (pnums == 4) + { + scan >> pnum1 >> ',' >> pnum2 >> ',' + >> pnum3;// >> ';'; + splines[i] = new CircleSeg(geompoints[pnum1-1], + geompoints[pnum2-1], + geompoints[pnum3-1]); + + } + + } + } + */ + + + + template + void SplineGeometry :: Load (const char * filename) + { + + ifstream infile; + Point x; + char buf[50]; + + + infile.open (filename); + + if ( ! infile.good() ) + throw NgException(string ("Input file '") + + string (filename) + + string ("' not available!")); + + TestComment ( infile ); + + infile >> buf; // file recognition + + tensormeshing.SetSize(0); + quadmeshing.SetSize(0); + + TestComment ( infile ); + if ( strcmp (buf, "splinecurves2dnew") == 0 ) + { + LoadDataNew ( infile ); + } + else if ( strcmp (buf, "splinecurves2dv2") == 0 ) + { + LoadDataV2 ( infile ); } else { - CopyEdgeMesh (splines[i]->copyfrom, i+1, mesh2d, searchtree); + LoadData(infile ); } -} + infile.close(); + } -template -void SplineGeometry :: CopyEdgeMesh (int from, int to, Mesh & mesh, Point3dTree & searchtree) -{ - int i; + template + void SplineGeometry :: LoadDataNew ( ifstream & infile ) + { - Array mappoints (mesh.GetNP()); - Array param (mesh.GetNP()); - mappoints = -1; - param = 0; + int nump, numseg, leftdom, rightdom; + Point x; + int hi1, hi2, hi3; + double hd; + char buf[50], ch; + int pointnr; - Point3d pmin, pmax; - mesh.GetBox (pmin, pmax); - double diam2 = Dist2(pmin, pmax); - if (printmessage_importance>0) - cout << "copy edge, from = " << from << " to " << to << endl; - - for (i = 1; i <= mesh.GetNSeg(); i++) - { - const Segment & seg = mesh.LineSegment(i); - if (seg.edgenr == from) - { - mappoints.Elem(seg[0]) = 1; - param.Elem(seg[0]) = seg.epgeominfo[0].dist; + TestComment ( infile ); + infile >> elto0; + TestComment ( infile ); + + infile >> nump; + geompoints.SetSize(nump); + + for (int i = 0; i < nump; i++) + { + TestComment ( infile ); + infile >> pointnr; + if ( pointnr > nump ) + { + throw NgException(string ("Point number greater than total number of points") ); + } + for(int j=0; j> x(j); - mappoints.Elem(seg[1]) = 1; - param.Elem(seg[1]) = seg.epgeominfo[1].dist; - } - } - bool mapped = false; - for (i = 1; i <= mappoints.Size(); i++) - { - if (mappoints.Get(i) != -1) - { - Point newp = splines.Get(to)->GetPoint (param.Get(i)); - Point<3> newp3; - for(int j=0; j> hd; + hd = 1; - mappoints.Elem(i) = npi; + Flags flags; - mesh.GetIdentifications().Add (i, npi, to); - mapped = true; - } - } - if(mapped) - mesh.GetIdentifications().SetType(to,Identifications::PERIODIC); - // copy segments - int oldnseg = mesh.GetNSeg(); - for (i = 1; i <= oldnseg; i++) - { - const Segment & seg = mesh.LineSegment(i); - if (seg.edgenr == from) - { - Segment nseg; - nseg.edgenr = to; - nseg.si = splines.Get(to)->bc; - nseg[0] = mappoints.Get(seg[0]); - nseg[1] = mappoints.Get(seg[1]); - nseg.domin = splines.Get(to)->leftdom; - nseg.domout = splines.Get(to)->rightdom; - - nseg.epgeominfo[0].edgenr = to; - nseg.epgeominfo[0].dist = param.Get(seg[0]); - nseg.epgeominfo[1].edgenr = to; - nseg.epgeominfo[1].dist = param.Get(seg[1]); - mesh.AddSegment (nseg); - } - } -} - + // get flags, + ch = 'a'; + // infile >> ch; + do + { -template -void SplineGeometry :: GetBoundingBox (Box & box) const -{ - if (!splines.Size()) - { - Point auxp = 0.; - box.Set (auxp); + infile.get (ch); + // if another int-value, set refinement flag to this value + // (corresponding to old files) + if ( int (ch) >= 48 && int(ch) <= 57 ) + { + infile.putback(ch); + infile >> hd; + infile.get(ch); + } + } + while (isspace(ch) && ch != '\n'); + while (ch == '-') + { + char flag[100]; + flag[0]='-'; + infile >> (flag+1); + flags.SetCommandLineFlag (flag); + ch = 'a'; + do { + infile.get (ch); + } while (isspace(ch) && ch != '\n'); + } + + if (infile.good()) + infile.putback (ch); + + if ( hd == 1 ) + hd = flags.GetNumFlag ( "ref", 1.0); + // geompoints.Append (GeomPoint(x, hd)); + geompoints[pointnr-1] = GeomPoint(x, hd); + geompoints[pointnr-1].hpref = flags.GetDefineFlag ("hpref"); + } + + TestComment ( infile ); + + infile >> numseg; + bcnames.SetSize(numseg); + for ( int i = 0; i < numseg; i++ ) + bcnames[i] = 0;//new"default"; + + SplineSeg * spline = 0; + for (int i = 0; i < numseg; i++) + { + TestComment ( infile ); + + infile >> leftdom >> rightdom; + + // cout << "add spline " << i << ", left = " << leftdom << endl; + + infile >> buf; + // type of spline segement + if (strcmp (buf, "2") == 0) + { // a line + infile >> hi1 >> hi2; + spline = new LineSeg (geompoints[hi1-1], + geompoints[hi2-1]); + } + else if (strcmp (buf, "3") == 0) + { // a rational spline + infile >> hi1 >> hi2 >> hi3; + spline = new SplineSeg3 (geompoints[hi1-1], + geompoints[hi2-1], + geompoints[hi3-1]); + } + else if (strcmp (buf, "4") == 0) + { // an arc + infile >> hi1 >> hi2 >> hi3; + spline = new CircleSeg (geompoints[hi1-1], + geompoints[hi2-1], + geompoints[hi3-1]); + // break; + } + else if (strcmp (buf, "discretepoints") == 0) + { + int npts; + infile >> npts; + Array< Point > pts(npts); + for (int j = 0; j < npts; j++) + for(int k=0; k> pts[j](k); + + spline = new DiscretePointsSeg (pts); + } + + // infile >> spline->reffak; + spline -> leftdom = leftdom; + spline -> rightdom = rightdom; + splines.Append (spline); + + // hd is now optional, default 1 + // infile >> hd; + hd = 1; + infile >> ch; + + // get refinement parameter, if it is there + // infile.get (ch); + // if another int-value, set refinement flag to this value + // (corresponding to old files) + if ( int (ch) >= 48 && int(ch) <= 57 ) + { + infile.putback(ch); + infile >> hd; + infile >> ch ; + } + + Flags flags; + while (ch == '-') + { + char flag[100]; + flag[0]='-'; + infile >> (flag+1); + flags.SetCommandLineFlag (flag); + ch = 'a'; + infile >> ch; + } + + if (infile.good()) + infile.putback (ch); + + splines.Last()->bc = int (flags.GetNumFlag ("bc", i+1)); + splines.Last()->hpref_left = int (flags.GetDefineFlag ("hpref")) || + int (flags.GetDefineFlag ("hprefleft")); + splines.Last()->hpref_right = int (flags.GetDefineFlag ("hpref")) || + int (flags.GetDefineFlag ("hprefright")); + splines.Last()->copyfrom = int (flags.GetNumFlag ("copy", -1)); + splines.Last()->reffak = flags.GetNumFlag ("ref", 1 ); + splines.Last()->hmax = flags.GetNumFlag ("maxh", 1e99 ); + + if ( flags.StringFlagDefined("bcname") ) + { + int mybc = splines.Last()->bc-1; + if ( bcnames[mybc] ) delete bcnames[mybc]; + bcnames[mybc] = new string (flags.GetStringFlag("bcname","") ); + } + + if ( hd != 1 ) + splines.Last()->reffak = hd; + } + if ( !infile.good() ) return; - } + TestComment ( infile ); + int numdomains; + int domainnr; + char material[100]; - Array > points; - for (int i = 0; i < splines.Size(); i++) - { - splines[i]->GetPoints (20, points); + if ( !infile.good() ) + return; - if (i == 0) box.Set(points[0]); - for (int j = 0; j < points.Size(); j++) - box.Add (points[j]); - } -} + infile >> numdomains; + materials.SetSize(numdomains) ; + maxh.SetSize ( numdomains ) ; + for ( int i = 0; i < numdomains; i++) + maxh[i] = 1000; -template -void SplineGeometry :: SetGrading (const double grading) -{ elto0 = grading;} + TestComment ( infile ); -/* -template -void SplineGeometry :: AppendPoint (const double x, const double y, const double reffac, const bool hpref) -{ - geompoints.Append (GeomPoint(x, y, reffac)); - geompoints.Last().hpref = hpref; -} -*/ + for ( int i=0; i -void SplineGeometry :: AppendPoint (const Point & p, const double reffac, const bool hpref) -{ - geompoints.Append (GeomPoint(p, reffac)); - geompoints.Last().hpref = hpref; -} + for ( int i=0; i> domainnr; + infile >> material; + strcpy(materials[domainnr-1], material); + + Flags flags; + ch = 'a'; + infile >> ch; + while (ch == '-') + { + char flag[100]; + flag[0]='-'; + infile >> (flag+1); + flags.SetCommandLineFlag (flag); + ch = 'a'; + infile >> ch; + } + + if (infile.good()) + infile.putback (ch); + + maxh[domainnr-1] = flags.GetNumFlag ( "maxh", 1000); + } + return; + } + template + void SplineGeometry :: LoadData ( ifstream & infile ) + { -template -void SplineGeometry :: AppendSegment(SplineSeg * spline, const int leftdomain, const int rightdomain, - const int bc, - const double reffac, const bool hprefleft, const bool hprefright, - const int copyfrom) -{ - spline -> leftdom = leftdomain; - spline -> rightdom = rightdomain; - spline -> bc = (bc >= 0) ? bc : (splines.Size()+1); - spline -> reffak = reffac; - spline -> hpref_left = hprefleft; - spline -> hpref_right = hprefright; - spline -> copyfrom = copyfrom; + int nump, numseg, leftdom, rightdom; + Point x; + int hi1, hi2, hi3; + double hd; + char buf[50], ch; + + materials.SetSize(0); + maxh.SetSize(0); + infile >> elto0; + + TestComment ( infile ); + + infile >> nump; + for (int i = 0; i < nump; i++) + { + TestComment ( infile ); + for(int j=0; j> x(j); + infile >> hd; + + Flags flags; + + ch = 'a'; + // infile >> ch; + do { + infile.get (ch); + } while (isspace(ch) && ch != '\n'); + while (ch == '-') + { + char flag[100]; + flag[0]='-'; + infile >> (flag+1); + flags.SetCommandLineFlag (flag); + ch = 'a'; + do { + infile.get (ch); + } while (isspace(ch) && ch != '\n'); + } + + if (infile.good()) + infile.putback (ch); + + geompoints.Append (GeomPoint(x, hd)); + geompoints.Last().hpref = flags.GetDefineFlag ("hpref"); + geompoints.Last().hmax = 1e99; + } + + PrintMessage (3, nump, " points loaded"); + TestComment ( infile ); + + infile >> numseg; + bcnames.SetSize(numseg); + for ( int i = 0; i < numseg; i++ ) + bcnames[i] = 0; // "default"; + + SplineSeg * spline = 0; + + PrintMessage (3, numseg, " segments loaded"); + for (int i = 0; i < numseg; i++) + { + TestComment ( infile ); + + infile >> leftdom >> rightdom; + + // cout << "add spline " << i << ", left = " << leftdom << ", right = " << rightdom << endl; + + infile >> buf; + // type of spline segement + if (strcmp (buf, "2") == 0) + { // a line + infile >> hi1 >> hi2; + spline = new LineSeg(geompoints[hi1-1], + geompoints[hi2-1]); + } + else if (strcmp (buf, "3") == 0) + { // a rational spline + infile >> hi1 >> hi2 >> hi3; + spline = new SplineSeg3 (geompoints[hi1-1], + geompoints[hi2-1], + geompoints[hi3-1]); + } + else if (strcmp (buf, "4") == 0) + { // an arc + infile >> hi1 >> hi2 >> hi3; + spline = new CircleSeg (geompoints[hi1-1], + geompoints[hi2-1], + geompoints[hi3-1]); + // break; + } + else if (strcmp (buf, "discretepoints") == 0) + { + int npts; + infile >> npts; + Array< Point > pts(npts); + for (int j = 0; j < npts; j++) + for(int k=0; k> pts[j](k); + + spline = new DiscretePointsSeg (pts); + } + + infile >> spline->reffak; + spline -> leftdom = leftdom; + spline -> rightdom = rightdom; + spline -> hmax = 1e99; + splines.Append (spline); + + + Flags flags; + ch = 'a'; + infile >> ch; + while (ch == '-') + { + char flag[100]; + flag[0]='-'; + infile >> (flag+1); + flags.SetCommandLineFlag (flag); + ch = 'a'; + infile >> ch; + } + + if (infile.good()) + infile.putback (ch); + + splines.Last()->bc = int (flags.GetNumFlag ("bc", i+1)); + splines.Last()->hpref_left = int (flags.GetDefineFlag ("hpref")) || + int (flags.GetDefineFlag ("hprefleft")); + splines.Last()->hpref_right = int (flags.GetDefineFlag ("hpref")) || + int (flags.GetDefineFlag ("hprefright")); + splines.Last()->copyfrom = int (flags.GetNumFlag ("copy", -1)); + if ( flags.StringFlagDefined("bcname") ) + { + int mybc = splines.Last()->bc-1; + if ( bcnames[mybc] ) delete bcnames[mybc]; + bcnames[mybc] = new string (flags.GetStringFlag("bcname","") ); + } + } + } + + + + void SplineGeometry2d :: PartitionBoundary (double h, Mesh & mesh2d) + { + enum { D = 2 }; + Box bbox; + GetBoundingBox (bbox); + double dist = Dist (bbox.PMin(), bbox.PMax()); + Point<3> pmin; + Point<3> pmax; - splines.Append(spline); -} + pmin(2) = -dist; pmax(2) = dist; + for(int j=0;j -void SplineGeometry :: AppendLineSegment (const int n1, const int n2, const int leftdomain, const int rightdomain, - const int bc, - const double reffac, const bool hprefleft, const bool hprefright, - const int copyfrom) -{ - SplineSeg * spline = new LineSeg(geompoints[n1],geompoints[n2]); - AppendSegment(spline,leftdomain,rightdomain,bc,reffac,hprefleft,hprefright,copyfrom); -} + Point3dTree searchtree (pmin, pmax); -template -void SplineGeometry :: AppendSplineSegment (const int n1, const int n2, const int n3, const int leftdomain, const int rightdomain, - const int bc, + for (int i = 0; i < splines.Size(); i++) + for (int side = 0; side <= 1; side++) + { + int dom = (side == 0) ? splines[i]->leftdom : splines[i]->rightdom; + if (dom != 0) splines[i] -> layer = GetDomainLayer (dom); + } + + for (int i = 0; i < splines.Size(); i++) + if (splines[i]->copyfrom == -1) + { + // astrid - set boundary meshsize to domain meshsize h + // if no domain mesh size is given, the max h value from the bounding box is used + double minimum = min2 ( GetDomainMaxh ( splines[i]->leftdom ), GetDomainMaxh ( splines[i]->rightdom ) ); + double maximum = max2 ( GetDomainMaxh ( splines[i]->leftdom ), GetDomainMaxh ( splines[i]->rightdom ) ); + minimum = min2 ( minimum, h ); + maximum = min2 ( maximum, h); + if ( minimum > 0 ) + splines[i]->Partition(minimum, elto0, mesh2d, searchtree, i+1); + else if ( maximum > 0 ) + splines[i]->Partition(maximum, elto0, mesh2d, searchtree, i+1); + else + splines[i]->Partition(h, elto0, mesh2d, searchtree, i+1); + } + else + { + CopyEdgeMesh (splines[i]->copyfrom, i+1, mesh2d, searchtree); + } + } + + + template + void SplineGeometry :: CopyEdgeMesh (int from, int to, Mesh & mesh, Point3dTree & searchtree) + { + int i; + + Array mappoints (mesh.GetNP()); + Array param (mesh.GetNP()); + mappoints = -1; + param = 0; + + Point3d pmin, pmax; + mesh.GetBox (pmin, pmax); + double diam2 = Dist2(pmin, pmax); + + if (printmessage_importance>0) + cout << "copy edge, from = " << from << " to " << to << endl; + + for (i = 1; i <= mesh.GetNSeg(); i++) + { + const Segment & seg = mesh.LineSegment(i); + if (seg.edgenr == from) + { + mappoints.Elem(seg[0]) = 1; + param.Elem(seg[0]) = seg.epgeominfo[0].dist; + + mappoints.Elem(seg[1]) = 1; + param.Elem(seg[1]) = seg.epgeominfo[1].dist; + } + } + + bool mapped = false; + for (i = 1; i <= mappoints.Size(); i++) + { + if (mappoints.Get(i) != -1) + { + Point newp = splines.Get(to)->GetPoint (param.Get(i)); + Point<3> newp3; + for(int j=0; jbc; + nseg[0] = mappoints.Get(seg[0]); + nseg[1] = mappoints.Get(seg[1]); + nseg.domin = splines.Get(to)->leftdom; + nseg.domout = splines.Get(to)->rightdom; + + nseg.epgeominfo[0].edgenr = to; + nseg.epgeominfo[0].dist = param.Get(seg[0]); + nseg.epgeominfo[1].edgenr = to; + nseg.epgeominfo[1].dist = param.Get(seg[1]); + mesh.AddSegment (nseg); + } + } + } + + + template + void SplineGeometry :: GetBoundingBox (Box & box) const + { + if (!splines.Size()) + { + Point auxp = 0.; + box.Set (auxp); + return; + } + + Array > points; + for (int i = 0; i < splines.Size(); i++) + { + splines[i]->GetPoints (20, points); + + if (i == 0) box.Set(points[0]); + for (int j = 0; j < points.Size(); j++) + box.Add (points[j]); + } + } + + template + void SplineGeometry :: SetGrading (const double grading) + { elto0 = grading;} + + /* + template + void SplineGeometry :: AppendPoint (const double x, const double y, const double reffac, const bool hpref) + { + geompoints.Append (GeomPoint(x, y, reffac)); + geompoints.Last().hpref = hpref; + } + */ + + template + void SplineGeometry :: AppendPoint (const Point & p, const double reffac, const bool hpref) + { + geompoints.Append (GeomPoint(p, reffac)); + geompoints.Last().hpref = hpref; + } + + + + + template + void SplineGeometry :: AppendSegment(SplineSeg * spline, const int leftdomain, const int rightdomain, + const int bc, + const double reffac, const bool hprefleft, const bool hprefright, + const int copyfrom) + { + spline -> leftdom = leftdomain; + spline -> rightdom = rightdomain; + spline -> bc = (bc >= 0) ? bc : (splines.Size()+1); + spline -> reffak = reffac; + spline -> hpref_left = hprefleft; + spline -> hpref_right = hprefright; + spline -> copyfrom = copyfrom; + + splines.Append(spline); + } + + template + void SplineGeometry :: AppendLineSegment (const int n1, const int n2, const int leftdomain, const int rightdomain, + const int bc, const double reffac, const bool hprefleft, const bool hprefright, const int copyfrom) -{ - SplineSeg * spline = new SplineSeg3(geompoints[n1],geompoints[n2],geompoints[n3]); - AppendSegment(spline,leftdomain,rightdomain,bc,reffac,hprefleft,hprefright,copyfrom); -} + { + SplineSeg * spline = new LineSeg(geompoints[n1],geompoints[n2]); + AppendSegment(spline,leftdomain,rightdomain,bc,reffac,hprefleft,hprefright,copyfrom); + } -template -void SplineGeometry :: AppendCircleSegment (const int n1, const int n2, const int n3, const int leftdomain, const int rightdomain, - const int bc, - const double reffac, const bool hprefleft, const bool hprefright, - const int copyfrom) -{ - SplineSeg * spline = new CircleSeg(geompoints[n1],geompoints[n2],geompoints[n3]); - AppendSegment(spline,leftdomain,rightdomain,bc,reffac,hprefleft,hprefright,copyfrom); -} + template + void SplineGeometry :: AppendSplineSegment (const int n1, const int n2, const int n3, const int leftdomain, const int rightdomain, + const int bc, + const double reffac, const bool hprefleft, const bool hprefright, + const int copyfrom) + { + SplineSeg * spline = new SplineSeg3(geompoints[n1],geompoints[n2],geompoints[n3]); + AppendSegment(spline,leftdomain,rightdomain,bc,reffac,hprefleft,hprefright,copyfrom); + } -template -void SplineGeometry :: AppendDiscretePointsSegment (const Array< Point > & points, const int leftdomain, const int rightdomain, - const int bc, - const double reffac, const bool hprefleft, const bool hprefright, - const int copyfrom) -{ - SplineSeg * spline = new DiscretePointsSeg(points); - AppendSegment(spline,leftdomain,rightdomain,bc,reffac,hprefleft,hprefright,copyfrom); -} + template + void SplineGeometry :: AppendCircleSegment (const int n1, const int n2, const int n3, const int leftdomain, const int rightdomain, + const int bc, + const double reffac, const bool hprefleft, const bool hprefright, + const int copyfrom) + { + SplineSeg * spline = new CircleSeg(geompoints[n1],geompoints[n2],geompoints[n3]); + AppendSegment(spline,leftdomain,rightdomain,bc,reffac,hprefleft,hprefright,copyfrom); + } + + template + void SplineGeometry :: AppendDiscretePointsSegment (const Array< Point > & points, const int leftdomain, const int rightdomain, + const int bc, + const double reffac, const bool hprefleft, const bool hprefright, + const int copyfrom) + { + SplineSeg * spline = new DiscretePointsSeg(points); + AppendSegment(spline,leftdomain,rightdomain,bc,reffac,hprefleft,hprefright,copyfrom); + } - template - void SplineGeometry :: GetMaterial( const int domnr, char* & material ) - { + template + void SplineGeometry :: GetMaterial( const int domnr, char* & material ) + { if ( materials.Size() >= domnr) material = materials[domnr-1]; else material = 0; - return; - } + return; + } - template - double SplineGeometry :: GetDomainMaxh( const int domnr ) - { - if ( maxh.Size() >= domnr && domnr > 0) - return maxh[domnr-1]; - else - return -1; - } + template + double SplineGeometry :: GetDomainMaxh( const int domnr ) + { + if ( maxh.Size() >= domnr && domnr > 0) + return maxh[domnr-1]; + else + return -1; + } - template - string SplineGeometry :: GetBCName( const int bcnr ) const - { - if ( bcnames.Size() >= bcnr) - if (bcnames[bcnr-1] ) - return *bcnames[bcnr-1]; - return "default"; - } + template + string SplineGeometry :: GetBCName( const int bcnr ) const + { + if ( bcnames.Size() >= bcnr) + if (bcnames[bcnr-1] ) + return *bcnames[bcnr-1]; + return "default"; + } -template -string * SplineGeometry :: BCNamePtr( const int bcnr ) -{ - if ( bcnr > bcnames.Size() ) + template + string * SplineGeometry :: BCNamePtr( const int bcnr ) + { + if ( bcnr > bcnames.Size() ) + return 0; + else + return bcnames[bcnr-1]; + } + + SplineGeometry2d :: ~SplineGeometry2d() + { + ; + } + + + extern void MeshFromSpline2D (SplineGeometry2d & geometry, + Mesh *& mesh, + MeshingParameters & mp); + + + int SplineGeometry2d :: GenerateMesh (Mesh*& mesh, MeshingParameters & mparam, + int perfstepsstart, int perfstepsend) + { + MeshFromSpline2D (*this, mesh, mparam); return 0; - else - return bcnames[bcnr-1]; -} - -SplineGeometry2d :: ~SplineGeometry2d() -{ - ; -} + } -extern void MeshFromSpline2D (SplineGeometry2d & geometry, - Mesh *& mesh, - MeshingParameters & mp); - - -int SplineGeometry2d :: GenerateMesh (Mesh*& mesh, MeshingParameters & mparam, - int perfstepsstart, int perfstepsend) -{ - cout << "SplineGeometry2d::GenerateMesh not only a dummy" << endl; - - MeshFromSpline2D (*this, mesh, mparam); - return 0; -} - - -Refinement & SplineGeometry2d :: GetRefinement () const -{ - return * new Refinement2d (*this); -} + Refinement & SplineGeometry2d :: GetRefinement () const + { + return * new Refinement2d (*this); + } diff --git a/libsrc/geom2d/splinegeometry.hpp b/libsrc/geom2d/splinegeometry.hpp index 11a2b82c..0a67ff89 100644 --- a/libsrc/geom2d/splinegeometry.hpp +++ b/libsrc/geom2d/splinegeometry.hpp @@ -15,7 +15,7 @@ in geom2d only 2D - Geometry classes (with material properties etc.) #ifndef _FILE_SPLINEGEOMETRY #define _FILE_SPLINEGEOMETRY -#include "../csg/csgparser.hpp" +// #include "../csg/csgparser.hpp" namespace netgen @@ -38,7 +38,8 @@ namespace netgen template < int D > class SplineGeometry { - protected: + // protected: + public: Array < GeomPoint > geompoints; Array < SplineSeg* > splines; double elto0; @@ -59,7 +60,7 @@ namespace netgen int Load (const Array & raw_data, const int startpos = 0); void Load (const char * filename); - void CSGLoad (CSGScanner & scan); + // void CSGLoad (CSGScanner & scan); void LoadData( ifstream & infile ); void LoadDataNew ( ifstream & infile ); diff --git a/libsrc/geom2d/vsgeom2d.cpp b/libsrc/geom2d/vsgeom2d.cpp index 71fabfb6..ac94f4db 100644 --- a/libsrc/geom2d/vsgeom2d.cpp +++ b/libsrc/geom2d/vsgeom2d.cpp @@ -3,9 +3,8 @@ #include #include -#include -#include +#include #include #include "vsgeom2d.hpp" @@ -14,9 +13,6 @@ namespace netgen { - - - /* *********************** Draw 2D Geometry **************** */ diff --git a/libsrc/occ/Makefile.am b/libsrc/occ/Makefile.am index 3df9f496..a3c7bd4b 100644 --- a/libsrc/occ/Makefile.am +++ b/libsrc/occ/Makefile.am @@ -20,5 +20,9 @@ libocc_la_SOURCES = Partition_Inter2d.cxx Partition_Inter3d.cxx \ Partition_Loop.cxx Partition_Loop2d.cxx Partition_Loop3d.cxx Partition_Spliter.cxx \ occconstruction.cpp occgenmesh.cpp occgeom.cpp occmeshsurf.cpp -liboccvis_la_SOURCES = occpkg.cpp vsocc.cpp +libocc_la_LIBADD = $(OCCLIBS) + +liboccvis_la_SOURCES = occpkg.cpp vsocc.cpp +liboccvis_la_LIBADD = libocc.la + diff --git a/libsrc/stlgeom/Makefile.am b/libsrc/stlgeom/Makefile.am index 6296628e..c8ea07aa 100644 --- a/libsrc/stlgeom/Makefile.am +++ b/libsrc/stlgeom/Makefile.am @@ -4,12 +4,12 @@ stltool.hpp stltopology.hpp vsstl.hpp AM_CPPFLAGS = -I$(top_srcdir)/libsrc/include $(TCL_INCLUDES) METASOURCES = AUTO -noinst_LTLIBRARIES = libstl.la libstlvis.la +lib_LTLIBRARIES = libstl.la libstlvis.la libstl_la_SOURCES = meshstlsurface.cpp stlgeom.cpp stlgeomchart.cpp \ stlgeommesh.cpp stlline.cpp stltool.cpp stltopology.cpp libstlvis_la_SOURCES = stlpkg.cpp vsstl.cpp - +libstlvis_la_LIBADD = libstl.la $(top_builddir)/libsrc/linalg/libla.la diff --git a/libsrc/stlgeom/stlpkg.cpp b/libsrc/stlgeom/stlpkg.cpp index f3c18598..79c416dc 100644 --- a/libsrc/stlgeom/stlpkg.cpp +++ b/libsrc/stlgeom/stlpkg.cpp @@ -595,8 +595,8 @@ namespace netgen using namespace netgen; -extern "C" int Ng_STL_Init (Tcl_Interp * interp); -int Ng_STL_Init (Tcl_Interp * interp) +extern "C" int Ng_stl_Init (Tcl_Interp * interp); +int Ng_stl_Init (Tcl_Interp * interp) { geometryregister.Append (new STLGeometryRegister); diff --git a/ng/Makefile.am b/ng/Makefile.am index 9971de67..eaaaec57 100644 --- a/ng/Makefile.am +++ b/ng/Makefile.am @@ -1,6 +1,6 @@ include_HEADERS = -AM_CPPFLAGS = -I$(top_srcdir)/libsrc/include -I$(top_srcdir)/libsrc/interface -DOPENGL -D$(TOGL_WINDOWINGSYSTEM) $(OCCFLAGS) $(TCL_INCLUDES) $(MPI_INCLUDES) $(FFMPEG_INCLUDES) $(JPEGLIB_INCLUDES) +AM_CPPFLAGS = -I$(top_srcdir)/libsrc/include -I$(top_srcdir)/libsrc/interface -DOPENGL -D$(TOGL_WINDOWINGSYSTEM) $(TCL_INCLUDES) $(MPI_INCLUDES) $(FFMPEG_INCLUDES) $(JPEGLIB_INCLUDES) bin_PROGRAMS = netgen netgen_SOURCES = demoview.cpp ngappinit.cpp ngpkg.cpp onetcl.cpp nginterface.cpp nginterface_v2.cpp parallelfunc.cpp parallelinterface.cpp demoview.hpp parallelfunc.hpp togl_1_7.h @@ -8,20 +8,19 @@ netgen_SOURCES = demoview.cpp ngappinit.cpp ngpkg.cpp onetcl.cpp nginterface.cpp netgen_LDADD = $(top_builddir)/libsrc/visualization/libvisual.a \ $(top_builddir)/libsrc/csg/libcsgvis.la \ - $(top_builddir)/libsrc/csg/libcsg.la \ $(top_builddir)/libsrc/geom2d/libgeom2dvis.la \ - $(top_builddir)/libsrc/geom2d/libgeom2d.la \ $(top_builddir)/libsrc/interface/libinterface.la \ - $(top_builddir)/libsrc/stlgeom/libstlvis.la \ - $(top_builddir)/libsrc/stlgeom/libstl.la \ $(top_builddir)/libsrc/meshing/libmesh.la \ $(top_builddir)/libsrc/gprim/libgprim.la \ $(top_builddir)/libsrc/linalg/libla.la \ $(top_builddir)/libsrc/general/libgen.la \ - $(OCCLIBS) -L$(TK_BIN_DIR)/Togl1.7 $(TOGLLIBDIR) -lTogl1.7 $(LIBGLU) $(TK_LIB_SPEC) $(TCL_LIB_SPEC) $(MPI_LIBS) $(FFMPEG_LIBS) $(JPEGLIB_LIBS) $(PKG_LIBS) + -L$(TK_BIN_DIR)/Togl1.7 $(TOGLLIBDIR) -lTogl1.7 $(LIBGLU) $(TK_LIB_SPEC) $(TCL_LIB_SPEC) $(MPI_LIBS) $(FFMPEG_LIBS) $(JPEGLIB_LIBS) $(PKG_LIBS) # # $(top_builddir)/libsrc/occ/liboccvis.la # $(top_builddir)/libsrc/occ/libocc.la +# $(top_builddir)/libsrc/stlgeom/libstlvis.la +# $(top_builddir)/libsrc/stlgeom/libstl.la +# $(top_builddir)/libsrc/geom2d/libgeom2d.la # add for static linkage of ngsolve: diff --git a/ng/dialog.tcl b/ng/dialog.tcl index 06ca0599..10b6f343 100644 --- a/ng/dialog.tcl +++ b/ng/dialog.tcl @@ -2494,50 +2494,6 @@ proc printlatestwarning { } { } -# for parallel visualization, overlapping meshes... -proc paralleldialog { } { - - set w .parallel_dlg - - if {[winfo exists .parallel_dlg] == 1} { - wm withdraw $w - wm deiconify $w - wm geometry $w =270x100 - - focus $w - } { - - toplevel $w - wm geometry $w =270x100 - -# frame $w.buttons -relief groove -borderwidth 3 -width 300 -# pack $w.buttons - set ww $w - - button $ww.visallb -text "View All" -width 20 -command\ - { Ng_VisualizeAll; } - pack $ww.visallb - - button $ww.visoneb -text "View One" -width 20 -command \ - { Ng_VisualizeOne; } - pack $ww.visoneb - - button $ww.overlap -text "overlap++" -width 20 -command \ - { Ng_IncrOverlap; } - - pack $ww.overlap - - wm withdraw $w - wm geom $w +100+100 - wm deiconify $w - wm title $w "Parallel Netgen" - focus .parallel_dlg - } -} - - -# paralleldialog -#wm withdraw $w proc runtestdialog { } { source $::ngdir/ngshell.tcl diff --git a/ng/drawing.tcl b/ng/drawing.tcl index 4d4ed9e4..a7219aaa 100644 --- a/ng/drawing.tcl +++ b/ng/drawing.tcl @@ -5,15 +5,16 @@ set oldmousex 0 set oldmousey 0 # -# if { 1 } { # use this one for Togl 2.0 +# if { 1 } { + + + # if {[catch {togl .ndraw -width 400 -height 300 -rgba true -double true -depth true -privatecmap false -stereo false -indirect true -create init -display draw -reshape reshape }] } { -# changed -indirect true/false !!! if {[catch {togl .ndraw -width 400 -height 300 -rgba true -double true -depth true -privatecmap false -stereo false -indirect true }] } { - puts "no OpenGL" } { # diff --git a/ng/menustat.tcl b/ng/menustat.tcl index d6083882..49b011f0 100644 --- a/ng/menustat.tcl +++ b/ng/menustat.tcl @@ -833,17 +833,6 @@ pack .bubar.exitb .bubar.surfm .bubar.stopm -side left #button .bubar.scan -text "Scan" \ # -command { Ng_ParseGeometry; set selectvisual geometry; Ng_SetVisParameters; redraw } -# fuer parallel - buttons :) -Ng_IsParallel; -if { $parallel_netgen } { -# catch{ -# source ${ngdir}/ngtcltk/parallel_dialog.tcl -# } - button .bubar.visallb -text "Parallel" -command \ - { paralleldialog; redraw } - pack .bubar.visallb -side left -} - button .bubar.zoomall -text "Zoom All" \ -command { Ng_ZoomAll; redraw } diff --git a/ng/ng.tcl b/ng/ng.tcl index 6c0d46d4..43904463 100644 --- a/ng/ng.tcl +++ b/ng/ng.tcl @@ -1,6 +1,6 @@ if {[catch {package require Tix } result ]} { - puts "cannot find package Tix" - puts "error : $result" + puts "cannot find package Tix" + puts "error : $result" } # if {[catch {package require Togl 2.0 } result ]} { diff --git a/ng/ngappinit.cpp b/ng/ngappinit.cpp index 5de522ce..a08f9ba1 100644 --- a/ng/ngappinit.cpp +++ b/ng/ngappinit.cpp @@ -101,7 +101,7 @@ int main(int argc, char ** argv) cout << "NETGEN-" << PACKAGE_VERSION << endl; cout << "Developed by Joachim Schoeberl at" << endl - << "2010-xxxx Vienna UT" << endl + << "2010-xxxx Vienna University of Technology" << endl << "2006-2010 RWTH Aachen University" << endl << "1996-2006 Johannes Kepler University Linz" << endl; diff --git a/ng/nginterface_v2.cpp b/ng/nginterface_v2.cpp index 4ae21158..d380bf63 100644 --- a/ng/nginterface_v2.cpp +++ b/ng/nginterface_v2.cpp @@ -1,20 +1,7 @@ #include - - #include -#include -#include -#include -#ifdef OCCGEOMETRY -#include -#endif - -#ifdef ACIS -#include -#endif - #ifdef SOCKETS #include "../sockets/sockets.hpp" #endif @@ -27,11 +14,6 @@ #include "nginterface.h" #include "nginterface_v2.hpp" -// #include - - -// #include - namespace netgen { @@ -43,16 +25,6 @@ namespace netgen extern Tcl_Interp * tcl_interp; #endif - extern AutoPtr geometry2d; - extern AutoPtr geometry; - extern STLGeometry * stlgeometry; - -#ifdef OCCGEOMETRY - extern OCCGeometry * occgeometry; -#endif -#ifdef ACIS - extern ACISGeometry * acisgeometry; -#endif #ifdef OPENGL extern VisualSceneSolution vssolution; @@ -217,23 +189,6 @@ namespace netgen double * dxdxi, size_t sdxdxi) { mesh->GetCurvedElements().CalcMultiPointSegmentTransformation<2> (elnr, npts, xi, sxi, x, sx, dxdxi, sdxdxi); - /* - for (int ip = 0; ip < npts; ip++) - { - Point<3> xg; - Vec<3> dx; - - mesh->GetCurvedElements().CalcSegmentTransformation (xi[ip*sxi], elnr, xg, dx); - - if (x) - for (int i = 0; i < 2; i++) - x[ip*sx+i] = xg(i); - - if (dxdxi) - for (int i=0; i<2; i++) - dxdxi[ip*sdxdxi+i] = dx(i); - } - */ } template <> @@ -257,7 +212,6 @@ namespace netgen return mesh->GetTopology().GetNFaces(); } - template <> DLL_HEADER Ng_Node<1> Ng_GetNode<1> (int nr) { Ng_Node<1> node; @@ -265,8 +219,6 @@ namespace netgen return node; } - - template <> DLL_HEADER Ng_Node<2> Ng_GetNode<2> (int nr) { Ng_Node<2> node; diff --git a/ng/ngpkg.cpp b/ng/ngpkg.cpp index 9688b295..88bfa90d 100644 --- a/ng/ngpkg.cpp +++ b/ng/ngpkg.cpp @@ -24,7 +24,6 @@ The interface between the GUI and the netgen library #include "../libsrc/sockets/socketmanager.hpp" #endif -// #include // to be sure to include the 'right' togl-version #include "togl_1_7.h" @@ -2012,10 +2011,10 @@ namespace netgen return TCL_ERROR; cout << "call Togl - load font (crash on my Linux64)" << endl; - // togl_font = Togl_LoadBitmapFont( togl, "Times"); // TOGL_BITMAP_8_BY_13 ); + togl_font = Togl_LoadBitmapFont( togl, "Times"); // TOGL_BITMAP_8_BY_13 ); // togl_font = Togl_LoadBitmapFont( togl, TOGL_BITMAP_8_BY_13 ); // togl_font = Togl_LoadBitmapFont( togl, NULL ); - // cout << "success" << endl; + cout << "success" << endl; glMatrixMode(GL_PROJECTION); glLoadIdentity(); @@ -2882,91 +2881,7 @@ namespace netgen return TCL_OK; } -#ifdef PARALLEL - int Ng_VisualizeAll (ClientData clientData, - Tcl_Interp * interp, - int argc, tcl_const char *argv[]) - { - int id, rc, ntasks; - MPI_Comm_size(MPI_COMM_WORLD, &ntasks); - MPI_Comm_rank(MPI_COMM_WORLD, &id); - string visualizationmode = Tcl_GetVar (interp, "::selectvisual", 0); - string scalfun = Tcl_GetVar (interp, "::visoptions.scalfunction", 0); - for ( int dest = 1; dest < ntasks; dest++) - { - MyMPI_Send ( "visualize", dest ); - MyMPI_Send ( visualizationmode, dest); - if ( visualizationmode == "solution" ) - MyMPI_Send ( scalfun, dest); - } - return TCL_OK; - } - - int Ng_VisualizeOne (ClientData clientData, - Tcl_Interp * interp, - int argc, tcl_const char *argv[]) - { - int id, rc, ntasks; - MPI_Comm_size(MPI_COMM_WORLD, &ntasks); - MPI_Comm_rank(MPI_COMM_WORLD, &id); - - string visualizationmode = Tcl_GetVar (interp, "::selectvisual", 0); - string scalfun = Tcl_GetVar (interp, "::visoptions.scalfunction", 0); - - MyMPI_Send ( "visualize", 1 ); - MyMPI_Send ( visualizationmode, 1); - - if ( visualizationmode == "solution" ) - MyMPI_Send ( scalfun, 1); - return TCL_OK; - } - - int Ng_IncrOverlap ( ClientData clientDate, - Tcl_Interp * interp, - int argc, tcl_const char * argv[] ) - { - int id, rc, ntasks; - MPI_Comm_size(MPI_COMM_WORLD, &ntasks); - MPI_Comm_rank(MPI_COMM_WORLD, &id); - - for ( int dest = 1; dest < ntasks; dest++) - { - MyMPI_Send ( "overlap++", dest ); - } - mesh->UpdateOverlap(); - return TCL_OK; - - } - - int Ng_SetSelectVisual ( ClientData clientDate, - Tcl_Interp * interp, - int argc, tcl_const char * argv[] ) - { - string visualizationmode; - MyMPI_Recv ( visualizationmode, 0); - Tcl_SetVar (interp, "::selectvisual", visualizationmode.c_str(), 0); - return TCL_OK; - } - - int Ng_SetScalarFunction ( ClientData clientDate, - Tcl_Interp * interp, - int argc, tcl_const char * argv[] ) - { - string visualizationmode; - string scalarfun; - visualizationmode = Tcl_GetVar (interp, "::selectvisual", 0); - - if ( visualizationmode == "solution" ) - { - MyMPI_Recv ( scalarfun, 0); - Tcl_SetVar (interp, "::visoptions.scalfunction", scalarfun.c_str(), 0); - } - return TCL_OK; - } - - -#endif int Ng_IsParallel (ClientData clientData, Tcl_Interp * interp, @@ -3073,7 +2988,7 @@ namespace netgen extern "C" int Ng_Init (Tcl_Interp * interp); extern "C" int Ng_CSG_Init (Tcl_Interp * interp); - extern "C" int Ng_STL_Init (Tcl_Interp * interp); + // extern "C" int Ng_stl_Init (Tcl_Interp * interp); #ifdef OCCGEOMETRY // extern "C" int Ng_occ_Init (Tcl_Interp * interp); @@ -3094,11 +3009,8 @@ namespace netgen #endif Ng_CSG_Init(interp); - Ng_STL_Init(interp); + // Ng_stl_Init(interp); -#ifdef OCCGEOMETRY - // Ng_occ_Init(interp); -#endif Ng_Geom2d_Init(interp); @@ -3358,29 +3270,6 @@ namespace netgen (ClientData)NULL, (Tcl_CmdDeleteProc*) NULL); - -#ifdef PARALLEL - Tcl_CreateCommand (interp, "Ng_VisualizeAll", Ng_VisualizeAll, - (ClientData)NULL, - (Tcl_CmdDeleteProc*) NULL); - - Tcl_CreateCommand (interp, "Ng_VisualizeOne", Ng_VisualizeOne, - (ClientData)NULL, - (Tcl_CmdDeleteProc*) NULL); - - Tcl_CreateCommand (interp, "Ng_IncrOverlap", Ng_IncrOverlap, - (ClientData)NULL, - (Tcl_CmdDeleteProc*) NULL); - - Tcl_CreateCommand (interp, "Ng_SetSelectVisual", Ng_SetSelectVisual, - (ClientData)NULL, - (Tcl_CmdDeleteProc*) NULL); - - Tcl_CreateCommand (interp, "Ng_SetScalarFunction", Ng_SetScalarFunction, - (ClientData)NULL, - (Tcl_CmdDeleteProc*) NULL); - -#endif Tcl_CreateCommand (interp, "Ng_IsParallel", Ng_IsParallel, (ClientData)NULL, (Tcl_CmdDeleteProc*) NULL); diff --git a/ng/occgeom.tcl b/ng/occgeom.tcl index 2bfcbaaf..efab5c33 100644 --- a/ng/occgeom.tcl +++ b/ng/occgeom.tcl @@ -1,6 +1,6 @@ if { [catch { load liboccvis[info sharedlibextension] Ng_OCC } result ] } { -# puts "cannot load occ" -# puts "error: $result" + puts "cannot load occ" + puts "error: $result" # dummy proc rebuildoccdialog { } { } diff --git a/ng/stlgeom.tcl b/ng/stlgeom.tcl index 664805f0..0f95572c 100644 --- a/ng/stlgeom.tcl +++ b/ng/stlgeom.tcl @@ -1,3 +1,9 @@ +if { [catch { load libstlvis[info sharedlibextension] Ng_STL } result ] } { + puts "cannot load stl" + puts "error: $result" +} + + .ngmenu.geometry add separator .ngmenu.geometry add command -label "STL Doctor..." \