Rozważmy przypadek, w którym końce pręta są swobodne. Zakładamy, że znana jest taka liczba NN, że rzeczywista liczba pęknięć nn spełnia nierówność nNn \leq N. Dodatkowo, zakładamy, że pęknięcia są „małe” w określonym sensie, co formalnie oznacza, że suma ich uogólnionych giętkości w przeliczeniu na odpowiednie wartości własne jest wystarczająco mała. W tym przypadku możliwa jest rekonstrukcja tych pęknięć na podstawie 2N2N wartości własnych układu z warunkami swoboda-swoboda.

Niech λk0\lambda_k^0 oznaczają wartości własne nieuszkodzonego pręta, a λk\lambda_k – wartości własne pręta z pęknięciami. Wprowadzając oznaczenia ξk=λk\xi_k = \sqrt{\lambda_k} oraz cˉi=EAci\bar{c}_i = EA c_i, przyjmujemy, że suma i=1ncˉi2λ2N0\sum_{i=1}^n |\bar{c}_i|^2 \ll \lambda_{2N}^0. To założenie prowadzi do reprezentacji macierzy QQ, której wyznacznik zeruje się dla ξk\xi_k będących pierwiastkami równania własnego z pęknięciami.

Analizując rozwinięcie Taylora funkcji detQ(cˉ,ξ)\det Q(\bar{c}, \xi) w otoczeniu punktu (0,,0,ξk0)(0, \ldots, 0, \xi_k^0), uzyskujemy relację, która pozwala związać zmiany wartości własnych z parametrami cˉi\bar{c}_i oraz lokalizacjami pęknięć xix_i. Kluczowe są tutaj pochodne wyznacznika względem cˉi\bar{c}_i i ξ\xi, które przyjmują postać wyrażeń trygonometrycznych zawierających funkcje sinus i cosinus w punktach lokalizacji pęknięć. Umożliwia to sformułowanie układu 2N2N równań względem 2N2N niewiadomych – elastyczności cic_i^* oraz lokalizacji xix_i.

Zauważalna jest niejednoznaczność rozwiązania – jedna sprężyna umieszczona w punkcie xix_i i posiadająca giętkość cic_i^* daje taki sam wkład do układu równań jak dwie sprężyny zlokalizowane w punktach xix_i oraz lxil - x_i, jeśli spełniają ci1+ci2=cic_i^{*1} + c_i^{*2} = c_i^*. W praktyce jednak możliwe jest jednoznaczne określenie wartości sin2(πxi/l)\sin^2(\pi x_i/l) oraz agregowanych giętkości cic_i^*.

Aby osiągnąć to jednoznaczne rozwiązanie, przekształca się układ do postaci trygonometrycznej, wykorzystując znane tożsamości sum sinusów podniesionych do potęg. W rezultacie uzyskuje się układ równań postaci:

i=1Ncizik=dk,zi=sin2(πxil),k=1,2,...,2N\sum_{i=1}^N c_i^* z_i^k = d_k, \quad z_i = \sin^2\left(\frac{\pi x_i}{l}\right), \quad k = 1, 2, ..., 2N

Taki układ pojawia się również w innych dziedzinach analizy funkcji meromorficznych. Rozwiązanie przebiega dwuetapowo. Najpierw konstruowany jest wielomian PN(z)P_N(z) o pierwiastkach ziz_i, na podstawie którego rozwiązuje się układ równań liniowych względem współczynników wielomianu. Po znalezieniu pierwiastków ziz_i, wraca się do układu podstawowego, w którym traktuje się ziz_i jako znane, a cic_i^* jako niewiadome.

Z uwagi na możliwość, że n<Nn < N, uzyskane pierwiastki zawierają zarówno rzeczywiste pozycje pęknięć, jak i pierwiastki „fałszywe”. Selekcja właściwych rozwiązań odbywa się poprzez kryteria:

  1. zi(0,1)z_i \in (0,1), czyli fizyczna lokalizacja musi być w obrębie pręta.

  2. ci>0c_i^* > 0, czyli fizyczna giętkość musi być dodatnia.

  3. Jeżeli wartość cic_i^* jest niewielka w porównaniu z maksymalną z uzyskanych wartości, uznaje się ją za nieistotną z punktu widzenia diagnostyki – dotyczy to albo pęknięć o znikomej wielkości, albo artefaktów obliczeniowych.

Powyższe podejście umożliwia identyfikację położenia najbardziej krytycznych uszkodzeń oraz ich względnych parametrów sztywnościowych, nawet w przypadku, gdy pełna jednoznaczność rozwiązania nie jest możliwa. W praktyce, skuteczność tej metody zależy od dokładności pomiaru częstotliwości drgań własnych i precyzji numerycznej rozwiązywanych układów równań. Literatura techniczna potwierdza skuteczność podejścia dla przypadków jednego lub dwóch pęknięć, jednak powyższe rozwinięcie pozwala na jego uogólnienie do dowolnej liczby małych pęknięć.

