FISH Contact Model

Introduction

Note

The project file for this example may be viewed/run in PFC.[1] The data files used are shown at the end of this example.

This example demonstrates basic usage of the FISH contact model. A model of 50 balls is settled under gravity in a box. The FISH contact model is used, where a FISH operator is defined to resolve the force-displacement law at each contact in a multi-threaded environment.

Numerical Model

The complete data file (3D) for this example is “FISHModel.dat (3D),” shown in its entirety at the end. Excerpted lines of this file are discussed below.

In this example, a simple linear model formulation is implemented, comprising normal and shear stiffness and a Coulomb slip criterion.

First, global variables are defined to formulate the model in terms of contact micro-deformability, mimicking the behavior of the Linear contact model deformability method.

; Define contact parameters (deformability and friction)
[emod = 1e9]
[krat = 2.5]
[fric = 0.2]

Contact property indices are also stored as global variables, and will be used by the FISH operator. Preferring indices to strings in multi-threaded computations results in substantial performance improvements.

; Store contact property indices to prevent using strings during cycling,
; for better performances
[forceInd = contact.model.prop.index("fish","force")]
[stifftInd = contact.model.prop.index("fish","stifft")]
[stiffaInd = contact.model.prop.index("fish","stiffa")]
[rgapInd = contact.model.prop.index("fish","rgap")]
[symbolInd = contact.model.prop.index("fish","symbol")]

The behavior of the desired model is then defined in a FISH operator. In this example, a simple linear force-displacement law with Coulomb slip is implemented, without viscous damping. The contact force is resolved in the contact local coordinate system, where the normal force corresponds to the \(x\)-component of the contact force, and the shear force corresponds to the \(y\)- and \(z\)-components of the contact force.

Note that the FISH operator takes a finite set of arguments that correspond to the contact state information, as detailed in the FISH Model description. Additional contact properties—such as the normal and shear stiffness in this case—can be stored in contact extra variables. In this example, the contact friction coefficient is assumed to be uniform for all contacts and is stored as global FISH variable. One could define different values for this coefficient as well, and also store it as a contact extra variable.

; Create a FISH operator to be used to resolve the Force-Displacement
; law for each contact during cycling.
; Specific properties (stiffness) are stored as contact extra variables
fish operator linearModel(cp,trans,ang,curv1,curv2,inertialMass ...
                          ,gap,canFail,activated) 
    ;do the force displacement law
    local kn = contact.extra(cp,1)
    local ks = contact.extra(cp,2)
    local overlap = contact.prop(cp,rgapInd) - contact.gap(cp)
    local lin_F_old = contact.prop(cp,forceInd)
    local force = lin_F_old
    force->x = overlap * kn
    force->x = math.max(force->x,0.0)
    local sforce = vector(0,0,0)
    sforce->y = force->y - trans->y * ks
    sforce->z = force->z - trans->z * ks
    if canFail == true
        ;check for sliding
        local crit = force->x * fric;
        local sfmag = math.mag(sforce)
        if (sfmag > crit) 
            local rat = crit / sfmag;
            sforce *= rat;
        endif
    endif
    force->y = sforce->y
    force->z = sforce->z
    ;set the force
    contact.prop(cp,forceInd) = force
end

Next, a FISH function is defined to be registered as a callback and executed anytime a contact is created in PFC. This will allow for the desired FISH operator to be set at the contacts, and for contact-specific properties to be computed and stored. The effective translational and rotational stiffnesses of the contact are also computed, which will be used for timestep computation as for any contact model.

; Use a FISH function to be registered as a callback and executed
; anytime a contact is created in PFC, to set the contact properties.
fish define catchContact(c)
    ; Compute kn and ks based on the specified emod and krat, 
    ; just like the deformability method for the linear model
    local rsum = 0.0
    local rsq = 0.0
    local rad1 = 0.0
    if type.pointer.id(c) == contact.typeid('ball-ball')
        rad1 = ball.rad(contact.end1(c))
        local rad2 = ball.rad(contact.end2(c))
        rsum = rad1 + rad2
        rsq = 1./math.min(rad1,rad2)        
    endif
    if type.pointer.id(c) == contact.typeid('ball-facet')
        rad1 = ball.rad(contact.end1(c))
        rsum = rad1
        rsq = 1.0/rad1
    endif
    local kn = math.pi() * emod / (rsq * rsq * rsum)
    local ks = 0.0
    if krat > 0
        ks = kn / krat
    endif
    contact.prop(c,stifftInd) = vector(kn,ks)
    contact.prop(c,stiffaInd) = vector(0,0,0)
    contact.prop(c,symbolInd) = "linearModel"
    contact.extra(c,1) = kn
    contact.extra(c,2) = ks
