~ubuntu-branches/debian/sid/yade/sid

« back to all changes in this revision

Viewing changes to examples/jointedCohesiveFrictionalPM/identifBis.py

  • Committer: Package Import Robot
  • Author(s): Anton Gladky
  • Date: 2014-01-13 20:13:14 UTC
  • mfrom: (1.1.9)
  • Revision ID: package-import@ubuntu.com-20140113201314-4x3stisjce4t0pdx
Tags: 1.07.0-1
* [d42c7de] Imported Upstream version 1.07.0
* [a421b10] Remove patch, applied by upstream.
* [d151f88] Set Standards-Version: 3.9.5. No changes.
* [a9616c6] Add upstream file.
* [c078bb4] Inject additional parameters for weak archs. 
            Thanks to Roland Stigge <stigge@antcom.de>. (Closes: #733152)
* [21a2430] Fix version definition for IPython>1.0.0.
* [8204093] Remove google-scripts from documentation.

Show diffs side-by-side

added added

removed removed

Lines of Context:
13
13
        InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=interactionRadius,label='is2aabb'),Bo1_Facet_Aabb()]),
14
14
        InteractionLoop(
15
15
                [Ig2_Sphere_Sphere_ScGeom(interactionDetectionFactor=interactionRadius,label='ss2d3dg'),Ig2_Facet_Sphere_ScGeom()],
16
 
                [Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(cohesiveTresholdIteration=1,alpha=1,tensileStrength=1e6,cohesion=1e6,jointNormalStiffness=1,jointShearStiffness=1,jointTensileStrength=1e6,jointCohesion=1e6,jointFrictionAngle=1,label='interactionPhys')],
 
16
                [Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(cohesiveTresholdIteration=1,label='interactionPhys')],
17
17
                [Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(smoothJoint=True,label='interactionLaw')]
18
18
        ),
19
19
        NewtonIntegrator(damping=1)
47
47
        nRef=normalRef/(normalRef.norm()) ## normalizes normalRef
48
48
        normalFacetSphere=i.geom.normal # geom.normal is oriented from id1 to id2 -> normalFacetSphere from facet (i.id1) to sphere (i.id2)
49
49
 
50
 
        if O.bodies[i.id2].mat.onJoint==False : ## particles has not yet been identified as belonging to a joint plane
51
 
            O.bodies[i.id2].mat.onJoint=True
52
 
            O.bodies[i.id2].mat.joint=1
 
50
        if O.bodies[i.id2].state.onJoint==False : ## particles has not yet been identified as belonging to a joint plane
 
51
            O.bodies[i.id2].state.onJoint=True
 
52
            O.bodies[i.id2].state.joint=1
53
53
            O.bodies[i.id2].shape.color=jointcolor1
54
54
            if nRef.dot(normalFacetSphere)>=0 :
55
 
                O.bodies[i.id2].mat.jointNormal1=nRef
 
55
                O.bodies[i.id2].state.jointNormal1=nRef
56
56
            elif nRef.dot(normalFacetSphere)<0 :
57
 
                O.bodies[i.id2].mat.jointNormal1=-nRef
58
 
        elif O.bodies[i.id2].mat.onJoint==True :  ## particles has already been identified as belonging to, at least, 1 facet
59
 
            if O.bodies[i.id2].mat.joint==1 and ((O.bodies[i.id2].mat.jointNormal1.cross(nRef)).norm()>0.05) : ## particles has already been identified as belonging to only 1 facet
60
 
                O.bodies[i.id2].mat.joint=2
 
57
                O.bodies[i.id2].state.jointNormal1=-nRef
 
58
        elif O.bodies[i.id2].state.onJoint==True :  ## particles has already been identified as belonging to, at least, 1 facet
 
59
            if O.bodies[i.id2].state.joint==1 and ((O.bodies[i.id2].state.jointNormal1.cross(nRef)).norm()>0.05) : ## particles has already been identified as belonging to only 1 facet
 
60
                O.bodies[i.id2].state.joint=2
61
61
                O.bodies[i.id2].shape.color=jointcolor2
62
62
                if nRef.dot(normalFacetSphere)>=0 :
63
 
                    O.bodies[i.id2].mat.jointNormal2=nRef
 
63
                    O.bodies[i.id2].state.jointNormal2=nRef
64
64
                elif nRef.dot(normalFacetSphere)<0 :
65
 
                    O.bodies[i.id2].mat.jointNormal2=-nRef
66
 
            elif O.bodies[i.id2].mat.joint==2 and ((O.bodies[i.id2].mat.jointNormal1.cross(nRef)).norm()>0.05) and ((O.bodies[i.id2].mat.jointNormal2.cross(nRef)).norm()>0.05): ## particles has already been identified as belonging to more than 1 facet
67
 
                O.bodies[i.id2].mat.joint=3
 
65
                    O.bodies[i.id2].state.jointNormal2=-nRef
 