Warto podkreślić, że identyfikacja pęknięć poprzez analizę spektrum własnego bazuje na teorii odwrotnego problemu drgań, która ma zastosowanie nie tylko w diagnostyce konstrukcji, lecz także w zagadnieniach analizy funkcji, optymalizacji oraz teorii sterowania.

W analizowanej metodzie zakłada się znajomość przybliżonych wartości własnych nieuszkodzonego układu. Jednakże w zastosowaniach praktycznych często dostępne są jedynie pomiary spektralne uszkodzonego pręta. Użycie dwóch różnych zestawów warunków brzegowych (np. swobodny-swobodny oraz zamocowany-swobodny) pozwala ominąć to ograniczenie, umożliwiając pełną rekonstrukcję zarówno liczby pęknięć, jak i ich parametrów – co zostało ujęte w twierdzeniu o jednoznacznej identyfikowalności pęknięć przy znajomości dwóch spektrów. Kluczowe w tym ujęciu jest przekształcenie równania różniczkowego z nieciągłościami (funkcjami delta Diraca) do postaci spektralnej, w której wpływ pęknięcia modelowany jest jako zaburzenie potencjału w równaniu Sturm–Liouville’a.

Jakie wyzwania niesie ze sobą odwrotne problemy w obrazowaniu fotoakustycznym i ich związki z nanocząstkami?

Problemy odwrotne w obrazowaniu fotoakustycznym, szczególnie w zastosowaniach z nanopartykulami jako środkami kontrastowymi, są coraz bardziej popularnym obszarem badań naukowych. Ich rozwiązanie pozwala na uzyskanie nowych, precyzyjnych metod diagnostycznych w medycynie i biologii. W niniejszym rozdziale omawiamy matematyczne podstawy tych problemów, ich podział na etapy oraz wyzwania związane z niestacjonarnym charakterem układów fotoakustycznych, w tym przypadki nieregularnych źródeł oraz zastosowanie nanocząsteczek jako kontrastujących agentów.

Podstawowym celem klasycznego modelu fotoakustycznego jest wykorzystanie zmierzonego ciśnienia akustycznego p(,)p(\cdot, \cdot) na granicy obszaru Ω×(0,T)\partial \Omega \times (0, T), wygenerowanego przez kilka iluminacji uiu_i, do rekonstrukcji współczynników akustycznych i optycznych w obrębie obszaru. Zwykle, problem ten dzieli się na dwa główne etapy: inwersję akustyczną oraz inwersję optyczną. Inwersja akustyczna polega na wykorzystaniu bocznych wartości ciśnienia akustycznego do rekonstrukcji funkcji źródłowej i współczynników akustycznych (np. masowej gęstości ρ\rho oraz prędkości cc). Z kolei inwersja optyczna opiera się na wcześniejszym wyznaczeniu funkcji źródłowej i wykorzystaniu jej do rekonstruowania współczynników optycznych, jak np. absorpcji (μab\mu_{ab}) czy rozpraszania (μsc\mu_{sc}).

W literaturze dotyczącej klasycznego obrazowania fotoakustycznego, przyjęto założenie, że gęstość masy ρ()\rho(\cdot) jest stała, a proces inwersji dzieli się na dwie główne części: inwersję akustyczną i optyczną. W przypadku sferycznego obszaru Ω\Omega oraz stałej prędkości c()c(\cdot), rekonstrukcja jest realizowana za pomocą transformacji Radona, jak to zostało opisane przez Natterera (2001). W przypadku ogólnych kształtów obszaru Ω\Omega oraz zmiennej prędkości c()c(\cdot), stosuje się rozkład spektralny, jak wykazali Agranovsky i in. (2009). Bardziej złożona sytuacja dotyczy przypadków, w których prędkość fali c()c(\cdot) jest zmienną funkcją, której wartość nie jest znana. W takich przypadkach jedynie niektóre wyniki dotyczące unikalności i stabilności rozwiązania zostały opracowane (Stefanov i Uhlmann, 2009; 2013).

Znajomość wyników inwersji akustycznej oraz optycznej pozwala na skuteczną rekonstrukcję współczynników akustycznych i optycznych z ograniczonej liczby pomiarów. Na przykład, Kirsch i Scherzer (2012) opracowali metodę rekonstrukcji dwóch współczynników (prędkości c()c(\cdot) oraz funkcji źródłowej) z pomiarów 6-wymiarowych. Takie podejście pozwala na precyzyjne wyodrębnienie parametrów akustycznych w przypadkach skomplikowanych geometrii.

