wall.append(O.bodies.append(pfacet) in triaxial test

I set Pfacet attribute wire=true to check body moving better, but it seems that the triax engine doesn‘t work. So is it workable to use Pfacet simulating triaxial test by set "wall.append(O.bodies.append(pfacet(nodes[0],nodes[1],nodes[2],wire=True,material='pFacetMat',color=color)))"? Thanks!!

Hi,
I would suggest to use boxes for the walls since they are completely flat. It would minimize the number of problems for a start.

Using pfacets should work but not the way you do it.
First problem: triax.wallIds is undefined, so triax engine doesn't know where to find walls .
Second problem: you can't impose the motion of a PFacet, it always follow the motion of the nodes. One exception is if you clump facet+nodes into a clump for each of the walls. Then the clumps' ids should be used to define wallIds.
I hope it helps
Bruno

Hi，Bruno.
> I would suggest to use boxes for the walls since they are completely flat. It would minimize the number of problems for a start.
I have tried boxes for the walls and it worked well except boxes can't interact with cylinder, and that's why I use Pfacet for walls.

> First problem: triax.wallIds is undefined, so triax engine doesn't know where to find walls .
So I should define a wallIDs list then append Pfacet-wall to it?

> Second problem: you can't impose the motion of a PFacet, it always follow the motion of the nodes. One exception is if you clump facet+nodes into a clump for each of the walls. Then the clumps' ids should be used to define wallIds.
One question, there are 8 nodes and I need to generate 12 Pfacets, so every node would be reused if I clump pfacet+nodes into a clump for each of the walls. Is that workable?

Besides, I find only Ig2_GridConnection_PFacet_ScGeom can generate contact between cylinder and other kind of wall-like body, and I have been working on fiber-reinforced sand so cylinder is the kind of body that I must use. Don't really know is there any other body that I can use to simulate wall?
Thanks for your help!

> -----原始邮件-----
> 发件人: "Bruno Chareyre" <email address hidden>
> 发送时间: 2021-01-27 02:41:01 (星期三)
> 收件人: <email address hidden>
> 抄送:
> 主题: Re: [Question #695192]: wall.append(O.bodies.append(pfacet) in triaxial test
>
> Your question #695192 on Yade changed:
> https://answers.launchpad.net/yade/+question/695192
>
> Status: Open => Answered
>
> Bruno Chareyre proposed the following answer:
> Hi,
> I would suggest to use boxes for the walls since they are completely flat. It would minimize the number of problems for a start.
>
> Using pfacets should work but not the way you do it.
> First problem: triax.wallIds is undefined, so triax engine doesn't know where to find walls .
> Second problem: you can't impose the motion of a PFacet, it always follow the motion of the nodes. One exception is if you clump facet+nodes into a clump for each of the walls. Then the clumps' ids should be used to define wallIds.
> I hope it helps
> Bruno
>
> --
> If this answers your question, please go to the following page to let us
> know that it is solved:
> https://answers.launchpad.net/yade/+question/695192/+confirm?answer_id=1
>
> If you still need help, you can reply to this email or go to the
> following page to enter your feedback:
> https://answers.launchpad.net/yade/+question/695192
>
> You received this question notification because you asked the question.

Hi，Bruno.
> I would suggest to use boxes for the walls since they are completely flat. It would minimize the number of problems for a start.
I have tried boxes for the walls and it worked well except boxes can't interact with cylinder, and that's why I use Pfacet for walls.

> First problem: triax.wallIds is undefined, so triax engine doesn't know where to find walls .
So I should define a wallIDs list then append Pfacet-wall to it?

> Second problem: you can't impose the motion of a PFacet, it always follow the motion of the nodes. One exception is if you clump facet+nodes into a clump for each of the walls. Then the clumps' ids should be used to define wallIds.
One question, there are 8 nodes and I need to generate 12 Pfacets, so every node would be reused if I clump pfacet+nodes into a clump for each of the walls. Is that workable?

Besides, I find only Ig2_GridConnection_PFacet_ScGeom can generate contact between cylinder and other kind of wall-like body, and I have been working on fiber-reinforced sand so cylinder is the kind of body that I must use. Don't really know is there any other body that I can use to simulate wall?
Or, can I just impose the motion of a Pfacet without triax engine?Is there any attribute like O.area.force(I make up it)?
My idea is that:
(1)Generate sample by Pfacets pushing spheres and cylinders till they achieve some certain positions(in a smaller box)
(2)triaxial test on sample.

Thanks for your help!

> I have tried boxes for the walls and it worked well except boxes can't interact with cylinder

Sure of that? Any example script?

> One question, there are 8 nodes and I need to generate 12 Pfacets, so every node would be reused if I clump pfacet+nodes into a clump for each of the walls.

I would expect 4 nodes / 2 facets for a rectangle, but whatever the numbers you can clump them all together in one body, hence the nodes won't be repeated.

Bruno

>Sure of that? Any example script?

Just replace PFacet code with:
walls=aabbWalls([mi,ma],thickness=0,material='walls')
wallIds=O.bodies.append(walls)
(try script below)

>I would expect 4 nodes / 2 facets for a rectangle, but whatever the numbers you can clump them all together in one body, hence the nodes won't be repeated.

4 nodes / 2 facets for a rectangle.
a box has 8 vertices and 6 sides(12 Pfacets).
hence 8 nodes would be repeated?

Can I just add force on Pfacets?

Thanks!
#################################
#run and you can see cylinders go through wall like there is not a wall

from yade.gridpfacet import *
from yade import pack, plot
#from random import random as r
import numpy as np
from numpy import *
import math
from yade import utils, qt
import random
from yade import polyhedra_utils, qt

#parameters

#mi,ma = (-100e-3,-100e-3,-200e-3),(100e-3,100e-3,200e-3)
mi,ma = (-10e-3,-10e-3,-20e-3),(10e-3,10e-3,20e-3)

color=[0.,1.,1.]

#### Parameter ####
young=4.0e6
poisson=3
density=1e3
stabilityThreshold=0.01
#=============================meterials========================================
O.materials.append(FrictMat(young=4.0e6,poisson=.3,frictionAngle=30,density=2630e0,label='sphereMat'))#for sphere

O.materials.append(FrictMat(young=3e9,poisson=.15,frictionAngle=20,density=910e+0,label='extcylMat'))#for sphere-cylinder
O.materials.append(CohFrictMat(young=3e9,poisson=.15,density=910e0,frictionAngle=20,normalCohesion=1e40,shearCohesion=1e40,momentRotationLaw=True,label='intcylMat'))#for cylinder-cylinder

O.materials.append( CohFrictMat( young=3e6,poisson=0.15,density=910e0,frictionAngle=20,normalCohesion=3e100,shearCohesion=3e100,momentRotationLaw=True,label='gridNodeMat' ) )#for gridNodes
#O.materials.append(CohFrictMat(young=3e9,poisson=.15,density=910e6,frictionAngle=20,normalCohesion=1e40,shearCohesion=1e40,momentRotationLaw=True,label='gridNodeMat'))#for gridNodes
O.materials.append(FrictMat(young=4e6,poisson=0.3,density=1000,frictionAngle=20,label='pFacetMat')) #for pfacet

O.materials.append(FrictMat(young=young,poisson=0.5,frictionAngle=0,density=2630,label='walls'))

triax=TriaxialStressController(

 maxMultiplier=1.+2e4/young, # spheres growing factor (fast growth)
 finalMaxMultiplier=1.+2e3/young, # spheres growing factor (slow growth)
 thickness = 0,

 stressMask = 0,
 max_vel = 3.7e-7,
 wall_front_activated=True,
 wall_back_activated=True,
 wall_top_activated=True,
 wall_bottom_activated=True,
 wall_left_activated=True,
 wall_right_activated=True,

 internalCompaction=False, # If true the confining pressure is generated by growing particles
)

#==============================Engines=========================================
O.engines=[
 ForceResetter(),
 InsertionSortCollider([
  Bo1_Sphere_Aabb(),
  Bo1_GridConnection_Aabb(),
  Bo1_PFacet_Aabb(),
  Bo1_Box_Aabb(),
  Bo1_Wall_Aabb()
 ]),
 InteractionLoop([
  Ig2_Sphere_Sphere_ScGeom(),
  Ig2_Box_Sphere_ScGeom(),
  Ig2_GridNode_GridNode_GridNodeGeom6D(),
  Ig2_Sphere_GridConnection_ScGridCoGeom(),
  Ig2_GridConnection_GridConnection_GridCoGridCoGeom(),
  Ig2_GridConnection_PFacet_ScGeom(),
  Ig2_PFacet_PFacet_ScGeom(),
  Ig2_Sphere_PFacet_ScGridCoGeom(),
  Ig2_Wall_PFacet_ScGeom(),
  Ig2_Wall_Sphere_ScGeom()
 ],
 [
  Ip2_CohFrictMat_CohFrictMat_CohFrictPhys(setCohesionNow=True,setCohesionOnNewContacts=False), # internal cylinder physics
  Ip2_FrictMat_FrictMat_FrictPhys(), # physics for external interactions, i.e., cylinder-cylinder, sphere-sphere, cylinder-sphere

  Ip2_FrictMat_CpmMat_FrictPhys(),
 ],
 [
  Law2_ScGeom_FrictPhys_CundallStrack(), # contact law for sphere-sphere
  Law2_ScGridCoGeom_FrictPhys_CundallStrack(), # contact law for cylinder-sphere
  Law2_ScGeom6D_CohFrictPhys_CohesionMoment(), # contact law for "internal" cylinder forces
  Law2_GridCoGridCoGeom_FrictPhys_CundallStrack(), # contact law for cylinder-cylinder interaction
  #Law2_GridCoGridCoGeom_CohFrictPhys_CundallStrack(),
  Law2_ScGeom_CpmPhys_Cpm(),
 ]
 ),
 GlobalStiffnessTimeStepper(timestepSafetyCoefficient=0.1,label='ts'),
 triax,
 NewtonIntegrator(gravity=(0,0,20),damping=0.5,label='newton'),
 ]

triax.goal1=triax.goal2=triax.goal3=-0.9
walls=aabbWalls([mi,ma],thickness=0,material='walls')
wallIds=O.bodies.append(walls)
#=======================================sphere=================================================
sp = yade.pack.SpherePack()
sp.makeCloud(mi,ma,porosity=0.66/1.66,psdSizes=[0.5e-3,0.75e-3,0.8e-3,0.9e-3,1.2e-3,2e-3],psdCumm=[0.,0.2,0.4,0.6,0.8,1.0],num=133)
spheres=sp.toSimulation(color=(0,0.3,0.7),material='sphereMat')

#=======================================fiber=================================================
length=12e-3
diameter=0.3e-3
Xw=0.25e-2 # fiber content
fiber_vol=(80e-3*(40e-3*0.5)**2*np.pi)*Xw
Vfiber=length*np.pi*(0.5*diameter)**2
numFiber=fiber_vol/Vfiber
print "the numer of fibres =",int(numFiber)
Ne=int(length/(0.9e-3))
print 'Ne=',Ne
nodesIds=[]
cylIds=[]
numFibre=0
target_num=int(numFiber)
target_num=int(36)
####fiberfilewrite
fiberwr=open('fibers.txt','w')
cylNodes=[]
time=0
HalleluYah=0
for n in range(1000):
    if numFibre<target_num:
  random.seed(n*5)
  z0=random.uniform(-20e-3, 20e-3 )
  random.seed( )
  sita0=random.uniform(-np.pi,np.pi)
  random.seed( n+7)
  l_y=random.uniform(-10e-3,10e-3)
  x0=l_y*sin(sita0)
  y0=l_y*cos(sita0)

  random.seed(n+8)
  sita=random.uniform(-0.5*np.pi, 0.5*np.pi)
  random.seed()
  phi=random.uniform(0, 2*np.pi)

  hz=length*sin(sita)
  hx=length*cos(sita)*cos(phi)
  hy=length*cos(sita)*sin(phi)

  x1=x0+hx
  y1=y0+hy
  z1=z0+hz

  if 0<abs(x1)<10e-3 and 0<abs(y1)<10e-3 and 0<abs(z1)<20e-3 and 0<abs(x0)<10e-3 and 0<abs(y0)<10e-3 and 0<abs(z0)<20e-3 :
   if len(cylNodes)<=0:
       fiberwr.write(str(numFibre)+'\t'+str(x0)+'\t'+str(y0)+'\t'+str(z0)+'\t'+str(x1)+'\t'+str(y1)+'\t'+str(z1)+'\n')
       cylNodes.append([x0,y0,z0,x1,y1,z1])
       numFibre+=1
       vertices=[]
       for i in range(0, Ne+1):
    px=float(i)*hx/float(Ne)+x0; py=float(i)*hy/float(Ne)+y0; pz=float(i)*hz/float(Ne)+z0;
    vertices.append([px,py,pz])
       cylinderConnection(vertices,0.15e-3,nodesIds,cylIds,color=[0.5,0.5,0],fixed=False,highlight=False,intMaterial='intcylMat',extMaterial='extcylMat')
   else :
    valid=True
    for c in cylNodes:
     if valid:
      A1=np.array([x0,y0,z0])
      B1=np.array([x1,y1,z1])
      r1=0.15e-3
      A2=np.array([c[0],c[1],c[2]])
      B2=np.array([c[3],c[4],c[5]])
      r2=0.15e-3
      l1=B1-A1
      l2=B2-A2
      l3=A2-A1

      delta1=(np.dot(l1,l2))/(np.linalg.norm(l1)*np.linalg.norm(l2))

      if abs(abs(delta1)-1)<0.0000000000000001:
           p=np.linalg.norm(np.cross(l2,l3))/np.linalg.norm(l2)
           if p<r1+r2 :
        print("PARALLEL INTERSECT!")
        valid=False
           else :
        print("PARALLEL DO NOT INTERSECT!")
        valid=True
      else :
       l4=np.cross(l1,l2)
       ll1=np.cross(l1,l4)
           ll2=np.cross(l2,l4)

           A11=np.matrix([ll1,[l2[1],-l2[0],0.0],[0.0,l2[2],-l2[1]]])
           x11=np.dot(ll1,A1)
           x12=l2[1]*A2[0]-l2[0]*A2[1]
           x13=l2[2]*A2[1]-l2[1]*A2[2]
           b1=np.matrix([[x11],[x12],[x13]])
           C1=np.linalg.solve(A11,b1)

           A21=np.matrix([ll2,[l1[1],-l1[0],0.0],[0.0,l1[2],-l1[1]]])
           x21=np.dot(ll2,A2)
           x22=l1[1]*A1[0]-l1[0]*A1[1]
           x23=l1[2]*A1[1]-l1[1]*A1[2]
           b2=np.matrix([[x21],[x22],[x23]])
           C2=np.linalg.solve(A21,b2)

           t1=(C1[0]-A1[0])/(B1[0]-A1[0])
           t2=(C2[0]-A2[0])/(B2[0]-A2[0])

           if t1<0 :
        E=A1
           elif t1>=0 and t1<=1 :
        E=C1
           else :
        E=B1

           if t2<0 :
        F=A2
           elif t2>=0 and t2<=1 :
        F=C2
           else :
        F=B2
           D=np.linalg.norm(E-F)

           if D<(r1 + r2) :
        print("INTERSECT!")
        valid=False

           else :
        valid=True

    if valid:
     print(n,valid)
     fiberwr.write(str(numFibre)+'\t'+str(x0)+'\t'+str(y0)+'\t'+str(z0)+'\t'+str(x1)+'\t'+str(y1)+'\t'+str(z1)+'\n')
         cylNodes.append([x0,y0,z0,x1,y1,z1])
         numFibre+=1
         vertices=[]
        for i in range(0, Ne+1):
      px=float(i)*hx/float(Ne)+x0; py=float(i)*hy/float(Ne)+y0; pz=float(i)*hz/float(Ne)+z0;
      vertices.append([px,py,pz])
        cylinderConnection(vertices,0.15e-3,nodesIds,cylIds,color=[0.5,0.5,0],fixed=False,highlight=False,intMaterial='intcylMat',extMaterial='extcylMat')
        #print(n,"A1",x0,y0,z0, "B1",x1,y1,z1)

fiberwr.close()

I got problems when clumping nodes and pfacets into a clump. There are several ways to creat Pfacets and I tried two of them, both did't work. Here is one way I tried. Please help me with it. Thanks!

####################################
nodes=[]
#up
aa= nodes.append(O.bodies.append( gridNode([-110e-3, 110e-3,-200e-3],RNode,wire=True,fixed=False,material='gridNodeMat') ))
bb=nodes.append(O.bodies.append( gridNode([ 110e-3, 110e-3,-200e-3],RNode,wire=True,fixed=False,material='gridNodeMat') ))
cc= nodes.append(O.bodies.append( gridNode([ 110e-3,-110e-3,-200e-3],RNode,wire=True,fixed=False,material='gridNodeMat') ))
dd= nodes.append(O.bodies.append( gridNode([-110e-3,-110e-3,-200e-3],RNode,wire=True,fixed=False,material='gridNodeMat') ))

ab= O.bodies.append(gridConnection(aa,bb,RNode,wire=True,material='gridNodeMat'))
ac= O.bodies.append(gridConnection(aa,cc,RNode,wire=True,material='gridNodeMat'))
ad= O.bodies.append(gridConnection(aa,dd,RNode,wire=True,material='gridNodeMat'))
bc= O.bodies.append(gridConnection(bb,cc,RNode,wire=True,material='gridNodeMat'))
bd= O.bodies.append(gridConnection(bb,dd,RNode,wire=True,material='gridNodeMat'))

abc= O.bodies.append(pfacet(aa,bb,cc,wire=True,material='pFacetMat',color=color))
acd= O.bodies.append(pfacet(aa,cc,dd,wire=True,material='pFacetMat',color=color))

clump1=O.bodies.appendClumped([\
aa, bb, cc, dd, ab, ac, ad, bc, bd, abc, acd
])

###################################################
Thanks!

And the error says:

  File "1-28.py", line 124, in <module>
    ab= O.bodies.append(gridConnection(aa,bb,RNode,wire=True,material='gridNodeMat'))
  File "/usr/lib/x86_64-linux-gnu/yade/py/yade/gridpfacet.py", line 109, in gridConnection
    sph1=O.bodies[id1] ; sph2=O.bodies[id2]
ArgumentError: Python argument types in
    BodyContainer.__getitem__(BodyContainer, NoneType)
did not match C++ signature:
    __getitem__(pyBodyContainer {lvalue}, int)

Thanks!

That's too many questions for me, sorry.
Simple answer: you did not show a problem with Box's, you showed a problem with Wall's (aabbWalls())
Walls are not Boxes. :)
This being said I would have expected this to work with walls too.

Could you please show the pfacet-box or pfacet-wall problem in a MWE [1]?

Bruno

[1] https://www.yade-dem.org/wiki/Howtoask

Thanks for your patience! I have read the article about Pfacet and I think it is what I want for wall as to act as soft membrane.And I tried to clump nodes and pfacs as you suggest but it didn't work. Could you help me with it?(I have commented on this question earlier) Thanks! zhanghanying 邮箱：<email address hidden> 签名由 网易邮箱大师 定制 On 01/28/2021 23:55, Bruno Chareyre wrote: Your question #695192 on Yade changed: https://answers.launchpad.net/yade/+question/695192    Status: Open => Needs information Bruno Chareyre requested more information: That's too many questions for me, sorry. Simple answer: you did not show a problem with Box's, you showed a problem with Wall's (aabbWalls()) Walls are not Boxes. :) This being said I would have expected this to work with walls too. Could you please show the pfacet-box or pfacet-wall problem in a MWE [1]? Bruno [1] https://www.yade-dem.org/wiki/Howtoask -- To answer this request for more information, you can either reply to this email or enter your reply at the following page: https://answers.launchpad.net/yade/+question/695192 You received this question notification because you asked the question.

Thanks Jérôme Duriez, that solved my question.