66
            elif O.bodies[i.id2].state.joint==2 and ((O.bodies[i.id2].state.jointNormal1.cross(nRef)).norm()>0.05) and ((O.bodies[i.id2].state.jointNormal2.cross(nRef)).norm()>0.05): ## particles has already been identified as belonging to more than 1 facet
 
67
                O.bodies[i.id2].state.joint=3
68
68
                O.bodies[i.id2].shape.color=jointcolor3
69
69
                if nRef.dot(normalFacetSphere)>=0 :
70
 
                    O.bodies[i.id2].mat.jointNormal3=nRef
 
70
                    O.bodies[i.id2].state.jointNormal3=nRef
71
71
                elif nRef.dot(normalFacetSphere)<0 :
72
 
                    O.bodies[i.id2].mat.jointNormal3=-nRef
73
 
            elif O.bodies[i.id2].mat.joint==3 and ((O.bodies[i.id2].mat.jointNormal1.cross(nRef)).norm()>0.05) and ((O.bodies[i.id2].mat.jointNormal2.cross(nRef)).norm()>0.05) and ((O.bodies[i.id2].mat.jointNormal3.cross(nRef)).norm()>0.05):
74
 
                O.bodies[i.id2].mat.joint=4
 
72
                    O.bodies[i.id2].state.jointNormal3=-nRef
 
73
            elif O.bodies[i.id2].state.joint==3 and ((O.bodies[i.id2].state.jointNormal1.cross(nRef)).norm()>0.05) and ((O.bodies[i.id2].state.jointNormal2.cross(nRef)).norm()>0.05) and ((O.bodies[i.id2].state.jointNormal3.cross(nRef)).norm()>0.05):
 
74
                O.bodies[i.id2].state.joint=4
75
75
                O.bodies[i.id2].shape.color=jointcolor5
76
76
 
77
77
#### second step -> find spheres interacting with spheres on facet (could be executed in the same timestep as step 1?)
82
82
        vertices=O.bodies[j.id1].shape.vertices
83
83
        normalRef=vertices[0].cross(vertices[1]) # defines the normal to the facet normalRef
84
84
        nRef=normalRef/(normalRef.norm()) ## normalizes normalRef
85
 
        if ((O.bodies[j.id2].mat.jointNormal1.cross(nRef)).norm()<0.05) :
86
 
            jointNormalRef=O.bodies[j.id2].mat.jointNormal1
87
 
        elif ((O.bodies[j.id2].mat.jointNormal2.cross(nRef)).norm()<0.05) :
88
 
            jointNormalRef=O.bodies[j.id2].mat.jointNormal2
89
 
        elif ((O.bodies[j.id2].mat.jointNormal3.cross(nRef)).norm()<0.05) :
90
 
            jointNormalRef=O.bodies[j.id2].mat.jointNormal3
 
85
        if ((O.bodies[j.id2].state.jointNormal1.cross(nRef)).norm()<0.05) :
 
86
            jointNormalRef=O.bodies[j.id2].state.jointNormal1
 
87
        elif ((O.bodies[j.id2].state.jointNormal2.cross(nRef)).norm()<0.05) :
 
88
            jointNormalRef=O.bodies[j.id2].state.jointNormal2
 
89
        elif ((O.bodies[j.id2].state.jointNormal3.cross(nRef)).norm()<0.05) :
 
90
            jointNormalRef=O.bodies[j.id2].state.jointNormal3
91
91
        else : continue
92
92
        facetCenter=O.bodies[j.id1].state.pos
93
93
        #### seek for each sphere interacting with the identified sphere j.id2
100
100
                    sphOnF=n.id2
101
101
                    othSph=n.id1
102
102
                facetSphereDir=(O.bodies[othSph].state.pos-facetCenter)
103
 
                if O.bodies[othSph].mat.onJoint==True :
104
 
                    if O.bodies[othSph].mat.joint==3 and ((O.bodies[othSph].mat.jointNormal1.cross(jointNormalRef)).norm()>0.05) and ((O.bodies[othSph].mat.jointNormal2.cross(jointNormalRef)).norm()>0.05) and ((O.bodies[othSph].mat.jointNormal3.cross(jointNormalRef)).norm()>0.05):
105
 
                        O.bodies[othSph].mat.joint=4
 
103
                if O.bodies[othSph].state.onJoint==True :
 
104
                    if O.bodies[othSph].state.joint==3 and ((O.bodies[othSph].state.jointNormal1.cross(jointNormalRef)).norm()>0.05) and ((O.bodies[othSph].state.jointNormal2.cross(jointNormalRef)).norm()>0.05) and ((O.bodies[othSph].state.jointNormal3.cross(jointNormalRef)).norm()>0.05):
 
105
                        O.bodies[othSph].state.joint=4
106
106
                        O.bodies[othSph].shape.color=jointcolor5
107
 
                    elif O.bodies[othSph].mat.joint==2 and ((O.bodies[othSph].mat.jointNormal1.cross(jointNormalRef)).norm()>0.05) and ((O.bodies[othSph].mat.jointNormal2.cross(jointNormalRef)).norm()>0.05):