W odniesieniu do inwersji optycznej, istotną rolę odgrywają wyniki dotyczące unikalności oraz stabilności rozwiązań. W badaniach Bal i Ren (2011) oraz Bal i Uhlmann (2010) wykazano, że pod pewnymi warunkami degeneracji rozwiązań, inwersja optyczna jest możliwa, a stabilność wyników może być zapewniona odpowiednim doborem danych wejściowych. Inne prace, takie jak Bonnetier i in. (2022), sugerują wykorzystanie punktowych źródeł światła do uzyskania oszacowań stabilności wyników rekonstrukcji.

Problem odwrotny w kontekście fotoakustycznego obrazowania przy użyciu nanocząsteczek jako agentów kontrastowych stanowi osobny obszar badań, który wymaga uwzględnienia unikalnych właściwości tych materiałów. Nanocząsteczki, dzięki swoim rozmiarom i specyficznym właściwościom optycznym, mogą pełnić rolę niezwykle skutecznych agentów kontrastujących, zwiększając precyzję uzyskiwanych obrazów. W tym kontekście, model fotoakustyczny z nanocząsteczkami opiera się na połączeniu trzech równań: dla ciśnienia akustycznego pp, pola elektrycznego EE, oraz temperatury TT. Źródło ciepła w tym modelu jest związane z wprowadzeniem nanocząsteczek, które poprzez oddziaływanie z polem elektromagnetycznym generują efekt fotoakustyczny.

Równania opisujące ten proces w formie matematycznej obejmują m.in. równanie dla ciśnienia akustycznego w przestrzeni R3\mathbb{R}^3 oraz dodatkowe warunki początkowe, w tym specyficzne warunki dotyczące pola elektrycznego i temperatury w układzie. W kontekście tego modelu, źródło wytwarzane jest przez całkowite pole elektryczne, a nie tylko jego intensywność, co sprawia, że model jest bardziej skomplikowany, ale również bardziej realistyczny w przypadku zastosowań medycznych. Wymaga to jednak uwzględnienia dodatkowych aspektów matematycznych, jak np. analiza osobliwości funkcji Green’a, która pozwala na uzyskanie dokładnych oszacowań dla pola ciśnienia generowanego przez źródło.

Równania dotyczące pola elektrycznego opierają się na układzie równań Maxwella, który jest podstawą dla wszystkich rozważań o elektromagnetyzmie w tym kontekście. Model zakłada również, że permittivity w obszarze nanocząsteczek jest opisana przez model Lorentza, a przewodnictwo magnetyczne jest zaniedbywane, ponieważ układ składa się z nie-magnetycznych tkanek i nanocząsteczek.

Aby uzyskać stabilność i precyzyjność rekonstrukcji, konieczne jest stosowanie szczegółowych rozszerzeń asymptotycznych i zaawansowanych metod analizy, takich jak te, które są wykorzystywane do analizy funkcji Green’a w kontekście fal akustycznych. Podejście to jest szczególnie przydatne w trudniejszych przypadkach, gdzie geometria układu nie jest regularna, a rozwiązania są nieliniowe.

Ważnym aspektem, który należy wziąć pod uwagę, jest odpowiednia analiza geometrii obszaru oraz właściwości medium. Kluczowym warunkiem dla poprawności rekonstruowanych parametrów jest założenie, że każde dwa punkty w obszarze Ω\Omega są połączone unikalną geodezyjną, co pozwala na wykorzystanie specyficznych metod matematycznych, takich jak rozwinięcia osobliwości.

Jak można oszacować rozmiar inkluzji na podstawie energii odkształcenia?

W badaniach nad odwrotnymi problemami w strukturach nanoskali, jednym z kluczowych zagadnień pozostaje szacowanie rozmiaru inkluzji materiałowej ukrytej w domenie, na podstawie danych brzegowych i energii mechanicznej układu. Podejście oparte na analizie energii odkształcenia umożliwia wyprowadzenie zarówno dolnych, jak i górnych oszacowań objętości inkluzji, przy spełnieniu precyzyjnie sformułowanych warunków regularności, symetrii i wypukłości tensora materiałowego.

Rozważmy domenę ograniczoną ΩR2\Omega \subset \mathbb{R}^2, której brzeg jest klasy C3,1C^{3,1} i spełnia warunek separacji od inkluzji Ω~Ω\tilde{\Omega} \Subset \Omega, tj. dist(Ω~,Ω)d0r0\text{dist}(\tilde{\Omega}, \partial\Omega) \geq d_0 r_0 dla pewnej stałej d0>0d_0 > 0. Materiały w domenie opisane są przez tensory sprężystości drugiego i trzeciego rzędu P,P~,Q,Q~P, \tilde{P}, Q, \tilde{Q}, spełniające warunki symetrii, silnej wypukłości oraz skokowe warunki przejścia na granicy inkluzji.

