Python plot.addData函数代码示例

def sim(angle):
	O.reset()
	O.dt = 1e-5
	# create two bricks with all blockedDOFs
	mortar = O.materials.append(MortarMat(young=young,poisson=GOverE,tensileStrength=tensileStrength,cohesion=cohesion,frictionAngle=frictionAngle,compressiveStrength=compressiveStrength,ellAspect=ellAspect))
	bs = b1,b2 = [polyhedron(((-1,-1,-1),(+1,-1,-1),(-1,+1,-1),(+1,+1,-1),(-1,-1,+1),(+1,-1,+1),(-1,+1,+1),(+1,+1,+1)),material=mortar) for i in (0,1)]
	b2.state.pos = (0,0,2)
	for b in bs:
		b.state.blockedDOFs = 'xyzXYZ'
	O.bodies.append(bs)
	#
	# factor to safely create interaction of just touching bricks
	factor=1.1
	O.engines=[
		ForceResetter(),
		InsertionSortCollider([Bo1_Polyhedra_Aabb(aabbEnlargeFactor=factor,label='bo1')]),
		InteractionLoop(
			[Ig2_Polyhedra_Polyhedra_ScGeom(interactionDetectionFactor=factor,label='ig2')], 
			[Ip2_MortarMat_MortarMat_MortarPhys()],
			[Law2_ScGeom_MortarPhys_Lourenco()]
		),
		NewtonIntegrator(),
	]
	O.step()
	ig2.interactionDetectionFactor = bo1.aabbEnlargeFactor = 1 # reset the interaction detection enlargement
	b2.state.vel = (sin(angle),0,cos(angle)) # sets velocity to produce desired 
	while len([i for i in O.interactions]) > 0: # run simulatinon until the interaction is broken
		sn,st = i.phys.sigmaN, i.phys.sigmaT.norm() # store last values
		O.step()
		if O.iter > 1e6:
			raise RuntimeError, "TODO" # not to run forever
	plot.addData(sn=sn,st=st) # after the interaction is broken, save stress to plot.data and return
	return

def dataCollector():
	zmax=hMax(2)
	zmin=hMin(2)
	V = S0*(zmax-zmin)
	#poro=porosity(V)	
	#F=O.forces.f(O.bodies[pfIds[45]].shape.node1.id)
	#print F
	S=pi*l**2
	Fnt=O.forces.f(topPlate)[2]
	sigmaN=Fnt/S0
	Fnb=abs(O.forces.f(bottomPlate)[2])
	sigmaNb=Fnb/S0
	pos=O.bodies[topPlate].state.pos[2]
	q=(sigmaNb-sigma)
	cui=(sigmaNb-sigma)/(sigmaNb+sigma)
	p=(sigmaNb+2*sigma)/2
	displ=O.bodies[topPlate].state.pos[2]-O.bodies[bottomPlate].state.pos[2]
	epsr=epsr0-(O.bodies[m0].state.pos[0]-O.bodies[mm].state.pos[0])
	epsa=epsa0-(O.bodies[topPlate].state.pos[2]-O.bodies[bottomPlate].state.pos[2])
	epsv=epsa+2*epsr
	if((displ0-displ)>0.008) and (load==True):
		O.pause()
		#print 'Real time = ', O.realtime
		print 'end of loading, O.realtime (min) = ', O.realtime/60.
		O.bodies[bottomPlate].state.vel=(0,0,0) 
		O.bodies[topPlate].state.vel=(0,0,0)
		saveData()
		#O.exitNoBacktrace()
	plot.addData(t=O.time,t2=O.time,displ=displ,d=(displ0-displ)/displ0,pos=pos,Fnt=Fnt,Fnb=Fnb,sigmaN=sigmaN,v=O.bodies[topPlate].state.vel[2],cui=cui,sigmaNb=sigmaNb,unbF=unbalancedForce(),p=p,q=q,p0=p0,sigma0=sigma0,epsv=epsv)

def history():
	plot.addData(e11=-triax.strain[0], e22=-triax.strain[1], e33=-triax.strain[2],
		ev=-triax.strain[0]-triax.strain[1]-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)

def addPlotData():
  fMove = Vector3(0,0,0)
  
  for i in idTop:
    fMove += O.forces.f(i)
  
  plot.addData(z=O.iter, pMove=fMove[2], pFest=fMove[2])

