On Sampling Methods for Inverse Biharmonic Scattering Problems in Supported Plates

On Sampling Metho ds for In v erse Biharmonic Scattering Problems in Supp orted Plates Carlos Borges ∗ Rafael Ceja-A y ala † P eter Nekrasov ‡ Marc h 24, 2026 Abstract W e study the inv erse problem of qualitativ ely recov ering a supp orted cavit y in a thin elastic plate go verned by the flexural (biharmonic) wa v e equation, using far-field pattern measuremen ts. W e deriv e a reciprocity principle and a factorization of the far- field op erator for the supp orted plate b oundary conditions, and we analyze its range prop erties to justify b oth the linear sampling metho d (LSM) and the direct sampling metho d (DSM). Numerical exp erimen ts assess the p erformance of LSM and DSM under noise, a limited amount of data, multiple scattering, and v ariations in the Poisson ’s ratio. The results show that both metho ds robustly recov er the obstacle’s location and con vex h ull, with DSM offering improv ed stabilit y and reduced computational cost. Keyw ords: Biharmonic scattering, flexural wa v es, supp orted plates, recipro city relations, Kirc hhoff–Lov e plate theory , in verse problems. 1 Intro duction In recen t years, biharmonic problems hav e garnered significant attention due to their appli- cations in the analysis of resonances in bridges and large structures [ 14 , 30 , 31 , 44 ], simulation of acoustic “black holes” [ 39 ], dev elopment of top ological wa ves in elastic systems [ 20 , 36 ], nondestructiv e testing of aerospace comp osites [ 34 ], design of platonic crystals and metama- terials [ 16 – 18 , 21 ], and the study of ice shelf flexure and ice-ocean interactions [ 3 , 4 , 37 , 38 ]. W e consider the inv erse problem of reconstructing the shape of a supp orted ca vity in a thin elastic plate from far-field scattering data. Letting the ca vity o ccup y a compact domain D ⊂ R 2 with smo oth b oundary ∂ D , the vertical displacement u of the plate satisfies the follo wing time-harmonic equation:          ∆ 2 u − k 4 u = 0 , in R 2 \ D , u = 0 , on ∂ D , ν ∆u + (1 − ν ) ∂ 2 n u = 0 , on ∂ D , (1.1) ∗ Sc ho ol of Data, Mathematical, and Statistical Sciences, Universit y of Central Florida, Orlando, FL 32816 Email:Carlos.Bor ges@ucf.e du † Sc ho ol of Mathematical and Statistical Sciences, Arizona State Univ ersity , T empe, AZ 85287 Email:r c ejaaya@asu.e du ‡ Committee on Computational and Applied Mathematics, Univ ersity of Chicago, Chicago, IL 60637 Email:pn3@uchic ago.e du 1 where ∆ denotes the t wo-dimensional Laplacian, ∂ n denotes the directional deriv ativ e with resp ect to the normal v ector n , ν ∈ [ − 1 , 0 . 5] is the Poisson’s ratio of the plate, and k > 0 is the w av en umber, which dep ends on the thickness of the plate and the time frequency of the inciden t field. The first b oundary condition corresp onds to the displacement of the plate, while the second b oundary condition corresp onds to the b ending moment. These b oundary conditions are frequen tly used to mo del plates that are supp orted by internal columns and ro ds [ 41 , 45 ] or for ice shelv es supp orted by land at the grounding line [ 35 ]. Previous w ork on biharmonic in verse scattering has fo cused on the recov ery of v olumet- ric p oten tials from b oundary and far-field data [ 12 , 27 , 28 , 33 , 42 , 43 ] or the reconstruction of obstacles with clamp e d plate b oundary conditions [ 9 , 10 , 23 , 24 ]. Many studies hav e fo cused on the global uniqueness problem for the recov ery of potentials or other material parameters asso ciated with biharmonic op erators [ 6 , 28 , 32 ]. Related w ork has also b een done to un- derstand the reconstruction of p oten tials and internal sources using passiv e measurements generated by unknown sources [ 11 , 22 ]. In the present study , we consider metho ds for determining the lo cation and shap e of a supp orted cavit y from far-field measurements of flexural w av es. In particular, we con- sider t wo sampling metho ds: the line ar sampling metho d (LSM) and the dir e ct sampling metho d (DSM). The LSM reconstructs the cavit y b oundary b y solving a far-field equation for synthetic inciden t fields, while the DSM offers a computationally efficient alternativ e via the application of the far-field op erator. Our analysis pro vides a rigorous justification of b oth the LSM and DSM by examining b oth the range and recipro cit y of the far-field op erator for the supp orted plate problem. Numerical simulations further demonstrate that b oth metho ds can accurately and robustly reco ver the shap e and p osition of the supported ca vity , ev en with a limited amount of data and noisy measurements. W e present sev eral new analytical and computational results for the in verse biharmonic scattering problem. First, w e present a recipro cit y relation that extends the clamp ed plate recipro cit y result of [ 25 ] to the supported plate setting. Next, we sho w that the w ell- established factorization of the far-field op erator holds for supp orted plates, pro viding a rigorous foundation for the LSM. W e also in tro duce a factorization of the data-to-pattern op erator based on the b oundary-in tegral form ulations of [ 19 , 38 ], allowing us to study the prop erties of the DSM. Finally , we p erform a comprehensive n umerical study of the fre- quency dep endence and robustness to noise and other parameters of the LSM and DSM. The pap er is organized as follo ws. In Section 2 , we form ulate the forw ard problem and pro vide definitions of the fundamental solution and far-field pattern. Section 3 fo cuses on the linear sampling metho d, establishing the recipro cit y principle and other results ab out the far-field op erator. In Section 4 , we fo cus on the direct sampling metho d b y deriving a factorization of the data-to-pattern op erator and describing the indicator function. Section 5 rep orts numerical exp eriments that illustrate the p erformance and h ighligh t the limitations of the metho ds across a range of scenarios. Finally , Section 6 offers concluding remarks and directions for future work. 2 2 Problem setup Let u b e the sum of a known inciden t field u i and an unknown scattered field u s . If u i satisfies the homogeneous equation in R 2 \ D , we ha ve that u s satisfies          ∆ 2 u s − k 4 u s = 0 , in R 2 \ D , u s = − u i , on ∂ D , M [ u s ] = − M [ u i ] , on ∂ D , (2.1) where the b ending momen t operator M is defined b y M [ u ] = ν ∆u + (1 − ν ) ∂ 2 n u, (2.2) and u s satisfies the Sommerfeld radiation condition at infinity: ∂ r u s − ik u s = o ( r − 1 / 2 ) , as r → ∞ . (2.3) Let Φ ( x, y ) denote the fundamental solution of the biharmonic wa ve equation satisfying ∆ 2 Φ ( x, y ) − k 4 Φ ( x, y ) = δ ( x − y ) , in R 2 , along with the Sommerfeld radiation condition ( 2.3 ). The formula for Φ is given b y: Φ ( x, y ) = i 8 k 2  H (1) 0 ( k | x − y | ) − H (1) 0 ( ik | x − y | )  , where H (1) 0 is the zeroth-order Hank el function of the first kind. Note that b ecause H (1) 0 ( ik | x − y | ) decays exp onen tially aw a y from any source(s), the asymptotic b eha vior of the funda- men tal solution is similar to that of H (1) 0 ( k | x − y | ) : Φ ( x, y ) = e i k | x | p | x |   Φ ∞ ( ˆ x, y ) + O 1 | x | !   as | x | → ∞ , where ˆ x : = x | x | ∈ S 1 and Φ ∞ ( ˆ x, y ) denotes the far-field pattern of the fundamen tal solution, giv en by: Φ ∞ ( ˆ x, y ) = − 1 2 k 2 e iπ / 4 √ 8 π k e − ik ˆ x · y , (2.4) The scattered field u s admits a similar expansion at infinit y: u s ( x ) = e i k | x | p | x |    u ∞ ( ˆ x ) + O 1 | x | !    as | x | − → ∞ , (2.5) where u ∞ ( ˆ x ) denotes the far-field pattern corresp onding to the scattered field u s . In this pap er, we will often use u s ( x, d ) to denote the solution of ( 2.1 ) for an inciden t plane wa ve in the direction d ∈ S 1 . Similarly , u ∞ ( ˆ x, d ) represen ts the far-field pattern of u s ( ˆ x, d ) . Finally , we define the far-field op erator F : L 2 ( S 1 ) → L 2 ( S 1 ) as: F [ g ]( ˆ x ) = Z S 1 u ∞ ( ˆ x, d ) g ( d ) d s ( d ) . (2.6) 3 This op erator will play a prominent role in b oth the linear sampling and direct sampling metho ds. In the next tw o sections, w e introduce and analyze the sampling metho ds considered in this pap er. W e start with the line ar sampling metho d (Section 3 ) follow ed b y the dir e ct sampling metho d (Section 4 ). The fundamen tal principle of b oth sampling metho ds is to determine if a giv en sampling p oin t z ∈ R 2 lies inside or outside the domain D by ev aluating a suitable indicator function. F or the linear sampling metho d, this indicator function is given b y the solution to a particular integral equation whic h b ecomes un b ounded for p oin ts not in D . Mean while, the direct sampling metho d uses an indicator function that comes from applying the far-field op erator directly to the fundamental solution. As we will see, this indicator function is more stable to ev aluate and leads to b etter numerical stabilit y . Using either of these indicator functions, one can infer an approximation of the obstacle for the supp orted plate problem. 3 Linea r sampling metho d The goal of the linear sampling metho d (LSM) is to find a function g z ∈ L 2 ( S 1 ) for each z ∈ R 2 suc h that F [ g z ]( ˆ x ) = ϕ z ( ˆ x ) , (3.1) where ϕ z ( ˆ x ) : = Φ ∞ ( ˆ x, z ) is the far-field pattern of the fundamental solution. F or the LSM, the L 2 norm of g z acts as an indicator function for the obstacle. In order to see this, we presen t the theoretical justification b ehind the linear sampling metho d (LSM). W e b egin by establishing a recipro city principle for the biharmonic scattering problem with supp orted plate b oundary conditions in Section 3.1 . Then, in Section 3.2 w e extend the factorization of the far-field op erator to the supp orted plate and pro ve that it is compact, injective, and has dense range, k ey comp onen ts in deriving the LSM for our problem. Finally , we presen t the result that s ho ws that the v alue of the indicator go es to infinity for sampling p oin ts outside the obstacle and remains b ounded for p oin ts inside. 3.1 Recip ro cit y p rinciple In this section, we derive a recipro cit y relation for the solution u s to the b oundary v alue problem ( 2.1 )-( 2.3 ). Similar recipro cit y principles are well kno wn for acoustic scattering [ 13 ] and for the clamp ed plate problem [ 25 ] but hav e not b een studied for the supp orted plate problem. These identities are used in in verse scattering to study the corresp onding far– field [ 24 ] or near–field [ 9 ] op erators whic h are used in qualitativ e metho ds for recov ering the obstacle. Theorem 3.1. L et u ∞ ( ˆ x, d ) denote the far-field p attern define d in ( 2.5 ) for the supp orte d plate pr oblem ( 2.1 ) with incident field u i = e ikx · d , wher e d ∈ S 1 . Then u ∞ ( ˆ x, d ) satisfies the r e cipr o city r elation: u ∞ ( ˆ x, d ) = u ∞ ( − d, − ˆ x ) . Pr o of. W e b egin by noting the standard integration by parts form ula [ 26 ]: Z R 2 \ D v ∆ 2 u − u∆ 2 v d x = Z ∂ D  v N [ u ] − uN [ v ] − ∂ n v M [ u ] + ∂ n uM [ v ]  d s (3.2) 4 where M is the b ending momen t op erator defined in ( 2.2 ) and N is the shear stress op erator giv en by: N [ u ] : = ∂ n ∆u + (1 − ν ) d d s ∂ n ∂ τ u, where ∂ τ denotes the tangential deriv ativ e and d d s denotes the arc-length deriv ativ e. The pro of pro ceeds as in [ 25 ]. Letting u = Φ ( x, y ) and v = u s ( x, d ) , and taking the limit as | x | → ∞ , one obtains the following expression for the far-field pattern: u ∞ ( ˆ x, d ) = Z ∂ D  u s ( y , d ) N [ u i ( · , − ˆ x )]( y ) − u i ( y , − ˆ x ) N [ u s ( · , d )]( y )  d s y − Z ∂ D  ∂ n y u s ( y , d ) M [ u i ( · , − ˆ x )]( y ) − ∂ n y u i ( y , − ˆ x ) M [ u s ( · , d )]( y )  d s y . where u i ( y , − ˆ x ) = e − ik ˆ x · y . A similar form ula holds for u ∞ ( − d, − ˆ x ) . Subtracting these t wo expressions and using the fact that u = u i + u s leads to the following: u ∞ ( ˆ x, d ) − u ∞ ( − d, − ˆ x ) = Z ∂ D  u ( y , − ˆ x ) N [ u ( · , d )]( y ) − u ( y , d ) N [ u ( · , − ˆ x )]( y )  d s y + Z ∂ D  − ∂ n y u ( y , − ˆ x ) M [ u ( · , d )]( y ) + ∂ n y u ( y , d ) M [ u ( · , − ˆ x )]( y )  d s y where u ( y , d ) and u ( y , − ˆ x ) are the solutions to ( 1.1 ) with incident directions d and − ˆ x ∈ S 1 , resp ectiv ely . Applying the b oundary conditions on ∂ D leads to the desired result. Note that this recipro cit y result also holds for the well-kno wn clamp e d and fr e e plate b oundary conditions. With this result, we are ready to analyze the op erator F giv en in ( 2.6 ) and obtain its factorization. 3.2 F acto rization of the fa r-field op erato r T o analyze the solv abilit y and stability of equation ( 3.1 ), w e examine a factorization of the far-field op erator F . First, w e define the Her glotz wave op er ator H : L 2 ( S 1 ) → H 3 / 2 ( ∂ D ) × H − 1 / 2 ( ∂ D ) : H [ g ] = − v g ν ∆v g + (1 − ν ) ∂ 2 n v g !      ∂ D , (3.3) where v g is the Herglotz wa ve function, defined as v g ( x ) = Z S 1 e ikx · d g ( d ) d s ( d ) . (3.4) W e also define the data-to-p attern op er ator G : H 3 / 2 ( ∂ D ) × H − 1 / 2 ( ∂ D ) → L 2 ( S 1 ) whic h maps supp orted plate b oundary data to the far-field pattern of the scattered field u s : G [ ( f 1 , f 2 ) T ]( ˆ x ) = u ∞ ( ˆ x ) , (3.5) In [ 8 ], solutions to the b oundary v alue problem ( 2.1 )-( 2.3 ) w ere shown to b e unique for u s ∈ H 2 loc ( R 2 \ D ) and for b oundary data ( f 1 , f 2 ) ∈ H 3 / 2 ( ∂ D ) × H − 1 / 2 ( ∂ D ) . Therefore, the data-to-pattern op erator G is well-defined. Then, by the linearit y of the b oundary v alue problem, the far-field op erator F can b e written as F = G H . (3.6) 5 This decomp osition allo ws us to describ e the existence of solutions to ( 3.1 ) in terms of the range of the op erator G . W e now present a lemma that establishes a necessary condition for the existence of a solution to ( 3.1 ): Lemma 3.1. L et G b e the data-to-p attern op er ator define d in ( 3.5 ) . Then ϕ z ∈ Range( G ) if and only if z ∈ D . The pro of of this lemma follo ws as in [ 25 ] and w e omit it here for succinctness. Due to the factorization of F , we hav e Range( F ) ⊂ Range( G ) . It follo ws that ϕ z ∈ Range( F ) only if z ∈ D . W e also ha ve the follo wing lemma for the Herglotz wa v e operator H : Lemma 3.2. The Her glotz wave op er ator H : L 2 ( S 1 ) → H 3 / 2 ( ∂ D ) × H − 1 / 2 ( ∂ D ) define d by ( 3.3 ) is c omp act and inje ctive. Pr o of. The proof follo ws as in [ 25 ] with minor mo difications. First, to sho w compact- ness, observe that the corresp onding Herglotz wa v e functions are smo oth solutions of the Helmholtz equation in R 2 . Hence, we ha ve v g ∈ H p loc ( R 2 ) for an y p > 0 and H [ g ] ∈ H p − 1 / 2 ( ∂ D ) × H p − 5 / 2 ( ∂ D ) for any p > 0 . The compactness of H in to H 3 / 2 ( ∂ D ) × H − 1 / 2 ( ∂ D ) then follo ws directly from standard Sob olev embedding theorems. T o prov e injectivit y , assume that H [ g ] = 0 . This implies that the asso ciated Herglotz w av e function satisfies v g = 0 and M [ v g ] = − ν k 2 v g + (1 − ν ) ∂ 2 n v g = (1 − ν ) ∂ 2 n v g = 0 on ∂ D . W e now mak e use of the standard decomp osition of the Laplacian in R 2 restricted to a smo oth b oundary , ∆ = ∂ 2 n + κ ∂ n + ∆ surf , (3.7) where κ denotes the curv ature of ∂ D and ∆ surf is the surface Laplacian; see [ 40 ]. Since v g v anishes identically on ∂ D , it follows that ∆ surf v g = 0 on ∂ D . Com bining this with the condition ∂ 2 n v g = 0 and the decomp osition ( 3.7 ), w e hav e ∂ n v g = 0 on ∂ D . By the unique con tinuation prop erty for solutions of the Helmholtz equation, it follo ws that v g = 0 in R \ D . Consequently , g = 0 , establishing the injectivity of H . Next, w e present a result concerning a key analytical prop ert y of the far-field op erator. Com bining the previous lemma with the reciprocity principle Theorem 3.1 , we arrive at the follo wing result: Theorem 3.2. A ssume that k > 0 is not a supp orte d plate tr ansmission eigenvalue, define d as k such that ther e exists a nontrivial ( p, q ) ∈ H 1 ( R 2 \ D ) × H 1 ( D ) satisfying              ∆p − k 2 p = 0 , in R 2 \ D , ∆q + k 2 q = 0 , in D , p + q = 0 , on ∂ D , M [ p + q ] = 0 , on ∂ D , (3.8) wher e p de c ays exp onential ly as | x | → ∞ . Then the far-field op er ator F : L 2 ( S 1 ) → L 2 ( S 1 ) is c omp act, inje ctive, and has dense r ange. Pr o of. The compactness of F follows directly from Lemma 3.2 . The pro of of injectivity follo ws as in [ 25 ] but we provide a sk etch of it here for completeness. Supp ose F [ g ] = 0 . Recall that F [ g ]( ˆ x ) is the far-field pattern of the solution u s for an incident field u i = v g . F or the far-field pattern to b e zero, u s m ust satisfy the mo dified Helmholtz equation in 6 R 2 \ D . Mean while, v g satisfies the Helmholtz equation in D . This implies that ( u s , v g ) solv es the supp orted plate transmission eigen v alue problem, which by assumption can only b e true if u s = v g = 0 . Therefore, g = 0 . Finally , to see that F has dense range, we can use the recipro cit y relation to write the adjoint of F as in [ 25 ]: F ∗ [ g ]( ˆ x ) = Z S 1 u ∞ ( − ˆ x, − d ) g ( d )d s ( d ) F ∗ m ust also be injectiv e b y the same argumen t as abov e, and therefore Range ( F ) = L 2 ( S 1 ) . Giv en the result in Theorem 3.2 , it b ecomes clear that due to the compactness of the far- field op erator F , problem ( 3.1 ) is ill-p osed. T o mitigate the ill-p osedness of the problem, w e need to apply some regularization metho d [ 15 ]. In th is pap er, we will use Tikhonov regularization where the problem b eing solved b ecomes the minimization of ∥F [ g z ] − ϕ z ∥ 2 L 2 ( S 1 ) + α ∥ g z ∥ 2 L 2 ( S 1 ) , (3.9) where α > 0 is the regularization parameter. In this case, we will define g α z as the minimizer of ( 3.9 ). The following result describes the behavior of the indicator function in the LSM as the regularization parameter α go es to zero. In essence, it shows that the approximate solution to the far-field equation ( 3.1 ) remains b ounded when the sampling p oin t lies within the scatterer D , which serves as a practical computational criterion for reconstructing the scatterer from the far-field data. Theorem 3.3. A ssume that the wavenumb er k is not an eigenvalue of the supp orte d tr ans- mission eigenvalue pr oblem given in ( 3.8 ) . Then, for any z ∈ R 2 \ D and a family of functions { g α z } ⊂ L 2 ( S 1 ) satisfying lim α → 0 ∥F [ g α z ] − ϕ z ∥ L 2 ( S 1 ) = 0 , (3.10) it fol lows that lim α → 0 ∥ g α z ∥ L 2 ( S 1 ) = ∞ . Pr o of. The pro of follo ws by an argumen t similar to [ 25 ]; we outline the argument here for completeness and readability . Because k is not a supp orted transmission eigen v alue, F is injectiv e, compact, and has dense range. In particular, w e know that for ϕ z ∈ L 2 ( S 1 ) suc h a family { g α z } satisfying ( 3.10 ) exists and can b e giv en as the family of functions that minimizes ( 3.9 ). No w, suppose b y contradiction that ∥ g α z ∥ L 2 ( S 1 ) remains b ounded as α → 0 for some z ∈ R 2 \ D . Then, there exists a subsequence { g α i z } ∞ i =0 suc h that g α i z ⇀ g z w eakly in L 2 ( S 1 ) as i → ∞ . Since F is compact, lim i →∞ F [ g α i z ] = F [ g z ] strongly in L 2 ( S 1 ) . This implies that there exists a g z ∈ L 2 ( S 1 ) such that F [ g z ] = ϕ z for z ∈ R 2 \ D . Ho wev er, this is a con tradiction of Lemma 3.1 , therefore g α z is unbounded. Finally , the LSM indicator function is defined as I LSM ( z ) = 1 ∥ g α z ∥ L 2 ( S 1 ) . (3.11) 7 By Theorem 3.3 , ∥ g α z ∥ L 2 ( S 1 ) blo ws up as α → 0 whenev er z / ∈ D , and therefore the indicator function tends to zero for sampling p oin ts z outside the ca vity . F or z ∈ D , approximate solutions with b ounded norm exist, but no uniform low er b ound for I LSM ( z ) can b e guar- an teed. Th us, the LSM characterizes the cavit y D through the growth b eha vior of the regularized solutions g α z . 4 Direct sampling metho d While the linear sampling metho d pro vides one approach for the qualitative reconstruction of obstacles, it requires solving a large n umber of ill-posed equations. T o address this limitation, w e in tro duce the direct sampling metho d, whic h ac hieves similar reconstruction capabilities with a reduced computational b ottlenec k. The goal of the direct sampling metho d is to apply the far-field op erator to the far-field pattern of a p oin t source to compute the following tw o indicator functions: I DSM - 1 ( z ) : =    ⟨ ϕ z , F [ ϕ z ] ⟩ L 2 ( S 1 )    ρ/ 2 , and I DSM - 2 ( z ) : = ∥F [ ϕ z ] ∥ ρ L 2 ( S 1 ) , (4.1) where ϕ z ( d ) = e − ikz · d is the far-field pattern of the fundamen tal solution (up to some rescaling) cen tered at a p oin t z and ρ > 0 is a regularization parameter used to refine the resolution of the imaging functionals. While direct sampling approaches hav e been widely explored for acoustic scattering problems (see, for instance, [ 29 ]), their applicability to biharmonic equations has only re- cen tly b een explored for the clamp ed plate problem [ 24 ], and our analysis follows closely . The primary theoretical results concern the factorization of the forw ard op erator G , which allo ws us to predict the decay of the indicator functions ab o ve. 4.1 Integral formulation fo r the supp o rted plate A dopting the form ulation for the supp orted plate problem app earing in [ 1 , 19 , 38 ], w e write the scattered field u s in terms of b oundary densities φ 1 and φ 2 : u s ( x ) = Z ∂ D  ∂ 3 n y Φ ( x, y ) + α 1 ∂ n y ∂ 2 τ y Φ ( x, y ) + α 2 κ ( y ) ∂ 2 n y Φ ( x, y ) + α 3 κ ′ ( y ) ∂ τ y Φ ( x, y )  φ 1 ( y ) d s y + Z ∂ D ∂ n y Φ ( x, y ) φ 2 ( y ) d s y , (4.2) where κ ( y ) is the curv ature at the p oin t y ∈ ∂ D , κ ′ ( y ) is the arc-length deriv ative of the curv ature, and ∂ n y and ∂ τ y represen t the normal and tangen tial deriv ativ es at y , resp ectiv ely . The co efficien ts α 1 , α 2 , and α 3 are given by: α 1 = 2 − ν , α 2 = ( − 1 + ν )(7 + ν ) 3 − ν , α 3 = (1 − ν )(3 + ν ) 1 + ν . Applying the supp orted plate b oundary conditions to u = u s + u i leads to the following system of b oundary in tegral equations deriv ed in [ 1 , 38 ]: − 1 2 I + K 11 K 12 c 0 κ 2 I + K 21 − 1 2 I + K 22 ! φ 1 φ 2 ! = − u i M [ u i ] ! , (4.3) 8 where c 0 = ( ν − 1)( ν +3)(2 ν − 1) 2(3 − ν ) and K 11 , K 12 , K 21 , and K 22 are boundary integral op erators with the following kernels: K 11 ( x, y ) = h ∂ 3 n y + α 1 ∂ n y ∂ 2 τ y + α 2 κ ( y ) ∂ 2 n y + α 3 κ ′ ( y ) ∂ τ y i Φ ( x, y ) , K 12 ( x, y ) = ∂ n y Φ ( x, y ) , K 21 ( x, y ) = h ∂ 2 n x ∂ 3 n y + α 1 ∂ 2 n x ∂ n y ∂ 2 τ y + α 2 κ ( y ) ∂ 2 n x ∂ 2 n y + α 3 κ ′ ( y ) ∂ 2 n x ∂ τ y + ν ∂ 2 τ x ∂ 3 n y + ν α 1 ∂ 2 τ x ∂ n y ∂ 2 τ y + ν α 2 κ ( y ) ∂ 2 τ x ∂ 2 n y + ν α 3 κ ′ ( y ) ∂ 2 τ x ∂ τ y i Φ ( x, y ) , K 22 ( x, y ) = h ∂ 2 n x ∂ n y + ν ∂ 2 τ x ∂ n y i Φ ( x, y ) , Next, we hav e the following statemen t ab out the integral equation ab o v e: Theorem 4.1. Equation ( 4.3 ) is F r e dholm se c ond kind on H 3 / 2 ( ∂ D ) × H − 1 / 2 ( ∂ D ) . Before proving the ab ov e theorem, we require the follo wing lemma. Lemma 4.1. The op er ator K 11 is c omp act fr om H 3 / 2 ( ∂ D ) → H 3 / 2 ( ∂ D ) , K 12 is c omp act fr om H − 1 / 2 ( ∂ D ) → H 3 / 2 ( ∂ D ) , K 21 is c omp act fr om H 3 / 2 ( ∂ D ) → H − 1 / 2 ( ∂ D ) , and K 22 is c omp act fr om H − 1 / 2 ( ∂ D ) → H − 1 / 2 ( ∂ D ) . Pr o of. Recall that the p o wer series for the fundamen tal solution Φ ( x, y ) cen tered at x = y is given by: Φ ( x, y ) = K S ( x, y ) + Φ B ( x, y ) + O ( | x − y | 6 ln | x − y | ) , (4.4) where K S ( x, y ) is a C ∞ function and the leading order singularity of Φ ( x, y ) is giv en b y the biharmonic Green’s function Φ B ( x, y ) : = 1 8 π | x − y | 2 ln | x − y | . Let K B ij denote the kernels ab ov e with deriv atives applied to Φ B ( x, y ) . In [ 38 ], it was sho wn that K B 21 and K B 22 are smooth for a smo oth curve and, therefore, K B 21 and K B 22 are compact. It w as also sho wn that h ∂ 3 n y + α 1 ∂ n y ∂ 2 τ y i Φ B ( x, y ) is smo oth when restricted to the curv e, therefore K B 11 − K B 11 is compact, where K B 11 is the remainder corresp onding to the terms h α 2 κ ( y ) ∂ 2 n y + α 3 κ ′ ( y ) ∂ τ y i Φ B ( x, y ) . T o sho w compactness of K B 11 , we note that since Φ B ( x, y ) is the fundamental solution of a fourth-order elliptic PDE, the theory of pseudo-differen tial op erators (e.g. Theorem 9.5.8 in [ 26 ]) guarantees that K B 11 is a contin uous op erator from H 3 / 2 ( ∂ D ) → H 5 / 2 ( ∂ D ) which is compactly em b edded in H 3 / 2 ( ∂ D ) . T o show that K B 12 is compact, we note that this k ernel can b e expanded on the surface as K B 12 ( γ ( s ) , γ (0)) = 1 8 π s 2 κ log( s ) + O ( s 3 log( s )) , where γ ( s ) is the arc-length parametrization of the curve. F or a smooth curve, this k ernel has the same leading-order singularit y as the k ernel of the biharmonic single lay er p oten tial R ∂ D G B ( x, y ) φ ( y ) d s y , which b y Theorem 9.5.8 in [ 26 ] is a con tinuous op erator from H − 1 / 2 ( ∂ D ) → H 5 / 2 ( ∂ D ) . Therefore, K B 12 is compact from H − 1 / 2 ( ∂ D ) → H 3 / 2 ( ∂ D ) . The next singular term in the p o w er series for Φ ( x, y ) is | x − y | 6 log( | x − y | ) , which is the fundamental solution to an eighth-order elliptic PDE. Since w e are taking at most five deriv ativ es in the kernels ab o ve, the rest of the terms in the p o wer series corresp ond to compact op erators. Theorem 4.1 follows immediately from the lemma ab o ve. F or conv enience, let Q denote the op erator app earing on the left-hand side of ( 4.3 ): Q : = − 1 2 I + K 11 K 12 c 0 κ 2 I + K 21 − 1 2 I + K 22 ! . 9 By analytic F redholm theory , the op erator Q is in vertible except at a discrete set of w av en umbers. Applying the b ounded inv erse theorem, there exists a Q − 1 : H 3 / 2 ( ∂ D ) × H − 1 / 2 ( ∂ D ) → H 3 / 2 ( ∂ D ) × H − 1 / 2 ( ∂ D ) for almost all w av en umbers such that: φ 1 φ 2 ! = −Q − 1 u i M [ u i ] ! . Then, from ( 2.4 ) and ( 4.2 ), the far-field pattern of the scattered field is giv en by: u ∞ ( ˆ x ) = Z ∂ D  ∂ 3 n y e − ik ˆ x · y + α 1 ∂ n y ∂ 2 τ y e − ik ˆ x · y + α 2 κ ( y ) ∂ 2 n y e − ik ˆ x · y + α 3 κ ′ ( y ) ∂ τ y e − ik ˆ x · y  φ 1 ( y ) d s y + Z ∂ D ∂ n y e − ik ˆ x · y φ 2 ( y ) d s y , (4.5) Therefore, we define the following op erator M ∞ : H 3 / 2 ( ∂ D ) × H − 1 / 2 ( ∂ D ) → L 2 ( S 1 ) which maps from b oundary densities to the far-field pattern: M ∞ [( φ 1 , φ 2 ) T ]( ˆ x ) : = Z ∂ D  ∂ 3 n y e − ik ˆ x · y + α 1 ∂ n y ∂ 2 τ y e − ik ˆ x · y + α 2 κ ( y ) ∂ 2 n y e − ik ˆ x · y + α 3 κ ′ ( y ) ∂ τ y e − ik ˆ x · y  φ 1 ( y ) d s y + Z ∂ D ∂ n y e − ik ˆ x · y φ 2 ( y ) d s y , Consequen tly , the far-field pattern for a supp orted ca vity D for an incident plane wa v e u i can b e represen ted as u ∞ ( ˆ x, d ) = M ∞ Q − 1 u i M [ u i ] ! , (4.6) allo wing us to write the far-field op erator as F = M ∞ Q − 1 H . (4.7) This factorization will b e useful for proving prop erties of the indicator functions used for the DSM in the follo wing section. 4.2 Prop erties of the imaging functions The next result pro vides an estimate on the rate of decay for the indicator function I DSM - 1 a wa y from the obstacle D . Theorem 4.2. L et I DSM - 1 b e define d as ab ove with r e gularization p ar ameter ρ > 0 . Then, for almost al l wavenumb ers, I DSM - 1 satisfies:    I DSM - 1 ( x ) | = O (dist( x, D ) − ρ/ 2 ) , as dist( x, D ) → ∞ , for some c omp act domain D , wher e dist( x, D ) : = inf {| x − y | | y ∈ D } . Pr o of. Given the factorization F = M ∞ Q − 1 H , the op erator Q − 1 is b ounded for almost all w av en umbers. Therefore, we can use Cauc hy-Sc hw arz to write:    ⟨ ϕ z , F [ ϕ z ] ⟩ L 2 ( S 1 )    ρ/ 2 =    ⟨ ( M ∞ ) ∗ [ ϕ z ] , Q − 1 H [ ϕ z ] ⟩    ρ/ 2 ≤ C ∥ ( M ∞ ) ∗ [ ϕ z ] ∥ ρ/ 2 H − 3 / 2 ( ∂ D ) × H 1 / 2 ( ∂ D ) ∥H [ ϕ z ] ∥ ρ/ 2 H 3 / 2 ( ∂ D ) × H − 1 / 2 ( ∂ D ) (4.8) 10 where ( M ∞ ) ∗ is given by ( M ∞ ) ∗ [ ϕ z ]( x ) = " ∂ 3 n x + α 1 ∂ n x ∂ 2 τ x + α 2 κ ( x ) ∂ 2 n x + α 3 κ ′ ( x ) ∂ τ x ∂ n x # Z S 1 e ikx · d ϕ z ( d ) d s ( d ) = 2 π " ∂ 3 n x + α 1 ∂ n x ∂ 2 τ x + α 2 κ ( x ) ∂ 2 n x + α 3 κ ′ ( x ) ∂ τ x ∂ n x # J 0 ( k | x − z | ) for x ∈ ∂ D where J 0 denotes the Bessel function of the first kind of order zero. Then, w e can use the trace theorem to write: ∥ ( M ∞ ) ∗ [ ϕ z ] ∥ H − 3 / 2 ( ∂ D ) × H 1 / 2 ( ∂ D ) ≤ C ∥ J 0 ( k | · − z | ) ∥ H 2 ( D ) The other term in inequalit y ( 4.8 ) can b e written in terms of J 0 as well, i.e. H [ ϕ z ]( x ) = 2 π  J 0 ( k | x − z | ) , M [ J 0 ( k | · − z | )]( x )  T , Again, we use the trace theorem to say that ∥H [ ϕ z ] ∥ H 3 / 2 ( ∂ D ) × H − 1 / 2 ( ∂ D ) ≤ C ∥ J 0 ( k | · − z | ) ∥ H 2 ( D ) Finally , in Lemma 3.2 of [ 24 ], the authors sho wed that for | α | ≤ 4 ∥ D α J 0 ( | x − · | ) ∥ L 2 ( D ) = O (dist( x, D ) − 1 ) as dist( x, D ) → ∞ , from whic h the result follows immediately . No w, we turn to the second imaging functional of interest, defined b y I DSM - 2 ( x ) = ∥F [ ϕ z ] ∥ ρ L 2 ( S 1 ) , ρ > 0 . Though the b eha vior of this second functional is similar to that of the first, w e include it here for completeness. It is p ossible to sho w that this functional has the same decay prop erties as the first functional, and indeed, this was recently shown for the clamp ed plate problem. W e state this result and refer the reader to [ 24 ] for the pro of: Theorem 4.3. L et I DSM - 2 b e define d as ab ove with r e gularization p ar ameter ρ > 0 . Then, for almost al l wavenumb ers, I DSM - 2 satisfies:    I DSM - 2 ( x ) | = O (dist( x, D ) − ρ/ 2 ) , as dist( x, D ) → ∞ , for some c omp act domain D , wher e dist( x, D ) : = inf {| x − y | | y ∈ D } . Using Theorem 4.3 , it is evident that the t wo imaging functionals under consideration are equiv alent up to multiplicativ e constan ts. Sp ecifically , there exist p ositiv e constants c 1 , c 2 > 0 such that c 1 |I DSM - 1 ( x ) | ≤ |I DSM - 2 ( x ) | ≤ c 2 |I DSM - 1 ( x ) | . for all x ∈ R 2 \ D . As a result, b oth imaging functionals are exp ected to provide comparable reconstructions of the obstacle D , and indeed, this is what is observ ed n umerically . 11 5 Numerical exp eriments In this section, w e first discuss the n umerical implementation of the metho ds b efore presen t- ing n umerical exp erimen ts for the LSM and DSM. In the first exp erimen t, we use simulated noise-free data at different frequencies to assess the qualitative p erformance of the tw o metho ds. In the second exp erimen t, we rep eat the same setup but introduce additive and m ultiplicative noise of v arying intensities. In the third exp erimen t, w e inv estigate the b e- ha vior of the metho ds for different v alues of the Poisson’s ratio. In the fourth exp erimen t, w e consider the reconstruction of m ultiple obstacles. Finally , w e study the effect of limited a v ailable data on the reconstruction of single and multiple obstacles. 5.1 Numerical implementation W e no w discuss the implemen tation of the metho ds. F or b oth approaches, we must discretize the far-field op erator. Since d = (cos( θ d ) , sin( θ d )) ∈ S 1 , we discretize the interv al [0 , 2 π ) uniformly to obtain d j = (cos(2 j π / N d ) , sin(2 j π / N d )) , j = 1 , . . . , N d , where N d denotes the n umber of inciden t directions used to prob e the domain. Using these p oints, we approximate the integral in the far-field op erator via the trap ezoidal rule: F [ g z ]( ˆ x ) ≈ 2 π N d N d X j =1 u ∞ ( ˆ x, d j ) g z ( d j ) . (5.1) F or the LSM, w e solve the regularized minimization problem ( 3.9 ), which is equiv alent to solving ( α I + F ∗ F )[ g z ]( ˆ x ) = F ∗ [ ϕ z ]( ˆ x ) , (5.2) where I is the identit y op erator. T aking ˆ x = (cos( θ r ) , sin( θ r )) ∈ S 1 , we can discretize the in terv al [0 , 2 π ) uniformly , to obtain ˆ x ℓ = (cos(2 ℓπ / N r ) , sin(2 ℓπ / N r )) , ℓ = 1 , . . . , N r , where N r is the num b er of receivers used to measure the scattered data. With these discretizations, we obtain the linear system ( α I + F ∗ F ) g z = F ∗ ϕ z , (5.3) where F is an N r × N d matrix with elements F ℓj = (2 π / N d ) u ∞ ( ˆ x ℓ , d j ) , F ∗ is its adjoint, I is the N d × N d iden tity matrix, g z is a v ector with N d comp onen ts g z ( d j ) b eing our v ariables, and ϕ z is a vector with N r comp onen ts ϕ z ( ˆ x ℓ ) . After computing g z at all sampling p oin ts z , we ev aluate the LSM indicator function ( 3.11 ). Regarding regularization, in the LSM we tested α = 10 − j , j = 1 , . . . , 6 , and in each exp erimen t w e present the reconstruction giving the b est qualitativ e result; reconstructions for other v alues of α are very similar. F or the DSM, since the scaling parameter ρ > 0 affects only visualization and do es not influence the lo calization prop erties of the indicators, we set ρ = 2 for I DSM - 1 and ρ = 1 for I DSM - 2 ; other c hoices lead to visually comparable results. Finally , in all figures, the graphs are rescaled so that the maximum v alue of the indicator function is one. W e display the rescaled DSM indicators directly , whereas for the LSM we plot the rescaled v alue of log ( I LSM ) to enhance contrast and impro ve visualization. T o obtain the far-field measuremen ts needed in the discretization of the far-field op er- ator F , we use the b oundary integral form ulation describ ed in Section 4.1 and originally app earing in [ 19 , 38 ]. The forward solv er is implemented using the MA TLAB soft ware pac k- age chunkIE [ 5 ]. In this package, the b oundary of the domain is divided into panels, and 12 the integral op erators are discretized using generalized Gaussian quadratures. The in tegral equation ( 4.3 ) is then solved using a Nyström metho d via GMRES. The calculation of the far-field pattern is done using the native smo oth quadrature. Because equation ( 5.3 ) must b e solv ed for ev ery sampling point z , it is efficien t to precompute a factorization of the left-hand side once (e.g., using an SVD of F ) and reuse it for eac h righ t-hand side. Computing such a factorization costs O  min( N d N 2 r , N r N 2 d )  , and applying it to each righ t-hand side costs O ( N 2 d ) . Th us, for a total of N g sampling p oin ts, the ov erall complexity of LSM is O  min( N d N 2 r , N r N 2 d )  + O ( N g N 2 d ) . If N d = N r , this simplifies to O ( N 3 d ) + O ( N g N 2 d ) . In contrast, the DSM implementation is muc h simpler. No ill-p osed equation must b e solved, and therefore no regularization is needed. The indicator functions ( 4.1 ) are simply ev aluated at each sampling p oin t. The dominant computational cost is computing the pro duct F ϕ z , where ( ϕ z ) j = ϕ z ( d j ) . In general, the cost of the DSM is O ( N g N r N d ) . When F is square with N d = N r , the ov erall DSM cost is O ( N g N 2 d ) . W e use three obstacle shapes in our exp erimen ts: a 5-arms star-shap ed domain param- eterized by γ A ( t ) = (1 + 0 . 3 cos(5 t ))(cos( t ) , sin( t )) , t ∈ [0 , 2 π ) , an 11-arms star-shap ed domain parameterized by γ B ( t ) = (1 + 0 . 5 cos(11 t ))(cos( t ) , sin( t )) , t ∈ [0 , 2 π ) , and a domain with a pronounced cavit y . W e will refer to those domains as 5-arms, 11-arms, and cavit y throughout the rest of this text. -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 (a) 5-arms -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 (b) 11-arms -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 (c) cavit y Figure 1: Ground-truth domains used in the simulations: a) 5-arms, b) 11-arms, and c) cavit y . In all exp erimen ts – except for Example 3, where we analyze the effect of the P oisson’s ratio – we use ν = 0 . 3 , which is t ypical for materials such as ice and steel. In Example 3, w e consider ν = 0 . 5 (close to rubb er), 0 (cork), and − 0 . 5 (auxetic metamaterials). W e test the n umerical accuracy of the forward solver for each of the w av enum b ers and eac h of the domains and rep ort it in T able 1 . The table lists the num b er of p oin ts N used on the b oundary for the discretization and the error in the calculation of the solution. The error was assessed using an analytic solution test b y placing a p oin t source inside the obstacle and solving for the field at a collection of 10 arbitrary p oints outside. Then, the ℓ ∞ norm of the error at those p oin ts w as calculated relativ e to the L 1 norm of the solution to the in tegral equation. Note that the forw ard solv er pro vides at least 7 digits of accuracy for all the discretizations used. 13 T able 1: Discretization details and accuracy of the computed field. Domain k N Relativ e error 5-arms 2 π 288 1 . 26 × 10 − 7 20 π 2976 2 . 46 × 10 − 12 11-arms 2 π 3776 2 . 42 × 10 − 7 20 π 9344 2 . 43 × 10 − 10 Ca vity 2 π 768 1 . 63 × 10 − 8 20 π 3584 3 . 84 × 10 − 14 Multiple 2 π 1248 2 . 22 × 10 − 9 10 π 5376 3 . 00 × 10 − 14 5.2 Example 1 – Resolution The main ob jectiv e of this example is to illustrate the p erformance of the method at differen t frequencies. In particular, we chose k = 2 π and 20 π so that the wa velength of the incoming w av e is λ = 1 and 0 . 1 , resp ectiv ely . The indicator function is calculated on a uniform grid o ver the domain [ − 3 , 3] with 300 p oints in eac h direction, totaling N g = 300 2 . W e use N d = N r = 128 for k = 2 π and N d = N r = 1024 for k = 20 π . F or the LSM, we set the regularization parameter to α = 10 − 4 based on empirical evidence that this choice of regularization parameter b est matc hes the original shap e. The reconstructions for k = 2 π are shown in Figure 2 , and those for k = 20 π in Figure 3 . In all cases, the general lo cation and what seems to be a substantial portion of the conv ex h ull of the obstacles are reco vered. F or the 5-arms domain, the reconstructions are very close to the ground truth; ho wev er, for the 11-arms domain, the metho ds fail to resolve the small ca vities along the b oundary . A similar effect o ccurs in the reconstruction of the ca vity-shaped domain, where most of the cavit y structure is lost. A t higher frequency , b oth metho ds yield more detailed reconstructions. Overall, b oth metho ds are effective for obtaining a coarse appro ximation of the domain, but they are unable to recov er fine-scale features, particularly ca vities. When extracting a lev el-set curve to represent the reconstructed obstacle, one should therefore expect a go od appro ximation of its conv ex h ull rather than its detailed b oundary structure. The b eha vior observed in this example mirrors observ ations in related Helmholtz scatter- ing problems; for example, [ 2 ] sho ws that cavities are difficult to recov er even with iterative metho ds in the Helmholtz setting. In [ 7 ], some cavities w ere reconstructed for p enetra- ble obstacles, and w e exp ect similar results could b e ac hieved in the biharmonic case for p enetrable inclusions as w ell, though this remains unclear for impenetrable ones. 14 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (a) LSM, k = 2 π , 5-arms -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (b) DSM–1, k = 2 π , 5-arms -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (c) DSM–2, k = 2 π , 5-arms -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (d) LSM, k = 2 π , 11-arms -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (e) DSM–1, k = 2 π , 11-arms -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (f ) DSM–2, k = 2 π , 11-arms -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (g) LSM, k = 2 π , cavit y -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (h) DSM–1, k = 2 π , cavit y -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (i) DSM–2, k = 2 π , cavit y Figure 2: Reconstructions using the LSM (left column), DSM–1 (middle column), and DSM–2 (righ t column) at wa ven um b er k = 2 π : the 5-arms domain (top ro w), the 11-arms domain (middle row), and the ca vity (b ottom ro w). 15 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.5 0.6 0.7 0.8 0.9 1 (a) LSM, k = 20 π , 5-arms -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (b) DSM–1, k = 20 π , 5-arms -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (c) DSM–2, k = 20 π , 5-arms -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.55 0.6 0.65 0.7 0.75 0.8 0.85 0.9 0.95 1 (d) LSM, k = 20 π , 11-arms -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (e) DSM–1, k = 20 π , 11-arms -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (f ) DSM–2, k = 20 π , 11-arms -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.5 0.6 0.7 0.8 0.9 1 (g) LSM, k = 20 π , cavit y -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (h) DSM–1, k = 20 π , cavit y -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (i) DSM–2, k = 20 π , cavit y Figure 3: Reconstructions using the LSM (left column), DSM–1 (middle column), and DSM–2 (righ t column) at w av enum b er k = 20 π for the 5-arms domain (top ro w), the 11-arms domain (middle row), and the cavit y (b ottom ro w). 16 5.3 Example 2 – Noisy data Next, we examine the b eha vior of the metho ds when the av ailable data con tain tw o types of noise: additive noise and m ultiplicative noise. W e use the 5-arms domain for this ex- p erimen t, and the data are generated at the wa ven umber k = 2 π . The n umber of incident directions and receivers is giv en by N d = N r = 128 . T o add noise, w e generate a complex random n umber δ ∈ C whose real and imaginary parts are independently sampled from a normal distribution with zero mean and unit v ari- ance. The measurements used in this exp eriment ˜ u ∞ are obtained by adding noise to the noiseless data u ∞ . F or the additiv e noise, the noise is added using the formula: ˜ u ∞ = u ∞ + c i δ | δ | | u ∞ | , and for the multiplicativ e noise, the noise is added using the formula: ˜ u ∞ = u ∞ + c i δ | δ | , where c i is the noise lev el. W e c ho ose c i = 0 . 05 , 0 . 5 , and 1 , corresp onding to 5% , 50% , and 100% , resp ectiv ely . F or the LSM reconstruction, the regularization parameter was chosen empirically as α = 10 − 1 for all noise lev els and noise types. The indicator function for all metho ds is computed on the same grid used in the previous example. The results are presen ted in Figures 4 and 5 , corresp onding to the additive and mul- tiplicativ e noisy data, resp ectiv ely . Note that the DSM is more robust to noise than the LSM, as it pro duces reconstructions that remain muc h closer to the original domain even as the noise level increases. In contrast, the LSM reconstruction deteriorates under b oth noise types, with an increased degradation in the presence of additiv e noise. 17 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.6 0.65 0.7 0.75 0.8 0.85 0.9 0.95 1 (a) LSM, additive, δ = 0 . 05 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (b) DSM–1, additive, δ = 0 . 05 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (c) DSM–2, additive, δ = 0 . 05 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.75 0.8 0.85 0.9 0.95 1 (d) LSM, additive, δ = 0 . 5 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (e) DSM–1, additive, δ = 0 . 5 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (f ) DSM–2, additive, δ = 0 . 5 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.8 0.85 0.9 0.95 1 (g) LSM, additive, δ = 1 . 0 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (h) DSM–1, additive, δ = 1 . 0 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (i) DSM–2, additive, δ = 1 . 0 Figure 4: Reconstructions using the LSM, DSM–1, and DSM–2 for the 5-arms domain at wa ven umber k = 2 π using data with 5% (top ro w), 50%(middle row), and 100%(bottom row) additiv e noise. 18 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.55 0.6 0.65 0.7 0.75 0.8 0.85 0.9 0.95 1 (a) LSM, mult., δ = 0 . 05 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (b) DSM–1, mult., δ = 0 . 05 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (c) DSM–2, mult., δ = 0 . 05 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.7 0.75 0.8 0.85 0.9 0.95 1 (d) LSM, mult., δ = 0 . 5 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (e) DSM–1, mult., δ = 0 . 5 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (f ) DSM–2, mult., δ = 0 . 5 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.75 0.8 0.85 0.9 0.95 1 (g) LSM, mult., δ = 1 . 0 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (h) DSM–1, mult., δ = 1 . 0 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (i) DSM–2, mult., δ = 1 . 0 Figure 5: Reconstructions using the LSM, DSM–1, and DSM–2 for the 5-arms domain at wa ven umber k = 2 π using data with 5% (top ro w), 50%(middle row), and 100%(bottom row) m ultiplicative noise. 19 5.4 Example 3 – Va riable P oisson’s Ratio In this exp erimen t, w e consider the effect of the P oisson’s ratio by testing the sampling metho ds on problems with differen t v alues of ν . W e examine the b eha vior of the methods for ν = − 0 . 5 , 0 , and 0 . 5 , which are approximately the P oisson’s ratio for auxetic materials, cork, and rubb er, resp ectiv ely . The other examples in this pap er use ν = 0 . 3 whic h corresp onds to ice and steel. W e use the 5-arms domain in this exp erimen t. W e consider t wo frequencies for the measuremen ts, k = 2 π and k = 20 π . F or k = 2 π , w e set N d = N r = 128 ; for k = 20 π , we use N d = N r = 1024 . T o a void inv erse crimes, w e insert 5% additiv e noise to the measured data. The indicator functions of the in verse metho ds are calculated on a uniform grid of 300 p oin ts in each direction o ver the domain [ − 3 , 3] . The regularization parameter for the LSM was chosen to b e α = 10 − 1 . The results are presented in Figures 6 for k = 2 π and 7 for 20 π . The DSM indicator functions do not app ear to v ary significantly with the Poisson’s ratio. In con trast, the LSM indicator function shows a sligh t dep endence on ν , though the results are similar. 20 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.6 0.65 0.7 0.75 0.8 0.85 0.9 0.95 1 (a) LSM, k = 2 π , ν = − 0 . 5 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (b) DSM–1, k = 2 π , ν = − 0 . 5 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (c) DSM–2, k = 2 π , ν = − 0 . 5 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.6 0.65 0.7 0.75 0.8 0.85 0.9 0.95 1 (d) LSM, k = 2 π , ν = 0 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (e) DSM–1, k = 2 π , ν = 0 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (f ) DSM–2, k = 2 π , ν = 0 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.6 0.65 0.7 0.75 0.8 0.85 0.9 0.95 1 (g) LSM, k = 2 π , ν = 0 . 5 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (h) DSM–1, k = 2 π , ν = 0 . 5 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (i) DSM–2, k = 2 π , ν = 0 . 5 Figure 6: Reconstructions using the LSM, DSM–1, and DSM–2 for the 5-arms domain at wa ven umber k = 2 π with Poisson’s ratio ν = − 0 . 5 (top ro w), 0 (middle row), and 0 . 5 (b ottom ro w). 21 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.75 0.8 0.85 0.9 0.95 1 (a) LSM, k = 20 π , ν = − 0 . 5 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (b) DSM–1, k = 20 π , ν = − 0 . 5 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (c) DSM–2, k = 20 π , ν = − 0 . 5 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.75 0.8 0.85 0.9 0.95 1 (d) LSM, k = 20 π , ν = 0 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (e) DSM–1, k = 20 π , ν = 0 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (f ) DSM–2, k = 20 π , ν = 0 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.75 0.8 0.85 0.9 0.95 1 (g) LSM, k = 20 π , ν = 0 . 5 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (h) DSM–1, k = 20 π , ν = 0 . 5 -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (i) DSM–2, k = 20 π , ν = 0 . 5 Figure 7: Reconstructions using the LSM, DSM–1, and DSM–2 for the 5-arms domain at w av enum b er k = 20 π with P oisson’s ratio ν = − 0 . 5 (top row), 0 (middle row), and 0 . 5 (b ottom row). 22 5.5 Example 4 – Multiple obstacles One of the ma jor adv antages of sampling metho ds is their ability to handle problems in- v olving multiple obs tacles without requiring any prior information ab out the underlying configuration. In this example, the domain consists of three 5-arms obstacles centered at (2 , 2 . 5) , (2 , − 2 . 5) , and ( − 2 , 0) . The incident wa ves ha ve wa ven umbers k = 2 π and k = 10 π . F or k = 2 π , we use N d = N r = 128 , and for k = 10 π , we use N d = N r = 512 . W e add 5% additive noise to the data to a void in verse crimes. The indicator functions for the reconstruction metho ds are calculated on a uniform grid ov er the domain [ − 10 , 10] with 500 uniformly spaced p oin ts along eac h axis. As in the other examples, for the LSM, we use the regularization parameter α = 10 − 1 . The results are shown in Figure 8 . All metho ds pro duce go od approximations of the obstacles, with the DSM yielding slightly more sharply defined reconstructions than the LSM. Note that recov ering fine details b etw een the obstacles is esp ecially challenging. -10 -5 0 5 10 -10 -5 0 5 10 0.6 0.65 0.7 0.75 0.8 0.85 0.9 0.95 1 (a) LSM, k = 2 π -10 -5 0 5 10 -10 -5 0 5 10 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (b) DSM–1, k = 2 π -10 -5 0 5 10 -10 -5 0 5 10 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (c) DSM–2, k = 2 π -10 -5 0 5 10 -10 -5 0 5 10 0.75 0.8 0.85 0.9 0.95 1 (d) LSM, k = 10 π -10 -5 0 5 10 -10 -5 0 5 10 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (e) DSM–1, k = 10 π -10 -5 0 5 10 -10 -5 0 5 10 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (f ) DSM–2, k = 10 π Figure 8: Reconstructions using the LSM, DSM–1, and DSM–2 for the multiple-obstacle example at w av enum b ers k = 2 π (top row) and k = 10 π (bottom row). 5.6 Example 5 – Limited Data Finally , w e rep ort ho w the amoun t of av ailable data affects the reconstructions. In our exp erimen ts, w e observ ed that when the data are insufficient–particularly for the LSM–the qualit y of the reconstructions deteriorates. W e consider t wo scenarios: a single 5-arms domain, and a multiple-obstacle configuration consisting of three 5-arms obstacles. F or the single obstacle case, w e take k = 20 π , while for the m ultiple-obstacle case w e use k = 10 π . In b oth settings, we use N d = N r = 128 . The indicator functions for the metho ds are calculated on the domain [ − 3 , 3] 2 with 300 p oin ts p er axis for the single-obstacle case, and on the domain [ − 10 , 10] 2 with 500 p oints p er axis for the m ultiple-obstacle case. The 23 regularization parameter is chosen empirically as 10 − 1 . A dditiv e noise of 5% is added to the measured data to a void inv erse crimes. The results are sho wn in Figure 9 . The upp er row presents the reconstructions for the single-obstacle exp erimen t, while the b ottom row displa ys the results for m ultiple obstacles. Note that the LSM reconstruction deteriorates more significan tly than those obtained with the DSM when the amoun t of data is limited. In particular, for the multiple-obstacle case, the LSM reconstruction b ecomes unclear to the p oin t where even the num b er of obstacles is difficult to identify . -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.84 0.86 0.88 0.9 0.92 0.94 0.96 0.98 1 (a) LSM, single obstacle -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (b) DSM–1, single obstacle -3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (c) DSM–2, single obstacle -10 -5 0 5 10 -10 -5 0 5 10 0.7 0.75 0.8 0.85 0.9 0.95 1 (d) LSM, multiple obstacle -10 -5 0 5 10 -10 -5 0 5 10 0.4 0.5 0.6 0.7 0.8 0.9 1 (e) DSM–1, multiple obstacle -10 -5 0 5 10 -10 -5 0 5 10 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 (f ) DSM–2, multiple obstacle Figure 9: Reconstructions using the LSM, DSM–1, and DSM–2 with a limited amount of data ( N d = N r = 128 ). T op: incident w av e k = 20 π and a single 5-arms domain. Bottom: incident wa ve k = 10 π and the three-obstacle configuration. 6 Conclusions W e hav e developed a theoretical foundation for solving the biharmonic inv erse scattering with supp orted plate b oundary conditions using the linear sampling metho d (LSM) and the direct sampling metho d (DSM). T o this end, we hav e established and c haracterized indicator functions that allow the recov ery of imp enetrable supp orted cavities. A comprehensiv e set of numerical exp erimen ts w ere used to examine the influence of frequency , noise, Poisson’s ratio, amount of data, and n umber of obstacles. The results sho w that b oth LSM and DSM reliably recov er the obstacle’s p osition and its con vex hull, with DSM offering sligh tly b etter robustness–esp ecially under noise or limited data av ailability . Consistent with related findings in Helmholtz, neither metho d resolves small cavities, reflecting intrinsic limitations of qualitativ e far-field approaches. LSM ad- ditionally requires careful regularization and is more sensitive to noise and data sparsity , though its p erformance remains comparable in well-resolv ed scenarios. F or multiple ob- 24 stacles, b oth metho ds correctly identify the configuration, with DSM again main taining stabilit y under reduced data. These observ ations highligh t the strengths and inherent limitations of sampling-based tec hniques for biharmonic scattering. F rom a theoretical standp oin t, we hop e to extend the existing in tegral equation form ulation for the supported plate to b e able to handle spurious resonances and establish the inv ertibility of the forward op erator for real frequencies. W e also hope to extend the theory of sampling metho ds to the free plate b oundary conditions. Lastly , we are in terested in using these methods to develop b oundary extraction heuristics and explore direct reconstruction methods based on neural netw orks as in [ 46 ]. A ckno wledgments: The authors thank Jeremy Hoskins, T ristan Go o dwill, Adriaan de Clercq, and Isaac Harris for generous and insightful discussions. Data a v ailability statemen t: Data and co de supp orting the findings of this study are a v ailable up on reasonable request from the authors. References [1] F. Ago cs, T. Go o dwill, J. G. Hoskins, and P . Nekraso v. Integral equations for flexural scattering problems with p eriodic b oundaries. arXiv pr eprint arXiv:2603.00366 , 2026. [2] T. Askham, C. Borges, J. Hoskins, and M. Rachh. Random walks in frequency and the reconstruction of obstacles with cavities from m ulti-frequency data. Journal of Scientific Computing , 98(1):15, 2024. [3] T. Askham, T. Go odwill, J. G. Hoskins, P . Nekrasov, and M. Rachh. Surface lay ers and linearized w ater wa ves: a b oundary integral equation framework. arXiv pr eprint arXiv:2510.22748 , 2025. [4] T. Askham, J. G. Hoskins, P . Nekrasov, and M. Rac hh. Integral equations for flexural- gra vity wa ves: analysis and numerical methods. arXiv pr eprint arXiv:2501.00887 , 2025. [5] T. Askham, M. Rac hh, M. O’Neil, J. Hoskins, D. F ortunato, S. Jiang, F. F ryklund, T. Go o dwill, H. Y ang, and H. Zhu. ch unkIE: a MA TLAB integral equation to olb o x. https://github.com/fastalgorithms/chunkie , 2024. [6] Y. M. Assylb ek ov and Y. Y ang. Determining the first order p erturbation of a p olyhar- monic operator on admissible manifolds. Journal of Differ ential Equations , 262(1):590– 614, 2017. [7] C. Borges, M. Rachh, and L. Greengard. On the robustness of in v erse scattering for p enetrable, homogeneous ob jects with complicated b oundary . Inverse Pr oblems , 39(3):035004, 2023. [8] L. Bourgeois and C. Hazard. On well–posedness of scattering problems in a Kirchhoff– Lo ve infinite plate. SIAM Journal on A pplie d Mathematics , 80(3):1546–1566, 2020. [9] L. Bourgeois and A. Reco quilla y . The linear sampling metho d for Kirchhoff-Lo ve infi- nite plates. Inverse Pr oblems and Imaging , 14(2):363–384, 2020. 25 [10] Y. Chang and Y. Guo. An optimization metho d for the inv erse scattering problem of the biharmonic wa v e. Communic ations on A nalysis & Computation (CA C) , 1(2), 2023. [11] Y. Chang, Y. Guo, T. Yin, and Y. Zhao. Analysis and computation of an in verse source problem for the biharmonic w av e equation. Inverse Pr oblems , 40(11):115011, 2024. [12] A. P . Choudhury and H. Heck. Stability of the inv erse b oundary v alue problem for the biharmonic operator: Logarithmic estimates. Journal of Inverse and Il l-p ose d Pr oblems , 25(2):251–263, 2017. [13] D. Colton and R. Kress. Inverse ac oustic and ele ctr omagnetic sc attering the ory , v ol- ume 93. Springer, 1998. [14] P . Dráb ek, G. Holub o vá, A. Matas, and P . Nečesal. Nonlinear mo dels of susp ension bridges: discussion of the results. Applic ations of Mathematics , 48(6):497–514, 2003. [15] H. W. Engl, M. Hank e, and A. Neubauer. R e gularization of inverse pr oblems , volume 375. Springer Science & Business Media, 1996. [16] M. F arhat, S. Guenneau, and S. Eno c h. Ultrabroadband elastic cloaking in thin plates. Physic al r eview letters , 103(2):024301, 2009. [17] M. F arhat, S. Guenneau, and S. Eno c h. Broadband cloaking of b ending w av es via homogenization of multiply p erforated radially symmetric and isotropic thin elastic plates. Physic al R eview B—Condense d Matter and Materials Physics , 85(2):020301, 2012. [18] M. F arhat, S. Guenneau, S. Eno c h, and A. B. Mo vc han. Cloaking bending w av es propagating in thin elastic plates. Physic al R eview B—Condense d Matter and Materials Physics , 79(3):033102, 2009. [19] P . F arkas. Mathematic al F oundations for F ast Metho ds for the Biharmonic Equation . PhD thesis, The Universit y of Chicago, 1989. (also found under T ec hnical Report 765, Departmen t of Computer Science, Y ale Universit y , 1990). [20] S. S. Ganti, T.-W. Liu, and F. Semp erlotti. T op ological edge states in phononic plates with embedded acoustic blac k holes. Journal of Sound and V ibr ation , 466:115060, 2020. [21] P . Gao, A. Climente, J. Sánchez-Dehesa, and L. W u. Theoretical study of platonic crystals with p erio dically structured n-b eam resonators. Journal of Applie d Physics , 123(9), 2018. [22] Y. Gao, H. Liu, and Y. Liu. On an inv erse problem for the plate equation with passive measuremen t. SIAM Journal on A pplie d Mathematics , 83(3):1196–1214, 2023. [23] J. Guo, Y. Long, Q. W u, and J. Li. On direct and inv erse obstacle scattering problems for biharmonic wa ves. Inverse Pr oblems , 40(12):125032, 2024. [24] I. Harris, H. Lee, and P . Li. Direct sampling for reco vering a clamp ed cavit y from the biharmonic far-field data. Inverse Pr oblems , 41(3):035013, 2025. [25] I. Harris, P . Li, and G. Ozo c hiaw aeze. Sampling metho ds for the inv erse cavit y scat- tering problem of biharmonic w av es. Inverse Pr oblems , 42(1):015002, 2026. 26 [26] G. C. Hsiao and W. L. W endland. Boundary in tegral equations. In Boundary Inte gr al Equations . Springer, 2021. [27] M. Ikehata. A sp ecial green’s function for the biharmonic operator and its application to an inv erse b oundary v alue problem. Computers & Mathematics w ith A pplic ations , 22(4-5):53–66, 1991. [28] K. Krup c hyk, M. Lassas, and G. Uhlmann. Determining a first order p erturbation of the biharmonic op erator b y partial b oundary measurements. Journal of F unctional A nalysis , 262(4):1781–1801, 2012. [29] J. Li and J. Zou. A direct sampling metho d for inv erse scattering using far-field data. arXiv pr eprint arXiv:1206.0727 , 2012. [30] A. E. Lindsay , B. Quaife, and L. W endelb erger. A boundary in tegral equation metho d for mo de elimination and vibration confinement in thin plates with clamp ed p oints. A dvanc es in Computational Mathematics , 44(4):1249–1273, 2018. [31] A. E. Lindsay , M. J. W ard, and T. K olokolnik ov. The transition to a p oin t constraint in a mixed biharmonic eigen v alue problem. SIAM Journal on Applie d Mathematics , 75(3):1193–1224, 2015. [32] B. Liu. Stability estimates in a partial data inv erse b oundary v alue problem for bi- harmonic op erators at high frequencies. Inverse Pr oblems and Imaging , 14(5):783–796, 2020. [33] B. Liu. Stabilit y estimates for an in verse b oundary v alue problem for biharmonic op er- ators with first order p erturbation from partial data. Inverse Pr oblems , 40(6):065001, 2024. [34] V. Memmolo, E. Monaco, N. D. Boffa, L. Maio, and F. Ricci. Guided w av e propagation and scattering for structural health monitoring of stiffened composites. Comp osite Structur es , 184:568–580, 2018. [35] M. H. Meylan, M. Ilyas, B. P . Lamic hhane, and L. G. Bennetts. Sw ell-induced flexural vibrations of a thick ening ice shelf ov er a shoaling seab ed. Pr o c e e dings of the R oyal So- ciety A: Mathematic al, Physic al and Engine ering Scienc es , 477(2254):20210173, 2021. [36] S. H. Mousavi, A. B. Khanikaev, and Z. W ang. T op ologically protected elastic w av es in phononic metamaterials. Natur e c ommunic ations , 6(1):8682, 2015. [37] P . Nekrasov and D. R. MacA y eal. Ocean wa ve blo c king b y p eriodic surface rolls fortifies arctic ice shelves. Journal of Glaciolo gy , 69(278):1740–1750, 2023. [38] P . Nekraso v, Z. Su, T. Askham, and J. G. Hoskins. Boundary integral form ulations for flexural wa v e scattering in thin plates. Journal of Computational Physics , 542:114355, 2025. [39] A. Pelat, F. Gautier, S. C. Conlon, and F. Semp erlotti. The acoustic black hole: A review of theory and applications. Journal of Sound and V ibr ation , 476:115316, 2020. [40] R. C. Reilly . Applications of the hessian op erator in a riemannian manifold. Journal of Differ ential Ge ometry , 12:465–482, 1977. 27 [41] P . Seide. Compressiv e Buckling of a Long Simply Supported Plate on an Elastic F oundation. Journal of the A er osp ac e Scienc es , 25(6):382–384, 1958. [42] T. Tyni and M. Harju. Inv erse backscattering problem for p erturbations of biharmonic op erator. Inverse Pr oblems , 33(10):105002, 2017. [43] Y. Y ang. Determining the first order p erturbation of a bi-harmonic op erator on b ounded and unbounded domains from partial data. Journal of Differ ential Equa- tions , 257(10):3607–3639, 2014. [44] D. Y u, Y. Liu, H. Zhao, G. W ang, and J. Qiu. Flexural vibration band gaps in Euler- Bernoulli b eams with lo cally resonant structures with tw o degrees of freedom. Phys. R ev. B , 73:064301, F eb 2006. [45] Y. Zhao, G. W ei, and Y. Xiang. Plate vibration under irregular internal supp orts. International Journal of Solids and Structur es , 39(5):1361–1383, 2002. [46] M. Zhou, J. Han, M. Rachh, and C. Borges. A neural netw ork warm-start approac h for the in verse acoustic obstacle scattering problem. Journal of Computational Physics , 490:112341, 2023. 28