21
Ana Sabino, João Tiago e Paula Antunes ANÁLISE DINÂMICA DIRECTA DO MOVIMENTO DA MARCHA Ana Sabino*, João Tiago** e Paula Antunes*** *Mestrado Integrado em Engenharia Biomédica, 4º ano Instituto Superior Técnico Av. Rovisco Pais, 1049-001 Lisboa e-mail: [email protected] ** Mestrado Integrado em Engenharia Biomédica, 4º ano Instituto Superior Técnico Av. Rovisco Pais, 1049-001 Lisboa e-mail: [email protected] *** Mestrado Integrado em Engenharia Biomédica, 4º ano Instituto Superior Técnico Av. Rovisco Pais, 1049-001 Lisboa e-mail: [email protected] Palavras-chave: Biomecânica, Marcha, Dinâmica, MATLAB, Modelo dos Corpos Múltiplos, Estabilização de Baumgarte. Resumo. Neste projecto pretendeu-se analisar a marcha humana recorrendo à análise dinâmica do modelo de corpos múltiplos implementada no software MATLAB. Para tal, recorreu-se ao laboratório de biomecânica para recolher dados antropométricos e posicionais. Posteriormente os dados obtidos foram filtrados por um filtro Butterworth de 2. ª ordem e a partir daqui calcularam-se guiamentos angulares e os comprimentos dos corpos rígidos considerados. O programa implementado é capaz de demonstrar graficamente a simulação da marcha humana no plano sagital ao longo de um ciclo de marcha, sendo que determina as posições, velocidades, acelerações, forças internas e momentos por integração das equações algébrico-diferenciais. Foi ainda implementado o método de Baumgarte de forma a estabilizar as equações de constrangimentos. Conclui-se que este programa apresenta resultados viáveis, conforme o esperado. 1. INTRODUÇÃO A palavra dinâmica tem origem no termo grego dynamike que significa forte. A dinâmica é um ramo da mecânica que estuda o movimento de um sistema tendo em consideração as forças exteriores aplicadas e as suas características inerciais [1]. Segundo Vaughan (1980), há dois tipos de análise de dinâmica de corpos rígidos, a directa e a inversa. Na análise dinâmica directa, as forças envolvidas num sistema mecânico são conhecidas, e o objectivo é determinar o movimento resultante da aplicação dessas forças. Na

Dynamic analysis of human gait

Embed Size (px)

Citation preview

Page 1: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

ANÁLISE DINÂMICA DIRECTA DO MOVIMENTO DA MARCHA

Ana Sabino*, João Tiago** e Paula Antunes***

*Mestrado Integrado em Engenharia Biomédica, 4º ano Instituto Superior Técnico

Av. Rovisco Pais, 1049-001 Lisboa e-mail: [email protected]

** Mestrado Integrado em Engenharia Biomédica, 4º ano

Instituto Superior Técnico Av. Rovisco Pais, 1049-001 Lisboa

e-mail: [email protected]

*** Mestrado Integrado em Engenharia Biomédica, 4º ano Instituto Superior Técnico

Av. Rovisco Pais, 1049-001 Lisboa e-mail: [email protected]

Palavras-chave: Biomecânica, Marcha, Dinâmica, MATLAB, Modelo dos Corpos Múltiplos, Estabilização de Baumgarte.

Resumo. Neste projecto pretendeu-se analisar a marcha humana recorrendo à análise dinâmica do modelo de corpos múltiplos implementada no software MATLAB. Para tal, recorreu-se ao laboratório de biomecânica para recolher dados antropométricos e posicionais. Posteriormente os dados obtidos foram filtrados por um filtro Butterworth de 2. ª ordem e a partir daqui calcularam-se guiamentos angulares e os comprimentos dos corpos rígidos considerados. O programa implementado é capaz de demonstrar graficamente a simulação da marcha humana no plano sagital ao longo de um ciclo de marcha, sendo que determina as posições, velocidades, acelerações, forças internas e momentos por integração das equações algébrico-diferenciais. Foi ainda implementado o método de Baumgarte de forma a estabilizar as equações de constrangimentos. Conclui-se que este programa apresenta resultados viáveis, conforme o esperado.

1. INTRODUÇÃO

A palavra dinâmica tem origem no termo grego dynamike que significa forte. A dinâmica é um ramo da mecânica que estuda o movimento de um sistema tendo em consideração as forças exteriores aplicadas e as suas características inerciais [1].

Segundo Vaughan (1980), há dois tipos de análise de dinâmica de corpos rígidos, a directa e a inversa. Na análise dinâmica directa, as forças envolvidas num sistema mecânico são conhecidas, e o objectivo é determinar o movimento resultante da aplicação dessas forças. Na

Page 2: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

análise dinâmica inversa, as variáveis cinemáticas do movimento são completamente conhecidas, sendo o objectivo encontrar as forças na origem dos movimentos [2].

