Skip to content
GitLab
Explore
Sign in
Register
Primary navigation
Search or go to…
Project
gmsh
Manage
Activity
Members
Labels
Plan
Issues
Issue boards
Milestones
Wiki
Code
Merge requests
Repository
Branches
Commits
Tags
Repository graph
Compare revisions
Snippets
Build
Pipelines
Jobs
Pipeline schedules
Artifacts
Deploy
Releases
Model registry
Operate
Environments
Monitor
Incidents
Analyze
Value stream analytics
Contributor analytics
CI/CD analytics
Repository analytics
Model experiments
Help
Help
Support
GitLab documentation
Compare GitLab plans
Community forum
Contribute to GitLab
Provide feedback
Terms and privacy
Keyboard shortcuts
?
Snippets
Groups
Projects
Show more breadcrumbs
Larry Price
gmsh
Commits
b3f71e99
Commit
b3f71e99
authored
13 years ago
by
Christophe Geuzaine
Browse files
Options
Downloads
Patches
Plain Diff
fix compile msvc
parent
9886cb4c
No related branches found
No related tags found
No related merge requests found
Changes
4
Expand all
Hide whitespace changes
Inline
Side-by-side
Showing
4 changed files
Geo/GModelFactory.cpp
+8
-7
8 additions, 7 deletions
Geo/GModelFactory.cpp
Mesh/CenterlineField.cpp
+85
-85
85 additions, 85 deletions
Mesh/CenterlineField.cpp
Mesh/meshGFaceDelaunayInsertion.cpp
+236
-236
236 additions, 236 deletions
Mesh/meshGFaceDelaunayInsertion.cpp
Post/PView.cpp
+1
-1
1 addition, 1 deletion
Post/PView.cpp
with
330 additions
and
329 deletions
Geo/GModelFactory.cpp
+
8
−
7
View file @
b3f71e99
...
...
@@ -3,6 +3,7 @@
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to <gmsh@geuz.org>.
#include
<stdlib.h>
#include
"GModelFactory.h"
#include
"ListUtils.h"
#include
"Context.h"
...
...
@@ -62,7 +63,7 @@ GEdge *GeoFactory::addLine(GModel *gm, GVertex *start, GVertex *end)
GFace
*
GeoFactory
::
addPlanarFace
(
GModel
*
gm
,
std
::
vector
<
std
::
vector
<
GEdge
*>
>
edges
)
{
//create line loops
int
nLoops
=
edges
.
size
();
std
::
vector
<
EdgeLoop
*>
vecLoops
;
...
...
@@ -77,8 +78,8 @@ GFace *GeoFactory::addPlanarFace(GModel *gm, std::vector< std::vector<GEdge *> >
if
(
!
c
){
GVertex
*
gvb
=
ge
->
getBeginVertex
();
GVertex
*
gve
=
ge
->
getEndVertex
();
Vertex
*
vertb
=
FindPoint
(
(
int
)
f
abs
(
gvb
->
tag
()));
Vertex
*
verte
=
FindPoint
(
(
int
)
f
abs
(
gve
->
tag
()));
Vertex
*
vertb
=
FindPoint
(
abs
(
gvb
->
tag
()));
Vertex
*
verte
=
FindPoint
(
abs
(
gve
->
tag
()));
if
(
!
vertb
){
vertb
=
Create_Vertex
(
gvb
->
tag
(),
gvb
->
x
(),
gvb
->
y
(),
gvb
->
z
(),
gvb
->
prescribedMeshSizeAtVertex
(),
1.0
);
...
...
@@ -121,7 +122,7 @@ GFace *GeoFactory::addPlanarFace(GModel *gm, std::vector< std::vector<GEdge *> >
List_Delete
(
temp
);
}
List_Add
(
temp
,
&
numEdge
);
}
}
int
num
=
gm
->
getMaxElementaryNumber
(
2
)
+
1
+
i
;
while
(
FindSurfaceLoop
(
num
)){
...
...
@@ -132,9 +133,9 @@ GFace *GeoFactory::addPlanarFace(GModel *gm, std::vector< std::vector<GEdge *> >
EdgeLoop
*
l
=
Create_EdgeLoop
(
num
,
temp
);
vecLoops
.
push_back
(
l
);
Tree_Add
(
gm
->
getGEOInternals
()
->
EdgeLoops
,
&
l
);
List_Delete
(
temp
);
}
List_Delete
(
temp
);
}
//create surface
int
numf
=
gm
->
getMaxElementaryNumber
(
2
)
+
1
;
Surface
*
s
=
Create_Surface
(
numf
,
MSH_SURF_PLAN
);
...
...
This diff is collapsed.
Click to expand it.
Mesh/CenterlineField.cpp
+
85
−
85
View file @
b3f71e99
...
...
@@ -48,7 +48,7 @@ double computeLength(std::vector<MLine*> lines){
double
length
=
0.0
;
for
(
int
i
=
0
;
i
<
lines
.
size
();
i
++
){
length
+=
lines
[
i
]
->
getLength
();
length
+=
lines
[
i
]
->
getLength
();
}
return
length
;
}
...
...
@@ -77,7 +77,7 @@ bool isClosed(std::set<MEdge, Less_Edge> &theCut){
// }
// // find a lonely MEdge
// for (std::list<MEdge>::iterator it = segments.begin();
// for (std::list<MEdge>::iterator it = segments.begin();
// it != segments.end(); ++it){
// MVertex *vL = it->getVertex(0);
// MVertex *vR = it->getVertex(1);
...
...
@@ -89,7 +89,7 @@ bool isClosed(std::set<MEdge, Less_Edge> &theCut){
// else boundv.erase(it2);
// }
// if (boundv.size() == 0) return true;
// if (boundv.size() == 0) return true;
// else return false;
}
...
...
@@ -107,7 +107,7 @@ void orderMLines(std::vector<MLine*> &lines, MVertex *vB, MVertex *vE){
}
// find a lonely MLine
for
(
std
::
list
<
MLine
*>::
iterator
it
=
segments
.
begin
();
for
(
std
::
list
<
MLine
*>::
iterator
it
=
segments
.
begin
();
it
!=
segments
.
end
();
++
it
){
MVertex
*
vL
=
(
*
it
)
->
getVertex
(
0
);
MVertex
*
vR
=
(
*
it
)
->
getVertex
(
1
);
...
...
@@ -129,7 +129,7 @@ void orderMLines(std::vector<MLine*> &lines, MVertex *vB, MVertex *vE){
if
(
v
==
vB
)
firstLine
=
(
boundv
.
rbegin
())
->
second
;
else
{
Msg
::
Error
(
"begin vertex not found for branch"
);
exit
(
1
);}
}
for
(
std
::
list
<
MLine
*>::
iterator
it
=
segments
.
begin
();
for
(
std
::
list
<
MLine
*>::
iterator
it
=
segments
.
begin
();
it
!=
segments
.
end
();
++
it
){
if
(
*
it
==
firstLine
){
segments
.
erase
(
it
);
...
...
@@ -149,7 +149,7 @@ void orderMLines(std::vector<MLine*> &lines, MVertex *vB, MVertex *vE){
while
(
first
!=
last
){
if
(
segments
.
empty
())
break
;
bool
found
=
false
;
for
(
std
::
list
<
MLine
*>::
iterator
it
=
segments
.
begin
();
for
(
std
::
list
<
MLine
*>::
iterator
it
=
segments
.
begin
();
it
!=
segments
.
end
();
++
it
){
MLine
*
e
=
*
it
;
if
(
e
->
getVertex
(
0
)
==
last
){
...
...
@@ -185,12 +185,12 @@ void orderMLines(std::vector<MLine*> &lines, MVertex *vB, MVertex *vE){
//lines is now a list of ordered MLines
lines
=
_m
;
//special case reverse orientation
if
(
lines
.
size
()
<
2
)
return
;
if
(
_orientation
[
0
]
&&
lines
[
0
]
->
getVertex
(
1
)
!=
lines
[
1
]
->
getVertex
(
1
)
&&
lines
[
0
]
->
getVertex
(
1
)
!=
lines
[
1
]
->
getVertex
(
0
)){
for
(
unsigned
int
i
=
0
;
i
<
lines
.
size
();
i
++
)
for
(
unsigned
int
i
=
0
;
i
<
lines
.
size
();
i
++
)
_orientation
[
i
]
=
!
_orientation
[
i
];
}
...
...
@@ -198,7 +198,7 @@ void orderMLines(std::vector<MLine*> &lines, MVertex *vB, MVertex *vE){
// else if (junctions.find(lines[0]->getVertex(1)) != junctions.end()) vB = lines[0]->getVertex(1);
// else printf("no vB junc found %d %d \n", lines[0]->getVertex(0)->getNum(), lines[0]->getVertex(1)->getNum());
// if (junctions.find(lines[lines.size()-1]->getVertex(0)) != junctions.end()) vE = lines[lines.size()-1]->getVertex(0);
// else if (junctions.find(lines[lines.size()-1]->getVertex(1)) != junctions.end()) vE = lines[lines.size()-1]->getVertex(1);
// else if (junctions.find(lines[lines.size()-1]->getVertex(1)) != junctions.end()) vE = lines[lines.size()-1]->getVertex(1);
// else printf("no vE junc found %d %d \n", lines[0]->getVertex(0)->getNum(), lines[0]->getVertex(1)->getNum());
// printf("in order vB= %d =%d \n", vB->getNum(), vE->getNum());
}
...
...
@@ -206,7 +206,7 @@ void orderMLines(std::vector<MLine*> &lines, MVertex *vB, MVertex *vE){
static
void
recurConnectByMEdge
(
const
MEdge
&
e
,
std
::
multimap
<
MEdge
,
MTriangle
*
,
Less_Edge
>
&
e2e
,
std
::
set
<
MTriangle
*>
&
group
,
std
::
set
<
MEdge
,
Less_Edge
>
&
touched
,
std
::
set
<
MEdge
,
Less_Edge
>
&
touched
,
std
::
set
<
MEdge
,
Less_Edge
>
&
theCut
)
{
if
(
touched
.
find
(
e
)
!=
touched
.
end
())
return
;
...
...
@@ -216,8 +216,8 @@ static void recurConnectByMEdge(const MEdge &e,
group
.
insert
(
it
->
second
);
for
(
int
i
=
0
;
i
<
it
->
second
->
getNumEdges
();
++
i
){
MEdge
me
=
it
->
second
->
getEdge
(
i
);
if
(
theCut
.
find
(
me
)
!=
theCut
.
end
()){
touched
.
insert
(
me
);
//break;
if
(
theCut
.
find
(
me
)
!=
theCut
.
end
()){
touched
.
insert
(
me
);
//break;
}
else
recurConnectByMEdge
(
me
,
e2e
,
group
,
touched
,
theCut
);
}
...
...
@@ -225,19 +225,19 @@ static void recurConnectByMEdge(const MEdge &e,
}
void
cutTriangle
(
MTriangle
*
tri
,
std
::
map
<
MEdge
,
MVertex
*
,
Less_Edge
>
&
cutEdges
,
void
cutTriangle
(
MTriangle
*
tri
,
std
::
map
<
MEdge
,
MVertex
*
,
Less_Edge
>
&
cutEdges
,
std
::
set
<
MVertex
*>
&
cutVertices
,
std
::
vector
<
MTriangle
*>
&
newTris
,
std
::
vector
<
MTriangle
*>
&
newTris
,
std
::
set
<
MEdge
,
Less_Edge
>
&
newCut
){
MVertex
*
c
[
3
]
=
{
0
,
0
,
0
};
for
(
int
j
=
0
;
j
<
3
;
j
++
){
MEdge
ed
=
tri
->
getEdge
(
j
);
MEdge
ed
=
tri
->
getEdge
(
j
);
std
::
map
<
MEdge
,
MVertex
*
,
Less_Edge
>
::
iterator
it
=
cutEdges
.
find
(
ed
);
if
(
it
!=
cutEdges
.
end
()){
c
[
j
]
=
it
->
second
;
}
}
MVertex
*
old_v0
=
tri
->
getVertex
(
0
);
...
...
@@ -316,7 +316,7 @@ void cutTriangle(MTriangle *tri,
}
Centerline
::
Centerline
(
std
::
string
fileName
)
:
kdtree
(
0
),
kdtreeR
(
0
),
nodes
(
0
),
nodesR
(
0
){
recombine
=
CTX
::
instance
()
->
mesh
.
recombineAll
;
index
=
new
ANNidx
[
1
];
...
...
@@ -378,14 +378,14 @@ void Centerline::importFile(std::string fileName){
}
if
(
triangles
.
empty
()){
Msg
::
Error
(
"Current GModel has no triangles ..."
);
Msg
::
Error
(
"Current GModel has no triangles ..."
);
exit
(
1
);
}
mod
=
new
GModel
();
mod
->
load
(
fileName
);
mod
->
removeDuplicateMeshVertices
(
1.e-8
);
current
->
setAsCurrent
();
current
->
setAsCurrent
();
int
maxN
=
0.0
;
std
::
vector
<
GEdge
*>
modEdges
=
mod
->
bindingsGetEdges
();
...
...
@@ -408,7 +408,7 @@ void Centerline::importFile(std::string fileName){
}
createBranches
(
maxN
);
}
...
...
@@ -435,7 +435,7 @@ void Centerline::createBranches(int maxN){
junctions
.
insert
(
*
it
);
}
}
//split edges
int
tag
=
0
;
for
(
unsigned
int
i
=
0
;
i
<
color_edges
.
size
();
++
i
){
...
...
@@ -446,7 +446,7 @@ void Centerline::createBranches(int maxN){
MVertex
*
vB
=
(
*
itl
)
->
getVertex
(
0
);
MVertex
*
vE
=
(
*
itl
)
->
getVertex
(
1
);
myLines
.
push_back
(
*
itl
);
erase
(
lines
,
*
itl
);
erase
(
lines
,
*
itl
);
itl
=
lines
.
begin
();
while
(
!
(
junctions
.
find
(
vE
)
!=
junctions
.
end
()
&&
junctions
.
find
(
vB
)
!=
junctions
.
end
())
)
{
...
...
@@ -481,14 +481,14 @@ void Centerline::createBranches(int maxN){
else
itl
++
;
}
if
(
vB
==
vE
)
{
Msg
::
Error
(
"Begin and end points branch are the same
\n
"
);
break
;}
orderMLines
(
myLines
,
vB
,
vE
);
orderMLines
(
myLines
,
vB
,
vE
);
std
::
vector
<
Branch
>
children
;
Branch
newBranch
=
{
tag
++
,
myLines
,
computeLength
(
myLines
),
vB
,
vE
,
children
,
1.e6
,
0.0
};
edges
.
push_back
(
newBranch
)
;
}
}
printf
(
"*** Centerline in/outlets =%d branches =%d
\n
"
,
(
int
)
color_edges
.
size
()
+
1
,
(
int
)
edges
.
size
());
//create children
for
(
unsigned
int
i
=
0
;
i
<
edges
.
size
();
++
i
)
{
MVertex
*
vE
=
edges
[
i
].
vE
;
...
...
@@ -513,7 +513,7 @@ void Centerline::distanceToSurface(){
Msg
::
Info
(
"Centerline: computing distance to surface mesh "
);
//COMPUTE WITH REVERSE ANN TREE (SURFACE POINTS IN TREE)
std
::
set
<
MVertex
*>
allVS
;
std
::
set
<
MVertex
*>
allVS
;
for
(
int
j
=
0
;
j
<
triangles
.
size
();
j
++
)
for
(
int
k
=
0
;
k
<
3
;
k
++
)
allVS
.
insert
(
triangles
[
j
]
->
getVertex
(
k
));
int
nbSNodes
=
allVS
.
size
();
...
...
@@ -522,13 +522,13 @@ void Centerline::distanceToSurface(){
std
::
set
<
MVertex
*>::
iterator
itp
=
allVS
.
begin
();
while
(
itp
!=
allVS
.
end
()){
MVertex
*
v
=
*
itp
;
nodesR
[
ind
][
0
]
=
v
->
x
();
nodesR
[
ind
][
1
]
=
v
->
y
();
nodesR
[
ind
][
2
]
=
v
->
z
();
nodesR
[
ind
][
0
]
=
v
->
x
();
nodesR
[
ind
][
1
]
=
v
->
y
();
nodesR
[
ind
][
2
]
=
v
->
z
();
itp
++
;
ind
++
;
}
kdtreeR
=
new
ANNkd_tree
(
nodesR
,
nbSNodes
,
3
);
for
(
unsigned
int
i
=
0
;
i
<
lines
.
size
();
i
++
){
MLine
*
l
=
lines
[
i
];
MVertex
*
v1
=
l
->
getVertex
(
0
);
...
...
@@ -536,9 +536,9 @@ void Centerline::distanceToSurface(){
double
midp
[
3
]
=
{
0.5
*
(
v1
->
x
()
+
v2
->
x
()),
0.5
*
(
v1
->
y
()
+
v1
->
y
()),
0.5
*
(
v1
->
z
()
+
v2
->
z
())};
kdtreeR
->
annkSearch
(
midp
,
1
,
index
,
dist
);
double
minRad
=
sqrt
(
dist
[
0
]);
radiusl
.
insert
(
std
::
make_pair
(
lines
[
i
],
minRad
));
radiusl
.
insert
(
std
::
make_pair
(
lines
[
i
],
minRad
));
}
}
void
Centerline
::
computeRadii
(){
...
...
@@ -552,7 +552,7 @@ void Centerline::computeRadii(){
edges
[
i
].
maxRad
=
std
::
max
(
itr
->
second
,
edges
[
i
].
maxRad
);
}
else
printf
(
"ARGG line not found
\n
"
);
}
}
}
}
...
...
@@ -570,12 +570,12 @@ void Centerline::buildKdTree(){
std
::
map
<
MVertex
*
,
int
>::
iterator
itp
=
colorp
.
begin
();
while
(
itp
!=
colorp
.
end
()){
MVertex
*
v
=
itp
->
first
;
nodes
[
ind
][
0
]
=
v
->
x
();
nodes
[
ind
][
1
]
=
v
->
y
();
nodes
[
ind
][
2
]
=
v
->
z
();
nodes
[
ind
][
0
]
=
v
->
x
();
nodes
[
ind
][
1
]
=
v
->
y
();
nodes
[
ind
][
2
]
=
v
->
z
();
itp
++
;
ind
++
;
}
for
(
unsigned
int
k
=
0
;
k
<
lines
.
size
();
++
k
){
for
(
unsigned
int
k
=
0
;
k
<
lines
.
size
();
++
k
){
MVertex
*
v0
=
lines
[
k
]
->
getVertex
(
0
);
MVertex
*
v1
=
lines
[
k
]
->
getVertex
(
1
);
SVector3
P0
(
v0
->
x
(),
v0
->
y
(),
v0
->
z
());
...
...
@@ -583,9 +583,9 @@ void Centerline::buildKdTree(){
for
(
int
j
=
1
;
j
<
nbPL
+
1
;
j
++
){
double
inc
=
(
double
)
j
/
(
double
)(
nbPL
+
1
);
SVector3
Pj
=
P0
+
inc
*
(
P1
-
P0
);
nodes
[
ind
][
0
]
=
Pj
.
x
();
nodes
[
ind
][
1
]
=
Pj
.
y
();
nodes
[
ind
][
2
]
=
Pj
.
z
();
nodes
[
ind
][
0
]
=
Pj
.
x
();
nodes
[
ind
][
1
]
=
Pj
.
y
();
nodes
[
ind
][
2
]
=
Pj
.
z
();
ind
++
;
}
}
...
...
@@ -619,23 +619,23 @@ void Centerline::createSplitCompounds(){
GEdge
*
pe
=
current
->
getEdgeByTag
(
i
+
1
);
//current edge
e_compound
.
push_back
(
pe
);
int
num_gec
=
NE
+
i
+
1
;
Msg
::
Info
(
"Parametrize Compound Line (%d) = %d discrete edge"
,
Msg
::
Info
(
"Parametrize Compound Line (%d) = %d discrete edge"
,
num_gec
,
pe
->
tag
());
GEdgeCompound
*
gec
=
new
GEdgeCompound
(
current
,
num_gec
,
e_compound
);
current
->
add
(
gec
);
}
// Parametrize Compound surfaces
std
::
list
<
GEdge
*>
U0
;
for
(
int
i
=
0
;
i
<
NF
;
i
++
){
std
::
list
<
GFace
*>
f_compound
;
GFace
*
pf
=
current
->
getFaceByTag
(
i
+
1
);
//current face
f_compound
.
push_back
(
pf
);
int
num_gfc
=
NF
+
i
+
1
;
GFace
*
pf
=
current
->
getFaceByTag
(
i
+
1
);
//current face
f_compound
.
push_back
(
pf
);
int
num_gfc
=
NF
+
i
+
1
;
Msg
::
Info
(
"Parametrize Compound Surface (%d) = %d discrete face"
,
num_gfc
,
pf
->
tag
());
GFaceCompound
::
typeOfMapping
typ
=
GFaceCompound
::
HARMONICPLANE
;
//GFaceCompound::typeOfMapping typ = GFaceCompound::CONFORMAL;
GFaceCompound
::
typeOfMapping
typ
=
GFaceCompound
::
HARMONICPLANE
;
//GFaceCompound::typeOfMapping typ = GFaceCompound::CONFORMAL;
GFaceCompound
*
gfc
=
new
GFaceCompound
(
current
,
num_gfc
,
f_compound
,
U0
,
typ
,
0
);
gfc
->
meshAttributes
.
recombine
=
recombine
;
...
...
@@ -663,8 +663,8 @@ void Centerline::cleanMesh(){
Msg
::
Info
(
"Writing new splitted mesh mySPLITMESH.msh"
);
current
->
writeMSH
(
"mySPLITMESH.msh"
,
2.2
,
false
,
false
);
std
::
set
<
MVertex
*>
allNod
;
discreteFace
*
mySplitMesh
;
std
::
set
<
MVertex
*>
allNod
;
discreteFace
*
mySplitMesh
;
std
::
vector
<
std
::
set
<
MVertex
*>
>
inOutNod
;
std
::
vector
<
discreteFace
*
>
inOutMesh
;
...
...
@@ -693,7 +693,7 @@ void Centerline::cleanMesh(){
mySplitMesh
->
quadrangles
.
push_back
(
new
MQuadrangle
(
v
[
0
],
v
[
1
],
v
[
2
],
v
[
3
]));
}
}
//Removing discrete Vertices - Edges - Faces
for
(
int
i
=
0
;
i
<
NV
;
i
++
){
GVertex
*
gv
=
current
->
getVertexByTag
(
i
+
1
);
...
...
@@ -702,7 +702,7 @@ void Centerline::cleanMesh(){
for
(
int
i
=
0
;
i
<
NE
;
i
++
){
GEdge
*
ge
=
current
->
getEdgeByTag
(
i
+
1
);
GEdge
*
gec
=
current
->
getEdgeByTag
(
NE
+
i
+
1
);
current
->
remove
(
ge
);
current
->
remove
(
ge
);
current
->
remove
(
gec
);
}
for
(
int
i
=
0
;
i
<
NF
;
i
++
){
...
...
@@ -711,19 +711,19 @@ void Centerline::cleanMesh(){
}
for
(
int
i
=
0
;
i
<
discFaces
.
size
();
i
++
){
GFace
*
gfc
=
current
->
getFaceByTag
(
NF
+
i
+
1
);
current
->
remove
(
gfc
);
current
->
remove
(
gfc
);
}
//Put new mesh in a new discreteFace
for
(
std
::
set
<
MVertex
*>::
iterator
it
=
allNod
.
begin
();
it
!=
allNod
.
end
();
++
it
)
mySplitMesh
->
mesh_vertices
.
push_back
(
*
it
);
mySplitMesh
->
meshStatistics
.
status
=
GFace
::
DONE
;
mySplitMesh
->
meshStatistics
.
status
=
GFace
::
DONE
;
current
->
createTopologyFromMesh
();
}
void
Centerline
::
createFaces
(){
std
::
vector
<
std
::
vector
<
MTriangle
*>
>
faces
;
std
::
multimap
<
MEdge
,
MTriangle
*
,
Less_Edge
>
e2e
;
...
...
@@ -799,14 +799,14 @@ void Centerline::createClosedVolume(){
GEdge
*
gec
=
current
->
getEdgeByTag
(
NE
+
boundEdges
[
i
]
->
tag
());
myEdges
.
push_back
(
gec
);
myEdgeLoops
.
push_back
(
myEdges
);
GFace
*
newFace
=
current
->
addPlanarFace
(
myEdgeLoops
);
GFace
*
newFace
=
current
->
addPlanarFace
(
myEdgeLoops
);
newFace
->
addPhysicalEntity
(
200
);
current
->
setPhysicalName
(
"in/out"
,
2
,
200
);
myFaces
.
push_back
(
newFace
);
}
}
Msg
::
Info
(
"Centerline action (closeVolume) has created %d in/out planar faces "
,
(
int
)
boundEdges
.
size
());
for
(
int
i
=
0
;
i
<
NF
;
i
++
){
GFace
*
gf
=
current
->
getFaceByTag
(
NF
+
i
+
1
);
myFaces
.
push_back
(
gf
);
...
...
@@ -829,7 +829,7 @@ void Centerline::closeVolume(){
}
void
Centerline
::
cutMesh
(){
is_cut
=
true
;
if
(
update_needed
){
...
...
@@ -854,16 +854,16 @@ void Centerline::cutMesh(){
// }
Msg
::
Info
(
"Splitting surface mesh (%d tris) with centerline %s "
,
triangles
.
size
(),
fileName
.
c_str
());
//splitMesh
for
(
unsigned
int
i
=
0
;
i
<
edges
.
size
();
i
++
){
std
::
vector
<
MLine
*>
lines
=
edges
[
i
].
lines
;
double
L
=
edges
[
i
].
length
;
double
D
=
(
edges
[
i
].
minRad
+
edges
[
i
].
maxRad
);
double
AR
=
L
/
D
;
printf
(
"*** Centerline branch %d (AR=%d)
\n
"
,
i
,
(
int
)
round
(
AR
));
printf
(
"*** Centerline branch %d (AR=%d)
\n
"
,
i
,
(
int
)
floor
(
AR
+
0.5
));
if
(
AR
>
4.0
){
int
nbSplit
=
(
int
)
round
(
AR
/
3.
);
int
nbSplit
=
(
int
)
floor
(
AR
/
3.
+
0.5
);
double
li
=
L
/
nbSplit
;
double
lc
=
0.0
;
for
(
int
j
=
0
;
j
<
lines
.
size
();
j
++
){
...
...
@@ -897,7 +897,7 @@ void Centerline::cutMesh(){
createFaces
();
current
->
createTopologyFromFaces
(
discFaces
);
current
->
exportDiscreteGEOInternals
();
//write
Msg
::
Info
(
"Writing splitted mesh 'myPARTS.msh'"
);
current
->
writeMSH
(
"myPARTS.msh"
,
2.2
,
false
,
false
);
...
...
@@ -911,25 +911,25 @@ void Centerline::cutMesh(){
}
void
Centerline
::
cutByDisk
(
SVector3
&
PT
,
SVector3
&
NORM
,
double
&
maxRad
){
double
a
=
NORM
.
x
();
double
b
=
NORM
.
y
();
double
c
=
NORM
.
z
();
double
d
=
-
a
*
PT
.
x
()
-
b
*
PT
.
y
()
-
c
*
PT
.
z
();
//printf("cut disk (R=%g)= %g %g %g %g \n", maxRad, a, b, c, d);
const
double
EPS
=
0.007
;
std
::
set
<
MEdge
,
Less_Edge
>
allEdges
;
for
(
unsigned
int
i
=
0
;
i
<
triangles
.
size
();
i
++
)
for
(
unsigned
int
i
=
0
;
i
<
triangles
.
size
();
i
++
)
for
(
unsigned
int
j
=
0
;
j
<
3
;
j
++
)
allEdges
.
insert
(
triangles
[
i
]
->
getEdge
(
j
));
allEdges
.
insert
(
triangles
[
i
]
->
getEdge
(
j
));
bool
closedCut
=
false
;
int
step
=
0
;
while
(
!
closedCut
&&
step
<
10
){
double
rad
=
1.2
*
maxRad
+
0.1
*
step
*
maxRad
;
std
::
map
<
MEdge
,
MVertex
*
,
Less_Edge
>
cutEdges
;
std
::
set
<
MVertex
*>
cutVertices
;
std
::
vector
<
MTriangle
*>
newTris
;
std
::
vector
<
MTriangle
*>
newTris
;
std
::
set
<
MEdge
,
Less_Edge
>
newCut
;
cutEdges
.
clear
();
cutVertices
.
clear
();
...
...
@@ -943,7 +943,7 @@ void Centerline::cutByDisk(SVector3 &PT, SVector3 &NORM, double &maxRad){
double
V1
=
a
*
P1
.
x
()
+
b
*
P1
.
y
()
+
c
*
P1
.
z
()
+
d
;
double
V2
=
a
*
P2
.
x
()
+
b
*
P2
.
y
()
+
c
*
P2
.
z
()
+
d
;
bool
inters
=
(
V1
*
V2
<=
0.0
)
?
true
:
false
;
bool
inDisk
=
((
norm
(
P1
-
PT
)
<
rad
)
||
(
norm
(
P2
-
PT
)
<
rad
))
?
true
:
false
;
bool
inDisk
=
((
norm
(
P1
-
PT
)
<
rad
)
||
(
norm
(
P2
-
PT
)
<
rad
))
?
true
:
false
;
double
rdist
=
-
V1
/
(
V2
-
V1
);
if
(
inters
&&
rdist
>
EPS
&&
rdist
<
1.
-
EPS
){
SVector3
PZ
=
P1
+
rdist
*
(
P2
-
P1
);
...
...
@@ -955,7 +955,7 @@ void Centerline::cutByDisk(SVector3 &PT, SVector3 &NORM, double &maxRad){
else
if
(
inters
&&
rdist
>=
1.
-
EPS
&&
inDisk
)
cutVertices
.
insert
(
me
.
getVertex
(
1
));
}
for
(
unsigned
int
i
=
0
;
i
<
triangles
.
size
();
i
++
){
for
(
unsigned
int
i
=
0
;
i
<
triangles
.
size
();
i
++
){
cutTriangle
(
triangles
[
i
],
cutEdges
,
cutVertices
,
newTris
,
newCut
);
}
if
(
isClosed
(
newCut
))
{
...
...
@@ -981,14 +981,14 @@ void Centerline::cutByDisk(SVector3 &PT, SVector3 &NORM, double &maxRad){
// l.getVertex(0)->x(), l.getVertex(0)->y(), l.getVertex(0)->z(),
// l.getVertex(1)->x(), l.getVertex(1)->y(), l.getVertex(1)->z(),
// 1.0,1.0);
// itp++;
// itp++;
// }
// fprintf(f2,"};\n");
// fclose(f2);
}
}
return
;
...
...
@@ -1004,19 +1004,19 @@ double Centerline::operator() (double x, double y, double z, GEntity *ge){
buildKdTree
();
update_needed
=
false
;
}
double
xyz
[
3
]
=
{
x
,
y
,
z
};
bool
isCompound
=
(
ge
->
dim
()
==
2
&&
ge
->
geomType
()
==
GEntity
::
CompoundSurface
)
?
true
:
false
;
std
::
list
<
GFace
*>
cFaces
;
std
::
list
<
GFace
*>
cFaces
;
if
(
isCompound
)
cFaces
=
((
GFaceCompound
*
)
ge
)
->
getCompounds
();
//take xyz = closest point on boundary in case we are on the planar in/out faces
if
(
ge
->
dim
()
==
3
||
(
ge
->
dim
()
==
2
&&
ge
->
geomType
()
==
GEntity
::
Plane
)
||
(
isCompound
&&
(
*
cFaces
.
begin
())
->
geomType
()
==
GEntity
::
Plane
)
){
int
num_neighbours
=
1
;
kdtreeR
->
annkSearch
(
xyz
,
num_neighbours
,
index
,
dist
);
xyz
[
0
]
=
nodesR
[
index
[
0
]][
0
];
xyz
[
0
]
=
nodesR
[
index
[
0
]][
0
];
xyz
[
1
]
=
nodesR
[
index
[
0
]][
1
];
xyz
[
2
]
=
nodesR
[
index
[
0
]][
2
];
}
...
...
@@ -1024,7 +1024,7 @@ double Centerline::operator() (double x, double y, double z, GEntity *ge){
int
num_neighbours
=
1
;
kdtree
->
annkSearch
(
xyz
,
num_neighbours
,
index
,
dist
);
double
d
=
sqrt
(
dist
[
0
]);
double
lc
=
2
*
M_PI
*
d
/
nbPoints
;
double
lc
=
2
*
M_PI
*
d
/
nbPoints
;
return
lc
;
}
...
...
@@ -1040,7 +1040,7 @@ void Centerline::operator() (double x, double y, double z, SMetric3 &metr, GEnt
buildKdTree
();
update_needed
=
false
;
}
//double lc = operator()(x,y,z,ge);
//metr = SMetric3(1./(lc*lc));
...
...
@@ -1050,20 +1050,20 @@ void Centerline::operator() (double x, double y, double z, SMetric3 &metr, GEnt
int
num_neighbours
=
2
;
kdtree
->
annkSearch
(
xyz
,
num_neighbours
,
index2
,
dist2
);
double
d
=
sqrt
(
dist2
[
0
]);
double
lc
=
2
*
M_PI
*
d
/
nbPoints
;
double
lc
=
2
*
M_PI
*
d
/
nbPoints
;
SVector3
p0
(
nodes
[
index2
[
0
]][
0
],
nodes
[
index2
[
0
]][
1
],
nodes
[
index2
[
0
]][
2
]);
SVector3
p1
(
nodes
[
index2
[
1
]][
0
],
nodes
[
index2
[
1
]][
1
],
nodes
[
index2
[
1
]][
2
]);
SVector3
dir0
=
p1
-
p0
;
SVector3
dir1
,
dir2
;
buildOrthoBasis
(
dir0
,
dir1
,
dir2
);
double
lcA
=
4.
*
lc
;
double
lam1
=
1.
/
(
lcA
*
lcA
);
double
lam2
=
1.
/
(
lc
*
lc
);
metr
=
SMetric3
(
lam1
,
lam2
,
lam2
,
dir0
,
dir1
,
dir2
);
delete
[]
index2
;
delete
[]
dist2
;
delete
[]
dist2
;
return
;
}
...
...
@@ -1093,7 +1093,7 @@ void Centerline::printSplit() const{
// fprintf(f3, "SP(%g,%g,%g){%g};\n",
// v->x(), v->y(), v->z(),
// (double)v->getNum());
// itj++;
// itj++;
// }
// fprintf(f3,"};\n");
// fclose(f3);
...
...
@@ -1110,7 +1110,7 @@ void Centerline::printSplit() const{
}
fprintf
(
f4
,
"};
\n
"
);
fclose
(
f4
);
}
#endif
This diff is collapsed.
Click to expand it.
Mesh/meshGFaceDelaunayInsertion.cpp
+
236
−
236
View file @
b3f71e99
This diff is collapsed.
Click to expand it.
Post/PView.cpp
+
1
−
1
View file @
b3f71e99
...
...
@@ -90,7 +90,7 @@ PView::PView(const std::string &xname, const std::string &yname,
_data
->
setFileName
(
yname
+
".pos"
);
_options
=
new
PViewOptions
(
PViewOptions
::
reference
);
_options
->
type
=
PViewOptions
::
Plot2D
;
_options
->
axes
=
2
;
_options
->
axes
=
3
;
//_options->lineWidth = 2.;
//_options->pointSize = 4.;
_options
->
axesLabel
[
0
]
=
xname
;
...
...
This diff is collapsed.
Click to expand it.
Preview
0%
Loading
Try again
or
attach a new file
.
Cancel
You are about to add
0
people
to the discussion. Proceed with caution.
Finish editing this message first!
Save comment
Cancel
Please
register
or
sign in
to comment