def measure():
	global qsMean,vxPartPY,phiPartPY
	#Evaluate the average depth profile of streamwise, spanwise and wall-normal particle velocity, particle volume fraction (and drag force for coupling with RANS fluid resolution), and store it in hydroEngine variables vxPart, phiPart, vyPart, vzPart, averageDrag.
	hydroEngine.averageProfile()
	#Extract the calculated vector. They can be saved and plotted afterwards. 
	vxPartPY = np.array(hydroEngine.vxPart)
	phiPartPY = np.array(hydroEngine.phiPart)

	#Evaluate the dimensionless sediment transport rate for information
	qsMean = sum(phiPartPY*vxPartPY)*dz/sqrt((densPart/densFluidPY - 1)*abs(gravityVector[2])*pow(diameterPart,3))
	plot.addData(SedimentRate = qsMean, time = O.time)	#Plot it during the simulation

	#Condition to stop the simulation after endTime seconds
	if O.time>=endTime:
		print('\n End of the simulation, simulated {0}s as required !\n '.format(endTime))
		O.pause()

	#Evaluate the Shields number from the maximum of the Reynolds stresses evaluated in the fluid resolution
	shieldsNumber = max(hydroEngine.ReynoldStresses)/((densPart-densFluidPY)*diameterPart*abs(gravityVector[2]))	
	print('Shields number', shieldsNumber)

	if saveData==1:	#Save data for postprocessing
		global fileNumber
		nameFile = scriptPath + '/sim'+ str(nbSim) + '/data/'+ str(fileNumber)+'.py'	# Name of the file that will be saved
		globalParam =  ['qsMean','phiPartPY','vxPartPY','vxFluidPY','zAxis']	# Variables to save
		Save(nameFile, globalParam)	#Save
		fileNumber+=1	#Increment the file number

	#Activate the fluid wall friction only at equilibrium. Not necessary for the transient.
	if O.time>10:
		hydroEngine.fluidWallFriction = True

def ploteo():
    a1=0
    a2=0
    a3=0
    a4=0
    a5=0
    a6=0
    b=0
    for i in range(5, n1+5):
        b=O.forces.f(i)
        a1=a1+math.sqrt(numpy.dot(b, b))
    for i in range(n1+5, n1+n2+5):
        b=O.forces.f(i)
        a2=a2+math.sqrt(numpy.dot(b, b))
    for i in range(n1+n2+5, n1+n2+n3+5):
        b=O.forces.f(i)
        a3=a3+math.sqrt(numpy.dot(b, b))
    for i in range(n1+n2+n3+5, n1+n2+n3+n4+5):
        b=O.forces.f(i)
        a4=a4+math.sqrt(numpy.dot(b, b))
    for i in range(n1+n2+n3+n4+5, n1+n2+n3+n4+n5+5):
        b=O.forces.f(i)
        a5=a5+math.sqrt(numpy.dot(b, b))
    for i in range(n1+n2+n3+n4+5, n1+n2+n3+n4+n5+5):
        b=O.forces.f(i)
        a5=a5+math.sqrt(numpy.dot(b, b))
    for i in range(n1+n2+n3+n4+n5+5, n1+n2+n3+n4+n5+n6+5):
        b=O.forces.f(i)
        a6=a6+math.sqrt(numpy.dot(b, b))

    fino=a2/n2+a4/n4+a6/n6
    gross=a1/n1+a3/n3+a5/n5
    G=fino/(fino+gross)
    time1=O.iter
    plot.addData(Gt=G,i=O.iter, i1=log10(time1))

def addPlotData():
	Fn = 0.
	Fnn = 0.
	if(axis==1):
		for i in posIds:
			Fn +=	O.forces.f(i)[axis]
			
		Fn=abs(Fn)	
		sigma=Fn/(1000*L)
		
		un = (O.bodies[posIds[0]].state.pos[axis] - O.bodies[posIds[0]].state.refPos[axis])
		eps=un/L
	
	else:
		for i in rightNodes:
			Fn += O.forces.f(i)[axis]
		Fn=abs(Fn)	
		sigma=Fn/(1000*L)
		
		un = (O.bodies[rightNodes[0]].state.pos[axis] - O.bodies[rightNodes[0]].state.refPos[axis])
		eps=un/L
	if un*1000 > 300:
		O.pause()
	#print eps	
	plot.addData( eps=eps*1000, un=un*1000, Fn=Fn,sigma=sigma )

