Optimal control theory : a method for the design of wind instruments

Optimal contr ol theory : a method f or the design of wind instruments G. Le V ey Abstract It has been asserted prev iously by the author that optimal control theory can be a v aluable framewo rk for theoretical studies about the shape that a wind instrument should hav e in order to sati sfy some optimization cri terion, inside a fairly general class. The purpose of the present work is to de velop this ne w approac h with a look at a specific criterion to be optimized. In this setting, the W ebster horn equation is regarded as a controlled dynamical equation in the space variable. Pressure i s t he state, the control being made of two parts : one variab le part, the inside diameter of the duct and one constant part, the weights of the elementary time-harmonic components of the velocity potential. T hen one looks for a control that optimizes a criterion related to the definition of an oscillation regime as the cooperation of se veral natural modes of vibration with the excitation, the playing frequency being t he one that maximizes the total generation of energy , as expo sed by A.H. Benade, following H . Bouasse. At the same time the relev ance of this criterion is questionned with the simulation results. I . I N T R O D U C T I O N Designing high quality musical instruments has b een the main concern of makers for cen turies, ac cumulating know-how from their prede cessors b y some kind of trial-an d-erro r p rocess, in volving musicians. On anoth er hand, sin ce Bernouilli and Lagrang e in the XVIII th century an d the important writing s of Helmholtz in the XIX th century a great d eal of research has also b een cond ucted in order to und erstand the pr inciples underly ing sound pro duction in tradition al o r mod ern musical instruments : see e.g. [2], [3], [6], [7], [13] for presentations of this v ast sub ject. W ind instruments su ch as br asses and woodwinds, ma king use o f an air column , must be designed in such a way as to properly arr ange the natur al fr equencies of this air column, in order that r egimes of o scillation can take p lace in conjunction with the nonlinea r flow-control device (reed or player’ s lips) [2]. As a conseq uence, harmonicity requirements for wind instruments between natural frequencies are desirable pro perties, d ue to the fact that inton ation, r esponsiveness an d tonal colou r can be attributed to well estab lished physical properties of the instrume nt [1]. I t is generally acc epted that these h armon icity require ments can lead to rather different shap es, which can b e (piecewis e) cones, cylinders or mor e complex, as can can be seen from the actual duct of real win d instruments. Nev ertheless, see [9] for interesting re sults about p iecewise c onic du cts, f ollowing an app roach based on characteristic impeda nces and transf er m atrices. Thu s there is r oom f or developping new method s that co uld h elp for the design. The present work is dedicated to such a task, with m ain accent on the metho dolog ical aspec ts. In [1 8], tra nsmission line mo delling of ho rns is used together with finite dimensional op timization te chnique s and a practical tool [ 19] was developed for instrument m akers, in or der to help improvin g existing instru ments as well as to design new o nes accordin g to a gi ven specification. Desired proper ties, specified in musical terms like intonation, response and pitch variability were compar ed with calculated values b ased o n an instrumen t’ s ac tual geometry . Similarly , in [22], a method is p resented for optimizing the shape o f a brass instru ment with r espect to its in tonation and impedan ce peak mag nitudes. The instrume nt is modelled using a one- dimension al tran smission line analog y with tr uncated cones. Th rough the use of an a pprop riate c hoice of design variables, the finite dim ensional optimization fin ds smooth horn profiles, that can also help in correctin g existing instruments, th e shape having been d esigned f rom an a prior i choice of a suc cession o f co nes. In [14], a frequen cy-domain method, using in verse quantum scatterin g for th e on e-dimen sional Klein-Gordon equatio n, allows to recover the ar ea functio n of a given aco ustical duc t in a noninv asi ve way , with out me asuring direc tly neither the inp ut imped ance nor the reflectance : this last one is m athematically derived fro m the wave radiated in r esponse to a high-imp edance so urce. In [16], parameter optimization tech niques have b een used to design shap es for b rass-trappin g Helm holtz reson ators that re sonate at a design set o f acoustic eig en values, taking into account physical and geometr ical constraints. In [ 15], a finite-elem ent eig enanalysis model of bars is used to compute optimal shap es fo r mallet percussion instrume nts such as vibrapho ne o r ma rimba-typ e bars. The objective function for the op timization p rocedu re is a targeted set of modal freq uencies. As one can see, o ptimization is at th e h eart of research o n the design question f or mu sical instrumen ts and th e eigenstru cture often plays a central role in the optimization criteria. The shape op timization of wind instruments, an in verse problem, appears to be a com plex on e, when co mpared to the analysis of what ph ysically o ccurs insid e a given in strument, a d irect, somewhat easier albeit not at all obvious, prob lem, as can be seen from the great amount o f researc h that it is th e subject of. The purpo se of the present work is a continuatio n of p reliminary ones by the autho r [20], [21], where the design qu estion for axisy mmetric wind instruments was posed as an o ptimal con trol problem fo r the so-called W ebster equatio n throug h a suitable refo rmulation as a dynamical equation in the spatial d imension, thereby using infinite dim ensional optimization of optimal control theory fo r continuous systems. This is in sharp contrast with the ab ove mention ned ref erences which use finite-dimen sional optimization schemes, i.e. with a fin ite number of paramete rs as un knowns in the theoretical f ormulatio n. IRCCyN UMR-CNRS 6597 and Ecole des Mines de Nantes - 4, rue A. Kastler - BP 20722 - 44307 NANTES Cedex 3, France. levey@emn. fr For mor e d etails abo ut the ho rn equation , see [12] o r th e more recen t [2 5] for a more mathe matical and updated exposition . Here, besides the expo sition o f how o ptimal control theory fits the wind instru ments design question , a special attention is giv en to the ch oice of one specific criterion, on the basis of p hysical consider ations, wh ereas in [2 0], only a fairly g eneral class o f cr iteria -a ffine in the contr ol- were c onsidered that led to a gener al result about the structure of the shape fo r the duct of a wind instrume nt. The cho sen criterion is in n o way un ique : many o ther ch oices could b e made , d epend ing on which characteristics (phy sical, perceptive and so on) are focused on. It a ppeared to be interesting to investigate from a musical acoustics po int of v iew as it comes from one definition of an o scillation r e gime inside wind instrumen ts, under stood as the c oopera tion of several natural mo des o f vibratio n, the ’playing freq uency’ b eing the one that maxim izes the to tal g eneration of energy [3], [2]. At this stage, it h as some arb itrariness in it : th e physical relev ance of this criterion is no t taken here for gr anted, as no othe r published material has b een found (than that of A.H. Ben ade) using this notion of o scillation regime. The reason for the choice of this criterion is twofold : first, to illustrate the g eneral methodolo gy that uses optim al control theory ; second, to question the rele vance of this cr iterion, that appeared to be appealing at first sight. The purpo se then is to question how interesting it is for desig ning win d instruments more than to assert such a relev ance a prio ri : n umerical results in section V show some q ualitative fea tures linked to th e chosen criterion. Nev ertheless, on the b asis o f this criterion , in th e co ntext of optim al contro l theory , necessary con ditions coming from the strong Pontryagin Max imum Principle [23], [8] are derived and some numerical simulations are con ducted, in or der to illustrate th e th eoretical predictions. In the presented approach, th e W ebster horn equation is regarded a s a contr olled dynamica l system in one spatial variable (the axis of an axisymmetric instrumen t). T he pressure (o r the velocity potential) is regarded as th e state and the contro l is ma de of two parts : one is the v ariable diameter (actually its deriv ati ve with respect to the axial space dim ension) of the duct. T he second part of the contr ol is made of the con stant weights of a time-harm onic decomposition of the velocity potential, which are unkn own, according to the ch osen criterion, the energy inside the du ct. Then one looks f or a con trol that ma ximizes the latter . This is anoth er difference with the approa ch in [1 8], [22], wh ere the measure is taken to be th e de viation from a measu red inpu t imped ance peaks location, mak ing use of a model-fitting oriented app roach. Wherea s these works can be d escribed as in verse pro blems, the present on e falls into the category of design pr oblems. But this d istinction is a m atter o f co n venience as bo th p roblems o nly differ fro m each other by th e ch osen criterio n, the unknown being the duct shape in each case. It is worth notin g in that respe ct that a similar criterion as in [22] c ould b e chosen here as well, while d efining the corresp onding o utput o f the model, leadin g to the same problem o f fitting a model to given data. This will be the subject of future research, with interesting com parisons to be made, the main point h ere being m ethodolo gical. T he appr oach f ollowed herea fter shows that, with the ch osen optimization criterion and the mod elling appro ach, the o btained shap e is piecewise continu ous. Moreover, when taking into account mor e and more time-harmon ic compo nents, is appea rs from the simulation s that the th e obtained duct sh ape is a cone : such a result tends to confirm that this criterion is likely no t to be sufficient to design high -quality instrum ents but th e realistic mode shapes that are obtained are en courag ing. It should be no ted that, to the best knowledge of th e author, this way of looking at the classical W ebster horn equ ation, as a co ntrolled dyn amical equ ation, appears to be new . One important poin t is that extensions of the optimal co ntrol methodology to other m ore r ealistic, n onlinear models of wind in struments and design criter ia can be co nsidered, althoug h Fourier analysis will generally be no mo re valid in such con texts : th e choice o f an adequate criterion will remain crucial. Numerical simulations show that th e proposed ap proach is ap pealing for design purpo ses. Nevertheless, the presented work has to b e considered as a p reliminary study that opens ne w perspectiv es for the design of wind instrumen ts : it is no t to be und erstood as an achievement in itself and has surely to be investigated more deeply for r eal application s. I n particula r , several aspects are not touch ed here such as taking accou nt o f tone holes in woodwinds (flutes, o boes, etc..), in corpo rating the nonlinear excitation mechan ism ( reed, lips for brasses) o r choo sing o ther criteria, in ord er to inclu de perceptive p arameters (such as given by a musician ), etc. ..Much remains to be don e. The pap er is organ ized as fo llows : section II, in ord er to fix the no tations, recalls the standard hor n equatio n and the definition of an oscillation regime as used h ere. Basic results from Optimal Control theo ry ar e briefly gi ven and lead to the formu lation of th e horn equatio n as a controlled dynamica l system. Th en a detailed exposition of how the oscillation regime is formu lated is g iv en in section III. Sec tion IV gives the math ematical formula tion o f the optimizatio n p roblem to solve, thr ough the necessary con ditions co ming from optim al co ntrol th eory . Nu merical results are given in sectio n V : they have to be considered as preliminary results as no thoroug h o r exhau sti ve investigation ha ve been conducted at this stage. Eventually , some con clusions and p erspectives are d rawn in section VI . An appendix provides a variational deriv ation of the horn wav e equation and reminds basic results from Sturm -Liouv ille theory . I I . P R E L I M I NA R I E S : H O R N W A V E E Q U A T I O N , O S C I L L AT I O N R E G I M E , O P T I M A L C O N T RO L T H E O RY A. W ave horn equation , oscillation re gime Consider an ax isymmetric horn with length L and sectio n diameter D ( x ) , a function o f the in depend ent variable, the space dimen sion x , acc ording to Fig. 1. Let ρ 0 be the mediu m m ass den sity , p ( x, t ) the acou stic pressure in the med ium, i.e. the deviation fro m th e atmosph eric pressur e, and v ( x, t ) the particular velocity inside the horn at abcissa x . Let k = ω /c be P S f r a g r e p l a c e m e n t s D D ( x ) x 0 L Fig. 1. Geometry of an axisymmetric horn the wa ve number, c b eing the sou nd velocity and ω the angular frequen cy in case of a har monic regime that is f ocused on here. For der iv a ti ves of any qu antity the index no tation will be u sed in th is section : for example the time ( resp. space) first deriv ati ve of φ will be denoted φ t (resp. φ x ). L et us first recall some usual mod elling hypo thesis : i) the diam eter D ( x ) is a slowly varying f unction o f x ; ii) acoustical quantities are fun ctions of x and t only ; iii) a velocity p otential φ exists, w hich means th at v = ∇ φ = φ x ; i v) one focuses on standing wave, time-harm onic solutions tha t write : φ ( x, t ) = ℜ ( ϕ ( x ) e iωt ) but on e will om it the ℜ part in the sequ el as it is related to the time depend ence that will no t be inv olved in the op timization process b elow ; v) in this work th ere is no dissipation . It can be shown tha t, in this situation, the pressure at x writes : p = − ρ 0 φ t . It is worth n oticing th at the first h ypothe sis is nec essary for th e 1D plane wave, W ebster equation to be valid [ 16] or [7], chap. 6. Anoth er limitation to th is equation is to lower frequen cy mod es. I n [7], the following appro ximate cond ition on the radiu s d eriv ati ve R ′ is given for the horn eq uation to be valid for plan e wa ves : 1 2 Z L 0 k R ′ 2 ( x ) dx ≪ 1 (1) This condition will b e used to fix bou nds o n th e control in th e simulations p resented in sectio n V. Wi th these conventions, the W ebster horn equation is easily obtained as the Eu ler equation of a variational pr oblem : 1 c 2 φ tt − φ xx − 2 D x D φ x = 0 (2) which reduces to : φ ′′ + 2 D ′ D φ ′ + k 2 φ = 0 (3) in case o f a harm onic regime of interest he re (see app endix a for details). Boundar y conditions are then g iv en for this equation, dependin g on the physical situation under study , mak ing it a Sturm-Liou ville p roblem (see appen dix b for a brief reminder of two useful results). On an other side, the notio n o f o scillation r e gime can be defin ed, following H. Bouasse [3 ] an d A.H. Benade [2]. It is “that state of the collective motion of a n air colu mn in which a nonlinear e xcitation mechan ism co llaborates with a set o f air column modes to mainta in a stead y o scillation contain ing several ha rmonically related fr equency comp onents, each with its own definite amplitud e. Then, the ’play ing frequency’ is th e one that maximizes the total gene ration of energy” [2]. Chosing su ch a definition as the basis of an optimization criterio n relies upon physical/aco ustical considerations of what happen s inside the instrumen t. Nevertheless, it has some arbitr ariness when c onsidering pe rceptual quality of sound, a s discussed in the introdu ction. But further investigations looked desirable. Other criteria will be considere d in f uture work so the present choice is only a first step to wards ph ysical/percep tual accounting in wind instrumen ts d esign. Notice also that the complex question of the nonlinear excitation mechanism, as the reed of an o boe or the lips of a br ass player e.g ., is d iscarded (equatio n (3) is h omogen eous), a lthough non linear effects can be imp ortant due to the nature of fluid dynamics, e .g. inside a reed or at the excitation localization in a flu te. The hypothe sis is made here that the n onlinear effects a re, in a first approx imation, en coded in the Fourier c oefficients o f the poten tial (see section III), c orrespon ding to an arbitrary and un known excitation . This is surely n ot satisfactory from a ph ysical p oint of view a nd will have to be modified carefully in the future . In the present context, where losses are a ssumed to v anish or to be inclu ded into the boun dary conditions, the energy inside the duct, note d E , is a c onserved quantity : E = T ( x, t ) + U ( x, t ) = π ρ 0 8 Z L 0 D 2 ( | φ x | 2 + | φ t c | 2 ) dx (4) with ∂ E ∂ t = 0 and T ( x, t ) , U ( x, t ) are compu ted in the appen dix. Following th e ab ove d efinition o f an oscillation regime, the energy can be taken as a design criter ion f or the duct shape. As it is a fu nction o f the diamete r , this o ne can be consider ed a s an unk nown fun ction, provided convenient data are given for the design. In that r espect, as m usical instrum ents are d esigned to o perate in harm onic regimes, the p otential, being in that case per iodic, is ame nable to a F ourier series analysis and ca n thus be written as a linear combinatio n of time- harmon ic elem entary componen ts, the co efficients of which will be other unknowns of the design problem s, while the frequ encies will be conside red as the data. This will be made more precise in section III. Th anks to the comp utations in the appen dix, the lin k is estab lished between an oscillation regime in the sense of A.H. Benade and the way the W ebster hor n is deri ved, as both rely up on the same definite integral, namely the lagr angian action, throu gh a Legendre transformation . As optimal con trol theory is variational to o, one has a coherent set of tools for the design method ology to be exposed. B. The wave ho rn equation as a contr ol system Suppose a dynam ical system is g i ven as : X ′ = f ( X , U ) (5) where X ( x ) ∈ R n is the state, U ( x ) ∈ [ U 1 , U 2 ] ⊂ R the contro l and a prime still denotes the deriv ati ve with respe ct to th e indepen dent variable x . The contr ol U is pur posedly restricted to be a scalar as th is suffices f or the need s of the present work (see [ 23], [8] for a mor e gen eral exp osition). Th e purpose of Optimal Control theory is to find a con trol U such that some function al (a scalar fu nction of the state and the contro l) : J ( X, U ) = Z x 1 x 0 g ( X , U ) dx (6) is maximiz ed or minimized , makin g it an op timization problem un der dynamic al constrain ts, nam ed for this reason infinite- dimensiona l. This v ariational-type p roblem is solved the f ollowing way : using a vector Lagrange multiplier µ , with compon ents µ i , i = 1 , . . . , n , d efine the Hamiltonian function of the pro blem, a scalar function, as : H ( X, U, µ ) = g ( X , U ) + µ T f ( X , U ) (7) One then obtain necessary cond itions fo r an optimu m as th e differential-algebraic system of equ ations ( D AE) :      X ′ = ∂ H ∂ µ − µ ′ = ∂ H ∂ X 0 = ∂ H ∂ U (8) together with suitable in itial/bound ary conditions. Notice first that, as H is a scalar, all of the thr ee above equations are vector o nes : ∂ H ∂ X e.g. is the vector with comp onents ∂ H ∂ X i , i = 1 , . . . , n . One has then a set of 2 n + 1 eq uations, the first 2 n ones being d ifferential equatio ns, th e last on e being of ord er z ero. Notice also th at these cond itions (8) are on ly necessary so that more in vestigations are needed to get sufficient cond itions, i.e. a complete and u nique solution . Also, the so-called Pontryag in maximu m princ iple [2 3], [8] allows to have mo re p recise results for the above wh en the control has e .g. bound constraints, as will be the case fo r the acoustic design problem here. Now , in o rder to fit the acoustic design problem at hand to this control theoretic setting, it is useful to reformulate the second or der equation (3) as a two dimensio nal first or der system, su itable fo r contr ol purp oses in view . T o this end , define the following variables : X 1 = D , Y = φ, Z = φ ′ (= v ) and W th e vector with co mpon ents X 1 , Y , Z . Also, as the sectio n diameter D of the horn is an u nknown to be determined as a f unction of x , it can b e co nsidered as a c ontrol variable to be designed. Actually , one can control either D itself o r the way it varies along the x axis, that is one can contro l its deriv ativ e D ′ with respect to x : as the modelling hypo thesis i) above is th at th e section is a slo wly varying function of x , on e n atural choice fo r the contro l is the deriv ati ve D ′ and bo unds will have to be imp osed to it in order to satisfy th e hypothesis given by in equality (1). As a consequen ce, d efining the control variable U as D ′ , con straints on it will be bou nd constraints : D 1 ≤ D ′ ≤ D 2 . T ogether w ith those on the c ontrol, constraints can also b e imposed to th e state through design/building constraints on the diameter D itself : an obvio us man datory state constraint is e.g. X 1 = D > 0 . After ha ving d ev elopped the second order deriv ativ es, and with these no tations, equation (3) rewrites as the first order d ifferential system :    X ′ 1 = U Y ′ = Z Z ′ = − 2 U X 1 Z − k 2 Y (9) which is a dyn amical control sy stem, affine in the co ntrol U : W ′ = f ( W , U ) = h 1 ( W ) + h 2 ( W ) U , with imm ediate definitions for h 1 and h 2 . On e has even a lin ear drift : h 1 ( W ) = AW , where matrix A is clearly defined . The singular ity in X 1 = 0 is not a real problem as it c orrespon ds to a vanishing d iameter, a hig hly u ninteresting situation (except po ssibly at one boun dary , th e apex of a com plete cone e.g.). Now , referring to theorem I , item 2 o f the app endix, for each eig en value λ n = k 2 n the eigenfunctio n ϕ n satisfies the wave ho rn eq uation with th e same bo undar y condition s. One can then re write for each eigenv alue λ n , the correspo nding first order d ifferential system of equ ations : ( X ′ 2 n = X 2 n +1 X ′ 2 n +1 = − 2 U X 1 X 2 n +1 − k 2 n X 2 n (10) with : X 2 n = ϕ n , X 2 n +1 = ϕ ′ n . Th is set of eq uations, i = 1 , . . . , to gether with the eq uation X ′ 1 = U constitutes the controlled dy namical model that w ill enter the optimizatio n p rocess as co nstraints. Remark th at on a practical side, only a finite numbe r of eigen values will be retain ed, as is detailed below in section III. One shou ld notice also that, whe reas in most contro l theoretic situations on e looks for a f eedback con trol, the searched after c ontrol for the ab ove ho rn equation is intrinsically op en-loop , as the du ct is designed on ce for all (at lea st in the pre sent state of techno logy . ..) so that optim al control th eory an d th e Pontryagin Maxim um Princip le [2 3], [8] are well suited in the present c ontext, althou gh giving only necessary cond itions, as alr eady mentionned . I I I . O S C I L L A T I O N R E G I M E A S A N O P T I M I Z A T I O N C R I T E R I O N Besides the ab ove simp le, albeit new , prob lem rewriting, the crux in the f ollowed a pproach of win d instrumen t design is in fo rmulating the oscillation regime, as defined in section II-A. It mu st be don e in a way suitab le for con ducting the optimization while ha ving the controlled dynamical model of section II-B above at hand . Considering the ab ove de finition for an o scillation r egime, the poten tial φ can be assumed to be a tim e-period ic functio n. For h armon icity requ irements of the sig nals in side th e d uct, the wav e numb ers k n s’ are fixed data given as m ultiples o f a fu ndamen tal fre quency : k n = nk 0 , k 0 = 2 π c f 0 . T o each of th e k n s’ correspon ds a solution ϕ n ( x ) of the W e bster equation with k n as wa ve num ber . The set { ϕ n ( x ) } n ∈ Z is a b asis of th e Hilber t space of squ are-integrab le fu nctions on [0 , L ] , thank s to Theo rem 1, item 3, append ix. Thus, thanks to Theor em 1 , item 4, append ix, th e overall potential φ can be written as : φ ( x, t ) = X n ∈ Z c n ϕ n e iω n t (11) which is seen to be the Fourier series deco mposition of φ ( x, t ) with respect to the time variable. Thus Parsev al theorem allows to write : | φ | 2 = X n ∈ Z c 2 n | ϕ n | 2 , | φ ′ | 2 = X n ∈ Z c 2 n | ϕ ′ n | 2 (12) The e nergy E can now be compu ted as a fun ction of c n , ϕ n , ϕ ′ n , simply sub stituting f or φ from (11) into (4), which gives, thanks to the relation s ( 12) : E = π ρ 0 8 Z L 0 D 2 X n ∈ Z c 2 n ( | ϕ ′ n | 2 + ω 2 n c 2 | ϕ n | 2 ) ! dx = π ρ 0 8 Z L 0 D 2 X n ∈ Z c 2 n ( | X 2 n +1 | 2 + k 2 n | X 2 n | 2 ) ! dx (13) the last equ ality stemming from th e d efinition of X 2 n , X 2 n +1 at the en d of section II. Now assume that the energy is conserved in side th e horn . At least, one can c onsider in a first app roxima tion th at the small par t of ene rgy which is r adiated outside the duct or d issipated at the inside boundar y is exactly co mpensated with the en ergy b rough t in by the excitation mechanism such as the ree d of the in strument. Such dissipation and o ther pheno mena sho uld surely be considered in future work, possibly throu gh the boun dary co nditions. The definition o f an oscillation r egime lead s then to m aximize E as d efined in equation (13). Th us, with the ω n s’ known an d fixed a s given d ata , th e optimization p roblem above has to find D ( x ) and the c n s’ . Observe ne vertheless that obtaining optimal c n s’ does not give any indication on some perceptive quality factor which is in deed importa nt for a hig h quality design . Thus it is likely at this stag e that optim izing the above-defin ed oscillation regime is unsufficient to that end. Section V will present simulatio ns for the above defined pro blem, after th e design prob lem itself h as been mathematically posed in the following sectio n. I V . O P T I M A L D E S I G N O F A H O R N S H A P E Thanks to the results of sectio n II an d II I, the d esign of axisym metric wind instru ments can be now fo rmulated as an optimal contr ol p roblem in the following way . Consider an oscillation regime as defined and developped in sectio ns I I-A, I II, each com ponen t bein g governed b y the W ebster horn equ ation with its own wa ve numb er k n and keepin g th e ap prox imation of the solution to a finite numb er, N , o f its first ter ms. The data ar e made of th e set of fixed multiples of a fundam ental frequen cy f 0 , leadin g to the set k n = j n k 0 , j n = 1 , . . . , N , with j n allowing to take into acco unt that a complete series of harmon ics o f f 0 or only the od d part of this series can be present. As φ is un known, the c n s’ of its time-ha rmonic decomp osition a re u nknown p arameters and co nsidered as constant contro ls. I n th e following, the c n s’ will be gather ed in a par ameter vector noted C = ( c 1 , c 2 , . . . , c N ) . Then the searche d after shape o f th e duct an d the c n s’ are o btained as the solution of the following pr oblem : ( P ) Max imize E as given by equatio n (1 3), with respect to D and C subject to : 1) The 2 N + 1 -dimension al dyna mical mo del ( n = 1 , . . . , N ) :    X ′ 1 = U X ′ 2 n = X 2 n +1 X ′ 2 n +1 = − 2 U X 1 X 2 n +1 − k 2 n X 2 n (14) each n th pair of the last 2 N eq uations being a W ebster equa tion with wa ve n umber k n . 2) Bound constraints on the variable c ontrol : D 1 ≤ U = D ′ ≤ D 2 . 3) State constraint : X 1 ≥ a > 0 , a a fixed real nu mber . It is apparent that, whenev er U > 0 , X 1 (0) > 0 , th e state constraint is satisfied. W ith in this co ntext, the stron g Pontr yagin Maximum Principle (PMP) [23], [8] is applicable. Necessary con ditions for an op timal control ( i.e. an optimal duct shape) are obtained in th e fo llowing way . Alo ng th e g eneral setting given by (5 ), (6), ( 7), ( 8), ad join the dyn amical constrain ts (14) to the criterion E , (13), through Lag range m ultipliers µ = ( µ n ) n =1 ,...,N , define the Hamilto nian of problem ( P ) as : H ( X, U, C, µ ) = π ρ 0 8 X 2 1 P N n =1 c 2 n ( X 2 2 n +1 + k 2 n X 2 2 n ) + µ 1 U + P N n =1 ( µ 2 n X 2 n +1 − µ 2 n +1 (2 U X 1 X 2 n +1 + k 2 n X 2 n )) (15) Then, following (8), µ is th e solution of the adjoin t d ifferential system ( n = 1 , . . . , N ) :          µ ′ 1 = − ∂ H ∂ X 1 = − π ρ 0 4 X 1 P N n =1 c 2 n ( X 2 2 n +1 + k 2 n X 2 2 n ) − 2 U X 2 1 P N n =1 µ 2 n +1 X 2 n +1 µ ′ 2 n = − ∂ H ∂ X 2 n = k 2 n ( µ 2 n +1 − π ρ 0 4 c 2 n X 2 1 X 2 n ) µ ′ 2 n +1 = − ∂ H ∂ X 2 n +1 = 2 U X 1 µ 2 n +1 − µ 2 n − π ρ 0 4 c 2 n X 2 1 X 2 n +1 (16) Thus th e o ptimization p roblem will b e co mplete when initial-b ound ary conditions are specified. Notice th at time- initial condition s are not rele vant here as the time variable has b een e liminated through th e hypothesis of time-h armonic regime and of en ergy conservation : the outp uts of th e optimiz ation pr ocess are sp ace-depe ndent, g iving e .g. the mo de shap es inside the d uct. As a starting point and in order to simply illustrate the metho d, witho ut g oing into deep physical co nsideration s of a specific instrume nt, that are po stponed to a fo rthcom ing pap er , th e end at x = 0 is assum ed to be c losed ( typically wh ere the excitation would be placed) so tha t the velocity and thus ϕ ′ vanish ther e : ∀ n = 1 , . . . , N , ϕ ′ 2 n (0) = X 2 n +1 (0) = 0 . Also the o ther end is assumed to be open (e.g. the oth er end o f a brass instrument) so that, a t the first order appro ximation , the pressure an d th us ϕ v anish : ∀ n = 1 , . . . , N , ϕ 2 n ( L ) = X 2 n ( L ) = 0 . These conditions lead to a T wo-Point Boundar y V alu e Problem (TPBVP) fo r the d ifferential system (14),(16) : one part of the state is fixed a t on e end and an other part is fixed at the other end. When some state comp onent X j is lef t u nspecified at one end , the correspond ing co state com ponent µ j must vanish there [23]. Thu s in th e pr esent situation : µ 2 n (0) = µ 2 n +1 ( L ) = 0 . Other cond itions can b e imp osed with more realistic considerations. For computing the op timal control, the Pontryag in max imum principle implies that, for ˆ X (resp. ˆ µ ) a solution of the state ( resp. ad joint state) eq uation, an o ptimal co ntrol ˆ U and o ptimal p arameter vector ˆ C are such that : ˆ H = H ( ˆ X , ˆ U , ˆ C , ˆ µ ) = max { D 1 ≤ U ≤ D 2 ] ,C = C st } H ( ˆ X , U, C, ˆ µ ) (17) But one can observe that H is affine with re spect to U so that : ∂ H ∂ U = µ 1 − 2 X 1 N X n =1 µ 2 n +1 X 2 n +1 = 0 (18) which is also named th e switching fun ction for this problem as its zeroes an d its deriv ati ve c an h elp determ ine the shap e of the duct. As equation (18) does not a llow to c ompute the contr ol U explicitely , one is faced with a problem with singula r extr emal ar c s [5], chap . 8. The stand ard method is to compute successive derivati ves of ∂ H /∂ U with r espect to x u ntil one is able to get an expression for U . In the present case, two such deri vati ves allo w this, after tedious but straightf orward computatio ns. In case U is bound ed, as it is the case her e, the m inimum of the hamiltonian H with respect to U , whe re ∂ H /∂ U = 0 , can occur at the b ounda ry of the dom ain ( see [5], [23]). Th is fact can b e observed in the simulatio ns : the control can b e o n its b ound alo ng several intervals of th e in tegration interval ( figure 4 ) o r can g o from o ne boun d to the o ther (figure 9), revealing a ban g-ban g type co ntrol [5]. Th e conc lusion is that on e g enerally ob tains quasi-con es. Con es could no t be obtained for exactly h armonic frequ encies [7], chap . 7 . This is compatible with the results obtain ed in [9]. N ow , thanks to the developments of this section , o ne has sufficient material to illustrate the de sign method through a fe w simu lations, on which som e q ualitative observations can b e mad e and co nfirmation of th e th eoretical resu lts is g iv en. Remem ber n ev ertheless that in th e following, only candidates for op timal sha pes are o btained at th is stage (necessary co nditions) b ecause, due to theorem 2, on ly q = 2 D ′ /D ca n be uniq uely determin ed fr om the given b ound ary spec tral d ata. Furth er data and theoretical in vestigations are necessary towards a really optimal shape . V . N U M E R I C A L R E S U LT S The numerical results presented here must be considered as qualitativ e illustratio ns of the above theoretical results. No interpretatio n in terms of mu sical quality will be attempted at this stage. Such results will be pur sued elsewhere mor e deeply in o rder to deri ve conclusion s o n realistic physical and pe rceptual basis. Ne vertheless, som e fea tures are worth noticing at this stage, when faced with some character istics of rea l instruments. A. Data In the following, data are as follows : the air mass density is taken to be ρ 0 = 1 , the sound velocity is that in free space at stan dard temperature : c = 340 m/s and on e fo cuses on a fun damental f requen cy for the note usually lab elled A 4 , i.e. f = 44 0 H z , lea ding to th e wa veleng th λ = c/f = 0 . 772 m . The ch osen wa ve number s are then indicated un der each figure : k 1 = 2 π /λ, k i = n i k 1 , i = 1 , . . . , n a nd n i an even or od d integer . At the narrow end, x = 0 , o ne has D (0) fixed with different v alues for each simu lation (see b elow) an d φ ′ (0) = 0 . At x = L , the co ndition φ (0) = 0 is satisfied o nly approx imately , by adjunction to the criter ion thr ough a simple penalty meth od. Th e du ct diam eter derivati ve is allowed to vary in the bound ed interval [ D ′ m , D ′ M ] , D ′ m , D ′ M being indicated on each figu re. They are chosen to satisfy inequality (1) with R = D/ 2 , for the maximu m value k of the k n s’, i.e. a re within the validity domain of the wave horn eq uation, accordin g to (1). The in terval of integration, i.e. the length of th e duct, is taken to be L = λ but could be taken as an unknown as well, leadin g to an an alogou s of the so- called time-op timal co ntrol. La st, all unspe cified initial values for the state variables and all in itial controls are taken to be rando m variables uniformly distributed on [0 , 1 ] . B. Numerical method The n umerical metho d used to solve the optim ization problem exposed in section IV is the f ollowing (see e.g. [5 ], chap. 7 for more details on several possible numerical methods used for optimal control pro blems) : an initial guess is giv en for the control vector and , using this control, the state equation is numerically integrated fro m X (0) , con sidering too the unsp ecified in itial condition s as co nstant contro ls, which are to be d etermined with the optimizatio n process. The integration was perfo rmed thank s to a predictor-correcto r scheme with an explicit Eu ler m ethod for the p rediction step an d a Crank-Nicholson schem e fo r the correction . Then the co state equ ation is integrated b ackwards -no tice that this eq uation is always linea r in the costate- using the obtained terminal values for X ( L ) . Th is allows to co mpute the gradien t of the objective function with respect to the control, thanks to ∂ H / ∂ U [8 ]. Then a standard optimiz ation rou tine is used, passing throug h the previous steps in a r ecursive way , in o rder to make this gr adient decrea se, until the specified term inal cond itions at x = L are -appro ximately here- achieved and the criter ion e v olves x no more, to a spec ified precision ( 10 − 5 ). For the time being, the criterion and the termina l constrain ts have been g athered in a sing le objecti ve function thr ough a simple penalty meth od but this could be imp roved. For th e op timization, a quasi-Ne wton method with projectio n has been used together with a BFGS meth od fo r a n estima tion of the hessian. Such an app roach is k nown to have, as usually Newton-type methods, the dr awback that it can give local min ima and to be sensiti ve to the initial guess for the solution [10]. Thu s, other methods su ch as direct solution of the TPBVP by shooting o r multiple shooting techn iques should b e interesting to in vestigate but this is defer red to futur e work because an extra work has to b e done in o rder to ha ve an explicit expre ssion for the co ntrol, as a f unction of the state and c ostate (see th e discussion in sectio n IV above). Nevertheless sev eral d ifferent initial co nditions for the state X we re tested : the results, not exposed here, did no t show to be that sensitiv e but this should be confirmed theoretically . C. S imulations The few simulations presente d here have been done while fixing, respectively , two, five an d ten co mpon ents for the overall potential inside the du ct, i.e. N = 2 , 5 , 10 in e quation ( 14). In each case, the duct shap e is sho wn first, follo wed by the correspo nding m odes shapes that hav e been norm alized to their m aximum value, at x = 0 . In addition fo r the case o f two compon ents only , the d iameter deriv ati ve h as been shown to illustrate the fact that it ca n be on ly p iecewise continu ous and that th e phen omeno n, men tionned in section IV, of singular extrem al arcs joining regular ar cs can app ear (see figure 4) : the control U = D ′ is o n its bo und on some sub intervals. For the fi ve compo nents cases, th e shape is also quasi-con ic but in a less obvious manner . But for the te n compon ent case, on e has an example of abrupt change in conicity and D ′ in figure 9 illustrates the p ossibility of ban g-bang type con trol mentionn ed in section I V. On each figure the impo sed lower ( D 1 ) and upper ( D 2 ) bound s are indicated. Notice that the imped ance, p/ u and the instantan eous p ower , pu ( p , the pressure, u th e volume velocity), for all x are e asily obtained from the outputs, ϕ ( x ) , ϕ ′ ( x ) and D ( x ) ∀ x ∈ [0 , L ] , of the op timization process. The simulations are as fo llows : 1) T wo com ponen ts : For this first simu lation, the fun damental frequen cy is A 4 and the seco nd compo nent is the second harmon ic (dou ble freq uency). The diameter a t the closed end is D (0 ) = 2 cm . The du ct shape is shown in figure 2, the two m odal shap es in figure 3. D ′ is shown in figure 4 as it makes appear the phen omeno n describe d in section IV, where th e control is at the boun d on two subinter vals, with the consequence tha t o ne ob tains a quasi-co ne for the duct, made of conic pieces joined by m ore complex but smooth parts. 2) Fiv e componen ts : In this second simulation (see fig ures 5,6), the data are the same as for two com ponen ts, except for D (0) = 1 cm . Th e fiv e co mpone nts have frequenc ies f i = ( i + 1) ∗ f 0 ; i = 0 , . . . , 4 . One has no t shown D ′ here but the result is here again a quasi-con e, although this is not apparent again on fig ure 5, because D ′ varies much less but still varies along the duct. 3) T en compon ents : In this th ird simulation (see figures 7,8,9), the data is the same as for fi ve c ompon ents, except for D (0) = 5 cm . T he ten com ponen ts have f requenc ies f i = ( i + 1) ∗ f 0 ; i = 0 , . . . , 9 . The results are noticeable as her e again, three co nic pieces are found so that th e deriv ati ve D ′ is shown in figure 9, making it here again a quasi-cone, with a bang- bang type behaviour f or D ′ . T he du ct is made o f conic pieces with different conicity joined tog ether . Only the five first modal shapes have been d isplayed, for better readability . D. Qua litative observations and commentaries The above simulatio ns sh ow o nly q ualitative albeit importa nt re sults at the p resent stage. The main one is that, as f oreseen by the th eoretical in vestigations, and fo r the specific ch osen criterion , the duct shapes are pie cewise continu ous, th is bein g particularly clear on the c ase of two and ten co mponen ts ; o bserve too in both cases the fact that the control is at th e bound on som e subin tervals, m eaning th at regular arcs and sin gular arcs can coexist. One can see also that qua si-cones are obtained, in an obvious way for th e first and second simulations. Recall nevertheless that, as mentionn ed, the numerical method can be sen siti ve to initial conditions in the optim ization process, so that it is difficult to interpre t these re sults in a precise mann er at the present stage. The qualitative obser vation is that the shap es are quasi-con es. At the same time, the modal shapes (i.e. the eigenfunction s of the Sturm-L iouville p roblem) show a beh aviour that is qualitatively in a greemen t with what is theoretically pr edicted fo r co nes ( see [7], chap. 7), i.e. one has m odal shap es that decrease as x increases to L . On o ne hand, th ese r esults are c oheren t with what is kn own for such geo metries [ 9]. On another han d, the fact that one gets q uasi-con es essentially would lend to conclude that the used n otion of oscillation regime does not allo w to captu re all the im portant featu res that are sou ght after fo r musical quality , a s flaring horns in brasses e.g . canno t be achieved with it : significant and impo rtant nonlin ear phenom ena [7] are no t taken into account with th e here r etained mo del and desig n criterion. Thus, using more realistic mod els an d refined criteria should be used instead . The computed modes have an extremum at the closed end, coherent with the ch osen bound ary con dition ( φ ′ (0) = 0 ). At the other ope n e nd, the modes do es no t vanish exactly : one reason for this is linked to th e way this end cond ition has b een taken into acco unt in the numer ical metho d : a simple pen alty method makes a compro mise between satisfying the criter ion and the bo undary con ditions there. The above suggested method of dire ctly solving the TPBVP thro ugh shoo ting techniq ues should giv e better results from this poin t of vie w . V I . C O N C L U S I O N A N D P E R S P E C T I V E S Thanks to a refor mulation within a control theory framework, the design of axisym metric wind instruments has been revisited. It has been shown that in order for an oscillation regime in the sense of A.H. Benade to take place inside the air column an d with the con sidered linear model, the shape tend s to b e c onic, as the nu mber of time-harmon ic componen ts grows. Th is ten ds to show that such a criterion, as it h as been mathe matically d escribed, is un sufficient to g rasp the necessary cond itions that lead to h igh-qu ality instrum ents. No other conclu sion f rom this point of v iew is given at this stage . Nev ertheless, and th is was o ne main purpo se o f this work, the approac h m akes it very flexible to de al with a g reat variety of design constraints, e ither on the co ntrol or o n the shap e itself. At th e present time, on ly q ualitative results have been giv en to illustrate the th eoretical results and th e field is wide open to numerical investigations as well as to experim ents on more physical an d musical premisses. I n that respect, some impo rtant issues can be straigh tforwardly put in p erspective : 1) A first possible development is to u se the meth od in a mod el fitting way , using an exper imentally m easured input imped ance fr om a real instrumen t an d taking a distance me asure with the corr espondin g model outp ut as cr iterion, in a similar way as in [18], [22]. Mor e realistic duct sh apes are to be expected and it will b e interesting to compare results fro m different approache s. 2) Refine the question of physical/perceptu al criteria and of the horn m odel : this relies upo n phy sical consider ations as well as on experiments with real existing instruments play ed by expert musicians. This is surely a central point to get a practical tool useful for real design. One consequence of modifying the criterion is that the structur e of the obtain ed sh ape will likely be different f rom that obtained in the p resent study . 3 ) De velop the method for woodwinds, thereby including toneholes in the design process. 4) T ake into account some imposed-shape d parts such as pieces of cylinders or cones fr om built-in comp onents, such as fo r brasses. This imp lies using path constraints in th e op timal co ntrol fr amework. 5 ) All the previous points depen d on the modelling question of how to include the nonlinea r excitatio n mech anism (ree d, lips...) in th e model as well as o ther nonlinear p henom ena, which were discarded in th e present stud y . 6) Last, a n interesting question is to try to retrieve, with the pr esented approa ch, the lattice o f sou nd tubes r esults fou nd in [ 9]. This will imply to modify th e mathematical framew ork o f th e controlled dy namical system, i.e. the func tion spaces where the state an d the con trol live, in order to be able to take into acco unt discon tinuous, and not only piece wise continuous, duct shapes. Th e mathematical theory of distributions co uld be an adequ ate framework for th is. T hese item s are but a small part of wh at sho uld be add ressed within the prop osed fra mew ork to achiev e a satisfyin g design method ology . V I I . A C K N O W L E D G E M E N T S The au thor is u ndebted to two ano nymous revie w ers an to the ed itor of the spe cial issue, J. Kergormard, for their num erous advices that helped greatly imp rove th e paper . R E F E R E N C E S [1] Musical acoustics network. Summ er meeting : wind instruments acoustics. http://www .music.ed.ac.uk/euch mi/man/mxhta.html , 2009. Edinbur gh. [2] A.H. Benade. Fundamental s of musical acoustics . Dov er , 1990. 2d edition. [3] H. Bouasse. Instrumen ts ` a vent, I, II . Delagra ve , 1929. New edition L ibrairi e A. Blanchard, 1986. [4] M. Bruneau. Manuel d’acoustique fondamentale . Etudes en m ´ ecani que des m at ´ eriaux et des structure s. Herm ` es, 1998. [5] A.E . Bryson and Y . C. Ho. Appli ed optimal contr ol . Hemisphere Pub . Corp., 1975. Re vised printing. [6] D.M. Campbell and C. Greated . The musician’ s guide to acoustics . Oxford Univ ersit y Press, 1987. [7] A. Chaigne and J. Ke rgomard. Acoustique des instruments de musique . Bel in, 2008. [8] F .C. Clarke. Optimization and nonsmooth analysis . W ile y , New Y ork, 1983. [9] J.P . Dalmont and J. Ker gomard. L attic es of sound tubes with harmonicall y related eigenfreq uencie s. Acta A custica , 2(5):421–4 30, 1994. [10] J.E . Dennis and R. B. Schnabel. Numerical methods for unconstr ained optimizat ion and nonline ar equation s . Prent ice Hall, Englew ood Cliffs, New Jersey , 1983. [11] J. Dieudonn ´ e. F oundations of modern analysis . Academic Press, 1960. [12] E. Eisner . T he complete solution of the ”webster” horn equation. J . Acoust. Soc. America , 41:1126–1 146, 1967. [13] N.H. Fletcher and T . Rossing. The physics of musical instruments . Sprin ger , 1991. [14] B.J. Forbes, E. Roy Pike, D.B. Sharp, and T . Aktosun. In ve rse potential scatterin g in duct acoustic s. J . Acoust. Soc. Am. , 119(1):65–73, 2006. [15] L. Henrique and J. Antunes. Optimal design and physical modelling of mallet percussion instruments. Acta A custica u. with Acustica , 89:948–962, 2003. [16] O. In ´ acio, L. Henrique, and J. Antunes. Design of duct cross section al areas in brass-trapp ing resonators for control rooms. Noise Contr ol Eng. J. , 55(2):172– 182, 2007. [17] A. Katchalo v , Y . Kuryle v , and M. Lassas. In verse boundary spectral prob lems . Cha pman & Hall/CRC, 2001. [18] W . Kausel. Opti mizatio n of brasswind instruments and its applic ation in bore reconstruct ion. J. New Music Researc h , 30(1):69–82, 2001. [19] W . Kausel, P . Anglmayer , and G. W idholm. A computer program for optimizat ion of brass instruments. In F orum Acusticum , Berlin, 1999. [20] G. Le V ey . Optimal control as a tool for the design of wind instruments. In Pr oc. International symposium on musical acoustic s, ISMA07 , Barcel ona, sept. 2007. [21] G. Le V ey . Optimal contro l theory and design of wind instruments. In Pr oc. IEEE Mediterrane an confer ence on contr ol, MED08 , Ajac cio, june 2008. inv ited conference . [22] D. Norela nd. Numerical T echnique s for A coustic Modelling and Design of Brass W ind Instruments . Acta uni ver sitatis upsaliensis, Uppsala Uni versit y , 2003. Comprehensi ve Summ aries of Uppsala Dissertati ons from the Faculty of Science and T echnol ogy , 862. [23] L. Pontriaguine, V . Boltianski, R. Gamkrelidze, and E. Mitchenk o. Th ´ eorie math ´ ematiqu e des proce ssus optimaux . Mir , Mosco w , 1978. [24] G.R. Putland. Model ling of horns and enclosur es of loudspeaker s . PhD thesis, Univ . Queensland, 1994. [25] S.W . Rienstra. W ebster’ s horn equati on re visited . SIAM J. Appl. Math. , 65(6):1981–2004, 2005. A P P E N D I X a) Th e horn wave equatio n : In stead of the usual, Newton-type mo delling appr oach, mak ing use of f orces an d m oments balance equations, a lagran gian o ne is fo llowed here, a s it fits better the v ariational co ntext for the design methods that form th e core of the p resent work. For stand ard deriv ation the r eader is referred to [24], [4] for phy sical co nsideration s and to [20] for contr ol the ory formulation . For a g eneral pressure field, the lagr angian a ction density at ea ch time in stant inside the duct is the difference between a kinetic term and a potential term, that are co mputed in a standar d way , as follows. Assumption ii) implies that the du ct can be consider ed as a c ontinuo us stack of cro ss-sections S ( x ) , p arametrized by the abcissa x . For each section S ( x ) located at x alo ng the h orn ax is, an actio n density is comp uted as the integral of den sities of the par ticles over the section. This leads to an expression propor tional to the cro ss-section ar ea, i.e. to D 2 ( x ) . Firstly , the kinetic term writes : T ( x, t ) = Z S ( x ) 1 2 ρ 0 v 2 dσ = π D 2 4 ρ 0 2 v 2 = π D 2 4 ρ 0 2 | φ x | 2 (19) Similarly , the poten tial e nergy term is g iv en as : U ( x, t ) = Z S ( x ) p 2 2 ρ 0 c 2 dσ = π D 2 4 p 2 2 ρ 0 c 2 = π D 2 4 ρ 0 2 | φ t c | 2 (20) Thus the lagran gian action d ensity writes : L ( x, t ) = T ( x, t ) − U ( x, t ) = π D 2 4 ρ 0 2 ( | φ x | 2 − | φ t c | 2 ) (21) Eventually , the lag rangian action inside the du ct writes : L = Z t 1 t 0 Z L 0 L ( x, t ) dxdt = π ρ 0 8 Z t 1 t 0 Z L 0 D 2 ( | φ x | 2 − | φ t c | 2 ) dxdt (22) According to Hamilton ’ s station ary action pr inciple, the dynamics in side the duct is obtained as the Eu ler e quation of th e above action L : ∂ L ∂ φ − ∂ ∂ x ( ∂ L ∂ φ x ) − ∂ ∂ t ( ∂ L ∂ φ t ) = 0 (23) But one can see that ∂ L ∂ φ = 0 , thus : − ∂ ∂ x (2 D 2 φ x ) + 2( D 2 /c 2 ) φ tt = 0 (24) and eventually , dividing by 2 D 2 : 1 c 2 φ tt − φ xx − 2 D x D φ x = 0 (25) which is recog nized to be the wa ve h orn eq uation or W ebster horn equation . The interest of the a bove variational de riv ation for this eq uation lies first in the fact th at the optim al control appro ach to the design is variational in nature too. Thu s, d ue to the harmo nic n ature of wa ves inside musical instru ments, the horn equation writes, for on e fixed value ω : φ xx + 2 D x D φ x + ω 2 c 2 φ = 0 (26) i.e. using primes fro m now on to deno te the deri vati ves with respect to the spatial dimension x along the axis, the only remaining indepe ndent variable : φ ′′ + 2 D ′ D φ ′ + k 2 φ = 0 (27) this equation being adju ncted a set of suitable b ound ary con ditions : a 1 φ (0) + b 1 φ ′ (0) = 0 , a 2 φ ( L ) + b 2 φ ′ ( L ) = 0 (28) which r epresents a resonator without acti ve co mpon ents and for wh ich lo sses can be taken into account in th e bo undary condition s that will b e precised for the simulations in section V. Thus one is faced with a homo geneou s Sturm-Lio uville problem . b) Ele ments of Sturm-Liouville theory : T wo theo rems from spectral theory of Sturm-Liouville problems are rec alled here for self-co ntainedn ess : one [11] for the so-called “ direct problems” and the second [17] for “inverse prob lems”. Let I = [ O , L ] and q ( x ) = 2 D ′ D . Then : Theorem 1 [11] : For every func tion q ( x ) contin uous in I : 1) The Sturm-Liouv ille pro blem has an in finite strictly increasing s equenc e of eigen values λ n ∈ R such that lim n →∞ λ n = + ∞ and the series P n 1 /λ 2 n conv erges. 2) For each eigenv alue λ n , the ho mogen eous Sturm- Liouville problem has a real-valued solution ϕ n ( x ) such that R b a ϕ 2 n ( x ) dx = 1 , which is u nique up to a multiplica ti ve real constan t. 3) The sequence ( ϕ n ) is an ortho norma l system in a convenient Hilbe rt space of functio ns. 4) Let w be co mplex-valued continuou s function d efined in I , the pr imitiv e of a r uled f unction w ′ such that : (i) w ′ is continuous in I , except possibly at a finite number of inter ior p oints. (ii) w ′ has a deri vati ve w ′′ continuo us in ev ery interval where w ′ is continu ous. (iii) w satisfies the bou ndary condition s in (28). Then, if c n = < w , ϕ n > = R b a w ( s ) ϕ n ( s ) ds , on e has : w ( x ) = P n c n ϕ n where the series converges unifor mly and abso lutely in I . On anoth er h and, the design question itself relies upon the following in verse bound ary spectr al theorem : Theorem 2 [17] : Assume that { λ 1 , λ 2 , . . . , ϕ ′ 1 (0) , ϕ ′ 2 (0) , . . . } a re th e bo undar y spec tral da ta of the Dirichlet-Sch r ¨ odinger operator, A 0 = − d 2 dx 2 + q , corr espond ing to the above Sturm-Liou ville problem, on an interval [0 , L ] . Then , these data determine L and q ( x ) uniquely . It is moreover worth no ting that only q can be determin ed u niquely from these da ta. This means here that only the ratio D ′ /D is so. T hus to determin e D itself needs supple mentary d ata. 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 −0.08 −0.06 −0.04 −0.02 0.00 0.02 0.04 0.06 0.08 P S f r a g r e p l a c e m e n t s x (in m) Duct shape D (in m) Fig. 2. Duct shape, 2 components, − 0 . 2 ≤ D ′ ≤ 0 . 2 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 −0.4 −0.2 0.0 0.2 0.4 0.6 0.8 1.0 P S f r a g r e p l a c e m e n t s ϕ 1 ( x ) , ϕ 2 ( x ) x (in m) ϕ n Fig. 3. Modal shapes, 2 components, − 0 . 2 ≤ D ′ ≤ 0 . 2 0 50 100 150 200 250 300 350 400 450 500 550 0.170 0.175 0.180 0.185 0.190 0.195 0.200 0.205 P S f r a g r e p l a c e m e n t s x (in m) D ′ ( x ) D ′ Fig. 4. Deri v ati ve of the duct diameter , 2 components, − 0 . 2 ≤ D ′ ≤ 0 . 2 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 −0.04 −0.03 −0.02 −0.01 0.00 0.01 0.02 0.03 0.04 P S f r a g r e p l a c e m e n t s x in m Duct shape D in m Fig. 5. Duct shape, 5 components, − 0 . 2 ≤ D ′ ≤ 0 . 2 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 −0.6 −0.4 −0.2 0.0 0.2 0.4 0.6 0.8 1.0 P S f r a g r e p l a c e m e n t s ϕ i ( x ) , i = 1 , . . . , 5 x (in m) ϕ n Fig. 6. Modal shapes, 5 components, − 0 . 2 ≤ D ′ ≤ 0 . 2 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 −0.08 −0.06 −0.04 −0.02 0.00 0.02 0.04 0.06 0.08 P S f r a g r e p l a c e m e n t s x in m Duct shape D in m Fig. 7. Duct shape, 10 components, − 0 . 2 ≤ D ′ ≤ 0 . 2 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 −0.8 −0.6 −0.4 −0.2 0.0 0.2 0.4 0.6 0.8 1.0 P S f r a g r e p l a c e m e n t s ϕ i ( x ) , i = 1 , . . . , 5 x (in m) ϕ n Fig. 8. Modal shapes, 10 components, − 0 . 2 ≤ D ′ ≤ 0 . 2 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 −0.25 −0.20 −0.15 −0.10 −0.05 0.00 0.05 0.10 0.15 0.20 0.25 P S f r a g r e p l a c e m e n t s x (in m) D ′ ( x ) D ′ Fig. 9. Deri v ati ve of the duct diamete r , 10 components, − 0 . 2 ≤ D ′ ≤ 0 . 2