Rafic Nader Lten

Andrey Smorodin
email: 1.andrey.v.smorodin@op.edu.ua

Natalia De

La Fournière
email: ndelafourniere@artefakt-ai.com

Yves Amouriq
email: yves.amouriq@univ-nantes.fr

Florent Autrusseau
email: florent.autrusseau@univ-nantes.fr

Automatic teeth segmentation on panoramic X-rays using deep neural networks

Keywords: Panoramic X-ray images, Teeth segmentation, U-Net, location prior, deep learning

In order to build an intelligent dental care process that both facilitates the treatment and improves the diagnosis, an accurate tooth segmentation and recognition on panoramic X-ray images might prove helpful. Although many studies have been conducted on teeth segmentation, few methods allow to perform tooth recognition and numbering at the same time. The existing methods allowing both those processes rely on instance segmentation architectures. To fill some gaps in the area of dental image segmentation, we propose a novel approach of automatic joint teeth segmentation and numbering using the pioneer U-Net model. We are first to employ the conventional U-Net model and show its limitations to provide accurate segmentation, being affected by noisy pixels outside the teeth region and by missing teeth in the X-ray images. To overcome this problem and reduce the misclassifications, we use a bounding box prior at the level of the skip connections. Such an approach helps guiding the network to better locate the teeth, and hence improves the segmentation. To validate the effectiveness of the method, we have conducted two experiments on the DNS Panoramic Dataset: a first one using manual bounding boxes and another one relying on a preliminary step of object detection. The implemented networks were evaluated using the Dice coefficient index and our results showed that consideration of location information onto the skip connections improves the performances of the semantic segmentation by 5% to 10% in average Dice accuracy depending on the quality of the bounding box labels.

I. INTRODUCTION

In order to formulate a diagnosis of oral diseases, panoramic X-rays are important tools for dentists to visualize the structure, shape and position of each tooth and confirm or discard a particular diagnosis such as a fracture, infection, tooth loss or simply to detect any previous dental treatment. The limitations of the image acquisition modality; its resolution, a bad contrast balance or the presence of high amplitude noise in panoramic X-rays can make the interpretation quite challenging and lead to some errors in formulating the diagnosis. So, providing an automatic analysis of panoramic X-ray images has become crucial in the dentistry community as it can increase the diagnosis accuracy while reducing the screening time and thus, the medical costs. One of such aids to dentists is medical image segmentation, which consists in the classification of each image pixel into an object of interest.

In the field of teeth segmentation, several unsupervised pixel-wise segmentation [START_REF] Silva | Automatic segmenting teeth in x-ray images: Trends, a novel data set, benchmarking and future perspectives[END_REF] approaches have been developed, mostly using intra-oral images (periapical or bitewing radiographs). Very few works have been conducted on panoramic X-ray images, including region-based [START_REF] Lurie | Recursive hierarchic segmentation analysis of bone mineral density changes on digital panoramic images[END_REF], threshold-based [START_REF] Amer | An efficient segmentation algorithm for panoramic dental images[END_REF], cluster-based [START_REF] Alsmadi | A hybrid fuzzy c-means and neutrosophic for jaw lesions segmentation[END_REF], or boundary-based [START_REF] Razali | Sobel and canny edges segmentations for the dental age assessment[END_REF] methods. Silva et al. [START_REF] Silva | Automatic segmenting teeth in x-ray images: Trends, a novel data set, benchmarking and future perspectives[END_REF], provides an interesting overview of 10 classical methods on dental imaging segmentation. They created a dataset of 1500 panoramic X-ray images, called UFBA-UESC Dental Images Dataset, on which they tested the previously cited methods. Unfortunately, they conclude that classical segmentation solutions failed to completely isolate the teeth from neighboring bone parts present inside the mouth.

Encouraged by their recent success on various computer vision tasks, they proposed Convolutional Neural Networks (CNN) for teeth segmentation. In [START_REF] Silva | Automatic segmenting teeth in x-ray images: Trends, a novel data set, benchmarking and future perspectives[END_REF], the authors propose separating the whole dental arch from the background using deep learning. This solution surpassed most of the traditional methods considerably. Jader et al. [START_REF] Jader | Deep instance segmentation of teeth in panoramic x-ray images[END_REF], investigated teeth segmentation and detection on a modified version of the previously mentioned dataset using a Mask R-CNN [START_REF] He | Mask R-CNN[END_REF] solution to carry out an instance segmentation task. However, all teeth were classified into a single category, and thus independent tooth recognition was not considered. Koch et al. [START_REF] Koch | Accurate segmentation of dental panoramic radiographs with U-Nets[END_REF] trained a semantic segmentation network using the pioneer U-Net [START_REF] Ronneberger | U-net: Convolutional networks for biomedical image segmentation[END_REF], with a patching scheme strategy and obtained good results at detecting the foreground class (presence of teeth). However, again, this approach does not allow an independent tooth extraction / recognition. Indeed, most of the recent teeth segmentation methods that can be found in the literature [START_REF] Zhao | Tsasnet: Tooth segmentation on dental panoramic x-ray images by twostage attention segmentation network[END_REF], [START_REF] Chen | Mslpnet: multi-scale location perception network for dental panoramic x-ray image segmentation[END_REF] are devoted to oneclass segmentation, i.e. all teeth are classified into one single category, thus ignoring both the morphological properties and It is quite challenging to come up with a computer vision strategy that will be able to precisely recognize and locate each of the 32 teeth on panoramic X-rays as such radiographs commonly suffer from a very large variability, due to a rather high amplitude noise or even low contrasts. Moreover, within each category (incisors, canines, premolars and molars), the teeth may sometimes be quite difficult to differentiate for non experts. Their shapes are very close, their grey level amplitudes are in the same range, and the position information can be misleading, as a missing tooth might lead its neighbor to shift and occupy its site. Aside from the teeth arrangement, which can vary significantly from person to person, two consecutive teeth can share some restorations and can be of very close proximity if not even overlapping on the X-ray image. This leads us to question CNN's ability to differentiate two distinct teeth sharing the same morphology and same pixel values. To automatically segment, detect and number the teeth, Silva et al. [START_REF] Silva | A study on tooth segmentation and numbering using end-to-end deep neural networks[END_REF] proposed instance segmentation methods using 4 different architectures, namely, PANet, HTC, ResNeSt and Mask R-CNN, and obtained very good performances.

In this work, we propose a different approach, where by using semantic segmentation we are able to accurately segment the image and locate the teeth positions, and thus, we can efficiently number the teeth in a panoramic X-ray using a U-Net based architecture. The latter is known for its ability to improve the effect of fine-grained segmentation compared to a Fully Convolutional Network (FCN) [START_REF] Long | Fully convolutional networks for semantic segmentation[END_REF] segmentation head, used in instance segmentation.

In order to increase the segmentation accuracy and reduce the misclassification of some teeth positions, we exploit the use of spatial prior onto U-Net training. For instance, recent works have studied the integration of prior knowledge, such as the shape or location of objects, used as constraints, in order to improve CNN's performance [START_REF] Zotti | Gridnet with automatic shape prior registration for automatic mri cardiac segmentation[END_REF]- [START_REF] Jurdi | BB-UNet: U-Net with bounding box prior[END_REF]. For example, Zotti et al. [START_REF] Zotti | Gridnet with automatic shape prior registration for automatic mri cardiac segmentation[END_REF] proposed extracting the cardiac center of mass by adding a regression model onto the bottleneck of U-Net. In [START_REF] Ghafoorian | Location sensitive deep convolutional neural networks for segmentation of white matter hyperintensities[END_REF], the authors injected multi-scale patches extracted from the image into their network to incorporate spatial information. In this context, we explore the effect of injecting the tooth location information by introducing a bounding box prior at the level of the skip connections [START_REF] Jurdi | BB-UNet: U-Net with bounding box prior[END_REF].

This paper is structured as follows: in section II, we provide some details on the modified U-Net architecture, show how it considers the teeth location for an optimized segmentation task. In section III we analyze improvements of the performances brought by this spatial localization addition onto the skip connections compared to an original U-Net model and finally, section IV concludes our work.

II. METHOD

In this section, we present the pipeline of our method. The architecture of the semantic segmentation model is described as well as the used prior knowledge supervision technique.

A. U-Net Architecture

As previously mentioned, the architecture of our method is based on the U-Net model, which is widely used for medical image semantic segmentation as it can achieve very good performances while relying on small sets of training data. As inputs, our model takes some images of size 512 × 1024 and at the output, produces a pixel-wise probability prediction for 33 classes; a separate class is dedicated to each tooth position (that is 32 classes), the last one being allocated to the background.

