Dr Johannes Karges

Dr Thirumal Yempala

Mickaël Tharaud

Dan Gibson
email: dang@ekmd.huji.ac.il

Dr Gilles Gasser
email: gilles.gasser@chimieparistech.psl.eu@www.gassergroup.com

Prof D Gibson

A Multiaction and Multitarget Ru(II)-Pt(IV) Conjugate Combining Cancer Activated Chemotherapy and Photodynamic Therapy to Overcome Drug Resistant Cancers

Cancer has emerged as one of the deadliest diseases worldwide. Pt(II) complexes are commonly used to treat this condition. To reduce their side effects and improve their pharmacological properties, Pt(IV) complexes are being developed as prodrug candidates that are activated by reduction in cancer cells. Concomitantly, photodynamic therapy (PDT) has received increasing attention over the last years. Among other compounds studied as photosensitizers (PSs), Ru(II) polypyridine complexes have gained much attention over the recent years due to their attractive characteristics. In this article, the first example of a novel Pt(IV)-Ru(II) conjugate, which combines cancer activated chemotherapy with photodynamic therapy, is presented. Upon entering the cancer cell, the Pt(IV) centre is reduced to Pt(II) and the axial ligands including the Ru(II) complex and phenylbutyrate are released. As each component has its individual target, the conjugate exerts a multitarget and multiaction effect with (photo-)cytotoxicity values upon irradiation up to clinically relevant 595 nm in the low nanomolar range in various (drug resistant) 2D monolayer cancer cells and 3D multicellular tumour spheroids.

Introduction

To date, the cure and treatment of cancer remains a major challenge in modern medicine. Cisplatin (cisdiamminedichloridoplatinum(II)) and its platinum based derivatives oxaliplatin and carboplatin are the clinically most frequently used anticancer drugs worldwide. Despite their enormous clinical success, the application of these compounds is limited due to severe side effects (e.g., kidney damage, nausea, vomiting and bone marrow suppression), low cancer cell selectivity and, more worryingly, an increasing number of platinum resistant tumours. [1] A known concept to improve the therapeutic outcome of chemotherapy is to administer cisplatin and its derivatives in combination with other biologically active molecules such as histone deacetylase inhibitors (Vorinostat) [2] , retinoid receptor regulators (vitamin D derivatives) [3] or microtubule network disturbing compounds (e.g., paclitaxel, vinblastine) [4] to generate a drug mixture with multiple targets and actions. Despite improvements in the treatment outcomes, the application of these drug mixtures is limited as the different components may reach the target after different circulation times leading to non-ideal drug doses. [5] To overcome these drawbacks, research efforts are devoted towards the development of novel Pt(IV) based compounds, which can be activated on site. As the compound enters the cell, the metal centre can be reduced from Pt(IV) to Pt(II), releasing its axial ligands. To improve the efficiency of the treatment, the axial ligands themselves can be biologically active molecules such as targeting moieties, carbon nanotubes, protein-binding moieties or receptor-binding peptides. [6] Recently, the use of phenylbutyrate (PhB), which is a histone deacetylase inhibitor that de-condenses chromatin to improve the interaction of cisplatin with DNA, was demonstrated. [7] As a consequence, Pt(IV) prodrugs with PhB as an axial ligand were shown to have a highly improved therapeutic profile, causing DNA damage and epigenetic effects. [8] Complementary to chemotherapy, cancer treatment using photodynamic therapy (PDT) has received increasing attention over the recent years. In PDT, a photosensitizer (PS) is either locally or systemically injected. Upon excitation at a specific wavelength, the PS acts as a photocatalyst to generate reactive oxygen species (ROS) and ultimately trigger cell death. [9] Worthy of note, Sadler and co-workers have recently reported the use of photoactivatable Pt(IV) complexes which can act by non-classical mechanism. [10] Among the investigated classes of compounds, metal complexes [11] and especially Ru(II) polypyridine complexes have received increasing attention due to their attractive chemical and photophysical properties including high water solubility, high ROS production, chemical stability and photostability. [12] Despite recent research efforts, the majority of studied Ru(II) polypyridine complexes are excited using blue or UV-A light. [13] Due to light scattering effects and absorption by biological compartments, the light penetration depth in the tissue is poor at these wavelengths, limiting the treatment of deep seated or large tumors. [14] To tackle this issue, some of us have recently reported the DFT-guided design of Ru(II) polypyridine complexes with a strong absorption red shift towards the biological spectral window (600-900 nm). In this pursuit, the complex [Ru(4,7-diphenyl-1,10-phenanthroline) 2 (4,4′-dimethyl-2,2′bipyridine)] 2+ , which showed a phototoxic effect from 480 nm up to clinically relevant 595 nm, was unveiled. [15] Herein, a novel Ru(II)-Pt(IV) conjugate (Ru-Pt), which combines cancer activated chemotherapy and PDT by acting on multiple targets inside the cancer cell with several biological mechanism to overcome drug resistances, is described. Upon entering cancer cells, the Pt(IV) metal centre is reduced to Pt(II) releasing PhB and the Ru(II) polypyridine complex. The Ru(II) complex was found to selectively accumulate in the Golgi apparatus where it is able to generate singlet oxygen ( 1 O 2 ) upon irradiation at various wavelength ranging from 480 up to clinically relevant 595 nm.

Results and Discussion

To ensure that the photophysical and biological properties of each component are not influenced, the moieties were conjugated with a long aliphatic linker. The linking unit between the Ru(II) and the Pt(IV) complex was synthesised starting from 4,4´-dimethyl-2,2´-bipyridine, which was asymmetrically oxidised to the carboxylic acid using SeO 2 . [16] The generated 4´-methyl-2,2´-bipyridinyl-4-carboxylic acid was treated with SOCl 2 to yield the acid chloride intermediate, which was used to attach tert-butyl-N-(6-aminohexyl)carbamate by peptide bond formation. [17] The generated 2,2´-bipyridine derivative was then coordinated to dichlorobis(4,7-diphenyl-1,10-phenanthroline)ruthenium(II), which was previously synthesised from 4,7-diphenyl-1,10-phenanthroline and RuCl 3 . Finally, the tert-butyloxycarbonyl protection group was removed using trifluoracetic acid to give the Ru(II) polypyridine complex Ru. To facilitate the conjugation of Ru to the Pt-PhB-OH the axial OH group of the Pt(IV) was activated to form Pt-PhB-MSC (MSC=monosuccinimidylcarbonate), whose synthesis was prepared in a similar manner to the structurally related derivatives. [8b] Both carboxylates (link to PhB) and carbamate (link to Ru) detach from the Pt IV following reduction. [8c] The Ru complex was then coupled with Pt-PhB-MSC in presence of diisopropylethylamine to yield the final Ru-Pt complex. The synthesis and characterisation of the complexes are described in the SI (Scheme S1, Figure S1-S8). The photophysical properties of Ru were compared with Ru-Pt to evaluate its potential as a PDT PS and to investigate whether the conjugation to the Pt(IV) complex influenced them (Table S1). As expected, the absorption and emission spectra of Ru and Ru-Pt (Figure S9-S10) showed no significant differences. Importantly, the spectra showed an absorption tail towards the biological spectral window potentially enabling the treatment of deep-seated and large tumours as the light penetration depth directly correlates with the used wavelengths. Worthy of note, previous studies have shown that PSs can show a phototoxic effect although their extinction coefficients are below 100 M -1 cm -1 . [12a,15] The conjugate was found to have a large Stokes shift, implying minimal interference between excitation and emission. The luminescence quantum yields of Ru an Ru-Pt were found to be similar (Φ em, Ru = 2.0 %, Φ em, Ru-Pt = 2.1 %), which is in the same range as other published Ru(II) polypyridine complexes. [18] Both compounds were found to have excited state lifetimes in the Table 1. IC50 values (μM) in the dark and upon irradiation at 480 (10 min, 3.1 J/cm 2 ) and 595 nm (120 min, 22.5 J/cm 2 ) for Pt-PhB-OAc, Ru and Ru-Pt in comparison to Protoporphyrin IX (PpIX) in cancerous in human cervical carcinoma (HeLa) and non-cancerous retinal pigment epithelium (RPE-1) cells. Average of three independent measurements. [19] The characteristic emission signal at 1270 nm in the infrared region using the direct method confirmed the production of 1 O 2 upon light exposure. Importantly, Ru and Ru-Pt are able to produce 1 O 2 upon irradiation at 450 and 540 nm. The measurements revealed high 1 O 2 quantum yields of 45-64 % in CH 3 CN and of 3-6 % in an aqueous solution (Table S2). Overall, these studies showed that Ru-Pt has attractive photophysical properties for applications as a PDT agent and the conjugation to the Pt(IV) centre did not significantly influence the photophysical properties of the Ru(II) complex.

The stability of the conjugate was then investigated as previous studies have shown that this could be problematic for metal complexes. [20] The conjugate was incubated in the dark in H 2 O at room temperature and physiological pH and the absorption recorded in time intervals up to 48 h (Figure S13). Additionally, the complex was analysed by HPLC after a 48 h incubation in H 2 O (Figure S14). As no significant differences were observed, the stability in an aqueous solution is indicated.

After a chemical and photophysical characterisation of the conjugate, its potential activity as a chemotherapeutic agent and PS for PDT was investigated. The cytotoxicity of cisplatin, Pt-PhB-OH, Pt-PhB-OAc, Ru and Ru-Pt towards noncancerous retinal pigment epithelium (RPE-1) and cancerous human cervical carcinoma (HeLa) cells in the dark as well as upon irradiation at various wavelengths between 480 -595 nm was determined (Table 1, Table S3). In both investigated cell lines the conjugate Ru-Pt (IC 50, HeLa = 3.26 ± 0.07 μM, IC 50, RPE-1 = 4.38 ± 0.57 μM) showed a higher cytotoxicity in the dark than Ru (IC 50, HeLa = 13.57 ± 0.28 μM, IC 50, RPE-1 = 4.22 ± 0.42 μM) or the platinum based complexes cisplatin, Pt-PhB-OH or Pt-PhB OAc (IC 50, HeLa = 5.9 -18.9 μM, IC 50, RPE-1 = 27.7 -81.4 μM) alone. All Pt(IV) prodrugs and importantly Ru-Pt showed a 1.3 -4.7 times higher cytotoxicity in the cancerous cell line in comparison to the noncancerous cell line. This effect was previously described in the literature and is caused by the preferential reduction of the Pt(IV) centre to Pt(II) in cancer cells, which are hypoxic. Strikingly, upon irradiation at various wavelengths between 480 -595 nm, the Ru-Pt conjugate is able to cause a phototoxic effect. Thanks to the additional PDT effect, the cytotoxicity improved drastically, with IC 50 values in the nanomolar range (IC 50, HeLa = 0.13 -0.39 μM, IC 50, RPE-1 = 0.45 -0.61 μM).