end

The initial model consists of a cloud of balls in a box:

; Setup the model
model domain extent -10.0 10.0 
wall generate box -5.0 5.0
model gravity 0 0 -10
ball generate number 50 radius 0.8 1.2 box -5.0 5.0
ball attribute density 2500 damp 0.7

The contact cmat default is used to select the FISH contact model as the default model to be installed whenever a new contact is created, and the callback function defined above is also registered to be executed as well.

; Set the CMAT to use the FISH contact model for all contacts 
contact cmat default model fish
; Register the FISH function with the contact_Create callback event
fish callback add catchContact event contact_create

Finally, the model is solved to the desired target criterion.

; Now solve the model
model history name 'aratio' mechanical ratio-average
model solve
model save 'final'

The final state of the model is shown in the figure below.

../../../../../../../_images/p3d-cmfish-final.png

Figure 1: The system after cycling.

Discussion

This tutorial introduces the main steps required to use the FISH contact model. This contact model can be used to facilitate contact model development and testing. For time sensitive models, implementation of contact models using a C++ Plugin is encouraged.

Data Files

FISHModel.dat (3D)

; fname: fishModel.dat (3D)
;
;  Exercise the FISH contact model
;=============================================================================
model new
model large-strain on
model random 10001

; Define contact parameters (deformability and friction)
[emod = 1e9]
[krat = 2.5]
[fric = 0.2]

; Store contact property indices to prevent using strings during cycling,
; for better performances
[forceInd = contact.model.prop.index("fish","force")]
[stifftInd = contact.model.prop.index("fish","stifft")]
[stiffaInd = contact.model.prop.index("fish","stiffa")]
[rgapInd = contact.model.prop.index("fish","rgap")]
[symbolInd = contact.model.prop.index("fish","symbol")]

; Create a FISH operator to be used to resolve the Force-Displacement
; law for each contact during cycling.
; Specific properties (stiffness) are stored as contact extra variables
fish operator linearModel(cp,trans,ang,curv1,curv2,inertialMass ...
                          ,gap,canFail,activated) 
    ;do the force displacement law
    local kn = contact.extra(cp,1)
    local ks = contact.extra(cp,2)
    local overlap = contact.prop(cp,rgapInd) - contact.gap(cp)
    local lin_F_old = contact.prop(cp,forceInd)
    local force = lin_F_old
    force->x = overlap * kn
    force->x = math.max(force->x,0.0)
    local sforce = vector(0,0,0)
    sforce->y = force->y - trans->y * ks
    sforce->z = force->z - trans->z * ks
    if canFail == true
        ;check for sliding
        local crit = force->x * fric;
        local sfmag = math.mag(sforce)
        if (sfmag > crit) 
            local rat = crit / sfmag;
            sforce *= rat;
        endif
    endif
    force->y = sforce->y
    force->z = sforce->z
    ;set the force
    contact.prop(cp,forceInd) = force
end

; Use a FISH function to be registered as a callback and executed
; anytime a contact is created in PFC, to set the contact properties.
fish define catchContact(c)
    ; Compute kn and ks based on the specified emod and krat, 
    ; just like the deformability method for the linear model
    local rsum = 0.0
    local rsq = 0.0
    local rad1 = 0.0
    if type.pointer.id(c) == contact.typeid('ball-ball')
        rad1 = ball.rad(contact.end1(c))
        local rad2 = ball.rad(contact.end2(c))
        rsum = rad1 + rad2
        rsq = 1./math.min(rad1,rad2)        
    endif
    if type.pointer.id(c) == contact.typeid('ball-facet')
        rad1 = ball.rad(contact.end1(c))
        rsum = rad1
        rsq = 1.0/rad1
    endif
    local kn = math.pi() * emod / (rsq * rsq * rsum)
    local ks = 0.0
    if krat > 0
        ks = kn / krat
    endif
    contact.prop(c,stifftInd) = vector(kn,ks)
    contact.prop(c,stiffaInd) = vector(0,0,0)
    contact.prop(c,symbolInd) = "linearModel"
    contact.extra(c,1) = kn
    contact.extra(c,2) = ks
end

; Setup the model
model domain extent -10.0 10.0 
wall generate box -5.0 5.0
model gravity 0 0 -10
ball generate number 50 radius 0.8 1.2 box -5.0 5.0
ball attribute density 2500 damp 0.7

; Set the CMAT to use the FISH contact model for all contacts 
contact cmat default model fish
; Register the FISH function with the contact_Create callback event
fish callback add catchContact event contact_create

; Now solve the model
model history name 'aratio' mechanical ratio-average
model solve
model save 'final'
program return

Endnote