mtvDisplayMath magenta!43!cyan!10!white

# Software Uncertainty in Integrated Environmental Modelling: the role of
Semantics and Open Science

Daniele de Rigo European Commission, Joint Research Centre, Institute for
Environment and Sustainability
Via E. Fermi 2749, I-21027 Ispra (VA), Italy Politecnico di Milano,
Dipartimento di Elettronica, Informazione e Bioingegneria
Via Ponzio 34/5, I-20133 Milano, Italy

Copyright © 2013 Daniele de Rigo. This work is licensed under a Creative
Commons Attribution 3.0 Unported License
(http://creativecommons.org/licenses/by/3.0/). See:
http://www.egu2013.eu/abstract_management/license_and_copyright.html This is
the author’s version of the work. The definitive version has been published in
the Vol. 15 of Geophysical Research Abstracts (ISSN 1607-7962) and presented
at the European Geosciences Union (EGU) General Assembly 2013, Vienna,
Austria, 07–12 April 2013 http://www.egu2013.eu/ Cite as: de Rigo, D., 2013.
Software Uncertainty in Integrated Environmental Modelling: the role of
Semantics and Open Science. Geophys Res Abstr 15, 13292+ Author’s version DOI:
10.6084/m9.figshare.155701 , arXiv: 1311.4762

Computational aspects increasingly shape environmental sciences [1]. Actually,
transdisciplinary modelling of complex and uncertain environmental systems is
challenging computational science (CS) and also the science-policy interface
[2, 3, 4, 5, 6, 7].

Large spatial-scale problems falling within this category – i.e. wide-scale
transdisciplinary modelling for environment (WSTMe) [8, 9, 10] – often deal
with factors (a) for which deep-uncertainty [2, 11, 7, 12] may prevent usual
statistical analysis of modelled quantities and need different ways for
providing policy-making with science-based support. Here, practical
recommendations are proposed for tempering a peculiar – not infrequently
underestimated – source of uncertainty. Software errors in complex WSTMe may
subtly affect the outcomes with possible consequences even on collective
environmental decision-making. Semantic transparency in CS [2, 8, 10, 13, 14]
and free software [15, 16] are discussed as possible mitigations (b).

Software uncertainty,

black-boxes and free software

Integrated natural resources modelling and management (INRMM) [17] frequently
exploits chains of nontrivial data-transformation models (D-TM), each of them
affected by uncertainties and errors.

Those D-TM chains may be packaged as monolithic specialized models, maybe only
accessible as black-box executables (if accessible at all) [18]. For end-
users, black-boxes merely transform inputs in the final outputs, relying on
classical peer-reviewed publications for describing the internal mechanism.

While software tautologically plays a vital role in CS, it is often neglected
in favour of more theoretical aspects.

(a) Complexity= { Transdisciplinary integration (e.g. systems of systems)
Environmental system(s) heterogeneity (e.g. geospatial fragmentation) Data
heterogeneity (formats, definitions, spatiotemporal density, …) Software
complexity (algorithms, dependencies, languages, interfaces, …) Uncertainty= {
Incomplete scientific knowledge (e.g. climate scenarios [19, 20, 21], tipping
points [22, 23, 24], … ) Modelling assumptions and simplifications [25, 26,
27] Uncertainty of measured/derived data Software uncertainty
Dynamicbehaviour= { Uncertainty propagation via: Propagation in the network of
interconnected WSTMe components [2, 14, 28, 17, 29, 30, 31, 32, 33] Iterations
within nonlinear optimization steps [5, 34, 35, 36, 37, 38, 39, 40] Data
fusion, harmonization, integration [9, 41, 42, 43, 44] Steps for computing and
aggregating criteria and indices [6, 11, 7, 45, 46, 47, 48]

This paradox has been provocatively described as “the invisibility of software
in published science. Almost all published papers required some coding, but
almost none mention software, let alone include or link to source code” [49].

Recently, this primacy of theory over reality [50, 51, 52] has been challenged
by new emerging hybrid approaches [53] and by the growing debate on open
science and scientific knowledge freedom [2, 54, 55, 56, 57].

In particular, the role of free software has been underlined within the
paradigm of reproducible research [18, 56, 57, 58]. In the spectrum of
reproducibility, the free availability of the source code is emphasized [56]
as the first step from non-reproducible research (only based on classic peer-
reviewed publications) toward reproducibility.

Y = f^*( X ) = f( θ^* , X ) Theoretic D-TM whose algorithm is typically
described in peer reviewed publications. The D-TM may e.g. implement a given
WSTMe as instance of a suitable family of functions $f$ by means of selected
parameters $\theta^{*}$. $\theta^{*}$ may be the result of an optimization
(regression, control problem, …). Y = f^$\zeta$ = f( θ^$\zeta$ , X , $\zeta$ )
Real D-TM where the software uncertainty $\zeta$ may affect both the function
family $f$ and the optimality of the selected parameters $\theta^{\zeta}$.
::|f( θ , X , $\zeta$ ) |::^^sem Semantically enhanced D-TM (e.g. SemAP).
The D-TM is subject to the semantic checks $sem$ as pre-, post-conditions and
invariants on inputs, outputs and the D-TM itself: Y = ::|f(θ,X,ζ )
|::^^sem⇔{Y=f(θ,X,ζ )□sem( Y, f, θ,X,ζ )(b)where  {X  is the input array of
data  X = { X_1, X_2, ⋯X_i ⋯X_n }X_i ∈C^N_i1 ×⋯×N_in_i  is a multi-dimensional
array (e.g. a two-dimensional  raster layer)Y  is analogously the output array
of data the modal/deontic logic operator  □p  means: it ought to be that  p .

Applying this paradigm to WSTMe, an alternative strategy to black-boxes would
suggest exposing not only final outputs but also key intermediate layers of
data and information along with the corresponding free software D-TM modules.

“Software errors in complex WSTMe may subtly affect the outcomes with possible
consequences even on collective environmental decision-making” “The chain of
free-software modules should be transparent” A concise, semantically-enhanced
modularization [13, 14] may help not only to see the code (as a very basic
prerequisite for semantic transparency) but also to understand – and correct –
it [59]. Semantically-enhanced, concise modularization is e.g. supported by
semantic array programming (SemAP) [13, 14] and its extension to geospatial
problems [8, 10].

Some WSTMe may surely be classified in the subset of software systems which
“are growing well past the ability of a small group of people to completely
understand the content”, while “data from these systems are often used for
critical decision making” [50]. In this context, the further uncertainty
arising from the unpredicted “(not to say unpredictable)” [51] behaviour of
software errors propagation in WSTMe should be explicitly considered as
software uncertainty [60, 61] (see b). The data and information flow of a
black-box D-TM is often a (hidden) composition of D-TM modules:

This chain of free-software D-TM modules (each of them semantically-enhanced)
should be transparent:

Semantics and design diversity

Silent faults [62] are a critical class of software errors altering
computation output without evident symptoms – such as computation premature
interruption (exceptions, error messages, …), obviously unrealistic results or
computation patterns (e.g. noticeably shorter/longer or endless computations).
As it has been underlined, “many scientific results are corrupted, perhaps
fatally so, by undiscovered mistakes in the software used to calculate and
present those results” [63].

“Semantic modularization might help to catch at least a subset of silent
faults, when misusing intermediate data outside the expected semantic context”
“Where the complexity and scale may lead to deep uncertainty, techniques such
as ensemble modelling may be recommendable” Despite the ubiquity of software
errors [60, 61, 62, 63, 64, 65, 66, 67, 68], the structural role of scientific
software uncertainty seems dramatically underestimated [2, 51]. Semantic D-TM
modularization might help to catch at least a subset of silent faults, when
misusing intermediate data outside the expected semantic context of a given
D-TM module (b). Where the complexity and scale of WSTMe may lead unavoidable
software-uncertainty to induce or worsen deep-uncertainty [2], techniques such
as ensemble modelling may be recommendable [11, 7, 12]. Adapting those
techniques for glancing at the software-uncertainty of a given WSTMe would
imply availability of multiple instances (implementations) of the same
abstract WSTMe. Independently re-implementing the same WSTMe (design diversity
[69]) might of course be extremely expensive. However, partly independent re-
implementations of critical D-TM modules may be more affordable and examples
of comparison between supposedly equivalent D-TM algorithms seem to
corroborate the interest of this research option [57, 70, 49].

## References

  * [1] Casagrandi, R., Guariso, G., 2009. Impact of ICT in environmental sciences: A citation analysis 1990-2007. _Environmental Modelling & Software 24_ (7), 865-871. DOI:10.1016/j.envsoft.2008.11.013 Google Scholar: 10214045160670186637
  * [2] de Rigo, D., (exp.) 2013. Behind the horizon of reproducible integrated environmental modelling at European scale: ethics and practice of scientific knowledge freedom. _F1000 Research_. Submitted
  * [3] Gomes, C. P., 2009. Computational sustainability: Computational methods for a sustainable environment, economy, and society. _The Bridge 39_ (4), 5-13. http://www.nae.edu/File.aspx?id=17673
  * [4] Easterbrook, S. M., Johns, T. C., 2009. Engineering the software for understanding climate change. _Computing in Science & Engineering 11_ (6), 65-74. DOI:10.1109/MCSE.2009.193 Google Scholar: 5281658010101603741
  * [5] Hamarat, C., Kwakkel, J. H., Pruyt, E., 2012. Adaptive robust design under deep uncertainty. Technological Forecasting and Social Change. DOI:10.1016/j.techfore.2012.10.004 Google Scholar: 1225653171705658912
  * [6] Bankes, S. C., 2002. Tools and techniques for developing policies for complex and uncertain systems. _Proceedings of the National Academy of Sciences of the United States of America 99_ (Suppl 3), 7263-7266. DOI:10.1073/pnas.092081399 Google Scholar: 2377027768488561757
  * [7] Kandlikar, M., Risbey, J., Dessai, S., 2005. Representing and communicating deep uncertainty in climate-change assessments. _Comptes Rendus Geoscience 337_ (4), 443-455. DOI:10.1016/j.crte.2004.10.010 Google Scholar: 8596177795590020821
  * [8] de Rigo, D., Corti, P., Caudullo, G., McInerney, D., Di Leo, M., San-Miguel-Ayanz, J., 2013. Toward Open Science at the European scale: Geospatial Semantic Array Programming for Integrated Environmental Modelling. _Geophysical Research Abstracts 15_ , 13245\+ DOI:10.6084/m9.figshare.155703 Google Scholar: 17118262245556811911
  * [9] Rodriguez Aseretto, D., Di Leo, M., de Rigo, D., Corti, P., McInerney, D., Camia, A., San Miguel-Ayanz, J., 2013. Free and Open Source Software underpinning the European Forest Data Centre. _Geophysical Research Abstracts 15_ , 12101\+ DOI:10.6084/m9.figshare.155700 Google Scholar: 14482024956822192435
  * [10] de Rigo, D., Corti, P., Caudullo, G., McInerney, D., Di Leo, M., San-Miguel-Ayanz, J., (exp.) 2013. Supporting Environmental Modelling and Science-Policy Interface at European Scale with Geospatial Semantic Array Programming. In prep.
  * [11] Lempert, R. J., 2002. A new decision sciences for complex systems. _Proceedings of the National Academy of Sciences of the United States of America 99_ (Suppl 3), 7309-7313. DOI:10.1073/pnas.082081699 Google Scholar: 3105059191542042306
  * [12] Gober, P., Kirkwood, C. W., 2010. Vulnerability assessment of climate-induced water shortage in Phoenix. _Proceedings of the National Academy of Sciences 107_ (50), 21295-21299. DOI:10.1073/pnas.0911113107 Google Scholar: 13404393835259418552
  * [13] de Rigo, D., 2012. Semantic Array Programming for Environmental Modelling: Application of the Mastrave library. In: Seppelt, R., Voinov, A. A., Lange, S., Bankamp, D. (Eds.), International Environmental Modelling and Software Society (iEMSs) 2012 International Congress on Environmental Modelling and Software. Managing Resources of a Limited Planet: Pathways and Visions under Uncertainty, Sixth Biennial Meeting. pp. 1167-1176. http://www.iemss.org/iemss2012/proceedings/D3_1_0715_deRigo.pdf Google Scholar: 6628751141895151391
  * [14] de Rigo, D., 2012. Semantic Array Programming with Mastrave - Introduction to Semantic Computational Modelling. http://mastrave.org/doc/MTV-1.012-1 Google Scholar: 6848554969929557252
  * [15] Free Software Foundation, 2012. What is free software? http://www.gnu.org/philosophy/free-sw.html (revision 1.118 archived at http://www.webcitation.org/6DXqCFAN3 ) Google Scholar: 7470103647360812109
  * [16] Stallman, R. M., 2009. Viewpoint: Why "open source" misses the point of free software. _Communications of the ACM 52_ (6), 31-33. DOI:10.1145/1516046.1516058 (free access version: http://www.gnu.org/philosophy/open-source-misses-the-point.html ) Google Scholar: 17751536887456926788
  * [17] de Rigo, D., 2012. Integrated Natural Resources Modelling and Management: minimal redefinition of a known challenge for environmental modelling. Excerpt from the _Call for a shared research agenda toward scientific knowledge freedom_ , Maieutike Research Initiative. http://www.citeulike.org/groupfunc/15400/home
  * [18] Morin, A., Urban, J., Adams, P. D., Foster, I., Sali, A., Baker, D., Sliz, P., 2012. Shining light into black boxes. _Science 336_ (6078), 159-160. DOI:10.1126/science.1218263 Google Scholar: 12575758499484368256
  * [19] Lempert, R., Schlesinger, M. E., 2001. Climate-change strategy needs to be robust. _Nature 412_ (6845), 375. DOI:10.1038/35086617 Google Scholar: 3059012525521592083
  * [20] Shell, K. M., 2012. Constraining cloud feedbacks. _Science 338_ (6108), 755-756. DOI:10.1126/science.1231083 Google Scholar: 17586687489186735372
  * [21] van der Sluijs, J. P., 2012. Uncertainty and dissent in climate risk assessment: A Post-Normal perspective. _Nature and Culture 7_ (2), 174-195. DOI:10.3167/nc.2012.070204 Google Scholar: 2453345655060253463
  * [22] Lenton, T. M., Held, H., Kriegler, E., Hall, J. W., Lucht, W., Rahmstorf, S., Schellnhuber, H. J., 2008. Tipping elements in the earth’s climate system. _Proceedings of the National Academy of Sciences 105_ (6), 1786-1793. DOI:10.1073/pnas.0705414105 Google Scholar: 9014054360036592667
  * [23] Hastings, A., Wysham, D. B., 2010. Regime shifts in ecological systems can occur with no warning. _Ecology Letters 13_ (4), 464-472. DOI:10.1111/j.1461-0248.2010.01439.x Google Scholar: 10306576888461248930
  * [24] Barnosky, A. D., Hadly, E. A., Bascompte, J., Berlow, E. L., Brown, J. H., Fortelius, M., Getz, W. M., Harte, J., Hastings, A., Marquet, P. A., Martinez, N. D., Mooers, A., Roopnarine, P., Vermeij, G., Williams, J. W., Gillespie, R., Kitzes, J., Marshall, C., Matzke, N., Mindell, D. P., Revilla, E., Smith, A. B., 2012. Approaching a state shift in earth’s biosphere. _Nature 486_ (7401), 52-58. DOI:10.1038/nature11018 Google Scholar: 17836601480741683208
  * [25] Milly, P. C. D., Betancourt, J., Falkenmark, M., Hirsch, R. M., Kundzewicz, Z. W., Lettenmaier, D. P., Stouffer, R. J., 2008. Stationarity is dead: Whither water management? _Science 319_ (5863), 573-574. DOI:10.1126/science.1151915 Google Scholar: 16810436092791821916
  * [26] Sloan, S., Pelletier, J., 2012. How accurately may we project tropical forest-cover change? A validation of a forward-looking baseline for REDD. _Global Environmental Change 22_ (2), 440-453. DOI:10.1016/j.gloenvcha.2012.02.001 Google Scholar: 10187012353585922701
  * [27] Nabuurs, G. J., van Putten, B., Knippers, T. S., Mohren, G. M. J., 2008. Comparison of uncertainties in carbon sequestration estimates for a tropical and a temperate forest. _Forest Ecology and Management 256_ (3), 237-245. DOI:10.1016/j.foreco.2008.04.010 Google Scholar: 6421873803139267676
  * [28] Green, D. G., Sadedin, S., 2005. Interactions matter – complexity in landscapes and ecosystems. _Ecological Complexity 2_ (2), 117-130. DOI:10.1016/j.ecocom.2004.11.006 Google Scholar: 558815488529198010
  * [29] Baker, R., Koch, F., Kriticos, D., Rafoss, T., Venette, R., van der Werf, W. (Eds.), 2012. Advancing risk assessment models for invasive alien species in the food chain: contending with climate change, economics and uncertainty. _Bioforsk FOKUS 7_. Bioforsk, Frederik A. Dahls vei 20, 1432 Ås, Norway. http://www.pestrisk.org/2012/BioforskFOKUS7-10_IPRMW-VI.pdf
  * [30] de Rigo, D., Caudullo, G., San-Miguel-Ayanz, J., Stancanelli, G., 2012. Mapping European forest tree species distribution to support pest risk assessment. In: Baker, R., Koch, F., Kriticos, D., Rafoss, T., Venette, R., van der Werf, W. (Eds.), Advancing risk assessment models for invasive alien species in the food chain: contending with climate change, economics and uncertainty. _Bioforsk FOKUS 7_. Bioforsk, Frederik A. Dahls vei 20, 1432 ås, Norway. http://www.pestrisk.org/2012/BioforskFOKUS7-10_IPRMW-VI.pdf Google Scholar: 6508055261897528514
  * [31] Thompson, I., Mackey, B., McNulty, S., Mosseler, A., 2009. Forest resilience, biodiversity, and climate change: a synthesis of the biodiversity/resilience/stability relationship in forest ecosystems. Vol. 43 of Technical Series. Secretariat of the Convention on Biological Diversity. ISBN: 9292251376 Google Scholar: 16541384134391848156
  * [32] Center for International Forestry Research., FAO Regional Office for Asia and the Pacific, 2005. Forests and floods: drowning in fiction or thriving on facts? Center for International Forestry Research; Food and Agriculture Organization of the United Nations, Regional Office for Asia and the Pacific. http://www.worldcat.org/isbn/9793361646
  * [33] Bonan, G. B., 2008. Forests and climate change: Forcings, feedbacks, and the climate benefits of forests. _Science 320_ (5882), 1444-1449. DOI:10.1126/science.1155121 Google Scholar: 12659549051072980947
  * [34] Ferreira, L., Constantino, M. F., Borges, J. G., Garcia-Gonzalo, J., 2012. A stochastic dynamic programming approach to optimize Short-Rotation coppice systems management scheduling: An application to eucalypt plantations under wildfire risk in portugal. _Forest Science_ , 353-365. DOI:10.5849/forsci.10-084 Google Scholar: 15823418833962927194
  * [35] de Rigo, D., Rizzoli, A. E., Soncini-Sessa, R., Weber, E., Zenesi, P., 2001. Neuro-dynamic programming for the efficient management of reservoir networks. In: Proceedings of MODSIM 2001, International Congress on Modelling and Simulation. Vol. 4. Modelling and Simulation Society of Australia and New Zealand, pp. 1949-1954. DOI:10.5281/zenodo.7481 Google Scholar: 16120008708786398621
  * [36] Bond, C. A., Champ, P., Meldrum, J., Schoettle, A., 2011. Investigating the optimality of proactive management of an invasive forest pest. In: Keane, R. E., Tomback, D. F., Murray, M. P., Smith, C. M. (Eds.), The future of high-elevation, five-needle white pines in Western North America: Proceedings of the High Five Symposium. U.S. Department of Agriculture, Forest Service, Rocky Mountain Research Station, pp. 295-302. http://www.treesearch.fs.fed.us/pubs/38241 Google Scholar: 15999737409186644552
  * [37] de Rigo, D., Castelletti, A., Rizzoli, A. E., Soncini-Sessa, R., Weber, E., 2005. A selective improvement technique for fastening neuro-dynamic programming in water resources network management. In: Zítek, P. (Ed.), Proceedings of the 16th IFAC World Congress. Vol. 16. International Federation of Automatic Control (IFAC), pp. 7-12. DOI:10.3182/20050703-6-CZ-1902.02172 Google Scholar: 3101827154043147997
  * [38] Phillis, Y. A., Kouikoglou, V. S., 2012. System-of-Systems hierarchy of biodiversity conservation problems. _Ecological Modelling 235-236_ , 36-48. DOI:10.1016/j.ecolmodel.2012.03.032 Google Scholar: 14693363263013674243
  * [39] Cavallo, A., Nardo, A., 2008. Optimal fuzzy management of reservoir based on genetic algorithm. In: Lowen, R., Verschoren, A. (Eds.), Foundations of Generic Optimization. Vol. 24 of Mathematical Modelling: Theory and Applications. Springer Netherlands, pp. 139-159. DOI:10.1007/978-1-4020-6668-9_2 Google Scholar: 596109054217603623
  * [40] Castelletti, A., de Rigo, D., Tepsich, L., Soncini-Sessa, R., Weber, E., 2008. On-Line design of water reservoir policies based on inflow prediction. In: Myung, C., Misra, P. (Eds.), Proceedings of the 17th IFAC World Congress. Vol. 17. International Federation of Automatic Control (IFAC), pp. 14540-14545. DOI:10.3182/20080706-5-KR-1001.02463 Google Scholar: 17611910683549006538
  * [41] Kempeneers, P., Sedano, F., Seebach, L. M., Strobl, P., San-Miguel-Ayanz, J., 2011. Data fusion of different spatial resolution remote sensing images applied to Forest-Type mapping. _IEEE Transactions on Geoscience and Remote Sensing 49_ (12), 4977-4986. DOI:10.1109/TGRS.2011.2158548 Google Scholar: 1622266323511963573
  * [42] Sedano, F., Kempeneers, P., Strobl, P., McInerney, D., San-Miguel-Ayanz, J., 2012. Increasing Spatial Detail of Burned Scar Maps Using IRS-AWiFS Data for Mediterranean Europe. _Remote Sensing 4_ (3), 726-744. DOI:10.3390/rs4030726 Google Scholar: 16654372463982241297
  * [43] de Rigo, D., Bosco, C., 2011. Architecture of a Pan-European Framework for Integrated Soil Water Erosion Assessment. _IFIP Advances in Information and Communication Technology 359_ , 310-318. DOI:10.1007/978-3-642-22285-6_34 Google Scholar: 3950024085016158193
  * [44] Voinov, A., Shugart, H. H., 2013. ’Integronsters’, integral and integrated modeling. _Environmental Modelling & Software 39_, 149-158. DOI:10.1016/j.envsoft.2012.05.014 Google Scholar: 17438791388360294950
  * [45] Mendoza, G. A., Martins, H., 2006. Multi-criteria decision analysis in natural resource management: A critical review of methods and new modelling paradigms. _Forest Ecology and Management 230_ (1-3), 1-22. DOI:10.1016/j.foreco.2006.03.023 Google Scholar: 17320958744182864860
  * [46] O’Farrell, P. J., Anderson, P. M. L., 2010. Sustainable multifunctional landscapes: a review to implementation. _Current Opinion in Environmental Sustainability 2_ (1-2), 59-65. DOI:10.1016/j.cosust.2010.02.005 Google Scholar: 11738195672718592653
  * [47] Dale, V. H., Beyeler, S. C., 2001. Challenges in the development and use of ecological indicators. _Ecological Indicators 1_ (1), 3-10. DOI:10.1016/S1470-160X(01)00003-6 Google Scholar: 6821471585672677713
  * [48] Gilbert, N., 2010. Balancing water supply and wildlife. _Nature_. DOI:10.1038/news.2010.505 Google Scholar: 4413680200627277742
  * [49] Barnes, N., Jones, D., 2011. Clear climate code: Rewriting legacy science software for clarity. _Software, IEEE 28_ (6), 36-42. DOI:10.1109/MS.2011.113 Google Scholar: 12280553201622082598
  * [50] Sanders, R., Kelly, D., 2008. Dealing with risk in scientific software development. _Software, IEEE 25_ (4), 21-28. DOI:10.1109/MS.2008.84 Google Scholar: 5979969088190020135
  * [51] Cerf, V. G., 2012. Where is the science in computer science? _Commun. ACM 55_ (10), 5. DOI:10.1145/2347736.2347737 Google Scholar: 15391984385653800704
  * [52] Pincas, U., 2011. Program verification and functioning of operative computing revisited: How about mathematics engineering? _Minds and Machines 21_ (2), 337-359. DOI:10.1007/s11023-011-9237-z Google Scholar: 8660462933112825884
  * [53] Sanders, P., 2009. Algorithm engineering – an attempt at a definition. In: Albers, S., Alt, H., Näher, S. (Eds.), Efficient Algorithms. Vol. 5760 of Lecture Notes in Computer Science. Springer Berlin Heidelberg, pp. 321-340. DOI:10.1007/978-3-642-03456-5_22 Google Scholar: 16310481551473795186
  * [54] Kleiner, K., 2011. Data on demand. _Nature Climate Change 1_ (1), 10-12. DOI:10.1038/nclimate1057 Google Scholar: 13711304819619162709
  * [55] Nature, 2011. Devil in the details. _Nature 470_ (7334), 305-306. DOI:10.1038/470305b
  * [56] Peng, R. D., 2011. Reproducible research in computational science. _Science 334_ (6060), 1226-1227. DOI:10.1126/science.1213847 Google Scholar: 905554772905069177
  * [57] Cai, Y., Judd, K. L., Lontzek, T. S., 2012. Open science is necessary. _Nature Climate Change 2_ (5), 299. DOI:10.1038/nclimate1509 Google Scholar: 10156741028963436768
  * [58] Ghisla, A., Rocchini, D., Neteler, M., Förster, M., Kleinschmit, B., 2012. Species distribution modelling and open source GIS: why are they still so loosely connected? In: Seppelt, R., Voinov, A. A., Lange, S., Bankamp, D. (Eds.), International Environmental Modelling and Software Society (iEMSs) 2012 International Congress on Environmental Modelling and Software. Managing Resources of a Limited Planet: Pathways and Visions under Uncertainty, Sixth Biennial Meeting. pp. 1481-1488. http://www.iemss.org/iemss2012/proceedings/D6_0897_Ghisla_et_al.pdf Google Scholar: 6265264564814017260
  * [59] Iverson, K. E., 1980. Notation as a tool of thought. _Communications of the ACM 23_ (8), 444-465. http://awards.acm.org/images/awards/140/articles/9147499.pdf Google Scholar: 15203139354397204728
  * [60] Lehman, M. M., 1989. Uncertainty in computer application and its control through the engineering of software. _J. Softw. Maint: Res. Pract. 1_ (1), 3-27. DOI:10.1002/smr.4360010103 Google Scholar: 11699944162229777276
  * [61] Lehman, M. M., Ramil, J. F., 2002. Software uncertainty. In: Bustard, D., Liu, W., Sterritt, R. (Eds.), Soft-Ware 2002: Computing in an Imperfect World. Vol. 2311 of Lecture Notes in Computer Science. Springer Berlin / Heidelberg, Ch. 14, pp. 477-514. DOI:10.1007/3-540-46019-5_14 Google Scholar: 15328955307076618013
  * [62] Hook, D., Kelly, D., 2009. Testing for trustworthiness in scientific software. In: Software Engineering for Computational Science and Engineering, 2009. SECSE ’09. ICSE Workshop on. IEEE, Washington, DC, USA, pp. 59-64. DOI:10.1109/SECSE.2009.5069163 Google Scholar: 9959381250960735693
  * [63] Hatton, L., 2007. The chimera of software quality. _Computer 40_ (8), 104-103. DOI:10.1109/MC.2007.292 Google Scholar: 691364900894680251
  * [64] Hatton, L., 1997. The t experiments: errors in scientific software. _Computational Science & Engineering, IEEE_ 4 (2), 27-38. DOI:10.1109/99.609829 Google Scholar: 8585904454769261330
  * [65] Hatton, L., 2012. Defects, scientific computation and the scientific method uncertainty quantification in scientific computing. _IFIP Advances in Information and Communication Technology 377_ , 123-138. DOI:10.1007/978-3-642-32677-6_8 Google Scholar: 1399353442026811576
  * [66] Lehman, M. M., 1996. Laws of software evolution revisited software process technology. In: Montangero, C. (Ed.), Software Process Technology. Vol. 1149 of Lecture Notes in Computer Science. Springer Berlin/Heidelberg, Ch. 12, pp. 108-124. DOI:10.1007/BFb0017737 Google Scholar: 3754503458653974527
  * [67] Oberkampf, W. L., DeLand, S. M., Rutherford, B. M., Diegert, K. V., Alvin, K. F., 2002. Error and uncertainty in modeling and simulation. _Reliability Engineering & System Safety 75_ (3), 333-357. DOI:10.1016/S0951-8320(01)00120-X Google Scholar: 10613684146262342334
  * [68] Wilson, G., 2006. Where’s the real bottleneck in scientific computing? _American Scientist 94_ (1), 5+. DOI:10.1511/2006.1.5 Google Scholar: 12828974620192082182
  * [69] Rebaudengo, M., Reorda, M., Violante, M., 2011. Software-Level Soft-Error mitigation techniques. In: Nicolaidis, M. (Ed.), Soft Errors in Modern Electronic Systems. Vol. 41 of Frontiers in Electronic Testing. Springer US, pp. 253-285. DOI:10.1007/978-1-4419-6993-4_9 Google Scholar: 15742657402676184752
  * [70] Beaudette, D., 2008. Simple comparison of two Least-Cost path approaches. In: Open Source Software Tools for Soil Scientists. http://casoilresource.lawr.ucdavis.edu/drupal/node/544 (archived at: http://www.webcitation.org/6D0LHBRXW ) Google Scholar: 14623451910171059303


# Dense Random Packing of Disks With a Power-Law Size Distribution in
Thermodynamic Limit:
Fractal-like Properties

Alexander Yu. Cherny  Joint Institute for Nuclear Research, Dubna 141980,
Russian Federation  Eugen M. Anitas  Joint Institute for Nuclear Research,
Dubna 141980, Russian Federation  Horia Hulubei, National Institute of Physics
and Nuclear Engineering, RO-077125 Bucharest-Magurele, Romania  Artem A.
Vladimirov  Joint Institute for Nuclear Research, Dubna 141980, Russian
Federation  Vladimir A. Osipov  Joint Institute for Nuclear Research, Dubna
141980, Russian Federation

###### Abstract

The correlation properties of a dense random system of disks with a power-law
size distribution are analyzed in momentum space in the thermodynamic limit.
The limit assumes that the total number the disks increases infinitely, while
the mean density of the disk centers and the range of the size distribution
are kept constant. The structure factor dependence on momentum transfer is
investigated for various numbers of disks and extrapolated to the
thermodynamic limit. The fractal power-law decay of the structure factor is
recovered in momentum space within the fractal range, which corresponds to the
range of the size distribution in real space. The fractal exponent coincides
with the exponent of the power-law size distribution as was shown previously
by the authors [Cherny, Anitas, Osipov, J. Chem. Phys. 158(4), 044114 (2023)].
The finite-size effects are also studied at very small momentum of order of
the inverse system size. We show that the structure factor is parabolic in
this region and calculate the prefactor analytically. The obtained results
reveal fractal-like properties of the packing. Our findings can be used to
analyze small-angle scattering from such systems.

dense random packing; power-law polydispersity; random sequential addition;
Delaunay triangulation;fractals; structure factor.

##  I Introduction

Dense random packings of granular particles or pores with sizes following a
power-law distribution are highly relevant for various applications [1]
including emulsion [2], soil [3], ceramics [4] and concrete technologies [5].
For hard materials this is due to the fact that power-laws give typically the
highest densities, as required in applications where very high resistance is
needed.

Besides the practical interest, there is a fundamental problem of
understanding the correlation properties in such systems. They exhibit complex
packings with fractal properties [6, 7, 8], which are controlled by two
parameters: the range of the power-law size distribution and its exponent.
Thus, this type of packing is characterized by a new kind of symmetry, termed
self-similarity [9].

The dense packing limit used in the previous paper [8] assumes that the total
number of the disks tends to infinity but the area occupied by the system
remains finite. In this limit, the lower border of the size distribution tends
to zero, which implies that maximum particle size ratio increases infinitely.
Here, we extend our results [8] to a more realistic case of a macroscopic
number of the disks and a fixed range of the size distribution. Then we study
the structure factor as a function of the momentum transfer in the
thermodynamic limit, where the area $A$ and the number of disks $N$ are
increased while $A/N=\mathrm{const}$.

To this aim, we propose a time-saving algorithm of dense random packing of
disks with a given size distribution, based on the Delaunay triangulation (DT)
[10]. The triangulation is commonly used for compact packing problems for
disks in 2D and spheres in 3D [11] and for analysis of pore space in a porous
media [12]. For the dense packing, the most time consuming task is to check
for possible intersections between a newly placed disk and a set of $N$
already placed disks. The DT helps us to optimize its computation time. This
algorithm needs only $O(\log N)$ repeated operations to select neighbours of
the newly placed disk and calculate distances to them. By contrast, without
the DT the same procedure requires $O(N)$ operations. For validation of DT-
based algorithm, the results from a modified version of random sequential
addition (RSA)-based algorithm used in Ref. [8] are also included.

In addition, we find that $S(q)\propto q^{2}$ for the small momentum transfer
$q\lesssim 2\pi/s$ (here $S(q)$ and $s$ are the structure factor and system
size, respectively). As is shown in Sec. III.1 below, the prefactor is
determined by the variance of the disk positions, which are randomly
distributed over the area.

##  II Algorithms of disk fractal generation

We consider a set of $N$ disks obeying the power-law distribution with the
exponent $D$ [6, 8]. The number of disks $dN(r)$ whose radii fall within the
range ($r,r+dr$) is proportional to $dr/r^{D+1}$. The radii vary from $a$ to
$R$. In two dimensions, the exponent satisfies the condition $0<D<2$.

###  II.1 Random-sequential-addition grid

The RSA-based algorithm was described in details in Ref. [8]. For a large
number of particles, it becomes too time-consuming, since the algorithm needs
about $O(N^{2})$ operations. For this reason, we use a simplified approach and
build the entire system from $k^{2}$ randomly generated cells. The system is a
square grid of $k$ by $k$ cells. Each cell is a trial, obtained by RSA [8] and
containing $M=N/k^{2}$ disks of radii $r_{i}=R\,{i}^{-1/D}$, $i=1,...,M$. Then
the thermodynamic limit is emulated by $k\to\infty$ with $M$ and $R/a$ being
constant. This scheme enables us to reduce the number of operations to
$O(M^{2}k^{2})=O(N/k^{2})$. It changes the spatial long-range correlations of
the disk positions at the distances bigger than $s/k$, but the fractal range
remains intact (see Sec. III.2 below).

Figure 1 (Left) shows a grid of $5\times 5$ cells with each cell containing
$M=5000$ disks with $D=1.5$ and $R/a\simeq 3.4\times 10^{-3}$. For these
control parameters, the packing fraction is about 0.95 and $s/R\simeq 15$.
Disks of the same radius are displayed in the same color: the largest disks
are red, and the smallest are purple.

Figure 1: (Left) A grid of $5\times 5$ cells; each one contains $M=5000$
disks, placed with the RSA algorithm. The total number of the disks
$N=125000$. (Right) A set of $N=125000$ disks placed by the DT algorithm. In
the both cases, the disks are compactly packed and randomly distributed, and
their radii follow the power-law distributions with the exponent $D=1.5$.
Disks with the same radius are displayed in the same color.

###  II.2 Delaunay triangulation algorithm

For the DT algorithm, we use a smoother realization of the same power-law
distribution as for the RSA grid

$r_{i}=R\left(\frac{m}{m+i-1}\right)^{1/D},\ \textrm{for}\ i=1,...,N.$ (1)

The smallest radius $a=R\big{[}m/(N+m-1)\big{]}^{1/D}$ slightly deviates from
$a=R({m}/{N})^{1/D}$, which is used in the RSA-based calculations at $m=k^{2}$
(Sec. II.1). The value of the biggest radius $R$ should be chosen to provide
the densest packing in the thermodynamic limit. It is realized when the area
of the embedding square is equal to the total area of the disks:

$s^{2}/m=\pi R^{2}m^{2/D-1}\sum_{i=m}^{\infty}i^{-2/D}=\pi
R^{2}m^{2/D-1}\zeta(2/D,m),$ (2)

where $s$ is the edge length of the square, and
$\zeta(t,a)=\sum_{k=0}^{\infty}(a+k)^{-t}$ is the Hurwitz zeta function. Then
the thermodynamic limit is $m\to\infty$ and $N/m=\text{const}$ and
$s^{2}/m\to\pi R^{2}D/(2-D)=\text{const}$ in full analogy with Sec. II.1 at
$m=k^{2}$. Note that the limit $\lim_{m\to\infty}s^{2}/m=\pi R^{2}D/(2-D)$ is
also obtained with the continuous power-law distribution as it should be.

Figure 2 compares the distributions of the radii used in the RSA grid and
given by Eq. (1). We observe that the latter distribution is smoother and the
both decay with the power-law $r^{-D}$ as it should be.

Figure 2: The cumulative power-law distribution in units of the largest radius
$R$ for the two structures shown in Figure 1. Vertical dashed line (blue)
denote the minimum radius $a$. Horizontal dashed line (green) denotes the
number of cells in the RSA grid.

To provide the compact packing of the disks within a finite area, the exponent
of the distribution should be restricted to $D<2$, since for $D\geqslant 2$
the total area of the infinite set of disks becomes infinite. Disks are
arranged into a square by placing them inside the square in the sequence from
largest to smallest. As the set of already placed disks grows, the
computational costs of checking for their collisions increase. Thus, we use
the DT to select neighbours of a newly placed disk and thus reduce the number
of the required calculations. The DT covers a convex hull of a given set of
points by a set of triangles in such a way that no point lies inside the
circumcircle of any triangle. There are very effective implementations of the
triangulation that need only $O(\log N)$ operations for both adding a new
point in the triangulation and finding neighbours of the point. We use the
code from the CGAL software package [13].

The algorithm consists of the following steps:

  1. 1.

Positions of the first six disks are searched randomly over the entire square.

  2. 2.

To place a new disk, we choose a random point for its center. The function
_vertices_on_conflict_zone_boundary_ of the CGAL package is used to obtain the
subset of already placed disks that are neighbours to the new disk and are
most likely to intersect with it. If the randomly selected point is too close
to their neighbours and collisions appear, the algorithm searches for an empty
space inside the polygon that includes the neighbours of the selected point.
To find an empty space, the new disk is bounced back with the shift vector
that amounts to the sum of the vectors of all found overlappings. If no
collisions appear then we go to the step 3. If a collision with the neighbours
still exists, a new random point is chosen once more and the step 2 is
repeated again. When the maximum number of unsuccessful attempts is reached,
the algorithm is restarted from the step 1.

  3. 3.

Before placing a new disk, algorithm does the final check for collisions with
all already placed disks. This is because the selection of the neighbors with
the DT might not work perfectly, although the probability of the failure is
very small.

  4. 4.

The center of the newly placed disk is added to the triangulation set and the
algorithm proceeds from the step 2.

Note that this algorithm does not place disks in a completely random way.
Although each new disc is initially dropped on the plane at random point,
after that the algorithm simulates the disc repelling from its neighbors and
the boundary of the square. In the case of $m=1$ we place the first disk at
the corner of the square, for $m>1$ the first disk is placed to a random
point. The larger the number of the disks and the packing fraction, the more
effective the suggested algorithm is in comparison to a simple random
placement of the disks such as RSA. The algorithm allows us to limit the
search for empty space to the neighborhood of any added disk, which makes it
especially effective at high packing fractions. The increase of performance is
approximately proportional to $1/(1-F)$ where $F$ is the current filling rate.

Figure 1(Right) shows one trial of a set of $N=125000$ disks generated using
the DT-based algorithm with the distribution (1) at the same control
parameters as described in Sec. II.1. The color coding is the same as in
Figure 1(Left). For the RSA grid, the spatial disk distribution is more
homogeneous by its construction. However, the spatial distributions of RSA and
DT are almost indistinguishable on scales of order $r\lesssim s/k$.

##  III The structure factor

We define the structure factor through the density _fluctuations_ over the
background of the density mean value

$S(q)\equiv\langle\delta\rho_{\bm{q}}\delta\rho_{\bm{-q}}\rangle_{\hat{q}}/N,$
(3)

where the Fourier components of the density fluctuations are given by

$\displaystyle\delta\rho_{\bm{q}}=\sum_{j=1}e^{-i\bm{q}\cdot\bm{r}_{j}}-\frac{N}{s^{2}}\int
d\bm{r}e^{-i\bm{q}\cdot\bm{r}}.$ (4)

Here the brackets $\langle\cdots\rangle_{\hat{q}}$ stand for the average over
all directions of unit vector $\hat{q}$ along $\bm{q}$, and the integral is
taken over the area of the embedding square. Note that in the previous papers
[6, 8], the structure factor was defined without subtracting the scattering
amplitude of the background. Such a definition leads to the appearance of a
spike of order of $N$ in the range $q\lesssim 2\pi/s$, which is not convenient
to study the thermodynamic limit. The definition (3) enables us to exclude the
spike.

Figure 3 shows the structure factor (3) for the dense packings depicted in
Figure 1. Each curve is calculated numerically with 20 trials. We observe that
the errors are almost negligible from the border of the fractal range
$q\gtrsim 2\pi/R$, see the discussion in Sec. III.2 below. This implies that
at the given control parameters $N$, $D$, and $R$, the structure factor is
independent of a specific packing configuration in the limit of high packing
fraction. However, when $q\lesssim 2\pi/R$ the errors are quite pronounced,
especially for the DT-based curves.

In order to illustrate graphically how the structure factor is obtained,
Figure 3 shows also the structure factor (of the densities) given by Eq. (3)
but without subtracting the scattering amplitude of the background square in
Eq. (4) (black curve). We observe that that for $q\lesssim 2\pi/s$ the
dominant contribution comes from the embedding square (green curve), while for
$q\gtrsim 2\pi/R$, the main contribution arises from the structure factor of a
single cell. The crossover point lies somewhere in between. The contribution
of a single cell was studied in detail in the previous paper [8].

A quantitative analysis of the structure factor is performed after smearing
the curves in Figure 3 with a log-normal distribution function, as described
in Ref. [8]. This procedure eliminates the maxima and minima by smoothing them
and provides a reliable way to estimates the power-law exponents by a linear
fit in a double-logarithmic plot. Moreover, the smoothing is a standard
procedure in analyzing small-angle scattering data, since there are three
contributions to the smearing of an ideal curve: the angular divergence of the
beam, the finite resolution of the detector, and the polychromatic nature of
the beam [14]. The results of smearing are shown in Figure 4 for various
number of disks $N$. For the RSA grid, the chosen values of $N$ correspond to
different grid sizes ($1\times 1$ at $N=5000$, $2\times 2$ at $N=20000$ etc).
In particular, in Figure 4 the violet curve at $N=125000$ is the smeared
version of the blue curve in Figure 3.

Figure 3:  The structure factor of the density fluctuations (blue) [i.e., Eq.
(3) with Fourier components (4)] in units of the inverse radius $1/R$ of the
largest disk (black). The structure factor of the densities (black) [i.e., Eq.
(3) with Fourier components (4) but without the second term]. The scattering
intensity from the embedding square is given in green. (a) RSA. (b) DT. In the
both cases, the red bars represent errors (the standard deviations for 20
trails). The number of disks $N=125000$ and $s/R\simeq 15$. The corresponding
structures are shown in Figure 1. The inset figures show the exact agreement
between a single trial (blue) and the analytical relation (5) (orange) at
small values of the wave vector $q\lesssim 2\pi/s$.

###  III.1 Small momenta: finite-size effects

We choose the origin of coordinates in the center of the square. Expanding Eq.
(4) at small momentum, we obtain:
$\delta\rho_{\bm{q}}=-iN\bm{q}\cdot\bm{r}_{\text{c}}+O(q^{2})$ with
$\bm{r}_{\text{c}}=\frac{1}{N}\sum_{j=1}^{N}\bm{r}_{j}$ being the center-of-
mass position of the total set of the points. By substituting this formula
into the definition (3) and using the identity
$\langle(\bm{q}\cdot\bm{r}_{\text{c}})^{2}\rangle_{\hat{q}}=q^{2}r_{\text{c}}^{2}/2$
in 2D, we are left with the structure factor at $q\lesssim 2\pi/s$

$S(q)\simeq\frac{N}{2}q^{2}r_{\text{c}}^{2}.$ (5)

The vector $\bm{r}_{\text{c}}$ is non-zero only due to the fluctuations of the
density. When the distribution of the points $\bm{r}_{j}$ over the square is
random and uniform then its variance
$\left\langle\left(\sum_{j=1}^{N}\bm{r}_{j}\right)\cdot\left(\sum_{j=1}^{N}\bm{r}_{j}\right)\right\rangle$
is proportional to $N$. It follows that the structure factor (5) is
independent of $N$ in the thermodynamic limit as it should be. Nevertheless,
in the thermodynamic limit, the range $q\lesssim 2\pi/s$ shrinks to zero, and,
hence, the formula (5) describes a finite-size effect anyway.

Figure 4 illustrates these features. For the RSA grid, the structure factor
remains practically unchanged with increasing $N$ when $q\lesssim 2\pi/s$. For
the DT based algorithm, the structure factor approaches the saturation curve
with increasing $N$. For $N>80000$, the corresponding scattering curves
(compare the cyan and violet ones) are practically indistinguishable in this
region. Note that all the curves fall within the error bars when $q\lesssim
2\pi/s$.

Figure 4 shows also that the power-law exponent at small momenta is exactly 2,
as predicted by Eq. (5). Moreover, the coefficient of $q^{2}$ in Eq. (5) is
exactly recovered from the corresponding curves (see the insets of Figure 3).

Figure 4: The smeared structure factor (3) vs. the wave vector (in units of
the inverse radius $1/R$ of the larges disk) at various numbers of particles
$N$. (a) RSA. (b) DT. In the both cases the red bars represent the errors. The
vertical dotted lines indicate the beginning of the fractal range.

###  III.2 High momenta: Fractal region and beyond

The fractal range in reciprocal space is given by [8]

$2\pi/R\lesssim q\ll 2\pi/a.$ (6)

The lower and upper borders of the fractal region are related to the largest,
and respectively smallest radii in the configuration. The left border of this
range is indicated by a vertical dotted line in Figure 4(a) and 4(b), as well
as in Figure 5(a) and 5(b). In the fractal range the structure factor decays
proportional to $q^{-\alpha}$. A linear fit in the range $22\lesssim
qR\lesssim 131$ gives the value $\alpha=1.508\pm 0.009$ for the power-law
exponent, which confirms numerically that $\alpha=D$.

Beyond the fractal region, when $q\gtrsim 2\pi/a$, the structure factor
attains its asymptotics with $S(q)\simeq 1$, as indicated by the horizontal
dotted line in Figure 4(a) and 4(b), as well as in Figure 3(a), 3(b), 5(a) and
5(b). The lengthy region where the structure factor is smaller than one
appears due to the finite values of smallest radius $a$ (see the detailed
discussion in Ref. [8]).

###  III.3 Random-sequential-addition grid vs. Delaunay triangulation
algorithm

A direct comparison of the structure factor (3) for the RSA grid and DT-based
algorithms is given in Figure 5 for the non-smeared (a) and smeared (b)
curves. The results show an excellent agreement between the two curves in the
fractal region and beyond. This indicates that the RSA and DT packing
algorithms yield the same fractal exponent, which coincides with the power-law
exponent of the size distribution. This is in contrast to jammed packings
generated using discrete element method simulations: it has been reported that
the structure-factor exponent is independent of the exponent of the power-law
size distribution in three dimensions and only slightly depends on it in two
dimensions [15].

At small momenta $q\gtrsim 2\pi/a$, the two curves have very distinct values
of the coefficients of $q^{2}$ in Eq. (5). This is related to the fluctuations
of the center-of-mass position of the disk centers. The fluctuations are
artificially suppressed by the method of construction of the RSA grid (see
Sec. II.1) and as a consequence, the coefficient is significantly reduced.

Figure 5: (a) The structure factors (3) vs. the wave vector (in units of the
inverse radius $1/R$ of the largest disk) for the RSA grid and DT at
$N=80000$. (b) The corresponding smeared structure factors.

##  IV Conclusions

In this work, the structure factor $S(q)$ [Eq. (3)] is used to study in the
thermodynamic limit the correlation properties of densely packed disks with a
power-law distribution of their sizes. The thermodynamic limit is simulated
with different methods of packings: the random sequential addition grid and
Delaunay triangulation-based algorithms. The spike in the structure factor at
small momenta, which is proportional to the total number of the disks, is
excluded by subtracting the scattering amplitude of the embedding square, see
Sec. III. Apart from the formal reasons, this definition can be quite
appropriate from the physical point of view, for instance, in studying the
internal structures by the small-angle neutron scattering with the contrast
variation method.

The main result of the paper is that in the thermodynamic limit, the fractal
range in $S(q)$ is still observable for the both algorithms: its borders (6)
are determined by ranges of the size distribution, and the fractal exponent
coincides with the exponent of the power-law size distribution [8].

The finite size effects are studied as well. It is shown (see Sec. III.1) that
in the range of small momenta, $S(q)\propto q^{2}$ and it is independent of
the number of particles, although this range shrinks in the thermodynamic
limit. We derive the expression for the coefficient of $q^{2}$, which is
related to the fluctuations of the center-of-mass of the disk positions.

The obtained results can be used for analyzing experimental small-angle
scattering data from densely packed objects, whose centers have a much higher
scattering length density as compared to its surrounding shell.

## References

  * [1] Torquato, S. Perspective: Basic understanding of condensed phases of matter via packing models. J. Chem. Phys. 2018, 149, Article ID: 020901. [CrossRef]
  * [2] Kwok, S.; Botet, R.; Shrapnack, L.; Cabane, B. Apollonian packing in polydisperse emulsions. Soft Matter 2020, 16, 2426–2430. [CrossRef]
  * [3] Iglauer, S.; Favretto, S.; Spinelli, G.; Schena, G.; Blunt, M. J. X-ray tomography measurements of power-law cluster size distributions for the nonwetting phase in sandstones. Phys. Rev. E 2010, 82, Article ID: 056315. [CrossRef]
  * [4] Bhattacharjee, A.; Hurley, R. C.; Graham-Bradi, L. Fragmentation and granular transition of ceramics for high rate loading. J. Am. Ceram. 2021, 105, 3062–3080. [CrossRef]
  * [5] Qifeng, L.; Yanwei, G.; Anguo, C.; Jie, J. Pore-network distribution laws of cementitious materials detected by high-energy X-ray computed tomography imaging system J. Phys.: Conf. Ser. 2021, 2011, Article ID: 012012. [CrossRef]
  * [6] Cherny, A. Yu.; Anitas, E. M.; Osipov, V. A.; Kuklin , A. I. Revised scattering exponents for a power-law distribution of surface and mass fractals Phys. Rev. E 2022, 106, Article ID: 024108. [CrossRef]
  * [7] Monti, J. M.; Clemmer, J. T.; Srivastava, I.; Silbert, L. E.; Grest, G. S. Large-scale frictionless jamming with power-law particle size distributions Phys. Rev. E 2022, 106, Article ID: 034901. [CrossRef]
  * [8] Cherny, A. Yu.; Anitas, E. M.; Osipov, V. A. Dense random packing with a power-law size distribution: The structure factor, mass–radius relation, and pair distribution function J. Chem. Phys. 2023, 158, Article ID: 044114. [CrossRef]
  * [9] Special Issue ”Symmetry and Fractals” in the journal Symmetry. Available online: [CrossRef] (accessed on ).
  * [10] Delaunay, B. Sur la sphère vide. A la mémoire de Georges Voronoï Bull. Acad. Sci. URSS 1934, VII(6), 793–800. [CrossRef]
  * [11] Reboul, N.; Vincens, E.; Cambou, B. Statistical analysis of void size distribution in a simulated narrowly graded packing of spheres Granul Matter. 2008, 10, 457–468. [CrossRef]
  * [12] Nguyen, N.-S.; Taha, H.; Marot, D. A new Delaunay triangulation based approach to characterize the pore network in granular materials Acta Geotech. 2021, 16, 2111–2129. [CrossRef]
  * [13] The CGAL Project CGAL User and Reference Manual, 5.5.2 ed.; Publisher: CGAL Editorial Board, 2023; [CrossRef]
  * [14] Wignall, G. D. Instrumental resolution effects in small-angle scattering J. Appl. Cryst. 1991, 24, 479–484. [CrossRef]
  * [15] Monti, J. M.; Srivastava, I; Silbert, L. E.; Lechman, J. B.; Grest, G. S. Fractal dimensions of jammed packings with power-law particle size distributions in two and three dimensions arXiv 2023, 2306.12499, 1–6. [CrossRef]