108
 
                        O.bodies[othSph].mat.joint=3
109
 
                        if facetSphereDir.dot(jointNormalRef)>=0:
110
 
                            O.bodies[othSph].mat.jointNormal3=jointNormalRef
111
 
                        elif facetSphereDir.dot(jointNormalRef)<0:
112
 
                            O.bodies[othSph].mat.jointNormal3=-jointNormalRef 
113
 
                    elif O.bodies[othSph].mat.joint==1 and ((O.bodies[othSph].mat.jointNormal1.cross(jointNormalRef)).norm()>0.05) :
114
 
                        O.bodies[othSph].mat.joint=2
115
 
                        if facetSphereDir.dot(jointNormalRef)>=0:
116
 
                            O.bodies[othSph].mat.jointNormal2=jointNormalRef
117
 
                        elif facetSphereDir.dot(jointNormalRef)<0:
118
 
                            O.bodies[othSph].mat.jointNormal2=-jointNormalRef
119
 
                elif  O.bodies[othSph].mat.onJoint==False :
120
 
                    O.bodies[othSph].mat.onJoint=True
121
 
                    O.bodies[othSph].mat.joint=1
 
107
                    elif O.bodies[othSph].state.joint==2 and ((O.bodies[othSph].state.jointNormal1.cross(jointNormalRef)).norm()>0.05) and ((O.bodies[othSph].state.jointNormal2.cross(jointNormalRef)).norm()>0.05):
 
108
                        O.bodies[othSph].state.joint=3
 
109
                        if facetSphereDir.dot(jointNormalRef)>=0:
 
110
                            O.bodies[othSph].state.jointNormal3=jointNormalRef
 
111
                        elif facetSphereDir.dot(jointNormalRef)<0:
 
112
                            O.bodies[othSph].state.jointNormal3=-jointNormalRef 
 
113
                    elif O.bodies[othSph].state.joint==1 and ((O.bodies[othSph].state.jointNormal1.cross(jointNormalRef)).norm()>0.05) :
 
114
                        O.bodies[othSph].state.joint=2
 
115
                        if facetSphereDir.dot(jointNormalRef)>=0:
 
116
                            O.bodies[othSph].state.jointNormal2=jointNormalRef
 
117
                        elif facetSphereDir.dot(jointNormalRef)<0:
 
118
                            O.bodies[othSph].state.jointNormal2=-jointNormalRef
 
119
                elif  O.bodies[othSph].state.onJoint==False :
 
120
                    O.bodies[othSph].state.onJoint=True
 
121
                    O.bodies[othSph].state.joint=1
122
122
                    O.bodies[othSph].shape.color=jointcolor4
123
123
                    if facetSphereDir.dot(jointNormalRef)>=0:
124
 
                        O.bodies[othSph].mat.jointNormal1=jointNormalRef
 
124
                        O.bodies[othSph].state.jointNormal1=jointNormalRef
125
125
                    elif facetSphereDir.dot(jointNormalRef)<0:
126
 
                        O.bodies[othSph].mat.jointNormal1=-jointNormalRef
 
126
                        O.bodies[othSph].state.jointNormal1=-jointNormalRef
127
127
 
128
128
#### for visualization:
129
129
#bj=0
130
130
#vert=(0.,1.,0.)
131
131
#hor=(0.,1.,1.)
132
132
#for o in O.bodies:
133
 
    #if o.mat.onJoint==True : # or o.shape.name=='Facet':
 
133
    #if o.state.onJoint==True : # or o.shape.name=='Facet':
134
134
        ##if  o.shape.name=='Facet':
135
135
            ##o.shape.wire=True
136
136
        ##o.state.pos+=(0,50,0)
137
137
        ##bj+=1
138
 
        #if o.mat.jointNormal1.dot(hor)>0 :
 
138
        #if o.state.jointNormal1.dot(hor)>0 :
139
139
            ##o.state.pos+=(0,50,0)
140
140
            #o.shape.color=jointcolor1
141
 
        #elif o.mat.jointNormal1.dot(hor)<0 :
 
141
        #elif o.state.jointNormal1.dot(hor)<0 :
142
142
            ##o.state.pos+=(0,55,0)
143
143
            #o.shape.color=jointcolor2
144
144
        #if o.mat.type>2 :
145
145
            #bj+=1
146
146
            #o.shape.color=jointcolor5
147
 
            ##print o.mat.jointNormal.dot(vert)
 
147
            ##print o.state.jointNormal.dot(vert)
148
148
 
149
149
 
150
150
##### to delete facets (should be OK to start the simulation after that!
154
154
 
155
155
O.resetTime()
156
156
O.interactions.clear()
157
 
print ''
158
 
print 'IdentificationSpheresOnJoint executed ! Spheres onJoint (and so on...) detected, facets deleted'
159
 
print ''
 
157
print '\n IdentificationSpheresOnJoint executed ! Spheres onJoint (and so on...) detected, facets deleted\n\n'
 
158
 
160
159
 
161
160