Overall, these results indicate the beneficial effect of combing the properties of the Pt(II) complex as a chemotherapeutic agent and the Ru(II) complex as a PS as a conjugate. Interestingly, Ru-Pt showed a slightly higher cytotoxicity in cancerous cells in comparison to non-cancerous cells.

Capitalising on these preliminary results, the potential of Ru-Pt to overcome resistance was further investigated in depth in human ovarian carcinoma (A2780) cells and its cisplatin resistant line (A2780 cis) as well as its doxorubicin resistant line (A2780 ADR). The uptake mechanism of the Ru-Pt conjugate was studied in A2780 cells by blocking different cellular pathways and determining the amount of Ru and Pt inside the cells (Figure 2a) by inductively coupled plasma mass spectrometry (ICP-MS). [21] The incubation with the cationic transporter inhibitor (tetraethylammonium chloride) did not significantly influence the cellular uptake, suggesting that this pathway does not majorly contribute to the internalisation of the compound. The incubation at lower temperature (4 °C) and with metabolic inhibitors (2-deoxy-D-glucose and oligomycin) decreased the uptake, indicating that the internalisation is energy dependent. In contrast, the incubation with endocytotic inhibitors (ammonium chloride or chloroquine) drastically decreased the uptake. Overall, these results indicate that Ru-Pt is internalised primarily by an energy dependent endocytotic mechanism. Worthy of note, the recently reported Ru-Pt metallacage was internalised by the same pathway. [START_REF] Zhou | Proc. Natl. Acad. Sci[END_REF] Importantly, similar levels of Ru and Pt were detected inside the cells indicating that the Ru-Pt is stable in the cell culture medium and was internalised as a single moiety. Following this, the localisation of the released Ru(II) complex and Pt(II) complex (cisplatin) were investigated by extraction of the major cellular organelles (nucleus, mitochondria, lysosomes, Golgi apparatus, endoplasmic reticulum) and determination of the amount of Ru and Pt inside of each by ICP-MS (Figure 2b). As expected, Ru and Pt were found in different cellular compartments, indicating that the Pt(IV) centre is reduced and the axial ligands (Ru and PhB) released and that the conjugate is acting simultaneously in different cellular organelles. As anticipated, the vast majority of Pt was found in the nucleus. This finding is in agreement with the in-depth investigation of the anticancer drug cisplatin, which is targeting the nucleus and binding to DNA. On the contrary, the majority of Ru was found in the Golgi apparatus. Additionally, small amounts of Pt and Ru were also found in the lysosomes and the endoplasmic reticulum. Taken together, these results indicate that Ru-Pt is internalised through an engulfing mechanism, in which the cell membrane draws the complex from the outside to the inside by formation of intracellular transport vesicles. These can then interact with the endomembrane system consisting of the lysosomes, Golgi apparatus and endoplasmic reticulum, where they release their content. As small amounts of Ru and Pt were both found in the lysosomes and the endoplasmic reticulum, it tends to demonstrate that the Ru-Pt conjugate is majorly released in these organelles. In a further process, the Pt(IV) centre is reduced and the Ru and Pt moieties localise in their respective target organelles (Golgi apparatus and nucleus).

To study the efficiency of the conjugate against drug resistant cell lines, the cytotoxicity of cisplatin, Pt-PhB-OH, Pt-PhB-OAc, Ru and Ru-Pt in the dark and upon irradiation at 480 (10 min, 3.1 J/cm 2 ) and 595 nm (60 min, 11.3 J/cm 2 ) against A2780, A2780 cis and A2780 ADR was investigated. All complexes (Table 2) were found to have lower cytotoxic values in all ovarian cancer cell lines in comparison to the HeLa and RPE-1 cell lines (Table 1). These results show the enormous potential of the conjugate Ru-Pt with cytotoxicity values in the dark in the low micromolar range (IC 50, dark = 0.98 -1.38 μM). It is worth mentioning that a PS should be ideally non-toxic in the dark. The irradiation drastically improved the cytotoxic profile (IC 50, 480nm = 0.08 -0.10 μM, IC 50, 595nm = 0.16 -0.19 μM) with PI values between 6.1 -13.8. While the cytotoxicity of cisplatin, Pt-PhB-OH and Pt-PhB-OAc is lower in the resistant cell lines A2780 cis and A2780 ADR, the Ru-Pt conjugate overcomes the resistance even in the dark and strikingly, by irradiation the cytotoxicity in all cell lines is further improved by a factor of ten. Comparing the cytotoxicity of Ru and Ru-Pt under irradiation we see that Ru-Pt is significantly more potent than Ru in all A2780 cell lines. Taken together, these results confirm the beneficial effect of the conjugation of the Ru(II) and Pt(IV) complex to the conjugate Ru-Pt.

Encouraged by these promising results, the mechanism of action of Ru-Pt was further examined by determining its cell death pathway. For this study, the cytotoxic profile of cisplatin, Ru in the dark, Ru-Pt in the dark and upon irradiation at 480 nm (10 min, 3.1 J/cm 2 ) at its IC 50 concentration in A2780 cells was measured in the presence of autophagy (3methyladenine), apoptosis (Z-VAD-FMK), paraptosis (cycloheximide) or necrosis (necrostatin-1) inhibitors (Figure 3a). The incubation with autophagy, paraptosis and necrosis inhibitors did not significantly influence the cell viability of cisplatin, Ru and Ru-Pt in the dark, indicating that these pathways are not primarily accessed for triggering of their cell death. In contrast, preincubation with an apoptosis inhibitor increased the survival of the cells treated with cisplatin, indicating that this mechanism is majorly contributing to the cell death. This result is in agreement with a study of cell death pathways of cisplatin in human ovarian carcinoma cells. [START_REF] Ormerod | [END_REF] In addition, elevated survival levels for the preincubation with an apoptotic inhibitor and treatment with Ru in the dark were detected, suggesting an apoptosis induced cell death mechanism. Therefore, it was not surprising to observe an apoptotic cell death pathway for Ru-Pt in the dark. Importantly, upon irradiation of Ru-Pt, the cell survival is highly increased with apoptotic and paraptotic inhibitors, indicating a cell death mechanism for the phototoxic effect involving both pathways. Worthy of note, several PDT agents were previously reported in the literature to cause cell death by apoptosis and paraptosis pathways. [24] For further examination of the cell death mechanism, its dependency on caspase 3/7 pathways was investigated using a Caspase-Glo 3/7 assay. These caspases are known executers of extrinsic and intrinsic apoptosis mechanism. [25] The caspase 3/7 activity was measured in A2780 cells for cisplatin, Ru and Ru-Pt in the dark and upon irradiation at 480 nm (10 min, 3.1 J/cm 2 ) using half of their IC 50 concentration (Figure 3b). As a positive control for caspase activity, the kinase inhibitor staurosporin was used. Cisplatin, Ru and Ru-Pt in the dark did not significantly increase caspase 3/7 activity, ruling out this pathway for cell death. The result for cisplatin is in agreement with a recent study in the same cell line that showed that apoptosis is induced by ERK/p53/PUMA activation. [26] In contrast to this, upon light irradiation and treatment with Ru and Ru-Pt, the caspase 3/7 activity is highly elevated, suggesting that the phototoxic effect is caused by these pathways. Worthy of note, other Ru(II) polypyridine complexes were also found to generate a cytotoxic effect by this pathway. [27] After evaluation of the (photo-)cytotoxicity on 2D monolayer cells, the effect of Ru-Pt on 3D multicellular tumour spheroids (MCTS) was investigated. This is of special interest as many anticancer agents have failed translation from 2D monolayer cells to in vivo models due to compromised drug delivery. MCTS can mimic the pathological conditions found in solid tumors such as proliferation gradients and its hypoxic centre and can therefore be used for the assessment of drug delivery. [28] A2780 MCTS with a diameter of about 400 μM were incubated with increasing concentrations of Ru-Pt for 12 h. After this time, the cytotoxic effect in the dark and upon irradiation at 480 nm (10 min, 3.1 J/cm 2 ) and 595 nm (120 min, 22.5 J/cm 2 ) was determined upon measurement of the ATP concentration of living cells inside the MCTS. Ru-Pt showed a strong cytotoxic effect (IC 50, dark = 7.32 ± 0.41 μM), indicating that the conjugate is able to completely penetrate the MCTS. Importantly, as the cytotoxicity is shifting to the nanomolar range (IC 50, 480nm = 0.49 ± 0.16 μM, IC 50, 595nm = 0.68 ± 0.27 μM) upon irradiation at 480 nm and 595 nm, the synergetic effect of the Pt(IV) complex as a chemotherapeutic agent and the Ru(II) complex as a PDT agent is demonstrated.

Conclusion

In summary, in this article, a novel conjugate (Ru-Pt) combining a Pt(IV) complex as a chemotherapy drug and a Ru(II) complex as a photosensitiser for long wavelength photodynamic therapy, which targets several cellular organelles and acts by various cytotoxic mechanism, is reported. The conjugate is able to enter cancerous cells by an energy depended endocytosis mechanism. While being stable in an aqueous solution, the Pt(IV) centre is reduced inside cancerous cells to Pt(II) and releases its axial ligands. Phenylbutyrate acts as a histone deacetylase inhibitor and can de-condense the chromatin, improving the binding of the Pt(II) complex (cisplatin) to DNA. The Ru(II) complex was found to accumulate at the Golgi apparatus. Next to the cytotoxic effect exerted by the anticancer properties of cisplatin, the Ru(II) complex can act as a photosensitizer to catalytically generate singlet oxygen from 480 up to clinically relevant 595 nm. Thanks to the multitargeting and multiaction mechanisms based on a combination of apoptosis and paraptosis pathways, (photo-)cytotoxicity values in the low nanomolar range were observed in various 2D monolayer cells as well as 3D MCTS. Importantly, the conjugate showed a synergetic effect and was able to overcome drug resistance, which is found in many clinically treated tumors. We strongly believe that this novel conjugate and the approach of combing a Pt(IV) complex as a chemotherapeutic drug and a Ru(II) complex as a photosensitizer for photodynamic therapy has a great potential for further development.