The U-Net architecture has a symmetric encoder/decoder structure with bottleneck and skip connections. We have chosen a 4 stages U-Net architecture. Each stage in the encoder consists of two successive convolutional blocks separated by max pooling layers in order to reduce the feature resolution, while each stage of the decoder consists of a transposed convolution followed by two convolutional blocks. Each convolutional block contains a 3×3 convolutional layer, a batch normalization layer and a Rectified Linear Unit (ReLU) activation. The last layer is a 1 × 1 convolution followed by a softmax activation function that produces the pixel-wise probability map. The encoder layers extract contextual information or features present in the image. On the other hand, the decoder layers are dedicated to determine the localization of the patterns and recover the image maps with their original input size. In order to combine both contextual and positional features, the skip connections are used between the encoder and the decoder.

Despite its success in incorporating global features via its symmetric architecture and the use of skip connections, the U-Net model has some limitations when incorporating location information, as we will see below. The way to address this issue is through the integration of prior information and more specifically, location prior, into the model architecture. In this

B. Modified U-Net architecture

In this work, we adopt a very recent U-Net adaptation, the BB-Conv layer introduced by Rosana et al. [START_REF] Jurdi | BB-UNet: U-Net with bounding box prior[END_REF]. A BB-Conv layer is composed of 2D max pooling layers followed by 2 convolutional layers. Bounding boxes relative to different teeth are fed independently through a multi-channel binary map onto the BB-Conv layer that gives, as an output, a feature map helping the network to improve the localization of a given tooth. The BB-conv layers are introduced onto each skip connection as shown in Fig. 2. At each stage, an element-wise multiplication is performed between the output of the BB-conv layer and the input features from the encoder before being concatenated with the features issued by the deconvolutional layers of the decoder. To distinguish this new network architecture from the original U-Net, in the following, it will be referred as the Modified U-Net.

C. Generating bounding boxes

In order to assess the effectiveness of the proposed model, two techniques are used to generate the bounding boxes. In the first case, the bounding boxes for each tooth were generated from the ground-truth as the smallest bounding box framing the tooth's mask. As we used manual bounding boxes for both our training and inference tasks, this setup supposedly provides the best results that can be reached, in the following it will be denoted as Optimal U-Net. Evidently, this approach only intends to prove the feasibility of the method and provide a glimpse of the best segmentation results one can reach. It isn't feasible in real world applications, as of course, the manual segmentation is not available. Ideally, automatic bounding boxes should be generated using an object detection framework for the training and inference steps. Our proposed system here is based on a region proposal approach, namely, the state-of-the-art Faster R-CNN model [START_REF] Ren | Faster R-CNN: towards realtime object detection with region proposal networks[END_REF]. A Faster R-CNN is a single unified network consisting of two stages: the Region Proposal Network (RPN) which proposes regions of interest (ROIs), in which an object can be localized, and the ROI classifier that identifies which category an object belongs to and refines the bounding box generated by the RPN. As an extension of the Faster R-CNN, the Mask R-CNN introduces a parallel branch (or the mask branch) to the network in order to carry out instance segmentation. For the sake of comparison between segmentation masks, we propose directly using bounding box coordinates predicted by the classification branch of the Mask R-CNN network to be fed onto the skip connections of our Modified U-Net, rather than training a separate Faster R-CNN network.

D. Pipeline of the Proposed method

We have developed a two-steps process for teeth segmentation using a Modified U-Net architecture. In the first step, we train a Mask R-CNN model to locate the teeth, and more precisely, to extract bounding boxes of each positional tooth. In the second step, we train a Modified U-Net network using the predicted bounding boxes from the first step.

III. EXPERIMENTAL RESULTS

Let us now give some details on the dataset composition, provide some information on the training step and show the performances our approach can achieve.

A. Dataset

The DNS Panoramic images from IvisionLAB [START_REF] Silva | A study on tooth segmentation and numbering using end-to-end deep neural networks[END_REF] were used. This dataset contains 543 images of size 1127 × 1191 annotated with the tooth numbering (FDI notation, i.e. position label) using the Coco format. The images can be divided into 8 categories regarding the presence or absence of all teeth in images, of restorations and appliances (please refer to [START_REF] Silva | A study on tooth segmentation and numbering using end-to-end deep neural networks[END_REF] for a detailed description of different categories). Mask R-CNN and U-Net models were trained using a 4-fold cross validation. In total, 111 images were retained to build the test set and the rest of the images were divided into 4 folds (108 images each) thus composing the train and validation data in a cross-validation fashion.

B. Object detection training procedure

Our experiments are based on the open source implementation of Mask R-CNN [START_REF] Abdulla | Mask R-CNN for object detection and instance segmentation on keras and tensorflow[END_REF]. We also benefit from transfer learning using pre-trained weights taken from the MS COCO dataset. The network's backbone was ResNet-101. We used the stochastic gradient descent optimizer with a learning rate of 0.001 and momentum of 0.9. For each combination of train/validation data, we trained a Mask R-CNN network using a batch size of 1 for 80 epochs (the validation loss didn't decrease after 80 epochs for any of the 4 networks). For training and inference, we set a detection score threshold at 0.05 to filter out all predictions with low confidence. Additional details of the several hyperparameters and parameters related to the implementation can be found in [START_REF] Abdulla | Mask R-CNN for object detection and instance segmentation on keras and tensorflow[END_REF]. The performances of the object detection task were evaluated using the mean Average Precision (mAP), i.e. the average of the Average Precisions calculated for all tooth classes with the traditional IoU (Intersection over Union) of 0.5 [START_REF] Everingham | The PASCAL Visual Object Classes Challenge[END_REF]. The mAPs of 4 trained neural networks (of 4-CV) on their respective validation data are 0.983, 0.963, 0.980 and 0.970 respectively. These results are comparable to the performances reached by previous works on teeth detection and numbering [START_REF] Silva | A study on tooth segmentation and numbering using end-to-end deep neural networks[END_REF], [START_REF] Tuzoff | Tooth detection and numbering in panoramic radiographs using convolutional neural networks[END_REF], [START_REF] Mahdi | Optimization technique combined with deep learning method for teeth recognition in dental panoramic radiographs[END_REF] and show that the object detection network performs very well, and most importantly, is consistent across the validation dataset.

After obtaining all of the bounding boxes and before training our Modified U-Net, for every image in our dataset, we choose one box to be selected per class at most, i.e., the box presenting the highest confidence score, based on the natural assumption that each tooth can occur in the image only once.

C. Loss function and model training

We use Dice loss as the criterion to optimize the U-Net model parameters. This loss is widely used in medical image segmentation and is defined as :

L dice = - 1 C C c=1 2 N i=1 p c (i).g c (i) N i=1 p c (i) 2 + N i=1 g c (i) 2
where N is the number of pixels in the image, C represents the number of class labels, p c (i) is the predicted probability of class c at pixel i and g c (i) ∈ {0, 1} represent the groundtruth label at position i. Note that C = 33 as we take into account the background class in the computation of the Dice loss; The idea is to improve the ability of the model to detect the boundaries of each tooth.

In order to study the effectiveness of the proposed method, we conducted three different experiments using three different model settings : a classical U-Net, a Modified U-Net with bounding boxes generated by Mask R-CNN and an Optimal U-Net described earlier. When building our dataset, we resized the images to 512 × 1024 pixels. For all three U-Net configurations, we used 64 feature maps in the highest stage. The number of feature maps double with each downsampling step until it reaches 512 within the fourth stage and within the bottleneck. Moreover, the Adam optimizer is used to train our model with an initial learning rate of 0.0001. We set the number of training epochs to be 60 with a batch size of 2. The learning rate was halved each 15 epochs if the validation performance did not improve.

D. Quantitative Results

Fig. 3 shows the overall results of the proposed method in terms of average Dice coefficient index (%) and its standard deviation for each tooth position from the cross validation step.

We observe that the Modified U-Net configuration consistently outperforms the original U-Net model for all teeth classes. Evidently, the Optimal U-Net model offers even better performances. The Optimal U-Net surpasses the other two models with an average dice index of 94.5% and with minimal standard deviations. This confirms the importance of properly integrating location knowledge onto the model structure. The performance gap between the Optimal and the Modified U-Net is due to the inaccuracies between the ground-truth bounding boxes and the ones detected by the Mask R-CNN network.

The quite large Dice score difference between the original U-Net and the Modified U-Net comes from the detection inaccuracies of the molar teeth. The inclusion of prior information via the bounding boxes onto the skip connections increases the dice accuracy by over 11% on tooth #46, 10% on teeth #28 and #36, 9% on tooth #16 and by 8% on teeth #18 and #38. This stems from the fact that the original U-Net model suffers from a misclassification issue on panoramic X-ray images with noisy pixels, degraded teeth and most importantly, missing teeth. The positional information carried by the bounding boxes guarantees that all teeth are located nearby their correct location. Finally, the Modified U-Net and U-Net architectures achieve an average Dice accuracy of 89.5% and 85% respectively.

