# Enumeration of Various Animals on the triangular lattice

Toufik Mansour  Department of Mathematics, University of Haifa, 3498838 Haifa,
Israel<EMAIL_ADDRESS>and  Reza Rastegar  Occidental Petroleum
Corporation, Houston, TX 77046 and Departments of Mathematics and Engineering,
University of Tulsa, OK 74104, USA - Adjunct Professor<EMAIL_ADDRESS>
###### Abstract.

In this paper, we consider various classes of polyiamonds that are animals
residing on the triangular lattice. By careful analyses through certain layer-
by-layer decompositions and cell pruning/growing arguments, we derive explicit
forms for the generating functions of the number of nonempty translation-
invariant baryiamonds (bargraphs in the triangular lattice), column-convex
polyiamonds, and convex polyiamonds with respect to their perimeter. In
particular, we show that the number of
(A) baryiamonds of perimeter $n$ is asymptotically

$\frac{(\xi+1)^{2}\sqrt{\xi^{4}+\xi^{3}-2\xi+1}}{2\sqrt{\pi
n^{3}}}\xi^{-n-2},$

where $\xi$ is a root of a certain explicit polynomial of degree 5.
(B) column-convex polyiamonds of perimeter $n$ is asymptotic to

$\frac{(17997809\sqrt{17}+3^{3}\cdot 13\cdot
175463)\sqrt{95\sqrt{17}-119}}{2^{7}\cdot 43^{2}\cdot 89^{2}\sqrt{6\pi
n^{3}}}\left(\frac{3+\sqrt{17}}{2}\right)^{n-1}.$

(C) convex polyiamonds of perimeter $n$ is asymptotic to

$\frac{1280}{441\sqrt{3\pi n^{3}}}3^{n}.$

###### Key words and phrases:

Triangular lattice; polyiamonds; generating functions

###### 2010 Mathematics Subject Classification:

82B41, 05B50, 05A16, 05A15

##  1\. Introduction

An animal, living on a two dimensional lattice, is an edge-connected set of
basic two-dimensional polygon-like cells on a two-dimension lattice, where the
connectivity of two cells is defined by having a common edge. Among several
classes of animals are squared polyominoes or simply polyominoes. These are
animals that are made by edge-connected squares on the squared lattice and are
well studied in statistical physics, combinatorics, and discrete geometry.
They specifically play a role in enumeration of graphs [20], modeling the
mechanics of macro-molecules [30] such as the collapse of branched polymers
[26], percolation processes [9], and cell growth processes [14, 22, 28], to
name a few. In combinatorics, the related fundamental problem is to enumerate
lattice animals of a specific size $n$; see [3, 4, 2, 19, 31, 11, 29] on this
problem and references therein. Even for its simplest instance of enumeration
of general animals, not much is known. Therefore, various constraints such as
directional convexity and/or directional growth have been deployed to reduce
the difficulty of the enumeration process in squared and hexagonal lattices.
See [2, 10, 6, 7, 8, 12, 13, 21] for a few examples of the squared lattice and
[2, 15, 16, 18, 32, 23] for the hexagonal lattice, where several important
subsets of animals including bargraphs, column-convex animals, and convex
animals are considered. Although, a natural extension to these works is to
explore the animals residing on the triangular lattice (see [2]), we were not
able to find much work in this direction. Therefore, in this paper, we will
strive to close some of the existing gaps concerning these objects.

Polyiamonds, which are animals residing on the triangular lattice, have been
the subject of a few studies in recent years. For instance, in one of the
first rigorous work on polyiamonds, Yang and Wilson [33] proved the minimum
perimeter of a polyiamond with $n$ cells is either $\lceil\sqrt{6n}\rceil$ or
$\lceil\sqrt{6n}\rceil+1$ depending on the parity of $n$. To prove this
result, they obtained a lower bound on the perimeter by considering maximal
polyiamonds and showed how to construct minimal polyiamonds that attain the
perimeter lower bounds. Extending in this direction, Malen and Roldán [24]
studied the maximum number of holes that a polyiamond with $n$ cells can
enclose. Also, they obtained the minimum number of cells required to construct
a polyiamond with $h$ holes. Very recently, in a series of works, by careful
analyses of specific composition and concatenation operators on polyiamonds,
Shahleh [29] and Barequet et al. [5] improved the lower and upper bounds on
the asymptotic growth rate with respect to the size to $2.8424$ and $3.6050$,
respectively. Here, the growth rate of a sequence $a_{n}$ with respect to $n$
is simply $\lim_{n}a_{n}^{1/n}$ provided that this limit exists.

Our main results in this paper are related to the enumeration of several
classes of polyiamonds; including baryiamonds, column-convex polyiamonds, and
convex polyiamonds with respect to their perimeters. See Table 1 for an
elaboration. This is mainly done by careful recursive analyses of these
objects through column-by-column decomposition and cell pruning/growing
arguments which lead to explicit forms of defined generating functions. We
remark that one major difference between the triangular lattice and the
squared and hexagonal lattices is the lack of some of the symmetrical
structures in the former, which introduces an increasing level of complexity
to the analysis. To that goal, we begin by providing a few definitions and
notations. Equip the two-dimensional plane with an $x-y$ system, where the
angle between $x$-axis and $y$-axis is $\frac{\pi}{3}$ counter clock-wise.
Mark each of these axes with integer points. We then partition this plane with
equilateral triangles; each with a horizontal edge parallel to $x$-axis and
the sides of length one. We refer to these triangular building blocks as
triangular cells or simply cells. The partition is uniquely defined such that
there is one cell in the first quarter with a vertex at the intersection of
the axes (i.e. origin) and two edges; coinciding with $x-y$ axes. We refer to
this lattice, by ${\mathcal{T}}.$ With some abuse of conventional notation,
let the line $x=j$ denote the line parallel with the $y$-axis, intersection
the $x$-axis at $j$. Define $y=j$ in a similar fashion. Considering the
orientation of cells we identify two types of cells: (1) down-cells where
their horizontal edges are at their bases, and (2) up-cells whose horizontal
sides are on the top.

A polyiamond is a finite edge-connected set $\nu$ of cells in ${\mathcal{T}}$,
where for each distinct pair of cells $\Delta_{1}$ and $\Delta_{2}$ in $\nu$,
there is a finite consecutive sequence of edge-connected cells in $\nu$
connecting $\Delta_{1}$ and $\Delta_{2}$. A cell in $\nu$ with at least one
edge in common with a cell in $\nu^{c}:={\mathcal{T}}\setminus\nu$ is referred
to as a boundary cell of $\nu$ and the common edge is referred to as a
boundary edge. The perimeter of a polyiamond $\nu$ is the total number of its
boundary edges. A column (resp. row) of a polyiamond $\nu$ is a non-empty
subset of cells in $\nu$, between lines $x=j$ and $x=j+1$ (resp. $y=j$ and
$y=j+1$) for some $j\in\mathbb{Z}$. A polyiamond is called column-convex
(resp. row-convex) if each of its columns (resp. rows) is a single contiguous
block of cells. A convex polyiamond is both column-convex and row-convex.
Depending on the orientation of the top and bottom cells of a convex column,
we identify four types of convex columns as seen in Figure 1.

Figure 1. From left to right, convex columns of type one, two, three, and four

A baryiamond is a column-convex polyiamond when either (a) it has only one
convex column or (b) it has more than one column and the bottom cells of all
consecutive columns are edge-connected. We use $\mathcal{B},$
$\mathcal{CCP},$, and $\mathcal{CP}$ to denote the sets of all unique
baryiamonds, column-convex polyiamonds, and convex polyiamonds up-to
translation in the triangular lattice. Clearly,
$\mathcal{B}\subset\mathcal{CCP}$, $\mathcal{CP}\subset\mathcal{CCP},$ and
$\mathcal{B}\setminus\mathcal{CP}$ and $\mathcal{CP}\setminus\mathcal{B}$ are
both non-empty. See Figure 2 for an example of elements in these sets.

Figure 2. An example of a (left) baryiamond; (middle) column-convex
polyiamond; (right) convex polyiamond

The rest of the paper is organized as follows. Section 2 is dedicated to the
analysis of $\mathcal{B}.$ In order to enumerate the number of baryiamonds of
perimeter $n$, we study the generating function

$\displaystyle B(p):=\sum_{\nu\in\mathcal{B}}p^{per(\nu)},$ (1.1)

where $per(\nu)$ denotes the perimeter of the polyiamond $\nu$. We obtain an
explicit form for (1.1). In particular, we prove

######  Theorem 1.1.

$B(p)=\frac{1+p-2p^{2}-4p^{3}-3p^{4}+p^{5}-(1+p)^{2}\sqrt{(1-p)(1-p-2p^{2}-2p^{3}+p^{4}-p^{5})}}{2p^{2}}.$

Moreover, the number of baryiamonds with perimeter $n$ is asymptotic to

$\frac{(\xi+1)^{2}\sqrt{\xi^{4}+\xi^{3}-2\xi+1}}{2\sqrt{\pi
n^{3}}}\xi^{-n-2},$

where $\xi=0.44617150675\cdots$ is the smallest root of the polynomial
$1-p-2p^{2}-2p^{3}+p^{4}-p^{5}$.

In Section 3, we focus on the column-convex polyiamonds, where we study the
generating function

$\displaystyle C(p):=\sum_{\nu\in\mathcal{CCP}}p^{per(\nu)}.$ (1.2)

The main result of this section is

######  Theorem 1.2.

$\displaystyle C(p^{2})$
$\displaystyle=p^{6}+\frac{p^{4}(1-u_{+}^{2})(1-u_{-}^{2})A}{B}$
$\displaystyle=2p^{6}+3p^{8}+6p^{10}+15p^{12}+40p^{14}+113p^{16}+330p^{18}+988p^{20}+O(p^{22}),$

where $u_{+}$ and $u_{-}$ are given in (3.9)-(3.10) and $A$ and $B$ are given
in Lemma 3.1. Moreover, the number of column-convex polyiamonds with perimeter
$n$ is asymptotic to

$\displaystyle\frac{(17997809\sqrt{17}+3^{3}\cdot 13\cdot
175463)\sqrt{95\sqrt{17}-119}}{2^{7}\cdot 43^{2}\cdot 89^{2}\sqrt{6\pi
n^{3}}}\left(\frac{3+\sqrt{17}}{2}\right)^{n-1}.$

Section 4 is finally devoted to the study of convex polyiamonds and the
generating function

$\displaystyle F(p):=\sum_{\nu\in\mathcal{CP}}p^{per(\nu)}.$ (1.3)

We find the exact expression for this generating function and prove that

######  Theorem 1.3.

$\displaystyle
F(p)=-\frac{(2p+1)(p-1)^{2}(p+1)^{3}p^{3}}{(4p^{2}+2p-1)(4p^{5}-7p^{3}+p^{2}+4p-1)}\sqrt{1-2p-3p^{2}}$
$\displaystyle+\frac{(2p+1)(p-1)(p+1)^{2}p^{5}}{2(2p^{4}-p^{3}-2p^{2}+p+1)(5p^{3}+4p^{2}-p-1)(2p-1)}\sqrt{1-4p^{2}}$
$\displaystyle+\frac{1}{8}p^{3}-\frac{67}{40}p^{2}-\frac{103}{50}p-\frac{15071}{8000}+\frac{30+35p}{352(4p^{2}+2p-1)}$
$\displaystyle+\frac{318p^{9}+2083p^{8}+275p^{7}-3632p^{6}-1050p^{5}+2721p^{4}+680p^{3}-1093p^{2}-164p+169}{88(4p^{10}-8p^{9}+p^{8}+19p^{7}-8p^{6}-20p^{5}+9p^{4}+10p^{3}-5p^{2}-2p+1)}$
$\displaystyle+\frac{-238p^{4}+3p^{3}+537p^{2}+214p-149}{704(4p^{5}-7p^{3}+p^{2}+4p-1)}$
$\displaystyle+\frac{119p^{2}+114p+19}{500(5p^{3}+4p^{2}-p-1)}+\frac{3p^{3}+2p^{2}-p-1}{8(2p^{4}-p^{3}-2p^{2}+p+1)}$
$\displaystyle=2p^{3}+3p^{4}+6p^{5}+15p^{6}+38p^{7}+102p^{8}+272p^{9}+739p^{10}+2006p^{11}+O(p^{12}).$

In addition, the number of convex polyiamonds with perimeter $n$ is asymptotic
to

$\frac{1280}{441\sqrt{3\pi n^{3}}}3^{n}.$

At last, we summarized our findings along with known results for other
lattices in Table 1 for the sake of comparison. We remark that most of the
calculations in this paper are completed by extensive use of Maple and we
chose to not include some of the details for the sake of brevity and space.

