Dump multispheres properties instead of atoms properties

Submitted by Matteo on Mon, 02/20/2017 - 12:25

Hello,

Is there any way to dump the multispheres properties so that I can visualise them correctly in Paraview, for their actual mass, volume or spherical equivalent radii?

Thanks

Best regards
Matteo

richti83's picture

richti83 | Mon, 02/20/2017 - 15:41

I use the following workflow:
compute and dump multisphere-information in LIGGGHTS Skript:

compute bid all rigid property id_multisphere #ID
compute btype all rigid property clumptype #Type
compute xcm all rigid property xcm #center of mass
compute fcm all rigid property fcm #force at center of mass
compute quat all rigid property quat #rotation arround COM
compute vel all rigid property vcm #velocity of COM
dump rigids all local 100 post/rigids*.dump c_bid c_btype c_xcm[1] c_xcm[2] c_xcm[3] c_quat[1] c_quat[2] c_quat[3] c_quat[4] c_vel[1] c_vel[2] c_vel[3] c_fcm[1] c_fcm[2] c_fcm[3]

Open rigids*.dump file in Paraview, using liggghts_rigid_reader [1]
Open STL File of scanned clump in Paraview,
Attach following python filter to both:

input1 = self.GetInputDataObject(0, 0)
input2 = self.GetInputDataObject(0, 1)
output = self.GetOutputDataObject(0)
newPoints = vtk.vtkPoints()
newcells = vtk.vtkCellArray()
P=input1.GetPointData()
B=input2.GetPointData()
BODIES=input1.GetNumberOfPoints()
EDGES=input2.GetNumberOfPoints()
SURFACES=input2.GetNumberOfCells()
outputarrayID = vtk.vtkDoubleArray()
outputarrayID.SetName("ID")
outputarrayID.SetNumberOfTuples(BODIES*SURFACES)
outputarrayT = vtk.vtkDoubleArray()
outputarrayT.SetName("type")
outputarrayT.SetNumberOfTuples(BODIES*SURFACES)
c=0
k=0
for body in range(0,BODIES):
r=input1.GetPoint(body)
CLUMP_ID=P.GetArray("id").GetValue(body);
CLUMP_TYPE=P.GetArray("type").GetValue(body);
M=P.GetAbstractArray('M')
T11=M.GetComponent(body,0)
T12=M.GetComponent(body,1)
T13=M.GetComponent(body,2)
T21=M.GetComponent(body,3)
T22=M.GetComponent(body,4)
T23=M.GetComponent(body,5)
T31=M.GetComponent(body,6)
T32=M.GetComponent(body,7)
T33=M.GetComponent(body,8)
for surface in range(0,SURFACES):
outputarrayID.SetValue(k,CLUMP_ID)
outputarrayT.SetValue(k,CLUMP_TYPE)
k+=1
cell = input2.GetCell(surface)
COUNT=cell.GetNumberOfPoints()
for j in range(0,COUNT): #all Points of this cell
id=cell.GetPointId(j)
#print id
coord = input2.GetPoint(id)
x, y, z = coord[:3]
#print x,y,z
#print r
xnew=T11*x+T12*y+T13*z+r[0]
ynew=T21*x+T22*y+T23*z+r[1]
znew=T31*x+T32*y+T33*z+r[2]
newPoints.InsertPoint(c, xnew, ynew, znew)
c+=1
newcells.InsertNextCell(3)
newcells.InsertCellPoint(c-1)
newcells.InsertCellPoint(c-2)
newcells.InsertCellPoint(c-3)
output.SetPoints(newPoints)
output.SetCells(5,newcells)
output.GetCellData().AddArray(outputarrayID)
output.GetCellData().AddArray(outputarrayT)

result: https://www.youtube.com/watch?v=iC5OZtnLc8Y

best, Christian.

[1] https://github.com/richti83/ParaView_Reader_for_LIGGGHTS/tree/master/rea...

I'm not an associate of DCS GmbH and not a core developer of LIGGGHTS®
ResearchGate | Contact

govind | Thu, 09/28/2017 - 16:09

If you don't mind I want to know how do u do this in paraview? I want to learn also these things. Can you given me some guidance?

Govind