Matias Martinez

Westley Weimer

Martin Monperrus

Do the Fix Ingredients Already Exist? An Empirical Inquiry into the Redundancy Assumptions of Program Repair Approaches

Keywords: D.2.5 [Software Engineering]: Testing and Debugging Automatic program fixing, automatic software repair, mining software repositories

Much initial research on automatic program repair has focused on experimental results to probe their potential to find patches and reduce development effort. Relatively less effort has been put into understanding the hows and whys of such approaches. For example, a critical assumption of the Gen-Prog technique is that certain bugs can be fixed by copying and re-arranging existing code. In other words, GenProg assumes that the fix ingredients already exist elsewhere in the code. In this paper, we formalize these assumptions around the concept of "temporal redundancy". A temporally redundant commit is only composed of what has already existed in previous commits. Our experiments show that a large proportion of commits that add existing code are temporally redundant. This validates the fundamental redundancy assumption of GenProg.

INTRODUCTION

To some extent, each program repair technique is based on an underlying assumption. GenProg's "secret sauce" [START_REF] Weimer | Automatically finding patches using genetic programming[END_REF][START_REF] Le Goues | A systematic study of automated program repair: Fixing 55 out of 105 bugs for $8 each[END_REF] is the assumption that large programs contain the seeds of their own repair and thus that rearrangements of existing statements can fix most bugs. This redundancy assumption is also behind four out of PAR's ten repair templates [START_REF] Kim | Automatic patch generation learned from human-written patches[END_REF]. By contrast, SemFix [START_REF] Nguyen | SemFix: Program Repair via Semantic Analysis[END_REF] makes different assumptions: it is based on the idea that some bugs can be repaired by changing only one variable assignment or if conditional expression. Although extensive empirical evaluations of the bug fixing powers and costs of these techniques have been carried out, less work has focused on exploring or validating such fundamental assumptions.

In this paper we focus on the assumptions made by techniques like PAR and GenProg, and ask whether it makes sense to re-arrange existing code or code changes to fix bugs. We formalize our investigation by defining a concept of "software temporal redundancy". A commit is said to be temporally redundant if it is only a re-arrangement of the code that has already existed in previous commits. We provide further insight into the problem by breaking down existing code along multiple dimensions, including the granularity at which it is considered (e.g., line-level vs. token-level) and the locality at which it is found (e.g., in the same file vs. anywhere in the program).

We measure temporal redundancy on a dataset of six open-source Java applications. For instance, at token level, we find that 31-52% of commits are temporally redundant. Our results show that a substantial fraction of software evolution takes place in a "closed world", with no creation of new material. To our knowledge, this paper is the first empirical quantification linking such evidence to automated program repair.

Our experiments enable us to understand the foundational assumptions of program repair approaches based on redundancy such as GenProg or PAR. First, we validate the redundancy assumption: it does makes sense to fix defects by rearranging existing statements. Second, we find that different assumptions about how patches should be formed are significant: depending on the way one picks code from elsewhere in the program, the search space is different in terms of size and potential for success. Those results are directly actionable for improving such automated program repair techniques. For example, our results indicate that only considering local redundancy at the line level would decrease the repair time without a similar reduction in success potential.

To sum up, our contributions are:

• the definition and formalization of software temporal redundancy;

• an experimental protocol to measure software temporal redundancy using versioning data;

• an empirical validation of the redundancy assumption of some program repair approaches on 16071 commits of open-source software applications;

• an empirical analysis on the relation between the way existing redundant code is selected for use in patches and the structure of the underlying search space.

SOFTWARE TEMPORAL REDUNDANCY

Research in software evolution uses time as a primary focus. For most software projects, software evolution artifacts are captured by version control systems (e.g. CVS, SVN or GIT). In this section, we present a new commit-level software evolution metric called "temporal redundancy".

Commits.

A commit in a version control system consists of new file versions. Conceptually, a commit can be viewed in two ways: as a set of file pairs (before and after the commit), or as a set of changes (applied to the version before the commit to obtain the version after it).

In this paper, we use this change-based view of commits. We consider that a commit adds and/or deletes new source code fragments. An update is considered as the combination of a deletion and an addition. We do not consider commits done on other artifacts than source code.

Fragment Redundancy.

We use fragment to denote a substring of source code. For instance, the source code line "for (int i=0;i<n;i++)" is a fragment. Fragments are always defined according to a level of granularity. In this paper we consider two different levels: lines (as separated by line breaks) and tokens (as separated by lexing rules).

A fragment F is snapshot redundant at time T if another instance of that same fragment F exists elsewhere in the program at time T . This is the redundancy studied by Gabel and Su [START_REF] Gabel | A study of the uniqueness of source code[END_REF] and used by GenProg [START_REF] Weimer | Automatically finding patches using genetic programming[END_REF]. In this paper we consider a richer notion of redundancy that includes historical context.

A fragment F is temporally redundant at time T if that same fragment F has already been seen during the history of the software under analysis (i.e., at time T ′ < T ). For instance, literal "42" might be added in version #1, be removed in version #2 and reused again in version #3. In the commit of version #3, the "42" fragment is temporally redundant. Thus, in this paper, once a fragment has appeared it is always subsequently viewed as a potential source of redundancy. We consider the first version of a program to be created by insertions from an empty initial program.

Commit Redundancy.

A source code commit is composed of added and/or deleted fragments. Using a cooking metaphor, the added fragments are the "ingredients" of the commit.

We define a temporally redundant commit as a commit for which all added fragments are (individually) temporally redundant. More formally, we define a commit Cj performed at the time Tj as a set Sj of added fragments and a set Rj of removed fragments. Let C be the set of all commits for a program. Then a temporally redundant commit Cj satisfies

∀f ∈ Sj| ∃Ci ∈ C| Ti < Tj ∧ f ∈ Si
For such commits, no new fragments are invented and no fresh material is introduced: the commit is only a re-arrangement of insertions that have already been seen in previous commits.

Scope of Temporal Redundancy.

A fragment is temporally redundant if that same fragment appeared in a previous commit. This kind of redundancy has a global scope: the location of the previous fragment instance does not matter.

We now define a more restricted local scope notion of temporal redundancy. A fragment is locally temporally redundant if that same fragment has been used in a previous commit to the same file.

For example, consider two files, F1 and F2, each containing 3 fragments: F1 = {a, b, c}, F2 = {d, e, f }. Suppose commit C1 adds fragment c to file F2. For that commit, fragment c is global temporally redundant (already available in F1), but not locally temporally redundant (never previously available in F2). Suppose commit C2 introduces another version of F2 replacing fragment e with d. In that commit, fragment d is locally redundant (since d was previously available in F2).

MEASURING REDUNDANCY

We now present an experimental design to measure the fraction of relevant commits to a program that are temporally redundant.

Experimental Protocol.

Given a level of granularity and a scoping level, our experimental protocol to measure the temporal redundancy of the evolution of a program consists of the following phases: a) Retrieving commits. All commits of the program under analysis are collected from the repository.

b) Filtering commit files. Only commits to executable code are retained; commits to test cases are discarded. c) Fragmenting commits. We split each relevant file into fragments at a given level of granularity (e.g., lines or tokens). This results in a before-commit and an after-commit sequence of fragments. We use the Myers differencing algorithm [START_REF] Myers | An o(nd) difference algorithm and its variations[END_REF] to compare both fragment sequences and obtain the added fragments of the commit.

d) Filtering fragments. We filter out whitespace and comments. In this paper we are only interested in the evolution of executable code and not in indentation or documentation.

e) Selecting acceptable commits. We select those commits that introduce at least one fragment after filtering. We call such commits acceptable. f) Indexing fragments. We consider each added fragment in each acceptable commit in ascending temporal order. If a fragment has not been encountered previously at the given scoping level (i.e., global or local), we index it with the date of its first introduction.

g) Measuring temporal redundancy. The temporal redundancy of the entire program's evolution is the fraction of acceptable commits that are temporally redundant (see Section 2). For example, a temporal redundancy of 0.1 means that 10% of acceptable commits are temporally redundant commits.

Dataset.

We use six open-source Java projects to measure the temporal redundancy. They are: Apache Log4j, JUnit, Picocontainer, Apache Commons Collections, Apache Commons Math, and Apache Commons Lang. The inclusion criterion is as follows: the Apache projects were used in previous research on automatic program repair [START_REF] Kim | Automatic patch generation learned from human-written patches[END_REF], while the remaining two are Java projects mentioned in previous research on software evolution [START_REF] Kagdi | A survey and taxonomy of approaches for mining software repositories in the context of software evolution[END_REF]. After applying the filters presented above on the 16071 commits of the dataset, we obtain 7076 acceptable commits.