Figure 1 .

 1 Figure 1. Chemical structures investigated in this work. Ru was isolated as a PF6 -salt whereas Ru-Pt was isolated as a TFA -salt.nanosecond range (FigureS11-S12) in an air saturated environment (147-183 ns) and in degassed environment (1003-1054 ns). Importantly, the presence of air had a drastic effect on the lifetime, indicating that the excited state of the complex can interact with molecular oxygen ( 3 O 2 ) to generate singlet oxygen ( 1 O 2 ). For quantification of the amount of 1 O 2 produced upon light irradiation, the singlet oxygen quantum yields were determined by two methods: 1) directly by measuring the phosphorescence signal of 1 O 2 , 2) indirectly by capturing 1 O 2 with a reporter molecule and monitoring its change by absorption spectroscopy.[19] The characteristic emission signal at 1270 nm in the infrared region using the direct method confirmed the production of 1 O 2 upon light exposure. Importantly, Ru and Ru-Pt are able to produce 1 O 2 upon irradiation at 450 and 540 nm. The measurements revealed high 1 O 2 quantum yields of 45-64 % in CH 3 CN and of 3-6 % in an aqueous solution (TableS2). Overall, these studies showed that Ru-Pt has attractive photophysical properties for applications as a PDT agent and the conjugation to the Pt(IV) centre did not significantly influence the photophysical properties of the Ru(II) complex. The stability of the conjugate was then investigated as previous studies have shown that this could be problematic for metal complexes.[20] The conjugate was incubated in the dark in H 2 O at room temperature and physiological pH and the absorption recorded in time intervals up to 48 h (FigureS13). Additionally, the complex was analysed by HPLC after a 48 h incubation in H 2 O (FigureS14). As no significant differences were observed, the stability in an aqueous solution is indicated. After a chemical and photophysical characterisation of the conjugate, its potential activity as a chemotherapeutic agent and PS for PDT was investigated. The cytotoxicity of cisplatin, Pt-PhB-OH, Pt-PhB-OAc, Ru and Ru-Pt towards noncancerous retinal pigment epithelium (RPE-1) and cancerous human cervical carcinoma (HeLa) cells in the dark as well as upon irradiation at various wavelengths between 480 -595 nm was determined (Table1, TableS3). In both investigated cell lines the conjugate Ru-Pt (IC 50, HeLa = 3.26 ± 0.07 μM, IC 50, RPE-1 = 4.38 ± 0.57 μM) showed a higher cytotoxicity in the dark than Ru (IC 50, HeLa = 13.57 ± 0.28 μM, IC 50, RPE-1 = 4.22 ± 0.42 μM) or the platinum based complexes cisplatin, Pt-PhB-OH or Pt-PhB OAc (IC 50, HeLa = 5.9 -18.9 μM, IC 50, RPE-1 = 27.7 -81.4 μM) alone. All Pt(IV) prodrugs and importantly Ru-Pt showed a 1.3 -4.7 times higher cytotoxicity in the cancerous cell line in comparison to the noncancerous cell line. This effect was previously described in the literature and is caused by the preferential reduction of the Pt(IV) centre to Pt(II) in cancer cells, which are hypoxic. Strikingly, upon irradiation at various wavelengths between 480 -595 nm, the Ru-Pt conjugate is able to cause a phototoxic effect. Thanks to the additional PDT effect, the cytotoxicity improved drastically, with IC 50 values in the nanomolar range (IC 50, HeLa = 0.13 -0.39 μM, IC 50, RPE-1 = 0.45 -0.61 μM). Overall, these results indicate the beneficial effect of combing the properties of the Pt(II) complex as a chemotherapeutic agent and the Ru(II) complex as a PS as a conjugate. Interestingly, Ru-Pt showed a slightly higher cytotoxicity in cancerous cells in comparison to non-cancerous cells. Capitalising on these preliminary results, the potential of Ru-Pt to overcome resistance was further investigated in depth in human ovarian carcinoma (A2780) cells and its cisplatin resistant line (A2780 cis) as well as its doxorubicin resistant line (A2780 ADR). The uptake mechanism of the Ru-Pt conjugate was studied in A2780 cells by blocking different cellular pathways and determining the amount of Ru and Pt inside the cells (Figure2a) by inductively coupled plasma mass spectrometry (ICP-MS).[21] The incubation with the cationic transporter inhibitor (tetraethylammonium chloride) did not significantly influence the cellular uptake, suggesting that this pathway does not majorly contribute to the internalisation of the compound. The incubation at lower temperature (4 °C) and with metabolic inhibitors (2-deoxy-D-glucose and oligomycin) decreased the uptake, indicating that the internalisation is energy dependent. In contrast, the incubation with endocytotic inhibitors (ammonium chloride or chloroquine) drastically decreased the uptake. Overall, these results indicate

Figure 2 .

 2 Figure 2. a) Cellular uptake mechanism study of Ru-Pt (10 μM) in A2780 cells in the presence of different inhibitors/conditions by determination of the amount of Ru and Pt inside the cells by ICP-MS. Endocytic inhibition: NH4Cl (50 mM) or chloroquine (100 μM), metabolic inhibition: 2-Deoxy-D-glucose (50 mM) and oligomycin (5 μM), cation transporter inhibition: Et4NCl, low temperature: incubation at 4°C, control: incubation at 37°C. b) Cellular localisation of Ru and Pt of the conjugate Ru-Pt (10 μM) in A2780 cells inside the major cellular organelles (endoplasmic reticulum, golgi apparatus, lysosomes, mitochondria, nucleus, whole cell) after 4 h incubation in the dark, extraction of their cellular organelles and determination of the amount of Ru and Pt inside each organelle by ICP-MS.

Figure 3 .

 3 Figure 3. Cell death mechanism study in A2780 cells in the dark or upon irradiation at 480 nm (10 min, 3.1 J/cm 2 ) a) upon treatment with the IC50 value of the corresponding compound and in the presence of different inhibitors by determination of the cell viability. Autophagy inhibitor: 3-Methyladenine (100 μM), apoptosis inhibitor: Z-VAD-FMK (20 μM), paraptosis inihibitor: Cycloheximide (0.1 μM), necrosis inihibitor: Necrostatin-1 (60 μM). b) upon treatment with a concentration of half of the IC50 value of the corresponding compound by determination of the caspase 3/7 activity.

Table 2 .

 2 IC50 values (μM) in the dark and upon irradiation at 480 (10 min, 3.1 J/cm 2 ) and 595 nm (60 min, 11.3 J/cm 2 ) for cisplatin,

Pt-PhB-OH, Pt-PhB OAc, Ru and Ru-Pt in

  comparison to Protoporphyrin IX (PpIX) in the human ovarian carcinoma (A2780), human cisplatin resistant ovarian carcinoma (A2780 cis) and human doxorubicin resistant ovarian carcinoma (A2780 ADR) cell lines. Average of three independent measurements.

				A2780					A2780 cis					A2780 ADR
		dark	480 nm	PI	595 nm	PI	dark	480 nm	PI	595 nm	PI	dark	480 nm	PI	595 nm	PI
	cisplatin	4.54 ±	-	-	-	-	19.53±	-	-	-	-	8.97 ±	-	-	-	-
		0.65					1.11					0.71
	Pt-PhB	1.39 ±	-	-	-	-	6.17 ±	-	-	-	-	3.69 ±	-	-	-	-
	OAc	0.22					0.36					0.32
	Pt-PhB-	1.07 ±	-	-	-	-	5.95 ±	-	-	-	-	2.56 ±	-	-	-	-
	OH	0.24					0.37					0.34
	Ru	15.35	2.86 ±	5.4	4.30 ±	3.6	17.26 ±	4.67 ±	3.7	7.10 ±	2.4	18.60 ±	3.41 ±	5.5	4.47 ±	4.2
		± 0.89	0.52		0.66		0.90	0.43		0.92		0.71	0.42		0.35
	Ru-Pt	0.98 ±	0.08 ±	12.3	0.16 ±	6.1	1.38 ±	0.10 ±	13.8	0.19 ±	7.1	1.13 ±	0.09 ±	12.6	0.17 ±	6.6
		0.22	0.04		0.05		0.29	0.03		0.06		0.28	0.03		0.04
	PpIX	>100	4.53 ±	22.1	8.19 ±	12.2	>100	6.95 ±	14.4	8.94 ±	11.2	>100	7.12 ±	14.0	10.38 ±	9.6
			0.61		1.13			1.46		1.53			1.96		1.62

Acknowledgements

We thank Dr. Philippe Goldner for access to state-of-the-art laser apparatus. This work was financially supported by an ERC Consolidator Grant PhotoMedMet to G.G. (GA 681679), has received support under the program "Investissements d' Avenir" launched by the French Government and implemented by the ANR with the reference ANR-10-IDEX-0001-02 PSL (G.G.). DG acknowledges the support of the Israel Science Foundation (grant 1002/18) and the Alex Grass Center for Drug Design and synthesis.


Jonathan Sarton
email: jonathan.sarton@univ-reims.fr.

Nicolas Courilleau
email: nicolas.courilleau@univ-reims.fr.

Yannick Remion
email: yannick.remion@univ-reims.fr.

Laurent Lucas
email: laurent.lucas@univ-reims.fr.

To appear in IEEE Transactions on Visualization and Computer Graphics Interactive Visualization and On-Demand Processing of Large Volume Data: A Fully GPU-Based Out-Of-Core Approach

Keywords: GPU, Caching system, Out-of-core data management, Large data, Interactive visualization, On-demand processing

In a wide range of scientific fields, 3D datasets production capabilities have widely evolved in recent years, especially with the rapid increase in their sizes. As a result, many large-scale applications, including visualization or processing, have become challenging to address. A solution to this issue lies in providing out-of-core algorithms specifically designed to handle datasets significantly larger than memory. In this article, we present a new approach that extends the broad interactive addressing principles already established in the field of out-of-core volume rendering on GPUs to allow on-demand processing during the visualization stage. We propose a pipeline designed to manage data as regular 3D grids regardless of the underlying application. It relies on a caching approach with a virtual memory addressing system coupled to an efficient parallel management on GPU to provide efficient access to data in interactive time. It allows any visualization or processing application to leverage the flexibility of its structure by managing multi-modality datasets. Furthermore, we show that our system delivers good performance on a single standard PC with low memory budget on the GPU.

INTRODUCTION

The needs for visualizing and/or processing large volume data are common today in different scientific fields and entertainment. Navigation inside such high-resolution volumes in real-time involves designing efficient out-of-core data management algorithms which address entire massive datasets from high-performance computing devices such as current GPUs.

The interactive addressing problem of any part of a data volume that exceeds the amount of GPU memory has been addressed by several methods. Nonetheless, these approaches have focused on visualization tools and more particularly volume rendering. Recent out-of-core proposed methods integrated into a fully visualization-driven pipeline are not adapted to manage efficiently on-demand processing of non-visible data during the visualization stage. In many contexts, however, it can be crucial to visualize and interact with data during visualization through processing. This involves providing out-of-core methods adapted to on-the-fly data modification by image processing algorithms applied during interactive visualization.

The proposed solution is a complete out-of-core pipeline from disk to GPU designed to access very large volumes exceeding GPU or CPU memory in interactive time. We base our work on modern methods already known in this field, such as the design of an output-sensitive algorithm with on-demand paging and data streaming, bricking and multi-resolution representation, the use of a brick pool as a cache on GPU texture memory and a virtual address translation mechanism.

This work aims to propose a solution which takes advantage of the many-core environment of the GPUs. This environment allows to carry out as many operations as possible in parallel while preserving enough GPU resources for an end-user application that could require significant computing resources (such as volume ray casting or convolution processing on the volume for instance).

