Emmanuel Pietriga
email: emmanuel.pietriga@inria.fr

Fernando Del Campo

Amanda Ibsen

Romain Primet

Caroline Appert

Olivier Chapuis

Maren Hempel

Roberto Muñoz

Susana Eyheramendy

Andres Jordan

Hervé Dole

Exploratory Visualization of Astronomical Data on Ultra-high-resolution Wall Displays

Keywords: Ultra-walls, Visualization, FITS Images, Large Datasets, Astronomical Database Queries

Ultra-high-resolution wall displays feature a very high pixel density over a large physical surface, which makes them well-suited to the collaborative, exploratory visualization of large datasets. We introduce FITS-OW, an application designed for such wall displays, that enables astronomers to navigate in large collections of FITS images, query astronomical databases, and display detailed, complementary data and documents about multiple sources simultaneously. We describe how astronomers interact with their data using both the wall's touchsensitive surface and handheld devices. We also report on the technical challenges we addressed in terms of distributed graphics rendering and data sharing over the computer clusters that drive wall displays.

INTRODUCTION

Ultra-high-resolution wall-sized displays, also called ultra-walls [START_REF] Nancel | Mid-air pointing on ultra-walls[END_REF] for short, feature a very high pixel density over a large physical surface. For instance, WILD, the first wall display we set up in our laboratory, has a total resolution of 20480 × 6400 = 131 megapixels for a surface area of 5.5m × 1.8m. Other, even larger, platforms include Stallion 2 (328 megapixels) and the Reality Deck [START_REF] Papadopoulos | The reality deck-an immersive gigapixel display[END_REF] (1,500 megapixels in an immersive configuration).

Ultra-walls have properties that make them well-suited to the visualization of very large datasets. They can represent the data with a high level of detail while at the same time retaining context: users can transition from an overview of the data to a detailed view simply by physically moving [START_REF] Ball | Move to improve: Promoting physical navigation to increase user performance with large displays[END_REF] in front of the display. Wall displays also offer good support for collaborative work, enabling multiple users to simultaneously visualize and interact with the displayed data. [START_REF] Beaudouin-Lafon | Multisurface Interaction in the WILD Room[END_REF][START_REF] Jakobsen | Up close and personal: Collaborative work on a high-resolution multitouch wall display[END_REF][START_REF] Liu | Shared interaction on a wall-sized display in a data manipulation task[END_REF] In this paper, we introduce FITS-OW, an application that enables astronomers to visualize and interact with very large FITS images and collections thereof. They can pan and zoom in images that are several hundred thousand pixels in both width and height, overlay the results of data analyses, fetch and display additional images of a specific object or region in the sky, showing observations in different ranges of the electromagnetic spectrum or made at different times. Additionally, FITS-OW lets astronomers query databases such as SIMBAD [START_REF] Wenger | The SIMBAD astronomical database -The CDS reference database for astronomical objects[END_REF] servers, and visualize the results of such queries in-place, right next to the corresponding source in the image. The very high pixel density of wall displays means that detailed information can be shown for multiple sources simultaneously, including multiple measurements as well as research articles retrieved dynamically through links found in relevant databases.

Figure 1. Zooming in GLIMPSE360, Spitzer's Infrared Milky Way image (540,000 × 15,000 pixels), on WILDER, a touch-sensitive ultra-wall at Inria Saclay, France, that features a resolution of 14,400 × 4,800 pixels.

We explain how we enable astronomers to perform these operations using interaction techniques that were designed specifically for wall displays, using direct manipulation and gestures performed on the wall's surface or on handheld tablets: adjusting the scale and color mapping used to render the raw, high-dynamic-range FITS data; overlaying and manually compositing multiple bands; measuring the brightness and colors of astronomical sources, displaying their light curve; making queries and adjusting object filters. We also explain how we addressed related technical challenges. Indeed, ultra-walls are often driven by clusters of computers (the abovementioned WILD platform uses 32+1 graphics processing units in 16+1 computers) and involve heterogeneous input devices, which causes problems of data sharing, graphics rendering, and handling multiple user input channels. We describe FITS-OW's architecture and our solution to the specific challenges that this application raised: the generation of FITS tile pyramids and their multi-scale rendering; the computation of sky coordinates; queries to sky catalogs; the dynamic adjustment of scale, color mapping and graphics compositing settings, and the framework for managing user interaction with various input sources.

BACKGROUND AND MOTIVATION

Multiple studies in the field of Human-Computer Interaction have been evaluating the potential benefits of very-high-resolution displays for visualizing very large datasets and very large imagery. One of the most recent examples to date is Ruddle et al.'s work 9 on the design and evaluation of high-end workstations for visualizing gigapixel images in the field of medical pathology. The display capacity of such dual-screen workstations remains necessarily limited, however, and researchers have also been investigating other types of interactive surfaces that feature a much higher display capacity: ultra-walls. Ultra-walls 1 are tiled displays, typically built by juxtaposing LCD or LED monitors, that feature a pixel density of 60-to-100 pixels per inch over a large physical surface, 14-to-20 feet in width and about 6 feet in height (Figure 1).

Early work on ultra-walls has mainly focused on the technical aspects of how to make such platforms: 10 how to display complex graphics, how to stream data across the nodes of the computer clusters that drive them. Now that the technical feasibility of such platforms has been established, researchers are now investigating other questions: evaluating the actual benefits of ultra-walls; [START_REF] Liu | Effects of display size and navigation type on a classification task[END_REF][START_REF] Shupp | Shaping the display of the future: The effects of display size and curvature on user performance and insights[END_REF][START_REF] Yost | Beyond visual acuity: The perceptual scalability of information visualizations for large displays[END_REF] studying co-located users' collaboration strategies, [START_REF] Jakobsen | Up close and personal: Collaborative work on a high-resolution multitouch wall display[END_REF][START_REF] Liu | Shared interaction on a wall-sized display in a data manipulation task[END_REF] or users' perception of visual variables based on their relative position in front of the display; [START_REF] Bezerianos | Perception of visual variables on tiled wall-sized displays for information visualization applications[END_REF] designing interaction techniques that support physical navigation, by coupling ultra-walls with input devices such as touch frames, motion-tracking systems and wireless multitouch devices, in order to enable multi-device and multiuser interaction with the displayed data. Examples of such techniques include high-precision remote pointing techniques to enable users to interact when they are not within arm's reach of the display, 1 or mid-air panning and zooming techniques for navigating maps, images and datasets that exceed the display capacity of ultrawalls. [START_REF] Nancel | Mid-air pan-and-zoom on wall-sized displays[END_REF] For instance, Figure 1 shows a small portion of GLIMPSE360, Spitzer's Infrared Milky Way image. [START_REF] Whitney | Glimpse360: Observing The Outback Of The Galaxy[END_REF] Even with the ultra-wall's display capacity, the viewport can still only accommodate 0.85% of the full image at maximum magnification, calling for multi-scale navigation techniques [START_REF] Pietriga | Pointing and beyond: An operationalization and preliminary evaluation of multi-scale searching[END_REF] to enable users to efficiently navigate in the data.

Application areas for ultra-walls are quite varied, ranging from the monitoring of complex infrastructures [START_REF] Prouzeau | Road traffic monitoring on a wall display[END_REF] and crisis management situations [START_REF] Chokshi | ePlan Multi-Surface: A Multi-Surface Environment for Emergency Response Planning Exercises[END_REF][START_REF] Onorati | WallTweet: A Knowledge Ecosystem for Supporting Situation Awareness[END_REF] to tools for exploratory data visualization, [START_REF] Tsandilas | Sketchsliders: Sketching widgets for visual exploration on wall displays[END_REF] visual analytics, 22, 23 medicine 24 and other scientific disciplines. [START_REF] Beaudouin-Lafon | Multisurface Interaction in the WILD Room[END_REF] In each of these application areas, subject-matter experts are faced with one common problem: the need to visualize and make sense of large amounts of data. The data can be fairly homogeneous but quite large, or they can be collections of more heterogeneous elements that need to be investigated together to perform analysis and sensemaking tasks.

Astronomers are actually faced with both cases: on one hand, they have to deal with the increasingly large observation datasets produced by modern instruments; and on the other hand, they often have to look at multiple, heterogeneous-but-complementary data simultaneously, including images of the same source observed at different wavelengths or taken during different observing campaigns (time domain astronomy), research articles, statistics and charts such as, e.g., light curves. A recent study by Meade et al. [START_REF] Meade | Are tiled display walls needed for astronomy?[END_REF] hints at the usefulness of ultra-walls for astronomers. A precursor to our work, the WeSpace 26 project provided astronomers with a shared workspace to conduct co-located collaborative work sessions. It consisted of a multi-touch table and a projected wall display, enabling several users to share their image data. While the astronomers already found WeSpace valuable to support their workflow, its wall display actually had a fairly limited display capacity (3072 × 1536 pixels) compared to the wall displays considered here. We believe that thanks to their much higher pixel density, ultra-walls can provide even more support for various tasks in an astronomer's workflow, including sharing and analyzing data collaboratively during group meetings, or getting an overview of results produced by data mining algorithms. We also believe that ultra-walls can be valuable tools for other activities such as education and outreach, as further discussed in Section 5. The following section gives an overview of FITS-OW's main features, describing how users can interact with the data displayed on the wall.