Following the same protocol, we computed the Dice metric using the test data set. Using unseen data, the Optimal U-Net offered the best performances with an average dice score of 94.49% followed by the Modified U-Net with an average Dice of 90%. Unsurprisingly, the U-Net model exhibits the worst Dice coefficient with an average of 86%. Fig. 4 summarizes the average and the boundaries (error bars) of the Dice ratio for each architecture with respect to the tooth type. As it can be observed from this Figure, the Original U-Net performances significantly decrease from incisors (0.89) to molars (0.82). This is due to teeth being segmented wrongly as one of their (missing) neighbors, which mostly occurs for premolars and molars. It is worth noting that 2.7% of incisors, 2.9% of canines, 3.8% of premolars and 7.9% of molars were missing in our test set (composed of 111 X-rays). The Modified U-Net presents a much smaller decrease in Dice ratio (0.91 for incisors, down to 0.88 for molars) explained by the rare miss-labeling of teeth by the object detection model. This can be clearly observed from Table I and molar teeth. However, the Optimal U-Net maintains consistent performances (0.94) across the different tooth types. These results confirm, once again, the importance of integrating some geographic location information into the segmentation task.

E. Qualitative results

In Fig. 6, we show some typical segmentation results produced by both the U-Net and our proposed Modified U-Net. Compared to U-Net, the Modified U-Net avoids segmenting surrounding tissues as teeth (the first row), most importantly, it avoids classifying pixels from one class (one tooth) into 2 or more classes (the second row). This mainly occurs when a tooth is degraded or missing. The third row shows an example of a false segmentation also occurring with the Modified U-Net due to a mislabeling of tooth position by the object detection network (tooth number 37 being labeled as number 38).

F. Comparison with instance segmentation

U-Net (or some of its derived forms) is among the gold standards for biomedical image segmentation, as it presents a high robustness despite a noisy image context. However, several works propose instance segmentation using the Mask R-CNN architecture on medical images [START_REF] Anantharaman | Utilizing mask r-cnn for detection and segmentation of oral diseases[END_REF]- [START_REF] Dogan | A two-phase approach using mask r-cnn and 3d u-net for highaccuracy automatic segmentation of pancreas in ct imaging[END_REF] as the latter offers interesting performances for some problems that require a sense of a whole object segmentation rather than focusing on individual pixels only. In fact, R-CNN (Region Based Convolutional Neural Networks) were initially developed for object detection purposes. The Mask R-CNN extends the Faster R-CNN by adding a branch, the purpose of which, is to predict the object mask, hence providing segmentation capabilities. However, at the pixel-level, the segmentation task of Mask R-CNN is performed by FCN, which is less accurate than U-Net especially for small datasets. To address this issue, the Modified U-Net combines the performance of U-Net with the advantage of an object detection step. For a thorough evaluation of our method and an exhaustive comparison with state-of-the-art approaches, we hereby compare our results with the segmentation masks obtained by the mask branch of Mask R-CNN. Table II shows that the Modified U-Net provides higher performance in terms of dice coefficient (2 percentage point) than Mask R-CNN for all teeth types. Finally, it is important to note that the Modified U-Net can be used easily with faster object detection approaches than region proposed methods such as the one stage detectors [START_REF] Liu | Ssd: Single shot multibox detector[END_REF]- [START_REF] Redmon | You only look once: Unified, real-time object detection[END_REF]. This method has the potential for practical and faster application for automated dental analysis, but further investigation is needed to assess its precision. On the other hand, one stage detectors are not quite suitable for instance segmentation as the mask branch incorporation is not straightforward.

IV. CONCLUSION

The vast majority of the research works on teeth detection and segmentation focuses on periapical or intraoral radiographs. Few works are devoted to the panoramic X-rays, as such radiographs unfortunately often suffer from a poor quality (high noise, unbalanced contrasts, etc.). The segmentation of panoramic X-rays is very challenging. In this work, we propose an automatic teeth segmentation and numbering using the Modified U-Net model that integrates location prior by connecting bounding boxes to the skip connections. We have demonstrated effectively that it is crucial to incorporate some geometric localization information into the U-Net, in order to improve its performances and avoid false segmentation occurring on panoramic images with noisy surrounding pixels, degraded or missing teeth. Experiments showed an increase of the Dice coefficient index for all teeth classes with an average increase of 5% and which can go up to 10% for an Optimal U-Net. Furthermore, the proposed segmentation scheme is superior to the Mask R-CNN segmentation in terms of the Dice coefficient. As the Modified U-Net depends on bounding boxes generated by an object detection approach, the future steps include using or developing an optimized object detection algorithm that minimizes the shift between the detected vs. ground-truth bounding boxes.

Fig. 1 .

 1 Fig. 1. FDI tooth notation

Fig. 2 .

 2 Fig. 2. Network architecture.

Fig. 3 .

 3 Fig. 3. Dice coefficients (a) and (c) and Standard Deviations (b) and (d) for U-Net, Modified U-Net, and Optimal U-Net configurations ((a)&(b): teeth from the upper jaw, (c)&(d): lower jaw).

Fig. 4 .

 4 Fig. 4. DICE coefficients for the three tested CNN architectures and for the four distinct teeth types. Results obtained from the test set (composed of 111 images).

Fig. 5 . 13 Fig. 5 .

 5135 Fig. 5. Sub-matrix of the confusion matrix for a 0.5 IoU detection threshold on the test dataset, corresponding to molar teeth.

Fig. 6 .

 6 Fig. 6. Examples of the segmentation results using the U-Net (c) and Modified U-Net (d) architectures.

ACKNOWLEDGMENT

This work was partially supported by the I-Site NExT fund (call for projects "Innovez"), the "Région des Pays de la Loire", as well as the "Fonds Européen de Développement Régional" (FEDER).


Oleksandr Ye Horniichuk

Karl Ridier

Gábor Molnár

Volodymyr O Kotsyubynsky

Sergiu Shova

Vladimir M Amirkhanov

L Il'ya Gural'skiy

Azzedine Salmon

Bousseksou

A Il'ya

Gural'skiy

Lionel Salmon

Azzedine Bousseksou

Solvatomorphism, polymorphism and spin

come

Introduction

Switchable molecular materials showing spin crossover (SCO) properties 1 attract attention of many research groups as they can be used for fabrication of functional devices, finding their application in information storage, molecular electronics, photonics, sensing and construction of mechanical actuators. 2- 6 In this context, the ability to fabricate high quality, pure and homogeneous thin films of SCO materials is crucial. To this aim, vacuum thermal deposition provides a powerful means. 7,8 However, there are only a handful of spin crossover complexes that could be deposited by this technique. [9][10][11][12][13][14][15][16][17][18][19][20][21][START_REF] Mahfoud | [END_REF][23] Among these, only a few sublimable SCO complexes undergo spin transition above room temperature, namely [Fe(pypyr(CF 3 ) 2 ) 2 (phen)] (pypyr = 2-(2′-pyridyl)pyrrolide, phen = 1,10-phenanthroline), 15 [Fe(HB(pz) 3 ) 2 ] (pz = pyrazolyl-1-yl) [START_REF] Mahfoud | [END_REF]23 and [Fe(HB(1,2,4-tz) 3 ) 2 ] (tz = triazol-1-yl), 21 the two latter belonging to the tris(1azolyl)borate family of compounds. [START_REF] Trofimenko | Scorpionates: The Coordination Chemistry of Polypyrazolylborate Ligands[END_REF] For Fe-pypyr thin films, the transition temperature T 1/2 is reduced from 390 K for the bulk compound to approximately 330 K, whereas the two other complexes display similar SCO properties in bulk and thin film forms (T 1/2  393 K and 335 K for the pyrazolyl and triazolyl derivatives, respectively). Notably, the high-quality thin films of [Fe(HB(1,2,4-triazol-1-yl) 3 ) 2 ] that were vacuum-deposited in the 20-200 nm thickness range, showed an abrupt SCO without any fatigue over more than 10 million thermal switching cycles, which could be exploited in photonic, electronic, sensor, actuator and temperature damping devices. 21,[START_REF] Manrique-Juarez | [END_REF][26][27][28][29] The synthesis of [Fe(HB(1,2,3-triazol-1-yl) 3 ) 2 ] 1 (Chart 1), a structural isomer of its 1,2,4-triazolyl counterpart, was described by Jenkins et al., 30 but there is no report on the SCO properties of this compound. In order to find out how the change of exodentate nitrogen atom position from 4 to 3 in the triazolyl ring impacts crystalline structure and SCO properties of the complex, as well as to analyze its solvato/polymorphism, we report herein the single crystal structure of 1 (both in solvated and desolvated forms) as well as its physical properties investigated by Mössbauer spectroscopy and variabletemperature magnetic susceptibility, Raman spectroscopy powder X-ray diffraction, calorimetry and thermogravimetry methods. Review on polymorphism in spin crossover complexes has been reported by Tao et al. 31 Such polymorphs are generally obtained serendipitously as mixture of crystals as illustrated by works of Ruben and collaborators. 32 On the other hand, numerous examples of solvatomorphs in SCO complexes have been also reported mainly in polymeric Fe(II) Hofmann type and in molecular Fe(II) and Fe(III) compounds. [33][34][35][36][37][38][39][40][41][42][43][44][45][46][47] These solvatomorphs result from a mixture of crystals obtained serendipitously, by a vapor-induced single-crystal-to-singlecrystal transformation or growing crystals by liquid-liquid diffusion methods. In the present case, solvatomorphs were Please do not adjust margins Please do not adjust margins obtained by changing the solvent during the recrystallization procedure.