Lattice | Animals sub-family |  The number of animals in sub-family with perimeter $n$ | Reference
---|---|---|---
Triangular | Bargraphs (Baryiamonds) | $\frac{(\xi+1)^{2}\sqrt{\xi^{4}+\xi^{3}-2\xi+1}}{2\sqrt{\pi n^{3}}}\xi^{-n-2}$ |  Theorem 1.1
| Column-Convex Animals | $\frac{(17997809\sqrt{17}+3^{3}\cdot 13\cdot 175463)\sqrt{95\sqrt{17}-119}}{2^{7}\cdot 43^{2}\cdot 89^{2}\sqrt{6\pi n^{3}}}\left(\frac{3+\sqrt{17}}{2}\right)^{n-1}$ |  Theorem 1.2
| Convex Animals | $\frac{1280}{441\sqrt{3\pi n^{3}}}3^{n}$ |  Theorem 1.3
Squared | Bargraphs | $\sqrt{\frac{2(1-\rho-\rho^{3})}{\pi\rho n^{3}}}\rho^{-n/2}$ |
|  |  with $\rho=\frac{1}{3}\left(-1-\frac{2^{8/3}}{(13+3\sqrt{3})^{1/3}}+2^{1/3}(13+3\sqrt{33})^{1/3}\right)$ | [27]
| Column-Convex Animals | $\frac{2\sqrt{5\sqrt{2}-7}((138\sqrt{2}+444)\sqrt{10-7\sqrt{2}}+589\sqrt{2}-410)}{2209}\cdot\frac{(1+\sqrt{2})^{n}}{n\sqrt{\pi n/2}}$ | [16]
| Convex Animals | $(n+3)2^{n-8}-4(n-7){n-8\choose n/2-4}$ | [13]
Hexagonal | Bargraphs | $\frac{\sqrt{3}(13-\sqrt{17})}{\sqrt{17\pi n^{3}}}\sqrt{2}^{\,n}$ | [25]
| Column-Convex Animals | $\frac{32(179\sqrt{7}-140)\sqrt{2}}{22743\sqrt{3\pi n^{3}}}\sqrt{3}^{n}$ |  g.f. reported in [16]
| Convex Animals | Unknown |  g.f. reported in [18]
Table 1.  Summary of enumeration of bargraphs, column-convex animals, and
convex animals on different lattices with respect to the perimeter. The
reported results for convex animals on the squared lattice are exact numbers,
whereas all the others are asymptotic results.

##  2\. Baryiamonds

In this section, we enumerate the baryiamonds with respect to their perimeter.
Recall the types of convex columns in Figure 1. Let $\mathcal{B}_{k}^{(i)}$ be
the set of all baryiamonds in $\mathcal{B}$ whose first columns are of type
$i$ and have $k$ up-cells. Set $\mathcal{B}^{(i)}=\cup_{k\geq
1}\mathcal{B}_{k}^{(i)}.$ We define $B^{(i)}_{k}(p)$ to be the generating
function with respect to perimeter over $\mathcal{B}_{k}^{(i)}$; that is,

$\displaystyle B^{(i)}_{k}(p)=\sum_{\nu\in\mathcal{B}_{k}^{(i)}}p^{per(\nu)},$

where $per(\nu)$ is the perimeter of $\nu.$ We similarly define the generating
function with respect to perimeter over $\mathcal{B}^{(i)}$ for $i=1,2,3,4$:

$\displaystyle B^{(i)}(p,u)=\sum_{k\geq
1}\sum_{\nu\in\mathcal{B}_{k}^{(i)}}u^{k}p^{per(\nu)},$

where $u$ marks the number of up-cells in the first column. It is clear that
$B(p)$, as defined in (1.1), is $\sum_{i=1}^{4}B^{(i)}(p,1).$ In the next few
paragraphs, we will obtain explicit expressions for $B^{(i)}(p,u)$ and prove
Theorem (1.1). The type of methods used to derive these expressions are
referred to cell growing/pruning methods throughout this paper. First, note
that each $\nu\in\mathcal{B}_{k}^{(1)}$ falls exclusively into one of the
cases described in Figure 3. Adding up the contribution of each of these
cases, we arrive at

$\displaystyle B_{k}^{(1)}(p)=p^{2k+2}$
$\displaystyle+\sum_{j=1}^{k}p^{2k+2-2j}(B^{(1)}_{j}(p)+B^{(2)}_{j-1}(p))$
$\displaystyle+\sum_{j\geq k+1}p^{2}(B^{(1)}_{j}(p)+B^{(2)}_{j-1}(p)).$ (2.1)

Figure 3. Various possibilities of baryiamonds in $\mathcal{B}_{k}^{(1)}$,
according to the type of their second column (if any)

Next, we describe $B^{(i)}_{k}(p)$, $i=2,3,4$, with respect to
$B^{(1)}_{k}(p)$. To that goal, note that

  * (i)

By removing the top down-cell of the first column in
$\nu\in\mathcal{B}_{k}^{(2)}$ with more than one cell, we get a bargraph in
$\mathcal{B}_{k}^{(1)}$; and hence for all $k\geq 1$,

$\displaystyle B^{(2)}_{k}(p)=pB^{(1)}_{k}(p);$ (2.2)

  * (ii)

By adding a down-cell to the bottom of the first column of
$\mathcal{B}_{k}^{(3)}$, we get a bargraph in $\mathcal{B}_{k}^{(1)}$; hence,
for all $k\geq 1$,

$\displaystyle B^{(3)}_{k}(p)=\frac{1}{p}B^{(1)}_{k}(p);$ (2.3)

  * (iii)

By removing the top down-cell of the first column and adding a down-cell to
the bottom of the same column of the $\nu\in\mathcal{B}_{k}^{(4)}$, we see
that, for all $k\geq 1,$

$\displaystyle B^{(4)}_{k}(p)=B^{(1)}_{k}(p).$ (2.4)

  * (iv)

The initial conditions are given by

$\displaystyle B_{0}^{(2)}(p)=p^{3},\quad B^{(i)}_{0}(p)=0\quad i=1,3,4.$
(2.5)

By multiplying each of equations in (2)-(2.5) by $u^{k}$ and adding up over
$k\geq 1$, we get

$\displaystyle
B^{(1)}(p,u)=\frac{p^{4}u}{1-p^{2}u}+\frac{p^{2}}{1-p^{2}u}(B^{(1)}(p,u)+uB^{(2)}(p,u))$
$\displaystyle\qquad\qquad\quad+\frac{p^{2}}{1-u}(uB^{(1)}(p,1)+uB^{(2)}(p,1)-B^{(1)}(p,u)-uB^{(2)}(p,u)),$
$\displaystyle B^{(2)}(p,u)=p^{3}+pB^{(1)}(p,u),\quad
B^{(3)}(p,u)=\frac{1}{p}B^{(1)}(p,u),$ $\displaystyle
B^{(4)}(p,u)=B^{(1)}(p,u).$ (2.6)

This implies

$\displaystyle B^{(1)}(p,u)$
$\displaystyle=\frac{p^{4}u}{1-p^{2}u}+\frac{p^{2}}{1-p^{2}u}((1+pu)B^{(1)}(p,u)+p^{3}u)$
$\displaystyle+\frac{p^{2}}{1-u}(u(1+p)B^{(1)}(p,1)-(1+pu)B^{(1)}(p,u)).$
(2.7)

This functional equation is then solved by the kernel method. See [1] for an
introduction to this method. To do so, we re-arrange the equation such that
all the terms containing $B^{(1)}(p,u)$ are on the left side and the rest of
terms are collected on the right side. Setting $u$ to

$u_{0}=\frac{1+p^{4}-(1+p)\sqrt{(1-p)(1-p-2p^{2}-2p^{3}+p^{4}-p^{5})}}{2p^{2}(1+p-p^{3})}$

will diminish the left hand side, and therefore, we obtain

$\displaystyle B^{(1)}(p,1)$
$\displaystyle=\frac{1-p-p^{2}-p^{3}-\sqrt{(1-p)(1-p-2p^{2}-2p^{3}+p^{4}-p^{5})}}{2p}.$
(2.8)

Then, by substituting (2.8) into the last three equations of (2) and adding
them up, we complete the proof of the first part of Theorem 1.1. The second
part is easily obtained by considering the singularity analysis (for example,
see Chapter 6 in [17]) applied to the generating function $B(p)$.

##  3\. Column-convex polyiamonds

This section is devoted to the enumeration of column-convex polyiamonds with
respect to their perimeter. Let $\mathcal{CCP}_{k}^{(i)}$ be the set of all
polyiamonds in $\mathcal{CCP}$ where their first columns are of type $i$ and
have $k$ up-cells. We define $C^{(i)}_{k}(p)$ to be the generating function
with respect to perimeter over $\mathcal{CCP}_{k}^{(i)}$; that is,

$\displaystyle
C^{(i)}_{k}(p)=\sum_{\nu\in\mathcal{CCP}_{k}^{(i)}}p^{per(\nu)}.$

Moreover, for $i=1,2,3,4$, we set

$C^{(i)}(p,u)=\sum_{k\geq 0}C_{k}^{(i)}(p)u^{k}.$

With these notations, $C(p)$ as defined by (1.2) is exactly
$\sum_{i=1}^{4}C^{(i)}(p,1).$ Therefore, we focus on obtaining $C^{(i)}(p,1)$
for $i=1,2,3,4$. First, let us review a few observations that relate
$C_{k}^{(i)}(p),$ $i=2,3,4,$ to $C_{k}^{(1)}(p).$

  * (i)

By removing the top down-cell of the first column of
$\nu\in\mathcal{CCP}_{k}^{(2)}$, we get another column-convex polyiamond in
$\mathcal{CCP}_{k}^{(1)}$; and hence for all $k\geq 1,$ we have

$\displaystyle C^{(2)}_{k}(p)=pC^{(1)}_{k}(p);$ (3.1)

  * (ii)

By adding a down-cell to the bottom of the first column of
$\nu\in\mathcal{CCP}_{k}^{(3)}$, we get another column-convex polyiamond in
$\mathcal{CCP}_{k}^{(1)}$; hence, for all $k\geq 1$

$\displaystyle C^{(3)}_{k}(p)=\frac{1}{p}C^{(1)}_{k}(p);$ (3.2)

  * (iii)

By removing the top down-cell of the first column and adding a down-cell to
the bottom of the same column of $\nu\in\mathcal{CCP}_{k}^{(4)}$, we see that,
for all $k\geq 1,$

$\displaystyle C^{(4)}_{k}(p)=C^{(1)}_{k}(p);$ (3.3)

  * (iv)

The initial conditions are given by

$\displaystyle C^{(2)}_{0}(p)=p^{3},\quad C^{(j)}_{0}(p)=0,\quad j=1,3,4.$
(3.4)

Multiplying the equations (3.1)-(3.4) by $u^{k}$ and adding all the terms up
over $k\geq 0$, we arrive at the system of equations

$\displaystyle C^{(2)}(p,u)=p^{3}+pC^{(1)}(p,u),\quad
C^{(3)}(p,u)=\frac{1}{p}C^{(1)}(p,u),\quad C^{(4)}(p,u)=C^{(1)}(p,u).$ (3.5)

Hence, to calculate $C(p),$ it is enough to find the explicit form of
$C^{(1)}(p,u).$ To that goal, we a write a recurrence for $C^{(1)}_{k}(p)$ by
investigating the types of first and second (if any) columns and their
relative positions to each other. See Figure 4 for a few examples of possible
cases. Let $\nu\in\mathcal{CCP}^{(1)}_{k}$ with $k\geq 2$.

Case (1)Case (2)Case (3)Case (4)Case
(5)Type(1),(2)(3),(4)Type(1),(2)(3),(4)Type(1),(2)(3),(4)Type(1),(2)(3),(4)
Figure 4. Various possibilities of polyiamonds in $\mathcal{CCP}_{k}^{(1)}$
with respect to the relative position of the first and second columns; Only
the cases whose second columns are of type $1$ are shown (See Figure 1). The
second column can be replaced with any other types mentioned in the picture.

Adding up the contributions of all cases in Figure 4, we obtain

$\displaystyle C_{k}^{(1)}(p)=p^{2k+2}$
$\displaystyle+\sum_{j=1}^{k}(k+1-j)p^{2k+2-2j}(C_{j}^{(1)}(p)+C_{j-1}^{(2)}(p)+C_{j+1}^{(3)}(p)+C_{j}^{(4)}(p))$
$\displaystyle+2\sum_{j=2}^{k}\frac{p^{2k+4-2j}-p^{2k+2}}{1-p^{2}}(C_{j}^{(1)}(p)+C_{j-1}^{(2)}(p)+C_{j+1}^{(3)}(p)+C_{j}^{(4)}(p))$
$\displaystyle+2\sum_{j\geq
k+1}\frac{p^{2}-p^{2k+2}}{1-p^{2}}(C_{j}^{(1)}(p)+C_{j-1}^{(2)}(p)+C_{j+1}^{(3)}(p)+C_{j}^{(4)}(p))$
$\displaystyle+\sum_{j\geq
k+2}(j-k-1)p^{2}(C_{j}^{(1)}(p)+C_{j-1}^{(2)}(p)+C_{j+1}^{(3)}(p)+C_{j}^{(4)}(p)),$

