yade-users team mailing list archive

[wangxiaoliang <question197661@xxxxxxxxxxxxxxxxxxxxx>]

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

    Status: Answered => Open

wangxiaoliang is still having a problem:
well, I have three parts,

in part 1 and part 2, I generate a sample with specific relative density
and then compress to  100kPa isotropicly, which generated a xml file to
be loaded in my triaxial test.

In part 3, I start my deviatoric load process to simulate the drained
triaxial test process. My code is below

nRead=utils.readParamsFromTable(
                knormal=170e6,
                alpha = 0.1,
                fricFinal = 30,
                beta = 1.8,
                eyta = 1.8,
                unknow=True
)
from yade.params import table

knormal = table.knormal
alpha = table.alpha
fricFinal = table.fricFinal
beta = table.beta
eyta = table.eyta

young = 6e6
poisson = 0.1
num_spheres = 10000
key = 'Ec'+str(young)+'poisson'+str(poisson)+'num'+str(num_spheres)
key_output = 'knormal_'+str(knormal)+'alpha_'+str(alpha)+'fric_'+str(fricFinal)+'beta_'+str(beta)+'eta_'+str(eyta)
rate = 0.02
damp=0.2 # damping coefficient
stabilityThreshold=0.01 

 O.load('consolidated_state_'+key+'_'+str(phi)+'.xml')

# we add new materials to our material container
        O.materials.append(MomentMat(eta=eyta,young=young,poisson=0.1,frictionAngle=radians(fricFinal),density=2600,label='spheres-moment'))
        O.materials.append(MomentMat(eta=0.0,young=young,poisson=0.1,frictionAngle=0,density=0,label='walls-moment'))

        # here we modify the material and interaction physics
        for i_body in O.bodies:
                if i_body.id < 6 :
                        i_body.material = O.materials[3]
                else:
                        i_body.material = O.materials[2]
        for i_interactions in O.interactions:
                if i_interactions.id1 > 5 and i_interactions.id2 > 5:
                        i_interactions.phys.Eta = eyta  # only interaction between particles are changed to eyta > 0
#here we move to modify the interaction loop
        Interaction = O.engines[2]
        Interaction =  InteractionLoop(
                [Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom()],
                [Ip2_MomentMat_MomentMat_MomentPhys(Knormal=knormal,Alpha = alpha,Beta=beta,useAlphaBeta=True,userInputStiffness=True)],
                [Law2_SCG_MomentPhys_CohesionlessMomentRotation()]
        )

        # here we redifine the triaxial engine
        triax=O.engines[4] #redifine 3d triaxial engine
        ##We move to deviatoric loading, let us turn internal compaction off to keep particles sizes constant
        triax.internalCompaction=False

        ## Change contact friction (remember that decreasing it would generate instantaneous instabilities)
        triax.setContactProperties(fricFinal)

        ##set independant stress control on each axis
        triax.stressControl_1=triax.stressControl_2=triax.stressControl_3=True
        ## We turn all these flags true, else boundaries will be fixed
        triax.wall_bottom_activated=True
        triax.wall_top_activated=True
        triax.wall_left_activated=True
        triax.wall_right_activated=True
        triax.wall_back_activated=True
        triax.wall_front_activated=True


        ##If we want a triaxial loading at imposed strain rate, let's assign srain rate instead of stress
        triax.stressControl_2=0 #we are tired of typing "True" and "False", we use implicit conversion from integer to boolean
        triax.strainRate2=rate
        triax.strainRate1=100*rate
        triax.strainRate3=100*rate

        ##we can change damping here. What is the effect in your opinion?
        #newton.damping=0.1

        ##Save temporary state in live memory. This state will be reloaded from the interface with the "reload" button.
        O.saveTmp()

#here we redefine the history recorder
        from yade import plot

### a function saving variables
        def history():
                plot.addData(e11=triax.strain[0], e22=triax.strain[1], e33=triax.strain[2],
                    s11=triax.stress(triax.wall_right_id)[0],
                    s22=triax.stress(triax.wall_top_id)[1],
                    s33=triax.stress(triax.wall_front_id)[2],
                    i=O.iter,CN = utils.avgNumInteractions())
        O.engines[4]=PyRunner(iterPeriod=1500,command='history()',label='recorder')

        while 1:
                O.run(1000,True)
                if O.iter%100 == 0:
                        print "iter = ", O.iter
                        print "Syy = ",triax.strain[1]
                if triax.strain[1] > 0.3:
                        break

        ### declare what is to plot. "None" is for separating y and y2 axis
        plot.plots={'i':('e11','e22','e33',None,'s11','s22','s33','CN')}

##display on the screen (doesn't work on VMware image it seems)
##plot.plot()

## In that case we can still save the data to a text file at the the end of the simulation, with: 
        O.save('Final'+key_output+'.xml')
        plot.saveDataTxt('shear'+key_output+'_'+str(phi)+'.txt')

-- 
You received this question notification because you are a member of
yade-users, which is an answer contact for Yade.