N S Ugwuanyi

X Kestelyn
email: xavier.kestelyn@ensam.eu

O Thomas

B Marinescu
email: bogdan.marinescu@ls2n.fr

B Wang
email: bin.wang@nrel.gov

A Normal Form-based Power System Out-of-step Protection

Keywords: Normal form, oscillation frequency drift, outof-step prediction, transient stability

This paper proposes a new system-level application for monitoring out-of-step (OOS) events in power systems. As already known, amplitude-dependent frequency shift is a nonlinear phenomenon of electromechanical oscillations under large disturbances. The frequency shift indicates the system's nearness to instability. This new tool utilizes the Normal Form method to identify the named phenomenon, leading to accelerated OOS detection. The proposed strategy is illustrated and compared to the equal-area criterion method in a single-machine-infinitebus power system. Extensive tests on IEEE 3-and IEEE 50machine power systems prove the efficacy and potential of the proposed method for online warnings of instability and ranking of vulnerable system modes.

I. INTRODUCTION

T HE he ultimate goal of any well-designed power system is to recover quickly following large disturbances such as faults, line switching, loss of large generation, e.t.c. Certain disturbances may cause the system to fail, resulting in a loss of synchronism. Out-of-step (OOS) or asynchronous operation refers to the loss of synchronism of one or more generators. OOS protection typically focuses on detecting an asynchronous operation as soon as it occurs. As a result, preventive measures such as islanding the affected generators, load-shedding, etc., are implemented.

Detecting OOS usually entails spotting significant changes in power system parameters at the start of an asynchronous action. These parameters include line impedance, rotor angle, phase, current, voltage, bus voltage frequency, e.t.c. Depending on the quantity tracked, different methods for detecting OOS exist in the literature. The authors in [START_REF] Alinezhad | Out-of-step protection based on equal area criterion[END_REF] extended the equal-area-criterion (EAC), which was previously restricted to two-machine systems, to multi-machine systems to detect OOS. This technique based on analysis of the power-angle characteristics is valuable for oscillations between two groups of generators or between an individual generator and the rest of the system. Ref. [START_REF] Paudyal | Application of equal area criterion conditions in the time domain for out-of-step protection[END_REF] employs it for OOS detection in protection applications based on power-time analysis. To find the cluster solution that gives the position leading to controlled islanding, [START_REF] Ayer | Online application of local oos protection and graph theory for controlled islanding[END_REF] uses graph theory combined with time-domain EAC. The method's strength lies in its use of local OOS protection to detect system-wide instability. Using Clarke's transformation and analytic geometry evaluation, [START_REF] Camarillo-Peñaranda | An approach for out-of-step protection based on swing center voltage estimation and analytic geometry parameters[END_REF] proposes an OOS protection method. It boasts a real-time implementation that combines time-domain analysis in Clarke's domain with the sudden-change-in-voltage (SCV) approach to enhance OOS protective feature action. Zhang et. al, [START_REF] Zhang | A novel out-of-step splitting protection based on the wide area information[END_REF] convert the bus voltage phase angles obtained from phasor measurement units into a new coordinate system in which all of the bus voltage phasor ends are stable during the out-of-step oscillation process, while only the common beginning moves. The out-of-step center position criterion and the bus grouping criterion are then proposed. This method, like the wide-area measurement (WAM) methods in [START_REF] Guo | Online identification of power system dynamic signature using pmu measurements and data mining[END_REF], [START_REF] Sauhats | Wide-area measurements-based out-of-step protection system[END_REF], depends on the accuracy of the measurement units and the fidelity of widearea communications to fully use it. Other existing approaches include variation in voltage frequency-based method [START_REF] Zhang | Characteristic analysis and calculation of frequencies of voltages in out-of-step oscillation power system and a frequency-based out-of-step protection[END_REF], line differential current-based method [START_REF] Hashemi | Current-based out-of-step detection method to enhance line differential protection[END_REF], instantaneous active power deviation-based method [START_REF] Rajabinasab | A real-time out-of-step protection strategy based on instantaneous active power deviation[END_REF], adaptive method [START_REF] Ariff | Adaptive protection and control in the power system for wide-area blackout prevention[END_REF], decision tree method [START_REF] Amraee | Transient instability prediction using decision tree technique[END_REF] and turbine generator behaviorbased method [START_REF] Abedini | Generator out-of-step prediction based on faster-than-realtime analysis: Concepts and applications[END_REF]. The majority of those references use EAC to supplement their approaches (e.g., [START_REF] Alinezhad | Out-of-step protection based on equal area criterion[END_REF]- [START_REF] Ayer | Online application of local oos protection and graph theory for controlled islanding[END_REF], [START_REF] Ariff | Adaptive protection and control in the power system for wide-area blackout prevention[END_REF]), so it appears to be the most commonly published.

Most of the previous methodologies concentrate on locating the OOS incident and separating the affected generators, which means that OOS must have occurred. As a result, most of them need the fault clearing time fed as an input to the relay. They include details in the form of a switch's-YES or NO-much like time-domain simulations. They do not give sufficient information on how close the system is to going down. Some OOS symptoms are hidden but dangerous, manifesting fully solely when the only available defense may be islanding. Furthermore, most of the methods examined are local, making them incapable of detecting system-wide instability. Any method that detects the power system's gradual descent into the OOS regime in a system-wide manner will undoubtedly improve its protection. Whether current, voltage, impedance or other parameters are monitored, they all lead to one thing: controlling the natural mode of the power system. The dynamics of a power system are related to its natural mode, in addition to the subjected disturbance.

In this paper, we propose to monitor the variations in the power system's electromechanical modal parameters (in this case, oscillation frequency) prior to the commencement of the asynchronous regime to detect and mitigate the OOS. Our attention is the electromechanical dynamics, which is the principal stability effect in many systems. There are practical reasons for this proposal. As examples, the Western Interconnection power grid experienced a massive breakup on August 10, 1996. According to measurement-based research, the 0.25-Hz oscillation, which was suggestive of an imminent breakup, resulted from a gradual decrease in oscillation frequency from 0.27 Hz to 0.25 Hz, as well as a decrease in oscillation damping ratio from 7.0% to 1.2% [START_REF] Hauer | A perspective on wams analysis tools for tracking of oscillatory dynamics[END_REF]. On June 4, 2003, there was another widespread oscillation resulting from weakened system conditions. Prony analysis determined the oscillation frequency to be 0.584 Hz and the mode shape to be an east-west interaction between the Colstrip and Kemano plants, although the dominant modes for these plants are typically near 0.78 Hz and 0.63 Hz, respectively [START_REF] Hauer | A perspective on wams analysis tools for tracking of oscillatory dynamics[END_REF]. The preceding real events indicate that the reported instabilities were accompanied by a decrease in modal frequency and damping, which can help detect the onset of an OOS. When the initial conditions of a power system change, the frequency of electromechanical oscillations (EOs) changes. In a linear case, the new frequency is constant and may be viewed as a new fundamental frequency at the new initial condition. In a nonlinear scenario, the EO's frequency does not keep constant but varies with the amplitude, a phenomenon often known as a frequency shift. In this paper, we use frequency shift for the two cases and view the former as a shift for a constant amplitude. This phenomenon has piqued the interest of researchers, notably the Pacific Northwest National Laboratory, which is currently developing tools for real-time monitoring of electromechanical modes' parameters, a technology Hauer et. al, [START_REF] Hauer | A perspective on wams analysis tools for tracking of oscillatory dynamics[END_REF] dubbed "Mode Meter." This paper employs the Normal Form (NF) method to detect the dominant mode's frequency shift towards the asynchronous regime. The efficient detection of oscillation frequency changes will supplement existing methods and increase the power system's reliability. The strength of the proposed method lies in its ability to detect system-wide instability speedily and easily. It is so because it is not based on any particular machine, but the prevailing system condition. Controlling the responsible generators is a second step, which can be done with any other method, e.g., participation factor analysis. While identifying the exact time for OOS and the responsible generators is crucial for islanding, the proposed method boasts the capability of indicating imminent OOS before it ever occurs. This added feature is key to preventive controls. In summary, this paper's contributions include (i) a new way to estimate and monitor the oscillation frequency without repeating the small-signal analysis (SSA), (ii) a new application of Normal Form method to power system protection, (iii) the proposed OOS detection algorithm is much fast, system-wide, and has potential for real-time implementation. Moreover, the proposed tool may be coordinated to work parallel with existing tools like WAM, EAC, e.t.c., to enhance system reliability. The rest of the paper is structured as follows: Section II discusses the theoretical foundations of the NF method. Sec-tion III provides the reader with clear motivational examples to help them understand the concept for the proposed technique. Consequently, the new method is presented. Section IV is devoted to in-depth simulations and discussions. Finally, Section V concludes the paper and suggests future research.