In the following, a state of the art of recent advances in external memory management is given in section 2. Section 3 gives an overview of the proposed pipeline. Our solution is presented and discussed regarding the data representation in section 4 and out-of-core management in section 5. Section 6 describes the impacts of on-demand processing during interactive visualization. In section 7, we introduce a compari-son of our approach with the closest previous ones. Our method will be evaluated in section 8 before discussing some specific points in section 9. Our perspective works are exposed in section 10.

Contributions. In this paper, we introduce an efficient out-of-core data management solution suitable for both interactive visualization and on-demand processing on visible or non-visible data of large 3D regular grids from GPU. Our method is based on the principle of virtual memory to address voxels through a GPU page table (PT) hierarchy introduced in [START_REF] Hadwiger | Interactive volume exploration of petascale microscopy data streams using a visualizationdriven virtual memory approach[END_REF]. In previous approaches, the cache miss management methods for accessing any part of a large volume were organized in the screen space and therefore directly bounded to the displayed data. This does not allow to handle efficiently data requests for non-visible or not currently visualized data. Moreover, these cache miss reports imply large memory transfers between the device and the host. They can become a very strong constraint in a context of visualization on very high-definition displays and even worse in a multi-view environment. In this article, we propose an approach that addresses these issues while providing an efficient parallel GPU implementation of the outof-core addressing data structure. Our method is designed to limit the communications between the CPU and the GPU to their strict minimum while ensuring not to overload the GPU occupancy. In addition, our system is created for a general-purpose context and provides an API available from the GPU for any applications such as standalone or combined visualization and processing of large volume data.

RELATED WORKS

External memory data management [START_REF] Vitter | Algorithms and data structures for external memory[END_REF] also called out-of-core data management, defines the class of techniques used to handle data that are too large to fit entirely into the main memory of the unit in charge of their processing. More specifically, this kind of approach can be used in interactive visualization of massive models [START_REF] Chiang | Outof-core algorithms for scientific visualization and computer graphics[END_REF][START_REF] Dietrich | Interactive massive model rendering[END_REF]. There is a vast body of literature for large datasets visualization on different data types, other than volume data [START_REF] Budge | Out-of-core Data Management for Path Tracing on Hybrid Resources[END_REF][START_REF] Goswami | An efficient multi-resolution framework for high quality interactive rendering of massive point clouds using multi-way kd-trees[END_REF][START_REF] Goswami | Scalable Parallel Out-of-core Terrain Rendering[END_REF][START_REF] Stone | Immersive Out-of-Core Visualization of Large-Size and Long-Timescale Molecular Dynamics Trajectories[END_REF]. The development of out-of-core methods more specifically for real-time visualization of regular voxel grids has been motivated by volume ray casting [START_REF] Kruger | Acceleration Techniques for GPU-based Volume Rendering[END_REF] of large datasets on GPU and has been widely used during the last decade.

In 2008, Gobbetti et al. [START_REF] Gobbetti | A single-pass GPU ray casting framework for interactive out-of-core rendering of massive volumetric datasets[END_REF] were the first to offer a complete, outof-core, multi-resolution volume renderer on GPU. Crassin et al. [START_REF] Crassin | Gigavoxels: Rayguided streaming for efficient and detailed voxel rendering[END_REF] proposed the year after, a more efficient system with Gigavoxel, a ray-guided streaming of opaque voxelized surfaces for entertainment purposes. In [START_REF] Engel | CERA-TVR: A framework for interactive high-quality teravoxel volume visualization on standard PCs[END_REF], Engel presented a framework for scientific visualization of tera-voxels, improving previous works by optimizing the GPU to CPU communications.

In all previously mentioned works, a tree structure is used to address out-of-core bricked data (an octree or a generalized N 3 -tree in Gigavoxel) with a kd-restart algorithm to go through the tree on the GPU. The basic principle is the use of a brick pool in GPU texture memory as a cache to store small bricks of voxels with, in the case of [START_REF] Crassin | Gigavoxels: Rayguided streaming for efficient and detailed voxel rendering[END_REF][START_REF] Engel | CERA-TVR: A framework for interactive high-quality teravoxel volume visualization on standard PCs[END_REF], a node pool to store the tree nodes. Cache misses are reported for the bricks that are not present in the video memory. These pools are updated at each frame to insert requested data by replacing unused ones if needed; all managed with a simple Least Recently Used (LRU) mechanism. While Gobbetti et al. used visibility information for culling, the others introduce a full ray-driven streaming that only loads visible data. Brix et al. [START_REF] Brix | Visualization of large volumetric multi-channel microscopy data streams on standard pcs[END_REF] approach relies on concepts described in [START_REF] Crassin | Gigavoxels: Rayguided streaming for efficient and detailed voxel rendering[END_REF] and adapts it to the specific needs of multi-channel microscopy on standard computers. More recently, Hoetzlein [START_REF] Hoetzlein | GVDB: Raytracing Sparse Voxel Database Structures on the GPU[END_REF] presented GVDB, a voxel database structure particularly suitable for large volumes with dynamic topology changes to render simulations with a ray tracer on a sparse hierarchy of grids. He also used a particular tree with their own short stack ray tracer. Moreover, several works for storage optimization or efficient construction have been proposed to improve tree structures in this area [START_REF] Baert | Out-of-Core Construction of Sparse Voxel Octrees[END_REF][START_REF] Kämpe | High Resolution Sparse Voxel DAGs[END_REF][START_REF] Laine | Efficient sparse voxel octrees[END_REF].

We can find many approaches based on a tree structure that are used to address out-of-core data (KD-tree [START_REF] Shih | Out-of-core visualization of time-varying hybrid-grid volume data[END_REF][START_REF] Zellmann | Rapid k-d Tree Construction for Sparse Volume Data[END_REF], binary space partitioning (BSP) tree [START_REF] Xue | Efficient volume rendering methods for out-of-Core datasets by semi-adaptive partitioning[END_REF] or octree [START_REF] Crassin | Gigavoxels: Rayguided streaming for efficient and detailed voxel rendering[END_REF][START_REF] Engel | CERA-TVR: A framework for interactive high-quality teravoxel volume visualization on standard PCs[END_REF][START_REF] Gobbetti | A single-pass GPU ray casting framework for interactive out-of-core rendering of massive volumetric datasets[END_REF]) by linking the data through the nodes of the tree. However, Hadwiger et al. [START_REF] Hadwiger | Interactive volume exploration of petascale microscopy data streams using a visualizationdriven virtual memory approach[END_REF] presented a new virtual memory approach to address several petabytes of biomedical data. They focused on electron microscopy volumes with a continuous stream of data. They compared their approach with tree traversal and observed that it scales better to extremely large volume sizes. Their work was extended to the visualization of segmented electron microscopy volumes [START_REF] Beyer | Exploring the Connectome: Petascale Volume Visualization of Microscopy Data Streams[END_REF]. In this extension, the data are segmented beforehand and then stored in an archive before being sent to the GPU on-demand. The segmentation is not performed during the interactive visualization stage on the data cached on the GPU. Jeong et al. [START_REF] Jeong | Scalable and Interactive Segmentation and Visualization of Neural Processes in EM Datasets[END_REF] proposed a volume renderer with on-the-fly processing for noise removal, but only on the currently visible part of the volume. More recently, a visibilitydriven method for on-the-fly filtering of 4D ultrasound data have been presented by Solteszova et al. [START_REF] Solteszova | Output-Sensitive Filtering of Streaming Volume Data: Output-Sensitive Filtering[END_REF][START_REF] Solteszova | Visibility-Driven Processing of Streaming Volume Data[END_REF].

Lastly, for a detailed analysis of out-of-core ray-guided volume rendering, one can refer to Fogal et al. [START_REF] Fogal | An analysis of scalable GPU-based ray-guided volume rendering[END_REF]. They present a study about the optimal brick subdivision, I/O disk access with or without compression and other characteristics analysis. A complete state of the art can be found in [START_REF] Beyer | State-of-the-Art in GPU-Based Large-Scale Volume Visualization[END_REF]. The above-mentioned works are presented with a complete description of different methods of data representation and storage or comparison of address translation approaches. In addition, some works focused on the use of out-of-core methods on clusters or supercomputers, to distribute the data on several nodes [START_REF] Fogal | Large Data Visualization on Distributed Memory multi-GPU Clusters[END_REF].

Although these different solutions address most of the out-of-core data management issues, they do not mention all the needs in a versatile way regardless of the end-user application. For instance, BSP and octree-based approaches may need intermediate structures to be built and thus updated for any change in the data stored in leaves. With animated scenes or scenes where content is unknown in advance, the whole hierarchy needs to be updated. Previous modern visualization-driven approaches proposed cache miss management that is not adapted to algorithms other than visualization. Hadwiger and Fogal, for example, use a set of hash tables to report missing elements organized in the screen space in a ray-guided context. This does not allow to generate cache misses for bricks that are not part of the working set required for visualization.

Besides the previously mentioned approaches, it is important to mention the NVidia's unified memory technology with CUDA [START_REF] Awan | Can Unified-Memory support on Pascal and Volta GPUs enable Out-of-Core DNN Training[END_REF]. This technology has been improved since the Pascal architecture and CUDA 8.0 and now allows to allocate buffers larger than the physical memory available on GPUs. By using this system with the address translation support (ATS) and the POSIX mmap, it is possible to address an entire volume stored on disk from the GPU. In that case, the page faults are managed automatically by the GPU and the OS. However and conversely to our method, this is strictly restricted to the use of the latest NVidia Tesla V100 cards with NVlink2, combined to a UNIX system with the latest kernel (4.16) on a Power 9 architecture [START_REF]Everything you need to know about unified memory[END_REF]. In addition, this technology does not allow to have a control on the page faults and does not allow to use the GPU texture memory (which is mandatory, for instance, in a volume ray casting renderer to take advantage of the interpolation during the sampling).

In the volume ray casting rendering context, Hadwiger et al. [START_REF] Hadwiger | Interactive volume exploration of petascale microscopy data streams using a visualizationdriven virtual memory approach[END_REF] introduced a data structure different from a tree to avoid maintaining and traversing it. Their choice was motivated by the major constraint of a constant flow of data coming from a microscope inducing the non-prior knowledge of the entire dataset. Moreover, their structure is interesting outside its application framework and scales better than octrees for extremely large volumes. The access time to the data is the same regardless of the resolution level, unlike a tree that requires a depth traversal to access the higher resolution data. Nonetheless, they use it only for visualization with a specific volume ray caster. Our method proposes to extend their approach in a general-purpose framework allowing to visualize and process volumes interactively on the GPU. In addition, Hadwiger proposed to manage the content of the caches present in its structure from the CPU, with per-frame read-back for brick usages and requests. To overcome this shortcoming, we take our inspiration from Cyril Crassin's work realized in another context and described in his PhD thesis [START_REF] Crassin | GigaVoxels: A Voxel-Based Rendering Pipeline For Efficient Exploration Of Large And Detailed Scenes[END_REF]. We adapt his method used to maintain an octree to propose a parallel version on the GPU of the management of the virtual addressing structure proposed by Hadwiger. Finally, we present our understanding of the implementation proposed by Hadwiger and give a detailed description of our own.

