Making a triaxial weave from three biaxial weaves

This notebook documents the path to implementation of triaxial weaves as the composition of three biaxial weaves generated using matrices. The functions in this notebook are early versions of those that are now implemented in triaxial-weave-units.R, biaxial-weave-units.R and render-weave-grids.R in the main ‘branch’.

Preamble

A triaxial weave can be composed from three independent biaxial weaves. This notebook shows how…

It uses functions already written in weave-grids.R. It should eventually also use functions from biaxial-weave-units.R to make up the component biaxial weave patterns.

library(dplyr)
library(sf)
library(pracma)
library(tmap)
library(stringr)

source("../render-weave-grids.R")
source("../biaxial-weave-units.R")

A triangular grid

Central to the approach is a triangular grid. grid_generator() returns a function, that given three integer coordinates returns x-y coordinates of the centroids of the corresponding triangular grid cell.

trigrid <- grid_generator(n_axes = 3)

The thing to understand about this grid coordinate system is that every triple of coordinates must sum to one of two consecutive integer values. By ‘default’ this would be 0 or 1, but any consecutive pair of integers will work (this flexibility is important later, when we want to match input biaxial coordinates to an output triaxial space).

We’re calling the lower of the two integer values the parity of the grid.

Note that the grid_generator function does no checking with regard to parity. If coordinates with different parity are supplied they are valid, but will return x-y coordinates coincident with the results of other coordinate inputs. This is because the unit vectors of the grid are non-orthogonal and the sum of any two is the inverse of the third.

Conversion of biaxial coordinates to triaxial

This constraint on coordinate sums means that given two coordinates and the desired parity, we can determine the missing coordinate (provided we know which one is missing). That is the basis for the following function which is used as a helper when populating an output combined triaxial grid with values from input biaxial grids.

Axis here is taken to mean: 1: layer A and B; 2: layers B and C; 3: layers C and A. The ordering of the biaxial combinations is important (cyclic).

# function to supply missing third coordinate in a triangular grid given two
# other coordinates, the parity, and the axis. The axis tells us which is missing:
# axis = 1 --> 3rd, axis = 2 --> 1st, axis = 3 --> 2nd
transform_ab_to_abc <- function(z1, z2, parity = 0, axis = 1) {
  # comments on first case, others are cyclic shifts - most easily read down the code
  # when axis = 1 a_coord and b_coord are preserved, c_coord is the paired residual values
  a_coord <- switch(axis, z1, c(parity + 1 - z1 - z2, parity - z1 - z2), z2)
  b_coord <- switch(axis, z2, z1, c(parity + 1 - z1 - z2, parity - z1 - z2))
  c_coord <- switch(axis, c(parity + 1 - z1 - z2, parity - z1 - z2), z2, z1)
  result <- switch( 
    axis, # return two triples that differ in the missing coordinate
    matrix(c(a_coord, b_coord, c_coord[1], 
             a_coord, b_coord, c_coord[2]), nrow = 2, ncol = 3, byrow = TRUE),
    matrix(c(a_coord[1], b_coord, c_coord, 
             a_coord[2], b_coord, c_coord), nrow = 2, ncol = 3, byrow = TRUE),
    matrix(c(a_coord, b_coord[1], c_coord, 
             a_coord, b_coord[2], c_coord), nrow = 2, ncol = 3, byrow = TRUE))
  return(result)
}

Helper functions to assemble the triaxial grid

Adding biaxial matrices to a triaxial grid at the right locations

First a function to take data from a biaxial pattern matrix and add them to a triaxial grid at the right location.

# function to add values to a supplied list tri, from a supplied matrix M, where
# values in the list will be indexed by triangular coordinates (converted to a
# string). This will be called 3 times, once for each axis. The calling context 
# should determine the parity, based on the sizes of the three matrices
add_biaxial_to_triaxial <- function(tri, M, axis = 1, parity = 0) {
  for (col in 1:ncol(M)) {
    for (row in 1:nrow(M)) {
      # get the target grid coordinates (there will be two)
      abc <- transform_ab_to_abc(col, row, parity = parity, axis = axis)
      for (par in 1:2) {
        # make a key (string) from each pair
        key <- abc[par, ] %>% str_c(collapse = ",")
        # if it is not in the target list, add a new thread order 
        if (is.null(tri[[key]])) {
          tri[[key]] <- list(decode_biaxial_to_order(M[row, col], axis = axis)) # c(M[row, col])
        } else { # if there is something there, then append the value
          tri[[key]] <- append(tri[[key]], list(decode_biaxial_to_order(M[row, col], axis = axis)))
          # tri[[key]] <- c(tri[[key]], M[row, col] )
        }
      }
    }
  }
  return(tri)
}

Encoding the state at each location in a biaxial weave

This function encodes the state at each site in a raw biaxial weave pattern, based on the pattern matrix and the warp and weft thread matrices. These are decoded by decode_biaxial_to_order() to axis-appropriate sequences of thread numbers.