def addData():
   stress=sum(normalShearStressTensors(),Matrix3.Zero)
   sigzz=stress[2,2]
   q1T=0.0
   q1N=0
   q2T=0.0
   q2N=0
   q3T=0.0
   q3N=0
   q4T=0.0
   q4N=0
   for i in O.bodies:
      if i.isClumpMember==True and i.state.pos[2]<0.25*O.cell.size[2]:
   	q1T+=i.state.temp
   	q1N+=1
      elif i.isClumpMember==True and i.state.pos[2]<0.5*O.cell.size[2]:
   	q2T+=i.state.temp
   	q2N+=1
      elif i.isClumpMember==True and i.state.pos[2]<0.75*O.cell.size[2]:
   	q3T+=i.state.temp
   	q3N+=1
      elif i.isClumpMember==True:
   	q4T+=i.state.temp
   	q4N+=1
   plot.addData(szz=sigzz,heat=totHeat,T1=q1T/q1N,T2=q2T/q2N,T3=q3T/q3N,T4=q4T/q4N,t=O.time)

def ploteo():
    a1=0.
    a2=0.
    a3=0.
    d1=0.
    d2=0.
    cant1=0
    cant2=0
    cant3=0
    actual=O.iter
    inicial=graf.firstIterRun
    velMaxFino=0.
    velMaxGros=0.
    velMeanFino=0.
    velMeanGros=0.
    time=O.engines[4].timeStepUpdateInterval*O.engines[4].previousDt
    for i in range(6, n1+6):
        b=O.forces.f(i)
        rii=O.bodies[i].shape.radius
        vel=O.bodies[i].state.vel
        velScalar=math.sqrt(numpy.dot(vel, vel))
        if 2*rii<D15:
         a2=a2+math.sqrt(numpy.dot(b, b))
         velFino[cant2]=velScalar
         d2=d2+velScalar
         cant2=cant2+1
        if 2*rii>D85:
         a1=a1+math.sqrt(numpy.dot(b, b))
         velGros[cant1]=velScalar
         d1=d1+velScalar
         cant1=cant1+1
        else: 
         a3=a3+math.sqrt(numpy.dot(b, b))
#         despRest[cant3]=deltaScalar
         cant3=cant3+1
    if actual<=inicial:
     __builtin__.time1=time
    else:     
     __builtin__.time1=time+__builtin__.time2
    velMaxFino=numpy.amax(velFino)
    velMeanFino=d2/cant2
    velMaxGros=numpy.amax(velGros)
    velMeanGros=d1/cant1
    gross=a1*1.
    fino=a2*1.
    resto=a3*1.
    Gfp=fino/(cant2*1.)/(fino/(cant2*1.)+gross/(cant1*1.)+resto/(cant3*1.))
    __builtin__.G=Gfp
    e=flow.porosity/(1-flow.porosity)
    Sf=0.15
    nf=e/Sf
    sigma=O.forces.f(5)[2]/(D*D)
    dz=O.bodies[5].bound.refPos[2] 
    __builtin__.icr=G/(rhoh*dz*g)*(sigma*tan(radians(angfric)))*1.+nf*rhos/rhoh
    it=flow.bndCondValue[4]/(O.bodies[5].bound.refPos[2]*g*rhoh)
    print "velFino",velMeanFino," Cant2",cant2
    plot.addData(Vmf=velMeanFino, VMaxf=velMaxFino,Vmg=velMeanGros, VMaxg=velMaxGros, i=__builtin__.time1, i1=__builtin__.time1 , i2=O.iter, ic=__builtin__.icr, ii=it)
    cant1=0
    cant2=0
    cant3=0

def defData():
	plot.addData(fy=O.forces.f(3)[1], # vertical component of the force sustained by the upper side of the shear box
		    fx=O.forces.f(3)[0], # horizontal component of the force sustained by the upper side of the shear box
		    step=O.iter,
		    gamma=O.bodies[3].state.pos[0] - length/2.0,
		    u=O.bodies[3].state.pos[1] - (height+thickness/2.0)
		    )