SYSTEM OVERVIEW AND CONTEXT

The system we propose is illustrated in Figure 1. It is composed of:

• Multi-resolution 3D mipmap pyramid subdivided into small bricks of voxels for each level and stored in a mass storage, obtained in a preprocessing step;

• Large and simple CPU RAM brick cache and a brick loader interface between the disk / CPU and the GPU;

• Virtual memory architecture: multi-level, multi-resolution PT hierarchy and a brick cache on GPU;

• Cache manager entirely on GPU to maintain the caches and offer an efficient cache miss management.

We propose to implement a multi-level, multi-resolution page table hierarchy to address a full volume using a virtual memory management scheme. This hierarchy is composed of several virtualization levels, with a page table for each level, to address very large volumes. The page table levels (except the root) as well as the bricks of voxels are cached in the GPU and managed individually by LRUs. We propose a full management of the GPU caches from the GPU in order to take advantage of the multi-threaded architecture and to limit communications and synchronization between CPU and GPU.

Caches management. This management can be summarized in two main tasks:

1. the updates of the usage information of cached bricks and the updates of LRU caches;

2. the cache miss reports.

Pipeline details. The main feature of our pipeline is to propose a GPU interface for any application that manipulates very large volumes represented as a regular grid of voxels. The navigation inside the volume is performed in a virtual normalized volume. At any time of the run-time application, a voxel can reside in several memory places (mass storage, system memory, or GPU memory). The voxel addressing is done as one uniform address space regardless of its physical storage location. The access to a voxel is determined by a pair (l, p) with l, the desired level of details (LOD) and p, the 3D normalized floating position in the virtual volume (p ≡ [x, y, z] ∈ [0, 1) 3 ). When the application requires a voxel (Fig. 1.(1)), the entire pipeline is triggered as follows: i) (Lookup) The first step is to check in the page table hierarchy if the brick containing the requested voxel is present or not in the GPU data cache (Sec. 5.1 & Fig. 1.( 2)). ii) (Hit) In the case where this brick is 1) is connected to the GPU interface of the pipeline and can require voxels (2) from the caching system. This system is entirely managed on the GPU to handle brick usage (3) and requests (2.1). When bricks are missing, a small request list is sent to the CPU (2.2) to load them asynchronously from the CPU cache manager (2.3) or from a mass storage device (2.4) where a bricked multi-resolution representation of the volume is stored. The bricks are sent back to the GPU (2.5) and the page table hierarchy is updated (2.6).

present in the brick cache, the cache manager can directly provide the requested voxel to the application. Then it reports the usage of that brick (Sec. 5.2). iii) (Miss) In the other case -when a cache miss occurs -a brick request is reported by the cache manager (Sec. 5.2).

Our implementation targets general usage. The optimal cache update strategy depends on the application and will surely be different for ray casting, volume data processing or slice-based volume visualization. On one hand, the cache update will trigger the LRU update from the previous use of bricks (Fig. 1.(3)) and, on the other hand, it will generate a complete request list with all the reported cache misses (Fig. 1.(2.1)). Then a small request list is sent to the CPU (Fig. 1.(2.2)) for an asynchronous processing. An asynchronous thread will request the bricks to the CPU cache manager (Fig. 1.(2.3)). If needed, the bricks are read from the mass storage and written into the CPU cache (Fig. 1.(2.4)). Then, the requested bricks will be written in a buffer accessible by the GPU (Sec. 5.3 & Fig. 1.(2.5)). When the bricks can be written in the GPU data cache, the page table hierarchy is updated accordingly (Sec. 5.6 & Fig. 1.(2.6)). The application efficiently manages the voxel requests, by only asking for the necessary bricks, so as never to load unnecessary ones (for instance, adopting a ray-guided method for volume ray casting).

This pipeline is relatively similar to that proposed by Hadwiger et al. [START_REF] Hadwiger | Interactive volume exploration of petascale microscopy data streams using a visualizationdriven virtual memory approach[END_REF]. However, its design differs in several aspects because the data accessed is not a constant stream output by a microscope. Therefore, we have an a priori knowledge of the dataset, and we store the bricks directly (e.g. 32 3 -256 3 voxels) rather than storing 2D tiles that require registrations and a step of 3D construction of bricks at run-time.

Communications. In our system, the communications from CPU to GPU are limited to the strict minimum. They consist of sending the bricks of voxels and a Boolean flag for each of them (for an empty brick information). The communications from GPU to CPU are restricted to a small request list containing (l, p) pairs of each requested bricks. All other mandatory actions of the out-of-core data management are performed inside the GPU, taking advantage of its computing power and limiting costly communications with the host.

DATA REPRESENTATION

The input data may come from a physical model acquisition (e.g. biomedical scans, MRI, etc.) as well as from the voxelization of any data as long as they can be represented as a regular 3D grid (Sec. 9).

Multi-resolution. Raw data volumes are preprocessed in order to create a 3D mipmap pyramid which represents the volume at different LODs. Data anisotropy in sampling, which often occurs in medical imaging, can be drastically reduced by applying different downsampling ratios along each axis in this multi-resolution representation.

Bricking. Each level of the multi-resolution pyramid is then divided into small independent blocks of voxels, called "bricks". All bricks have the same size, regardless of their resolution level. This object space decomposition method allows for the manipulation of small amounts of data rather than the whole large volume. The brick shape does not necessarily have to be cubic; different edge sizes can be used. In any case, when the volume is not sufficient to complete all bricks, the last bricks are simply completed with empty voxels. Our approach avoids generating unnecessary empty bricks like most approaches that manipulate only power of two volume sizes. In addition, an overlapping area can be created, between neighbouring bricks, to ensure coherence of future treatments on the voxels (interpolation, convolution, etc.).

Our system is designed for the use of volumes whose entire dataset is already defined and known. It means that treatments could be performed on the data before using them in the system. As we are dealing with large volumes of data, both steps mipmapping and bricking take a substantial time to be done; therefore it is interesting to perform these steps outside the stream.

Compression. Finally, the bricks are stored on a large storage device in raw data or compressed format. Using a compressed brick format has two advantages: the volumes will take less memory space on the storage; and it will reduce the loading time of the bricks from the mass storage to the CPU (Sec. 8.4). Nevertheless, we are targeting a real-time application using scientific images; this means the compression algorithm needs to be lossless and has to provide a fast time decompression. To this end, we have opted for an LZ41 compression scheme.

OUT-OF-CORE DATA MANAGEMENT

We detail in this section the out-of-core addressing data structure used to access all the dataset. We present its implementation on GPU and the parallelization strategy of its associated cache manager. We also propose our implementation version and technical details for the virtual hierarchy management and update.

Memory virtualization

In order to address the bricks of a large multi-resolution volume from the GPU, we adopt an approach of multi-level, multi-resolution page table hierarchy. This structure is based on a large page table organized on several levels, which consists of mapping the virtual address of each voxel with its physical address. This can be seen as a pyramidal structure where each level virtualizes a set of entries (e.g. 16 3 -64 3 ) organized into 3D blocks in the level below (Fig. 2). The top of the hierarchy is composed of the multi-resolution page directory (MRPD) which corresponds to the starting point of the virtual addressing of any voxels, at any resolution level. This is the only level which is not virtualized and entirely present on the GPU with a low memory impact (Sec. 8.4). For very large volumes of data, the MRPD itself is too large to be held in GPU memory. The virtualization principle can therefore be applied to the page table itself. Hence, below the MRPD, there may be N intermediate PT levels with, for each one, a cache containing PT blocks. Adding intermediate levels increases the amount of virtualization and reduces the size of the MRPD. Finally, the data cache containing the bricks is located at the bottom of the hierarchy. The theoretical complexity of the algorithm to access data, given its tree-based data structure, is O(log(b)), where b is the number of bricks at leaf level. However, the constant of this complexity is very small, and in practice, with two intermediate levels, petabytes of data can be addressed. In addition, the access time does not increase with LOD (unlike an octree).

Our implementation of this hierarchy has been planned generically in order to dynamically create as much intermediate PT level as necessary to address any volume size. This allows the use of smaller bricks and therefore a finer data addressing.

Note, it is important not to confuse the multi-resolution representation volume and its virtualization used for the address indexing on the GPU. While, the former refers to the resolution levels l (LODs), the latter refers to the virtualization levels.

Implementation

Each hierarchy level is implemented with the CUDA surface API [START_REF]Nvidia CUDA programming guide 8.0[END_REF] in order to read and write elements into 3D texture memory with hardware accelerated 3D indexing. In the case of a graphic end-user application, the spatial coherent access patterns of the surface memory have a significant benefit. All the caches are managed by an LRU implemented as a simple device vector with the Thrust parallel template library [START_REF]Nvidia CUDA programming guide 8.0[END_REF].

The data cache implementation is templated in order to store any type of voxels from GRAY8 to RGBA32 and others (Sec. 9). Conversely, a page table entry is represented with four 32-bit integers to store the 3D coordinates c of the block virtualized by this entry in the next hierarchy level and a flag representing the state of this block. This flag can be either empty, mapped or unmapped (not present on the GPU).

The MRPD is represented in memory as a buffer containing one 3D subgrid per LOD added next to each other on the x-axis. This representation is illustrated in figure 4. The size of the MRPD is equal to the sum of the x-size of each subgrid and the size of the largest subgrid in the yand zaxes. This representation implies that a part of the buffer will not be used.

A cache manager is used to maintain the bricks and page tables usage as well as to manage the brick requests. These operations are done in parallel on the GPU with data parallel primitive stream compaction operations.

LRU update. Updating the used bricks is necessary for caches management. This is assured by an LRU algorithm on a list for each cache. All the GPU threads report the use of elements by marking the corresponding elements, in a shared 3D buffer (called usage buffer) on global GPU memory, with a global 32-bit integer timestamp. This timestamp is incremented after each update; it is not necessary to use atomic writes to handle concurrent threads access because each thread will write the same timestamp integer in the corresponding entry. There is one such buffer for each cache and it contains as many entries as elements present in its corresponding cache (one to one 3D mapping). Then, a mask is created and filled with a Boolean value that indicates if the corresponding usage buffer entry contains the current timestamp or a previous one. The stream compaction ensures the sorting of the lists by moving the most recently used elements to the beginning while maintaining the order of the others (Fig. 3.(a)).

Cache misses. Figure 3.(b) illustrates the creation of the requested brick list. When a requested brick is not present in the GPU (flag set to unmapped), a cache miss occurs and triggers a request for this brick. In the same way as the element usage information, the requested bricks are maintained in a list through a stream compaction operation from a buffer (called request buffer) shared by all the CUDA threads. This buffer contains as many entries as bricks in the whole multiresolution volume, and each thread will tag the corresponding entries to the required bricks with the global integer timestamp. The resulting request list is limited to a small size in order to limit the number of bricks loading at each update. This list contains the requested bricks IDs. These IDs are guaranteed to be unique and are constructed from the spatial position in the volume and the resolution level of the bricks. The proposed cache miss management strategy is entirely generic. It makes it possible to request any brick at any resolution level and at any time without being oriented towards a specific application type.