Neste projecto implementou-se um programa de análise dinâmica directa para determinar as forças internas e as reacções articulares num modelo biomecânico bidimensional durante a marcha.

Num problema de carácter dinâmico tem-se o número de variáveis dependentes maior que o número de equações de constrangimento que guiam o sistema. Para se estabelecer o equilíbrio cria-se um sistema de equações do movimento. Estas equações algébrico-diferenciais são complexas e, por isso, procede-se à sua integração pelo método directo. Para estabilizar as equações de constrangimento recorre-se ao método de Baumgarte.

A análise dinâmica apresenta várias aplicações na biomecânica para além do contexto do nosso trabalho, por exemplo, uma delas é no estudo da mudança do padrão da marcha em pessoas após a amputação total ou parcial de um membro inferior [3], outra é na afectação da marcha de pessoas portadoras de paralisia cerebral [4]. Fora do contexto da Biomecânica este tipo de análise também é aplicável, como por exemplo, em sistemas de tubagens [5]. 2. NOVA ABORDAGEM DE ANÁLISE EM BIOMECÂNICA

A análise dinâmica de sistemas multicorpo é um método de análise que tem vindo a sofrer desenvolvimentos. Neste sentido têm surgido estudos com o intuito de optimizar este método, um dos quais é a utilização de modelos multicorpo fléxiveis aliada a métodos de elementos finitos. Nestas simulações são geradas equações que modelam o comportamento dinâmico do sistema, cargas e reacções. A flexibilidade deve ser incluída nos modelos computacionais de modo a que haja uma maior aproximação do modelo ao corpo humano, garantindo maior precisão dos resultados. Para melhorar este modelo propõe-se a inclusão de uma estrutura flexível criada num software de análise de sistemas multicorpos obtidos por um modelo de elementos finitos do corpo. Utiliza-se por exemplo um software DADS® para modelação de sistemas multicorpos e o ANSYS® para elementos finitos. Os resultados de estudos recentes demonstram a importância de se avaliar os efeitos da flexibilidade estrutural nos esforços desenvolvidos durante a marcha [6,7]. 3. RECOLHA DE DADOS

A recolha de dados foi realizada no laboratório de mecânica equipado com uma plataforma

de pressão e 4 câmaras digitais Qualisys-ProREFLEX. Os trâmites iniciais passaram por recolher dados do paciente, tais como, o sexo, a idade, peso e altura. Colocaram-se nos pacientes 17 marcadores reflectores em 17 locais específicos do corpo (em cada metatarso, em cada tornozelo, em cada calcanhar, na lateral de cada joelho, na lateral de cada anca, em cada pulso, na lateral de cada cotovelo, em cada ombro e na zona superior da cabeça). Para minimizar o erro, os marcadores foram colocados em pontos cuja variação relativa da posição fosse mínima. Em seguida fizeram-se as medições dos comprimentos de cada segmento,

Page 3: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

medindo-se a partir de cada ponto médio de um marcador para outro. Assim, inicou-se a actividade experimental propriamente dita. As 4 câmaras de vídeo

digitais Qualisys-ProREFLEX forneceram dados da trajectória dos marcadores reflectores que após tratados no computador podem ser reproduzidos tridimensionalmente.

É do notar que se devem recalibrar as câmaras, caso seja necessário.

4. DESCRIÇÃO DA METODOLOGIA DINÂMICA DE SISTEMAS MULTICORPO

4.1. Princípio das Potências Virtuais

Utilizando este princípio é possível formular as equações do movimento de um sistema de

coordenadas naturais n dependentes. Estabelece-se que P* representa a potência virtual de valor nulo:

* * * *

1

0 0cn

Ti i

i

P f q P=

= = ⇔ = =∑ q f&& (1)

A incógnota *q& representa o vector de velocidades virtuais dependentes. Já f é o vector de

forças virtuais caracterizado por um conjunto de velocidades imaginárias consistentes com a forma homogénea das equações de velocidade dos constrangimentos cinemáticos:

* =qΦ q 0 (2)

O vector das forças virtuais caracteriza-se pelo conjunto de todas as forças que produzem potência virtual, ou seja, é constituído por todas as forças exteriores:

= −f Mq g&& (3)

onde M é a matriz de massa do sistema, q&& o vector das acelerações naturais e g o vector das forças exteriores.

Figura 1. Esquema de montagem durante a aquisição de dados no laboratorio.

Page 4: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

O vector das forças virtuais, f , não é constituído pelas forças internas, uma vez que estas se anulam (são pares acção-reacção) e portanto, não contribuem como potência virtual. Assim, as forças internas obtêm-se pelo método dos multiplicadores de Lagrange:

Ti = qg Φ λ (4)

onde ig é o vector das forças internas do sistema, qΦ a matriz jacobiana dos constrangimentos e λ o vector dos multiplicadores de Lagrange, que representam a magnitude das forças internas.

Considerando (3) e (4), a potência virtual pode escrever-se:

( )* * 0T TP = − + =qq Mq g Φ λ& && (5)