EMPIRICAL RESULTS

We now present our empirical results on the temporal redundancy of software (as defined in Section 2) following the experimental design presented in Section 3.

Line-Level Temporal Redundancy

Research question: What is the amount of line-based temporal redundancy?

For each application of our dataset, we measured the total number of acceptable commits within the analysis timespan (the complete data is available at http://goo.gl/k0rZWc) and the global-scope line-level temporal redundancy. Columns #2 and #3 of Table 1 report the results.

For instance, for log4j (the first row), there are 1687 commits which add at least one executable line. Only 9% of those 1681 are temporally redundant commits.

Overall, at the level of lines, 3-17% of the accepted commits are temporally redundant commits. Their basic ingredients are only previously-inserted code.

This has additional implications for automatic repair and code synthesis: for synthesizing those commits, the searchspace has a finite number of atomic building blocks (previouslyobserved line-level fragments). Theoretically, a redundancybased approach should be able to synthesize all the temporallyredundant commits.

The interval 3-17% is large and the reasons behind this variation are not obvious. One reason could be that the commit conventions used by the developers of a project are different. For instance, some projects prefer to have small and atomic commits (one bug fix or feature per commit). Other projects are less restrictive on this point. This has a direct impact on the redundancy: small and atomic commits are more likely to be redundant. A deeper analysis of this question is future work.

Token-Level Temporal Redundancy

Research question: Is there a difference between line-based and token-based temporal redundancy?

Before answering this question, we note that, analytically, all temporally redundant commits at the level of lines are necessarily temporally redundant commits at the level of tokens. Furthermore, a unique new line might be exclusively composed of existing tokens. Consequently, the token-based temporal redundancy must be equal to or greater than linebased temporal redundancy. We now measure the temporal redundancy at the line and token level.

Table 1 reports global scope line-level (column #3) and token-level temporal redundancy (column #7). For instance, in log4j, there is a line-level temporal redundancy of 9% but a token-level redundancy of 39%. Overall, at the token level, 29-52% of commits are temporally redundant.

For all projects, token-based temporal redundancy exceeds line-based. This follows from the analytical argument above and gives confidence in the experiment's construct validity.

Overall, token-level temporal redundancy (between 29% and 52%) is much higher than the line-level temporal redundancy. This is very good news. For automated repair and code synthesis, a high temporal redundancy implies a smaller search space. This holds for both the line and token level of granularity. Our token-level temporal redundancy measurements imply that for between 29% and 52% of accepted commits, synthesis and repair need never invent a new token. For instance, repair or synthesis of arithmetic code need only consider recombining existing literals and operations for one-third to one-half of commits.

We thus propose that the repair search space is composed of two components: the search space of fragments (the atomic building blocks) and the search space of their combinations. Our experiment enables us to precisely measure the former: we can count the total number of fragments seen up to a given point in time. We call this a fragment pool. There is one fragment pool per level of granularity. This is conceptually similar to the pool of four DNA bases or the pool of amino acids in biology.

Table 1 gives the size of the global scope fragment pool for line-level (column #4) and token-level (column #8) analyses. For instance, in the considered slice of history of log4j, there are 43313 different lines and 14294 different tokens that are involved in the software evolution. For all applications under study, the token pool at the point in time of the last commit is much smaller than the line pool.

For automated repair or code synthesis, there is a tension between working with the line pool or the token pool. To some extent, the temporally redundant commits correspond to the number of commits that can be synthesized. With the line pool, the combination of lines is much smaller (the combination space is smaller) but fewer commits can be synthesized (∼10%). With the token pool, more commits can be synthesized (∼40%), but at the price of exploring a much bigger combination space.

Redundancy Scope Experiment

Research question: Do file scope restrictions impact the amount of temporal redundancy?

We now measure the temporal redundancy available at the local scope in the same file (as defined in Section 2).

Table 1 reports global and local scope temporal redundancy in our dataset. For each granularity, there is one column "Global" and one column "Local" corresponding to the different scope. For instance, at the line level, column #3 is the temporal redundancy at the global scope and column #5 gives it when considering a local scope.

As discussed in Section 4.1, at the line granularity, there are between 3% and 17% of temporally redundant commits at the global scope. At the local scope, there are between 2% and 16%.

The temporal redundancy of both scopes is of the same order of magnitude. In all projects, more than half of the temporally redundant commits actually have local temporal redundancy. Consequently, at line granularity, most of the temporal redundancy is localized in the same file.

At token-level granularity, the results are similar: we find a large amount of token-level local-scope temporal redundancy. Tokens are likely to be reused in commits impacting the same file. This further indicates that the fragment locality matters during software evolution. We note that the difference of redundancy between global and local scopes is slightly higher at the token level (col. #7 vs #9) than at the line level (col. #3 vs #5). This is again a promising result with respect to the search space of automatic repair or synthesis. First, at the line level, the local scope pool is able to seed the same order of magnitude of commits as the global one. In other words, it is almost as fertile as the global pool. Second, when one considers the local scope pool, the search space is much smaller. For instance, for log4j, the median local pool size at the line level is 57 lines, compared to 43313 at the global scope level. Restricting attention to the local scope reduces the search space greatly while still enabling the synthesis of a large number of commits. Those results are directly actionable for improving GenProg and other redundancy-based approaches: our results indicate that only considering local redundancy would decrease the repair time while keeping a high repair success potential. This is in line with our body of previous experiments.

RELATED WORK

Some work on code clone devetion also focuses on software redundancy, studying the differences between line-and token-level granularity [START_REF] Kamiya | Ccfinder: a multilinguistic token-based code clone detection system for large scale source code[END_REF]. Kim et al. [START_REF] Kim | An empirical study of code clone genealogies[END_REF] considered code clones via a temporal perspective, linking clones together across versions. By contrast, we study whether changes only contain existing code. Most importantly, a change composed of many redundant fragments may not be a clone in itself: each redundant fragment could come from different, unrelated, locations or appear in a rearranged order.

In previous work [START_REF] Martinez | Mining software repair models for reasoning on the search space of automated program fixing[END_REF], we discussed the topology of the repair space and efficient manners to navigate it. In this paper, we focus on the relation between commits and the repair space, which has not been previously addressed.

Gabel and Su [START_REF] Gabel | A study of the uniqueness of source code[END_REF] studied the uniqueness of source code. Their "syntactic redundancy" calculates the degree to which portions of software applications are redundant. Hindle and colleagues [START_REF] Hindle | On the naturalness of software[END_REF] studied the repetitiveness and predictability of code. Both consider software redundancy from a spatial viewpoint. By contrast, we study temporal redundancy. Our measures are not present in, and cannot be inferred from, their experiments.

Nguyen et al. [START_REF] Nguyen | A study of repetitiveness of code changes in software evolution[END_REF] presented a study of repetitiveness of code changes in software evolution. Where they measure the repetitiveness of code changes by abstracting over commits and literals, we measure the redundancy of code changes at the commit level considering real, non-abstracted code (e.g., we care about the actual values of literals).

CONCLUSION

In this paper, we define a notion of software temporal redundancy and empirically measure its presence in source code. For example, as many as 52% of commits are composed entirely of previously-existing tokens. Our results have direct application to automated program repair. For example, we find that searches for line-level repairs can focus on a single file, dramatically reducing the search space without substantially reducing repair success potential.

Our future work will focus on measuring temporal redundancy on bug fix commits and on studying different granularities beyond lines and tokens (e.g. at the AST level, focusing on some kind of source code elements).

  To appear in the Proceedings of the 2014 International Conference on Software Engineering, track New Ideas and Emerging Results (ICSE NIER 2014)

Table 1 :

 1 The temporal redundancy of six open-source applications.

			Line granularity		Token granularity
	Program	Acceptable Commits	Global Temporal	Pool	Local Temporal	Pool	Global Temporal	Pool	Local Temporal	Pool
			redundancy	Size	redundancy	Size	redundancy	Size	redundancy	Size
	log4	1687	9% 43313	6%	57	39% 14294	19%	71
	junit	713	17%	8855	16%	18	43%	3256	29%	72.5
	pico	157	3% 16911	2% 22.5	31%	6273	8%	46
	collections	1019	7% 25406	4%	35	52%	4163	23%	85.5
	math	2210	6% 69943	4%	37	45% 20742	18% 100.5
	lang	1290	8% 22330	6%	63	50%	6692	29%	98
	C 1	C 2	C 3	C 4	C 5	C 6	C 7	C 8	C 9	C 10

ACKNOWLEDGMENTS

This research is done with support from EU Project Diversify FP7-ICT-2011-9 #600654 and the NSF (SHF-0905236, CCF-1116289, CCF-0954024), AFOSR (FA9550-07-1-0532, FA9550-10-1-0277), and DARPA (P-1070-113237). We thank Maxence G. de Montauzan and Sylvain Magnier for their participation to the experimentation.


Isabel Guillamon

Hermann Suderow

Sebastian Vieira

I Guillamón

P Rodière

Scanning tunneling spectroscopy with superconducting tips of Al

Keywords: Superconductivity, Scanning Tunneling Microscopy and Spectroscopy, 2H-NbSe 2

published or not. The documents may come

Scanning tunneling spectroscopy with superconducting tips of Al

Introduction

The development of new experimental methods for scanning tunneling microscopy and spectroscopy (STM/S) at very low temperatures has opened important paths to obtain fundamental information at the nanoscale. Studies with unprecedented detail can be made thanks to the construction of low temperature STM/S apparatus. For instance, the temperature stability of STS experiments in dilution refrigerators is excellent and allows for long term studies of the electronic properties of surfaces with a resolution in energy which is orders of magnitude below characteristic superconducting energies. Not surprisingly, detailed 1 * Corresponding author: hermann.suderow@uam.es measurements in different kinds of superconducting materials have unraveled the presence of local size inhomogeneities of the superconducting state [1,2]. These inhomogeneities reveal a number of different effects either of intrinsic nature or due to the presence of defects or impurities. Intrinsic effects are related to multiband properties of the Fermi surface in 2H-NbSe 2 [2], and to the pseudo gap of the High T c superconductors [1]. Impurities or defects produce electron scattering near the surface, whose energy and position dependence permits obtaining important information about local pair breaking effects or the bandstructure of the High T c superconductors [3,4].

One of the most promising future developments is the use of superconducting tips [5,6,7,8]. To see their advantages, let us write down the expression for the tunneling current vs. bias voltage I(V ) that is typically used to explain tunneling spectroscopy experiments [1]. As it is well known,

I(r 0 , V ) ∝ dE[f (E -eV ) -f (E)] ×N T (r 0 , E -eV )N S (r 0 , E)T (r 0 , E, eV )
where r 0 is the position of the tip apex on the sample's surface, f the Fermi function, N T and N S the densities of states of tip and sample respectively and T the tunneling matrix element. When the tip is normal the whole position and voltage dependence is usually assumed to be due to the superconducting density of states of the sample N S , and it is easy to find the generally accepted result of a tunneling conductance G(V ) = dI/dV that is directly proportional to the convolution between N S (r 0 , E) and the derivative of the Fermi function. When the tip is superconducting, its density of states strongly varies in the energy range of interest. For example, a junction at zero temperature between two isotropic BCS superconductors (note that in STM/S the insulating barrier is vacuum), will have zero current until the bias voltage reaches the sum of gaps of tip ∆ T and sample ∆ S (until V = ∆ T + ∆ S ). In the case of having a sample with a density of states differing from BCS theory, it is easy to de-convolute N S (E) directly from G(V ), provided N T (E) is well known. In addition, the high quasiparticle peaks at N T (∆ T ) may in some cases enhance small features of N S (E), making their measurement easier [9]. On the other hand, eventual local changes in the superconducting densities of states of the tip N T , produced, for instance, by a varying local magnetic field b(r 0 ) at the sample's surface, will produce variations in N T (r 0 , E) that can be used to measure b(r 0 ). Moreover, at finite temperature, the rounding of the Fermi edge produces additional peaks in G(V ) at V ≈ ∆ T -∆ S , whose height exponentially depends on temperature. Obviously, the Josephson effect, absent in studies made with normal tips, can be readily measured using STM with superconducting tips. Although the experiment is not straightforward, it opens many new possibilities [10,11,12].

Up to now, tips of Nb and Pb have been used to obtain important new insight into several superconducting materials [5,6,7,8,13]. Here we present first successful attempts in using superconducting tips of Al.

Experimental

We use a STM set-up in a dilution refrigerator which allows for scanning ranges up to two microns and has a resolution in energy of 15 µV [10,14]. The set-up has been tested previously with several measurements on different superconducting systems. The sample holder is located below the tip and can be moved in-situ by a mechanical pulling mechanism. We glue to the sample holder a sample of the same material as the tip (Al), a sample of Au to test the tunneling characteristics of the tip and a sample of 2H-NbSe 2 . The magnetic field is always applied parallel to the tip and perpendicular to the sample surface.

The procedure for obtaining the tips is described in Ref. [7]. In summary, we position a tip of Al on top of the Al sample, cut the feedback loop, and bring in a controlled way the tip into contact by moving the Z piezo. The surface oxide layer around the tip is easily removed by a repeated indentation procedure, in which the tip is indented up to one micron into the sample. At the same time, additional smaller sized (up to a hundred nm) pushing and pulling ramps are added to the Z piezo signal. In this way, the contacting region is mechanically annealed and shaped [15,16,6]. During the final pulling process, the formation of atomic size contacts is followed in detail. These contacts lead to a characteristic step like conductance vs. distance G(z) curves, shown in Fig. 1a. In the tunneling regime, when the contact is fully broken, G(z) has a clear exponential dependence (inset of Fig. 1a), from which a work function of several eV is extracted. Thousands of G(z) curves are recorded to make histograms that count the values of the conductance found during the experiment. The result is shown in Fig. 1b, and gives, as found previously [17], a bell shaped curve with a maximum at about 2/3 of the quantum of conductance G 0 =h/2e, and a strong decrease at G 0 . Subsequent smaller peaks are found at values close, but not equal, to multiples of G 0 . All these features, as well as the form of the curves shown in Fig. 1a, are related to the atomic orbital properties of the Al atoms [17,18]. Tips with different shapes can be built and moved in-situ to a normal sample for its characterization and then further to the superconducting sample under study. G(V ) curves obtained using this procedure on Au have been discussed previously, and are clean, with a negligible density of states close to the Fermi level, and fully explained by isotropic BCS theory using superconducting gaps that range between 0.175 mV and 0.2 mV [7,10].

Fig. 2. The magnetic field behavior of the tunneling conductance, normalized to the high bias value, of a tip of Al with a critical field of about 0.1 T. Data have been taken at zero field (triangles), and then from 0.02 T up to 0.1 T, in steps of 0.01 T.

Typically, Al tips remain superconducting up to critical fields H c,tip of some hundreds of Gauss. In Fig. 2 we show the tunneling characteristics obtained after having moved in-situ a tip to the normal sample of Au. Superconducting features are observed here up to H c,tip = 0.1 T, although the critical field of Al is of about 0.01 T. The same kind of enhancement of the critical field, has been observed previously in Pb [7]. Indeed, application of a magnetic field on superconducting tips leads to a curious system in which the superconducting properties strongly vary as a function of the position with respect to the tip apex [19]. Superconductivity is destroyed in the bulk by the applied field, but survives at the apex up to magnetic fields that can be two orders of magnitude higher than the bulk critical field [6,20]. As the lateral dimensions of the tip close to the apex are easily well below the London penetration depth, superconductivity remains there up to high fields, provided the tip is long and sharp enough so that the proximity effect of the bulk does not destroy superconductivity. H c,tip is found to strongly depend on the geometry of the tip [7,6], and highest values are obtained for the sharpest tips. Ginzburg-Landau calculations have carefully considered the way the order parameter varies when approaching the apex [21,22]. In Fig. 3 we show the tunneling characteristics obtained in 2H-NbSe 2 at zero field using a tip of Al. The tunneling conductance vs. bias voltage G(V) is at first sight rather similar to the result obtained using normal tips [23,24,25]. A closer look shows however that the quasiparticle peaks are slightly steeper, and located at somewhat higher bias voltages. The density of states of 2H-NbSe 2 is very easily de-convoluted from the experimental tunneling conductance data using a BCS expression with a gap of 0.2 mV for the tip (lower right inset of Fig. 3). When increasing temperature the curves remain qualitatively the same, with some additional broadening. Actually, the features expected in tunnel junctions between two different superconductors at ∆ T -∆ S produced by thermal broadening are not observed in our experiment, because superconductivity in Al is lost already at 1.2 K. This is a significant difference with respect to the results obtained in the same material using Pb tips in Ref. [9], where the T c 's of tip and sample are of the same order.

In the mixed state, at a field below H c,tip , we clearly observe vortices, as shown in Fig. 4, where we compare STS results taken at 0.03 T with a superconducting Al tip (left, its H c,tip is of 0.05 T), and results found when using a normal tip (Au, right). The G(V ) curves along a line into the vortex core are presented in Fig. 4a. Let us recall that, as shown previously in Refs. [23,24,25,26], the electronic density of states of a vortex in 2H-NbSe 2 is very peculiar, and reflects the localized electronic states that appear within vortex cores as well as the intrinsic anisotropic superconducting proper-ties of this material. A peak develops at the Fermi level density of states at the center of the vortex core that moves to higher voltages when leaving the core (right panel of Fig. 4a). Additional peaks appear at different voltage values and positions, and the STS images show a strong dependence on bias voltage. At the Fermi level (right panel of Fig. 4b), a characteristic star is observed, whose rays split into two around 0.3 mV (right panel of Fig. 4c) and lead to another star shape structure that is turned by about 30 • at 0.45 mV (right panel of Fig. 4d). When the measurement is with a tip of Al, the superconducting density of states of Al is superposed to the local density of states of 2H-NbSe 2 . Therefore, there is a sharp downturn of the tunneling conductance close to the Fermi level and inside the vortex core (left panel of Fig. 4a). STS images constructed with the tip of Al show qualitatively the expected features of the local density of states of 2H-NbSe 2 , namely the star shape, ray splitting and turning of the star (left panels of Fig. 4b, c and d respectively), at, however, somewhat higher bias voltages due to the gap in the density of states of the tip.

The variation of the magnetic field at the surface b(r 0 ) in and around vortices (along distances of the order of several tens of nm) at fields of several tens of mT is relatively small, of only around 10% -20%. Detection of the effect of such small changes in b(r 0 ) on the tunneling conductance of the tip is rather cumbersome, and requires the use of tips whose superconducting density of states N T (E) strongly depends on the magnetic field. Previous careful measurement of N T (E) on a normal sample is obviously crucial. In Fig. 5, we show data taken along a line into a vortex core at 0.04 T with an Al tip whose superconducting density of states shows a stronger field dependence than in the previous experiment. At the center, we observe the zero bias peak expected with a normal tip (top curves of Fig. 5). However, when leaving the center, the superconducting gap of Al appears around the Fermi level (bottom curves of Fig. 5 are not flat around zero bias, as expected when using a normal tip). The local magnetic field sensed by Fig. 5. Tunneling conductance as a function of the position taken with a tip of Al at 0.04 T. At the vortex core, we find that there is no gap structure corresponding to the tip (top curves of the figure), and that the gap of the tip slowly opens when leaving the vortex core (bottom curves of the figure). The data have been taken along a ray of the star shape vortex, as in previous Fig. 4a,i.e. along the arrow schematically shown in the bottom left inset. The tip acts as a local probe of the magnetic field, whose position dependence is shown in the right inset.

the tip is then easily computed by de-convoluting the Fermi level density of states N T (r 0 , E = 0) of Al from G(r 0 , V ) (this removes effects due to the density of states of 2H-NbSe 2 ) and using the previously measured dependence of N T (E = 0) as a function of the magnetic field. The result (inset of Fig. 5) has a large associated absolute uncertainty, which we estimate to be about 50% in the absolute value of the magnetic field. However, it gives first data of the magnetic field close to the vortex cores, and demonstrates the feasibility of measuring the local magnetic field b(r 0 ) with a superconducting Al tip as a probe. In addition, the overall order of magnitude of the variation agrees well with the expected variation of the magnetic field using reasonable values for the London penetration depth in this compound [27].

Summary and outlook

In summary, we have shown STM/S experiments on the surface of 2H-NbSe 2 , which demonstrate that Al tips are very valuable local probes. The vortex lattice is viewed with great precision, and the superconducting properties of the tip allow for the obtention of both its electronic and magnetic properties. Future developments should also lead to scanning spectroscopy studies of the local Josephson effect. Its observation on Al with STM was indeed already reported in Ref. [20], although further work is needed to exploit its possibilities.

On more general grounds, it is clear that increasing the palette of elements from which tips can be made is highly desirable. As a matter of fact, it is well known that tunneling conductance critically depends on the tip orbital states [28,29]. Therefore, different and complementary results should be obtained at the atomic level by using tips with significantly different available orbital states. The use of light elements, as Al, which only has s and p orbitals, as opposed to the more extended use of transition metal elements, could reveal unexpected insight into the properties of the surfaces under study.

Fig. 1 .

 1 Fig.1. In a. we show typical conductance vs. distance curves close to the atomic contact regime obtained between tip and sample of Al (line is for indentation and dashed line for retraction of the tip). The conductance is normalized to the quantum of conductance G 0 =h/2e. In the inset we show the conductance vs. distance obtained in the tunneling limit. In b. we show conductance histograms obtained from thousands of curves similar to the ones shown in a.

Fig. 3 .

 3 Fig. 3. In the upper left inset we show atomic resolution images made with a tip of Al at zero field. The atomic Se surface lattice of 2H-NbSe 2 is clearly viewed, as well as the charge density wave order characteristic of this material. In the figure we show the tunneling G(V ) (normalized to its high bias value) obtained in a sample of 2H-NbSe 2 with an Al tip at zero magnetic field. Convolution of the density of states shown in the lower right inset, which corresponds to the density of states of 2H-NbSe 2 typically measured with normal tips, with a simple BCS density of states with ∆=0.2 mV leads to the line in the figure, which closely follows experiment.

3 .

 3 STM/S in 2H-NbSe 2 with superconducting tips of Al.

Fig. 4 .

 4 Fig. 4. STS vortex imaging with a tip of Al (left), compared to imaging with a normal tip of Au (right), at 0.1 K. In a we show tunneling conductance curves obtained along a line of about 60 nm that goes into the vortex core, following one of the rays of the star shown in b (top to bottom is from the center to out of the vortex core). Below, we show images of a single vortex at different bias voltages, marked at the upper left corner of each panel. Color scale corresponds to the conductance normalized at high bias. In b, the characteristic star shape of the local density of states of a vortex is visualized. When increasing the bias voltage, the rays of the star split (c) and the star turns by 30 • (d).

Acknowledgement We specially acknowledge support from NES program, and discussions with V. Crespo and J.G. Rodrigo. This work was also supported by the Spanish MEC (Consolider and FIS-2004-02897 programs), by the Comunidad de Madrid through program "Science and Technology in the Millikelvin", and by ECOM. The Laboratorio de Bajas Temperaturas is associated to the ICMM of the CSIC.


Ilya Sheikin

Pierre Rodiere

Rikio Settai

Yoshichika Ōnuki

High-Field de Haas-van Alphen Effect in non-centrosymmetric CeCoGe 3 and LaCoGe 3

Keywords: CeCoGe3, LaCoGe3, non-centrosymmetric, de Haas-van Alphen effect

We report on de Haas-van Alphen effect measurements in the non-centrosymmetric systems CeCoGe3 and LaCoGe3 in magnetic field up to 28 Tesla. In both compounds, two new high frequencies were observed in high fields. The frequencies were not detected in previous lower field measurements. The frequencies do not originate from magnetic breakdown, and, therefore, are likely to be intrinsic features of the compounds. In CeCoGe3, the corresponding effective masses are strongly enhanced, being of the order of 30 bare electron masses.

Following the discovery of superconductivity in the noncentrosymmetric heavy-fermion compound CePt 3 Si, 1) materials without inversion symmetry in their crystal structure have attracted a lot of experimental and theoretical interest. The interest is mainly due to a fascinating theoretical prediction 2,3) that superconductive pairing in such systems requires an admixture of a spin-singlet with a spintriplet state. Subsequent observation of superconductivity in other non-centrosymmetric compounds CeRhSi 3 , 4) CeIrSi 3 , 5) CeCoGe 3 , 6,7) and CeIrGe 3 8) has further stimulated the research efforts. Very unusual superconducting properties were indeed observed in some systems, such as enormously high upper critical field in CeIrSi 3 9) and CeRhSi 3 . 10) The absence of inversion symmetry in the crystal lattice of a metal brings about a strong spin-orbit coupling, which in turn leads to the splitting of the electronic energy bands. The Fermi-surface of such a metal is, therefore, also split into two surfaces characterized by different chirality. This naturally results in the appearance of two distinct frequencies in the spectra of de Haas-van Alphen (dHvA) oscillations. It was demonstrated theoretically, 11) that the analysis of the oscillatory spectra in such materials provides a direct measure of the strength of the spin-orbit coupling. It is thus very important to obtain detailed and precise dHvA frequencies and their angular dependence in materials without inversion center. In most non-centrosymmetric compounds this coupling is extremely strong being of the order of 1000 K.

