Using Python with PFC
This section gives an example of interacting with PFC through a Python program.
To begin, go to the IPython pane by selecting it from the i menu. All the lines in this example may be entered at the prompt of this pane. Also, see Example File below for further information about running this example.
First, import the c itasca module with the Python
import
statement. The as
in the import
statement shorthands the
name of the itasca
module to it
. The itasca
module defines the
interaction between Python and PFC.
import itasca as it
Next some PFC particles are created. The it.command
function
(itasca.command
) is used to issue a set of PFC commands. A cubic
packing of 8,000 spherical particles is created. Triple quotes are used in
Python to define a multi-line string.
it.command("""
model new
model large-strain on
model domain extent -5e-2 6e-2 -6e-2 5e-2 -5e-2 5e-2
contact cmat default model linear property kn 1e1 dp_nratio 0.2
ball generate cubic box -0.02375 0.02375 rad 1.25e-3
ball attr dens 2600
""")
Confirm that 8000 balls were created using the it.ball.count
(itasca.ball.count
) function.
it.ball.count()
- output:
8000
The c itasca module defines functions and classes for interacting with the PFC model. Some of these functions return objects.
ball = it.ball.find(1)
print(ball)
- output:
<itasca.ball.Ball object at 0x000001686BD452B0, ID : 1>
Here it.ball.find(1)
returns a c Ball object
for the ball with id 1. Most PFC model items are objects (instances of a
class) like this. The c Ball object has many
methods.
ball.radius()
- output:
0.00125
In the above, the variable ball
is the c Ball
object created in the preceding snippet. radius
is a ball method that
returns the ball radius.
The for
statement is used for iteration.
radius_sum = 0.0
for b in it.ball.list():
radius_sum += b.radius()
In the above, the for
statement is used to loop over all the PFC
balls. During each loop the variable b
is a different
c Ball object. The result of the iteration can be
checked to see that the sum of the ball radii is as expected.
print(radius_sum)
print(ball.radius() * it.ball.count())
- output:
10.000000000000503 10.0
There are two places where one may interact with PFC balls:
the functions that are part of the c itasca.ball module, and
the methods of an individual c Ball object.
In Python terminology, ball
is a submodule of the itasca
module (hence
its name c itasca.ball). In the snippets above,
it.ball.count
, it.ball.list
, and it.ball.find
are functions
defined in the it.ball
module (that is, the c itasca.ball
module: remember, it
was defined at the outset as a shorthand substitute
for itasca
). On the other hand, b.radius
is a method defined as part
of the c Ball class.
See the c itasca.ball module reference for the full list of module functions.
See the c Ball class description for a complete list of the
Ball
object methods.
PFC contacts are linked to Python in the same way as with ball. This is seen by cycling PFC for one step so contacts are created.
it.command("model cycle 1")
Use near
to find the ball nearest the origin,
b = it.ball.near((0,0,0))
and confirm its position with the pos
method.
b.pos()
- output:
vec3(( 1.250000e-03, 1.250000e-03, 1.250000e-03 ))
Note
The pos
method returns a length three vector represented as a vec3
object. See the c vec module for more information about
vec3
objects and working with vectors in Python.
The number of this ball’s contacts is obtained with the contacts
method.
len(b.contacts())
- output:
6
The built-in function len
is used in Python to determine the
length of a sequence. The following code loops over all the contacts
that the ball b
has and prints the contact id and contact position.
for c in b.contacts():
print("contact with id: {} at {}".format(c.id(), c.pos()))
- output:
contact with id: 11063 at vec3(( 1.250000e-03, 1.250000e-03, -1.951564e-18 )) contact with id: 12163 at vec3(( 1.250000e-03, -1.951564e-18, 1.250000e-03 )) contact with id: 12218 at vec3(( -1.951564e-18, 1.250000e-03, 1.250000e-03 )) contact with id: 12221 at vec3(( 2.500000e-03, 1.250000e-03, 1.250000e-03 )) contact with id: 12222 at vec3(( 1.250000e-03, 2.500000e-03, 1.250000e-03 )) contact with id: 12223 at vec3(( 1.250000e-03, 1.250000e-03, 2.500000e-03 ))
Note
The string object format
method is used to substitute values
into a string. The {} parts of the string are replaced by the
arguments to the format method — in this case to show the contact id and
position.
Similar inquiries may be made of a contact object.
c = b.contacts()[0]
The brackets ([ ]) are used to get an item from the tuple returned by the
contacts
method. (Note in Python indices start a 0, so this method call
returns the first contact.) Contact force (in the global frame) is obtained
with force_global
.
c.force_global()
- output:
vec3(( 0.000000e+00, 0.000000e+00, 0.000000e+00 ))
The properties defined on this contact are obtained with props
.
c.props()
- output:
{'kn': 10.0, 'ks': 0.0, 'fric': 0.0, 'lin_force': vec3(( 0.000000e+00, 0.000000e+00, 0.000000e+00 )), 'lin_slip': False, 'lin_mode': 0, 'rgap': 0.0, 'emod': 5092.958178940651, 'kratio': 0.0, 'dp_nratio': 0.2, 'dp_sratio': 0.0, 'dp_mode': 0, 'dp_force': vec3(( 0.000000e+00, -0.000000e+00, -0.000000e+00 )), 'user_area': 0.0}
The props
method returns a Python dictionary object.
To access an individual part of a dictionary, brackets ([ ]) are used to
specify the desired key name.
c.props()['fric']
- output:
0.0
Alternatively, the prop
method can be used.
c.prop('fric')
- output:
0.0
Contact (or ball) properties can be set with the set_prop
method.
c.set_prop('fric', 0.5)
print(c.prop('fric'))
- output:
0.5
The bodies on either side of the contact can be accessed with the
end1
and end2
methods.
print(c.end1())
print(c.end2())
- output:
<itasca.ball.Ball object at 0x000001686BD45430, ID : 3811> <itasca.ball.Ball object at 0x000001686BD45430, ID : 4211>
Test which end of the contact corresponds to the Ball object b
.
print(c.end1() == b)
print(c.end2() == b)
- output:
False True
In the following example a Python list comprehension is used to create
a list of all the balls contacting the ball represented by the object
b
.
neighbor_list = [c.end1() if c.end2() == b else c.end2() for c in b.contacts()]
print("central ball id: {}, position: {}".format(b.id(), b.pos()))
print
for i, neighbor in enumerate(neighbor_list):
print("neighbor ball {} id: {}, position: {}".format(i,
neighbor.id(),
neighbor.pos()))
- output:
central ball id: 4211, position: vec3(( 1.250000e-03, 1.250000e-03, 1.250000e-03 )) neighbor ball 0 id: 3811, position: vec3(( 1.250000e-03, 1.250000e-03, -1.250000e-03 )) neighbor ball 1 id: 4191, position: vec3(( 1.250000e-03, -1.250000e-03, 1.250000e-03 )) neighbor ball 2 id: 4210, position: vec3(( -1.250000e-03, 1.250000e-03, 1.250000e-03 )) neighbor ball 3 id: 4212, position: vec3(( 3.750000e-03, 1.250000e-03, 1.250000e-03 )) neighbor ball 4 id: 4231, position: vec3(( 1.250000e-03, 3.750000e-03, 1.250000e-03 )) neighbor ball 5 id: 4611, position: vec3(( 1.250000e-03, 1.250000e-03, 3.750000e-03 ))
As expected, this ball has 6 neighbors because of the cubic packing.
Python Type System
In the next example, in a new model three balls and a wall are created to demonstrate the Python type system.
it.command("""
model new
model large-strain on
model domain extent -1 1 -1 1 -1 1
contact cmat default model linear property kn 1e1 dp_nratio 0.2
""")
from vec import vec
origin = vec((0.0, 0.0, 0.0))
rad = 0.1
eps = 0.001
b1 = it.ball.create(rad, origin)
b2 = it.ball.create(rad, origin + (rad-eps, 0, 0))
# create a third ball close enough to generate a virtual contact
b3 = it.ball.create(rad, origin + (rad*3 + eps, 0, 0))
it.command("""
ball attribute density 1200
wall create vertices ...
-{rad} -{rad} -{rad} ...
{rad} -{rad} -{rad} ...
{rad} {rad} -{rad}
model cycle 1
""".format(rad=rad))
Executing this code and plotting the balls and walls should look like the figure below.
There are three active contacts and one inactive contact in this model. By default the it.contact.list
function only returns active contacts.
for c in it.contact.list():
print(c)
- output:
<itasca.BallBallContact object at 0x000001687388E030, ID : 1> <itasca.BallFacetContact object at 0x000001686BD45210, ID : 1> <itasca.BallFacetContact object at 0x000001687388E030, ID : 2>
The inactive contacts can be included in the list by adding the
all=True
keyword argument.
for c in it.contact.list(all=True):
print(c)
- output:
<itasca.BallBallContact object at 0x000001686BD45210, ID : 1> <itasca.BallBallContact object at 0x000001687388E030, ID : 2> <itasca.BallFacetContact object at 0x000001686BD45210, ID : 1> <itasca.BallFacetContact object at 0x000001687388E030, ID : 2>
Different contact types have different Python types. Use with the
built-in type
function to determine the type of an object.
c1, c2, c3, c4 = tuple(it.contact.list(all=True))
print(type(c1) is it.BallBallContact)
print(type(c1) is it.BallFacetContact)
print(type(c3) is it.BallFacetContact)
- output:
True False True
Here it.BallBallContact
is a Python type object. If only the Ball-Ball
contacts are wanted, pass this type object to the it.contact.list
function.
for c in it.contact.list(type=it.BallBallContact, all=True):
print(c)
- output:
<itasca.BallBallContact object at 0x000001686BD44FF0, ID : 1> <itasca.BallBallContact object at 0x000001687388E030, ID : 2>
Python Callback Functions
Python functions can be run at specific points during the PFC
calculation sequence. The itasca.set_callback
function registers
Python functions to be called during cycling while the
itasca.remove_callback
removes these function callbacks.
i=0
def my_callback(*args):
global i
i += 1
print("in Python callback function.")
it.set_callback("my_callback", -1)
it.command("model cycle 5")
print("The Python callback function was called {} times".format(i))
it.remove_callback("my_callback",-1)
i=0
it.command("model cycle 5")
print("The Python callback function was called {} times".format(i))
- output:
in Python callback function. in Python callback function. in Python callback function. in Python callback function. in Python callback function. The Python callback function was called 5 times The Python callback function was called 0 times
Note
Callback function registration is cleared when a model new
command
is given. Be careful not to register the same Python function to
be called more than once per cycle. Use the fish list callback
command to see which functions have already been registered.
Example File
The source code for the example above is python_pfc.py. This file can be loaded into the program’s i Editor pane and run using the “Execute” button.
Was this helpful? ... | Itasca Software © 2024, Itasca | Updated: Aug 13, 2024 |