4.2. Equações do movimento

Da potência virtual tiram-se as equações de movimento para sistemas de corpos múltiplos constrangidos:

T+ =qMq Φ λ g&& (6)

Adicionam-se mais nh equações, resultando no seguinte sistema de equações:

T + =

=

q

q

Mq Φ λ g

Φ q γ

&&

&& (7)

Convertendo a equação (7) na forma matricial obtém-se:

T

=

q

q

M Φ q g

Φ 0 λ γ

&&

(8)

sendo q&& e λ incógnitas. O passo seguinte, neste projecto, é aplicar o método da integração directa para resolver (8).

4.3. Dinâmica directa

Esta análise consiste na determinação do movimento de um sistema multicorpo que ocorre devido à aplicação de forças exteriores. A partir destas, é então possível calcular a resposta dinâmica do modelo biomecânico ( q , q& , q&& ), bem como as forças internas λ .

O esquema seguinte representa a resolução do Problema de Valor Inicial:

Page 5: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

4.4. Estabilizaçao de Baumgarte

Através deste método de estabilização converte-se um sistema de equações algébrico-diferenciais num sistema de equações diferenciais ordinárias (ODE). Assim garante-se assim uma melhor aproximação da solução analítica à real (convergência do método).

Equações do movimento resultam em equações diferenciável instável, por isso utiliza-se de forma a obter uma equação diferenciável mais estável usa-se a equação:

22α β+ + =Φ Φ Φ 0&& & (9)

Sabendo a primeira e segunda derivadas das equações dos constrangimentos o sistema de equações (8) fica então estabilizado:

( ) 22

T

α β

= − − −

q

qq

gM Φ qγ Φ q ν ΦΦ 0 λ

&&

& (10)

Onde α e β são parâmetros de estabilização. .

4.5. Matriz de massa

A matriz de massa M é generalizada para corpos planares desde que o elemento tenha dois pontos a partir dos quais ela possa ser formulada. Esta matriz simétrica positiva e constante no tempo contém informações de todos os momentos de inércia e massas para cada segmento corporal considerado.

Para corpos planares, é possível definir coordenadas naturais globais (x,y) e locais ( ,x y ) de um corpo rígido definido por i e j:

Determinar q&& e λ através das

Equações do Movimento 0

00

0 t t=

=

qY

q&

t

t

=

qY

q

&&

&&

t t

tdt

+∆= ∫Y Y&

Nota: ODE 45 t t

t t

+∆

+∆

=

qY

q& t ← t+∆t

Figura 2. Representação do Problema do Valor Inicial.

Page 6: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

Figura 3. Representação do referencial local e global.

A matriz de massa M é então dada por:

T dρΩ

= Ω∫M C C (11)

sendo que ρ representa a densidade mássica. Obtém-se então a matriz mássica:

G G G2 2 2

G G G2 2

G G2 2

G G2 2

20

20

0

0

i i i

ij ij ij ij ij ij

i i

ij ij ij ij ij

i i

ij ij ij ij

i i

ij ij ij ij

mx I mx I my Im

L L L L L L

mx I my mx Im

L L L L L

mx I my IL L L L

my mx I IL L L L

− + − − −

− + − =

− −

M (12)

Os coeficientes c1 e c2 de (Error! Bookmark not defined.) são representados por:

1

2

1 010 1ij

c x

c yL

=

(13)

Repare-se que Lij dá-nos o comprimento do segmento selecionado, a massa do elemento é

representada por m, Ii é o momento polar de inércia em relação ao ponto i e ( )G G,x y são as coordenadas locais do centro de massa.

CM

ξ

x

y

j

i

P

r rj

ri

rg

η

Page 7: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

4.6. Forças de reacção

As forças internas são obtidas através da matriz de constrangimentos e dos multiplicadores de Lagrange. Pode-se representar o vector das forças internas do sistema por:

J

PI PE J PI PI PE PE JJxT T T

i i i iy

= + + + = + + +

q q q

λg g g g Φ λ Φ λ Φ

λL L (14)

Relativamente aos constrangimentos de junta, tem-se:

J

J JxJ JJ J

J

1 0

1 1

0 1

x x

x x y yTi

y y xy x

y y

R

R

R

R

λλ λ λλλ λ λλ

λ

− − − − = = = − = =

q

Ig Φ

I (15)

Quanto aos constrangimentos de produto interno:

( )( )( )( )( )( )( )( )

klPI PI PI PI PIr

-

l k

l k

l kkl

l kT

i j iij

j iij

j i

j i

x xy y

x x

y y

x x

y y

x x

y y

λ λ λ

− − − − −

− −

= = = − − − −

− −

q

r

g Φr

r

(16)

Para determinar os momentos articulares iM tem-se que:

( )( ) ( )( )PU1 sin sini ij kl ijM f L t L L tθ λ θ= = (17)

( )( ) ( )( )PU2 sin sini kl kl ijM f L t L L tθ λ θ= = (18)