FITS-OW OVERVIEW

FITS-OW, which stands for FITS images On Walls, is an application that enables visualizing FITS images, [START_REF] Wells | FITS -a Flexible Image Transport System[END_REF] and other data about the sources in those images, on ultra-walls driven by clusters of computers. Figure 2 shows the application running on the WILDER ultra-wall. It is based on earlier work by our group that explored the technical feasibility of displaying large FITS images by slicing and sharing them over the nodes of a computer cluster. [START_REF] Beaudouin-Lafon | Multisurface Interaction in the WILD Room[END_REF] At the time of writing, some features are still in the process of being transferred from this earlier project * to the FITS-OW codebase. †

Features

FITS-OW can load two main types of FITS images. First, users can display multiple regular 2D FITS images, either accessed through the network or stored locally. M31 (Figure 2-a) and the three images of the Eagle nebula (Figure 2-b) are instances of images of this type. Second, users can also load very large FITS images (up to several hundred thousand pixels in both width and height) that have been pre-processed by FITS-OW to generate the equivalent of mipmaps 28 (tiled pyramids) for FITS images. Further information about how these multi-scale FITS images are generated is given in Section 4. The background image in Figure 2-c is an instance of such a very large, pre-processed multi-scale FITS image.

Multi-scale FITS images often far exceed ultra-walls' display capacity and thus span the entire display surface, even at low magnification. Several such background images can be loaded simultaneously, but in most cases, these images will correspond to the same on-sky region with different filters or at different points in time ‡ . Such sets of multi-scale, spatially-aligned FITS images will be superimposed on different layers whose visibility and opacity can be interactively adjusted. Regular FITS images are displayed on top of these larger background images, and can be freely moved and resized. They do not necessarily have to correspond to a subregion of the background multi-scale images.

For each individual image, either regular or multi-scale, astronomers thus have access to the whole dynamic range and can look for sources at the very limits of the noise-dominated regime. They can select an appropriate color lookup table using the color map selector (Figures 2-f and3-a) and an appropriate scaling function (linear, square root, log, histogram equalization) as in FITS image viewers such as, e.g., SAOImage DS9. [START_REF] Mandel | New features of SAOImage DS9[END_REF] A histogram widget lets them specify the corresponding domain using a range slider. If those parameters are adjusted for a regular FITS image, the changes are applied dynamically to the image, as soon as a new color table or scaling function is selected. If they are adjusted for a multi-scale image, FITS-OW shows the resulting rendering in two preview windows only (Figure 3-b,c): one corresponds to a thumbnail view of the entire image, the other to a 1:1 close-up on a region of the image. The changes are applied to the entire image only after confirmation by the user. The purpose of this preview feature is to keep the interface responsive, as applying a new scaling function or color lookup table to a large multi-scale FITS image cannot be achieved at interactive frame rates.

Individual FITS images can be moved and resized freely. The multi-scale background FITS image, whose dimensions often exceed the wall's display capacity, can also be panned and zoomed at will. Users can also instantiate DragMags, [START_REF] Ware | The DragMag Image Magnifier[END_REF] that provide magnified, more detailed views on user-selected regions of the multi-scale FITS image. As illustrated in Figure 4, multiple dragmags can coexist, enabling users to view multiple sources featured in the background image (whose field of view will typically be relatively large) at the same time. Figure 5 ‡ For example, one image we refer to in Section 4 comes from the VISTA Variable in the Via Lactea survey 29 (VVV). The 6-year-long survey has collected 80 TB of reduced data, covering a 560 deg 2 field on the sky in 5 optical/near-infrared filters. The same fields were observed up to 100 times to follow luminosity changes of the stars within the survey area. illustrates a complementary scenario, where a source is rendered with different color lookup tables and scaling functions.

Images of the same region of the sky corresponding to observations at different wavelengths can be freely juxtaposed or superimposed. When superimposed, the opacity of those images can be adjusted, enabling users to perform manual composites. A typical scenario would be to take FITS images of the same source observed at different wavelengths, assign them different color lookup tables such as the red, green and blue color ramps, overlay them and adjust their opacity interactively. As mentioned earlier, this can also be achieved on the different layers of the multi-scale background FITS image.

Users can instantiate multiple cursors that can be moved seamlessly across the display surface. Those cursors can be controlled by touching the display, or remotely, using a handheld device such as a tablet (Figure 6-a). Each individual user can control one or more cursors. The displayed coordinates are obtained using Astropy's [START_REF] Robitaille | Astropy: A community Python package for astronomy[END_REF] WCS module to manage FITS image to World Coordinate System (WCS) transformations.

Astropy's WCS module is also used to convert user selections in an image into RA-Dec coordinates when querying astronomical catalogs. So far, we have only implemented support for SIMBAD queries. However, it should be relatively straightforward to add support for queries to other databases. As Figures 7 and8 illustrate, users can select a region in any of the FITS images displayed and perform queries, specifying filters on the object type and other parameters, and specifying which measurements they are interested in. In addition to superimposing sources from the result-set on top of the FITS image on which the query was made, FITS-OW also displays a sorted list of all items in that result-set. Selecting a source updates the info card, which gives detailed information about it, including the type of object, and measurements such as distance or radial velocity.

An additional feature related to SIMBAD queries we are currently working on is the possibility for users to display the research articles relevant to some particular source directly on the wall. The underlying UI toolkit already enables us to display PDF documents on ultra-walls by seamlessly sharing them dynamically with the cluster's nodes (Figure 2-e), and queries to SIMBAD can be configured to return URLs to research papers related to a given source. We are now working on an interface to enable users to select relevant papers to display, and 9. The contextual pie menu can be invoked on any image. After invoking the menu, the user selects the Scale item by moving leftward, which pops-up a submenu with four options corresponding to the available scaling functions. She then moves rightward, continuing the same gesture, selecting the log scaling function, which is immediately previewed. Releasing the fingers from the surface while inside this menu item confirms the selection of this new scaling function.

the underlying mechanism to automatically retrieve the PDF document from the URL provided by SIMBAD, which is usually an HTML page about the article in a digital library rather than a direct link to the PDF itself.

Interaction

While many of the features described above can be found in FITS image viewers such as DS9 [START_REF] Mandel | New features of SAOImage DS9[END_REF] or JSky, [START_REF] Brighton | JSky: Building a library of reusable Java components for Astronomy[END_REF] the latter have been developed to run on desktop workstations. All features are designed to be accessible through conventional input devices such as keyboard and mouse. But as discussed earlier in Section 2, these input devices are not adapted to interact with wall-sized displays. In order to promote physical navigation and to effectively enable collaborative work in front of ultra-walls, [START_REF] Beaudouin-Lafon | Multisurface Interaction in the WILD Room[END_REF] we started working on interaction techniques that better enable users to navigate in the data, adjust rendering parameters, and perform queries.

The two main types of input considered are handheld devices such as tablets, which enable users to interact with the display from any location, and direct touch on the wall (when available), which requires users to be within arm's reach of the display. Interaction is based on a combination of widgets such as pie menus, lists and buttons, multi-touch gestures, and text input.

As mentioned earlier, users can control multiple cursors remotely using tablets or smartphones. Those cursors can be used to obtain the sky coordinates of a particular image pixel, to make measurements, to bound queries to particular regions, and to adjust parameters such as the chosen color table. Multiple users can interact simultaneously, each controlling one or more cursors with their own handheld device. [START_REF] Chapuis | Smarties: An input system for wall display development[END_REF] They also have the possibility to take control of other users' cursors.

A contextual hierarchical pie menu [START_REF] Callahan | An empirical comparison of pie vs. linear menus[END_REF] provides access to some of the main features (Figure 9): initiating a query to SIMBAD, getting a global view of the workspace, changing an image's scaling function, its color lookup table. This menu can be invoked from the tablets, or by touching the wall with four or more fingers. Text input is only supported from handheld devices for now.

Multi-touch gestures complement the above interactions. The usual flick and pinch gestures used to pan and zoom can be performed on the wall. In addition, a set of gestures drawn from Olafsdottir & Appert's design space for multi-touch gestures for discrete and continuous control 36 are mapped to commands as described in Figure 10, enabling eyes-free interaction from a distance on the tablets.

IMPLEMENTATION

