A2. Form Finding

Objective

We will learn how to form find a cable net with the force density method using the compas function fd_numpy. For doing so we must create the data structure required and compile it to the numeric data.

Procedure

• o. setup the paths and proxy

• a. setup the datastructure for the cable net mesh

  • a1. create the mesh

  • a2. update its default attributes

  • a3. set boundary conditions

• b. form find it into its equilibrium state

  • b1. compile the numeric data

  • b2. compute equilibrium

  • b3. update the datastructure

• c. visualise the form found cable net

• d. export the updated mesh datastructure

o. Paths and Proxy

We will build upon the mesh data structure that we imported and saved in the JSON file in the data folder in A1. Import section. We will output the resulting mesh in a new JSON file.

# ==============================================================================
# Paths
# ==============================================================================

HERE = os.path.dirname(__file__)
DATA = os.path.abspath(os.path.join(HERE, '..', 'data'))
FILE_I = os.path.join(DATA, 'cablenet_import.json')
FILE_O = os.path.join(DATA, 'cablenet_fofin.json')

In order to run the numerical compas function fd_numpy that bases on NumPy in Rhino we must create a Proxy.

# ==============================================================================
# Proxy
# ==============================================================================

proxy = Proxy('compas.numerical')
fd = proxy.fd_numpy

a. Datastructure

We want to create the COMPAS mesh data structure, then we want to extend it with attributes and methods related to form finding and analysis of flexible cable net formwork for concrete shells, and we want to set its boundary conditions.

imported mesh

a1. Create the Mesh

We start by creating a COMPAS mesh from the JSON file that we generated by the import from Rhino.

# ==============================================================================
# Cablenet mesh datastructure
# ==============================================================================

# create the mesh from imported geometry
mesh = Mesh.from_json(FILE_I)

We use a mesh instead of a network because it provides functionalities for the concrete weight and fabric shuttering etc.

a2. Update Default Vertex and Edge Attributes

The necessary default attributes for the vertices are (similar to the nodes in section A0):

• 'rx': 0.0, 'ry': 0.0, 'rz': 0.0 – component of an residual force

• 'px': 0.0, 'py': 0.0, 'pz': 0.0 – component of an externally applied load

• 'is_anchor': False, – to indicate that a vertex is anchored and can take reaction forces in XYZ.

The necessary attributes for the edges are:

• 'q': 1.0 – force densities of an edge.

• 'f': 0.0 – force in an edge.

• 'l': 0.0 – stressed length of an edge.

• 'l0': 0.0 –unstressed length of an edge.

The meaning of the unstressed length l0 will be explained later in section C1. materialise cables.

# a2. set default vertex and edge attributes
dva = {
    'rx': 0.0,            # X-component of an residual force.
    'ry': 0.0,            # Y-component of an residual force.
    'rz': 0.0,            # Z-component of an residual force.
    'px': 0.0,            # X-component of an externally applied load.
    'py': 0.0,            # Y-component of an externally applied load.
    'pz': 0.0,            # Z-component of an externally applied load.
    'is_anchor': False,   # Indicate that a vertex is anchored and can take reaction forces in XYZ.
}
mesh.update_default_vertex_attributes(dva)

dea = {
    'q': 1.0,             # Force densities of an edge.
    'f': 0.0,             # Force in an edge.
    'l': 0.0,             # Stressed Length of an edge.
    'l0': 0.0,            # Unstressed Length of an edge.
}
mesh.update_default_edge_attributes(dea)

all of this is already implemented within compas_fofin with the class compas_fofin.datastructures.Cablenet.

a3. Set Boundary Conditions

We identify the anchors on the boundary and update the 'is_anchor' attribute of these vertices:

# a3. Boundary conditions
boundary = mesh.vertices_on_boundaries()[0]
mesh.vertices_attribute('is_anchor', True, keys=boundary)

We visualise the cable net with anchored vertices in red:

# ==============================================================================
# Visualize
# ==============================================================================

baselayer = "CSD2::A2::FormFinding"

artist = MeshArtist(mesh, layer=baselayer+"::InputMesh")
artist.clear_layer()

artist.draw_vertices(color={vertex: (255, 0, 0) for vertex in mesh.vertices_where({'is_anchor': True})})
artist.draw_edges()
artist.draw_faces()

b. Form Finding

Now we want to form find the cable net into its equilibrium state under pure prestress. For this, we will use the numeric COMPAS function fd_numpy. The function does not 'understand' the datastructure, it requires the input data in lists. Thus, we provide it with lists and it form finds the equilibrium numerically and it outputs arrays, but it has no idea of the network itself, so we must update the mesh attributes ourselves afterwards.

