Skip to content

Commit 9c70de1

Browse files
committed
2 parents 90f64d3 + 05a852c commit 9c70de1

File tree

6 files changed

+213
-44
lines changed

6 files changed

+213
-44
lines changed

app/main.cpp

Lines changed: 70 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -54,7 +54,7 @@ void case_check()
5454
#endif
5555
}
5656

57-
# define TIDBG
57+
//# define TIDBG
5858

5959

6060
#ifdef TIGHT_INCLUSION_RUN_EXAMPLES
@@ -100,7 +100,10 @@ void run_rational_data_single_method(
100100
{
101101
for (int ff = 0; ff < 2; ff++)
102102
{
103-
103+
if (folders[fnbr] == "erleben-spike-hole"&&ff == 1)
104+
{
105+
continue;
106+
}
104107
all_V = read_rational_csv(
105108
#ifdef TIDBG
106109
"D:\\vs\\collision\\interval\\Tight-Inclusion\\build\\edge-edge-0474.csv",
@@ -205,10 +208,12 @@ void run_rational_data_single_method(
205208
}
206209
}
207210
}
211+
#ifdef TIDBG
208212
std::cout << "fps " << new_false_positives << " fns " << new_false_negatives << " total queries " << total_number + 1
209213
<<" total positives "<< total_positives << std::endl;
210214

211215
exit(0);
216+
#endif
212217
}
213218
}
214219

@@ -228,20 +233,78 @@ void run_code()
228233
{
229234
double ms = 0.0;
230235
double tolerance = 1e-6;
231-
std::cout << "\nRunning Vertex-Face data" << std::endl;
236+
237+
std::cout << "\nRunning simulation Vertex-Face data" << std::endl;
238+
run_rational_data_single_method(
239+
/*is_edge_edge*/ false, /*is_simulation*/ true, ms, tolerance);
240+
std::cout << "finish simulation Vertex-Face data" << std::endl;
241+
std::cout << "\nRunning simulation Edge-Edge data" << std::endl;
242+
run_rational_data_single_method(
243+
/*is_edge_edge*/ true, /*is_simulation*/ true, ms, tolerance);
244+
std::cout << "finish simulation Edge-Edge data" << std::endl;
245+
246+
std::cout << "\nRunning handcrafted Vertex-Face data" << std::endl;
232247
run_rational_data_single_method(
233248
/*is_edge_edge*/ false, /*is_simulation*/ false, ms, tolerance);
234-
std::cout << "finish Vertex-Face data" << std::endl;
235-
std::cout << "\nRunning Edge-Edge data" << std::endl;
249+
std::cout << "finish handcrafted Vertex-Face data" << std::endl;
250+
std::cout << "\nRunning handcrafted Edge-Edge data" << std::endl;
236251
run_rational_data_single_method(
237252
/*is_edge_edge*/ true, /*is_simulation*/ false, ms, tolerance);
238-
std::cout << "finish Edge-Edge data" << std::endl;
253+
std::cout << "finish handcrafted Edge-Edge data" << std::endl;
254+
255+
239256
}
240257
#endif
241258