CeCoGe 3 and its non-4f analog LaCoGe 3 crystallize in the tetragonal BaNiSn 3 -type crystal structure. CeCoGe 3 is a moderate heavy-fermion system with a Sommerfeld coefficient of the specific heat γ = 0.032 J/molK 2 . It undergoes an untiferromagnetic transition at T N 1 = 21 K, followed by two more transitions at T N 2 = 12 K and T N 3 = 8 K. 12) Three consecutive metamagnetic transitions were observed for magnetic field applied along the easy magnetic c-axis. The high pressure phase diagram of CeCoGe 3 is quite complicated demonstrating six different phases. 13) Magnetic order vanishes around 55 kbar, and superconductivity is observed in the range of 54-75 kbar. * E-mail address: ilya.sheikin@grenoble.cnrs.fr

The dHvA effect in CeCoGe 3 and LaCoGe 3 was previously investigated in magnet fields up to 17 Tesla. 14) In both compounds, four fundamental frequencies were identified, all of them split due to spin-orbit interaction. The frequencies themselves and their angular dependencies are similar in both compounds and are in a rather good agreement with the results of theoretical band structure calculations performed for LaCoGe 3 . As compared to other non-centrosymmetric compounds, the splitting of dHvA frequencies in CeCoGe 3 and LaCoGe 3 is relatively small, implying a moderate spin-orbit coupling of about 100 K. In LaCoGe 3 , all the dHvA frequencies were doubled implying either a slightly unperfect quality of the sample or a field dependence of the frequencies. The highest effective mass observed in CeCoGe 3 is 12 bare electron masses corresponding to the β-branch. Finally, only one frequency originating from the α-branch representing the biggest Fermi surface was initially observed in CeCoGe 3 . 15) We present here the results of the dHvA effect measurements in non-centrosymmetric CeCoGe 3 and LaCoGe 3 in magnetic fields up to 28 Tesla produced by a resistive magnet in GHMFL. The measurements were performed on single crystals similar to those studied by Thamizhavel et al. 14) The details of crystals preparation and characterization are given elsewhere. 12) The measurements were performed using a torque cantilever magnetometer. The magnetometer was mounted in a top-loading dilution refrigerator equipped with a low-temperature rotation stage.