for all $k\geq 2$.

By multiplying by $u^{k}$ and summing over $k\geq 2$, we then obtain

$\displaystyle C^{(1)}(p,u)$ $\displaystyle=\frac{p^{4}u}{1-p^{2}u}$
$\displaystyle+\frac{p^{2}}{(1-p^{2}u)^{2}}(C^{(1)}(p,u)+uC^{(2)}(p,u)+\frac{1}{u}(C^{(3)}(p,u)-uC^{(3)}_{1}(p))+C^{(4)}(p,u))$
$\displaystyle-\frac{2p^{2}}{(1-p^{2}u)(1-u)}(C^{(1)}(p,u)+uC^{(2)}(p,u)+\frac{1}{u}(C^{(3)}(p,u)-uC^{(3)}_{1}(p))+C^{(4)}(p,u))$
$\displaystyle+\frac{2p^{2}u}{(1-u)(1-p^{2}u)}(C^{(1)}(p,1)+C^{(2)}(p,1)+C^{(3)}(p,1)-C^{(3)}_{1}(p)+C^{(4)}(p,1))$
$\displaystyle+\frac{p^{2}}{(1-u)^{2}}(C^{(1)}(p,u)+uC^{(2)}(p,u)+\frac{1}{u}(C^{(3)}(p,u)-uC^{(3)}_{1}(p))+C^{(4)}(p,u))$
$\displaystyle-\frac{p^{2}u(2-u)}{(1-u)^{2}}(C^{(1)}(p,1)+C^{(2)}(p,1)+C^{(3)}(p,1)-C^{(3)}_{1}(p)+C^{(4)}(p,1))$
$\displaystyle+\frac{p^{2}u}{1-u}\frac{\partial}{\partial
u}(C^{(1)}(p,u)+uC^{(2)}(p,u)+\frac{1}{u}(C^{(3)}(p,u)-uC^{(3)}_{1}(p))+C^{(4)}(p,u))\mid_{u=1}.$

To compute $C^{(1)}(p^{2},1),$ we next define
$D^{(i)}(p,u):=C^{(i)}(p^{2},u^{2})$. By (3.5), the last equation can be
written as

$\displaystyle\left(1-\frac{p^{2}u^{2}(1+p^{2}u^{2})^{2}(1-p^{4})^{2}}{(1-p^{4}u^{2})^{2}(1-u^{2})^{2}}\right)D^{(1)}(p,u)$
$\displaystyle=\frac{p^{2}u^{2}(1+p^{2})(p^{6}u^{2}-2p^{4}u^{4}+3p^{4}u^{2}-2p^{2}u^{2}+p^{2}-1)}{(1-p^{4}u^{2})(1-u^{2})^{2}}D^{(1)}(p,1)$
(3.6)
$\displaystyle-\frac{p^{10}u^{4}}{(1-p^{4}u^{2})^{2}}C^{(1)}_{1}(p^{2})+\frac{p^{2}u^{2}(1+p^{2})^{2}}{1-u^{2}}\frac{\partial}{\partial
u}C^{(1)}(p^{2},u)\mid_{u=1}-\frac{(p^{4}u^{2}-p^{2}-1)p^{8}u^{2}}{(1-p^{4}u^{2})^{2}}.$

Recall that $D^{(1)}(p,0)=C^{(1)}(p^{2},0)=0.$ We take the derivative at
$u=0$, and solve for $C^{(1)}_{1}(p^{2})$. This leads to

$\displaystyle
C^{(1)}_{1}(p^{2})=(1+p)p^{4}+p(p^{2}-1)D^{(1)}(p,1)+p(p+1)^{2}\frac{\partial}{\partial
u}C^{(1)}(p^{2},u)\mid_{u=1}.$ (3.7)

Then, substituting (3.7) into (3), we arrive at

$\displaystyle\left(1-\frac{p^{2}u^{2}(1+p^{2}u^{2})^{2}(1-p^{4})^{2}}{(1-p^{4}u^{2})^{2}(1-u^{2})^{2}}\right)D^{(1)}(p,u)$
$\displaystyle=\biggl{(}\frac{p^{2}u^{2}(1+p^{2})(1+(2u^{2}-1)p^{2}+2u^{2}(u^{2}-2)p^{4}-2u^{4}p^{6})}{(1-p^{4}u^{2})^{2}(1-u^{2})^{2}}$
$\displaystyle\qquad+\frac{p^{10}u^{4}(1+p^{2})(u^{2}(3-2u^{2})-(u^{4}-3u^{2}+1)p^{2}+(1-u^{2})^{2}p^{4})}{(1-p^{4}u^{2})^{2}(1-u^{2})^{2}}\biggr{)}D^{(1)}(p,1)$
(3.8)
$\displaystyle+\frac{p^{2}u^{2}(1+p^{2})^{2}(1-2u^{2}p^{4}+u^{4}p^{8}-u^{2}(1-u^{2})p^{10})}{(1-p^{4}u^{2})^{2}(1-u^{2})}\frac{\partial}{\partial
u}C^{(1)}(p^{2},u)\mid_{u=1}$
$\displaystyle+\frac{p^{8}u^{2}(1+p^{2}-u^{2}p^{4}-u^{2}p^{10}-u^{2}p^{12})}{(1-p^{4}u^{2})^{2}}.$

Again, we deploy the kernel method to solve this functional equation. To that
end, let
$K(u)=1-\frac{p^{2}u^{2}(1+p^{2}u^{2})^{2}(1-p^{4})^{2}}{(1-p^{4}u^{2})^{2}(1-u^{2})^{2}}$
and set

$\displaystyle
u_{+}^{2}=u_{+}^{2}(p)=\frac{p^{10}-2p^{6}+2p^{4}+p^{2}+2}{4p^{4}}-\frac{(1-p^{4})(1+p^{2})}{4p^{3}}\sqrt{1+p^{2}}\sqrt{4-3p^{2}+p^{4}}$
(3.9)
$\displaystyle-\frac{1-p^{4}}{2p^{4}}\sqrt{1+p^{2}(1+p^{2})^{2}+\frac{p(p^{10}-2p^{6}+2p^{4}+p^{2}+2)}{2(1-p^{2})}(p-p^{3}-\sqrt{1+p^{2}}\sqrt{4-3p^{2}+p^{4}})}$
$\displaystyle=1+p+\frac{1}{2}p^{2}+\frac{9}{8}p^{3}+2p^{4}+\frac{239}{128}p^{5}+O(p^{6}),$
$\displaystyle u_{-}^{2}=u_{-}^{2}(p)=u_{+}^{2}(-p),$ (3.10)

where $u_{\pm}$ are the roots of $K(u)=0$. Therefore,

$\displaystyle\biggl{(}1+(2u_{\pm}^{2}-1)p^{2}+2u_{\pm}^{2}(u_{\pm}^{2}-2)p^{4}-2u_{\pm}^{4}p^{6}$
$\displaystyle\qquad+p^{8}u_{\pm}^{2}(u_{\pm}^{2}(3-2u_{\pm}^{2})-(u_{\pm}^{4}-3u_{\pm}^{2}+1)p^{2}+(1-u_{\pm}^{2})^{2}p^{4})\biggr{)}D^{(1)}(p,1)$
$\displaystyle=(1+p^{2})(1-u_{\pm}^{2})(1-2u_{\pm}^{2}p^{4}+u_{\pm}^{4}p^{8}-u_{\pm}^{2}(1-u_{\pm}^{2})p^{10})\frac{\partial}{\partial
u_{\pm}}C^{(1)}(p^{2},u_{\pm})\mid_{u_{\pm}=1}$
$\displaystyle+\frac{p^{6}(1-u_{\pm}^{2})^{2}(1+p^{2}-u_{\pm}^{2}p^{4}-u_{\pm}^{2}p^{10}-u_{\pm}^{2}p^{12})}{1+p^{2}}.$

Solving this system of equations, we obtain an explicit formula for
$C^{(1)}(p^{2},1)$. More specifically,

######  Lemma 3.1.

The generating function $C^{(1)}(p^{2},1)$ is given by

$\displaystyle\frac{p^{6}(1-u_{+}^{2})(1-u_{-}^{2})A}{(1+p^{2})^{2}B},$

where

$\displaystyle A$
$\displaystyle=2p^{6}+p^{4}-p^{2}-1-p^{4}(p^{6}+p^{4}-1)(u_{+}^{2}+u_{-}^{2})$
$\displaystyle\qquad\qquad+p^{8}(p^{12}+p^{10}-2p^{8}-2p^{6}+2p^{4}+2p^{2}-1)u_{+}^{2}u_{-}^{2},$
$\displaystyle B$
$\displaystyle=1-2p^{4}+2p^{8}(u_{+}^{2}+u_{-}^{2})-p^{8}(u_{+}^{4}+u_{-}^{4})-2p^{4}(p^{14}-2p^{10}+p^{8}+p^{6}+1)u_{+}^{2}u_{-}^{2}$
$\displaystyle\qquad\qquad+p^{8}(p^{10}-2p^{6}+p^{2}+2)u_{+}^{2}u_{-}^{2}(u_{+}^{2}+u_{-}^{2})+p^{12}(p^{4}-2)u_{+}^{4}u_{-}^{4}.$

This lemma along with (3.5) completes the proof of the first part in Theorem
1.2. To show the second part of Theorem 1.2, we apply the singularity analysis
to the generating function $C(p^{2})$. Note that $C(p^{2})$ can be further
rewritten as

$\frac{a_{00}+a_{10}\sqrt{\alpha}+a_{01}\sqrt{\beta}+a_{11}\sqrt{\alpha\beta}}{2(4p^{22}-16p^{20}+20p^{18}+32p^{16}-75p^{14}-76p^{12}+182p^{10}+152p^{8}-155p^{6}-164p^{4}+48p^{2}+64)},$

where

$\displaystyle a_{00}$
$\displaystyle=(1-p^{2})\bigl{(}14p^{24}-34p^{22}+3p^{20}+143p^{18}-21p^{16}-377p^{14}-57p^{12}+449p^{10}$
$\displaystyle+181p^{8}-245p^{6}-112p^{4}+72p^{2}+32\bigr{)},$ $\displaystyle
a_{10}$
$\displaystyle=\frac{\sqrt{2}}{4}(p^{4}-1)\bigl{(}4p^{24}-12p^{22}+2p^{20}+54p^{18}-19p^{16}-138p^{14}+54p^{12}+250p^{10}$
$\displaystyle-3p^{8}-234p^{6}-62p^{4}+80p^{2}+32\bigr{)}$
$\displaystyle+\frac{\sqrt{2}p}{4}(p^{2}+1)(p^{2}-1)^{2}\bigl{(}4p^{18}-4p^{16}-10p^{14}+30p^{12}+63p^{10}-22p^{8}-99p^{6}$
$\displaystyle-28p^{4}+56p^{2}+32\bigr{)}\sqrt{p^{2}(p^{2}-1)^{2}+4},$
$\displaystyle a_{01}$
$\displaystyle=\frac{\sqrt{2}}{4}(p^{4}-1)\bigl{(}4p^{24}-12p^{22}+2p^{20}+54p^{18}-19p^{16}-138p^{14}+54p^{12}+250p^{10}$
$\displaystyle-3p^{8}-234p^{6}-62p^{4}+80p^{2}+32\bigr{)}$
$\displaystyle-\frac{\sqrt{2}p}{4}(p^{2}+1)(p^{2}-1)^{2}\bigl{(}4p^{18}-4p^{16}-10p^{14}+30p^{12}+63p^{10}-22p^{8}-99p^{6}$
$\displaystyle-28p^{4}+56p^{2}+32\bigr{)}\sqrt{p^{2}(p^{2}-1)^{2}+4},$
$\displaystyle a_{11}$
$\displaystyle=p^{2}(p^{2}+1)(p^{2}-1)^{2}(2p^{12}+p^{10}+4p^{6}+5p^{4}-3p^{2}-4),$
$\displaystyle\alpha$
$\displaystyle=p^{12}-2p^{8}+4p^{6}+5p^{4}+4p^{2}+2+\frac{p(p^{10}-2p^{6}+2p^{4}+p^{2}+2)}{1-p^{2}}\sqrt{p^{2}(p^{2}-1)^{2}+4},$
$\displaystyle\beta$
$\displaystyle=p^{12}-2p^{8}+4p^{6}+5p^{4}+4p^{2}+2-\frac{p(p^{10}-2p^{6}+2p^{4}+p^{2}+2)}{1-p^{2}}\sqrt{p^{2}(p^{2}-1)^{2}+4}.$