Comparison with an octree. With an octree approach, the PT entries of the virtual memory management are the octree nodes themselves (stored in a node pool). This implies to transfer these nodes from the CPU to the GPU memory by raising PT entry requests (node requests). In the proposed method, the PT entries are updated directly in the GPU, thus avoiding transfers to the device. The single communications between the central memory and the GPU texture memory are data (bricks) transfers. However, our request buffer has to be sized to the number of total bricks. For instance, let us consider a large volume of 16384 3 voxels RGBA (which is about 17.6 TB) with a ratio of two between two consecutive LODs and a brick size of 64 3 voxels. Since the request buffer is represented in memory in the same way as the MRPD, it will require, with nine LODs entries the z-dimension. With a timestamp on 4 bytes, the total size of the request buffer will be around 134 MB. With a data cache of 4 GB (which is correct for modern GPUs), it represents less than 3.5% overhead in memory.

Data fetch

The brick requests are asynchronously handled by the CPU to fetch the data. As shown by Crassin [START_REF] Crassin | GigaVoxels: A Voxel-Based Rendering Pipeline For Efficient Exploration Of Large And Detailed Scenes[END_REF], it is more efficient to use CUDA zero-copy [START_REF]Nvidia CUDA programming guide 8.0[END_REF] in order to load the bricks into GPU memory. It ensures optimal uses of the hardware rather than transferring those manually using copy operations triggered by the CPU. To achieve this, we use a requested buffer allocated in CPU memory with CUDA; this buffer is pinned on a physical address region for the PCI Express controller to use directly. A dedicated CPU thread takes care of fetching bricks from the CPU cache or the mass storage and signals the GPU that bricks are available in this buffer. The writing of the requested bricks from the requested buffer to the GPU data cache is performed with a single pass in a CUDA kernel with one thread per voxel per brick.

Empty bricks. Empty brick information makes possible to distinguish the bricks that do not participate in the running algorithm because they are considered as empty, transparent or without any important information. It is important to never load such bricks on the GPU to avoid cluttering the data cache. The access to a brick whose flag is set to empty in the virtual addressing structure does not generate a cache miss. However, a brick flagged to unmapped, which generates a cache miss, may potentially be empty. Figure 3.(c) shows the adopted strategy to efficiently deal with this behaviour. For each request of brick, we test on the CPU whether it is empty or not (regardless of which criteria). Only non-empty bricks are added in the requested buffer. A Boolean vector indicating for each requested brick whether it is empty or not is transferred to the GPU. With this information, it is possible for the GPU to identify the bricks to transfer and to write inside the data cache. Thus, a dual stream compaction operation is used on the initial brick request list from the empty information mask. This step generates a sorted list to discriminate non-empty from empty bricks. The former will be added to the data cache and referenced as mapped in the hierarchy. The latter will be referenced with the empty flag.

Voxel access

This section is a reminder of the method presented by Hadwiger et al. [START_REF] Hadwiger | Interactive volume exploration of petascale microscopy data streams using a visualizationdriven virtual memory approach[END_REF]. For the sake of reproducibility, we adopt a formal way of presenting the involved concepts. In particular, we present the virtual volume representation in section 5.5 and the update of the virtualization hierarchy in section 5.6.

Using figure 4, a voxel is addressed from a LOD l and a position p ∈ [0, 1) 3 , with the following two equations:

PDentry(l, p) = PageDirBase[l] + p v[l] 0 ∈ N 3
(1)

Centry n (l, p) = Cbase n (l, p) + p v[l] n mod b n ∈ N 3 (2)
with

Cbase n (l, p) = MRPD[PDentry(l, p)].c if n = 1 PTCache n-1 [Centry n-1 (l, p)].c if n > 1 (3)
Equation ( 1) returns the 3D address to look at in the MRPD. In this formula, PageDirBase[l] corresponds to the 3D position of the beginning of the subgrid for the resolution level l in the MRPD and v[l] 0 its number of entries (3D size). The operator is the elementwise vector floor, i.e. For instance, the access to the highlighted voxel (voluntarily restricted to the plane z = 0) of figure 4.(a), using the configuration of figure 4.(b), will be:

PDentry(0, 1 8 , 1 4 , 0 1 ) =      x = 0 + 1 8 × 4 = 0 y = 0 + 1 4 × 2 = 0 z = 0 + 0 1 × 1 = 0 It is important to note that
even if p is composed of real values, the results are given as integer values which is required to match with the positions in the caches. The entry located at PDentry(0, [ 1 8 , 1 4 , 0 1 ]) = [0, 0, 0] in the MRPD contains the beginning of the entry (here a brick) in the data cache and a flag set to mapped. In our case the entry will be [x, y, z] = [2, 0, 0]. We can now calculate the coordinate of the voxel in the data cache:

Centry 1 (0, 1 8 , 1 4 , 0 1 ) =      x = 2 + 1 8 × 8 mod 2 = 3 y = 0 + 1 4 × 4 mod 2 = 1 z = 0 + 0 1 × 1 mod 1 = 0
By accessing the coordinates [x, y, z] = [3, 1, 0] in the data cache, we can read the values of the requested voxel. The behaviour will be the same independently of the number of intermediate caches (Fig. 4.(c)) with equation ( 1) for the MRPD and equation ( 2) for all other caches.

Real vs. virtual volume representations

The virtualization hierarchy may introduce changes in the size of the volume. These changes are induced by the successive block sizes b n between the different virtualization levels n. Hence, a volume of a given size may be virtually larger than its real size. To express these differences we note the size of the real volume r[l] n and the size of the virtual volume v[l] n . The processing of the real and the virtual sizes (r[l] n and v[l] n ) is performed in two steps using the following equations:

r[l] n = r[l] n+1 b n+1 (4)
and

v[l] n = r[l] n if n = 0 v[l] n-1 b n if n ≥ 1 (5)
The first step uses equation ( 4 Position. The previously mentioned vector p ∈ [0, 1) 3 represents the 3D normalized coordinates of a voxel in the virtual volume representation. Thereby, a position p = [0.5, 0.5, 0.5] will refer to the center of the virtual volume which is not necessary the center of the real volume.

To correct this, one needs to change its coordinate system as follow:

p real = p virtual (r[l] n max v[l] n max )
Block size. The choices of the different b n are at the user's discretion. We chose small values b n for this paper in order to simplify the understanding. In real situations the b n are commonly sized from 16 3 to 64 3 . The greater these values, the lower the MRPD, however, larger blocks will require more memory space. A balance must be struck between the b n , the size of the caches and the size of the MRPD.

Hierarchy updates

Figure 5 illustrates the update of the PT hierarchy according to the configuration of figure 4. This update consists to load new bricks into the data cache and to reference them. It follows three sequential steps: i) removing the references of the bricks that will be removed from the GPU data cache (if the cache is already full), ii) writing the new bricks in the GPU data cache and updating its LRU, and iii) referencing these new bricks in the PT hierarchy.

For instance, let us suppose we request the brick containing the values of the voxel positioned at l = 0 and p = [ 2 8 , 0 4 , 0 1 ] from figure 4.(a). To add this new brick in the data cache illustrated in figure 5.(c), it is required to dereference the oldest brick given by the LRU (the one positioned at coordinates [0, 0, 0] in the data cache). Using equation [START_REF]Nvidia CUDA programming guide 8.0[END_REF] with the coordinates of the oldest brick, it is possible to set the flag of its corresponding entry in the MRPD to unmapped (Fig. 5.(a)). The brick is then considered as not present in the cache. Once this first step is done, the new brick is loaded in place of the oldest one for the second step (Fig. 5.(c)). Finally, knowing the coordinates of the new brick and using equation ( 1), we set in the corresponding MPRD entry the coordinates of this brick in the data cache and its flag to mapped (Fig. 5.(a)).

In case of intermediate PT caches level (Fig. 5.(b)), the behaviour is the same. However, if we want to load the brick containing the voxel l = 0 and p = [ 2 8 , 0 4 , 0 1 ], it will be necessary to add a new PT block in the PT cache in order to be able to reference the added brick. To add a new PT block, the oldest one (given by the LRU) is dereferenced from the level just above. Then a new one is generated with all flags set to unmapped. Finally the oldest one is replaced by this new one which is then reference in the level just above.

ON-DEMAND PROCESSING DURING VISUALIZATION

The proposed approach allows the introduction to on-demand data processing algorithms on visible or not currently visible data during the visualization stage. A processing stage is initiated by an action of the user during the interactive visualization. It can be launched in a different CUDA stream in order to be able to run in parallel (kernel overlapping) and to be independent of the visualization algorithm. Both algorithms (visualization and processing) can request bricks at the same time using the shared request buffer without concurrency problem. Regardless of which step it is generated, these requests will be handled in the same way by the cache manager. The processing step can impact the data in different ways. It can provide global information that does not modify the data itself (e.g. number or size of an element). It can indirectly modify a voxel by providing information about it (whether or not a voxel is present in a set, e.g. classification), which will be then interpreted by the visualization step to possibly treat this voxel differently (e.g. applying a certain color). Finally, it can modify the value of a voxel (e.g. filtering). This modification can be permanent: the new value of the voxel is written back in the cache and then propagated with a CPU read-back to be written on the storage; or not permanent: the treatment should be reapplied to this voxel if it is removed and put back later into the cache. If it is necessary to propagate the modification of the value of a voxel to all the LODs, it is up to the application to support this behaviour by requesting the corresponding brick in each LOD present in the multi-resolution volume representation. It is necessary to build as many addresses (l, p) as resolution levels l, keeping the same p (the normalized voxel position in the virtual volume).

Working set. Figure 6 shows the different configurations of the possible working set. An approach that only allows visualization or processing is represented by the first configuration (Fig. 6.(a)). The working set can be defined by one of the other representations in the context of on-demand processing during visualization. A processing could operate only on the data that are visualized during its initialization (Fig. 6.(c)). In that case, the data involved in the processing algorithm are already present in the GPU data cache (or already requested by the visualization stage). Conversely, if the processing algorithm sets up data that are not part of the visualization-specific working set (Fig. 6.(b, d)), this can lead to new brick requests. This may increase the risk of having a working set that is too large to be stored in the cache. In practice, the last configuration would probably be the least frequent. However, it may be found in a distributed environment with a subset of the resources dedicated to the visualization and possibly another subset in charge of processing a part of the volume that is not visualized.

COMPARISON WITH OTHER SYSTEMS

By contrast to above-mentioned approaches (Sec. 2), our proposal is the first system that uses a virtual memory architecture based on a multi-level, multi-resolution page table hierarchy with a complete management on the GPU. In this way, in addition to taking advantage of its multi-threaded environment, we designed a system with communications between CPU and GPU that are reduced to their bare minimum.