Figure 1 shows the oscillatory torque and the corresponding Fourier transform in LaCoGe 3 in field from 20 to 28 Tesla. All the previous dHvA frequencies 14) can be clearly identified, and their values are very close. The most remarkable difference with lower field measurements is the observation of two new fundamental frequencies denoted as A and B in figure 1. The new frequencies, F A = 14.1 kT and F B = 21 kT, are considerably higher than the highest frequency, F α = 9.15 kT, reported for lower fields.

The dHvA oscillations observed in CeCoGe 3 between 20 and 28 T are shown in the upper panel of figure 2. Like in LaCoGe 3 , all the fundamental frequencies previously observed in CeCoGe 3 are still present at high field (lower panel of figure 2). In addition, both components of the α-branch are now clearly resolved. With magnetic field applied at 9 • from the crystallographic c-axis, the two frequencies are close to each other, F α + = 10.17 kT and F α -= 10.47 kT. The most significant result, however, is the presence of the new frequencies, A and B, in the Fourier spectrum of CeCoGe 3 . The frequencies are somewhat lower than in LaCoGe 3 being F A = 12.31 kT and F B = 19.42 kT. Like in LaCoGe 3 , both new frequencies are fundamental. Figure 3 shows the temperature dependence of the dHvA amplitudes of the new frequencies A and B as well as α and β-branches in CeCoGe 3 for magnetic field applied at 9 • from the crystallographic c-axis, the same orientation as in figure 2. These data allow one to determine effective masses by fitting the experimental points to the temperature-dependent part of the Lifshits-Kosevich formula. 16) The best fits to the formula along with the extracted effective masses are also shown in figure 3. It was previously reported 14) that α and β-branches possess the highest effective masses of 8 and 12 bare electron masses respectively. These values are very close to the current results if only the α -frequency is considered and taking into account that only a single frequency from the α-branch was initially observed in previous measurements. Interestingly, the affective masses of the two frequencies originating from the α-branch differ by a factor of more than two, being 14.9 and 7.3 bare electron masses for α + and α -frequencies respec- tively. This is in contrast with all the other branches where the effective masses of the two components are quite close to each other. The most surprising result, however, is that the effective masses of the new frequencies are strongly enhanced being 27 and 31 bare electron masses for A and B respectively. These values by far exceed the masses of the other previously observed frequencies.