II. NORMAL FORM THEORIES

The Poincaré and Poincé-Dulac NF theory [START_REF] Poincaré | Les mthodes nouvelles de la mcanique cleste[END_REF] uses a polynomial transformation to simplify nonlinear differential equations. The NF method requires, to begin with, the representation of the nonlinear system with a polynomial of appropriate degree followed by a transformation that removes or simplifies the nonlinear couplings. Ref. [START_REF] Wang | Power system transient stability analysis using high-order taylor expansion systems[END_REF] conducted detailed studies that revealed that a polynomial of degree 3 substantially retains the power system's nonlinear behavior in some vicinity of the stable equilibrium point. The NF method is particularly useful in power systems for studying the effects of modal interactions in control designs [START_REF] Lomei | A New Approach to Reduce the Non-Linear Characteristics of a Stressed Power System by Using the Normal Form Technique in the Control Design of the Excitation System[END_REF], [START_REF] Liu | Assessing placement of controllers and nonlinear behavior using normal form analysis[END_REF] and power system stability [START_REF] Wang | Nonlinear Modal Decoupling of Multi-Oscillator Systems With Applications to Power Systems[END_REF]- [START_REF] Amano | Proposal of nonlinear stability indices of power swing oscillation in a multimachine power system[END_REF]. This technique has also been used in vibration mechanics to investigate the frequency shift phenomenon of a natural mode subjected to a displacement [START_REF] Touzé | Hardening/softening behaviour in non-linear oscillations of structural systems using non-linear normal modes[END_REF]. The frequency shift phenomenon occurs in power systems as well, but it is not given the attention it deserves. These critical applications sparked a slew of studies aimed at easing NF's computational burden [START_REF] Ugwuanyi | Power System Nonlinear Modal Analysis Using Computationally Reduced Normal Form Method[END_REF]- [START_REF] Hernandez | Recursive linearization of higherorder for power system models[END_REF]. In the following subsections, the NF theory for first and second-order systems is reviewed.

A. Normal Form of First-order System

Consider a dynamical system:

ẋ = f (x), (1)
where x is the N -dimensional state vector, f is a smooth vector field. The system in (1) is assumed to have a stable equilibrium point at the origin, i.e. f (0) = 0, (otherwise, a coordinate transformation can move it to the origin) and can be approximated by Taylor series around the equilibrium point. Therefore, system (1) can be written as:

ẋ = Āx + Γ(x), (2)
where Ā is the Jacobian of f , Γ(x) represents nonlinear terms expanded in polynomial series and Γ(0) = 0.

Assuming that Ā is diagonalizable, a linear transformation x = Φy decouples the linear part of (2) as

ẏ = Λy + f nl (y), (3)
where Λ is a diagonal matrix of the eigenvalues of Ā and f nl (y) = Φ T Γ(Φy) gathers the nonlinearities in the new coordinates. Φ is the right eigenvector of Ā, normalized such that Φ T Φ = 1. In power systems, Λ has complex elements. NF theorem states that if no resonance exist among the eigenvalues of Ā, there is a nonlinear transformation

y = z + h(z), (4)
which simplifies the system (3) to

ż = Λz. (5)
h(z) is a higher order polynomial. In other words, in [START_REF] Zhang | A novel out-of-step splitting protection based on the wide area information[END_REF], the counterpart of the nonlinear terms in the system (3) is removed. In the case of resonances, all mononomial terms in f nl (y) that are resonant cannot be removed. NF finds h to eliminate the removable nonlinear terms. The NF for the first-order system has a conceptual disadvantage in that the physical meanings of the state variables may be unclear, and using second-order systems can relieve this nuisance.

B. Normal Form of Second-order Systems

The NF derivations for second-order mechanical systems have been extensively discussed in the literature and will not be replicated in this article. For details, see [START_REF] Touzé | Hardening/softening behaviour in non-linear oscillations of structural systems using non-linear normal modes[END_REF] and references therein, however, the basic steps required to understand the method are presented here. This paper will present a novel application of previous NF findings to the study of power system oscillations, inspired by the relationships between second-order mechanical and electrical models.

1) Power System Model: The equation of motion for an N -machine classical model power system can be represented by a 2nd-order differential equation [START_REF] Guo | Online identification of power system dynamic signature using pmu measurements and data mining[END_REF].

M i δi + C i δi + P ei = P mi , (6)
P ei = E 2 i G ii + N k=1,k =i E i E k [G ik cos δ ik + B ik sin δ ik ].
M i , C i , P ei , P mi , δ i are respectively inertia constant, damping constant, electrical power, mechanical power and absolute rotor angle of generator i, and δ ik = δ i -δ k . G ik and B ik are respectively the conductance and suceptance between generator i and k while E i is the constant emf behind the transient reactance of generator i.

2) Taylor Expansion: The 3rd-order Taylor expansion of (6) for a constant input reads:

Mij δj + Cij δj + N j=1 Kijδj + N j=1 N k=1 F 2 i jk δjδ k + . . . N j=1 N k=1 N l=1 F 3 i jkl δjδ k δ l = 0, (7)
where

K ij = ∂Pe i ∂δj | δ=δ0 , F 2 i jk = 1 2 ∂ 2 Pe i ∂δj ∂δ k | δ=δ0 , F 3 i jkl = 1 6 ∂ 3 Pe i ∂δj ∂δ k ∂δ l | δ=δ0 .
The eigenvalues of K are all real. In the conventional NF, the Taylor expansion of (6) has to be re-written in first-order, an operation that doubles the size of the system and leads to the computation of complex quantities.

3) Linear Transformation: Let Ω be a vector of the eigenvalues, obtained by solving

(K -Ω 2 i M)Φ i = 0. Let the i-th element be Ω 2 i .
If uniform damping is assumed, matrices M, C and K can be diagonalized simultaneously, leading to:

     Φ T i M Φj = σij ∀i, j Φ T i KΦj = Ω 2 i σij ∀i, j Φ T i CΦj = 2ζiΩiσij ∀i, j, (8)
where ζ is the damping factor, and σ ij = 1 ∀i = j, 0 ∀i = j. If Φ ip denotes the i th element of Φ p and using one generator as a reference, a linear transformation

δ = Φy =⇒ δ i = Φ ip y p (9) transforms (7), ∀p = 1, . . . N , to ÿp = -Ω 2 p yp - N q=1 N r=1 G p qr yqyr - N q=1 N r=1 N s=1
H p qrs yqyrys, [START_REF] Rajabinasab | A real-time out-of-step protection strategy based on instantaneous active power deviation[END_REF] where G p qr = F 2 i jk Φ ip Φ jq Φ kr , H p qrs = F 3 i jkl Φ ip Φ jq Φ kr Φ ls . y is the modal variable. The damping term has been canceled in [START_REF] Rajabinasab | A real-time out-of-step protection strategy based on instantaneous active power deviation[END_REF] by subtracting the equation corresponding to the reference generator from all other equations.

To avoid complex quantities, ( 10) is put to first-order using the velocity ȳp = ẏp as auxiliary variable and then, each oscillator represented by a block diagonal matrix

0 1 -Ω 2 p 0 .
The auxiliary variable allows [START_REF] Rajabinasab | A real-time out-of-step protection strategy based on instantaneous active power deviation[END_REF] to be written in the form of ( 3) and can always be canceled to recover 2nd-orderlike oscillator equations. The first difference compared to the traditional NF method is that the 2nd-order system model is expanded by the Taylor series before putting it to first-order. The frequency of the undamped linear mode is unchanged while admitting simpler computation.

4) Nonlinear Transformation: To decouple ( 10), the nonlinear change of variables defined in ( 11) is used.