a força f exercida (constrangimentos de corpo rígido) é:

PI2 ijf Lλ= (19)

5. IMPLEMENTAÇÃO COMPUTACIONAL

Para a concretização dos objectivos propostos começou-se por avaliar o desempenho do programa desenvolvido no trabalho anterior [1] aquando da introdução dos dados adquiridos

Page 8: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

em laboratório. Registaram-se e analisaram-se os ângulos, as posições e as velocidades dos pontos, do nosso modelo biomecânico, mais importantes para a análise da marcha humana (secção lkj). Após a validação do programa de cinemática optou-se por recolher os dados necessários ao início de uma análise dinâmica, também de sistema multicorpo, da marcha humana.

O programa implementado permite uma análise dinâmica directa de dados recolhidos do mesmo modo que os utilizados neste trabalho, seja qual for o movimento realizado. As posições e velocidades iniciais, de todos os pontos que constituem o nosso modelo, são o ponto de partida para o estudo das propriedades dinâmicas do sistema e análise das características das forças exteriores que lhe são aplicadas ao longo da simulação. No final são representados vários parâmetros cinemáticos e dinâmicos dos pontos relevantes do sistema de forma a facilitar a interpretação qualitativa e quantitativa dos resultados obtidos.

A componente computacional deste projecto evolui ao longo de 4 etapas: tratamento das coordenadas das posições de pontos de interesse ao longo de um ciclo da marcha, análise cinemática, tratamento de dados adquiridos, em simultâneo, num tapete de pressão, desenvolvimento do simulador de dinâmica e apresentação de resultados.

5.1 Tratamento de dados cinemáticos

Os dados obtidos no laboratório consistem em valores de coordenadas de pontos, considerados de interesse, ao longo do intervalo de tempo equivalente à duração de um ciclo da marcha, do aluno observado (João Tiago, 22 anos, 1.69 m e 60 kg). Os marcadores reflectores utilizados foram colocados ao longo do corpo do aluno segundo o esquema da figura 4.

Figura 4. Modelo multicorpo 3D adoptado pelo software de aquisição utilizado no laboratório

Page 9: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

O processamento destes dados permitiu a obtenção dos guiamentos dos vários ângulos entres os segmentos do modelo 2D (Figura 5), adoptado por nós, assim como dos comprimentos médios dos mesmos.

Antes de extrair qualquer informação filtraram-se os dados com um filtro Butterworth,

passa baixo, de segunda ordem. Existe a necessidade de definir a frequência de amostragem, que depende do procedimento de aquisição de sinal no laboratório, e a frequência de corte do filtro, que é mais elevada para pontos que possuam movimentos mais rápidos. Devido à falta de estabilidade que o filtro apresenta no início e no final de cada conjunto de dados filtrado foi necessário adicionar dois pontos em cada extremo desses mesmos conjuntos de modo a evitar essa situação. Após a filtragem estes são desprezados.

Uma vez que com os dados filtrados há que ter atenção com o tipo de dados que possuímos, 3D, e com o tipo de dados que necessitamos para conseguir o nosso modelo cinemático, em 2D. Para a determinação dos comprimentos médios de cada segmento utilizamos a definição de norma com três coordenadas, excepto para o tronco, que de todos os segmentos é o único que difere nos dois modelos (Figura 4 e 5). Como tal o comprimento deste foi definido a partir do ponto médio das ancas e do ponto médio dos ombros. Para o cálculo dos guiamentos angulares a terceira dimensão foi desprezada, e o seno e o coseno de cada ângulo, ao longo do tempo, foram conseguidos através do produto externo 2D e do produto interno 2D, respectivamente.

Após este tratamento, os dados foram interpolados e introduzidos na função cinemática já

Figura 4.1. Esquema do modelo implementado

Page 10: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

pormenorizadamente explicada em [11]. Os resultados desta simulação encontram-se discutidos mais è frente neste relatório.

Para além das matrizes de constrangimento de corpo rígido, cc, de juntas de revolução, jr, de guiamentos angulares, ang, de translação global do sistema, tr, e de rotação global do sistema, rot, descritas também em [11], são necessárias as matrizes pos, vel e acel que contêm, respectivamente, as posições, velocidades e acelerações de cada ponto do modelo ao longo do tempo de simulação. Todos estes dados são guardados no ficheiro dadoscinematica.mat que será carregado no início do simulador de dinâmica directa.

5.2. Tratamento de dados dinâmicos

Para a análise dinâmica da marcha foi necessário conhecer a massa (20), o centro de massa (21) e o momento de inércia (22) de cada segmento. A massa ponderada foi calcula através de:

(20)

em que K, tabelado em A1 em anexo, é a fracção contribuída pelo corpo em questão e M é a massa total do corpo (em quilogramas). No final da distribuição do peso pelas variadas fracções (respeitantes a cada segmento) a soma de todas as contribuições terá de ser igual ao peso total do indivíduo. O cálculo da massa da mão foi considerado em conjunto com o antebraço, utilizando o K correspondente a esta situação. Os centros de massa de cada segmento foram calculados em relação ao ponto proximal de cada segmento e a partir da seguinte expressão:

(21)

sendo R proximal o valor tabelado em A1 em anexo e Lij o comprimento do mesmo (Tabela I). Os centros de massa foram determinados em termos de coordenadas locais. Por último, os momentos de inércia I foram conseguidos através da expressão:

(22)

em que K próximal é o valor do raio de giração proximal de cada segmento (Tabela A1 em anexo). Todos os resultados deste tratamento de dados estão disponíveis na tabela II.

Segmento Corpo Rígido Comprimento (mm)

Calcanhar direito - Tornozelo direito 1 102

Calcanhar direito – Metatarso direito 2 205

Page 11: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

Metatarso direito - Tornozelo direito 3 142

Tornozelo direito - Joelho direito 4 390

Joelho direito - Anca direito 5 382

Joelho esquerdo - Anca esquerdo 6 395

Tornozelo esquerdo - Joelho esquerdo 7 410

Metatarso esquerdo - Tornozelo esquerdo 8 142

Metatarso esquerdo - Calcanhar esquerdo 9 200

Calcanhar esquerdo - Tornozelo esquerdo 10 105

Mão direita - Cotovelo direito 11 255

Cotovelo direito - Ombro direito 12 280

Cotovelo esquerdo - Ombro esquerdo 13 290

Mão esquerda - Cotovelo esquerdo 14 250

Pescoço + Cabeça 15 350

Tronco 16 460

Tabela I – Dados antropométricos do modelo biomecânico implementado

Segmento Corporal Corpo Rígido

Comprimento do segmento (mm)

Gy (mm) Massa do Segmento

(Kg)

Raio de giração (mm)

Momento de inércia

(Kg.m2)

Pé esquerdo 4_6 128,84 64,4195 0,87 61,198525 0,015889

Pé direito 1_3 127,70 63,8515 0,87 60,658925 0,015832

Coxa esquerda 13_14 380,63 164,812357 6 122,943167 0,959961

Coxa direita 11_12 387,44 167,760654 6 125,142474 0,994613

Perna esquerda 9_10 393,73 170,483791 2,79 118,905554 0,471955

Perna direita 7_8 412,88 178,778339 2,79 124,690666 0,518996

Antebraço esquerdo 16_19 269,74 117,607512 1,32 81,731826 0,104862

Antebraço direito 16_17 268,33 116,991008 1,32 81,303384 0,103765

Braço esquerdo 19_20 237,78 162,164596 1,68 76,564516 0,104833

Page 12: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

Braço Direito 17_18 240,55 164,053054 1,68 77,456134 0,107289

Cabeça + Pescoço 16_21 373,14 209,420726

6 4,86

103,6627615

0,728915

Tronco 15_16 466,85 233,425 29,82 - 0,924455

Tabela II. Valores de massa, centro de massa, raio de giração e momento de inércia para o modelo biomecânico em estudo.

Para se proceder à análise dinâmica foram fornecidos dois ficheiros com dados obtidos no laboratório: um com informação temporal do apoio das diferentes partes de cada pé na placa de pressão e outro com o valor da força de reacção normal da placa e as respectivas coordenadas (locais) de aplicação em cada pé, e ao longo da marcha.

Tal como na cinemática, procedeu-se primeiramente à filtragem, de modo idêntico ao da cinemática, mas com uma frequência de corte de 8 Hz. Posteriormente, foi preciso sincronizar as forças relativas ao pé direito pois o momento de entrada deste na placa de pressão não coincidia com o momento inicial da marcha, captada pelas câmaras. Sabendo que a marcha é um movimento ciclico, localizou-se o instante em que a força atinge o seu valor máximo, e imediatamente a seguir começa a diminuir, e admitiu-se ser esse o instante em que o pé direito começa a deixar o tapete. Esta parte final foi retirada e movida para o início do vector da força, colmatando assim a falta de valores nesta parte. Deste modo garantiu-se a presença de força de reacção no pé direito, ao longo de todo o ciclo da marcha.

O início da aplicação das forças foi situado na abcissa inicial do metatarso direito e na abcissa inicial do calcanhar esquerdo. Obtidas as referências para as posições globais de aplicação das forças, estas foram processadas para coordenadas locais, como indicado no início desta secção.

Estes dados foram armazenados na matriz fmiR (2x2) em que a primeira linha corresponde ao pé direito e a segunda ao pé esquerdo. Cada linha está estruturada do seguinte modo:

fmiR = índice do corpo rígido [110 x 4] (23)

sendo as quatro colunas da matriz [110 x 4] preenchidas pelas posições globais e pela força

normal.

5.3. Simulador de dinâmica directa

O primeiro ficheiro a ser corrido é o ensaio.m. onde é carregado o ficheiro dadoscinematicos.mat que contem a informação já explicitada anteriormente. Aqui é criada a matriz mi que possui a seguinte estrutura:

mi = [índice do corpo rígido mi xi yi Ii] (24)

sendo xi e yi as coordenadas do centro de massa, que coincidem com o ponto de

aplicação de cada força em cada segmento e Ii o momento de inércia de cada segmento, já

Page 13: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

calculado na Tabela II. Esta em conjunto com a matriz fmiR é utilizada para calcular as forças concentradas em cada um dos segmentos bem como o ponto onde cada uma é aplicada em cada um deles. Tanto as forças como o ponto de aplicação não variam ao longo do tempo e ambos são armazenados na matriz fmi.

De seguida procedemos à introdução de todos estes dados na função dinâmica. Esta é definida como:

[pos2,vel2,acel2,lam,Fh,Fu,Mg]=dinamica[qi,qpi,h,u,ipreg,ipretr,iprerot,cr,iprefcr,alfa,beta,fa,kmax] (25)

Onde os vectores qi e qpi correspondem às aproximações iniciais para os vectores de

posições e de velocidades generalizadas do sistema, h à matriz de constrangimento de segmento de corpo rígido, u à matriz de constrangimento de junta explícita entre dois pontos, ipreg, ipretr e iprerot às matrizes de células de guiamento de ângulos relativos, de translação e de rotação global, cr à matriz que, em cada linha, introduz um corpo rígido definido por um segmento do mesmo tipo de mi, iprefcr é uma matriz de células que em cada linha introduz uma força externa aplicada a um corpo rígido do tipo de fmi, alfa e beta são parâmetros do método de Baumgarte, fa é a frequência de amostragem e kmax o número de instantes da análise. Os índices pi, pj, etc. devem ser especificados em concordância com o modelo 2D e coerentes com a ordem das coordenadas introduzidas nos vectores qi e qpi .

Tal como na cinemática os dados de entrada carecem de interpolação para cálculo das velocidades e acelerações, com base nos guiamentos passados à função. Para tal construiu-se a função interpola21.m que, através das funções spline e mmspder do Matlab, cria o polinómio interpolador das variáveis de entrada e suas primeira e segunda derivadas. Para avaliar os valores saídos da interpola21 criou-se a função interpola22 que retorna os valores interpolados para o instante pretendido através da função ppval.

A interpolação é também aplicada à matriz das forças externas. Após a interpolação dos dados de entrada o programa constrói a matriz de massa global

concatenando as matrizes de massa calculadas para cada corpo rígido. De seguida são calculadas as posições e as velocidades por integração do sistema de

equações do movimento. Esta integração é feita pela função ode45 que, integrando a eqdemov reconstrói o comportamento da solução. À eqdemov são apenas passadas as formas interpoladas dos guiamentos e das forças e estas são processadas no instante específico dentro da eqdemov. A ode45 integra tudo de uma só vez, para todos os instantes do intervalo de interesse, utilizando iterativamente a equação do movimento.

As acelerações e os multiplicadores de Lagrange são obtidos recorrendo à equação do movimento eqdemov22 à qual são fornecidas as posições e as velocidades resultantes da interpolação da eqdemov21. Dentro de cada uma destas equações de movimento encontra-se um conjunto de funções iguais às usadas na análise cinemática (calcpi, calcpe, clacjr, etc.) de modo a calcular o vector de constrangimentos e a sua jacobiana.

As tensões nos segmentos, as forças de reacção nas juntas explícitas e os momentos aplicados nas juntas são obtidos através do tratamento dos multiplicadores de Lagrange.

Page 14: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

Esta função devolve quatro matrizes: uma com as posições, outra com as velocidades, a terceira com as acelerações e a última com os multiplicadores de Lagrange, ao longo do movimento. Todas elas têm uma primeira linha com tantas colunas quantos os instantes de tempo considerados, seguida de tantas linhas quanto o número de coordenadas generalizadas. Apesar de não serem devolvidas para o exterior da função, são também calculadas as forças de tensão nos diferentes segmentos rígidos Fh, as forças de reacção nas juntas Fu e os momentos de binário nas juntas Mg.

Para visualização gráfica das reacções em cada pé ao longo da marcha foi criada a função gráficos2 que consiste num ciclo principal que devolve a imagem em cada instante de tempo. Para além do plot dos pontos do vector de posições são representadas também as forças de reacção.

Após correr o ficheiro ensaio deverá correr abrir o ficheiro resultados e corrê-lo também. Os gráficos obtidos foram construídos a partir das variáveis criadas na função dinâmica.

6. APRESENTAÇÃO E DISCUSSÃO DOS RESULTADOS

De modo a analisar os resultados da simulação da marcha foram recolhidos os gráficos respeitantes à posição, velocidade, aceleração e ângulo de pontos considerados de interesse.

Figura 5. Gráfico da variação das posições a) da cabeça, b) dos ombros e c) da anca ao longo de um ciclo da