Figure 4 shows the angular dependence of the dHvA frequencies in CeCoGe 3 for the magnetic field rotation from [001] towards [100] direction. The frequencies were determined from the Fourier transform over the field range from 20 to 28 Tesla. The angular dependence of the frequencies observed at lower fields is very similar to the previously reported one. 14) The highest of the new frequencies, B, exists only over a small angular range at about 10 • from the c-axis. The other new frequency, A, however, is observed over an angular range of about 15 • and is centered around the c-axis.

It is not perfectly clear whether the new dHvA frequencies A and B observed both in LaCoGe 3 and CeCoGe 3 emerge only at high field above 17 T or exist at lower fields as well, but can not be detected due to experimental limitations. Very small amplitudes of the corresponding oscillations as well as strongly enhanced effective masses found in CeCoGe 3 make the latter scenario quite possible. In LaCoGe 3 , however, while the amplitudes of A and B oscillations are still relatively small, the effective masses must be small too. The new frequencies, if present, should therefore be possible to observe at lower field. An appealing possibility is a decrease of the Dingle temperature at high field. Since the oscillatory amplitude depends exponentially on the Dingle temperature, the oscillations might be too small to detect at low field, but become "visible" at higher field where the Dingle temperature is lower. Such scenario is indeed realized in CeCoIn 5 , [START_REF] Sheikin | [END_REF] where some of the dHvA frequencies are detected only at high field above an abrupt change of the Dingle temperature. Nonetheless, if the new frequencies exist, even if experimentally undetectable, already at low field, it would be difficult to explain their absence in theoretical band structure calculations at least for LaCoGe 3 , where neither magnetic order nor 4felectrons complicate the picture. Moreover, the Fermi-surface cross-sections corresponding to both frequencies do not exceed the area of the Brillouin zone perpendicular to the applied magnetic field. On the other hand, if the emergence of the new frequencies is intrinsically related to high magnetic fields, this would easily account for their absence in the band structure calculations performed for zero magnetic field. If this indeed is the case, the new frequencies certainly do not originate from a magnetic breakdown as there are no two fundamental frequencies that would sum up to yield any of them. Furthermore, the high effective masses observed in CeCoGe 3 also rule out such a possibility. Thus, the new frequencies are likely to be an intrinsic property of the compounds and seem to originate from a field-induced modification of the Fermisurface.

In conclusion, two new high dHvA frequencies have been observed both in LaCoGe 3 and CeCoGe 3 in high magnetic field up to 28 T. The frequencies are similar in both compounds. Neither of the new frequencies was detected in the previous lower-field measurements or revealed by theoretical band structure calculations. In CeCoGe 3 , the frequencies correspond to strongly enhanced effective masses implying that they represent thermodynamically important parts of the Fermi-surface. While it is not certain if the frequencies appear in high magnetic field only or simply experimentally undetectable at lower field, they certainly do not originate from magnetic breakdown. They are therefore likely to be intrinsic for both compounds and are possibly due to a field-induced modification of the Fermi-surface.

Part of this work has been supported by the EuroMagNET II under the EU contract number 228043.

Fig. 1 .

 1 Fig. 1. High field (20-28 T) dHvA oscillations (upper panel) and its Fourier spectrum (lower panel) in LaCoGe 3 for magnetic field applied at 2 • from the c-axis at 60 mK. The frequencies that were observed in the previous lower-field measurements are denoted by Greek letters, while the new frequencies observed at high field only are denoted A and B.

3 Fig. 2 .

 32 Fig. 2. dHvA oscillatory signal (upper panel) and its Fourier spectrum (lower panel) observed in CeCoGe 3 with magnetic field (20-28 T) applied at 9 • from [100] to [110] at T = 60 mK. The frequencies denoted θ, ǫ, β and α were observed in the previous lower field measurements (with the exception of α + -branch). The frequencies denoted A and B are observed only at high field and were not detected in the previous measurements.

