← Back to team overview

yade-users team mailing list archive

Re: [Question #688203]: Simulation blocked after O.run in a conditional translation motion

 

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

    Status: Needs information => Open

Rioual gave more information on the question:
...Here is the code:

#Start

import time
from yade import qt 
from yade import ymport

#Import packing in the cylinder from prepare-packing.py
O.load('init_Final_state_packing.yade')

###########################################
## PhysicalParameters 
Density=2400
frictionAngle=radians(35)
tc = 0.001
en = 0.05
et = 0.05

## Import wall's geometry
facetMat=O.materials.append(ViscElMat(frictionAngle=frictionAngle,tc=tc, en=en, et=et)) # **params sets kn, cn, ks, cs
sphereMat=O.materials.append(ViscElMat(density=Density,frictionAngle=frictionAngle,tc=tc,en=en,et=et))
from yade import ymport
#******************************************
global TransEngload2
#******************************************
# Importation de la roue Warren-Spring

global fctIdsWS
fctIdsWS=O.bodies.append(ymport.stl('WS-FR-1.stl',color=(1,0,0),material=facetMat))


#fctIdscylinder=O.bodies.append(ymport.stl('PotN.stl',color=(1,0,0),material=facetMat))


## Timestep 
O.dt=.2*tc

## Engines
O.engines=[
## Resets forces and momenta the act on bodies
ForceResetter(),
## Associates bounding volume to each body.
#BoundDispatcher([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()]),
## Using bounding boxes find possible body collisions.
#InsertionSortCollider(),
InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()]),
## Interactions
InteractionLoop(
		# the loading plate is a facet, we need to handle sphere+sphere, sphere+facet
		[Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom(),Ig2_Wall_Sphere_ScGeom()],
		[Ip2_FrictMat_FrictMat_FrictPhys()],
		[Law2_ScGeom_FrictPhys_CundallStrack()]
		),  
## Apply gravity
GravityEngine(gravity=[0,0,-9.81]),
## Cundall damping must been disabled!
NewtonIntegrator(damping=0),


## Apply kinematics to wheel
PyRunner(command='kinematics_WS()',realPeriod=1,label='kine'),
]

from yade import qt
qt.View()


def kinematics_WS():


	h_WS = calc_h()[0]
	hmaxSpheres = calc_h()[1]


	print 'test0','h_WS=',h_WS,'hmaxSpheres=',hmaxSpheres

        TransEngload2 =
TranslationEngine(ids=fctIdsWS,translationAxis=[0,-1,0],velocity=10)

        O.engines=O.engines+[TransEngload2]


	while h_WS > hmaxSpheres: 

		print 'test01'
		
		TransEngload2.dead = False

		print 'test02'
		
		O.run(1,True)

#               O.engines=O.engines+[PyRunner(command='calc_h()')]


		print 'test03'

		h_WS = calc_h()[0]
		hmaxSpheres = calc_h()[1]

		print 'test1','h_WS=',h_WS,'hmaxSpheres=',hmaxSpheres
	
		TransEngload2.dead = True

        else:

                amTOT = sum((O.bodies[facetid].state.angMom)[1] for
facetid in fctIdsWS)

	
		while (amTOT< 1e-10):

			O.engines=O.engines+[PyRunner(command='addTorque()',iterPeriod=1)] 
			
			amTOT = sum((O.bodies[facetid].state.angMom)[1] for facetid in fctIdsWS)

                else:

# Stop simulation et measurement of the shear torque and cohesion

			C_WS = 3* (imposed_T)/(2*Pi*(R0^3-R1^3))
			print 'C_WS=',C_WS

			print '********End of the calculation for Cohesion***********' #
			O.pause()

#*********************************************************************************************************************************
#This function calculate the height of the boundary h_WS and the maximum height of the packing hmax_Spheres
def calc_h():

###########################################
#Calculate h_WS 

	minY = 1e99
	maxY = -1e99

	for facet in fctIdsWS:
#		print 'facet=', facet
   		facet= O.bodies[facet]	
		vs = [facet.state.pos + facet.state.ori * v for v in facet.shape.vertices] # vertices in global coord system
#		print 'vs=',vs 

		minY = min(minY,min(v[1] for v in vs))
		maxY = max(maxY,max(v[1] for v in vs)) ###				

#       print 'maxY=',maxY,'minY=',minY

	h_WS = maxY
        
###########################################
#Calculate hmax_Spheres

	idHMax=0				# on definit une variable pour identifier quel corps a la hauteur la plus haute
	hMax=0.0				# initialisation de la hauteur a zero
	for i in O.bodies:			# on parcours tout les corps du systeme
		h=i.state.pos[1]		# on extrait la position selon l'axe y
		if (type(i.shape).__name__ == 'Sphere' and h>hMax):	#si le solide est une sphere et sa position est plus haut que hmax
			hMax=h			# on enregistre sa hauteur 
			idHMax=i.id		# on garde en memoire de quel corps il s'agit
#	if (idHMax == 0): return rMax		# valeur de retour par default, rMAX
#	else: return (hMax+O.bodies[idHMax].shape.radius)
	
	hmaxSpheres = hMax+O.bodies[idHMax].shape.radius
	
        return h_WS,hmaxSpheres

#*****************************************************************************************************************************
#*****************************************************************************************************************************
# Add incremental torque on the WS boundary

def addTorque():


	increment_T= 0.001
	
	for Idfacets in fctIdsWS:
		O.forces.addT(Idfacets,(0,increment_T,0))
		imposed_T = O.forces.permT(Idfacets)
	return imposed_T

#*****************************************************************************************************************************

...So this is the complete script, maybe the problem is in relation with
the problem of dealing with an aditional engine in a while loop (?) . I
am not sure of myself on this....

Thanks for your feedback,

best

Vincent,

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