Skip to content
GitLab
Menu
Projects
Groups
Snippets
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in
Toggle navigation
Menu
Open sidebar
iwr
meshconv
Commits
3d6dc36b
Commit
3d6dc36b
authored
Oct 10, 2019
by
Stenger, Florian
Browse files
meshconv v3.16
parent
6425d4a3
Changes
11
Expand all
Hide whitespace changes
Inline
Side-by-side
examples/decimate_surface.mcs
View file @
3d6dc36b
...
...
@@ -11,6 +11,10 @@ if[{target_elements} undefined] {target_elements} = 60000 end
if[{feature_angle} undefined] {feature_angle} = 140.0 end
//minimum inner angle of triangles in the output mesh:
if[{inner_angle} undefined] {inner_angle} = 25.0 end
//factor controlling the impact of curvature on mesh-resolution:
if[{curvature_factor} undefined] {curvature_factor} = 1 end
//cutoff-factor for curvature-values:
if[{cutoff_factor} undefined] {cutoff_factor} = 0.1 end
{input_file} = examples/data/ch_surface.vtu
{output_file} = ch_surface_decimated.vtu
...
...
@@ -29,7 +33,7 @@ if[{split_elements} > {target_elements}]
split_edges_element_count[m, {split_elements} {vectors}]
end
decimate_element_count[m, {target_elements}, {feature_angle} {vectors}]
decimate_element_count[m, {target_elements}, {feature_angle}
, {curvature_factor}, {cutoff_factor}
{vectors}]
delete_degenerate_elements[m]
delete_duplicate_elements[m]
...
...
manual.pdf
View file @
3d6dc36b
No preview for this file type
src/backgrid.cc
View file @
3d6dc36b
This diff is collapsed.
Click to expand it.
src/commands.cc
View file @
3d6dc36b
...
...
@@ -192,6 +192,32 @@ void call_arh_type(vector<string> &args, bool test){
}
}
//--------------------< call_average_curvature_vector >
void
call_average_curvature_vector
(
vector
<
string
>
&
args
,
bool
test
){
string
task_description
=
"calculating average-curvature-vector"
;
int
prompt_length
=
task_description
.
size
()
+
5
;
if
(
!
test
){
if
(
verbosity
>=
1
){
cout
<<
task_description
<<
fill_dots
(
prompt_length
,
task_description_length
)
<<
flush
;
}
if
(
undefined_variable_warning
!=
""
)
cout
<<
undefined_variable_warning
<<
flush
;
check_argument_count
(
args
,
2
,
3
);
Mesh
*
m
=
lookup_mesh
(
args
[
0
]);
if
(
m
->
dim_of_elements
!=
2
)
error
(
"average_curvature_vector supports triangles only."
);
DataVector
*
vec
=
lookup_vector
(
args
[
1
],
/*may_create=*/
true
);
int
levels
=
3
;
if
(
args
.
size
()
==
3
)
levels
=
ival
(
args
[
2
]);
if
(
levels
<
1
)
error
(
"level-argument must be positive."
);
average_curvature_vector
(
*
m
,
*
vec
,
levels
);
if
(
verbosity
>=
1
)
cout
<<
"ok"
<<
endl
;
}
else
{
if
(
task_description_length
<
prompt_length
)
task_description_length
=
prompt_length
;
}
}
//--------------------< call_backgrid_search_distance >
void
call_backgrid_search_distance
(
vector
<
string
>
&
args
,
bool
test
){
if
(
!
test
){
...
...
@@ -865,8 +891,24 @@ void call_decimate_edge_length(vector<string> &args, bool test){
double
min_edge_length
=
dval
(
args
[
1
]);
double
max_feature_angle
=
0.0
;
if
(
args
.
size
()
>=
3
)
max_feature_angle
=
dval
(
args
[
2
]);
int
vec_stop
=
args
.
size
()
-
1
;
int
vec_stop
=
args
.
size
()
-
1
,
vec_start
=
3
;
double
curvature_factor
=
1.0
;
double
cutoff_factor
=
0.1
;
int
min_feature_net_size
=
10
;
if
(
args
.
size
()
>=
4
){
char
curv_factor_first
=
args
[
3
][
0
];
if
(
curv_factor_first
>=
'0'
&&
curv_factor_first
<=
'9'
){
++
vec_start
;
curvature_factor
=
dval
(
args
[
3
]);
}
}
if
(
args
.
size
()
>=
5
){
char
cutoff_factor_first
=
args
[
4
][
0
];
if
(
cutoff_factor_first
>=
'0'
&&
cutoff_factor_first
<=
'9'
){
++
vec_start
;
cutoff_factor
=
dval
(
args
[
4
]);
}
}
if
(
args
.
size
()
>=
4
){
char
last_first
=
args
[
args
.
size
()
-
1
][
0
];
if
(
last_first
>=
'0'
&&
last_first
<=
'9'
){
...
...
@@ -875,14 +917,14 @@ void call_decimate_edge_length(vector<string> &args, bool test){
}
}
vector
<
DataVector
*>
vecs
,
el_vecs
;
for
(
int
i
=
3
;
i
<=
vec_stop
;
++
i
){
for
(
int
i
=
vec_start
;
i
<=
vec_stop
;
++
i
){
if
(
args
[
i
][
0
]
==
'^'
)
el_vecs
.
push_back
(
lookup_vector
(
args
[
i
]));
else
vecs
.
push_back
(
lookup_vector
(
args
[
i
]));
}
if
(
m
->
dim_of_elements
==
2
){
decimate_3d
(
*
m
,
min_edge_length
,
INT_MAX
,
0.0
,
max_feature_angle
,
min_feature_net_size
,
vecs
,
el_vecs
);
curvature_factor
,
cutoff_factor
,
vecs
,
el_vecs
);
}
else
if
(
m
->
dim_of_world
==
2
||
m
->
dim_of_elements
==
1
){
decimate_2d
(
*
m
,
min_edge_length
,
INT_MAX
);
}
else
if
(
m
->
dim_of_world
==
3
||
m
->
dim_of_elements
==
1
){
...
...
@@ -909,7 +951,7 @@ void call_decimate_inner_angle(vector<string> &args, bool test){
if
(
undefined_variable_warning
!=
""
)
cout
<<
undefined_variable_warning
<<
flush
;
check_argument_count_at_least
(
args
,
2
);
Mesh
*
m
=
lookup_mesh
(
args
[
0
]);
double
m
ax
_inner_angle
=
dval
(
args
[
1
]);
double
m
in
_inner_angle
=
dval
(
args
[
1
]);
double
max_feature_angle
=
0.0
;
if
(
args
.
size
()
>=
3
)
max_feature_angle
=
dval
(
args
[
2
]);
int
vec_stop
=
args
.
size
()
-
1
;
...
...
@@ -928,10 +970,10 @@ void call_decimate_inner_angle(vector<string> &args, bool test){
}
if
(
m
->
dim_of_elements
==
2
){
decimate_3d
(
*
m
,
0.0
,
INT_MAX
,
m
ax
_inner_angle
,
max_feature_angle
,
min_feature_net_size
,
vecs
,
el_vecs
);
decimate_3d
(
*
m
,
0.0
,
INT_MAX
,
m
in
_inner_angle
,
max_feature_angle
,
min_feature_net_size
,
/*curvature_factor=*/
1.0
,
/*cutoff_factor=*/
0.1
,
vecs
,
el_vecs
);
}
else
{
error
(
"decimate_inner_angle is only implemented for triang
ular
meshes."
);
error
(
"decimate_inner_angle is only implemented for triang
le-
meshes."
);
}
if
(
verbosity
>=
1
)
cout
<<
"ok"
<<
endl
;
...
...
@@ -955,8 +997,24 @@ void call_decimate_element_count(vector<string> &args, bool test){
int
max_elements
=
ival
(
args
[
1
]);
double
max_feature_angle
=
0.0
;
if
(
args
.
size
()
>=
3
)
max_feature_angle
=
dval
(
args
[
2
]);
int
vec_stop
=
args
.
size
()
-
1
;
int
vec_stop
=
args
.
size
()
-
1
,
vec_start
=
3
;
double
curvature_factor
=
1.0
;
double
cutoff_factor
=
0.1
;
int
min_feature_net_size
=
10
;
if
(
args
.
size
()
>=
4
){
char
curv_factor_first
=
args
[
3
][
0
];
if
(
curv_factor_first
>=
'0'
&&
curv_factor_first
<=
'9'
){
++
vec_start
;
curvature_factor
=
dval
(
args
[
3
]);
}
}
if
(
args
.
size
()
>=
5
){
char
cutoff_factor_first
=
args
[
4
][
0
];
if
(
cutoff_factor_first
>=
'0'
&&
cutoff_factor_first
<=
'9'
){
++
vec_start
;
cutoff_factor
=
dval
(
args
[
4
]);
}
}
if
(
args
.
size
()
>=
4
){
char
last_first
=
args
[
args
.
size
()
-
1
][
0
];
if
(
last_first
>=
'0'
&&
last_first
<=
'9'
){
...
...
@@ -965,14 +1023,14 @@ void call_decimate_element_count(vector<string> &args, bool test){
}
}
vector
<
DataVector
*>
vecs
,
el_vecs
;
for
(
int
i
=
3
;
i
<=
vec_stop
;
++
i
){
for
(
int
i
=
vec_start
;
i
<=
vec_stop
;
++
i
){
if
(
args
[
i
][
0
]
==
'^'
)
el_vecs
.
push_back
(
lookup_vector
(
args
[
i
]));
else
vecs
.
push_back
(
lookup_vector
(
args
[
i
]));
}
if
(
m
->
dim_of_elements
==
2
){
decimate_3d
(
*
m
,
-
1.0
,
max_elements
,
0.0
,
max_feature_angle
,
min_feature_net_size
,
vecs
,
el_vecs
);
curvature_factor
,
cutoff_factor
,
vecs
,
el_vecs
);
}
else
if
(
m
->
dim_of_world
==
2
||
m
->
dim_of_elements
==
1
){
decimate_2d
(
*
m
,
-
1.0
,
max_elements
);
}
else
if
(
m
->
dim_of_world
==
3
||
m
->
dim_of_elements
==
1
){
...
...
@@ -1428,7 +1486,7 @@ void call_eliminate_flat_triangles(vector<string> &args, bool test){
if
(
undefined_variable_warning
!=
""
)
cout
<<
undefined_variable_warning
<<
flush
;
check_argument_count_at_least
(
args
,
1
);
Mesh
*
m
=
lookup_mesh
(
args
[
0
]);
if
(
m
->
dim_of_elements
!=
2
)
error
(
"Input has to be a triang
ular
mesh."
);
if
(
m
->
dim_of_elements
!=
2
)
error
(
"Input has to be a triang
le-
mesh."
);
double
angle_threshold
=
0.1
;
int
vec_stop
=
args
.
size
()
-
1
;
if
(
args
.
size
()
>
1
){
...
...
@@ -2203,7 +2261,7 @@ void call_leole(vector<string> &args, bool test){
if
(
undefined_variable_warning
!=
""
)
cout
<<
undefined_variable_warning
<<
flush
;
check_argument_count
(
args
,
1
);
m
=
lookup_mesh
(
args
[
0
]);
if
(
m
->
dim_of_elements
!=
2
)
error
(
"Only triang
ular
elements supported."
);
if
(
m
->
dim_of_elements
!=
2
)
error
(
"Only triang
le-
elements supported."
);
delete_refinement_manager
(
args
[
0
]);
i
=
break_elements_6
(
*
m
,
m_be6
);
...
...
@@ -2233,7 +2291,7 @@ void call_levelset(vector<string> &args, bool test){
check_argument_count
(
args
,
3
);
Mesh
m_levelset
(
"$TEMPORARY"
);
Mesh
*
m
=
lookup_mesh
(
args
[
0
]);
if
(
m
->
dim_of_elements
!=
3
)
error
(
"Only tetrahedra
l
elements supported."
);
if
(
m
->
dim_of_elements
!=
3
)
error
(
"Only tetrahedra
-
elements supported."
);
delete_refinement_manager
(
args
[
0
]);
DataVector
*
vec
=
lookup_vector
(
args
[
1
]);
double
threshold1
=
dval
(
args
[
2
]);
...
...
@@ -4973,6 +5031,8 @@ void call_write(vector<string> &args, bool test){
if
(
suffix
==
"3d"
||
suffix
==
"2d"
){
write_amdis_file
(
args
[
0
],
*
m
);
}
else
if
(
suffix
==
"msh"
){
write_msh_file
(
args
[
0
],
*
m
);
}
else
if
(
suffix
==
"arh"
){
if
(
arh_type
==
arh1_parallel
){
#ifdef ARH1_SUPPORT
...
...
@@ -5322,6 +5382,7 @@ void execute_next_command(string &cmdline, bool test){
else
if
(
srcmd
==
"arh_info"
)
call_arh_info
(
cmdargs
,
test
);
else
if
(
srcmd
==
"arh_type"
)
call_arh_type
(
cmdargs
,
test
);
else
if
(
srcmd
==
"assign"
)
call_assign
(
cmdargs
,
test
);
else
if
(
srcmd
==
"average_curvature_vector"
)
call_average_curvature_vector
(
cmdargs
,
test
);
else
if
(
srcmd
==
"backgrid_search_distance"
)
call_backgrid_search_distance
(
cmdargs
,
test
);
else
if
(
srcmd
==
"backgrid_segmentation"
)
call_backgrid_segmentation
(
cmdargs
,
test
);
else
if
(
srcmd
==
"blur_vector"
)
call_blur_vector
(
cmdargs
,
test
);
...
...
src/filereaders.cc
View file @
3d6dc36b
...
...
@@ -360,7 +360,7 @@ struct VtuValueDecompressor : public VtuValueReader{
};
//--------------------< read_vtu_file >
//triang
ular
meshes in 2d- and 3d-worlds and tetrahedra
l
meshes in 3d worlds supported at the moment
//triang
le-
meshes in 2d- and 3d-worlds and tetrahedra
-
meshes in 3d worlds supported at the moment
void
read_vtu_file
(
string
file
,
Mesh
&
m
,
vector
<
DataVector
*>
&
vecs
,
vector
<
DataVector
*>
&
element_vecs
,
int
&
dim_of_world
,
bool
read_mesh
,
bool
add_mesh
,
bool
really_3d
){
...
...
@@ -1071,7 +1071,7 @@ void read_pvtu_file(string file, string mp, vector<string> &vps, bool unload, in
}
//--------------------< read_amdis_file (from file) >
//triang
ular
meshes in 2d- or 3d-worlds and tetrahedra
l
meshes in 3d worlds supported at the moment.
//triang
le-
meshes in 2d- or 3d-worlds and tetrahedra
-
meshes in 3d worlds supported at the moment.
void
read_amdis_file
(
string
file
,
Mesh
&
m
,
bool
add_mesh
,
uint64_t
byte_offset
){
//open input file
ifstream
in
(
file
.
c_str
(),
ios_base
::
in
|
ios_base
::
binary
);
...
...
@@ -1083,7 +1083,7 @@ void read_amdis_file(string file, Mesh &m, bool add_mesh, uint64_t byte_offset){
}
//--------------------< read_amdis_file (from ifstream) >
//triang
ular
meshes in 2d- or 3d-worlds and tetrahedra
l
meshes in 3d worlds supported at the moment.
//triang
le-
meshes in 2d- or 3d-worlds and tetrahedra
-
meshes in 3d worlds supported at the moment.
void
read_amdis_file
(
ifstream
&
in
,
Mesh
&
m
,
bool
add_mesh
,
uint64_t
byte_offset
){
string
buf
;
int
tmp_numelements
=
0
,
vertex_offset
=
0
,
boundary_position
=
-
1
;
...
...
@@ -1210,7 +1210,7 @@ void read_amdis_file(ifstream &in, Mesh &m, bool add_mesh, uint64_t byte_offset)
}
//--------------------< read_vrml_file >
//triang
ular
meshes in 3d-worlds supported at the moment. This seems to be a very
//triang
le-
meshes in 3d-worlds supported at the moment. This seems to be a very
//special vrml-kind...
void
read_vrml_file
(
string
file
,
Mesh
&
m
){
string
buf
;
...
...
@@ -1486,7 +1486,7 @@ void read_wsxm_file(string file, Mesh &m, bool mirror, double bottom){
}
//--------------------< read_dae_file >
//only triang
ular
meshes in 3d-worlds supported at the moment
//only triang
le-
meshes in 3d-worlds supported at the moment
void
read_dae_file
(
string
file
,
Mesh
&
m
){
string
buf
[
3
],
current_component
=
""
,
tmp
;
int
i
,
p
,
v
[
3
];
...
...
@@ -1639,7 +1639,7 @@ void read_dae_file(string file, Mesh &m){
}
//--------------------< read_rwx_file >
//only triang
ular
meshes in 3d-worlds supported at the moment
//only triang
le-
meshes in 3d-worlds supported at the moment
void
read_rwx_file
(
string
file
,
Mesh
&
m
,
int
skip_objects
){
string
buf
[
3
];
int
i
=
0
,
p
,
v
[
3
],
skips
=
0
;
...
...
@@ -1716,7 +1716,7 @@ void read_rwx_file(string file, Mesh &m, int skip_objects){
}
//--------------------< read_off_file >
//only triang
ular
meshes in 2d- or 3d-worlds supported at the moment
//only triang
le-
meshes in 2d- or 3d-worlds supported at the moment
void
read_off_file
(
string
file
,
Mesh
&
m
,
int
dim_of_world
){
string
buf
[
4
];
int
i
,
p
,
v
[
4
],
tmp_numelements
;
...
...
@@ -1770,7 +1770,7 @@ void read_off_file(string file, Mesh &m, int dim_of_world){
}
//--------------------< read_cdt_file >
//only triang
ular
meshes in 2d- or 3d-worlds supported at the moment
//only triang
le-
meshes in 2d- or 3d-worlds supported at the moment
void
read_cdt_file
(
string
file
,
Mesh
&
m
,
int
dim_of_world
){
string
buf
[
3
];
int
i
,
p
,
v
[
3
],
tmp_numelements
,
tmp2_numelements
;
...
...
@@ -1821,7 +1821,7 @@ void read_cdt_file(string file, Mesh &m, int dim_of_world){
}
//--------------------< read_msh_file >
//only triang
ular
meshes in 2d- or 3d-worlds supported at the moment
//only triang
le-
meshes in 2d- or 3d-worlds supported at the moment
void
read_msh_file
(
string
file
,
Mesh
&
m
,
int
dim_of_world
,
bool
store_boundary_field
){
string
buf
[
4
];
int
p
,
t
,
v
[
4
],
tmp_numelements
=
0
,
tmp2_numelements
=
0
;
...
...
src/filewriters.cc
View file @
3d6dc36b
...
...
@@ -1062,17 +1062,14 @@ void write_pvtu_file(string file, MeshPartitioning &mp, const vector<VectorParti
//--------------------< write_matlab_file >
void
write_matlab_file
(
string
file
,
Mesh
&
m
,
const
vector
<
DataVector
*>
&
vecs
){
ofstream
out
;
int
i
,
p
;
//prepare output file
o
ut
.
open
(
file
.
c_str
(),
ios_base
::
out
|
ios_base
::
trunc
|
ios_base
::
binary
);
o
fstream
out
(
file
.
c_str
(),
ios_base
::
out
|
ios_base
::
trunc
|
ios_base
::
binary
);
if
(
!
out
.
good
())
error
(
"Could not write file
\"
"
+
file
+
"
\"
."
);
out
.
precision
(
output_precision
);
//iterate through vertices
for
(
p
=
0
;
p
<
m
.
dim_of_world
;
++
p
){
for
(
i
=
0
;
i
<
m
.
numvertices
;
++
i
){
for
(
int
p
=
0
;
p
<
m
.
dim_of_world
;
++
p
){
for
(
i
nt
i
=
0
;
i
<
m
.
numvertices
;
++
i
){
#ifdef DEBUG_REDUCED_PRECISION_IN_FILEWRITERS
out
<<
reduce_precision
(
m
.
vertices
[
i
][
p
])
<<
" "
;
#else
...
...
@@ -1081,8 +1078,8 @@ void write_matlab_file(string file, Mesh &m, const vector<DataVector*> &vecs){
}
out
<<
"
\n
"
;
}
for
(
p
=
0
;
p
<
int
(
vecs
.
size
());
++
p
){
for
(
i
=
0
;
i
<
m
.
numvertices
;
++
i
){
for
(
int
p
=
0
;
p
<
int
(
vecs
.
size
());
++
p
){
for
(
i
nt
i
=
0
;
i
<
m
.
numvertices
;
++
i
){
#ifdef DEBUG_REDUCED_PRECISION_IN_FILEWRITERS
out
<<
reduce_precision
((
*
vecs
[
p
])[
i
])
<<
" "
;
#else
...
...
@@ -1095,6 +1092,84 @@ void write_matlab_file(string file, Mesh &m, const vector<DataVector*> &vecs){
out
.
close
();
}
//--------------------< write_msh_file >
//only line-segment meshes in 2D Worlds or triangle-meshes in 3D worlds supported.
void
write_msh_file
(
string
file
,
Mesh
&
m
){
if
(
m
.
dim_of_elements
!=
1
&&
m
.
dim_of_elements
!=
2
){
error
(
"Only line-segment-meshes and triangle-meshes supported."
);
}
//prepare output file
ofstream
out
(
file
.
c_str
(),
ios_base
::
out
|
ios_base
::
trunc
|
ios_base
::
binary
);
if
(
!
out
.
good
())
error
(
"Could not open file for writing."
);
out
.
precision
(
output_precision
);
if
(
!
m
.
boundary_checked
)
determine_boundary
(
m
);
//format
out
<<
"$MeshFormat
\n
"
;
out
<<
"2.2 0 8
\n
"
;
out
<<
"$EndMeshFormat
\n
"
;
//vertices
out
<<
"$Nodes
\n
"
;
out
<<
m
.
numvertices
<<
"
\n
"
;
for
(
int
i
=
0
;
i
<
m
.
numvertices
;
++
i
){
Vertex
&
v
=
m
.
vertices
[
i
];
out
<<
i
<<
" "
<<
v
[
0
]
<<
" "
<<
v
[
1
]
<<
" "
<<
v
[
2
]
<<
"
\n
"
;
}
out
<<
"$EndNodes
\n
"
;
//collect edges with boundary-information (for triangles only)
map
<
pair
<
int
,
int
>
,
int
>
edges
;
if
(
m
.
dim_of_elements
==
2
){
for
(
list
<
Element
>::
iterator
it
=
m
.
elements
.
begin
(),
it_end
=
m
.
elements
.
end
();
it
!=
it_end
;
++
it
){
Element
&
e
=
*
it
;
for
(
int
b
=
0
;
b
<
3
;
++
b
){
if
(
e
.
b
[
b
]
<
0
){
pair
<
int
,
int
>
edge
(
e
.
v
[(
b
+
1
)
%
3
],
e
.
v
[(
b
+
2
)
%
3
]);
if
(
edge
.
first
>
edge
.
second
){
int
tmp
=
edge
.
first
;
edge
.
first
=
edge
.
second
;
edge
.
second
=
tmp
;
}
edges
.
insert
(
pair
<
pair
<
int
,
int
>
,
int
>
(
edge
,
e
.
b
[
b
]));
}
}
}
}
//elements
out
<<
"$Elements
\n
"
;
out
<<
m
.
numelements
+
edges
.
size
()
<<
"
\n
"
;
int
i
=
0
;
if
(
m
.
dim_of_elements
==
1
){
//line-segments
for
(
list
<
Element
>::
iterator
it
=
m
.
elements
.
begin
(),
it_end
=
m
.
elements
.
end
();
it
!=
it_end
;
++
it
){
Element
&
e
=
*
it
;
out
<<
++
i
<<
" 1 0 "
<<
e
.
v
[
0
]
<<
" "
<<
e
.
v
[
1
]
<<
"
\n
"
;
}
}
else
if
(
m
.
dim_of_elements
==
2
){
//triangles
//boundary-edges
for
(
map
<
pair
<
int
,
int
>
,
int
>::
iterator
it
=
edges
.
begin
(),
it_end
=
edges
.
end
();
it
!=
it_end
;
++
it
){
pair
<
pair
<
int
,
int
>
,
int
>
edge
=
*
it
;
out
<<
++
i
<<
" 1 1 "
<<
wafb
(
edge
.
second
)
<<
" "
<<
edge
.
first
.
first
<<
" "
<<
edge
.
first
.
second
<<
"
\n
"
;
}
//triangles
for
(
list
<
Element
>::
iterator
it
=
m
.
elements
.
begin
(),
it_end
=
m
.
elements
.
end
();
it
!=
it_end
;
++
it
){
Element
&
e
=
*
it
;
out
<<
++
i
<<
" 2 0 "
<<
e
.
v
[
0
]
<<
" "
<<
e
.
v
[
1
]
<<
" "
<<
e
.
v
[
2
]
<<
"
\n
"
;
}
}
out
<<
"$EndElements
\n
"
;
out
.
close
();
}
//--------------------< write_geo_file_2d >
//1D-line-segment-meshes in 2D worlds only!
void
write_geo_file_2d
(
string
file
,
Mesh
&
m
,
DataVector
*
vec
,
int
step
){
...
...
@@ -1305,7 +1380,7 @@ void write_geo_file_3d(string file, Mesh &m, DataVector* vec){
vector
<
int
>
cindex_map
;
if
(
m
.
dim_of_world
!=
3
||
m
.
dim_of_elements
!=
2
){
error
(
"Only triang
ular
meshes in 3D worlds supported."
);
error
(
"Only triang
le-
meshes in 3D worlds supported."
);
}
line_loop_offset
=
m
.
numvertices
+
3
*
m
.
numelements
;
...
...
@@ -1456,7 +1531,7 @@ void write_geo_file_3d(string file, Mesh &m, DataVector* vec){
}
//--------------------< write_pos_file >
//triang
ular
and tetrahedra
l
meshes only!
//triang
le-
and tetrahedra
-
meshes only!
void
write_pos_file
(
string
file
,
Mesh
&
m
,
DataVector
&
vec
){
ofstream
out
(
file
.
c_str
(),
ios_base
::
out
|
ios_base
::
trunc
|
ios_base
::
binary
);
if
(
!
out
.
good
())
error
(
"Could not write file
\"
"
+
file
+
"
\"
."
);
...
...
@@ -3162,9 +3237,9 @@ void write_periodic_file(string file, Mesh &m, vector<vector<vector<int> > > &as
for
(
a
=
1
;
a
<=
m
.
dim_of_elements
;
++
a
){
out
<<
ft
(
assoc
[
i
][
e
][
a
],
2
,
' '
)
<<
" "
;
}
out
<<
" "
<<
ft
(
assoc
[
i
][
e
+
1
][
0
],
8
,
' '
)
<<
" "
;
out
<<
" "
<<
ft
(
assoc
[
i
][
e
+
1
][
0
],
8
,
' '
)
<<
" "
;
for
(
a
=
1
;
a
<=
m
.
dim_of_elements
;
++
a
){
out
<<
ft
(
assoc
[
i
][
e
+
1
][
a
],
2
,
' '
)
<<
" "
;
out
<<
ft
(
assoc
[
i
][
e
+
1
][
a
],
2
,
' '
)
<<
" "
;
}
out
<<
"
\n
"
;
}
...
...
src/geometry.cc
View file @
3d6dc36b
...
...
@@ -251,7 +251,7 @@ double distance_point_line_2d(const Vertex &point, const Vertex &lin0, const Ver
//calculates the distance between a given point and the line-segment between 'lin0' and 'lin1'.
//This version also returns the intersection-point of the perpendicular with the given line-segment.
double
distance_point_line_with_intersection_2d
(
const
Vertex
&
point
,
const
Vertex
&
lin0
,
const
Vertex
&
lin1
,
Vertex
&
intersection
){
const
Vertex
&
lin1
,
Vertex
&
intersection
){
double
gx
,
gy
,
dx
,
dy
,
t
;
//the (complete) line through 'lin0' and 'lin1' is given by 'l = lin0 + t*g' with g as follows:
...
...
@@ -264,17 +264,15 @@ double distance_point_line_with_intersection_2d(const Vertex &point, const Verte
//if the orthogonal projection of 'point' onto l is not on the given line segment then one of the
//vertices 'lin0' or 'lin1' is the nearest point to 'point'
if
(
t
<
0
){
if
(
t
<
0
){
//'lin0' is nearest
intersection
[
0
]
=
lin0
[
0
];
intersection
[
1
]
=
lin0
[
1
];
}
else
if
(
t
>
1
){
}
else
if
(
t
>
1
){
//'lin1' is nearest
intersection
[
0
]
=
lin1
[
0
];
intersection
[
1
]
=
lin1
[
1
];
}
else
{
}
else
{
//orthogonal projection is nearest
intersection
[
0
]
=
lin0
[
0
]
+
t
*
gx
;
intersection
[
1
]
=
lin0
[
1
]
+
t
*
gy
;
...
...
@@ -323,9 +321,10 @@ double angle_between_vectors_3d(const Vertex &n1, const Vertex &n2){
//--------------------< normal_vector_3d >
void
normal_vector_3d
(
const
Vertex
&
tri0
,
const
Vertex
&
tri1
,
const
Vertex
&
tri2
,
Vertex
&
normal
){
Vertex
p
,
q
;
for
(
int
i
=
0
;
i
<
3
;
++
i
)
p
[
i
]
=
tri0
[
i
]
-
tri1
[
i
];
for
(
int
i
=
0
;
i
<
3
;
++
i
)
q
[
i
]
=
tri0
[
i
]
-
tri2
[
i
];
for
(
int
i
=
0
;
i
<
3
;
++
i
){
p
[
i
]
=
tri0
[
i
]
-
tri1
[
i
];
q
[
i
]
=
tri0
[
i
]
-
tri2
[
i
];
}
normal
[
0
]
=
p
[
1
]
*
q
[
2
]
-
p
[
2
]
*
q
[
1
];
normal
[
1
]
=
p
[
2
]
*
q
[
0
]
-
p
[
0
]
*
q
[
2
];
normal
[
2
]
=
p
[
0
]
*
q
[
1
]
-
p
[
1
]
*
q
[
0
];
...
...
@@ -334,24 +333,25 @@ void normal_vector_3d(const Vertex &tri0, const Vertex &tri1, const Vertex &tri2
//--------------------< unit_normal_vector_3d >
void
unit_normal_vector_3d
(
const
Vertex
&
tri0
,
const
Vertex
&
tri1
,
const
Vertex
&
tri2
,
Vertex
&
normal
){
double
b
;
Vertex
p
,
q
;
for
(
int
i
=
0
;
i
<
3
;
++
i
)
p
[
i
]
=
tri0
[
i
]
-
tri1
[
i
];
for
(
int
i
=
0
;
i
<
3
;
++
i
)
q
[
i
]
=
tri0
[
i
]
-
tri2
[
i
];
for
(
int
i
=
0
;
i
<
3
;
++
i
){
p
[
i
]
=
tri0
[
i
]
-
tri1
[
i
];
q
[
i
]
=
tri0
[
i
]
-
tri2
[
i
];
}
normal
[
0
]
=
p
[
1
]
*
q
[
2
]
-
p
[
2
]
*
q
[
1
];
normal
[
1
]
=
p
[
2
]
*
q
[
0
]
-
p
[
0
]
*
q
[
2
];
normal
[
2
]
=
p
[
0
]
*
q
[
1
]
-
p
[
1
]
*
q
[
0
];
b
=
sqrt
(
normal
[
0
]
*
normal
[
0
]
+
normal
[
1
]
*
normal
[
1
]
+
normal
[
2
]
*
normal
[
2
]);
double
b
=
sqrt
(
normal
[
0
]
*
normal
[
0
]
+
normal
[
1
]
*
normal
[
1
]
+
normal
[
2
]
*
normal
[
2
]);
for
(
int
i
=
0
;
i
<
3
;
++
i
)
normal
[
i
]
/=
b
;
}
//--------------------< triangle_area_3d >
double
triangle_area_3d
(
const
Vertex
&
tri0
,
const
Vertex
&
tri1
,
const
Vertex
&
tri2
){
Vertex
p
,
q
,
n
;
for
(
int
i
=
0
;
i
<
3
;
++
i
)
p
[
i
]
=
tri0
[
i
]
-
tri1
[
i
];
for
(
int
i
=
0
;
i
<
3
;
++
i
)
q
[
i
]
=
tri0
[
i
]
-
tri2
[
i
];
for
(
int
i
=
0
;
i
<
3
;
++
i
){
p
[
i
]
=
tri0
[
i
]
-
tri1
[
i
];
q
[
i
]
=
tri0
[
i
]
-
tri2
[
i
];
}
n
[
0
]
=
p
[
1
]
*
q
[
2
]
-
p
[
2
]
*
q
[
1
];
n
[
1
]
=
p
[
2
]
*
q
[
0
]
-
p
[
0
]
*
q
[
2
];
n
[
2
]
=
p
[
0
]
*
q
[
1
]
-
p
[
1
]
*
q
[
0
];
...
...
@@ -636,8 +636,8 @@ bool intersection_line_line_with_intersection_3d(const Vertex &p1, const Vertex
}
if
(
ci
!=
-
1
){
//determine best intersection-point
//priority is p3
(*)
> p4 > p3 > p2 (one of these three points must be part of both line
segments
//
(*):
p3 is only priorized over p4 if it is not a vertex of the first line segment
//priority is p3 > p4 > p3 > p2 (one of these three points must be part of both line
//
segments.
p3 is only priorized over p4 if it is not a vertex of the first line segment
if
(
mua
>
epsilon
&&
mua
-
1
<
-
epsilon
){
intersection
[
0
]
=
p3
[
0
];
intersection
[
1
]
=
p3
[
1
];
...
...
@@ -1049,7 +1049,7 @@ double distance_point_line_3d(const Vertex &point, const Vertex &lin0, const Ver
dy
=
point
[
1
]
-
fy
;
dz
=
point
[
2
]
-
fz
;
return
sqrt
(
dx
*
dx
+
dy
*
dy
+
dz
*
dz
);
return
sqrt
(
dx
*
dx
+
dy
*
dy
+
dz
*
dz
);
}
//--------------------< distance_point_triangle_3d >
...
...
@@ -1090,22 +1090,22 @@ double distance_point_triangle_3d(const Vertex &point,
}
//now test whether the orthogonal projection of 'point' onto p lies inside the triangle
if
(
u
>=
0
&&
v
>=
0
&&
u
+
v
<=
1
){
if
(
u
>=
0
&&
v
>=
0
&&
u
+
v
<=
1
){
//yes, inside. The length of the normal-vector is the desired distance
fx
=
point
[
0
]
-
tri0
[
0
]
-
u
*
gx
-
v
*
hx
;
fy
=
point
[
1
]
-
tri0
[
1
]
-
u
*
gy
-
v
*
hy
;
fz
=
point
[
2
]
-
tri0
[
2
]
-
u
*
gz
-
v
*
hz
;
distance
=
sqrt
(
fx
*
fx
+
fy
*
fy
+
fz
*
fz
);
fx
=
point
[
0
]
-
tri0
[
0
]
-
u
*
gx
-
v
*
hx
;
fy
=
point
[
1
]
-
tri0
[
1
]
-
u
*
gy
-
v
*
hy
;
fz
=
point
[
2
]
-
tri0
[
2
]
-
u
*
gz
-
v
*
hz
;
distance
=
sqrt
(
fx
*
fx
+
fy
*
fy
+
fz
*
fz
);
}
else
{
//no, not inside. The desired distance is the distance to one of the triangle-edges
distance
=
9e100
;
tmp
=
distance_point_line_3d
(
point
,
tri0
,
tri1
);
if
(
tmp
<
distance
)
distance
=
tmp
;
if
(
tmp
<
distance
)
distance
=
tmp
;
tmp
=
distance_point_line_3d
(
point
,
tri1
,
tri2
);
if
(
tmp
<
distance
)
distance
=
tmp
;
if
(
tmp
<
distance
)
distance
=
tmp
;
tmp
=
distance_point_line_3d
(
point
,
tri2
,
tri0
);
if
(
tmp
<
distance
)
distance
=
tmp
;
if
(
tmp
<
distance
)
distance
=
tmp
;
}
return
distance
;
}
...
...
src/interface.cc
View file @
3d6dc36b
...
...
@@ -216,7 +216,9 @@ void import_vector(DataVector &v, string name){
//--------------------< check_argument_count >
void
check_argument_count
(
vector
<
string
>
&
args
,
int
min_count
,
int
max_count
){
if
(
max_count
==
-
1
)
max_count
=
min_count
;
if
(
int
(
args
.
size
())
<
min_count
||
int
(
args
.
size
())
>
max_count
)
error
(
"Invalid number of arguments."
);
if
(
int
(
args
.
size
())
<
min_count
||
int
(
args
.
size
())
>
max_count
){
error
(
"Invalid number of arguments."
);
}
}
//--------------------< check_argument_count_at_least >
...
...
src/mesh_operations.cc
View file @
3d6dc36b
...
...
@@ -494,7 +494,7 @@ void combine_meshes(Mesh &m1, const Mesh &m2, vector<DataVector*> &vecs){
int
num_vecs
=
vecs
.
size
();
int
num_new_vertices
=
m2
.
numvertices
;
//mark all vertices which are in
both
"m2" but not in "m1" as obsolete in "new_vertices".
//mark all vertices which are in "m2" but not in "m1" as obsolete in "new_vertices".
for
(
int
i
=
0
;
i
<
m1
.
numvertices
;
++
i
){
vector
<
pair
<
Vertex
,
int
>*>
find_result
=
new_vertices
.
find
(
m1
.
vertices
[
i
]);
if
(
!
find_result
.
empty
()){
...
...
@@ -503,7 +503,7 @@ void combine_meshes(Mesh &m1, const Mesh &m2, vector<DataVector*> &vecs){
}
}
//make room for the new vertices in "v
1
".
//make room for the new vertices in "v
ecs
".
for
(
int
e
=
0
;
e
<
num_vecs
;
++
e
){
vecs
[
e
]
->
resize
(
vecs
[
e
]
->
size
()
+
num_new_vertices
);
}
...
...
@@ -852,8 +852,8 @@ void repair_arh_vertex_order(Mesh &m, vector<DataVector*> &vecs){
}
//--------------------< levelset_3d >
//Calculates the
triangular
contour-surface-mesh of a
tetrahedral
volume-mesh at a given threshold
//value
for value
-vector 'vec'.
//Calculates the contour-surface-mesh of a volume-mesh at a given threshold
value for
//value-vector 'vec'.