FITS-OW is implemented mainly in Java using zvtm-cluster, 37 a module of our zvtm toolkit [START_REF] Pietriga | A Toolkit for Addressing HCI Issues in Visual Language Environments[END_REF][START_REF] Pietriga | ZVTM -Zoomable Visual Transformation Machine[END_REF] (which was also used to develop several of ALMA's operations monitoring and control interactive visualizations [START_REF] Pietriga | Interaction Design Challenges and Solutions for ALMA Operations Monitoring and Control[END_REF] ). zvtm-cluster uses aspect-oriented programming (AspectJ 41 ) to enable zvtm-based applications to run almost unmodified on cluster-driven wall displays 10 as the one depicted in Figure 11, taking care of sharing the scene graph across all cluster nodes driving the display tiles. FITS-OW runs on the master computer (MSTR), with generic zvtm-cluster SlaveApp clients running on the display nodes. Those generic clients can render any zvtm graphics primitives, loading the corresponding external resources (such as FITS image tiles or PDF document pages) when necessary.

The configuration of zvtm-cluster is flexible, enabling applications such as FITS-OW to run on different hardware setups. Figure 11 The four directions are mapped to the scaling functions according to their position in the Scale pie menu (Figure 9).

1CP/Free/ 2CP/Free/ 3CP/Anchored/ 3CP/Free/ 3CP/Free/ External/Linear External/Linear External/Circular External/Circular External/Linear Linux Fedora) drives two display tiles, with two graphics cards per node (one for each tile). In this case, we typically choose to run two instances of SlaveApp, one per tile. A similar configuration is used on WILD (20,480 × 6,400 pixels), that depends on a cluster of 16 Mac Pro computers to drive 32 tiles 30" in diagonal.

The configuration used on WILDER (Figures 1 and2) is different: in this case, a cluster of 10 computers drives 15 × 5 = 75 square tiles (960 × 960 pixels each). Those tiles are grouped in 2 columns, 8 and 7 tiles in width respectively. Each cluster node (running Linux Ubuntu) drives one row of one of these columns, the Xorg server being configured so that the 7 or 8 tiles are seen as one single display workspace. The SlaveApp configuration that yields the best frame rates in this case consists of 2 instances per node, that each handle half of the node's logical display, whose dimensions are 7680 × 960 pixels.

As mentioned in Section 3.1, very large FITS images (up to 10 5 by 10 5 pixels) need to be pre-processed in order to generate multi-scale representations of them. These representations are basically mipmap [START_REF] Williams | Pyramidal parametrics[END_REF] versions of the original FITS images, consisting of square FITS tiles stored as individual files and arranged into a quadtree described in an XML hierarchy. The tiles themselves are still FITS files, and not pre-rendered JPEG or PNG images, as users need to be able to change the scaling function or color lookup table. The XML file describes the quadtree, and also holds the original image's min and max values. Indeed, this information is required to render the tiles across the cluster in a homogeneous manner, applying the scaling function to all tiles with the same parameters. Otherwise, the scaling function would be applied to individual tiles locally, based on each tile's own histogram, resulting in a very inhomogeneous rendering. This FITS pre-processor is implemented as a Python script that uses Astropy. [START_REF] Robitaille | Astropy: A community Python package for astronomy[END_REF] Computing mipmaps for these very large FITS images enables a more efficient sharing and rendering of the data on the cluster's individual nodes: each SlaveApp node has a synchronized copy of multi-scale images' quadtree, and only loads the tiles that intersect the current viewport (position and level of detail), thus limiting memory consumption and only rendering tiles on the nodes by which they are actually displayed. The largest FITS pyramid we tested so far is a 3-filter image that corresponds to a small part of the VVV survey. [START_REF] Minniti | VISTA Variables in the Via Lactea (VVV): The public ESO near-IR variability survey of the Milky Way[END_REF] The three images, 86,499 × 13,474 pixels each, are put on superimposed layers. We also tested JPEG pyramids up to 540,000 × 15,000 pixels (GLIMPSE360 [START_REF] Whitney | Glimpse360: Observing The Outback Of The Galaxy[END_REF] ), which gives an indication that FITS-OW can support multi-scale FITS images several tens of gigapixels in size, given that FITS tile rendering performance is similar to that of JPEG tiles thanks to the pre-rendering strategy described next.

The parsing and initial rendering of FITS images in Java is handled by JSky 33 v3.0. JSky FITSImage objects are then wrapped in a Glyph compatible with zvtm's graphical object model, [START_REF] Pietriga | A Toolkit for Addressing HCI Issues in Visual Language Environments[END_REF] which draws a pre-rendered AWT BufferedImage version of the FITS image based on the current scaling function and color table settings in order to speed-up rendering. A similar mechanism is implemented for wrapping PDF document pages, which are also pre-rendered using ICEpdf 42 at an appropriate level of detail and also stored in a BufferedImage. External resources such as FITS images and PDF documents can be served to the cluster's nodes in two ways. The first solution consists in copying the data on all cluster nodes and having SlaveApp instances access them locally. The second solution consists in serving them from a separate computer, which can be the one running FITS-OW (MSTR in Figure 11) or a dedicated server. The first solution works best for large but relatively static resources such as multi-scale FITS images (which typically amount to several gigabytes). On the ANDES and WILDER platforms, we store them on the individual nodes' solid state drive (SSD) to minimize access times. The second solution is less efficient but much more flexible, enabling the dynamic sharing of resources in a more transparent manner from the user's perspective. The latter is thus best suited to the sharing of FITS images and PDF documents retrieved dynamically from the Internet through, e.g., SIMBAD or Aladin during a work session. In this case, FITS-OW downloads these resources, caches them, and serves them to the cluster nodes through a dedicated HTTP server running locally.

As mentioned in Section 3.1, transformations between image pixels and sky coordinates (RA-Dec) are handled via Astropy. [START_REF] Robitaille | Astropy: A community Python package for astronomy[END_REF] We originally used JSky's utilities, but eventually realized that JSky may return incorrect results (this was confirmed by comparing the results returned by DS9, Astropy, and JSky). We then decided to use Astropy, which is already part of our framework. This adds some complexity, however, as Astropy depends on numpy, which only runs with cpython § . A first prototype used RabbitMQ to exchange messages between the Java and Python processes. We then decided to instantiate a cpython interpreter from the JVM running FITS-OW, which is made possible by JEP [START_REF] Johnson | JEP: Embed Python in Java[END_REF] (Java Embedded Python). Queries to astronomical databases use this coordinate conversion utility to generate the queries based on the user's selections in the FITS images on one hand, and to superimpose the sources from the result-set on the FITS images. For now, as discussed earlier, only SIMBAD queries are supported. The query strings/scripts are generated dynamically based on the parameters (selected filters and measurements) set by the user through the user interface (Figure 7), and sent as HTTP requests to the SIMBAD server.

Finally, regarding user input management, FITS-OW depends essentially on two libraries. Smarties 34 enthe rapid prototyping of interfaces for handheld devices such as tablets and smartphones, facilitating the development of collaborative interfaces for ultra-walls. We use the Java+Android implementation of Smarties, and couple it with the ILDAGestureManager (developed as an internal project for now). The latter allows to recognize gestures by analyzing input data streams coming from heterogeneous input devices (e.g., tactile screens, trackpads, motion tracking systems). The library enables the incremental recognition of elaborate gestures such as those defined in Olafsdottir & Appert's design space. [START_REF] Olafsdottir | Multi-touch gestures for discrete and continuous control[END_REF]

DISCUSSION AND FUTURE WORK

As present, FITS-OW is a proof of concept project and lacks many of the features that would make it a suitable research tool for astronomers. But rather than trying to replicate all the features present in mature software such as DS9, our goal is to investigate the potential benefits and limitations of ultra-walls using real-world usage scenarios, considering specific application areas such as, in this case, astronomy. FITS-OW is part of a larger research effort 5 to go beyond the design and evaluation of interaction and visualization techniques considered in isolation. We want to take a more holistic approach, investigating how these techniques can be integrated, from both the interaction design and software engineering perspectives, to create actual interactive systems that support subject-matter experts in their large-data analysis tasks.

Astronomy has been one of our privileged application areas since we started researching ultra-walls in 2009. Astronomy is a very good candidate for our purposes, because the data that astronomers have to handle as part of their workflow are both very large, varied, and distributed. Of course, astronomers are not simply looking at images to perform their analyses. Given how large datasets have become (e.g., all-sky surveys such as GLIMPSE, 2MASS, SDSS or wide-field surveys such as VVV and DES), they now rather rely on algorithms, including data mining techniques, to identify relevant sources and phenomenons. But although data analyses can be partially automatized using sophisticated software, data exploration and quality assessment still relies on the interactive visualization of images and other data on screen, a task that astronomers perform on their own computer. As observed by Meade et al.: "As valuable as automatic analyses of these enormous datasets are, astronomy still relies heavily on visual inspection". [START_REF] Meade | Are tiled display walls needed for astronomy?[END_REF] We believe that ultra-walls can be especially interesting for such tasks, especially those that involve astronomical surveys, as ultra-walls make it possible to visualize a significantly larger field of view at the same level of detail than a desktop workstation. For instance, the results of data mining algorithms can be visualized on the wall, which will make visual inspection (finding the most promising sources, removing false positives, etc.) more efficient, both because as discussed earlier in Section 2 it will promote physical navigation, [START_REF] Ball | Move to improve: Promoting physical navigation to increase user performance with large displays[END_REF] and because it will enable users to display multiple complementary data in situ, that will help make decisions.