Fig. 3 . 3 Fig. 4 .

 334 Fig. 3. (Color online) Temperature dependence of the dHvA amplitude is shown for α and β-branches as well as the new frequencies A and B of CeCoGe 3 with magnetic field applied at 9 • from c to a-axis. Lines are the fits to the temperature dependent part of the Lifshits-Kosevich formula A(T ) ∝ αm * T /B sinh(αm * T /B) , where α ≈ 14.69 T/K. The effective masses, m * , obtained from the fits are also shown. m 0 is the bare electron mass.


Ken Arroyo Ohori

Guillaume Damiand

Hugo Ledoux

Constructing an n-dimensional cell complex from a soup of (n-1)-dimensional faces

There is substantial value in the use of higher-dimensional (>3D) digital objects in GIS that are built from complex real-world data. This use is however hampered by the difficulty of constructing such objects. In this paper, we present a dimension independent algorithm to build an n-dimensional cellular complex with linear geometries from its isolated (n -1)-dimensional faces represented as combinatorial maps. It does so by efficiently finding the common (n -2)-cells (ridges) along which they need to be linked. This process can then be iteratively applied in increasing dimension to construct objects of any dimension. We briefly describe combinatorial maps, present our algorithm using them as a base, and show an example using 2D, 3D and 4D objects which was verified to be correct, both manually and using automated methods.

Introduction

Higher-dimensional digital objects represent well-defined extents of space in arbitrary dimensions. In geographic information systems (GIS), these can be generated when 2D/3D space, time [START_REF] Peuquet | Representations of Space and Time[END_REF], scale [START_REF] Peter Van Oosterom | Towards a true vario-scale structure supporting smooth-zoom[END_REF], semantics [START_REF] Baglatzi | On the formulation of conceptual spaces for land cover classification systems[END_REF], or others are all treated as independent axes of a coordinate system. An example of such an object is presented in Fig. 1. This form of representation offers interesting advantages compared to traditional ones where multiple representations of the same object are stored separately and linked in an ad hoc fashion, such as conceptual simplicity, immediate access to all existing topological relationships, the possibility to represent complex events like motion and the ease of maintaining consistency between objects [START_REF] Stoter | 5D modeling -applications and advantages[END_REF].

However, creating higher dimensional digital objects from real-world data is inherently difficult on multiple levels. Since we are usually only familiar with up to 3D physical space, describing such objects might be straightforward mathematically, but it is nonetheless unintuitive. Higher dimensional data models are also complex, and thus realising even very simple objects requires a large number of operations on abstract elements. Finally, manipulating the related data structures while ensuring that all its required references are correctly kept is already non-trivial in 3D [START_REF] Mäntylä | An introduction to solid modeling[END_REF], and increasingly difficult in higher dimensions.

Nonetheless, it is possible to use the concepts behind boundary representation to significantly reduce the difficulty of the problem. In boundary representation (b-rep or BREP), an n-dimensional closed object can be unambiguously described by the (n -1)-dimensional boundary that encloses it, as originally defined by Baumgart [START_REF] Bruce | A polyhedron representation for computer vision[END_REF] and Braid [START_REF] Braid | The synthesis of solids bounded by many faces[END_REF] for a (3D) solid with a (2D) surface boundary composed of flat polygonal patches. This concept is valid in higher dimensions as well, and is related to the concept of a cell complex [START_REF] Hatcher | Algebraic Topology[END_REF] in topology, where an n-dimensional cell (n-cell) in the complex has a number of (n -1)-cells (faces) as its boundary, and these faces are also part of the complex. An n-cell is an abstract object which is considered to be homeomorphic to an n-ball (e.g. point, segment, disk, ball, etc.). In this paper we use combinatorial maps to represent a cell complex, which are described in Sec. 2.

The (n-1)-dimensional boundary of an n-cell is much easier to conceive than the original n-cell, since the (n-1)-cells that it is composed of can be themselves described individually. However, this requires the existence of an algorithm that is able to connect these separate (n -1)-cells to form the n-cell in an efficient manner, abstracting such issues as incompatible orientations in the model, the handling of duplicate cells and the identification of common boundaries. This operation, fully dimension independent, presented in Sec. 3 and the focus of the present paper, can be performed recursively in increasing dimension to generate arbitrary cell complexes in any dimension, and thus we refer to it as incremental construction.

Another possible use of incremental construction is to generate the topological information, i.e. incidence and adjacency, between a set of existing objects. This is in fact simply a subset of what the problem incremental construction is meant to solve-the identification of common boundaries-, and fits very well within the frame of GIS, where data models tend to contain very limited topological information but topological queries are of great importance [START_REF] Egenhofer | Point-set topological spatial relations[END_REF]. Some GIS models, like the OGC Simple Features Specification [START_REF]OpenGIS Implementation Specification for Geographic Information -Simple Feature Access -Part 1: Common Architecture[END_REF] have no topology in their structure, even repeating the coordinates of individual points when these appear in multiple line segments or polygons. Others, such as CityGML [START_REF]OGC City Geography Markup Language (CityGML) Encoding Standard[END_REF], only have implicit topological information (e.g. these surfaces should form a closed shell) which is often unenforced in their geometry.

We have implemented our algorithm based on the CGAL Combinatorial Maps and Linear Cell Complex packages, which are described in Sec. 4 together with the details of our implementation, including a discussion of the computational complexity of our approach. In Sec. 5, we present an example in 4D that shows this approach in practice, comparing its output with correct results which have been generated by the Linear Cell Complex package and verifying it analytically. We finish with conclusions and discussion in Sec. 6.

Combinatorial maps

Combinatorial maps (or simply maps) are an ordered topological model originally proposed by Edmonds [START_REF] Edmonds | A combinatorial representation of polyhedral surfaces[END_REF] to describe the 2D surfaces of 3D objects. Their extension to arbitrary dimensions is described by Lienhardt [START_REF] Lienhardt | N-dimensional generalized combinatorial maps and cellular quasi-manifolds[END_REF] for objects without boundaries and extended by Poudret et al. [START_REF] Poudret | Cartes combinatoires ouverts[END_REF] to objects with boundaries. They are able to describe subdivisions of orientable quasi-manifolds 3 .

Intuitively, a combinatorial map is composed of two elements: darts and relations between them (β). The precise definition of a dart is related to an underlying simplicial decomposition of the object, each dart being equivalent to a simplex in it. However, intuitively it can be seen as an oriented edge on the boundary of a facet, which itself is on the boundary of a solid, which itself is on the boundary of a 4-cell, and so on. It is therefore equivalent to the half-edge data structure in 2D, but extends naturally to higher dimensions. Meanwhile, the relations are functions connecting darts that are related along a certain dimension. In this manner, β 1 joins consecutive oriented edges within a facet forming a loop, β 2 joins adjacent facets within a solid, β 3 joins adjacent solids within a 4-cell, and so on. As in other models based on directed elements, β ijoined darts for i > 1 have opposite orientations.

More formally, an n-dimensional combinatorial map (or n-map) is defined by an (n + 1)-tuple M = (D, β 1 , . . . , β n ) where D is a finite and non-empty set of darts, β 1 is a partial permutation on D (a function f :

D ∪ {∅} → D ∪ {∅} such that ∀d 1 ∈ D, ∀d 2 = d 1 ∈ D, f (d 1 ) = ∅ and f (d 2 ) = ∅ ⇒ f (d 1 ) = f (d 2 )), ∀2 ≤ i ≤ n, β i is a partial involution (a partial permutation f such that ∀d ∈ D, f (d) = ∅ ⇒ f (f (d)) = d), and ∀1 ≤ i < i + 2 ≤ j ≤ n, β i • β j is
also a partial involution (see an example in Fig. 2).

Here, ∅ is a special value used to indicate that a given dart d has no other dart in relation by a given β i . In such a case we have β i (d) = ∅ and that means that d belongs to the i-boundary of the described object (it belongs to only one i-cell). A dart d is said to be i-free when β i (d) = ∅. Otherwise it is i-sewn with a second dart d ′ and we have

β i (d) = d ′ = ∅.
In order to traverse an n-map, the orbit operator < A > (d) =< β a1 , . . . , β ak > (d) obtains all the darts that can be reached from dart d by successive applications of the links β a1 , . . . , β ak ∈ A. Certain orbits are particularly interesting: for any 1

≤ i ≤ n, < ✁ ✁ β i > (d) =< β 1 , . . . , β i-1 , β i+1 , . . . , β n > (d) contains all the darts in the i-cell of d, while < {β i • β j |∀i, j 1 ≤ i < j ≤ n} > (d) contains
all the darts in the 0-cell (vertex) of d. As two darts linked by a β i have opposite directions (for 2 ≤ i ≤ n), they belong to the two vertices of a same edge. Thus by combining two β, we obtain a dart of the same vertex than d.