In the work of Hadwiger et al. and Fogal et al. [START_REF] Fogal | An analysis of scalable GPU-based ray-guided volume rendering[END_REF][START_REF] Hadwiger | Interactive volume exploration of petascale microscopy data streams using a visualizationdriven virtual memory approach[END_REF] the use of the same type of structure could be found but with a cache management carried out on the CPU. Hadwiger et al. attached a use-bit to each brick present in the cache to report brick usages. To report cache misses in their parallel volume ray casting application, they hosted on the GPU, one hash table per image tile to avoid contention and synchronization.

For an image tile size of 64 × 64 with a full HD frame buffer, this sums up to 510 hash tables with expensive atomic operations to report brick requests. Those hash tables are then sent back to the CPU to handle the requests. It is possible to quantify the transfer of these hash tables. For this, it is necessary to note that they used a 32-bit or a 64-bit ID to report each cache misses, depending on the size of the volume. They also limit the number of requests per ray to a small number M in order to distribute cache misses over multiple frames. This represents a message set of about 1920 × 1080 × 4 × 4 ≈ 33 MB (in the worst-case scenario) for an HD viewport with M = 4 and a 32-bit ID, at each frame. It reaches 3840 × 2160 × 4 × 4 ≈ 133 MB for a 4K viewport with the same configuration. These values may have to be multiplied by N in a context of visualization with N views (stereo or multiscopic) [START_REF] Sarton | Virtual Review of Large Scale Image Stack on 3D Display[END_REF]. For instance, a 4K screen with 11 views it represents 133 × 11 ≈ 1.4 GB. These quantities are smaller in the case of Fogal et al. since they use a compressed format of lock-free hash tables. We propose to limit globally the request number instead of a per-ray limitation. This allows us to work in a general-purpose context and to be more scalable on the cache miss transfers. As already mentioned, we draw our inspiration from Crassin [START_REF] Crassin | GigaVoxels: A Voxel-Based Rendering Pipeline For Efficient Exploration Of Large And Detailed Scenes[END_REF] for a comprehensive cache management on the GPU. For the brick usage, we proceed in the same way. However, one notable difference compared to his work is the use of a tree structure with nodes as page tables to address the data bricks. In his system, the octree nodes (equivalent to our page table entries) are requested by the GPU and induce communications between the central memory and the video memory. Moreover, in contrast to previous systems that are particularly focused on volume ray casting, our approach takes place in a general-purpose context.

EVALUATION

In order to test our GPU cache system, we developed two kinds of visualization applications and a processing stage. The first is a 2D rendering virtual microscope with a navigation through an image stack [START_REF] Sarton | Virtual Review of Large Scale Image Stack on 3D Display[END_REF] that provides easy access to high-resolution images. The second is a GPU-based volume ray casting [START_REF] Battin | A Premixed autostereoscopic OptiX-based Volume Rendering[END_REF]. Furthermore, to test our system with on-demand processing, we offer a convolutional image processing algorithm.

To illustrate the results, we used three datasets (Fig. 7):

• (a) A 114 histological slices stack of a mouse brain with a resolution of 64000 × 50000 RGBA pixels (1.459 TB);

• (b) A 3D Mandelbulb fractal of 4352 3 RGBA voxels (330 GB) generated with Mandelbulb3D;

• (c) A 2160 × 2560 × 1072 volume with grayscale 16 bits voxels (11 GB) of a primate hypocampus from a light sheet microscope.

All the tests were carried out on an NVidia GeForce Titan X with 6 GB of VRAM, a CPU Intel i7 4790K 4 GHz, 32 GB of RAM. We used CUDA 8.0 and OpenGL 4.5 interoperability to render on a full HD viewport 1920 × 1080.

Use case 1 -Virtual microscope

The first developed application allows us to navigate (pan and zoom) in a stack of images with 2D multi-resolution rendering. The interest is to simulate the behaviour of a microscope to visualize and navigate into large slices. The principle is to compose a 2D texture from the bricks intersecting all or part of the plane perpendicular to the z-axis in the volume considering the 3D camera position in the volume. We opt for a strategy that promotes the quality of the visual feedback for the user. In this sense, when a request is lifted because a brick of resolution level l is missing, the cache manager provides our renderer with a lower resolution brick (if there is one in the cache) until the brick of level l arrives in the data cache.

Use case 2 -OptiX™ based volume rendering

In order to maximize the GPU load (cache and end-user application included), a ray casting volume renderer module based on the NVidia OptiX [START_REF] Parker | OptiX: A General Purpose Ray Tracing Engine[END_REF] engine was developed. This module is based on the proposed GPU volume ray casting of Kruger and Westermann [START_REF] Kruger | Acceleration Techniques for GPU-based Volume Rendering[END_REF] and uses a modern ray-guided system to ensure a visibility-driven rendering of multi-resolution volumes. The sample integration scheme simulates an emission-absorption model based on the equation presented by Max [START_REF] Max | Optical models for direct volume rendering[END_REF].

On-demand image processing

In order to test the system with on-demand processing, we have developed a convolution image processing algorithm. This is an edge detection algorithm that is integrated into the main rendering loop. At any time, the user can initiate it by interacting with the visualized data.

It is not relevant here to provide details about this algorithm which can be seen as a black box. This algorithm was used to validate the correct behaviour of an on-demand processing step during the interactive visualization phase, especially to test the sharing of out-of-core data request management of such a process with the rest of the application.

Performance analysis

Figure 8, obtained by benchmarking the application described in section 8.1, presents the comparison of the time spent on three operations: the rendering and the caches management, both performed on the GPU, and the bricks loading processed asynchronously the CPU. The caches management operations cover: the LRUs updates and the data request management done at each iteration, the writing of the bricks in the data cache and the hierarchy updates done when new bricks are sent to the GPU. These times are calculated on the datasets (b) by performing the following navigation scenario: zooming from the lowest to the highest resolution level, navigating through the slices of the image stack along the z-axis, pan navigating within a slice (x-,y-axes), and an edge detection convolutional algorithm for image processing. A low amount of time without any interaction is deliberately maintained between each step in order to illustrate the null bricks loading time during these periods. Throughout the entire scenario, the rendering time is constant and very low (∼ 0.2 ms). Moreover, we note that our application allows the rendering to be carried out with a rate of about several hundreds of frames per second. Such high frame-rates can be achieved because the proposed application does not require significant resources at the scale of the data management. In addition, this time is constant because the loading of bricks is performed asynchronously and does not interfere with rendering. The caches update times are less than 2 ms with additional peaks according to the brick loading. These update times are reasonable and allow interactive navigation in the scene. Furthermore, in the case of an application that requires more computational resources, they become negligible.

This scenario allows us to validate an on-demand image processing step within an interactive visualization application. Thus, we can see in figure 8 some of the brick loading peaks that emanate from requests caused by this processing step. Once again, we can note that the rendering time is not affected by this step. As a consequence, these two steps can cohabit by using our out-of-core management system, both at the same time.

Despite the fact that the entire management of the caches and the virtualization hierarchy are realized on the GPU, its occupancy is around 5%. The bottleneck occurs at the data loading, between the mass storage and the CPU. The computation footprint on the GPU is small enough to allow heavy rendering workloads.

Bricks loading. We measured i) the average loading time for one brick, from its request by the GPU, to its rendering (or processing), and ii) the average time to render a complete correct view (with all the needed bricks addressed in the GPU cache). These two measurements, shown in the table 2, were performed from empty GPU caches. We studied two scenarii: firstly, a worst-case scenario where the bricks are only present on the mass storage. Secondly, when all the bricks are present in our intermediate CPU cache. In the first scenario, the given times indicate: 1) the bricks requests management (request list creation and sending to the CPU), 2) the disk access time, 3) the LZ4 decompression, 4) the CPU to GPU transfer and 5) the update of the addressing hierarchy (page table referencing and LRU update). The second scenario only includes the steps 1, 4 and 5. Note that like Hadwiger [START_REF] Hadwiger | Interactive volume exploration of petascale microscopy data streams using a visualizationdriven virtual memory approach[END_REF], our worst-case scenario takes few seconds to generate a complete view.

Ray casting performance. The second part details the GPU memory footprint for all the resources needed in the out-of-core data management. The third part gives the average brick loading times for i) one brick and ii) all the bricks needed to a complete full HD view. The last row gives the frame rates for the volume ray casting rendering. The FPS measurement environment is described in section 8.4 and the virtualization parameters used to measure the amount of memory, depending on the datasets, are described in section 8.5.

age frame rates of the volume rendering using the application described in section 8.2. These measurements used a linear transfer function and a full screen rendering of the highly zoomed volumes. These frame rates can be considered at least as efficient as the previous approaches similar to ours. Therefore, we do not lose in rendering performance using a volume ray caster while proposing a more general method allowing different types of visualization and / or processing.

Memory occupancy

The second part of the table 2 shows the required memory on the GPU for all the elements of our structure and its management (without the caches themselves). For (a) and (b), we used a PT cache storing 500 PT blocks in the second virtualization level and a data cache storing 2860 bricks (∼ 3 GB). For (c), we use a data cache storing 5720 bricks (∼ 3 GB). The memory occupancy is about 12.65 MB for the first dataset (a) with two virtualization levels, using a brick size of 64 3 and a PT block size of 32 3 voxels; 2.65 MB for the second (b) with two virtualization levels, using a brick and a PT block size of 64 3 voxels; 1.11 MB for the last (c) with only one virtualization level and a brick size of 64 3 voxels. The low costs of these values maintain a limited pressure rate between the memory allocated for the caches management and those necessary for a well-functioning state of the end-user application.

DISCUSSION

The proposed approach offers benefits for a wide range of applications, but suffers a few limitations detailed here.

Limitations

Race conditions. This system implies some limitations regarding race conditions. Indeed, the thread in charge of loading and preparing the bricks cannot load them in the cache and update the virtual hierarchy without inducing possible collision with the GPU rendering (or processing) thread.

Page table hierarchy updates. There are two main limitations related to the update of the PT hierarchy. The first occurs when an entry is removed from a PT cache. In figure 5.(b) when the new PT block overwrites the oldest one, we lose the link to all the bricks that the last one was referencing. In essence, some data may still be present in the caches but not referenced (garbage that is not taken out). However, these orphan blocks/bricks will converge to the oldest part of the LRUs and will be replaced in future iterations. The second limitation is related to the way to reference, in the PT, the bricks recently added to the data cache. Although we propose an implementation of this step on the GPU, it is done sequentially for each brick. For instance, let us assume we want to add, at the same time, two spatially close bricks, both being referenced by a single PT block. A parallel system would add two PT blocks in the cache that would both reference a single brick while only one block may be required to reference these two new bricks. To avoid such a scenario, each brick needs to be referenced one by one.

