RESUMO
Nesta contribuição analisaremos o problema inverso de identificação do coeficiente de rigidez em vigas modeladas pela equação de Euler-Bernoulli, a partir de medidas da deflexão. Apresentaremos o problema na forma de uma equação de operador parâmetro-para-medidas, para o qual, provaremos propriedades importantes, como compacidade e continuidade. Mostraremos ainda que o operador parâmetropara-medidas é Fréchet diferenciável e que satisfaz a condição do cone tangente. Essas propriedades são suficientes para recuperarmos de forma estável e convergente (método de regularização) o coeficiente de rigidez através de métodos iterativos como Landweber e Steepest descent. Por fim, apresentamos os efeitos de estabilidade das soluções aproximadas com relação as medidas com diferentes níveis de ruídos, através de exemplos numéricos.
Palavras-chave: identificação do coeficiente de rigidez; métodos iterativos; regularização; problemas inversos
ABSTRACT
In this contribution, we analyze the identifiability of the stiffness coefficient inverse problem in beams modeled by the Euler-Bernoulli equation, from measurements of the beam deflection. We present the problem in the form of a parameter-to-measure operator equation, for which we prove important properties, such as compactness and continuity. We also show that the parameter-to-measure operator is Fréchet differentiable and that it satisfies the tangential cone condition. These properties are sufficient to recover in a stable and convergent way (regularization method) the stiffness coefficient through iterative methods such as Landweber and Steepest descent. Finally, we present the stability effects of the approximate solutions concerning measurements with different noise levels through numerical examples.
Keywords: beam stiffness coefficient identification; iterative methods; regularization; inverse problems
1 INTRODUÇÃO
As vigas são um dos elementos mais básicos em sistemas estruturais, as quais são as responsáveis pela sustentação em projetos de engenharia 9. Manter esses sistemas seguros continua a ser um tema muito importante da pesquisa. Como os danos estão relacionados às alterações dos materiais internos da estrutura (propriedades dos coeficientes do modelo), conhecer essas propriedades de forma não-destrutiva possibilitam avaliar as condições de risco 3), (11), (14), (19.
Nesta contribuição, assumiremos que a deflexão u(x) de uma viga estática unidimensional, de comprimento L, devido a força transversal f(x)∈L2[0,L], possa ser modelada pela equação de Euler-Bernoulli 4), (9
onde, a(x)=E(x)I(x) é coeficiente de rigidez da viga, resultante do produto do módulo de elasticidade E e o momento de inércia I.
Utilizando o balanceamento das força e dos momentos, a equação (1.1) pode ser dividida em um sistema de equações diferenciais ordinárias de segunda ordem 15, dado por
Existe na literatura duas abordagens associadas à teoria de vigas dadas pela equação de EulerBernoulli (1.1) (equivalentemente1 (1.2)). A primeira delas diz respeito ao problema tradicional na mecânica, a qual consiste em avaliar a resposta de uma viga, quando conhecemos a força aplicada e os parâmetros da estrutura 17. Matematicamente, significa determinar a deflexão u = u(x) satisfazendo a equação (1.1), quando o coeficiente de rigidez a=a(x)>0 e a força atuante f(x)≥0 são dados, conjuntamente com quatro condições de contorno apropriadas, e.g., 9. A este problema denominamos problema direto 11), (17), (19. Trataremos da boa colocação do problema direto na Seção 2.
A segunda abordagem está relacionada aos problemas inversos: estes consistem em um problema de identificação da fonte f(x), quando u, a e as condições de contorno são conhecidas ou um problema de identificação do coeficiente de rigidez a(x), quando u, f e as condições de contorno são conhecidas 14), (15), (18), (19. É importante realçar o fato de que ambos os problemas inversos relacionados com a equação de Euler-Bernoulli são mal-postos no sentido de Hadamard 1), (7, pois, na prática, ambos os problemas dependem da diferenciação de medidas da deflexão u, que inclui ruídos, e, por isso, instável com relação as medidas, e.g. 6), (10), (12), (15 e referências.
O problema de identificação da fonte f foi investigado em 8 e referências. O problema de identificação do coeficiente de rigidez a(x) foi investigado anteriormente em 14), (15), (18), (19 para a equação de Euler-Bernoulli dada por (1.1). Em 19, um coeficiente de rigidez constante por partes é analisado. Uma fórmula explícita para obter tal coeficiente é apresentado. No entanto, a técnica utilizada por 19 exige o conhecimento prévio de pontos da viga onde os coeficientes mudam de valor. Tal hipótese não tem sentido prático, pois ter esse conhecimento implica em saber onde a viga sofreu rupturas. Além disso, a equação obtida depende de um quociente que envolve a segunda derivada dos dados (da deformação da viga). Já é bem conhecido que o problema de derivação nos dados é um problema instável 6), (12. Em 14 uma fórmula explícita também é proposta para obter o coeficiente de rigidez, mas utilizando as funções de Green. Contudo, a equação obtida também depende da derivada dos dados. Então para tratar o problema da diferenciação nos dados, 14 propõe a utilização de métodos de regularização como Molificação e Landweber afim de estabilizar o problema de diferenciação nos dados. Já 15 propõe como método de regularização, o método bem conhecido na literatura denominado método de Tikhonov para determinar uma aproximação estável com relação as medidas. Resultados de unicidade de identificação do coeficiente a(x) foram derivados em 15 e com condições menos restritivas na suavidade de a(x) em 20.
A identificação do coeficiente de rigidez a(x), além de ser instável com relação as medidas, ainda impõem restrições extras com relação a identificabilidade. De fato, note que, se u*=u(a*) é uma solução da equação (1.1), então, para qualquer função linear ζ, temos que u ∗ satisfaz a equação (1.1), com coeficiente de rigidez ˆa=a*+ζ(d2u*dx2)-1, sempre que d2u*dx2≠0. Portanto, uma condição necessária para que a equação (1.1) possua uma solução associada a um único coeficiente de rigidez a(x)>0, implica que o conjunto
seja não vazio. Por outro lado, se J é um subconjunto de [0, L] com medida não-nula, então o coeficiente de rigidez a(x) não é mais identificável em J (veja a equação (1.2). As condições que permitam a identificabilidade de a(x) em [0, L] e que J≠∅, incluem u'(0)=u'(L) ou M(0)M(L)=0.
Desta forma, dentre todas as possibilidades para as condições de contorno do tipo Dirichlet e/ou Neumann2, nesta abordagem consideraremos somente as condições de contorno que garantam: a) condições de contorno associadas com vigas estudadas na teoria de elasticidade, b) que J≠∅, c) que a viga não seja estaticamente indeterminada 9. Assim, assumiremos que a equação (1.2) está associada com as seguintes condição de contorno
Nesta contribuição, daremos ênfase ao problema de identificação estável do coeficiente de rigidez a(x). Este problema é importante, pois o conhecimento de parâmetros como coeficiente de rigidez a(x) em (1.2) permite avaliar imperfeições ou falhas nessas estruturas 5),(11),(16. No entanto, não existe uma maneira de obter informações sobre esse parâmetro de maneira direta 10), (15. A única informação disponível são um conjunto de medidas indiretas u δ , sujeitas a erros de nível δ≥0, que dizem respeito a medições da deflexão u(x) da viga, quando conhecida a força atuante f(x) do sistema. Dessa forma, para obtermos um coeficiente a(x) estável e convergente, com relação as medidas u δ , é necessário utilizar alguma estratégia de regularização 6), (10), (12. Um dos métodos mais tradicionais na teoria de problemas inversos é a regularização de Tikhonov 6), (10), (12. Neste trabalho, mostraremos propriedades suficientes do operador parâmetro-paramedidas (veja Seção 2) de forma que métodos iterativos, como Landweber e Steepest descent 1), (10, sejam caracterizados como estratégias de regularização para o problema de identificação do coeficiente de rigidez a(x).
As principais contribuições deste trabalho podem ser resumidas em: na Seção 2 formularemos o problema de identificação do coeficiente a(x) na forma de uma equação de operadores (operador parâmetro-para-medidas). Provaremos que tal operador é contínuo, compacto, Fréchet diferenciável e que satisfaz a condição de cone tangente. Na Seção 3, discutiremos como tais propriedades são suficientes para provar que métodos iterativos, como Landweber e Steepest descent, forneçam coeficientes iterados aδk*, estáveis e convergentes com relação ao nível de ruídos δ nas medidas da deflexão da viga, se a estratégia de parada para as iterações pelo princípio de discrepância são aplicados. Por fim, na Seção 4, apresentaremos alguns resultados numéricos para o problema de identificação do coeficiente a(x), com níveis de ruídos distintos, mostrando a estabilidade das soluções aproximadas pelos métodos iterativos propostos.
2 PROBLEMA COMO UMA EQUAÇÃO DE OPERADORES
Nesta seção, apresentaremos o problema de identificação do coeficiente de rigidez da viga a(x), como uma equação de operadores parâmetros-para-medidas. Provaremos, primeiramente, que tal operador parâmetros-para-medidas está bem definido (veja definição em (2.5) abaixo). Em seguida, provaremos propriedades deste operador, como continuidade e compacidade, caracterizando assim o problema de identificação do coeficiente de rigidez a(x) como um problema malposto no sentido de Hadamard, e.g. 6), (10), (12. Por fim, provaremos que o operador parâmetrospara-medidas é Fréchet diferenciável e que satisfazer a condição do cone tangente. Tais propriedades do operador parâmetros-para-medida formam um conjunto condições necessárias e suficientes para caracterizar alguns métodos iterativos, como os que apresentaremos na Seção 3, como métodos de regularização 6), (10), (12, para problema de identificação do coeficiente de rigidez da viga. Aproveitaremos ainda esta seção para obter resultados preliminares que servirão como o suporte teórico para a implementação numérica na Seção 4.
Para tal, consideraremos as seguintes hipóteses gerais:
H1) O coeficiente de rigidez a∈Ad:={amensurável em[0,L]:ˉa≥a(x)≥a¯>0, com ˉa,a¯,a(0),a(L)conhecidos}. Ainda, assumiremos que ∃a*∈Ad solução de (2.5). A distribuição de forças f∈C[0,L] e satisfaz f(x)>0 ao longo da viga.
H1’) O coeficiente de rigidez a(x) além de satisfazer H1) é tal que ║a ' ║ é uniformemente limitada, ou seja, o coeficiente pertence ao conjunto A={a∈L∞∣λ1≤a≤λ2,∥a'∥≤Q}, onde ║·║ denota a norma L2(0,L) com 0<λ1≤λ2<∞ e Q>0.
H2) As medidas uδ∈L2[0,L] diferem dos dados exatos u∈L2[0,L], por um nível de ruídos δ>0,i.e.,‖u-uδ‖L2[0,L]≤δ.
O primeiro resultado auxiliar para definirmos o operador parâmetro-para-medidas é o seguinte:
Lema 2.1. Considere a hipótese H1) satisfeita. Então, existe uma única solução “fraca” u∈H10[0,L] , i.e., uma função u que satisfaz
para o problema de valor de contorno (1.3)-(1.2).
Proof. Como f∈C[0,L], segue da teoria elíptica clássica 2), (13 que existe uma única M∈C2[0,L] solução3 da primeira equação em (1.2) com condições de contorno (1.3).
Vamos definir a forma bilinear A:H1[0,L]×H10[0,L]→ℝ dada por
e o seguinte funcional linear l:H10[0,L]→ℝ dado por
Notemos que u é solução fraca da segunda equação do sistema (1.2) se e somente se u satisfaz
O próximo passo é mostrar que a forma bilinear (2.2) associada a formulação variacional (2.1) é contínua e coerciva em H1[0,L]. Ou seja, existem constantes ˉC,C¯, independentes de ϕ∈H1[0,L], tais que
e
respectivamente.
De fato, segue da definição da forma bilinear, da desigualdade de Cauchy-Schwarz e da hipótese H1) que
o que prova a continuidade da forma bilinear. A coercividade segue de
onde usamos a desigualdade de Poincaré 2), (13 para obtermos a última estimativa.
O próximo passo consiste em provarmos que o funcional linear l definido em (2.3) é contínuo em H1[0,L]. De fato, segue da definição do funcional l em (2.3) e da desigualdade de Cauchy-Schwarz que
A imersão contínua de H1[0,L] em L2[0,L] finaliza a prova.
Portanto, segue do Lema de Lax-Milgram 2 a existência de uma única função u∈H10[0,L] satisfazendo (2.1). □
Do Lema 2.1 segue que a equação de operadores (operador parâmetro-para-medidas),
onde u(a) é a única solução fraca de (2.1), está bem definida. Assim, podemos considerar o problema de identificação do coeficiente de rigidez a(x), como o problema de determinar a∈Ad na Equação (2.5) para dados medidos uδ∈L2[0,L] que satisfazem a hipótese H2).
2.1 Propriedades do operador
Nesta subseção, provaremos propriedades do operador F em (2.5). O teorema a seguir trata da continuidade e compacidade do operador F na topologia de L2[0,L].
Teorema 2.1. Assuma a hipótese H1) satisfeita. Então, o operador F definido em (2.5) é contínuo e compacto.
Proof. Seja F(an)=u(an):=un e F(a)=u(a):=u as respectivas soluções de (1.3)-(1.2), correspondentes a a n ,a∈Ad. Para provar a continuidade, basta mostrar que F é sequencialmente contínuo visto que L2[0,L] é um espaço separável 2. Dessa forma, consideremos (an)∈Ad com an→a quando n→∞ na topologia de L2[0,L]. Como Ad é fechado, segue que a∈Ad.
Pela linearidade de (1.2), un-u satisfaz
com condições de contorno homogêneas, haja visto que ambas as condições de contorno são as mesmas. Considere a função teste vn:=un-u∈H10[0,L]. Para concluirmos o teorema, basta mostrar que vnL2⟶0 quando anL2⟶a Para demonstrar este fato, consideraremos vn∈H10[0,L], como função teste para o problema (2.6). Portanto,
Integrando por partes o lado esquerdo de (2.7) obtemos
Como vn∈H10[0,L], segue que o termo de fronteira se anula. Assim, da desigualdade de Hölder, segue que
onde ˜C=Ka¯2 pois sup[0,L]|M(x)|≤K, haja visto que M é função contínua em [0, L].
Pela desigualdade de Poincaré 2), (13 aplicado a vn∈H10[0,L] obtemos
onde C é uma constante que depende somente de ˜C e da constante de Poincaré.
Assim, obtemos de (2.9) que se anL2⟶a então vnL2⟶0. Provando assim que o operador F é contínuo de L2[0,L] para L2[0,L].
Para demonstrar a compacidade do operador, considere a sequência (an)∈Ad convergindo fracamente para a. Como Ad é convexo e fechado, é, portanto, fracamente fechado 2. Disto segue que a∈Ad. Vamos mostrar que un→u, onde un=u(an) e u=u(a) são as respectivas soluções fracas de (1.3)-(1.2).
Considerando φ∈C∞0 na formulação variacional, obtemos
Segue da regularidade de M20 e de φ que aMφ∈L2[0,L]. Como por hipótese de an⇀a em L2[0,L] segue que o lado direito de (2.10) converge para zero, isto é,
Com isso, o lado esquerdo de (2.10) converge para zero.
Por outro lado, da Hipótese H1), do teorema de convergência limitada e da densidade de C∞0[0,L] em L2[0,L], temos que un-u⇀z em L2[0,L]. Segue de (2.10) que z satisfaz, no sentido das distribuições, o problema az''=0 com condições de contorno homogêneas, haja visto que u n e u possuem condições de contorno homogêneas como em (1.3). Como a∈Ad, e pelo Princípio do Máximo 2), (21, concluímos que z=0. Portanto, un⇀u. Além disso, tomando φ=uk-u como função teste em (2.10), segue que ‖uk-u‖L2→0. Assim, uk-u convergem fraco e em norma para zero. O que completa a demonstração. □
É bem conhecido da teoria de problemas inversos 10),(12 que a compacidade do operador implica na necessidade de utilizarmos alguma estratégia de regularização, a fim de obtermos soluções estáveis e convergentes para dados ruidosos. Para garantir condições suficientes para o operador (2.5), de modo que possamos utilizar algum método iterativo de regularização, necessitamos do seguinte resultado auxiliar.
Proposição 2.1. Assuma as hipóteses H1) H2) satisfeitas. Seja a∈Ad e u:=u(a) a única solução fraca de (1.2) (dada pelo Lema 2.1). Defina
Então:
-
O operador A(a) é linear, limitado e bijetivo.
-
O operador inverso de A(a), denotado por A−1(a), está bem definido, é linear e limitado.
-
O operador A(a) possui um operador adjunto que satisfazA*(a)v=-(av)'' >, o qual é inversível. O operador adjunto inverso(A*(a))-1:=(A-1(a))*satisfaz(A*(a))-1w=v >, onde v é a única solução de-(av)''=w, com condições de contorno homogêneas.
-
O operador F(a) é Fréchet diferenciável. O operador derivada Fréchet (denotamos por F' (a)) é linear e limitado e satisfaz
O operador adjunto de F' (a) também é um operador linear e limitado e satisfaz
Proof. Que o operador A(a) é linear, segue diretamente da definição. Vamos então provar que A(a) é uma bijeção. Considere w∈L2[0,L]. Como a∈Ad, segue do Teorema de Lax-Milgram 2 que existe uma única solução u∈H10[0,L] de -au''=w com condições de contorno (1.3). Além disso, da regularidade elíptica 2 segue que u∈H2[0,L], o que conclui a afirmativa. Portanto, existe A −1(a) (o operador inverso de A(a)) e também é um operador linear.
Para concluir a prova de i) e ii) nos falta provar que A(a) é limitado. Observe que, como a∈Ad e usando que ‖u‖H2={‖u''‖L2+‖u'‖L2+‖u‖L2}, obtemos
provando que A(a) é, de fato, limitado. O Teorema da Aplicação Aberta 13 garante que A −1(a) é limitado.
A demonstração de iii), segue diretamente da teoria de operadores lineares limitados 13 e de i) e ii). O que temos que verificar é como o operador adjunto atua. Pela definição de adjunto, temos que
para qualquer φ∈C∞0[0,L]. Por densidade de C∞0[0,L]⊂L2[0,L], segue que A*(a)v=-(av)''.
Para provarmos iv), seja h tal que ˜a:=a+h∈Ad. Denote por F(a+h)=u(a+h):=˜u,F(a)=u(a):=u as respectivas soluções de (1.2). Pela linearidade de (1.2), segue que a diferença F(a+h)-F(a)=˜u-u:=z satisfaz
com condições do contorno homogêneas. Como, por hipótese, ũ é solução de (1.2), temos que ˜u''=˜a-1M∈L2[0,L]. Assim, com argumentos similares aos usados no Lema 2.1, garantimos a existência de uma única solução z=z(h)∈H2[0,L]⊂L2[0,L] de (2.16). Segue de (2.16) que esta é linear em h. Defina F'(a)h=:z única solução de (2.16). Assim, F' (a) é um operador linear em h, sempre que a+h∈Ad.
Como consequência da regularidade elíptica 2, obtemos ainda que F'(a)h=:z (única solução de (2.16), satisfaz
Portanto, F' (a) pode ser estendido como um operador linear e limitado (que denotaremos ainda por F' (a)) em L2[0,L]. Este satisfaz ||||F'(a)||||≤a¯-1||||M||||L2[0,L]:=γ-1.
De i) e ii) segue que F'(a)h=z=z(h)=A(a)-1(hu'').
Falta mostrar que o operador F' (a) definido acima é de fato a derivada de Fréchet do operador F, definido em (2.5). Para isso, vamos definir w=F(a+h)-F(a)-F'(a)h=˜u-u-z, para qualquer h∈L2[0,L] com a+h∈Ad. Segue de (1.2) e (2.16) que w satisfaz
com condições de contorno homogêneas. Pelo princípio de continuação única 2, segue que w(x)=0,∀x∈[0,L] é a única solução da equação acima. Logo, o operador linear e limitado F' (a) definido é a derivada de Fréchet de F.
Além do mais, existe o operador adjunto de F' (a), denotado por F' (a)∗, que é linear e limitado. Para acabar a demonstração do item iv) temos que verificar (2.14). Seja h∈C∞0 Portanto,
para qualquer que seja h∈C∞0. Como C∞0[0,L] o qual é denso em L2[0,L], o resultado enunciado segue. Isso acaba a demonstração. □
O próximo resultado é um resultado auxiliar para provarmos que o operador F(a) satisfaz a condição de cone tangente.
Lema 2.2. Assuma as hipóteses H1’) e H2) satisfeitas. Sejam a, ˜a∈A. DefinaF(a)=u(a):=u,F(˜a)=˜u(˜a):=˜ue seja F' (a) a derivada de Fréchet do operador F (o qual existe pela Proposição 2.1). Então
para qualquerw∈C∞0, onde(A-1(a))*representa o operador adjunto do operador A−1(a) (o qual existe pela Proposição 2.1).Proof. Segue das hipóteses que
Multiplicando o termo (˜u-u) na equação (2.19) pelo operador identidade I=A-1(a)A(a), pela linearidade e continuidade de A −1(a) e como A −1(a) é linear e limitado, definido entre espaços de Hilbert, sabemos que existe um operador adjunto, o qual foi denotado por (A-1(a))*.
Fazendo as manipulações apropriadas e colocando em evidência (˜a-a) obtemos de (2.19) que
concluindo a demonstração. □
Para garantirmos a condição do cone tangente precisamos fazer algumas estimativas que serão apresentadas nos seguintes lemas:
Lema 2.3. Assumindo as hipóteses do Lema 2.2 satisfeitas. Então existe uma constante C1>0 , tal que
Proof. Segue do Lema 2.1 que ‖(a(A-1(a))*w)''‖L2=‖(av)''‖L2=‖w‖L2.
Agora a estimativa segue do fato de a∈A⊂Ad, com C1=a¯-1/2. □
Lema 2.4. Assumindo as hipóteses do Lema 2.2 satisfeitas. Suponha ainda que a, ˜a∈A . Então existe uma constante C2>0 tal que
Proof. Utilizando a desigualdade de Cauchy-Schwarz temos
e como w∈C∞0, segue da desigualdade de Poincaré 2 e dos mesmos argumentos do Lema 2.3 que
Por outro lado, como ((˜a-a)a-1)'=a'(˜a-a)+a(˜a-a)a2, temos da desigualdade triangular e do fato de a∈A, que
de onde o resultado enunciado segue. □
Lema 2.5. Assumindo as hipóteses do Lema 2.2, e do Lema 2.4 satisfeitas. Suponha ainda para ζ=0 ou ζ=L tenhamos a(ζ)=˜a(ζ) e a'(ζ)=˜a'(ζ) . Então existe uma constante C3>0 tal que
Proof. De acordo com as hipóteses do Lema e de uma integração por partes temos
Agora a estimativa segue com os argumentos do Lema 2.4. □
Vamos provar agora o último resultado antes de verificarmos que a condição do cone tangente é válida para o problema de identificação do coeficiente.
Proposição 2.2. Assumindo as hipóteses dos Lemas 2.2 - 2.5 sejam satisfeitas. Então, existe uma constante K (que depende somente de L, ā, a), tal que
Proof. Primeiramente, da hipótese de a∈A, da integração por partes duas vezes e algumas manipulações, obtemos
onde usamos novamente a notação de que F(˜a)=˜u e F(a)=u.
Derivando [(˜a-a)a-1a((A-1(a))*w]'' e usando a desigualdade triangular, obtemos que
Agora, o resultado segue diretamente das estimativas acima e dos Lemas 2.3, 2.4 e 2.5, respectivamente. □
O próximo teorema garante a condição de cone tangente para o operador F(a).
Teorema 2.2. Assuma as hipóteses da Proposição 2.2 satisfeitas. Seja K como na Proposição 2.2. Então, a condição do cone tangente
é satisfeita, para η:=K‖˜a-a‖H1<1/2 , sempre que ã, a∈B2ρ(a*) , para ρ suficientemente pequeno.
Proof. Segue do Lema 2.2 e da Proposição 2.2 que existe uma constante K>0 (independente de a, ã) tal que
Tomando o supremo ‖w‖C∞0≤1 em ambos os lados da equação acima, juntamente com a densidade de C∞0 em L 2, obtemos
Disto, a demonstração segue da escolha de ρ como na hipótese.
3 MÉTODOS ITERATIVOS PARA IDENTIFICAÇÃO DO COEFICIENTE DE RIGIDEZ
Existe uma gama bastante variada de métodos iterativos de regularização para problemas inversos. Um apanhado geral de tais métodos podem ser consultados em 1), (10 e referências.
Neste trabalho, para a identificação do coeficiente de rigidez a(x), focaremos em métodos cuja iteração é descrita por
onde, F(·) e F ' (·)∗ são os operadores parâmetros-para-medidas definido em (2.5) e o operador adjunto da derivada de Frechet como na Proposição 2.1item iv). O sufixo δ representa a iteração para dados com ruídos u δ . A iteração para dados exatos δ=0 é definida da mesma forma, sem o sufixo δ.
Nesta contribuição consideraremos somente os casos em que ωδk é definido por
cuja escolha define os métodos iterativos de steepest descent ou Landweber, respectivamente. Veja 1), (10 para maiores detalhes da derivação de tais métodos iterativos.
O próximo teorema resume os resultados de convergência e estabilidade dos métodos iterativos que serão utilizados para os testes numéricos deste trabalho.
Teorema 3.3. [Convergência e estabilidade] Assuma as hipóteses H1) H2) satisfeitas. Além disso que a0=aδ0∈Bρ(a*) . Então a iteração (3.1), com qualquer escolha de ωδk dada por (3.2) satisfaz:
-
[Convergência] limk→∞ak=a* quando δ=0 .
-
[Estabilidade] Seja δ>0 . Assuma que a iteração (3.1) é parada de acordo com o princípio da discrepância, isto é, no primeiro kδ* tal que ||uδ-F(aδk*)||≤τδ<||uδ-F(aδk)|| para algum τ>1 Entao a iteraçao (3.1) produz uma sequência (aδk*) que é estavel e convergente, com relação ao nível de ruídos δ, quando δ→0 .
Proof. A demonstração segue da teoria geral de convergência de métodos iterativos em [10, Cap 2], haja visto que o operador F(a) satisfaz os resultados de continuidade (Teorema 2.1), da Frechet diferenciabilidade (Proposição 2.1) e da condição do cone tangente (Teorema 2.2), respectivamente. □
4 RESULTADOS NUMÉRICOS
Nesta seção, apresentamos algumas simulações numéricas para a identificação do coeficiente a(x) utilizando os métodos iterativos discutidos na Seção 3.
Neste trabalho, os dados são obtidos de forma sintética, a partir da solução analítica u=u(a*) para (1.2), para o coeficiente a ∗(x) (solução do problema inverso) conhecido, para x∈[0,1]. Os dados com ruídos u δ são obtidos a partir da solução u em (1.2) com um ruído de magnitude δ, obtido aditivamente, gerado por uma variável aleatória com distribuição uniforme Γ, tomando valores em [−1, 1]. Em outras palavras, consideramos que os dados são tal que uδ=u+δΓ, satisfazendo ‖u-uδ‖≤δ. Os valores de δ serão especificados em cada experimento numérico.
A Tabela 1 descreve, esquematicamente, os passos a fim de obtermos o coeficiente de forma iterada4.
Em todos os exemplos numéricos que seguem, as EDO’s a serem resolvidas nos passos do algoritmo descrito na Tabela 1 são obtidas pelo método de diferenças finitas, para x∈[0,1], nos pontos xi=i/N, onde i = 0, 1, 2, . . . , N, para N = 50 pontos. O custo computacional está associado, basicamente, a solução de duas EDO’s de segunda ordem nos Passos 3 e 5 do algoritmo apresentado na Tabela 1. Isto é, uma EDO, para determinar u(a k ), solução do problema direto (assumindo que M é conhecida) e, em seguida, outra EDO para determinar F'(ak)rk (direção do passo no método iterativo).
Apresentaremos dois exemplos numéricos, o primeiro para o coeficiente a*(x)=1 e o outro para a*(x)=12-x5. Para cada exemplo iremos mostrar as iterações dos métodos para dois níveis de ruıdos diferentes, obtendo os coeficientes aproximados em cada um deles. Em todos os exemplos apresentados abaixo, o chute inicial foi a 0(x) = 0, 75, para x∈[0,1]. Esta escolha para a 0 foi tomada para os experimentos numéricos por representar uma escolha neutra, em termos de informações a-priori de ambos os experimentos numéricos aprestados. Como o problema de identificação do coeficiente de rigidez é não-linear, outras escolhas para a 0, satisfazendo as hipóteses do Teorema 3.3, podem, eventualmente, levar a resultados qualitativos e quantitativos distintos. Mas esse não é o foco de nossa abordagem.
Exemplo 1: Neste exemplo, supomos que estamos aplicando uma força f(x)≡1 em uma viga que satisfaça a equação de Euler-Bernoulli (1.2), com condições de contorno (1.3), com coeficiente a*(x)=1. A solução de (1.2) para este caso é u(x)=x424-x36+x24.
As Figuras 1-(a)-(f) mostram os coeficientes recuperados aδk* pelos métodos iterativos apresentados na Seção 3, implementados de acordo com o algoritmo apresentado na Tabela 1, para diferentes níveis de ruídos δ. A Figura 1-(a)-(b) mostram o coeficiente identificado aδk* para o nível de ruídos δ = 0,01. A iteração de Landweber e de steepest descent satisfizeram o princípio da discrepância para k ∗ = 118 e k ∗ = 85, espectivamente, determinando a regra de parada de cada método. Na Figura 1-(c)-(d) apresentamos a iteração aδk*, quando diminuímos o nível de ruído para δ = 0,001. Para este caso, os critérios de parada pelo principio de discrepância foram atingidos em k ∗ = 651 para a iteração de Landweber e em k ∗ = 132 para a iteração de steepest descent. Por fim, a Figura 1-(e)-(f) apresentamos a iteração aδk*, para o nível de ruído para δ = 0,0001. Para este caso, os critérios de parada pelo principio de discrepância foram atingidos em k ∗ = 905 para a iteração de Landweber e em k ∗ = 161 para a iteração de steepest descent.
É possível verificar nos exemplos numéricos apresentados na Figura 1 que o critério de parada k*=k*(δ), dado pelo princípio da discrepância, é um função monótona do nível de ruídos δ, como enunciado no Teorema 3.3-ii).
Na Figura 2, apresentamos o comportamento do resíduo para a iteração de Landweber (Figura 2-a)) e para a iteração de steepest-descent (Figura 2-b)), respectivamente, para a recuperação do coeficiente a do Exemplo 1, com nível de ruídos δ = 0,001, até o critério de parada dado pelo princípio da discrepância (como no Teorema 3.3) ser atingido.
Exemplo 2: Para este exemplo, suponha que uma força f(x)≡1 é aplicada em uma viga satisfazendo a equação de Euler-Bernoulli (1.2), com condições de contorno (1.3), cujo coeficiente de rigidez é dado por a*(x)=12-x5. A solução de (1.2)-(1.3) para este caso é u(x)=-x9144+x856-x784+x412-x33+x22.
Nas Figuras 3-(a)-(d) mostram os coeficientes recuperados aδk* pelos métodos iterativos apresentados na Seção 3, implementados de acordo com o algoritmo apresentado na Tabela 1, para diferentes níveis de ruídos δ, para o Exemplo 2. As Figuras 3-(a)-(b) mostram o coeficiente identificado aδk* para o nível de ruídos δ = 0,0001. A iteração de Landweber e de steepest descent satisfizeram o princípio da discrepância para k ∗ = 6075 e k ∗ = 173, respectivamente, determinando a regra de parada de cada método. Na Figura 3-(c)-(d) apresentamos a iteração aδk*, quando aumentamos o nível de ruído para δ = 0,001. Para este caso, os critérios de parada pelo principio de discrepância forma atingidos em k ∗ = 5601 para a iteração de Landweber e em k ∗ = 152 para a iteração de steeped descent.
Podemos perceber pelos testes numéricos apresentados que os métodos iterativos recuperam o coeficiente de rigidez aδk* de forma estável e convergente comprovando o que discutimos na Seção 3. Além disso, os resultados obtidos nos exemplos quando comparados os métodos, em cada nível de ruído, são muito semelhantes. Contudo, a disparidade está no número de iterações calculadas até que cada método atinja o critério de parada. Isso nos mostra que para o problema de identificação do coeficiente estudado a iteração de steepest descent converge mais rápido que a iteração dada por Landweber. De fato, esse resultado já era esperado pois, enquanto que em Landweber a escolha de ω k = 1 é fixa, na iteração de steepest descent a escolha de ωδk depende de cada nova iteração. Esta escolha adaptativa de ω k na iteração de steepest descent favorece o algoritmo ainda mais na recuperação de um coeficiente menos trivial, como o apresentado no Exemplo 2. Ainda, o fato do critério de parada nos resultados numéricos serem uma função monótona do nível de ruídos δ confirmam os resultado teórico previstos no Teorema 3.3.
Os resultados de reconstrução apresentados nesta contribuição são comparáveis, em termos qualitativos, aos obtidos em 14 usando regularização por molificação e em 15 usando o método de regularização de Tikhonov.
5 CONCLUSÕES
Neste trabalho provamos condições suficientes com relação ao operador parâmetro-para-medidas F(a), num problema de vigas modelado pela equação de Euler-Bernoulli, de forma a garantir que métodos iterativos de regularização, como os de Landweber e steepest descent forneçam soluções estáveis e convergentes para a reconstrução do coeficiente de rigidez a(x), com relação ao nível de ruídos nos dados das deflexões da viga u δ . Tal abordagem teórica se mostra necessária, haja visto que o problema de identificação do coeficiente de rigidez a(x) é mal posto no sentido de Hadamard. De fato, este problema é instável com relação as medidas, uma vez que, qualquer inversa para o operador F(a) (mesmo que linearizado) deve ser um operador ilimitado, pois F(a) é compacto (veja Teorema 2.1). Apresentamos ainda uma série de exemplos numéricos que suportam os resultados teóricos de convergência e estabilidade das soluções aproximadas para o problema de identificação do coeficiente de rigidez, através dos métodos iterativos de Landweber e steepest descent. De fato, os resultados numéricos mostraram que tais métodos recuperam o coeficiente aproximado de forma estável com diferentes níveis de ruídos nas medidas. Ainda, em termos qualitativos, é possível perceber que as iterações de Landweber e steepest descent proporcionam resultados bem semelhantes. No entanto, estes diferem sensivelmente quanto ao número de iterações necessárias até que seja atingido o critério de parada. Haja visto que, em cada iteração são necessárias a solução de duas EDO’s, a iteração de steepest descent é muito mais rápida e acumula menos erros (devido a discretização) que a iteração de Landweber.
Outros métodos iterativos de regularização, em especial, métodos de segunda ordem como métodos do tipo-Newton 10, serão fruto de investigações futuras.
REFERÊNCIAS
- 1 J. Baumeister. “Stable Solution of Inverse Problems”. Vieweg, 1 ed. (1987).
- 2 H. Brezis. “Functional analysis, Sobolev spaces and partial differential equations”. Springer Science & Business Media (2010).
- 3 S. Caddemi & A. Morassi. Multi-cracked Euler-Bernoulli beams: Mathematical modeling and exact solutions. International Journal of Solids and Structures, 50(6) (2013), 944-956.
- 4 J. GERE. “Mecânica dos Materiais”. Editora Cengage (2009).
- 5 F. Ghrib, L. Li & P. Wilbur. Damage identification of euler-bernoulli beams using static responses. Journal of engineering mechanics, 138(5) (2011), 405-415.
- 6 W. Groetch. “Inverse Problems in the Mathematical Sciences”. Vieweg, 1 ed. (1993).
- 7 J. Hadamard. “Lectures on the Cauchy Problem in Linear Partial Differencial Equations”. Yale University Press, 1 ed. (1923).
- 8 A. Hasanov & O. Baysal. Identification of an unknown spatial load distribution in a vibrating cantilevered beam from final overdetermination. J. Inverse Ill-Posed Probl., 23(1) (2015), 85-102.
- 9 R. Hibbeler. “Resistência dos Materiais”. Pearson Prentice Hall, 7 ed. (2006).
- 10 B. Kaltenbacher, A. Neubauer & O. Scherzer. “Iterative regularization methods for nonlinear ill-posed problems”, volume 6. Walter de Gruyter (2008).
- 11 A. Kawano. Uniqueness in the determination of unknown coefficients of an Euler-Bernoulli beam equation with observation in an arbitrary small interval of time. Journal of Mathematical Analysis and Applications, 452(1) (2017), 351-360.
- 12 A. Kirsch. “An introduction to the mathematical theory of inverse problems”, volume 120. Springer Science & Business Media (2011).
- 13 E. Kreyszig. “Introductory Functional Analysis with Applications”. Canadá, 1 ed. (1978).
- 14 A. Lerma & D. Hinestroza. Coefficient identification in the Euler-Bernoulli equation using regularization methods. Applied Mathematical Modelling, 41 (2017), 223-235.
- 15 D. Lesnic, L. Elliott & D. Ingham. Analysis of coefficient identification problems associated to the inverse Euler-Bernoulli beam theory. IMA Journal of Applied Mathematics, 62(2) (1999), 101-116.
- 16 G. Liu & S. Chen. A novel technique for inverse identification of distributed stiffness factor in structures. Journal of Sound and Vibration, 254(5) (2002), 823-835.
- 17 E. Luchinetti & E. Stüssi. Measuring the flexural rigidity in non-uniform beams using an inverse problem approach. Inverse problems, 18(3) (2002), 837.
- 18 T. Marinov, R. MAarinova & A. Vatsala. Coefficient Identification in Euler-Bernoulli Equation from Over-posed Data. Neural, Parallel, and Scientific Computations, 23 (2015), 193-218.
- 19 T. Marinov & A. Vatsala . Inverse problem for coefficient identification in the Euler-Bernoulli equation. Computers & Mathematics with Applications, 56(2) (2008), 400-410.
- 20 E.F. Medeiros. “Identificação do coeficiente de rigidez no modelo de EULER-BERNOULLI para vigas”. Master’s thesis (2019).
- 21 M. Protter & H. Weinberger. “Maximum principles in differential equations”. Springer Science & Business Media (2012).
-
1
Vale a pena ressaltar que a equivalência entre os problemas apresentados da forma (1.1) e (1.2) são uma consequência da comparação dos mesmos e da unicidade de soluções de (1.1). No entanto, para o sistema (1.2) com certas condições de contorno, o número de reações excede o número de pontos de equilíbrio independentes, tornando a viga estaticamente indeterminada 9. Para evitar esse percalço, nesta abordagem, consideraremos as condições de contorno (1.3, descritas abaixo.
-
2
Que somam 70 no total, se considerarmos as condições de contorno que são invariantes pela mudança de variável x→(L−x) .
-
3
Vale ressaltar aqui que M(x) solução da primeira equação em (1.2) é independente de a(x).
-
4
O algoritmo apresentado na Tabela 1 foi implementado, passo a passo pelos autores, utilizado a linguagem Phyton.