When the relations in an orbit are applied in a well-defined order, these orbits are generated in a consistent manner. This makes it possible to generate a canonical representation of (a subset of) the darts in a map, which combined with labelled darts can be used to test combinatorial map isomorphism in quadratic time, as shown by Gosselin et al. [START_REF] Gosselin | Efficient search of combinatorial maps using signatures[END_REF] and demonstrated by searching for patterns in images.

The incremental construction of cells in a combinatorial map is based on the sewing operator, which joins two i-cells along an (i -1)-cell which after the operation lies in their common boundary. Intuitively, given two i-free darts d 1 and d 2 , the i-sew operation between d 1 and d 2 will pairwise link the orbits of these two darts by β i so that we will obtain

β i (d 1 ) = d 2 .
We can see in Fig. 3 an example of the 3-sew operation. Starting from a 3D combinatorial map describing two isolated tetrahedra, we identify two faces by using the 3-sew operations on darts d and d ′ . This operation puts in relation all the darts of the two initial faces by pairs so that we obtain a valid combinatorial map (i.e. the constraint that β 1 • β 3 is a partial involution is still satisfied).

A combinatorial map only describes the topological part of an object in term of a cell complex, i.e. the set of cells in all dimensions and all the incidence and adjacency relations. Applications often require adding information associated with specific cells (e.g. to associate a colour to each vertex, or a normal to each face). This is possible thanks to the attribute notion. An i-attribute is the information associated with i-cells. As cells are implicitly represented by sets of darts in combinatorial maps, links between i-attributes and i-cells are done through the darts of the i-cells: all the darts belonging to a same i-cell are linked to the same i-attribute.

These attributes are very useful as they allow to associate any information to any cell, and given a dart, we have a direct access to all of its associated attributes. Moreover, these attributes can be used to describe the geometry of the objects. Indeed, we can associate to each vertex of a combinatorial map a point in R d2 by using 0-attributes. A combinatorial map with this type of embedding gives a linear cell complex. Indeed, in this case, the geometry of each edge is a segment, the geometry of each face is a planar polygon, and so on. d2 is the dimension of the geometry, i.e. the dimension of the ambient space. Generally we have d2 ≥ d (d being the combinatorial dimension). For example we can use d = d2 = 2 to describe a planar graph embedded in a plane or d = 2 and d2 = 3 to describe a polyhedral mesh embedded in R 3 . Now given two i-cells in a linear cell complex, we can test if these two cells are identical, i.e. if they have both the same topology and the same embedding information. To test if they have the same topology, we use the technique introduced above for combinatorial map isomorphism, by considering only the two i-cells instead of two whole combinatorial maps. Testing if the two cells have the same embedding can be done during the isomorphism test, simply by testing if two darts d and d ′ considered by the isomorphism function are linked with two points with the same coordinates. Note that we have to consider the two possible orientations of one of the two i-cells in order to detect also if the two cells are identical but with reverse orientations. We make use of this technique to efficiently test for the identity of two cells in Sec. 3.

Incremental Construction

Since an n-cell in a cell complex can be described by the (n -1)-cells on its boundary, the same thinking can be applied in reverse: a yet unbuilt n-cell can be constructed based on a set of (n-1)-cells which are known to form its complete (closed) boundary.

The algorithm is applied object by object in increasing dimension, constructing isolated 0-cells first, and continuing through 2-cells, 3-cells and further. Our explanation follows the construction of single cell of a certain dimension, starting with unique isolated vertices in Sec. 3.1. 1-cells are skipped since 2-cells can be easily described as a succession of 0-cells. 2-cells are then built from an ordered sequence of 0-cells, as covered in Sec. 3.2. Finally, for n-cells (n > 2), the method receives an unordered set (soup) of (n -1)-cells, the geometries of which are used as faces of the finished n-cell that is returned, as explained in Sec. 3.3. These individual cell creation algorithms are applied object by object in order of increasing dimension, so that the lower dimensional cells generated as the output of earlier stages can be used as the input of latter stages.

The incremental construction algorithm as a whole keeps track of already built cells by maintaining a reference to one of their darts, making sure that no identical cells (with equal geometry and topology) are ever created. For efficiency reasons, it is convenient to use darts which are embedded into the lexicographically smallest points of each cell, which combined with smallest point indices for all existing n-and (n -1)-cells, and the (n -2)-cells on the boundary of the (n -1)-cells to be used, greatly accelerates this process.

Vertices (0D)

A single point, as defined by a unique tuple of coordinates, can be present in multiple higher dimensional cells and thus described multiple times. However, to ensure the correct generation of topology and enforce the orientable quasimanifold criterion, it should only be created once. When doing so, a new point embedding at that location should be created, and a new free dart representing the 0-cell should be created and linked to it. Every new instance of the same point should then link to this original dart.

The output of this step is thus a map from each input point to a 0-cell (dart embedded into a point) at that location (see an example in Fig. 4(a)).

Facets (2D)

In order to create a 2-cell from a sequence of 0-cells, three general steps are needed:

1. The unique 0-cells, as resulting from the evaluation of the 0-cells' point embeddings in the map obtained when processing all 0-cells, are obtained. The result of this evaluation might be the same 0-cell as provided, or an already existing 0-cell at that location. Each of these 0-cells might be a single 1-free dart, in which case it can be used directly (see an example in Fig. 4(b)), or it can be a non 1-free group of darts, in which case it has already been used as part of a different 1-cell (and possibly other higher dimensional cells). These darts used in 1-cells can be reused only when they would become part of the 2-cell that will be built and are not part of any 2cell. If the darts cannot be reused, a copy of them with opposite orientation has to be created, linking it to the same embedding as the original (see examples in Fig. 4(c) to (h). In all these cases, there is at least one dart which is duplicated). The opposite orientation ensures that the two (old and new) can then be 2-sewn when constructing a 3-cell in a subsequent step. 2. The darts obtained in the previous step are 1-sewn sequentially, and the last is 1-sewn to the first, forming a closed loop. 3. Just as in the creation of 0-cells, a 2-cell can be on the boundary of multiple higher dimensional cells (or simply described multiple times in the input), and as such, it is necessary to check if the 2-cell has already been created.

If any comparison returns that the 2-cell already exists, the newly created darts are deleted 4 and the existing one is returned, otherwise the new 2-cell is returned.

3-cells and higher (dimension independent)

The method to create n-cells from their (n -1)-cell boundaries is identical for all n > 2, allowing a dimension-independent function to be created. As with the creation of 2-cells from 0-cells, it consists of three general steps:

1. First of all, whether each (n -1)-cell has already been created beforehand is checked. This is meant so that multiple identical (n -1)-cells are never created in the final cell complex, even when they are given as input. If an (n-1)-cell already exists and it is (n-1)-free, it is reused as part of the n-cell (see an example in Fig. 5(a)). If it exists but is already part of a different ncell, it is duplicated with reverse orientation (see an example in Fig. 5(b) for face labeled d). As in the case of 2-cells, this is done so that it has the same geometric embeddings and attributes, but its opposite orientation ensures that the two (old and new) can be directly n-sewn together. 2. The (n -1)-cells (faces) are (n -1)-sewn along their common (n -2)dimensional boundaries (ridges). If two groups of connected (n -1)-cells with incompatible orientations would be joined by this operation (i.e. the pair of two corresponding ridges have the same orientation), the orientation of one of the groups is reversed before the link is created. If more than one match for a ridge is found, the object being represented is not a quasi-cellular manifold, and thus cannot be represented using combinatorial maps. 3. The newly constructed n-cell is finally compared to other n-cells to check if it already exists. If an n-cell is found to exist, the algorithm should delete the darts that are part of the newly created n-cell and instead return the existing cell. This ensures that only a single instance of an n-cell is created.

Implementation and Complexity

We have implemented the incremental construction algorithm in C++ using the Combinatorial Maps and Linear Cell Complex packages in CGAL 5 . In order to improve the performance of the incremental construction algorithm, we use some indices that map the lexicographically smallest point embedding of some cells of a given dimension to a dart embedded at that location. These indices are implemented as C++ Standard Library 6 maps with point embeddings as keys and lists of darts as values, using a custom compare function so that the points are internally sorted in lexicographical order. Because std::map is normally implemented as a self-balancing binary search tree, O(log n) search, insertion and deletion times and O(n) space can be expected. Since we create objects dimension by dimension, it is not necessary to maintain indices for all the cells of all dimensions at the same time. The only ones used are: all n-and (n -1)-cells, and the (n -2)-cells on the boundary of the (n -1)-cells for that step. Most of these can be built incrementally, adding new cells as they are created in O(log c), with c the number of cells of that dimension, assuming that the smallest vertex and a dart embedded there are kept during its construction. The complexity of building any index of cells of any dimension is thus O(c log c) and it uses O(c) space. Note that this also gives the computational complexity of creating a map of all unique 0-cells in the cell complex.

