FreeFOAM The Cross-Platform CFD Toolkit
slidingInterface.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "slidingInterface.H"
28 #include <OpenFOAM/polyMesh.H>
31 #include <OpenFOAM/plane.H>
32 
33 // Index of debug signs:
34 // p - adjusting a projection point
35 // * - adjusting edge intersection
36 
37 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
38 
39 namespace Foam
40 {
41  defineTypeNameAndDebug(slidingInterface, 0);
43  (
44  polyMeshModifier,
45  slidingInterface,
46  dictionary
47  );
48 }
49 
50 
51 template<>
53 {
54  "integral",
55  "partial"
56 };
57 
58 
61 
62 
63 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
64 
65 void Foam::slidingInterface::checkDefinition()
66 {
67  const polyMesh& mesh = topoChanger().mesh();
68 
69  if
70  (
71  !masterFaceZoneID_.active()
72  || !slaveFaceZoneID_.active()
73  || !cutPointZoneID_.active()
74  || !cutFaceZoneID_.active()
75  || !masterPatchID_.active()
76  || !slavePatchID_.active()
77  )
78  {
80  (
81  "void slidingInterface::checkDefinition()"
82  ) << "Not all zones and patches needed in the definition "
83  << "have been found. Please check your mesh definition."
84  << abort(FatalError);
85  }
86 
87  // Check the sizes and set up state
88  if
89  (
90  mesh.faceZones()[masterFaceZoneID_.index()].empty()
91  || mesh.faceZones()[slaveFaceZoneID_.index()].empty()
92  )
93  {
94  FatalErrorIn("void slidingInterface::checkDefinition()")
95  << "Master or slave face zone contain no faces. "
96  << "Please check your mesh definition."
97  << abort(FatalError);
98  }
99 
100  if (debug)
101  {
102  Pout<< "Sliding interface object " << name() << " :" << nl
103  << " master face zone: " << masterFaceZoneID_.index() << nl
104  << " slave face zone: " << slaveFaceZoneID_.index() << endl;
105  }
106 }
107 
108 
109 void Foam::slidingInterface::clearOut() const
110 {
111  clearPointProjection();
112  clearAttachedAddressing();
113  clearAddressing();
114 }
115 
116 
117 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
118 
119 
120 // Construct from components
121 Foam::slidingInterface::slidingInterface
122 (
123  const word& name,
124  const label index,
125  const polyTopoChanger& mme,
126  const word& masterFaceZoneName,
127  const word& slaveFaceZoneName,
128  const word& cutPointZoneName,
129  const word& cutFaceZoneName,
130  const word& masterPatchName,
131  const word& slavePatchName,
132  const typeOfMatch tom,
133  const bool coupleDecouple,
134  const intersection::algorithm algo
135 )
136 :
137  polyMeshModifier(name, index, mme, true),
138  masterFaceZoneID_
139  (
140  masterFaceZoneName,
141  mme.mesh().faceZones()
142  ),
143  slaveFaceZoneID_
144  (
145  slaveFaceZoneName,
146  mme.mesh().faceZones()
147  ),
148  cutPointZoneID_
149  (
150  cutPointZoneName,
151  mme.mesh().pointZones()
152  ),
153  cutFaceZoneID_
154  (
155  cutFaceZoneName,
156  mme.mesh().faceZones()
157  ),
158  masterPatchID_
159  (
160  masterPatchName,
161  mme.mesh().boundaryMesh()
162  ),
163  slavePatchID_
164  (
165  slavePatchName,
166  mme.mesh().boundaryMesh()
167  ),
168  matchType_(tom),
169  coupleDecouple_(coupleDecouple),
170  attached_(false),
171  projectionAlgo_(algo),
172  trigger_(false),
173  pointMergeTol_(pointMergeTolDefault_),
174  edgeMergeTol_(edgeMergeTolDefault_),
175  nFacesPerSlaveEdge_(nFacesPerSlaveEdgeDefault_),
176  edgeFaceEscapeLimit_(edgeFaceEscapeLimitDefault_),
177  integralAdjTol_(integralAdjTolDefault_),
178  edgeMasterCatchFraction_(edgeMasterCatchFractionDefault_),
179  edgeCoPlanarTol_(edgeCoPlanarTolDefault_),
180  edgeEndCutoffTol_(edgeEndCutoffTolDefault_),
181  cutFaceMasterPtr_(NULL),
182  cutFaceSlavePtr_(NULL),
183  masterFaceCellsPtr_(NULL),
184  slaveFaceCellsPtr_(NULL),
185  masterStickOutFacesPtr_(NULL),
186  slaveStickOutFacesPtr_(NULL),
187  retiredPointMapPtr_(NULL),
188  cutPointEdgePairMapPtr_(NULL),
189  slavePointPointHitsPtr_(NULL),
190  slavePointEdgeHitsPtr_(NULL),
191  slavePointFaceHitsPtr_(NULL),
192  masterPointEdgeHitsPtr_(NULL),
193  projectedSlavePointsPtr_(NULL)
194 {
195  checkDefinition();
196 
197  if (attached_)
198  {
200  (
201  "Foam::slidingInterface::slidingInterface\n"
202  "(\n"
203  " const word& name,\n"
204  " const label index,\n"
205  " const polyTopoChanger& mme,\n"
206  " const word& masterFaceZoneName,\n"
207  " const word& slaveFaceZoneName,\n"
208  " const word& cutFaceZoneName,\n"
209  " const word& cutPointZoneName,\n"
210  " const word& masterPatchName,\n"
211  " const word& slavePatchName,\n"
212  " const typeOfMatch tom,\n"
213  " const bool coupleDecouple\n"
214  ")"
215  ) << "Creation of a sliding interface from components "
216  << "in attached state not supported."
217  << abort(FatalError);
218  }
219  else
220  {
221  calcAttachedAddressing();
222  }
223 }
224 
225 
226 // Construct from components
227 Foam::slidingInterface::slidingInterface
228 (
229  const word& name,
230  const dictionary& dict,
231  const label index,
232  const polyTopoChanger& mme
233 )
234 :
235  polyMeshModifier(name, index, mme, Switch(dict.lookup("active"))),
236  masterFaceZoneID_
237  (
238  dict.lookup("masterFaceZoneName"),
239  mme.mesh().faceZones()
240  ),
241  slaveFaceZoneID_
242  (
243  dict.lookup("slaveFaceZoneName"),
244  mme.mesh().faceZones()
245  ),
246  cutPointZoneID_
247  (
248  dict.lookup("cutPointZoneName"),
249  mme.mesh().pointZones()
250  ),
251  cutFaceZoneID_
252  (
253  dict.lookup("cutFaceZoneName"),
254  mme.mesh().faceZones()
255  ),
256  masterPatchID_
257  (
258  dict.lookup("masterPatchName"),
259  mme.mesh().boundaryMesh()
260  ),
261  slavePatchID_
262  (
263  dict.lookup("slavePatchName"),
264  mme.mesh().boundaryMesh()
265  ),
266  matchType_(typeOfMatchNames_.read((dict.lookup("typeOfMatch")))),
267  coupleDecouple_(dict.lookup("coupleDecouple")),
268  attached_(dict.lookup("attached")),
269  projectionAlgo_
270  (
271  intersection::algorithmNames_.read(dict.lookup("projection"))
272  ),
273  trigger_(false),
274  cutFaceMasterPtr_(NULL),
275  cutFaceSlavePtr_(NULL),
276  masterFaceCellsPtr_(NULL),
277  slaveFaceCellsPtr_(NULL),
278  masterStickOutFacesPtr_(NULL),
279  slaveStickOutFacesPtr_(NULL),
280  retiredPointMapPtr_(NULL),
281  cutPointEdgePairMapPtr_(NULL),
282  slavePointPointHitsPtr_(NULL),
283  slavePointEdgeHitsPtr_(NULL),
284  slavePointFaceHitsPtr_(NULL),
285  masterPointEdgeHitsPtr_(NULL),
286  projectedSlavePointsPtr_(NULL)
287 {
288  // Optionally default tolerances from dictionary
289  setTolerances(dict);
290 
291  checkDefinition();
292 
293  // If the interface is attached, the master and slave face zone addressing
294  // needs to be read in; otherwise it will be created
295  if (attached_)
296  {
297  if (debug)
298  {
299  Pout<< "slidingInterface::slidingInterface(...) "
300  << " for object " << name << " : "
301  << "Interface attached. Reading master and slave face zones "
302  << "and retired point lookup." << endl;
303  }
304 
305  // The face zone addressing is written out in the definition dictionary
306  masterFaceCellsPtr_ = new labelList(dict.lookup("masterFaceCells"));
307  slaveFaceCellsPtr_ = new labelList(dict.lookup("slaveFaceCells"));
308 
309  masterStickOutFacesPtr_ =
310  new labelList(dict.lookup("masterStickOutFaces"));
311  slaveStickOutFacesPtr_ =
312  new labelList(dict.lookup("slaveStickOutFaces"));
313 
314  retiredPointMapPtr_ = new Map<label>(dict.lookup("retiredPointMap"));
315  cutPointEdgePairMapPtr_ =
316  new Map<Pair<edge> >(dict.lookup("cutPointEdgePairMap"));
317  }
318  else
319  {
320  calcAttachedAddressing();
321  }
322 }
323 
324 
325 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
326 
328 {
329  clearOut();
330 }
331 
332 
333 void Foam::slidingInterface::clearAddressing() const
334 {
335  deleteDemandDrivenData(cutFaceMasterPtr_);
336  deleteDemandDrivenData(cutFaceSlavePtr_);
337 }
338 
339 
340 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
341 
343 {
344  return masterFaceZoneID_;
345 }
346 
347 
349 {
350  return slaveFaceZoneID_;
351 }
352 
353 
355 {
356  if (coupleDecouple_)
357  {
358  // Always changes. If not attached, project points
359  if (debug)
360  {
361  Pout<< "bool slidingInterface::changeTopology() const "
362  << "for object " << name() << " : "
363  << "Couple-decouple mode." << endl;
364  }
365 
366  if (!attached_)
367  {
368  projectPoints();
369  }
370  else
371  {
372  }
373 
374  return true;
375  }
376 
377  if
378  (
379  attached_
380  && !topoChanger().mesh().changing()
381  )
382  {
383  // If the mesh is not moving or morphing and the interface is
384  // already attached, the topology will not change
385  return false;
386  }
387  else
388  {
389  // Check if the motion changes point projection
390  return projectPoints();
391  }
392 }
393 
394 
396 {
397  if (coupleDecouple_)
398  {
399  if (attached_)
400  {
401  // Attached, coupling
402  decoupleInterface(ref);
403  }
404  else
405  {
406  // Detached, coupling
407  coupleInterface(ref);
408  }
409 
410  return;
411  }
412 
413  if (trigger_)
414  {
415  if (attached_)
416  {
417  // Clear old coupling data
418  clearCouple(ref);
419  }
420 
421  coupleInterface(ref);
422 
423  trigger_ = false;
424  }
425 }
426 
427 
429 {
430  if (debug)
431  {
432  Pout<< "void slidingInterface::modifyMotionPoints("
433  << "pointField& motionPoints) const for object " << name() << " : "
434  << "Adjusting motion points." << endl;
435  }
436 
437  const polyMesh& mesh = topoChanger().mesh();
438 
439  // Get point from the cut zone
440  const labelList& cutPoints = mesh.pointZones()[cutPointZoneID_.index()];
441 
442  if (cutPoints.size() && !projectedSlavePointsPtr_)
443  {
444  return;
445  }
446  else
447  {
448  const pointField& projectedSlavePoints = *projectedSlavePointsPtr_;
449 
450  const Map<label>& rpm = retiredPointMap();
451 
452  const Map<Pair<edge> >& cpepm = cutPointEdgePairMap();
453 
454  const Map<label>& slaveZonePointMap =
455  mesh.faceZones()[slaveFaceZoneID_.index()]().meshPointMap();
456 
457  const primitiveFacePatch& masterPatch =
458  mesh.faceZones()[masterFaceZoneID_.index()]();
459  const edgeList& masterEdges = masterPatch.edges();
460  const pointField& masterLocalPoints = masterPatch.localPoints();
461 
462  const primitiveFacePatch& slavePatch =
463  mesh.faceZones()[slaveFaceZoneID_.index()]();
464  const edgeList& slaveEdges = slavePatch.edges();
465  const pointField& slaveLocalPoints = slavePatch.localPoints();
466  const vectorField& slavePointNormals = slavePatch.pointNormals();
467 
468  forAll (cutPoints, pointI)
469  {
470  // Try to find the cut point in retired points
471  Map<label>::const_iterator rpmIter = rpm.find(cutPoints[pointI]);
472 
473  if (rpmIter != rpm.end())
474  {
475  if (debug)
476  {
477  Pout << "p";
478  }
479 
480  // Cut point is a retired point
481  motionPoints[cutPoints[pointI]] =
482  projectedSlavePoints[slaveZonePointMap.find(rpmIter())()];
483  }
484  else
485  {
486  // A cut point is not a projected slave point. Therefore, it
487  // must be an edge-to-edge intersection.
488 
489  Map<Pair<edge> >::const_iterator cpepmIter =
490  cpepm.find(cutPoints[pointI]);
491 
492  if (cpepmIter != cpepm.end())
493  {
494 // Pout << "Need to re-create hit for point " << cutPoints[pointI] << " lookup: " << cpepmIter() << endl;
495 
496  // Note.
497  // The edge cutting code is repeated in
498  // slidingInterface::coupleInterface. This is done for
499  // efficiency reasons and avoids multiple creation of
500  // cutting planes. Please update both simultaneously.
501  //
502  const edge& globalMasterEdge = cpepmIter().first();
503 
504  const label curMasterEdgeIndex =
505  masterPatch.whichEdge
506  (
507  edge
508  (
509  masterPatch.whichPoint
510  (
511  globalMasterEdge.start()
512  ),
513  masterPatch.whichPoint
514  (
515  globalMasterEdge.end()
516  )
517  )
518  );
519 
520  const edge& cme = masterEdges[curMasterEdgeIndex];
521 // Pout << "curMasterEdgeIndex: " << curMasterEdgeIndex << " cme: " << cme << endl;
522  const edge& globalSlaveEdge = cpepmIter().second();
523 
524  const label curSlaveEdgeIndex =
525  slavePatch.whichEdge
526  (
527  edge
528  (
529  slavePatch.whichPoint
530  (
531  globalSlaveEdge.start()
532  ),
533  slavePatch.whichPoint
534  (
535  globalSlaveEdge.end()
536  )
537  )
538  );
539 
540  const edge& curSlaveEdge = slaveEdges[curSlaveEdgeIndex];
541 // Pout << "curSlaveEdgeIndex: " << curSlaveEdgeIndex << " curSlaveEdge: " << curSlaveEdge << endl;
542  const point& a = projectedSlavePoints[curSlaveEdge.start()];
543  const point& b = projectedSlavePoints[curSlaveEdge.end()];
544 
545  point c =
546  0.5*
547  (
548  slaveLocalPoints[curSlaveEdge.start()]
549  + slavePointNormals[curSlaveEdge.start()]
550  + slaveLocalPoints[curSlaveEdge.end()]
551  + slavePointNormals[curSlaveEdge.end()]
552  );
553 
554  // Create the plane
555  plane cutPlane(a, b, c);
556 
557  linePointRef curSlaveLine =
558  curSlaveEdge.line(slaveLocalPoints);
559  const scalar curSlaveLineMag = curSlaveLine.mag();
560 
561  scalar cutOnMaster =
562  cutPlane.lineIntersect
563  (
564  cme.line(masterLocalPoints)
565  );
566 
567  if
568  (
569  cutOnMaster > edgeEndCutoffTol_
570  && cutOnMaster < 1.0 - edgeEndCutoffTol_
571  )
572  {
573  // Master is cut, check the slave
574  point masterCutPoint =
575  masterLocalPoints[cme.start()]
576  + cutOnMaster*cme.vec(masterLocalPoints);
577 
578  pointHit slaveCut =
579  curSlaveLine.nearestDist(masterCutPoint);
580 
581  if (slaveCut.hit())
582  {
583  // Strict checking of slave cut to avoid capturing
584  // end points.
585  scalar cutOnSlave =
586  (
587  (
588  slaveCut.hitPoint()
589  - curSlaveLine.start()
590  ) & curSlaveLine.vec()
591  )/sqr(curSlaveLineMag);
592 
593  // Calculate merge tolerance from the
594  // target edge length
595  scalar mergeTol =
596  edgeCoPlanarTol_*mag(b - a);
597 
598  if
599  (
600  cutOnSlave > edgeEndCutoffTol_
601  && cutOnSlave < 1.0 - edgeEndCutoffTol_
602  && slaveCut.distance() < mergeTol
603  )
604  {
605  // Cut both master and slave.
606  motionPoints[cutPoints[pointI]] =
607  masterCutPoint;
608  }
609  }
610  else
611  {
612  Pout<< "Missed slave edge!!! This is an error. "
613  << "Master edge: "
614  << cme.line(masterLocalPoints)
615  << " slave edge: " << curSlaveLine
616  << " point: " << masterCutPoint
617  << " weight: " <<
618  (
619  (
620  slaveCut.missPoint()
621  - curSlaveLine.start()
622  ) & curSlaveLine.vec()
623  )/sqr(curSlaveLineMag)
624  << endl;
625  }
626  }
627  else
628  {
629  Pout<< "Missed master edge!!! This is an error"
630  << endl;
631  }
632  }
633  else
634  {
636  (
637  "void slidingInterface::modifyMotionPoints"
638  "(pointField&) const"
639  ) << "Cut point " << cutPoints[pointI]
640  << " not recognised as either the projected "
641  << "or as intersection point. Error in point "
642  << "projection or data mapping"
643  << abort(FatalError);
644  }
645  }
646  }
647  if (debug)
648  {
649  Pout << endl;
650  }
651  }
652 }
653 
654 
656 {
657  if (debug)
658  {
659  Pout<< "void slidingInterface::updateMesh(const mapPolyMesh& m)"
660  << " const for object " << name() << " : "
661  << "Updating topology." << endl;
662  }
663 
664  // Mesh has changed topologically. Update local topological data
665  const polyMesh& mesh = topoChanger().mesh();
666 
667  masterFaceZoneID_.update(mesh.faceZones());
668  slaveFaceZoneID_.update(mesh.faceZones());
669  cutPointZoneID_.update(mesh.pointZones());
670  cutFaceZoneID_.update(mesh.faceZones());
671 
672  masterPatchID_.update(mesh.boundaryMesh());
673  slavePatchID_.update(mesh.boundaryMesh());
674 
675 //MJ:Disabled updating
676 // if (!attached())
677 // {
678 // calcAttachedAddressing();
679 // }
680 // else
681 // {
682 // renumberAttachedAddressing(m);
683 // }
684 }
685 
686 
688 {
689  if (!projectedSlavePointsPtr_)
690  {
691  projectPoints();
692  }
693 
694  return *projectedSlavePointsPtr_;
695 }
696 
698 {
699  pointMergeTol_ = dict.lookupOrDefault<scalar>
700  (
701  "pointMergeTol",
702  pointMergeTol_
703  );
704  edgeMergeTol_ = dict.lookupOrDefault<scalar>
705  (
706  "edgeMergeTol",
707  edgeMergeTol_
708  );
709  nFacesPerSlaveEdge_ = dict.lookupOrDefault<label>
710  (
711  "nFacesPerSlaveEdge",
712  nFacesPerSlaveEdge_
713  );
714  edgeFaceEscapeLimit_ = dict.lookupOrDefault<label>
715  (
716  "edgeFaceEscapeLimit",
717  edgeFaceEscapeLimit_
718  );
719  integralAdjTol_ = dict.lookupOrDefault<scalar>
720  (
721  "integralAdjTol",
722  integralAdjTol_
723  );
724  edgeMasterCatchFraction_ = dict.lookupOrDefault<scalar>
725  (
726  "edgeMasterCatchFraction",
727  edgeMasterCatchFraction_
728  );
729  edgeCoPlanarTol_ = dict.lookupOrDefault<scalar>
730  (
731  "edgeCoPlanarTol",
732  edgeCoPlanarTol_
733  );
734  edgeEndCutoffTol_ = dict.lookupOrDefault<scalar>
735  (
736  "edgeEndCutoffTol",
737  edgeEndCutoffTol_
738  );
739 
740  if (report)
741  {
742  Info<< "Sliding interface parameters:" << nl
743  << "pointMergeTol : " << pointMergeTol_ << nl
744  << "edgeMergeTol : " << edgeMergeTol_ << nl
745  << "nFacesPerSlaveEdge : " << nFacesPerSlaveEdge_ << nl
746  << "edgeFaceEscapeLimit : " << edgeFaceEscapeLimit_ << nl
747  << "integralAdjTol : " << integralAdjTol_ << nl
748  << "edgeMasterCatchFraction : " << edgeMasterCatchFraction_ << nl
749  << "edgeCoPlanarTol : " << edgeCoPlanarTol_ << nl
750  << "edgeEndCutoffTol : " << edgeEndCutoffTol_ << endl;
751  }
752 }
753 
754 
756 {
757  os << nl << type() << nl
758  << name()<< nl
759  << masterFaceZoneID_.name() << nl
760  << slaveFaceZoneID_.name() << nl
761  << cutPointZoneID_.name() << nl
762  << cutFaceZoneID_.name() << nl
763  << masterPatchID_.name() << nl
764  << slavePatchID_.name() << nl
765  << typeOfMatchNames_[matchType_] << nl
766  << coupleDecouple_ << nl
767  << attached_ << endl;
768 }
769 
770 
771 // To write out all those tolerances
772 #define WRITE_NON_DEFAULT(name) \
773  if( name ## _ != name ## Default_ )\
774  { \
775  os << " " #name " " << name ## _ << token::END_STATEMENT << nl; \
776  }
777 
778 
780 {
781  os << nl << name() << nl << token::BEGIN_BLOCK << nl
782  << " type " << type() << token::END_STATEMENT << nl
783  << " masterFaceZoneName " << masterFaceZoneID_.name()
785  << " slaveFaceZoneName " << slaveFaceZoneID_.name()
787  << " cutPointZoneName " << cutPointZoneID_.name()
789  << " cutFaceZoneName " << cutFaceZoneID_.name()
791  << " masterPatchName " << masterPatchID_.name()
793  << " slavePatchName " << slavePatchID_.name()
795  << " typeOfMatch " << typeOfMatchNames_[matchType_]
797  << " coupleDecouple " << coupleDecouple_
799  << " projection " << intersection::algorithmNames_[projectionAlgo_]
801  << " attached " << attached_
803  << " active " << active()
804  << token::END_STATEMENT << nl;
805 
806  if (attached_)
807  {
808  masterFaceCellsPtr_->writeEntry("masterFaceCells", os);
809  slaveFaceCellsPtr_->writeEntry("slaveFaceCells", os);
810  masterStickOutFacesPtr_->writeEntry("masterStickOutFaces", os);
811  slaveStickOutFacesPtr_->writeEntry("slaveStickOutFaces", os);
812 
813  os << " retiredPointMap " << retiredPointMap()
814  << token::END_STATEMENT << nl
815  << " cutPointEdgePairMap " << cutPointEdgePairMap()
816  << token::END_STATEMENT << nl;
817  }
818 
819  WRITE_NON_DEFAULT(pointMergeTol)
820  WRITE_NON_DEFAULT(edgeMergeTol)
821  WRITE_NON_DEFAULT(nFacesPerSlaveEdge)
822  WRITE_NON_DEFAULT(edgeFaceEscapeLimit)
823  WRITE_NON_DEFAULT(integralAdjTol)
824  WRITE_NON_DEFAULT(edgeMasterCatchFraction)
825  WRITE_NON_DEFAULT(edgeCoPlanarTol)
826  WRITE_NON_DEFAULT(edgeEndCutoffTol)
827 
828  os << token::END_BLOCK << endl;
829 }
830 
831 
832 // ************************ vim: set sw=4 sts=4 et: ************************ //