← Back to team overview

yade-users team mailing list archive

Re: [Question #692998]: principal axis of clump

 

Question #692998 on Yade changed:
https://answers.launchpad.net/yade/+question/692998

    Status: Answered => Open

jsonscript is still having a problem:
hello:
       very thanks for your advice, I choose this simple artificial clump consists of three particles.  I want to compute principle axes of  minimum moment of inertia, I think the result should be (0.57735,0.57735,0.57735), but I always get (-0.80473785, -0.50587936,0.31061722)..could you help check it?


Here are the codes:

from __future__ import print_function
import numpy as np
from builtins import range
from yade import pack,export,qt
from numpy import linalg as LA
import math


#define material for all bodies:
id_Mat=O.materials.append(FrictMat(young=1e6,poisson=0.3,density=1000,frictionAngle=1))
Mat=O.materials[id_Mat]

#define engines:
O.engines=[
	ForceResetter(),
	InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb()]),
	InteractionLoop(
		[Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom()],
		[Ip2_FrictMat_FrictMat_FrictPhys()],
		[Law2_ScGeom_FrictPhys_CundallStrack()]
	),
	NewtonIntegrator(damping=0.7,gravity=[0,0,-10])
]


#### show how to use appendClumped():


#create 3 clumps:
clump1=O.bodies.appendClumped([\
sphere([0,0,0], material=Mat, radius=1),\
sphere([1,1,1], material=Mat, radius=1),\
sphere([2,2,2], material=Mat, radius=1)\
])


#definition for getting informations from all clumps:
def getClumpInfo():
	for b in O.bodies:
		if b.isClump:
			print('Clump ',b.id,' has following members:')
			keys = list(b.shape.members.keys())
			for ii in range(0,len(keys)):
				print('- Body ',keys[ii])
			r = b.state.ori.toRotationMatrix()
			inertiaTensor = r*b.state.inertia.asDiagonal()*r.transpose()	
  			w,v = LA.eig(inertiaTensor)
  			emax=max(abs(w[0]),abs(w[1]),abs(w[2]))
			print(inertiaTensor)
			print(w)
  			print(v)
  			emin=min(abs(w[0]),abs(w[1]),abs(w[2]))
  			o=np.zeros(3)
			for j in range(0,3):
			  if abs(w[j])==emin:
			    o[0]=v.item((0,j))
			    print(v.item((0,j)))
			    o[1]=v.item((1,j))
			    o[2]=v.item((2,j))			
			print(o)
                        print('inertia:',b.state.inertia)
			print('mass:',b.state.mass,'\n')
			print('ori:',b.state.ori,'\n')

getClumpInfo()
renderer = qt.Renderer()
qt.View()

Thanks

-- 
You received this question notification because your team yade-users is
an answer contact for Yade.