encode_biaxial_weave <- function(pattern, warp, weft) {
  pattern[which(pattern == 1)] <- 5        # warp present and on top
  pattern[which(pattern == 0)] <- 4        # weft present and on top
  pattern[which(warp < 0)] <- 1            # warp absent
  pattern[which(weft < 0)] <- 2            # weft absent
  pattern[which(weft < 0 & warp < 0)] <- 3 # both absent
  return(pattern)
}

Converting biaxial encoding to a strand order

The data it adds are strand orderings at each location, where a pair of values indicate that the first strand is on top of the second, a single value indicates only that strand is present in the biaxial weave at that location, and a NULL means no strand is present from the biaxial at that location. Which strands are which depends on the axis.

decode_biaxial_to_order <- function(code, axis = 1) {
  if (axis == 1) {
    return(switch(
      code, 2, 1, NULL, 2:1, 1:2
    ))
  }
  if (axis == 2) {
    return(switch(
      code, 3, 2, NULL, 3:2, 2:3
    ))
  }
  if (axis == 3) {
    return(switch(
      code, 1, 3, NULL, c(1,3), c(3,1)
    ))
  }
}

Combining orderings from biaxial weaves

This function takes a list of orderings on the values 1:n, and combines them into a single consistent ordering. If a value is not present in any of the members of the orderings list it will not appear in the output. If all the orderings are empty, NULL is returned.

# combines a set of orderings on the values
# the orderings are a list of vectors (which may be empty)
# for example list(c(1, 2), c(2, 3), c(1, 3))
combine_orderings <- function(orderings, values = 1:3, verbose = FALSE) {
  # assemble a matrix of the positions of entries
  # in each ordering among the values
  ranks <- matrix(0, length(orderings), length(values))
  for (i in seq_along(orderings)) {
    ranks[i, ] <- match(values, orderings[[i]])
  }
  # replace any missing matches with a high score
  ranks[which(is.na(ranks))] <- 100
  # sum the ranks of each value
  scores <- colSums(ranks)
  max_score <- length(orderings) * 100
  number_present <- sum(scores < max_score)
  if (number_present == 0) {
    result <- NULL
  } else {
    result <- values[order(scores)[1:number_present]]
  }
  if (verbose) {
    return(list(result = result, ranks = ranks, scores = scores)) 
  } else {
    return(result)
  }
}

Doing the weave…

k <- 2
type = "cube"

Make up the paired biaxial weaves

Setup three tie-up matrices… and for now just use I for the treadling and threading matrices. Note that we have a cyclic arrangement of biaxial pairs A > B > C > A, so we treat the first of these as the ‘warp’ in each pairing.

# three tie-up matrices - no idea if these are right
# cube weave
tu_AB <- make_twill_matrix(c(1,2)) %>% repmat(k)
tu_BC <- make_twill_matrix(c(1,2)) %>% repmat(k)
tu_CA <- make_twill_matrix(c(1,2)) %>% repmat(k)

# open hex basic weave
if (type == "hex") {
  tu_AB <- ones(3) %>% repmat(k)
  tu_BC <- ones(3) %>% repmat(k)
  tu_CA <- ones(3) %>% repmat(k)
}
tr_AB <- diag(nrow(tu_AB))
th_AB <- diag(ncol(tu_AB))
tr_BC <- diag(nrow(tu_BC))
th_BC <- diag(ncol(tu_BC))
tr_CA <- diag(nrow(tu_CA))
th_CA <- diag(ncol(tu_CA))

Make up strand matrices… later these can include more than one colour… for now we’re just getting it to work. (See modify_for_missing_threads() in biaxial_weave_units.R)

# thread patterns, these would be coded x > 0 for a thread, -1 for missing
# cube weave
AB_A <- ones(3) %>% repmat(k)
AB_B <- t(AB_A)
BC_B <- ones(3) %>% repmat(k) 
BC_C <- t(BC_B)
CA_C <- ones(3) %>% repmat(k) 
CA_A <- t(CA_C)

#open hex basic weave
if (type == "hex") {
  AB_A <- matrix(c(-1, -1, 1), 3, 3, byrow = TRUE) %>% repmat(k)
  AB_B <- t(AB_A)
  BC_B <- matrix(c(-1, -1, 1), 3, 3, byrow = TRUE) %>% repmat(k)
  BC_C <- t(BC_B)
  CA_C <- matrix(c(-1, -1, 1), 3, 3, byrow = TRUE) %>% repmat(k)
  CA_A <- t(CA_C)
}
# for example
AB_A

##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    1    1    1    1    1    1
## [2,]    1    1    1    1    1    1
## [3,]    1    1    1    1    1    1
## [4,]    1    1    1    1    1    1
## [5,]    1    1    1    1    1    1
## [6,]    1    1    1    1    1    1

Do the matrix multiplication and encode the outcome as follows depending on the state of the warp and weft thread matrices (here AB_A, AB_B etc.). The encode_biaxial_weave() function is in biaxial-weave-units.R.