mesh after form finding

to see what this function does exactly check out thefd_numpyfunction.

b1. Compile to Numeric Data

This requires the following input parameters and returns:

Parameters

• vertices (list) – XYZ coordinates of the vertices of the network

• edges (list) – Edges between vertices represented by

• fixed (list) – Indices of fixed vertices.

• q (list) – Force density of edges.

• loads (list) – XYZ components of the loads on the vertices.

Returns

• xyz (array) – XYZ coordinates of the equilibrium geometry.

• q (array) – Force densities in the edges.

• f (array) – Forces in the edges.

• l (array) – Lengths of the edges

• r (array) – Residual forces.

the numerical function inputs are lists that are ordered in the index and not the keys (vertex ids), so we must convert for the input the keys to indices and for the returns the indices to keys!

For this, we use the Mesh.key_index() method that returns a dictionary that maps vertex dictionary keys to the corresponding index in a vertex list or array.

# ==============================================================================
# Compile numerical data as fofin input
# ==============================================================================

# to map vertex dictionary keys to the corresponding index in a vertex list
vertex_index = mesh.key_index()

# complile lists (ordered in indices) of geometry, loads and force densities
X = mesh.vertices_attributes('xyz')
P = mesh.vertices_attributes(['px', 'py', 'pz'])
Q = mesh.edges_attribute('q')

# translate anchored vertex keys to index list of fixed vertices
fixed = [vertex_index[vertex] for vertex in mesh.vertices_where({'is_anchor': True})]
# translate edges from vertex keys to indexes
edges = [(vertex_index[u], vertex_index[v]) for u, v in mesh.edges()]

b2. Compute Equilibrium

Now we hand this over to the fd function that we made available in Rhino with the Proxy (see o.Paths and Proxy):

# ==============================================================================
# Compute equilibrium
# ==============================================================================

X, Q, F, L, R = fd(X, edges, fixed, Q, P)

If you were to run it like this with the visualisation you would see that the geometry actually hasn't changed at all yet.

b3. Update the Cablenet Mesh Datastructure

After the numerical calculation, we must update the mesh datastructure with its new equilibrium xyz-coordinates and the remaining residual forces by setting its vertex attributes and the new force density, force and length by setting its edge attributes.

# ==============================================================================
# Update cablenet mesh
# ==============================================================================

for vertex in mesh.vertices():
    # translate back from indices to keys
    index = vertex_index[vertex]
    mesh.vertex_attributes(vertex, 'xyz', X[index])
    mesh.vertex_attributes(vertex, ['rx', 'ry', 'rz'], R[index])

for index, edge in enumerate(mesh.edges()):
    # translate back from indices to keys
    mesh.edge_attribute(edge, 'q', Q[index][0])
    mesh.edge_attribute(edge, 'f', F[index][0])
    mesh.edge_attribute(edge, 'l', L[index][0])

all of this is already implemented within compas_fofin with the compas_fofin.fofin.update_xyz_numpy method that already takes care of the update of the attributes as it has access to them.

The simple visualisation remains the same as in a. Datastructures.

c. Visualization

If we want to visualise the reaction forces, residuals and edge forces to better understand what's going on, we can do so as follows:

mesh with visualized reaction forces (green) and edge forces (red)

from compas_rhino.artists import MeshArtist
from A3_visualization import draw_reactions
from A3_visualization import draw_residuals
from A3_visualization import draw_forces
from A3_visualization import draw_loads

...

# ==============================================================================
# Visualize
# ==============================================================================

baselayer = "CSD2::A2::FormFinding"

artist = MeshArtist(mesh, layer=baselayer+"::Mesh")
artist.clear_layer()

artist.draw_vertices(color={vertex: (255, 0, 0) for vertex in mesh.vertices_where({'is_anchor': True})})
artist.draw_edges()
artist.draw_faces()

draw_reactions(mesh, baselayer=baselayer)
draw_residuals(mesh, baselayer=baselayer)
draw_forces(mesh, baselayer=baselayer)
draw_loads(mesh, baselayer=baselayer)

These functions are not just available from somewhere in space, we import them from a file that we implement in the next section A3. Visualisation.

d. Export

At the very end, we export the equilibrium cable net by serialising the mesh data structure into a JSON file.

# ==============================================================================
# Export
# ==============================================================================

mesh.to_json(FILE_O)