Methods

Infrared spectrum was recorded on a Perkin Elmer 1725 spectrometer. Elemental analysis (for C, H and N) was performed using a Perkin Elmer 2400 series II Instrument. Single-crystal X-ray diffraction analysis was performed on an Oxford-Diffraction XCALIBUR E CCD diffractometer with graphite-monochromated Mo-Kα radiation. The unit cell determination and data integration were carried out using the CrysAlisPro package from Oxford Diffraction. Multi-scan correction for absorption was applied. All the structures were solved using an intrinsic phasing method using SHELXT and refined by full-matrix least-squares on F2 with SHELXL 48,49 using the graphical interface of Olex2. 50 Non-hydrogen atoms were refined anisotropically. Aromatic hydrogen atoms were geometrically fixed and refined using a riding model. Variable temperature powder X-ray diffraction patterns were recorded using a Panalytical X'Pert diffractometer equipped with a Cu Xray tube, a Ge(111) incident beam monochromator (λ = 1.5406 Å) and an X'Celerator detector. Variable-temperature Raman spectra were acquired by means of an Xplora (Horiba) Raman microspectrometer, using the 532 nm line of a Nd-YAG laser (0.1 mW). The sample temperature was controlled using a Linkam THMS600 liquid nitrogen cryostat. Differential thermal analysis and thermogravimetric (DTA-TG) data were acquired simultaneously using a Perkin-Elmer Diamond thermal analyzer under a 100 ml min -1 nitrogen flow at 10 K min -1 . DSC analysis was carried out on a Netsch DSC 3500 Sirius instrument under nitrogen purging gas (20 ml min -1 ) at a heating/cooling rate of 10 K min -1 . Variable-temperature magnetic susceptibility data were obtained at cooling and heating rates of 2 K min -1 under a field of 1 kOe using a Quantum Design MPMS magnetometer. 57 Fe-Mössbauer spectra were recorded using a constant acceleration Mössbauer spectrometer with a 57 Co source. Fitting of the experimental data was performed with Recoil software. 51

Materials and synthesis

Commercial reactants and solvents were used directly without further purification. The tris(1,2,3-triazolyl)borate ligand (NaL) and its neutral [Fe(HB(1,2,3-tz) 3 ) 2 ] iron (II) complex (hereafter FeL 2 or 1) were synthesized according to experimental procedures reported in the literature 30 . The coordination compound isolated from the mother aqueous solution appears as dodecahydrate. The anhydrous powder form can be obtained via heating the hydrate above 90℃. FeL 2

Results and discussion

X-ray crystallography

Complex 1 is soluble in polar solvents, such as water, methanol, dimethylsulfoxide (DMSO), dimethylformamide (DMF) and poorly soluble in much less polar chloroform. The complex tends to form solvates, which can be easily crystallized from the above-mentioned solvents. The number of the solvent molecules per one molecule of the complex depends on the nature of the solvent. For instance, in the case of DMSO, the formation of two solvates 1 3 space group. DMSO molecules in the latter crystals are disordered over six and over two positions, respectively. Upon long-term exposure to ambient air (~hours to days), crystals of solvatomorphs tend to lose their solvents. This process is followed by disintegration of the crystals along with structural changes and can be promoted by heating. Remarkably, crystals of 1•2CHCl 3 and 1•2MeOH undergo single-crystal-to-singlecrystal transformations upon heating. When heated at 313 K for one day, solvatomorph 1•2CHCl 3 completely loses its solvent to form stable crystals of 1. The same transformation occurs in case of its methanol counterpart, but the latter requires heating at 323 K for two days. Incomplete solvent loss can be observed in 1•7DMSO: after two days of evaporation, the single crystal degrades to powder, the diffraction pattern of which is consistent with the simulated one for 1•2DMSO (Fig. S7). The desolvated compound 1 crystallizes in the monoclinic P2 1 /c space group with 2 molecules of FeL 2 per unit cell (Fig. 1). It should be noted that the crystals used for refining the structure were obtained via a solid-state crystal-to-crystal transformation. This huge structural rearrangement is related to a chloroform solvent loss, namely c.a. 25% of the volumic contraction occurs (1412 to 1063 A 3 ) that leads to low quality multi-twin crystals. Three domains were considered upon the refinement. The room temperature (low spin state) unit cell Please do not adjust margins Please do not adjust margins volume per molecule in the crystals of 1 (531.8 Å 3 ) is greater when compared to its 1,2,4-triazolyl isomer (516.6 Å 3 ). 52 While the structure of [Fe(HB(1,2,3-tz) 3 ) 2 ] is less dense, the mean intermolecular Fe•••Fe separation (8.49 Å) is slightly shorter in 1 than in the 1,2,4-triazolyl isomer (8.87 Å). The coordination spheres of iron(II) ions in the solvatomorphs (as well as in 1) have slightly distorted [FeN 6 ] octahedral geometries (see Table S2). The Fe-N distances lie in the range from 1.952(3) to 1.972(2) Å, while the mean value of coordination octahedra volume (V p ) is 10.05 Å 3 , as expected for the low spin iron(II). [START_REF] Guionneau | Spin Crossover in Transition Metal Compounds II[END_REF] Considering all the solvatomorphs, the average value of the octahedral distortion parameter equals ∑|90°-𝜃 𝑁 -𝐹𝑒 -𝑁 | 21.66˚. Unfortunately, the high spin structure of 1 could not be assessed because crystals crash at the corresponding high temperatures. Indeed, single crystals of 1, obtained via crystalto-crystal transformation through thermally induced solvent loss, appear as low quality, highly cracked at room temperature, while further heating is a trigger for crystal degradation. In the crystals of most of the solvatomorphs, as well as in 1, π-π stacking interactions play a key role in the packing of FeL 2 molecules. The interactions occur between aromatic 1,2,3triazolyl rings of adjacent complex molecules (Fig. 1, S1-S5). It is noteworthy, that these interactions are spread along one direction in 1•2DMSO, however they are bidirectional in 1 (Fig. S1) and in the rest of the solvatomorphs where they are observed. The centroid-centroid distances range from 3.560(3) to 3.6959( 17) Å (Table 1). Another type of interactions that hold molecules of the complex together in the crystal lattices of 1 and 1•Solv is represented by rather weak C-H•••N interactions. The exodentate nitrogen atoms and C-H moieties of 1,2,3triazolyl rings of neighboring FeL 2 molecules are involved in such intermolecular interactions (Fig. 2, S2-S5). Likely to π-π stacking interactions, C-H•••N contacts propagate along a and b crystallographic directions (Fig. S2), however the structure of 1•2CHCl 3 reveals only a single propagation direction (Fig. S4). The parameters of the respective C-H•••N contacts are gathered in Table 1. In 1•2DMSO and 1•2DMF•H 2 O, the oxygen-containing polar solvents are bound to molecules of 1 via C-H•••O contacts (Fig. 2,S5a). In 1•7DMSO, the molecules of FeL 2 do not display any direct contacts with each other, but interact through solvent molecules by means of C-H•••O contacts (Fig. S6). In 1•2MeOH, a hydrogen bonding occurs involving the proton of the OH moiety of the solvent and the exodentate nitrogen atom of the triazolyl ring of the complex (d(O•••N) = 2.867(4) Å, Fig. S3a). Chloroform molecules in 1•2CHCl 3 reveal C-H•••N and C-H•••Cl interactions with FeL 2 molecules (Fig. S4). Compound 1 undergoes a virtually complete, gradual spin crossover between the low spin (S = 0) and the high spin (S = 2) states, with a transition temperature T 1/2 = 373 K. We must note that the SCO is not fully completed at 400 K, which is the highest attainable temperature by our magnetometer, but variabletemperature Raman spectroscopic measurements conducted in the 303-433 K range reveal the flattening of the SCO curve slightly above 400 K (see Fig S9b). On the other hand, room temperature Mössbauer spectra recorded both for 1•12H 2 O and for the anhydrous powder 1 evidence only a low spin ferrous species (see Fig. S10). The variable temperature magnetic data could be conveniently fitted by the conventional Slichter-Drickamer (SD) model 37 (see Fig. 3), providing estimation of the enthalpy change (H = 14 kJ mol -1 ), entropy change (S = 38 J K -1 mol -1 ) and cooperativity ( = 5387 ± 52 J mol -1 ) associated with the SCO phenomenon. One can remark here the relatively low entropy change and high cooperativity of the spin crossover in 1, similar to the values reported previously for the 1,2,4-triazolyl isomer (S = 48 J K -1 mol -1 ,  = 5700 ± 50 J/mol). [START_REF] Slichter | [END_REF] At first glance, one may ascribe the different magnetic behavior on the first heating to solvent loss. Indeed, TG-DTA data acquired for compound 1•12H 2 O reveals dehydration in the temperature range between ca. 313 and 363 K (Fig. 4). However, the DSC curve obtained for the hydrated powder displays several intense and broad endothermic peaks during the first heating (Fig. 5), which cannot be solely ascribed to the water loss, but denote also other (overlapped) physical phenomena, including obviously the SCO, but also various structural changes (vide infra). Remarkably, the DSC trace on cooling and on successive thermal cycling display a single exothermic (endothermic) peak at 366 K (373 K) on cooling (heating). The peak temperature as well as the associated enthalpy change of ca 14 kJ mol -1 agrees with the SD fit of the magnetic susceptibility data, allowing us to assign this peak to the sole SCO phenomenon.