Loss of virtual pages. The proposed approach of virtual memory architecture does not guarantee the meta-information perennity of pages inside the page table. For the possible intermediate virtualization levels (PT caches), there is a problem of information preservation. The pages are created on-the-fly on the GPU as needed and then cached. When a page is removed from the cache, all its contained information is lost. A typical example is the empty brick state indicated by the empty flag which will not be preserved. In contrast, this problem does not occur with the use of tree structures (like octrees) since the nodes of the tree are actually used as pages. These nodes are written and saved and then transferred to the GPU if necessary. Processing. Our system can support any data processing if it is used alone. However, if it is initiated during an interactive visualization stage, it needs to be restricted. Indeed, an algorithm manipulating the whole volume will generate many cache flushes, which may severely impact the possibility to maintain the interactive visualization working set stored in the data cache. Thus, we must limit the processing to a lower working set size corresponding to a set of smaller local treatments.

Strengths

With this data structure, the system gets access to the data in O(log(b)) where b is the number of bricks. However, the complexity constant is such that in practice it dwarfs the access time to a very small value (two or three levels in practice), as detailed in section 5.1.

Communications. The communications between the CPU and the GPU are restricted in order to reduce to data flow. The information transferred is simply the requested brick coordinates for the GPU to CPU communications and the voxels themselves, organized in bricks, for the CPU to GPU communications. Restricting the communications with these two kinds of data prevents bottlenecks.

Data Type. With the proposed approach, it is possible to virtualize any kind of data as long as these data can be represented as a regular 3D grid. Our approach uses the CUDA texture memory that can store up to 4 channels with a maximum of 32 bits each. However the approach also applies to general data configuration with the use of global memory instead of CUDA surfaces to extend to a more general purpose context; the same algorithms and the logic still apply.

Scalability. This system may be easily scaled in a high performance computing environment: cache managers deployed on each node are independent, each one handling a subworking set. The main difference would lie in the cache misses, where the nodes could communicate together in order to speed up the brick loading time.

CONCLUSION & PERSPECTIVES

We propose an out-of-core caching system fully managed on GPU to address very large volumes of data. Our pipeline can be used by any type of application that requires manipulating volumes represented as a regular grid of voxels, with the ability to manipulate multi-modal data. In particular, this allows us to introduce the possibility of on-demand processing of visible or not currently visible data during an interactive visualization stage. We demonstrate the validity of our implementation using a 2D high-resolution slice visualization tool with interactive navigation in volumes exceeding the amount of GPU and CPU memory with a very high rendering frame-rate and a low memory footprint on the GPU. We also show that our method performs at least as effectively as the previous ones in a context of volume ray casting application.

Hierarchy updates. In the current approach, the PT hierarchy updates are performed on the GPU but in a sequential context (one thread kernel). It prevents duplications in the load of bricks or race conditions in the PT entry updates but it does not take advantages of the parallelism of the GPUs. It could be interesting to consider a parallel strategy to improve performance of this update step.

Transfer times. The main bottleneck in the system is the brick loading time from the storage device to the GPU. Data compression addresses part of this problem, but it could be more interesting to be able to decompress the bricks directly on the GPU. This could take advantage of the multi-threaded environment to implement an efficient decompression algorithm and reduce the amount of transferred data to the GPU. As many bricks are loaded, there is a natural parallelism that can be exploited. Another solution would be to bring this system to a high performance computing environment. By properly sharing the work, the different nodes of such an environment could then distribute the transfers and thus reduce this bottleneck.

Cache misses. We could consider a strategy to improve the cache miss management system by prioritizing requests. For instance, with a counter we could determine which bricks have been the most requested. Such bricks would have a higher priority for the loading step. This could potentially improve visual quality of the rendering for the same overall performance level of the system.

Fig. 1 .

 1 Fig.1. Overall out-of-core pipeline. An end-user application (1) is connected to the GPU interface of the pipeline and can require voxels (2) from the caching system. This system is entirely managed on the GPU to handle brick usage (3) and requests (2.1). When bricks are missing, a small request list is sent to the CPU (2.2) to load them asynchronously from the CPU cache manager (2.3) or from a mass storage device (2.4) where a bricked multi-resolution representation of the volume is stored. The bricks are sent back to the GPU (2.5) and the page table hierarchy is updated (2.6).

Fig. 2 .

 2 Fig. 2. Multi-resolution, multi-level page table hierarchy virtualization principle. A 2D representation of our 3D addressing structure. In this example, only one intermediate virtualization level is used (PT1). One entry of the MRPD addresses a PT block of 3 × 3 entries in the PT1 cache and one entry in the PT1 cache addresses a brick of 4 × 4 voxels in the data cache.

Fig. 3 .

 3 Fig. 3. Parallel cache manager mechanisms on the GPU (a) Cache usage update: a single pass to update LRU with double stream compaction (Cartesian products) from a shared buffer to report bricks usage. (b) Bricks request list creation: a single stream compaction from a shared request buffer to create a small request list. (c) Empty bricks management strategy at loading.

  [x, y, z] = [ x , y , z ]. The operator is the element-wise vector multiplication, i.e. [x, y, z] [a, b, c] = [xa, yb, zc].If the flag contained in the entry at the coordinates PDentry(l, p) of the MRPD is set to mapped, we can continue to equation[START_REF]Everything you need to know about unified memory[END_REF]. This second gives the 3D address to look at in all different caches (PT and data caches). The operator modulo is defined on vectors as [x, y, z] mod [a, b, c] = [x mod a, y mod b, z mod c]. The value n corresponds to the virtualization level from n = 0 (MRPD) to n = n max (data volume) and passing through intermediate cache levels (0 < n < n max ). Since equation (1) computes the case n = 0, n goes from n = 1 (the first cache after the MRPD) to n = n max -1 (the data cache) in equations (2) and (3). Hence, as shown in the equation (3), Cbase n (l, p) is a vector containing 3D cache coordinates c. It is stored in the MRPD for n = 1 or in the cache just above otherwise. This vector will indicate the 3D position of the beginning of the block in the cache below. v[l] n and b n are respectively: the virtualized volume size of LOD l according to the virtualization level n and the size of a PT block, for the PT cache of level n; or the volume size in voxels and the size of a brick, for the data cache (Sec. 5.5).

Fig. 4 .

 4 Fig. 4. 2D representation of voxel addressing in the virtual hierarchy. (a) Volume used as an example, represented with three LODs cut in 2 × 2 voxel bricks. It contains 8 × 4 voxels at the highest resolution level (L0), 4 × 2 voxels at the intermediate LOD (L1) and 2 × 1 voxels at the coarsest LOD (L2) filled with empty voxels (black) to get one entire brick. (b) Data structure to address the voxels using one virtualization level. One entry of the MRPD addresses one brick. (c) Data structure to address the voxels using two virtualization levels. One entry in the MRPD addresses one PT block of 1 × 2 entries in the intermediate PT cache. One entry in this PT cache addresses one brick in the data cache. In both cases (b) and (c), we try to get the value of the voxel positioned at l = 0 and p = [ 1 8 , 1 4 , 0]. In the MRPD and the PT cache, white means the flag is set to mapped and black means it is set to unmapped. The circles represent the PT entries and the squares represent the voxels.

  ) and computes the real dimensions of the volume for each virtualization level. Starting from the dimension of the volume in voxels r[l] n max , this allows us to obtain the dimension of the MRPD r[l] 0 according to the block size at each virtualization level n. The operator is the element-wise vector division, i.e. [x, y, z] [a, b, c] = [ x a , y b , z c ]. The operator is the element-wise vector ceil, i.e. [x, y, z] = [ x , y , z ]. The second step, using equation (5), computes the virtual sizes v[l] n from the MRPD v[l] 0 = r[l] 0 to the volume v[l] n max in voxels according to the different block size b n . Tables 1.(a, b) indicate the real and the virtual dimensions of the volume shown in section 5.4 and figure 4. The table 1.(b), especiallyfor LOD L1 and L2, illustrates the fact that virtual sizes can be larger than real sizes, due to the ceiling operator in equation (4).

Table 1 .

 1 Real and virtual sizes of the volume illustrated on figure 4.(a). Table (a) corresponds to the configuration of figure 4.(b). Table (b) corresponds to the configuration of figure 4.(c). The r[l] n and the v[l] n are the virtual representation of the volume at the virtualization level n. For a PT cache level, v[l] n corresponds to the set of PT entries mapped on the volume representation at the virtualization level n.

Fig. 5 .

 5 Fig. 5. 2D representation of the hierarchy update. Update of the PT hierarchy when adding the brick containing the voxel positioned at l = 0 and p = [ 2 8 , 0 4 , 0 1 ] in the configuration of figure 4.(a). The status before the update is shown on the left of the arrows and the status after the update on their right. (a) The update of the MRPD in the case of one virtualization level (according to figure 4.(b)). (b) The update of the MRPD and the intermediate PT cache (and its LRU) in the case of two virtualization levels (according to figure 4.(c)). (c) The update of the data cache (and its LRU) in both cases.

Fig. 6 .

 6 Fig. 6. Working set configurations. The blue set represents the whole multi-resolution volume data, the green and the red ones represent the visualization and/or processing working sets. (a) a unique application (visualization or processing), (b, c, d) on-demand processing during interactive visualization.

Fig. 7 .

 7 Fig. 7. Datasets used for system validation and performance analysis. (a) A 1.459 TB histological slices stack render with our high-resolution 2D slices visualizer, (b, c) a 330 GB Mandelbulb fractal and a 11 GB primate hypocampus from a light sheet microscope, render with our own OptiX-based DVR solution [27].

Fig. 8 .

 8 Fig. 8. GPU usage and bricks loading time. Comparison of the time spent (in ms) for the rendering, the caches management, and the bricks loading on a 25-second exploration sequence inside the Mandelbulb dataset (b), with the following scenario: a zoom-in, a z-axis navigation, a pan navigation and an edge detection algorithm as image processing. This timeline was obtained with the 2D visualizer described in section 8.1. The five images (a) to (e) are key snapshots of the previously described scenario.

Table 2 .

 2 Volumes description, memory occupancy and ray casting frame rates. The first part of this table lists the characteristics of the 3 datasets used.

	-

The last row of table

2

shows the aver

http://lz4.github.io/lz4/

ACKNOWLEDGMENTS

This work is supported by the French national funds (PIA2'program "Intensive Computing and Numerical Simulation" call) under contract No. P112331-3422142 (3DNeuroSecure project). The purpose of this project is to propose a collaborative solution to process and interactively visualize massive multi-scale data from ultra-high resolution 3D imaging. This secure solution also aims to break therapeutic innovation by allowing the exploitation of 3D images and complex data of large dimensions as part of applications framework linked to neurodegenerative diseases like Alzheimer's. We would like to thank all the partners of the consortium led by Neoxia, the three French clusters (Cap Digital, Systematic and Medicen), Florent Duguet for his participation and advice on the problems we encountered, Thierry Delzescaux and the Mircen team (CEA, France) for the datasets (a) & (c) Fig. 7 as well as NVidia for all their advice.