yade-users team mailing list archive

[Daria <question689000@xxxxxxxxxxxxxxxxxxxxx>]

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

    Status: Answered => Open

Daria is still having a problem:
scriptYade.py:

from __future__ import print_function
import sys
from yadeimport import *
from yade.utils import *

initMPI()                           #Initialize the mpi environment, always required.
fluidCoupling = yade.FoamCoupling();     #Initialize the engine
fluidCoupling.getRank();            #part of Initialization.


#example of spheres in shear flow : two-way point force coupling
class simulation():

	def __init__(self):
		epsilon = 1e-08
		minval = epsilon
		maxval = 1-epsilon
		length = maxval - minval
		halflength = length / 2
		radius = 0.01
		radiusEpsilon = 0.001

#		O.periodic = True
#		O.cell.setBox(maxval, maxval, maxval)

		numspheres=1000
		young = 1e6#ice#5e6#1
		poisson = 0.2#0.4
		density = 1000

		O.materials.append(FrictMat(young=young,poisson=poisson,frictionAngle=radians(15),density=density,label='spheremat'))
		O.materials.append(FrictMat(young=200e9,poisson=0.001,frictionAngle=0,density=0,label='wallmat'))

		#wall coords, use facets for wall BC:
		v0 = Vector3(minval,minval,minval)
		v1 = Vector3(minval,minval,maxval)
		v2 = Vector3(maxval,minval,minval)
		v3 = Vector3(maxval,minval,maxval)

		v4 = Vector3(minval,maxval,minval)
		v5 = Vector3(minval,maxval,maxval)
		v6 = Vector3(maxval,maxval,minval)
		v7 = Vector3(maxval,maxval,maxval)

		df0 = facet([v0,v1,v2],material='wallmat')
		O.bodies.append(df0)
		df1 = facet([v1,v3,v2],material='wallmat')
		O.bodies.append(df1)

		uf0 = facet([v4,v5,v6],material='wallmat')
		O.bodies.append(uf0)
		uf1 = facet([v5,v7,v6],material='wallmat')
		O.bodies.append(uf1)

		lf0 = facet([v0,v1,v5],material='wallmat')
		O.bodies.append(lf0)
		lf1 = facet([v1,v5,v4],material='wallmat')
		O.bodies.append(lf1)

		rf0 = facet([v2,v3,v7],material='wallmat')
		O.bodies.append(rf0)
		rf1 = facet([v2,v7,v6],material='wallmat')
		O.bodies.append(rf1)

		ff0 = facet([v1,v5,v4],material='wallmat')
		O.bodies.append(ff0)
		ff1 = facet([v5,v7,v4],material='wallmat')
		O.bodies.append(ff1)

		bf0 = facet([v0,v4,v2],material='wallmat')
		O.bodies.append(bf0)
		bf1 = facet([v4,v6,v2],material='wallmat')
		O.bodies.append(bf1)

		#spheres
		mn, mx= Vector3(minval+2*radius, minval+2*radius, minval+2*radius), Vector3(maxval-2*radius, maxval-2*radius, maxval-2*radius)

		sp = pack.SpherePack()
		sp.makeCloud(mn,mx,rMean=radius,rRelFuzz=radiusEpsilon, num=numspheres)
		O.bodies.append([sphere(center,rad,material='spheremat') for center,rad in sp])

                sphereIDs = [b.id for b in O.bodies if
type(b.shape)==Sphere]

                #coupling engine settings

		fluidCoupling.setNumParticles(len(sphereIDs))
		fluidCoupling.setIdList(sphereIDs)
		fluidCoupling.isGaussianInterp=False  #use pimpleFoamYade for gaussianInterp

		# Integrator
		newton=NewtonIntegrator(damping=0.0, gravity = (0.0 ,-1.81, 0.0))
		# add small damping in case of stability issues.. ~ 0.1 max, also note : If gravity is needed, set it in constant/g dir.

#		O.dt=1e-5
#		O.dynDt=False

		O.engines=[
			ForceResetter(),
			InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()], allowBiggerThanPeriod=True),
			InteractionLoop(
				[Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom()],
				[Ip2_FrictMat_FrictMat_FrictPhys()],
				[Law2_ScGeom_FrictPhys_CundallStrack()]
			),
			GlobalStiffnessTimeStepper(timestepSafetyCoefficient=0.5, label = "ts"),
		    fluidCoupling, #to be called after timestepper
		    PyRunner(command='sim.printMessage()', iterPeriod= 1, label='outputMessage'),
			newton,
		    VTKRecorder(fileName='yadep/3d-vtk-',recorders=['spheres'],virtPeriod=0.02)#iterPeriod=1000)
		]

	def printMessage(self):
		print("********************************YADE-ITER = " + str(O.iter) +" **********************************")
#		if O.iter == 4000:
#			maxVel = 0.05
#			for b in O.bodies:
#				if type(b.shape)==Sphere:
#					bodyVel = abs(b.state.vel.norm())
#					if bodyVel > maxVel:
#						raise ValueError("Body velocity exceeds imposed shear velocity by ", abs(bodyVel-maxVel))

	def irun(self,num):
		O.run(num,1)

if __name__=="__main__":
	sim = simulation()
	sim.irun(100000)
	# print("body id = ", O.bodies[34].id)
	fluidCoupling.killMPI()

import builtins
builtins.sim=sim

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