259+
//void run_dbg() {
260+
// std::string file = "D:\\vs\\collision\\interval\\Tight-Inclusion\\build\\Release\\ee1simu1.csv";
261+
// std::vector<double> tois;
262+
// std::cout << "before read" << std::endl;
263+
// Eigen::MatrixXd all_V = read_csv(file, tois);
264+
// std::cout << "readed" << std::endl;
265+
// if (all_V.rows() != tois.size() || all_V.rows() % 8 != 0) {
266+
// std::cout << "wrong sizes, " << all_V.rows() << " " << tois.size() << std::endl;
267+
// }
268+
//
269+
// int v_size = all_V.rows() / 8;
270+
// int counter = 0;
271+
// double tolerance = 1e-6;
272+
// double minimum_seperation = 0;
273+
// int max_itr = 1e6;
274+
// int total_positives = 0;
275+
// int no_zero = 0;
276+
// for (int i = 0; i < v_size; i++) {
277+
// double toi_float = tois[i * 8];
278+
// if (toi_float > 0) continue;
279+
// Eigen::Matrix<double, 8, 3> V = all_V.middleRows<8>(8 * i);
280+
// const std::array<double, 3> err = { {-1, -1, -1} };
281+
//
282+
// double toi;
283+
// const double t_max = 1;
284+
//
285+
// double output_tolerance = tolerance;
286+
//
287+
// int CCD_TYPE = 1;
288+
// bool new_result = edgeEdgeCCD_double(
289+
// V.row(0), V.row(1), V.row(2), V.row(3), V.row(4),
290+
// V.row(5), V.row(6), V.row(7), err, minimum_seperation,
291+
// toi, tolerance, t_max, max_itr, output_tolerance,
292+
// CCD_TYPE);
293+
// if (new_result) {
294+
// total_positives += 1;
295+
// }
296+
// if (toi > 0) {
297+
// no_zero += 1;
298+
// }
299+
//
300+
// }
301+
// std::cout << "total positives " << total_positives << std::endl;
302+
// std::cout << "no zero " << no_zero << std::endl;
303+
//}
304+
242305
int main(int argc, char *argv[])
243306
{
244-
// inclusion_ccd::Rational a;
307+
//run_dbg();
245308
#ifdef TIGHT_INCLUSION_RUN_EXAMPLES
246309
run_code();
247310
#else

app/read_rational_csv.cpp

Lines changed: 79 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -10,6 +10,85 @@
1010
namespace inclusion_ccd
1111
{
1212

13+
// for debug
14+
Eigen::MatrixXd read_csv(const std::string& inputFileName, std::vector<double>& results)
15+
{
16+
// be careful, there are n lines which means there are n/8 queries, but has
17+
// n results, which means results are duplicated
18+
results.clear();
19+
std::vector<std::array<double, 3>> vs;
20+
vs.clear();
21+
std::ifstream infile;
22+
infile.open(inputFileName);
23+
std::array<double, 3> v;
24+
if (!infile.is_open()) {
25+
std::cout << "Path Wrong!!!!" << std::endl;
26+
std::cout << "path, " << inputFileName << std::endl;
27+
return Eigen::MatrixXd(1, 1);
28+
}
29+
30+
int l = 0;
31+
while (infile) // there is input overload classfile
32+
{
33+
l++;
34+
std::string s;
35+
if (!getline(infile, s))
36+
break;
37+
if (l == 1) continue;// skip the first row
38+
if (s[0] != '#') {
39+
std::istringstream ss(s);
40+
std::array<std::string, 5> record;
41+
int c = 0;
42+
while (ss) {
43+
std::string line;
44+
if (!getline(ss, line, ','))
45+
break;
46+
try {
47+
48+
record[c] = line;
49+
c++;
50+
51+
}
52+
catch (const std::invalid_argument e) {
53+
std::cout << "NaN found in file " << inputFileName
54+
<< " line " << l << std::endl;
55+
e.what();
56+
}
57+
}
58+
59+
double x = std::stod(record[0]);
60+
double y = std::stod(record[1]);
61+
double z = std::stod(record[2]);
62+
63+
v[0] = x;
64+
v[1] = y;
65+
v[2] = z;
66+
vs.push_back(v);
67+
// if (record[6] != "1" && record[6] != "0") {
68+
// std::cout
69+
// << "ERROR:result position should be 1 or 0, but it is "
70+
// << record[6] << std::endl;
71+
// }
72+
73+
results.push_back(std::stod(record[4]));
74+
75+
76+
}
77+
}
78+
79+
Eigen::MatrixXd all_v(vs.size(), 3);
80+
for (int i = 0; i < vs.size(); i++) {
81+
all_v(i, 0) = vs[i][0];
82+
all_v(i, 1) = vs[i][1];
83+
all_v(i, 2) = vs[i][2];
84+
}
85+
if (!infile.eof()) {
86+
std::cerr << "Could not read file " << inputFileName << "\n";
87+
}
88+
89+
return all_v;
90+
}
91+
1392
Eigen::MatrixXd read_rational_csv(
1493
const std::string &inputFileName, std::vector<bool> &results)
1594
{

app/read_rational_csv.hpp

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -9,5 +9,6 @@ namespace inclusion_ccd
99

1010
Eigen::MatrixXd read_rational_csv(
1111
const std::string &inputFileName, std::vector<bool> &results);
12+
Eigen::MatrixXd read_csv(const std::string& inputFileName, std::vector<double>& results);
1213

1314
} // namespace inclusion_ccd

tight_inclusion/avx.cpp

Lines changed: 7 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -45,7 +45,13 @@ namespace inclusion_ccd
4545
// }
4646

4747
// calculate a*(2^b)
48-
long power(const long a, const int b) { return a << b; }
48+
long power(const long a, const int b)
49+
{
50+
// The fast bit shifting power trick only works if b is not too larger.
51+
assert(b < 8 * sizeof(long) - 1);
52+
// WARNING: Technically this can still fail with `b < 8 * sizeof(long) - 1` if `a > 1`.
53+
return a << b;
54+
}
4955
std::array<Scalar, 8> function_ee(
5056
const Scalar &a0s,
5157
const Scalar &a1s,

tight_inclusion/inclusion_ccd.cpp

Lines changed: 55 additions & 36 deletions
Original file line numberDiff line numberDiff line change
@@ -298,11 +298,12 @@ namespace inclusion_ccd
298298
// unsigned so can be larger than MAX_NO_ZERO_TOI_ITER
299299
unsigned int no_zero_toi_iter = 0;
300300
bool is_impacting, tmp_is_impacting;
301-
301+
Scalar tolerance_in = tolerance;
302+
Scalar ms_in = ms;
302303
do
303304
{
304305
Vector3d tol = compute_edge_edge_tolerance_new(
305-
a0s, a1s, b0s, b1s, a0e, a1e, b0e, b1e, tolerance);
306+
a0s, a1s, b0s, b1s, a0e, a1e, b0e, b1e, tolerance_in);
306307

307308
// this should be the error of the whole mesh
308309
std::array<Scalar, 3> err1;
@@ -335,15 +336,16 @@ namespace inclusion_ccd
335336
if (CCD_TYPE == 0)
336337
{
337338
tmp_is_impacting = interval_root_finder_double_normalCCD(
338-
tol, toi, false, err1, ms, a0s, a1s, b0s, b1s, a0e, a1e,
339+
tol, toi, false, err1, ms_in, a0s, a1s, b0s, b1s, a0e, a1e,
339340
b0e, b1e);
341+
return tmp_is_impacting;
340342
}
341343
else
342344
{
343345
assert(CCD_TYPE == 1);
344346
assert(t_max >= 0 && t_max <= 1);
345347
tmp_is_impacting = interval_root_finder_double_horizontal_tree(
346-
tol, tolerance, toi, false, err1, ms, a0s, a1s, b0s, b1s,
348+
tol, tolerance_in, toi, false, err1, ms_in, a0s, a1s, b0s, b1s,
347349
a0e, a1e, b0e, b1e, t_max, max_itr, output_tolerance);
348350
}
349351
assert(!tmp_is_impacting || toi >= 0);
@@ -360,26 +362,34 @@ namespace inclusion_ccd
360362
toi = tmp_is_impacting ? toi : t_max;
361363
}
362364

363-
if (tmp_is_impacting && toi == 0)
364-
{
365-
// This modification is for CCD-filtered line-search (e.g., IPC)
366-
// we rebuild the time interval since tol is conservative:
367-
// this is the new time range
368-
t_max = std::min(tol[0] * 10, Scalar(0.1));
369-
370-
// if early terminated, use tolerance; otherwise, use smaller tolerance
371-
// althouth time resolution and tolerance is not the same thing, but decrease
372-
// tolerance will be helpful
373-
tolerance =
374-
output_tolerance > tolerance ? tolerance : 0.1 * tolerance;
375-
}
365+
// This modification is for CCD-filtered line-search (e.g., IPC)
366+
// strategies for dealing with toi = 0:
367+
// 1. shrink t_max (when reaches max_itr),
368+
// 2. shrink tolerance (when not reach max_itr and tolerance is big) or
369+
// ms (when tolerance is too small comparing with ms)
370+
if (tmp_is_impacting && toi == 0 && no_zero_toi) {
371+
372+
// meaning reaches max_itr, need to shrink the t_max to return a more accurate result to reach target tolerance.
373+
if (output_tolerance > tolerance_in) {
374+
t_max *= 0.9;
375+
}
376+
else {// meaning the given tolerance or ms is too large. need to shrink them,
377+
if (10 * tolerance_in < ms_in) {// ms is too large, shrink it
378+
ms_in *= 0.5;
379+
}
380+
else {// tolerance is too large, shrink it
381+
tolerance_in *= 0.5;
382+
}
383+
}
384+
}
385+
376386

377387
no_zero_toi_iter++;
378388

379389
// Only perform a second iteration if toi == 0.
380390
// WARNING: This option assumes the initial distance is not zero.
381391
} while (no_zero_toi && no_zero_toi_iter < MAX_NO_ZERO_TOI_ITER
382-
&& tmp_is_impacting && toi == 0 && CCD_TYPE == 1);
392+
&& tmp_is_impacting && toi == 0);
383393
assert(!no_zero_toi || !is_impacting || toi != 0);
384394

