e-mail : thierry.paillard@univ-pau.fr

Résumé

Avec l'avancée en âge, la diminution de la densité minérale osseuse est inéluctable. L'exercice physique régulier constitue un moyen physiologique susceptible d'atténuer les effets de la déminéralisation osseuse normale (i.e. non pathologique) particulièrement chez la personne âgée. En effet, l'exercice physique induit des contraintes mécaniques générant une déformation de l'os qui stimulent l'ostéogénèse et favorisent le remodelage osseux. La fonction ostéogénique des entraînements aérobies (e.g. marche, course) n'est effective que si l'intensité de l'exercice (i.e. les impacts au sol et donc les déformations osseuses) est de haute intensité et celle des entraînements de renforcement musculaire est effective que si les contractions musculaires effectuées sont dynamiques et de forte intensité. Les activités physiques réalisées en condition de décharge corporelle (e.g. natation, vélo) ou en condition statique (e.g. étirements, exercices d'équilibration) sont en revanche des activités parmi les moins ostéogéniques. La durée optimale d'un programme d'exercices significatif en termes de remodelage osseux se situe entre 4 et 6 mois. Les effets ostéogéniques de l'exercice régulier débute à partir de 2-3 séances hebdomadaires. L'activation de la fonction ostéogénique par le biais de l'exercice physique est plus difficile pour la femme vieillissante que pour l'homme vieillissant à cause de ses facteurs hormonaux moins favorables.

Abstract

With age advancement, the decrease of bone mineral density is ineluctable. Physical exercise constitutes a physiological approach likely to attenuate or limit the effects of normal bone demineralization (i.e. not pathological) particularly in elderly subjects. Indeed, physical exercise induces mechanical constraints generating bone deformation which stimulates osteogenesis and favors bone remodelage. Physical activities achieved in condition of body discharge (e.g. swimming, cycling) or in static condition (e.g. stretching, balance) do not stimulate (or very weakly) osteogenesis. The osteogenic function of aerobic training (e.g. walking, running) is effective only if the intensity of exercise is high (i.e. the impacts on the ground and thus the bone deformation) and that of strength training is effective only if the completed muscular contractions are dynamic and carried out with heavy loads. The calcium concentration increase is greater on the concave side than on the convex side for the bones which undergo strong mechanical pressures during exercise. Hence, it is advisable to vary the directions of mechanical constraints during physical activity to strengthen the resistance of the bone in all the plans. In order to obtain significant effects in terms of bone remodelage, the optimal duration of training programs should last at least 4 to 6 months. The osteogenic effects of regular exercise begin from 2-3 weekly sessions. The activation of osteogenesis by means of physical exercise is more difficult in aging women than in aging men because of hormonal factors that are not favorable in aging women. At last, regular exercise is fundamental not only to maintain bone mineral density but also to reduce the risk of bone fracture since there is a relationship between the bone mineral density and the risk of bone fracture.

Key words: elderly subjects, physical exercise, bone mineral density, amount of physical activity, risk of fracture Avec l'avancée en âge, la diminution de la densité minérale osseuse est inéluctable. Elle se traduit par une diminution de la résistance mécanique de l'os à toute sorte de force (e.g. choc, pression, contrainte mécanique) qui s'exerce sur lui. Le risque de chute lié à l'altération de la fonction d'équilibration augmentant également sous l'effet de l'avancée en âge, le risque de fracture croît en conséquence avec l'âge [START_REF] Nitz | Health-related predictors of falls and fractures in women over 40[END_REF].

L'exercice physique régulier constitue un moyen physiologique (avec l'alimentation) susceptible d'atténuer les effets de la déminéralisation osseuse normale (i.e. non pathologique) particulièrement chez la personne âgée. En effet, l'exercice physique induit des contraintes mécaniques générant une déformation de l'os qui stimulent l'ostéogénèse et favorisent le remodelage osseux chez les personnes vieillissantes [START_REF] Srinivasan | Bone mechanotransduction may require augmentation in order to strengthen the senescent skeleton[END_REF].

A ce titre, les investigations observationnelles et interventionnelles ayant pour objectif d'étudier les effets préventifs de l'exercice régulier contre la déminéralisation osseuse ont été relativement conséquentes au cours des deux dernières décennies. Toutefois, même s'il est actuellement possible d'identifier les activités physiques favorables à l'ostéogenèse, la question des effets quantitatifs et qualitatifs de l'exercice considéré comme ostéogénique en termes de paramètres d'intensité, de durée et de fréquence sur la prévention de la déminéralisation osseuse chez les sujets âgés reste encore ouverte. C'est pourquoi l'objectif de cette courte revue consiste à présenter les principaux résultats des études traitant des effets préventifs des différents paramètres (cités ci-dessus) des exercices pratiqués chez la personne âgée contre l'altération de son appareil osseux.

Vieillissement du tissu osseux

La matière osseuse est un tissu vivant en perpétuel remaniement. Des cellules osseuses disparaissent et se créent quotidiennement. Deux voies antagonistes interagissent l'une sur l'autre (voie Wnt ou wingless-related integration site et voie RANK/RANKL/OPG ou receptor activator of nuclear factor-kappa B/receptor activator of nuclear factor-kappa B ligand/osteoprotegerin) pour assurer l'équilibre entre la formation et la résorption osseuses [START_REF] Mitchell | Clinical impact of recent genetic discoveries in osteoporosis[END_REF]. Ces deux voies régulent l'activité des ostéoblastes (responsables de la formation du tissu osseux) et des ostéoclastes (responsables de la résorption du tissu osseux). La prévalence de l'ostéoporose augmente avec l'avancée en âge. Les fractures ostéoporotiques affectent particulièrement l'extrémité proximale du fémur (e.g. le col du fémur), les corps vertébraux et l'extrémité distale du radius [START_REF] Khosla | Osteoporosis in men[END_REF]. L'ostéoporose se manifeste après 40 ans et se montre plus sévère chez les femmes que chez les hommes surtout après la ménopause [START_REF] Guggenbuhl | Osteoporosis in males and females: Is there really a difference?[END_REF]. Après 65 ans, le risque fracturaire serait plus de 2 fois plus élevé chez la femme que chez l'homme [START_REF] Mitchell | Clinical impact of recent genetic discoveries in osteoporosis[END_REF][START_REF] Khosla | Osteoporosis in men[END_REF]. En revanche, les fractures de la hanche chez l'homme sont associés à une plus grande morbidité (et mortalité) que chez la femme [START_REF] Khosla | Osteoporosis in men[END_REF][START_REF] Guggenbuhl | Osteoporosis in males and females: Is there really a difference?[END_REF]. Outre le facteur sexe, d'autres facteurs tels que l'inactivité physique, la déficience hormonale (en oestrogène chez la femme ménopausée), l'insuffisance nutritionnelle, les facteurs génétiques et le style de vie sont susceptibles d'accentuer ou d'atténuer le processus ostéoporotique [START_REF] Stransky | Nutrition as prevention and treatment of osteoporosis[END_REF]. A partir d'un certain seuil de déminéralisation osseuse, les os concernés par cet état biologique perdent de leur solidité et sont particulièrement exposés au risque fracturaire (Figure 1).

Effets préventifs de l'exercice

Le processus d'involution du tissu osseux peut être ralenti voire éventuellement inversé dans certaines circonstances très favorables chez les personnes âgées. Outre les conditions nutritionnelles et environnementales propices pour lutter contre ce processus, les sujets âgés doivent soumettre leur squelette à des contraintes mécaniques en pratiquant des exercices musculaires. En effet, les contraintes mécaniques induisant un certain seuil de déformation de l'os stimulent l'ostéogénèse et favorisent le remodelage osseux chez les personnes vieillissantes [START_REF] Srinivasan | Bone mechanotransduction may require augmentation in order to strengthen the senescent skeleton[END_REF]. L'exercice physique stimule l'activité des ostéoblastes et facilite l'activité des ostéocytes [START_REF] Ozcivici | Mechanical signals as anabolic agents in bone[END_REF] et l'absorption intestinale du calcium [START_REF] Stransky | Nutrition as prevention and treatment of osteoporosis[END_REF]. Il impacterait négativement la production des ostéocytes en sclérostine dont la fonction est d'inhiber les ostéoblastes et (donc) la formation osseuse et induirait une augmentation de la production de facteurs anaboliques liée au métabolisme osseux [START_REF] Bosković | Adapted physical activity in the prevention and therapy of osteoporosis[END_REF]. L'activité physique agirait favorablement sur la masse et l'architecture osseuses chez la personne âgée.

Cependant, les effets de l'activité physique sur la DMO (exprimée en g/cm2) ne sont pas univoques. Pour certains auteurs, un programme d'activité physique peut retarder ou même inverser le processus de déminéralisation osseuse. Pour d'autres auteurs, les effets d'un programme d'entraînement seraient moins significatifs sur la DMO. En fait, le processus de minéralisation osseuse peut être plus ou moins influencé par la nature, l'intensité, la durée et la fréquence des activités physiques entreprises. Chez les sujets âgés, le facteur sexe peut également influencer les capacités d'adaptation de leurs appareils osseux sous l'effet de l'exercice musculaire.