Consider $r:=\frac{\sqrt{\sqrt{17}-3}}{2}$. Then, $\alpha$ and $\beta$ can be
written as

$\displaystyle\alpha^{2}$
$\displaystyle=\frac{95\sqrt{17}-119}{6}r^{4}(1+p/r)+O((1+p/r)^{2}),$
$\displaystyle\beta^{2}$
$\displaystyle=\frac{95\sqrt{17}-119}{6}r^{4}(1-p/r)+O((1-p/r)^{2}).$

Consequently, when $p^{2}$ is near the dominant singularity of $C(p^{2})$,
$C(p^{2})$ can be approximated by

$\displaystyle-\frac{(17997809\sqrt{17}+61587513)\sqrt{2}}{937339456}\frac{\sqrt{95\sqrt{17}-119}r^{2}}{\sqrt{6}}\left(\sqrt{1+p/r}+\sqrt{1-p/r}\right),$

from which the second part of Theorem 1.2 follows by the singularity analysis.

##  4\. Convex polyiamonds

In this section, we enumerate the number of convex polyiamonds with respect to
their perimeter. To that end, we first define a few notations. We say that the
bottom (resp. top) cell of the $i$-th column of a polyiamond $\nu$ is at the
position $\ell$ if its lowest (resp. highest) point coincides with the line
$y=\ell$. Let $b(\nu,i)$ and $u(\nu,i)$ denote the position of the bottom and
top cells of the $i$-th column in the polyiamond $\nu$. Denote the number of
columns of $\nu$ by $col(\nu)$. We define $\mathcal{CP}_{k}^{bu(i)}$ to be the
set of all polyiamonds $\nu$ in $\mathcal{CP}$ where ($1$) the first column of
$\nu$ is of type $i$ and has $k$ up-cells, and ($2$) $b(\nu,j+1)>b(\nu,j)$ and
$u(\nu,j+1)<u(\nu,j)$ for all $j=1,2,\ldots,col(\nu)-1$. Also, we set
$\mathcal{CP}_{k}^{u(i)}$ to be the set of polyiamonds where the condition
($2$) is replaced by the condition ($2^{\prime}$) $u(\nu,j+1)<u(\nu,j)$ for
all $j=1,2,\ldots,col(\nu)-1$. Similarly, we define $\mathcal{CP}_{k}^{b(i)}$
so that, ($2$) is replaced by the condition; ($2^{\prime\prime}$)
$b(\nu,j+1)>b(\nu,j)$ for all $j=1,2,\ldots,col(\nu)-1$. Finally, let
$\mathcal{CP}_{k}^{(i)}$ to be the set of convex polyiamonds $\nu$ satisfying
only ($1$).

The enumeration of the convex polyiamonds with respect to perimeter is
conducted as follows. The first step is to perform the counting in
$\mathcal{CP}^{bu}:=\cup_{i=1}^{4}\cup_{k\geq 1}\mathcal{CP}_{k}^{bu(i)}$ in a
relatively straightforward manner. Then, we extend the result from
$\mathcal{CP}^{bu}$ to $\mathcal{CP}^{u}:=\cup_{i=1}^{4}\cup_{k\geq
1}\mathcal{CP}_{k}^{u(i)}$ and $\mathcal{CP}^{b}:=\cup_{i=1}^{4}\cup_{k\geq
1}\mathcal{CP}_{k}^{b(i)}$. This can be done since the generating functions
over these subsets can be written recursively in terms of the corresponding
generating functions in $\mathcal{CP}^{bu}$. Finally, for the last step, we
derive the results for $\mathcal{CP}$ by lifting up the result obtained for
$\mathcal{CP}^{u}$, $\mathcal{CP}^{b}$, and $\mathcal{CP}^{bu}$.

###  4.1. Enumeration over $\mathcal{CP}^{bu}$

Suppose $F_{k}^{bu(i)}(p)$ is the generating function of polyiamonds in
$\mathcal{CP}_{k}^{bu(i)}$ with respect to their perimeter; that is

$\displaystyle
F_{k}^{bu(i)}(p):=\sum_{\nu\in\mathcal{CP}_{k}^{bu(i)}}p^{per(\nu)}.$

Moreover, define $F^{bu(i)}(p,u)=\sum_{k\geq 0}F^{bu(i)}_{k}(p)u^{k}$ for all
$i=1,2,3,4$, where $u$ marks the number of up-cells in the first column.

Note that any $\nu\in\mathcal{CP}_{k}^{bu(1)}$, depending on the type of the
first and the second columns (if any) and their relative position, falls into
one of cases described in Figure 5.

Case (1)Case (2)Case (3) Figure 5. Decomposition of a polyiamond in
$\mathcal{CP}_{k}^{bu(1)}$ with respect the first and second columns relative
position

Adding up the contributions of these terms, we arrive at

$\displaystyle F^{bu(1)}_{k}(p)$
$\displaystyle=p^{2k+2}+\sum_{j=1}^{k}(k+1-j)p^{2k+2-2j}(F^{bu(1)}_{j}(p)+F^{bu(2)}_{j-1}(p)),$
(4.1)

for all $k\geq 1$. In addition, by a cell growing/pruning based argument,
similar to the one used in the previous sections, we obtain

$\displaystyle F^{bu(2)}_{k}(p)=pF^{bu(1)}_{k}(p),\quad
F^{bu(3)}_{k}(p)=\frac{1}{p}F^{bu(1)}_{k}(p),\quad k\geq 1,$ $\displaystyle
F^{bu(4)}_{k}(p)=F^{bu(1)}_{k}(p)\quad F^{bu(2)}_{0}(p)=p^{3},\quad
F^{bu(i)}_{0}(p)=0,\ k\geq 1,\,i=1,3,4.$ (4.2)

Multiplying by $u^{k}$ and adding up over $k\geq 1$, the recursions (4.1) and
(4.1) can be written as the system of equations

$\displaystyle
F^{bu(1)}(p,u)=\frac{p^{4}u}{(1-p^{2}u)^{2}}+\frac{p^{2}}{(1-p^{2}u)^{2}}(F^{bu(1)}(p,u)+uF^{bu(2)}(p,u)),$
$\displaystyle F^{bu(2)}(p,u)=p^{3}+pF^{bu(1)}(p,u),$ $\displaystyle
F^{bu(3)}(p,u)=\frac{1}{p}F^{bu(1)}(p,u),\quad F^{bu(4)}(p,u)=F^{bu(1)}(p,u).$
(4.3)

By solving this system, we obtain the following result.

######  Lemma 4.1.

$F^{bu(1)}(p,u)=\frac{p^{4}u(1+p-p^{2}u)}{1-p^{2}-2p^{2}u-p^{3}u+p^{4}u^{2}}.$

###  4.2. Enumeration over $\mathcal{CP}^{u}$

Suppose $F_{k}^{u(i)}(p)$ is the generating function of polyiamonds in
$\mathcal{CP}_{k}^{u(i)}$ with respect to their perimeter; that is

$\displaystyle
F_{k}^{u(i)}(p):=\sum_{\nu\in\mathcal{CP}_{k}^{u(i)}}p^{per(\nu)}.$

Moreover, define $F^{u(i)}(p,u)=\sum_{k\geq 0}F^{u(i)}_{k}(p)u^{k}$ for all
$i=1,2,3,4$, where $u$ marks the number of up-cells in the first column.

Note that, inspecting the type of the first and second columns (if any) and
their relative position, any $\nu\in\mathcal{CP}_{k}^{u(1)}$ belongs to one of
cases described in Figure 6.

Case (1)Case (2)Case (3)Case (4)Type(1),(2)Type(1),(2)Type(1),(2)(3),(4)
Figure 6. Decomposition of a polyiamond in $\mathcal{CP}_{k}^{u(1)}$ with
respect to the relative position of the first and second columns. Only the
cases whose second column is of type $1$ are shown. The rest of cases are
obtained by simply considering other types of second columns stated in the
picture.

Adding up the contributions of these cases, we get

$\displaystyle
F^{u(1)}_{k}(p)=p^{2k+2}+\sum_{j=1}^{k-1}(k-j)p^{2k+2-2j}(F^{bu(1)}_{j}(p)+F^{bu(2)}_{j-1}(p))$
$\displaystyle+\sum_{j=1}^{k}p^{2k+2-2j}(F^{u(1)}_{j}(p)+F^{u(2)}_{j-1}(p))$
(4.4)
$\displaystyle+\sum_{j=1}^{k}(p^{2k+2-2}+p^{2k+2-4}+\cdots+p^{2k+2-2(j-1)})(F^{u(1)}_{j}(p)+F^{u(2)}_{j-1}(p))$
$\displaystyle+\sum_{j=1}^{k}(p^{2k+2-2}+p^{2k+2-4}+\cdots+p^{2k+2-2j})(F^{u(3)}_{j+1}(p)+F^{u(4)}_{j}(p))$
$\displaystyle+\sum_{j\geq
k+1}(p^{2k+2-2}+p^{2k+2-4}+\cdots+p^{2})(F^{u(1)}_{j}(p)+F^{u(2)}_{j-1}(p)+F^{u(3)}_{j+1}(p)+F^{u(4)}_{j}(p)),$

for $k\geq 1$. Once again, we apply the cell growing/pruning based argument
and get

$\displaystyle F^{u(2)}_{k}(p)=pF^{u(1)}_{k}(p),\quad
F^{u(3)}_{k}(p)=\frac{1}{p}F^{u(1)}_{k}(p),\quad k\geq 1,$ $\displaystyle
F^{u(4)}_{k}(p)=F^{u(1)}_{k}(p),\quad F^{u(2)}_{0}(p)=p^{3},\quad
F^{u(i)}_{0}(p)=0,\quad k\geq 1,\,i=1,3,4.$ (4.5)

Multiplying by $u^{k}$ and summing up over $k\geq 1$, equations (4.2) and
(4.2) can be written as

$\displaystyle
F^{u(1)}(p,u)=\frac{p^{4}u}{1-p^{2}u}+\frac{p^{4}u}{(1-p^{2}u)^{2}}(F^{bu(1)}(p,u)+uF^{bu(2)}(p,u))$
$\displaystyle\quad-\frac{p^{2}u}{(1-u)(1-p^{2}u)}(F^{u(1)}(p,u)+uF^{u(2)}(p,u)+\frac{1}{u}F^{u(3)}(p,u)+F^{u(4)}(p,u))$
$\displaystyle\quad+\frac{p^{2}u}{(1-u)(1-p^{2}u)}(F^{u(1)}(p,1)+F^{u(2)}(p,1)+F^{u(3)}(p,1)+F^{u(4)}(p,1))$
(4.6)

Then by (4.1), we get

$\displaystyle F^{u(1)}(p,u)$
$\displaystyle=\frac{p^{4}(1+p)u}{1-p^{2}u}+\frac{p^{4}u}{(1-p^{2}u)^{2}}(F^{bu(1)}(p,u)+uF^{bu(2)}(p,u))$
(4.7)
$\displaystyle-\frac{p}{(1-u)(1-p^{2}u)}((1+pu)^{2}F^{u(1)}(p,u)-u(1+p)^{2}F^{u(1)}(p,1)).$

Again, we apply the kernel method to (4.7). To that goal, let
$u=u_{0}=\frac{1-p-\sqrt{1-2p-3p^{2}}}{2p^{2}}$ (which is the generating
function for Motzkin numbers) and use Lemma 4.1 and (4.1). Then (4.7) gives

$\displaystyle
F^{u(1)}(p,1)=\frac{p^{3}(3p^{2}+p-1+(p+1)\sqrt{1-2p-3p^{2}})}{(1-2p-4p^{2})(1+p)}.$
(4.8)

By substituting (4.8) into (4.7) and using Lemma 4.1 and (4.1), we obtain the
following.

######  Lemma 4.2.

$\displaystyle F^{u(1)}(p,u)$
$\displaystyle=\frac{p^{4}(p+1)(u\sqrt{1-2p-3p^{2}}+(1-p)u-2)}{(1-2p-4p^{2})(p^{2}u^{2}+pu-u+1)}-p^{2}$
$\displaystyle+\frac{p^{2}((u+2)p^{5}u-(3u+2)p^{4}-(2u+4)p^{3}+(u+3)p^{2}+2p-1)}{(p^{4}u^{2}-p^{3}u-2p^{2}u-p^{2}+1)(4p^{2}+2p-1)}$

This lemma and (4.2) imply that the generating function
$\sum_{i=1}^{4}F^{u(i)}(p,1)$ is:

$\displaystyle\frac{p^{2}(-1+p+2p^{2}-p^{3}+(p+1)\sqrt{1-2p-3p^{2}})}{1-2p-4p^{2}}.$