385395
return is_impacting;
@@ -408,13 +418,14 @@ namespace inclusion_ccd
408418
// unsigned so can be larger than MAX_NO_ZERO_TOI_ITER
409419
unsigned int no_zero_toi_iter = 0;
410420
bool is_impacting, tmp_is_impacting;
411-
421+
Scalar tolerance_in = tolerance;
422+
Scalar ms_in = ms;
412423
do
413424
{
414425
Vector3d tol = compute_face_vertex_tolerance_3d_new(
415426
vertex_start, face_vertex0_start, face_vertex1_start,
416427
face_vertex2_start, vertex_end, face_vertex0_end,
417-
face_vertex1_end, face_vertex2_end, tolerance);
428+
face_vertex1_end, face_vertex2_end, tolerance_in);
418429

419430
//////////////////////////////////////////////////////////
420431
// this is the error of the whole mesh
@@ -445,16 +456,17 @@ namespace inclusion_ccd
445456
if (CCD_TYPE == 0)
446457
{
447458
tmp_is_impacting = interval_root_finder_double_normalCCD(
448-
tol, toi, true, err1, ms, vertex_start, face_vertex0_start,
459+
tol, toi, true, err1, ms_in, vertex_start, face_vertex0_start,
449460
face_vertex1_start, face_vertex2_start, vertex_end,
450461
face_vertex0_end, face_vertex1_end, face_vertex2_end);
462+
return tmp_is_impacting;
451463
}
452464
else
453465
{
454466
assert(CCD_TYPE == 1);
455467
assert(t_max >= 0 && t_max <= 1);
456468
tmp_is_impacting = interval_root_finder_double_horizontal_tree(
457-
tol, tolerance, toi, true, err1, ms, vertex_start,
469+
tol, tolerance_in, toi, true, err1, ms_in, vertex_start,
458470
face_vertex0_start, face_vertex1_start, face_vertex2_start,
459471
vertex_end, face_vertex0_end, face_vertex1_end,
460472
face_vertex2_end, t_max, max_itr, output_tolerance);
@@ -473,26 +485,33 @@ namespace inclusion_ccd
473485
toi = tmp_is_impacting ? toi : t_max;
474486
}
475487