def measure():
	global qsMean
	global vxPartPY
	global phiPartPY
	global zAxis
	#Evaluate the average depth profile of streamwise, spanwise and wall-normal particle velocity, particle volume fraction (and drag force for coupling with RANS fluid resolution), and store it in hydroEngine variables vxPart, vyPart, vzPart, phiPart, averageDrag
	hydroEngine.averageProfile()
	#Extract the calculated vector. They can be saved and plotted afterwards. 
	vxPartPY = np.array(hydroEngine.vxPart)
	vyPartPY = np.array(hydroEngine.vyPart)
	vzPartPY = np.array(hydroEngine.vzPart)
	phiPartPY = np.array(hydroEngine.phiPart)
	averageDragPY = np.array(hydroEngine.averageDrag)

	#Evaluate the dimensionless sediment transport rate for information
	qsMean = sum(phiPartPY*vxPartPY)*dz/sqrt((densPart/densFluidPY - 1)*abs(gravityVector[2])*pow(diameterPart,3))
	plot.addData(SedimentRate = qsMean, time = O.time)	#Plot it during the simulation

	#Condition to stop the simulation after endTime seconds
	if O.time>=endTime:
		print('\n End of the simulation, simulated {0}s as required !\n '.format(endTime))
		O.pause()

		#Z scale used for the possible plot at the end
		global zAxis
		for i in range(0,ndimz):
			zAxis[i] = i*dz/diameterPart

def defData():
	i=O.interactions[1,0]
	vecFn=i.phys.normalForce
	vecDist=upperSphere.state.pos-lowerSphere.state.pos
	plot.addData(normFn=vecFn.norm(),normFnBis=vecFn.norm(),fnY=vecFn[1],step=O.iter,
	  unPerso=lowerSphere.shape.radius+upperSphere.shape.radius-vecDist.norm(),unTrue=i.geom.penetrationDepth,
	  gamma=upperSphere.state.pos[0]-lowerSphere.state.pos[0],fx=O.forces.f(0)[0],torque=O.forces.t(1)[2])

def addPlotData():
	plot.addData(unbalanced=utils.unbalancedForce(),i=O.iter,
		sxx=triax.stress[0],syy=triax.stress[1],szz=triax.stress[2],
		exx=triax.strain[0],eyy=triax.strain[1],ezz=triax.strain[2],
		# save all available energy data
		Etot=O.energy.total(),**O.energy
	)

def plotAddData():
	i = O.interactions[0,1]
	if i.phys:
		plot.addData(
			fn = i.phys.normalForce.norm(),
			dspl = O.bodies[1].state.displ().norm(),
		)

def measure():
	global qsMean,vxPartPY,phiPartPY

	#Evaluate the average depth profile of streamwise, spanwise and wall-normal particle velocity, particle volume fraction (and drag force for coupling with RANS fluid resolution), and store it in hydroEngine variables vxPart, vyPart, vzPart, phiPart, averageDrag
	hydroEngine.averageProfile()
	#Extract the calculated vector. They can be saved and plotted afterwards. 
	vxPartPY = np.array(hydroEngine.vxPart)
	vyPartPY = np.array(hydroEngine.vyPart)
	vzPartPY = np.array(hydroEngine.vzPart)
	phiPartPY = np.array(hydroEngine.phiPart)
	averageDragPY = np.array(hydroEngine.averageDrag)

	#Evaluate the dimensionless sediment transport rate for information
	qsMean = sum(phiPartPY*vxPartPY)*dz/sqrt((densPart/densFluidPY - 1)*abs(gravityVector[2])*pow(diameterPart,3))
	plot.addData(SedimentRate = qsMean, time = O.time)	#Plot it during the simulation

	#Condition to stop the simulation after endTime seconds
	if O.time>=endTime:
		print('\n End of the simulation, simulated {0}s as required !\n '.format(endTime))
		O.pause()

		#Z scale used for the possible plot at the end
		global zAxis
		for i in range(0,ndimz):
			zAxis[i] = i*dz/diameterPart


	if saveData==1:	#Save data for postprocessing
		global fileNumber
		nameFile = scriptPath + '/data/'+ str(fileNumber)+'.py'	# Name of the file that will be saved
		globalParam =  ['qsMean','phiPartPY','vxPartPY','vxFluidPY','zAxis']	# Variables to save
		Save(nameFile, globalParam)	#Save
		fileNumber+=1	#Increment the file number