marcha e d) a variação da velocidade dos mesmos. Na figura 5 observa-se que as variações da posição da cabeça, da anca e dos ombros são

muito semelhantes. Esses pontos não influenciam a marcha, dado que a postura do executante se mantém durante um ciclo. Relativamente às velocidades, observa-se no gráfico 5 d) que estas são sobreponíveis para a cabeça, anca e ombros, pois a variação do deslocamento ao longo do tempo é uniforme para estes três pontos, sendo que, como referido anteriormente estes não são pontos com influência directa na marcha.

Page 15: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

Figura 6. Gráfico a) da variação da posição dos calcanhares esquerdo e direito; b) da variação da velocidade dos calcanhares; c) durante uma passada na posição vertical do calcanhar que ataca o solo no inicio do ciclo: d) durante a passada da velocidade vertical e horizontal do calcanhar que ataca o solo no início do ciclo.

Comparando a figura 6.a) com a figura 6.c), mais precisamente o gráfico das posições do calcanhar esquerdo (gráfico verde) com o gráfico retirado da literatura [8], verifica-se a semelhança tanto da fase como na forma. Nesta mesma curva verde do gráfico 7 a) observa-se um período entre o início e os 0.4segundos que corresponde à fase de apoio, a partir deste período a variação da posição é muito maior atingindo-se um pico que representa a elevação máxima da respectiva perna durante a marcha.

O gráfico de velocidades 6.b) obtido é muito semelhante ao da referência bibliográfica [8], sendo que na figura 6.d) na fase de apoio a velocidade é exactamente zero e no gráfico b) apesar de não ser zero, aproxima-se muito deste valor.

Figura 7. Gráfico correspondente à variação dos ângulos das juntas dos ombros esquerdo (verde) e direito (azul).

No gráfico da figura 7 verifica-se uma assimetria entre as duas curvas, isto deve-se à maior

amplitude que o braço esquerdo apresenta durante a marcha comparativamente ao braço direito. Este facto, não se deve a um erro de programação mas indica uma característica física do executante.

Page 16: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

Nos gráficos das velocidades dos joelhos na figura 8, verifica-se que a velocidade segundo x é maior na fase de balanço, como seria de esperar. A velocidade segundo z atinge o seu máximo nesta mesma fase.

No gráfico relativo à variação de posição do joelho segundo z (figura 8 a)) observam-se oscilações acentuadas quando comparadas às da cabeça, ombros e anca, já que os joelhos intervêm directamente na marcha.

A variação angular dos joelhos está de acordo com o esperado. A figura 8 b) está em fase com a curva cinzenta da figura 9 e apresenta a mesma forma.

As forças de reacção são medidas por uma plataforma de forças colocadas no percurso do

executante da marcha tal como referido anteriormente. O resultado da aquisição da placa tal encontra-se na figura 10 e da sua observação depreende-se que o resultado obtido pela simulação é satisfatório e consistente com o gráfico presente na literatura [8] (também presente na figura 10).

Figura 8. Gráfico a) da variação da posição dos joelhos esquerdo e direito; b) da variação de velocidade dos joelhos segundo o eixo zz; c) da variação angular dos joelhos; d) da variação de velocidade dos joelhos segundo o eixo xx.

Figura 9. Variação padrão do ângulo do joelho segundo a referencia [9].

Page 17: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

Figura 10. Reacções adquiridas pela placa e relação com o disponibilizado em [8].

De uma observação mais cuidada concluímos que as forças de reacção equivalentes a cada

pé se mantém durante cerca de 60% da passada e existe a sobreposição das reacções aquando se encontram os dois pés apoiados, o que acontece em cerca de 10% da passada (caso se procedesse a uma análise da corrida seria de prever que tal não acontece-se visto que durante o movimento da corrida apenas é apoiado um pé no chão de cada vez).

Na figura 11 encontra-se a comparação entre os resultados para a posição, velocidade e aceleração do calcanhar ao longo do tempo para as simulações cinemática e dinâmica. Procedeu-se a esta comparação para verificar se se garantia a fiabilidade da simulação. Como se pode observar na figura ambas as simulações representam satisfatoriamente o sucedido com o calcanhar durante um ciclo da marcha e são concordantes entre si.

Figura 11 - Posições, velocidades e acelerações do calcanhar ao longo do tempo de simulação cinemática (a azul)

e dinâmica (a vermelho).

Ao utilizar-se parâmetros de Baumgarte menos rigorosos, os resultados acabam por tender

Page 18: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

a divergir, o que é uma situação espectável, visto que a solução apresentada tem maior liberdade para estender afastamentos em relação à solução real.

Após a conclusão da validação do simulador a nível das variáveis cinemáticas procedeu-se ao cálculo dos multiplicadores de Lagrange para os vários instantes.

Figura 12. Convenção de Sinais no que diz respeito aos momentos medidos para as articulações ao membro

inferior. Existem diversas patologias ou lesões que são caracterizadas pela alteração em relação aos