In order to further investigate the possible interplay between solvent loss, structural changes and SCO in 1, variable temperature powder XRD measurements were carried out in the 333-493 K temperature range following heating of 1•12H 2 O at 328 K for 1 h. The powder patterns recorded at 333, 363 and 393 K are substantially different -denoting successive dehydration/structural changes (Fig. 6). Of particular interest is that the experimental powder pattern recorded at 363 K shows the best match with the diffraction pattern simulated using the single crystal structural data of 1 (see Fig. S8). Then, further heating up to 493 K and cooling back to 303 K does not reveal any substantial change in the diffractograms (apart small shifts and broadenings of the peak positions -depending on the temperature). This finding is in line with the calorimetric and Please do not adjust margins Please do not adjust margins magnetic measurements, pointing to the occurrence of a very gradual and reversible SCO following a first heating of the sample up to ca. 400 K.

Conclusions

In the present work, we report the crystal structures of solvated and desolvated forms of bis[hydrotris(1,2,3-triazol-1yl)borate]iron(II) as well as their spin crossover properties. Similar to other unsubstituted tris(1-azolyl)borate iron(II) complexes, the complex exhibits reversible spin crossover above room temperature. The SCO in this complex is gradual and is centred at 373 K. Interestingly, the first heating of the hydrated complex is accompanied by a series of structural phase transitions resulting in different polymorphic phases of the compound. Crystal structures of the anhydrous complex and most of the solvatomorphs share the same types of interactions between the molecules of complex, including π-π stackings and weak C-H•••N hydrogen bond interactions. In the course of this work, we have also evidenced the possibility of sublimating complex 1, which opens attractive perspectives for the deposition of thin films and will be the matter of our further research.

  Fig. S1-S4 crystal structures of some solvatomorphs with intermolecular interactions displayed, Fig. S5 packing diagram of 1•7DMSO, Fig. S6 X-Ray diffraction patterns for DMSO solvatomorph, Fig. S7 TGA and DTA curves for 1•12H 2 O, Fig. S8 Simulated and experimental PXRD patterns of 1, Fig. S9 Variable temperature Raman spectroscopy data for 1, Fig. S10 Mössbauer spectra of 1•12H 2 O and polymorphs of 1. CCDC 2161993-2161998. For ESI and crystallographic data in CIF or other electronic format see DOI: 10.1039/x0xx00000x

  1•2DMSO and 1•7DMSO was observed serendipitously from attempts of recrystallization of 1 in the same condition. The methanol-or chloroform-based solvates appear only as 1•2MeOH and 1•2CHCl 3 , respectively. The dimethylformamide solvate tends to crystallize with inclusion of one water molecule to form 1•2DMF•H 2 O. The solvatomorphs 1•2MeOH and 1•2CHCl 3 crystallize in monoclinic space groups P2 1 /n and P2 1 /c, respectively. 1•2DMF•H 2 O crystallizes in the orthorhombic Pbcn space group. In the case of DMSO-based solvatomorphs, 1•2DMSO crystallizes in the triclinic P space group, while 1•7DMSO adopts a trigonal R c

Chart 1

 1 Molecular structure of [Fe(HB(1,2,3-triazlol-1-yl) 3 ) 2 ] (1).

Fig. 1 .

 1 Fig. 1. Crystal structure of 1•2DMSO showing intermolecular π-π stacking interactions developing along the a axis. Hydrogen atoms are omitted for clarity. Displacement elipsoids are drawn at 50% probability level.

Fig. 2

 2 Fig. 2 Crystal structure of 1•2DMSO showing intermolecular C-H•••N and C-H•••O contacts. Parts of the complex molecules are omitted for clarity. Displacement elipsoids are drawn at 50% probability level. Hydrogen atoms in the solvent molecules are omitted for clarity

Fig. 3 χ

 3 Fig. 3 χ M T vs T plot, recorded over two thermal cycles starting with the powder sample 1•12H 2 O. Slichter-Drickamer fit was performed for the second heating curve.

Fig. 4 Fig. 5

 45 Fig. 4 TGA and DTA curves for 1•12H 2 O.

•12H 2 O (1•12H 2 O)

  . C 12 H 38 B 2 FeN 18 O 12 (M = 704.07). Anal. calc.: C, 20.47; H, 5.44; N, 35.81%. Found: C, 20.55; H, 5.42; N, 35.55%. FeL 2 (1). C 12 H 14 B 2 N 18 Fe (M = 487.88) Anal. calc.: C, 29.55; H, 2.89; N, 51.68%. Found: C, 29.51; H, 2.86; N, 51.57%.

Table 1 .

 1 Intermolecular FeL 2 -FeL 2 interaction parameters.

| J. Name., 2012, 00, 1-3

J. Name., 2013, 00, 1-3 | 3

J. Name., 2013, 00, 1-3 | 5

Acknowledgements

Authors acknowledge the financial support from the Ministry of Education and Science of Ukraine (grant 22BF037-09) and the courage of the Armed Forces of Ukraine that made the submission of this manuscript possible. We acknowledge Laure Vendier and Jean-Francois Meunier (LCC) for performing variable temperature PXRD and magnetic measurements, respectively. OYH acknowledges the French Embassy in Ukraine and Campus France for a BGF scholarship.

Conflicts of interest

There are no conflicts to declare.


David Schmoranzer

James Butterworth

Sébastien Triqueneaux

Eddy Collin

Andrew Fefferman

Design Evaluation of Serial and Parallel sub-mK Continuous Nuclear Demagnetization Refrigerators

Keywords: ultra low temperatures, adiabatic nuclear demagnetization, continuous refrigeration techniques, PrNi

We present the evaluation of two different design configurations of a two-stage PrNi 5 continuous nuclear demagnetization refrigerator. Serial and parallel configurations of the two stages are considered, with emphasis on the attainable cooling power at sub-mK temperatures and the impact of the design choices on the operation of the refrigerator. Numerical simulations of heat transfer in the setup are used to evaluate the performance of the refrigerator as well as the technological requirements for the essential thermal links. In accord with similar findings for adiabatic demagnetization refrigerators [Shirron, Cryogenics 62, 2014], our simulations show that the performance of both configurations improves as the thermal links improve, and that the parallel configuration yields a higher cooling power than the series design for a given thermal link resistance and sample temperature.

Introduction