yp = Rp + N i=1 N j≥i (a p ij RiRj + b p ij SiSj) + N i=1 N j≥i N k≥j r p ijk RiRjR k + N i=1 N j≥i N k≥j u p ijk RiSjS k . ( 11a
) ȳp = Sp + N i=1 N j=1 γ p ij RiSj + N i=1 N j≥i N k≥j µ p ijk SiSjS k + N i=1 N j≥i N k≥j v p ijk SiRjR k . (11b)
R p and S p are new variables, respectively corresponding to the angle and speed. The unknown coefficients a p ij , b p ij , γ p ij , u p ijk , r p ijk , v p ijk , µ p ijk (equivalence of h-coefficients in firstorder case), are determined by just solving sets of linear equations formulated by equating like terms. Details of their computations can be found in [START_REF] Touzé | Hardening/softening behaviour in non-linear oscillations of structural systems using non-linear normal modes[END_REF], [START_REF] Touzé | A normal form approach for nonlinear normal modes[END_REF]Appendix C]. Unlike the h-coefficients in (4), these coefficients are all real.

By substituting ( 11) into (10), followed by simplifications, the NF for the p-th mode is given by [START_REF] Touzé | Hardening/softening behaviour in non-linear oscillations of structural systems using non-linear normal modes[END_REF]:

Rp = -Ω 2 p R p -(A p ppp + H p ppp )R 3 p -B p ppp R p Ṙ2 p , (12)
where the coefficients A p ijk and B p ijk are the third-order residual terms expressed as [START_REF] Amraee | Transient instability prediction using decision tree technique[END_REF] assumes that there is no strong modal resonance, thereby neglecting the effects of cross-coupling nonlinearities. A reverse transformation is possible via [START_REF] Ariff | Adaptive protection and control in the power system for wide-area blackout prevention[END_REF] and [START_REF] Hashemi | Current-based out-of-step detection method to enhance line differential protection[END_REF].

A p ijk = N l≥i G p il a l jk + N l≤i G p il a l jk , B p ijk = N l≥i G p il b l jk + N l≤i G p il b l jk . Equation

III. PROPOSED METHOD FOR OUT-OF-STEP TRACKING

The proposed method for OOS protection is described in this section. To begin, the reader is given clear motivational examples to help them understand the key concept that follows.

A. Motivations for the Proposed Technique

Firstly, consider a single-machine-infinite-bus (SMIB) power system whose data and model are described in Section IV-A and equation ( 6) respectively. It can easily be shown that if [START_REF] Guo | Online identification of power system dynamic signature using pmu measurements and data mining[END_REF] is written in first-order and C is zero, its eigenvalues are given by:

λ 1,2 = ω 2 s M ω 2 0 P eo ± j ω 2 s M ω 2 0 P eo 2 -4 ω 2 s M ω0 EV Xs cos δ 0 2 , ( 13
)
where P eo is the electrical power output of the generator, δ 0 , the initial rotor angle, ω s , the synchronous speed, and ω 0 , the initial rotor speed. Note that the simplifying assumption ω s = ω 0 was not used in [START_REF] Abedini | Generator out-of-step prediction based on faster-than-realtime analysis: Concepts and applications[END_REF]. Equation [START_REF] Abedini | Generator out-of-step prediction based on faster-than-realtime analysis: Concepts and applications[END_REF] shows the effects of the generator loading on the oscillation's damping and frequency. As reported in [START_REF] Wang | Formulation and Characterization of Power System Electromechanical Oscillations[END_REF], if the generator's output increases, the frequency and damping decrease, as shown in Fig. 1a. To show this quantitatively, let the equilibrium point of the SMIB system be 15 • . The frequency of the mode is 7.82 rad/s (1.24 Hz). Using the EAC, the stability bound for this operating point is 165 • (2.88 rad). Now, consider an initial condition (δ = 164.90 • , ω = 0), i.e., 99.94% of the stability bound. The nonlinear response of this condition is shown in One can, therefore, conclude that in a stable mode, the frequency decreases as the generator output or fault level increases. This frequency's descent reaches a critical point below which the mode loses stability, and the probably affected generators are those participating significantly in that oscillation mode. This phenomenon motivates the mode's stability monitoring for OOS protection.

B. Out-of-step Detection Criterion

The OOS detection criterion in this paper is quite simple and based on careful exploitation of [START_REF] Amraee | Transient instability prediction using decision tree technique[END_REF]. As seen in Fig. 1b, a nonlinear phenomenon of EOs under large disturbances is amplitude-dependent frequency drift. With the first-order representation, an analytical expression of the frequency shift (or nonlinear frequency) can be found for (12) using the perturbation methods, e.g. the multiple-scale method, as:

Ω p nl = Ω p 1 + 3(A p ppp + H p ppp ) + Ω 2 p B p ppp 8Ω 2 p a 2 , ( 14
)
where Ω p nl is a frequency, which depends on the amplitude a considered, while Ω p is the fundamental frequency. The term in parenthesis is a correction to the fundamental frequency as the amplitude changes. Equation ( 14), also known as the backbone curve in mechanics, has already been derived and used to investigate the hardening and softening behavior of nonlinear structural vibration [START_REF] Touzé | Hardening/softening behaviour in non-linear oscillations of structural systems using non-linear normal modes[END_REF]. To the best of the authors' knowledge, the application of such a concept to the study of multi-machine power system oscillation is original.

Estimating the critical frequency shift necessitates an evaluation of the critical amplitude of the oscillation. Since the rotor angle δ is transformed to NF, the maximum displacement of R = R cr can be obtained from ( 12) by setting R = Ṙ = 0 as:

R pcr = R p = -Ω 2 p (A p ppp + H p ppp ) , (15)
where R pcr stands for the critical displacement of p-th NF variable, R. Since a is the amplitude of the p-th mode, the critical displacement R pcr has to be projected to the modal coordinate using the modal model (11a) as follows:

a pcr = y p (R pcr ) = R pcr + a p pp R 2 pcr + r p ppp R 3 pcr . (16)
Therefore, the critical frequency shift is computed using ( 14) with a = a pcr . With the critical frequency shift computed, the stability is evaluated anytime as follows:

1) Obtain the system condition δ sc , where sc means system condition. 2) Obtain the modal amplitude of the system condition using (9) as a psc = y psc . 3) With the amplitude determined, compute the frequency shift, Ω p nlsc using (14). 4) If Ω p nlsc > Ω p nlcr , the system condition is considered stable; otherwise, it is unstable. The following benefits of the proposed method are noteworthy:

• Since it is based on the oscillation mode, a system-wide parameter, it detects system-wide instability. • It does not need the fault clearing time fed as an input to the relay, so it has a similar advantage advertised in [START_REF] Alinezhad | Out-of-step protection based on equal area criterion[END_REF]. • It requires only computations related to the mode of interest, making it fast and easy, thanks to the selective computation technique promoted in [START_REF] Ugwuanyi | A New Fast Track to Nonlinear Modal Analysis of Power System Using Normal Form[END_REF]. • It does not require computing the initial conditions in the NF space, so this challenging step is avoided. • It is unnecessary to transform each system initial condition to NF space as in [START_REF] Wang | Nonlinear Modal Decoupling of Multi-Oscillator Systems With Applications to Power Systems[END_REF], [START_REF] Wang | Nonlinear modal decoupling based power system transient stability analysis[END_REF]. • The impact of generator loading, load loss/addition, faults, e.t.c, on the oscillation frequency may be investigated for warnings of impending OOS. The proposed method may incur errors due to modal interactions besides the inherent truncation error caused by Taylor expansion. Simplification to NF [START_REF] Amraee | Transient instability prediction using decision tree technique[END_REF] assumes negligible effects of inter-modal interactions. This assumption is more valid for the SMIB, which has one mode or if the initial conditions for k = p are zero for the multi-machine system. When a multi-machine power system experiences a large disturbance, the initial condition for k = p is not zero in general, and the effect of modal interaction increases. Close to the stability boundary, cross-coupling interactions, as well as 4 th and higher-order interactions, may be significant. As a result, Ω p nl could be slightly farther/closer, leading to a lightly conservative or optimistic outcome. A normalizing factor is added to steer the cumulative effects of all modes to the particular mode of interest, which improves the challenge of inter-modal interactions. The normalizing factor is discussed in the next paragraph.