Zasadniczym narzędziem do uzyskania oszacowań jest tzw. lemat energetyczny, który wiąże różnicę energii rozwiązania niezakłóconego u0u_0 i zakłóconego uu – odpowiednio W0W_0 i WW – z integralnym udziałem energii odkształcenia w inkluzji. Dokładnie:

ηξ0t3Ω~D2u02+t2D3u02W0W(δ1)ξ1t3Ω~D2u02+t2D3u02,\eta^* \xi_0^* t^3 \int_{\tilde{\Omega}} |D^2 u_0|^2 + t^2 |D^3 u_0|^2 \leq W_0 - W \leq (\delta^* - 1) \xi_1^* t^3 \int_{\tilde{\Omega}} |D^2 u_0|^2 + t^2 |D^3 u_0|^2,

gdzie η,δ,ξ0,ξ1\eta^*, \delta^*, \xi_0^*, \xi_1^* są stałymi zależnymi od materiałowych parametrów inkluzji i ośrodka macierzystego. Warunki te zapewniają, że energia zmagazynowana w obszarze inkluzji zależy w sposób jednoznaczny i kontrolowany od samego pola odkształceń w nieobecności inkluzji.

Wynik ten umożliwia wyprowadzenie dolnego ograniczenia objętości inkluzji:

Ω~W0WC1r02W0,|\tilde{\Omega}| \geq \frac{W_0 - W}{C_1 r_0^{ -2} W_0},

gdzie C1>0C_1 > 0 zależy wyłącznie od geometrii i parametrów materiałowych. Z kolei oszacowanie górne wymaga zastosowania bardziej zaawansowanych technik kontynuacji jednoznacznej, w szczególności tzw. nierówności trzech kul oraz nierówności podwajania dla hessianu pola przemieszczenia D2u0D^2 u_0. Istotne jest tu uwzględnienie możliwości zanikania pochodnych drugiego i trzeciego rzędu w głębi domeny, co wymaga precyzyjnego ujęcia propagacji małości – znanego jako Lipschitz propagation of smallness.

Rezultaty te umożliwiają ostatecznie zapisanie nierówności:

Ω~(W0WC2r02W0)p,|\tilde{\Omega}| \leq \left( \frac{W_0 - W}{C_2 r_0^{ -2} W_0} \right)^p,

gdzie C2>0C_2 > 0, a p>1p > 1 zależą od zestawu parametrów opisujących materiał (moduły Lamégo, parametry izotropii), geometrię domeny oraz dane brzegowe wyrażone poprzez normy Soboleva rozmaitych rzędów.

Konstrukcja tych oszacowań opiera się na teorii regularności dla równań różniczkowych cząstkowych eliptycznych wyższego rzędu, zanurzeniach Soboleva, oszacowaniach typu Agmon-Douglis-Nirenberg oraz fundamentalnych własnościach rozwiązań słabych, w tym unikalności i ciągłości względem danych wejściowych. Kluczowe jest również zastosowanie pokryć lokalnych oraz szacowań wnętrzowych dla funkcji klasy H6H^6, umożliwiających kontrolę norm LL^\infty poprzez globalną energię układu.

Nie mniej istotne pozostają warunki kompatybilności danych brzegowych, zapewniające istnienie rozwiązania zadania Neumanna i jego unikalność po odpowiednim znormalizowaniu. Dane te muszą należeć do przestrzeni Soboleva ujemnych rzędów: H5/2(Ω),H3/2(Ω),H1/2(Ω)H^{ -5/2}(\partial\Omega), H^{ -3/2}(\partial\Omega), H^{ -1/2}(\partial\Omega), co uwzględnia realia modelowania odwrotnego problemu w rzeczywistych eksperymentach fizycznych lub pomiarach numerycznych.

Istotne jest zrozumienie, że skuteczność tej metody zależy od jakości danych wejściowych oraz poprawnej identyfikacji parametrów materiałowych, takich jak moduły Lamégo μ,λ\mu, \lambda i ich odpowiedniki w inkluzji μ~,λ~\tilde{\mu}, \tilde{\lambda}. Skokowe warunki przejścia muszą być sformułowane w sposób spójny, np. poprzez różnice κ~κ\tilde{\kappa} - \kappa, gdzie κ=2μ(2μ+3λ)2μ+λ\kappa = \frac{2\mu(2\mu+3\lambda)}{2\mu + \lambda}, co pozwala wyrazić przejścia między ośrodkiem a inkluzją w sposób materiałowo sensowny.

W ramach szacowań nie można pominąć wpływu nieregularności geometrycznych oraz lokalnych anomalii w danych brzegowych. Te mogą istotnie zaburzyć propagację informacji z brzegu do wnętrza domeny, co bezpośrednio wpływa na dokładność oszacowania objętości inkluzji. Z tego względu, teoria propagacji małości i kontynuacji jednoznacznej odgrywa fundamentalną rolę w pełnym ujęciu problemu.