476-
if (tmp_is_impacting && toi == 0)
477-
{
478-
// This modification is for CCD-filtered line-search (e.g., IPC)
479-
// we rebuild the time interval since tol is conservative:
480-
// this is the new time range
481-
t_max = std::min(tol[0] * 10, Scalar(0.1));
482-
483-
// if early terminated, use tolerance; otherwise, use smaller tolerance
484-
// althouth time resolution and tolerance is not the same thing, but decrease
485-
// tolerance will be helpful
486-
tolerance =
487-
output_tolerance > tolerance ? tolerance : 0.1 * tolerance;
488-
}
488+
// This modification is for CCD-filtered line-search (e.g., IPC)
489+
// strategies for dealing with toi = 0:
490+
// 1. shrink t_max (when reaches max_itr),
491+
// 2. shrink tolerance (when not reach max_itr and tolerance is big) or
492+
// ms (when tolerance is too small comparing with ms)
493+
if (tmp_is_impacting && toi == 0&& no_zero_toi) {
494+
495+
// meaning reaches max_itr, need to shrink the t_max to return a more accurate result to reach target tolerance.
496+
if (output_tolerance > tolerance_in) {
497+
t_max *= 0.9;
498+
}
499+
else {// meaning the given tolerance or ms is too large. need to shrink them,
500+
if (10 * tolerance_in < ms_in) {// ms is too large, shrink it
501+
ms_in *= 0.5;
502+
}
503+
else {// tolerance is too large, shrink it
504+
tolerance_in *= 0.5;
505+
}
506+
}
507+
}
489508

490509
no_zero_toi_iter++;
491510

492511
// Only perform a second iteration if toi == 0.
493512
// WARNING: This option assumes the initial distance is not zero.
494513
} while (no_zero_toi && no_zero_toi_iter < MAX_NO_ZERO_TOI_ITER
495-
&& tmp_is_impacting && toi == 0 && CCD_TYPE == 1);
514+
&& tmp_is_impacting && toi == 0);
496515
assert(!no_zero_toi || !is_impacting || toi != 0);
497516

498517
return is_impacting;

0 commit comments

Comments
 (0)