Recently, significant effort was devoted to the development of a novel ultra-low temperature refrigerator [START_REF] Toda | Design and expected performance of a compact and continuous nuclear demagnetization refrigerator for sub-mK applications[END_REF][START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF] which would be able to overcome the limitations of traditional nuclear demagnetization setups [START_REF] Lounasmaa | Experimental Principles and Methods below 1 K[END_REF][START_REF] Pobell | Matter and Methods at Low Temperatures[END_REF] and continually maintain sub-mK temperatures. Such a continuous nuclear demagnetization refrigerator (CNDR) would find use in ultra-low temperature research of numerous physical phenomena such as superfluidity of 3 He [START_REF] Mkinen | Half-quantum vortices and walls bounded by strings in the polar-distorted phases of topological superfluid 3He[END_REF][START_REF] Mkinen | Mutual friction in superfluid 3He-B in the lowtemperature regime[END_REF][START_REF] Tsepelin | Visualization of quantum turbulence in superfluid 3He-B: Combined numerical and experimental study of Andreev reflection[END_REF], quantum states of macroscopic objects [START_REF] Wollman | Quantum squeezing of motion in a mechanical resonator[END_REF][START_REF] Ockeloen-Korppi | Stabilized entanglement of massive mechanical oscillators[END_REF], or dissipation in amorphous matter [START_REF] Liu | Hydrogen-free amorphous silicon with no tunneling states[END_REF][START_REF] Rogge | Nonlinear dielectric response of glasses at low temperature[END_REF][START_REF] Fefferman | Acoustic properties of amorphous silica between 1 and 500 mK[END_REF][START_REF] Schmoranzer | Observations on Thermal Coupling of Silicon Oscillators in Cryogen-Free Dilution Refrigerators[END_REF] or modern materials such as graphene [START_REF] Liu | Shear modulus of monolayer graphene prepared by chemical vapor deposition[END_REF]. The main factor limiting the use of traditional nuclear demagnetization setups is the fact that at sub-mK temperatures, the duration of the experimental window is often less than the thermal relaxation times of the materials of interest [START_REF] Pobell | Matter and Methods at Low Temperatures[END_REF] (spin-lattice relaxation in non-conducting materials, relaxation processes in amorphous materials, heat exchange between metal walls and liquid helium), which precludes proper sample thermalization. This issue becomes crucial if the demagnetization setup is mounted on a dry dilution refrigerator, as these typically have higher intrinsic heat leaks due to vibrations from the pulse tube assembly [START_REF] Riabzev | Vibration generation in a pulse tube refrigerator[END_REF], but even on "wet" systems, the thermal cycling might be disruptive for precise measurements. Additionally, the ability to maintain sub-mK temperature indefinitely may provide an important stepping stone for the development of further cooling techniques and eventually lead to the experimental realization of novel states of matter, such as the dual Bose-/Fermi-type superfluidity in 3 He- 4 He mixtures [START_REF] Tuoriniemi | Physics at its coolest[END_REF].

Previously, similar designs of a continuous (electronic) adiabatic demagnetization refrigerator (CADR) were discussed in detail in the literature [START_REF] Shirron | A Multi-Stage Continuous-Duty Adiabatic Demagnetization Refrigerator[END_REF][START_REF] Shirron | Development of a cryogen-free continuous ADR for the constellation-X mission[END_REF], and a review may be found in Ref. [START_REF] Shirron | Development of a cryogen-free continuous ADR for the constellation-X mission[END_REF]. The design presented below builds on the above-mentioned work and extends the temperature range of the refrigerator to sub-mK temperatures by replacing paramagnetic salt pills with a suitable nuclear magnetic refrigerant, in our case PrNi 5 . The fundamental principle of operation of the CNDR is the same as for CADR, except that it is the entropy of nuclear spins rather than of outer-shell electronic spins that is manipulated using the external magnetic field.

The construction of the CNDR requires the use of at least two nuclear demagnetization stages with separate superconducting solenoids providing the magnetic field. The stages have to be connected to each other and/or to the sample space via thermal links consisting of metal wires/blocks and superconducting heat switches. Gas-gap heat switches are unavailable in the relavant temperature range (below 10 mK) as they become inefficient at temperatures below ≈250 mK [START_REF] Shirron | Development of a cryogen-free continuous ADR for the constellation-X mission[END_REF]. In previous work on CNDR [START_REF] Toda | Design and expected performance of a compact and continuous nuclear demagnetization refrigerator for sub-mK applications[END_REF][START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF] it was shown that thermal links are indeed the crucial elements of the entire setup, requiring the use of highest purity metals prepared using specialized treatments and contacting techniques. Specifically, it was demonstrated by means of numerical simulations [START_REF] Toda | Design and expected performance of a compact and continuous nuclear demagnetization refrigerator for sub-mK applications[END_REF][START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF] that the equivalent electrical resistance (from the Wiedemann-Franz law) of the thermal link between the two demagnetization stages in a series configuration should be comparable to 150 nΩ or less (corresponding to a thermal conductance of 1.6 × 10 -4 WK -1 at 1 mK) to attain a cooling power of 20 nW at 1 mK, see also Fig. 6 further below.

Practically, to obtain such a low resistance, one would need to use 5N or 6N purity aluminium blocks (RRR over 5000 was obtained in CNRS Grenoble on 6N aluminium), 5N or better copper wires annealed under oxygen (RRR over 10000 was obtained) or silver wires of equal purity as well as specialized contacting techniques (stable Al to Cu contacts were made with contact resistances of 30 nΩ or less on the area of several cm 2 ). Another crucial issue is the contact to PrNi 5 and the purity of the nuclear refrigerant (RRR values near 30 have been reported [START_REF] Folle | Nuclear refrigeration properties of PrNi5[END_REF]). The technical developments and the details of processes used in the construction of the CNDR setup represent a separate topic and will be published elsewhere. The necessary material properties are mentioned solely to justify the values of the model parameters and link them to the ongoing development of the refrigerator.

In this manuscript, we use numerical simulations to compare two fundamentally different design approaches to the construction of the CNDR. Namely, the already discussed series configuration [START_REF] Toda | Design and expected performance of a compact and continuous nuclear demagnetization refrigerator for sub-mK applications[END_REF][START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF] is contrasted to the parallel configuration of the two demagnetization stages. Design considerations and the thermal models for both configurations are presented and, subsequently, numerical simulations similar to the calculations in Ref. [START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF] are used to provide quantitative comparison of the performance of the CNDR in these two configurations, with special attention devoted to the quality of the thermal links used.

Comparing CNDR and CADR

As the CADR setup is the starting point for the development of CNDR, it might be useful at this point to highlight some of the practical differences between CADR and CNDR setups due to the different operating range of temperatures. Besides the choice of the nuclear refrigerant and heat switches mentioned above, additional differences arise due to the significantly suppressed thermal conductivities of all materials at mK or sub-mK temperatures. Electronic thermal conductivity in metals is directly proportional to temperature, and hence will be reduced by three orders of magnitude between 1 K and 1 mK. Furthermore, the thermal conductivity of the heat switches in the superconducting state is suppressed by an even greater factor upon cooling, with a switching ratio as great as 10 6 observed at 50 mK [START_REF] Pobell | Matter and Methods at Low Temperatures[END_REF].

This leads to several important differences on the practical level. While serial CADR setups are often limited by the high heat leak through the final heat switch in the superconducting state (off-state), for CNDR, the opposite is true. While the off-state conductivities are negligibly low, the heat switches, along with their leads and contacts represent a significant thermal resistance even in the normal state (on-state), limiting the heat transfer rate through any such element. As a result, for the serial CNDR, the heat transfer rate between the two demagnetization stages may become comparable to the external heat leak to the sample space (both typically in the nW range), which further affects the overall performance of the refrigerator and to some extent complicates its analysis.

CNDR Design and Operation

Both principal design configurations have already been discussed for CADR in Ref. [START_REF] Shirron | Development of a cryogen-free continuous ADR for the constellation-X mission[END_REF] and are illustrated for CNDR in Fig. 1. In the series configuration (S-CNDR), the first of the nuclear demagnetization stages (NDS-1) is connected to the mixing chamber plate (MXC) of a dilution refrigerator via a superconducting heat switch (HS-1). The second stage (NDS-2) is directly linked to the sample space (load) and is separated from the first NDS via another heat switch unit (HS-2). The stages may be designed asymmetrically, with NDS-2 containing a lower quantity of the nuclear refrigerant than NDS-1, or with a different nuclear refrigerant altogether.

The operation of the S-CNDR consists of three main steps and can be described as follows (see Fig. 2). Assuming steady-state operation, we start with NDS-2 at the sample temperature and at a moderate magnetic field (of order 60 mT), and thermally decoupled from NDS-1, which is at low temperature and minimum magnetic field (≈10 mT).