The design choice that consists in treating very large FITS images as zoomable backgrounds was made with astronomical surveys in mind. Again, such images are not aimed at being displayed alone for free-form, purely visual, exploration. Still, displaying them as background images gives a frame of reference on top of which additional information can be displayed, much like a base map does in geographical information systems. [START_REF] Lobo | An evaluation of interactive map comparison techniques[END_REF] The information displayed on top of these background images can be anything relevant, including, as detailed earlier: portals showing a magnified view of a given region or source in the background image (Figure 4), other portals showing the same source observed at different points in time or at different wavelengths, using different scaling functions and color maps (Figure 5), pages of research articles referencing a given source featured in the images (Figure 2), data coming from catalogs such as measurements and statistics: distance, flux, charts showing, e.g., spectra, or light curves for variable stars. As catalog cross-matching algorithms 45 make it increasingly easy to find data about a given source, the capacity of ultra-walls to display many such data without having to constantly manipulate windows to reveal hidden pieces of information -as typically happens on desktop workstationsshould prove valuable to astronomers. A background image is not limited to a single layer at a time. As discussed in Section 3.1, multiple (spatiallyaligned) images can be superimposed and composited. Individual layers can correspond to different filters from the same survey, or they can come from different surveys. For instance, the VVV survey 29 consists of a multiwavelength and multi-epoch imaging of the sky, giving astronomers access to more than a hundred different images of the same object. One of its goals is to discover stars that change their intrinsic brightness within days. The large FITS image mentioned earlier in Section 4 comes from this survey and actually corresponds to only a fraction of it. That specific image, centered on Baade's Window -a region in the Milky Way bulge with reduced dust extinction -could be combined with, e.g., OGLE, [START_REF] Udalski | The optical gravitational lensing experiment. final reductions of the OGLE-III Data[END_REF] an optical survey also studying variable stars, that also covers this region of the sky.