Following the definitions in [START_REF] Matsushita | Modal Analysis of Multi-Degree-of-Freedom Systems[END_REF] for mechanical systems, the modal mass m * p = Φ T p MΦ p is equal to the sum of (mass) × squared (mode displacement) for each mode (i.e., m * p = N q=1 m q Φ 2 pq ). Similarly, the modal stiffness k * p = Φ T p KΦ p is the sum of the strain energy stored in each oscillator, and the total strain energy indicates that k * p = m p Ω 2 p . Thus, we define T se p = m p Ω 2 p , T se p nl = m p Ω 2 p nl , where T se p and T se p nl are the total modal strain energies of the p th mode at linear and nonlinear frequencies, respectively. The normalizing factor is thus defined for p-th mode as

S.F = T se p nl T se p . (17)
So, for any system condition, the frequency shift is tuned as:

Ωnlsc = Ω nlsc × S.F. (18)
The results will later show that the proposed approach incurs tolerable errors for the considered model. In the following sections, numerical simulations are presented to test the method.

IV. POWER SYSTEM CASE STUDIES

In this section, the proposed method is first illustrated with an SMIB power system and then applied on the IEEE 3-and the IEEE 50-machine power systems.

A. Case Study 1 -SMIB Power System

Consider the previously mentioned SMIB power system, which is shown in Fig. 2.

The proposed method is applied following the steps below: • Therefore, the NF model [START_REF] Amraee | Transient instability prediction using decision tree technique[END_REF] becomes: R = -61.1866R + (0.7322 + 10.1978)R 3 -0.0239R Ṙ2 .

• With R = Ṙ = 0, evaluate the critical amplitude in R from [START_REF] Poincaré | Les mthodes nouvelles de la mcanique cleste[END_REF] as R cr = 2.366, and substitute into [START_REF] Wang | Power system transient stability analysis using high-order taylor expansion systems[END_REF] to get the critical modal amplitude as: a cr = 2.366 + 0.0447 × (2.366) 2 + 0 = 2.616. • With the critical modal amplitude a cr determined, use [START_REF] Hauer | A perspective on wams analysis tools for tracking of oscillatory dynamics[END_REF] to calculate the critical frequency shift as: Ω nlcr = 7.8222× 1 -3×10.93-61.1866×0.0239 8×61.1866 × 2.616 2 = 4.0764 rad/s.

• Therefore, the oscillation frequency should not go below 4.0764 rad/s to maintain synchronism. • Since the SMIB has one mode only, the normalizing factor S.F = 1. Thus, Ωnlsc = Ω nlsc for any disturbance.

Observe that the critical frequency shift we just calculated is nearly identical to the one estimated in [Fig. 1b, Section III-A].

Because it only has one mode, it may not be obvious what the NF did in the preceding example, but the job of the NF is to decouple the multi-oscillator system into separate oscillators of the type exemplified, which will be shown in next subsection.

The SMIB was subjected to various disturbances by varying the rotor angle initial condition to put the criterion to the test. The frequency shift is computed for each initial condition, and the nonlinear solution is evaluated by integrating the nonlinear model in MATLAB ® . In this and subsequent sections, the nonlinear solution obtained with MATLAB ® is regarded as the "true" solution. The outcomes for the various conditions are presented in Table I. For comparison, the result of the equal-area criterion is also shown. We assume the reader is familiar with the EAC method, so we leave out the specifics. The Table shows that when the rotor angle amplitude reaches 2.8798 rad, the proposed method classifies the system as unstable because the Ωnlsc = 4.0705 rad/s, which is less than the Ω nlcr = 4.0764 rad/s. The EAC and the true solution are stable up to 2.8815 rad, implying that the proposed method is conservative with an error of only 0.06% in angle in this case. The response of the reconstructed NF, linear, and true solutions is shown in Fig. 3, confirming that the proposed method becomes unstable almost at the same time as the true solution. The NF's accuracy originates from the frequency shift against the constant frequency during a disturbance, as the linear solution assumes. It is noteworthy that if the swinging of a generator against the rest is of interest, as is the concept of multi-machine EAC [START_REF] Alinezhad | Out-of-step protection based on equal area criterion[END_REF], the proposed method will work very similarly to the EAC. Such limitations, however, are unnecessary. In the sections that follow, the proposed method will be deployed for multi-machine applications. The IEEE 3-machine power system [START_REF] Sauer | Power system dynamics and stability[END_REF] shown in Fig. 4 has two electromechanical modes of frequencies 13.72 rad/s (Mode 1) and 8.82 rad/s (Mode 2). According to participation factor analysis, Mode 1 is most associated with generator G3, whereas Mode 2 is associated with G2. Mode 2 also has inter-area effects because it is excited regardless of the fault location. A three-phase fault was added to Bus 4, and it to Table II, the proposed method detects instability at 0.25s, while the true instability occurs at 0.24s. In this case, the results also show that for the chosen fault conditions, Mode 2 is less stable than Mode 1. This can be seen in the table, where the frequency of Mode 2 decreases faster as the stress increases, whereas the frequency of Mode 1 decreases slowly.

A second case is considered by adding a three-phase fault to Bus 9 and clearing the fault at increasing time intervals until instability occurs. Mode 1 is dominant in this case. Table III lists the values of Ωnlsc after each clearing time. According to the table, the proposed method detected instability at 0.26s, which is earlier than the true solution (i.e., 0.29s). The table also shows that Mode 1 is driving the instability in this case. The numbers in the 3rd column, i.e., Ω2 nlsc , are the same, which appears quite ideal based on the participation factors. Notice that by observing a significant decrease in the oscillation frequency, one can already tell if there is an intended OOS. The affected generators can be controlled to reinforce the system, preventing the OOS and potential blackout. Such insight is in contrast to traditional tools' YES or NO feedback.

C. Case Study 3 -IEEE 50-Machine Power System

The proposed method is tested using the IEEE 50-machine power system to show the application to large systems. A section of the test system is shown in Fig. 5. The power Fig. 5. Portion of one-line diagram of the IEEE 50-machine system [START_REF] Xu | A decompositionbased practical approach to transient stability-constrained unit commitment[END_REF] flow data was extracted from [START_REF] Zimmerman | MATPOWER (Version 7.0)[END_REF], [START_REF] Zimmerman | Matpower: Steady-state operations, planning, and analysis tools for power systems research and education[END_REF] while the dynamic data was taken from [START_REF] Liu | Assessing placement of controllers and nonlinear behavior of electrical power system using normal form information[END_REF]. We know from previous studies that a fault at Bus 7 causes very high amplitude oscillation, which is primarily associated with local modes and leads to early instability. Besides, a fault at Bus 34 activates both local and inter-area modes, resulting in the simultaneous instability of several generators. Thus, two scenarios are thoroughly examined: 1) a fault on Bus 7 and 2) a fault on Bus 34.

The results of mode stability monitoring for fault at Bus 7 using the proposed method are shown in Table IV. Because mode 8 is the first to exhibit instability, only its stability is displayed in Table IV. It is discovered that Mode 8 becomes unstable at a time 0.26s with Ω8 nlsc less than its critical value 9.05 rad/s. The true instability occurs at 0.27s when generator 7 (G7) loses synchronism with the rest of the system (see Fig. 6a). Mode 8 is local, and is associated most with G7 as confirmed by its shape (see Fig. 6b).

-0.5 0 0.5 In the second case of fault at Bus 34, we observed Mode 18 and Mode 28 becoming unstable almost simultaneously, resulting in widespread instability. Table V shows the results for monitoring these two modes. The results show that syn- chronism is lost at 0.44s with our method because Ω18 nlsc falls below its critical value of 8.21 rad/s (see column 2). At the same time, Ω28 nlsc is about 1% away from its critical value, as shown in column 3. As shown in Fig. 7a, the true instant of instability is 0.40s. Mode 18 is associated with a single generator, whereas Mode 28 involves multiple generators. The shape of Mode 28, shown in Fig. 7b, confirms the separation of G1, G3, and G4 from the rest, as observed in Fig. 7a.

Table VI provides a detailed comparison of the proposed method's decisions with the true decisions for different cases. Our observation after other tests is that the method tends to be more accurate when several modes are not simultaneously excited. The maximum error in the 5 th column is 10%. However, the cases are not exhaustive, and higher errors may exist because NF results are only reliable within the approximation's validity boundary. Modern power systems also include additional equipment such as FACTS, HVDC, IBR, etc. As a result, the proposed method should be tested on such a system in the future. We emphasize that the method's strength is in the warning for impending OOS. As a result, even when the precise time eludes the proposed method, it remains faithful in indicating impending OOS. Thus, the proposed method should be used for mode instability warning signals.