Le processus ostéogénique dépend en partie de la nature de l'activité physique pratiquée et de l'impact des contraintes mécaniques qu'elle induit sur le squelette. Les entraînements aérobies de haute intensité et le travail de renforcement musculaire sont de nature à activer l'ostéogénèse. Un entraînement aérobie de 18 mois composé de marche, de montées d'escalier, d'efforts sur ergocycle et de jogging à raison de 3 fois par semaine, pendant 50 minutes, à intensité soutenue (55 à 75% du VO2 max), favorise le maintien de la DMO au niveau du col du fémur chez des femmes ménopausées âgées de 52-53 ans [START_REF] Heinonen | Effect of two training regimens on bone mineral density in healthy perimenopausal women: a randomized controlled trial[END_REF]. Globalement, les activités qui impliquent les plus fortes contraintes mécaniques sont celles qui engendrent les effets positifs les plus importants sur le tissu osseux. En effet, des activités comme la danse (rythmique), la gymnastique (dynamique) et le step qui impliquent la réalisation d'impulsions et de réceptions corporelles générant des impacts au sol importants, sont des activités parmi les plus ostéogéniques susceptibles d'être pratiquées par les personnes âgées alors que des activités comme la natation et le vélo qui se pratiquent en condition de décharge corporelle n'engendrant pas d'impacts au sol sont des activités parmi les moins ostéogéniques. Enfin, pour une activité physique donnée telle que le renforcement musculaire au pouvoir ostéogénique important, certains mouvements sont davantage ostéogéniques que d'autres. Les actions musculaires de tirage seraient plus efficaces que les actions de portées car elles entraînent une déformation (contrainte mécanique) osseuse plus importante [START_REF] Stransky | Nutrition as prevention and treatment of osteoporosis[END_REF].

