/[escript]/trunk/pycad/py_src/primitives.py
ViewVC logotype

Contents of /trunk/pycad/py_src/primitives.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 916 - (show annotations)
Fri Dec 15 08:25:57 2006 UTC (12 years, 9 months ago) by gross
File MIME type: text/x-python
File size: 23858 byte(s)
curves implemenred and tested.
1 # $Id:$
2
3 """
4 Geometrical Primitives
5
6 the concept is inspired by gmsh and very much focused on the fact that
7 the classes are used to wrk with gmsh.
8
9 @var __author__: name of author
10 @var __copyright__: copyrights
11 @var __license__: licence agreement
12 @var __url__: url entry point on documentation
13 @var __version__: version
14 @var __date__: date of the version
15 """
16
17
18 __author__="Lutz Gross, l.gross@uq.edu.au"
19 __copyright__=""" Copyright (c) 2006 by ACcESS MNRF
20 http://www.access.edu.au
21 Primary Business: Queensland, Australia"""
22 __license__="""Licensed under the Open Software License version 3.0
23 http://www.opensource.org/licenses/osl-3.0.php"""
24 __url__="http://www.iservo.edu.au/esys/escript"
25 __version__="$Revision:$"
26 __date__="$Date:$"
27
28 import numarray
29 from transformations import _TYPE, Translation, Dilation, Transformation
30
31
32 def resetGlobalPrimitiveIdCounter():
33 global global_primitive_id_counter
34 global_primitive_id_counter=1
35
36 def setToleranceForColocation(tol=1.e-11):
37 global global_tolerance_for_colocation
38 global_tolerance_for_colocation=tol
39
40 def getToleranceForColocation():
41 return global_tolerance_for_colocation
42
43 resetGlobalPrimitiveIdCounter()
44 setToleranceForColocation()
45
46 class Primitive(object):
47 """
48 template for elementary geometrical object
49 """
50 def __init__(self):
51 """
52
53 """
54 global global_primitive_id_counter
55 self.__ID=global_primitive_id_counter
56 global_primitive_id_counter+=1
57
58 def getID(self):
59 return self.__ID
60
61 def __repr__(self):
62 return "%s(%s)"%(self.__class__.__name__,self.getID())
63 def __cmp__(self,other):
64 return cmp(self.getID(),other.getID())
65
66 def getConstructionPoints(self):
67 """
68 returns the points used to construct the primitive
69 """
70 out=set()
71 for i in self.getPrimitives():
72 if isinstance(i,Point): out.add(i)
73 return list(out)
74
75 def getPrimitives(self):
76 """
77 returns primitives used to construct the primitive
78 """
79 return []
80
81 def copy(self):
82 """
83 returns a deep copy of the object
84 """
85 return Primitive()
86
87
88 def modifyBy(self,transformation):
89 """
90 modifies the coordinates by applying a transformation
91 """
92 for p in self.getConstructionPoints(): p.modifyBy(transformation)
93
94
95 def __add__(self,other):
96 """
97 returns a new object shifted by other
98 """
99 return self.apply(Translation(numarray.array(other,_TYPE)))
100
101 def __sub__(self,other):
102 """
103 returns a new object shifted by other
104 """
105 return self.apply(Translation(-numarray.array(other,_TYPE)))
106
107 def __iadd__(self,other):
108 """
109 shifts the point by other
110 """
111 self.modifyBy(Translation(numarray.array(other,_TYPE)))
112 return self
113
114 def __isub__(self,other):
115 """
116 shifts the point by -other
117 """
118 self.modifyBy(Translation(-numarray.array(other,_TYPE)))
119 return self
120
121 def __imul__(self,other):
122 """
123 modifies object by applying L{Transformation} other. If other is not a L{Transformation} it will try convert it.
124 """
125 if isinstance(other,int) or isinstance(other,float):
126 trafo=Dilation(other)
127 elif isinstance(other,numarray.NumArray):
128 trafo=Translation(other)
129 elif isinstance(other,Transformation):
130 trafo=other
131 else:
132 raise TypeError, "cannot convert argument to Trnsformation class object."
133 self.modifyBy(trafo)
134 return self
135
136 def __rmul__(self,other):
137 """
138 applies L{Transformation} other to object. If other is not a L{Transformation} it will try convert it.
139 """
140 if isinstance(other,int) or isinstance(other,float):
141 trafo=Dilation(other)
142 elif isinstance(other,numarray.NumArray):
143 trafo=Translation(other)
144 elif isinstance(other,Transformation):
145 trafo=other
146 else:
147 raise TypeError, "cannot convert argument to Transformation class object."
148 return self.apply(trafo)
149
150 def __neg__(self):
151 return ReversedPrimitive(self)
152
153 def setLocalScale(self,factor=1.):
154 """
155 sets the local refinement factor
156 """
157 for p in self.getConstructionPoints(): p.setLocalScale(factor)
158
159 def getGmshCommand(self, local_scaling_factor=1.):
160 """
161 returns the Gmsh command(s) to create the primitive
162 """
163 raise NotImplementedError("getGmshCommand is not implemented for this class %s."%self.__class__.__name__)
164
165 def apply(self,transformation):
166 """
167 returns a new object by applying the transformation
168 """
169 raise NotImplementedError("apply is not implemented for this class %s."%self.__class__.__name__)
170
171 def isColocated(self,primitive):
172 """
173 returns True is the two primitives are located at the smae position
174 """
175 raise NotImplementedError("isColocated is not implemented for this class %s."%self.__class__.__name__)
176
177 class Point(Primitive):
178 """
179 a three dimensional point
180 """
181 def __init__(self,x=0.,y=0.,z=0.,local_scale=1.):
182 """
183 creates a point with coorinates x,y,z with the local refinement factor local_scale
184 """
185 super(Point, self).__init__()
186 self.setCoordinates(numarray.array([x,y,z],_TYPE))
187 self.setLocalScale(local_scale)
188
189 def setLocalScale(self,factor=1.):
190 """
191 sets the local refinement factor
192 """
193 if factor<=0.:
194 raise ValueError("scaling factor must be positive.")
195 self.__local_scale=factor
196 def getLocalScale(self):
197 """
198 returns the local refinement factor
199 """
200 return self.__local_scale
201 def getCoordinates(self):
202 """
203 returns the coodinates of the point as L{numarray.NumArray} object
204 """
205 return self._x
206 def setCoordinates(self,x):
207 """
208 returns the coodinates of the point as L{numarray.NumArray} object
209 """
210 if not isinstance(x, numarray.NumArray):
211 self._x=numarray.array(x,_TYPE)
212 else:
213 self._x=x
214
215 def getPrimitives(self):
216 """
217 returns primitives used to construct the primitive
218 """
219 return [self]
220
221 def isColocated(self,primitive):
222 """
223 returns True if L{Point} primitive is colocation (same coordinates)
224 that means if |self-primitive| <= tol * max(|self|,|primitive|)
225 """
226 if isinstance(primitive,Point):
227 primitive=primitive.getCoordinates()
228 c=self.getCoordinates()
229 d=c-primitive
230 return numarray.dot(d,d)<=getToleranceForColocation()**2*max(numarray.dot(c,c),numarray.dot(primitive,primitive))
231 else:
232 return False
233
234 def copy(self):
235 """
236 returns a deep copy of the point
237 """
238 c=self.getCoordinates()
239 return Point(c[0],c[1],c[2],local_scale=self.getLocalScale())
240
241 def modifyBy(self,transformation):
242 """
243 modifies the coordinates by applying a transformation
244 """
245 self.setCoordinates(transformation(self.getCoordinates()))
246
247 def apply(self,transformation):
248 """
249 returns a new L{Point} by applying the transformation
250 """
251 new_p=self.copy()
252 new_p.modifyBy(transformation)
253 return new_p
254
255 def getGmshCommand(self, local_scaling_factor=1.):
256 """
257 returns the Gmsh command(s) to create the primitive
258 """
259 c=self.getCoordinates()
260 return "Point(%s) = {%s , %s, %s , %s };"%(self.getID(),c[0],c[1],c[2], self.getLocalScale()*local_scaling_factor)
261
262 class Primitive1D(Primitive):
263 """
264 general one-dimensional primitive
265 """
266 def __init__(self,*args):
267 """
268 create a one-dimensional primitive
269 """
270 super(Primitive1D, self).__init__()
271
272 class Curve(Primitive1D):
273 """
274 a curve defined through a list of control points.
275 """
276 def __init__(self,*points):
277 """
278 defines a curve form control points
279 """
280 if len(points)<2:
281 raise TypeError("Curve needs at least two points")
282 super(Curve, self).__init__()
283 i=0
284 for p in points:
285 i+=1
286 if not isinstance(p,Point): raise TypeError("%s-th argument is not a Point object."%i)
287 self.__points=points
288
289 def __len__(self):
290 """
291 returns the number of control points
292 """
293 return len(self.__points)
294
295 def getStartPoint(self):
296 """
297 returns start point
298 """
299 return self.__points[0]
300
301 def getEndPoint(self):
302 """
303 returns end point
304 """
305 return self.__points[-1]
306
307 def getControlPoints(self):
308 """
309 returns a list of the points
310 """
311 return self.__points
312
313 def getPrimitives(self):
314 """
315 returns primitives used to construct the Curve
316 """
317 out=set()
318 for p in self.getControlPoints(): out|=set(p.getPrimitives())
319 out.add(self)
320 return list(out)
321
322 def copy(self):
323 """
324 returns a deep copy
325 """
326 new_p=[]
327 for p in self.getControlPoints(): new_p.append(p.copy())
328 return self.__class__(*tuple(new_p))
329
330
331 def apply(self,transformation):
332 """
333 applies transformation
334 """
335 new_p=[]
336 for p in self.getControlPoints(): new_p.append(p.apply(transformation))
337 return self.__class__(*tuple(new_p))
338
339 def isColocated(self,primitive):
340 """
341 returns True curves are on the same position
342 """
343 if isinstance(primitive,self.__class__):
344 if len(primitive) == len(self):
345 cp0=self.getControlPoints()
346 cp1=primitive.getControlPoints()
347 for i in range(len(cp0)):
348 if not cp0[i].isColocated(cp1[i]):
349 return False
350 return True
351 else:
352 return False
353 else:
354 return False
355
356 class Spline(Curve):
357 """
358 a spline curve defined through a list of control points.
359 """
360 def getGmshCommand(self):
361 """
362 returns the Gmsh command(s) to create the Curve
363 """
364 out=""
365 for i in self.getControlPoints():
366 if len(out)>0:
367 out+=", %s"%i.getID()
368 else:
369 out="%s"%i.getID()
370 return "Spline(%s) = {%s};"%(self.getID(),out)
371
372
373 class BezierCurve(Curve):
374 """
375 a Bezier curve
376 """
377 def getGmshCommand(self):
378 """
379 returns the Gmsh command(s) to create the Curve
380 """
381 out=""
382 for i in self.getControlPoints():
383 if len(out)>0:
384 out+=", %s"%i.getID()
385 else:
386 out="%s"%i.getID()
387 return "Bezier(%s) = {%s};"%(self.getID(),out)
388
389 class BSpline(Curve):
390 """
391 a BSpline curve. Control points may be repeated.
392 """
393 def getGmshCommand(self):
394 """
395 returns the Gmsh command(s) to create the Curve
396 """
397 out=""
398 for i in self.getControlPoints():
399 if len(out)>0:
400 out+=", %s"%i.getID()
401 else:
402 out="%s"%i.getID()
403 return "BSpline(%s) = {%s};"%(self.getID(),out)
404
405 class Line(Curve):
406 """
407 a line is defined by two L{Point}s
408 """
409 def __init__(self,*points):
410 """
411 defines a line with start and end point
412 """
413 if len(points)!=2:
414 raise TypeError("Line needs two points")
415 super(Line, self).__init__(*points)
416 def getGmshCommand(self):
417 """
418 returns the Gmsh command(s) to create the Curve
419 """
420 return "Line(%s) = {%s, %s};"%(self.getID(),self.getStartPoint().getID(),self.getEndPoint().getID())
421
422
423 class Arc(Primitive1D):
424 """
425 defines an arc
426 """
427 def __init__(self,center,start,end):
428 """
429 creates an arc by the start point, end point and center
430 """
431 if isinstance(center,Point): raise TypeError("center needs to be a Point object.")
432 if isinstance(end,Point): raise TypeError("end needs to be a Point object.")
433 if isinstance(start,Point): raise TypeError("start needs to be a Point object.")
434 super(Arc, self).__init__()
435 self.__center=center
436 self.__start=start
437 self.__end=end
438
439 def getStartPoint(self):
440 """
441 returns start point
442 """
443 return self.__start
444
445 def getEndPoint(self):
446 """
447 returns end point
448 """
449 return self.__end
450
451 def getCenterPoint(self):
452 """
453 returns center
454 """
455 return self.__center
456
457 def getPrimitives(self):
458 """
459 returns C{set} of primitive used to construct the Curve
460 """
461 out=set()
462 out|=set(self.getStartPoint().getPrimitives())
463 out|=set(self.getEndPoint().getPrimitives())
464 out|=set(self.getCenterPoint().getPrimitives())
465 out.add(self)
466 return list(out)
467
468 def getGmshCommand(self):
469 """
470 returns the Gmsh command(s) to create the primitive
471 """
472 return "Circle(%s) = {%s, %s, %s};"%(self.getID(),self.getStartPoint().getID(),self.getCenterPoint().getID(),self.getEndPoint().getID())
473
474 def copy(self):
475 """
476 returns a deep copy
477 """
478 return Arc(self.getCenterPoint().copy(),self.getStartPoint().copy(),self.getEndPoint().copy())
479
480 def apply(self,transformation):
481 """
482 applies transformation
483 """
484 return Arc(self.getCenterPoint().apply(transformation),self.getStartPoint().apply(transformation),self.getEndPoint().apply(transformation))
485
486 def isColocated(self,primitive):
487 """
488 returns True curves are on the same position
489 """
490 if isinstance(primitive,Arc):
491 if not self.getCenterPoint().isColocated(primitive.getCenterPoint()): return False
492 if not self.getEndPoint().isColocated(primitive.getEndPoint()): return False
493 if not self.getStartPoint().isColocated(primitive.getStartPoint()): return False
494 return True
495 else:
496 return False
497
498
499 #=================================================================================================================================
500 class CurveLoop(Primitive):
501 """
502 An oriented loop of curves.
503
504 The loop must be closed and the L{Curves}s should be oriented consistently.
505 """
506 def __init__(self,*curves):
507 """
508 creates a polygon from a list of line curves. The curves must form a closed loop.
509 """
510 super(CurveLoop, self).__init__()
511 self.__curves=[]
512 self.addCurve(*curves)
513 def addCurve(self,*curves):
514 for i in range(len(curves)):
515 if not curves[i].isCurve():
516 raise TypeError("%s-th argument is not a Curve object."%i)
517 self.__curves+=curves
518
519 def isCurveLoop(self):
520 return True
521 def getCurves(self):
522 return self.__curves
523 def __add__(self,other):
524 return CurveLoop(*tuple([c+other for c in self.getCurves()[::-1]]))
525 def __len__(self):
526 return len(self.__curves)
527 def getPrimitives(self):
528 out=set([self])
529 for i in self.getCurves(): out|=i.getPrimitives()
530 return out
531 def getConstructionPoints(self):
532 out=set()
533 for i in self.getCurves(): out|=i.getConstructionPoints()
534 return out
535 def getGmshCommand(self):
536 out=""
537 for i in self.getCurves():
538 if len(out)>0:
539 out+=", %s"%i.getID()
540 else:
541 out="%s"%i.getID()
542 return "Line Loop(%s) = {%s};"%(self.getID(),out)
543
544 class Surface(Primitive):
545 """
546 a surface
547 """
548 def __init__(self,loop):
549 """
550 creates a surface with boundary loop
551
552 @param loop: L{CurveLoop} defining the boundary of the surface
553 """
554 super(Surface, self).__init__()
555 if not loop.isCurveLoop():
556 raise TypeError("argument loop needs to be a CurveLoop object.")
557 self.__loop=loop
558 def isSurface(self):
559 return True
560 def getBoundaryLoop(self):
561 return self.__loop
562 def __add__(self,other):
563 return Surface(self.getBoundaryLoop()+other)
564 def getPrimitives(self):
565 out=set([self]) | self.getBoundaryLoop().getPrimitives()
566 return out
567 def getConstructionPoints(self):
568 return self.getBoundaryLoop().getConstructionPoints()
569 def getGmshCommand(self):
570 return "Ruled Surface(%s) = {%s};"%(self.getID(),self.getBoundaryLoop().getID())
571
572 class PlaneSurface(Surface):
573 """
574 a plane surface with holes
575 """
576 def __init__(self,loop,holes=[]):
577 """
578 creates a plane surface.
579
580 @param loop: L{CurveLoop} defining the boundary of the surface
581 @param holes: list of L{CurveLoop} defining holes in the surface.
582 @note: A CurveLoop defining a hole should not have any lines in common with the exterior CurveLoop.
583 A CurveLoop defining a hole should not have any lines in common with another CurveLoop defining a hole in the same surface.
584 """
585 super(PlaneSurface, self).__init__(loop)
586 for i in range(len(holes)):
587 if not holes[i].inCurveLoop():
588 raise TypeError("%i th hole needs to be a CurveLoop object.")
589 self.__holes=holes
590 def getHoles(self):
591 return self.__holes
592 def __add__(self,other):
593 return PlaneSurface(self.getBoundaryLoop()+other, holes=[h+other for h in self.getHoles()])
594 def getPrimitives(self):
595 out=set([self]) | self.getBoundaryLoop().getPrimitives()
596 for i in self.getHoles(): out|=i.getPrimitives()
597 return out
598 def getConstructionPoints(self):
599 out=self.getBoundaryLoop().getConstructionPoints()
600 for i in self.getHoles(): out|=i.getConstructionPoints()
601 return out
602 def getGmshCommand(self):
603 out=""
604 for i in self.getHoles():
605 if len(out)>0:
606 out+=", %s"%i.getID()
607 else:
608 out="%s"%i.getID()
609 if len(out)>0:
610 return "Plane Surface(%s) = {%s, %s};"%(self.getID(),self.getBoundaryLoop().getID(), out)
611 else:
612 return "Plane Surface(%s) = {%s};"%(self.getID(),self.getBoundaryLoop().getID())
613
614 class RuledSurface(Surface):
615 """
616 A ruled surface, i.e., a surface that can be interpolated using transfinite interpolation
617 """
618 def __init__(self,loop):
619 """
620 creates a ruled surface from a
621
622 @param loop: L{CurveLoop} defining the boundary of the surface. There is a restriction of composed of either three or four L{Curve} objects.
623 """
624 if not loop.isCurveLoop():
625 raise TypeError("argument loop needs to be a CurveLoop object.")
626 if len(loop)<3:
627 raise TypeError("the loop must contain at least three Curves.")
628 super(RuledSurface, self).__init__(loop)
629 def __add__(self,other):
630 return RuledSurface(self.getBoundaryLoop()+other)
631 def getGmshCommand(self):
632 return "Ruled Surface(%s) = {%s};"%(self.getID(),self.getBoundaryLoop().getID())
633
634 class SurfaceLoop(Primitive):
635 """
636 a surface loop. It defines the shell of a volume.
637
638 The loop must represent a closed shell, and the L{Surface}s should be oriented consistently.
639 """
640 def __init__(self,*surfaces):
641 """
642 creates a surface loop
643 """
644 super(SurfaceLoop, self).__init__()
645 self.__surfaces=[]
646 self.addSurface(*surfaces)
647 def addSurface(self,*surfaces):
648 for i in range(len(surfaces)):
649 if not surfaces[i].isSurface():
650 raise TypeError("%s-th argument is not a Surface object."%i)
651 self.__surfaces+=surfaces
652
653 def isSurfaceLoop(self):
654 return True
655 def getSurfaces(self):
656 return self.__surfaces
657 def __add__(self,other):
658 return SurfaceLoop([c+other for c in self.getSurfaces])
659 def __len__(self):
660 return len(self.__surfaces)
661 def getPrimitives(self):
662 out=set([self])
663 for i in self.getSurfaces(): out|=i.getPrimitives()
664 return out
665 def getConstructionPoints(self):
666 out=set()
667 for i in self.getSurfaces(): out|=i.getConstructionPoints()
668 return out
669 def getGmshCommand(self):
670 out=""
671 for i in self.getSurfaces():
672 if len(out)>0:
673 out+=", %s"%i.getID()
674 else:
675 out="%s"%i.getID()
676 return "Surface Loop(%s) = {%s};"%(self.getID(),out)
677
678 class Volume(Primitive):
679 """
680 a volume with holes.
681 """
682 def __init__(self,loop,holes=[]):
683 """
684 creates a volume
685
686 @param loop: L{SurfaceLoop} defining the boundary of the surface
687 @param holes: list of L{SurfaceLoop} defining holes in the surface.
688 @note: A SurfaceLoop defining a hole should not have any surfaces in common with the exterior SurfaceLoop.
689 A SurfaceLoop defining a hole should not have any surfaces in common with another SurfaceLoop defining a hole in the same volume.
690 """
691 super(Volume, self).__init__()
692 if not loop.isSurfaceLoop():
693 raise TypeError("argument loop needs to be a SurfaceLoop object.")
694 for i in range(len(holes)):
695 if not holes[i].isSurfaceLoop():
696 raise TypeError("%i th hole needs to be a SurfaceLoop object.")
697 self.__loop=loop
698 self.__holes=holes
699 def getHoles(self):
700 return self.__holes
701 def getSurfaceLoop(self):
702 return self.__loop
703 def __add__(self,other):
704 return Volume(self.getSurfaceLoop()+other, holes=[h+other for h in self.getHoles()])
705 def getPrimitives(self):
706 out=set([self]) | self.getSurfaceLoop().getPrimitives()
707 for i in self.getHoles(): out|=i.getPrimitives()
708 return out
709 def getConstructionPoints(self):
710 out=self.getSurfaceLoop().getConstructionPoints()
711 for i in self.getHoles(): out|=i.Points()
712 return out
713 def getGmshCommand(self):
714 out=""
715 for i in self.getHoles():
716 if len(out)>0:
717 out+=", %s"%i.getID()
718 else:
719 out="%s"%i.getID()
720 if len(out)>0:
721 return "Volume(%s) = {%s, %s};"%(self.getID(),self.getSurfaceLoop().getID(), out)
722 else:
723 return "Volume(%s) = {%s};"%(self.getID(),self.getSurfaceLoop().getID())
724
725 class ReversedPrimitive(object):
726 def __init__(self,prim):
727 self.__prim=prim
728 def __getattr__(self,name):
729 if name == "getID":
730 return self.getReverseID
731 else:
732 return getattr(self.__prim,name)
733 def getReverseID(self):
734 return -self.__prim.getID()
735
736 class PropertySet(Primitive):
737 """
738 defines a group L{Primitive} objects.
739 """
740 def __init__(self,tag=None,*items):
741 super(PropertySet, self).__init__()
742 self.__items=items
743 self.__tag=tag
744 def getPrimitives(self):
745 out=set([self, self.getBoundaryLoop().getPrimitives()])
746 for i in self.getHoles(): out|=i.getPrimitives()
747 return out
748
749 class PrimitiveStack(object):
750 def __init__(self,*items):
751 self.__prims=set()
752 for i in items:
753 self.__prims|=i.getPrimitives()
754 self.__prims=list(self.__prims)
755 self.__prims.sort()
756
757 def getGmshCommands(self):
758 out=""
759 for i in self.__prims:
760 out+=i.getGmshCommand()+"\n"
761 return out

  ViewVC Help
Powered by ViewVC 1.1.26