D. Computational Complexity and Scalability

The NF method is generally considered difficult because: (i) it needs nonlinear iterations to find the initial condition in the NF space for any considered operating point, which can be time-consuming; and (ii) it requires Taylor and modal expansions to construct its model, which is also time-consuming when considering large systems. Fortunately, the proposed method avoids these primary challenges. It does not require the computation of the initial conditions in the NF space. The NF model is constructed once and may merely be repeated if the power flow is updated. The solution to attempting a very large system (already used in [START_REF] Vittal | Determination of generator groupings for an islanding scheme in the manitoba hydro system using the method of normal forms[END_REF] for systems with over 300 generators) is to reduce the network and concentrate on a specific region of the network when constructing the model. Remarkably, only the computations related to the modes of interest are necessary with our method. The computation can drastically reduce if the dominant modes are identified, and then the analysis can focus on them.

Once the NF model is constructed, the proposed tool is very speedy. It takes about 3.5 milli seconds (4 th column, Table VII) to check all the modes in the IEEE 50-machine power system with an Intel(R) Core (TM)i7-7700HQ, 2.8GHz personal computer. Construction of the NF model is the only step that takes a little longer. The total time to develop the NF model in the case of the IEEE 50-machine power system is approximately 5.13 minutes (2.21 minutes for ( 6)-( 7) and 2.92 minutes for ( 11)-( 16) ). As earlier noted, the NF is built once and may be updated if the power flow updates. This runtime for developing NF is not prohibitively expensive and can be performed regularly in power system operation. The proposed method is therefore useful for online applications. For very large systems, the NF model can be built offline and then used online. It should be noted that the runtime for developing NF can reduce drastically if (i) more powerful computing resources are utilized, e.g., HPC, or (ii) only selected modes are analyzed, which deserve further detailed investigations.

A computational comparison with similar works [START_REF] Wang | Nonlinear Modal Decoupling of Multi-Oscillator Systems With Applications to Power Systems[END_REF], [START_REF] Wang | Nonlinear modal decoupling based power system transient stability analysis[END_REF], where the modal model (i.e, ( 6)-( 7)) is built in first order is shown in Table VII. Indeed, avoiding first-order representation saves enormous time, especially for large systems. This action, however, updates all prior computations and may be an issue for online applications in terms of computation. Instead of the traditional power flow, static state estimation can be utilized to improve measurement errors. The first block in the right section, which regularly scans the system condition (i.e., checks the rotor angles), may be realized by real-time measurement tools. We have assumed that the rotor angles will be estimated correctly.

E. Practical Significance

The proposed method has the potential for online implementation, and power system operators can use it to estimate the proximity of modes to instability quickly and roughly. The proposed approach, like many others, inherits the shortcomings of classical modeling. However, the nature of the EO is the same regardless of the model used, and detailed modeling does not significantly change the frequency of the oscillation. The inability to detect the exact time of OOS, unlike many previous schemes, does not rule the proposed tool out because the new tool can provide useful insights for controls by e.g., ranking vulnerable modes and suggesting generators to regulate based on their participation in those modes. The proposed method may also be utilized as an extra tool by power system planners to plan system operation and expansion, taking OOS mitigation measures into account. Given the current state of this study, we advocate the proposed approach as a backup tool. For practical application, the new tool can be run routinely (say every 10 minutes) in the power system operation for warning signals. Once it raises a warning flag, one can use other more robust approaches to analyze the concerned generators in detail for appropriate control actions.

F. Compatibility with Existing Tools

As stated at the outset of this article, the proposed tool may be used in conjunction with existing ones to improve system reliability. The algorithm for the majority of previous methods is as follows: Is there an OOS? If yes, isolate the responsible generators. If no, no control actions are needed. The new scheme adds more decision questions: If there is no OOS, will it occur soon? If so, control the potential generators (e.g., loadshedding). The new tool may supplement WAM techniques, in which the proposed method determines the critical frequency shift while WAM estimates the modal frequency in real-time. Since our scheme requires only the rotor angle as an input, one may use it parallel with EAC-based methods. One step of the EAC-based approach reported in [START_REF] Alinezhad | Out-of-step protection based on equal area criterion[END_REF], for example, is a fast estimation of the rotor angle, based on local measurements of voltages and currents at the generator terminal. Both methods can use the same input estimated by [START_REF] Alinezhad | Out-of-step protection based on equal area criterion[END_REF]. However, the coordination between the proposed method and the existing OOS detection scheme may be a concern that merits further investigation.

V. CONCLUSIONS

This paper expanded the real-valued Normal Form proposed in [START_REF] Touzé | Hardening/softening behaviour in non-linear oscillations of structural systems using non-linear normal modes[END_REF] to investigate practical power system oscillations. Employing second-order modeling of the power system, the NF allows multi-machine power systems decoupled as separate oscillators. A novel method for monitoring out-of-step occurrence is proposed using the NF that results. The new approach was demonstrated with an SMIB power system, which showed favorable competition compared to the well-known equal area criterion, and then applied successfully to the IEEE 3-and IEEE 50-machine power systems with promising results.

Compared to time-domain simulation, the proposed tool has a maximum error of ±10% and computes quickly for tested cases. Thus, it has potential for online applications, and power system operators can use it to make a quick and rough estimation of OOS proximity and rank vulnerable modes. Therefore, the tool may provide important information leading to efficient corrective/preventive control actions. The paper considered only short-circuit scenarios, but the idea is general, as illustrated in Section III-A, and can extend to investigate other contingencies such as loss of load or generation.

Future work will compare the method with other OOS protection approaches. Strategies for pre-selecting modes to check for any disturbance should be investigated; to further reduce the computation time. The frequency shift phenomenon is general, and we envisage the extension of the proposed method to detailed models. However, with detailed generator models, putting the entire system model into the 2nd-order form can be challenging and should be re-visited.

Fig. 1 .

 1 Fig.1. Power output and fault affect the frequency and damping peaks indicated in the figure is 0.65 Hz (4.0841 rad/s), which is 48% less than the linear frequency. It implies that the two principal causes of OOS (loading and faults) manifest in decreasing the oscillation frequency.One can, therefore, conclude that in a stable mode, the frequency decreases as the generator output or fault level increases. This frequency's descent reaches a critical point below which the mode loses stability, and the probably affected generators are those participating significantly in that oscillation mode. This phenomenon motivates the mode's stability monitoring for OOS protection.

Fig. 2 .

 2 Fig. 2. Single-line diagram of a SMIB Power System • The electrical output power P e = EV Xs sin δ. • Evaluate the Taylor expansion terms in (7) as: K = EV Xs cos δ 0 = 1.1361, F 2 = -EV 2Xs sin δ 0 = -0.1522, F 3 = -EV 6Xs cos δ 0 = -0.1894. • Solve (K -Ω 2 M )Φ = 0 to get Ω = 7.8222 and Φ = 1. • Use (8) and (9) to obtain the modal model (10) as: ÿ + 61.1866y -8.1974y 2 -10.1978y 3 = 0. • Following [23], compute the NF coefficients in (11) as: a = 0.0447, b = 0.0015, γ = -0.0893, r = u = µ = v = 0, and in (12) as A = -0.7322, B = -0.0239.

Fig. 3 .

 3 Fig. 3. Response of the system during unstable fault

Fig. 4 .

 4 Fig. 4. One line representation IEEE 3-machine power system was cleared over time by removing the fault until instability occurred. Mode 2 is dominant in this case. The values of Ωnlsc after each clearing time are listed in TableII. According

Fig. 6 .

 6 Fig. 6. (a) Trajectory at 0.27s for fault at Bus 7, (b) Shape of Mode 8

Fig. 7 .

 7 Fig. 7. (a) Trajectory at 0.40s for fault at Bus 34, (b) Shape of Mode 28

Fig. 8

 8 Fig.8is a pictorial representation of the necessary steps of the proposed method. The left half of the figure represents the construction of the NF model, which is a one-time process, while the right portion represents OOS monitoring and admits online operation. Because the NF method is local and based on an asymptotic expansion in the neighborhood of an equilibrium point, the NF model can be reinitialized if the power flow changes to reset the critical frequency shift. This action, however, updates all prior computations and may be an issue for online applications in terms of computation. Instead of the traditional power flow, static state estimation can be utilized to improve measurement errors. The first block in the right section, which regularly scans the system condition (i.e., checks the rotor angles), may be realized by real-time measurement tools. We have assumed that the rotor angles will be estimated correctly.