Case (1)Case (2)Case (3)Case (4)Type(1),(2)Type(1),(2)Type(1),(2) Figure 7.
Decomposition of a polyiamond in $\mathcal{CP}_{k}^{b(1)}$ with respect to the
relative position of the first and second columns. Only the cases whose second
column is of type $1$ are shown. The rest of cases are obtained by simply
considering other types of second columns stated in the picture.

###  4.3. Enumeration over $\mathcal{CP}^{b}$

Suppose $F_{k}^{b(i)}(p)$ is the generating function of polyiamonds in
$\mathcal{CP}_{k}^{b(i)}$ with respect to their perimeter; that is

$\displaystyle
F_{k}^{b(i)}(p):=\sum_{\nu\in\mathcal{CP}_{k}^{b(i)}}p^{per(\nu)}.$

Moreover, define $F^{b(i)}(p,u)=\sum_{k\geq 0}F^{b(i)}_{k}(p)u^{k}$ for all
$i=1,2,3,4$, where $u$ marks the number of up-cells in the first column.

Note that any $\nu\in\mathcal{CP}_{k}^{b(1)}$, depending on the type of the
first and second columns (if any) and their relative position, falls into one
of cases described in Figure 7. Adding up the contributions of these terms, we
get

$\displaystyle
F^{b(1)}_{k}(p)=p^{2k+2}+\sum_{j=1}^{k-1}(k-j)p^{2k+2-2j}(F^{bu(1)}_{j}(p)+F^{bu(2)}_{j-1}(p))$
$\displaystyle\quad+\sum_{j=1}^{k}p^{2k+2-2j}(F^{b(1)}_{j}(p)+F^{bu(2)}_{j-1}(p))$
$\displaystyle\quad+\sum_{j=1}^{k}(p^{2k+2-2}+p^{2k+3-4}+\cdots+p^{2k+2-2(j-1)})(F^{b(1)}_{j}(p)+F^{bu(2)}_{j-1}(p))$
$\displaystyle\quad+\sum_{j\geq
k+1}(p^{2k+2-2}+p^{2k+2-4}+\cdots+p^{2})(F^{b(1)}_{j}(p)+F^{bu(2)}_{j-1}(p)),$

for all $k\geq 1$. An application of the cell growing/pruning based argument
gives

$\displaystyle F^{b(2)}_{k}(p)=F^{bu(2)}_{k}(p),\quad
F^{b(3)}_{k}(p)=\frac{1}{p}F^{b(1)}_{k}(p),\quad k\geq 1,$ $\displaystyle
F^{b(4)}_{k}(p)=F^{bu(4)}_{k}(p),\quad F^{b(2)}_{0}(p)=p^{3},\quad
F^{b(i)}_{0}(p)=0,\quad k\geq 1,\,i=1,3,4.$

By multiplying by $u^{k}$ and summing up over $k\geq 1$, all these recursions
can be written as

$\displaystyle
F^{b(1)}(p,u)=\frac{p^{4}u}{1-p^{2}u}+\frac{p^{4}u}{(1-p^{2}u)^{2}}(F^{bu(1)}(p,u)+uF^{bu(2)}(p,u))$
$\displaystyle+\frac{p^{2}u}{(1-u)(1-p^{2}u)}(F^{b(1)}(p,1)+F^{bu(2)}(p,1)-F^{b(1)}(p,u)-uF^{bu(2)}(p,u))$
(4.9) $\displaystyle F^{b(2)}(p,u)=F^{bu(2)}(p,u),\
F^{b(3)}(p,u)=\frac{1}{p}F^{b(1)}(p,u),\ F^{b(4)}(p,u)=F^{bu(4)}(p,u).$

This implies

$\displaystyle
F^{b(1)}(p,u)=\frac{p^{4}u}{1-p^{2}u}+\frac{p^{4}u}{(1-p^{2}u)^{2}}(F^{bu(1)}(p,u)+uF^{bu(2)}(p,u))$
$\displaystyle\quad+\frac{p^{2}u}{(1-u)(1-p^{2}u)}(F^{b(1)}(p,1)+F^{bu(2)}(p,1)-F^{b(1)}(p,u)-uF^{bu(2)}(p,u)).$
(4.10)

We let $u=u_{0}=\frac{1-\sqrt{1-4p^{2}}}{2p^{2}}$ (the generating function for
the Catalan numbers) in (4.10) and use Lemma 4.1 and equation (4.1) to get

$\displaystyle F^{b(1)}(p,1)$ (4.11)
$\displaystyle=\frac{p^{3}((2p^{2}-1)(p^{5}+6p^{3}+5p^{2}-p-1)-(2p+1)(p^{3}-2p^{2}-p+1)(p+1)^{2}\sqrt{1-4p^{2}})}{2(p+1)(5p^{3}+4p^{2}-p-1)(p^{3}-2p^{2}-p+1)}.$

By substituting (4.11) into (4.10) and using Lemma 4.1, we get

######  Lemma 4.3.

$\displaystyle
F^{b(1)}(p,u)=\frac{-u(p+1)(2p+1)p^{5}}{2(p^{2}u^{2}-u+1)(5p^{3}+4p^{2}-p-1)}\sqrt{1-4p^{2}}-p^{2}$
$\displaystyle+\frac{p^{2}(up^{8}+2(u+1)p^{7}+(u+3)p^{6}-(5u+6)p^{5}-(4u+7)p^{4}+(u+5)(1+p)p^{2}-p-1)}{(p^{4}u^{2}-p^{3}u-2p^{2}u-p^{2}+1)(5p^{3}+4p^{2}-p-1)}$
$\displaystyle-\frac{(u-2)(2p+1)(p+1)p^{5}}{(p^{2}u^{2}-u+1)(10p^{3}+8p^{2}-2p-2)},$

An implication of this Lemma and (4.3) is that

$\displaystyle\sum_{i=1}^{4}F^{b(i)}(p,1)=\frac{p^{2}(2p^{4}-6p^{3}-7p^{2}+1-(2p+1)(p+1)^{2}\sqrt{1-4p^{2}})}{2(5p^{3}+4p^{2}-p-1)}.$

###  4.4. Enumeration over $\mathcal{CP}$

Suppose $F_{k}^{(i)}(p)$ is the generating function of polyiamonds in
$\mathcal{CP}_{k}^{(i)}$ with respect to their perimeter; that is

$\displaystyle
F_{k}^{(i)}(p):=\sum_{\nu\in\mathcal{CP}_{k}^{(i)}}p^{per(\nu)}.$

Moreover, define $F^{(i)}(p,u)=\sum_{k\geq 0}F^{(i)}_{k}(p)u^{k}$ for all
$i=1,2,3,4$, where $u$ marks the number of up-cells in the first column.

It is easy to see that each $\nu\in\mathcal{CP}_{k}^{(1)}$ belongs to one of
cases described in Figure 8.

Case (1)Case (2)Case (3)Case (4)Case
(5)Type(1),(2)Type(1),(2)Type(1),(2)Type(1),(2)Case (6)Case (7)Case
(8)Type(1),(2)(3),(4)Type(1),(2)Type(1),(2)(3),(4) Figure 8. Decomposition of
a polyiamond in $\mathcal{CP}_{k}^{(1)}$ with respect to the relative position
of the first and second columns. Only the cases whose second column is of type
$1$ are shown. The rest of cases are obtained by simply considering other
types of second columns stated in the picture.

Adding up the contributions of these terms, we get

$\displaystyle F^{(1)}_{k}(p)=p^{2k+2}+p^{2}(F^{(1)}_{k}(p)+F^{(2)}_{k-1}(p))$
$\displaystyle+\sum_{j=1}^{k-2}(k-1-j)p^{2k+2-2j}(F^{bu(1)}_{j}(p)+F^{bu(2)}_{j-1}(p))$
$\displaystyle+\sum_{j=1}^{k-1}p^{2k+2-2j}(F^{u(1)}_{j}(p)+F^{b(1)}_{j}(p)+F^{u(2)}_{j-1}(p)+F^{bu(2)}_{j-1}(p))$
$\displaystyle+\sum_{j=1}^{k}(p^{2k+2-2}+p^{2k+2-4}+\cdots+p^{2k+2-2(j-1)})(F^{u(1)}_{j}(p)+F^{u(2)}_{j-1}(p))$
$\displaystyle+\sum_{j=1}^{k}(p^{2k+2-2}+p^{2k+2-4}+\cdots+p^{2k+2-2j})(F^{u(3)}_{j+1}(p)+F^{u(4)}_{j}(p))$
$\displaystyle+\sum_{j\geq
k+1}(p^{2k+2-2}+p^{2k+2-4}+\cdots+p^{2})(F^{u(1)}_{j}(p)+F^{u(2)}_{j-1}(p))$
$\displaystyle+\sum_{j\geq
k+1}(p^{2k+2-2}+p^{2k+2-4}+\cdots+p^{2})(F^{u(3)}_{j+1}(p)+F^{u(4)}_{j}(p))$
$\displaystyle+\sum_{j=1}^{k}(p^{2k+2-2}+p^{2k+2-4}+\cdots+p^{2k+2-2(j-1)})(F^{b(1)}_{j}(p)+F^{bu(2)}_{j-1}(p))$
$\displaystyle+\sum_{j\geq
k+1}(p^{2k+2-2}+p^{2k+2-4}+\cdots+p^{2})(F^{b(1)}_{j}(p)+F^{bu(2)}_{j-1}(p))$
$\displaystyle+\sum_{j\geq
k+1}(j-1-k)p^{2}(F^{(1)}_{j}(p)+F^{u(2)}_{j-1}(p))+\sum_{j\geq
k+1}(j-k)p^{2}(F^{(3)}_{j+1}(p)+F^{u(4)}_{j}(p)),$

for all $k\geq 1$. As before, we apply the cell growing/pruning based argument
to get

$\displaystyle F^{(2)}_{k}(p)=pF^{u(1)}_{k}(p),\quad
F^{(3)}_{k}(p)=\frac{1}{p}F^{(1)}_{k}(p),\quad k\geq 1,$ $\displaystyle
F^{(4)}_{k}(p)=F^{u(1)}_{k}(p),\quad F^{(2)}_{0}(p)=p^{3},\quad
F^{(i)}_{0}(p)=0,\quad k\geq 1,\,i=1,3,4.$ (4.12)

By multiplying equations in (4.4) by $u^{k}$ and summing up over $k\geq 1$, we
obtain

$\displaystyle F^{(2)}(p,u)=p^{3}+pF^{u(1)}(p,u),\quad
F^{(3)}(p,u)=\frac{1}{p}F^{(1)}(p,u),\quad F^{(4)}(p,u)=F^{u(1)}(p,u),$ (4.13)

and

$\displaystyle(1+p)(1-p-\frac{p}{(1-u)^{2}})F^{(1)}(p,u)=\frac{p^{4}u}{1-p^{2}u}+\frac{p^{6}u^{2}}{(1-p^{2}u)^{2}}(F^{bu(1)}(p,u)+uF^{bu(2)}(p,u))$
$\displaystyle-\frac{(p^{2}u^{2}-p^{2}u+1)p^{2}}{(1-u)(1-p^{2}u)}(F^{u(1)}(p,u)+uF^{u(2)}(p,u)+F^{b(1)}(p,u)+uF^{bu(2)}(p,u))$
$\displaystyle+p^{2}uF^{(2)}(p,u)-\frac{p^{2}u}{(1-u)(1-p^{2}u)}(\frac{1}{u}F^{u(3)}(p,u)+F^{u(4)}(p,u))$
$\displaystyle+\frac{p^{2}u}{(1-u)(1-p^{2}u)}(F^{u(1)}(p,1)+F^{u(2)}(p,1)+F^{u(3)}(p,1)+F^{u(4)}(p,1))$
$\displaystyle+\frac{p^{2}u}{(1-u)(1-p^{2}u)}(F^{b(1)}(p,1)+F^{bu(2)}(p,1))$
$\displaystyle-\frac{p^{2}u(2-u)}{(1-u)^{2}}F^{(1)}(p,1)+\frac{2p^{2}u}{1-u}\frac{\partial}{\partial
u}F^{(1)}(p,u)\mid_{u=1}+\frac{p^{2}u}{1-u}\frac{\partial}{\partial
u}F^{u(4)}(p,u)\mid_{u=1}$
$\displaystyle+\frac{p^{2}u}{(1-u)^{2}}F^{(2)}(p,u)-\frac{p^{2}u(2-u)}{(1-u)^{2}}F^{(2)}(p,1)+\frac{p^{2}u}{1-u}\frac{\partial}{\partial
u}F^{(2)}(p,u)\mid_{u=1}$
$\displaystyle+\frac{p^{2}u}{(1-u)^{2}}F^{u(4)}(p,u)-\frac{pu}{(1-u)^{2}}F^{(1)}(p,1)-\frac{p^{2}u}{(1-u)^{2}}F^{u(4)}(p,1).$