La nature des contraintes mécaniques demeure fondamentale au processus de minéralisation osseuse. Les contraintes dynamiques sont davantage ostéogéniques que les contraintes statiques. En outre, l'orientation des contraintes est également déterminante. L'idéal étant de pratiquer une activité physique dynamique exerçant des pressions mécaniques variées dans leurs directions sur un grand nombre de segments. La concentration calcique est plus importante du côté concave que du côté convexe de la déformation osseuse induite par la pression mécanique [START_REF] Paillard | Vieillissement et Condition Physique[END_REF]. D'où l'importance de varier la direction des contraintes mécaniques au cours de l'activité physique pour renforcer la résistance de l'os dans les 3 plans (quand c'est possible !).

Une méta analyse de Layne et Nelson [START_REF] Layne | The effects of progressive resistance training on bone density: a review[END_REF] a fait le point sur l'optimisation de l'intensité de l'entraînement à prescrire afin de favoriser la minéralisation osseuse. Cette étude atteste qu'une grande partie des travaux effectués sur ce thème abonde dans le sens où un programme de renforcement musculaire de haute intensité agit de façon bénéfique sur la DMO. Taaffe et al. [START_REF] Taaffe | Comparative effects of high-and lowintensity resistance training on thigh muscle strength, fiber area, and tissue composition in elderly women[END_REF] ont effectué une étude d'une durée de 52 semaines chez des femmes âgées de 65 à 69 ans pratiquant 3 fois par semaine des exercices de renforcement musculaire. Ils ont constitué 3 groupes, dont un a servi de groupe témoin. Les deux autres groupes ont réalisé des programmes de musculation de haute (80% 1 RM i.e. répétition maximale) ou de basse intensité (40% 1 RM). Seul, le groupe ayant pratiqué à une haute intensité a pu maintenir sa valeur initiale de DMO au niveau des fémurs. Le groupe témoin et le groupe ayant pratiqué à 40% d'intensité ont vu leur DMO régresser. Par ailleurs, selon Pruitt et al. [START_REF] Pruitt | Effects of a one -year high -intensity versus lowintensity resistance training program on bone mineral density in older women[END_REF], ce même programme n'engendre pas de différence de DMO entre les 3 groupes au niveau de la colonne vertébrale. A partir de là, la DMO se stabiliserait seulement sur les sites osseux particulièrement stimulés par l'exercice musculaire. Deux revues de la littérature concluent que le renforcement musculaire de haute intensité peut même augmenter la DMO sur les sites anatomiques qui subissent les plus fortes pressions mécaniques chez des femmes ménopausées [START_REF] Zehnacker | Effect of weighted exercises on bone mineral density in post menopausal women. A systematic review[END_REF][START_REF] Gutin | Can vigourous exercise play a role in osteoporosis prevention. A review[END_REF]. Une activité physique stimulant particulièrement les quadriceps engendre une augmentation de la DMO du col du fémur, alors qu'une stimulation de l'iliopsoas génère une élévation de la DMO de la colonne vertébrale [START_REF] Paillard | Vieillissement et Condition Physique[END_REF].

La durée d'un programme d'entraînement influe sur le processus de minéralisation osseuse.

Un programme de renforcement musculaire relativement intense (monter 300 marches d'escaliers par séance, 4 séances hebdomadaires) d'une durée de 6 semaines chez des femmes âgées de 62 à 75 ans n'a pas suffi à augmenter la DMO du col du fémur et de la colonne lombaire (L1-L4) [START_REF] Paillard | Comparison between three strength development methods on body composition in healthy elderly women[END_REF]. La durée optimale d'un programme d'exercices significatif en termes de remodelage osseux se situe entre 4 et 6 mois [START_REF] Bonaiuti | The conclusions of the systematic revision with the results from the major clinical and experimental studies[END_REF]. Cependant, une durée plus longue ne constitue pas un gage d'amélioration systématique puisqu'un programme de musculation de forte intensité étalé sur deux ans (2 périodes de 42 semaines séparées par 10 semaines de repos et de tests), n'a pas induit d'augmentation de la DMO de la colonne lombaire et du corps entier chez les 57 sujets du groupe expérimental par rapport aux 56 sujets du groupe contrôle (tous âgés entre 60 et 80 ans) [START_REF] Mccartney | A longitudinal trial of weight training in the elderly: continual improvements in years 2[END_REF]. Par ailleurs, le facteur sexe influence également les effets induits par l'activité physique sur la DMO. En effet, une activité comme la marche qui induit certaines contraintes mécaniques et demeure donc relativement ostéogénique, peut s'avérer inefficace pour lutter contre la diminution de la DMO chez les femmes ménopausées [START_REF] Gutin | Can vigourous exercise play a role in osteoporosis prevention. A review[END_REF][START_REF] Kerr | Exercise effects on bone mass in postmenopausal women are site-specific and load-dependent[END_REF]. En revanche, la pratique de la marche n'engendre certes, pas d'augmentation de la DMO chez des sujets masculins ayant sensiblement le même âge mais elle la stabilise [START_REF] Paillard | Effects of brisk walking on postural balance, locomotion, body composition and aerobic capacity in healthy ageing men[END_REF]. Les contraintes mécaniques même fréquentes induites par la marche sont probablement insuffisantes pour amplifier le processus de minéralisation osseuse chez les sujets âgés tandis qu'elles peuvent vraisemblablement, au mieux, stabiliser la DMO surtout chez les femmes âgées. Seule la pratique de la marche intense et pratiquée à un rythme très soutenu aurait un rôle ostéogénique significatif. Hatori et al. [START_REF] Hatori | The effects of walking at the anaerobic threshold level on vertebral bone loss in postmenoposal women[END_REF] ont en effet comparé 3 groupes composés de femmes ménopausées, deux ont réalisé un programme de marche (30 minutes, 3 fois par semaine, pendant 7 mois) dont un à une intensité correspondant à une valeur de fréquence cardiaque supérieure au seuil anaérobie (n=12) et l'autre à une intensité inférieure à ce seuil (n=9), le troisième était le groupe contrôle (n=12). Seul le groupe ayant pratiqué à une intensité supérieure au seuil anaérobie a élevé sa DMO lombaire (+1,1%) évaluée par absorptiométrie biphotonique à rayon X par comparaison au groupe ayant pratiqué à une intensité inférieure au seuil anaérobie (-1,0%) et au groupe contrôle (-1,7%) dont la DMO lombaire a diminué [START_REF] Hatori | The effects of walking at the anaerobic threshold level on vertebral bone loss in postmenoposal women[END_REF]. Par ailleurs, dans le cadre d'un programme d'activités physiques aérobies combinées (ergocycle, marche, stretching) de 14 mois chez des sujets âgés de plus de 60 ans, la DMO a augmenté chez les hommes mais pas chez les femmes [START_REF] Blumenthal | Effects of exercise training on bone density in older men and women[END_REF]. Globalement, il est probable que l'exercice physique chez les femmes ménopausées puisse ralentir ou stabiliser le processus de déminéralisation osseuse mais pas l'inverser ou seulement dans le cadre d'activités fortement ostéogéniques (à impacts mécaniques élevés). En effet, dans le cadre d'un programme de renforcement musculaire, des femmes âgées (n=20) soumises à 2 entraînements hebdomadaires pendant 1 an des muscles ilio-psoas, rotateurs latéraux (ou pelvi-trochantériens) et érecteurs du rachis lombaire ont vu leur DMO fémorale et lombaire s'accroître de 0,9% et 1,0% tandis que les femmes du groupe témoin (n=19) de ce protocole expérimental ont vu leur DMO diminuer de -2,5% et -1,8% respectivement (p=0,02 et p=0,04 respectivement) [START_REF] Nelson | Effects of high-intensity strength training on multiple risk factors for osteoporotic fractures. A randomized controlled trial[END_REF].

Lorsque les conditions environnementales, nutritionnelles et physiologiques sont favorables, les activités de renforcement musculaire intenses pourraient augmenter modérément la DMO chez les sujets âgés, tandis que les activités de type endurance (e.g. la marche) pratiquées isolément ne peuvent probablement que la maintenir (ou l'augmenter faiblement). Les activités pratiquées en décharge corporelle (vélo, natation…) s'avèrent, quant à elles, manifestement inopérantes dans une optique d'augmentation ou même de stabilisation de la DMO. Par ailleurs, l'efficacité d'un exercice physique relève plus de son intensité et de sa fréquence que de sa durée.

Risque de fracture

Le taux d'incidence de fractures de hanches et de fractures vertébrales demeure plus faible chez des sujets âgés physiquement actifs comparativement à des sujets âgés non actifs [START_REF] Marks | Hip fractures among the elderly: causes, consequences and control[END_REF][START_REF] Sinaki | Stronger back muscles reduce the incidence of vertebral fractures: a prospective 10 year follow-up of postmenopausal women[END_REF][START_REF] Cumming | Epidemiology of hip fractures[END_REF]. Le suivi de 9704 femmes de plus de 65 ans durant environ 8 ans (7,6 années) indique que les personnes appartenant aux quintile et quartile les plus élevés en termes de niveau (volume) d'activité physique présentent un taux de fractures de hanche 36% (RR=0.64 ; IC 95% [0.45-0.891]) plus faible que celles appartenant au quintile le plus bas [START_REF] Gregg | Physical activity and osteoporotic fracture risk in older women. Study of osteoporotic fractures research group[END_REF]. Ces auteurs précisent que les femmes actives qui pratiquent à une intensité allant de modérée à vigoureuse diminuent de 42% et 33% le risque de fracture de hanche et de vertèbre respectivement en comparaison avec les femmes inactives. Dans le même sens, il existerait une relation entre la DMO et le risque fracturaire ostéoporotique [START_REF] Sakai | Shorter unipedal standing time and lower bone mineral density in women with distal radius fractures[END_REF][START_REF] Lane | Epidemiology, etiology, and diagnosis of osteoporosis[END_REF][START_REF] Dargent-Molina | Use of clinical risk factors in elderly women with low bone mineral density to identify women at higher risk of hip fracture: The EPIDOS prospective study[END_REF]. Néanmoins, outre la valeur de la DMO, d'autres facteurs environnementaux et génétiques ainsi que la probabilité de chute accentuent ou atténuent ce risque fracturaire [START_REF] Kaptoge | Low BMD is less predictive than reported falls for future limb fractures in women across Europe: results from the European Prospective Osteoporosis Study[END_REF]. Le temps de maintien en appui monopodal demeure particulièrement discriminant en termes de risque fracturaire du radius chez les sujets de plus de 65 ans, particulièrement s'il est inférieur à 15 secondes [START_REF] Sakai | Shorter unipedal standing time and lower bone mineral density in women with distal radius fractures[END_REF]. En pratique, à titre préventif, il conviendrait de soumettre les sites particulièrement fragiles au regard du risque fracturaire à des contraintes mécaniques fréquentes et relativement intenses avec des activités particulièrement adaptées au statut physiologique de l'individu concerné. Il s'avère par ailleurs que plus un individu est actif, plus son risque de chute et donc de fracture augmente, surtout s'il pratique en milieu naturel comportant certains risques aggravants (e.g. terrain glissant, bosselé).

Conclusion

  Exercice et densité minérale osseuse chez les sujets âgés : implications théoriques et pratiques Exercise and bone mineral density in old subjects: theorical and practical implications Thierry Paillard Laboratoire Activité Physique Performance et Santé (UPRES EA 4445), Université de Pau & Pays de l'Adour, Département STAPS de Tarbes, ZA Bastillac Sud, 65000 Tarbes, France. Thierry Paillard Université de Pau et des Pays de l

  Mots clés : sujet âgé, exercice physique, densité minérale osseuse, volume d'activité physique, risque de fracture Points clés :  L'exercice physique dynamique (i.e. des activités physiques générant des pressions (impacts) mécaniques sur les os versus des activités n'en générant pas telles que celles réalisées en condition de décharge corporelle e.g. natation, vélo) régulier stimule l'activité des ostéoblastes et facilite l'activité des ostéocytes et l'absorption intestinale du calcium.  Il impacterait négativement la production des ostéocytes en sclérostine dont la fonction est d'inhiber les ostéoblastes et (donc) la formation osseuse et induirait une augmentation de la production de facteurs anaboliques liée au métabolisme osseux.  Il augmenterait la concentration calcique de façon plus importante du côté concave que du côté convexe de la déformation osseuse induite par la pression mécanique ; d'où l'importance de varier la direction des contraintes mécaniques au cours de l'activité physique pour renforcer la résistance de l'os dans tous les plans.  Plus l'intensité et la fréquence de l'exercice sont élevées, meilleure demeure la réponse ostéogénique tandis que la régularité de l'exercice représente un autre paramètre fondamental dans le maintien de la densité minérale osseuse.  Chez la personne âgée, il existerait une relation entre le volume d'activité physique et la densité minérale osseuse (sur les sites anatomiques spécifiquement soumis à des contraintes mécaniques intenses et répétées) ainsi qu'une relation entre la densité minérale osseuse et le risque de fracture osseuse.

  Chez l'enfant et l'adolescent, la formation de la matière osseuse est excédentaire, c'est-à-dire que l'activité des ostéoblastes est supérieure à celle des ostéoclastes. Cet excédent permet la croissance osseuse. Chez le jeune adulte, les processus de résorption et de formation sont normalement étroitement liés de façon à ce que la quantité de résorption osseuse soit compensée par une quantité sensiblement égale de nouvelle formation osseuse ; la masse osseuse régresse néanmoins légèrement et linéairement entre 20 et 50 ans. Ces processus sont régulés au niveau local par une multitude de cytokines et de facteurs de croissance [4]. En revanche, chez le sujet âgé, l'activité des ostéoclastes dépasse celle des ostéoblastes. Le nombre et l'activité des ostéoblastes diminuent et certaines trabécules individuelles disparaissent ou subissent une perforation [5]. Ceci engendre d'une part, une régression de l'espace nécessaire à la formation de nouvelles cellules osseuses, et d'autre part, une désorganisation du réseau trabéculaire. Ainsi, le vieillissement se caractérise par une réduction de la masse et un changement de l'architecture du tissu osseux. Ces modifications histologiques de l'os se définissent par le terme d'ostéopénie et sont également liées à un ralentissement de la vitesse du processus de remodelage osseux. Selon l'organisation mondiale de la santé (OMS), une valeur de la densité minérale osseuse (DMO) comprise entre 1 et 2,5 déviations standard au-dessous de la densité normale moyenne des adultes jeunes correspond à l'ostéopénie et une valeur inférieure à 2,5 déviations standard correspond à l'ostéoporose.

  La fréquence des séances hebdomadaires d'activité physique demeure également un critère déterminant pour lutter efficacement contre le processus de déminéralisation osseuse. Deux séances d'activité physique hebdomadaires impliquant de fortes contraintes mécaniques semblent être le minimum pour favoriser la conservation de la DMO chez le sujet âgé. Dans une optique d'augmentation de la DMO, une fréquence supérieure paraît nécessaire mais ne constitue toutefois pas une garantie de résultat positif selon le statut physiologique et les conditions nutritionnelles et environnementales de l'individu. Chez des sujets âgés fragiles présentant de faibles capacités d'équilibration et de force musculaire ainsi qu'une faible DMO, le risque de chute et/ou de fractures spontanées demeure important. Pour ce type de sujets, le dilemme consiste à proposer une activité susceptible d'améliorer les paramètres physiques et physiologiques précédemment déclinés tout en évitant de mettre en péril leur intégrité physique à travers une activité physique stimulante et dynamique. La technique de vibration entière du corps (le sujet placé en position érigée sur une plate-forme vibrante) constitue à ce titre une certaine opportunité. En effet, certaines études montrent que la pratique régulière de la vibration améliore la DMO (ainsi que la force musculaire et l'équilibre postural) [e.g. 21, 22]. Toutefois, les conclusions d'une très récente méta-analyse [23] divergent notablement à ce propos. Leung et al. [23] suggèrent que les différents paramètres de stimulation (e.g. fréquence de vibration corporelle, inclinaison de la forme-forme vibrante, durée des séances, nombre de séances hebdomadaires, durée du programme, sites anatomiques osseux évalués) peuvent influer sur les résultats. Ainsi, d'autres travaux liés à la technique de vibration corporelle faisant varier ces paramètres méritent d'être effectués pour maitriser davantage ses effets sur la fonction ostéogénique des sujets âgés.

  Chez le sujet âgé, il existerait une relation entre son volume d'activité physique et sa DMO (sur les sites anatomiques spécifiquement soumis à des contraintes mécaniques intenses et répétées) ainsi qu'une une relation entre sa DMO et son risque de fracture osseuse. Dans une optique de maintien de la DMO, la régularité de l'exercice demeure fondamentale (durant toute la vie). Plus l'intensité de l'exercice est élevée meilleure demeure la réponse ostéogénique. Une fréquence d'au moins 2-3 séances d'activité physique hebdomadaires est nécessaire pour optimiser les effets préventifs sur la déminéralisation osseuse.

Figure :

 : Figure :Représentation schématique du processus de vieillissement osseux et de ses conséquences (adaptée de Paillard[START_REF] Paillard | Vieillissement et Condition Physique[END_REF]).

Figure :

 : Figure: Schematic representation of the process of bone alteration related to age and its consequences.

  trame protéique Diminution de l'action des ostéoblastes Facteurs combinés au vieillissement osseux normal Facteurs qui limitent ou accentuent l'ampleur et la vitesse du processus de vieillissement osseux Sexe (homme ou femme) Style de vie (consommation ou non d'alcool, de tabac et autres produits) toxiques) Facteurs génétiques (favorables ou défavorables)


The paper presents the recent developments of the acceleration techniques for the Method Of Characteristics (MOC) in the code APOLLO2. The main contribution concerns the introduction of a multidomain DP N technique where all regions belonging to the same macro-domain are coupled by an integral operator that is strictly equivalent to the MOC. Different macro-domains are coupled via currents that are defined with the DP N formalism. This new Integral DP N (IDP N ) operator is build by using transmission and escape probabilities factors that respect symmetries/anti-symmetries and complementary properties that are enforced to preserve the physics of the problem and to save computational effort. These factors are computed using the numerical tracking of the MOC operator. The paper presents results on realistic assembly calculations that demonstrate the effectiveness of the IDP N operator as a synthetic acceleration tool.

INTRODUCTION

In this paper we introduce an evolution of the DP N operator currently used in the two dimensional MOC module of the APOLLO2 code [START_REF] Sanchez | APOLLO II: A user-oriented, portable, modular code for multigroup transport assembly calculations[END_REF] (called TDT) to accelerate the free iterations for the MOC. [START_REF] Halsall | CACTUS, a Characteristics Solution to the Neutron Transport Equation in Complicated Geometries[END_REF][START_REF] Smith | CASMO-4 characteristic methods for two-dimensional PWR and BWR core calculations[END_REF][START_REF] Dahmani | An efficient preconditioning technique using Krylov subspace methods for 3D characteristics solvers[END_REF][START_REF] Hong | CRX: a code for rectangular and hexagonal lattices based one the method of characteristics[END_REF][START_REF] Zim | Fundamentals and recent developments of reactor physics methods[END_REF] Much work has been dedicated in recent years to the theme of acceleration techniques for the MOC. The reason for all this interest is the capability of the MOC to solve unstructured geometry problems without using any homogenization procedure. As with every method that uses the source iteration algorithm, the MOC can be very slow for collision dominated problems. Therefore acceleration techniques are essential in order to permit a wider practical application of the MOC. The DP N method, originally introduced by Sanchez and Chetaine [START_REF] Sanchez | A Synthetic acceleration for a two dimensional characteristic method in unstructured meshes[END_REF], is the operator to be used in the framework of the synthetic acceleration procedure for TDT. Generalizations of this technique have been recently introduced in [START_REF] Santandrea | Acceleration Techniques for the characteristics method in unstructured meshes[END_REF][START_REF] Santandrea | Analysis and improvements of the DP N acceleration technique for the method of characteristics in unstructured meshes[END_REF] to improve the physical approximations of the initial formulation and to speed up the solution procedure. Further improvements are presented here and results are shown to illustrate the speed up factors of the TDT code (which is currently used for industrial calculations as a flux solver of the APOLLO2 [START_REF] Sanchez | APOLLO II: A user-oriented, portable, modular code for multigroup transport assembly calculations[END_REF] code system).

The principal contribution of this paper is a new way of using the DP N acceleration. Much alike the classical Interface Current (IC) [START_REF] Sanchez | New development in the Interface-Current Module TDT[END_REF] treatment of the Collision Probability Methods (CPM), we show that the MOC admits a similar Integral DP N (IDP N ) treatment. In this approach one subdivides the system into macro-domains in each of which one applies an integral formalism that does not introduce any supplementary approximation with respect to the MOC. The coupling between regions belonging to different macro-domains is done via interface DP N currents that result from some simplifying assumptions on the interface flux anisotropy. These assumptions are the only approximations of the method with respect to the MOC. One can expect this operator to behave better than the classic DP N operator, where the DP N approximations were made at the interfaces between every region and not only at the interfaces between macro-domains. The main differences with respect to the standard IC method are that our method is not based on the CPM but on the MOC, that the interface approximations are given by the DP N method and that the resulting equations have a response-matrix form that is different from that of the IC one. Since the numerical derivation of the IDP N method is consistent with the MOC, one can expect that this method will be a good acceleration tool. We have not investigated, as it has been done for the IC method, the behavior of the IDP N method as a stand-alone operator. Only its capability as a means to accelerate the MOC free iterations is considered here. Moreover, even if the IDP N method theory is valid also for three dimensional unstructured geometries, our practical implementations are limited to two dimensional systems, due to the actual capabilities of TDT.

The paper is organized as follows. In section 2 we recall the basis of the MOC and introduce some notation for the following sections. In section 3 the IDP N synthetic operator is introduced and its detailed properties are discussed together with the computational algorithm. In section 4 the external synthetic acceleration is revised and discussed. Section 5 is devoted to the analysis and discussion of results of calculations for a realistic assembly. Conclusions follow. Appendix A and B are dedicated to the description of improvements for the treatment of the equations for the boundary conditions, and to show many symmetries/anti-symmetries and reciprocity relations of the IDP N . We also demonstrate that this operator is robust in the limit of nearly voided media.

THE METHOD OF CHARACTERISTICS

The DP N acceleration equations are used as a low-order operator for the acceleration of the MOC. [START_REF] Sanchez | A Synthetic acceleration for a two dimensional characteristic method in unstructured meshes[END_REF][START_REF] Santandrea | Acceleration Techniques for the characteristics method in unstructured meshes[END_REF][START_REF] Santandrea | Analysis and improvements of the DP N acceleration technique for the method of characteristics in unstructured meshes[END_REF] The equations are obtained by the same projective technique that is used to derive the transport MOC equations but over different representation functions. We present here the transport MOC equations and in the next section the low-order DP N equations. Details can be found in [START_REF] Santandrea | Qualification Physique et amélioration du schéma numérique d'une méthode des caractéristiques pour l'équation de transport des neutrons[END_REF][START_REF] Sanchez | Treatment of boundary condition in trajectory based deterministic transport methods[END_REF].

The method of characteristics provides an iterative solution for the discrete ordinate formulation of the transport equation in a domain composed of unstructured homogeneous regions. The MOC equations are based on two approximations. The first consists of introducing a region-based flat spatial representation for the source term :

( ) ( ) ,
,

s i i i q A S φ ⎡ ⎤ ⋅ Σ + ∈ ⎣ ⎦ r ∼ Ω Ω , c e l l , i r (1) where ( ) ( ) { } A A ρ = Ω Ω is a vector of real spherical harmonics, { } , , s i s i diag ρ Σ = Σ
contains for each region the scattering coefficients associated to each spherical harmonic and i S denotes the external source. The components of vector i φ are the angular flux moments ( ) ( )

(4 ) 1 , 4 i d A π φ π = ∫ Ω Ω Ω i ψ (2)
where is the region averaged angular flux. All angular integrals are calculated with the discreteordinates quadrature formula.

( )

i ψ Ω
The second approximation is to use an expansion for the angular flux on region boundaries. For unstructured meshes one uses a set of parallel trajectories for each discrete ordinate direction Ω. The angular flux in direction Ω is then assumed to be constant over the cross sectional area associated with each trajectory,

( ) ( ) ( ) , , ,
t t t z z ψ θ ψ ⊥ Ω = + ∑ r r r r ∼ Ω Ω ⊥ Ω , (3)
and it is computed analytically along the trajectory. In this formula z and ⊥ r are, respectively, the coordinate in direction Ω and the coordinates over a plane orthogonal to Ω and ( )

t θ ⊥
r is the characteristic function associated to the cross sectional area of trajectory t.

The MOC equations comprise a balance equation and a transmission equation. Only the latter concerns us here. The transmission equation, which results from the analytical integration along each trajectory, gives the angular flux exiting the trajectory in terms of the incoming value and of the internal source. The homogeneity of the region as well as the flat source approximation are used to obtain:

( ) ( ) ( ) ( ) ( 1) , , ( ) 1 (
) ( ) ( ), i i i i i i R t R t n n n i i i i E R e t e t q ψ ψ Ω -Σ -Σ - + - Σ - = + Σ ( 4
)
where n is the iteration index, R i (t) is the length of trajectory t within region i and Ω its direction.

A transport iteration consists of a sweep over all the trajectories crossing the domain. This provides the fluxes leaving the domain and, via use of the region balance equations, the updated values of the region averaged angular fluxes. The latter are then used in equation ( 2) to update the region averaged flux moments that are compared to the previous ones to check convergence. If the iteration has not converged, the sources are updated via Eq. ( 1) and the transport iteration is repeated. Details of the full MOC equations can be found in [START_REF] Halsall | CACTUS, a Characteristics Solution to the Neutron Transport Equation in Complicated Geometries[END_REF][START_REF] Sanchez | Treatment of boundary condition in trajectory based deterministic transport methods[END_REF].

THE INTEGRAL DP N SYNTHETIC OPERATOR

A technique similar to the MOC is used to derive the DP N acceleration equations but with approximations (1) and ( 3) replaced by low-order ones. The source region expansion is as in (1) but with a smaller number of angular moments. The second approximation is more severe. The region boundary is decomposed into surfaces, a surface α being typically the common interface with a neighboring region, and a spatially constant angular approximation is assumed for the boundary fluxes:

( )

( , ) , , s A s S α α ψ ψ Ω Ω ± ± ∈ r r ∼ (5)
where S α stands for a generic surface and α ψ ± indicates exiting (+) or entering (-) surface averaged angular moments. The vector contains all the spherical harmonic of polar order k ( )

A Ω N ≤ . The index "N" indicates that a "Double P N " spherical harmonics expansion is used in (5) (one per side), this gives its name to the acceleration. Together with equation ( 5) the following currents are defined:

(2 ) ( ) ( , ) , s S J dS d n A A α ρ ρ α υ π ψ Ω Ω Ω Ω ± ± ± ± ± = ρ υ υ α α ψ = ∑ ∫ ∫ r i (6)
where coefficient A are defined below in [START_REF] Santandrea | Analysis and improvements of the DP N acceleration technique for the method of characteristics in unstructured meshes[END_REF]. We note that for the lowest moment this quantity is exactly the neutron current.

The first equation type of the DP N synthetic method expresses the current exiting a surface α in function of the interface fluxes of neighboring regions. This is done by using the same set of tracking trajectories used for the MOC. In the standard DP N , method approximation (5) is applied to the surfaces of all regions for all entering fluxes. This allows writing a transmission equation that couples all ,α ψ ± quantities. In this work we consider a partition of the system into macro-domains composed of homogeneous regions and we apply [START_REF] Hong | CRX: a code for rectangular and hexagonal lattices based one the method of characteristics[END_REF] to only the surfaces between macro-domains. We use the integral transport equation over a given trajectory "t". The angular flux exiting a region contained in macro-domain I is:

( ) ( , ) , 1 , j ij j i e x i j I j i j q e e e τ τ τ ψ ψ - - + - ∂ - ∈ < = + - Σ ∑ (7)
where τ(i,∂) is the optical length between the boundary of the macro-domain I and the entering boundary of region i along trajectory t, τ j the optical length over a region j, Σ j the total cross section in region j and τ ij the optical distance between the entering point in region i and the exiting point of region j. Here we consider a point as "exiting" or "entering" with respect to the given trajectory. Figure 1 gives a sketch of the situation.

Using the projection defined in (6) over an arbitrary surface, and applying the DP N approximation (1) to the macro-domain interface entering fluxes, gives: , , ,

,

i i I i I J T E q T ρ ρ υ υ ρ υ υ α α β β α υ β υ ψ ψ + + - - + - ∈∂ ∈ = + = ∑ ∑ Eq + (8)
where the source q υ i has been defined implicitly in [START_REF] Sanchez | APOLLO II: A user-oriented, portable, modular code for multigroup transport assembly calculations[END_REF] and E and T are escape and transmission coefficients respectively. The sum over β is done on region surfaces that belong to the surface ∂I of the macro-domain I, and the sign "+/-" indicates that one has to consider the entering/exiting current. Equation [START_REF] Santandrea | Acceleration Techniques for the characteristics method in unstructured meshes[END_REF] gives us an expression for all currents appearing in the system, regardless of their position with respect to macrointerfaces. The only approximation in ( 8) is to assume a reduced DP N flux anisotropy in the entering flux. In order to allow coupling between macro-domains, another DP N expansion is assumed for the currents exiting a macro-domain interface, that is to say:

, J A ρ ρ υ α α α υ ψ + + + = ∑ υ . ( 8a
)
The definitions and properties of the coefficients A,E and T are quite similar to those of the standard DP N method but, in order to stress some important differences, we recall them here. The matrix coefficients A α appearing in ( 6) and (8a), the transmission and escape coefficients in equation ( 8), are:

( ) ( ) (2 ) 1 , 4 A d S d nA A ρυ ρ υ α α π π Ω Ω Ω Ω + + = ⋅ ∫ ∫ ( ) ( ) ( ) ( ) ,( ) ( ) , ( ) 1 , 4 1 ( ), 4
where (β→α) (respectively (i→α)) indicates that the integration domain has to include only the trajectories going from surface β to surface α (respectively from region i to surface α), τ αβ is the optical distance between α and β, τ α ,(i) is the optical distance along direction Ω between the exiting side from region "i" and the surface α, and, finally, E(R i ) is defined in (4) has been used. It is important to note that for consistency all these coefficients are computed using the MOC trajectories for the integral over surface α. Moreover it is easy to see that the coefficients T and E obey symmetries/anti-symmetries and complementary relations. By "complementary" relation we mean a relation that links the two set of coefficients (E and T). The symmetry/anti-symmetry relation for the T coefficients are:

, T T s s T ρυ υρ ρυ υ ρ α β α β β α + - + - + - = = (10)
where s ρ stands for the parity of the harmonic ρ. It is important to stress how (10) can be used to greatly reduce the computational effort necessary to obtain the transmission coefficients. Actually, instead of computing the coefficients T by using [START_REF] Santandrea | Analysis and improvements of the DP N acceleration technique for the method of characteristics in unstructured meshes[END_REF] and sweeping over all angles, we strongly reduce the integration domain thanks to [START_REF] Sanchez | New development in the Interface-Current Module TDT[END_REF], and compute, sweeping the track, the following integrals:

( ) ( ) [ ] [ ] ( ) 0, 0, / 2 1 ˆ, 4 T d S d nA A e ρ αβ τ ρυ υ α β α β α ρ υ π π π + - - + - → ≤ ∈ × = ⋅ ∫ ∫ Ω Ω Ω Ω Ω (11)
where the conditions listed below the integral sign mean that we are considering only the lowest diagonal angular harmonics and the angular directions belonging to the first two octants. Thanks to [START_REF] Sanchez | New development in the Interface-Current Module TDT[END_REF] it is easy to see that:

ˆˆ, T T s s T ρυ ρυ ρυ υ ρ α β α β β α + - + - + - = + (11a)
for ρ<υ, while for ρ>υ the transmission terms are obtained directly from [START_REF] Sanchez | New development in the Interface-Current Module TDT[END_REF]. Calculations made with [START_REF] Santandrea | Qualification Physique et amélioration du schéma numérique d'une méthode des caractéristiques pour l'équation de transport des neutrons[END_REF] imply the use of less than half of the trajectories that would be needed were we to use [START_REF] Santandrea | Analysis and improvements of the DP N acceleration technique for the method of characteristics in unstructured meshes[END_REF]. A set of complementary relations are also valid for the escape coefficients:

( ) , , [ i I I i i E T T T A ρυ ρυ ρυ ρυ ρυ β β αβ α α β α β α β β δ δ δ Σ + + + + - + - ∉∂ ∈∂ ∈ = - - + ∑ ] . α + (12)
Equation ( 12) expresses the escape coefficients from a region "i" through a surface α as a sum over all the surfaces of the region. Three terms inside the sum must be distinguished. The first one, in which appears the symbol I β δ ∉∂ , comes from the surfaces that are completely internal to the macro-region I. The second one has to be included when the boundary surface β belongs to the border of the macro-region I. The last one has to be considered when β coincides with α. Expression [START_REF] Sanchez | Treatment of boundary condition in trajectory based deterministic transport methods[END_REF] gives the escape coefficients only in terms of the transmission ones. Therefore, the transmission coefficients are the only quantities directly computed from tracking.

Equation ( 8) must be completed to include the boundary conditions, that is to say, by expressing the entering boundary currents in terms of the exiting ones. A detailed description of these equations can be found in [START_REF] Santandrea | Qualification Physique et amélioration du schéma numérique d'une méthode des caractéristiques pour l'équation de transport des neutrons[END_REF]. In this paper we introduce an important improvement for the boundary currents equations that is detailed in appendix A. This treatment is common to both the standard DP N and IDP N methods.

Since the source q contains flux moments, supplementary equations are needed to close the problem.

These are obtained from the DP N region balance relations [START_REF] Santandrea | Analysis and improvements of the DP N acceleration technique for the method of characteristics in unstructured meshes[END_REF]:

( ) , . i i i i i i J J V V B υ υ α α q ρ ρ υ υ α υ υ φ + - - + = Σ Σ ∑ ∑ (13)
Here B is an operator that takes into account the lack of numerical orthogonality for the spherical harmonic functions and is defined as:

( ) ( ) ( ) 4 1 Ω . 4 B d A A ρυ ρ υ π π = Ω ∫ Ω (13a)
Equation ( 13) becomes ill-conditioned for voided or weakly absorbent regions and is therefore cast into a more stable form. The algebraic manipulations are parallel to the ones used in the standard DP N method [START_REF] Santandrea | Analysis and improvements of the DP N acceleration technique for the method of characteristics in unstructured meshes[END_REF], but many substantial differences appear due to the full coupling of the fluxes for regions belonging to the same macro-region. Another benefit of the procedure is to reduce the size of the problem to be solved iteratively. In fact, the iteration problem then deals only with the surface currents.

Substituting the expression (8) of the currents into formula ( 13), one obtains the following system per macro-region I :

,

I I I I I e x t M I Gq φ ψ - = + (14)
where the expressions for matrices I, G, M are given in appendix B. It is important to note that the coefficients of the "incoming" matrix I can be directly obtained from the escape coefficients E given by ( 12), as detailed in equation (B2). This relation allows further computational savings. The matrix G of ( 14) is defined in (B3a) as

, BV G= -E Σ
where the operator E can be defined as a function of the transmission factors as showed in (B5), and B is the matrix (13a). An analytical definition of the operator G is given in appendix B (formulas (B7) and (B7a)), and there we show also that this operator is well defined also for nearly voided media. Moreover we also prove that:

ij ij ji G G s G s ρυ υ ρ ρ ρ υ = = υ . ( 15
)
Equation ( 15) enables us to reduce the computational effort by adopting a strategy similar to the one sketched in [START_REF] Santandrea | Qualification Physique et amélioration du schéma numérique d'une méthode des caractéristiques pour l'équation de transport des neutrons[END_REF] and (11a) for the transmission factors.

The numerical counterpart of (B7) and (B7a) in the framework of the MOC is straightforward, and it is the way we have chosen to compute the operator E. The other possibility to compute E I , the escape matrix for a given macro-region I, is described in appendix B, and makes use of reciprocity relationships between E I and the transmission factors between all surfaces of I. We have opted for the first way since the storage required by E I is by far lower than the one required by the transmission factors. System ( 14) is then used to express the flux moments as:

ˆ, ext I Gq

φ ψ - = + (14a)
where

1 Î M I - =
and . These expressions are inserted back into the transmission equation [START_REF] Santandrea | Acceleration Techniques for the characteristics method in unstructured meshes[END_REF] through the source term (1). This gives:

1 Ĝ M G - = ( ) ( ) 1 1 ŝ s d ext A T E I A E G I q ψ ψ + - - - = + Σ + Σ + , ( 16
)
where

{ } , , s i s i diag ρ Σ =
Σ and I d is the identity matrix. It is worth noting that the matrix T+EΣ s Î, that appears in ( 16), is sign-symmetric (as shown in appendix B), a property that is exploited to reduce the computational effort. Equation ( 16) is the final iterative form used to solve the problem and, since it is written for currents between macro-regions, it has a reduced dimension with respect to the standard DP N method. This means that problem (16) should be faster to solve than its equivalent in the standard DP N method. Once the currents have been obtained by solving (16), volume fluxes are obtained using (14a).

All transmission and escape coefficients can be shown to be bounded. For the transmission coefficients, this fact is clear from their definition [START_REF] Santandrea | Analysis and improvements of the DP N acceleration technique for the method of characteristics in unstructured meshes[END_REF]. On the contrary, to show boundedness for the escape coefficients, it is necessary to develop and implement numerically an asymptotic expansion of equation [START_REF] Sanchez | Treatment of boundary condition in trajectory based deterministic transport methods[END_REF].

In fact, it is clear that expression [START_REF] Sanchez | Treatment of boundary condition in trajectory based deterministic transport methods[END_REF] is greatly ill-conditioned for nearly voided systems. This development is given in appendix B.

Problem (16) has the same response-matrix form as the standard DP N operator. In [START_REF] Santandrea | Analysis and improvements of the DP N acceleration technique for the method of characteristics in unstructured meshes[END_REF] it has been shown that a numerous set of preconditioning and optimization techniques can be used for this problem. This is an important point for our acceleration technique since all previous work, done in the framework of the standard DP N operator, can be re-used for the IDP N acceleration.

All the geometry partitions that we have used in this article have been obtained by the Spectral Bisection Algorithm (SBA) [START_REF] Saad | Iterative Methods for Sparse Linear Systems[END_REF], that is a standard way to compute domain decompositions.

EXTERNAL SYNTHETIC ACCELERATIONS

In this work we have also applied the new IDP N acceleration as a multigroup acceleration technique to improve the rate of convergence in multigroup calculations. Until now, most of the work that has been done to accelerate the external iterations used some non-linear strategy, [START_REF] Zim | Fundamentals and recent developments of reactor physics methods[END_REF] some classical linear extrapolations (like the Chebyshev one) or some Krylov subspace techniques. [START_REF] Warsa | Krylov Subspace Iterations for Deterministic k-Eigenvalue Calculations[END_REF] Recently, Suslov [START_REF] Suslov | An algebraic collapsing acceleration method for the acceleration of the inner (scattering) iterations in long characteristics transport theory[END_REF] proposed a linear algorithm using a mono-group synthetic operator in the external framework. We implemented this idea in APOLLO2 with slight modifications, and used the standard DP N and the new IDP N operator to test its efficiency. In this section we describe the algorithm and discuss implementation details.

The multigroup external equations are driven by fission and, written with the iteration index, read :

( ) ( ) ( ) ( , ,1 ( 1)
1

p j p j p d u p L H H F ψ ψ ψ λ - - - = + ) 1 -, ( 17
)
where L is the one-group transport operator, H d is the group-to-group and slowing down scattering operator, H u is the up-scattering operator and F is the fission operator. The upper index p denotes the external iteration index, while j is the thermal iteration index. Denoting by J tot the number of thermal iterations per external sweep, it is assumed in (17) that ( ) ( )

, tot p pJ ψ ψ =
for every external iteration p. One can then write the converged solution as

( ) ( ) ( ) , p j p j f ψ ψ ∞ = + , , (18)
where

( ) ( ) ( ) F L H ψ ψ λ ∞ - = ∞ . ( 19
)
We substitute formula ( 18) in ( 19), taking into account [START_REF] Zmijarevic | [END_REF], to obtain:

( ) ( ) ( ) ( ) ( ) ( ) ( ) , 1 , , , 1
, given by (20). This last problem is solved using once again the power iteration solution strategy, that is to say by fixing the terms coming from fission and compute a multigroup solution with the help of "j" thermal iterations. The new estimated eigenvalue for the next iteration is given by the usual condition:

p j p p j p j p j s u u p F L H f H f f H F ψ ψ ψ λ λ λ - - ⎛ ⎞ - = + + ∆ + - ⎜ ⎜ ⎝ ⎠ ⎟ ⎟ (20)
( ) ( ) ( ) ( ) ( ) ( ) ( ) , , , , , ,
p j p j k k p j p j k F f F f ψ λ ψ - + == + 1 ,
where the index k refers to the external iterations necessary to solve the problem (20).

RESULTS AND COMPARISONS

In this section we will prove the effectiveness of the IDP N acceleration, presenting some results related to the Atrium Benchmark.

[16] Here we will compute a version of this benchmark that makes use of a sixteen group P 0 transport corrected cross section library that has been obtained by a previous self-shielding phase in APOLLO2. [START_REF] Zmijarevic | [END_REF] The Atrium benchmark is a severe realistic assembly calculation, that has strong transport effects and slow convergence rate. The spatial computational mesh used for our calculation is shown in figure 2, and it consists of 3560 regions. A full description of the case can be found in [16], and we here only insist on some basic features. The Atrium benchmark represents a BWR assembly. The fuel and gadolinium rods are surrounded by an external as well as an internal box, where homogeneous water flows.

As it is typical for these reactors, liquid water concentration varies inside the assembly. In fact, for the considered version of the benchmark, we have a 0% vapor fraction inside both the internal and the external box of the assembly, while we have a 40% vapor fraction for the water in direct contact with fuel rods. Fuel rods are characterized by three different enrichments, and strong flux heterogeneities appear next to gadolinium rods, that are distinguishable in figure 2 by the very fine mesh refinement necessary to represent flux gradients. We used the transport code APOLLO2 where all the developments we present here are now implemented. The tracking data used to compute the benchmark are shown in table 1. We used a product angular formula with sixteen uniform azimutal angles, and two optimal Bickley polar angles. [START_REF] Leonard | Optimal Polar Angles and Weights for the Characteristics Method[END_REF] The transversal spatial step is 0.025 cm. This results in a total of 643,617 horizontal intersections with computational regions.

A first comparison deals with two calculations that have used a DP 1 or IDP 1 multigroup initialization and a thermal synthetic acceleration, a procedure that is also discussed in [START_REF] Santandrea | Analysis and improvements of the DP N acceleration technique for the method of characteristics in unstructured meshes[END_REF]. Table 2 gives the performances of the old DP 1 operator and of the new IDP 1 one. Free iterations are also shown in order to illustrate the efficiency of the acceleration methods. All the results in the table, with the exception of the free iterations ones, were obtained using a mesh grid that was constructed using the recursive SBA. [START_REF] Saad | Iterative Methods for Sparse Linear Systems[END_REF] For the standard DP 1 method, the chosen grid is used to make a two-level multigrid acceleration method [START_REF] Santandrea | Analysis and improvements of the DP N acceleration technique for the method of characteristics in unstructured meshes[END_REF], while for the IDP 1 operator we use this same grid to compute macro-domain integral probabilities. We have always used the adapted ILU0 preconditioning technique [START_REF] Santandrea | Analysis and improvements of the DP N acceleration technique for the method of characteristics in unstructured meshes[END_REF] to speed up the DP N and IDP N equations. Due to their great impact over the total computational time the building times are also shown. All computing times are expressed in seconds on a 2.8Ghz 512Kb cache PC.

Building times strongly depend on the preliminary manipulations required: algebraic elimination process for the standard DP 1 method, and computation of the transmission probabilities for the IDP 1 method.

As the number of regions per macro-region increases, the building time also increases, but the burden of the accelerated iterative resolution strongly decreases. Therefore the appropriate number of regions to be included in each macro-region should be optimized. This phenomenon was already described in [START_REF] Santandrea | Analysis and improvements of the DP N acceleration technique for the method of characteristics in unstructured meshes[END_REF] and we

will not insist on it here. Our numerical tests show that a good number of regions per macro one for the Atrium benchmark is 70, so we always retained the grid mesh resulting from the recursive SBA with this fixed maximum number of regions. It is important to stress that, while for the standard DP N method the grid mesh has no physical impact, for the IDP N it should be beneficial to take advantage of a "physical" grid, that is a grid with reduced flux anisotropy on grid interfaces. This is the case for example for the classical IC method where the grid is given by the region decomposition of the assembly. In fact one can suppose that flux anisotropy is reduced at region interfaces. In this work we will always use the "numerical" grid given by the SBA, even if it certainly does not guarantee the best result.

As one can see from table 2, the best time with consistent tracking, 26.4 sec., is obtained for the IDP 1 operator. The standard DP 1 operator spends a lot of time in the initial building phase, and requires 151 internal transport iterations. In spite of the non-"physical" SBA grid mesh, the spectral behavior of the IDP 1 operator is quite better since it converges with only 102 iterations. The thermal synthetic acceleration, used with the IDP 1 method, reduces the number of total transport iterations by almost a factor of 40 with respect to free iterations. Unfortunately, the IDP 1 acceleration (as well as the DP 1 one) requires an expensive building phase that accounts for half of the total time; this explains why the accelerated calculation is only 15 times faster than the free iterations.

In order to reduce the expensive building phase, we introduced a low-order tracking described at the third row of table 1. The low-order tracking is now used to compute a synthetic IDP 1 operator: this choice entails a lower building cost but can cause a poorer convergence rate. If one compares the data for the reference and low-order tracking, one can see that the number of tracks that are used to build the IDP 1

operator is now decreased by a factor of four with respect to the reference case. We will call this approach "inconsistent" since it is generally required that the acceleration operator be obtained in a "consistent" way with respect to transport, in order to avoid problems of poor convergence or even of divergence of the accelerated algorithm. Using an inconsistent acceleration can therefore lead to lack of robustness. At this time we do not suggest a method to prevent these drawbacks. In the results we show, and more generally on the basis of our experience, we can state that the behavior of the acceleration is not degraded too much if one maintains the low-order tracking formula at an acceptable level.

The new IDP 1 results are given in fourth row of table 2. The use of a low-order tracking should entail poorer IDP 1 performances. The results in table 2 show that this effect is very limited, since the total number of iterations goes from 102 to 107. On the contrary, thanks to the faster building time, the total computational time decreases from 26.4 to 23.6 seconds. Note that while most of the IDP 1 building time is directly proportional to tracking, the standard DP 1 building time strongly depends on the geometrical structure of the problem, and weakly on tracking. The use of a low-order tracking has practically no effect on the standard DP 1 operator, and therefore we do not present inconsistent DP 1 results in table 2.

A last batch of results (table 3) deals with the impact of the external synthetic acceleration. We used this algorithm in conjunction with the standard DP N operator or the new IDP N operator. Both have been computed inconsistently with respect to transport: we used the low-order tracking described by the third row of table 1. In table 3, one can see that the external acceleration is very effective for the synthetic IDP 1 operator, since it reduces the total number of internal transport iterations from 106 to 48. Unfortunately this yields only a slight time reduction, the total computational time goes from 23 to 22.3 seconds. This is due to the fact that the synthetic external acceleration "transfers" a great amount of work from the transport to the synthetic operator. And in fact the total number of synthetic operator problems to be solved is multiplied by more than a factor two, with respect to the cases presented in table 2.

The only possibility to further reduce the computational time is to solve the IDP N problem faster, without increasing building times. We suppose that this could be possible by adopting a full multilevel algorithm (and not a purely two-level one as in this paper). We will address this issue in future works.

CONCLUSIONS

In this paper a new DP N synthetic operator has been introduced. Its main approximation consists in supposing reduced flux anisotropy over a limited number of surfaces. In fact, after having subdivided the system into some macro-regions, the IDP N operator supposes limited flux anisotropy only over those surfaces that divide macro-regions. Thanks to this preliminary hypothesis, and by making use of the basic MOC equations, one can find an integral operator whose unknowns are the flux components over macrointerfaces. Many symmetries/anti-symmetries and complementary relations help to minimize the computational effort, and taking into account important physical properties in the numerical operator.

Calculations made with APOLLO2 on the Atrium Benchmark show that the IDP 1 method is more efficient than its standard counterpart, and that it can be retained as an efficient acceleration tool. The preliminary time necessary to build the operator is reduced and the spectral performance is improved. The behavior of the new IDP N operator in conjunction with the external synthetic acceleration approach has also been evaluated. Even if the overall computational time is reduced, further improvements of the IDP N operator are necessary to make this reduction significant. In particular, we suggest adopting a true multilevel algorithm that will help to reduce building and solution times. We will address this issue in future works.

APPENDIX A: BOUNDARY CONDITIONS TREATMENT

In this section we give a detailed description of the boundary current treatment. The complete form of the boundary current equations has been given elsewhere. [START_REF] Santandrea | Acceleration Techniques for the characteristics method in unstructured meshes[END_REF] What we describe here is an improvement in the treatment of these equations for the DP N problem. We shall show that for all local boundary conditions (i.e. all boundary conditions where the neutrons exiting the system re-enter through the same surface), exiting currents can be eliminated from the iterative part of the problem with almost no computational and memory cost. Let us write the transmission problem (8) in the following manner:

i i T T E S ψ ψ ψ ψ βψ + - - ∂ ∂ - + ∂ ∂ ⎧ = + + ⎪ ⎨ = ⎪ ⎩ , ( A1
)
where we have distinguished the transmission contributions coming from internal currents "ψ -" and from external currents "ψ - ∂ " . Each entering external current can be written as

( ) t α α ψ βψ -= + S . (A2)
For local boundary conditions t(α)=α, but in general surfaces do not coincide. For both cases equation (A1) can be re-written as

i i T T E ψ ψ ψ + - + ∂ ∂ = + + (A3)
where T T β ∂ ∂ =

.

The system of equations (A3) is a reduced form of the initial system (A1): it contains less unknowns and equations but shares the same profile for the remaining equations. This means that contrarily to standard elimination processes, where fill-in always appears, here no elements are added to the initial operator. In the case of no local boundary conditions, the further elimination of the boundary exiting currents appearing in (A3) would cause fill-in and, for this reason, we decided not to eliminate them. For local boundary conditions the situation is different. In fact one can extract the boundary external currents from (A3), and write them as

( ) ( ) ( , ) ( ) ( ) , i i ij j i j i T T E , q ψ ψ ψ ∂ + - ∂ ∂ + ≠ = + + ∑ , (A4) where ( ) ,
i ψ ∂ + is the outgoing external modes related to the i-th external surface, T (i) is the transmission matrix giving the contribution of internal currents to the i-th external boundary current, is the transmission matrix describing the link between the j-th external surface and the i-th one and describes the escape.

It is clear that the only case when is different from zero is when the i-th and the j-th surface share the same region. For convenience we will suppress from now on the tilde notation over the T ( , )

i j T ∂ ( ) i E q ( , ) i j T ∂ ∂ symbol since
each external transmission term will always be post-multiplied by the albedo operator. The objective is to show that building the formal inversion

,

i j i j I T - ∂ ⎛ ⎞ - ⎜ ⎝ ⎠
∑ ⎟ is not computationally expensive and can be done without any fill-in in the LU decomposition necessary to compute this inverse. To see this, let us eliminate the first external current from the system (A4), to

i i i i i j ij j i i j i j D I T T T T T T T T E T E q ψ ψ ψ ∂ ∂ ∂+ ∂ - ∂ ∂ ∂ ∂+ ∂ ≠ > ⎡ ⎤ ⎡ ⎤ ⎡ - = + + + + + ⎣ ⎦ ⎣ ⎦ ⎣ ⎤ ⎦ ∑ , obtain ( ) 1 ( ,1) (1, ) ( ) ( ) ( ,1) (1) ( ,1) (1, ) ( , ) ( ) ( ) ( ,1) (1) , , 1 i i
which is valid for all i>1. Rearranging these equations to a unit diagonal gives:

( ) ( ) ( , ) ( ) ( ) , 1 , 1 , 1 1
,

i i i j j i j i j T T E q 1 i ψ ψ ψ ∂ + - ∂ ∂ + ≠ > = + + ∀ ∑ > , (A6) where , ( ) ( ) ( ) ( ) ( ) ( ) 1 ,1 1 1 1
i i i i T D T T T - ∂ = + ( ) ( ) ( ) ( ) ( ) ( ) 1 , , 1 1 , 1 ,1 i
j i j i j i T D T T T - ∂ ∂ ∂ = + , ∂ and
( ) ( ) ( ) ( ) ( ) ( ) 1 ,1 1 1 1 i i i i E D E T E - ∂ = + .
System (A6) is formally similar to system (A4), with the exception that the number of unknowns is reduced by one in all equations but the first, and that the profile of the operator has been correspondingly reduced since each coupling term due to the first external current has been eliminated from all equations. A similar reasoning can be applied to all the other currents, so that at the final step one obtains an upper diagonal operator. Before the l-th step one has to eliminate the l-th current from the equations greater than l.

After the elimination of this current from the i-th equation of the system (where clearly i>l) this equation reads ( )

( ) ( , ) ( ) ( ) , , ,
i i i j j l l j i j l T T E i l q ψ ψ ψ ∂ + - ∂ ∂ + ≠ > = + + ∑ (A7)
where, after defining

( ) ( ) , , i l l i l i D I T T ∂ ∂ = - , one has ( ) ( ) ( ) ( ) ( ) ( 1 , 1 , 1 1
i i i l l i l l l T D T T T - - ∂ - - = + ) l l , 1 l ∂ - , (A8a) ( ) ( ) 1 , ( , )
,

i j i j l l i T D T - ∂ = . (A8b)
Equations (A8a-b) allow the calculation of all matrices

( ) i l
T from all matrices at the preceding step ( ( )

1 i l T -).
As before, the elimination process gives rise only to the cancellation of boundary terms, without any fill-in. This back-substitution process can be easily repeated, until each external current is expressed in function only of internal currents and collision sources. In fact when treating the j-th equation in the backward procedure, one can use the following situation:

( ) ( ) ( , ) ( ) ( ) , , , ( ) ( ) ( ) ,
, ,

j j jk k j j i k j i i i i back i back T T E j q for i j T E q ψ ψ ψ ψ ψ ∂ + - ∂ ∂ + > ∂ + - = + + > = + ∑ (A9)
where one can note that all equations beyond the j-th one have already been diagonalised. The backsubstitution of these equations into the j-th one gives:

( ) ( ) , , ( ) ( ) ( , ) ( ) ( ) ( , ) ( ) , , , , , j
j j back j back j j j i i j j i i j i i b a c k j i i b a c k i j i j T E T T T E T E Vq ψ ψ ∂ + ∂ - ∂ > > ⎡ ⎤ ⎡ = + + + ⎢ ⎥ ⎢ ⎣ ⎦ ⎣ ∑ ∑ ⎤ ⎥ ⎦
, which gives the recursive definition of the back matrices. Finally the system (A4) is reduced to a system of the form:

( ) ( ) ( ) , i i i back back T E q ψ ψ ∂ + - = + , ( A10
)
where all external elements have disappeared. Eqs (A10) are substituted into (A4) to completely eliminate the boundary currents from the problem and finally obtain:

, ( ( ) ( i back back T T T T E E ) ) q ψ ψ + +∂ - + ∂ ∂ - ∂ ∂ ∂ = + + + ∑ ∑ ,
where the sums over is done on boundary currents. The balance equations have also to be consistently modified. In fact the original balance equations (15a), written without the hat sign, can be split in the following manner:

∂ H GVq H H GV φ ψ βψ ψ - + - ∂ ∂ = + = + + q ) , ( A11
)
where H ∂ stands for the boundary part of the balance operator H. By making use of (A10) one obtains:

( ) ( ˆĤ T H G H E Vq φ β ψ β - ∂ ∂ = + + + , if , if . j j j I j I j j I j I T T E E q J J T A E q I ρυ ρυ υ ρυ ρυ υ α β α β β α α β υ υ ρ ρ α α ρυ υ ρυ υ ρυ υ α β β α α α β υ υ υ ψ α ψ ψ α + - - - - + - + - + - - - - + ∈∂ ∈ ∈∂ ∈ ⎧ - + - ⎪ - = ⎨ - + ∈∂ ⎪ ⎩ ∑ ∑ ∑ ∑ ∑ I ∉∂ (B1)
From equation (B1) we derive the following expression for the leakage term of the balance equation:

( ) ( ) ( ) ( ) ( ) , , , , , , , , , , , 1 1
I i I ij I I i I i i i I I I j j j i i I i E J J T T T A E E E q ρυ β ρυ ρ ρ α α ρυ ρ υ ρ υ ρ υ α α α β υ α β α β α β α β α β υ α ρυ ρυ ρυ α α α α α α δ δ δ δ δ ψ + - + - - + - + + - ∉∂ ∈∂ ∈ ∈ ∂ ∈ ∉∂ ∈∂ ∈ ∈ - ⎧ ⎫ ⎡ ⎤ = - + - ⎨ ⎬ ⎣ ⎦ Σ Σ ⎩ ⎭ ⎧ ⎫ ⎡ ⎤ + - + ⎨ ⎬ ⎣ ⎦ Σ ⎩ ⎭ ∑ ∑ ∑ ∑ , . j j I υ υ ∈ ∑ - (B2)
Equation (B2) implicitly defines the matrix I of equation [START_REF] Suslov | An algebraic collapsing acceleration method for the acceleration of the inner (scattering) iterations in long characteristics transport theory[END_REF], and another auxiliary matrix E I used to construct the matrices G and M of equation [START_REF] Warsa | Krylov Subspace Iterations for Deterministic k-Eigenvalue Calculations[END_REF]. The incoming matrix I can be better re-written in function of the escape coefficients [START_REF] Sanchez | Treatment of boundary condition in trajectory based deterministic transport methods[END_REF]:

( ) ( ) ( ) ( ) , , , , , , , 1 1 ,
I I i i i I I i i i I T T T A s s T T T A s s E ρυ ρυ ρυ ρυ ρυ α α α β β α β α β α β α α , ρ υ ρυ ρυ ρυ ρυ ρ υ ρυ α α α β β α β α β α α β α δ δ δ δ δ δ - + - - - + - - + - + + + - + + ∉∂ ∈∂ ∈ ∉∂ ∈∂ ∈ ⎡ ⎤ = - - + - = ⎣ ⎦ Σ ⎡ ⎤ - - + - = ⎣ ⎦ Σ ∑ ∑ - (B3)
where [START_REF] Sanchez | New development in the Interface-Current Module TDT[END_REF] has been used together with the fact that A s s A 13), one obtains the following expression for the matrix M I of equation ( 15) :

( ) ( ) , [ j ] I i ij I ij s j ij ij i B V M V E ρυ ρυ ρυ υ ρυ δ δ δ = + -Σ Σ , (B4)
where i,jOEI. Finally, the matrix G that appears in ( 15) is given by: ,

I I i BV G = - Σ E (B4a)
where E I is defined in (B1) and B in [START_REF] Saad | Iterative Methods for Sparse Linear Systems[END_REF]. Equations (B3) and ( 10) allow us to show that the matrix T+EΣ s Î appearing in ( 16) is sign-symmetric that is to say that it satisfies the relation ( ) ( )

ˆt s s T E I s T E I s + Σ = + Σ .
Indeed, first note that (B3) and ( 10) can be written in matrix form as t and . Moreover let us suppose that (and then also

T sTs = t E sI = -s 1 I I M - Σ 1 I I M - Σ ) is sign-symmetric. This implies that ( ) ( 1 t s I s I T E M I s T E M I s - + Σ = + Σ ) 1 - , (B5)
which is the expected result. We still need to demonstrate Interchanging i with j in the left hand side amounts to interchanging α with β in the right hand side.

Using [START_REF] Sanchez | New development in the Interface-Current Module TDT[END_REF], one then sees that each term in (B6) is transformed either into itself multiplied by s ρ s υ , as it is the case for the terms into the third line of (B6), either in other terms of the same row multiplied by s ρ s υ , as it is the case for the terms in the first two lines. We will find this result in another way in formula (B8). Equations (B4) and (B8) can be used to reduce the computational effort as well as the memory storage, since they show that only halves of the matrices E I and T+EΣ s Î have to be computed and stored.

Thanks to (B6) an analytical definition of operator G can be given:

( ) ( ) ( , ) ( ) ( , ) 1
( ) ( ) 4 lm s i i s j ij ij s l m l T r V j i V m T r V G E d dr n A A E e E i τ ρυ ρυ ρ υ τ τ π - ∈ Ω ∩ → ∂ ∈ Ω ∩ = - = Ω Ω ⋅ Ω Ω ≠ ∑ ∫ ∫ j.(B7)
In equation (B7) the symbol (j→i) stands for the angular directions going from region j to region i, the symbol T(r,Ω) stands for the set of all the chords intersected by the trajectories that pass through the point r with angle Ω, so that T(r,Ω)∩V i (T(r,Ω)∩V j ) stands for the subset of chords that belong to region i (or region j) and symbol E is defined in (4). Therefore the sum in (B7) takes into accounts the case of the intersections with non-convex regions. The optical chord lengths τ l ,τ m and τ lm are related respectively to the chord in region i, to the chord in region j and to the distance between the exiting point of segment l and the entering point of segment m. An analogous expression can be written for the diagonal terms:

( ) ( ) (2 ) ( , ) , ( , ) 1 4 {
[ ( )]/ ( ) ( ) }, i lm s i s i i ii ii s i V l i l m l T r V l m T r V l m B V G E d r d nA A l E E e E ρυ ρυ ρυ ρ υ π τ π τ τ + ∂ - ∈ Ω ∩ ∈ Ω ∩ ≠ = - = ΩΩ⋅ Ω Ω Σ - Σ + ∫ ∫ ∑ ∑ i τ (B7a)
where the second sum inside parentheses comes again from non-convex regions. Inspecting (B7) and (B7a), one can readily realize the symmetric role of indexes i and j as well as of the angular indexes. This assures us

that ij ij ji G G s G s ρυ υ ρ ρ ρ υ = = υ . (B8)
Next, we show that the escape coefficient for the source in Eq. ( 9) remains bounded at the vacuum limit. From the expression for i E ρυ α we obtain where we used ( )

2 2 0 1 lim 2 i i i R i i i i R e R O -Σ Σ → Σ - → - + Σ Σ
i . Formula (B9) shows that the escape coefficients remain bounded, but unfortunately does not give any complementary relation. This means that, in case of nearly voided media a new integral, given by (B9), must be computed.

Another term that must be treated specifically at the vacuum limit is the diagonal term of operator G in equation (B3a). One can show that this diagonal term has the following limit that is to be computed explicitly for the nearly voided cases.

Figure 1: Sketch for the transport equation [START_REF] Sanchez | A Synthetic acceleration for a two dimensional characteristic method in unstructured meshes[END_REF]. The macro-region "I" contains many micro regions "r". For a trajectory along direction Ω, many intersections have to be considered from the entering to the exiting point. In order to express each exiting angular flux along this trajectory, the exponential of the optical length between these regions must be computed. For example, to express thanks to equation [START_REF] Sanchez | A Synthetic acceleration for a two dimensional characteristic method in unstructured meshes[END_REF], a sum over the preceding collision sources (weighted with and ) and over the entering flux must be done. Table 3: Results for the Atrium assembly calculation using the standard DP 1 operator and the new integral IDP 1 operator with double tracking. All cases are DP N -initialized eigenvalue calculations and external iterations are accelerated with the external synthetic algorithm. The computational times are expressed in seconds. The reference k-effective is 1.13811.

Figure 1: Sketch for the transport equation [START_REF] Sanchez | A Synthetic acceleration for a two dimensional characteristic method in unstructured meshes[END_REF]. Macro-region "I" contains many micro regions "r". For a trajectory along direction Ω many intersections have to be considered from the entering to the exiting point.

In order to express each exiting angular flux along this trajectory the exponential of the optical length between these regions must be computed. For example to express thank to equation ( 7) a sum over the preceding collision sources, weighted with and , and over the entering flux must be done.

  between two consecutive transport iterates. Moreover, λ is the new eigenvalue obtained as the result of the multigroup correction problem (20). Problem (20) is an eigenvalue problem but with the presence of an external inhomogeneous term. An accelerated external iteration consists of a free MOC external iteration followed by one synthetic multigroup solution correction

  The remaining step is to complete the inversion of the original by substituting back, from the last until the first, the new computed equations. Let us call N bd the total number of boundary equations in (A4). After the N bd steps described above, one has an upper diagonal operator where the last equation gives the N bd -th external current in function only of internal unknowns and collision sources. The backsubstitution of this equation into the (N bd -1)-th one allows again to express the (N bd -1)-th external current in function only of internal currents and collision sources.

  coefficients of the incoming matrix I and escape coefficients E. Re-inserting (B1) into the balance equation (

FIGURE 2

 2 FIGURE 2

Figure 2 :

 2 Figure 2: Computational mesh for the BWR atrium assembly.

Figure 2 :

 2 Figure 2: Computational mesh for the BWR Atrium assembly.

TABLE 1

 1

	trajectories	tracks	azimuthal angles	polar angles	Transversal
					Integration step
	162	643,617	16 (uniform)	2 (Bickley)	0.025cm
	143	162,610	10 (uniform)	2 (Bickley)	0.05 cm

Table 1 :

 1 Tracking parameters for the calculation of the Atrium assembly with cyclic trajectories. The number of tracks is the total number of intersections of trajectories with regions. The second and the third rows gives, respectively, the data of the reference tracking and of the DP N low-order one.

			TABLE 2
		building time Total time Nb. External it. Nb. Internal it.
	DP 1 Consistent tracking	22.7	41.4	7	151
	IDP 1 Consistent tracking	12.8	26.4	5	102
	IDP 1 Inconsistent tracking	9.44	23.6	5	107
	Free Iterations		412	13	3848

Table 2 :

 2 Results for the Atrium assembly calculation using the standard DP 1 operator and the new integral IDP 1 operator. All cases are DP N -initialized eigenvalue calculations. The computational times are expressed in seconds. The reference k-effective is 1.13811.

TABLE 3

 3

		building time Total time Nb. External it. Nb. Internal it.
	DP 1	19.5 sec	40.3	6	96
	IDP 1	9.4 sec	22.3	3	48
	Free Iterations		412	13	3848

APPENDIX B: MATRIX COEFFICIENTS OF THE INTEGRAL DP N METHOD AND ASYMPTOTIC LIMITS OF ESCAPE PROBABILITIES

In this appendix we give some numerical details about the DP N equations. First of all we detail the expression of matrices M, G and H appearing in [START_REF] Warsa | Krylov Subspace Iterations for Deterministic k-Eigenvalue Calculations[END_REF]. Next we consider the vacuum limit for the escape coefficients in equation [START_REF] Santandrea | Analysis and improvements of the DP N acceleration technique for the method of characteristics in unstructured meshes[END_REF]. Substituting ( 8) into (13) yields: