rsvs3D  0.0.0
Codes for the c++ implementation of the 3D RSVS
snakstruct_test.cpp
1 #include <cmath>
2 #include <cstdlib>
3 #include <ctime>
4 #include <fstream>
5 #include <iostream>
6 #include <sstream>
7 
8 #include "RSVSalgorithm.hpp"
9 #include "RSVScalc.hpp"
10 #include "RSVSintegration.hpp"
11 #include "RSVSmath.hpp"
12 #include "main.hpp"
13 #include "matrixtools.hpp"
14 #include "meshrefinement.hpp"
15 #include "postprocessing.hpp"
16 #include "snake.hpp"
17 #include "snakeengine.hpp"
18 #include "tetgenrsvs.hpp"
19 #include "triangulate.hpp"
20 #include "warning.hpp"
21 
22 using namespace std;
23 
24 // Implementation in snakstruct.cpp
25 
26 void Test_mathRSVS_FD(snake &testSnake, triangulation &RSVStri, vector<double> &dt, vector<int> &isImpact,
27  RSVScalc &calcObj, tecplotfile &outSnake2, double totT);
28 
29 void PrintRSVSSnake2D(tecplotfile &outSnake, snake &testSnake, double totT, triangulation &testTriangle, mesh &triMesh,
30  triangulation &triRSVS, mesh &voluMesh, int nVoluZone, int ii)
31 {
32  vector<int> vertList;
33  int jj;
34  if (testSnake.snaxs.size() > 0)
35  {
36  // testSnake.snakeconn.TightenConnectivity();
37  outSnake.PrintMesh(testSnake.snakeconn, 1, totT);
38 
39  testSnake.snakeconn.PrepareForUse();
40  testTriangle.stattri.clear();
41  testTriangle.trivert.clear();
42  testTriangle.PrepareForUse();
43  triRSVS.PrepareForUse();
44  TriangulateMesh(testSnake.snakeconn, testTriangle);
45  MeshTriangulation(triMesh, testSnake.snakeconn, testTriangle.stattri, testTriangle.trivert);
46  outSnake.PrintMesh(triMesh, 2, totT);
47  MeshTriangulation(triMesh, testSnake.snakeconn, triRSVS.dynatri, triRSVS.trivert);
48  outSnake.PrintMesh(triMesh, 3, totT);
49  outSnake.PrintTriangulation(triRSVS, &triangulation::dynatri, 4, totT);
50  if (ii == 0)
51  {
52  outSnake.PrintTriangulation(triRSVS, &triangulation::dynatri, 5, totT, rsvs3d::constants::tecplot::line);
53  outSnake.PrintTriangulation(triRSVS, &triangulation::dynatri, 6, totT, rsvs3d::constants::tecplot::line);
54  outSnake.PrintTriangulation(triRSVS, &triangulation::dynatri, 7, totT, rsvs3d::constants::tecplot::line);
55  }
56  outSnake.PrintTriangulation(triRSVS, &triangulation::intertri, 5, totT, rsvs3d::constants::tecplot::line);
57  outSnake.PrintTriangulation(triRSVS, &triangulation::trisurf, 6, totT, rsvs3d::constants::tecplot::line);
58  if (int(triRSVS.acttri.size()) > 0)
59  {
60  outSnake.PrintTriangulation(triRSVS, &triangulation::stattri, 7, totT, rsvs3d::constants::tecplot::line,
61  triRSVS.acttri);
62  }
63 
64  vertList.clear();
65  for (jj = 0; jj < int(testSnake.isMeshVertIn.size()); ++jj)
66  {
67  if (testSnake.isMeshVertIn[jj])
68  {
69  vertList.push_back(testSnake.snakemesh()->verts(jj)->index);
70  }
71  }
72  if (int(testSnake.isMeshVertIn.size()) == 0)
73  {
74  vertList.push_back(testSnake.snakemesh()->verts(0)->index);
75  }
76  outSnake.PrintMesh(*(testSnake.snakemesh()), 8, totT, rsvs3d::constants::tecplot::point, vertList);
77  outSnake.PrintVolumeDat(voluMesh, nVoluZone, 9, totT);
78  outSnake.PrintSnake(testSnake, 10, totT);
79  }
80 }
81 
82 void PrintRSVSSnake(tecplotfile &outSnake, snake &testSnake, double totT, triangulation &testTriangle, mesh &triMesh,
83  triangulation &triRSVS, mesh &voluMesh, int nVoluZone, int ii)
84 {
85  vector<int> vertList;
86  int jj;
87  if (testSnake.snaxs.size() > 0)
88  {
89  // testSnake.snakeconn.TightenConnectivity();
90  outSnake.PrintMesh(testSnake.snakeconn, 1, totT);
91 
92  testSnake.snakeconn.PrepareForUse();
93  testTriangle.stattri.clear();
94  testTriangle.trivert.clear();
95  testTriangle.PrepareForUse();
96  TriangulateMesh(testSnake.snakeconn, testTriangle);
97  MeshTriangulation(triMesh, testSnake.snakeconn, testTriangle.stattri, testTriangle.trivert);
98  outSnake.PrintMesh(triMesh, 2, totT);
99  MeshTriangulation(triMesh, testSnake.snakeconn, triRSVS.dynatri, triRSVS.trivert);
100  outSnake.PrintMesh(triMesh, 3, totT);
101  outSnake.PrintTriangulation(triRSVS, &triangulation::dynatri, 4, totT);
102  if (ii == 0)
103  {
104  outSnake.PrintTriangulation(triRSVS, &triangulation::dynatri, 5, totT, rsvs3d::constants::tecplot::line);
105  outSnake.PrintTriangulation(triRSVS, &triangulation::dynatri, 6, totT, rsvs3d::constants::tecplot::line);
106  outSnake.PrintTriangulation(triRSVS, &triangulation::dynatri, 7, totT, rsvs3d::constants::tecplot::line);
107  }
108  outSnake.PrintTriangulation(triRSVS, &triangulation::intertri, 5, totT, rsvs3d::constants::tecplot::line);
109  outSnake.PrintTriangulation(triRSVS, &triangulation::trisurf, 6, totT, rsvs3d::constants::tecplot::line);
110  if (int(triRSVS.acttri.size()) > 0)
111  {
112  outSnake.PrintTriangulation(triRSVS, &triangulation::stattri, 7, totT, rsvs3d::constants::tecplot::line,
113  triRSVS.acttri);
114  }
115 
116  vertList.clear();
117  for (jj = 0; jj < int(testSnake.isMeshVertIn.size()); ++jj)
118  {
119  if (testSnake.isMeshVertIn[jj])
120  {
121  vertList.push_back(testSnake.snakemesh()->verts(jj)->index);
122  }
123  }
124  if (int(testSnake.isMeshVertIn.size()) == 0)
125  {
126  vertList.push_back(testSnake.snakemesh()->verts(0)->index);
127  }
128  outSnake.PrintMesh(*(testSnake.snakemesh()), 8, totT, rsvs3d::constants::tecplot::point, vertList);
129  outSnake.PrintVolumeDat(voluMesh, nVoluZone, 9, totT);
130  outSnake.PrintSnake(testSnake, 10, totT);
131  }
132 }
133 
134 void PrepareMultiLvlSnake(mesh &snakeMesh, mesh &voluMesh, snake &testSnake, vector<int> &dims, triangulation &triRSVS)
135 {
136  vector<int> elmMapping;
137  RSVScalc calcVolus;
138  int ii;
139 
140  snakeMesh.PrepareForUse();
141  snakeMesh.OrientFaces();
143  testSnake.SetSnakeMesh(&snakeMesh);
144  // testSnake.disp();
145  if (snakeMesh.WhatDim() == 3)
146  {
147  for (ii = 0; ii < snakeMesh.volus.size(); ++ii)
148  {
149  elmMapping.push_back(1);
150  }
151  }
152  else if (snakeMesh.WhatDim() == 2)
153  {
154  for (ii = 0; ii < snakeMesh.surfs.size(); ++ii)
155  {
156  elmMapping.push_back(1);
157  }
158  }
159  else
160  {
161  RSVS3D_ERROR_ARGUMENT("Incorrect dimension");
162  }
163  CartesianMapping(snakeMesh, elmMapping, dims);
164  CoarsenMesh(snakeMesh, voluMesh, elmMapping);
165  snakeMesh.AddParent(&voluMesh, elmMapping);
166 
167  sort(elmMapping);
168  unique(elmMapping);
169  DisplayVector(elmMapping);
170  for (ii = 0; ii < voluMesh.volus.size(); ++ii)
171  {
172  voluMesh.volus[ii].target = (double(rand() % 1001) / 1000.0);
173  }
174  voluMesh.PrepareForUse();
175  voluMesh.OrientFaces();
176 
177  calcVolus.CalculateMesh(voluMesh);
178  calcVolus.ReturnConstrToMesh(voluMesh, &volu::volume);
179 
180  triRSVS.stattri.clear();
181  triRSVS.trivert.clear();
182  triRSVS.PrepareForUse();
183  TriangulateMesh(snakeMesh, triRSVS);
184 
185  testSnake.PrepareForUse();
186 }
187 
188 void PrepareMultiLvlSnakeNoVoluGen(mesh &snakeMesh, mesh &voluMesh, snake &testSnake, triangulation &triRSVS)
189 {
190 
191  RSVScalc calcVolus;
192 
193  snakeMesh.PrepareForUse();
194  snakeMesh.OrientFaces();
196  testSnake.SetSnakeMesh(&snakeMesh);
197  // testSnake.disp();
198 
199  voluMesh.PrepareForUse();
200  voluMesh.OrientFaces();
201 
202  calcVolus.CalculateMesh(voluMesh);
203  calcVolus.ReturnConstrToMesh(voluMesh, &volu::volume);
204 
205  triRSVS.stattri.clear();
206  triRSVS.trivert.clear();
207  triRSVS.PrepareForUse();
208  TriangulateMesh(snakeMesh, triRSVS);
209 
210  testSnake.PrepareForUse();
211 }
212 
213 void Test_randvelstep(snake &testSnake, vector<double> dt, vector<int> &isImpact)
214 {
215  CalculateSnakeVelRand(testSnake);
216  testSnake.CalculateTimeStep(dt, 0.25);
217  testSnake.UpdateDistance(dt);
218  testSnake.PrepareForUse();
219  testSnake.UpdateCoord();
220  testSnake.PrepareForUse();
221  Test_stepalgo(testSnake, isImpact);
222 }
223 
224 void Test_randvelstep_mc(snake &testSnake, vector<double> dt, vector<int> &isImpact)
225 {
226  CalculateSnakeVelRand(testSnake);
227  testSnake.CalculateTimeStep(dt, 0.25);
228  testSnake.UpdateDistance(dt);
229  testSnake.PrepareForUse();
230  testSnake.UpdateCoord();
231  testSnake.PrepareForUse();
232  Test_stepalgo_mergeclean(testSnake, isImpact);
233 }
234 
235 void Test_unitvelstep(snake &testSnake, vector<double> dt, vector<int> &isImpact, bool reflectVel = true)
236 {
237  if (reflectVel)
238  {
239  CalculateSnakeVelUnitReflect(testSnake);
240  }
241  else
242  {
243  CalculateSnakeVelUnit(testSnake);
244  }
245  testSnake.CalculateTimeStep(dt, 0.25);
246  testSnake.UpdateDistance(dt);
247  testSnake.PrepareForUse();
248  testSnake.UpdateCoord();
249  testSnake.PrepareForUse();
250  Test_stepalgo(testSnake, isImpact);
251 }
252 
253 // -------------------------------------------------------------------------------------------
254 // TEST CODE
255 // -------------------------------------------------------------------------------------------
256 
257 int Test_SnakeStructures()
258 {
259  // Test the functionality provided by arraystructures
260 
261  int errFlag, errTest;
262 
263  errFlag = 0;
264 
265  cout << "--------------------------------------------" << endl;
266  cout << " testing coordvec" << endl;
267  cout << "--------------------------------------------" << endl;
268  errTest = Test_coordvec();
269  errFlag = errFlag | (errTest != 0);
270 
271  cout << "--------------------------------------------" << endl;
272  cout << " testing snax" << endl;
273  cout << "--------------------------------------------" << endl;
274  errTest = Test_snax();
275  errFlag = errFlag | (errTest != 0);
276 
277  cout << "--------------------------------------------" << endl;
278  cout << " testing snaxedge" << endl;
279  cout << "--------------------------------------------" << endl;
280  errTest = Test_snaxedge();
281  errFlag = errFlag | (errTest != 0);
282 
283  cout << "--------------------------------------------" << endl;
284  cout << " testing Snake" << endl;
285  cout << "--------------------------------------------" << endl;
286  errTest = Test_snake();
287  errFlag = errFlag | (errTest != 0);
288 
289  return (errFlag);
290 }
291 
292 int Test_coordvec()
293 {
294 
295  coordvec testCoord, unitCoord;
296  try
297  {
298  testCoord.assign(1.0, 2.0, 3.0);
299 
300  cout << "base vector: ";
301  testCoord.disp();
302 
303  unitCoord = testCoord.Unit();
304  cout << "unit vector: ";
305  unitCoord.disp();
306 
307  cout << "unit access: ";
308  cout << "coord vec [" << testCoord.Unit(0) << "," << testCoord.Unit(1) << "," << testCoord.Unit(2) << "] norm 1"
309  << endl;
310 
311  cout << "base oper(): ";
312  cout << "coord vec [" << testCoord(0) << "," << testCoord(1) << "," << testCoord(2) << "] " << endl;
313  cout << "base oper(): ";
314  testCoord.disp();
315  cout << "base oper[]: ";
316  cout << "coord vec [" << testCoord[0] << "," << testCoord[1] << "," << testCoord[2] << "] " << endl;
317  cout << "base oper[]: ";
318  testCoord.disp();
319 
320  cout << "base ope()=: {compile error}";
321  // testCoord(0)=0;
322  testCoord.disp();
323  cout << "base ope[]=: ";
324  testCoord[0] = 0;
325  testCoord.disp();
326  testCoord.GetNorm();
327  cout << "base getnor: ";
328  testCoord.disp();
329  }
330  catch (exception const &ex)
331  {
332  cerr << "Exception: " << ex.what() << endl;
333  return -1;
334  }
335  return (0);
336 }
337 
338 int Test_snax()
339 {
340  int i;
341 
342  snaxarray snaxStack, snaxStack2; // stack of ints
343  snax singleSnax;
344 
345  try
346  {
347  singleSnax.disp();
348  snaxStack.assign(5, singleSnax);
349  snaxStack.disp();
350  snaxStack.PopulateIndices();
351  cout << endl;
352  snaxStack[0].index = 10;
353  snaxStack.disp();
354 
355  snaxStack.PrepareForUse();
356  i = snaxStack.find(10);
357  cout << "found position " << i << " Index " << snaxStack(i)->index << endl;
358 
359  snaxStack2 = snaxStack;
360  snaxStack2.disp();
361  snaxStack[0] = singleSnax;
362  snaxStack.disp();
363 
364  cout << "Are the Same " << CompareDisp(snaxStack, snaxStack2) << endl;
365  }
366  catch (exception const &ex)
367  {
368  cerr << "Exception: " << ex.what() << endl;
369  return -1;
370  }
371  return (0);
372 }
373 
374 int Test_snaxedge()
375 {
376  snaxedgearray snaxStack, snaxStack2; // stack of ints
377  snaxedge singleSnax;
378  try
379  {
380  singleSnax.normvector[1] = 2;
381  snaxStack.Init(4);
382  snaxStack.PopulateIndices();
383  snaxStack[1] = singleSnax;
384 
385  snaxStack.PrepareForUse();
386 
387  snaxStack.disp();
388  }
389  catch (exception const &ex)
390  {
391  cerr << "Exception: " << ex.what() << endl;
392  return -1;
393  }
394  return (0);
395 }
396 
397 int Test_snake()
398 {
399  snake testSnake, testSnake2, testSnake3;
400  mesh snakeMesh, snakeMesh2;
401 
402  bool errFlag;
403  int errTest = 0;
404  // int nVert,nSurf,nEdge,nVolu;
405  try
406  {
407  snakeMesh.Init(8, 12, 6, 1);
408  snakeMesh2.Init(8, 12, 6, 1);
409 
410  testSnake.Init(&snakeMesh, 8, 12, 6, 1);
411  testSnake2 = testSnake;
412 
413  errFlag = CompareDisp(testSnake, testSnake2);
414  cout << "Compare snakes after = assignement: 1=" << errFlag << endl;
415  errTest = errTest + int(!errFlag);
416  testSnake.ChangeIndices(1, 2, 3, 4);
417  cout << "Succesfully changed indices (ChangeIndices)" << endl;
418  testSnake.ChangeIndicesSnakeMesh(5, 6, 7, 8);
419  cout << "Succesfully changed indices (ChangeIndicesSnakeMesh)" << endl;
420 
421  testSnake.PrepareForUse();
422 
423  cout << "-----------------------testSnake1-----------------------" << endl;
424  // testSnake.disp();
425  testSnake.displight();
426 
427  testSnake3 = testSnake.MakeCompatible(testSnake2);
428  cout << "-----------------------testSnake2-----------------------" << endl;
429  // testSnake2.disp();
430  cout << "-----------------------testSnake3-----------------------" << endl;
431  // testSnake3.disp();
432  testSnake.Concatenate(testSnake3);
433  testSnake.PrepareForUse();
434 
435  testSnake3.Init(&snakeMesh2, 8, 12, 6, 1);
436  testSnake.MakeCompatible_inplace(testSnake3);
437 
438  try
439  {
440  testSnake.Concatenate(testSnake3);
441  cerr << "Error : Concatenation between snakes on different meshes was allowed" << endl;
442  }
443  catch (exception const &ex)
444  {
445  cout << "Succesfully generated failure" << endl;
446  }
447  cout << "-----------------------testSnake-----------------------" << endl;
448  testSnake.displight();
449  }
450  catch (exception const &ex)
451  {
452  cerr << "Exception: " << ex.what() << endl;
453  return -1;
454  }
455  return (errTest);
456 }
457 
458 int Test_snakeinit_velselect(int velCase = 0, int nStep = 300)
459 {
460  snake testSnake;
461  mesh snakeMesh, triMesh;
462  triangulation testTriangle;
463  std::string fileToOpen;
464  tecplotfile outSnake;
465  double totT = 0.0;
466  vector<double> dt;
467  vector<int> isImpact;
468  int start_s, stop_s, ii;
469  // bool errFlag;
470  int errTest = 0;
471 
472  try
473  {
474  fileToOpen = "../TESTOUT/TestSnake";
475  switch (velCase)
476  {
477  case 0: // random
478  fileToOpen += "_rand";
479  break;
480  case 1: // unit with reflection
481  fileToOpen += "_unitref";
482  break;
483  case 2: // unit step with no reflection
484  fileToOpen += "_unitnoref";
485  break;
486  default:
487  RSVS3D_ERROR_ARGUMENT("velCase option not known");
488  break;
489  }
490  fileToOpen += std::string("_") + std::to_string(nStep) + ".plt";
491 
492  outSnake.OpenFile(fileToOpen.c_str());
493  errTest += snakeMesh.read("../TESTOUT/mesh203010.dat");
494  snakeMesh.PrepareForUse();
495  testSnake.SetSnakeMesh(&snakeMesh);
496  outSnake.PrintMesh(*(testSnake.snakemesh()));
497 
498  snakeMesh.OrientFaces();
499  start_s = clock();
500  testSnake.PrepareForUse();
501 
502  SpawnAtVertex(testSnake, 1022);
503  SpawnAtVertex(testSnake, 674);
504  SpawnAtVertex(testSnake, 675);
505  SpawnAtVertex(testSnake, 728);
506  SpawnAtVertex(testSnake, 729);
507  SpawnAtVertex(testSnake, 731);
508  testSnake.displight();
509  stop_s = clock();
510  cout << "time: " << (stop_s - start_s) / double(CLOCKS_PER_SEC) * 1000 << "ms" << endl;
511 
512  start_s = clock();
513  testSnake.PrepareForUse();
514  for (ii = 0; ii < nStep; ++ii)
515  {
516  cout << ii << " ";
517 
518  if (testSnake.snaxs.size() > 0)
519  {
520  // testSnake.snakeconn.TightenConnectivity();
521  outSnake.PrintMesh(testSnake.snakeconn, 1, totT);
522 
523  testSnake.snakeconn.PrepareForUse();
524  testTriangle.stattri.clear();
525  testTriangle.trivert.clear();
526  testTriangle.PrepareForUse();
527  TriangulateMesh(testSnake.snakeconn, testTriangle);
528  testTriangle.CalcTriVertPos();
529  MeshTriangulation(triMesh, testSnake.snakeconn, testTriangle.stattri, testTriangle.trivert);
530  outSnake.PrintMesh(triMesh, 2, totT);
531  }
532  switch (velCase)
533  {
534  case 0: // random
535  Test_randvelstep(testSnake, dt, isImpact);
536  break;
537  case 1: // unit with reflection
538  Test_unitvelstep(testSnake, dt, isImpact, true);
539  break;
540  case 2: // unit step with no reflection
541  Test_unitvelstep(testSnake, dt, isImpact, false);
542  break;
543  default:
544  RSVS3D_ERROR_ARGUMENT("velCase option not known");
545  break;
546  }
547  cout << endl;
548 
549  totT = totT + 1;
550  }
551 
552  if (testSnake.snaxs.size() > 0)
553  {
554  outSnake.PrintMesh(testSnake.snakeconn, 1, totT);
555  testTriangle.stattri.clear();
556  testTriangle.trivert.clear();
557  testTriangle.PrepareForUse();
558  TriangulateMesh(testSnake.snakeconn, testTriangle);
559  testTriangle.CalcTriVertPos();
560 
561  MeshTriangulation(triMesh, testSnake.snakeconn, testTriangle.stattri, testTriangle.trivert);
562  outSnake.PrintMesh(triMesh, 2, totT);
563  }
564  stop_s = clock();
565  cout << "time: " << (stop_s - start_s) / double(CLOCKS_PER_SEC) * 1000 << "ms" << endl;
566  testSnake.displight();
567  // the code you wish to time goes here
568 
569  // testSnake.disp();
570  }
571  catch (exception const &ex)
572  {
573  cerr << "Exception: " << ex.what() << endl;
574  return -1;
575  }
576  return (errTest);
577 }
578 
579 int Test_snakeinit_random()
580 {
581  // rand
582  return (Test_snakeinit_velselect(0));
583 }
584 int Test_snakeinit_unit()
585 {
586  // unit reflect
587  return (Test_snakeinit_velselect(1));
588 }
589 int Test_snakeinit_unitnoreflect()
590 {
591  // unit no reflect
592  return (Test_snakeinit_velselect(2));
593 }
594 
595 int Test_snakeinit_random_short()
596 {
597  // rand
598  return (Test_snakeinit_velselect(0, 20));
599 }
600 int Test_snakeinit_unit_short()
601 {
602  // unit reflect
603  return (Test_snakeinit_velselect(1, 20));
604 }
605 int Test_snakeinit_unitnoreflect_short()
606 {
607  // unit no reflect
608  return (Test_snakeinit_velselect(2, 20));
609 }
610 
611 int Test_snakeinit_MC()
612 {
613  snake testSnake;
614  mesh snakeMesh, triMesh;
615  triangulation testTriangle;
616  const char *fileToOpen;
617  tecplotfile outSnake;
618  double totT = 0.0;
619  vector<double> dt;
620  vector<int> isImpact;
621  int start_s, stop_s, ii;
622  // bool errFlag;
623  int errTest = 0;
624 
625  try
626  {
627  fileToOpen = "../TESTOUT/TestSnake.plt";
628 
629  outSnake.OpenFile(fileToOpen);
630  errTest += snakeMesh.read("../TESTOUT/mesh203010.dat");
631  snakeMesh.PrepareForUse();
632  testSnake.SetSnakeMesh(&snakeMesh);
633  outSnake.PrintMesh(*(testSnake.snakemesh()));
634 
635  snakeMesh.OrientFaces();
636  start_s = clock();
637  testSnake.PrepareForUse();
638 
639  SpawnAtVertex(testSnake, 1022);
640  SpawnAtVertex(testSnake, 674);
641  SpawnAtVertex(testSnake, 675);
642  SpawnAtVertex(testSnake, 728);
643  SpawnAtVertex(testSnake, 729);
644  SpawnAtVertex(testSnake, 731);
645  testSnake.displight();
646  stop_s = clock();
647  cout << "time: " << (stop_s - start_s) / double(CLOCKS_PER_SEC) * 1000 << "ms" << endl;
648 
649  start_s = clock();
650  testSnake.PrepareForUse();
651  for (ii = 0; ii < 50; ++ii)
652  {
653  cout << ii << " ";
654 
655  if (testSnake.snaxs.size() > 0)
656  {
657  // testSnake.snakeconn.TightenConnectivity();
658  outSnake.PrintMesh(testSnake.snakeconn, 1, totT);
659 
660  testSnake.snakeconn.PrepareForUse();
661  testTriangle.stattri.clear();
662  testTriangle.trivert.clear();
663  testTriangle.PrepareForUse();
664  TriangulateMesh(testSnake.snakeconn, testTriangle);
665  testTriangle.CalcTriVertPos();
666  MeshTriangulation(triMesh, testSnake.snakeconn, testTriangle.stattri, testTriangle.trivert);
667  outSnake.PrintMesh(triMesh, 2, totT);
668  }
669  if (ii == 30)
670  {
671  cout << "break here" << endl;
672  }
673  Test_randvelstep_mc(testSnake, dt, isImpact);
674  cout << endl;
675 
676  totT = totT + 1;
677  }
678 
679  if (testSnake.snaxs.size() > 0)
680  {
681  outSnake.PrintMesh(testSnake.snakeconn, 1, totT);
682  testTriangle.stattri.clear();
683  testTriangle.trivert.clear();
684  testTriangle.PrepareForUse();
685  TriangulateMesh(testSnake.snakeconn, testTriangle);
686  testTriangle.CalcTriVertPos();
687 
688  MeshTriangulation(triMesh, testSnake.snakeconn, testTriangle.stattri, testTriangle.trivert);
689  outSnake.PrintMesh(triMesh, 2, totT);
690  }
691  stop_s = clock();
692  cout << "time: " << (stop_s - start_s) / double(CLOCKS_PER_SEC) * 1000 << "ms" << endl;
693  testSnake.displight();
694  // the code you wish to time goes here
695 
696  // testSnake.disp();
697  }
698  catch (exception const &ex)
699  {
700  cerr << "Exception: " << ex.what() << endl;
701  return -1;
702  }
703  return (errTest);
704 }
705 
706 int Test_snakeinitflat()
707 {
708  snake testSnake;
709  mesh snakeMesh;
710  const char *fileToOpen;
711  tecplotfile outSnake;
712  int start_s, stop_s, ii;
713  vector<double> dt;
714  vector<int> isImpact;
715  double totT = 0.0;
716  // bool errFlag;
717  int errTest = 0;
718 
719  try
720  {
721  fileToOpen = "../TESTOUT/TestSnake2.plt";
722 
723  outSnake.OpenFile(fileToOpen);
724  errTest += snakeMesh.read("../TESTOUT/tecout100100.dat");
725  snakeMesh.PrepareForUse();
726  snakeMesh.SetBorders();
727  snakeMesh.PrepareForUse();
728 
729  testSnake.SetSnakeMesh(&snakeMesh);
730  outSnake.PrintMesh(*(testSnake.snakemesh()));
731 
732  start_s = clock();
733  testSnake.PrepareForUse();
734  SpawnAtVertex(testSnake, 103);
735  SpawnAtVertex(testSnake, 4020);
736  SpawnAtVertex(testSnake, 780);
737  testSnake.displight();
738 
739  testSnake.PrepareForUse();
740  for (ii = 0; ii < 50; ++ii)
741  {
742  cout << ii << " ";
743  if (testSnake.snaxs.size() > 0)
744  {
745  outSnake.PrintMesh(testSnake.snakeconn, 1, totT);
746  }
747 
748  Test_stepalgo(testSnake, isImpact);
749  cout << endl;
750  totT = totT + 1;
751  }
752  // testSnake.disp();
753  // the code you wish to time goes here
754  stop_s = clock();
755  cout << "time: " << (stop_s - start_s) / double(CLOCKS_PER_SEC) * 1000 << "ms" << endl;
756 
757  // testSnake.disp();
758  }
759  catch (exception const &ex)
760  {
761  cerr << "Exception: " << ex.what() << endl;
762  return -1;
763  }
764  return (errTest);
765 }
766 
767 void Test_stepalgo(snake &testSnake, vector<int> &isImpact)
768 {
769  vector<double> dt;
770  CalculateSnakeVelRand(testSnake);
771  testSnake.CalculateTimeStep(dt, 0.9);
772  testSnake.UpdateDistance(dt, 0.9, true);
773  testSnake.PrepareForUse();
774  testSnake.UpdateCoord();
775  SnakeConnectivityUpdate(testSnake, isImpact);
776 }
777 
778 void Test_stepalgo_mergeclean(snake &testSnake, vector<int> &isImpact)
779 {
780 
781  int start_s;
782 
783  start_s = clock();
784 
785  start_s = rsvs3d::TimeStamp("position: ", start_s);
786 
787  testSnake.PrepareForUse();
788  MergeCleanSnake(testSnake, isImpact);
789  testSnake.PrepareForUse();
790  start_s = rsvs3d::TimeStamp("MergeClean: ", start_s);
791 
792  testSnake.SnaxImpactDetection(isImpact);
793  SpawnArrivedSnaxels(testSnake, isImpact);
794  start_s = rsvs3d::TimeStamp("Spawn: ", start_s);
795 
796  testSnake.PrepareForUse();
797  MergeCleanSnake(testSnake, isImpact);
798  testSnake.PrepareForUse();
799 
800  testSnake.OrientFaces();
801  start_s = rsvs3d::TimeStamp("Clean: ", start_s);
802 }
803 
804 int Test_snakeOrderEdges()
805 {
806 
807  mesh snakeMesh;
808  const char *fileToOpen;
809  tecplotfile outSnake;
810 
811  // bool errFlag;
812  int errTest = 0;
813 
814  try
815  {
816  fileToOpen = "../TESTOUT/TestOrderEdges.plt";
817 
818  outSnake.OpenFile(fileToOpen);
819  errTest += snakeMesh.read("../TESTOUT/mesh234.dat");
820  snakeMesh.HashArray();
821  snakeMesh.SetMaxIndex();
822 
823  outSnake.PrintMesh(snakeMesh);
824 
825  snakeMesh.OrderEdges();
826  outSnake.PrintMesh(snakeMesh);
827  // testSnake.disp();
828  }
829  catch (exception const &ex)
830  {
831  cerr << "Exception: " << ex.what() << endl;
832  return -1;
833  }
834  return (errTest);
835 }
836 
837 int Test_MeshRefinement()
838 {
839 
840  mesh snakeMesh, voluMesh, triMesh;
841  const char *fileToOpen;
842  tecplotfile outSnake;
843  vector<int> elmMapping, dims;
844  triangulation testTriangle;
845  RSVScalc calcObj, calcObj2;
846  int ii;
847 
848  // bool errFlag;
849  int errTest = 0;
850 
851  try
852  {
853  fileToOpen = "../TESTOUT/Test_Multimesh.plt";
854 
855  outSnake.OpenFile(fileToOpen);
856  errTest += snakeMesh.read("../TESTOUT/mesh203010.dat");
857  snakeMesh.PrepareForUse();
858  outSnake.PrintMesh(snakeMesh);
859 
860  snakeMesh.OrderEdges();
861  snakeMesh.OrientFaces();
862  outSnake.PrintMesh(snakeMesh);
863  // testSnake.disp();
864  for (ii = 0; ii < snakeMesh.volus.size(); ++ii)
865  {
866  elmMapping.push_back(1);
867  }
868  dims.assign(3, 0);
869  dims[0] = 2;
870  dims[1] = 3;
871  dims[2] = 1;
872  CartesianMapping(snakeMesh, elmMapping, dims);
873  CoarsenMesh(snakeMesh, voluMesh, elmMapping);
874  snakeMesh.AddParent(&voluMesh, elmMapping);
875 
876  sort(elmMapping);
877  unique(elmMapping);
878  DisplayVector(elmMapping);
879  for (ii = 0; ii < voluMesh.volus.size(); ++ii)
880  {
881  voluMesh.volus[ii].fill = (double(rand() % 1001) / 1000.0);
882  voluMesh.volus[ii].target = voluMesh.volus[ii].fill;
883  voluMesh.volus[ii].error = voluMesh.volus[ii].fill;
884  }
885  voluMesh.PrepareForUse();
886  voluMesh.TightenConnectivity();
887  voluMesh.OrderEdges();
888  snakeMesh.OrientFaces();
889  voluMesh.OrientFaces();
890  triMesh.OrientFaces();
891 
892  testTriangle.PrepareForUse();
893  TriangulateMesh(snakeMesh, testTriangle);
894  testTriangle.PrepareForUse();
895  testTriangle.CalcTriVertPos();
896 
897  testTriangle.acttri.clear();
898  testTriangle.acttri.reserve(testTriangle.stattri.size());
899  for (ii = 0; ii < testTriangle.stattri.size(); ++ii)
900  {
901  testTriangle.acttri.push_back(testTriangle.stattri(ii)->index);
902  }
903 
904  testTriangle.PrepareForUse();
905  calcObj.CalculateTriangulation(testTriangle);
906  calcObj.ReturnConstrToMesh(testTriangle);
907  calcObj.Print2Screen();
908 
909  calcObj2.CalculateMesh(voluMesh);
910  calcObj2.ReturnConstrToMesh(voluMesh, &volu::error);
911  calcObj2.Print2Screen();
912 
913  MeshTriangulation(triMesh, voluMesh, testTriangle.stattri, testTriangle.trivert);
914 
915  outSnake.PrintMesh(voluMesh);
916  outSnake.PrintMesh(triMesh);
917  }
918  catch (exception const &ex)
919  {
920  cerr << "Exception: " << ex.what() << endl;
921  return -1;
922  }
923  return (errTest);
924 }
925 
926 int Test_surfcentre()
927 {
928  // int ii,n;
929  // vector<int> vertind;
930  grid::coordlist veccoord;
931  SurfCentroid tempCalc;
932  vector<double> v1 = {0.0, 0.0, 0.0};
933  vector<double> v2 = {1.0, 1.0, 0.0};
934  vector<double> v3 = {0.0, 1.0, 1.0};
935  vector<double> v4 = {1.0, 0.0, 0.0};
936  vector<double> v5 = {1.0, 0.0, 0.0};
937  // ArrayVec<double> tempCoord,jac,hes;
938 
939  // coord.assign(0,0,0);
940  // n=int(surfin.edgeind.size());
941 
942  // veccoord.reserve(n);
943  // ConnVertFromConnEdge(meshin, surfin.edgeind,vertind);
944 
945  // for(ii=0; ii<n; ++ii){
946  // veccoord.push_back(&(meshin.verts.isearch(vertind[ii])->coord));
947  // }
948  veccoord.push_back(&v1);
949  veccoord.push_back(&v2);
950  veccoord.push_back(&v3);
951  veccoord.push_back(&v4);
952  veccoord.push_back(&v5);
953  veccoord.push_back(&v5);
954  veccoord.push_back(&v5);
955  veccoord.push_back(&v5);
956  tempCalc.assign(veccoord);
957  tempCalc.Calc();
958 
959  // tempCalc.ReturnDat(tempCoord,jac,hes);
960  // coord.assign(tempCoord[0][0],tempCoord[1][0],tempCoord[2][0]);
961 
962  return (0);
963 }
964 
965 int Test_RSVSalgo_init()
966 {
967  // int nVoluZone, ii;
968 
969  snake testSnake, testSnake2;
970  mesh snakeMesh, voluMesh, voluMesh2;
971  // mesh triMesh;
972  triangulation testTriangle, triRSVS, triRSVS2;
973  vector<int> dims;
974  const char *fileToOpen;
975  tecplotfile outSnake;
976  // double totT=0.0;
977  // vector<double> dt;
978  // vector<int> isImpact;
979  int start_s, stop_s;
980  // bool errFlag;
981  int errTest = 0;
982 
983  dims.assign(3, 0);
984  dims[0] = 2;
985  dims[1] = 3;
986  dims[2] = 1;
987  try
988  {
989  fileToOpen = "../TESTOUT/TestAlgoRSVS.plt";
990 
991  outSnake.OpenFile(fileToOpen);
992  errTest += snakeMesh.read("../TESTOUT/mesh203010.dat");
993 
994  PrepareMultiLvlSnake(snakeMesh, voluMesh, testSnake, dims, triRSVS);
995  PrepareMultiLvlSnake(snakeMesh, voluMesh2, testSnake2, dims, triRSVS2);
996  voluMesh.volus[0].target = 0.0;
997  voluMesh.volus[3].target = 0.0;
998  voluMesh.volus[4].target = 1.0;
999  voluMesh.PrepareForUse();
1000  outSnake.PrintMesh(*(testSnake.snakemesh()));
1001  outSnake.PrintMesh(voluMesh);
1002  // nVoluZone=outSnake.ZoneNum();
1003 
1004  start_s = clock();
1005  SpawnRSVS(testSnake2, 0);
1006  SpawnRSVS(testSnake, 1);
1007  testSnake.PrepareForUse();
1008  testSnake2.PrepareForUse();
1009 
1010  stop_s = clock();
1011  cout << "time: " << (stop_s - start_s) / double(CLOCKS_PER_SEC) * 1000 << "ms" << endl;
1012  testSnake.displight();
1013  outSnake.PrintMesh(testSnake.snakeconn);
1014  outSnake.PrintSnakeInternalPts(testSnake);
1015  outSnake.PrintMesh(testSnake2.snakeconn);
1016  outSnake.PrintSnakeInternalPts(testSnake2);
1017 
1018  if (testSnake.snaxs.size() == 0 && testSnake2.snaxs.size() == 0)
1019  {
1020  errTest = 3;
1021  RSVS3D_ERROR_NOTHROW("Both test snakes had no snaxels.");
1022  }
1023  else if (testSnake2.snaxs.size() == 0)
1024  {
1025  errTest = 2;
1026  RSVS3D_ERROR_NOTHROW("test snake 2 had no snaxels.");
1027  }
1028  else if (testSnake.snaxs.size() == 0)
1029  {
1030  errTest = 1;
1031  RSVS3D_ERROR_NOTHROW("test snake 1 had no snaxels.");
1032  }
1033  }
1034  catch (exception const &ex)
1035  {
1036  cerr << "Exception: " << ex.what() << endl;
1037  return -1;
1038  }
1039  return (errTest);
1040 }
1041 
1042 int Test_RSVSvoro_init()
1043 {
1044  // int nVoluZone, ii;
1045 
1046  snake testSnake, testSnake2;
1047  mesh snakeMesh, voluMesh;
1048  // mesh triMesh;
1049  triangulation testTriangle, triRSVS, triRSVS2;
1050  vector<int> dims;
1051  vector<double> vecPts;
1052  const char *fileToOpen;
1053  tecplotfile outSnake;
1054  tetgen::apiparam voroparam;
1055  // double totT=0.0;
1056  // vector<double> dt;
1057  // vector<int> isImpact;
1058  int start_s, stop_s;
1059  // bool errFlag;
1060  int errTest = 0;
1061 
1062  dims.assign(3, 0);
1063  dims[0] = 2;
1064  dims[1] = 3;
1065  dims[2] = 1;
1066  try
1067  {
1068  fileToOpen = "../TESTOUT/TestVoroRSVS.plt";
1069  outSnake.OpenFile(fileToOpen);
1070 
1071  voroparam.edgelengths = {0.0, 3.0};
1072  voroparam.distanceTol = 0.5;
1073  vecPts.assign(4, 0.66);
1074  vecPts.insert(vecPts.end(), 4, 0.33);
1075  tetgen::RSVSVoronoiMesh(vecPts, voluMesh, snakeMesh, voroparam);
1076  PrepareMultiLvlSnakeNoVoluGen(snakeMesh, voluMesh, testSnake, triRSVS);
1077  PrepareMultiLvlSnakeNoVoluGen(snakeMesh, voluMesh, testSnake2, triRSVS2);
1078  // for (int i = 0; i < voluMesh.volus.size(); ++i)
1079  // {
1080  // voluMesh.volus[i].target = 0.5;
1081  // }
1082  int i1 = 1;
1083  for (int i = 0; i < 3; ++i)
1084  {
1085  voluMesh.volus[i1++].target = 1.0;
1086  }
1087  // voluMesh.volus[3].target=0.0;
1088  voluMesh.PrepareForUse();
1089  outSnake.PrintMesh(*(testSnake.snakemesh()));
1090  outSnake.PrintMesh(voluMesh);
1091  // nVoluZone=outSnake.ZoneNum();
1092 
1093  start_s = clock();
1094  SpawnRSVS(testSnake2, 1);
1095  SpawnRSVS(testSnake, 0);
1096  testSnake.PrepareForUse();
1097 
1098  stop_s = clock();
1099  cout << "time: " << (stop_s - start_s) / double(CLOCKS_PER_SEC) * 1000 << "ms" << endl;
1100  testSnake.displight();
1101  outSnake.PrintMesh(testSnake.snakeconn);
1102  outSnake.PrintSnakeInternalPts(testSnake);
1103  outSnake.PrintMesh(testSnake2.snakeconn);
1104  outSnake.PrintSnakeInternalPts(testSnake2);
1105 
1106  if (testSnake.snaxs.size() == 0 && testSnake2.snaxs.size() == 0)
1107  {
1108  errTest = 3;
1109  RSVS3D_ERROR_NOTHROW("Both test snakes had no snaxels.");
1110  }
1111  else if (testSnake2.snaxs.size() == 0)
1112  {
1113  errTest = 2;
1114  RSVS3D_ERROR_NOTHROW("test snake 2 had no snaxels.");
1115  }
1116  else if (testSnake.snaxs.size() == 0)
1117  {
1118  errTest = 1;
1119  RSVS3D_ERROR_NOTHROW("test snake 1 had no snaxels.");
1120  }
1121  }
1122  catch (exception const &ex)
1123  {
1124  cerr << "Exception: " << ex.what() << endl;
1125  return -1;
1126  }
1127  return (errTest);
1128 }
1129 
1130 int Test_RSVSalgo()
1131 {
1132  // int nVoluZone, ii;
1133 
1134  snake testSnake;
1135  mesh snakeMesh, voluMesh, triMesh;
1136  // mesh triMesh;
1137  triangulation testTriangle, triRSVS;
1138  vector<int> dims;
1139  const char *fileToOpen;
1140  tecplotfile outSnake, outSnake2;
1141  // double totT=0.0;
1142  // vector<double> dt;
1143  // vector<int> isImpact;
1144  int start_s, stop_s;
1145  // bool errFlag;
1146  int errTest = 0;
1147  int nVoluZone;
1148  RSVScalc calcObj;
1149 
1150  dims.assign(3, 0);
1151  dims[0] = 2;
1152  dims[1] = 3;
1153  dims[2] = 1;
1154  try
1155  {
1156  fileToOpen = "../TESTOUT/TestAlgoRSVSstep.plt";
1157  outSnake.OpenFile(fileToOpen);
1158  fileToOpen = "../TESTOUT/TestAlgoRSVSstep_snake.plt";
1159  outSnake2.OpenFile(fileToOpen);
1160  errTest += snakeMesh.read("../TESTOUT/mesh203010.dat");
1161 
1162  PrepareMultiLvlSnake(snakeMesh, voluMesh, testSnake, dims, triRSVS);
1163  calcObj.BuildMathArrays(1, 1);
1164  voluMesh.volus[0].target = 0.0;
1165  voluMesh.volus[3].target = 0.0;
1166  voluMesh.volus[4].target = 1.0;
1167  voluMesh.PrepareForUse();
1168  outSnake.PrintMesh(*(testSnake.snakemesh()));
1169  outSnake.PrintMesh(voluMesh);
1170  nVoluZone = outSnake.ZoneNum();
1171 
1172  start_s = clock();
1173 
1174  SpawnRSVS(testSnake, 1);
1175  testSnake.PrepareForUse();
1176 
1177  triRSVS.PrepareForUse();
1178 
1179  TriangulateSnake(testSnake, triRSVS);
1180  triRSVS.PrepareForUse();
1181  triRSVS.CalcTriVertPos();
1182  int ii;
1183  double totT = 0.0;
1184  vector<double> dt;
1185  vector<int> isImpact;
1186  MaintainTriangulateSnake(triRSVS);
1187 
1188  for (ii = 0; ii < 50; ++ii)
1189  {
1190  cout << ii << " ";
1191  PrintRSVSSnake(outSnake, testSnake, totT, testTriangle, triMesh, triRSVS, voluMesh, nVoluZone, ii);
1192  stop_s = clock();
1193  Test_stepalgoRSVS(testSnake, triRSVS, dt, isImpact, calcObj, outSnake2, totT);
1194  stop_s = rsvs3d::TimeStamp("Total: ", stop_s);
1195  cout << endl;
1196  totT = totT + 1;
1197  }
1198 
1199  PrintRSVSSnake(outSnake, testSnake, totT, testTriangle, triMesh, triRSVS, voluMesh, nVoluZone, ii);
1200 
1201  stop_s = clock();
1202  cout << "time: " << (stop_s - start_s) / double(CLOCKS_PER_SEC) * 1000 << "ms" << endl;
1203  testSnake.displight();
1204  outSnake.PrintMesh(testSnake.snakeconn);
1205  outSnake.PrintSnakeInternalPts(testSnake);
1206  }
1207  catch (exception const &ex)
1208  {
1209  cerr << "Exception: " << ex.what() << endl;
1210  return -1;
1211  }
1212  return (errTest);
1213 }
1214 
1215 int Test_RSVSalgoflat()
1216 {
1217  // int nVoluZone, ii;
1218 
1219  snake testSnake;
1220  mesh snakeMesh, voluMesh, triMesh;
1221  // mesh triMesh;
1222  triangulation testTriangle, triRSVS;
1223  vector<int> dims;
1224  const char *fileToOpen;
1225  tecplotfile outSnake, outSnake2;
1226  // double totT=0.0;
1227  // vector<double> dt;
1228  // vector<int> isImpact;
1229  int start_s, stop_s;
1230  // bool errFlag;
1231  int errTest = 0;
1232  RSVScalc calcObj;
1233 
1234  dims.assign(3, 0);
1235  dims[0] = 2;
1236  dims[1] = 2;
1237  dims[2] = 0;
1238  try
1239  {
1240  fileToOpen = "../TESTOUT/TestAlgoRSVS2Dstep.plt";
1241  outSnake.OpenFile(fileToOpen);
1242  fileToOpen = "../TESTOUT/TestAlgoRSVS2Dstep_snake.plt";
1243  outSnake2.OpenFile(fileToOpen);
1244  errTest += snakeMesh.read("../TESTOUT/tecout100100.dat");
1245 
1246  PrepareMultiLvlSnake(snakeMesh, voluMesh, testSnake, dims, triRSVS);
1247 
1248  voluMesh.surfs[0].target = 1.0;
1249  voluMesh.surfs[1].target = 1.0;
1250  voluMesh.surfs[3].target = 1.0;
1251  voluMesh.surfs[2].target = 1.0;
1252 
1253  calcObj.BuildMathArrays(1, 1);
1254  voluMesh.PrepareForUse();
1255  outSnake.PrintMesh(*(testSnake.snakemesh()));
1256  outSnake.PrintMesh(voluMesh);
1257  int nVoluZone;
1258  nVoluZone = outSnake.ZoneNum();
1259 
1260  start_s = clock();
1261 
1262  SpawnRSVS(testSnake, 1);
1263  testSnake.PrepareForUse();
1264 
1265  triRSVS.PrepareForUse();
1266 
1267  TriangulateSnake(testSnake, triRSVS);
1268  triRSVS.PrepareForUse();
1269  triRSVS.CalcTriVertPos();
1270  int ii;
1271  double totT = 0.0;
1272  vector<double> dt;
1273  vector<int> isImpact;
1274  MaintainTriangulateSnake(triRSVS);
1275  for (ii = 0; ii < 20; ++ii)
1276  {
1277  // testSnake.displight();
1278  cout << ii << " ";
1279  PrintRSVSSnake2D(outSnake, testSnake, totT, testTriangle, triMesh, triRSVS, voluMesh, nVoluZone, ii);
1280  stop_s = clock();
1281 
1282  Test_stepalgoRSVS(testSnake, triRSVS, dt, isImpact, calcObj, outSnake2, totT);
1283  // Test_stepalgo(testSnake, isImpact);
1284  // MaintainTriangulateSnake(triRSVS);
1285  stop_s = rsvs3d::TimeStamp("Total: ", stop_s);
1286  cout << endl;
1287  totT = totT + 1;
1288  }
1289 
1290  PrintRSVSSnake2D(outSnake, testSnake, totT, testTriangle, triMesh, triRSVS, voluMesh, nVoluZone, ii);
1291 
1292  stop_s = clock();
1293  cout << "time: " << (stop_s - start_s) / double(CLOCKS_PER_SEC) * 1000 << "ms" << endl;
1294  testSnake.displight();
1295  outSnake.PrintMesh(testSnake.snakeconn);
1296  outSnake.PrintSnakeInternalPts(testSnake);
1297  }
1298  catch (exception const &ex)
1299  {
1300  cerr << "Exception: " << ex.what() << endl;
1301  return -1;
1302  }
1303  return (errTest);
1304 }
1305 
1306 int Test_snakeRSVS()
1307 {
1308  int nVoluZone;
1309  snake testSnake;
1310  mesh snakeMesh, triMesh, voluMesh;
1311  triangulation testTriangle, triRSVS;
1312  vector<int> dims;
1313  const char *fileToOpen;
1314  tecplotfile outSnake, outSnake2;
1315  double totT = 0.0;
1316  vector<double> dt;
1317  vector<int> isImpact;
1318  int start_s, stop_s, ii;
1319  // bool errFlag;
1320  int errTest = 0;
1321  RSVScalc calcObj;
1322 
1323  dims.assign(3, 0);
1324  dims[0] = 2;
1325  dims[1] = 3;
1326  dims[2] = 1;
1327  try
1328  {
1329  fileToOpen = "../TESTOUT/TestSnakeRSVS.plt";
1330  outSnake.OpenFile(fileToOpen);
1331  fileToOpen = "../TESTOUT/TestSnakeRSVS_snake.plt";
1332  outSnake2.OpenFile(fileToOpen);
1333 
1334  errTest += snakeMesh.read("../TESTOUT/mesh203010.dat");
1335  PrepareMultiLvlSnake(snakeMesh, voluMesh, testSnake, dims, triRSVS);
1336 
1337  outSnake.PrintMesh(*(testSnake.snakemesh()));
1338  outSnake.PrintMesh(voluMesh);
1339  nVoluZone = outSnake.ZoneNum();
1340 
1341  // SpawnAtVertex(testSnake,1022);
1342  // SpawnAtVertex(testSnake,674);
1343  // SpawnAtVertex(testSnake,675);
1344  // SpawnAtVertex(testSnake,728);
1345  // SpawnAtVertex(testSnake,729);
1346  SpawnAtVertex(testSnake, 731);
1347  testSnake.displight();
1348 
1349  start_s = clock();
1350  testSnake.PrepareForUse();
1351  triRSVS.PrepareForUse();
1352 
1353  TriangulateSnake(testSnake, triRSVS);
1354  triRSVS.PrepareForUse();
1355  triRSVS.CalcTriVertPos();
1356  MaintainTriangulateSnake(triRSVS);
1357 
1358  for (ii = 0; ii < 20; ++ii)
1359  {
1360  cout << ii << " ";
1361  PrintRSVSSnake(outSnake, testSnake, totT, testTriangle, triMesh, triRSVS, voluMesh, nVoluZone, ii);
1362  stop_s = clock();
1363  Test_stepalgoRSVS(testSnake, triRSVS, dt, isImpact, calcObj, outSnake2, totT);
1364  stop_s = rsvs3d::TimeStamp("Total: ", stop_s);
1365  cout << endl;
1366  totT = totT + 1;
1367  }
1368 
1369  PrintRSVSSnake(outSnake, testSnake, totT, testTriangle, triMesh, triRSVS, voluMesh, nVoluZone, ii);
1370 
1371  stop_s = clock();
1372  cout << "time: " << (stop_s - start_s) / double(CLOCKS_PER_SEC) * 1000 << "ms" << endl;
1373  testSnake.displight();
1374  }
1375  catch (exception const &ex)
1376  {
1377  cerr << "Exception: " << ex.what() << endl;
1378  return -1;
1379  }
1380  return (errTest);
1381 }
1382 
1383 void repositionatquarteredge(snake &snakein)
1384 {
1385  /*Assumes there are only 2 snakes per edge*/
1386 
1387  for (int ii = 0; ii < int(snakein.snaxs.size()); ++ii)
1388  {
1389  if (snakein.snaxs(ii)->orderedge > 1)
1390  {
1391  snakein.snaxs[ii].d = 0.75;
1392  }
1393  else
1394  {
1395  snakein.snaxs[ii].d = 0.25;
1396  }
1397  }
1398 
1399  snakein.PrepareForUse();
1400  snakein.UpdateCoord();
1401 }
1402 
1403 int Test_RSVSalgo_singlevol(int sparseCutOff)
1404 {
1405  // int nVoluZone, ii;
1406 
1407  snake testSnake;
1408  mesh snakeMesh, voluMesh, triMesh, triMesh2, triMesh3, triMesh4;
1409  // mesh triMesh;
1410  triangulation testTriangle, triRSVS;
1411  vector<int> dims;
1412  const char *fileToOpen;
1413  tecplotfile outSnake, outSnake2;
1414  // double totT=0.0;
1415  // vector<double> dt;
1416  // vector<int> isImpact;
1417  int start_s, stop_s;
1418  // bool errFlag;
1419  int errTest = 0;
1420  RSVScalc calcObj;
1421  int ii;
1422  double totT = 0.0;
1423  vector<double> dt;
1424  vector<int> isImpact;
1425 
1426  dims.assign(3, 0);
1427  dims[0] = 1;
1428  dims[1] = 1;
1429  dims[2] = 3;
1430  try
1431  {
1432  fileToOpen = "../TESTOUT/TestAlgoRSVSstep.plt";
1433 
1434  outSnake.OpenFile(fileToOpen);
1435 
1436  fileToOpen = "../TESTOUT/TestAlgoRSVSstep_snake.plt";
1437  outSnake2.OpenFile(fileToOpen);
1438  errTest += snakeMesh.read("../TESTOUT/mesh203010.dat");
1439 
1440  PrepareMultiLvlSnake(snakeMesh, voluMesh, testSnake, dims, triRSVS);
1441  for (ii = 0; ii < voluMesh.volus.size(); ++ii)
1442  {
1443  voluMesh.volus[ii].target = 0.8;
1444  }
1445  voluMesh.volus[1].target = 0.1;
1446  voluMesh.volus[2].target = 0.3;
1447  // voluMesh.volus[0].target=0.001;//0.05;//0.0001;
1448  // voluMesh.volus[1].target=0.05;//0.05;
1449  // voluMesh.volus[2].target=0.001;//0.05;//0.0001;
1450  // voluMesh.volus[3].target=0.3;//0.05;//0.0001;
1451  // voluMesh.volus[4].target=0.3;//0.05;//0.0001;
1452  // voluMesh.volus[3].target=1;//0.05;//0.0001;
1453  // voluMesh.volus[4].target=1;//0.05;
1454  // voluMesh.volus[5].target=1;//0.05;//0.0001;
1455  // voluMesh.volus[3].target=1;//0.05;//0.0001;
1456  voluMesh.PrepareForUse();
1457  outSnake.PrintMesh(*(testSnake.snakemesh()));
1458  outSnake.PrintMesh(voluMesh);
1459  int nVoluZone;
1460  nVoluZone = outSnake.ZoneNum();
1461 
1462  start_s = clock();
1463 
1464  SpawnRSVS(testSnake, 1);
1465  testSnake.PrepareForUse();
1466 
1467  triRSVS.PrepareForUse();
1468 
1469  TriangulateSnake(testSnake, triRSVS);
1470  triRSVS.PrepareForUse();
1471  triRSVS.CalcTriVertPos();
1472  MaintainTriangulateSnake(triRSVS);
1473  calcObj.SetSparseDVcutoff(sparseCutOff);
1474  for (ii = 0; ii < 2; ++ii)
1475  {
1476  cout << ii << " ";
1477  // if (ii%2==0){
1478  PrintRSVSSnake(outSnake, testSnake, totT, testTriangle, triMesh, triRSVS, voluMesh, nVoluZone, ii);
1479  // }
1480  stop_s = clock();
1481  Test_stepalgoRSVS(testSnake, triRSVS, dt, isImpact, calcObj, outSnake2, totT);
1482  stop_s = rsvs3d::TimeStamp("Total: ", stop_s);
1483  cout << endl;
1484  totT = totT + 1;
1485  }
1486 
1487  PrintRSVSSnake(outSnake, testSnake, totT, testTriangle, triMesh, triRSVS, voluMesh, nVoluZone, ii);
1488 
1489  stop_s = clock();
1490  cout << "time: " << (stop_s - start_s) / double(CLOCKS_PER_SEC) * 1000 << "ms" << endl;
1491  testSnake.displight();
1492  testSnake.PrepareForUse();
1493  // outSnake.PrintMesh(testSnake.snakeconn);
1494  // outSnake.PrintSnakeInternalPts(testSnake);
1495 
1496  repositionatquarteredge(testSnake);
1497  triRSVS.CalcTriVertPos();
1498  triRSVS.PrepareForUse();
1499  FILE *fid;
1500  try
1501  {
1502  triMesh2 = TriarrayToMesh(triRSVS, triRSVS.intertri);
1503  fid = fopen("../TESTOUT/triintertestbig.dat", "w");
1504  triMesh2.PrepareForUse();
1505  triMesh2.write(fid); //(triMesh,3,totT);
1506  fclose(fid);
1507  }
1508  catch (exception const &extriMesh2)
1509  {
1510  RSVS3D_ERROR_NOTHROW("intertri output fails");
1511  errTest++;
1512  }
1513 
1514  try
1515  {
1516  triMesh3 = TriarrayToMesh(triRSVS, triRSVS.dynatri);
1517  fid = fopen("../TESTOUT/snakeconnouttri.dat", "w");
1518  triMesh3.PrepareForUse();
1519  triMesh3.write(fid); //(triMesh,3,totT);
1520  fclose(fid);
1521  }
1522  catch (exception const &extriMesh2)
1523  {
1524  RSVS3D_ERROR_NOTHROW("dynatri output fails");
1525  errTest++;
1526  }
1527 
1528  try
1529  {
1530  triMesh4 = TriarrayToMesh(triRSVS, triRSVS.stattri);
1531  fid = fopen("../TESTOUT/stattri.dat", "w");
1532  triMesh4.PrepareForUse();
1533  triMesh4.write(fid); //(triMesh,3,totT);
1534  fclose(fid);
1535  }
1536  catch (exception const &extriMesh2)
1537  {
1538  RSVS3D_ERROR_NOTHROW("stattri output fails");
1539  errTest++;
1540  }
1541  }
1542  catch (exception const &ex)
1543  {
1544  cerr << "Exception: " << ex.what() << endl;
1545  return -1;
1546  }
1547  return (errTest);
1548 }
1549 
1550 int Test_RSVSalgo_singlevol_fullmath()
1551 {
1552  return (Test_RSVSalgo_singlevol(10000000));
1553 }
1554 
1555 int Test_RSVSalgo_singlevol_sparse()
1556 {
1557  return (Test_RSVSalgo_singlevol(0));
1558 }
1559 
1560 int Test_snakeRSVS_singlevol()
1561 {
1562  int nVoluZone;
1563  snake testSnake;
1564  mesh snakeMesh, triMesh, voluMesh;
1565  triangulation testTriangle, triRSVS;
1566  vector<int> dims;
1567  const char *fileToOpen;
1568  tecplotfile outSnake, outSnake2;
1569  double totT = 0.0;
1570  vector<double> dt;
1571  vector<int> isImpact;
1572  int start_s, stop_s, ii;
1573  // bool errFlag;
1574  int errTest = 0;
1575  RSVScalc calcObj;
1576 
1577  dims.assign(3, 0);
1578  dims[0] = 1;
1579  dims[1] = 1;
1580  dims[2] = 1;
1581  try
1582  {
1583  fileToOpen = "../TESTOUT/TestSnakeRSVS.plt";
1584  outSnake.OpenFile(fileToOpen);
1585  fileToOpen = "../TESTOUT/TestSnakeRSVS_snake.plt";
1586  outSnake2.OpenFile(fileToOpen);
1587 
1588  errTest += snakeMesh.read("../TESTOUT/mesh203010.dat");
1589 
1590  PrepareMultiLvlSnake(snakeMesh, voluMesh, testSnake, dims, triRSVS);
1591  voluMesh.volus[0].target = 0.01;
1592  voluMesh.volus.PrepareForUse();
1593  outSnake.PrintMesh(*(testSnake.snakemesh()));
1594  outSnake.PrintMesh(voluMesh);
1595  nVoluZone = outSnake.ZoneNum();
1596 
1597  // SpawnAtVertex(testSnake,1022);
1598  // SpawnAtVertex(testSnake,674);
1599  // SpawnAtVertex(testSnake,675);
1600  // SpawnAtVertex(testSnake,728);
1601  // SpawnAtVertex(testSnake,729);
1602  SpawnAtVertex(testSnake, 731);
1603  testSnake.displight();
1604 
1605  start_s = clock();
1606  testSnake.PrepareForUse();
1607  triRSVS.PrepareForUse();
1608 
1609  TriangulateSnake(testSnake, triRSVS);
1610  triRSVS.PrepareForUse();
1611  triRSVS.CalcTriVertPos();
1612  MaintainTriangulateSnake(triRSVS);
1613  for (ii = 0; ii < 20; ++ii)
1614  {
1615  cout << ii << " ";
1616  PrintRSVSSnake(outSnake, testSnake, totT, testTriangle, triMesh, triRSVS, voluMesh, nVoluZone, ii);
1617  stop_s = clock();
1618  Test_stepalgoRSVS(testSnake, triRSVS, dt, isImpact, calcObj, outSnake2, totT);
1619  stop_s = rsvs3d::TimeStamp("Total: ", stop_s);
1620  cout << endl;
1621  totT = totT + 1;
1622  }
1623 
1624  PrintRSVSSnake(outSnake, testSnake, totT, testTriangle, triMesh, triRSVS, voluMesh, nVoluZone, ii);
1625 
1626  stop_s = clock();
1627  cout << "time: " << (stop_s - start_s) / double(CLOCKS_PER_SEC) * 1000 << "ms" << endl;
1628  testSnake.displight();
1629  }
1630  catch (exception const &ex)
1631  {
1632  cerr << "Exception: " << ex.what() << endl;
1633  return -1;
1634  }
1635  return (errTest);
1636 }
1637 
1638 int Test_MeshOrient()
1639 {
1640 
1641  // int nVoluZone, ii;
1642 
1643  snake testSnake;
1644  mesh snakeMesh, voluMesh;
1645  // mesh triMesh;
1646  triangulation testTriangle, triRSVS;
1647  vector<int> dims;
1648  const char *fileToOpen;
1649  tecplotfile outSnake;
1650  // double totT=0.0;
1651  // vector<double> dt;
1652  // vector<int> isImpact;
1653 
1654  // bool errFlag;
1655  int errTest = 0;
1656  RSVScalc calcObj;
1657  FILE *fid;
1658 
1659  dims.assign(3, 0);
1660  dims[0] = 1;
1661  dims[1] = 3;
1662  dims[2] = 1;
1663 
1664  fileToOpen = "../TESTOUT/TestAlgoRSVSstep.plt";
1665 
1666  outSnake.OpenFile(fileToOpen);
1667  errTest += snakeMesh.read("../TESTOUT/mesh6612.dat");
1668 
1669  PrepareMultiLvlSnake(snakeMesh, voluMesh, testSnake, dims, triRSVS);
1670  voluMesh.volus[0].target = 1; // 0.05;//0.0001;
1671  voluMesh.volus[1].target = 1; // 0.05;
1672  voluMesh.volus[2].target = 1; // 0.05;//0.0001;
1673  // voluMesh.volus[3].target=1;//0.05;//0.0001;
1674  voluMesh.PrepareForUse();
1675  outSnake.PrintMesh(*(testSnake.snakemesh()));
1676  outSnake.PrintMesh(voluMesh);
1677 
1678  fileToOpen = "../TESTOUT/volumesh6612.dat";
1679  fid = fopen(fileToOpen, "w");
1680  if (fid != NULL)
1681  {
1682  voluMesh.PrepareForUse();
1683  voluMesh.SetBorders();
1684  voluMesh.OrientFaces();
1685  voluMesh.write(fid);
1686  }
1687  return (errTest);
1688 }
1689 
1690 void Test_stepalgoRSVS(snake &testSnake, triangulation &RSVStri, vector<double> &dt, vector<int> &isImpact,
1691  RSVScalc &calcObj, tecplotfile &outSnake2, double totT)
1692 {
1693  int start_s;
1694 
1695  start_s = clock();
1696 
1697  // calcObj.Print2Screen(1);
1698  // calcObj.Print2Screen(2);
1699  // CalculateSnakeVel(testSnake);
1700  // Small step away from edge without crossover.
1701  // Need to develop that.
1702  // Check if impact detect crossovers
1703  start_s = rsvs3d::TimeStamp(" triangulate:", start_s);
1704  // calcObj.limLag=10000.0;
1705  calcObj.CalculateTriangulation(RSVStri);
1706  calcObj.ReturnConstrToMesh(RSVStri);
1707  start_s = rsvs3d::TimeStamp(" deriv:", start_s);
1708  calcObj.CheckAndCompute(2);
1709  calcObj.ReturnVelocities(RSVStri);
1710  start_s = rsvs3d::TimeStamp(" solve:", start_s);
1711 
1712  calcObj.Print2Screen();
1713  // calcObj.Print2Screen(2);
1714  CalculateNoNanSnakeVel(testSnake);
1715  outSnake2.PrintSnake(testSnake, 1, totT);
1716  testSnake.CalculateTimeStep(dt, 0.4);
1717  testSnake.UpdateDistance(dt, 0.9, true);
1718  testSnake.PrepareForUse();
1719  testSnake.UpdateCoord();
1720 
1721  SnakeConnectivityUpdate(testSnake, isImpact);
1722  MaintainTriangulateSnake(RSVStri);
1723 }
1724 
1725 void Test_mathRSVS_FD(snake &testSnake, triangulation &RSVStri, vector<double> &dt, vector<int> &isImpact,
1726  RSVScalc &calcObj, tecplotfile &outSnake2, double totT)
1727 {
1728  int start_s;
1729 
1730  RSVScalc calcObj2;
1731  double fdStep = 1e-3;
1732  start_s = clock();
1733 
1734  // calcObj.Print2Screen(1);
1735  // calcObj.Print2Screen(2);
1736  // CalculateSnakeVel(testSnake);
1737  // Small step away from edge without crossover.
1738  // Need to develop that.
1739  // Check if impact detect crossovers
1740  start_s = rsvs3d::TimeStamp(" triangulate:", start_s);
1741  // calcObj.limLag=10000.0;
1742  // for(int ii=0; ii<testSnake.snakemesh()->verts.size();++ii){
1743  // for(int jj=0; jj<3;++jj){
1744  // testSnake.snakemesh()->verts[ii].coord[jj] +=1;
1745  // }
1746  // }
1747  // testSnake.PrepareForUse();
1748  // testSnake.UpdateCoord();
1749  // testSnake.PrepareForUse();
1750  calcObj.CalculateTriangulation(RSVStri, 2);
1751  testSnake.snaxs[137].d += fdStep;
1752  testSnake.PrepareForUse();
1753  testSnake.UpdateCoord();
1754  testSnake.PrepareForUse();
1755  calcObj2.CalculateTriangulation(RSVStri);
1756 
1757  cout << endl;
1758  for (int ii = 0; ii < calcObj.numConstr(); ++ii)
1759  {
1760  cout << calcObj2.constr[ii] - calcObj.constr[ii] << " " << (calcObj2.constr[ii] - calcObj.constr[ii]) / fdStep
1761  << endl;
1762  }
1763 
1764  calcObj.ReturnConstrToMesh(RSVStri);
1765  start_s = rsvs3d::TimeStamp(" deriv:", start_s);
1766  calcObj.CheckAndCompute(3);
1767  calcObj.ReturnVelocities(RSVStri);
1768  start_s = rsvs3d::TimeStamp(" solve:", start_s);
1769  outSnake2.PrintSnake(testSnake, 1, totT);
1770 
1771  calcObj.Print2Screen();
1772  // calcObj.Print2Screen(2);
1773  CalculateNoNanSnakeVel(testSnake);
1774  testSnake.CalculateTimeStep(dt, 0.5);
1775  testSnake.UpdateDistance(dt, 0.34, true);
1776  testSnake.PrepareForUse();
1777  testSnake.UpdateCoord();
1778 
1779  SnakeConnectivityUpdate(testSnake, isImpact);
1780  MaintainTriangulateSnake(RSVStri);
1781 }
1782 
1783 int Test_SurfCentreDerivatives()
1784 {
1785 
1786  std::vector<std::vector<double>> coords;
1787  double pi = 3.14159265358979323846;
1788  ofstream outfile;
1789  double stepLength;
1790 
1791  for (auto nCoords : {4, 5, 6, 7})
1792  {
1793 
1794  SurfCentroid surfCentre(nCoords);
1795  std::ostringstream filename;
1796  filename << "../TESTOUT/derivatives/surfcentre_jac" << nCoords << ".csv";
1797  outfile.open(filename.str());
1798  outfile.precision(16);
1799 
1800  coords.assign(nCoords, {0.0, 0.0, 0.0});
1801  for (int i = 0; i < nCoords; ++i)
1802  {
1803  surfCentre.assign(i, coords[i]);
1804  }
1805  for (int ord = 0; ord < 8; ++ord)
1806  {
1807  stepLength = pow(10.0, -ord);
1808  for (int i = 0; i < nCoords; ++i)
1809  {
1810  double x = sin(2.0 * pi / nCoords * i);
1811  double y = cos(2.0 * pi / nCoords * i);
1812  double z = -(x);
1813  coords[i][0] = x;
1814  coords[i][1] = y;
1815  coords[i][2] = z;
1816  }
1817 
1818  for (int i = 0; i < 3; ++i)
1819  {
1820  coords[0][i] = coords[1][i] + (coords[0][i] - coords[1][i]) * stepLength;
1821  }
1822  for (int i = 0; i < nCoords; ++i)
1823  {
1824  std::cout << "coord vector " << i << ": ";
1825  DisplayVector(coords[i]);
1826  std::cout << std::endl;
1827  }
1828  for (int i = 0; i < nCoords; ++i)
1829  {
1830  for (int j = 0; j < 3; ++j)
1831  {
1832  outfile << coords[i][j] << ", ";
1833  }
1834  outfile << std::endl;
1835  }
1836  surfCentre.Calc();
1837  surfCentre.jac_ptr().write(outfile);
1838  outfile << std::endl;
1839  }
1840  outfile.close();
1841  }
1842  return 0;
1843 }
1844 
1845 int Test_Matrix3D()
1846 {
1847  Eigen::MatrixXd array3d(2, 4), vec(1, 2), returnVec(2, 2), expected(2, 2);
1848 
1849  array3d << 0, 1, 4, 5, 2, 3, 6, 7;
1850  vec << 1, 2;
1851  expected << 8, 11, 14, 17;
1852 
1853  std::cout << "3d:" << std::endl << array3d << std::endl;
1854  std::cout << "vec:" << std::endl << vec << std::endl;
1855  std::cout << "expected:" << std::endl << expected << std::endl;
1856 
1857  VecBy3DimArray(vec, array3d, returnVec);
1858 
1859  std::cout << "returnVec:" << std::endl << returnVec << std::endl;
1860 
1861  return (0);
1862 }
RSVSalgorithm.hpp
Functions which are part of the RSVS algorithm but not core to the snaking process.
RSVScalc.hpp
Provides the infrastructure for calculation of the RSVS equations.
RSVSintegration.hpp
Integration into the full 3 dimensional Restricted Snake Volume of Solid method.
RSVSmath.hpp
Performs Volume and Area calculations for the RSVS process.
RSVScalc
Class to handle the RSVS calculation.
Definition: RSVScalc.hpp:130
RSVScalc::CalculateMesh
void CalculateMesh(mesh &meshin)
Calculates the mesh volumes.
Definition: RSVScalc.cpp:294
RSVScalc::constr
Eigen::VectorXd constr
Constraint value vector, size: [nConstr, 1].
Definition: RSVScalc.hpp:177
RSVScalc::ReturnConstrToMesh
void ReturnConstrToMesh(triangulation &triRSVS) const
Returns a constraint to the triangulation::meshDep.
Definition: RSVScalc.cpp:400
RSVScalc::CheckAndCompute
void CheckAndCompute(int calcMethod=0, bool sensCalc=false)
Prepare the active arrays for SQP calculation and calculate the SQP step.
Definition: RSVScalc_SQP.cpp:264
RSVScalc::CalculateTriangulation
void CalculateTriangulation(const triangulation &triRSVS, int derivMethod=0)
Calculates the triangulation volume and area derivatives.
Definition: RSVScalc.cpp:95
RSVScalc::Print2Screen
void Print2Screen(int outType=0) const
Prints different amounts of RSVScalc owned data to the screen.
Definition: RSVScalc.cpp:334
RSVScalc::numConstr
int numConstr() const
Getter for the number of constraints.
Definition: RSVScalc.hpp:455
RSVScalc::ReturnVelocities
void ReturnVelocities(triangulation &triRSVS)
Returns velocities to the snaxels.
Definition: RSVScalc.cpp:143
RSVScalc::BuildMathArrays
void BuildMathArrays(int nDv, int nConstr)
Builds mathematics arrays.
Definition: RSVScalc.cpp:443
coordvec
Handles the use and norm of a vector for which the norm and the unit value might be needed.
Definition: mesh.hpp:114
mesh
Class for mesh handling.
Definition: mesh.hpp:592
mesh::OrientFaces
void OrientFaces()
Orients either surfaces or edges depending on the dimensionality of the object.
Definition: mesh.cpp:4948
snake
Definition: snake.hpp:83
snax
Definition: snake.hpp:192
tetgen::apiparam::distanceTol
double distanceTol
Distance tolerance.
Definition: tetgenrsvs.hpp:147
tetgen::apiparam::edgelengths
std::vector< double > edgelengths
Controls the edgelengths at regular intervals.
Definition: tetgenrsvs.hpp:145
main.hpp
file containing the main functions and the command line parser.
matrixtools.hpp
Tools to support conversion, display and derivatives of Eigen matrices.
grid::coordlist
std::vector< const std::vector< double > * > coordlist
Defines a list of coordinates.
Definition: mesh.hpp:87
meshrefinement.hpp
Tools for the refinement and coarsening of meshes.
CoarsenMesh
void CoarsenMesh(const mesh &meshchild, mesh &newparent, const std::vector< int > &elmMapping)
Definition: meshrefinement.cpp:11
postprocessing.hpp
Provide tecplot file formating for mesh and snake outputs.
snake.hpp
Provides the core restricted surface snake container.
snakeengine.hpp
Functions needed to evolve the r-surface snake.
tetgenrsvs.hpp
Interface between the RSVS project and tetgen.
tetgen::RSVSVoronoiMesh
std::vector< int > RSVSVoronoiMesh(const std::vector< double > &vecPts, mesh &vosMesh, mesh &snakMesh, tetgen::apiparam &inparam)
Genrates Snaking and VOS RSVS meshes from the voronoi diagram of a set of points.
Definition: tetgenrsvs.cpp:876
triangulate.hpp
Provides a triangulation for snake and meshes.
warning.hpp
Provides the error and warning system used by the RSVS3D project.
RSVS3D_ERROR_NOTHROW
#define RSVS3D_ERROR_NOTHROW(M)
Generic rsvs warning.
Definition: warning.hpp:120
RSVS3D_ERROR_ARGUMENT
#define RSVS3D_ERROR_ARGUMENT(M)
Throw a invalid_argument.
Definition: warning.hpp:148