Meade et al. further observe that "[some] phenomena are not even predicted, making automatic analysis meaningless. It is often a case of not knowing what you are looking for until you see it", [START_REF] Meade | Are tiled display walls needed for astronomy?[END_REF] which resonates with the famous description of Data Visualization as "a way to answer questions you didn't know you had". [START_REF] Plaisant | The challenge of information visualization evaluation[END_REF] Beyond such serendipitous forms of discovery, ultra-walls can also be useful for classification tasks, [START_REF] Liu | Effects of display size and navigation type on a classification task[END_REF] which can also be performed collaboratively. [START_REF] Jakobsen | Up close and personal: Collaborative work on a high-resolution multitouch wall display[END_REF][START_REF] Liu | Shared interaction on a wall-sized display in a data manipulation task[END_REF] Better support for collaborative work was actually identified by Wigdor et al. [START_REF] Wigdor | WeSpace: The Design Development and Deployment of a Walk-up and Share Multi-surface Visual Collaboration System[END_REF] as something that would be highly valuable in some of the typical workflows of astronomers, including data analysis as well as the preparation of observation proposals. Their WeSpace platform, designed to facilitate face-to-face collaborative research sessions (which were not normally part of the astronomers' workflows), was informed by interviews and observations of astronomer group meetings, who expressed a "desire for synchronous co-located collaboration at various phases of their process". [START_REF] Wigdor | WeSpace: The Design Development and Deployment of a Walk-up and Share Multi-surface Visual Collaboration System[END_REF] Ultra-wall rooms provide similar environments, with the added benefit that the wall display features a much higher pixel density, enabling astronomers to share and render even more data and documents. One important feature in that respect, which is not yet supported by FITS-OW but enabled by the underlying UI software toolkits, [START_REF] Pietriga | Rapid development of user interfaces on cluster-driven wall displays with jBricks[END_REF] will be to enable users to display instances of applications running on their laptop on the wall. Such type of sharing, which enables only a relatively crude form of visual comparison and correlation, is still useful because it does not require sharing the actual data with other users, but just a presentation of it. This will usually make the process easier from a technical perspective, and can also be convenient if the data is still under embargo.

Finally, another avenue for future work we plan to investigate is to enable astronomers to perform some simple forms of on-the-fly data analyses directly on the data loaded on the wall. For now, all data displayed on the wall come from local files or from queries to data catalogs, i.e., data that were processed offline, prior to the work session in the ultra-wall room. Many large-scale analyses will only be available this way, especially for surveys, as the computing power and time required to perform them is not compatible with the time frame of the typical collaborative work sessions we envision. However, our integration of Astropy in FITS-OW through JEP, which was originally implemented as a workaround for WCS transformations (Section 4), is not limited to the WCS module. It could for instance be used to perform some simple on-the-fly analyses for several sources of interest displayed on the wall, possibly taking advantage of the cluster's computing and memory power; knowing that this power is sometimes under-utilized on ultra-wall visualization clusters depending on the actual hardware configuration.

Ultra-walls can also be of great benefit for education and outreach activities in astronomy. Here again, the very high display capacity means that a lot of complementary information about sources of interest can be shown to the audience at the same time, instead of being presented sequentially, as in typical slideshows. As both computer and display technologies keep advancing, we can expect ultra-walls to leave the confines of HCI and visualization laboratories to appear in classrooms and public spaces in the maybe-not-so-distant future. In the longer term, ultra-walls could also appear in scientists' offices, replacing, e.g., a bulkhead. Coupled with their laptop or workstation, they would then enable an even tighter integration of this novel generation of interactive display surfaces in their research workflow.

Figure 2 .

 2 Figure 2. FITS-OW running on the WILDER platform, showing: multiple FITS images, (a) M31 on the left side, (b) three juxtaposed images that show observations of the Eagle nebula at different wavelengths, and (c) a much larger FITS image (86,499 × 13,474 pixels) used as a zoomable background over the entire wall; (d) the result-set of a SIMBAD query restricted to observations about galaxies; (e) basic measurements for galaxy M31; (e) a page of a research paper (PDF) discussing that particular galaxy; (f) the color map selector.

Figure 3 .

 3 Figure 3. Color lookup table selector (a) with previews of the new selection applied to a thumbnail of the entire multi-scale FITS image (b) and to a 1:1 close-up on a region of the image (c).

Figure 4 .

 4 Figure 4. Multiple DragMags magnifying different sources of interest from the same multi-scale FITS image on the ANDES ultra-wall (11,520 × 4,320 pixels) set up at Inria Chile in Santiago.

Figure 5 .

 5 Figure 5. The same source rendered on the ANDES wall display using different color lookup tables and scaling functions.

Figure 6 .

 6 Figure 6. (a) Remote control of cursors using a handheld device. (b) Distance measurements using multiple cursors.

Figure 7 .

 7 Figure 7. SIMBAD Query Interface.

Figure 8 .

 8 Figure 8. Results of a SIMBAD query superimposed on the corresponding FITS image, along with a sorted list of all items in the result-set. Selecting an element in this list updates the detailed info in the lower right window and highlights the source in the image. All windows can be freely repositioned on the wall.

Figure

  Figure9. The contextual pie menu can be invoked on any image. After invoking the menu, the user selects the Scale item by moving leftward, which pops-up a submenu with four options corresponding to the available scaling functions. She then moves rightward, continuing the same gesture, selecting the log scaling function, which is immediately previewed. Releasing the fingers from the surface while inside this menu item confirms the selection of this new scaling function.

Figure 10 .Figure 11 .

 1011 Figure 10. Mapping of tablet multi-touch gestures to commands in FITS-OW.

  illustrates one such configuration, where each cluster node (Dell workstations running # fingers Gesture Type[START_REF] Olafsdottir | Multi-touch gestures for discrete and continuous control[END_REF]

		Action
	1CP Free/External/Linear	Move selected cursor: absolute repositioning if initiated outside sym-
		bol, relative repositioning if initiated inside it
	2CP Free/External/Linear	Pan
	3CP Anchored/External/Circular Smoothly zoom in (CW)/out (CCW)
	3CP Free/External/Circular	Switch to the previous (CCW)/next (CW) color table
	3CP Dwell (500 milliseconds)	Switch to the corresponding scaling function ( * )
	+ Free/External/Linear
	5CP Dwell (1 second)	Zoom-out to a global view of the workspace.
	CP = Contact Point
	( * )

§ cpython is the default implementation of the Python programming language, implemented in C. Jython is an alternative implementation of Python in Java, that enables better integration of Python and Java code.

ACKNOWLEDGMENTS

This work has been partially supported by: • Chile's Ministry for the Economy, Development, and Tourism's Programa Iniciativa Científica Milenio through grant IC 120009, awarded to the Millennium Institute of Astrophysics (MAS); • The Chilean Economic Development Agency (CORFO) and Inria through their support of Inria Chile as an International Center of Excellence in Chile; • France's ANR EquipEx DIGISCOPE -Investissements d'avenir ANR-10-EQPX-26-01.


Azeddine Cheref
email: cheref@laas.fr

Alessandro Agnetis
email: agnetis@dii.unisi.it

Christian Artigues
email: artigues@laas.fr

Jean-Charles Billaut
email: jean-charles.billaut@univ-tours.fr

Fixed sequence integrated production and routing problems

Introduction

The problem considered in this paper can be described as follows. A set of n jobs must be produced by a single machine and then delivered by a single vehicle.. Each job j requires a certain processing time p j on a the machine, and must be delivered after its completion to the location j. We denote by t ij the transportation time from destination i to destination j. We use M to denote the depot (manufacturer) and we assume that transportation times are symmetric and satisfy the triangle inequality. The set of jobs delivered during a single round trip constitutes a batch. The vehicle has a capacity c which is the maximum number of jobs it can load and hence deliver in a round trip. These models for coordinating production and delivery schedules have been largely analyzed and reviewed by Chen [START_REF] Chen | Integrated production and outbound distribution scheduling: review and extensions[END_REF], who proposed a detailed classification scheme. In this paper the production sequence is fixed and the jobs must be delivered in the order in which they are released, hence a production sequence also specifies the sequence in which the customers have to be reached. Since the production sequence is given, the problem consists in determining a partition of all jobs into batches, and each batch will then be routed according to the manufacturing sequence.

The performance measures we consider in this paper (denoted f ) is the total delivery time, i.e., f = n j=1 D j , where D j is the time at which the job j = 1, . . . , n is delivered at its destination.

Li et al [START_REF] Li | Machine scheduling with deliveries to multiple customer locations[END_REF] proved the NP-hardness of the general problem in which the sequence is not fixed and has to be decided. In this paper we show that the problem (denoted P ) is already NP-hard when the sequence is fixed and we deal with the special case in which all distances t ij are identical. For this special case, we show that the problem can be efficiently solved. Finally, we briefly enounce additional results. Detailed proofs are given in [START_REF] Agnetis | Integrated production and routing problems: algorithms and complexity[END_REF].

For our purposes, we introduce the Even-Odd Partition (EOP) problem. Given a set of n pairs of positive integers (a

1 , b 1 ), (a 2 , b 2 ), . . . , (a n , b n ), in which, for each i, a i > b i . Letting K = n i=1 (a i + b i )
, is there a partition (S, S) of the index set I = {1, 2, . . . , n} such that i∈S a i + i∈ S b i = K/2 ? EOP is proved NP-hard in the ordinary sense by Garey et al [START_REF] Garey | One processor scheduling with earliess and tardiness penalities[END_REF]. We will actually use a slightly modified version of the problem which remains equivalent to (EOP).

Modified Even-Odd Partition (MEOP). A set of n pairs of positive integers

(a 1 , b 1 ), . . . , (a n , b n ) is given, in which, for each i, a i > b i . Letting Q = n i=1 (a i -b i ), is there a partition (S, S) of the index set I = {1, 2, . . . , n} such that i∈S (a i -b i ) = Q/2 ? Theorem 1. P ( n j=1 D j ) is NP-hard. Proof (Sketch).
Given an instance of MEOP, we build an instance of P as follows. There are 3n + 3 jobs. The processing times of the jobs are defined as follows:

p 1 = p 2 = p 3 = p 3n+2 = p 3n+3 = 0, p 3n+1 = 4x n + b n + Q/2
and for all i = 1, . . . , n -1, p 3i+1 = 4x i + b i -2, p 3i+2 = 1 and p 3i+3 = 1. Where the x i are defined as

x i = (3a i -2b i + 3(n -i)(a i -b i ))
/2 for all i = 1, . . . , n and x n+1 = 0.

In the following, we refer to the set of jobs (3(i -1) + 1, 3(i -1) + 2, 3i), i = 1, . . . , n, as the triple T i . For what concerns the travel times, we let:

t M,3(i-1)+1 = t 3(i-1)+1,M = t M,3(i-1)+2 = t 3(i-1)+2,M = t M,3i = t 3i,M = x i ∀i = 1, . . . , n t 3(i-1)+1,3(i-1)+2 = a i , t 3(i-1)+2,3i = b i , t 3i,3i+1 = x i + x i+1 ∀i = 1, . . . , n t M,3n+1 = t 3n+1,M = t M,3n+2 = t 3n+2,M = t M,3n+3 = t 3n+1,3n+2 = t 3n+2,3n+3 = 0
Finally, vehicle capacity is c = 2. The problem consists in determining whether a solution exists such that the total delivery time does not exceed f * .

f * = n i=1 (3C 3i + 7x i + b i ) + C 3n+1 + C 3n+2 + C 3n+3 -Q/2.
(1)

Then we use the following Lemma.

Lemma 2. If a schedule satisfying (1) exists, then there exists one satisfying the following property: for all i = 1, . . . , n + 1, jobs 3i and 3i + 1 are NOT in the same batch.

Proof (omitted) Since the schedule satisfies Lemma 2 and c = 2, for each triple T i , i = 1, . . . , n, there are exactly two batches, and only two possibilities (see Figure 1), namely:

• option A: the first batch is {3(i -1) + 1, 3(i -1) + 2} and the second is {3i}.

• option B : the first batch is {3(i -1) + 1} and the second is From there, one can see that a schedule of value f * exists if and only if EOP is a yes-instance.

{3(i -1) + 2, 3i} .... .... Ti+1 3n+1, 3n+2, 3n+3 T2 1 2 3 b 1 3(n-1)+1 3(n-1)+2 3n bn Q/2 3(i-1)+1 3(i-1)+2 3i b i .... .... Ti+1 3n+1, 3n+2, 3n+3 i j=1 aj -bj 3(i-1)+1 3(i-1)+2 3i a i 3(n-1)+1 3(n-1)+2 3n an Q/2 a1 -b1 T2 1 2 3 a 1

A special case: constant travel times

In this section we address the special case in which all travel times are identical. We start by analyzing some properties of the optimal solution then we give polynomial solution algorithm.

Clearly, every time the vehicle is back at the depot, it can either (i) restart immediately with a new batch consisting of jobs already completed, or (ii) it can wait for the completion of some jobs to be delivered. Suppose that a vehicle that departed at time t returns at the depot, and starts again at a certain time t . Let J be the set of jobs released between t and t (extremes included). The round trip starting at time t is called maximal if either (i) the batch contains c jobs of J , or (ii) it contains all jobs of J . The next proposition gives a key feature of the optimal solutions. Proposition 3. There exists an optimal solution in which all round trips are maximal.

Proof (omitted)

Theorem 4. Problem P with constant travel times can be solved in polynomial time.

Proof (Sketch). Following Li et al [START_REF] Li | Machine scheduling with deliveries to multiple customer locations[END_REF], we call N SS (non stop shipment) a sequence of consecutive round trips during which the vehicle is never waiting at the depot, followed by a waiting time. We denote by [i, j] an NSS starting at time C i (hence, i is the last job of the first round trip of the NSS) and ending before C j (when another NSS will start).

Starting from the proposition 3, we can construct in polynomial time our [i, j] for all i < j. Note that, the last trip of [i, j] contains a maximum number of jobs in order to finish the trip before C j . We denote by f (j) the following value of the optimal solution of the problem restricted to the first j jobs. Then, we have that recursive formula: f (j) = min ∀i<j {f (i) + F ij }, where F ij = j k=i+1 D k and [i, j] an NSS. We also prove that the preemptive case in NP-hard. We study also the case where there is a fixed number of different locations. We propose pseudo-polynomial algorithms for the general case.

Figure 1 :

 1 Figure 1: Round trips with options A or B.

Acknowledgment

This work was supported by the ANR project no. ANR-13-BS02-0006-01 named Athena.


Azeddine Cheref
email: azeddine.cheref@univ-tours.fr

Christian Artigues
email: artigues@laas.fr

Jean-Charles Billaut
email: jean-charles.billaut@univ-tours.fr

A new robust approach for a production scheduling and delivery routing problem

Keywords: Integrated production scheduling and delivery routing, robust optimization, recoverable robustness, groups of permutable jobs, tabu search

In a supply chain, two essential elements are the production and the distribution. These two problems used to be solved separately, but it it well known that considering them jointly may lead to a global optimization of the supply chain performances. In this paper, we incorporate the delivery plan of a vehicle into a single machine scheduling problem, representing a single manufacturing facility. We consider the simplified case where there is only one vehicle available to serve the customers, with infinite capacity. However, we assume that the data are known with uncertainty. The objective is to find a schedule and a delivery plan so that a robustness criteria is minimized, under a scenario-based uncertainty modeling. In contrast with standard robust optimization approaches for scheduling, we do not propose a single complete solution to the problem that has to be feasible w.r.t. all scenarios and minimizes a worstcase criterion over the scenarios. Instead, we adopt the recoverable robustness framework that considers first-stage decisions and second-stage recovery options. We propose for the first-stage a set of solutions by using the concept of groups of permutable jobs. At second stage, a greedy and online recovery algorithm exploits the revealed information about the jobs available to be scheduled or to be delivered at decision time. We propose two tabu search algorithms, one based on the standard robust optimization scheme and one based on the new approach. We compare the two robust heuristics on a set of randomly generated problem instances.

INTRODUCTION

More and more research is devoted to integrated production scheduling and delivery problems, see for instance [START_REF] Chen | Integrated production and outbound distribution scheduling: review and extensions[END_REF]. However most of the work that has been done in the literature consider deterministic problems. However such integrated problems are highly subject to uncertainty. Data such as processing times or release dates are generally not known with an absolute certitude, and for delivery problems, transportation durations are time dependent and subject to fluctuations due to traffic. In this context, it is generally important to find a solution to the problem which is not only a good solution, but also a solution that keeps a good performance in presence of uncertainties. Robustness considerations in scheduling are frequently treated in the literature and it is not possible to make here an exhaustive review of the state-of-theart. We voluntarily restrain our work to robust scheduling approaches where no probability distribution is available for the uncertain data. Instead there is a set of possible scenarios for the input data. In particular we consider, as a reference method, the standard robust discrete opti-mization framework proposed by [START_REF] Kouvelis | Robust discrete optimisation and its applications[END_REF] (see also [START_REF] Bertsimas | Robust discrete optimization and network flows[END_REF]) . In this framework, a complete sequence of jobs on each machine is the output of the robust scheduling method. The uncertainty is represented by a discrete finite set of scenarios, each scenario corresponding to a deterministic value of the problem parameters. On a particular scenario, the performance of the sequence is generally measured by computing the earliest start schedule compatible with the prescribed sequence and the realized scenario. Hence, the worst case performance of a sequence on the scenario set can be computed a priori. The well-known problem with this framework is its conservativeness, as mentioned by [START_REF] Bertsimas | The price of robustness[END_REF], due to the fact that issuing a complete sequence lacks flexibility w.r.t. the realized scenario.

As an alternative we propose to extend these framework by integrating two other families of approaches for scheduling under uncertainty. The first family originates from production scheduling, and was described in [START_REF] Billaut | A new method for workshop real time scheduling[END_REF]. It aims at computing on a machine a sequence of groups of permutable jobs such that any permutation of the job inside each group yields a feasible job sequence. The final sequence is then selected in real time so as to better react to disruptions. In [START_REF] Wu | A Graph-Theoretic Decomposition of the Job Shop Scheduling Problem to Achieve Scheduling Robustness[END_REF], the author have shown through simulation of various processing time disturbances that computing a static way such a solution structure and then allowing the remaining decisions to be taken dynamically via dispatching rules yields promising results in a job-shop scheduling setting. The second family of approaches is based on recoverable robustness. This framework, defined in [START_REF] Liebchen | The Concept of Recoverable Robustness, Linear Programming Recovery, and Railway Applications[END_REF], aims at reducing the conservativeness of robust optimization. Recoverable robustness considers two sets of decision variables and consists in providing an assignment for the first-stage decision variables and a family of recovery algorithms such that, for each scenario, there exists one recovery algorithm in the family that computes feasible assignments for the second-stage decision variables. Recently the recoverable robustness framework has been successfully applied as such to a scheduling problem by [START_REF] Caprara | Delay-Robust Event Scheduling[END_REF] and to vehicle routing problem by [START_REF] Agra | The robust vehicle routing problem with time windows[END_REF] that used the close concept of adjustable robustness as proposed by [START_REF] Ben-Tal | Adjustable robust solutions of uncertain linear programs Mathematical Programming[END_REF].

To integrate these frameworks, we consider as in the standard robust scheduling approach a set of finite scenarios such that a scenario corresponds to a realization of uncertain parameters. We address a single machine scheduling problem, representing a single manufacturing facility in which a set of jobs with release dates is scheduled non preemptively. In addition, there is one vehicle available to deliver the jobs, with infinite capacity which has to be routed according to a travel time matrix between job locations. The vehicle takes a batch of completed jobs and delivers the customers at each job location preferably before the job delivery time and then return to the manufacturing facility to select another batch and, so on until all jobs are delivered. All parameters are assumed to be uncertain, i.e. job release dates, job processing times, job delivery dates and vehicle travel times The objective is to find a schedule and a delivery plan so that the maximum lateness over all jobs and all scenarios is minimized.

In the recoverable robustness framework, we propose to compute, as the output of the first decision level, a set of groups of permutable jobs on the machine, a set of batches and, inside each batch, a set of permutable job deliveries. To evaluate this first set of decision, a second level decision set adapts to the realized scenario in a greedy and online manner : by scheduling the production jobs inside each production group according to the earliest release date criterion (jobs are schedule as soon as they are available) and by scheduling the delivery jobs according to the nearest neighbor criterion (the vehicle always delivers the nearest job). This limited recovery algorithms simulates a real time decision making process which implies that a decision has to de made as soon as a portion of the scenario is revealed in an online framework. This means that when a group of production or delivery jobs is completed, the recovery algorithms wait for the next job to be released in the new group in the realized scenario for the production part or for the nearest job to visit according to the realized scenario for the delivery part.

In this work we formally present the standard robust and the online recovery robust integrated scheduling an delivery routing problem in Section 2. Then we present two tabu search heuristic, one implementing the standard robust scheme and one implementing the online recovery scheme in Section 3. In Section 4 we provide computational experiments to validate our proposal. Concluding remarks and directions for further research are pointed out in Section 5.

STANDARD ROBUST VS ONLINE RECOVERABLE ROBUST APPROACH FOR INTEGRATED PRODUCTION SCHEDULING AND DELIVERY ROUTING

Following the classical three-field notation of scheduling problems ?, Z-L. Chen introduces in Chen (2010) a fivefield notation α|β|π|δ|γ for the problems where scheduling and routing are integrated at an operational level.

In this notation, α, β and γ correspond to the machine configuration, scheduling constraints and objective function respectively, similarly to the notation of scheduling problems. Fields π and δ specify the characteristics of the delivery process and the number of customers, respectively.

In the following, we focus on a single machine environment (α = 1), with routing as a delivery problem (i.e. solving a vehicle routing problem is part of the decision), we have a single vehicle with unbounded capacity and n customers to deliver

(π = {routing, V (1, ∞)}, δ = n).
We consider a set J = {J 1 , ..., J n } of n jobs to schedule on a single machine representing a manufacturing facility (plant). We denote by p j the processing time of J j , r j the release date of J j that corresponds to the earliest possible start time of production and d j the delivery due date of J j , i.e. the date at which J j is supposed to be delivered to the customer. It is assumed that the storage capacity of jobs after production is not limited. However, we also consider that the jobs have strict due dates, which is equivalent to say that the maximum lateness is bounded. Each job J j is associated with a customer location j, and 0 denotes the location of the plant. t 0,j denotes the travel time between the plant and the customer location j and t i,j is the travel time between customer locations i and j (∀i, j, 1 ≤ i, j, ≤ n). Variables C j denote the production completion time of J j at the plant (∀j, 1 ≤ j ≤ n). Variables D j denote the delivery completion time of J j , i.e. the date of delivery at the customer location (∀j, 1 ≤ j ≤ n).

Remember that jobs are delivered by a single vehicle that may return to the plant several times (as in the multitrip traveling salesman problem) with unlimited capacity (capacity K = n). Without uncertainty considerations, the problem can be denoted by 1|r j |routing, V (1, ∞)|n|γ where γ is the objective function. In the following, we have γ = L max where L max denotes the maximum delivery lateness, defined by L max = max 1≤j≤n (D j -d j ). The deterministic problem is to find a production sequence of jobs (equivalently, to determine jobs production completion times), a batching of jobs and a route for each batch, so that the objective function is minimized. A batch can only start when all the jobs of the batch are completed.

We assume that some uncertainty affect the input data and that our objective is to find a robust solution. To represent the uncertainty, we use a scenario based approach similar to the one introduced in [START_REF] Kouvelis | Robust discrete optimisation and its applications[END_REF]. The data associated to the jobs and to the transportation times between sites vary according to predefined scenarios. We denote by S the set of possible scenarios over the planning horizon. Now, we say that to each job J j is associated a release date r s j , a processing time p s j and a due date of delivery d s j for each scenario s ∈ S. The matrix of travel times T s = (t s i,j ) 0≤i≤n,0≤j≤n gives the transportation time from any site i to any site j under scenario s ∈ S.

We denote by C s j the production completion time of job J j under scenario s and by D s j its delivery completion time under scenario s. The maximum lateness of J j under scenario s is denoted by L s j = D s j -d s j . Under uncertainty, the maximum lateness is L max = max s∈S (max 1≤j≤n L s j ). We now formally define the concept of online recoverable robustness. Let us consider the standard robust optimization problem that we can denote as min

x∈∪ s∈S Xs max s∈S f (x, s)
where S is the set of scenarios, X s is the set of solutions feasible for scenario s and f (x, s) is some performance indicator of solution x ∈ X s under scenario s ∈ S. This definition implies that x has to be feasible for any scenario, inducing a high degree of conservativeness [START_REF] Bertsimas | The price of robustness[END_REF]. The concept of recoverable robustness [START_REF] Liebchen | The Concept of Recoverable Robustness, Linear Programming Recovery, and Railway Applications[END_REF] is often considered in the literature as a two-stage extension to the robust optimization problem, in which the (first-stage) solution x can be modified in a limited way in the second-stage to fit with the realized scenario to obtain solution y = A(x, s) where A(x, s) is the application-dependent recovery algorithm of solution x under scenario s. Hence the recoverable robust optimization problem becomes:

min x∈X max s∈S f (A(x, s), s)
where X is the scenario-independent set of first-stage solutions, which is a more realistic definition in practical applications. Note that x does not necessarily assign all the decision variables of the problem. As in stochastic programming x can be only the subset of first-stage decision variables, while algorithms A(x, s) returns a complete solution, giving values to both first-stage and second stage decision variables.

Instantiating this concept to our integrated scheduling and vehicle routing problem, we define two different robust problems by using different definition of the first stage decisions x and of the second-stage online algorithm A.

The standard robust integrating scheduling and delivery routing problems defines x = (π, B, Σ) where π a complete of the jobs on the machines B a batching (i.e. a partition of the jobs) and Σ = {σ(B), B ∈ B} the complete ojob delivery sequence inside each batch B ∈ B). Algorithm A(x, s) is a simple completion time adjustment algorithm that sets feasible completion times C s j and D s j according to the realized parameters following the prescribed production and delivery sequences.

We illustrate the robust optimization scheme with the following example.

The matrices (t 1 i,j ) and (t 2 i,j ) are the following:

s = 1 J 1 J 2 J 3 J 4 J 5 r 1 j 0 7 3 4 3 p 1 j 3 4 1 2 4 d 1 j 11 18 9 10 17 s = 2 J 1 J 2 J 3 J 4 J 5 r 2 j 3 3 0 1 7 p 2 j 2 5 1 3 3 d 2 j 9 17 10 11 18 (t 1 i,j ) =        0 2 4 3 2 2 2 0 3 2 1 3 3 3 0 3 2 2 3 2 3 0 1 3 2 1 2 1 0 2 2 3 2 3 2 0        (t 2 i,j ) =        0 3 2 2 2 3 3 0 4 1 2 3 2 4 0 4 3 2 2 1 4 0 2 4 2 2 3 2 0 4 3 3 2 4 4 0       
Let consider the schedule (3, 4, 1, 5, 2) for the production phase, two batches composed by jobs {1, 3, 4} and {2, 5} with routings (4, 3, 1) and (5, 2), respectively. The production and transportation phases are illustrated in Fig. 1 for each scenario. t r denotes the departure time of tour r.

- (4) - (0) r 1 (1) - - (0) (3) (2) r 2 -- (5) - (0)
0 2 4 6 8 10 12 14 16 18 20 22 24 26

L max = 3 - Scenario s = 2 J 1 J 2 J 3 J 4 J 5
0 2 4 6 8 10 12 14 16 18 20 22 24 26

- (3) (2) r 2 (4) - (0) r 1 (1) - - (0) (0) -- (5)
-

L max = 3 J 1 J 2 J 3 J 4 J 5 Scenario s = 1 - Fig. 1 (0)

. Robust solution for two scenarios

The online recovery robust integrating scheduling and delivery routing problems defines x = (π , B , Σ ) where π a sequence of groups of permutable jobs of the jobs on the machines B a batching (i.e. a partition of the jobs) and Σ = {σ(B), B ∈ B} the sequence of groups of permutable delivery jobs inside each batch B ∈ B). Algorithm A (x , s) is a greedy algorithms that takes the jobs inside a group of permutable production jobs in the earliest release date first order and the jobs inside a group of permutable delivery jobs in the nearest neighbor order, to simulate an online reaction to the revealed scenario.

We consider here the same example as for the standard robust approach.

Let consider the sequence of groups of permutable jobs

G J 1 = {J 1 , J 3 , J 4 } and G J 2 = {J 2 , J 5 }.
For the transportation phase, one batch contains jobs {J 1 , J 3 , J 4 } and the other batch contains the jobs {J 2 , J 5 } (notice that the composition of batches may be different than the definition of groups). In the first batch, there are two groups of permutable deliveries G D 1 = {3} and G D 2 = {1, 4} and in the second batch there is only one group G D 3 = {2, 5}. Transportation phases are illustrated in Fig. 2 for each scenario and global solution that includes the production and transportation phases are illustrated in Fig. 3 for each scenario. For this example the maximum lateness found is L max = 0. This illustrates the potential gain of applying the online recoverable robustness scenario. Indeed, the recoverable robustness framework has two merits: it is more realistic as such repairing elements are often present in practice and it is able in theory to reach better worst-case performances on the scenario set. However the question remains whether under limited CPU time, recoverable robustness algorithms are able to actually reach better performance than robust algorithms for reasonable problem instance size. To make such a comparison we propose two heuristics in Section 3 and we compare them under a fair setting in Section 4. batch 1 batch 2

- m 0 m 2 m 5 Y - y - m 0 m J 3 m J 1 m J 4 1 y - m 0 m 3 m 4 m 1 ? ? i m 0 m 2 m 5 1 Scenario 2 q 9 - m 0 m 3 m 4 m 1 6 6 q m 0 m 2 m 5 ) Scenario 1
- (3) r 1 (4) - (0) - - (1) (0) J 3 J 4 J 1 J 2 J 5 L max = 0 Scenario s = 2 (5) r 2 - (2) - (0) L max = 0 Scenario s = 1 J 1 J 3 J 4 J 5 J 2 0 2 4 6
8 10 12 14 16 18 20 22 24 26