valores padrão de um dos tipos de forças internas mais importantes denominadas binários nas juntas. Desta forma procedeu-se à análise dos binários na junta do tornozelo, no joelho e na anca (figura 13) visto que são as juntas com maior relevância para o estudo da marcha e patologias associadas. Quando comparados os resultados obtidos com os disponibilizados na literatura [8] verifica-se que estão bastante aproximados o que revela que o individuo executante do movimento de marcha aparenta não ter lesões nem patologias associadas a alterações significativas neste parâmetro.

Figura 13. Relação entre os resultados para o momento do tornozelo, joelho e anca obtidos através da

implementação computacional e os resultados da referência [8] para os momentos da anca joelho e tornozelo respectivamente.

É de notar que os momentos nas juntas de revolução não são uma medida estanque e que

Page 19: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

existem desvios previstos e estudados para certas situações como é o caso ilustrado na figura 14 retirado da referência [8] e que mostra os desvios nos momentos das três articulações supra-referidas para o caso da variação da velocidade de execução do movimento em estudo.

Figura 14. Desvios em relação aos momentos do tornozelo, joelho e anca causados pela execução mais rápida ou mais lenta da marcha.

Na figura 14 está a análise das forças de reacção nas juntas correspondentes a articulações

de ambos membros inferiores (nomeadamente o tornozelo, joelho e anca) em relação ao eixo vertical.

Figura 15. Gráfico representativo da força de reacção propagada através do corpo. A curva azul representa o

tornozelo, a verde o joelho e a magenta a anca.

Pode facilmente observar-se que o gráfico apresenta a mesma forma que o presente na figura 10 correspondente às curvas de reacção na placa de pressão. Isto acontece devido à propagação dessas forças para o resto do corpo. Há a dizer que esta transmissão de forças ao longo da estrutura óssea é de extrema importância para a manutenção da vitalidade óssea através de um fenómeno denominado por lei de Wolff [10]. Esta lei enuncia que para uma remodelação óssea correcta é vital que ocorra solicitação mecânica do osso.

CONCLUSÃO Tendo em conta os resultados obtidos conclui-se que o procedimento computacional foi

implementado com êxito pois, comparados com a literatura, estes são bastante aceitáveis.

Page 20: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

Quanto à análise cinemática multicorpo implementada conseguiu-se um bom estudo da marcha humana. Os valores finais conseguidos foram introduzidos na análise dinâmica directa, com sucesso. No entanto, existem valores discrepantes e sem muito significado prático, ao longo do algoritmo.

Os conceitos de análise dinâmica adquiridos nas aulas, foram todos aplicados seguindo a estruturação proposta nas mesmas.

O tempo de execução é bastante maior que o da cinemática, sendo no entanto compreensível. Um modo de melhoramento dos resultados seria, por exemplo: uma aquisição laboratorial com mais câmaras ou uma maneira mais eficiente e precisa de colocação dos marcadores.

Em relação ao trabalho anterior, esta é uma análise muito mais completa e eficaz para estudo de todo o tipo de actividades e/ou patologias relacionadas com o movimento.

REFERÊNCIAS

[1] – http://pt.wikipedia.org/wiki/Dinâmica [2] – Amádio, A.; Barbanti, V. “A biodinâmica do movimento humano e suas relações

interdisciplinares”, in Estação Liberdade, São Paulo, p. 45-70, 2000. [3] – Farahmand, F.; et al., “Kinematic and Dynamic Analysis of the Gait Cycle of Above-

Knee Amputees, in Sharif University of Technology, July 2006 [4] – Filippin, N.; Bonamigo, E., “Implicações terapéuticas da análise dinámica da marcha

na paralisia cerebral – um estudo de caso, in Revista de Fisioterapia da Universidade de Cruz Alta, Julho de 2003

[5] – Nieves, V., “Static and Dynamic Analysis of a Piping System”, University of Puerto Rico, December 2004

[6] - Schiehlen, W., “Recent developments in multibody dynamics” in Journal of Mech. Sci.

and Tech. Vol. 19, Sup. 1. Springer, January 2005. [7] – Zaeh, M.; Siedl, D., “A nex method for simulation of machining performance by

integrating finite element and multi-body simulation for machine tools”, University of München, June 2007

[8] - Winter, David A., The Biomechanics and Motor Control of Human Gait: Normal, Elderly.

[9] – Bronzino, J.; Peterson D., “Biomechanics -Principles and Applications”. CRC Press, 2008

[10] – Wolff, J., "The Law of Bone Remodeling". Berlin Heidelberg New York: Springer, 1986 (translation of the German 1892 edition)

[11] – Sabino, A.; Tiago, J.; Antunes, P.; “Análise Cinemática do Movimento Humano”, Portugal, 2009

Page 21: Dynamic analysis of human gait

Ana Sabino, João Tiago e Paula Antunes

ANEXOS

Tabela A1. Massa, localização do centro de massa e raio de giração de diversos segmentos