We set $u^{\pm}=v_{1}=1\pm\sqrt{\frac{p}{1-p}}$ and use Lemmas 4.1, 4.2 and
4.3, and obtain a system of two equations with two variables $F^{(1)}(p,1)$
and $\frac{\partial}{\partial u}F^{(1)}(p,u)\mid_{u=1}$. Solving this system,
we get the following result. We omit the lengthy calculations for the sake of
brevity.

######  Lemma 4.4.

The generating function $F^{(1)}(p,1)$ is given by

$\displaystyle\frac{p^{4}(-2p^{5}+p^{4}+3p^{3}-3p^{2}-2p+2)}{(4p^{2}+2p-1)(4p^{5}-7p^{3}+p^{2}+4p-1)}(p-\sqrt{1-2p-3p^{2}})$
$\displaystyle+\frac{p^{6}(p^{2}-1)(2p+1)}{2(2p^{4}-p^{3}-2p^{2}+p+1)(5p^{3}+4p^{2}-p-1)(2p-1)}\sqrt{1-4p^{2}}$
$\displaystyle-\frac{9}{5}p^{2}-\frac{13}{50}p-\frac{1577}{1000}+\frac{1}{p+1}$
$\displaystyle+\frac{570p^{9}+1009p^{8}-671p^{7}-1764p^{6}+210p^{5}+1251p^{4}-4p^{3}-459p^{2}-20p+63}{88(4p^{10}-8p^{9}+p^{8}+19p^{7}-8p^{6}-20p^{5}+9p^{4}+10p^{3}-5p^{2}-2p+1)}$
$\displaystyle+\frac{24p^{4}+7p^{3}-p^{2}+8p-3)}{88(4p^{5}-7p^{3}+p^{2}+4p-1)}+\frac{-p^{2}+19p+24}{500(5p^{3}+4p^{2}-p-1)}+\frac{p^{3}+2p^{2}-p-1}{8(2p^{4}-p^{3}-2p^{2}+p+1)}$
$\displaystyle=p^{4}+p^{5}+3p^{6}+8p^{7}+{\bf
20}p^{8}+58p^{9}+152p^{10}+427p^{11}+1155p^{12}+3211p^{13}+O(p^{14}).$

Note that, for instance, the coefficient of $p^{8}$ is $20$, indicating that
there are $20$ polyiamonds in $\mathcal{CP}^{(1)}$ of perimeter $8$ as seen in
Figure 9.

Figure 9. Convex polyiamonds in $\mathcal{CP}^{(1)}$ with perimeter 8

By Lemma 4.4, (4.8) and (4.13), we obtain the first part of Theorem 1.3.
Finally, when $p$ is near the dominant singularity of the generating function
$F^{(1)}(p,1)$ (see Lemma 4.4), $F^{(1)}(p,1)$ can be approximated by

$\displaystyle\frac{2114724967}{3292828497}-\frac{80}{15219}\sqrt{5}-\frac{2560}{441\sqrt{3}}\sqrt{1-3p}+O((1-3p)),$

which, by the singularity analysis, gives the second part of Theorem 1.3.

Acknowledgements: The authors would like to thank the referee for pointing out
the mistake in the statement of Theorem 1.2 in the original draft.

## References

  * [1] C. Banderier, M. Bousquet-Mélou, A. Denise, P. Flajolet, D. Gardy, and D. Gouyou-Beauchamps, Generating functions for generating trees, Discrete Math., 246:1-3 (2000), 29–55.
  * [2] E. Barcucci, F. Bertoli, A. Del Lungo, and R. Pinzani, The Average Height of Directed Column-Convex Polyominoes Having Square, Hexagonal and Triangular Cells, Mathl. Comput. Modelling, 26:8-10 (1997), 27–36.
  * [3] G. Barequet and M. Shalah, Improved upper bounds on the growth constants of polyominoes and polycubes, Available at https://arxiv.org/abs/1906.11447.
  * [4] G. Barequet and M. Shalah, $\lambda>4$: an improved lower bound on the growth constant of polyominoes, Communications of the ACM, 59:7 (2016), 88–95.
  * [5] G. Barequet, M. Shalah, and Y. Zheng, An improved lower bound on the growth constant of polyiamonds, In Proceedings of the 22nd International Computing and Combinatorics Conference, (2017), 50–61. (Appeared also in Journal of Combinatorial Optimization).
  * [6] E. Barcucci, A. Frosini, and S. Rinaldi, Direct-convex polyominoes: ECO method and bijective results, in: R. Brak, O. Foda, C. Greenhill, T. Guttmann, A. Owczarek (Eds.). Proc. of Formal Power Series and Algebraic Combinatorics 2002, Melbourne, 2002.
  * [7] E. A. Bender, Convex $n$-ominoes, Discrete Math., 8 (1974), 219–226.
  * [8] M. Bousquet-Mélou, A method for the enumeration of various classes of column-convex polygons, Discrete Math., 154 (1996), 1–25.
  * [9] S. R. Broadbent and J. M. Hammersley, Percolation processes: I. crystals and mazes, Math. Proc. of the Cambridge Philosophical Society, 53:3 (1957), 629–641.
  * [10] K. Buchin, M.-K. Chiu, S. Felsner, G. Rote, and A. Schulz, The number of convex polyominoes with given height and width, Available at https://arxiv.org/abs/1903.01095
  * [11] A. R. Conway and A.J. Guttmann, On two-dimensional percolation, J. of Phys. A: Mathematical and General, 28:4 (1995), 891.
  * [12] M. P. Delest and S. Dulucq, Enumeration of directed column–convex animals with given perimeter and area, Croatica Chem. Acta, 66 (1993), 59–80.
  * [13] M. P. Delest and G. Viennot, Algebraic languages and polyominoes enumeration, Theoretical Computer Science, 34 (1984), 169–206.
  * [14] M. Eden, A two-dimensional growth process, Dynamics of fractal surfaces, 4 (1961), 223–239.
  * [15] S. Feretić, A perimeter enumeration of column-convex polyominoes, Discrete Mathematics and Theoretical Computer Science, 7 (1997), 57–84.
  * [16] S. Feretić, A perimeter enumeration of column-convex polyominoes, Discrete Mathematics and Theoretical Computer Science, 9 (2007) 57–84.
  * [17] P. Flajolet and R. Sedgewick, Analytic Combinatorics, Cambridge University Press, 2009.
  * [18] D. Gouyou-Beauchamps and P. Leroux, Enumeration of symmetry classes of convex polyominoes on the honeycomb lattice, Theoretical Computer Science, 346:2–3 (2005) 307–334.
  * [19] A.J. Guttmann, Polygons, polyominoes and polycubes, Lecture Notes in Physics, Springer, 775 (2009).
  * [20] F. Harary, Unsolved problems in the enumeration of graphs, Publications of the Mathematical Institute of the Hungarian Academy of Sciences, 5 (1960), 1–20.
  * [21] D. Kim, The number of convex polyominos with given perimeter, Discrete Math., 70:1 (1988), 47–51.
  * [22] D. A. Klarner, Cell growth problems, Canad. J. Math., 19 (1967), 851–863.
  * [23] K. Y. Lin and S. J. Chang, Rigorous results for the number of convex polygons on the square and honeycomb lattices, J. Phys. A: Math. Gen., 21 (1988), 2635–2642.
  * [24] G. Malen and É. Roldán, Polyiamonds attaining extremal topological properties, Part I, available at https://arxiv.org/abs/1906.08447
  * [25] T. Mansour and R. Rastegar, On enumeration of hexobars, work in progress.
  * [26] P. J. Peard and D. S. Gaunt, 1/d-expansions for the free energy of lattice animal models of a self-interacting branched polymer, J. Phys. A: Math. Gen., 28:21 (1995), 6109.
  * [27] J.-G. Penaud, Animaux dirigés diagonalement convexes et arbres ternaires, Technical Report 90-62, LaBRI, Université Bordeaux 1, 1990.
  * [28] R. C. Read, Contributions to the cell growth problem, Canad. J. Math., 14 (1962), 1–20.
  * [29] M. Shalah, Formulae and growth rates of animals on cubical and triangular lattices, PhD Dissertation, Technion, 2017.
  * [30] H. N. V. Temperley, Combinatorial problems suggested by the statistical mechanics of domains and of rubber-like molecules, Phys. Rev., 103 (1956), 1–16.
  * [31] X.G. Viennot, A survey of polyominoes enumeration, 4th FPSAC Proc., Publications du LACIM, 11 (1992), 399–420.
  * [32] M. Vöge and A. J. Guttmann, On the number of hexagonal polyominoes, Theoretical Computer Science, 307 (2003), 433–453.
  * [33] W. C. Yang and R. R. Meyer, Maximal and minimal polyiamonds, Technical Report, University of Wisconsin-Madison, 2002.


L PREPRINT - To Appear at ASP-DAC 20234.5cm,1cm1.001.8

# Hardware Security Primitives using Passive RRAM Crossbar Array: Novel TRNG
and PUF Designs