- (0) r 1 - (4) (3) - (1) - - (0)
(2)

r 2 -- (5) - (0)

TABU SEARCH HEURISTICS

We propose in this section two tabu search algorithms for the integrated scheduling and vehicle routing problem. The first one is called SRT (Scheduling, Routing, Tabu search) to solve the standard robust approach and, a second called OSRT (Online, Scheduling, Routing, Tabu search) to solve the online recoverable robust approach. In order to compare them, the two algorithms proposed are relatively the same and can be presented together. Note that to simplify the procedure we assume that OSRT makes only one delivery group by delivery batch. This means that the routing of deliveries of each batch is decided by the greedy NN algorithm on the realized scenario. We consider that Algorithm 1. Tabu search algorithm

(1) Generate an initial solution (π, B, Σ) for SRT and (π , B , Σ ) for OSRT (2) Perform tabu search on the production scheduling part with fixed batches and routing with NN (on the reference scenario for SRT and on the realized scenario for OSRT).

(3) Perform tabu search on the batching part with fixed production schedule and routing with NN (on the reference scenario for SRT and on the realized scenario for OSRT). ( 4) (For SRT only) Post optimize the routing part with the 2-opt heuristic. (5) Return to step 1 until a stopping criterion (global time limit) is met among the scenario there is a particular scenario s 0 corresponding to the average case, or to a reference scenario (see Section 4). The routing part of SRT is also simplified in a first step by routing the jobs according to NN on the reference scenario s 0 , but the routing is post-optimized in a second phase. The Tabu algorithm implemented is quite standard. This metaheuristic starts from an initial solution and then searches the most appropriate solution among neighbors. A solution can be chosen if it is not tabu (i.e. not in the tabu list) or if it satisfies the aspiration criterion. The chosen neighbor immediately becomes tabu and the search process is then repeated using this chosen neighbor as a new basic solution. Algorithm 3 summarized the common process, with an additional step top postoptimize the routing for SRT. This post-optimization is used to compensate the fact that SRT is not allowed to react to the realized scenario.

We detail the different elements hereafter.

Inital solution

We first describe how we obtain a feasible solution for SRT.

We first find a sequence of production jobs by sorting the jobs according to their increasing due dates on scenario s 0 (Earliest Due Date rule). From this sequence, we build batch by selecting the best solution in terms of maximum lateness over all scenarios between creating a new batch with the current job or inserting this job into the previous batch. The routing is performed by the NN heuristic on the average scenario.

For OSRT, we also first sort the jobs according to their increasing due dates on scenario s 0 but we then build groups of permutable production jobs and batches simultaneously. More precisely, there are four options for each job, the best option is selected at each iteration:

• Insert the job in the current group; insert the job in the current batch For OSRT, each option is evaluated by ERD on the realized scenario for the groups of permutable production jobs and by NN on the realized scenario for the routing of each delivery batch.