Step 1: a slow demagnetization is started on NDS-2 (Fig. 2:a) to keep the sample temperature stable and offset any external heat leaks, while NDS-1 is magnetized to full magnetic field (of order 1 T) and warms up (Fig. 2:A 1 ). Eventually, it is thermally coupled to the MXC and allowed to exchange heat (Fig. 2:A 2 ), transferring its excess entropy to the MXC. Step 2: once the final magnetic field and a temperature sufficiently close to that of the MXC are reached, NDS-1 is decoupled from the MXC and its demagnetization begins (Fig. 2:B), aiming for a temperature just below that of the sample space. Step 3: once the temperature of NDS-1 drops below that of NDS-2 (this should coincide with the end of the slow demagnetization of NDS-2), the two stages are thermally linked via the HS and heat exchange between them takes place (Fig. 2:C), leading to transfer of entropy from NDS-2 to NDS-1. During this time, the demagnetization of NDS-1 may continue for some time until the minimum field (≈10 mT) is reached, while NDS-2 undergoes a magnetization to its starting value of the magnetic field (Fig. 2:b) accompanied with a thermal relaxation to NDS-1 temperature (Fig. 2:c). We note that the heat transfer rate needed at this point is given jointly by the amount of heat produced in NDS-2 during its magnetization and the sample space heat leak. As the performance-limiting factor of the S-CNDR is the heat transfer rate between the two stages, it might be advantageous to perform the magnetization of NDS-2 in the initial phase of this step, and thus increase the total heat transferred by inducing a higher temperature difference between NDS-2 and NDS-1, albeit at the cost of reduced temperature stability of the sample space.

On the other hand, the parallel configuration (P-CNDR) requires each of the two NDS to be linked to both the MXC and the sample space via separate HS units. In this case, both NDS should be designed symmetrically, with balanced heat capacities. The extra heat switches necessary for the P-CNDR configuration represent some design complications as each HS unit requires its own small superconducting solenoid which must be mounted on the dilution refrigerator and operated independently. When choosing between S-CNDR and P-CNDR, one must therefore weigh any benefits found in the performance of the P-CNDR against the increased complexity of the setup and its larger footprint on the dilution refrigerator. The operation of the P-CNDR has two main steps (see Fig. 3) during which either of the NDS is in thermal contact with the sample space and providing the cooling via a slow demagnetization (Fig. 3:A), while the other NDS is being "regenerated" at the MXC temperature. The regeneration procedure involves decoupling from the sample and magnetization to the maximum field (Fig. 3:B 1 ), heat exchange with the MXC (Fig. 3:B 2 ), followed by thermal decoupling and demagnetization to or below the sample temperature (Fig. 3:B 3 ), where the newly regenerated NDS will be switched for its counterpart. We note that the P-CNDR configuration to some extent mitigates the requirements on the quality of the thermal links, as the necessary heat transfer rate through HS-2 or HS-3 is given only by the sample space heat leak, whereas for S-CNDR, the heat produced by the magnetization of the second stage needs to be transferred through HS-2 as well.

Based on the comparison with CADR setups [START_REF] Shirron | Development of a cryogen-free continuous ADR for the constellation-X mission[END_REF], one would expect a higher cooling power for the P-CNDR setup than for S-CNDR, as with the two stages in parallel, heat rejection at MXC occurs during a longer part of the operating cycle. Nevertheless, it remains to be seen if other complications, specific to ultra-low temperatures arise due to the increased number of heat switches and thermal links, potentially influencing the overall heat leak in the sample space or introducing additional dissipative processes into the operation of the refrigerator, such as heat switch manipulation.

Thermal Models

The thermal models corresponding to both configurations are shown in Fig. 4. In the S-CNDR setup, the thermal link between NDS-1 and NDS-2 (including the heat switch HS-2) is the crucial component of the entire setup [START_REF] Toda | Design and expected performance of a compact and continuous nuclear demagnetization refrigerator for sub-mK applications[END_REF][START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF], as Fig. 2 The steady-state operation cycle of the series CNDR. In this design, NDS-1 provides precooling for NDS-2 during its up-magnetization, and rejects entropy at the MXC temperature. Ideally, the loop of NDS-2 should be as narrow as possible to minimize losses and provide a stable temperature, but with realistic durations of each step, a finite span of temperatures is necessary. The span of entropies covered by NDS-2 is determined by the choice of its maximum and minimum magnetic field values.

it determines the heat exchange rate at which NDS-2 can be precooled from NDS-1 during its magnetization. The duration of this step, representing the longest part of the refrigerator cycle, together with the heat leaks involved determine the efficiency of the S-CNDR setup and hence the final temperature it can attain.

The necessary heat exchange between NDS-1 and NDS-2 could be, in principle, shortened if a lower capacity stage was used for NDS-2, or if the range of applied magnetic fields is reduced. Depending on the intended application and the required stability of sample space temperature, it may be possible tune the heat capacity of the second stage, together with the duration of the refrigerator cycle. An optimum choice exists, leading to peak performance. For the purposes of a direct comparison with P-CNDR, we assume a symmetrically designed S-CNDR system, with the "active" heat capacity of the second stage determined by the limiting values of applied magnetic field, even though this need not be the optimum case in terms of S-CNDR performance.

In the P-CNDR setup, properties of the thermal links between the NDS stages and the sample space (Fig. 4: wires 1d and 2d) determine the rate of maximum available cooling power at a given temperature. As mentioned in the previous Section, less heat needs to be transferred at low temperature via these links than for the S-CNDR setup and hence the requirements on the equivalent electrical resistance of the thermal links may be expected to be less stringent.

Fig. 3 The steady-state operation cycle of the parallel CNDR. While the initial cooldown and the first few cycles would differ slightly for both demagnetization stages, once a steadystate is reached, NDS-1 and NDS-2 follow the same operation cycle (out of phase with each other) and both reject entropy at the MXC temperature. As the rate of entropy rejection per cycle is effectively doubled compared to S-CNDR, a higher cooling power is expected for the P-CNDR configuration.

Nevertheless, thermal decoupling due to a finite thermal resistance between the sample and the stages will necessarily reduce the refrigeration efficiency, as any heat removed from the sample will generate additional entropy, in excess of the amount that would be generated if the stage and the sample were at the same temperature. This is fully taken into account in the numerical simulations presented below.

To directly compare the two configurations, we consider in each case two identical demagnetization stages, each containing 0.2 mol of PrNi 5 . The thermal behaviour of the linking elements (wires) is modeled based on an equivalent electrical resistance from the Wiedemann-Franz law, assuming the wires would be made from a high-purity metal such as silver or copper. The same electrical resistance is then used for all the wires.

For most simulations, the heat capacity of the sample space is chosen so to approximate that of 10 g of high-purity Cu at 1 mK and zero magnetic field. Specific heat capacity of 1.15×10 -5 J kg -1 is used, in close agreement with the electronic contribution to the heat capacity of Cu given in Ref. [START_REF] Pobell | Matter and Methods at Low Temperatures[END_REF]). We note that with the same heat load and thermal resistance, the simulation results do not change appreciably if a mass of 100 g or 1000 g is used for the sample space instead. Additional simulations were performed with a sample configuration modeling 100 g of Cu in a magnetic field of 100 mT, including the much larger heat capacity of the nuclear magnetic moments [START_REF] Pobell | Matter and Methods at Low Temperatures[END_REF]. As expected, 4 Thermal models of series (A) and parallel (B) CNDR setups. Note that in the discussion of the serial CNDR setup in Ref. [START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF], the sample link is assumed to have negligible thermal resistance. This can be partly justified by the final design of this element, if a larger piece of high-conductivity metal (copper or silver) is used to cover most of the distance from the stage to the sample space. However, thinner wires would still have to be used to connect to the PrNi 5 inside its solenoid. While this would not greatly impact the results in Ref. [START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF], as NDS-2 is permanently connected to the sample space and has a much higher heat capacity, the results nonetheless represent an ideal case. On the other hand, in the simulations of P-CNDR discussed here, no such idealized thermal links are used.

the extra heat capacity had no effect on the final temperatures obtained, but resulted in a longer settling time, in this case comparable to the duration of one cycle of the CNDR.

The thermal models of individual components of the setup are detailed in Ref. [START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF]. The thermal resistance, R(T ), of any linking element consisting of an Al heat switch and connecting leads is given as the sum of the thermal resistances of all its parts and of boundary (contact) resistances. For the Al heat switches as such, different relations were employed for their normal and superconducting state, describing heat conduction by electrons (∝ T ) and phonons (∝ T 3 ). A residual resistivity ratio (RRR) of 11000 was used to describe the normal state, resulting in heat conductivity κ N = 1.01 × 10 4 [Wm -1 K -2 ]T, while the heat conductivity κ S = 0.68[Wm -1 K -4 ]T 3 was used for the superconducting state, yielding a switching ratio κ S /κ N = 6.73[K -2 ]T 2 , in agreement with the low temperature limit of Al heat switch properties given in Ref. [START_REF] Pobell | Matter and Methods at Low Temperatures[END_REF]. For convenience, the contact thermal resistances were modeled using a temperature-independent equivalent electrical resistance and the Wiedemann-Franz law. A similar approach has been used to model the wires in thermal links. For the nuclear refrigerant, an empirical model of its thermo-magnetic properties was used [START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF] that interpolates between the paramagnetic behaviour at higher temperatures and the ferromagnetic ordering that sets in near 0.4 mK. Agreement with experimental data of Ref. [START_REF] Folle | Nuclear refrigeration properties of PrNi5[END_REF] was the driving criterion in the devising of this model (see Ref. [START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF]).