Pattern	Warp	Weft	Overall state	Code
1	>0	>0	Warp on top, both present	5
0	>0	>0	Weft on top, both present	4
*	<0	<0	Both absent	3
*	>0	<0	Weft absent, warp present	2
*	>0	<0	Warp absent, weft present	1

# three weave-patterns 1 for the warp on top 0 for the weft
pat_AB <- (((tr_AB %*% tu_AB %*% th_AB) > 0) * 1) %>%
  encode_biaxial_weave(AB_A, AB_B)
pat_BC <- (((tr_BC %*% tu_BC %*% th_BC) > 0) * 1) %>%
  encode_biaxial_weave(BC_B, BC_C)
pat_CA <- (((tr_CA %*% tu_CA %*% th_CA) > 0) * 1) %>%
  encode_biaxial_weave(CA_C, CA_A)

pat_AB

##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    5    4    4    5    4    4
## [2,]    4    5    4    4    5    4
## [3,]    4    4    5    4    4    5
## [4,]    5    4    4    5    4    4
## [5,]    4    5    4    4    5    4
## [6,]    4    4    5    4    4    5

Combine biaxial weaves into a triaxial ‘loom’

Now combine the results into a single loom. This could be a 4D array (3 x 3D arrays), but it would be very sparse, so it’s more convenient to use a list keyed by coordinate triples as strings. This might be slower than a 4D array, but the code is likely clearer this way, and we are not dealing with very large grids (a few hundred sites at most).

The parity calculation here makes a triangular grid that accomodates the weft strands (rows) in all three matrices. We need to consider how to handle unequal sized matrices. At present two should work:

All three matrices same size will yield a hexagonal unit. Note that even-sized matrices are probably required, else the parity value will be odd and this means the orientation of cells is reversed
If one matrix is twice as large as the (equal-sized) other two, then we will get a diamond

Other combinations likely (untested) make skewed hexagons, and it’s not clear what to do with those (for now).

This operation should probably be wrapped in a function.

par <- (3 + ncol(pat_AB) + ncol(pat_BC) + ncol(pat_CA)) %/% 2
# par <- ifelse(par %% 2 == 0, par, par + 1)
loom <- list() %>%
  add_biaxial_to_triaxial(pat_AB, axis = 1, parity = par) %>% 
  add_biaxial_to_triaxial(pat_BC, axis = 2, parity = par) %>% 
  add_biaxial_to_triaxial(pat_CA, axis = 3, parity = par)

# only locations where we have complete information should be used
loom <- loom[which(lengths(loom) == 3)] 
loom <- lapply(loom, combine_orderings)

Drawing the weave…

This uses get_all_cell_strands() from render-weave-grids.R. The weave appearance here is produced simply by overdrawing in the appropriate order. That will need to be changed to draw only the visible portions of strands, but differencing the shapes in the correct sequence.

weave <- list()
colours <- c()
for (key in names(loom)) {
  # recover the coordinates from the key string, and get the x-y offset
  coords <- str_split(key, ",")[[1]] %>% as.numeric()
  xy <- trigrid(coords)
  # get the strand order
  strand_order <- loom[[key]]
  if (is.null(strand_order)) { next }
  # orientations from the coordinate parity and the strand axis
  o <- 0:2 * 120
  w <- 0.8
  weave <- add_shapes_to_list(weave, 
                              get_visible_cell_strands(n = 3, width = w,
                              parity = sum(coords), strand_order = strand_order,
                              orientations = o) + xy)
  colours <- c(colours, strand_order)
}
tris <- weave %>% st_sfc()
unit <- st_sf(data.frame(strand = as.factor(colours)), geometry = tris)
unit %>% plot(border = NA)

Superceded code

Determining the layering order based on input layer orders

Each biaxial is coded TRUE (==1) where the warp is on top. So we have three truth values. These are combined by the function below to produce a permutation order on the layers.

# determines the permutation order based on three T/F values for pairs of layers
# the truth table is. Values provided as integers, but listed below as boolean
#
# 1 > 2 | 2 > 3 | 3 > 1 | permutation
# ------|-------|-------|-------------
# TRUE  | TRUE  | FALSE | (1, 2, 3)
# TRUE  | FALSE | FALSE | (1, 3, 2)
# FALSE | TRUE  | FALSE | (2, 1, 3)
# FALSE | TRUE  | TRUE  | (2, 3, 1)
# TRUE  | FALSE | TRUE  | (3, 1, 2)
# FALSE | FALSE | TRUE  | (3, 2, 1)
#
# Note that TTT and FFF are invalid and will return an NA
pat_to_order <- function(x) {
  if (sum(x) == 3 || sum(x) == 0) {return(NA)}
  if (x[1]) {
    if (x[2]) { return(1:3) }           # 110 -> 123
    else {
      if (!x[3]) { return(c(1, 3, 2)) } # 100 -> 132
      else { return(c(3, 1, 2)) }       # 101 -> 312
    }
  }
  if (x[2]) {
    if (x[3]) { return(c(2, 3, 1)) }    # 011 -> 231
    else { return(c(2, 1, 3)) }         # 010 -> 213
  }
  return(3:1)                           # 001 -> 321
}