Solution encoding

The problem of routing the jobs of a batch is systematically solved by the nearest neighbor algorithm (NN) in the STR algorithm as well as in the OSTR algorithm. In fact, assuming that there is only one group of permutable delivery for each batch fact that for each batch and each scenario the nearest neighbor search is used in OSTR algorithm. Therefore, it suffices to have the scheduling part and the composition of the batches to encode a solution. A solution of the problem S is encoded directly by a vector of size 2n. The first part of size n contains the scheduling part which is not the same for both algorithms STR and OSTR, whereas, the second part of size n contains the batching part witch is similar in both algorithms. For STR the first vector gives the position of each production job in the sequence while in OSR, the first vector gives the group number of the job.

Neighbor definition

Starting from this encoding, several neighborhoods have been defined. The neighborhood generation of the scheduling sub-problem are first presented for the two algorithms, then we present those of the batching sub-problem whose are similar in the two algorithms. We use a perimeter δ as the search area from a solution and fix it to 0.4. In each case, some details are given for more accuracy.

In the neighborhood generation of the scheduling subproblem of STR, two types of neighborhood generation methods are used: (1) Position sequence exchange: this neighborhood performs the exchange of two values v i and v j in the first part of the encoding. The perimeter δ is taken into consideration in the sense that |v i -v j | ≤ δn.

Then, the best move is retained over all possible changes.

(2) Insertion at position in the sequence: this neighborhood performs the exchange of the position v i of the job J i by inserting the job J i at a position v k and where, |v iv k | ≤ δn. The best move is retained and the jobs at positions v i+1 , .., v k are advanced in the sequence if

v i > v k or, jobs at positions v k , .., v i-1 are delayed if v i < v k .
In a neighborhood generation of the scheduling subproblem of the OSTR, neighbors concerns groups of permutable jobs. We denote by |G| the number of groups and have the following four types of neighborhood generation methods: (1) Groups exchange: this neighborhood performs the exchange of two values in the first part of the encoding, i.e., the swap of the assignment of two jobs to groups with taking into consideration the perimeter δ. Let v j the group of the job J i and v i the group of the job J j , the exchange is possible if and only if

v i = v j and v i -δ|G| ≤ v j ≤ v i + δ|G|.
(2) Insertion group at position: this neighborhood states that a job J i of the group v j can be inserted in an existing group

v k with v i - δ|G| ≤ v k ≤ v i + δ|G|.
(3) Group split: this neighborhood states that one group can be split into two groups, which delays the next groups. Since the sequence into the groups depends to the scenario, the number of possible cases to cut a group is exponential. Therefore, we order the jobs of a group according to the average due dates of the jobs over the scenarios, which reduces the number of cases for the sizes of the groups. (4) Group merging: this neighborhood performs the grouping of two consecutive groups into one, which advances the next groups.

Since the second sub-problem is the constitution of the batches, witch is the same for SRT and tOSRT, the neighborhood generation is then the same for the two problems. One can see also that the types of neighbors used for the scheduling part of OSTR algorithm where the goal was to form groups can be used to form batches. Then, we have the previously four types of neighborhood generation methods which are applied to the batches.

Other elements : tabu list, aspiration criterion and stop condition

To discourage the tabu algorithm to select the solutions that have been visited recently, a tabu list is integrated and it contains solutions can not be selected. These solutions are called Tabu. In our case, since the scheduling part and the routing part are treated separately, we opted for two tabu list with a fixed size equal to the number of jobs. The first list contains the part of the vector corresponding to the scheduling, so that the second list contains the part corresponding to the batches.The fact of having two tabu list can ensure that a solution is considered tabu only if her scheduling part or her routing part is tabu. For this reason, a solution is aspired if it is better than the best solution ever found. We have a fixed time limit, the general algorithm stops when the time limit is reached. However, for the scheduling sub-problem and the batching sub-problem, there is no time limit but the algorithm stops when the solution is not improved.

COMPUTATIONAL EXPERIMENTS

Data generation

To validate our models, we generate random data sets.

For the generation of scenarios, we call scenario s = 1 the "reference scenario" and we first generate this instance (r 1 j , p 1 j , d 1 j ). Given that r 1 j and d 1 j depend on the processing times p 1 j , we first generate the p 1 j for all the jobs. We take p 1 j ∈ [1, 50] and we denote by P = n j=1 p 1 j . The release dates r 1 j are generated in the interval [0, γP ] where γ is a given parameter. Dues dates d 1 j are generated in the interval [α -β 2 P, α + β 2 P ], where α and β are given parameters. For the generation of the matrix of distances between sites, we generate randomly the coordinates (X 1 j , Y 1 j ) of each site j in the interval [0, 50] and the distance between two sites is the classical euclidian distance. For modeling the uncertainty on the data, scenario s (s ≥ 2) is generated by a random perturbation of the reference scenario. For a given parameter v s j , we perturb the reference scenario [(1-

ω v )v 1 j , (1 + ω v v 1 j )]
where ω v represents the perturbation magnitude. To have a fair comparison, the same time limit was used for SRT and OSRT.

Computational Results

In this section, we compare the tabu search algorithms SRT and OSRT, as well as the initials solutions for both cases, on the instances presented above. The number of jobs varies from 10 to 100, the parameters α, β and γ are fixed to 1. There are different classes of problems depending on the number of jobs n and the coefficient of variation ω between the scenarios. An instance is then represented by his type called (n, ω) and for each n and ω ∈ {20, 40, 60}, 10 instances are generated with a number of scenarios |S| ∈ 2, 5, 10. The size of the tabu lists is equal to the number of jobs n and the perimeter δ to 0.4. For each (n, S, ω), the algorithms are tested for the ten instances and the displayed results represent an average for the ten instances. In Table 1, column ∆ (1) = (GSR -GOSR)/GSR is the average gap between the initial solutions GSR and GOSR and ∆ (2) = (SRT -OSRT )/SRT is the average gap between the solutions SRT and OSRT.

CONCLUSION

The results shown in table 1 show that the recoverable robustness framework based on groups of permutable jobs generally allows to find most robust solutions than the standard robust optimization framework, under the same experimental framework. This is especially true for an important number of jobs and the largest number of scenarios (up to near 18% improvement). The improvement is even drastic for the comparison of the simple greedy initial solution procedures. In the future we will compare the proposed methods on a larger scenario set. It would be indeed interesting to compare the robust framework and the recoverable robust framework on a test bed of scenarios that are generated according to the same random perturbation scheme on the reference scenario, but which was not used to compute the robust solutions. More precisely, we would have a training set of scenarios of reduced cardinality to compute the robust solution, while these solutions will be evaluated on a larger set, which would be closer to the real life context.

Fig

  Fig. 2. Tours {(3), (1, 4)} and {(2, 5)} under the two possible scenarios

Fig. 3 .

 3 Fig. 3. Production and transportation phases under the two possible scenarios

•

  Insert the job in the current group; create a new batch and insert the job • Create a new group insert the job; insert the job in the current batch • Create a new group and insert the job; create a new batch and insert the job

Table 1 .

 1 Comparison of methods on the instances with |S| ∈ {2, 5, 10}

				SRT Vs OSRT
	n	S	ω	∆ 1 (%)	∆ 2 (%)
	10	2	20	38.31	5.58
			40	36.08	1.95
			60	28.75	-4.32
		5	20	35.77	1.30
			40	29.67	0.54
			60	27.46	10.24
		10	20	38.26	2.40
			40	31.19	-0.24
			60	26.85	-0.83
	20	2	20	56.38	2.44
			40	47.68	1.32
			60	51.43	3.78
		5	20	54.60	3.92
			40	47.08	4.72
			60	47.79	2.00
		10	20	48.98	4.88
			40	45.03	5.62
			60	44.12	5.01
	30	2	20	55.82	0.88
			40	55.30	4.50
			60	57.20	1.44
		5	20	54.56	3.59
			40	52.19	4.22
			60	49.20	3.43
		10	20	54.50	5.01
			40	48.80	2.49
			60	45.54	5.46
	50	2	20	64.22	2.31
			40	59.17	3.10
			60	59.16	4.59
		5	20	59.45	4.27
			40	54.49	5.26
			60	51.00	5.95
		10	20	59.24	6.76
			40	49.48	3.19
			60	48.37	6.44
	100	2	20	62.70	1.30
			40	60.54	5.41
			60	61.09	4.90
		5	20	61.21	4.39
			40	52.08	3.13
			60	49.11	3.77
		10	20	58.37	16.15
			40	51.64	9.28
			60	47.90	17.82

ACKNOWLEDGEMENTS

This work was supported by the financial support of the ANR ATHENA project, grant ANR-13-BS02-0006 of the French Agence Nationale de la Recherche.

This work was supported by the financial support of the ANR ATHENA project, grant ANR-13-BS02-0006 of the French Agence Nationale de la Recherche.