Additional heat leaks (vibrational heating Qvibr = 10 -8 [WT -1 ] |B|, eddy current heating Qeddy = 0.03[WT -2 s 2 ] Ḃ2 , plus an additional constant heat leak of 2 nW) are considered for each nuclear demagnetization stage in the same way as in Ref. [START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF]. The dependencies used and the values of the prefactors are adapted from Ref. [START_REF] Parpia | Optimization procedure for the cooling of liquid 3He by adiabatic demagnetization of praseodymium nickel[END_REF], where such behaviour was observed experimentally on a nuclear demagnetization refrigerator in the relevant range of magnetic fields. We note that other works have reported different dependencies of the vibrational heating on the magnetic field [START_REF] Todoshchenko | Dry demagnetization cryostat for sub-millikelvin helium experiments: Refrigeration and thermometry[END_REF], and that in general, vibrational heating will depend on the details of the experimental setup. The numerical computations are again described in detail in Ref. [START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF], as well as the model of entropy of PrNi 5 used throughout this work.

Results and Discussion

For the comparison of the performance of the two CNDR configurations, the data on final temperature vs. heat leak from Ref. [START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF] will be used for the S-CNDR setup, while new results will be presented here for the P-CNDR setup.

First, a direct comparison of sample time-traces of temperatures and magnetic fields for the S-CNDR and P-CNDR setups is given in Fig. 5 for a 150 nΩ equivalent electrical resistance of the thermal links (modeling two bundles of copper wires and an aluminium heat switch in normal state) and 5 nW sample space heat leak (representing an optimistic estimate, see, e.g., Ref. [START_REF] Parpia | Optimization procedure for the cooling of liquid 3He by adiabatic demagnetization of praseodymium nickel[END_REF] for comparison). The small jumps in sample temperature for the S-CNDR setup are due to the relatively rapid magnetization of the second stage. While these could be avoided by magnetizing more slowly, the overall performance would be degraded due to a lower net heat transfer between the stages, resulting in a higher cycle-averaged temperature. Regardless, the S-CNDR setup fails to reach the 0.5 mK mark (highlighted in Fig. 5) even in the cycle minimum temperature, while the P-CNDR setup, once cooled down, stays below this mark consistently.

Similar numerical simulations as those shown in Fig. 5 were performed for the P-CNDR setup operating with different values of the equivalent electrical resistance of the thermal links and different values of sample space heat leak. For each case, steady-state operation was established and the maximum sample temperature within a cycle was obtained. These results are summarized in Fig. 6 and compared to the data of Ref. [START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF]. For the S-CNDR setup, due to considerable temperature variations present, cycle-averaged temperatures were 1. In all cases, the P-CNDR setup reaches lower temperatures than the S-CNDR setup. Moreover, the results indicate that the rather stringent requirements for the operation of the S-CNDR setup (thermal links with equivalent electrical resistance comparable to 150 nΩ) are to some extent alleviated in the P-CNDR configuration. Temperatures close to 1 mK should be attainable with P-CNDR even with significantly worse thermal links of 500 nΩ equivalent resistance under a heat load of 20 nW.

Conclusions

While one must bear in mind that the simulations do not necessarily describe all experimental facts accurately, our results clearly show that, given the same quality of thermal links, the P-CNDR configuration has superior performance to the S-CNDR design in terms of the final temperature reached at a given heat leak, or in terms of the cooling power attained at a given temperature.

We note, however, that the superior cooling power of the P-CNDR comes at the cost of additional system complexity. Additionally, this design may lead to reduced stability of sample space temperature not captured in the presented simulations, as switching the sample connection between the two stages may in practice induce significant temperature variations due to the required manipulation of two heat switches. Indeed, one drawback of P-CNDR is the necessity to operate a total of four heat switches together with their solenoids instead of just two such units in the S-CNDR design, leading to additional heat leaks and delays within the refrigerator cycle along with the setup footprint and complexity. Nevertheless, aside from special cases that may be better suited to the S-CNDR design, we expect the P-CNDR to perform better than the S-CNDR for a given thermal link conductance, and the P-CNDR design clearly warrants experimental study. (bottom) setups, as given by the numerical simulations for a 150 nΩ equivalent electrical resistance of the thermal links and 5 nW sample space heat leak. Parabolic profiles for magnetization/demagnetization are used. After a stabilization time of approximately two cycles, a steady state is reached and the maximum sample temperature during a cycle can be extracted. The durations of (de-)magnetization steps and values of magnetic fields are tuned by hand to achieve near-optimum performance. The starting conditions are chosen with the MXC temperature of ≈8 mK and the maximum magnetic field used on the demagnetization stages is 1.2 T, as in Ref. [START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF]. The (gray) dashed line marks the temperature of 0.5 mK. The data in the left panel are adapted from Ref. [START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF]. It is clearly shown that the P-CNDR setup can achieve lower temperatures with the same quality of thermal links under the same heat load. Conversely, this implies that at a given temperature (in the mK range), the P-CNDR setup indeed has a higher cycle-averaged cooling power than the S-CNDR setup, as expected. We note that if cycle-averaged temperatures were used instead of the maximal ones, this conclusion would hold, see text and Table 1.

Fig. 1

 1 Fig.1The schematic design of serial (A) and parallel (B) configuration of the CNDR.

Fig. 5

 5 Fig.5Time traces of magnetic fields and temperatures for the S-CNDR (top) and P-CNDR (bottom) setups, as given by the numerical simulations for a 150 nΩ equivalent electrical resistance of the thermal links and 5 nW sample space heat leak. Parabolic profiles for magnetization/demagnetization are used. After a stabilization time of approximately two cycles, a steady state is reached and the maximum sample temperature during a cycle can be extracted. The durations of (de-)magnetization steps and values of magnetic fields are tuned by hand to achieve near-optimum performance. The starting conditions are chosen with the MXC temperature of ≈8 mK and the maximum magnetic field used on the demagnetization stages is 1.2 T, as in Ref.[START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF]. The (gray) dashed line marks the temperature of 0.5 mK.

Fig. 6

 6 Fig.6Comparison of the final temperature (maximum sample space temperature during a cycle) vs. sample space heat leak for the serial (left) and parallel (right) CNDR setups. The data in the left panel are adapted from Ref.[START_REF] Schmoranzer | Development of a sub-mK Continuous Nuclear Demagnetization Refrigerator[END_REF]. It is clearly shown that the P-CNDR setup can achieve lower temperatures with the same quality of thermal links under the same heat load. Conversely, this implies that at a given temperature (in the mK range), the P-CNDR setup indeed has a higher cycle-averaged cooling power than the S-CNDR setup, as expected. We note that if cycle-averaged temperatures were used instead of the maximal ones, this conclusion would hold, see text and Table1.

Table 1

 1 Final temperatures obtained from numerical simulations of the S-CNDR and P-CNDR setups. For the S-CNDR setup, both maximal and cycle-averaged temperatures are shown, for P-CNDR, only the maximal temperatures are given. See also Fig.6.

	Link resistance	Heat leak	S-CNDR max.	S-CNDR avg. P-CNDR max.
	(nΩ)	(nW)	(mK)	(mK)	(mK)
	50	5	0.52	0.48	-
	50	10	0.61	0.57	-
	50	20	0.78	0.70	-
	150	5	0.62	0.57	0.47
	150	10	0.75	0.70	0.60
	150	20	1.03	0.94	0.76
	500	5	0.99	0.94	0.62
	500	10	1.33	1.29	0.82
	500	20	1.95	1.90	1.13
	1000	5	-	-	0.83
	1000	10	-	-	1.13
	1000	20	-	-

obtained as well. All values of temperature are listed in Table

Acknowledgements We acknowledge support from the ERC StG grant UNIGLASS No. 714692 and ERC CoG grant ULT-NEMS No. 647917. The research leading to these results has received funding from the European Unions Horizon 2020 Research and Innovation Programme, under Grant Agreement no 824109. DS acknowlegdes support from the Czech Science Foundation (GACR) under project No. 20-13001Y.