error using saveVtk() to get tesselation geometry
Hi There,
I try to use saveVtk() to take a look at the tesselation geometry.
I have added the following when defining the engines :
FlowEngine(
GlobalStiffness
and
try
flow.saveVtk()
However I get the error: Triangulation does not exist. Sorry.
I have tried adding the tesselationwrapper as follows
TW=Tesselation
TW.triangulate()
TW.computeVolu
flow.saveVtk()
However I still get the same error. Is there another step missing which computes the triangulation and makes it available for saveVtk?
many thanks
Jesse
Question information
- Language:
- English Edit question
- Status:
- Invalid
- For:
- Yade Edit question
- Assignee:
- No assignee Edit question
- Solved by:
- Bruno Chareyre
- Solved:
- Last query:
- Last reply:
Revision history for this message
|
#1 |
Hi,
FlowEngine contains a triangulation/
Yes, definitely, saveVtk will only output the result of some fluid solving. If the fluid is not solved (flow.dead=True or no timestepping) the output is empty.
You may try flow.emulateAct
I actually doubt outputing a vtk is what you need, but it might be another question.
HTH
Bruno
Revision history for this message
|
#2 |
Hi Bruno,
Thank you for the quick reply. Ive added flow.emulateAct
flow.dead=0
flow.defTolera
flow.meshUpdat
flow.useSolver=3
flow.permeabil
flow.viscosity=10
flow.bndCondIs
flow.bndCondVa
flow.boundaryU
O.dt=0.1e-3
O.dynDt=False
flow.emulateAc
flow.saveVtk()
The console output however from saveVtk() is:
0 : Vh==NULL!! id=0 Point=1.92862 1.65482 4.78545 rad=0.738095
1 : Vh==NULL!! id=1 Point=3.64074 2.16205 4.14654 rad=0.714286
2 : Vh==NULL!! id=2 Point=2.03866 3.15345 0.0887077 rad=0.690476
3 : Vh==NULL!! id=3 Point=3.56674 0.0664596 3.07013 rad=0.666667
4 : Vh==NULL!! id=4 Point=0.54744 0.471346 2.35463 rad=0.642857
5 : Vh==NULL!! id=5 Point=0.469151 1.40331 4.06784 rad=0.619048
6 : Vh==NULL!! id=6 Point=3.80346 1.90056 1.84533 rad=0.595238
7 : Vh==NULL!! id=7 Point=0.882771 2.24981 0.484466 rad=0.571429
8 : Vh==NULL!! id=8 Point=4.90384 3.31689 4.64228 rad=0.547619
9 : Vh==NULL!! id=9 Point=2.64444 4.79931 2.32225 rad=0.52381
10 : Vh==NULL!! id=10 Point=2.22052 1.92293 1.89341 rad=0.5
11 : Vh==NULL!! id=11 Point=1.71494 1.92762 3.09753 rad=0.47619
12 : Vh==NULL!! id=12 Point=0.988165 4.71568 0.649083 rad=0.452381
13 : Vh==NULL!! id=13 Point=3.12047 3.11159 4.68382 rad=0.428571
14 : Vh==NULL!! id=14 Point=2.58058 2.02467 2.95446 rad=0.404762
15 : Vh==NULL!! id=15 Point=4.8716 2.04733 1.46188 rad=0.380952
16 : Vh==NULL!! id=16 Point=1.43367 4.63815 3.07543 rad=0.357143
17 : Vh==NULL!! id=17 Point=4.30153 4.29709 4.99369 rad=0.333333
18 : Vh==NULL!! id=18 Point=3.33361 1.48195 0.298009 rad=0.309524
19 : Vh==NULL!! id=19 Point=1.96329 0.208186 4.32673 rad=0.285714
CHOLMOD error: invalid xtype
CHOLMOD error: argument missing
Segmentation fault (core dumped)
What does Vh=NULL mean and is this the cause of the CHOLMOD and segmentation fault?
thank you
Jesse
Revision history for this message
|
#3 |
Vh=NULL is when a particle is refused by thriangulation algorithm. I don't think it is the result of saveVtk since it is not triangulating, it must happen before.
Main reason for this problem is excessive overlap of some particles (center of one sphere inside another sphere), are you sure you have a consistent input?
B
Revision history for this message
|
#4 |
Hi Bruno,
I am trying with several different packings, indeed one of the testing packings has overlapping spheres. However I have also just used makeClound to generate a packing of non overlapping spheres and I have the same issue.
Ive put a MWE at https:/
To test with the other packing would be to uncomment
sp = ymport.
n = max(int(a[0]) for a in attrs)+1
colors = [randomColor() for _ in xrange(n)]
agglomerates = [[] for _ in xrange(n)]
for s,a in zip(sp,attrs):
aa = int(a[0])
s.agglomerate = aa
s.shape.color = colors[aa]
agglomerate
for a in agglomerates:
O.bodies.
and the packing https:/
hope this helps to clarify
many thanks
Jesse
Revision history for this message
|
#5 |
Hi,
In your MWE the error does not come from saveVtk() (comment line 75-> still the same error).
And in this script also there are large overlaps between the spheres and the walls (and the walls are just very large spheres, i.e. same issue), that's because you generated a periodic cloud in a finite box.
Bruno
Revision history for this message
|
#6 |
Hi Bruno,
thank you for the reply. I try with other packings but cant get around the Vh=NULL error
for example the code below with the following input file https:/
I understand that the error isnt necessarily coming from saveVTK, should I open a new question instead? Is there another way to access the triangulation data? My goal is actually to be able to read the edges and vertices of the regular triangulation.
code:
O.reset()
sphs = ymport.text('clump' + str(clumpNum) + '.txt')
O.bodies.
print('Number of particles: ', len(sphs))
##########
raw_input('Before FlowEngine saveVtk')
O.engines = [
ForceResett
InsertionSo
#Interactio
#[Ig2_
# [Ig2_Sphere_
# [Ip2_FrictMat_
# [Law2_ScGeom_
#),
FlowEngine(
GlobalStiff
#PeriTriaxC
#NewtonInte
NewtonInteg
]
factor = .5
O.dt = factor * PWaveTimeStep()
#B. Activate flow engine and set boundary conditions in order to get permeability
#TW=Tesselatio
#TW.triangulate()
#TW.computeVol
flow.dead=0
flow.defTolera
flow.meshUpdat
flow.useSolver=3
flow.permeabil
flow.viscosity=10
flow.bndCondIs
flow.bndCondVa
flow.boundaryU
O.dt=0.1e-3
O.dynDt=False
flow.emulateAc
flow.saveVtk()
raw_input(
Revision history for this message
|
#7 |
Hi, better reproduce the error with makeCloud to make a true MWE please.
Bruno
p.s. Why are you using FlowEngine instead of TesselationWrapper to get list of edges? Just curious.
Revision history for this message
|
#8 |
Hi Bruno,
I looked at the tesselationwrapper, but didnt see a way to access the edges as a list pertaining to each voronoi cell. Is this possible?
I saw saveVtk can save out the voronoi diagram and seemed to be the easier way to access the geometry...
I should really say (taking yet another step back and looking at the bigger picture), what I really want to do is obtain a pore size distribution for my packing. (The pore size I am defining to be the radius of an inscribed sphere that fits in the space between particles.)
So to calculate something that looks like this https:/
many thanks
Jesse
Revision history for this message
|
#9 |
TesselationWrapper also has some vtk-export functions, e.g. defToVtk(), but it doesn't seem to be the question in your case.
What you want is access to internals of the triangulation and going through vtk export just for this is a bit awkward. It would be probably more convenient to manipulate the data directly in python.
FlowEngine received more attention recently and it's python interface gives more access to internals. It is thus interesting to hear your feedback on what would me missing in TesselationWrapper.
But in the end, if what you want is the set of inscribed spheres my impression is that you can get it directly from TwoPhaseFlowEngine (I don't have the precise procedure in mind right now unfortunately).
Bruno
Revision history for this message
|
#10 |
I've tried using the TwoPhaseFlowEngine to obtain the inscribed spheres. However I had issues where the pore sizes calculated were all very very small. I guess this lead me to try to vizualise the pore structure by exporting it to vtk. Perhaps these are all the same issue. For example I have added :
TwoPhase=
TwoPhase.
to the script I posted on 2018-06-19 and I get the following error:
negative volume for an ordinary pore (temp warning, should still be safe)
*** stack smashing detected ***: /usr/bin/python terminated
Aborted (core dumped)
how do you suggest getting to the bottom of the triangulation/
many thanks
Jesse
Revision history for this message
|
#11 |
I can't understand this triangulation failure unless you prove it with a MWE (using makeCloud). It seems you are using some weird input, that's my guess at the moment. It is completely unrelated to the other questions like getting inscribed radius.
Bruno
Revision history for this message
|
#12 |
Hi Bruno,
thank you again for the reply. I started making a MWE using makeCloud and it ends up that both saveVtk and the TwoPhaseFlowEngine work fine (although when I load the vtk in paraview I dont know how to see the triangulation, is there documentation for this?).
However, when I try loading my geometry the tesselation doesnt seem work any more. Does this mean if I would like to twophaseflowEngine or the tesselation in general the geometry needs to be regenerated with makeCloud? Could O.bodies.
many thanks
Jesse
Revision history for this message
|
#13 |
Hi, loading the vtk file in paraview should already display a mesh in principle. Not the case?
Both methods work: loading a geometry or makeCloud, they both go through O.bodies.append at some point.
The problem is not the method, it is the data. There must be some inconsistency in your specific way or reloading something resulting in unexpected geometry (did you inspect visually at least? for answer #5 it was enough to open the 3D view to understand the problem).
Revision history for this message
|
#14 |
Hi,
loading the vtk file in paraview I see a convex mesh. I dont see how to see the triangulated mesh. http://
I've made another MWE where I load clump0.txt reset the scene and use O.bodies.
(I also dont see anything wrong when I look at the spheres in the 3d viewport. )
When I run the twophaseflowengine I print the cell id and pore size when TwoPhaseFlowEng
cout << endl << "Test3 Cell: " << cell->info().id << ": Pore: " << Rin;
The output is a long list of pores all the same size ending with *** stack smashing detected ***: /usr/bin/python terminated
.........
Test3 Cell: 914: Pore: 4.92e-10
Test3 Cell: 915: Pore: 4.92e-10
Test3 Cell: 916: Pore: 4.92e-10
Test3 Cell: 917: Pore: 4.92e-10
Test3 Cell: 918: Pore: 4.92e-10
Test3 Cell: 919: Pore: 4.92e-10
Test3 Cell: 920: Pore: 4.92e-10negative volume for an ordinary pore (temp warning, should still be safe)
negative volume for an ordinary pore (temp warning, should still be safe)
negative volume for an ordinary pore (temp warning, should still be safe)
*** stack smashing detected ***: /usr/bin/python terminated
Aborted (core dumped)
When I try to run saveVtk the error is as follows.
......
148 : Vh==NULL!! id=148 Point=7.42392e-07 6.92963e-07 8.02883e-07 rad=1.23e-08
149 : Vh==NULL!! id=149 Point=8.03578e-07 7.18035e-07 8.15755e-07 rad=1.23e-08
150 : Vh==NULL!! id=150 Point=7.43404e-07 7.82494e-07 7.64931e-07 rad=1.23e-08
151 : Vh==NULL!! id=151 Point=6.77428e-07 7.66925e-07 7.91695e-07 rad=1.23e-08
CHOLMOD error: invalid xtype
CHOLMOD error: argument missing
Segmentation fault (core dumped)
clump0.txt I have pasted here: https:/
def tryvtk():
O.engines = [
ForceResett
InsertionSo
#Interactio
#[Ig2_
# [Ig2_Sphere_
# [Ip2_FrictMat_
# [Law2_ScGeom_
#),
FlowEngine(
GlobalStiff
#PeriTriaxC
#NewtonInte
NewtonInteg
]
factor = .5
O.dt = factor * PWaveTimeStep()
flow.dead=0
flow.defTolera
flow.meshUpdat
flow.useSolver=3
flow.permeabil
flow.viscosity=10
flow.bndCondIs
flow.bndCondVa
flow.boundaryU
O.dt=0.1e-3
O.dynDt=False
flow.emulateAc
flow.saveVtk()
raw_input(
def main():
from yade import export,ymport
import random
random.seed(1)
dim = Vector3(
O.materials.
O.periodic=True
O.cell.
attrs = []
#sp = ymport.
sp = ymport.
spheres = []
for b in O.bodies:
if b.isClump == False:
spheres.
raw_input(
O.resetThisScene() #####!!!!!!!
raw_input(
for pos,r in spheres:
O.bodies.
#print(pos,r)
raw_input(
#quit()
TwoPhase=
TwoPhase.
raw_input(
tryvtk()
raw_input(
#quit()
print('Number of particles = ' , len(O.bodies))
print('Total particles inc. periodic', len(O.bodies))
try:
from yade import qt
qt.View()
except:
pass
if __name_
main()
Revision history for this message
|
#15 |
You Vtk data file is ok, if you want to see internal mesh you need to "clip" or "slice" (and display "surface with edges" instead of wireframe).
For the second part: I guess you have overlapping spheres in the clump, and triangulation is possible only for limited overlaps.
Bruno
Revision history for this message
|
#16 |
Thanks for the quick reply. With surface with edges and a clip I can see the internal mesh.
I made sure this time to avoid overlapping spheres. The clump is generated by using makeCloud, compressing the volume with triaxcontroller and saving only spheres which are inside a predicate.
Revision history for this message
|
#17 |
Not sure what happens with clumps. If you can summarize the problem and reproduce in a simple way with a self-contained script (sorry but I'm tired of downloading extra files on external servers today...) I may have a look.
Bruno
Revision history for this message
|
#18 |
Thank you for the reply.
Ive made a new MWE which also generates clump0.txt.
https:/
Revision history for this message
|
#19 |
WARNING! rmin>rmax. rmin=1e-11 ,rmax=1.66533e-12
WARNING! rmin>rmax. rmin=1e-11 ,rmax=1.66533e-12
WARNING! rmin>rmax. rmin=1e-11 ,rmax=3.93018e-12
WARNING! rmin>rmax. rmin=1e-11 ,rmax=3.93018e-12
Segmentation fault
In addition, as I was seeing random lines from your last script I found:
def createAggregate
O.reset()
.....
then:
for idx, radAgg in enumerate(
createAggrega
Which basically means that you keep creating new clumps, adding them to O.bodies, and finally erasing everything (including O.engines) before creating the next clump...
Another exotic thing is at the end:
if __name_
main()
No clue what it's doing. Maybe it is not wrong, but for sure it is not what we expect in a MWE.
This was supposed to be about saveVtk() but it is much too confused.
I think I'll give up on this question, please really read [1] again if it can help.
B
Revision history for this message
|
#21 |
My apologies that the script is convoluted. Ive had to merge many scripts into one to make a single MWE that doesnt rely on any external files. However your description is correct. I create a clump, add it to the simulation save it. Then delete everything again to before moving on to create a new clump.
Once I have created a clump I want to have information from the tesselation. When I just create a particle cloud with makeCloud I dont have an issue running the tesselationwrapper or twophaseflowengine. However when I try to run it on the clump I get those segmentations faults.
The if __name_
doesnt have an effect on the script it just lets me the python file and use the functions in other scripts. Sorry for including it
I really would appreciate your continued help on this question.
Revision history for this message
|
#22 |
> delete everything again to before moving on to create a new clump
If it crashes on clump n°1 then there is no need for the MWE to erase and restart.
If it crashes on clump n°X then send a MWE using only clump n°X.
> when I try to run it on the clump I get those segmentations faults.
May I point out that the warnings+segfault which popped up in #19 had never been mentioned before and I have no reason to believe it is related to saveVtk() or anything else discussed previously?
> The if __name_
Why is it in the script then?
> I really would appreciate your continued help on this question.
No sorry. Please keep this question closed. When you have things cleared out maybe try another question which will follow strictly every bit of [1] (especially point 3).
Good luck
Bruno