Fig. 8 .

 8 Fig. 8. Summary of the proposed method

TABLE I COMPARING

 I PROPOSED METHOD WITH EAC ON SMIB POWER SYSTEM

	Post-fault Ω nlcr = 4.0764 rad/s Proposed	EAC	True
	(rad)	Ωnlsc (rad/s)	Decision	Decision Decision
	1.7453	6.6175	Stable	Stable	Stable
	2.4435	5.2168	Stable	Stable	Stable
	2.6180	4.7833	Stable	Stable	Stable
	2.8798	4.0705	Unstable	Stable	Stable
	2.8815	-	Unstable	Unstable Unstable

TABLE II CASE

 II STUDY 2: FAULT AT BUS 4 MODE 1 (Ω 1 nlcr = 8.55), MODE 2 (Ω 2 nlcr = 5.44 )

	Fault	Ω1 nlsc	Ω2 nlsc	Proposed	True
	duration (s)	(rad/s)	(rad/s)	Decision	Decision
	0.15	13.65	8.32	Stable	Stable
	0.22	13.41	6.67	Stable	Stable
	0.24	13.30	5.89	Stable	Unstable
	0.25	13.22	5.36	Unstable	Unstable

TABLE III CASE

 III STUDY 2: FAULT AT BUS 9 MODE 1 (Ω 1 nlcr = 8.55), MODE 2 (Ω 2 nlcr = 5.44)

	Fault	Ω1 nlsc	Ω2 nlsc	Proposed	True
	duration (s)	(rad/s)	(rad/s)	Decision	Decision
	0.15	13.00	8.81	Stable	Stable
	0.22	10.67	8.81	Stable	Stable
	0.25	9.00	8.81	Stable	Stable
	0.26	8.30	8.81	Unstable	Stable
	0.29	-	8.81	Unstable	Unstable

TABLE IV CASE

 IV STUDY 3: FAULT AT BUS 7; MODE 8 (Ω 8 nlcr = 9.05)

	Fault	Ω8 nlsc	Proposed	True
	duration (s)	(rad/s)	Decision	Decision
	0.15	15.27	Unstable	Stable
	0.22	11.78	Unstable	Stable
	0.25	9.27	Unstable	Stable
	0.26	8.33	Unstable	Stable
	0.27	-	Unstable	Unstable

TABLE V CASE

 V STUDY 3: FAULT AT BUS 34 MODE 18 (Ω 18 nlcr = 8.21), MODE 28 (Ω 28 nlcr = 5.76)

	Fault	Ω18 nlsc	Ω28 nlsc	Proposed	True
	duration (s)	(rad/s)	(rad/s)	Decision	Decision
	0.15	13.26	9.16	Stable	Stable
	0.25	12.83	8.75	Stable	Stable
	0.35	11.44	7.67	Stable	Stable
	0.40	9.25	6.41	Stable	Unstable
	0.44	8.13	5.83	Unstable	Unstable

ACKNOWLEDGMENT

This work was authored in part by the National Renewable Energy Laboratory, operated by Alliance for Sustainable Energy, LLC, for the U.S. DOE under Contract No. DE-AC36-08GO28308. The views expressed in the article do not necessarily represent the views of the DOE or the U.S. Government. The U.S. Government retains and the publisher, by accepting the article for publication, acknowledges that the U.S. Government retains a nonexclusive, paid-up, irrevocable, worldwide license to publish or reproduce the published form of this work, or allow others to do so, for U.S. Government purposes.

This study is a part of a PhD project supported by TETFUND Nigeria from 2017 to 2020. The funder's only role was to provide the funds.

N.S.


J.-M Malherbe
email: jean-marie.malherbe@obspm.fr

The quadruple spectroheliograph of Meudon observatory (1909-1959)

Keywords: Spectroheliograph, spectroscopy, monochromatic imaging, Sun, photosphere, chromosphere

The spectroheliograph was invented independently by Henri Deslandres (France) and George Hale (USA) in 1892, following the spectroscopic method suggested by Jules Janssen in 1869. This instrument is dedicated to the production of monochromatic images of the Sun in order to reveal the structures of the photosphere and the chromosphere at various altitudes. Sporadic observations started in Paris, but Deslandres moved soon to Meudon and designed, with Lucien d'Azambuja, an universal and powerful instrument, the quadruple spectroheliograph. It was devoted to systematic observations of the Sun (the long-term activity survey since 1908) and scientific research in solar physics. This paper describes the instrument and presents some original observations made with the high dispersion 7-metre spectrograph. It was dismantled in the sixties, but the solar patrol continued with the 3-metre chambers with Hα and CaII K lines, and is still working today with the numerical version of the spectroheliograph.

-INTRODUCTION