Checking whether a given cell already exists in the cell complex is more complex. Finding a list of cells that have a certain smallest vertex is done in O(log c). Theoretically, all existing cells in the complex could have the same smallest vertex, leading to up to c quadratic time cell-to-cell comparisons just to find whether one cell exists. However, every dart is only part of a single cell of any given dimension, so while every dart could conceivably be a starting point for the identity comparison, a single dart cannot be used as a starting point in more than one comparison, and thus a maximum of d complex identity comparisons will be made for all cells, with d complex the total number of darts in the cell complex. From these d complex darts, two identity comparisons are started, one assuming that the two cells (new and existing) have the same orientation, and one assuming opposite orientations. Each of these involves a number of dart-to-dart comparisons in the canonical representations that cannot be higher than the number of darts in the smallest of the two cells. The number of darts in the existing cell is unknown, but starting from the number of darts in the newly created cell (d cell ), it is safe to say that no more than d cell dart-to-dart comparisons will be made in each identity test, leading to a worst-case time complexity of O(d complex d cell ). Note that this is similar to an isomorphism test starting at every dart of the complex.

Finally, creating an n-cell from a set of (n -1)-cells on its boundary is more expensive, since the (n -2)-cell (ridge) index needs to be computed for every ncell. Following the same reasoning as above, it can be created in O(r log r) with r the number of ridges in the n-cell, and uses O(r) space. Checking whether a single ridge has a corresponding match in the index is done in O(d cell d ridge ), with d cell the number of darts in the n-cell and d ridge the number of darts in the ridge to be tested. Since this is done for all the ridges in an n-cell, the total complexity of this step, which dominates the running time of the algorithm, is The analyses given above give an indication of the computational and space complexity of the incremental algorithm as a whole. However, it is worth noting that in realistic cases the algorithm fares far better than in these worst-case scenarios: the number of cells that have a certain smallest vertex is normally far lower than the total number of cells in the complex, most of their darts are not embedded at the smallest vertex, and from these darts most identity comparisons will fail long before reaching the end of their canonical representation.

Finally, one more nuance can affect the performance of this approach. We have discussed that when two groups of darts with incompatible orientations have to be joined, the orientation of one of these has to be reversed. This is easily done by obtaining all the connected darts of one of the groups, preferably the one that is expected to be smaller, and reversing their orientation 2-cell by 2-cell. Every dart d in a 2-cell is then 1-sewn to the previous dart in the polygonal curve of the 2-cell (β -1 1 (d)). A group of n darts can then have its orientation reversed in O(n) time. This is not a problem in practice since GIS datasets generally store nearby objects close together, but if a cell complex is incrementally constructed in the worst possible way, i.e. creating as many disconnected groups as possible, this could have to be repeated for every cell of every dimension.

Fig. 1 .

 1 Fig. 1. The Aula Congress Centre in the TU Delft campus represented as extruded 3D and 4D models.

Fig. 2 .

 2 Fig. 2. (Left) A 3D cellular complex composed with two tetrahedra sharing a common face. There are 2 3-cells, 7 2-cells, 9 1-cells and 5 0-cells. (Right) The corresponding 3D combinatorial map having 48 darts.

Fig. 3 .

 3 Fig. 3. (Left) A 3D combinatorial map describing 2 isolated tetrahedra. (Right) The 3D combinatorial map obtained after 3-sewing darts d and d ′ . All the darts of the two initial faces are 3-sewed by pairs.

Fig. 4 .

 4 Fig. 4. Illustration of the different steps of the reconstruction of 2-cells. (a) Initial configuration: one dart per vertex. (b) After a = make 2 cell(1, 2, 3). (c) After b = make 2 cell(2, 4, 3). (d) After d = make 2 cell(1, 4, 3). (e) After c = make 2 cell(1, 4, 2). (f) After e = make 2 cell(1, 3, 5). (g) After f = make 2 cell(5, 3, 4). (h) Final result, after g = make 2 cell(4, 5, 1).

Fig. 5 .

 5 Fig. 5. Illustration of the two steps of the reconstruction of 3-cells. We start from the combinatorial map given in Fig. 4(h) which is the result of the reconstruction of 2-cells. (a) After make 3 cell(a, b, c, d). (b) Final result, after make 3 cell(d, e, f, g).

  ridges O(d cell d ridge ) = O(d 2 cell ).

A specific combinatorial interpretation of the concept of a manifold. See Lienhardt[START_REF] Lienhardt | N-dimensional generalized combinatorial maps and cellular quasi-manifolds[END_REF] for details.

Since the 2-cell already exists, all the used darts are copies of existing ones. This deletion thus does not erase any unique instance of a 0-or 1-cell.

The Computational Geometry Algorithms Library: http://www.cgal.org

For instance, the GNU Standard C++ Library: http://gcc.gnu.org/libstdc++/

Acknowledgments

This research is supported by the Dutch Technology Foundation STW, which is part of the Netherlands Organisation for Scientific Research (NWO), and which is partly funded by the Ministry of Economic Affairs (Project code: 11300).

Example

The CGAL Linear Cell Complex package provides functions to generate a series of primitives which are known to be created with correct geometry and topology, and can then be sewn together to generate more complex models. We have therefore created various 2D, 3D and 4D cell complexes using both these functions and our approach. In this manner it was possible to test the validity of our models using the identity comparison described in Sec. 2, as well as manually verifying all β-links.

In the following we show an interesting case, how a tesseract (see Fig. 6) can be generated using our approach. A tesseract is the 4D analogue of a cube, and is a 4-cell bounded by 8 cubical 3-cells, each of which is bounded by 6 square 2-cells. It thus consists of one 4-cell, 8 3-cells, 24 2-cells, 32 1-cells and 16 0-cells.

Using our approach, an empty 0-cell index is first created. Then, the 16 vertices of the tesseract, each vertex p i described by a tuple of coordinates (x i , y i , z i , w i , can be created as p i = make 0 cell(x i , y i , z i , w i ), which returns a unique dart embedded at each location, and added to the 0-cell index. At this point, the algorithm would have built an unconnected cell complex consisting solely of 16 completely free darts.

An empty index of 2-cells is then created. Each of its 24 square facets can be built based on its vertices as f i = make 2 cell(p j , p k , p l , p m ), which 1-sews (copies of) these darts in a loop and returns the dart embedded at the smallest vertex of the facet. These are added to the index of 2-cells. Since every vertex is used in 6 different 2-cells, each dart would be copied 5 times. The cell complex at this point thus consists of 24 disconnected groups of 4 darts each.

Next, an empty index of 3-cells is created and the index of 0-cells can be deleted. For each of the 8 cubical 3-cells, a function call of the form v i = make 3 cell(f j , f k , f l , f m , f n , f o ) is made. At this point, an index of the 1D ridges of each face is built, which is used to find the 12 pairs of corresponding ridges that are then be 2-sewn together. When a 3-cell is created, it is added to the index. Since every facet bounds two 3-cells, each dart is duplicated once again, resulting in a cell complex of 8 disconnected groups of 24 darts each.

Finally, the tesseract is created with the function t = make 4 cell(v 1 , v 2 , . . . , v 8 ). This can use the index of 2-cells to find the 24 corresponding pairs of facets that are then 3-sewn to generate the final cell complex.

We tested the validity of this object by performing a series of tests on the structure (complete and symmetric sewing), testing whether each cube was identical to the expected outcome, and manually verified the β-links of its 192 darts.

Conclusions and future work

We have shown that it is possible to apply the fundamental concept of boundary representation, describing an n-cell by its (n-1)-dimensional faces, to incrementally construct cell complexes of any dimension. To the best of our knowledge, this technique is the only one that has been described and/or implemented for 4D cell complexes or higher. Using a variety of indices is efficient, generating an n-cell in O(d 2 ) in the worst case, with d the total number of darts in the cell, and our algorithm should fare markedly better in realistic datasets.

We intend to use this approach, supplemented with techniques under development that generate the required (n -1)-dimensional faces, to generate objects for higher dimensional geographic information systems, as well as other applications. This will allow us to take real-world 3D city models and incorporate additional dimensions to them, such as time and scale, to create 4D/5D objects to which higher dimensional analyses can be performed.