Simranjeet Singh1§, Furqan Zahoor 2§, Gokulnath Rajendran2§, Sachin Patkar1,
Anupam Chattopadhyay2, Farhad Merchant3 {simranjeet<EMAIL_ADDRESS>{anupam@, furqan.zahoor<EMAIL_ADDRESS><EMAIL_ADDRESS>aachen.de  1Indian Institute of Technology, Bombay, 2Nanyang Technological
University, Singapore, 3RWTH Aachen University, Germany

###### Abstract

With rapid advancements in electronic gadgets, the security and privacy
aspects of these devices are significant. For the design of secure systems,
physical unclonable function (PUF) and true random number generator (TRNG) are
critical hardware security primitives for security applications. This paper
proposes novel implementations of PUF and TRNGs on the RRAM crossbar
structure. Firstly, two techniques to implement the TRNG in the RRAM crossbar
are presented based on write-back and 50% switching probability pulse. The
randomness of the proposed TRNGs is evaluated using the NIST test suite. Next,
an architecture to implement the PUF in the RRAM crossbar is presented. The
initial entropy source for the PUF is used from TRNGs, and challenge-response
pairs (CRPs) are collected. The proposed PUF exploits the device variations
and sneak-path current to produce unique CRPs. We demonstrate, through
extensive experiments, reliability of 100%, uniqueness of 47.78%, uniformity
of 49.79%, and bit-aliasing of 48.57% without any post-processing techniques.
Finally, the design is compared with the literature to evaluate its
implementation efficiency, which is clearly found to be superior to the state-
of-the-art.

###### Index Terms:

PUF, TRNG, RRAM, Memristors, Hardware Security

§§footnotetext: Equal contribution

##  I Introduction

Physical unclonable functions (PUFs) are being considered the most reliable
hardware security primitives especially in IoT era for device authentication
and key generation [1]. PUF usually exploits the variations in the
manufacturing process and physical randomness inherent to the device, making
them highly reliable and unclonable [2]. PUF utilizes the idea that a device-
specific unique key can be derived each time rather than storing a set of keys
in the memory. PUF works on the principle of computing a unique output
response for the applied input challenge, thus creating a challenge-response
pair (CRP), which is crucial for designing security protocols, ranging from
device attestation to data encryption.

Various CMOS-based PUF desingns such as Arbiter PUF [3, 4] and SRAM PUF [5, 6]
are employed, which depend on the randomness introduced during their
fabrication process for their operation. However, they suffer from various
demerits, as they are area inefficient [7], sensitive to temperature
variations [8] and require data post-processing [9]. Also, carbon nanotube-
based PUF have been proposed, which utilize the random nature of the bits
generated from the unpredictable connections in the CNT devices [10]. This
approach showed some reliability. However, it is not usually employed due to
the complexity involved.

Figure 1: RRAM based PUF design

In recent years, resistive random access memory (RRAM) has been investigated
widely for use as the main memory and in-memory computing. A major drawback to
RRAMs being widely adopted as next-generation memories is their stochastic
nature and the intrinsic variation in switching parameters [11]. This
uncertain behavior of RRAM is critical for implementing hardware security
primitives, including PUF design [12]. The schematic structure of RRAM
consists of a metal-insulator-metal stack as depicted in Fig. 1b. The RRAM can
be considered a variable resistance, where memory states are stored as a high-
resistance state (HRS) and a low-resistance state (LRS). The RRAM state is
switched by applying a suitable voltage pulse. The devices significantly
impact the resistance state of the device due to the device-to-device (D2D)
and cycle-to-cycle (C2C) variations [13].

This stochasticity of RRAM is desirable for implementing PUF and other
hardware security primitives. It utilizes variations during the manufacturing
and fabrication process for authentication and secret key storage purposes
[14]. Lately, multiple research efforts have been focused for utilizing
emerging memories such as RRAM [15] and STT-MRAM (Spin Transfer Torque
Magnetic RAM) [16] for various security applications. RRAM has emerged as the
most feasible option for designs of PUFs since it allows extraction of
intrinsic hardware secrets from process variations, similar to that of the
traditional silicon-based PUFs. For RRAM, the key is dependent only on the CF
variations in the device during LRS and HRS which are unique for each
switching cycle. The secret is independent of the configuration that is
selected by a challenge [17]. These primitives are usually employed for
generating unique keys that connect inseparably to the hardware. A typical
RRAM-based PUF design is depicted in Fig. 1. A typical PUF schematic based on
RRAM crossbars is depicted along with the crossbar structure for collecting
the CRP.

For RRAM, the key is dependent only on the CF variations in the device during
LRS and HRS which are unique for each switching cycle.

In addition to PUFs, TRNGs are another important element for hardware
encryption modules. TRNGs typically exploit the randomness in physical
processes to generate a stream of random numbers [18]. Although various CMOS-
based designs have been implemented, they have limitations as they only
provide limited security-specific properties. The emerging technologies,
especially RRAM, have shown great potential due to high density, low power
operation, and new randomness sources [19]. The major contributions of this
work are as follows:

  * •

Two novel TRNG constructions using RRAM crossbar are proposed. The randomness
of designed TRNGs is tested using NIST tests.

  * •

Investigations on the impact of D2D and C2C variations on the crossbar and how
these variations can be exploited to design the PUF are studied. As passive
crossbar is used to design the PUF, we exploited sneak-path current as well to
construct the PUF.

  * •

Lastly, the performance of the proposed 16-bit PUF is evaluated and compared
with other RRAM-based PUF designs.

The remainder of the paper is organized as follows: Section II details the
methods employed to realize the TRNG and PUF designs based on RRAM. Section
III presents the results of the RRAM-based security primitives and discusses
some critical performance metrics such as uniformity, bit aliasing, and
reliability. Section IV concludes the paper.

##  II Methodology

Firstly, we demonstrate and validate a novel TRNG based on RRAM that relies on
the variations in the material stack of a device in the resistive switching
layer. We utilize the inherent stochasticity of the RRAM as the randomness
source to achieve the TRNG design. In this work, we propose two variants of
TRNG based on RRAM as design elements, as depicted in Fig. 2. The first design
utilizes the write-back technique for generating random numbers in the
crossbar, whereas the second design produces random numbers by applying a 50%
switching probability pulse [20]. The proposed RRAM-based TRNG designs have
lower circuit complexity advantages as randomness is generated and harvested
using simple elements. The proposed TRNG designs can be incorporated into
memory subsystems, enhancing security and area efficiency. The random bit
streams generated using both designs are depicted in Fig. 2a and Fig. 2b,
respectively.

Fig. 2a shows the flow to implement the TRNG using the write-back technique.
Initially, all the devices are formed into the LRS, followed by resetting all
devices back to HRS. Next, a pair of devices is selected, and a READ voltage
(0.2V) is applied across the rows of selected devices. As the variations on
the crossbar are normally distributed, one device will provide more current
than another. The current following through these devices are compared, and
the fast device in the selected pair is explicitly switched to LRS. This will
ensure an equal probability of 0’s and 1’s in each row of the crossbar, and
similar steps keep sliding over the crossbar till the last device.

Figure 2: Algorithm to implement TRNG in passive RRAM crossbar using the
write-back technique in (a) and 50% probability pulse in (b).

The second technique characterizes the switching behavior of devices appearing
due to D2D and C2C variations. Fig. 2b shows the flow of this technique.
Experimentally, a 50% switching probability pulse is analyzed and applied to
the crossbar, which results in 50% switching in the crossbar from HRS to LRS.
Initially, all the devices are in HRS, and Due to the Gaussian distribution of
variation on a crossbar, each device will have a different switching
threshold. When a 50% switching probability pulse is applied, random devices
in the crossbar switch their state (HRS to LRS), and other remains in the
previous state. The switching of devices is a direct function of variation in
the crossbar. The bitmap generated using this technique is shown in Fig. 2b.
Additionally, to evaluate the randomness of the generated bit streams, the
bits generated by the RRAM-based TRNG are passed through the NIST test suite,
and the results confirm the random nature of the obtained bit streams.
Further, the design TRNGs are used to construct the hardware security
primitives.

Next, we propose an architecture to implement the PUF on the passive RRAM
crossbar array depicted in Fig. 3. The crossbar is first initialized to an
entropy source based on the earlier TRNGs algorithm. The resposes of the PUF
are collected by applying challenges to the crossbar’s rows. These challenges
are mapped to a read voltage pulse (1 $\rightarrow$ 0.2V and 0$\rightarrow$
0V) before applying it to the crossbar. Based on the input challenge and
device variations, Kirchoff’s current law collects the current following
through each device at the column lines. In addition, the sneak-path current
in the crossbar contributes to the current at each column, which is entirely
random. The current sense amplifier (CSA) converts the analog current values
to boolean response bits at the output. Further, many CRPs are collected to
check the different PUF properties. The proposed architecture shows an
implementation of TRNG and PUF on a single RRAM crossbar structure.

Figure 3: Passive RRAM-based PUF used for challenge-response generation and
sneak-path current is marked in red. Where ’n’ and ’m’ represent the number of
rows and columns in the crossbar, respectively.

##  III Experimental Results

This section evaluates the different properties of the proposed PUF on the
RRAM crossbar, such as uniformity, uniqueness, bit-aliasing, and reliability.
First, we collect the CRPs from the proposed design (refer to Fig. 1), and all
mentioned properties are evaluated on this data. We also examine the
implementations of different TRNGs in the RRAM crossbar and their validation
using the NIST test. Later, we study the impact of D2D and C2C variations on
the crossbar and exploitation of the same to construct the TRNG and PUF.

###  III-A RRAM cell Characterization

The RRAM cell used for this study consists of a Pt/Ti/TiOx/HfO2/Pt material
stack, which has been proven to have higher stability in terms of
electroforming voltage stability and thermal stability [21, 22]. Fig. 1b shows
the equivalent device stack of one RRAM cell. The device has two states HRS
$\rightarrow$ logic ‘0’ and LRS $\rightarrow$ logic ‘1’. Fig. 1c shows the I-V
characteristics of a device where SET (Programming to LRS) and RESET
(Programming to HRS) paths are marked with arrows. A voltage pulse with
specific duration and amplitude is applied to program the device into HRS or
LRS. The device provides resistance between $60-100K\Omega$ in HRS and
$1.5-1.6K\Omega$ in LRS. The devices in the crossbar are used without any
selector (Passive) in series. Generally, the passive crossbar faces the
problem of sneak-path current, which can be exploited to design the TRNGs and
PUF. Next, we examined the impact of process variations on RRAM cells and how
these variations can be exploited to design novel hardware security
primitives.

###  III-B Variation Analysis

In this section, we study the impact of manufacturing process variations on
RRAM and how we can exploit these variations to design TRNGs and PUF. There
are majorly two types of variations known as D2D and C2C which can affect the
switching behavior of RRAM devices. The following are the parameters affected
during the manufacturing process of RRAM used in this study:

  * •

Oxygen vacancy concentration in the disc, maximum ($N_{Disc,max}$) and minimum
($N_{Disc,min}$).

  * •

Radius of the disc ($r_{disc}$).

  * •

Length of the disc ($l_{disc}$).

Table I gives the ranges of variable parameters those are used for this study.
These variations are experimentally verified and that’s provides a validation
to our design choices [21].

#### C2C variations

C2C variations in the devices are analyzed by changing the variable parameters
before every SET and RESET operation. The change of parameters in the next
cycle is based on the previous state parameters and random walk of parameters
in +ve and -ve directions with an equal probability. The variability
parameters are exported in a CSV file using MATLAB scripts. The final
parameter file contains variation values and the simulation time stamp at
which these parameters will be adopted. This file is then loaded into the
_Cadence Spectre_ simulation as a parameter list.

Fig. 4a shows the distribution of HRS and LRS for 1000 cycles. Each cycle
represents a voltage pulse with a 20ns pulse duration, 2.5ns rise time, and
2.5ns fall time. The new set of values are applied after the 100ns period.
Within these 100ns, the device is RESET with 25ns pulse of -2V and SET with 2V
pulse with the same duration. The state of the device is read by applying a
READ pulse of 0.2V after every RESET/SET operation. In the next phase (after
100ns), new values are fetched from the parameter file and applied in the
simulation. It has been observed that there is a $\pm$5% change in HRS and a
$\approx$1% change in LRS during C2C variations.

TABLE I: Parameter ranges for variation study RSD = 0.5
---
- | low / mean / up
$N_{min,var}$ | 4 / 8 / 16
$N_{max,var}$ | 18 / 20 / 22
$r_{var}$ | 40.5 / 45 / 49.5
$l_{var}$ | 0.36 / 0.4 /0.44

#### D2D variations

To simulate D2D variations, the random set of values for each device’s variate
parameters are drawn from the experimentally verified Gaussian distribution
[21]. These variations were then independently applied to the available
devices in the crossbar. Fig. 4b shows the distribution of HRS and LRS states
in a 10x10 crossbar (100 devices) without a transistor. To measure the HRS
distributions, all the devices are first switched to LRS by applying a 2V SET
pulse with 150ns pulse width along with 50ns rise and fall time. In the next
cycle, the devices are switched to the HRS state by applying RESET voltage of
-2V with a similar pulse duration used during SET operation. The resistance
state of the device is read at 0.2 V. The same procedure is used to find LRS
distribution. It has been observed that HRS has more impact of process
variations compared to LRS. Due to variations, the HRSvaries from 31$K\Omega$
to 155$K\Omega$ with an average of 65.56$K\Omega$ and LRS from 1.55$K\Omega$
\- 1.67$K\Omega$ with an average of 1.64$K\Omega$. However, HRS has a wide
distribution; all HRSs are distinguishable from LRS. A carefully chosen pulse
to RESET the devices will lead to switching the devices randomly from LRS to
HRS or vice-versa, which can eventually be exploited to design the TRNGs.

Figure 4: LRS and HRS distribution for 1000 cycles under C2C variations in (a)
D2D variation on 10x10 crossbar in (b).

#### Pulse width study

This study gives the pulse width required to switch a device from LRS to HRS
and vice-versa. We analyze the switching characteristic of devices on
different pulse duration with a fixed voltage. Table II shows the switching
probability of crossbar with different pulse duration. To conduct this
experiment, all the devices are considered initially into HRS, and a fixed
voltage of -0.8V with 1ns rise & fall time for pulse duration less than 100ns
and 10ns for other pulse duration is applied. A total of 100 devices in a
10x10 crossbar are connected to an individual voltage source. A READ voltage
of 0.2V is applied to read the state of the devices. The crossbar contains D2D
and C2C variations as discussed in Section III-B.

Fig. 5 shows the switching behavior of devices with different pulse widths.
The devices’ states are read at 0.2V and sorted according to their resistance
value. Initially, all the devices were in HRS. As the pulse width increases,
more devices start switching into LRS from the HRS state. This analysis gives
a distribution of HRS and LRS, which can further be used to choose the 50%
probability pulse. Finding a pulse width that can randomly switch 50% of
devices in the crossbar is exploited to design the TRNG. A 3us pulse is used
for these particular devices to switch the 50% of devices in the crossbar, as
given in Table II.

Figure 5: Switching behavior and distribution of LRS and HRS of 100 devices in crossbar under different pulse frequencies. TABLE II: Pulse width study to find out 50% pulse to generate random numbers in RRAM crossbar. SET voltage = -0.8V, READ voltage = 0.2V, Rise and Fall time = 1ns. Pulse duration | SET probability | LRS ($K\Omega$) | HRS ($K\Omega$)
---|---|---|---
10 ns | 0% | 1.67 - 1.7 | 23-140.7
100 ns | 11% | 1.62-1.67 | 38-140.8
1 us | 36% | 1.6-1.675 | 32.2-140.8
3 us | 50% | 1.58 - 1.675 | 47- 140.8
10 us | 64% | 1.585-1.67 | 47-140
1 ms | 94% | 1.58-1.67 | 94-128
10 ms | 100% | 1.58-1.68 | -

###  III-C PUF properties

For evaluating a PUF design, a number of properties need to be analyzed, such
as uniformity, uniqueness, bit-aliasing, and reliability. In the following
subsections, we provide the analysis of the proposed design on these
properties in detail.

The first step is calculating the intra Hamming distance (intra-HD), which
shows the difference in the output response bits corresponding to different
input challenges. Fig. 6(a) shows the histogram for intra-HD in responses of
16-bit RRAM-based PUF, which is calculated by the Equation 1. Ideally, PUF
should result in a 50% intra-HD in the response bits. In this case, with a
particular entropy and variations, the proposed PUF shows 39.68% intra-HD.

$\mbox{\rm Intra-HD}=\frac{1}{n}\sum_{i=1}^{n}HD(R_{i},R_{i+1})\times 100$ (1)

Where ‘n’ is the total number of challenges given to the PUF, $R_{i}$ and
$R_{i+1}$ are the responses to the challenges $C_{i}$ and $C_{i+1}$
respectively. Note that $R_{i}$ is a bit vector structured from each response
bit position at the columns from $r_{i}$ to $r_{m}$ for $i_{th}$ challenge.

(a) Intra-HD

(b) Uniformity

(c) Uniqueness

Figure 6: The different properties of 16-bit RRAM-based PUF

#### Uniformity

The uniformity property shows how the response bit from each device is
distributed in 0’s and 1’s. It is expected from the PUF for an equal number of
switching between 0 and 1 for each response bit. The ideal value for
uniformity is 50%, calculated by taking the average of all the response bits’
distribution. Fig. 6(b) shows the distribution of 0’s (green) and 1’s (red) in
each response bit position. Even though some response bits are more inclined
towards 0’s, and some are towards 1’s, the proposed PUF shows 49.79%
uniformity. The distributions of 0’s and 1’s are a function of the entropy
source and will change its inclination on another entropy set for other
crossbars.

#### Uniqueness

Uniqueness characterizes the ability to distinguish between different devices
based on their response to the same challenge. As the applied challenges are
identical, the difference between the responses is entirely based on the
process variations. It is expected that the PUF should respond with 50% HD in
the responses to the same applied challenges for different devices. The
uniqueness experiment is conducted using different crossbar variations
randomly drawn from the Gaussian distribution. A total of 50K CRPs are
collected from each device. Further, Fig. 6(c) shows the histogram of
uniqueness, which is calculated as per Equation 2. The proposed PUF gives
47.78% uniqueness, which is very close to the ideal value.

$\mbox{\rm Uniqueness
}=\frac{2}{k(k-1)}\sum_{i=1}^{k-1}\sum_{j=i+1}^{k}\frac{HD(R_{i},R_{j})}{n}\times
100$ (2)

Where, $i$ and $j$ are two different devices. $R_{i}$ and $R_{j}$ are the
responses from $device_{i}$ and $device_{j}$ for the challenge C respectively.
$k$ is the number of devices and $n$ is total number of responses.

#### Bit-aliasing

Bit-aliasing investigates the switching of each response bit across the number
of available crossbars and can be calculated as per Equation 3. It is
essential to ensure that no bit positions of responses are strongly biased
towards 0 or 1. To calculate bit-aliasing for the proposed design, CRPs are
collected from two different crossbars, and then each bit position of response
is analyzed. Ideally, bit-aliasing should be 50% for a well-balanced design,
and the proposed design shows a bit-aliasing of 48.57%.

$\mbox{\rm Bit-aliasing}=\frac{1}{n}\sum_{i=1}^{n}(R_{i,p}==1)\times 100$ (3)

where $R_{i,p}$ is the response bit at the $pth$ cell location in the $ith$
chip and ‘n’ is the total number of responses including response from all
devices.

#### Reliability

Reliability captures how efficient a PUF is in repeating the same response for
the same challenge. The higher the reliability, the less cost for the error-
correcting techniques. Reliability is calculated by intra-HD between responses
during the same applied challenges. Ideally, intra-HD for the same response
will be zero, and reliability is given as $(100\%-\mbox{\rm(intra-HD)})$. To
check the reliability of the proposed design, we applied the same challenge of
10K times and collected the responses. The proposed design shows a reliability
of 100% without any error-correcting techniques.

###  III-D Machine Learning-based Modelling Attacks

The machine learning (ML) models can be trained to predict the responses of a
given PUF. It has been shown that nearly all variants of CMOS-based PUFs can
be efficiently modelled using ML [23] and correspondingly attack resistant PUF
structures are also being proposed [24]. These attacks characterize the
silicon behavior from CRPs. However, in RRAM-based PUF, the underlined devices
have different characteristics and physics compared to the CMOS-based design,
which requires a deeper analysis. We reserve the study of detailed analysis of
modeling attacks for the future.

TABLE III: Comparison of the performance metrics of the proposed PUF Parameters |  Proposed PUF |  Nano PUF [25] |  RRAM PUF [26] |  Memristor PUF [27] |  MR-PUF [28]
---|---|---|---|---|---
Crossbar | 16x16 | 1 cell | 64x8 | 2x10x10 (3D) | 1 cell
Uniqueness | 47.78 | 47 | 49.85 | 50 | 49.3
Uniformity | 49.79 | 47 | 47.28 | 49.5 | 49.9
Reliability | 100 | 90 | 98.67 | 97 | 100
Bit-aliasing | 48.57 | - | 49.86 |  $\approx$ 50 | 49.6

###  III-E Randomness tests and Comparison

The process variations in RRAM devices results in different switching behavior
of these devices and are exploited to design the TRNG. We proposed two methods
(50% probability pulse and write-back) to generate TRNGs in the 16x16 RRAM
crossbar, which are further used to construct the PUF. The national institute
of standards and technology (NIST) randomness tests are conducted to analyze
the randomness. The method with a 50% probability pulse passes 11 tests, and
the write-back technique passes 08 out of 15 tests. NIST tests are conducted
on just 16x16 crossbar data, and it passes all critical tests. The NIST tests
can further be improved by increasing the size of the crossbar.

Table III compares the proposed PUF with other RRAM-based PUFs. Table III
shows that the proposed design on a 16x16 crossbar outperforms other RRAM-
based designs, even without the post-processing. MR-PUF [28] analyzes the
single-cell behavior (not in crossbar) and presents the result after post-
processing. Our technique is implemented on the passive crossbar, which
exploits the process variations along with the sneak-path current.
Construction using a passive crossbar leads to less overall area.

##  IV Conclusions

This work presents novel TRNG and PUF designs based on RRAM technology in the
passive crossbar. The inherent variability of the switching mechanism of RRAMs
is exploited to generate random bits. The TRNG output was evaluated using the
NIST test suite, confirming the generated bits’ random nature. Next, the
exploitation of device variations and sneak-path currents in the passive RRAM
crossbar is shown to implement the PUF. The results are evaluated with
extensive SPICE simulations, which verify the robust nature of the proposed
PUF as 100% reliability is achieved without the need for post-processing. The
proposed PUF design also demonstrates good performance for various metrics
such as uniformity, uniqueness, and bit-aliasing. The results are encouraging
and demonstrate the potential for using RRAM-based designs in hardware
security applications. We will explore the prototyping of the proposed designs
and a detailed study of resilience against attacks in the future.

## Acknowledgment

This work was partially funded by Deutsche Forschungsgemeinschaft (DFG –
German Research Foundation) under the priority programme SPP 2253 and was
partially funded by the grant NRF-CRP21-2018-0003.

## References

  * [1] R. Liu _et al._ , “A highly reliable and tamper-resistant rram puf: Design and experimental validation,” in _IEEE HOST_ , 2016, pp. 13–18.
  * [2] R. Zhang _et al._ , “Nanoscale diffusive memristor crossbars as physical unclonable functions,” _Nanoscale_ , vol. 10, no. 6, pp. 2721–2726, 2018\.
  * [3] D. Ganta _et al._ , “Easy-to-build arbiter physical unclonable function with enhanced challenge/response set,” in _IEEE International Symposium on Quality Electronic Design (ISQED)_ , 2013, pp. 733–738.
  * [4] D. Lim _et al._ , “Extracting secret keys from integrated circuits,” _IEEE Transactions on Very Large Scale Integration (VLSI) Systems_ , vol. 13, no. 10, pp. 1200–1205, 2005.
  * [5] D. E. Holcomb _et al._ , “Power-up sram state as an identifying fingerprint and source of true random numbers,” _IEEE Transactions on Computers_ , vol. 58, no. 9, pp. 1198–1210, 2008.
  * [6] J. Guajardo _et al._ , “Fpga intrinsic pufs and their use for ip protection,” in _International workshop on cryptographic hardware and embedded systems_. Springer, 2007, pp. 63–80.
  * [7] F. Armknecht _et al._ , “Memory leakage-resilient encryption based on physically unclonable functions,” in _Towards Hardware-Intrinsic Security_. Springer, 2010, pp. 135–164.
  * [8] C. Herder _et al._ , “Physical unclonable functions and applications: A tutorial,” _Proceedings of the IEEE_ , vol. 102, no. 8, pp. 1126–1141, 2014\.
  * [9] J. Rajendran _et al._ , “Nano-ppuf: A memristor-based security primitive,” in _IEEE Computer Society Annual Symposium on VLSI_ , 2012, pp. 84–87.
  * [10] B. Srinivasu _et al._ , “Cycle puf: A cycle operator based puf in carbon nanotube fet technology,” in _IEEE 21st International Conference on Nanotechnology (NANO)_ , 2021, pp. 13–16.
  * [11] F. Zahoor _et al._ , “Resistive random access memory (rram): an overview of materials, switching mechanism, performance, multilevel cell (mlc) storage, modeling, and applications,” _Nanoscale research letters_ , vol. 15, no. 1, pp. 1–26, 2020.
  * [12] Y. Pang _et al._ , “Optimization of rram-based physical unclonable function with a novel differential read-out method,” _IEEE Electron Device Letters_ , vol. 38, no. 2, pp. 168–171, 2017.
  * [13] M. Gonzalez _et al._ , “Analysis of the switching variability in ni/hfo2-based rram devices,” _IEEE Transactions on Device and Materials Reliability_ , vol. 14, no. 2, pp. 769–771, 2014.
  * [14] S. Balatti _et al._ , “Physical unbiased generation of random numbers with coupled resistive switching devices,” _IEEE Transactions on Electron Devices_ , vol. 63, no. 5, pp. 2029–2035, 2016.
  * [15] B. Gao _et al._ , “A unified puf and trng design based on 40-nm rram with high entropy and robustness for iot security,” _IEEE Transactions on Electron Devices_ , vol. 69, no. 2, pp. 536–542, 2022.
  * [16] L. Zhang _et al._ , “Highly reliable memory-based physical unclonable function using spin-transfer torque mram,” in _IEEE International Symposium on Circuits and Systems (ISCAS)_ , 2014, pp. 2169–2172.
  * [17] Y. Gao _et al._ , “Physical unclonable functions,” _Nature Electronics_ , vol. 3, no. 2, pp. 81–91, 2020.
  * [18] R. Govindaraj _et al._ , “Csro-based reconfigurable true random number generator using rram,” _IEEE Transactions on Very Large Scale Integration (VLSI) Systems_ , vol. 26, no. 12, pp. 2661–2670, 2018.
  * [19] B. Yang _et al._ , “Rram random number generator based on train of pulses,” _Electronics_ , vol. 10, no. 15, p. 1831, 2021.
  * [20] R. John _et al._ , “Halide perovskite memristors as flexible and reconfigurable physical unclonable functions,” _Nature Communications_ , vol. 12, 06 2021.
  * [21] C. Bengel _et al._ , “Variability-aware modeling of filamentary oxide-based bipolar resistive switching cells using SPICE level compact models,” _IEEE TSCAS I_ , vol. 67, no. 12, pp. 4618–4630, 2020.
  * [22] A. Hardtdegen _et al._ , “Improved switching stability and the effect of an internal series resistor in HfO2/TiOx bilayer ReRAM cells,” _IEEE Trans on Electron Devices_ , vol. 65, no. 8, pp. 3229–3236, 2018.
  * [23] U. Rührmair _et al._ , “Puf modeling attacks on simulated and silicon data,” _IEEE Transactions on Information Forensics and Security_ , vol. 8, no. 11, pp. 1876–1891, 2013.
  * [24] N. Shah _et al._ , “Introducing recurrence in strong pufs for enhanced machine learning attack resistance,” _IEEE Journal on Emerging and Selected Topics in Circuits and Systems_ , vol. 11, no. 2, pp. 319–332, 2021.
  * [25] J. Rajendran _et al._ , “Nano meets security: Exploring nanoelectronic devices for security applications,” _Proceedings of the IEEE_ , vol. 103, no. 5, pp. 829–849, 2015.
  * [26] J. Kim _et al._ , “A physical unclonable function with redox-based nanoionic resistive memory,” _IEEE Transactions on Information Forensics and Security_ , vol. 13, no. 2, pp. 437–448, 2018.
  * [27] H. Nili _et al._ , “Hardware-intrinsic security primitives enabled by analogue state and nonlinear conductance variations in integrated memristors,” _Nature Electronics_ , vol. 1, 03 2018.
  * [28] H. Ibrahim _et al._ , “Memristor-based puf for lightweight cryptographic randomness,” _Scientific Reports_ , vol. 12, 05 2022.