The principle of imaging spectroscopy of the solar surface was early introduced by Janssen (1869) after the discovery of how to observe prominences at any time, outside eclipses, using spectroscopic means: "with a spectroscope rotating around its axis, the entrance slit scans the Sun; a second slit in the spectrum selects a spectral line of interest and, using an ocular, the retinal persistence forms a monochromatic image". George Hale (USA) and Henri Deslandres (France) invented independently the spectroheliograph in 1892 [START_REF] Malherbe | 130 years of spectroheliograms at Paris-Meudon observatories (1893-2023)[END_REF] on this basis. Deslandres (figure 1, see the biographies of d 'Azambuja, 1948 and[START_REF] Lequeux | Henri Deslandres: from molecular to solar physics[END_REF], produced in Paris the first spectroheliograms of CaII K intensities using a thin slit in the spectrum integrating K2v, K3 and K2r (figure 2) and the first "section" spectroheliograms using a large slit (figure 3), which are made of the juxtaposition of CaII K spectra (cross sections by spatial steps of 30") for the measurement of radial velocities (see [START_REF] Deslandres | Recherches sur l'atmosphère solaire. Photographie des couches gazeuses supérieures. Etude simultanée des formes et des mouvements des vapeurs[END_REF].

Figure 2 : the first CaII K spectroheliograms of the disk ("chromosphère du disque", top) and of the limb ("chromosphère du bord", bottom) with a mask upon the disk. After [START_REF] Deslandres | Recherches sur l'atmosphère solaire. Photographie des couches gazeuses supérieures. Etude simultanée des formes et des mouvements des vapeurs[END_REF]. Courtesy CNUM/CNAM. Deslandres moved to Meudon in 1898 and hired there an assistant, Lucien d'Azambuja (figure 1, see the biographies of [START_REF] Rösch | Lucien d'Azambuja (1884-1970)[END_REF][START_REF] Martres | Lucien d'Azambuja[END_REF]. They together built a much more powerful instrument than the one from Paris, the large quadruple spectroheliograph (section 2), which was completed in 1908. A removable spectrohelioscope was added later (section 3); systematic observations of the photosphere and chromosphere (section 4) started immediately, together with scientific research (section 5) leading to d'Azambuja's thesis (1930) and to the big memoir of d' Azambuga & d'Azambuja (1948) about the physics of solar filaments and prominences.

-THE LARGE QUADRUPLE SPECTROHELIOGRAPH OF MEUDON

The quadruple spectroheliograph was built by [START_REF] Deslandres | Recherches sur l'atmosphère solaire. Photographie des couches gazeuses supérieures. Etude simultanée des formes et des mouvements des vapeurs[END_REF] and d'Azambuja in a new house (1907) and was fed by a two mirror coelostat (figure 4). Four chambers were installed in a 20-metre long laboratory. Chambers n° I and n° II had moderate dispersion and were devoted to systematic observations, respectively with CaII K and Hα lines. Chambers n° III and n° IV were dedicated to scientific research, for multi-line or high dispersion observations. They could not work simultaneously. The entrance objective (250 mm diameter, 4 m focal length), entrance slit, collimator, grating and prisms were common to all combinations, which are detailed in figures 5. ,4,5,6 (4 possible combinations) Removable parts: t,m of chamber III, F4 (spectrograph IV), P' (disperser IVbis), tubes of chamber I and reducer IVbis Motorized translators: t,m for entrance objective Op and photographic plate carriers Manual translators: R or P, M1 (for chamber IV) or O1 (for chamber III) Figure 5 : optical parts of the Meudon quadruple spectroheliograph (after [START_REF] Deslandres | Recherches sur l'atmosphère solaire. Photographie des couches gazeuses supérieures. Etude simultanée des formes et des mouvements des vapeurs[END_REF]d'Azambuja, 1930). Some of the different combinations were not simultaneously available (courtesy OP).

The chamber n° II was used for systematic observations of Hα with a plane grating (figure 6), while chamber n° I was dedicated to CaII K with a set of three prisms (figure 6). Both chambers were compatible with chamber n° III which was there for research purpose, allowing multi-line observations, with similar dispersion provided by the same focal length (spectra of 86 mm height, formed by the O1, O2, and O'1 lenses, respectively for chambers n° I, II, III). The beam was directed towards the spectrograph n° III by the flat mirror M1. The use of the high dispersion spectrograph n° IV was more complicated, because some parts had to be introduced, and other removed. As the focal length of the concave mirror M2 was 7 m, the height of the spectrum at the slit F4 was 205 mm. An associated afocal reducer IVbis, composed of the combination of a second 7 m concave mirror M3 and one of the four lenses O3, O4, O5, or O6, allowed to reduce the image size without any alteration of the spectral resolution determined by the width of the slit F4. P' was a 30° prism, it was there to eliminate parasitic light (or orders) from the grating. The support of F4 and the tube of the reducer IVbis had to be dismounted to permit daily observations with the CaII K Removable chamber n° II. The dispersion provided by the four spectrographs is given by figure 7. It was possible to attain about 1 Å/mm for scientific research with chamber n° IV from the near Ultra Violet to the near Infra Red. , after d'Azambuja (1930).

Figure 8 shows the quadruple spectroheliograph configured for the use of the high dispersion chamber n° IV (7 m focal length, spectrum of 205 mm height). The associated afocal reducer n° IVbis was used to reduce the height of the spectrum in order to be more compatible with classical photographic plates (13 x 18 cm² max). The reduction factor was equal to the ratio between the focal length of the concave mirror M3 (f = 7 m) and the focal length of one of the available lenses O3, O4, O5 or O6. The reducer included a low dispersion 30° prism to eject parasitic light. The dispersion of the chamber n° IV was convenient to produce monochromatic images with a very narrow bandpass (as low as 0.05 Å). For that reason, it was employed for scientific purpose, not for daily observations. Alternatively, Figure 9 shows the quadruple spectroheliograph configured for systematic observations of the photosphere and the chromosphere with chambers n° I (CaII K) and n° II (Hα). For that purpose, the slit F4 of chamber n° IV had to be displaced, together with several parts of the associated reducer n° IVbis which was used in combination with it. CaII K (K1v for the photosphere and K3 for the chromosphere) was observed with the set of three prisms, while Hα was formed by the plane grating. The instrument and observing program are described by d'Azambuja (1920 a and 1920 b)

-THE VISUAL SPECTROHELIOSCOPE

Figure 10 shows the spectroheliosocope which was installed later (1940) by d'Azambuja. It was a visual device, using two Anderson's square prisms (figure 11). The first prism, in front of the entrance slit F of the spectrograph, rotated in order to scan active regions of the Sun. The second prism, rotating at the same speed, was located after the slit F2 which was isolating the Hα line in the spectrum. The size of the square prisms were in the ratio (2.3) of the focal lengths of the collimator (1.30 m) and the chamber (3.0 m). The observer used an ocular (figure 11) to see the monochromatic image formed by the output slit F2 and the retinal persistence. This system was very fast to install and remove (to allow photographic observations), and was of great interest for visual inspection of solar activity, so that the observers used it daily.

-SYSTEMATIC OBSERVATIONS WITH CHAMBERS n° I AND n° II

Systematic observations were performed with the configuration of figure 9: spectrograph n° I for CaII K using 3 flint prisms of 60° at minimum deviation, providing a spectral resolution of 0.15 Å, and spectrograph n° II for Hα using a Rowland plane grating, 568 groves/mm, providing in the first order a spectral resolution of about 0.40 Å. The magnification of the spectrographs was 2.31 (3.0 m chambers, 1.30 m collimator), so that the size of the Sun on the photographic glass plates was 86 mm. The width of the entrance slit was 0.03-0.04

Chamber II, Hα Rotating prism

Rotating prism

F

Motor

F2 ocular

Mechanical transmission

mm (1.5"-2") well adapted to the usual seeing of 2" in Meudon; while the second slit in the spectrum was 0.075 mm wide. From 1908 to WW1 (1914), observations consisted of daily spectroheliograms in CaII K3 (line centre, chromosphere), CaII K1v (violet line wing, photosphere) and Hα centre (chromosphere). Observations were totally interrupted during WW1. Figure 12 shows typical observations performed daily from 1919 to 1939 until the start of WW2. During this period, CaII K3, CaII K1v and Hα continued, and a "section" spectroheliogram with an enlarged slit in the spectrum (1 mm = 2 Å) was added for Dopplershift measurements of the CaII K line; for that purpose, full spectra of 86 cross sections of the solar surface were recorded by spatial steps of about 1 mm = 22". Contrarily to WW1, observations were not stopped by WW2, but the lack of manpower was an important constraint, so that "section" spectroheliograms were abandoned, and never restarted after WW2. Typical systematic observations done after WW2 include the initial choice of 1909 (CaII K3, K1v and Hα), plus a long exposure CaII K3 for limb features and prominences, with a disk attenuator (because prominences are 5-10 times fainter than the disk).

Systematic observations continued under the direction of Lucien d'Azambuja until his retirement in 1954 (he was 70 years old). At this date, Mrs d'Azambuja became in charge of observations, but his husband was still working as a volunteer; she retired in 1959. Gualtiero Olivieri (1928-2018) and Marie-Josèphe Martres (1924Martres ( -2014) ) became respectively responsible of observations and of the synoptic map program initiated in 1919 by L. d'Azambuja (these are synthetic drawings of solar activity for each solar rotation, based on daily spectroheliograms, and showing the average position of filaments, sunspots and facular regions).

-SPECIAL SPECTROHELIOGRAMS WITH THE HIGH DISPERSION CHAMBER n° IV

During his thesis work, that he defended in 1930, Lucien d'Azambuja explored many spectral lines with the large 7 m spectroheliograph devoted to research (figure 8). High dispersion is required to study thin photospheric lines such as those chosen by d'Azambuja: MgI 3838 Å, SrII 4078 Å, FeI 4202 Å, CaI 4227 Å, FeI 4384 Å, MgI 5184 Å, HeI 5876 Å, NaI 5890 Å, as well as chromospheric infrared lines of CaII 8498 Å, CaII 8542 Å, HeI 10830 Å. This large spectroheliograph (dismantled in the sixties) was able to offer the dispersion of 1 mm/ Å, providing excellent spectral resolutions on line profiles of about 0.05 Å. The set of 3 prisms was used in the blue part of the spectrum, and the plane grating for longer wavelengths.

Figure 13 displays three images, respectively in the centre of Hα line and in the blue wing. Filaments are still visible at -0.25 Å, but vanish at -0.50 Å (near the inflexion points of the line profile) and the photosphere (below the chromosphere) begins to appear at -0.50 Å, except in the case of Dopplershifts (0.5 Å = 25 km/s). Such observations are of interest, because they were never done daily, contrarily to the Hα centre. , K2v, K1v). 29 July 1927. After d'Azambuja (1930). Courtesy OP.

Figure 14 displays three images, respectively in the centre of CaII K line (K3) and in the blue wing. Filaments and bright plages are perfectly visible in K3, but this is not the case of sunspots. Filaments are fairly observable at -0.25 Å (K2v) and sunspots begin to appear. But filaments completely vanish at -1.27 Å (K1v) while the photosphere, below the chromosphere, becomes clearly visible, with sunspots and bright faculae (which are associated to chromospheric plages). Observations of K2v are particularly rare in the collection, contrarily to K3 and K1v which are systematic.

Figure 15 shows an exceptional observation of the Balmer series, Hα, Hβ, Hγ, Hδ and Hε. The contrast of filaments and plages (chromosphere) vanishes with decreasing wavelength, because lines form deeper. On the contrary, sunspots appear progressively along the series, which reveals that Hα is the best choice for chromospheric structures, while Hδ or Hε are not better than CaII K1v for sunspots and faculae (however Hε is contaminated by the red wing H1r of the CaII H line). This observational series is the only one of the Meudon collection. Helium lines where also studied by Mr and Mrs d'Azambuja. Figure 16 shows an exceptional observation in HeI D3 5876 Å of an active region producing the huge flare of 26 July 1946, while figure 17 shows the first world-wide observation done in the infrared HeI 10830 Å line. These observations are absolutely unique in the collection (d'Azambuja & d'Azambuja, 1938). Several decades later, Kitt Peak National Observatory started systematic observations in this infrared line of neutral Helium.

Infrared lines of CaII where investigated by d' Azambuja (1930) in his thesis work. Figure 18 shows images in CaII 8498 Å and 8542 Å, on 4 September 1928, in comparison with the usual CaII K3 3934 Å. These observations are extremely rare in Meudon collection and showed that chromospheric structures are better seen with CaII K3 (figure 19). Several decades later, Kitt Peak National Observatory started systematic observations in CaII 8542 Å in order to measure chromospheric magnetic fields (the Zeeman sensitivity is much better for this infrared line than for violet lines, such as CaII K or CaII H).

Figure 19 shows also a detail of figure 18, for comparison of HeI 10830 Å and CaII K. While filaments appear dark in both lines, this is not the case of plages, which are bright in CaII K but dark in HeI 10830 Å. , compared to CaII K3 3934 Å (top). 29 June 1928. After d'Azambuja (1930). HeI 10830 Å, compared to CaII K3 3934 Å (bottom). 17 September 1938. Courtesy OP.

CaI 4227 Å line, as well as many other photospheric lines, was studied by d' Azambuja (1930) during the course of his thesis. Figure 20 shows images in the line core and wing, in comparison with the usual CaII K1v spectroheliogram (the blue wing of CaII K line) for the photosphere. It sounded clearly that this Calcium line (and many thin lines) provided spectroheliograms that are comparable to CaII K1v, so that observations in these explored wavelengths were never reproduced later.

-CONCLUSION

The spectrographs n° III, n° IV and associated reducer n° IVbis remained unique during many decades and were dismantled after the retirement of Lucien and Marguerite d'Azambuja. However, parts of the 7 m spectrograph n° IV were transformed in the same building after 1960 for the measurements of magnetic fields through the Zeeman effect, using a new telescope (a 40 cm Newtonian telescope fed by a 80 cm Foucault siderostat, independent of systematic observations). The Zeeman analyser was the Hale's grid, composed of alternatively quarter (λ/4) and three quarter (3λ/4) wavelength strips (of a few arc seconds wide) and a polaroid (transmitting σ+ and σ-Zeeman components). λ/4 and 3λ/4 measurements were simultaneous, but not cospatial. This was the starting point of many developments of solar magnetometry techniques which will be related in another paper. Systematic observations in CaII K, CaII H and Hα are still active, with a modernized version of chamber n° II, but now use a fast scientific CMOS detector.

-REFERENCES

L. d'Azambuja, "Le spectrohéliographe et les photographies quotidiennes de l'atmosphère du soleil. Part 1", 1920a, L' Astronomie,34,[433][434][435][436][437][438][439][440][441][442][443][444][445][446] L. d'Azambuja, "Le spectrohéliographe et les photographies quotidiennes de l'atmosphère du soleil. Part 2", 1920b, L'Astronomie, 34, 487-498. L. d'Azambuja, "Recherches sur la structure de la chromosphère solaire", thesis, 1930, Annales de l'Observatoire de Paris, section d'astrophysique, à Meudon, volume 8, fascicule 2, Gauthier-Villars editor, Paris. On-line at Gallica/BNF. https://gallica.bnf.fr/ark:/12148/bpt6k6537953n L. d'Azambuja and M. d'Azambuja, "Images Nouvelles de la chromosphère solaire obtenues au spectrohéliographe avec la raie infrarouge de l'hélium 10830 et avec la raie Hε de l'hydrogène ", 1938, Bulletin Astronomique, 11, 349-359.

Figure 1 :

 1 Figure 1 : Lucien d'Azambuja and Henri Deslandres (left), 1913. Lucien d'Azambuja's jubilee (50 years of astronomy) occured in 1949 and a ceremony was organized in 1950; among participants, (1) Bernard Lyot, (2) Marguerite d'Azambuja, (3) Lucien d'Azambuja, (4) André Danjon, director (courtesy OP).

Figure 3 :

 3 Figure 3 : the first CaII K "cross section" spectroheliograms of the disk (top) and of the limb (bottom) with a mask upon the disk, showing CaII K spectra by spatial steps of 30". After Deslandres (1897). Courtesy CNUM/CNAM.

Figure 4 :

 4 Figure4: the two-mirror coelostat feeding the quadruple spectroheliograph (courtesy OP).

Figure 6 :

 6 Figure 6 : the plane grating (left), the set of three prisms (centre) and the motorized photographic plate carrier (right). Courtesy OP.

Figure 7 :

 7 Figure 7 : the dispersion of the spectrographs in Å/mm. Chamber n° I (prisms) and n° II (grating) were used for systematic observations of CaII K and Hα. The spectrograph n° IV, with 7 m focal length instead of 3 m, provided the best dispersion and was dedicated to solar research. Courtesy OP, after d'Azambuja (1930).

Figure 8 :

 8 Figure 8 : Meudon quadruple spectroheliograph configured for the use of the high spectral resolution chamber n° IV (7 m focal length) and associated focal reducer n° IVbis (the letters refer to figure 5). The moderate resolution chamber II for systematic observations of Hα remained available (courtesy OP).

  .

Figure 9 :

 9 Figure 9 : Meudon quadruple spectroheliograph configured for systematic observations of CaII K and Hα with moderate resolution chambers n° I and n° II (courtesy OP). The chamber n° IV and associated afocal reducer n° IVbis were partly dismounted for daily observations (the letters refer to figure 5, courtesy OP).

Figure 10 :

 10 Figure 10 : Meudon visual spectrohelioscope for Hα using two high speed rotating prisms and the moderate resolution chamber n° II (courtesy OP).

Figure 11 :

 11 Figure 11 : the rotating prism in front of the entrance slit F of the spectrograph (left), the ocular (centre) and the rotating prism at the output slit F2 (right). Courtesy OP.

Figure 12 .

 12 Figure12.Typical spectroheliograms of the period 1919-1939 (example of 11 July 1937). CaII K1v for the photosphere (sunspots, faculae); CaII K3 for the chromosphere (filaments, bright plages); CaII K line profiles (cross sections of the Sun by 22" steps for Dopplershift measurements); Hα for the chromosphere (filaments, plages). Courtesy OP.

Figure 13 .

 13 Figure 13. Observations at three positions in the Hα line. Filaments vanish in the line wing, except in the case of Dopplershifts, while photospheric structures begin to appear. 5 June 1930. After d'Azambuja (1930).

Figure 14 .

 14 Figure 14. Observations at three positions in the CaII K line (K3, K2v, K1v). 29 July 1927. After d'Azambuja (1930). Courtesy OP.

Figure 15 .

 15 Figure 15. The Balmer series Hα, Hβ, Hγ, Hδ and Hε of Hydrogen, 14 August 1947, and the CaII K3 spectroheliogram for comparison. Hα is chromospheric, but Hδ and Hε are rather photospheric. Courtesy OP.

Figure 16 .

 16 Figure 16. The flare of 25 July 1946 observed in HeI D3 at 5876 Å from 16:58 to 17:03. Courtesy OP.

Figure 17 .

 17 Figure 17. The first observation of the infrared HeI 10830 Å line (right), where both filaments and plages appear darker than the quiet Sun. CaII K3 (left) is given for comparison. 17 September 1938. Courtesy OP.

Figure 18 .

 18 Figure18.CaII 8498 Å and CaII 8542 Å, compared to CaII K3 3934 Å (right). 4 September 1928. After d'Azambuja (1930) and courtesy OP.

Figure 19 .

 19 Figure19.CaII 8542 Å, compared to CaII K3 3934 Å (top). 29 June 1928. After d'Azambuja (1930).HeI 10830 Å, compared to CaII K3 3934 Å (bottom). 17 September 1938. Courtesy OP.

Figure 20 .

 20 Figure 20. Spectroheliograms of the CaI 4227 Å line, compared to CaII K1v for the photosphere (at right, for reference, the chromospheric Hα line). After d'Azambuja (1930) and courtesy OP.