← Back to team overview

yade-users team mailing list archive

[Question #678789]: FEM-DEM coupling (multiscale)

 

New question #678789 on Yade:
https://answers.launchpad.net/yade/+question/678789

Hi:
I did a dem-fem multiscale couplering base on the updated code of Jan Stranský multi-example and modified.
1.I slowed down the loading speed in fem. I changed the loading time from 0.2s to 200s. 
2.I changed the dt from 2e-4 to 2e-2
3.In order to make sure the dem is stable, I changed the density in dem from 1e3 to 1e7 and set the gravity(0,0,0.0)

Here are the codes of dem:
#from libyade import yade
from yade import *
from yade import pack, Vector3, Vector3

nREVs = 16

def randomPeriPack(radius,initSize,seed):
	O.switchScene(); O.resetThisScene()
	sp=pack.SpherePack()
	O.periodic=True
	O.cell.setBox(initSize)
	sp.makeCloud(Vector3().Zero,O.cell.refSize,radius,0.,-1,True,seed=seed)
	O.engines=[ForceResetter(),InsertionSortCollider([Bo1_Sphere_Aabb()],verletDist=.05*radius),InteractionLoop([Ig2_Sphere_Sphere_ScGeom()],[Ip2_FrictMat_FrictMat_FrictPhys()],[Law2_ScGeom_FrictPhys_CundallStrack()]),PeriIsoCompressor(charLen=2*radius,stresses=[-100e9,-1e8],maxUnbalanced=1e-2,doneHook='O.pause();',globalUpdateInt=20,keepProportions=True),NewtonIntegrator(damping=.8)]
	O.materials.append(FrictMat(young=30e9,frictionAngle=.6,poisson=.3,density=1e7))
	for s in sp: O.bodies.append(utils.sphere(s[0],s[1]))
	O.dt=utils.PWaveTimeStep()
	O.timingEnabled=True
	O.run(); O.wait()
	for b in O.bodies: b.state.pos = O.cell.wrap(b.state.pos)
	ret=pack.SpherePack()
	ret.fromSimulation()
	O.switchScene()
	return ret

for i in xrange(nREVs):
	if i != 0:
		newScene = O.addScene()
		O.switchToScene(newScene)
	young = 4e6 if i<8 else 2e6
	O.materials.append(CpmMat(young=young,poisson=.2,epsCrackOnset=1e100,sigmaT=1e100,relDuctility=2))
	sp = randomPeriPack(.01,.1*Vector3.Ones,i+10)
	sp.toSimulation()
	print('scene',i,'created')
	O.engines = [
		ForceResetter(),
		InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=1.5,label='is2aabb')],allowBiggerThanPeriod=True),
		InteractionLoop(
			[Ig2_Sphere_Sphere_ScGeom(interactionDetectionFactor=1.5,label='ss2d3dg')],
			[Ip2_CpmMat_CpmMat_CpmPhys()],
			[Law2_ScGeom_CpmPhys_Cpm()]),
		NewtonIntegrator(damping=.2,gravity=(0,0,0.0)),
	]
	O.dt = 0.
	O.step()
	is2aabb.aabbEnlargeFactor = 1.
	ss2d3dg.interactionDetectionFactor = 1.

def vtkExport(i):
	base = '/tmp/multi1_yade_{:02d}'
	from yade import export
	for rvei in xrange(nREVs):
		O.switchToScene(rvei)
		name = base.format(rvei)
		export.VTKExporter(name,i).exportSpheres(useRef=True,what=[('dspl','b.state.displ()')])

However, the stress got extremely large after several steps.

Solving [Step number        0, Time    2.000000e-02]
('stress', Matrix3(-2.7232542208268342e-11,2.1607057838358735e-11,-7.565900182489134e-12, 2.1607057838358735e-11,-3.5733805174291325e-11,5.910176084244153e-12, -7.565900182489134e-12,5.910176084244153e-12,2.6839670851857137e-12))
('strain', Matrix3(0,0,0, 0,0,0, 0,0,0))
v:
Matrix3(0,0,0, 0,0,0, 0,0,0)

Solving [Step number        1, Time    4.000000e-02]
('stress', Matrix3(-1.083433491135522e-8,-1.3503776707434368e-8,1.7486952135573553e-8, -1.3503776707434368e-8,3.6347064127015305e-8,-2.926186112822691e-8, 1.7486952135573553e-8,-2.926186112822691e-8,8.99554151810844e-10))
('strain', Matrix3(0,0,0, 0,0,0, 0,0,0))
v:
Matrix3(0,0,0, 0,0,0, 0,0,0)

Solving [Step number        2, Time    6.000000e-02]
('stress', Matrix3(0.000015525172297319974,0.0000795526387192835,-0.00009283124932038486, 0.0000795526387192835,-0.0003630625371313257,0.0002590613138648092, -0.00009283124932038486,0.0002590613138648092,-0.00005720234018106285))
('strain', Matrix3(0.00005392650492956297,-1.7164899182052643e-18,-9.869235370759615e-6, -1.7164899182052643e-18,0.000053926504929562556,-9.869235370764533e-6, -9.869235370759615e-6,-9.869235370764533e-6,-0.00013675213675213405))
v:
Matrix3(0,1.3597608349867324e-16,2.353049132322896e-17, 1.3597608349867324e-16,0,6.850607918032826e-17, 2.353049132322896e-17,6.850607918032826e-17,0)

Solving [Step number        3, Time    8.000000e-02]
('stress', Matrix3(-11.3107674295793,1.7129703607246294,-2.787768422239802, 1.7129703607246294,-9.609769538283032,-0.23642769428808477, -2.787768422239802,-0.23642769428808477,-12.643860844491378))
('strain', Matrix3(0.00006091697779284176,4.488899208845239e-15,-0.000011148580698508899, 4.488899208845239e-15,0.00006091697779162887,-0.00001114858069311587, -0.000011148580698508899,-0.00001114858069311587,-0.0003103221564777195))
v:
Matrix3(1.2212453270876722e-13,1.1814561132037645e-13,1.6167065989061103e-13, 1.1814561132037645e-13,-1.6653345369377348e-14,-1.1915983548294906e-13, 1.6167065989061103e-13,-1.1915983548294906e-13,5.551115123125783e-14)

Solving [Step number        4, Time    1.000000e-01]
('stress', Matrix3(-7.308025222980424e8,9.503633740807763e7,-2.2945586023590547e8, 9.503633740807763e7,-9.180035574511294e8,1.3301196286263055e8, -2.2945586023590547e8,1.3301196286263055e8,-9.263253771518726e8))
('strain', Matrix3(0.000014233498104035896,-2.3639298424766302e-11,-2.6049070587392296e-6, -2.3639298424766302e-11,0.000014233512498104343,-2.6049467336185456e-6, -2.6049070587392296e-6,-2.6049467336185456e-6,-0.00045189774370015776))
v:
Matrix3(-5.176525874617255e-10,-1.0401586226850383e-9,-9.210391726013966e-10, -1.0401586226850383e-9,-1.3655188091377113e-10,5.564194956105476e-10, -9.210391726013966e-10,5.564194956105476e-10,-5.627942556429844e-10)

Why did the stress get so large? The strain seems stable.
Yours, Joe

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