UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Modelling the swarming behaviour of army ants Watmough, James 1992

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata


831-ubc_1992_fall_watmough_james.pdf [ 1.07MB ]
JSON: 831-1.0079673.json
JSON-LD: 831-1.0079673-ld.json
RDF/XML (Pretty): 831-1.0079673-rdf.xml
RDF/JSON: 831-1.0079673-rdf.json
Turtle: 831-1.0079673-turtle.txt
N-Triples: 831-1.0079673-rdf-ntriples.txt
Original Record: 831-1.0079673-source.json
Full Text

Full Text

MODELLING THE SWARMING BEHAVIOUR OF ARMY ANTSByJames WatmoughB.A.Sc (Engineering Physics) University of British ColumbiaA THESIS SUBMITTED IN PARTIAL FULFILLMENT OFTHE REQUIREMENTS FOR THE DEGREE OFMASTER OF SCIENCEinTHE FACULTY OF GRADUATE STUDIESDEPARTMENT OF MATHEMATICSANDINSTITUTE OF APPLIED MATHEMATICSWe accept this thesis as conformingto the required standard£ uUMANovember 1992© James Watmough, 1992In presenting this thesis in partial fulfilment of the requirements for an advanceddegree at the University of British Columbia, I agree that the Library shall make itfreely available for reference and study. I further agree that permission for extensivecopying of this thesis for scholarly purposes may be granted by the head of mydepartment or by his or her representatives. It is understood that copying orpublication of this thesis for financial gain shall not be allowed without my writtenpermission.Department of____________The University of British ColumbiaVancouver, CanadaDate OQ-’r. C5DE-6 (2/88)abstractA model of the collective motion of a group of social insects interacting through a trail networkis developed around a set of simple rules. The model assumes that individuals interact withthe existing trail network both by following the trails, and by adding to and reenforcing trails.The main hypothesis of the model is that the mechanics of the motion of the individuals aresufficient to determine the collective motion of the group, or equivalently, that the parametersgoverning the behaviour of the individuals are sufficient to determine the patterns of the trailnetworks produced. Both a simulation and a mathematical model are derived from this sameset of rules. The mathematical model is derived for a single spatial dimension to study thepropagation of the raiding column. The results are compared both to the swarming behaviourof army ants observed in the field, and to the results of simulations of other authors whichwere based on different assumptions. We conclude that it is possible that many aspects of theswarming pattern are controlled by the parameters governing the behaviour of the individual.11Table of Contentsabstract iiList of Tables vList of Figures vi1 Introduction 11.1 The Swarm 11.2 Organization in Social Insects 51.3 Communication Through a Trail Network 51.4 Using Simulations to Model Trail Formation 71.5 Overview 8ft. 2 Mechanics of Motion 92 1 Rules of Motion 92.2 Results of the Computer Simulation 113 A Model in One Spatial Dimension 183.1 Derivation of the Model 183.2 Reduction to the Travelling Wave Form 233.3 Analysis of the State Space 283.3.1 Local Analysis of the Steady States 283.3.2 Global Analysis 303.4 Stability of the Travelling Wave as a Solution to the PDEs 373.5 Implications and Generalizations 37‘UDiscussion 43Bibliography 48A Dimensional Analysis 50B Elgenvalues in 3 Dimensions 52C The Location ofD Details of the Computer Simulation 58List of Notations 60ivList of TablesB.1 Real eigenvalues of a 3x3 matrix 54B.2 Imaginary eigenvalues of a 3x3 matrix 54VList of Figures1.1 Development of the swarm 31.2 Typical raiding patterns of army ants 41.3 Mechanism of trail detection 72.1 Comparison of simulation to field observations 132.2 Dependence of trail network on forward persistence 142.3 Dependence of trail network on the affinity of the ants to the trails 152.4 Dependence of trail network on the rate of trail deposition 162.5 The formation of a trunk trail 173.1 The State Space 323.2 Intersection of W with T° 343.3 A schematic showing the boundaries of A— and their intersections with W? . 353.4 A section showing the intersection of W and W? with a plane of constant T 363.5 Stable travelling wave solution: E .1, c = .5, and -y = 1 383.6 Stable wave resulting from oscillating input 393.7 Simulations of Deneubotirg et al 46C.1 Locating the eigenspace of u0 57viChapter 1Introduction1.1 The SwarmAspects of trail following and trail formation are utilized to varying degrees by many socialorganisms, from the repeated use of the same trail by an individual, to the formation of vasttrail networks by large numbers of insects. Many means of communication and interactionwithin the trail network are used. Individuals may use familiar landmarks as roadsigns; or theymay orient themselves with respect to the sun, moon, stars, or even the earth’s magnetic field.Able [1] gives many good examples. One of the most striking examples of trail formation arethe raiding columns formed by army ants. Here, up to 20 million’ individuals will march as asingle entity. Although these ant colonies display an extensive caste system, it is known [21]that the networks form without the use of any preset pattern, and without direction from anyant, or group of ants, which may be considered an authority. It is possible that the ants maytake some clues from their surroundings [22]; however, in this paper we study the formation ofa trail network independent of these effects.To model the collective motion in the raiding column we take the hypothesis that there isno preset or ‘inherited’ pattern followed by the ants, and further, that there are no ‘leader’ants to be given the task of deciding the direction of the day’s pillage. In fact, an amusingcharacteristic of the ‘front line’ of the raiding column is that any ant who discovers herself to beon unexplored territory will immediately ‘about face’ and attempt to reenter the swarm. Thecolumn itself is maintained not by any central intelligence on behalf of the upper caste, but bythe determination of individuals to push forward and follow their reluctant leaders [21].‘The number of ants in a. swarm is species dependent, this figure is for Dorylus (Anomrna) wilverthis. (seeSchneirla [29])1A typical scenario has the ants forming a bivouac at night, linking hundreds of thousands ofworkers around the queen and her brood. The bivouac begins to break apart towards morningas the light levels increase [21]. Pressure from this dissolving mass of ants then forces the outerants away from the nest site. Observations of Schneirla [29] on the species Eciton burchelli,a new world ant, suggest that these ants begin a radial expansion which grows fastest in thedirection offering the least resistant terrain. After about 30 minutes a trunk trail will havedeveloped as the ants attempt to remain together (see figure 1.1). Here we have the firstindication that the motion is not simply due to diffusion, since the regions in which travel isdifficult, or prey is scarce, are abandoned entirely. If the motion were simply diffusive, theseareas would be explored, albeit on a slower time scale. After an hour a typical swarm patternhas formed (figure 1.2), and a few of the ants begin to return along the newly formed trailsladen with pieces of their favourite prey.The fully developed E. burchelli swarm is composed of a loose network of trails 10 to 15meters wide, and about one meter deep. It will progress at a fairly steady rate of 10 to 12meters per hour. Behind the swarm, the trail network gradually converges to a single trailtrunk. The swarm front expands quickest into areas which are rich in booty, and those whichpresent level, unobstructed terrain. Thus, the swarm does not simply sweep forward, but willveer left and right into these preferred areas as it proceeds. As pressure from the unfoldingbivouac increases, the swarm front will sprout subswarms from the flanks. Some of these willgrow, and some will perish, depending on the terrain and the abundance of prey.For the species Eciton hcimaturn, the swarm front will be several orders of magnitude smaller,and the raid will propagate as a narrow column rather than the broad fan of Eciton burchelli [29].Ants moving along the trail trunk have been clocked at speeds varying erratically between 3and 5 cm/s while heading away from the nest. Those returning to the nest with food proceedat a steadier pace usually in the range of 4 to 5 cm/s. Note that this is an order of magnitudelarger than the aforementioned column speed (10 meters/hour 3 mm/s).2PsIFigure 1.1: The development of the swarm: (a) initial radial expansion as the ants begin toleave the bivouac; (b) the formation of a truck trail; (c) the fully developed swarm. Taken fromSchneirla [29].2343Ec/fon homo turn Column Raid Ec/fol? burche/// Swarm Raidraid front‘ ..:.cS F/11 (7’ booty cache(..rb ivouoc\bose column5melergFigure 1.2: Typical raiding patterns of army ants: (left) E. hamaturn during the nomadic phase;(right) E. burchelli. Taken from Rettenmeyer [28].swarmbooty cache5 meterscolumnbivouac41.2 Organization in Social InsectsCamazine [7] and Deneubourg et al. [10] have both presented arguments suggesting that organizational phenomena can be broken into two broad classes. As one extreme the society employsa central scrutinizer to oversee the entire organizational process, and at the other they relyon a large number of simple automata in a parallel arrangement. In the first case, commonlyreferred to as the Blueprint Hypothesis2,there is a requirement for a complex algorithm, anda higher level of information processing. Unfortunately, such systems are in general not robustto changes in the environment, and are known to be inflexible to changes in the needs of thesociety.The second case has long been known a.s Self-Organization [7, 8, 11]. Under the selforganization hypothesis there is no need to relay all information to a central processing unitfor decision making. This in turn leads to a system that is more robust to changes in the information supplied to the system. Since all the processing is done on the level of the individual,and in a parallel arrangement, the algorithms ca.n be far simpler. in addition the informationrequired by the system can be both localized and unrefined. If any mistakes occur in the algorithm of a single individual, this will have little effect on the operation of the entire system. Infact, if there is a. small amount of randomness in the behaviour of the individuals, such systemswill lead to an optimal solution even in a changing environment [10]. Several examples of thepossible use of this system by ants are given in [2, 5, 8, 9, 10, 11], and for honey bees in [7].1.3 Communication Through a Trail NetworkThe dominant form of communication and orientation along the trails is through the use ofchemicals. The chemicals used, and the complexity of the communication varies a great dealbetween species. No less tha.n ten anatomical structures have been noted as being used in20n a philosophical note, it is possible that the term ‘blueprint’ does not convey the correct meaning. Weare not begging the question of what gives rise to the pattern, we are distinguishing between two distinct levelsof intelligence possible in the system. In the blueprint hypothesis, decisions are made by a well defined centralintelligence; whereas iii the self-organizational hypothesis, the apparent intelligent behaviour of the system hasno tangible source, but arises out of the coordinated efforts of a large number of individuals.5depositing trail indicators. In addition, as many 14 chemicals can be used as indicators in asingle species [21]. Fortunately, in the large colonies typical of the army ants the means ofcommunication appear to be simpler. The species Eciton hamatum use only excretions fromthe hindgut [21]. Although many species of ants will deposit pheromone only as a periodic trailmarker, using several different combinations in accordance with some measure of worth of thebooty they have fouiid along a particular trail, army ants are known to lay pheromone in acontinuous track both when foraging and when returning to the nest [3, 12]. Although it hasalso been noted that many differences in the collective motion of a group can be attributedto its size [14], as well as the type of food foraged [12], we will concern ourselves herein onlywith simple communication through the use of chemical trail markers. In this manner weattempt to determine the dependence of the collective motion on the parameters describingthis communication.According to Calenbuhr [5], ants use a bi-receptor mechanism to determine their motionrelative to the pheromone trail. This allows them to detect not only the local concentrations ofpheromone, but to determine the local gradients in this concentration as well. As a result, theants move in a continual weaving fashion in and out of the trail (see figure 1.3). There is evidencethat the resolution degrades at higher concentrations due to saturation of the antennae [6]. Thismeans that two trails, each marked with a large level of pheromone, will appear to the ant tobe equal in concentration; whereas, if the two trails are marked with pheromone in the samerelative proportions but at lower concentrations, the ant will be able to distinguish betweenthem. The ability of the ants to detect and follow trails is not perfect; thus, they will followa given trail for a finite distance, and this distance will be chosen from a random distributionwith a certain mean and variance. It has been noted by several authors (see the discussionin Calenbuhr [5]) that the mean distance that an ant will follow a trail increases with theconcentration of pheromone on that trail.61’ilFigure 1.3: Adapted from Schneirla [29]. The ant will use it’s antenna to detect the stronger oftwo trails. As the ant attempts to follow a single trail it will weave in and out of the pheromonetrack.1.4 Using Simulations to Model Trail FormationWith problems involving simple interactions between many similar units, it is useful to studythe formation of the pattern using a simulation. In our case we are dealing with the interactionsof large numbers of ants through a chemical trail network. Although a mathematical modelof the motion has been derived in two dimensions, it is difficult to analyze except for specialcases. Using a simulation gives us the advantage of being able to explore the development of thenetwork in two dimensions without being forced to make specific assumptions. Furthermore,the results of the simulation are displayed in ‘real time’ on the computer screen. These visualimages provide a excellent means to both study the problem, and convey the results to a wideraudience. We may expand on the initial hypothesis of the model and view the results of variousassumptions on the formation of the trail network. This can sometimes be done more easily,7and often far faster by using a simulation than by using a more rigorous mathematical analysis.The insight provided by this experimentation can also be used to direct a further, more detailed,analysis of the problem. Such an analysis will still be necessary as it is not possible, using thesimulation, to place definite bounds on region of parameter space in which different resultsoccur.1.5 OverviewIn section 2.1 we construct a set of simple rules governing the response of individuals both toeach other, and to the trail network. In section 2.2, these rules are incorporated into a computersimulation.In section 3.1 the rules mentioned above are used to derive a system of partial differentialequations (P DEs). However, due to the complexity of the analysis of a model with two spatialdimensions, we limit our analysis to a model with only one spatial dimension. This modelrepresents the propagation of the leading edge of the trail network. In section 3.2, usingvarious techniques of applied analysis, we simplify the PDEs to a system of ordinary differentialequations (ODEs) whose solutions represent travelling wave solutions of the original PDEs. Insection 3.3, these ODEs are then studied analytically, and we are able to show how the speed ofpropagation of the raiding column depends on the parameters of the model. Section 3.5 dealswith generalizations on this model which provide us with a method for relating the motionresulting from our rules to similar motions caused by other mechanisms.8Chapter 2Mechanics of Motion2.1 Rules of MotionTo model the behaviour of the swarm using a system of automatons, we develop a set ofrules or algorithms to be followed by each individual. As noted in section 1.2, one advantageof self-organization is that it is flexible, and requires little processing on the part of eachindividual. Hence, our rules are simple, and yet applicable over a wide range of settings. Also,although it is clear that external conditions such as food distribution and terrain are importantin determining the pattern of the swarm, we focus on the dependence of the pattern on thebehaviour of the individual ants. Hence, our rules will depend not on external constraints,but on the parameters governing the behaviour of the individuals. The computer simulationis an ideal environment for this type of experimentation. In any real system, it is generallydifficult to separate the subject from its environment. Many physiological experiments can beperformed on an individual; however, it becomes increasingly difficult to design experiments ifone’s interest is in the interactions between many individuals. With the above points in mind,and using the observations of section 1.3 as a guide, we propose the following rules of motion:1. The automatons (ants) move at a fixed speed.2. As each ant moves it deposits a trail pheromone in its path. This deposition will occurregardless of whether the ant was following a trail, or searching randomly.3. The trail pheromone evaporates at a steady rate.4. If an ant is following a trail there is a probability that it will lose the trail. When the antloses the trail it becomes ‘lost’ (see rule 7).95. If an ant following a trail comes to a point at which two or more trails intersect, there isa probability that it will turn and follow the new trail. This will be modelled as follows:i) If there is a trail pheromone directly in front of the ant, it will continue moving hiits current direction.ii) If there is no trail pheromone ahead of the ant, it will move either to the left or tothe right depending on which side has the stronger trail.iii) If both the above tests fail, or if the trails to the left and right are of equal strength,then the ant will move a.s if it were lost (rule 7).6. If a lost ant encounters a trail, there is a probability that it will turn and follow it. Thisprobability will depend on the density of pheromone on the trail, and the angle the trailmakes with the ant’s current direction of motioii. If the ant does not follow the trail itwill remain ‘lost’ (see rule 7).7. If the ant is lost, it will move randomly. However, this random motion will have a certaindirectional persistence (see below).These rules are meant only as a. starting point for the simulations. They are flexible, and canbe generalized to include many types of interactions. For example, the rate of trail depositionneed not be constant. Also, the probability that a lost ant turns to follow a trail may dependnot only on the concentration of the trail pheromone, but also on the density of ants followingthe trail. We have purposefully kept the behaviour of the ants simple in order to minimize theparameters involved in the initial model. As the model is studied, further parameters may beadded, and a more general model developed.The above rules were used to design an algorithm for the computer simulation. The simulation begins by releasing ants one at a time from a central nest location. Each ant is given aposition and a velocity. This velocity is composed of a speed, which is assumed to be constant,and a direction. At the beginning of each time interval each ant takes a survey of its local environment, and chooses a new direction depending on the local concentration and concentration10gradient (computed by discrete differences) of trail pheromone. If there are no trails at the ant’sposition, or if the ant fails to detect a trail, then it will choose this new direction randomly.However, as mentioned in rule 7, the random turn will have a certain directional persistence.Each ant then deposits a pheromone droplet at its current position, moves one step in its newdirection, and assumes a new position. At the end of each time interval the concentration oftrail pheromone at each position is decremented by one unit.The behaviour of each individual is captured by the following parameters: the affinity ofthe ants to the trails, the directional persistence of the random turning, and the rate of traildeposition. If the ants have a higher affinity to the trails, they will tend to remain on a traillonger. A higher directional persistence indicates that the ant will tend to continue in a singledirection for a longer period of time. Likewise, a lower persistence indicates that turns willoccur more frequently. These parameters are explained in more detail in appendix D. Severalother parameters have been included in the simulation, but have not as yet been fully explored.We will discuss these in chapter 3.52.2 Results of the Computer SimulationSeveral runs of the simulation were performed for various values of the parameters. The keyresults of the runs can be summarized as follows:• If there was a strong directional persistence to the movement of the lost ants, then atendency for a few dominant trails to form was observed. If, on the other hand, the antshad lower directional persistence, then a denser network of trails generally formed (seefigure 2.2).• The ratio of followers to lost ants did not depend on the directional persistence of the lostants.• The ratio of followers to lost ants decreased as the affinity of the ants to the trails decreased.11• There was a slight increase in the number of trails in the network, and a drop in thestrength of the trails formed as the affinity of the ants to the trails was decreased. Thisled to a denser network of weaker trails (see figure 2.3).• The density of the network increased as the rate of trail deposition decreased (see figure 2.4).Each of the runs of figures 2.2 to 2.4 was of the same length.In figure 2.1 we compare the results of two runs with the raiding patterns of figure 1.2.This shows that by varying the three parameters in concert, we can shift from the observedraiding pattern of E. burcheili to that of F. hamaturn. That is, the pattern shifts from a highlybranched, area. covering fan to a more tree like structure. It is possible, then, to achieve therange of patterns observed on the field by varying the parameters governing the behaviour ofthe individual, without explicitly including a dependence of the behaviour on precise locationsof food, or other spatial information. Deneubourg [12] suggests that variations in the fowidistribution might be responsible for differences in the raiding pattern. Our results do notin fact disagree with this proposal. We will discuss in the chapter 3.5 how variations in thesurroundings ma affect the parameters governing the motion of the individuals. For example,a homogeneous food distribution may lead to a high degree of random motion, whereas a patchydistribution may lead to a more directed searching behavior.Figure 2.5 shows the development of a stable raiding column. This pattern was observedoccasionally for a high affinity, a strong directional persistence, and a narrow range of depositionrates. In general, the results of the simulation show that a network of trails can form for anyset of parameter values. By varying individual parameters we observe a shift in the resultingpattern from trails which are strongly marked and carry a high proportion of traffic, to weakertrails, and from a dense trail network to a loose one.12(a) (b)(c) (d)Figure 2.1: A comparison of the results of the simulation to field observations: (a) observationsof an E. burchelli swarm; (b) results of a simulation with‘4d = .749, = 255, and r = 12; (c)observations of a E. harnatum swarm; (d) results of a simulation with = .925, = 255, andr = 6. The pattern has shifted from a highly branched network in which the central trails leadto a large area (b), to a network with fewer dominant trails which concentrate activity overseveral smaller areas (d).-313Figure 2.2: The dependence of the pattern of the trail network on the forward persistence of thelost ants: (top) the results of a simulation with a high forward persistence (d = .925, = 255,and T = 12); (bottom) the results of a simulation with a low forward persistence (d = .749,= 255, and T = 12). An increase in the density of the network is observed as the forwardpersistence is decreased...14am*Figure 2.3: The dependence of trail network on the affinity of the ants to the trails. (top) theresults of a simulation with a high affinity (Wd = .749, = 255, and r = 12); (bottom) theresults of a simulation with a low affinity (1d = .749, = 245, and r = 12).4..ISSIs.15(a)(b)Figure 2.4: The dependence of trail network on the rate of trail deposition: a) the results of asimulation with a low rate of deposition (d = .749, = 255, and r = 3; b) the results of asimulation with a medium rate of deposition (‘I’d = .749, = 255, and r = 6); c) the results ofa simulation with a high rate of deposition (‘I’d = .749, = 250, and T = 12).aaa(c)164.4I time 395 %Figure 2.5: Simulation results showing the formation of a trunk trail: ‘Pd .238, = 250,and r = 3. The positions of the ants have been included for clarity. Such trunk trails occuroccasionally in a. small region of parameter space. In a. larger region of parameter space, a trunktrail may form, but not persist.4,vV4.time 248I•S444.SS S- 4• 4.SS4.IS4.Stime 3284.V4.—44.4.S*S4.I17Chapter 3A Model in One Spatial DimensionTo gain a more complete understanding of the effects of the behaviour of the individual on theresulting pattern, we will derive and study a simple mathematical model of the group’s motion.The derivation is based on the velocity jump process of Othmer et al. [25]. In sections 3.2and 3.3, using various techniques of applied analysis we examine special solutions of thismodel known as travelling waves. These can be studied qualitatively using the state spaceof an associated dynamical system. The stability of the travelling wave solution is examinednumerically in section 3.4. Finally in section 3.5 we discuss possible extensions to the model,and examine -the general class of equations to which our model belongs.3.1 Derivation of the ModelConsider again the trails formed in figure 1.2, and their development in figure Li. We maydevelop equations for the full model in two spatial dimensions; however, this proves to be adifficult model to analyze. We instead consider how the parameters of the model may effect thepropagation of the trail network in a single dimension. Let x E (0, oc) represent a distance fromthe nest, measured along the dominant direction of the motion of the swarm. Note that this isnot equivalent to assuming the motion is one dimensional, but merely integrates the variationsin density across the dominant direction of motion. We then study the advance of the swarminto an unexplored region.The ants are divided into two groups. Let F(x, t) represent the total density of ants that arecurrently following a trail, and let L(x, t) represent the total density of ants that are lost andlaying a new trail. Let T(x, L) dx be the total length of trails between x and x+dx. Although the18trails are not oriented, the ants may move along them in either direction. Assume that the antsmove at a constant speed; we denote this speed s for the lost ants, and v for the followers. Todistinguish between the two directions of motion, let F+(x, i) and L+(x, t) denote the densitiesants moving away from the nest (with velocities +s, and +v), and F(x,t),L(x,t) denotethe densities of ants moving towards the nest (with velocities —s, and —v). Pfistner [26] uses asimilar basis for a model of swarming in myxobacteria.As stated in rule 2 of section 2.2, each of these ants continuously secretes trail pheromoneat a constant rate. Since this rate could be different for followers and lost ants, we define Tf asthe deposition rate for the followers, and TI as that for the lost ants. This deposition producesa trail which evaporates at a constant rate -y. We may write the governing equation for theevolution of the trails asajT(x t) = —7T(x, t) + TfF(x, t) +r1L(x, t). (3.1)Since the ants are not perfect followers, they have some probability per unit time of hecominglost ants (see rule 4 of section 2.1). To model this we introduce the parameter E as the rateat which followers become lost. Those ants which remain following the trails may reverse theirdirection. Let p+ be the rate at which followers moving away from the nest turn towards thenest, and let p be the rate at which followers moving towards the nest turn and begin movingaway. The rate at which lost ants find trails (see rule 6 of section 2.1) will increase with thelocal density of the trail network. Let c be the rate at which lost ants encountering trails turnand follow these trails. As an ant encountering a trail can follow it in either direction, let bethe probability that an ant ‘choosing’ to follow a trail reverses its direction. The rate of changeof the density of followers at any given point (x, t) is then governed by the equationsOF F-+ = —cF + c(1 — /3)LT + a/3LT — pF + pF (3.2)8F- 0F— = —cF + u/3LT + c(1 — /3)LT + pF — pF (3.3)If we assume that the ants following a trail have some sense of which direction leads to the19nest and that they prefer to travel away from the nest, then we can assume that p is muchlarger than the other parameters in the system. Equation 3.3 then implies that F(x,t) decaysexponentially on the order of F(x,0)e’. The assumption may not be valid in the laterstages of the swarm; however, for the initia.l stages, when there is a large outward flux of antsfrom the unfolding bivouac, it. will hold. Adding equations (3.2) and(3.3) produces8F a+ - F) = -cF + aLT (3.4)If we assume F to be exponentially small, then F+— F F+ F, and equation 3.4 becomesaF 8F-—- + v— = -eF + aLT. (3.5)at axTo model the motion of the lost ants, assume that the they reverse their direction at randomtimes chosen from a Poisson distribution with mean \. That is, their movements will becomposea of runs at a constant spend either towards,. or away from the nest, se-arate-d atrandom intervals by 180° turns. As with the lost ants which become followers, those followerswhich become lost may either continue moving in the same direction, or reverse their directionwhen they lose the trail . Assume that the same constant governs this probability. We cannow write the following equations for the mean rate of change of the lost ants with time:a aL(x,t) + s—L(x,t) = e(l — /3)F(x,t) + e/3F(x,t) (3.6)— aLT — AL(x,t) + AL(x,t)a a— s—L(x,t) = c/3F(x,t) + e(1 — /3)F(x,t) (3.7)— aLT + AL(x,t) — L(x,t)Adding equations (3.6) and (3.7) we obtaina a+ s—(L — L) = eF — aLT, (3.8)20and subtracting them,0 8— L) + sj_L = — 2i3)(F — F-) — — LiT — 2A(L — Li. (3.9)To eliminate the terms involving L± and F± in the last two equations we first differentiate(3.8) with respect to t, and then multiply (3.9) by s and differentiate the result with respect tox. Taking the difference of the resulting equations yields02L 92L 8 8— =— (F—nLT)— 2As—(L — L) (3.10)0+ sb— [E(1 — 2/3)(F — F-) — a(L — L)T].Although this equation is hyperbolic, if we assume that s2 >> 1 and that A >> 1 then theequation can be approximated by the parabolic equation02L 0(32/2x)= s—(L — 17). (3.11)This is equivalent to observing the motion of the ants on a time scale which is larger than therate of reversal of the lost ants, and on a spatial scale which is much smaller than the distancea lost ant would have travelled if it had moved at a constant speed s without reversing. Thisis known as the diffusion limit and is discussed in more detail in Othmer et al. [25].Finally, using (3.8) we replace the flux term on the r.h.s of the above equation by quantitiesinvolving only the total densities F,L, and T. This yields the following interaction-diffusionequation for L:OL 2 02L=- eF + aLT. (3.12)Collectively, equations (3.1),(3.5), and (3.12) represent the evolution of the trail networkin the case where the frequency of the random turns of the lost ants (A), and the rate atwhich ants following trails towards the nest turn away from the nest (pj are large relativeto the remaining parameters of the system. After deriving these equations from a simple21random process (microscopic), it is useful to discuss the physical and biological significance(macroscopic) of each term. For convenience we repeat the complete system below.= T1L + Tf F — (3.13)IJt(i) (ii)a--= —----(vF) —€F + aLT} (3.14)(iv)(ni)DL 62L= ---+eF-aLT. (3.15)(v)T represents the total length of trails per unit length of the strip (0, oo), and F and Lrepresent the total densities of ants moving left (towards the nest) and right (away from thenest) respectively. The parameter j = is known as a motility coefficient, and representsthe diffusivity of the lost ants due to random motion. The terms over brace (i) depict thereinforcement of trails by each ant, and term (ii) the decay due to evaporation of the trails.Term (iii) of equation (3.14) represents the directed motion, or convection of followers alongthe trails at the velocity v. The terms grouped under (iv) which appear in equations (3.14) and(3.15) represent the rate, E, at which followers lose the trail, and the rate, a, per unit length oftrail, at which the lost ants begin to follow an existing trail. Finally, the exploratory motionof the lost ants is represented by the diffusion term (v) with a motility coefficient ii. A furtherdiscussion of the general form of these equations appears in section 3.5For the remainder of this chapter we will be interested in finding solutions to equations (3.13-3.15) which represent the propagation of the swarm. One distinguishing trait of any swarm isthat it maintains a semi-rigid profile as it moves. Okubo [24] gives an excellent account of thedifferences between coherent swarming motion and simple diffusion, and describes the forcesthat are necessary to maintain the swarm shape. Solutions representing the movement of a22fixed profile are known as travelling waves. For examples of such waves and their applicationsto biology and chemistry see [15, 23].3.2 Reduction to the Travelling Wave FormIn the following analysis we assume that each of the coefficients in the system (3.13 - 3.15) isa constant. Deviations from these assumptions will be discussed at the end of this chapter. Afirst step in the analysis is a reduction of the equations to a dimensionless form. The details ofthis are left to appendix A. Here we merely note that upon introducing the rescalings:au2 T — —T - = —v ‘= —t,II EILF*= —i-F, —, (3.16)V -‘IV V2= —x,/1 aILTj TfL * = —L, r” = —.-yv2 Tiand dropping the *‘s, we obtain the dimensionless equations:T = -y(L + TF — T). (3.17)F = —F—EF+LT, (3.18)= L + eF — LT. (3.19)We are interested in examining solutions to this system which represent waves of a population propagating into an unexplored area. These solutions are characterized by a fixed profilemoving at a constant speed c. Such solutions are known as travelling waves, and are studiedby transforming the system to the moving coordinates z = x — ct, and t’ = t. The steady statesolutions of the PDEs in these new coordinates, found by setting the derivatives of T, F, andL with respect to the new time t’ to zero, correspond to waves with a fixed profile moving at aconstant speed c in the original coordinate system. These steady state solutions will satisfy anautonomous system of ODEs in the variable z.23To represent the system (3.13-3.15) in the moving coordinates we use the chain rule ofcalculus:OF OFOz OFOt’—(3.20)OF OF 9z OF at’=Oz at’and — are found by again applying the chain rule to the transformations z(x, t) = x— etax Oxand t’(x,t) = t:az az azdz = —dx + —dt dx—cdt —=1, —=—c;ax at ax at (3.21)at’ at’ at’ at’dt’ = -—--dx + —dt = dt — = 0, — = 1.ax at ax atIf we substitute these results back into (3.20) we find that:OF OFa2 — ax‘ (3.22)OF OF OFaTUsing these results we transform (3.17-3.19) to the moving coordinates, and set-b—- =aF OL= 0, and = 0. This produces the following system of ordinary differential equations:—cT’ = -y(L + rF — T); (3.23)(1 — c)F’ = —cF + LT; (3.24)—cL’ = L” + F — LT. (3.25)Here, the ‘ indicates differentiation with respect to the wave variable z.These equations can be further simplified by the following operations. First, add equations(3.24) and (3.25) to produce the equation(1— c)F’— cL’ = L”. (3.26)24This equation can be integrated once from —oo to z to give(1— c)F—cL = L’+k. (3.27)The constant of integration, k, is given byk = cL0 — (1 — c)Fo + (L’)o, (3.28)where the subscripts’ indicate evaluation at the point z = —.Replacing equation 3.25 with equation 3.27 leaves us with an autonomous system of firstorder ODEs. Such a system is commonly referred to as a dynamical system. The analysis ofthese systems is introduced for applications to travelling waves in [15, 23], and in a mathematicalsetting in [20, 31]. Since the system of ODEs is autonomous, its solutions can be representedas curves in the three dimensional state space (T, F, L). In addition, at each point in this statespace there is an associated vector (T’, F’, L’) which is tangent to the solution curves passingthough that point. A fixed point is defined to be a point where the tangent vector (T’, F’, L’)vanishes. A homogeneous spatial distribution which is a steady state solution of the equationsin the original coordinate system, (3.13-3.15), appears in the state space as a single fixed point.We can study the shape of solutions to the PDEs by examining the trajectories of thesolution curves in the T-F-L state space. If these trajectories are to represent a wave profiletraveffing through a population, then they must satisfy the following criteria:• Populations must remain bounded. Hence, the trajectories are restricted to homodinicorbits (closed curves passing through a single fixed point), heterodinic orbits (curvesconnecting two fixed points), or limit cycles (closed curves which do not pass through anyfixed points).• A population density must remain either positive or zero. Hence the trajectories must becontained in the positive octant of the state space.1The use of zero as a subscript will become apparent soon.25• The waves represent a swarm propagating into an empty region of (physical) space. Hencethe trajectory must end at the origin (T, F, L) = 0. Further, the origin represents ahomogeneous spatial distribution (no ants, and no trails) and must therefore be a fixedpoint. With regards to the variable z, the solution will asymptotically approach (0, 0, 0)as z cc. This fixed point will be referred to as the trivial steady state, and denoted byU’.• Finally, we expect that some fixed density is established behind the wave; hence thetrajectory must originate at a fixed point in the state space. With regards to the variablez, the solution vi1l asymptoticajly approach this fixed point as z— —cc. This fixed pointwill he referred to as the populated steady state, and denoted by u0Thus the travelling wave of interest will he represented in the state space by a curve connectinga point in the first octant to the origin.An inspection of (3.27) shows that the origin u1, will be a fixed point if and only if k = 0.In addition, if it0 is a fixed point then by definition (L’)o = 0. With these simplifications (3.28)reduces to cL0 (1— c)Fo. Thus, given any initial steady state (To, F0,L0)T, there will be aunique wave speedF0= (L + F0) (3.29)for which the final state will be unexplored territory. We have now reduced the original systemof PDEs to the following system of first order ODEs:7TI = —[T— L— TF]; (3.30)F’= (l_C)[+LT]; (3.31)L’ = (1— c)F — cL. (3.32)The solution of interest satisfies the boundary conditions26T T0F = F0 =uo,L L0(3.33)T T1 0F = F1 = 0 ti,L L1 0+and the restriction that T, F, L he positive for —00 < z < +00. Any solution to this system willbe a solution to the original system of PDEs which represents a travelling wave propagatingfrom the populated steady state u0 into the trivial steady state u1 with a velocity c.Setting tile i.h.s of equations (3.30-3.32) to zero yields the system of equations for tile fixedpoint s:T L+TF; (3.34)= LT; (3.35)(1 — c)F cL. (3.36)There are only two points satisfying these constraints. One of these is the origin, and the other,which must be u0 if the travelling wave is to exist, is given by(C1—c= (C (337)(1— c)(i— c + rc)(C1 — c + rcNote that since u0 lies in the positive octant of tile state space, we mllst restrict tile values ofc so that 0 < c < 1. Recalling that we are dealing with a dimensionless system, tills 15 simplya statement that the wave propagate forward, and that tile speed of propagation of the swarmfront be 110 faster that the velocity v of the followers in the original coordinates.273.3 Analysis of the State SpaceThe problem of finding travel]ing wave solutions to the initial system of PDEs has now beenreduced to that of finding bounded, positive solutions to the system of three first order ODEs.To show that such solutions exist will require a further examination of the state space of thesystem. First, we perform a local analysis of the fixed points u0 and u. This will provide uswith information on the behaviour of the dynamical system near these two points. Next, weperform a global analysis to determine the nature of the solutions in regions containing boththe fixed points. In this manner, we establish conditions for the existence of solution curvesconnecting these two points.3.3.1 Local Analysis of the Steady StatesLinearizing equations (3.30) thru (3.32) about a point (T,F,T) yields the system:L T-TF’_____—__________F—F (3.38)(i-c) (1-c) (1-c)L—L0 1—c —cSubstituting (T, F, L) = (0, 0,0), it follows that at the fixed point u1, the linearized system hasthe eigenvalues:= 7/c; >0= —E/(1—c); < 0 (3.39)A3 , <028and eigenvectors:(1—c)2+r(c(1—c)—€) 71____________—i—C+Ec/7 C= 0 — C 0 (3.40)i—c0 1—c—+cCSince 7, c, c > 0, we have that A > 0, and A2, A3 < 0 over the entire parameter space. Thisindicates that u is the intersection of a two dimensional stable manifold Wj9, and a onedimensional unstable manifold 1’V.At the populated steady state u0 the eigenvalues of 3.38 are roots of the cubic equationA3+AA2B +C=O, (3.41)whereA = c + 1 — c 7/c, (3.42)I’€B=— I I + I (3.43)c(1— C + TC)(7= i—c (3.44)Following the technique of Derrick [13] outlined in appendix B, we can use the fact that B < 0,and C > 0, for all values of the parameters to show that two of the eigenvalues will havereal parts greater than zero, and the third will be real and negative. Thus, although it is nota simple task to determine the exact elgenvalues and eigenvectors of the populated state forgeneral parameter values, we can assert that there will always be a two dimensional unstablemanifold W’ and a one dimensional stable manifold W intersecting at u0. Further, a detailedanalysis (see Appendix C) of the local dynamics near this fixed point shows that both of theunstable eigenvectors must point into the region of phase space where each of the derivativesT’, F’, and L’ are negative, and that each of their components must be negative.29We now know the behaviour of the vector field, and thus the dynamics of (3.30-3.32) nearthe fixed points of the system We also know that the system has only two such fixed pointsIf there is to be a trajectoiy tiavelling fiom u0 to u1 it must be contained in both the unstablemanifold of Uo, and the stable manifold of tt. To show that such a trajectory exists will requirea knowledge of the global behaviour of the dynamical system.3.3.2 Global AnalysisFor systems of two ODEs, a qualitative analysis of the global structure of the state space canoften be performed quickly and easily. In contrast, our system consists of three ODEs, andsuch analysis proves more difficult. By a careful examination of the geometry of the statespace and the flow of the associated vector field, and by using some numerical experiments weare able determine conditions under which the existence of the desired solution is guaranteed.Furthermore, we can determine the shape of the wave, and the effect of the parameter valueson that shape.-To begin, we note that the heterodinic orbit must be contained in both Wj9, the stablemanifold of u1, and T’V’, the unstable manifold of ‘u0. Conversely, since there is a uniquesolution passing through each point in state space which is not a fixed point, if W intersectswith W at any point then this is sufficient to guarantee the existence of the heterodinic orbit.We show that W will fail to intersect W only if there is a separatrix2 that bounds W.However, numerical experiments gave convincing evidence of a structurally stable heterodinicorbit for each set of parameter values tested, and so we suspect that such a separatrix does notarise in this system.Let T°, F°, and L° denote the nuhisurfaces of the vector space. These are the surfaces wherethe derivatives of the state variables vanish, and are given by the following equations:T°= { (T, F, L) I T — L — TF = 0 }; (3.45)2A solution that is topologically different from neighbouring solutions. In JR3 this may be a family of suchsolutions (see [4, 20])30F°= { (T, F, L) I LT — eF = 0 }; (3.46)L° = {(T,F,L) (1— c)F— cL = 0}. (3.47)These nulisurfaces are shown in figure 3.1. Let T+ and T be the regions in which T’ is positiveand negative respectively. These are the sets defined by:T = {(T,F,L)IT—L—rF>0};T = {(T,F,L)IT—L—TF<0}.The sets F±, L± are defined in a similar fashion. Using the directions ‘left’, ‘right’, ‘up’, ‘down’,‘ahead’, and ‘behind’ defined in figure 3.1, T lies above T°, F— lies below F°, and L lies tothe right of L°. Finally, we define the following sets (as illustrated in figures 3.1 and 3.2):Q1 = {(T,F,L)IT,F,L>0};111= {(T,F,L)IF=0, T,L>0};112= {(T,F,L)IT=0, F,L>0};= F+flL; (3.48)FflL;B T°flA;E T°flW.The region Qi is the first octant of the state space. 111 is the back wall of the first octant,and 112 the floor. is the region in which both F and L are increasing, and is located aboveF° and to the left of L°. A is the region in which both F and L are decreasing. It is locatedbelow F° and to the right of L°. The regions B and E lie in the plane T° as shown in figure 3.2.To analyze the flow of the vector field in the state space, the following observations arenecessary:1. The flow through the surface F°flT flL (the portion of F° to the right of L° andbelow T°), is oriented from F+ into F—. That is, it passes through F° from above.31T0downnulisurface intersecUonsF-L.:.:::.:... T-FT-L (hidden)Figure 3.1: The State Space (Octant I). showing the nulisurfaces T°, F°, and L°, and theregions A±, 1 and 2•backforwardleftT-axsF°rightF-axis— L-axisI —322. The flow through 11’ is into Qi (forward), since F’> 0 on 11’.3. The flow through 112 is out of Q (down), since T’ > 0 on 112.4. The flow through L° is from L into L (left to right) below the F° surface, and fromL— into L+ (right to left) above it. To see this note that, since L° is a vertical plane, thedirection of flow across this plane will be determined by the sign of F’.From these points we can conclude that any trajectory initiated in A— flQi flT—, that is,the portion of the first octant below T°, below F°, and to the right of L°, must exit that regionthrough either B, or 112. this follows since the flow on all other boundaries of this region isdirected into this region. Further, a trajectory cannot remain in this region since at each point,the flow velocity is bounded away from zero, and each component of the flow is decreasing.Thus, not oniy are there no fixed points, but there are no closed curves (or limit sets) containedentirely in this region.The unstable manifold W’ must pass through u0 tangent to the unstable eigenspace ofthe linearized system near u0 (see Wiggins [31]). By the analysis of Appendix C, the unstableeigenvectors of this system are directed into the sets A— and A+. If we consider the local unstablemanifold as a disc tangent to W’ spanned by these unstable eigenvectors, then this disc mustbe imbedded in the state space near u0, and intersect both A+ and A—. This immediatelyshows that the intersection of T’V’ with T°, which is the intersection E, will be non-empty.Figure 3.2 shows the intersection of T° with Q as viewed from below. The intersection Ewill be a continuous curve on this plane which must originate at u0. By the observations ofthe preceding paragraph, E must also originate in B. Also, a trajectory which intersects Acannot exit this region through F° unless it first passes through B. Hence the curve E cannotcross the line F° fl T°. It can also be shown that trajectories in W0U which pass through frombelow must do so from within A—. Thus it follows that E does not cross the line L° fl T°. Thisargument involves first showing that the portion of W’ forward of u0 must also lie to the rightof u0, and then that a trajectory in this region must pass through L° into A— before passing33T° ntersect L0Figure 3.2: A view of T° from below, showing the sets B, E, and the intersection of W withT°.through T°. These restrictions on the curve E imply that it is contained entirely in B. Thus,it appears that unless l’V’ is bounded in some manner, the curve E will extend to u. Thiswould imply the existence of a heteroclinic orbit connecting UQ to u1.Now, we examine the stable manifold 1’V. From the eigenvectors given by 3.40, it followsthat l47j intersects A— for all parameter values. Also, a similar analysis to that used to examineW in the above paragraph, it follows that WjS intersects F°, and L°, but cannot intersect 112or B. Further, the intersections of T47j with L° and F° cannot rise above T°, nor can they fallbelow 112 (see figure 3.3). Since the flow in A— below T° is bounded and oriented downwards,we can conclude that Wj5 must either approach arbitrarily close to u0, or pass below it.If W? passes below u0, and does not intersect W’, then there must be a third manifoldpositioned between u0 and l19. This manifold will intersect B, and prevent E from reaching u.If fVj does not intersect u, it is still possible that E extends to u1. In this case the heterocinicwill exist, and may not be unique. 1-lowever, it will not be a structurally stable feature of thevector field. Finally, if W19 approaches arbitrarily close to u0, then the previous observationsabout the location of l’V imply tl1at the two manifolds must intersect. This is an intersectionof two two-dimensional manifolds, and will be structurally stable to small perturbations of thevector field. This means that such an orbit occurs at a single point in parameter space, it mustT0 intersect F°U134T0 intersectT0 intersect F°U0Figure 3.3: A schematic showing the boundaries of A— and their intersections with Wj. TheF° surface has been distorted for clarity.also occur at all parameter values in a neighbourhood of that point.We can now use numerical techniques to show the existence of the homoclinic orbit forspecific parameter values. To do this we show that the intersection of 14’ and Wj9 containsat least one point. As mentioned above, since the tangent vector at each point is unique, therecan be oniy one solution passing through each point in the state space. By definition, thesolution passing through any point in the unstable manifold of the fixed point u0 tends to u0as z — —cc. Also, the solution of any point in the stable manifold of the fixed point u1 tendsto u as z — cc. Thus, if a point x lies in both of these manifolds, then there must be aunique solution passing through x which connects u0 to u1. Numerically we can compute thelocal stable and unstable manifolds of u1 and u0 respectively (that is, the eigenspaces of thelinearized system), and solve the system forwards in time from u0, and backwards in time fromu1. To determine if the two manifolds intersect, we compute their intersection with a plane ofconstant T. This intersection is shown in figure 3.4. The points represent the intersections ofF°Intersections of the stable manifold of u1with F0. and L0.35the manifolds with a plane of constant T. The fact that these intersections cross implies thatSS.— .S.• 4’b•SFigure 3.4: A section showing the intersection of W0’ and W? with a plane of constant T:c = 1,7 = 4, c = .25, T = .17. The crossing of the intersections implies that W and Wj havea common point, and thus a common trajectory. The figure was created using PhasePlane [18]there is a point on the plane which is common to both W’ and Wf. In this manner we candetermine the existence of a heteroclinic orbit for any choice of parameter values. Further, theintersection will persist for small perturbations from these values, and by the above arguments,the heterodinic orbit will persist as well. That is, the heterodinic orbIt is stable at least tosmall perturbations of parameter space. These experiments were performed on various pointsin parameter space, with similar results in each case. Thus we provide a strong case for the36existence of the travelling wave solution, if not for all parameter values, then over a large portionof the parameter space.3.4 Stability of the Travelling Wave as a Solution to the PDEs.Note that c is not a parameter of the original system of PDEs. This implies that there is aredundancy in the solutions to the ODEs. That is, we only require the existence of the heterocinic orbit in the dynamical system for a single speed c given any point in 7, E, T-space. Todetermine the stability of the travelling wave as a solution to the original system of PDEs, wesolve the full system of PDEs (equations 3.13 to 3.15) for various initial conditions. The solutions were found using the method of lines with a fourth order implicit Runge-Ktitta algorithm.This method consists of discretizing the spatial domain of the PDEs, yielding a coupled systemof ODEs which can then be solved using the Runge-Kutta algorithm. An upwind scheme wasused for the discretization of the parabolic equation for F. The travelling wave was found tobe stable for each case tested. Further, an oscillating wave was never observed. The results areshown in the figures 3.5 and Implications and GeneralizationsIn summary, although we have been unable to show the existence of the travelling wave solutionfor every parameter value, we have shown that the wave does exist at certain parameter values.Also, we have shown that it persists for small perturbations from these parameter values.Thus, the heteroclinic orbit which represents the travelling wave in the (T, F, L)-state space isstructurally stable to small perturbations of the parameters. Using a qualitative global analysisof the dynamical system, we were also able to determine the shape of the travelling wave, andthe effect of the parameters on that shape. These effects can be summarized as follows:1. We noted that the speed of the wave would depend on the ratio of the density of thefollowers to the total density of ants at the nest. That is, a higher proportion of followersto lost ants at the nest would result in a faster propagation of the trail network.37LegendFollowers - time = 100-time= 80-time= 60-time= 40Lostants -time= 100-time=80-time=60-time=40— . —.— Trails- time = 100-time= 80-time= 60-time= 40Lost ants - initial dataFollowers - initial dataSLOPED !WP1J\ N \\\. \\ \.\,“.\‘.•\.C000. 20 40 60 80 100distance from nestFigure 3.5: Stable travelling wave solution: = .1, c = .5, and = 138C0D0.00.60.4- 0.20OSC!LLA11WG !t4P’J1’ - ?NtT1AW( LOS1100Figure 3.6: Stable wave resulting from oscillating input: e = .1, c = .5, and y = 1. All antswere initially lost.Legendfollowers - time = 20______- time = 40-time =60- time = 80-time = 100-time=20-time=40-time=60-time=80-time=100• initial datalost antslost ants0 20 40distance from nest80392. The travelling wave solution represents a realistic profile of a population density; particularly in that it is bounded, and everywhere positive. This follows from the fact that ifthe heteroclinic orbit exists, it must lie in the first octant of the state space.3. If the eigenvalues of the linearized system near u0 are real, then the solution will bemonotonic decreasing in each of the variables T, F, and L. If they are complex, then thesolution will oscillate as it approaches u0. The oscillating wave form was not observed asa stable solution of the PDEs.4. The monotonic travelling wave was observed to be a stable solution to the original systemof PDEs.The model also shows that the density of the lost ants attains its largest value close to the nest.That is, L(z) < L0 at any point z along the wave. It would be desirable to develop a model iiiwhich the lost ants obtain a maximum near the front of the wave. To examine this possibilityfurther we study the class of equations intO which our model falls.The equations are part of a larger class given by the following system:= G(T, F, L)— T; (3.49)F =_- (vF) — R(T,F,L); (3.50)02L= ii----+R(T,F,L). (3.51)Equation 3.51 represents the growth and decay of the network. The network decays exponentially, and its growth is governed by the two interacting species F and L according to thefunction G. The first species, F, is transported along the network at a velocity v; whereas thesecond, L, diffuses with a motility coefficient of t. The the interactions of the two species arerepresented in the function R. Thus, the behaviour of the lost individuals is analogous to thatof a diffusing substance, and the behaviour of the followers is similar to that of a substancewhich is being either transported or convected through a network. These two substaices can40be considered to interact both through and with the trail network. The network itself can becompared to a reactive medium.Using tile techniques of section 3.2 we may determine the associated dynamical system forthe travelling waves admitted by the PDEs (3.49-3.51):—cT’ = G — (3.52)=— (vF)’ — (3.53)L’ = (v — c)F — cL. (3.54)(3.55)From this we see that the shape of the nullsurface L° is determined by the parameter v, andthat changes in the reaction term R do not affect the shape of this uullsiirface. In fact theequation R(T, F, L) = 0 determines the nuilsurface F°.One method of generalizing the original model is to allow the parameters to be functions ofthe position in the state space. In fact, the experiments cited by Calenhuhr [6] indicate that€ is an increasing function of the trail strength. Further experimentation via both numericalanalysis, and the simulation, is ullder way to determine the effects of altering the parameters inthis fashion. The results of this study will be presented in a later paper. Note that the only wayof forcillg the density of the lost ants to attain its maximum nearer the wave front is to changethe nullsurface L°. From an examination of the generalized equations, we see that this can beaccomplished only by introducing a dependence of v on the state variables T, F, and L. TheT2 + v*F2dependarice= T2 + F2was introduced in a numerical analysis of the PDEs; however, nosignificant difference was noted in the results. The shape of the stable wave solution remainedunchanged.When generalizing equations in this fashion, one must be careful with the interpretation ofthe results since we have not maintained a rigorous connection between the behaviour of theindividuals and the parameters in the model. However, we may now note the similarities of theform of these equations with other models; two interesting examples are the motion of gliding41myxobacteria [26], and the transport of chemicals within cells [27].42DiscussionThe main goals of this thesis were to study the collective behaviour of a group of interactingindividuals, and to determine how the patterns of motion of the group are affected by variationsin the parameters governing the behaviour of the individuals. We focused on the specificexample of the propagation of a swarm of army ants. This problem was modeled by a system ofautomata interacting with a trail network. The behaviour of the each automaton was governedby a set of rules which were dependent only on the local density of the trail network. Theserules did not explicitly refer to any external variables such as the distribution of food sourcesor variations in the terrain. Rather we suggest that such variations affect the parametersgoverning the behaviour of the individuals and do not directly influence the collective motion ofthe population. The model was studied using both a computer simulation, and a mathematicalmodel.To visualize the motion of the swarm in two dimensions, a computer simulation of a systemof automata. was developed. Each automaton was equipped with a set of simple rules whichgoverned its behaviour, and its interactions with the network of trails. These interactionswere represented by the following parameters: directional persistence of the ants not followinga trail, affinity of ants to the trails, and the rate of deposition of new trail material. Thepatterns were studied to qualitatively determine the density of the trail network, the tortuosityof the pattern, and the distribution of lost ants relative to followers. At first, each parameterwas varied individually to determine its effect on the resulting patteri. By varying all threeparameters in concert, we were able to obtain patterns characteristic of the swarms observed innature. Hence, we have shown that the parameters governing the behaviour of the individualare sufficient to model the problem. External influences may then be responsible for reinforcingthese behavioural responses, and maintaining a stable pattern.43In addition to the simulation, a mathematical model was derived from the rules governingthe behaviour of the individual automata. To simplify the analysis, the mathematical modelexamines the behaviour of the swarm in only a single spatial dimension. That is, we examinethe density of both the ants and the trail network as a function of the distance from the nest.The mathematical model provides a means of determining the relationships between certainaspects of the swarm pattern and the parameters of the model. Specifically, we have shownthat the swarm will propagate as a travelling wave for any of the parameter values, and thatthe speed of this wave will increase as the affinity of the ants for the trail increases. Througha further qualitative analysis we have also determined the shape of the travelling wave for anyparameter value. Numerical experiments show the wave to be stable. This model was alsogeneralized to consider the class of equations that could be used to model similar behaviour.In this way, it is possible to link the phenomenon in question with other problems.The process used to generate the equations for the one dimensional model can also beextended to derive the equation of motion in two dimensions. We have extended the analysisof Othmer et al. [25] to model the motion of automata interacting with a trail network intwo dimensions. Although the full system in two dimensions is difficult to analyze, severalsimplifying assumptions can be made. For example, Edelstein-Keshet [16] has developed amodel with no spatial dependence to examine the orientational aspects of the pattern, and todetermine the conditions necessary for a transition from a loose network of trails to a networkwith a strong directional order. By deriving this and other simple models from the more generalmodel in two dimensions it is possible to accurately integrate the results obtained with each ofthese models.Having developed both the mathematical model and the simulation, we may now use themto study the effects of slight alterations in the behaviour of the individual on the formationand stability of the trail network. It is obvious tl1at the full behaviour of the swarm cannot becaptured in any simple mathematical model. However, by using simple models we may studythe dependance of the resulting patterns on the behaviour of the individuals. Once these simple44models have been analyzed, they may be extended to include a more complex behaviour. Forexample, Calenbuhr [61 has shown that the mean distance an ant follows a trail will increaseas the strength of the trail increases. This implies that the probability that a follower loses thetrail at any point will decrease as the trail strength increases. This has been implemented inrecent versions of the simulation, but studies of the model are not as yet complete. Anotherpossibility is that the followers may reverse their direction of travel along the trails if they donot encounter other ants following the same trail. Gordon [19] has recently experimented withsuch behaviour.We are not the first to use a simulation to model the swarming behaviour of army ants.Deneubourg et al [12] have used a simulation to study the motion of ants swarming into a regioncovered with various food distributions. His automata (ants) initially swarmed forward froma nest leaving a small amount of trail pheromone in their path. At each step forward, eachant’s choice of a new direction depended on the relative distribution of pheromone ahead leftand ahead right of its current position. If an ant found food, it returned to the nest leaving alarger pheromone deposition in its wake. Deneubourg observed that as the swarm progressedforward a network of trails formed immediately behind the front, and a single trail dominatednearer the nest. By varying the food distribution, various patterns of trail networks arose.Three possible patterns are shown in figure 3.7. Case (a) shows a raiding pattern typical of E.harnatum, obtained for a patchy food distribution. Case (c) was obtained using a homogeneousfood distribution, and is more characteristic of the broad swarm of E. burchelli.To examine the formation and stability of these trail networks Ermentrout and EdelsteinKeshet [17] studied a system of automatons (ants) interacting on a torus through the use oftrail pheromones. They began with an initially random distribution of ants and found thatthree distinct final patterns arose: random milling; the formation of strong trails, but no stablepattern; and the formation of a small number of very strong, long lived trails. In addition tothe simulation, they also analyzed the formation and stability of strong and weak trails usinga mathematical model. Their results indicated that no stable trails would form if the density45Step 900Figure 3.7: Simulations of Deneubourg et al. [12] showing the dependence of the raiding patternon the food distribution. The inserts are sketches of actual swarms taken from Schneirla [28].a mathematical model. Their results indicated that no stable trails would form if the densityof ants was below a certain threshold; strong stable trails would develop and persist at largerdensities; and if the density was increased still further, both weak and strong trails wouldpersist.The work of Deneubourg et. al .suggests that there is a significant dependence of the finalpattern on the distribution of food in the environment. They concluded tha.t a single mechanismis implemented by different species, and that the environment is responsible for any differencesin behaviour. We have shown tha.t it is possible to obtain a wide range of patterns by varying theparameters governing the behaviour of the individuals without referring explicitly to changesin the environment. This is not in contradiction with the results of the previous simulationsand experiments of Deneubourg and Franks. We suggest that variations in the environment toE. burchelliE. hamatumE. raPaxStep 1100 1 Step 130046not directly influence the motion of the ants, but rather influence tile parameters governing thebehaviour of the ants. That is, a change in the environment effects a change in the collectivemotion of tile swarm through the chain3: environment — individual behaviour — individualmotion — collective motion. For example, an increase in the local food density may lead to adecrease in the directional persistence of the individuals, which then leads to a collective motionwhich will cover a given area more effectively. Our simulation then, studies the effect of thebehaviour of the individual on the collective motion of the group.We have developed a model for tile collective motion of army ant raiding swarms. This modelis based on parameters governing the motion of individual ants. We have shown, using both asimulation and a mathematical analysis, that the model exhibits similar swarming behaviourboth to previous simulations based on different assumptions, and to patterns observed in nature.31t is interesting to note that this is not a deterministic chain of cause and effect; there is a possibility ofrandom fluctuations in the effects at each link in the chain.47Bibliography[1] K. Able. Mechanisms of orientation, navigation, and homing. In S. Gauthreaux, editor,Animal Migration, Orientation, and Navigation, pages 283—373. Academic Press, N.Y.,1980.[2] S. Aron, J. M. Pasteels, and J. L. Deiioubourg. Trail-laying behavior during exploratoryrecruitment in the argentine ant, iridomyrmex humilis (Mayr). Biology of Behaviour,14(3):207—2i7, 1989.[3] R. Beckers, J. L. Deneubourg, S. 4oss, and J. M. Pasteels. Collective decision makingthrough food recruitment. Insectes ,Sociaux (Paris), 37(3):258—267, 1990.[4] G. W. Bluman and S. Kumei. Symmetries and Differential Equations. Applied Mathematical Sciences 81. Springer-Verlag, Berlin, 1989.[5] V. Calenbuhr and J. Deneubourg. A model for trail following in ants: individual andcollective behaviour. In W. Alt and G. Hoffmann, editors, Biological Motion, Proceedings, Königswinter, 1989), Lecture iVotes in Biomathernatics, 89, pages 453—469. SpringerVerlag, Berlin, 1989.-[6] V. Calenbuhr and J. L. Deneubourg. A model for osmotropotactic orientation. Journal ofTheoretical Biology, submitted.[7] S. Carnazine. Self-organizing pattern formation on the combs of honeybee colonies. Behavioral Ecology and Sociobiology, 28:61—76, 1991.[8] J. Deneubourg, S. Aron, S. Goss, and J. Pasteels. The self-organizing exploratory patternof the Argentine ant .Jo’urna.i of Insect Behavior, 3(2):159—168, 1990.[9] J. Denetibourg, S. Coss, Sendva-Franks, C. Detrain, and L. Crétian. The dynamics ofcollective sorting; robot-like ants and ant-like robots. In Simulation of Animal Behaviour;From Animals to Animuts. MIT Press, Cambridge.[10] J. L. Deneubourg, Goss, Beckers, and Sandini. Collectively self-solving problems. InBabloyantz, editor, Self organization, Emergent Properties, and Learning. Plenum, 1991.[11] J. L. Deneubourg and S. Goss. Collective patterns and decision making. Ethology, Ecologyand Evolution, 1, 1989.[12] J. L. Deneuhourg, S. Goss, Franks, and J. M. Pasteels. The blind leading the blind: Chemically mediated morphogenesis and army ant raid patterns. Journal of Insect Behavior,2:719—725, 1989.48[13] B. Derrick, personal communication.[14] C. Detrain, J. Deneuborg, S. Goss, and Y. Quinet. Dynamics of collective exploration inthe ant Pheidole pallidula. Psyche, in press.[15] L. Edelstein-Keshet. Mathematical Models in Biology. Random House, 1988.[16] L. Edelstein-Keshet. Trail Following as an Adaptable Mechanism for Population Behaviour. preprint, 1992.[17] L. E.delstein-Keshet and B. Ermentrout. Trail following in social insects: Patterns ofmotion and behaviour of individuals. preprint, 1991.[18] B. Ermentrout. PhasePlane. 1990.[191 D. Gordon. personal communication.[20] J. Gnckenheimer and P. Holmes. Nonlinear Oscillations, Dynamical Systems, and Bifurcations of Vector Fields. Applied Mathematical Sciences 42. Springer-Verlag, Berlin,1983.[21] B. Hölldobler aud E. 0. Wilson. The Ant.s. Harvard University Press, Cambridge, 1990.[22] K. Jaffe, C. R.amos, C. Lagalla, and L. Parra. Orientation cues used by ants. InsectesSociaux (Paris), 37(2):101—1L5, 1990.-[23] J. Murray. Mathematical Biology. Springer Verlag, Berlin, 1989.[24] A. Okubo. Diffusion and Ecological Problems: Mathematical Models. Springer Verlag,Berlin, 1986.[25] H. G. Othmer, S. Dunbar, and W. Alt. Models of dispersal in biological systems. Journalof Mathematical Biology, 26:263—298, 1988.[26] B. Pfistner. A one-dimensional model of the swarming behavior of myxobacteria. In W. Altaud G. Hoffmann, editors, Biological Motion, Proceedings, Königswinter, (1989), LectureNotes in Biomathematics, 89, pages 556—565. Springer Verlag, Berlin, 1989.[27] M. Reed and J. Blum. A model for fast axonal transport. Journal of Cell Motility, 5:507—527, 1985.[28] C. W. Rettenmeyer. Behavioral studies of army ants. University of Kansas ScientificBulletin, 44:281—465, 1963.[29] T. C. Schneirla. Army Ants, a study in social organization. W. H. Freeman & Co., 1971.[30] G. Strang. Linear Algebra and its Applications. Academic Press, INC., Orlando, Florida,1976.[31] S. Wiggins. Introduction to Applied Nonlinear Dynamical Systems and Chaos. Springer—Verlag, Berlin, 1990.49Appendix ADimensional AnalysisThe ideas of dimensional analysis were introduced by the engineering physicist Buckingham,and are used extensively in physics to reduce the number of independent parameters appearingin a problem. These ideas are part of a larger theory of using symmetries in the reduction ofdifferential equatiois (see Bluman and Kurnei [4]). We present a simple method of using anarbitrary scaling of a differential equation to reduce the number of independent parameters.Equations (3.13 - 3.15) have two independent, and three dependent variables, each of thesevariables can be rescaled independently without changing the nature of the solutions. To rescalethe system we introduce the new coordinates t*,x*,T*,F*,L* EJ defined byt x T = TT, F = FF*, L = LL*, Le B?. (A.1)Substituting these values into the PDE system yields.TaT*= rlLL*+r/.F*_7TT*,- (v7) : - EF* + cJJTL*T*,jJ____= 8(xj2 +EF*— aJJTL*T*.These equations may be reorganized to the form1OT*L* + F*-T*,\.yTJDF 1v’\ ÔF* - oLTT=— I I — EtF* + — L*T*,F50/ —\ 2r* —/ ‘-‘- Ft= t—l +=_F*_o1TL*T*.Dt \2) D(x*)2 LSince all solutions to the original equations remain as solutions to the rescaled equations regardless of the scalings we may choose so as to eliminate as manycoefficients as possible. For example, choosing:a’p.2 T-—-,v v2= = (A.2)Vx = —,IL cx TI2’7Vleaves us with the system:i DT—-b--— L+?*F*_T*,DES DF=OL* L12*— +c*F*_L*T*Dt* —This rescaled system contains only the three parameters:7* =([/ES =V*TfT = —.TiThe Buckingham Pi theorem provides a more rigorous method to ensure that the largestpossible reduction in the number of parameters has been achieved. This will be the case providedthe new parameters and coordinates are all dimensionless.51Appendix BEigenvalues in 3 DimensionsIn this appendix we develop the relations of Derrick [13] between the properties of a 3x3 matrix,and its eigenvalues. Consider the system= Ax, (BJ)near the fixed point x=O. We will examine the case where A is a 3x3 matrix with real elements.The eigenvalues of the system are given by the roots of the equation+ A + B\ + C = 0, (B.2)whereA = —trace A, (B.3)B = A11 + A22 + A33, (B.4)C —det A. (B.5)are the cofactors of A [29]. For the case where each of these roots are real, we will namethem a,,8. and‘,such that < ‘y. We iow have thatA = —a—/3—-y, (B.6)B = cq3+/3’y+-ya’, (B.7)C’ = —a, (B.8)and we can methodically construct table B.1. Tins table shows the relationship between theproperties of the matrix A, and its possible real eigenvalues.52Similarly, for the case of one real root and two complex roots, let the roots be given by a,± iw. We then have thatA = —(a+23), (B.9)B = 2a/9 + 32 + w2, (B.lo)C = —n’(f32+w), (B.11)and we construct table B.2.One immediate consequence is that a necessary condition for a single eigenvalue to passthrough zero (a bifurcation) is that C = 0. Also, necessary conditions for a Hopf bifurcation(pure imaginary eigenvalue pair) are:B > 0, (B.12).4B = C. (B.13)53sgnA sgnB sgnC sgna’ sgn/3 sgw’ Conditions+ + + ——— none- + + /37>IaI(/3+7), IaH/3+7.+ +—— + a/3>(a+/3)7, (a+/3)>7.+— +— + + Iaj(/3+y)>/37, IaI>!3+7.+— — —— + 7k+/3I>a, Ic+/3I>7.— + + — + + /3+7>Ial, /37>ja/3+7).- + - + + + none—— + 7>Ia+/31, a!3>71a+/31.—— + — + + Iaj(/3+7)>/37, IaI>/3+7.————— + 7Ia+/3a/3, 7>Ia+/3l.Table B.1: Real eigenvalues of a 3x3 matrix, a < /3 < 7, a,,B,7 E JRsgnA sgnB sgnC sgna sgn/3 Conditions+ +— + IaI>2/3, /32+>2Ia/3—— none+ +— +— 2I/3 > a /32 + w2 > 2a1/31.+— +— + IaI>2/3, 21a1/3>/3+w.+—— +—2j/3>a, 2a1/31>/3+— + +— + 2/3>IaI, />2IaI/3.— +— + + none+— a>21/31, /32+’>2a/3I.—— +— + 2/3> al, 21a1/3 > /32 + 2,——— +— a>21/31, 2a1/31>/3+w.Table B.2: Imaginary eigenvalues of a 3x3 matrix. A = a, A2,3= /3 + icy’54Appendix CThe Location of WUo,locHerin we show that the eigenvectors of the linearized system near u0 must be directed into thenegative octant of state space relative to u0, and into the region where all three components ofthe tangent vectors are negative.To analyze the eigenvectors near the fixed point u0 we divide the local state space intoeight regions using the nuliplanes of T, F, and L. Over these regions we superimpose theeight octants of the local state space. The positive octant being the region where all threecomponents of u0 — (T, F, L) are greater than zero, and the negative octant the region whereall three components of u0 — (T F, L) are less than zero. We will take the convention that ‘up’and ‘down’ refer to the positive and negative T directions respectively; ‘left’/’right’ to the Ldirections; and ‘forward’/’back’ to the F directions (see figure 3.1).Figure C.1 shows the intersections of the F and L nullsurfaces with planes of constant T. Ifthe eigenvectors are to point into the positive octant relative to u0, then they must also pointinto the region where T’, F’, and L’ are positive. To do otherwise would lead immediately toa contradiction. Similar arguments lead to restrictions on the possibility of the eigenvectorspointing into the other local octants. Imposing these restrictions on the L components of theeigenvectors, it follows that they must point into the shaded regions of figure C.1 (b). Figure C.1(a) shows that F’ > 0 immediately above u0, and thus the eigenvectors cannot point into theshaded region behind and above u0. Thus, this region is no longer shaded in figure C.1 (a).Similarly, from figure C.1 (c) it follows that the eigenvectors cannot point ahead and below u.Thus this region is no longer shaded in figure C.1 (c). If the eigenvectors are to point into aregion above u0, then they must also point into T+. Thus the eigeiivectors cannot point into55the region ahead of u0 and below T°. Similarly, they cannot point into the region behind u0and above T°. This restricts the eigenvectors to point into the negative octant of state spacerelative to u0 intersected with the region where all three components of the tangent vectors arenegative.56Figure C.1: Sections of constant T showing the location of the eigenspace of u0. The sectionsare for values of T where: a) T > T0, b) T = T0, and c) T < T0.The shaded regions indicatethe possible locations of the eigenvectors uo.L LF0FF°(a) T > T0L0FOF FL0(b) T=T0 (C) T <T057Appendix DDetails of the Computer SimulationThe computer simulation was developed to simulate the motion of the ants on a rectangulargrid of arbitrary size. Data stacks were used to keep track of the positions and velocities of theants, and of the positions and strengths of the trails. The velocity of the ants is composed of aconstant speed, and a direction which is an integer multiple of 450• The program begins at thebottom of the ant stack and performs the following algorithm for each ant:1. The affinity’ is computed according to the following formula:STRENGTH *low+ SATURATION (D.1)where STRENGTH is the strength of the trail at the ant’s location, and SATURATION,low’ and are the parameters governing the detection mechanism.2. A random number in the range [0,255] is generated and compared to the parameter .If the random number is larger than , then the ant is assumed to be lost (step 4).Otherwise, the ant attempts to follow the trail (step 3).3. The point directly ahea.d of the ant is tested for the presence of trail pheromone. If thetest is successful, the ant will move to that location. If the test is negative (no trail) theant will check for trails to the left and right, and eitheri) turn (either incrementing of decrementing2 the velocity by one) and move to thestronger trail if the strengths differ,ii) or, move randomly (step 4) if the strengths are equal.‘An integer in the range [0,255]2These calculations are performed modulo eight584. A random number is generated, and compared to the TURNING KERNEL enteredat the start of the run. This will yield a random velocity which is added to that of theants (addition is performed modulo eight). The ant then moves one step in this directionand its new position and velocity are stored in the ant stack.5. After the ant is moved, the trail stack is searched for the new position. If there is a trailthere, the pheromone level is increased by r, otherwise, a new trial is added to the stack,and given a strength T.After processing all the ants in the ant stack, the program decrements the strength of each trailin the trail stack by one, removing any trails of zero strength.The turning kernel governs the relative bias that the lost ant shows towards turning a certainangle. The kernel is assumed to have no external bias, hence it depends oniy on the differencesbetween the ants current direction and the new direction. Further, the kernel is assumed to besymmetric with respect to left and right turns. The forward persistence, is an indkatiom ofthe forward bias of the kernel, and is computed by= jcosB()d. (D.2)B(i,b) is the relative probability of turning through an angle ‘, and is subject to the constraintsB(—’b) = B(b), (D.3)= 1. (D.4)59List of NotationsThe rate at which lost ants follow trails per encounter 19Eigenvalue of A 50/3 The probability of an ant reversing when it moves on or off of a trail ... 19/3 Eigenvalue of A 507 The rate at which the trail pheremone decays 19nondimensionalized version 23Eigenvalue of A 50The rate at which followers lose trails 19nondimensionalized version 23A The mean rate of reversal of the lost ants -. 20Eigenvalues of the fixed point u1 29The region of state space where F’ < 0, and L’ < 0 31A+ The region of state space where F’ > 0, and L’ > 0 31The motility coefficient of the lost ants 22111 The ‘floor’ of the first octant of state space 31112 The ‘back wall’ of the first octant of state space 31p The rate at which followers turn away from the nest 19The rate at which followers turn towards the nest 19r The ratio of follower deposition rate to that of the lost ants 23Tj The rate at which followers deposit trails 19Ti The rate at which lost ants deposit trails 19The affinity of the ants to the trails 58low Lower bound of 5860Turn angle of lost ant 58Forward persistence of the lost ants 58Eigenvectors of the fixed point u1 29w Complex portion of eigenvalue of A 50A 3x3 matrix 50Cofactors of A 50A—trace A 50B Intersection of A with T° 31B ZA, 50B(I’) The turning kernel governing the random motion of the lost ants 58C—detA 50c The speed of the travelling wave 23E The intersection of W’ with B 31F(x,t) Total density of followers 18F(x,t) Density of followers moving towards the nest 18F+(x, t) Density of followers moving away from the nest 18F° F’ nulisurface 30F Region of state space where F’ < 0 30Region of state space where F’ > 0 30G(T, F, L) Growth function of the generalized system 40L(x, t) Total density of lost ants 18L(x,t) Density of lost ants moving towards the nest 18L+(x, t) Density of lost ants moving away from the nest 18rO /L nuilsurface 30L Region of state space where L’ < 0 30L+ Region of state space where L’> 0 30The first octant of state space 3161R(T, F, L) Reaction function of the generalized system 40s Speed of lost ants 18T(x, t) Total length of trails per unit length of network 18time 18time in the moving coordinate system 23T’ nulisurface 30T Region of state space where T’ < 0 30T+ Region of state space where T’ > 0 30uo The populated steady state 26The trivial steady state 26v Speed of followers 18W5’ The stable manifold of the populated steady state 29W0” The unstable manifold of the populated steady state 29Wj9 The stable manifold of the trivial steady state 29W The unstable manifold of the trivial steady state 29x Distance away from the nest 18z Distance from the wavefront 2362


Citation Scheme:


Citations by CSL (citeproc-js)

Usage Statistics



Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            async >
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:


Related Items