system dynamics ( sd ) is a computer simulation problem - solving approach with a foundation in concept of system feedbacks with the purpose of gaining insight into real - world system behavior ( forrester , 1969 ) .
this approach was first introduced by jay forrester in the mid-1950s as a framework for the conceptual representation as well as the quantitative modeling of economic systems ( radzicki and taylor , 1997 ) . 
the paradigm of system dynamics itself assumes that things are interconnected in complex patterns ; that the world is made up of stocks , flows , and feedback loops ; that information flows are intrinsically different from physical flows ; that nonlinear processes and delays are important elements in systems ; and that behavior arises out of system structure  ( meadows , 1989 , p16 ) . besides economical and business management problems ,
a relatively large number of sd models have already been applied to study environmental processes , especially in the field of hydrology and water resource management ( e.g. ahmad and simonovic , 2000 ; gastlum et  al . , 2010 ; khan et  al . , 2009 ) .
the main purpose of these models is to provide insight into non - linear system behavior for the ultimate purpose of assisting informed decision making by stakeholders and policymakers .
however , sd was originally developed as an approach for modeling non - spatial systems .
efforts to integrate spatial modeling capacities into sd resulted in the development of software as for instance spatial modeling environment ( sme ) pioneered by maxwell and costanza ( 1997 ) or simile presented in muetzelfeldt and massheder ( 2003 ) .
mazzoleni et  al . , 2006 ) were proposed to address shortcomings of sd in modeling spatial processes .
concepts as well as the availability of mature software products enabled the application of ssd to a number of case studies in different fields such as ecosystem modeling , hydrology or invasive species control ( e.g. ahmad and simonovic , 2004 ; bendor and metcalf , 2006 ; voinov et  al . , 2004 ) .
the general principle of this model type is that elements spread from areas with high concentration to areas with lower concentrations .
in contrast , the simulation of evolving spatial structures has been widely neglected so far in ssd . the term structure may be defined as  the arrangement of and relations between the parts or elements of something complex  ( oxford dictionaries , 2014 ) .
a structural change model deals with the temporal evolution of spatial structures and the processes driving this evolution .
an example of this type of model is a land cover model ;  land cover is defined by the attributes of the earth 's land surface and immediate subsurface  ( lambin et  al .
natural and anthropogenic processes modify the spatial arrangement of these attributes over time . the spatial arrangement in turn
so far , this type of feedback system has been implemented by few modelers in the area of ssd research .
furthermore , models which are currently available either represent structures as relational space ( e.g. bendor and kaza , 2012 )  which omits absolute positions of objects  or neglects feedback between processes and structures ( e.g. lauf et  al . , 2012 ) .
therefore , the objectives of this article are : ( i ) to introduce an approach which addresses feedbacks between processes and dynamic spatial structures in ssd models ; ( ii ) to initiate a discussion on data exchange between raster or vector - based spatial objects and traditional sd models ; ( iii ) to describe the implementation of a general - purpose python - based software coupling tool ; and ( iv ) to illustrate the impact of process - structure feedbacks based off of a well - established theoretical sd model called daisyworld .
daisyworld is a fictional system introduced by watson and lovelock ( 1983 ) , which serves as a simple parable to show effects of interactions between life and its environment ( wood et  al . , 2008 ) . in this example , the close linkage of the biota to the atmospheric environment causes self - regulation . in order to determine the effects of structural feedbacks on the system 's behavior , simulation results of two different daisyworld versions with and without spatially
a more detailed explanation of the daisyworld model and its basic idea is given in section 4 . in the following sections established sd and ssd concepts are explained and a selection of current ssd models in literature is discussed .
this is followed by a methodological discussion on the technical implementation and the nature of structural change models in ssd .
subsequently , the spatial and non - spatial daisyworld simulations are explained in detail and the simulation results are presented and discussed .
system dynamics ( sd ) combines mathematics and computer simulation to explore the behavior of real - world systems , relationships and processes over time .
the sd approach is effective for formalizing system structure , providing a better understanding of what drives system behavior and examining future dynamics based on a given set of assumptions . in this way ,
modeling system behavior over time can provide insight into significant relationships , reveal patterns , expose sources of undesirable system behavior , and help avoid unforeseen consequences of future policy implementations .
system structures can be conceptualized by means of a causal loop diagram ( see fig . 
arrows indicate relationships between system elements and can expose feedback loops inherent in the system .
the relationship between variables can be characterized by assigning a polarity ( positive or negative ) based on the direction of the change .
for example , a negative polarity indicates that the change of a variable in one direction causes the second variable to change in opposite direction . a positive polarity ,
the behavior of a particular loop is driven by system structure and the polarities of variables in a particular loop .
goal seeking  or  stabilizing  and conversely , positive loops are often considered  diverging from equilibrium  or  destabilizing  .
the conceptual causal loop diagram serves as a blueprint for the implementation of an executable stock and flow model ( see fig . 
this type of model is comprised of stocks ( accumulations or depletions ) , flows ( rates of change ) , delays ( time lags ) and feedbacks ( reciprocal relationships ) . in order to simulate behavior
, the sd model uses mathematics to describe and relate model stocks and variables to each other .
stocks are mathematically represented as integrals , which provide a form of memory to the system .
flows are rates of change expressed as a set of ordinary first - order differential equations .
there is a need for spatial and temporal dynamic simulations , as the non - linear behavior of highly complex systems can not be intuitively forecasted .
real - world systems are highly dynamic in both time and space and unanticipated , even counter - intuitive behavior can often arise in complex systems ( cf .
forrester , 1971 ) , because oftentimes  cause and effect are distant in time and space  ( sanders and sanders , 2004 , p3 ) .
therefore , in order to more accurately simulate behavior of real - world systems , it is necessary to consider both the short and long - term dynamics of system behavior including the impacts and interactions between time and space .
this ultimately assists in assessing the effectiveness of policy decisions and helps identify those actions which would help achieve desired goals .
although sd research is trending in this direction ,  the spatial dimension has not received a great deal of attention in system dynamics modeling .
an intensive literature review showed that there are only a number of articles dealing with this subject 
this is not due to a limitation in the method , but rather due to the nature of sd tools ( kelly et  al . , 2013 ) .
the question of how to represent physical space in a model and how to link this representation to a dynamic sd model is crucial in the field of ssd modeling .
raster or vector data models with geographical information systems ( gis ) are commonly used to represent spatial environments .
the linkage of sd variables to a spatial gis data model includes both data association as well as semantic association ( zhang , 2008 ) .
data association refers to technical coupling and the way of managing the mutual exchange of data between sd and gis .
( 2007 ) distinguish between loose , moderate and tight coupling approaches for the connection of two standalone computer applications by data transfer . in a loosely coupled ssd model ,
, the gis may be used to preprocess spatial inputs , which are then passed to a sd model ; thereafter , execution results are returned to the gis for visualization and additional analysis .
a moderately coupled ssd model enables an indirect communication of systems by means of shared access to a database .
finally , tight coupled applications are characterized by the synchronous operation of systems , allowing direct inter - system communication during program execution .
this is achieved by invoking commands from both sd and gis in a single script .
( 2005 ) suggest the term  embedded  to classify solutions which fully integrate modeling functionality in one independent piece of software .
since variables exhibit different manifestations within sd and gis ( zhang , 2008 ) , a semantic association also needs to be defined for the ssd model .
for instance , raster cells or vector features ( point , line or area objects ) of the gis data model may be linked to sd stocks . in this way , dynamically changing values of the stock element in sd can be transferred to the spatially decomposed simulation environment in gis and vice versa . due to structural similarities , two - dimensional stock arrays are often associated to gis raster ( e.g. bendor and metcalf , 2006 ; ahmad and simonovic , 2004 ) .
one significant advantage is the spatial topology ( relative position of array elements ) inherent in this data type .
ssd applications based on vector data models are usually lacking information on the relative position of object features .
in addition to structural peculiarities , the type of spatial representations used in ssd models strongly depends on the processes being modeled .
for instance , the continuous raster view tends to work best in describing physical quantities ( goodchild , 2005 ) .
huggett ( 1993 ) , for instance , introduced a ssd model of nitrogen accumulation based on terrain data represented by means of a spreadsheet raster . in this way
, the spreadsheet provides the sd model with terrain information being used for the derivation of nitrogen flow directions .
the relative positions of raster elements ( neighborhood conditions ) , which is required to determine nitrogen flow directions , is defined in the stock array .
the link to the spreadsheet is provided as a built - in function of the sd software and allows for spatial results to be displayed in two - dimensions .
a similar model was implemented by bendor and metcalf ( 2006 ) for the simulation of invasive species spread using the software sme . in sme models
are implemented as a three - part view - modelbase - driver architecture ( maxwell and costanza , 1997 ) .
models are graphically constructed in the view component using a graphical modeling environment ( e.g. stella or vensim ) and archived in the modelbase as reusable modules .
the driver is the actual simulation environment which links modules as a set of distributed c++ objects . in this way
sd models are replicated to create a spatial grid covering the study area ( costanza and voinov , 2004 ) .
in contrast to the model introduced by huggett ( 1993 ) , more complex neighborhood conditions are considered .
the propagation of bark beetles is modeled as a function of beetle density in an origin raster cell and the distance to a destination raster cell .
moreover , the migration of beetles to adjacent cells is restricted by a maximum migration distance which is defined prior to the simulation .
conventional gis is only used in this model to initially parameterize spatial arrays and to visualize spatial data in a raster format .
in contrast to ssd models with spatial extensions , ssd models based on gis coupling benefit from a comprehensive set of spatial modeling and data management functionality.ahmad and simonovic ( 2004 ) , for instance , coupled sd to gis for the simulation of damages caused by floods in the red river basin ( canada ) .
similar to the models presented so far , diffusion processes were modeled based on raster data .
sd captures spatio - temporal dynamics of discharge depths , which are sent to a corresponding gis raster .
damages to the infrastructure affected by the modeled flood is reported and visualized in gis .
the data association was achieved by passing data from sd to gis through a spreadsheet , using the dynamic data exchange ( dde ) protocol .
however , the authors concluded that this coupling approach does n't enable a fully automated update of results in gis .
zhang ( 2008 ) overcame this problem by the integration of sd and gis functionality in a universal development environment ( ude ) which serves as middleware to tightly couple these systems .
this model was applied to the spatio - temporal simulation of contaminant propagation in the shoghua river , china .
the process is modeled on a raster - basis in sd and visualized in a customized component gis .
in contrast to ahmad and simonovic 's flood model , a one - dimensional approach was selected for simulation of the pollutant concentration along the river course .
compared to raster - based simulations introduced so far , vector approaches are underrepresented in the field of ssd modeling.https://docs.google.com/a/stud.sbg.ac.at/document/d/1rnccnprajvzg5ytt6xchs359lsf-rbtutdw2dutmx6a/edit 
 lowry and taylor ( 2009 ) , for instance , proposed a technical solution for coupling sd to google earth , which enables the update of kml ( vector format ) color attributes based on sd variables .
similar to the model on pollution propagation by zhang , the systems were tightly coupled by invoking commands in a separate script .
the definition of relative position of vector polygons was neglected , since no interaction between vector polygons was foreseen .
ssd applications are predominantly raster - based approaches which apply a one - to - one raster cell - to - stock ratio . moreover , ssd models typically fail to capture two - way interactions between sd and gis .
communication between sd and gis software is usually unidirectional ; with gis conventionally used as a visualization tool for sd outputs .
an application of ssd for modeling evolving spatial structures requires concepts and techniques which partially differ from those used in diffusion models .
the following sections constitute a first attempt to define the structural change model type and to describe a general framework for an implementation .
a change of spatial structures concerns a qualitative change of spatial objects as , for example , a change from one land cover type to another .
processes interact with space and lead to observable patterns and structures in space ( cf .
structures determine how processes come into effect as , for instance , studied in landscape ecology ( turner , 1989 ) .
this implies that there exists a feedback loop between processes and structures which has to be considered in structural change models . in sd ,
feedback elements are non - spatial and expressed numerically . in the case of the structural change model however , landscape topologies are involved as a spatial feedback element .
the characteristics introduced can be illustrated by means of a simple grassland farming example . in this example , cultivation activities such as harvesting , transport and storage of forage constitute processes .
the intensity of the ongoing farming process can be quantified by the forage gained per area .
this intensity affects the infrastructure , since maintenance of farm lanes depends on their usefulness .
the quality of infrastructure , however , also affects the intensity of farming as it can be assumed that well exploited areas are farmed more intensively and vice versa .
this system can be formalized in a simple non - spatial sd model ( see fig .  2 ) .
causal dependencies are represented by means of proportionality factors ( factors p1 , p2 and p3 ) . in this example
, the quality of infrastructure deteriorates proportional to itself , whereas rate of maintenance is proportional to the level of farming intensity .
farming intensity is dependent on the level of infrastructure quality as a function of the rate of intensity change . in the presented realization of this system ,
the infrastructure is rated based on a single quantitative measure . in a spatial model ,
however , infrastructure may actually be a collection of various different infrastructure types ( e.g. roads , bridges , buildings , other ) .
in addition , the infrastructure may be further described by defining characteristics ( e.g. if the infrastructure is roads , they may be further categorized as paved or unpaved ) .
this enables a derivation of travel times from an infrastructure map , which can be returned as variable to update farming processes .
this two - way interaction provides improved realism of structural system feedbacks and simulation results .
the way structures and processes are integrated in this model type is explained in section 3.2 .
semantic association in ssd concerns the way non - spatial sd variables are linked to the spatial environment being modeled . in the presented approach ,
these structural elements are related to the stock values and updated based on spatial allocation rules .
for many applications , the growth of an object is validly described by assuming least cost or gravitational rules .
for instance , a reduction in object size may be based on the principle of cost optimization .
the example mentioned in section 3.1 related to agricultural infrastructures similarly aims at coverage which most efficiently services the most profitable areas .
raster cells are assigned to stocks and ownership is exchanged between them throughout the simulation .
thus , the proposed method breaks the common one - to - one raster cell - to - stock ratio of established ssd applications ( see also section 2.2 ) .
, a stock could also be related to multiple , non - adjacent zones like land use patches of the same type ( see fig . 
the relationship of regions or zones of different type and their interactions are defined in the stock and flow diagram .
alternatively , clustered multi - raster objects could be replaced by single vector objects . in this way , it is possible to link a more sophisticated view of spatial structures to sd models .
once the stock variable has a spatial counterpart , spatial implications such as the resistance to further growth can be derived from georeferenced data layers .
dynamically changing spatial properties are analyzed and transferred back as variables to the sd model .
this link between evolving objects at different points in time acts as a spatio - temporal ( structural ) feedback loop in the system .
the technical implementation of the bidirectional interaction between process and structure is described in the following section . in order to lay the foundation for managing spatial complexity in sd ,
the following section discusses the approach taken to develop a custom coupling program and the conversion of time - continuous quantities to spatial representations and vice versa .
for instance , the stellar application developed by naimi and voinov ( 2012 ) converts stella sd models into r scripts , which may be useful to modelers who desire more sophisticated statistical and spatial analysis for their sd models .
in contrast to stellar , the proposed approach aims at establishing a link between sd and conventional gis by accessing software libraries and pre - existing program functions available from both sd and gis software .
vensim software provides the sd modeling and simulation capabilities and arcgis software introduces the spatial element to the model .
despite the fact that both vensim and arcgis are proprietary software packages , they are also mainstream programs used widely in the sd and gis communities .
the coupling tool was developed using arcgis version 10.1 ; vensim dss version 5.10 ; and python version 2.7.2 .
the custom coupling tool uses dlls ( .dll ) , python libraries ( .lib ) , packaged vensim models ( .vpm ) , and raster - based gis databases ( .mdb ) to create an environment for dynamically linking vensim model simulation and gis .
vensim dynamic library links ( dlls ) are used to access vensim models and execute commands via the vensim dll to give control over vensim model settings , variables and simulation parameters .
the coupling tool operates differently than a regular vensim simulation ; it does not just execute a pre - determined simulation sequence ; the coupling tool actually uses vensim 's game command to provide control over model simulation .
first , the coupling program initializes libraries , dlls , save drives , and other necessary start functions ; loads a published model file ( .vpm ) and loads a gis raster file ( .mdb ) .
the user is then asked to define simulation settings , following which , pertinent spatial data is extracted from raster files and sent to variables in the vensim simulation model .
simulation then begins by executing a single time step in the sd simulation model . upon completion of a time step
, the simulation data is used as input into the spatial analysis component of the program .
spatial analysis is completed for a single time step and results of spatial analysis are sent back to the sd model .
the same simulation process continues until the specified final time is reached and the simulation ends . during a simulation , the coupling program tracks simulation progress and calls the vensim dll to sketch temporal results as graphs in real - time and uses arcgis to generate a semi - continuous sequence of maps .
the python language is one that is freely available and quite popular in the programming community .
the logic behind selecting python as the sole programming language is for code reusability , consistency and program transferability . despite the fact that vensim software has java , c++ and visual basic support tools , python was selected to develop the custom coupling program for reasons of software compatibility , cross platform functionality , open source development and intended future use .
the python packages imported into the custom coupling program include : numpy ( for scientific computations and array operations ) ; arcpy ( for getting access to arcgis functionality ) ; matplotlib ( for visualizing simulation output ) ; and ctypes ( a foreign function library for compatibility with c data types ) .
the custom coupling tool has been developed using the numpy package to replace the equivalent spatial analysis functions of arcpy .
thus , arcpy provides spatial references as well as the interface to the arcgis database whereas the spatial analysis is conducted in numpy . in the further course of the paper the spatial component of the model
one of the unique elements of the approach and specific application in this paper is the synchronous data flow between sd and gis .
to effectively have sd simulate system processes and gis to perform spatial analysis requires bidirectional synchronous operation between sd and gis software .
is preserved between time steps while spatial analysis is performed . in effect , generating snapshots of changes in spatial structure over time .
this differs from the simulation of diffusion models and diffusion processes , where no interruption of the sd model simulation is required , because diffusion models do not require data to be returned from gis to the sd simulation .
this coupling approach could be used for any number of potential ssd applications which involve feedback between process and space .
the program gives access to a comprehensive set of spatial analysis and dynamic simulation functions and coordinates exchange between them .
moreover , any packaged vensim model can be loaded by simply adding a model directory to the code .
the code structure serves as a framework ( see fig .  5 ) , which can be retained , while code lines referring to variable names , spatial data or spatial analysis functions may be customized to the needs of the model .
an additional strength of the proposed coupling tool is that the sd model is not converted into code .
this enables users to modify the stock and flow structure after loading the model into the middleware .
furthermore , non - spatial results can be visualized in vensim as tables and diagrams . at the same time
the embedding of sd models in a high level and freely available programming language facilitates an extension of model functionality by linking to additional libraries .
hence , even though a graphical user interface would have improved the usability , the fact that such systems are closed becomes restricting in case of more complex , spatially explicit models ( naimi and voinov , 2012 ) . in the following section , daisyworld is described and expanded to serve as first use case for this approach .
daisyworld is a fictitious environment introduced by watson and lovelock ( 1983 ) to describe biological homeostasis of the global environment .
the term homeostasis originated in the field of physiology and refers to  the ability of living beings to maintain their own stability 
the intention of the original daisyworld model was to illustrate that homeostasis could also be observed in environmental systems as described by the gaia hypothesis ( lovelock and margulis , 1974 ) . in the daisyworld model
a fictitious planet accommodates three land cover types classified as black daisies , white daisies and fertile soil ( known as  bare ground  in the original daisyworld model ) .
a distant sun emits luminosity and warms up the planet 's surface . due to the difference in albedo between black and white daisies ,
black daisies cause a local heating , whereas temperature next to areas vegetated by white daisies falls below average planetary temperatures .
since both types of daisies reach a growth optimum of 100% at a temperature of 22.5  c , localized heating leads to beneficial conditions for black daisies in a daisyworld with average planetary temperature below the optimum temperature .
conversely , local cooling favors white daisies in an environment where temperature is above this threshold .
this mechanism counteracts changes of the sun 's luminosity by means of two negative feedback loops .
increasing luminosity leads to higher temperature which favors the growth of white daisies ; which , in turn , decreases the temperature in daisyworld due to a lower planetary albedo .
the same principle applies to decreasing luminosity which causes a spread of black daisies . in this way
, type of land cover acts as a stabilizing component within the system and can be considered a driving factor of homeostatic control in daisyworld .
the albedo of the planet is calculated as an area - weighted average of the albedos from the three surfaces ( see table  1 , equation ( 1 ) ) .
planetary temperatures are based on the stefan boltzmann law , which determines temperatures as a fraction of the absorbed energy input . for the sake of simplicity ,
the planet is treated as a planar surface with uniform energy input ( see table  1 , equation ( 2 ) ) .
local temperatures in the vicinity of black and white daisies are expressed by the difference in daisy and planetary albedos ( see table  1 , equation ( 3 ) ) .
daisies grow as a parabolic function respective to their local temperature ( see table  1 , equation ( 4 ) ) .
moreover , growth is affected by the globally available fertile soil patches not covered by daisies .
in addition , it is assumed that a fraction of daisies die off each time step ( see table  1 , equation ( 5 ) ) .
ford ( 2010 ) implemented daisyworld as a non - spatial , dynamic model in sd .
the proposed ssd daisyworld model fulfills the aforementioned criteria of process - structure interaction defined for structural change models ( see section 3.1 ) and may therefore be considered an appropriate example in demonstrating spatio - temporal effects for this model type . in order to consider a spatial daisyworld model , there were two significant modifications made to the original model : ( i ) spatialization of daisyworld ; and ( ii ) the introduction of  barren land  as a fourth land cover type .
the main purpose of introducing this additional land cover type is to emphasize spatial modeling mechanics and its potential implications for sd modeling .
the barren land cover type is considered to be immutable and unproductive land in the daisyworld simulation environment .
to accommodate the new cover type , a constant was introduced into the daisyworld sd model with the assumption that barren land has an albedo value of 0.5 .
the initial values and settings used in the spatial daisyworld model are presented in table  2 , according to the numeration seen in fig . 
the area of different land cover types such as the area of black daisies ( see fig .  6 ; no .
11 ) are represented by stocks in the sd model . in order to initialize the model , stock values
this is implemented in such a way that the areas of respective land cover types in the spatial raster correspond to the initial stock values .
the growth of daisy populations in sd requires updating in the spatial representation based on a set of allocation rules .
therefore , in a first step fertile soil cells sharing a border to a daisy cell of the respective species are selected . from this pool of candidates
this also implies that vegetated areas grow laterally , provided that sufficient fertile soil is available .
alternatively , this simple allocation rule may be extended in the future by enabling growth within a distance buffer to daisies of the same species .
in addition to the temperature dependent growth of daisies , there is a reduction in daisy area as they die .
this acts as a negative feedback loop in the sd model which triggers an increase in decay as the daisy area increases ( see fig .  6 , no . 12 and 14 ) .
in order to spatially allocate daisy decay , a timestamp is provided to every raster cell that indicates the age of daisy growth .
the age of the daisies is used as criteria in determining which raster cell subsequently experiences daisy death .
the raster cell with the oldest daisies are sequentially selected and updated to reduce the amount of daisy growth .
when a raster cell 's daisies die , the cell returns to fertile soil type .
due to the dynamic evolution of the landscape ( see fig .  7 ) , the opportunities for a species to spread changes as a function of soil availability .
the growth reduction control in the original version assumes that the remaining soil is equally shared between the two species .
the scarcity of soil is imposed as a globally defined factor by reducing the growth rate of both species as a linear function of fertile soil ( see table  1 , equation ( 5 ) ) .
for instance , if 50% of daisyworld is fertile soil , growth rates of both species are multiplied by 0.5 . in a spatial daisyworld
however , the distribution of soil resources is a matter of locally available fertile soil patches .
thus , the original model was supplemented with an individual growth reduction multiplier for each type of daisy ( see fig .  6 ; no .
the edge length between landscape elements is important for the growth of species ( turner , 1989 ) .
this is considered by counting the number of shared borders between daisy cells and fertile soil cells . in theory
the proportion of actual adjacencies to the maximum number of adjacencies is used as a proxy of available expansion areas for each daisy species . in order to sequentially search the raster for adjacencies , a 3  x  3 moving window is used .
this procedure is repeated each time the daisyworld landscape changes to update the sd model with new information . in order to take these modifications into consideration , equations ( 1 ) and ( 5 ) ( see table  1 )
were adapted as follows:(6)a=gag+bab+waw+barrabarr(7)b , w/t = xb , w(db , wb , w ) the planetary albedo was calculated including the area and albedo of barren land .
furthermore , the area of fertile soil in equation ( 5 ) was replaced by the growth reduction multiplier for each species ( b , w ) calculated by(8)db , w = gb , w / cb , w8where gb , w is the number of fertile soil raster cells next to black or white daisies and cb , w is the number of raster cells occupied by black or white daisies .
a spatial and a non - spatial version of daisyworld are simulated and compared for various scenarios to verify the impact of structural feedback on system behavior .
additionally , the performance and functionality of the proposed coupling approach is illustrated . in order to make the two versions comparable ,
the non - spatial model is also supplemented with  barren land  as an additional land type .
for all scenarios , a change in luminosity was assumed which is balanced out by means of land cover adaptation in the system . since the initial composition of the landscape effects the simulation results ,
results generated for these landscapes were compared amongst each other and to the non - spatial simulation results .
furthermore , a stress test scenario was conducted by assuming disproportional black and white daisy areas .
this scenario was intended to compare the resilience of spatial and non - spatial models under unfavorable conditions .
in addition , luminosity was gradually increased over a longer time period in order to observe the continuous adaptation of both systems to a changing environment .
both spatial and non - spatial daisyworld simulations were run for 100 time steps ( see fig . 
the land cover is composed of 100  ha of each black and white daisies , 500  ha of fertile soil and 300  ha of barren land . due to the initial random distribution of land cover types in the spatial model , multiple simulations were completed to identify effects of land cover distributions on the simulation results .
stable conditions are achieved after approximately 15 time steps since the onset of initial cooling .
this cooling phase is much more pronounced in the non - spatial simulation than in the spatial simulation .
whereas the non - spatial system is able to raise temperatures in this low radiation environment , the spatial system fails to completely resolve the outside disturbance to the system .
the variation in the results for five randomly generated landscapes , however , reveals considerable dependencies on land type distribution .
nevertheless , the range of spatial results clearly differs from the non - spatial simulation .
the main reason for these differences is a result of the quick response of daisy growth to changing luminosity and temperature in the non - spatial simulation which causes more pronounced temperature compensations .
the scenario assumes a land composition of 190  ha of black daisies , 10  ha of white daisies , 500  ha of fertile soil and 300  ha of barren land .
it can be assumed that these initial conditions cause an increase in black daisies and a decay of white daisies as a result of decreasing luminosity .
since white daisies are initialized with such a small area , the response of daisy growth to changing environmental conditions endangers the remaining areas .
initially , growth rates of white daisies are significantly higher in the non - spatial simulation than in the spatial simulation ( see fig . 
9 ) . however , the negative feedback between areas of white daisies and daisy decay causes higher decay rates as a result of decreasing luminosity after ten time steps in the non - spatial simulation .
due to the slower decay of white daisies , this land cover can be sustained on a low level in the spatial daisyworld .
therefore , although the non - spatial system is able to raise temperatures to a higher level , the spatial system is more resilient in the stress test scenario .
the main reason for this behavior may be attributed to the different daisy growth rates , which , in turn , are a result of the way available space for daisy propagation is defined in each of the respective models . in the non - spatial daisyworld simulation ,
propagation of daisies does n't require spatial proximity to fertile soil , and accessibility of fertile soil is explicitly considered in the spatial model by means of a growth reduction multiplier .
this modification hinders daisy growth and crucially affects the response of the system to changing luminosity .
for instance , gradually increasing luminosity results in a period of cooling in the non - spatial system , since growth rates of white daisies unravel and even reverse effects of luminosity on temperatures ( see fig . 
the spatial simulation however , reveals less of an impact of land cover adaptation on temperature reduction .
this is due to the restriction of suitable soils available for white daisies in the spatial simulation .
the predefined growth restriction , allowing growth only next to daisies of the same species , leads to a local clustering of daisies .
fertile soil which may be available elsewhere on daisyworld can not be exploited if spatial barriers prevent daisy clusters from further expansion .
land cover adaptation has a diminishing effect on temperatures , but there is no effective control on temperature when daisyworld is simulated using a more realistic , spatial algorithm which considers structural feedbacks .
the raster resolution in the spatial model was varied in order to illustrate effects on the model output .
an increase in spatial resolution can be associated with decreasing adaptive capacities of the system ( see fig . 
thus , results do not only depend on the composition of the landscape structure as shown in fig .  8 , but also on the way the landscape is represented in the data model .
the models ' sensitivity to the spatial resolution is a consequence of results from the raster data model and model calculus .
consider a cluster of daisies to be composed of multiple daisy - type raster cells .
there are daisy cells in the center of the daisy cluster , which are surrounded by other daisy cells .
there are also daisy cells at the edge of the cluster that have a fraction of cell neighbors which are daisies and a fraction of cell neighbors which are fertile soil .
a change in spatial resolution will modify the proportion of daisy raster cells on the edge of the cluster and the number of daisy cells at the center of the cluster , which biases the modeling results . as it turns out
, if spatial resolution is increased there tends to be a more pronounced daisy growth reduction .
the sensitivity analysis was also used to test the computational performance of the coupling program .
an increase in the spatial resolution partially leads to an overproportional rise of the processing time .
the bisection of the resolution from 100 to 50  m  almost doubles the processing time from less than 10  s to about 17  s. if higher resolutions are used , the increase in processing times exhibits more of an exponential character .
a simulation based on a 25  m grid takes about 1  min ; more than 12  min  are required to run the simulation with a resolution of 12.5  m. the increase in processing times results from spatial analyses which are performed on a raster basis . increasing the time horizon (
i.e. the number of time steps ) of a simulation , reflects an increase in the number of interactions required between sd and gis . however
consequently , shortcomings in performance are not caused by the interaction between sd and gis , but instead depend on numpy raster analysis functionality invoked on the gis - side of the program .
furthermore , problems associated with the relatively limited number of stock variables in most sd products , as stated in section 2.2 , are circumvented by the proposed approach . due to the aggregation of raster cells
this makes the structural change model type favorable for larger study areas . even though the number of links is smaller for this model type , interactions are bidirectional .
information exchange and update is not only for the purpose of visualizing spatial results , but introduces spatial ( structural ) feedback to the model .
therefore , modeling structural change in ssd makes higher claims against the quality and performance of links between sd and gis .
thus , a link based on data in a script is preferable over an interaction via files . the usage of georeferenced raster data types provided with arcpy in combination with multi - dimensional numpy data types proved suitable for this task .
next to coupling requirements , the need for a synchronous operation of functions ( see fig . 
the ability of the vensim dll to pause and continue the simulation at predefined points in time enabled an automation of this otherwise quite cumbersome task .
the comparison of the performance of simulations with and without vensim interruptions for 100 time steps revealed negligible differences in processing times ( less than one second ) .
this implies that the context switch does n't produce a large overhead and indicates an overall reasonable performance of the proposed synchronization of operations .
synchronization is one approach to incorporate dynamic process - structure relations in ssd by integrating iterative spatial analysis into the simulation process .
the approach presented in this paper especially benefits ssd modeling research and projects whose requirements for using gis go well beyond visualization purposes .
practitioners may use this framework for modeling systems which require that spatial structures are dynamic and coupled to feedback processes .
the application of ssd to the simulation of structural change requires a different semantic association of sd variables with the gis data model .
the commonly used one - to - one raster cell - to - stock association may be replaced by a greater than one - to - one multiple raster cells - to - stock ratio .
a stock value represents the area of a structural element such as a specific land cover type .
the area exchanged between these stocks is assigned to the raster data model by a set of predefined allocation rules .
processes which modify spatial structures and structures which in turn affect processes define an important feedback component of this system . as a consequence of this mutual interdependence , interaction between sd and gis needs to be bidirectional .
the tight coupling of the sd software vensim to arcgis via a python middleware enables efficient system interactions during program execution .
furthermore , operations performed in gis and sd need to be synchronized for the structural change model type .
the computational overhead associated with these specific coupling features is virtually negligible , as demonstrated in the example of a fictitious system called daisyworld . in daisyworld , a dynamically adapting landscape of black and white daisies balances out temperature changes caused by varying luminosities .
the comparison of a spatial daisyworld simulation to a non - spatial simulation emphasizes the significance of incorporating structural feedbacks inherent in spatial models .
the initial composition and land type distribution of the daisyworld landscape has a consequence on the assignment of future land cover types and simulation results . moreover , simulation results suggest daisies die to point of extinction in the non - spatial model .
this is likely due to the more pronounced response of land cover in the lumped simulation , which does n't take local growth restrictions ( as an expression of spatial structures ) into account .
locally restricted daisy growth leads to increased temperatures in the case where luminosity is increasing , whereas the pronounced land cover adaptation in the non - spatial version actually causes decreasing temperatures .
thus , the spatial system exhibits higher robustness , whereas the non - spatial version has higher adaptive capacities .
however , results of the spatial model are not only influenced by the initial landscape structure , but also depend strongly on raster resolution .
an increase in spatial resolution goes hand in hand with decreasing adaptive capacities and vice versa .
this is due to changes in the raster topology ( shared borders between daisy and fertile soil cells ) .
nevertheless , results support the idea of utilizing ssd for the simulation of structural change . compared to diffusion process models , less data links between sd and gis are required .
the generic design of the coupling tool lends itself for easy implementation of more realistic ssd models .
in the presented approach objects are generated by assigning consistent numeric identifiers to a discrete number of raster cells .
objects assembled in this way are visually interpreted as objects ; however they are not fully recognized as objects in gis .
therefore , the next step towards a more efficient approach is to partially vectorize objects and assign them distinct attributes .
this modification of the approach would also allow establishing object hierarchies for multi - scale modeling .
hierarchical organization enables simultaneous simulation processes on multiple scales as well as interactions between them .
moreover , biasing effects of different raster resolutions may be prevented by using a vector model instead .
the full integration of gis functions facilitates the visualization , analysis and management of non - linear spatial arrangements in ssd models .
therefore , the fusion of gis and sd functionality , as presented in this paper , opens up new applications for ssd in fields where regular grids and one - directional stock - to - raster associations are not conducive to representing reality .
possible future applications of the ssd approach include simulating land ownership changes and simulating emergency response during disasters .
more specifically , the coupling program is be used for modeling structural changes in grassland agriculture in austria .
this includes simulating competition between different farms to assess ownership changes . in this example , processes at the farm level are implemented as part of a vensim model .
transportation costs are calculated iteratively , as the agricultural areas for each farm changes over time . therefore ,
additionally , an application of the structural change approach is planned for the dynamic simulation of serviceable emergency response areas in response to natural disasters .
street networks are governed by changing traffic during hazards and access to critical infrastructure may become restricted .
thus , emergency response is functioning as part of a dynamic temporal and spatial landscape .
the coupling program aspires to capture these changes and improve disaster planning , response , and recovery efforts .
the proposed approach may be especially suitable for applications in the fields of land use and land cover modeling .
this is due to the clear distinction between processes and structures as well as the restrictions to qualitative changes in space .
the structure change model type may be used as an alternative to commonly used ca and abm approaches in this field .
the open design of the coupling tool facilitates an integration of other approaches to address the various requirements of integrated environmental assessment ( cf .
the current version of middleware was designed as a first prototype , which is generic - form rudimentary coding that invokes commands from both arcgis and vensim .
so far , front end features such as a gui have not been implemented , but may be included in future versions .
 system dynamics ( sd ) is an effective approach for helping reveal the temporal behavior of complex systems . 
 although there have been recent developments in expanding sd to include systems spatial dependencies , most applications have been restricted to the simulation of diffusion processes ; this is especially true for models on structural change ( e.g. lulc modeling ) . to address this shortcoming , 
 a python program is proposed to tightly couple sd software to a geographic information system ( gis ) . 
 the approach provides the required capacities for handling bidirectional and synchronized interactions of operations between sd and gis . in order to illustrate the concept and the techniques proposed for simulating structural changes , 
 a fictitious environment called daisyworld has been recreated in a spatial system dynamics ( ssd ) environment . 
 the comparison of spatial and non - spatial simulations emphasizes the importance of considering spatio - temporal feedbacks . 
 finally , practical applications of structural change models in agriculture and disaster management are proposed . 
Introduction
Spatial system dynamics
Methodology
Daisyworld
Discussion of results
Conclusions
Potential extensions