( PDF ) Estudo teórico e implementação computacional do método de elementos de contorno aplicado a problemas elastostáticos lineares bidimensionais

UNIVERSIDADE FEDERAL DO PARANÁ

PÓS-GRADUAÇÃO EM MÉTODOS NUMÉRICOS EM ENGENHARIA

SETOR DE TECNOLOGIA

ESTUDO TEÓRICO E IMPLEMENTAÇÃO COMPUTACIONAL

DO MÉTODO DE ELEMENTOS DE CONTORNO

APLICADO A PROBLEMAS ELASTOSTÁTICOS LINEARES

BIDIMENSIONAIS

CURITIBA

2007

ads:

Livros Grátis

http://www.livrosgratis.com.br

Milhares de livros grátis para download.

RAPHAEL FERNANDO SCUCIATO

ESTUDO TEÓRICO E IMPLEMENTAÇÃO COMPUTACIONAL

DO MÉTODO DE ELEMENTOS DE CONTORNO

APLICADO A PROBLEMAS ELASTOSTÁTICOS LINEARES

BIDIMENSIONAIS

Dissertação apresentada como requisito parcial à

obtenção do grau de Mestre em Métodos Numéricos em

Engenharia, Curso de Pós-Graduação em Métodos

Numéricos em Engenharia, Setor de Tecnologia,

Universidade Federal do Paraná.

Orientador: Prof. Dr. Luiz Alkimin de Lacerda

Co-Orientador: Prof. Dr. José Antônio Marques Carrer

CURITIBA

2007

TERMO DE APROVAÇÃO

RAPHAEL FERNANDO SCUCIATO

ESTUDO TEÓRICO E IMPLEMENTAÇÃO COMPUTACIONAL DO

MÉTODO DE ELEMENTOS DE CONTORNO APLICADO A PROBLEMAS

ELASTOSTÁTICOS LINEARES BIDIMENSIONAIS

Dissertação aprovada como requisito parcial à obtenção do grau de Mestre em Métodos Numéricos

em Engenharia, Curso de Pós-Graduação em Métodos Numéricos em Engenharia, Setor de

Tecnologia, Universidade Federal do Paraná, pela seguinte banca examinadora:

Orientador:

_________________________________________________

Prof. Luiz Alkimin de Lacerda, D. Sc.

Pós-Graduação em Métodos Numéricos em Engenharia, UFPR

Co-Orientador:

_________________________________________________

Prof. José Antonio Marques Carrer, D. Sc.

Departamento de Matemática, UFPR

_________________________________________________

Prof. Webe João Mansur, Ph. D.

Programa de Engenharia Civil, COPPE/UFRJ

_________________________________________________

Prof. Luiz Antonio Soares de Souza, D. Sc.

Departamento de Estruturas, UEL

Curitiba, 02 de março de 2007

Dedico esta obra ao meu Pai,

que me ensinou os preceitos mais

importantes da Engenharia: capricho,

precisão, criatividade, dedicação e respeito.

AGRADECIMENTOS

Primeiramente, desejo externar meu mais profundo agradecimento a Deus,

que inundou a minha vida e o meu ser com a sua maravilhosa graça; ao Senhor

Jesus Cristo, meu Salvador, que morreu por mim para me salvar dos meus próprios

pecados; ao Espírito Santo, meu Consolador, que me ensinou o significado da

verdadeira adoração; os quais deram-me força, paciência, persistência,

determinação, conforto, paz e, inquestionavelmente, inteligência e discernimento

para concretizar esta obra.

Aos meus pais, Mario e Eledir, e a minha irmã, Lisiane, pela paciência,

compreensão e incentivo durante os momentos mais difíceis e sem os quais, com

certeza, eu não teria alcançado sucesso.

A minha amada Patricia, a maior bênção que o Senhor concedeu-me depois

de minha salvação eterna, instrumento nas mãos do Pai para mostrar-me a direção

certa nos momentos mais confusos; que o nosso amor transcenda à eternidade.

Aos mestres, Prof. Alkimin e Prof. Carrer, pela orientação inestimável, pela

confiança depositada em mim e pela paciência em todos os momentos de dúvidas;

espero que este trabalho deixe-os tão orgulhosos quanto eu estou.

Ao colega e amigo Roberto, pelo companheirismo, espirituosidade e bom

humor, características que sempre me ajudaram a perceber que não estava só em

meio às adversidades.

A CAPES (Coordenação de Aperfeiçoamento de Pessoal de Nível Superior),

pela ajuda financeira que foi tremendamente importante para a realização de alguns

sonhos pessoais.

À Maristela, pela alegria, bom humor e, claro, pelo café e chá sempre

quentinhos que, inegavelmente, ajudaram a criar valiosos momentos de integração.

Aos professores e colegas do CESEC que, direta ou indiretamente,

contribuíram para a realização deste trabalho.

Ao Pr. Roberto, Pra. Gema, Pr. Rechier e Gislaine, pelas orações constantes

que certamente empreenderam em meu favor e que, em retribuição, somente posso

rogar a Deus que os recompense grandemente na Glória.

Ao colega, amigo e irmão MC Palu, pela amizade verdadeira.

"O temor do SENHOR é o princípio da sabedoria..."

Salmo 111:10a

vii

SUMÁRIO

LISTA DE FIGURAS................................................................................................xiii

LISTA DE TABELAS .............................................................................................xviii

LISTA DE GRÁFICOS ..............................................................................................xx

RESUMO................................................................................................................xxiii

ABSTRACT............................................................................................................xxiv

1 INTRODUÇÃO .........................................................................................................1

1.1 Revisão Bibliográfica .................................................................................3

1.2 Objetivos......................................................................................................7

1.3 Metodologia.................................................................................................7

1.4 Organização ................................................................................................8

2 ELASTOSTÁTICA ...................................................................................................9

2.1 Preliminares Matemáticas..........................................................................9

viii

2.1.1 Notação Indicial Cartesiana...............................................................9

2.1.2 Delta de Kronecker..........................................................................11

2.1.3 Função Delta de Dirac.....................................................................12

2.1.4 Teorema da Divergência .................................................................12

2.1.5 Série de Taylor ................................................................................13

2.1.6 Quadratura Gaussiana ....................................................................13

2.1.7 Requisitos de Continuidade.............................................................14

2.1.7.1 Condição de Continuidade de Hölder ..................................15

2.1.8 Valor Principal de Cauchy ...............................................................16

2.2 Equações Básicas da Elasticidade Linear..............................................18

2.2.1 Equações de Equilíbrio....................................................................18

2.2.2 Estado Plano de Tensões ...............................................................20

2.2.3 Tensões Principais ..........................................................................25

2.2.4 Deformações Específicas................................................................29

2.2.5 Equações de Compatibilidade.........................................................30

2.2.6 Relações Tensão-Deformação........................................................32

2.2.7 Equações Governantes ...................................................................34

2.2.8 Condições para Estado Plano de Tensões .....................................36

2.2.9 Condições para Estado Plano de Deformações..............................37

3 FORMULAÇÃO INTEGRAL DE CONTORNO ......................................................39

3.1 Teorema do Trabalho Recíproco de Betti...............................................39

3.2 Equações Integrais de Contorno.............................................................42

3.3 Soluções Fundamentais...........................................................................44

3.3.1 Problemas Tridimensionais .............................................................45

ix

3.3.2 Problemas Bidimensionais ..............................................................46

3.4 Equação Integral de Contorno para Deslocamentos no Contorno ......47

3.5 Regiões Infinitas .......................................................................................52

3.6 Discretização Numérica ...........................................................................55

3.7 Elementos Constantes .............................................................................56

3.7.1 Formação do Sistema de Equações................................................56

3.7.2 Integração .......................................................................................61

3.7.2.1 Submatrizes

cc

H ..................................................................62

3.7.2.2 Submatrizes

cc

G ..................................................................63

3.7.2.3 Submatrizes

cn

G quando nc

≠

............................................64

3.7.2.4 Submatrizes

cn

H quando nc

≠

...........................................65

3.7.3 Deslocamentos e Tensões nos Pontos Internos .............................67

3.8 Elementos Isoparamétricos .....................................................................69

3.8.1 Formação do Sistema de Equações................................................73

3.8.1.1 Elementos Lineares .............................................................73

3.8.1.2 Elementos Quadráticos........................................................77

3.8.2 Integração .......................................................................................80

3.8.2.1 Elementos Lineares – Coeficientes

α

cn

ij

Q .............................80

3.8.2.2 Elementos Lineares – Coeficientes

α

cn

ij

P .............................85

3.8.2.3 Elementos Quadráticos – Coeficientes

α

cn

ij

Q .......................87

3.8.2.4 Elementos Quadráticos – Coeficientes

α

cn

ij

P .......................96

3.8.3 Deslocamentos e Tensões nos Pontos Internos .............................98

4 CÓDIGOS COMPUTACIONAIS ............................................................................99

x

4.1 Elementos Constantes .............................................................................99

4.1.1 Organização e Funcionamento do Código ......................................99

4.1.2 Modelo de Arquivo de Entrada de Dados......................................101

4.2 Elementos Lineares ................................................................................102

4.2.1 Organização e Funcionamento do Código ....................................103

4.2.2 Modelo de Arquivo de Entrada de Dados......................................104

4.3 Elementos Quadráticos..........................................................................106

4.3.1 Organização e Funcionamento do Código ....................................107

4.3.2 Modelo de Arquivo de Entrada de Dados......................................108

4.4 Múltiplos Contornos ...............................................................................110

4.4.1 Detalhes da Implementação Computacional .................................111

4.4.2 Código em Fortran 95....................................................................112

4.4.3 Organização e Funcionamento do Código ....................................112

4.4.4 Modelo de Arquivo de Entrada de Dados......................................112

4.5 Considerações de Simetria....................................................................114

4.5.1 Problema Exemplo para Análise de Simetria ................................115

4.5.1.1 Simetria em Relação ao Eixo Horizontal............................118

4.5.1.2 Simetria em Relação ao Eixo Vertical ................................122

4.5.1.3 Simetria em Relação ao Eixo Horizontal e Vertical ............126

4.5.1.4 Deslocamentos e Tensões nos Pontos Internos................131

4.5.2 Detalhes da Implementação Computacional .................................132

4.5.3 Código em Fortran 95....................................................................132

4.5.4 Organização e Funcionamento do Código ....................................133

4.5.5 Modelo de Arquivo de Entrada de Dados......................................133

4.6 Peculiaridades da Modelagem de Problemas Utilizando MEC ...........135

xi

5 PROBLEMAS RESOLVIDOS..............................................................................137

5.1 Cavidade Pressurizada...........................................................................137

5.1.1 Considerações Gerais ...................................................................137

5.1.2 Solução Analítica...........................................................................138

5.1.3 Dados do Problema.......................................................................139

5.1.4 Malhas de Elementos Constantes.................................................139

5.1.4.1 Malha CAV12C .....................................................................140

5.1.4.2 Malha CAV24C .....................................................................140

5.1.4.3 Malha CAV48C .....................................................................141

5.1.5 Malhas de Elementos Lineares .....................................................141

5.1.5.1 Malha CAV12L .....................................................................142

5.1.5.2 Malha CAV24L .....................................................................142

5.1.5.3 Malha CAV48L .....................................................................143

5.1.6 Malhas de Elementos Quadráticos................................................143

5.1.6.1 Malha CAV6Q .......................................................................144

5.1.6.2 Malha CAV12Q .....................................................................144

5.1.6.3 Malha CAV24Q .....................................................................145

5.1.7 Resultados.....................................................................................145

5.1.7.1 Tensões Radiais ................................................................145

5.1.7.2 Tensões Angulares ............................................................146

5.1.7.3 Deslocamentos Radiais .....................................................147

5.1.7.4 Análise dos Resultados......................................................148

5.2 Cilindro Pressurizado.............................................................................159

5.2.1 Considerações Gerais ...................................................................159

5.2.2 Solução Analítica...........................................................................159

5.2.3 Dados do Problema.......................................................................160

5.2.4 Malhas de Elementos Quadráticos................................................161

xii

5.2.4.1 Malha CIL18Q .....................................................................161

5.2.4.2 Malha CIL30Q .....................................................................162

5.2.4.3 Malha CIL36Q .....................................................................163

5.2.5 Resultados.....................................................................................164

5.2.5.1 Tensões Radiais ................................................................164

5.2.5.2 Tensões Angulares ............................................................164

5.2.5.3 Deslocamentos Radiais .....................................................165

5.2.5.4 Análise dos Resultados......................................................165

5.3 Placa com Orifício...................................................................................172

5.3.1 Considerações Gerais ...................................................................172

5.3.2 Solução Analítica...........................................................................172

5.3.3 Dados do Problema.......................................................................173

5.3.4 Malhas de Elementos Quadráticos................................................174

5.3.4.1 Malha PL68Q .......................................................................175

5.3.4.2 Malha PL34QX .....................................................................176

5.3.4.3 Malha PL34QY .....................................................................177

5.3.4.4 Malha PL17Q .......................................................................178

5.3.5 Resultados.....................................................................................179

5.3.6 Análise dos Resultados .................................................................183

6 CONCLUSÃO ......................................................................................................186

6.1 Sugestões para Trabalhos Futuros.......................................................187

REFERÊNCIAS.......................................................................................................188

xiii

LISTA DE FIGURAS

FIGURA 2.1 – Componentes de tensões atuando nas faces de um volume elementar

representando um ponto no interior de um corpo tridimensional ...........................................

19

FIGURA 2.2 – Forças de superfície e forças internas atuando em um elemento localizado no

contorno de um corpo bidimensional de espessura unitária constante .................................

21

FIGURA 2.3 – Componentes ou projeções do vetor unitário normal ao contorno ...................... 22

FIGURA 2.4 – As tensões atuando em planos inclinados podem ser obtidas a partir de um

elemento triangular extraído de um elemento localizado no interior do corpo ......................

24

FIGURA 2.5 – Elemento triangular resultante e as tensões que atuam sobre o mesmo............. 24

FIGURA 2.6 – Solução geométrica para determinar os ângulos que definem os Planos

Principais nos quais atuam as Tensões Principais....................................................................

27

FIGURA 2.7 – Solução geométrica para determinar os ângulos que definem os planos nos

quais atuam as tensões tangenciais máxima e mínima ............................................................

28

FIGURA 2.8 – Componentes de deformações específicas atuando em um volume elementar

representando um ponto no interior de um corpo tridimensional: a) Volume elementar sem

deformação; b) Volume elementar deformado ...........................................................................

30

FIGURA 3.1 – Uma região arbitrária com domínio

*

Ω e contorno

*

Γ envolvendo o problema

com domínio

Ω e contorno Γ .................................................................................................... 40

xiv

FIGURA 3.2 – Processo limite quando o ponto fonte está localizado no contorno de um

problema bidimensional................................................................................................................

48

FIGURA 3.3 – Processo limite quando o ponto fonte está localizado no contorno de um

problema tridimensional ...............................................................................................................

48

FIGURA 3.4 – Região infinita envolvendo uma cavidade................................................................ 53

FIGURA 3.5 – Definição da normal e do sentido de integração: a) Problemas com domínio

infinito; b) Problemas com domínio finito...................................................................................

55

FIGURA 3.6 – Malha modelo formada por elementos de contorno constantes ........................... 57

FIGURA 3.7 – Sistema algébrico tGuH

=

destacando os valores das condições de contorno

prescritas para um problema arbitrário....................................................................................... 60

FIGURA 3.8 – Sistema algébrico FYA

=

resultante da reordenação das colunas das matrizes

H e G de modo a agrupar todas as incógnitas do problema no lado direito da igualdade e

os valores das condições de contorno prescritas no lado esquerdo da igualdade...............

61

FIGURA 3.9 – Representação geométrica dos termos necessários para o cálculo dos

coeficientes das submatrizes

cc

G ............................................................................................... 64

FIGURA 3.10 – Representação geométrica dos termos necessários para o cálculo dos

coeficientes das submatrizes

cn

G quando nc

≠

..................................................................... 65

FIGURA 3.11 – Representação geométrica dos termos necessários para o cálculo dos

coeficientes das submatrizes

cn

H

quando nc

≠

..................................................................... 66

FIGURA 3.12 – Elementos de contorno bidimensionais definidos em termos de um sistema de

coordenadas adimensional...........................................................................................................

70

FIGURA 3.13 – Malha modelo formada por elementos de contorno lineares............................... 74

FIGURA 3.14 – Malha modelo formada por elementos de contorno quadráticos........................ 78

xv

FIGURA 3.15 – Relação entre os sistemas de coordenadas adimensionais para elementos de

contorno lineares conforme a posição do ponto fonte dentro do elemento...........................

82

FIGURA 3.16 – Relação entre os sistemas de coordenadas adimensionais para elementos de

contorno quadráticos conforme a posição do ponte fonte dentro do elemento ....................

89

FIGURA 4.1 – Modelo de arquivo de entrada de dados do código MEC_ELASTIC_CONST ........ 101

FIGURA 4.2 – Ilustração da malha do modelo de arquivo de entrada de dados do código

MEC_ELASTIC_CONST ..................................................................................................................

102

FIGURA 4.3 – Modelo de arquivo de entrada de dados do código MEC_ELASTIC_LIN ............ 105

FIGURA 4.4 – Ilustração da malha do modelo de arquivo de entrada de dados do código

MEC_ELASTIC_LIN.......................................................................................................................

106

FIGURA 4.5 – Modelo de arquivo de entrada de dados do código MEC_ELASTIC_QUAD .......... 109

FIGURA 4.6 – Ilustração da malha do modelo de arquivo de entrada de dados do código

MEC_ELASTIC_QUAD.....................................................................................................................

110

FIGURA 4.7 – Exemplo de problema que apresenta múltiplos contornos.................................. 111

FIGURA 4.8 – Modelo de arquivo de entrada de dados do código MEC_ELASTIC_QUAD_MC.... 113

FIGURA 4.9 – Ilustração da malha do modelo de arquivo de entrada de dados do código

MEC_ELASTIC_QUAD_MC ..............................................................................................................

114

FIGURA 4.10 – Exemplo de problema que apresenta simetria..................................................... 115

FIGURA 4.11 – Malha modelo para análise de problemas que apresentam simetria ................ 117

FIGURA 4.12 – Malha modelo considerando apenas simetria em relação ao eixo horizontal.. 119

FIGURA 4.13 – Técnica de reflexão do ponto fonte considerando apenas simetria em relação ao

eixo horizontal..............................................................................................................................

122

xvi

FIGURA 4.14 – Malha modelo considerando apenas simetria em relação ao eixo vertical ...... 123

FIGURA 4.15 – Técnica de reflexão do ponto fonte considerando apenas simetria em relação ao

eixo vertical...................................................................................................................................

126

FIGURA 4.16 – Malha modelo considerando simetria em relação ao eixo horizontal e vertical

simultaneamente..........................................................................................................................

127

FIGURA 4.17 – Técnica de reflexão do ponto fonte considerando apenas simetria em relação ao

eixo vertical para um problema que apresenta simetria em relação ao eixo horizontal e

vertical simultaneamente ............................................................................................................

130

FIGURA 4.18 – Técnica de reflexão do ponto fonte considerando apenas simetria em relação ao

eixo horizontal e vertical para um problema que apresenta simetria em relação ao eixo

horizontal e vertical simultaneamente.......................................................................................

130

FIGURA 4.19 – Técnica de reflexão do ponto fonte considerando apenas simetria em relação ao

eixo horizontal para um problema que apresenta simetria em relação ao eixo horizontal e

vertical simultaneamente ............................................................................................................

131

FIGURA 4.20 – Modelo de arquivo de entrada de dados do código MEC_ELASTIC_QUAD_MC_SIM

........................................................................................................................................................

134

FIGURA 4.21 – Ilustração da malha do modelo de arquivo de entrada de dados do código

MEC_ELASTIC_QUAD_MC_SIM .....................................................................................................

134

FIGURA 5.1 – Cavidade em um meio infinito sob pressão uniforme........................................... 138

FIGURA 5.2 – Representação geométrica da malha CAV12C........................................................ 140

FIGURA 5.3 – Representação geométrica da malha CAV24C........................................................ 140

FIGURA 5.4 – Representação geométrica da malha CAV48C........................................................ 141

FIGURA 5.5 – Representação geométrica da malha CAV12L........................................................ 142

FIGURA 5.6 – Representação geométrica da malha CAV24L........................................................ 142

xvii

FIGURA 5.7 – Representação geométrica da malha CAV48L........................................................ 143

FIGURA 5.8 – Representação geométrica da malha CAV6Q.......................................................... 144

FIGURA 5.9 – Representação geométrica da malha CAV12Q........................................................ 144

FIGURA 5.10 – Representação geométrica da malha CAV24Q...................................................... 145

FIGURA 5.11 – Cilindro sob pressão uniforme .............................................................................. 159

FIGURA 5.12 – Representação geométrica da malha CIL18Q...................................................... 161

FIGURA 5.13 – Representação geométrica da malha CIL30Q...................................................... 162

FIGURA 5.14 – Representação geométrica da malha CIL36Q...................................................... 163

FIGURA 5.15 – Placa retangular com orifício central sob tração uniforme ................................ 172

FIGURA 5.16 – Representação geométrica da malha PL68Q........................................................ 175

FIGURA 5.17 – Representação geométrica da malha PL34QX...................................................... 176

FIGURA 5.18 – Representação geométrica da malha PL34QY...................................................... 177

FIGURA 5.19 – Representação geométrica da malha PL17Q........................................................ 178

xviii

LISTA DE TABELAS

Tabela 5.1 – Dados utilizados na análise do problema de cavidade pressurizada.................... 139

Tabela 5.2 – Tensões radiais obtidas com o código MEC_ELASTIC_CONST ............................... 145

Tabela 5.3 – Tensões radiais obtidas com o código MEC_ELASTIC_LIN.................................... 146

Tabela 5.4 – Tensões radiais obtidas com o código MEC_ELASTIC_QUAD ................................. 146

Tabela 5.5 – Tensões angulares obtidas com o código MEC_ELASTIC_CONST.......................... 146

Tabela 5.6 – Tensões angulares obtidas com o código MEC_ELASTIC_LIN .............................. 147

Tabela 5.7 – Tensões angulares obtidas com o código MEC_ELASTIC_QUAD............................ 147

Tabela 5.8 – Deslocamentos radiais obtidos com o código MEC_ELASTIC_CONST................... 147

Tabela 5.9 – Deslocamentos radiais obtidos com o código MEC_ELASTIC_LIN ....................... 148

Tabela 5.10 – Deslocamentos radiais obtidos com o código MEC_ELASTIC_QUAD ................... 148

Tabela 5.11 – Dados utilizados na análise do problema de cilindro pressurizado .................... 160

Tabela 5.12 – Tensões radiais obtidas com o código MEC_ELASTIC_QUAD_MC......................... 164

Tabela 5.13 – Tensões angulares obtidas com o código MEC_ELASTIC_QUAD_MC ................... 164

xix

Tabela 5.14 – Deslocamentos radiais obtidos com o código MEC_ELASTIC_QUAD_MC............. 165

Tabela 5.15 – Dados utilizados na análise do problema de placa com orifício.......................... 173

Tabela 5.16 – Tensões angulares obtidas com o código MEC_ELASTIC_QUAD_MC_SIM para a

malha PL68Q .................................................................................................................................

179

Tabela 5.17 – Tensões angulares obtidas com o código MEC_ELASTIC_QUAD_MC_SIM para a

malha PL34QX...............................................................................................................................

180

Tabela 5.18 – Tensões angulares obtidas com o código MEC_ELASTIC_QUAD_MC_SIM para a

malha PL34QY...............................................................................................................................

181

Tabela 5.19 – Tensões angulares obtidas com o código MEC_ELASTIC_QUAD_MC_SIM para a

malha PL17Q .................................................................................................................................

182

xx

LISTA DE GRÁFICOS

Gráfico 5.1 – Solução analítica para as tensões radiais versus raio versus malhas CAV12C,

CAV24C e CAV48C .........................................................................................................................

150

Gráfico 5.2 – Solução analítica para as tensões radiais versus raio versus malhas CAV12L,

CAV24L e CAV48L .........................................................................................................................

150

Gráfico 5.3 – Solução analítica para as tensões radiais versus raio versus malhas CAV6Q,

CAV12Q e CAV24Q .........................................................................................................................

151

Gráfico 5.4 – Erro no cálculo das tensões radiais versus raio versus malhas CAV12C, CAV24C e

CAV48C ..........................................................................................................................................

151

Gráfico 5.5 – Erro no cálculo das tensões radiais versus raio versus malhas CAV12L, CAV24L e

CAV48L ..........................................................................................................................................

152

Gráfico 5.6 – Erro no cálculo das tensões radiais versus raio versus malha CAV6Q, CAV12Q e

CAV24Q ..........................................................................................................................................

152

Gráfico 5.7 – Solução analítica para as tensões angulares versus raio versus malhas CAV12C,

CAV24C e CAV48C .........................................................................................................................

153

Gráfico 5.8 – Solução analítica para as tensões angulares versus raio versus malhas CAV12L,

CAV24L e CAV48L .........................................................................................................................

153

xxi

Gráfico 5.9 – Solução analítica para as tensões angulares versus raio versus malhas CAV6Q,

CAV12Q e CAV24Q .........................................................................................................................

154

Gráfico 5.10 – Erro no cálculo das tensões angulares versus raio versus malhas CAV12C,

CAV24C e CAV48C .........................................................................................................................

154

Gráfico 5.11 – Erro no cálculo das tensões angulares versus raio versus malhas CAV12L,

CAV24L e CAV48L .........................................................................................................................

155

Gráfico 5.12 – Erro no cálculo das tensões angulares versus raio versus malha CAV6Q, CAV12Q

e CAV24Q .......................................................................................................................................

155

Gráfico 5.13 – Solução analítica para os deslocamentos radiais versus raio versus malhas

CAV12C, CAV24C e CAV48C ..........................................................................................................

156

Gráfico 5.14 – Solução analítica para os deslocamentos radiais versus raio versus malhas

CAV12L, CAV24L e CAV48L ..........................................................................................................

156

Gráfico 5.15 – Solução analítica para os deslocamentos radiais versus raio versus malhas

CAV6Q, CAV12Q e CAV24Q ............................................................................................................

157

Gráfico 5.16 – Erro no cálculo dos deslocamentos radiais versus raio versus malhas CAV12C,

CAV24C e CAV48C .........................................................................................................................

157

Gráfico 5.17 – Erro no cálculo dos deslocamentos radiais versus raio versus malhas CAV12L,

CAV24L e CAV48L .........................................................................................................................

158

Gráfico 5.18 – Erro no cálculo dos deslocamentos radiais versus raio versus malha CAV6Q,

CAV12Q e CAV24Q .........................................................................................................................

158

Gráfico 5.19 – Solução analítica para as tensões radiais versus raio versus malha CIL18Q... 166

Gráfico 5.20 – Solução analítica para as tensões radiais versus raio versus malha CIL30Q... 166

Gráfico 5.21 – Solução analítica para as tensões radiais versus raio versus malha CIL36Q... 167

Gráfico 5.22 – Erro no cálculo das tensões radiais versus raio versus malha CIL18Q, CIL30Q e

CIL36Q ..........................................................................................................................................

167

xxii

Gráfico 5.23 – Solução analítica para as tensões angulares versus raio versus malha CIL18Q

........................................................................................................................................................

168

Gráfico 5.24 – Solução analítica para as tensões angulares versus raio versus malha CIL30Q

........................................................................................................................................................

168

Gráfico 5.25 – Solução analítica para as tensões angulares versus raio versus malha CIL36Q

........................................................................................................................................................

169

Gráfico 5.26 – Erro no cálculo das tensões angulares versus raio versus malha CIL18Q,

CIL30Q e CIL36Q .........................................................................................................................

169

Gráfico 5.27 – Solução analítica para os deslocamentos radiais versus raio versus malha

CIL18Q ..........................................................................................................................................

170

Gráfico 5.28 – Solução analítica para os deslocamentos radiais versus raio versus malha

CIL30Q ..........................................................................................................................................

170

Gráfico 5.29 – Solução analítica para os deslocamentos radiais versus raio versus malha

CIL36Q ..........................................................................................................................................

171

Gráfico 5.30 – Erro no cálculo dos deslocamentos radiais versus raio versus malha CIL18Q,

CIL30Q e CIL36Q .........................................................................................................................

171

Gráfico 5.31 – Solução analítica para as tensões angulares versus raio ................................... 183

Gráfico 5.32 – Erro no cálculo das tensões angulares versus raio versus malha PL68Q ......... 184

Gráfico 5.33 – Erro no cálculo das tensões angulares versus raio versus malha PL34QX....... 184

Gráfico 5.34 – Erro no cálculo das tensões angulares versus raio versus malha PL34QY....... 185

Gráfico 5.35 – Erro no cálculo das tensões angulares versus raio versus malha PL17Q ......... 185

xxiii

RESUMO

Este trabalho tem por objetivo, primeiramente, realizar um estudo teórico dos

aspectos matemáticos e computacionais relacionados à formulação do Método de

Elementos de Contorno aplicado a problemas elastostáticos lineares bidimensionais.

Posteriormente, são desenvolvidos códigos em linguagem de programação Fortran

95 capazes de analisar problemas de Estado Plano de Tensões e Estado Plano de

Deformações, utilizando elementos de contorno constantes, lineares e quadráticos,

múltiplos contornos e considerações de simetria.

Palavras-chave: Método de Elementos de Contorno, Elastostática, Elasticidade

Linear Bidimensional, Fortran 95.

xxiv

ABSTRACT

The first purpose of this work is to accomplish a theoretical study of the mathematical

and computational aspects related to the Boundary Element Method formulation

applied to bidimensional linear elastostatic problems. Later, computational codes are

developed in Fortran 95 language, which are capable of solving problems using

constant, linear and quadratic boundary elements, taking multiboundaries and

symmetry conditions into consideration.

Key-words: Boundary Element Method, Elastostatic, Two-dimensional Linear

Elasticity, Fortran 95.

1

1 INTRODUÇÃO

Problemas de engenharia são frequentemente descritos por leis físicas as

quais são comumente expressas por equações diferenciais parciais. Em muitos

casos, uma representação matemática alternativa e equivalente do problema é

encontrada em termos de equações integrais de contorno. Com o avanço nas

técnicas de modelagem numérica e o incremento na capacidade de processamento

dos computadores, métodos de modelagem baseados em equações integrais de

contorno podem ser agora usados para a simulação de muitos problemas práticos

de engenharia. A mais geral e efetiva técnica numérica para solução de equações

integrais de contorno é o Método de Elementos de Contorno.

1

A redução das dimensões do problema analisado constitui uma das razões

pela qual o Método de Elementos de Contorno é tão atrativo: em problemas

bidimensionais, apenas o contorno unidimensional do domínio necessita ser

discretizado em elementos; em problemas tridimensionais, apenas as superfícies do

contorno necessitam ser discretizadas. Isso significa que, quando comparada ao

Método de Elementos Finitos e a outras técnicas de análise de domínio, uma análise

1

ALIABADI, M. H. The Boundary Element Method, Vol. 2, Applications in Solids and

Structures. Chippenham: John Wiley & Sons, 2002. p. xiii.

2

de contorno resulta em uma substancial redução no esforço de modelagem

computacional. Além do mais, como os parâmetros incógnitos (tais como

deslocamentos e forças de superfície) são aproximados nos contornos discretizados,

um sistema de equações muito menor do que o gerado pelo Método de Elementos

Finitos é obtido.

2

Uma peculiaridade do Método de Elementos de Contorno é que ele provê

um modelo contínuo do domínio, uma vez que nenhuma discretização do mesmo é

requerida. Isso conduz a uma ótima resolução dos parâmetros incógnitos no domínio

do problema. As soluções nos pontos internos são calculadas após as incógnitas de

contorno terem sido calculadas, de maneira semelhante a um pós-processamento. A

densidade, distribuição e localização dos pontos internos não interferem na malha

de contorno, tampouco nos valores das incógnitas de contorno.

3

Outra característica do Método de Elementos de Contorno é a satisfação

automática de condições de contorno para domínios infinitos e semi-infinitos, o que

elimina a necessidade de discretização numérica em contornos remotos.

4

Como dito anteriormente, o Método de Elementos de Contorno constitui-se

numa técnica numérica, baseada num procedimento de discretização, para a

solução das equações oriundas das denominadas formulações integrais de

contorno. Essas formulações integrais de contorno são divididas em duas categorias

diferentes. A primeira, e talvez a mais popular, é a chamada formulação direta, na

qual as funções incógnitas presentes na formulação são variáveis físicas do

problema. Em elasticidade, por exemplo, essas funções incógnitas são os campos

de deslocamentos e forças de superfície. A outra abordagem é chamada de

formulação indireta, na qual as funções incógnitas são representadas por fontes de

densidade fictícias. Uma vez que essas fontes de densidade são calculadas, os

valores dos parâmetros físicos podem ser obtidos através de simples integrações.

5

Neste trabalho, o Método de Elementos de Contorno será desenvolvido

usando a formulação direta.

2

Ibid., p. 1.

3

Ibid., p. 2,3.

4

Ibid., p. 3.

5

Ibid., p. 15.

3

1.1 Revisão Bibliográfica

Apesar do Método de Elementos de Contorno ser uma técnica relativamente

nova em análises de problemas de engenharia, seus fundamentos podem ser

relacionados às formulações matemáticas clássicas de FREDHOLM

6

e MIKHILIN

7

,

em teoria de potencial, e BETTI

8

, SOMIGLIANA

9

e KUPRADZE

10

, em elasticidade.

O desenvolvimento das formulações no contexto de equações integrais de contorno

é devido à JASWON

11

, HESS e SMITH

12

, MASSONNET

13

, RIZZO

14

e CRUSE

15

. O

trabalho de LACHAT

16

e LACHAT e WATSON

17

são provavelmente as contribuições

iniciais mais significativas na direção de tornar o Método de Elementos de Contorno

uma técnica numérica efetiva.

18

Com relação à teoria da elasticidade, mais especificamente, os fundamentos

da equação integral de contorno podem ser relacionados às formulações

matemáticas clássicas de BETTI

19

, SOMIGLIANA

20

, MUSKHELISHVILI

21

e

6

FREDHOLM, I. Sur une classe d'equations fonctionelles. Acta Mathemactica, n. 27, p.

365-390,1903.

7

MIKHILIN, S. G. Integral Equations. London: Pergamon Press, 1957.

8

BETTI, E. Teoria dell elasticita. Il Nuovo Ciemento, p. 7-10, 1872.

9

SOMIGLIANA, C. Sopra l'equilibrio di un corpo elastico isotrope. Il Nuovo Ciemento, n.

20, p. 181-185, 1886.

10

KUPRADZE, V. D. Potential methods in the theory of elasticity. Israel Programme for

Scientific Translations, 1965.

11

JAWSON, M. A. Integral equation method in potential theory. Proceedings of the Royal

Society of London, n. 275, p. 23-32, 1963.

12

HESS, J. L., SMITH, A. M. O. Calculation of potential flows about arbitrary bodies.

Progress in Aeronautical Sciences, n. 8, 1967.

13

MASSONNET, C. E. Numerical use of integral procedures. Stress Analysis, cap. 10, p.

198-235,1965.

14

RIZZO, F. J. An integral approach to boundary-value problems of classical elastostatics.

Quaterly Journal of Applied Mathematics, n. 25, p. 83-95, 1967.

15

CRUSE, T. A. Numerical solutions in three-dimensional elastostatics. International

Journal of Solids and Structures, n. 5, p. 1259-1274, 1969.

16

LACHAT, J. C. A Further Development of the Boundary Integral Technique for

Elastostatics. PhD Thesis, University of Southampton, 1975.

17

LACHAT, J. C., WATSON, J. O. Effective numerical treatment of boundary integral

equations: A formulation for three-dimensional elastostatics. International Journal for Numerical

Methods in Engineering, n. 10, p. 991-1005, 1976.

18

ALIABADI, op. cit., p. 1.

19

BETTI, op. cit.

20

SOMIGLIANA, op. cit.

4

KUPRADZE

22

. O desenvolvimento dessas formulações, no contexto de métodos

numéricos, tomou forma em meados da década de 60 através dos trabalhos de

MASSONNET

23

, RIZZO

24

e CRUSE

25

, sendo que esse último introduziu a primeira

formulação para a elasticidade tridimensional. Progressos posteriores, durante a

década de 70, agregaram a representação paramétrica à formulação, de maneira

similar àquela usada no Método de Elementos Finitos. RICCARDELLA

26

e CRUSE

27

melhoraram a implementação numérica do Método de Elementos de Contorno por

meio da introdução de variação linear sobre os elementos. No entanto, como

mencionado anteriormente, a contribuição mais significativa para a construção do

Método de Elementos de Contorno resultou do trabalho de LACHAT

28

e LACHAT e

WATSON

29

, os quais desenvolveram uma formulação isoparamétrica para o

método.

30

Importantes aplicações da elasticidade tridimensional são as estruturas

formadas por placas e cascas, muito comuns em projetos de engenharia (por

exemplo: asas de aviões e painéis de fuselagem, tanques de armazenamento, lajes

de edifícios, entre outros). Há duas teorias amplamente utilizadas para placas e

cascas: a primeira foi desenvolvida por KIRCHHOFF

31

, e é comumente referenciada

como a "teoria clássica", e a outra foi desenvolvida por REISSNER

32

, e é conhecida

como teoria de deformação por cisalhamento.

33

21

MUSKHELISHVILI, N. I. Some Basic Problems of the Mathematical Theory of

Elasticity. Noordhoff: Leyden, 1953.

22

KUPRADZE, op. cit.

23

MASSONNET, op. cit.

24

RIZZO, op. cit.

25

CRUSE, op. cit.

26

RICCARDELLA, P. C. An Implementation of the Boundary Integral Technique for

Planar Problems of Elasticity and Elastoplasticity. PhD Thesis, Carnegie-Mellon University, 1973.

27

CRUSE, T. A., SWEDLOW, J. L. Interctive program for analysis and design problems in

advanced composite technology. Report AFML-TR-71-268, Carnegie-Mellon University, 1971.

28

LACHAT, op. cit.

29

LACHAT, J. C., WATSON, J. O., op. cit.

30

ALIABADI, op. cit., p. 15.

31

KIRCHHOFF, G. Uber das gleichgewicht und die bewegung einer elastischen scheibe. J.

Rein Angew Math., n. 40, p. 51-88, 1850.

32

REISSNER, E. On bending of elastic plates. Quartely Journal of Applied Mathematics,

n. 5, p. 55-68, 1947.

33

ALIABADI, op. cit., p. 103.

5

A aplicação do Método de Elementos de Contorno para a teoria clássica de

placas foi primeiramente proposta por JASWON, MAITI e SYMM

34

. O

desenvolvimento da formulação direta do Método de Elementos de Contorno para a

teoria de placas de Kirchhoff foi introduzido por FORBES e ROBINSON

35

, seguido

pelos trabalhos de BÉZINE

36

, STERN

37

e TOTTENHAM

38

. HARTMANN e

ZORMANTEL

39

ampliaram o desenvolvimento do método e implementaram-no para

placas com geometrias, carregamentos e apoios relativamente complexos. Seguindo

os trabalhos anteriores, outros (ABDEL-AKHER e HARTLEY

40

, KARAMI,

ZARRINCHANG e FOROUGHI

41

, TANAKA e MIYAZAKI

42

, VENTURINI e PAIVA

43

e

PAIVA e ALIABADI

44

) desenvolveram e aplicaram com maior profundidade a

formulação integral de contorno para a teoria clássica de flexão de placas.

45

Por sua vez, o desenvolvimento da formulação direta do Método de

Elementos de Contorno para a teoria de placas de Reissner foi apresentado por

VANDER WEEËN

46

. Mais tarde, KARAM e TELLES

47

estenderam a formulação para

34

JASWON, M. A., MAITI, M., SYMM, G. T. Numerical biharmonic analysis and some

applications. International Journal of Solids and Structures, n. 3, p. 309-332, 1967.

35

FORBES, D. J., ROBINSON, A. R. Numerical analysis of elastic plates and shalow shells

by an integral equation method. Structural Research Series Report, n. 345, University of Illinois,

Urbana, 1969.

36

BÉZINE, G. Boundary integral formulation for plate flexure with arbitrary boundary

conditions. Mech. Research Communications, n. 5(4), p. 197-206, 1978.

37

STERN, M. A general boundary integral formulation for plate bending problems.

International Journal of Solids and Structures, n. 15, p. 769-782, 1979.

38

TOTTENHAM, H., The boundary element method for plates and shells, in Developments

in Boundary Element Methods I, Applied Science, London, 173-205, 1979.

39

HARTMANN, F., ZORMANTEL, R. The direct boundary element method in plate bending.

International Journal for Numerical Methods in Engineering, n. 23, p. 2049-2069, 1986.

40

ABDEL-AKHER, A., HARTLEY, G. A. Evaluation of boundary integrals for plate bending.

International Journal for Numerical Methods in Engineering, n. 28, p. 75-93, 1989.

41

KARAMI, G., ZARRINCHANG, J., FOROUGHI, B. Analitycal treatment of boundary

integrals in direct boundary element analysis of plate bending problems. International Journal for

Numerical Methods in Engineering, n. 37, p. 2409-2427, 1994.

42

TANAKA, M., MIYAZAKI, K. A direct BEM for elastic plate-structures subjected to arbitrary

loadings. Boundary Elements VII, p. 4-16, 1985.

43

VENTURINI, W. S., PAIVA, J. B. Plate bending analysis by the boundary element method

considering zoned thickness domain. Software for Engineering Workstations, n. 4, p. 183-185,

1988.

44

PAIVA, J. B., ALIABADI, M. H. Boundary element analysis of zoned plates in bending.

Computacional Mechanics, n. 25, p. 560-566, 2000.

45

ALIABADI, op. cit., p. 103.

46

VANDER WEEËN, F. Application of the boundary integral equation method to Reissner’s

plate model. International Journal for Numerical Methods in Engineering, n. 18, p. 1-10, 1982.

6

englobar regiões infinitas e também reportaram que o modelo de placa de Reissner

é adequado tanto para placas delgadas quanto espessas. BARCELLOS e SILVA

48

apresentaram uma formulação similar àquela de Vander Weeën para o modelo de

placa de Mindlin. A formulação deles difere da formulação de Reissner no fator de

cisalhamento. WESTPHAL e BARCELLOS

49

investigaram a importância dos termos

negligenciados nas soluções fundamentais derivadas por Vander Weeën. Eles

concluíram que esses termos não têm efeito sobre os resultados.

50

Formulações não-lineares para o Método de Elementos de Contorno podem

ser encontradas em TELLES e KARAM

51

, para problemas elastoplásticos, e XIAO-

YAN, MAO-KUANG

e XIUXI

52

, para problemas geometricamente não-lineares. Mais

recentemente, WEN, ALIABADI e YOUNG

53

apresentaram uma formulação de

elementos de contorno para a teoria de placas de Kirchhoff no domínio da

freqüência. Em ALIABADI

54

, avanços recentes na análise de problemas de flexão de

placas com elementos de contorno podem ser encontrados.

55

Trabalhos mais recentes têm-se concentrado em desenvolver métodos mais

eficientes para o tratamento matemático das integrais singulares e hiper-singulares

presentes nas formulações integrais de contorno. Alguns exemplos de tais trabalhos

seriam GUIGGIANI e CASALINI

56

, GUIGGIANI

57

, HUANG e CRUSE

58

e FLEURY,

MANSUR e AZEVEDO

59

.

47

KARAM, V. J., TELLES, J. C. F. On boundary elements for Reissner’s plate theory.

Engineering Analysis, n. 5, p. 21-27, 1988.

48

BARCELLOS, C. A., SILVA, L. H. M. A boundary element formulation for the Midlin’s plate

model. Boundary Element Technology, p. 123-130, 1989.

49

WESTPHAL, T. Jr., BARCELLOS, C. A.. Application of the boundary element method to

Reissner’s and Midlin’s plate model. Boundary Elements XII, n. 1, p. 467-477, 1990.

50

ALIABADI, op. cit., p. 103,104.

51

TELLES, J. F. C., KARAM, V. J. Nonlinear material analysis of Reissner's plates. Plate

Bending with Boundary Elements, p. 127-164, 1998.

52

XIAO-YAN, L., MAO-KUANG, H., XIUXI, W. Geometrically nonlinear analysis of Reissner

type plate by the boundary element method. Computer and Structures, n. 37, p. 911-916, 1990.

53

WEN, P. H., ALIABADI, M. H., YOUNG, A. Boundary element method for dynamic plate

bending problems. International Journal of Solids and Structures, n. 37, p. 5177-5188, 2000.

54

ALIABADI, M. H. Plate Bending Analysis with Boundary Elements. Southampton:

Computational Mechanics Publications, 1998.

55

ALIABADI, op. cit., p. 104.

56

GUIGGIANI, M., CASALINI, P. Direct computation of Cauchy principal value integrals in

advanced boundary elements. International Journal for Numerical Methods in Engineering, vol.

24, p. 1711-1720, 1987.

7

1.2 Objetivos

Os principais objetivos deste trabalho são:

a) revisão da formulação matemática da teoria da elasticidade linear;

b) apresentação da formulação integral de contorno para a elastostática

linear bidimensional;

c) discussão dos aspectos computacionais relacionados ao processo de

implementação numérica do Método de Elementos de Contorno aplicado

à elastostática linear bidimensional;

d) criação de códigos em linguagem de programação Fortran 95 capazes de

analisar problemas elastostáticos lineares bidimensionais de Estado

Plano de Tensões e Estado Plano de Deformações.

1.3 Metodologia

Com vistas a alcançar os objetivos estabelecidos anteriormente, a seguinte

metodologia é adotada neste trabalho: em primeiro lugar, realiza-se um estudo

teórico-matemático dos tópicos afins; em segundo lugar, descrevem-se as

implementações computacionais (características, organização e funcionamento) e,

finalmente, faz-se a validação dos códigos desenvolvidos através da análise de

problemas que possuem solução analítica conhecida.

57

GUIGGIANI, M. The evaluation of Cauchy principal value integrals in the Boundary

Element Method – A review. Mathematical Computational Modelling, vol. 15, n. 3-5, p. 175-184,

1991.

58

HUANG, Q., CRUSE, T. A. Some notes on singular integral techniques in boundary

element analysis. International Journal for Numerical Methods in Engineering, vol. 36, p. 2643-

2659, 1993.

59

FLEURY, P., MANSUR, W. J., AZEVEDO, J. P. S. Um novo enfoque para formulação

hipersingular do MEC para a solução de problemas potenciais bidimensionais. Anais do XVIII

CILAMCE, vol. 1, p. 261-268, 1997.

8

1.4 Organização

O restante deste trabalho está organizado do seguinte modo: a primeira

parte do Capítulo 2 apresenta, inicialmente, uma revisão de alguns conceitos

matemáticos básicos utilizados ao longo dos desenvolvimentos posteriores. A

segunda parte traz uma revisão da formulação matemática da teoria da elasticidade

linear, deduzindo as principais relações que serão necessárias para a apresentação

do Método de Elementos de Contorno e suas aplicações a problemas de Estado

Plano de Tensões e Estado Plano de Deformações.

O Capítulo 3, por sua vez, apresenta a formulação do Método de Elementos

de Contorno aplicado a problemas elastostáticos lineares bidimensionais. Os

aspectos computacionais relacionados à implementação numérica do método são

discutidos. São considerados os casos de elementos de contorno constantes e

isoparamétricos (lineares e quadráticos).

A descrição da organização e funcionamento dos códigos elaborados com

base nas formulações desenvolvidas é apresentada no Capítulo 4. Particularidades

da elaboração de códigos que utilizam malhas com múltiplos contornos e

considerações de simetria, são discutidas nas últimas seções desse capítulo.

A validação dos códigos desenvolvidos, apresentada no Capítulo 5, é feita

por meio da análise de problemas que possuem solução analítica conhecida. Cada

seção desse capítulo trata de um problema específico, estando organizadas da

seguinte maneira: características e dados gerais do problema, solução analítica,

malhas de elementos de contorno utilizadas nas análises numéricas, tabelas

comparativas de resultados e discussões pertinentes.

Conclusões e considerações finais acerca dos resultados globais e

contribuições deste trabalho são abordadas no Capítulo 6.

9

2 ELASTOSTÁTICA

2.1 Preliminares Matemáticas

Nesta seção, alguns conceitos matemáticos básicos, que serão utilizados ao

longo do desenvolvimento de todo este trabalho, serão brevemente revisados.

2.1.1 Notação Indicial Cartesiana

A Notação Indicial Cartesiana faz uso dos índices subscritos 1, 2, 3 para

representar as direções Cartesianas x, y, z, e dispensa símbolos de somatório

quando ocorrem índices subscritos repetidos uma vez num mesmo termo da

expressão. Por exemplo:

∑

=

=++=++=

3

1i

ii

2

3

2

1332211ii

aaaaaaaaaaaaa

e

10

∑

=

=++=

3

1i

ii332211ii

aaaaa

Essa convenção de somatório é conhecida como Convenção de Somatório

de Einstein. É importante ressaltar que expressões tais como a

i

b

i

x

i

não estão

definidas nessa convenção, ou seja, um índice nunca deve ser repetido mais de uma

vez quando a Convenção de Somatório de Einstein for usada. Consequentemente,

uma expressão da forma

∑

=

3

1i

iii

xba

deve sempre reter seu sinal de somatório.

Com essa convenção, pode-se simplificar um sistema algébrico de N por N

expressões

∑

=

N

1j

ijij

yxA i = 1, 2, ..., N

da seguinte maneira:

ijij

yxA

=

i, j = 1, 2, ..., N

Os índices repetidos na expressão anterior implicam somatório. Vale

ressaltar que A

ij

≠ A

ji

a menos que a matriz A seja simétrica. Outra característica

dessa mesma expressão é a presença do índice i uma única vez em cada termo.

Um índice que aparece uma única vez em cada termo de uma equação é

chamado de índice livre e assume os valores 1, 2, ..., N por vez. Como exemplo,

deve-se expandir a expressão

jiji

nt

σ

=

i, j = 1, 2, 3

primeiramente, nas seguintes expressões:

jj11

nt

σ

=

jj22

nt

σ

=

jj33

nt

σ

=

Finalmente, a Convenção de Somatório de Einstein completa o sistema algébrico

equivalente:

3132121111

nnnt

σ

+

=

3232221212

nnnt

σ

+

=

11

3332321313

nnnt

σ

+

=

Se há dois índices livres aparecendo numa equação, como em

jmimij

AAT

=

i, j = 1, 2, 3

então essa equação é uma contração de nove equações, cada uma com um

somatório de três termos no lado direito. De fato:

131312121111m1m111

AAAAAAAAT

+

=

231322122111m2m112

AAAAAAAAT

+

=

331332123111m3m113

AAAAAAAAT

+

=

..................................................................

333332323131m3m333

AAAAAAAAT

+

=

Utilizando-se a Notação Indicial Cartesiana, a distância r entre dois pontos x

e x' pode ser escrita na forma

(

)

2/1

ii

rrr = (2.1)

com

'xxr

iii

−

=

(2.2)

As derivadas de r em relação aos pontos x e x' são dadas por

r

,r

x

r

i

==

∂

r

,r

'x

r

i

−=−=

∂

(2.3)

A derivada normal de r é dada pela expressão

ii

n,r

n

r

=

∂

(2.4)

onde n

i

denota as componentes do vetor unitário normal à superfície de contorno

direcionado para fora do corpo.

2.1.2 Delta de Kronecker

O Delta de Kronecker, simbolizado por

δ

ij

, é definido como

⎩

⎨

⎧

≠

=

jise0

jise1

ij

δ

(2.5)

12

2.1.3 Função Delta de Dirac

A representação matemática de excitações concentradas, tais como forças

pontuais em mecânica dos sólidos e dos fluidos, cargas elétricas pontuais em

eletrostática ou forças de impulso em acústica, pode ser alcançada com a adoção da

função Delta de Dirac.

A função Delta de Dirac, simbolizada por

Δ(X',X), é definida como

()

⎩

⎨

⎧

=∞

=

contráriocaso0

para

,

X'X

XX'

Δ

(2.6)

A principal característica da função Delta de Dirac é possuir o valor zero em

todos os pontos X, exceto no ponto X = X', denominado ponto fonte, onde a função

torna-se infinita. Assim, essa função representa singularidades pontuais no ponto X'.

A função Delta de Dirac possui a seguinte propriedade:

()() ()

X'XXX'X fdf,

b

a

=

∫

Δ

(2.7)

onde – ∞ ≤ a, b ≤ + ∞ e a ≤ X' ≤ b.

2.1.4 Teorema da Divergência

O Teorema da Divergência é uma identidade fundamental que relaciona uma

integral de volume com uma integral de superfície, ou seja,

∫∫

=

S

ii

V

i,i

dSnfdVf (2.8)

onde f

i,i

= ∂f

i

/ ∂x

i

e n

i

denota as componentes do vetor unitário normal à superfície S

direcionado para fora do corpo.

13

2.1.5 Série de Taylor

A Série de Taylor de uma função f (real ou complexa) infinitamente

diferenciável, definida em um intervalo aberto (a – r , a + r), é a série de potências

dada por

()

(

)

(

)

()

∑

∞

=

−=

0n

n

ax

!n

af

xT (2.9)

onde n! é o fatorial de n e f

(n)

(a) denota a n-ésima derivada de f no ponto a. Se a = 0,

a série também é chamada Série de MacLaurin.

Se f não for indefinidamente diferenciável em a, isto é, se f só admitir n

derivadas em a, então vale a fórmula de Taylor

() () ( ) ()

(

)

()

(

)

()

++

−

+

−

+−+= ...af

!3

ax

af

!2

ax

afaxafxf

3

2

1

(

)

()

() ()

xraf

!n

ax

...

n

+

−

++ (2.10)

onde r

n

(x) é uma função de x tal que

(

)

()

0

ax

xr

lim

n

ax

=

−

→

(2.11)

2.1.6 Quadratura Gaussiana

Muitas das integrais que necessitam ser calculadas no âmbito da aplicação

do Método de Elementos de Contorno não são triviais, ou seja, ou a primitiva da

função integrando não existe explicitamente, ou é demasiado complicado obtê-la,

inviabilizando, desse modo, a sua utilização prática, principalmente em códigos

computacionais. Por esse motivo é essencial recorrer a técnicas de integração

numérica, que também recebem a designação de regras de quadratura.

Em análise numérica, uma regra de quadratura é uma aproximação de uma

integral definida de uma função, usualmente definida como uma soma ponderada

dos valores da função em pontos específicos do domínio de integração.

14

A regra de quadratura de Gauss, ou Gaussiana, é uma regra de quadratura

de n pontos construída de modo a obter resultados exatos para polinômios de grau

2n-1 através da escolha adequada de n pontos ξ

i

e n pesos w

i

. O domínio de

integração para tal regra é normalmente tomado como sendo o intervalo fechado

[-1,1], assim a quadratura de Gauss é definida como

() () ()

∑∑

∫

==

+

−

≅+=

n

1i

iin

n

1i

ii

1

fwEfwdf

ξξξξ

(2.12)

onde n é o número de pontos de integração, ξ

i

é a coordenada do i-ésimo ponto de

integração, w

i

é o peso associado ao i-ésimo ponto de integração ξ

i

, e E

n

é o erro ou

resíduo, dado pela expressão

(

)

()()

[]

(

)

n2

3

4

1n2

n

d

fd

!n21n2

!n2

E

ξ

+

=

+

(2.13)

Os pontos de integração ξ

i

correspondem às raízes do Polinômio de

Legendre, definido como

()

⎥

⎦

⎤

⎢

⎣

⎡

−=

n

2

n

1

d

!n2

1

P

ξ

(2.14)

e os pesos w

i

são calculados através da relação

()

2

i

n

2

i

d

dP

1

2

w

ξξ

ξ

=

⎥

⎦

⎤

⎢

⎣

⎡

−

= (2.15)

É muito comum encontrar os valores dos pontos de integração ξ

i

e dos

pesos

w

i

já calculados e tabelados para uma determinada faixa de valores de ξ

i

,

geralmente para valores de 2 a 10 (ver, por exemplo, BREBBIA e DOMINGUEZ

1

).

2.1.7 Requisitos de Continuidade

Na matemática, muitas vezes faz-se necessário assegurar que as funções

ou as classes de funções com as quais se está trabalhando satisfaçam certos

1

BREBBIA, C. A., DOMINGUEZ, J. Boundary Elements, An Introductory Course. Boston:

Computational Mechanics Publications, 1989, p. 268.

15

requisitos de continuidade, de modo a garantir a existência, validade ou mesmo

unicidade de certas propriedades ou características.

Intuitivamente, uma função contínua em um ponto p de seu domínio é uma

função cujo gráfico não apresenta "salto" em p.

2

Uma definição mais rigorosa para o

conceito de continuidade é obtida do Cálculo Diferencial e Integral

3

, e pode ser

enunciada na forma

() ( )

⎪

⎩

⎪

⎨

⎧

<−⇒<−

∈∀>∃>∀

⇔

εδ

δε

pfxfpx

,Dx,quetal0,dado0

pemcontínuaf

f

(2.16)

Diz-se que f é continua em

f

DA ⊂ se f for contínua em todo Ap ∈ . Diz-se,

simplesmente, que f é uma função contínua se f for contínua em todo p de seu

domínio. A notação C(Ω) ou C

0

(Ω) é usada para denotar o conjunto de todas as

funções contínuas com domínio Ω. De maneira similar, C

1

(Ω) é usado para denotar o

conjunto das funções diferenciáveis com domínio Ω cujas primeiras derivadas são

contínuas, C

2

(Ω) é usado para denotar o conjunto das funções duplamente

diferenciáveis com domínio Ω cujas segundas derivadas são contínuas, e assim por

diante.

2.1.7.1 Condição de Continuidade de Hölder

Uma função satisfaz a condição de continuidade de Hölder se existem

constantes k > 0 e 0 <

α

≤ 1 tal que

() ()

α

sxksfxf −≤− (2.17)

para x suficientemente próximo de s. A continuidade de Hölder é "mais restritiva" que

a continuidade normal, de modo que apenas funções "patológicas" são contínuas

sem satisfazer também a condição de Hölder.

4

2

GUIDORIZZI, H. L. Um Curso de Calculo, Vol. 1. 5. ed. Rio de Janeiro: LTC, 2001 (reimp.

2004), p. 54.

3

Ibid., p. 62.

4

GUIGGIANI, M. The evaluation of Cauchy principal value integrals in the Boundary

Element Method – A review. Mathematical Computational Modelling, vol. 15, n. 3-5, p. 176, 1991.

16

2.1.8 Valor Principal de Cauchy

Integrais de funções que se tornam infinitas em um ponto do intervalo de

integração são denominadas integrais impróprias. Como será visto no próximo

capítulo, o desenvolvimento matemático do Método de Elementos de Contorno

conduz a várias integrais dessa natureza, sendo oportuno, portanto, realizar uma

breve descrição do tratamento matemático empregado em tais situações.

Dada a integral imprópria

()

∫

=

b

a

dxxfI

que apresenta uma assíntota vertical em x = c (ou seja, uma descontinuidade infinita

em x = c), sendo a < c < b, então I pode ser calculada como

() ()

∫∫

+

→

−

→

+=

b

c

0

c

a

0

dxxflimdxxflimI

δ

ε

Se os dois limites existirem, a integral converge ou é convergente. Se, por

outro lado

()

±∞=

∫

−

→

ε

c

a

0

dxxflim

e/ou

()

±∞=

∫

+

→

b

c

0

dxxflim

δ

então a integral diverge ou é divergente ou é não convergente.

Fazendo δ = ε, a integral imprópria não convergente pode existir no sentido

do Valor Principal de Cauchy, possuindo um valor finito, definido como

() () ()

⎥

⎦

⎤

⎢

⎣

⎡

+=

∫∫∫

+

−

→

b

c

a

0

b

a

dxxfdxxflimdxxfVPC

ε

(2.18)

embora

()

±∞=

∫

−

→

ε

c

a

0

dxxflim

e/ou

17

()

±∞=

∫

+

→

b

c

0

dxxflim

ε

Para ilustrar esse conceito considere-se a função

()

α

x

1

xf =

para a qual x = 0 é uma assíntota vertical. A integral

()

∫∫

+

−

+

−

==

1

dx

x

1

dxxfI

α

deve ser calculada como

=

−

+

−

=+=

+

−

→

−

→

+

→

−

→

∫∫

1

0

1

0

1

0

1

0

1

x

lim

1

x

limdx

x

1

limdx

x

1

limI

δ

α

δ

ε

α

ε

δ

α

δ

ε

α

ε

αα

() ()

[]

[

]

⎭

⎬

⎫

⎩

⎨

⎧

−

+

⎭

⎬

⎫

⎩

⎨

⎧

−−−

−

=

−

→

−−

→

α

δ

αα

ε

δ

α

ε

α

1

0

11

0

1

lim1

1

lim

Se

α

< 1, então k = 1-

α

> 0 e a integral imprópria converge, pois

(

)

⎪

⎩

⎪

⎨

⎧

=

=−

→

0lim

k

0

k

0

δ

ε

δ

ε

No entanto, se

α

= 3, por exemplo, então

()

3

x

1

xf =

e

() ()

()

⎪

⎩

⎪

⎨

⎧

+∞===

+∞=

−

=−=−

→

−

→→

→

−

→→

2

0

2

0

k

0

2

0

2

0

k

0

1

limlimlim

1

limlimlim

δ

δδ

ε

εε

δδδ

εεε

ou seja, a integral imprópria diverge.

Calculando o Valor Principal de Cauchy para a função

()

α

x

1

xf =

obtém-se

() ()

[

]

α

αα

ε

α

εε

α

−

−−

→

+

−

−+−−−

−

=

∫

1

11

0

1

11lim

1

dx

x

1

VPC

Quando

α

= 3 encontra-se

() ()

[

]

=−+−−−

−

=

−

−−

→

+

−

∫

2

22

0

1

3

11lim

31

1

dx

x

1

VPC

εε

ε

18

() ()

[]

011lim

2

1

11

lim

2

1

22

0

2

22

0

=−+−−=

⎥

⎦

⎤

⎢

⎣

⎡

−+

−

−=

−−

→

−

→

εεε

ε

εε

2.2 Equações Básicas da Elasticidade Linear

Nesta seção serão apresentadas as principais relações matemáticas,

concernentes à teoria da elasticidade linear, que serão necessárias para o

desenvolvimento do Método de Elementos de Contorno e suas aplicações a

problemas de Estado Plano de Tensões e Estado Plano de Deformações. Para uma

explanação mais detalhada sobre a formulação matemática da teoria da elasticidade

linear podem ser consultadas as obras de TIMOSHENKO

5

, HIBBELER

6

,

CRANDALL

7

, entre outros.

2.2.1 Equações de Equilíbrio

O estado de tensões em um ponto de um corpo tridimensional pode ser

descrito pelas componentes de tensões atuando em um volume elementar de lados

dx, dy, dz, como apresentado na Figura 2.1. O primeiro índice subscrito de cada

componente de tensão representa a direção da normal ao plano no qual a tensão

atua, e o segundo, a direção na qual a tensão propriamente dita atua.

Das nove componentes de tensões, apenas seis são essenciais para

descrever o estado de tensões de um ponto no interior de um corpo tridimensional.

Isso porque nem todas as componentes são independentes, mas estão relacionadas

5

TIMOSHENKO, S. P., GOODIER, J. N. Theory of Elasticity. 3. ed. Tokyo: McGraw-Hill,

1970.

6

HIBBELER, R. C. Resistência dos Materiais. 5. ed. São Paulo: Pearson Prentice Hall,

2004.

7

CRANDALL, Stephen H. et al. An Introduction to the Mechanics of Solids. 2. ed.

Singapore: McGraw-Hill, 1978.

19

através de equações de equilíbrio. Aplicando a equação de equilíbrio de momentos,

que é escrita tomando os momentos das componentes de tensão em relação a um

determinado ponto do volume elementar, produzem-se as chamadas relações

complementares de cisalhamento, isto é,

yxxy

σ

=

zxxz

σ

=

zyyz

σ

=

ou então

jiij

σ

=

i, j = 1, 2, 3 (2.19)

FIGURA 2.1 – Componentes de tensões atuando nas faces de um volume elementar

representando um ponto no interior de um corpo tridimensional

As equações de equilíbrio, para um sistema tridimensional sujeito a forças

externas e forças de corpo b

x

, b

y

e b

z

, são dadas por

0b

zyx

x

xz

xy

xx

=+

∂

+

∂

+

∂

σ

0b

zyx

y

yzyyyx

=+

∂

+

∂

+

∂

σ

0b

zyx

z

zz

zy

zx

=+

∂

+

∂

+

∂

σ

as quais podem ser reescritas como

20

0b

xxx

i

k

ik

j

ij

i

ii

=+

∂

+

∂

+

∂

σ

i ≠ j ≠ k i, j, k = 1, 2, 3

ou ainda, de forma mais compacta,

0b

ij,ij

=

+

σ

i, j = 1, 2, 3 (2.20)

Essas equações devem ser satisfeitas em todo o interior do corpo e, por

isso, são denominadas condições de equilíbrio no domínio Ω do corpo.

2.2.2 Estado Plano de Tensões

As equações de equilíbrio podem ser simplificadas para a condição

bidimensional de Estado Plano de Tensões. Nessa condição, a qual é aplicável a

placas delgadas, admite-se que as tensões na direção da espessura da placa são

negligenciáveis. Se o eixo z estiver orientado segundo a espessura da placa, então

as componentes σ

zz

, σ

xz

e σ

yz

serão todas nulas. As componentes não nulas σ

xx

, σ

yy

e

σ

xy

serão avaliadas sobre a espessura da placa e serão consideradas independentes

de z. As equações de equilíbrio simplificam para

0b

yx

x

xy

xx

=+

∂

+

∂

σ

0b

yx

y

yyyx

=+

∂

+

∂

σ

ou

0b

ij,ij

=

+

σ

i, j = 1, 2 (2.21)

Essas equações satisfazem os requisitos de equilíbrio em todos os pontos

do domínio Ω do corpo. No entanto, o equilíbrio também deve ser satisfeito em todos

os pontos do contorno Γ do corpo, onde as componentes de forças de superfície por

unidade de área são denotadas por t

x

e t

y

.

Para o elemento triangular da Figura 2.2, localizado no contorno de um

corpo bidimensional de espessura unitária constante, as forças de superfície na face

AB do contorno devem estar em equilíbrio com as forças internas nas faces internas

AC e CB.

21

FIGURA 2.2 – Forças de superfície e forças internas atuando em um elemento localizado no

contorno de um corpo bidimensional de espessura unitária constante

Somando todas as forças na direção x e isolando para t

x

encontra-se

0dxdydst

yxxxx

=

−

σ

ds

dx

ds

dy

t

yxxxx

σσ

+=

Da Figura 2.3 é possível obter as relações

ds

dy

lcosn

x

===

θ

ds

dx

mcosn

y

===

γ

(2.22)

As derivadas dy/ds e dx/ds são os cossenos diretores l e m dos ângulos que o

vetor unitário normal à face AB forma com os eixos coordenados x e y,

respectivamente. Desse modo,

mlt

yxxxx

σ

+

=

e, de maneira análoga, para a direção y,

mlt

yyxyy

σ

+

=

As derivadas dy/ds e dx/ds também podem ser consideradas como sendo

iguais às componentes ou projeções do vetor unitário normal à face

AB segundo os

22

eixos coordenados x e y, respectivamente. Assim, as equações anteriores são

reescritas como

yyxxxxx

nnt

σ

+

=

yyyxxyy

nnt

σ

+

=

ou, de forma mais compacta,

jiji

nt

σ

=

i, j = 1, 2 (2.23)

FIGURA 2.3 – Componentes ou projeções do vetor unitário normal ao contorno

A extensão para sistemas tridimensionais é relativamente simples, a saber:

nmlt

zxyxxxx

σ

+

=

nmlt

zyyyxyy

σ

+

=

nmlt

zzyzxzz

σ

+

=

onde l, m e n são os cossenos diretores dos ângulos que o vetor unitário normal à

superfície do contorno forma com os eixos coordenados

x, y e z, respectivamente.

Logo, as equações anteriores podem ser reescritas na forma

zzxyyxxxxx

nnnt

σ

+

=

zzyyyyxxyy

nnnt

σ

+

=

23

zzzyyzxxzz

nnnt

σ

+

=

resultando em

jiji

nt

σ

=

i, j = 1, 2, 3 (2.24)

Essas equações devem ser satisfeitas em todo o contorno do corpo e, por

isso, são denominadas equações ou condições de equilíbrio no contorno Γ do corpo.

As expressões para tensões em planos inclinados podem ser obtidas tendo

em vista um elemento triangular extraído de um elemento localizado no interior do

corpo, através de um corte inclinado de um ângulo θ em relação ao eixo y (Figura

2.4). Para que o elemento resultante ABC (Figura 2.5) esteja em equilíbrio, as

tensões atuando nas faces AB, BC e CA necessitam estar em equilíbrio. Como o

elemento possui espessura unitária constante, a área das faces é proporcional aos

comprimentos dos lados do triângulo.

Considerando, primeiramente, que o elemento triangular esteja sujeito

apenas a tensões normais e calculando o somatório das forças que atuam segundo

a direção normal à face AB, obtém-se

0sinBCcosACAB

yyxxn

=

−

θ

σ

θ

σ

Dividindo todos os termos por AB,

θσθσσ

sin

AB

BC

cos

AB

AC

yyxxn

+=

utilizando a Figura 2.5,

θσθσσ

2

yy

2

xxn

sincos +=

e aplicando as identidades trigonométricas pertinentes, resulta

(

)

(

)

θσσσσσ

2cos

2

1

2

1

yyxxyyxxn

−++= (2.25)

De maneira similar, mas desta vez calculando o somatório das forças que

atuam segundo a direção paralela à face AB, resulta

0cosBCsinACAB

yyxxt

=

−

+

θ

σ

θ

σ

θσθσσ

cos

AB

BC

sin

AB

AC

yyxxt

+−=

θ

σ

θ

σ

cossinsincos

yyxxt

+

−

=

(

)

θσσσ

2sin

2

1

yyxxt

−−= (2.26)

24

FIGURA 2.4 – As tensões atuando em planos inclinados podem ser obtidas a partir de um

elemento triangular extraído de um elemento localizado no interior do corpo

FIGURA 2.5 – Elemento triangular resultante e as tensões que atuam sobre o mesmo

25

Repetindo o processo anterior, mas considerando, desta vez, que o

elemento triangular esteja sujeito apenas a tensões tangenciais, obtém-se

0cosBCsinACAB

yxxyn

=

−

θ

σ

θ

σ

θσθσσ

cos

AB

BC

sin

AB

AC

yxxyn

+=

θ

σ

θ

σ

cossinsincos

yxxyn

+

=

(

)

θ

σ

cossin2

xyn

=

θ

σ

2sin

xyn

=

(2.27)

e

0sinBCcosACAB

yxxyt

=

+

−

θ

σ

θ

σ

θσθσσ

sin

AB

BC

cos

AB

AC

yxxyt

+=

θσθσσ

2

yx

2

xyt

sincos −=

(

)

θθσσ

22

xyt

sincos −=

⎟

⎠

⎞

⎜

⎝

⎛

+−+=

θθσσ

2cos

2

1

2

1

2cos

2

1

2

1

xyt

θ

σ

2cos

xyt

=

(2.28)

Finalmente, para um estado de tensões bidimensionais geral, as expressões

das tensões atuando num plano inclinado são escritas como combinações dos

resultados anteriormente obtidos, ou seja,

()

(

)

θσθσσσσσ

2sin2cos

2

1

2

1

xyyyxxyyxxn

+−++=

(2.29)

e

(

)

θσθσσσ

2cos2sin

2

1

xyyyxxt

+−−= (2.30)

2.2.3 Tensões Principais

Dado que σ

n

varia com o ângulo θ, ou seja, é função de θ, então σ

n

terá um

valor máximo ou mínimo quando dσ

n

/ dθ = 0. Da expressão (2.29) tem-se

26

()

02cos22sin

d

xyyyxx

n

=+−−=

θσθσσ

θ

σ

(2.31)

yyxx

xy

2

2tan

2cos

2sin

σσ

σ

θ

−

==

(2.32)

A equação (2.32) apresenta duas soluções distintas: os ângulos 2θ e 2θ + π,

que diferem de 180º, ou ainda os ângulos θ e θ + π/2, que diferem de 90º. A Figura

2.6 ilustra esses ângulos graficamente. Os planos orientados segundo esses

ângulos, nos quais as tensões normais são um máximo ou um mínimo, são

chamados Planos Principais, e as tensões normais que atuam nesses planos são

denominadas Tensões Principais. Tomando a expressão (2.31) e dividindo todos os

termos por dois resulta

()

02cos2sin

2

1

txyyyxx

==+−−

σθσθσσ

isto é, a tensão tangencial é nula nos Planos Principais.

Reportando-se novamente à Figura 2.6 verificam-se as seguintes relações:

()

2

xy

2

yyxx

xy

4

2

2sin

σσσ

σ

θ

+−

=

(2.33)

()

2

xy

2

yyxx

4

2cos

σσσ

σ

θ

+−

−

= (2.34)

()

2

xy

2

yyxx

xy

4

2

2sin

σσσ

σ

π

θ

+−

−

=

⎟

⎠

⎞

⎜

⎝

⎛

+ (2.35)

(

)

()

2

xy

2

yyxx

4

2

2cos

σσσ

σ

π

θ

+−

−

=

⎟

⎠

⎞

⎜

⎝

⎛

+ (2.36)

Substituindo (2.33) e (2.34) em (2.29) resulta

(

)

(

)

⎥

⎦

⎤

⎢

⎣

⎡

+−++=

2

xy

2

yyxxyyxxn

4

2

1

σσσσσσ

Substituindo (2.35) e (2.36) em (2.29) resulta

(

)

(

)

⎥

⎦

⎤

⎢

⎣

⎡

+−−+=

2

xy

2

yyxxyyxxn

4

2

1

σσσσσσ

Desse modo, é possível concluir que as Tensões Principais são dadas por

()

(

)

(

)

⎥

⎦

⎤

⎢

⎣

⎡

+−±+=

2

xy

2

yyxxyyxx

mín

máx

n

4

2

1

σσσσσσ

(2.37)

27

FIGURA 2.6 – Solução geométrica para determinar os ângulos que definem os Planos

Principais nos quais atuam as Tensões Principais

As tensões tangenciais máxima e mínima podem ser encontradas de

maneira semelhante, considerando agora dσ

t

/ dθ = 0. Da expressão (2.30) obtém-se

()

02sin22cos

d

xyyyxx

t

=−−−=

θσθσσ

θ

σ

(2.38)

(

)

xy

yyxx

2

2tan

2cos

2sin

σ

θ

−

==

(2.39)

A Figura 2.7 ilustra graficamente os ângulos que são soluções da equação

(2.39). Conclui-se que os planos nos quais atuam as tensões tangenciais máxima e

mínima são defasados 45º dos Planos Principais.

Da referida Figura 2.7 verificam-se as seguintes relações:

(

)

()

2

xy

2

yyxx

4

2sin

σσσ

σ

θ

+−

−

= (2.40)

()

2

xy

2

yyxx

xy

4

2

2cos

σσσ

σ

θ

+−

= (2.41)

()

2

xy

2

yyxx

4

2

2sin

σσσ

σ

π

θ

+−

−

=

⎟

⎠

⎞

⎜

⎝

⎛

+ (2.42)

28

()

2

xy

2

yyxx

xy

4

2

2cos

σσσ

σ

π

θ

+−

−

=

⎟

⎠

⎞

⎜

⎝

⎛

+ (2.43)

FIGURA 2.7 – Solução geométrica para determinar os ângulos que definem os planos nos

quais atuam as tensões tangenciais máxima e mínima

Substituindo (2.40) e (2.41) em (2.30) resulta

(

)

2

xy

2

yyxxt

4

2

1

σσσσ

+−=

Substituindo (2.42) e (2.43) em (2.30) resulta

(

)

2

xy

2

yyxxt

4

2

1

σσσσ

+−−=

Desse modo, é possível concluir que as tensões tangenciais máxima e

mínima são dadas por

()

(

)

2

xy

2

yyxx

mín

máx

t

4

2

1

σσσσ

+−±= (2.44)

Substituindo (2.40) e (2.41) ou (2.42) e (2.43) em (2.29) resulta

(

)

yyxxn

2

1

σσσ

+= (2.45)

ou seja, os planos nos quais atuam as tensões tangenciais máxima e mínima

possuem uma tensão normal associada, a qual é fornecida pela expressão (2.45).

29

Outra relação que pode ser obtida a partir das expressões (2.37) e (2.44) é:

() () ()

[]

mín

n

máx

n

máx

t

2

1

σσσ

−= (2.46)

Os valores das tensões normais máxima e mínima e das tensões

tangenciais máxima e mínima, bem como os valores dos ângulos de orientação dos

planos nos quais atuam essas tensões, podem ser obtidos utilizando-se o processo

gráfico do Círculo de Mohr que, em geral, é fácil de usar e lembrar. Além disso, essa

abordagem permite visualizar como as tensões normais e tangenciais variam

conforme o plano em que atuam esteja orientado segundo diferentes ângulos.

Detalhes sobre a construção do Círculo de Mohr e exemplos de sua aplicação

podem ser encontrados em HIBBELER

8

.

2.2.4 Deformações Específicas

O estado de deformações em um ponto de um corpo tridimensional pode ser

descrito pelas componentes de deformações específicas atuando em um volume

elementar de lados dx, dy, dz, como apresentado na Figura 2.8. As componentes ε

xx

,

ε

yy

e ε

zz

são denominadas deformações específicas longitudinais e as componentes

ε

xy

, ε

xz

e ε

yz

deformações específicas transversais. Observa-se que, em particular,

deformações específicas longitudinais provocam mudança de volume do elemento

infinitesimal, enquanto deformações específicas transversais provocam mudança na

sua geometria.

As seis componentes de deformações específicas, correspondentes às

componentes de tensões, podem ser escritas, em notação indicial, como

jiij

ε

=

i, j = 1, 2, 3 (2.47)

Representando os deslocamentos por u

x

, u

y

e u

z

(u

i

e i = 1, 2, 3 em notação

indicial), as deformações específicas estão relacionadas aos mesmos por intermédio

das expressões

8

HIBBELER, R. C. Resistência dos Materiais. 5. ed. São Paulo: Pearson Prentice Hall,

2004, p. 360-367.

30

x

u

x

xx

∂

=

ε

y

u

y

yy

∂

=

ε

z

u

z

zz

∂

=

ε

(2.48)

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

=

y

u

x

u

2

1

x

y

xy

ε

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

=

z

u

x

u

2

1

x

z

xz

ε

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

=

z

u

y

u

2

1

y

z

yz

ε

(2.49)

Alternativamente, em notação indicial, é possível reescrever essas relações como

(

)

i,jj,iij

uu

2

1

+=

ε

i, j = 1, 2, 3 (2.50)

FIGURA 2.8 – Componentes de deformações específicas atuando em um volume elementar

representando um ponto no interior de um corpo tridimensional: a) Volume elementar sem

deformação; b) Volume elementar deformado

2.2.5 Equações de Compatibilidade

Diferenciando a primeira expressão em

(2.49) em relação a x e y,

⎟

⎠

⎞

⎜

⎝

⎛

∂

∂∂

∂

+

⎟

⎠

⎞

⎜

⎝

⎛

∂

∂∂

∂

=

∂∂

∂

y

u

yxx

u

yxyx

2

x

2

y

2

xy

2

ε

reordenando os termos e a ordem das derivadas parciais,

31

0

x

u

yyy

u

xxyx

2

x

2

y

2

xy

2

=

⎟

⎠

⎞

⎜

⎝

⎛

∂

∂∂

∂

+

⎟

⎠

⎞

⎜

⎝

⎛

∂

∂∂

∂

+

∂∂

∂

−

ε

e utilizando as relações em (2.48), obtém-se

0

yx

2

xx

2

yy

2

xy

2

=

∂

+

∂

+

∂∂

∂

−

ε

εε

(2.51)

De maneira similar, para as demais expressões em (2.49), deduz-se

0

zx

2

xx

2

zz

2

xz

2

=

∂

+

∂

+

∂∂

∂

−

ε

(2.52)

0

yz

zy

2

zz

2

yy

2

yz

2

=

∂

+

∂

+

∂∂

∂

−

ε

εε

(2.53)

Diferenciando cada expressão em (2.49) do seguinte modo

⎟

⎠

⎞

⎜

⎝

⎛

∂

−

∂

=

⎟

⎠

⎞

⎜

⎝

⎛

∂

⇒

⎥

⎦

⎤

⎢

⎣

⎡

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

⎟

⎠

⎞

⎜

⎝

⎛

∂

=

∂

y

u

zz

2

x

u

zy

u

zx

u

z2

1

z

x

xyy

x

yxy

εε

⎟

⎠

⎞

⎜

⎝

⎛

∂

−

∂

=

⎟

⎠

⎞

⎜

⎝

⎛

∂

⇒

⎥

⎦

⎤

⎢

⎣

⎡

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

⎟

⎠

⎞

⎜

⎝

⎛

∂

=

∂

x

u

yy

2

z

u

yz

u

yx

u

y2

1

y

z

xzxx

z

xz

εε

⎟

⎠

⎞

⎜

⎝

⎛

∂

−

∂

=

⎟

⎠

⎞

⎜

⎝

⎛

∂

⇒

⎥

⎦

⎤

⎢

⎣

⎡

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

⎟

⎠

⎞

⎜

⎝

⎛

∂

=

∂

z

u

xx

2

y

u

xz

u

xy

u

x2

1

x

yyz

z

y

z

yz

εε

e diferenciando a primeira expressão em (2.48) em relação a y e z, reordenando as

derivadas parciais e substituindo progressivamente cada um dos termos isolados

anteriormente, obtém-se

⎟

⎠

⎞

⎜

⎝

⎛

⎟

⎠

⎞

⎜

⎝

⎛

∂

−

∂

=

⎟

⎠

⎞

⎜

⎝

⎛

⎟

⎠

⎞

⎜

⎝

⎛

∂

=

⎟

⎠

⎞

⎜

⎝

⎛

⎟

⎠

⎞

⎜

⎝

⎛

∂

=

∂∂

∂

x

u

yy

2

xz

u

yxx

u

zyzy

z

xzxxxx

2

εε

⎟

⎠

⎞

⎜

⎝

⎛

⎟

⎠

⎞

⎜

⎝

⎛

∂

−

∂

+

∂

−

∂

=

⎟

⎠

⎞

⎜

⎝

⎛

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

−

∂

=

∂∂

∂

y

u

zz

2

x

2

y

2

xz

u

xx

2

y

2

xzy

x

xyyz

xz

yyz

xzxx

2

εε

ε

εε

⎟

⎠

⎞

⎜

⎝

⎛

⎟

⎠

⎞

⎜

⎝

⎛

∂

−

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

+

∂

−

∂

=

∂∂

∂

y

u

zxzyxx

2

zy

x

xy

xz

yz

xx

2

ε

zyzyxx

2

zy

xx

2

xy

xz

yz

xx

2

∂∂

∂

−

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

+

∂

−

∂

=

∂∂

∂

ε

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

+

∂

−

∂

=

∂∂

∂

zyxxzy

xy

xz

yz

xx

2

ε

(2.54)

De modo análogo, para as demais expressões em (2.48), encontra-se

32

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

−

∂

=

∂∂

∂

zyxyzx

xy

xz

yzyy

2

ε

εε

(2.55)

⎟

⎠

⎞

⎜

⎝

⎛

∂

−

∂

+

∂

=

∂∂

∂

zyxzyx

xy

xz

yz

zz

2

ε

(2.56)

As equações (2.51) a (2.56) são denominadas equações de compatibilidade

das deformações específicas e podem ser reescritas em notação indicial como

0

xx

2

i

jj

2

j

ii

2

ji

ij

2

=

∂

+

∂

+

∂∂

∂

−

ε

(2.57)

0

xxxxxx

k

ij

j

ik

i

jk

ikj

ii

2

=

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

+

∂

−

∂

−

∂∂

∂

ε

(2.58)

onde i, j, k = 1, 2, 3 e i ≠ j ≠ k.

2.2.6 Relações Tensão-Deformação

Para uma barra de seção transversal constante, sujeita a um carregamento

axial uniforme (direcionado segundo o eixo x), as seguintes relações são válidas

para pontos localizados suficientemente longe das extremidades da barra:

E

xx

σ

ε

=

E

xx

yy

σ

νε

−=

E

xx

zz

σ

νε

−= (2.59)

onde E é o módulo de elasticidade longitudinal ou módulo de Young e ν é o

coeficiente de Poisson. As equações em

(2.59) são expressões da Lei de Hooke.

Para um corpo sujeito a um carregamento geral, as deformações específicas são

dadas por

(

)

[

]

zzyyxxxx

E

1

σσνσε

+−=

(2.60)

()

[

]

zzxxyyyy

E

1

σσνσε

+−= (2.61)

(

)

[

]

yyxxzzzz

E

1

σσνσε

+−= (2.62)

xyxy

E

1

σ

ν

ε

+

=

xzxz

E

1

σ

ν

ε

+

=

yzyz

E

1

σ

ν

ε

+

= (2.63)

33

Tomando as equações (2.60) a (2.62), é possível expressar as tensões

normais em função das deformações específicas. Primeiramente, essas equações

devem ser reescritas da seguinte maneira:

E

xxzzyyxx

ε

ν

σ

ν

σ

=

−

E

yyzzyyxx

ε

ν

σ

ν

σ

=

−

+

−

E

zzzzyyxx

ε

σ

ν

σ

ν

σ

=

+

−

Colocando o sistema de equações em forma matricial

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

=

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

−−

E

1

zz

yy

xx

zz

yy

xx

ε

σ

νν

e invertendo a matriz de coeficientes, obtém-se

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

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+−

−

+−+−

−

+−

+−+−+−

−

=

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

E

121

1

121121

1

121

121121121

1

zz

yy

xx

zz

yy

xx

ε

νν

ν

νν

ν

νν

ν

νν

ν

νν

ν

νν

ν

νν

ν

νν

ν

νν

ν

σ

Resolvendo para σ

xx

resulta

()()

νν

νε

ε

σ

+−

+

−

=

121

EEEE

zzyyxxxx

xx

ou, somando e subtraindo νε

xx

E

()()

νν

νε

ε

σ

+−

+

−

=

121

EEEE2E

zzyyxxxxxx

xx

Finalmente

(

)

()()

(

)

()()

νν

νε

νν

ε

ν

σ

+−

−

+

+−

+

=

121

21E

121

E

xx

zzyyxx

xx

Definindo

()

ν

μ

+

=

12

E

(2.64)

()

ν

νμ

λ

21

2

−

= (2.65)

zzyyxx

e

ε

+

=

(2.66)

34

respectivamente denominados módulo de elasticidade transversal, constante de

Lamé e deformação específica volumétrica, tem-se

xxxx

2e

με

λ

σ

+

=

(2.67)

e, para as demais tensões normais,

yyyy

2e

με

λ

σ

+

=

(2.68)

zzzz

2e

με

λ

σ

+

=

(2.69)

Para as tensões tangenciais o processo é direto. Das relações (2.63) resulta

xyxyxy

2

1

E

μεε

ν

σ

=

+

= (2.70)

xzxzxz

2

1

E

μεε

ν

σ

=

+

= (2.71)

yzyzyz

2

1

E

μεε

ν

σ

=

+

=

(2.72)

As equações (2.67) a (2.72) são denominadas Equações Constitutivas da

Elasticidade, e podem ser reescritas em notação indicial como

ijkkijij

2

με

ε

λδ

σ

+

= i, j, k = 1, 2, 3 (2.73)

2.2.7 Equações Governantes

As Equações Governantes da Elasticidade, conhecidas como Equações de

Navier, são as condições de equilíbrio expressas em termos dos deslocamentos.

Primeiramente, deve-se obter expressões para as tensões em função das derivadas

dos deslocamentos substituindo as relações (2.48) e (2.49) nas equações (2.67) a

(2.72), do que resulta

x

u

2e

x

xx

∂

+=

μλσ

y

u

2e

y

yy

∂

+=

μλσ

z

u

2e

z

zz

∂

+=

μλσ

35

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

=

y

u

x

u

x

y

xy

μσ

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

=

z

u

x

u

x

z

xz

μσ

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

=

z

u

y

u

y

z

yz

μσ

Substituindo essas relações na primeira equação (2.20) e utilizando as

relações (2.48) e (2.66) obtém-se

0b

zyx

x

xz

xy

xx

=+

∂

+

∂

+

∂

σ

0b

z

u

xz

u

y

u

xy

u

x

u

2

x

e

x

2

x

2

z

2

x

2

y

2

x

2

=+

∂

+

∂∂

∂

+

∂

+

∂∂

∂

+

∂

+

∂

μμμμμλ

0b

z

u

xz

u

y

u

xy

u

x

u

x

u

x

e

x

2

x

2

z

2

x

2

y

2

x

2

x

2

=+

∂

+

∂∂

∂

+

∂

+

∂∂

∂

+

∂

+

∂

+

∂

μμμμμμλ

0b

z

u

x

y

u

xx

x

u

x

e

x

2

x

2

zz

2

x

2

yy

xx

2

x

2

=+

∂

+

∂

+

∂

+

∂

+

∂

+

∂

+

∂

μ

ε

μμ

ε

μ

ε

μμλ

()

0b

x

z

u

y

u

x

u

x

e

xzzyyxx

2

x

2

x

2

x

2

=+++

∂

+

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

+

∂

+

∂

εεεμμλ

0b

x

e

u

x

e

xx

2

=+

∂

+∇+

∂

μμλ

()

0b

x

e

u

xx

2

=+

∂

++∇

μλμ

(2.74)

onde

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

+

∂

=∇

2

zyx

é o Operador Laplaciano. De maneira equivalente para

as demais equações

(2.20) encontra-se

()

0b

y

e

u

yy

2

=+

∂

++∇

μλμ

(2.75)

()

0b

z

e

u

zz

2

=+

∂

++∇

μλμ

(2.76)

As equações (2.74) a (2.76) podem ser reescritas em notação indicial na

forma

(

)

0buu

iji,jjj,i

=

+

μ

λ

μ

i, j = 1, 2, 3 (2.77)

36

sendo essa última equação a condição de equilíbrio no domínio Ω do corpo.

Utilizando um processo semelhante, mas desta vez substituindo as relações

apresentadas no início desta seção nas equações (2.24), deduz-se

zzxyyxxxxx

nnnt

σ

+

=

z

x

z

y

x

y

x

n

z

u

x

u

n

y

u

x

u

n

x

u

2et

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

+

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

+

⎟

⎠

⎞

⎜

⎝

⎛

∂

+=

μμμλ

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

+

∂

+

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

+

∂

+=

z

y

x

z

x

y

x

xx

n

x

u

n

x

u

n

x

u

n

z

u

n

y

u

n

x

u

ent

μμλ

(i)

zzyyyyxxyy

nnnt

σ

+

=

z

y

z

y

x

y

n

z

u

y

u

n

y

u

2en

y

u

x

u

t

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

+

⎟

⎠

⎞

⎜

⎝

⎛

∂

++

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

=

μμλμ

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

+

∂

+

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

+

∂

+=

z

y

x

z

y

x

y

yy

n

y

u

n

y

u

n

y

u

n

z

u

n

y

u

n

x

u

ent

μμλ

(ii)

zzzyyzxxzz

nnnt

σ

+

=

z

y

z

x

z

n

z

u

2en

z

u

y

u

n

z

u

x

u

t

⎟

⎠

⎞

⎜

⎝

⎛

∂

++

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

+

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

=

μλμμ

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

+

∂

+

⎟

⎠

⎞

⎜

⎝

⎛

∂

+

∂

+

∂

+=

z

y

x

z

y

z

x

z

zz

n

z

u

n

z

u

n

z

u

n

z

u

n

y

u

n

x

u

ent

μμλ

(iii)

É possível reescrever essas equações em notação indicial na forma

ji,jjj,iij,ji

nununut

μ

λ

+

= i, j = 1, 2, 3 (2.78)

ou ainda

(

)

i,jj,ijij,ji

uunnut

+

=

μ

λ

i, j = 1, 2, 3 (2.79)

sendo essa última equação a condição de equilíbrio no contorno Γ do corpo.

2.2.8 Condições para Estado Plano de Tensões

Para placas delgadas, se nenhuma carga é aplicada perpendicurlamente à

superfície da placa, as tensões na direção da espessura da placa podem ser

consideradas negligenciáveis. Partindo do pressuposto que o eixo z esteja orientado

37

segundo a espessura da placa, as tensões σ

zz

, σ

xz

e σ

yz

serão todas nulas, e as

tensões não nulas σ

xx

, σ

yy

e σ

xy

serão avaliadas sobre a espessura da placa e serão

admitidas independentes da direção z. Essas considerações podem ser resumidas

da seguinte maneira:

⎪

⎩

⎪

⎨

⎧

==

=

0

zyyz

zxxz

zz

σσ

σ

⎪

⎩

⎪

⎨

⎧

≠=

≠

0

yxxy

yy

xx

σσ

σ

Aplicando essas relações nas equações (2.60) a (2.63) obtém-se

(

)

yyxxxx

E

1

σνσε

−=

(

)

xxyyyy

E

1

σνσε

−=

()

yyxxzz

E

σσ

ν

ε

+−=

(

)

yyxxzz

1

εε

ν

ε

+

−

−=

xyyxxy

E

1

σ

ν

εε

+

==

0

zxxz

=

ε

0

zyyz

=

ε

2.2.9 Condições para Estado Plano de Deformações

As condições para Estado Plano de Deformações são admitidas para corpos

nos quais a geometria e o carregamento não variam significativamente na direção do

eixo longitudinal do corpo, que pode ser adotado como sendo o eixo z. Nesses

casos, as demais variáveis são funções apenas das coordenadas x e y. O

deslocamento na direção do eixo

z é nulo em todos os pontos da seção transversal

e, desse modo, as deformações específicas ε

zz

, ε

xz

e ε

yz

também se anulam. Essas

considerações podem ser resumidas da seguinte maneira:

⎪

⎩

⎪

⎨

⎧

==

=

0

zyyz

zxxz

zz

εε

ε

⎪

⎩

⎪

⎨

⎧

≠=

≠

0

yxxy

yy

xx

εε

ε

Aplicando essas relações nas equações (2.60) a (2.63) obtém-se

(

)

yyxxzz

σ

ν

σ

+

=

38

0

zxxz

=

σ

0

zyyz

=

σ

⎟

⎠

⎞

⎜

⎝

⎛

−

=

yyxx

2

xx

1E

1

σ

ν

σ

ν

ε

⎟

⎠

⎞

⎜

⎝

⎛

−

=

xxyy

2

yy

1E

1

σ

ν

σ

ν

ε

xyyxxy

E

1

σ

ν

εε

+

==

Das considerações anteriores verifica-se que a única equação de

compatibilidade de deformações específicas para a condição de Estado Plano de

Deformações é

0

yx

2

xx

2

yy

2

xy

2

=

∂

+

∂

+

∂∂

∂

−

ε

εε

As expressões para a condição de Estado Plano de Deformações podem

simular problemas de Estado Plano de Tensões quando se adota o coeficiente de

Poisson ν' e o módulo de Young E' modificados, definidos como

ν

+

=

1

'

⎥

⎦

⎤

⎢

⎣

⎡

⎟

⎠

⎞

⎜

⎝

⎛

+

−=

2

1

1E'E

ν

(2.80)

39

3 FORMULAÇÃO INTEGRAL DE CONTORNO

Existem diferentes métodos para a derivação do Método de Elementos de

Contorno, entre eles o teorema da reciprocidade e o Método de Resíduos

Ponderados. Neste trabalho, o Método de Elementos de Contorno será derivado

usando o teorema do trabalho recíproco de Betti.

3.1 Teorema do Trabalho Recíproco de Betti

A formulação direta do Método de Elementos de Contorno para elastostática

pode ser derivada do teorema do trabalho recíproco de Betti para dois estados auto-

equilibrados (Figura 3.1) representados por (u

i

, t

i

, b

i

) e (u

i

*, t

i

*, b

i

*); u

i

e u

i

* são

deslocamentos; t

i

e t

i

* são forças de superfície; e b

i

e b

i

* são forças de corpo.

Das equações de equilíbrio

0b

ij,ij

=

+

σ

i, j = 1, 2, 3

é possível escrever a seguinte relação:

(

)

0dub

*

iij,ij

=+

∫

Ω

σ

ou

40

0dubdu

*

ii

*

ij,ij

=+

∫∫

ΩΩ

σ

(3.1)

onde Ω denota o domínio, com contorno Γ, do problema. As tensões, forças de

corpo e de superfície e deslocamentos são funções de X

∈

Ω [X ≡ (x, y, z) para um

corpo tridimensional ou X ≡ (x, y) para um corpo bidimensional].

FIGURA 3.1 – Uma região arbitrária com domínio

*

Ω e contorno

*

Γ envolvendo

o problema com domínio

Ω e contorno Γ

Esse primeiro passo corresponde a um processo de ponderação da equação

diferencial governante, processo esse idêntico ao que é utilizado no Método de

Resíduos Ponderados. Nesse método a função de ponderação, normalmente

designada por

w (correspondendo a u

i

* neste caso), é arbitrária, exceto pelo fato de

que deve ser contínua no domínio Ω do problema e deve possuir derivadas

contínuas até um determinado grau, sendo que o grau de continuidade requerido

varia conforme o problema em questão. (Para uma explicação mais detalhada a

41

respeito da derivação do Método de Elementos de Contorno a partir do Método de

Resíduos Ponderados, ver BREBBIA e DOMINGUEZ

1

).

A primeira integral em (3.1) pode ser reescrita, usando a regra da derivada

de um produto de funções, como

(

)

∫∫∫

−=

ΩΩΩ

ΩΩΩ dududu

*

j,iij

j,

*

iij

*

ij,ij

σσσ

Observando que

()

(

)

*

j,iij

*

i,jji

*

j,iij

*

i,jij

*

j,iij

*

ijij

uuu

2

1

uu

2

1

σσσσσεσ

=+=+=

tem-se

(

)

∫∫∫

−=

ΩΩΩ

ΩΩΩ ddudu

*

ijij

j,

*

iij

*

ij,ij

εσσσ

e, usando o Teorema da Divergência na primeira integral à direita,

()

∫∫∫

==

ΓΓΩ

ΓΓΩ dutdundu

*

ii

*

ijij

j,

*

iij

σσ

resulta

∫∫∫

−=

ΩΓΩ

ΩΓΩ ddutdu

*

ijij

*

ii

*

ij,ij

εσσ

Substituindo esse último resultado na equação (3.1) deriva-se a expressão

∫∫∫

=+

ΩΩΓ

ΩΩΓ ddubdut

*

ijij

*

ii

*

ii

εσ

(3.2)

Obs.: essa equação dá origem à designação "trabalho", presente no nome

do teorema que se está demonstrando, uma vez que o integrando à direita do sinal

de igualdade pode ser entendido sob esse ponto de vista da seguinte maneira:

trabalho é igual à força (ou tensão) multiplicada por deslocamento (ou deformação

específica).

Utilizando as equações constitutivas da elasticidade, definidas anteriormente

como

ijkkijij

2

με

ε

λδ

σ

+

=

i, j, k = 1, 2, 3 (3.3)

a integral à direita em (3.2) pode ser escrita da seguinte maneira:

(

)

∫∫

+=

ΩΩ

ΩΩ d2d

*

ijijkkij

*

ijij

εμεελδεσ

1

BREBBIA, C. A., DOMINGUEZ, J. Boundary Elements, An Introductory Course. Boston:

Computational Mechanics Publications, 1989.

42

ou

(

)

∫∫

+=

ΩΩ

ΩΩ d2d

*

ijijkk

*

ijij

*

ijij

εμεεελδεσ

Desde que

*

kk

*

ijij

εεδ

=

e

ijijkk

ε

δ

ε

=

tem-se

(

)

∫∫

+=

ΩΩ

ΩΩ d2d

*

ijijijij

*

kk

*

ijij

εμεεδλεεσ

ou

(

)

∫∫

+=

ΩΩ

ΩΩ d2d

ij

*

ij

*

kkij

*

ijij

εμεελδεσ

e, novamente usando as equações constitutivas, tem-se

∫∫

=

ΩΩ

ΩΩ dd

ij

*

ij

*

ijij

εσεσ

(3.4)

Utilizando (3.2) no lado esquerdo de (3.4) e, por analogia de notação,

também no lado direito, a expressão resultante

∫∫∫∫

+=+

ΩΓΩΓ

ΩΓΩΓ dubdutdubdut

i

*

ii

*

i

*

ii

*

ii

(3.5)

é denominada teorema do trabalho recíproco de Betti.

3.2 Equações Integrais de Contorno

A equação integral de contorno para problemas de elastostática pode ser

derivada do teorema do trabalho recíproco de Betti tomando as forças de corpo

como sendo forças pontuais num meio infinito, representadas pela função Delta de

Dirac da seguinte maneira:

(

)

i

*

i

e,b XX'

Δ

=

onde o vetor unitário e

i

corresponde a um força unitária positiva na direção i aplicada

no ponto X', denominado ponto fonte. Em problemas bidimensionais, e

i

é uma força

por unidade de espessura e, em problemas tridimensionais, uma força concentrada.

43

Por sua vez, os deslocamentos e forças de superfície correspondentes à solução do

problema para o caso particular de uma força pontual são escritos como

(

)

jij

*

i

e,Uu XX'=

(

)

jij

*

i

e,Tt XX'=

onde U

ij

(X',X) e T

ij

(X',X) representam, respectivamente, as componentes de

deslocamento e forças de superfície na direção

j no ponto X devido a uma força

pontual unitária atuando na direção

i no ponto fonte X'.

Utilizando essas soluções na equação (3.5) obtém-se

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

∫∫∫

+−=

ΩΓΓ

ΩXXX'ΓxxX'ΓxxX'X' db,Udu,Tdt,Uu

jijjijjiji

(3.6)

onde X' e X ∈ Ω e x ∈ Γ, e os termos U

ij

e T

ij

são denominados soluções

fundamentais. A equação acima é conhecida como identidade de Somigliana para

deslocamentos. Ela relaciona o valor dos deslocamentos em um ponto do domínio Ω

com os valores de forças de superfície e de deslocamentos no contorno Γ. Em forma

matricial tem-se

∫∫∫

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

−

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

⎭

⎬

⎫

⎩

⎨

⎧

ΩΓΓ

ΩΓΓ d

b

UU

d

u

TT

d

t

UU

u

2

1

2221

1211

2

1

2221

1211

2

1

2221

1211

2

1

para problemas bidimensionais, e

+

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

−

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

∫∫

ΓΓ

ΓΓ d

u

TTT

d

t

UUU

u

3

2

1

333231

232221

131211

3

2

1

333231

232221

131211

3

2

1

∫

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

Ω

Ωd

b

UUU

3

2

1

333231

232221

131211

para problemas tridimensionais.

As deformações específicas num ponto qualquer do domínio Ω podem ser

calculadas por diferenciação da equação de deslocamentos (3.6) em relação ao

ponto fonte X' e substituição em (2.50), resultando

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

∫∫∫

+−=

ΩΓΓ

ΩXXX'ΓxxX'ΓxxX'X' db,Udu,Tdt,Uu

jk,ijjk,ijjk,ijk,i

(3.7)

onde U

ij,k

e T

ij,k

são derivadas das soluções fundamentais.

44

Finalmente, a identidade de Somigliana para tensões nos pontos do domínio

é encontrada substituindo a equação (3.7) nas equações constitutivas (3.3),

resultando em

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

∫∫∫

+−=

ΩΓΓ

ΩXXX'ΓxxX'ΓxxX'X' db,Ddu,Sdt,D

kkijkkijkkijij

σ

(3.8)

onde D

kij

e S

kij

são obtidas de U

ij,k

e T

ij,k

e da aplicação da lei de Hooke.

3.3 Soluções Fundamentais

A existência de soluções para forças pontuais desempenha um importante

papel na formulação do Método de Elementos de Contorno. É possível reescrever as

equações de Navier para uma força pontual unitária, aplicada ao corpo no ponto

fonte X', como

()

0e,u

21

u

i

*

ji,j

*

jj,i

=+

−

+ XX'

Δ

ν

μ

As soluções das equações governantes, devidas a uma força pontual, são

comumente chamadas de "soluções fundamentais". Existem várias maneiras para

obter-se essas soluções. A técnica mais popular consiste no uso do vetor de

Galerkin

2

; no entanto, métodos alternativos baseados em funções complexas de

tensões também são eficazes

3

. As soluções para uma força pontual num meio

infinito foram originalmente derivadas por Lorde Kelvin, e são conhecidas como

soluções fundamentais de Kelvin.

2

ALIABADI, M. H. The Boundary Element Method, Vol. 2, Applications in Solids and

Structures. Chippenham: John Wiley & Sons, 2002. p. 31-33.

3

Ibid., p. 525-529.

45

3.3.1 Problemas Tridimensionais

Para problemas tridimensionais, as soluções fundamentais são dadas por

4

()

[

]

jiijij

,r,r43

r116

1

,U +−

−

=

δν

νπμ

xX' (3.9)

onde U

ij

(X',x) representa o deslocamento na direção j no ponto x devido a uma força

pontual unitária atuando na direção i no ponto fonte X', e também por

()

[

]

()

⎭

⎬

⎫

⎩

⎨

⎧

−−++−

∂

−

=

ijjijiij

2

ij

,rn,rn21,r,r321

n

r

r18

1

,T

νδν

νπ

xX' (3.10)

onde T

ij

(X',x) representa a força de superfície na direção j no ponto x devido a uma

força pontual unitária atuando na direção i no ponto fonte X'. Diferenciando (3.9) e

(3.10) em relação a X' tem-se

()

[

]

jikijkkjikij

2

k,ij

,r,r,r,r,r3,r43

r116

1

,U

δδδν

νπμ

−−+−

−

=xX'

(3.11)

()

[

]

{

++−

−

=

kjikij

3

k,ij

n,r,r3,r21

r18

1

,T

δν

νπ

xX'

()

(

)

(

)

(

)

+

−

+

−

+

jkiikikjjk

n,r,r321n,r,r321

δ

ν

δ

ν

()

[

]

}

kllkjikijjikijk

,,n,r,r,r,r15,r213,r3,r3 −−−++

δνδδ

(3.12)

Usando as relações tensão-deformação, resulta

()

(

)

[

]

kjikijijkjik

2

kij

,r,r,r3,r,r,r21

r18

1

,D +−+−

−

=

δδδν

νπ

xX'

(3.13)

()

[

]

+

⎩

⎨

⎧

−++−

∂

−

=

kjiijkjikkij

3

kij

,r,r,r5,r,r,r21

n

r

3

r14

,S

δδνδν

νπ

μ

xX'

(

)

(

)

(

)

(

)

}

kijijkjikjikkijkji

n41nn,r,rn321,r,rn,r,rn3

δ

ν

δ

ν

−

+

−

++ (3.14)

4

Ibid., p. 31,32.

46

3.3.2 Problemas Bidimensionais

Para problemas bidimensionais, as soluções fundamentais são dadas por

5

()

⎥

⎦

⎤

⎢

⎣

⎡

+

⎟

⎠

⎞

⎜

⎝

⎛

−

=

jiijij

,r,r

r

1

ln43

18

1

,U

δν

νπμ

xX' (3.15)

onde U

ij

(X',x) representa o deslocamento na direção j no ponto x devido a uma força

pontual unitária atuando na direção i no ponto fonte X', e também por

()

[

]

()

⎭

⎬

⎫

⎩

⎨

⎧

−−++−

∂

−

=

ijjijiijij

,rn,rn21,r,r221

n

r

r14

1

,T

νδν

νπ

xX' (3.16)

onde T

ij

(X',x) representa a força de superfície na direção j no ponto x devido a uma

força pontual unitária atuando na direção i ponto fonte X'. Diferenciando (3.15) e

(3.16) em relação a X' tem-se

()

[

]

jikijkkjikijk,ij

,r,r,r,r,r2,r43

rE14

1

,U

δδδν

νπ

ν

−−+−

−

+−

=

xX' (3.17)

()

[

]

+

⎩

⎨

⎧

−−−+

∂

−

=

kjikijijkjik

2

k,ij

,r,r,r4,r21,r,r

n

r

2

r14

1

,T

δνδδ

νπ

xX'

()

(

)

(

)

(

)

+

−

+

jkiikikjjk

n,r,r221n,r,r221

δ

ν

δ

ν

()

[]

⎭

⎬

⎫

+−+

kjiij

n,r,r221

δν

(3.18)

Usando as relações tensão-deformação, resulta

()

(

)

[

]

kjikijijkjikkij

,r,r,r2,r,r,r21

r14

1

,D +−+−

−

=

δδδν

νπ

xX'

(3.19)

()

[

]

+

⎩

⎨

⎧

−++−

∂

−

=

kjiijkjikkij

2

kij

,r,r,r4,r,r,r21

n

r

2

r12

,S

δδνδν

νπ

μ

xX'

(

)

(

)

(

)

()

}

kijijkjikjikkijkji

n41nn,r,rn221,r,rn,r,rn2

δ

ν

δ

ν

−

+

−

+

++ (3.20)

As soluções fundamentais, anteriormente apresentadas, são válidas para a

condição de Estado Plano de Deformações. As soluções para a condição de Estado

Plano de Tensões podem ser obtidas a partir dessas mesmas soluções através da

5

Ibid., p. 33.

47

introdução do coeficiente de Poisson ν' e do módulo de Young E' modificados

definidos como

ν

+

=

1

'

⎥

⎦

⎤

⎢

⎣

⎡

⎟

⎠

⎞

⎜

⎝

⎛

+

−=

2

1

1E'E

ν

(3.21)

3.4 Equação Integral de Contorno para Deslocamentos no Contorno

A equação integral de contorno (3.6) é válida para qualquer ponto fonte X'

localizado dentro do domínio Ω do corpo; porém, para que possa ser utilizada,

necessita dos valores dos deslocamentos e forças de superfície nos pontos

pertencentes ao contorno. É possível calcular tais valores avaliando-se o limite

dessa expressão quando X' ∈ Ω → x'

∈

Γ.

Para tal, considere-se que o domínio do problema tenha sido aumentado,

em torno do ponto fonte x', de uma região circular, no caso de problemas

bidimensionais, ou de uma região esférica, no caso de problemas tridimensionais,

com contorno

+

ε

Γ , raio ε e centro em x', como ilustrado nas Figuras 3.2 e 3.3,

respectivamente.

O contorno aumentado

A

Γ é representado por

(

)

+

+−=

εε

ΓΓΓΓ

A

(3.22)

onde Γ

ε

é a porção do contorno original que foi removida. O comportamento das

soluções fundamentais pode ser estudado tomando o limite quando ε → 0 e,

portanto,

Γ

A

→ Γ.

A primeira integral de contorno no lado direito de

(3.6) é reescrita como

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

∫∫∫

+

→

−

→

+=

ε

ΓΓΓΓ

Γxxx'Γxxx'ΓxxX' dt,Ulimdt,Ulimdt,U

jij

0

jij

0

jij

(3.23)

onde ambas as expressões no lado direito contém integrandos singulares de ordem

O(ln 1/r) em duas dimensões e O(

1

r

−

) em três dimensões.

48

FIGURA 3.2 – Processo limite quando o ponto fonte está localizado no contorno de um

problema bidimensional

FIGURA 3.3 – Processo limite quando o ponto fonte está localizado no contorno de um

problema tridimensional

49

A primeira integral no lado direito de (3.23) é integrável na forma de uma

integral imprópria, pois o tipo de singularidade presente é o que comumente

denomina-se "singularidade fraca". Como será mostrado mais adiante neste

trabalho, para problemas bidimensionais, ou seja, problemas com singularidade de

ordem O(ln 1/r), é possível dividir essa integral em duas partes: uma singular e outra

não singular. A parte não singular é integrada utilizando a quadratura de Gauss

convencional, apresentada no capítulo anterior. A parte singular é integrada

utilizando uma versão modificada da referida quadratura, denominada quadratura de

Gauss logarítmica (para valores calculados e tabelados de pontos e pesos de Gauss

para a quadratura Gaussiana logarítmica ver, por exemplo, BREBBIA e

DOMINGUEZ

6

). Já a segunda integral no lado direito de (3.23) tende a zero à

medida que ε → 0. Esse resultado pode ser obtido empregando coordenadas

polares, para o caso bidimensional, e coordenadas esféricas, para o caso

tridimensional, ao avaliar-se a integral sobre

+

ε

Γ à medida que ε → 0.

A segunda integral de contorno no lado direito de (3.6) pode ser reescrita

como

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

∫∫∫

+

→

−

→

+=

ε

ΓΓΓΓ

Γxxx'Γxxx'ΓxxX' du,Tlimdu,Tlimdu,T

jij

0

jij

0

jij

(3.24)

onde ambas as expressões no lado direito contém integrandos singulares de ordem

O(

1

r

−

) em duas dimensões e O(

2

r

−

) em três dimensões.

A primeira integral é calculada no sentido do Valor Principal de Cauchy

7

. A

segunda integral é regularizada pelo primeiro termo da expansão dos

deslocamentos em série de Taylor

8

em torno do ponto fonte, processo que conduz a

()()

=

∫

+

→

ε

Γ

Γxxx' du,Tlim

jij

0

6

BREBBIA, C. A., DOMINGUEZ, J. Boundary Elements, An Introductory Course. Boston:

Computational Mechanics Publications, 1989, p. 271.

7

GUIGGIANI, M. The evaluation of Cauchy principal value integrals in the Boundary

Element Method – A review. Mathematical Computational Modelling, vol. 15, n. 3-5, p. 175-184,

1991 e GUIGGIANI, M., CASALINI, P. Direct computation of Cauchy principal value integrals in

advanced boundary elements. International Journal for Numerical Methods in Engineering, vol.

24, p. 1711-1720, 1987.

8

Os deslocamentos são assumidos como satisfazendo a condição de continuidade de

Hölder, ou seja, u

∈ C

0,

α

se existem constantes k > 0 e 0 <

α

≤ 1 tal que | u(x') – u(x) | ≤ kr

α

.

50

()() ()

[

]

() ( )

∫∫

+

→

+

→

+

⎭

⎬

⎫

⎩

⎨

⎧

−=

ε

ΓΓ

Γxx'x'Γx'xxx' d,Tlimuduu,Tlim

ij

0

jjjij

0

(3.25)

onde a primeira integral no lado direito deve ser igual a zero pelo requisito de

continuidade dos deslocamentos. A segunda integral conduz a um coeficiente

denotado por

() ( ) () ()

x'x'Γxx'x'

Γ

jijij

0

j

ud,Tlimu

α

ε

=

∫

+

→

(3.26)

e, desse modo,

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

x'x'Γxxx'ΓxxX'

ΓΓ

jijjijjij

udu,TVPCdu,T

α

+=

∫∫

(3.27)

onde

()

∫

Γ

Γd...VPC

indica uma integral avaliada no sentido do Valor Principal de

Cauchy.

Finalmente, a equação integral de contorno para deslocamentos no contorno

reduz-se a

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

+=+

∫∫

ΓΓ

Γxxx'Γxxx'x'x' dt,Udu,TVPCuC

jijjijjij

()()

∫

+

Ω

ΩXXx' db,U

jij

(3.28)

onde C

ij

(x')u

j

(x') = [δ

ij

+ α

ij

(x')]u

j

(x') e é denominado termo livre.

Para contornos bidimensionais "suaves", como o apresentado na Figura 3.2,

a avaliação dos coeficientes α

ij

(x') pode ser obtida analiticamente. Denotando por r a

distância entre o ponto fonte x' e um ponto x no contorno

+

ε

Γ , e usando coordenadas

polares, resulta

jijir

θ

ε

θ

ε

θ

sincossinrcosr

+

=

+

=

θεθ

ε

ddrd ==

+

Γ

1

r

rn

r

===

•

=

•

=

∂

ε

nnnr

θ

cosn

1

=

θ

sinn

2

=

θ

cos,r

1

=

θ

sin,r

2

=

Utilizando essas relações, as constantes α

ij

(x') são avaliadas como

() ( )

()

[

]

=+−

−

==

∫∫

→

+

→

π

ε

θεθν

ενπ

α

0

2

0

11

0

11

dcos221

1

lim

14

1

d,Tlim

Γ

Γxx'x'

51

()

[

]

()

=

⎥

⎦

⎤

⎢

⎣

⎡

+−

−

=+−

−

=

∫∫∫

πππ

θθθν

νπ

θθν

νπ

0

2

00

2

dcos2d21

14

1

dcos221

14

1

()

=

⎥

⎦

⎤

⎢

⎣

⎡

⎟

⎠

⎞

⎜

⎝

⎛

++−

−

=

∫∫

ππ

θθθν

νπ

00

d2cos

2

1

2

1

2d21

14

1

()

=

⎥

⎦

⎤

⎢

⎣

⎡

++−

−

=

∫∫∫

πππ

θθθθν

νπ

000

d2cosdd21

14

1

()

[]

(

)

()

(

)

()

2

1

14

12

14

2

021

14

1

−=

−

−=

−

+

−

−=++−

−

=

νπ

ν

π

νπ

π

ν

π

ππν

νπ

(i)

() ( )

()

=

−

==

∫∫

→

+

→

π

ε

θεθθ

ενπ

α

0

12

0

12

dsincos2

1

lim

14

1

d,Tlim

Γ

Γxx'x'

() ()

=

−

=

−

=

∫∫

ππ

θθθ

νπ

θθθ

νπ

00

dsincos

14

2

dsincos2

14

1

() ()

0d2sin

14

1

d2sin

2

1

12

1

00

=

−

=

−

=

∫∫

ππ

θθ

νπ

θθ

νπ

(ii)

() ( )

()

=

−

==

∫∫

→

+

→

π

ε

θεθθ

ενπ

α

0

21

0

21

dcossin2

1

lim

14

1

d,Tlim

Γ

Γxx'x'

() ()

=

−

=

−

=

∫∫

ππ

θθθ

νπ

θθθ

νπ

00

dcossin

14

2

dcossin2

14

1

() ()

0d2sin

14

1

d2sin

2

1

12

1

00

=

−

=

−

=

∫∫

ππ

θθ

νπ

θθ

νπ

(iii)

() ( )

()

[

]

=+−

−

==

∫∫

→

+

→

π

ε

θεθν

ενπ

α

0

2

0

22

0

22

dsin221

1

lim

14

1

d,Tlim

Γ

Γxx'x'

()

[

]

()

=

⎥

⎦

⎤

⎢

⎣

⎡

+−

−

=+−

−

=

∫∫∫

πππ

θθθν

νπ

θθν

νπ

0

2

00

2

dsin2d21

14

1

dsin221

14

1

()

=

⎥

⎦

⎤

⎢

⎣

⎡

⎟

⎠

⎞

⎜

⎝

⎛

−+−

−

=

∫∫

ππ

θθθν

νπ

00

d2cos

2

1

2

1

2d21

14

1

()

=

⎥

⎦

⎤

⎢

⎣

⎡

−+−

−

=

∫∫∫

πππ

θθθθν

νπ

000

d2cosdd21

14

1

()

[]

(

)

()

(

)

()

2

1

14

12

14

2

021

14

1

−=

−

−=

−

+

−

−=−+−

−

=

νπ

ν

π

νπ

π

ν

π

ππν

νπ

(iv)

Conclusão:

()

2

ij

δ

α

−=x' (3.29)

e

52

() ()

2

C

ij

ijijij

δ

αδ

=+= x'x' (3.30)

Para contornos tridimensionais "suaves", como o apresentado na Figura 3.3,

o processo é análogo (utilizando, no entanto, coordenadas esféricas em lugar de

coordenadas polares) e o resultado é idêntico ao encontrado para contornos

bidimensionais.

Para contornos "não-suaves" os termos livres podem ser avaliados

indiretamente a partir de uma consideração de movimento de corpo rígido, que será

explicada posteriormente, em momento oportuno.

3.5 Regiões Infinitas

A extensão da equação (3.28) para regiões infinitas deve ser feita levando-se

em consideração hipóteses adicionais relativas às funções envolvidas. Essas

hipóteses estão associadas ao comportamento das funções em uma superfície

infinitamente distante do ponto fonte e são denominadas condições de regularidade.

Seja ρ o raio de uma esfera com contorno Γ

ρ

, centrada no ponto fonte x', que

envolve a cavidade, ou as cavidades, do problema, como ilustrado na Figura 3.4.

A equação

(3.28) é reescrita para a região entre Γ e Γ

ρ

(considerando as

forças de corpo nulas, para simplificar) como

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

=++

∫∫

ρ

ΓΓ

Γxxx'Γxxx'x'x' du,Tdu,TuC

jijjijjij

( ) () ( ) ()

∫∫

+=

ρ

ΓΓ

Γxxx'Γxxx' dt,Udt,U

jijjij

Tomando o limite quando ρ → ∞, a equação anterior pode ser escrita apenas

em termos de integrais sobre o contorno Γ se a condição de regularidade seguinte

for satisfeita:

( ) () ( ) ()

[]

0dt,Uu,Tlim

jijjij

=

⎭

⎬

⎫

⎩

⎨

⎧

−

∫

∞→

ρ

Γ

Γxxx'xxx' (3.31)

53

FIGURA 3.4 – Região infinita envolvendo uma cavidade

Para problemas tridimensionais, tem-se

(

)

(

)

1

ij

O,U

−

=

ρ

xx'

(

)

(

)

2

ij

O,T

−

=

ρ

xx'

θϕ

ρ

ddd JΓ =

(

)

2

O

ρ

=J

portanto, se na pior hipótese,

(

)

(

)

1

j

Ou

−

=

ρ

x

e

(

)

(

)

2

j

Ot

−

=

ρ

x

no infinito, a condição de regularidade (3.31) é satisfeita. Deve-se observar que se a

carga total aplicada na superfície Γ não for auto-equilibrada, o princípio de Saint-

Venant

9

mostra que u

j

(x) e t

j

(x) terão o mesmo comportamento que a solução

fundamental correspondente a uma carga concentrada na direção da resultante.

Desse modo,

9

MALVERN, L. E. Introduction to the Mechanics of a Continuous Medium. Englewood

Cliffs, N. J.: Prentice-Hall, 1969.

54

(

)

(

)

1

j

Ou

−

=

ρ

x

e

(

)

(

)

2

j

Ot

−

=

ρ

x

são obtidos e cada termo da condição de regularidade se anula separadamente.

Para problemas bidimensionais, tem-se

()

(

)

(

)

()

⎩

⎨

⎧

≠

=+

=

jise1O

jise1lnO

,U

ij

ρ

xx'

(

)

(

)

1

ij

O,T

−

=

ρ

xx'

θ

ρ

dd JΓ =

(

)

ρ

O=J

portanto, para garantir que cada termo da condição de regularidade (3.31) se anule

separadamente, deve-se ter

(

)

(

)

1

j

Ou

−

=

ρ

x

e

(

)

(

)

2

j

Ot

−

=

ρ

x

como antes. Esse caso, no entanto, não corresponde ao comportamento da solução

fundamental no infinito. Baseado no mesmo argumento do caso tridimensional,

substituindo u

j

(x) e t

j

(x) pelos tensores correspondentes à solução fundamental

bidimensional, verifica-se que a condição de regularidade é satisfeita. A única

diferença é que, agora, os termos não se cancelam separadamente, mas se

cancelam quando ρ → ∞.

A discussão anterior permite concluir que a condição de regularidade (3.31) é

sempre satisfeita se u

j

(x) e t

j

(x) comportam-se, na pior hipótese, como as

correspondentes soluções fundamentais no infinito. Desse modo, problemas de

cavidades em meios infinitos podem ser representados pela equação (3.28), desde

que a normal ao contorno Γ e o sentido de integração sejam definidos conforme a

Figura 3.5.

55

FIGURA 3.5 – Definição da normal e do sentido de integração: a) Problemas com domínio

infinito; b) Problemas com domínio finito

3.6 Discretização Numérica

As equações integrais de contorno, apresentadas até agora, só podem ser

resolvidas analiticamente para problemas muito simples. Para problemas

constituídos de geometrias e carregamentos mais complexos, somente a adoção de

métodos numéricos permite solucioná-las. O Método de Elementos de Contorno é

um método numérico para solução dessas equações integrais de contorno partindo

de um procedimento de discretização.

O primeiro passo nesse procedimento de discretização é dividir, ou

aproximar, o contorno Γ em Ne segmentos ou elementos. Desse modo, na ausência

de forças de corpo, para simplificar, a equação (3.28) pode ser reescrita na forma

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

∑

∫

∑

∫

==

=+

Ne

1n

n

njij

Ne

1n

n

njijjij

dt,Udu,TuC

ΓΓ

Γxxx'Γxxx'x'x'

(3.32)

onde

∑

=

Ne

1n

n

ΓΓ .

56

3.7 Elementos Constantes

Nesse trabalho, a geometria do contorno é dividida em Ne elementos retos, e

os valores das incógnitas e das condições de contorno u

j

(x) e t

j

(x) são considerados

constantes ao longo do comprimento de cada elemento e iguais aos valores dessas

mesmas incógnitas nos pontos centrais (também denominados pontos nodais ou

nós funcionais) dos respectivos elementos. Partindo dessas considerações,

reescreve-se a equação (3.32) como

∑

∫

∑

∫

==

=+

Ne

1n

n

c

ij

n

j

Ne

1n

n

c

ij

n

j

c

j

c

ij

dUtdTuuC

ΓΓ

ΓΓ (3.33)

onde

n

j

u e

n

j

t são os valores dos deslocamentos e forças de superfície nos pontos

centrais dos elementos de contorno;

c

j

u ,

c

j

c

ij

uC ,

c

ij

T e

c

ij

U são, respectivamente, os

valores dos deslocamentos, termos livres e soluções fundamentais para o ponto

fonte x' = x

c

.

Fazendo os pontos fonte percorrerem, sucessivamente, os pontos centrais

dos elementos de contorno, isto é, fazendo c = 1, 2, ..., Ne, são geradas Ne equações

para Ne valores de deslocamentos e Ne valores de forças de superfície. Aplicando as

condições de contorno do problema, ou seja, deslocamentos e forças de superfície

prescritos para os pontos nodais, produz-se um sistema algébrico que, ao ser

resolvido, resulta nos valores das incógnitas do problema. Outra observação

pertinente: todo ponto fonte x

c

está localizado em um contorno "suave", portanto, os

coeficientes associados aos termos livres podem ser calculados através da

expressão

(3.30).

3.7.1 Formação do Sistema de Equações

Para ilustrar a aplicação da equação

(3.33) far-se-á uso da malha modelo

apresentada na Figura 3.6. A seguinte convenção é estabelecida quando tratar-se

da representação geométrica de malhas de elementos de contorno: os elementos de

57

contorno são representados por linhas azuis delimitadas por pequenos traços

vermelhos e sua numeração correspondente é aquela dentro dos círculos verdes; os

pontos azuis, sobre os elementos de contorno, são os nós de contorno e sua

numeração é aquela localizada próxima aos mesmos; caso existam pontos no

domínio do problema (também chamados de pontos internos), os mesmos

apresentarão cor vermelha e numeração independente daquela dos elementos ou

nós de contorno.

FIGURA 3.6 – Malha modelo formada por elementos de contorno constantes

Aplicando a equação (3.33), quando o ponto fonte é o Nó 1, resulta

++++++

∫∫∫∫

2

1

12

2

1

2

1

11

1

2

1

12

1

11

1

2

1

12

1

11

udTudTudTudTuCuC

ΓΓΓΓ

=++++

∫∫∫∫

4

2

4

1

12

4

1

4

1

11

3

2

3

1

12

3

1

3

1

11

udTudTudTudT

ΓΓΓΓ

++++=

∫∫∫∫

2

1

12

2

1

2

1

11

1

2

1

12

1

11

tdUtdUtdUtdU

ΓΓΓΓ

4

2

4

1

12

4

1

4

1

11

3

2

3

1

12

3

1

3

1

11

tdUtdUtdUtdU

∫∫∫∫

++++

ΓΓΓΓ

ΓΓΓΓ (i)

++++++

∫∫∫∫

2

1

22

2

1

2

1

21

1

2

1

22

1

21

1

2

1

22

1

21

udTudTudTudTuCuC

ΓΓΓΓ

58

=++++

∫∫∫∫

4

2

4

1

22

4

1

4

1

21

3

2

3

1

22

3

1

3

1

21

udTudTudTudT

ΓΓΓΓ

++++=

∫∫∫∫

2

1

22

2

1

2

1

21

1

2

1

22

1

21

tdUtdUtdUtdU

ΓΓΓΓ

4

2

4

1

22

4

1

4

1

21

3

2

3

1

22

3

1

3

1

21

tdUtdUtdUtdU

∫∫∫∫

++++

ΓΓΓΓ

ΓΓΓΓ (ii)

A expressão (i) corresponde a uma força unitária atuando no ponto fonte na

direção

1; a expressão (ii) corresponde a uma força unitária atuando no ponto fonte

na direção

2.

Em forma matricial, essas expressões são agrupadas do seguinte modo:

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

∫∫

2

1

2

1

22

1

21

1

12

1

11

1

2

1

22

1

21

1

12

1

11

1

2

1

22

1

21

1

12

1

11

u

d

TT

u

d

TT

u

CC

ΓΓ

=

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

∫∫

4

2

4

1

4

1

22

1

21

1

12

1

11

3

2

3

1

3

1

22

1

21

1

12

1

11

u

d

TT

u

d

TT

ΓΓ

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

∫∫

2

1

2

1

22

1

21

1

12

1

11

1

2

1

22

1

21

1

12

1

11

t

d

UU

t

d

UU

ΓΓ

∫∫

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

4

2

4

1

4

1

22

1

21

1

12

1

11

3

2

3

1

3

1

22

1

21

1

12

1

11

t

d

UU

t

d

UU

ΓΓ

ou, de forma mais compacta,

∑

∫

∑

∫

=

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

4Ne

1n

n

2

n

1

n

1

22

1

21

1

12

1

11

4Ne

1n

n

2

n

1

n

1

22

1

21

1

12

1

11

1

2

1

22

1

21

1

12

1

11

t

d

UU

u

d

TT

u

CC

ΓΓ

De maneira análoga, para o Nó 2,

∑

∫

∑

∫

=

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

4Ne

1n

n

2

n

1

n

2

22

2

21

2

12

2

11

4Ne

1n

n

2

n

1

n

2

22

2

21

2

12

2

11

2

1

2

22

2

21

2

12

2

11

t

d

UU

u

d

TT

u

CC

ΓΓ

para o Nó 3,

∑

∫

∑

∫

=

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

4Ne

1n

n

2

n

1

n

3

22

3

21

3

12

3

11

4Ne

1n

n

2

n

1

n

3

22

3

21

3

12

3

11

3

2

3

1

3

22

3

21

3

12

3

11

t

d

UU

u

d

TT

u

CC

ΓΓ

e para o Nó 4,

∑

∫

∑

∫

=

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

4Ne

1n

n

2

n

1

n

4

22

4

21

4

12

4

11

4Ne

1n

n

2

n

1

n

4

22

4

21

4

12

4

11

4

2

4

1

4

22

4

21

4

12

4

11

t

d

UU

u

d

TT

u

CC

ΓΓ

Ao adotar-se a notação

59

⎥

⎦

⎤

⎢

⎣

⎡

=

n

22

n

21

n

12

n

11

n

CC

C

⎭

⎬

⎫

⎩

⎨

⎧

=

n

2

n

1

n

u

⎭

⎬

⎫

⎩

⎨

⎧

=

n

2

n

1

n

t

⎥

⎦

⎤

⎢

⎣

⎡

=

n

22

n

21

n

12

n

11

n

UU

U

⎥

⎦

⎤

⎢

⎣

⎡

=

n

22

n

21

n

12

n

11

n

TT

T

as equações anteriores passam a ser reescritas como

∑

∫

∑

∫

=

=+

4Ne

1n

n

1

4Ne

1n

n

111

dd

ΓΓ

tΓUuΓTuC

∑

∫

∑

∫

=

=+

4Ne

1n

n

2

4Ne

1n

n

222

dd

ΓΓ

tΓUuΓTuC

∑

∫

∑

∫

=

=+

4Ne

1n

n

3

4Ne

1n

n

333

dd

ΓΓ

tΓUuΓTuC

∑

∫

∑

∫

=

=+

4Ne

1n

n

4

4Ne

1n

n

444

dd

ΓΓ

tΓUuΓTuC

ou, de maneira equivalente,

∑

∫

∑

∫

=

=+

4Ne

1n

n

c

4Ne

1n

n

ccc

dd

ΓΓ

tΓUuΓTuC c = 1, 2, ..., Ne = 4

Denotando as integrais das soluções fundamentais por

∫

=

n

ccn

d

Γ

ΓTH

∫

=

n

ccn

d

Γ

ΓUG

tem-se

∑∑

=

=+

4Ne

1n

ncn

4Ne

1n

ncncc

tGuHuC c = 1, 2, ..., Ne = 4

Definindo

⎪

⎩

⎪

⎨

⎧

=+=

≠=

ncquando

ccncn

cncn

CHH

HH

resulta

∑∑

=

4Ne

1n

ncn

4Ne

1n

ncn

tGuH c = 1, 2, ..., Ne = 4

A equação anterior, reescrita em forma matricial, produz o sistema algébrico

tGuH

=

o qual, quando reordenado através da aplicação das condições de contorno

prescritas para o problema, produz o sistema algébrico final

FYA

=

60

onde o vetor Y contém todas as incógnitas (deslocamentos ou forças de superfície)

do problema; a matriz de coeficientes A é formada de colunas de H e G; o vetor F é

formado de colunas de H e G multiplicadas pelas condições de contorno prescritas.

Uma vez que esse sistema esteja resolvido, todos os valores das incógnitas de

contorno serão conhecidos.

O processo de introdução das condições de contorno no sistema algébrico

Hu = Gt, transformando-o, dessa maneira, no sistema algébrico final AY = F, será

ilustrado através das Figuras 3.7 e 3.8.

A Figura 3.7 apresenta o sistema algébrico Hu = Gt resultante para a malha

modelo da Figura 3.6, destacando, na cor azul, os valores dos deslocamentos

prescritos e, na cor vermelha, os valores das forças de superfície prescritas para um

problema arbitrário. Também estão destacadas as colunas das matrizes H e G que

multiplicam os valores das condições de contorno prescritas.

FIGURA 3.7 – Sistema algébrico tGuH

=

destacando os valores das condições de contorno

prescritas para um problema arbitrário

Para agrupar todas as incógnitas do problema no lado esquerdo da

igualdade e todos os valores das condições de contorno prescritas no lado direito da

igualdade, basta trocar as colunas das matrizes H e G da forma indicada na Figura

61

3.8 e inverter o sinal dos coeficientes. Agora, o lado direito da igualdade possui

apenas valores conhecidos, o que dá origem ao vetor F, enquanto que o lado

esquerdo da igualdade possui a matriz de coeficientes A que multiplica o vetor de

incógnitas Y, isto é, o sistema algébrico AY = F.

FIGURA 3.8 – Sistema algébrico FYA

=

resultante da reordenação das colunas das matrizes

H

e G de modo a agrupar todas as incógnitas do problema no lado direito da igualdade e os

valores das condições de contorno prescritas no lado esquerdo da igualdade

3.7.2 Integração

Todas as integrais, nas expressões anteriores, podem ser avaliadas

numericamente. Para o caso específico de elementos constantes, entretanto, é mais

simples e mais exato calcular algumas dessas integrais analiticamente,

particularmente aquelas sobre o elemento que apresenta singularidade, isto é,

quando

c = n. Nos demais elementos, os valores das integrais em

cn

H e

cn

G são

computados utilizando-se o processo de quadratura de Gauss.

62

3.7.2.1 Submatrizes

cc

H

As submatrizes

cc

H , que formam a diagonal de submatrizes de H, são

dadas por

⎥

⎦

⎤

⎢

⎣

⎡

=

⎥

⎦

⎤

⎢

⎣

⎡

=

2/10

02/1

HH

cc

22

cc

21

cc

12

cc

11

cc

H

Esse resultado é obtido considerando-se a relação

10

ρ

1

r

nr

lim

0r

=

∂

→

onde ρ denota o raio de curvatura do elemento. Para elementos com geometria reta,

ρ = ∞, logo

0

11

r

nr

lim

0r

=

∞

==

∂

→

ρ

Quando essa última relação é utilizada na avaliação da integral da solução

fundamental apresentada em (3.16), em torno do ponto fonte x

c

e ao longo do

elemento que contém o ponto fonte, verifica-se que essa integral é nula. Desse

modo, as submatrizes

cc

H são dadas por

⎥

⎦

⎤

⎢

⎣

⎡

+

⎥

⎦

⎤

⎢

⎣

⎡

=

⎥

⎦

⎤

⎢

⎣

⎡

=

c

22

c

21

c

12

c

11

cc

22

cc

21

cc

12

cc

11

cc

CC

00

HH

H

Como mencionado anteriormente, em elementos constantes com geometria

reta, todo ponto fonte

x

c

está localizado em um contorno "suave", portanto, os

coeficientes

c

ij

C podem ser calculados através da expressão (3.30), conduzindo ao

resultado apresentado inicialmente.

10

MANSUR, W. J., TELLES, J. C. F., PRODANOFF, J. H. A., FRAUCHES, E. On BEM

singular integrals for two-dimensional potencial applications. Engeneering Analysis with Boundary

Elements, vol. 9, p. 185-187, 1992.

63

3.7.2.2 Submatrizes

cc

G

As submatrizes

cc

G , que formam a diagonal de submatrizes de G,

calculadas analiticamente a partir da solução fundamental apresentada em (3.15),

em torno do ponto fonte x

c

e ao longo do elemento que contém o ponto fonte,

resultam em

11

⎥

⎦

⎤

⎢

⎣

⎡

=

cc

22

cc

21

cc

12

cc

11

cc

GG

G

sendo que

()

()( )

(

)

⎥

⎦

⎤

⎢

⎣

⎡

+−−

−

=

2

1

cc

11

R4

r

Rln143

18

R2

G

ν

νπμ

()

⎥

⎦

⎤

⎢

⎣

⎡

−

==

2

21

cc

21

cc

12

R4

rr

18

R2

GG

νπμ

()

()( )

(

)

⎥

⎦

⎤

⎢

⎣

⎡

+−−

−

=

2

cc

22

R4

r

Rln143

18

R2

G

ν

νπμ

onde R, r

1

e r

2

são definidos, com a ajuda da Figura 3.9, como

1

2

111

xxrr −==

r

1

2

222

xxrr −==

r

() ()

2

1

rr

2

1

RR +==

r

11

BREBBIA, C. A., DOMINGUEZ, J. Boundary Elements, An Introductory Course. Boston:

Computational Mechanics Publications, 1989, p. 199-202.

64

FIGURA 3.9 – Representação geométrica dos termos necessários para o cálculo dos

coeficientes das submatrizes

cc

G

3.7.2.3 Submatrizes

cn

G

quando nc

≠

As submatrizes

cn

G quando c ≠ n (as submatrizes de G posicionadas fora da

diagonal de submatrizes) são calculadas através da integração numérica da solução

fundamental apresentada em (3.15). Para tanto é preciso definir os termos r, r,

1

e r,

2

.

Usando a Figura 3.10, onde

g denota um ponto de integração de Gauss, conclui-se

que

c

1

g

1

11

xxrr −==

r

c

2

g

2

22

xxrr −==

r

() ()

(

)

(

)

2

c

2

g

2

c

1

g

1

2

1

xxxxrrrr −+−=+==

r

()()

2

c

2

g

2

c

1

g

1

c

1

g

1

xxxx

xx

r

,r

−+−

−

==

r

65

()()

2

c

2

g

2

c

1

g

1

c

2

g

2

xxxx

xx

r

,r

−+−

−

==

r

FIGURA 3.10 – Representação geométrica dos termos necessários para o cálculo dos

coeficientes das submatrizes

cn

G quando nc

≠

3.7.2.4 Submatrizes

cn

H quando nc

≠

As submatrizes

cn

H quando c ≠ n (as submatrizes de H posicionadas fora da

diagonal de submatrizes) são calculadas através da integração numérica da solução

fundamental apresentada em

(3.16). Para tanto é preciso definir, além dos termos já

definidos na seção anterior, os termos ∂r/∂n, n

1

e n

2

. Usando a Figura 3.11 conclui-se

que

(

)

θθθθαα

sinº90sinsinº90coscosº90coscoscosnnn

11

=+=−====

r

(

)

θθθθαα

cosº90sincosº90cossinº90sinsinsinnnn

22

−=−=−====

r

66

FIGURA 3.11 – Representação geométrica dos termos necessários para o cálculo dos

coeficientes das submatrizes

cn

H quando nc

≠

O ângulo θ representa a inclinação do elemento em relação ao eixo x

1

, isto é,

()()

2

1

2

1

2

1

2

1

xxxx

xx

cos

−+−

−

=

θ

()()

2

1

2

1

2

1

2

xxxx

xx

sin

−+−

−

=

θ

e assim

()()

2

1

2

1

2

1

2

1

xxxx

xx

n

−+−

−

=

()()

2

1

2

1

2

1

2

1

2

xxxx

xx

n

−+−

−

=

Finalmente,

2211

n,rn,r

n

r

+=

∂

()()

+

−+−

−

−+−

−

=

∂

2

1

2

1

2

1

2

c

2

g

2

c

1

g

1

c

1

g

1

xxxx

xx

xxxx

xx

n

r

67

()()

2

1

2

1

2

1

2

1

2

c

2

g

2

c

1

g

1

c

2

g

2

xxxx

xx

xxxx

xx

−+−

−

−+−

−

+

3.7.3 Deslocamentos e Tensões nos Pontos Internos

Os deslocamentos e tensões nos pontos pertencentes ao domínio Ω do

problema podem, uma vez que os valores das incógnitas de contorno tenham sido

obtidos, ser calculados utilizando-se as identidades de Somigliana (3.6) e (3.8),

respectivamente.

Para X

d

∈

Ω, a expressão (3.6) será reescrita, na ausência de forças de

corpo, como

()

(

)

()

(

)

()

∫∫

−=

ΓΓ

ΓxxXΓxxXX du,Tdt,Uu

j

d

ijj

d

ij

d

i

e, de modo equivalente, para o contorno discretizado, na forma

∑

∫

∑

∫

==

−=

Ne

1n

n

d

ij

n

j

Ne

1n

n

d

ij

n

j

d

i

dTudUtu

ΓΓ

Desenvolvendo a notação indicial e reagrupando matricialmente, de maneira

análoga ao que foi feito anteriormente, obtém-se

∑

∫

∑

∫

==

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

−

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

⎭

⎬

⎫

⎩

⎨

⎧

Ne

1n

n

2

n

1

n

d

22

d

21

d

12

d

11

Ne

1n

n

2

n

1

n

d

22

d

21

d

12

d

11

d

2

d

1

u

d

TT

t

d

UU

u

ΓΓ

Adotando as mesmas notações já utilizadas, resulta

∑

∫

∑

∫

==

−=

Ne

1n

n

d

Ne

1n

n

dd

ΓΓ

uΓTtΓUu

ou

∑∑

==

−=

Ne

1n

ndn

Ne

1n

ndnd

uHtGu

Essa equação relaciona os deslocamentos em um ponto qualquer do

domínio com os valores dos deslocamentos e forças de superfície de todo o

contorno, e pode ser facilmente avaliada dado que as soluções fundamentais não

apresentam singularidades.

68

Considerando o mesmo ponto X

d

∈

Ω, a expressão (3.8) será reescrita, na

ausência de forças de corpo, para simplificar, como

()

(

)

()

(

)

()

∫∫

−=

ΓΓ

ΓxxXΓxxXX du,Sdt,D

k

d

kijk

d

kij

d

ij

σ

e, de modo equivalente, para o contorno discretizado, na forma

∑

∫

∑

∫

==

−=

Ne

1n

n

d

kij

n

k

Ne

1n

n

d

kij

n

k

d

ij

dSudDt

ΓΓ

σ

Desenvolvendo a notação indicial e reagrupando matricialmente, de maneira

análoga ao que foi feito anteriormente, obtém-se

∑

∫

∑

∫

==

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

−

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

Ne

1n

n

2

n

1

n

d

222

d

122

d

221

d

121

d

212

d

112

d

211

d

111

Ne

1n

n

2

n

1

n

d

222

d

122

d

221

d

121

d

212

d

112

d

211

d

111

d

22

d

21

d

12

d

11

u

d

SS

t

d

DD

ΓΓ

σ

Entretanto, levando em conta que σ

12

= σ

21

, a expressão anterior simplifica para

∑

∫

∑

∫

==

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

−

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

Ne

1n

n

2

n

1

n

d

222

d

122

d

212

d

112

d

211

d

111

Ne

1n

n

2

n

1

n

d

222

d

122

d

212

d

112

d

211

d

111

d

22

d

12

d

11

u

d

SS

t

d

DD

ΓΓ

σ

Adotando as mesmas notações já utilizadas, acrescidas de

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

=

d

22

d

12

d

11

d

σ

⎥

⎦

⎤

⎢

⎣

⎡

=

d

222

d

122

d

212

d

112

d

211

d

111

d

DD

D

⎥

⎦

⎤

⎢

⎣

⎡

=

d

222

d

122

d

212

d

112

d

211

d

111

d

SS

S

resulta

∑

∫

∑

∫

==

−=

Ne

1n

n

d

Ne

1n

n

dd

ΓΓ

uΓStΓDσ

Denotando as integrais por

∫

=

n

ddn

d

Γ

ΓDD

∫

=

n

ddn

d

Γ

ΓSS

obtém-se

∑∑

==

−=

Ne

1n

ndn

Ne

1n

ndnd

uStDσ

Essa equação relaciona as tensões em um ponto qualquer do domínio com

os valores dos deslocamentos e forças de superfície de todo o contorno. Assim

69

como no cálculo dos deslocamentos nos pontos internos, pode ser facilmente

avaliada dado que as soluções fundamentais não apresentam singularidades.

3.8 Elementos Isoparamétricos

Uma contribuição significativa na elaboração do Método de Elementos de

Contorno foi a introdução de uma representação paramétrica para a geometria e

para as funções incógnitas, similar à formulação isoparamétrica do Método de

Elementos Finitos

12

. Nesse tipo de formulação, o parâmetro geométrico x, as

incógnitas de deslocamento u

j

(x) e as incógnitas de forças de superfície t

j

(x) são

aproximadas utilizando funções de interpolação, da seguinte maneira:

()

⎪

⎩

⎪

⎨

⎧

=

∑

=

m

1

jj

m

1

jj

m

1

jj

tNt

uNu

xNx

α

η

(3.34)

onde N

α

são polinômios de grau m-1 e apresentam, como característica fundamental,

o fato de terem valor unitário no ponto

α e valor nulo nos demais pontos. Os termos

α

j

x ,

α

j

u e

α

j

t são, respectivamente, as coordenadas geométricas, os deslocamentos

e as forças de superfície no ponto

α. Essas funções de interpolação são definidas

em termos de coordenadas adimensionais

η, sendo que -1 ≤ η ≤ 1, como mostra a

Figura 3.12.

12

LACHAT, J. C. A Further Development of the Boundary Integral Technique for

Elastostatics. PhD Thesis, University of Southampton, 1975 e LACHAT, J. C., WATSON, J. O.

Effective numerical treatment of boundary integral equations: A formulation for three-dimensional

elastostatics. International Journal for Numerical Methods in Engineering, n. 10, p. 991-1005,

1976.

70

FIGURA 3.12 – Elementos de contorno bidimensionais definidos em termos de um sistema de

coordenadas adimensional

Para elementos lineares, ou seja, m = 2, tem-se

()

η

−= 1

2

1

N

1

()

η

+= 1

2

1

N

2

(3.35)

e para elementos quadráticos, onde m = 3,

()

1

2

1

N

1

−=

ηη

2

1N

η

−=

()

1

2

1

N

2

+=

ηη

(3.36)

Em geral, as funções de interpolação podem ser derivadas dos polinômios

Lagrangianos, os quais são definidos, para o grau

m-1, como

()

∏

≠=

−

=

m

i0i

i

N

α

ηη

η

(3.37)

É possível observar que N

α

(η) é dado pelo produto de m fatores lineares. Os

polinômios de Lagrange possuem as seguintes propriedades: para um ponto

qualquer β, N

α

(η

β

) = δ

αβ

, onde η

0

≤ η ≤ η

m

; a soma dos polinômios é igual à unidade,

ou seja,

()

1N

m

1

=

∑

=

α

η

; a soma das derivadas primeiras dos polinômios é igual a zero.

71

Uma formulação discretizada de elementos de contorno é alcançada através

da substituição das expressões definidas em (3.34) na equação integral de contorno

(3.32), resultando em

() () () ()

∑∑∑∑

====

=+

Ne

1n

m

1

n

j

n

ij

Ne

1n

m

1

n

j

n

ijjij

tQuPuC

α

αα

α

αα

x'x'x'x'

(3.38)

onde os coeficientes

()

x'

α

n

ij

P e

(

)

x'

α

n

ij

Q são definidos em termos das integrais das

soluções fundamentais sobre os segmentos de contorno Γ

n

(com dΓ

n

= J

n

(η)dη, onde

J

n

(η) é o Jacobiano da transformação de coordenadas) da seguinte maneira:

( ) () ()

[]

()

∫

+

−

=

1

nij

n

ij

dJ,TNP

ηηηη

α

xx'x'

(3.39)

( ) () ()

[]

()

∫

+

−

=

1

nij

n

ij

dJ,UNQ

ηηηη

α

xx'x' (3.40)

O Jacobiano da transformação de coordenadas é dado por

13

()

2

m

1

n

2

m

1

n

1

2

1

n

x

d

dN

x

d

dN

d

dx

d

dx

J

⎟

⎠

⎞

⎜

⎝

⎛

+

⎟

⎠

⎞

⎜

⎝

⎛

=

⎟

⎠

⎞

⎜

⎝

⎛

+

⎟

⎠

⎞

⎜

⎝

⎛

=

∑∑

==

α

ηηηη

η

(3.41)

O método mais direto de solução da equação (3.38) é o método de

colocação pontual, ou seja, essa equação é avaliada nos pontos nodais x

c

, para c =

1, 2, ..., Nn

, onde Nn corresponde ao número total de nós de contorno, gerando,

desse modo, um sistema algébrico de equações que pode ser facilmente resolvido,

após a aplicação das condições de contorno prescritas do problema, para

determinar os valores das incógnitas de deslocamentos e forças de superfície.

Assim:

∑∑∑∑

====

=+

Ne

1n

m

1

n

j

cn

ij

Ne

1n

m

1

n

j

cn

ij

c

j

c

ij

tQuPuC

α

αα

α

αα

c = 1, 2, ..., Nn (3.42)

O duplo somatório no lado esquerdo de (3.42) deve ser avaliado mantendo-

se em mente que alguns nós são compartilhados entre elementos. Uma vez que os

valores de deslocamentos são definidos de maneira única nesses nós, esse

somatório pode, então, ser combinado em um único somatório sobre todos os nós

de contorno, ou seja,

13

ALIABADI, M. H. The Boundary Element Method, Vol. 2, Applications in Solids and

Structures. Chippenham: John Wiley & Sons, 2002. p. 45.

72

∑∑∑

===

=+

Ne

1n

m

1

n

j

cn

ij

Nn

1

j

c

ij

c

j

c

ij

tGuHuC

α

αα

γ

γγ

c = 1, 2, ..., Nn (3.43)

onde

γ

c

ij

H é formado da combinação de coeficientes

α

cn

ij

P e

α

cn

ij

G é igual aos

coeficientes

α

cn

ij

Q .

Definindo

⎪

⎩

⎪

⎨

⎧

=+=

≠=

γ

γγ

cquandoCHH

cquandoHH

c

ij

c

ij

c

ij

c

ij

c

ij

resulta

∑∑∑

===

=

Ne

1n

m

1

n

j

cn

ij

Nn

1

j

c

ij

tGuH

α

αα

γ

γγ

c = 1, 2, ..., Nn (3.44)

Essa equação pode ser escrita em forma matricial como

tGuH

=

onde H é uma matriz de dimensões 2Nn x 2Nn; G é uma matriz de dimensões 2Nn x

2(Ne)(m); u é um vetor com 2Nn componentes e t é um vetor com 2(Ne)(m)

componentes. Após a substituição das condições de contorno prescritas para o

problema, o sistema de equações anterior produz o sistema algébrico

FYA

=

onde o vetor Y contém todas as incógnitas (deslocamentos ou forças de superfície)

do problema; a matriz de coeficientes A é formada de colunas de H e G; o vetor F é

formado de colunas de

H e G multiplicadas pelas condições de contorno prescritas.

Uma vez que esse sistema esteja resolvido, todos os valores das incógnitas de

contorno serão conhecidos.

Da expressão (3.44) é possível deduzir que as submatrizes diagonais de H,

as quais contêm os coeficientes

c

ij

C , podem ser avaliadas indiretamente a partir de

uma consideração de movimento de corpo rígido, ou seja, se for imposto ao

problema somente condições de contorno de deslocamentos, tem-se

0uH

=

e, desse modo,

∑

≠=

−=

Nn

c1

c

ij

cc

ij

HH

γγ

γ

(3.45)

73

A expressão (3.45) só é válida para domínios finitos. Para domínios infinitos

as condições de regularidade no infinito não são mais satisfeitas, isso porque um

movimento de corpo rígido assume

(

)

(

)

1Ou

j

=

x

Nessa situação, a condição de regularidade torna-se

()() ()()

[]

() ( ) ()

xΓxx'xΓxxx'xxx'

ΓΓ

jij

0

jjijjij

0

ud,Tlimudt,Uu,Tlim ==

⎭

⎬

⎫

⎩

⎨

⎧

−

∫∫

→→

ρ

e a seguinte relação é obtida

∑

≠=

−=

Nn

c1

c

ijij

cc

ij

HH

γγ

γ

δ

(3.46)

3.8.1 Formação do Sistema de Equações

Para ilustrar a aplicação da equação (3.42) far-se-á uso das malhas modelo

apresentadas na Figura 3.13 e na Figura 3.14, para os casos de elementos de

contorno lineares e quadráticos, respectivamente.

3.8.1.1 Elementos Lineares

Aplicando a equação

(3.42), quando o ponto fonte é o Nó 1, resulta

()

=+++++

∑

=

4Ne

1n

2n

2

2n1

12

2n

1

2n1

11

1n

2

1n1

12

1n

1

1n1

11

1

2

1

12

1

11

uPuPuPuPuCuC

()

∑

=

+++=

4Ne

1n

2n

2

2n1

12

2n

1

2n1

11

1n

2

1n1

12

1n

1

1n1

11

tQtQtQtQ (i)

()

=+++++

∑

=

4Ne

1n

2n

2

2n1

22

2n

1

2n1

21

1n

2

1n1

22

1n

1

1n1

21

1

2

1

22

1

21

uPuPuPuPuCuC

()

∑

=

+++=

4Ne

1n

2n

2

2n1

22

2n

1

2n1

21

1n

2

1n1

22

1n

1

1n1

21

tQtQtQtQ (ii)

74

FIGURA 3.13 – Malha modelo formada por elementos de contorno lineares

Em forma matricial, essas equações passam a ser agrupadas do seguinte

modo:

=

⎟

⎠

⎞

⎜

⎝

⎛

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

∑

=

4Ne

1n

2n

2

2n

1

1n

2

1n

1

2n1

22

2n1

21

1n1

22

1n1

21

2n1

12

2n1

11

1n1

12

1n1

11

1

2

1

22

1

21

1

12

1

11

u

PPPP

u

CC

∑

=

⎟

⎠

⎞

⎜

⎝

⎛

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

4Ne

1n

2n

2

2n

1

1n

2

1n

1

2n1

22

2n1

21

1n1

22

1n1

21

2n1

12

2n1

11

1n1

12

1n1

11

t

QQQQ

De maneira análoga, para o Nó 2,

=

⎟

⎠

⎞

⎜

⎝

⎛

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

∑

=

4Ne

1n

2n

2

2n

1

1n

2

1n

1

2n2

22

2n2

21

1n2

22

1n2

21

2n2

12

2n2

11

1n2

12

1n2

11

2

1

2

22

2

21

2

12

2

11

u

PPPP

u

CC

75

∑

=

⎟

⎠

⎞

⎜

⎝

⎛

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

4Ne

1n

2n

2

2n

1

1n

2

1n

1

2n2

22

2n2

21

1n2

22

1n2

21

2n2

12

2n2

11

1n2

12

1n2

11

t

QQQQ

para o Nó 3,

=

⎟

⎠

⎞

⎜

⎝

⎛

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

∑

=

4Ne

1n

2n

2

2n

1

1n

2

1n

1

2n3

22

2n3

21

1n3

22

1n3

21

2n3

12

2n3

11

1n3

12

1n3

11

3

2

3

1

3

22

3

21

3

12

3

11

u

PPPP

u

CC

∑

=

⎟

⎠

⎞

⎜

⎝

⎛

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

4Ne

1n

2n

2

2n

1

1n

2

1n

1

2n3

22

2n3

21

1n3

22

1n3

21

2n3

12

2n3

11

1n3

12

1n3

11

t

QQQQ

e para o Nó 4,

=

⎟

⎠

⎞

⎜

⎝

⎛

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

∑

=

4Ne

1n

2n

2

2n

1

1n

2

1n

1

2n4

22

2n4

21

1n4

22

1n4

21

2n4

12

2n4

11

1n4

12

1n4

11

4

2

4

1

4

22

4

21

4

12

4

11

u

PPPP

u

CC

∑

=

⎟

⎠

⎞

⎜

⎝

⎛

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

4Ne

1n

2n

2

2n

1

1n

2

1n

1

2n4

22

2n4

21

1n4

22

1n4

21

2n4

12

2n4

11

1n4

12

1n4

11

t

QQQQ

ou, de maneira equivalente,

=

⎟

⎠

⎞

⎜

⎝

⎛

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

∑

=

4Ne

1n

2n

2

2n

1

1n

2

1n

1

2cn

22

2cn

21

1cn

22

1cn

21

2cn

12

2cn

11

1cn

12

1cn

11

c

2

c

1

c

22

c

21

c

12

c

11

u

PPPP

u

CC

∑

=

⎟

⎠

⎞

⎜

⎝

⎛

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

4Ne

1n

2n

2

2n

1

1n

2

1n

1

2cn

22

2cn

21

1cn

22

1cn

21

2cn

12

2cn

11

1cn

12

1cn

11

t

QQQQ

c = 1, 2, ..., Nn = 4

76

Tendo em mente que alguns nós são compartilhados entre elementos e que

os valores dos deslocamentos são definidos de maneira única nesses nós, o

somatório à esquerda na expressão anterior pode ser avaliado da seguinte maneira:

∑

=

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

4Ne

1n

2n

2

2n

1

2cn

22

2cn

21

2cn

12

2cn

11

1n

2

1n

1

1cn

22

1cn

21

1cn

12

1cn

11

u

PP

u

PP

Expandindo o somatório obtém-se

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

11

2

11

1

11c

22

11c

21

11c

12

11c

11

u

PP

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

12

2

12

1

12c

22

12c

21

12c

12

12c

11

u

PP

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

21

2

21

1

21c

22

21c

21

21c

12

21c

11

u

PP

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

22

2

22

1

22c

22

22c

21

22c

12

22c

11

u

PP

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

31

2

31

1

31c

22

31c

21

31c

12

31c

11

u

PP

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

32

2

32

1

32c

22

32c

21

32c

12

32c

11

u

PP

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

41

2

41

1

41c

22

41c

21

41c

12

41c

11

u

PP

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

42

2

42

1

42c

22

42c

21

42c

12

42c

11

u

PP

As parcelas que apresentam a mesma cor devem ser combinadas, pois

possuem nós comuns a um mesmo elemento, isto é: o nó 2 do elemento 1 é o nó 1

do elemento 2; o nó 2 do elemento 2 é o nó 1 elemento 3; o nó 2 do elemento 3 é o

nó 1 do elemento 4; o nó 2 do elemento 4 é o nó 1 do elemento 1. Adotando uma

notação do tipo nó

γ

= 1, 2, ..., Nn, ao invés da notação "nó

α

do elemento n", é

possível agrupar esses termos do seguinte modo:

∑

=

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

4Nn

1

2

1

c

22

c

21

c

12

c

11

u

HH

γ

γγ

A equação geral pode, agora, ser reescrita como

=

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

∑

=

4Nn

1

2

1

c

22

c

21

c

12

c

11

c

2

c

1

c

22

c

21

c

12

c

11

u

HH

u

CC

γ

γγ

∑

=

⎟

⎠

⎞

⎜

⎝

⎛

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

4Ne

1n

2n

2

2n

1

1n

2

1n

1

2cn

22

2cn

21

1cn

22

1cn

21

2cn

12

2cn

11

1cn

12

1cn

11

t

GGGG

c = 1, 2, ..., Nn = 4

onde

αα

cn

ij

cn

ij

QG = .

77

Agrupando todos os deslocamentos incógnitos, resulta

∑∑

=

⎟

⎠

⎞

⎜

⎝

⎛

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

4Ne

1n

2n

2

2n

1

1n

2

1n

1

2cn

22

2cn

21

1cn

22

1cn

21

2cn

12

2cn

11

1cn

12

1cn

11

4Nn

1

2

1

c

22

c

21

c

12

c

11

t

GGGG

u

HH

γ

γγ

c = 1, 2, ..., Nn = 4

A equação anterior pode ser reescrita em forma matricial como

tGuH

=

onde H é uma matriz de dimensões 2Nn x 2Nn; G é uma matriz de dimensões 2Nn x

4Ne; u é um vetor com 2Nn componentes e t é um vetor com 4Ne componentes.

3.8.1.2 Elementos Quadráticos

Aplicando a equação (3.42), quando o ponto fonte é o Nó 1, resulta

()

=+++++++

∑

=

4Ne

1n

3n

2

3n1

12

3n

1

3n1

11

2n

2

2n1

12

2n

1

2n1

11

1n

2

1n1

12

1n

1

1n1

11

1

2

1

12

1

11

uPuPuPuPuPuPuCuC

()

∑

=

+++++=

4Ne

1n

3n

2

3n1

12

3n

1

3n1

11

2n

2

2n1

12

2n

1

2n1

11

1n

2

1n1

12

1n

1

1n1

11

tQtQtQtQtQtQ (i)

()

=+++++++

∑

=

4Ne

1n

3n

2

3n1

22

3n

1

3n1

21

2n

2

2n1

22

2n

1

2n1

21

1n

2

1n1

22

1n

1

1n1

21

1

2

1

22

1

21

uPuPuPuPuPuPuCuC

()

∑

=

+++++=

4Ne

1n

3n

2

3n1

22

3n

1

3n1

21

2n

2

2n1

22

2n

1

2n1

21

1n

2

1n1

22

1n

1

1n1

21

tQtQtQtQtQtQ (ii)

Procedendo de maneira análoga para os demais nós de contorno, é possível

resumir as equações resultantes numa única equação equivalente, dada por

=

⎟

⎠

⎞

⎜

⎝

⎛

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

∑

=

4Ne

1n

3n

2

3n

1

2n

2

2n

1

1n

2

1n

1

3cn

22

3cn

21

2cn

22

2cn

21

1cn

22

1cn

21

3cn

12

3cn

11

2cn

12

2cn

11

1cn

12

1cn

11

c

2

c

1

c

22

c

21

c

12

c

11

u

PPPPPP

u

CC

78

FIGURA 3.14 – Malha modelo formada por elementos de contorno quadráticos

∑

=

⎟

⎠

⎞

⎜

⎝

⎛

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

4Ne

1n

3n

2

3n

1

2n

2

2n

1

1n

2

1n

1

3cn

22

3cn

21

2cn

22

2cn

21

1cn

22

1cn

21

3cn

12

3cn

11

2cn

12

2cn

11

1cn

12

1cn

11

t

QQQQQQ

c = 1, 2, ..., Nn = 8

Tendo em mente que alguns nós são compartilhados entre elementos e que

os valores dos deslocamentos são definidos de maneira única nesses nós, o

somatório a esquerda na expressão anterior pode ser avaliado da seguinte maneira:

∑

=

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

4Ne

1n

3n

2

3n

1

3cn

22

3cn

21

3cn

12

3cn

11

2n

2

2n

1

2cn

22

2cn

21

2cn

12

2cn

11

1n

2

1n

1

1cn

22

1cn

21

1cn

12

1cn

11

u

PP

u

PP

u

PP

Expandindo o somatório obtém-se

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

11

2

11

1

11c

22

11c

21

11c

12

11c

11

u

PP

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

12

2

12

1

12c

22

12c

21

12c

12

12c

11

u

PP

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

13

2

13

1

13c

22

13c

21

13c

12

13c

11

u

PP

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

21

2

21

1

21c

22

21c

21

21c

12

21c

11

u

PP

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

22

2

22

1

22c

22

22c

21

22c

12

22c

11

u

PP

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

23

2

23

1

23c

22

23c

21

23c

12

23c

11

u

PP

+

79

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

31

2

31

1

31c

22

31c

21

31c

12

31c

11

u

PP

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

32

2

32

1

32c

22

32c

21

32c

12

32c

11

u

PP

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

33

2

33

1

33c

22

33c

21

33c

12

33c

11

u

PP

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

41

2

41

1

41c

22

41c

21

41c

12

41c

11

u

PP

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

42

2

42

1

42c

22

42c

21

42c

12

42c

11

u

PP

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

43

2

43

1

43c

22

43c

21

43c

12

43c

11

u

PP

As parcelas que apresentam a mesma cor devem ser combinadas, pois

possuem nós comuns a um mesmo elemento, isto é: o nó 3 do elemento 1 é o nó 1

do elemento 2; o nó 3 do elemento 2 é o nó 1 elemento 3; o nó 3 do elemento 3 é o

nó 1 do elemento 4; o nó 3 do elemento 4 é o nó 1 do elemento 1. Adotando uma

notação do tipo nó

γ

= 1, 2, ..., Nn, ao invés da notação "nó

α

do elemento n", é

possível agrupar esses termos do seguinte modo:

∑

=

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

8Nn

1

2

1

c

22

c

21

c

12

c

11

u

HH

γ

γγ

A equação geral pode, agora, ser reescrita como

=

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

+

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

∑

=

8Nn

1n

2

1

c

22

c

21

c

12

c

11

c

2

c

1

c

22

c

21

c

12

c

11

u

HH

u

CC

γ

γγ

∑

=

⎟

⎠

⎞

⎜

⎝

⎛

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

4Ne

1n

3n

2

3n

1

2n

2

2n

1

1n

2

1n

1

3cn

22

3cn

21

2cn

22

2cn

21

1cn

22

1cn

21

3cn

12

3cn

11

2cn

12

2cn

11

1cn

12

1cn

11

t

QQQQQQ

c = 1, 2, ..., Nn = 8

onde

αα

cn

ij

cn

ij

QG = .

Agrupando todos os deslocamentos incógnitos, resulta

∑∑

=

⎟

⎠

⎞

⎜

⎝

⎛

⎪

⎭

⎪

⎬

⎫

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

=

⎟

⎠

⎞

⎜

⎝

⎛

⎭

⎬

⎫

⎩

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

4Ne

1n

3n

2

3n

1

2n

2

2n

1

1n

2

1n

1

3cn

22

3cn

21

2cn

22

2cn

21

1cn

22

1cn

21

3cn

12

3cn

11

2cn

12

2cn

11

1cn

12

1cn

11

8Nn

1

2

1

c

22

c

21

c

12

c

11

t

QQQQQQ

u

HH

γ

γγ

c = 1, 2, ..., Nn = 8

A equação anterior pode ser reescrita em forma matricial como

tGuH

=

80

onde H é uma matriz de dimensões 2Nn x 2Nn; G é uma matriz de dimensões 2Nn x

6Ne; u é um vetor com 2Nn componentes e t é um vetor com 6Ne componentes.

3.8.2 Integração

A avaliação dos coeficientes

α

cn

ij

P e

α

cn

ij

Q , dados pelas expressões (3.39) e

(3.40), respectivamente, quando o ponto fonte c não coincide com nenhum dos nós

do elemento n, não apresenta maiores dificuldades e pode ser alcançada através do

processo de quadratura de Gauss convencional. No entanto, quando o ponto fonte c

é um dos nós do elemento n, as expressões (3.39) e (3.40) apresentarão

singularidades e precisam ser avaliadas de maneira ligeiramente diferente. As

próximas seções descrevem o processo matemático necessário para contornar essa

dificuldade, tanto para elementos lineares quanto para elementos quadráticos.

3.8.2.1 Elementos Lineares – Coeficientes

α

cn

ij

Q

A expressão (3.40) pode ser reescrita como

() ()

[]

()

∫

+

−

=

1

n

cn

ij

cn

ij

dJrUNQ

ηηηη

α

(3.47)

Para problemas lineares, as funções de interpolação são dadas por

() ()

ηη

−= 1

2

1

N

1

() ()

ηη

+= 1

2

1

N

2

(3.48)

e a geometria do contorno, por sua vez, é dada por

()

∑

=

2

1

jj

xNx

α

η

j = 1, 2 (3.49)

Desse modo, as coordenadas de um ponto qualquer no contorno são escritas como

(

)

(

)

() ()

⎪

⎩

⎪

⎨

⎧

+=

2

22

1

21

g

2

12

1

11

g

1

xNxNx

ηη

(3.50)

Substituindo (3.48) em (3.50) resulta

81

(

)

(

)

()()

⎪

⎩

⎪

⎨

⎧

++−=

1

2

1

2

g

2

1

2

1

2

1

g

1

xx

2

1

xx

2

1

x

xx

2

1

xx

2

1

x

η

(3.51)

Utilizando (3.51) em (3.41) é possível determinar J

n

(η), logo

()

()()

L

2

1

xx

2

1

xx

2

1

d

dx

d

dx

J

2

1

2

1

2

1

2

1

n

=

⎥

⎦

⎤

⎢

⎣

⎡

−+

⎥

⎦

⎤

⎢

⎣

⎡

−=

⎟

⎠

⎞

⎜

⎝

⎛

+

⎟

⎠

⎞

⎜

⎝

⎛

=

ηη

η

(3.52)

onde L é o comprimento do elemento n que está sendo integrado.

Voltando à expressão (3.47), a solução fundamental é expressa por

()

[]

()

() ()

⎥

⎦

⎤

⎢

⎣

⎡

+

⎟

⎠

⎞

⎜

⎝

⎛

−

=

ηηδ

η

ν

νπμ

η

jiij

cn

ij

,r,r

r

1

ln43

18

1

rU (3.53)

ou, então,

()

[]

()

() ()

ηη

νπμ

δ

ηνπμ

ν

η

jiij

cn

ij

,r,r

18

1

r

1

ln

18

43

rU

−

+

⎟

⎠

⎞

⎜

⎝

⎛

−

=

(3.54)

Definindo

(

)

()

νπμ

ν

−

=

18

43

1CTE

()

νπμ

−

=

18

1

2CTE (3.55)

a solução fundamental torna-se

()

[]

()

() ()

ηηδ

η

jiij

cn

ij

,r,r2CTE

r

1

ln1CTErU +

⎟

⎠

⎞

⎜

⎝

⎛

=

(3.56)

e, desse modo, a expressão (3.47) é reescrita na forma

()

() () ()

2

L

d,r,rN2CTE

2

L

d

r

1

lnN1CTEQ

1

ji

1

ij

cn

ij

∫∫

+

−

+

−

+

⎟

⎠

⎞

⎜

⎝

⎛

=

ηηηηηδ

η

αα

α

(3.57)

Na determinação dos termos

(

)

η

r

e

(

)

(

)

η

ji

,r,r , dois casos devem ser

considerados, de acordo com a posição do ponto fonte c dentro do elemento n,

conforme ilustra a Figura 3.15: o Caso 1 corresponde ao ponto fonte localizado no nó

inicial do elemento; o Caso 2 corresponde ao ponto fonte localizado no nó final do

elemento.

82

FIGURA 3.15 – Relação entre os sistemas de coordenadas adimensionais para elementos de

contorno lineares conforme a posição do ponto fonte dentro do elemento

Caso 1:

()

(

)

(

)

2

1

2

g

2

1

g

1

xxxxr −+−=

η

(

)

()

(

)

()

⎪

⎩

⎪

⎨

⎧

−+=−

1

2

1

2

g

2

1

2

1

g

1

xx1

2

1

xx

xx1

2

1

xx

η

() ()

(

)

(

)

()

ηηη

+=−+−+= 1L

2

1

xxxx1

2

1

r

2

1

2

1

2

1

(3.58)

() ()

(

)

()

(

)

()

(

)

2

1

2

1

g

1

g

1

11

L

xx

r

xx

r

xx

,r,r

−

=

−−

=

ηη

(3.59)

() ()

(

)

()

(

)

()

(

)

2

1

2

1

2

g

2

1

2

g

2

22

L

xx

r

xx

r

xx

,r,r

−

=

−−

=

ηη

(3.60)

() ()

(

)

()

(

)

()

(

)

(

)

2

1

2

1

2

1

2

g

2

1

g

1

21

L

xxxx

r

xx

r

xx

,r,r

−−

=

−−

=

ηη

(3.61)

83

Caso 2:

()

(

)

(

)

2

g

2

1

g

1

xxxxr −+−=

η

(

)

()

(

)

()

⎪

⎩

⎪

⎨

⎧

−−=−

2

1

2

g

2

1

2

1

g

1

xx1

2

1

xx

xx1

2

1

xx

η

() ()

(

)

(

)

()

ηηη

−=−+−−= 1L

2

1

xxxx1

2

1

r

2

1

2

1

(3.62)

() ()

(

)

()

(

)

()

(

)

2

1

2

1

g

1

2

1

g

1

11

L

xx

r

xx

r

xx

,r,r

−

=

−−

=

ηη

(3.63)

() ()

(

)

()

(

)

()

(

)

2

1

2

g

2

g

2

22

L

xx

r

xx

r

xx

,r,r

−

=

−−

=

ηη

(3.64)

() ()

(

)

()

(

)

()

(

)

(

)

2

1

2

1

2

g

2

1

g

1

21

L

xxxx

r

xx

r

xx

,r,r

−−

=

−−

=

ηη

(3.65)

Substituindo esses resultados na equação (3.57) resulta

Caso 1:

() () ()

2

L

dN,r,r2CTE

2

L

d1L

2

1

lnN1CTEQ

1

ji

1

ij

cn

ij

∫∫

+

−

+

−

+

⎟

⎠

⎞

⎜

⎝

⎛

+−=

ηηηδηη

αα

α

() () ()

+

⎥

⎦

⎤

⎢

⎣

⎡

⎟

⎠

⎞

⎜

⎝

⎛

+

+−=

∫∫

+

−

+

−

2

L

d

2

1

lnNdLlnN1CTEQ

1

ij

1

ij

cn

ij

ηδ

η

ηηδη

αα

α

()

2

L

dN,r,r2CTE

1

ji

∫

+

−

+

ηη

α

() () ()

+−

⎟

⎠

⎞

⎜

⎝

⎛

+

−=

∫∫

+

−

+

−

2

L

dNLln1CTE

2

L

d

2

1

lnN1CTEQ

1

ij

1

ij

cn

ij

ηδηηδ

η

αα

α

()

2

L

dN,r,r2CTE

1

ji

∫

+

−

+

ηη

α

(3.66)

Caso 2:

() () ()

2

L

dN,r,r2CTE

2

L

d1L

2

1

lnN1CTEQ

1

ji

1

ij

cn

ij

∫∫

+

−

+

−

+

⎟

⎠

⎞

⎜

⎝

⎛

−−=

ηηηδηη

αα

α

() () ()

+

⎥

⎦

⎤

⎢

⎣

⎡

⎟

⎠

⎞

⎜

⎝

⎛

−

+−=

∫∫

+

−

+

−

2

L

d

2

1

lnNdLlnN1CTEQ

1

ij

1

ij

cn

ij

ηδ

η

ηηδη

αα

α

()

2

L

dN,r,r2CTE

1

ji

∫

+

−

+

ηη

α

84

() () ()

+−

⎟

⎠

⎞

⎜

⎝

⎛

−

−=

∫∫

+

−

+

−

2

L

dNLln1CTE

2

L

d

2

1

lnN1CTEQ

1

ij

1

ij

cn

ij

ηδηηδ

η

αα

α

()

2

L

dN,r,r2CTE

1

ji

∫

+

−

+

ηη

α

(3.67)

As singularidades na primeira integral em (3.66) e na primeira integral em

(3.67) podem ser contornadas através das transformações de coordenadas

apresentadas na Figura 3.15, ou seja,

Caso 1:

2

1

1c

η

ξ

+

=⇒=

(3.68)

Caso 2:

2

1

2c

η

ξ

−

=⇒= (3.69)

Consequentemente, novas funções de interpolação serão obtidas para cada

caso da posição do ponto fonte, isto é,

Caso 1:

(

)

(

)

ξ

−

=

1N

1

(

)

ξ

=

2

N (3.70)

Caso 2:

(

)

ξ

=

1

N

(

)

(

)

ξ

−

=

1N

2

(3.71)

Também são encontradas novas relações (3.51) e (3.52) para cada caso da

posição do ponto fonte.

Caso 1:

(

)

()

⎪

⎩

⎪

⎨

⎧

+−=

1

2

1

2

g

2

1

2

1

g

1

xxxx

ξ

(3.72)

()

()()

Lxxxx

d

dx

d

dx

J

2

1

2

1

2

1

2

1

n

=−+−=

⎟

⎠

⎞

⎜

⎝

⎛

+

⎟

⎠

⎞

⎜

⎝

⎛

=

ξξ

ξ

(3.73)

Caso 2:

(

)

()

⎪

⎩

⎪

⎨

⎧

+−−=

2

1

2

g

2

1

2

1

g

1

xxxx

ξ

(3.74)

()

()()

Lxxxx

d

dx

d

dx

J

2

1

2

1

2

1

2

1

n

=−+−=

⎟

⎠

⎞

⎜

⎝

⎛

+

⎟

⎠

⎞

⎜

⎝

⎛

=

ξξ

ξ

(3.75)

Substituindo (3.68) e (3.69) em (3.66) e (3.67), respectivamente, obtém-se

85

Caso 1:

() () () ()

+−−=

∫∫

+

−

+

2

L

dNLln1CTE

2

L

d2lnN1CTEQ

1

ij

1

0

ij

cn

ij

ηδηξδξξ

αα

α

()

2

L

dN,r,r2CTE

1

ji

∫

+

−

+

ηη

α

() () ()

+−

⎟

⎠

⎞

⎜

⎝

⎛

=

∫∫

+

−

+

2

L

dNLln1CTELd

1

lnN1CTEQ

1

ij

1

0

ij

cn

ij

ηδηξδ

ξ

αα

α

()

2

L

dN,r,r2CTE

1

ji

∫

+

−

+

ηη

α

(3.76)

Caso 2:

() () ( ) () ()

+−−−=

∫∫

+

−+

2

L

dNLln1CTE

2

L

d2lnN1CTEQ

1

ij

0

1

ij

cn

ij

ηδηξδξξ

αα

α

()

2

L

dN,r,r2CTE

1

ji

∫

+

−

+

ηη

α

() () ()

+−

⎟

⎠

⎞

⎜

⎝

⎛

=

∫∫

+

−

+

2

L

dNLln1CTELd

1

lnN1CTEQ

1

ij

1

0

ij

cn

ij

ηδηξδ

ξ

αα

α

()

2

L

dN,r,r2CTE

1

ji

∫

+

−

+

ηη

α

(3.77)

A primeira integral em (3.76) e a primeira integral em (3.77) são avaliadas por

meio de uma quadratura de Gauss modificada, denominada quadratura de Gauss

logarítmica (para valores calculados e tabelados de pontos e pesos de Gauss para a

quadratura Gaussiana logarítmica ver, por exemplo, BREBBIA e DOMINGUEZ

14

).

Essas integrais poderiam, alternativamente, ser calculadas analiticamente.

3.8.2.2 Elementos Lineares – Coeficientes

α

cn

ij

P

A expressão (3.39) pode ser reescrita como

() ()

[]

()

∫

+

−

=

1

n

cn

ij

cn

ij

dJrTNP

ηηηη

α

(3.78)

sendo que

14

BREBBIA, C. A., DOMINGUEZ, J. Boundary Elements, An Introductory Course. Boston:

Computational Mechanics Publications, 1989, p. 271.

86

()

[]

()()

(

)

()

[

]

()

⎭

⎬

⎫

⎩

⎨

⎧

−−++−

∂

−

=

ijjijiij

cn

ij

,rn,rn21,r,r221

n

r

r14

1

rT

νδν

η

ηνπ

η

ou

()

[]

()

[]

η

rf

r

1

rT

cn

ij

=

onde f [r (η)] é, nitidamente, não singular.

Verificando o comportamento da expressão (3.78), para cada caso da

posição do ponto fonte, observa-se:

Caso 1:

()

[]

() ()

()

[]

=

⎟

⎠

⎞

⎜

⎝

⎛

⎟

⎠

⎞

⎜

⎝

⎛

+

⎟

⎠

⎞

⎜

⎝

⎛

−==

∫∫

+

−

+

−

1

n1

1cn

ij

dL

2

1

rf

1L

2

1

2

1

dJrf

r

1

NP

ηη

η

ηηηη

η

()

[]

∫

+

−

+

−

=

1

drf

1

2

1

ηη

η

()

[]

() ()

()

[]

=

⎟

⎠

⎞

⎜

⎝

⎛

⎟

⎠

⎞

⎜

⎝

⎛

+

⎟

⎠

⎞

⎜

⎝

⎛

+==

∫∫

+

−

+

−

1

n2

2cn

ij

dL

2

1

rf

1L

2

1

2

1

dJrf

r

1

NP

ηη

η

ηηηη

η

()

[]

∫

+

−

=

1

drf

2

1

ηη

(Conclusão: apenas a primeira expressão apresenta singularidade).

Caso 2:

()

[]

() ()

()

[]

=

⎟

⎠

⎞

⎜

⎝

⎛

⎟

⎠

⎞

⎜

⎝

⎛

−

⎟

⎠

⎞

⎜

⎝

⎛

−==

∫∫

+

−

+

−

1

n1

1cn

ij

dL

2

1

rf

1L

2

1

2

1

dJrf

r

1

NP

ηη

η

ηηηη

η

()

[]

∫

+

−

=

1

drf

2

1

ηη

()

[]

() ()

()

[]

=

⎟

⎠

⎞

⎜

⎝

⎛

⎟

⎠

⎞

⎜

⎝

⎛

−

⎟

⎠

⎞

⎜

⎝

⎛

+==

∫∫

+

−

+

−

1

n2

2cn

ij

dL

2

1

rf

1L

2

1

2

1

dJrf

r

1

NP

ηη

η

ηηηη

η

()

[]

∫

+

−

+

=

1

drf

1

2

1

ηη

η

(Conclusão: apenas a segunda expressão apresenta singularidade).

Os coeficientes

α

cn

ij

P singulares são exatamente os coeficientes que

compõem as submatrizes diagonais de H, as quais são calculadas indiretamente

87

através de consideração de movimento de corpo rígido, como discutido

anteriormente. A conclusão é que, para efeitos de implementação computacional, a

avaliação dos coeficientes

α

cn

ij

P , quando c é um dos nós do elemento n, pode ser

feita através da quadratura de Gauss convencional, pois os coeficientes singulares,

calculados imprecisamente através desse processo, são desprezados e

recalculados em função dos demais coeficientes, resultantes de integrações não

singulares e, portanto, precisamente avaliados.

3.8.2.3 Elementos Quadráticos – Coeficientes

α

cn

ij

Q

A expressão (3.40) pode ser reescrita como

() ()

[]

()

∫

+

−

=

1

n

cn

ij

cn

ij

dJrUNQ

ηηηη

α

(3.79)

Para problemas quadráticos, as funções de interpolação são dadas por

()

1

2

1

N

1

−=

ηη

2

1N

η

−=

()

1

2

1

N

2

+=

ηη

(3.80)

e a geometria do contorno, por sua vez, é dada por

()

∑

=

3

1

jj

xNx

α

η

j = 1, 2 (3.81)

Desse modo, as coordenadas de um ponto qualquer no contorno são escritas como

(

)

(

)

(

)

() () ()

⎪

⎩

⎪

⎨

⎧

++=

3

23

2

22

1

21

g

2

3

13

2

12

1

11

g

1

xNxNxNx

ηηη

(3.82)

Substituindo (3.80) em (3.82) resulta

(

)

(

)

()()

⎪

⎩

⎪

⎨

⎧

+−++−=

2

1

2

3

2

3

2

1

2

g

2

1

3

1

3

1

2

1

2

g

1

xxx

2

1

xx2x

2

1

x

xxx

2

1

xx2x

2

1

x

ηη

(3.83)

Utilizando (3.83) em (3.41) é possível determinar J

n

(η), logo

()

=

⎟

⎠

⎞

⎜

⎝

⎛

+

⎟

⎠

⎞

⎜

⎝

⎛

=

2

1

n

d

dx

d

dx

J

ηη

η

()()( )()

2

1

2

3

2

3

2

1

2

1

3

1

3

1

2

1

xx

2

1

xx2xxx

2

1

xx2x

⎥

⎦

⎤

⎢

⎣

⎡

−++−+

⎥

⎦

⎤

⎢

⎣

⎡

−++−=

ηη

(3.84)

88

Voltando à expressão (3.79), a solução fundamental é expressa por

()

[]

()

() ()

⎥

⎦

⎤

⎢

⎣

⎡

+

⎟

⎠

⎞

⎜

⎝

⎛

−

=

ηηδ

η

ν

νπμ

η

jiij

cn

ij

,r,r

r

1

ln43

18

1

rU

(3.85)

ou, então,

()

[]

()

() ()

ηη

νπμ

δ

ηνπμ

ν

η

jiij

cn

ij

,r,r

18

1

r

1

ln

18

43

rU

−

+

⎟

⎠

⎞

⎜

⎝

⎛

−

=

(3.86)

Definindo novamente

(

)

()

νπμ

ν

−

=

18

43

1CTE

()

νπμ

−

=

18

1

2CTE

(3.87)

a solução fundamental torna-se

()

[]

()

() ()

ηηδ

η

jiij

cn

ij

,r,r2CTE

r

1

ln1CTErU +

⎟

⎠

⎞

⎜

⎝

⎛

=

(3.88)

e, desse modo, a expressão (3.79) é reescrita na forma

()

+

⎟

⎠

⎞

⎜

⎝

⎛

=

∫

+

−

1

nij

cn

ij

dJ

r

1

lnN1CTEQ

ηηδ

η

α

() () () ()

∫

+

−

+

1

nji

dJ,r,rN2CTE

ηηηηη

α

(3.89)

Na determinação dos termos

(

)

η

r e

(

)

(

)

η

ji

,r,r , três casos devem ser

considerados, de acordo com a posição do ponto fonte c dentro do elemento n,

conforme ilustra a Figura 3.16: o Caso 1 corresponde ao ponto fonte localizado no nó

inicial do elemento; o Caso 2 corresponde ao ponto fonte localizado no nó

intermediário do elemento; o

Caso 3 corresponde ao ponto fonte localizado no nó

final do elemento.

Caso 1:

()

(

)

(

)

2

1

2

g

2

1

g

1

xxxxr −+−=

η

()

()()

()

()()

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

−++−+=−

⎥

⎦

⎤

⎢

⎣

⎡

−++−+=−

1

2

3

2

1

2

1

2

g

2

1

2

1

3

1

2

1

g

1

xxxx2x

2

1

1xx

xxxx2x

2

1

1xx

ηη

89

()

()()

2

1

2

3

2

1

2

1

2

1

3

1

2

1

xxxx2x

2

1

xxxx2x

2

1

1r

⎥

⎦

⎤

⎢

⎣

⎡

−++−+

+

⎥

⎦

⎤

⎢

⎣

⎡

−++−

+=

η

ηη

(3.90)

FIGURA 3.16 – Relação entre os sistemas de coordenadas adimensionais para elementos de

contorno quadráticos conforme a posição do ponte fonte dentro do elemento

() ()

(

)

()

(

)

()

=

−−

=

ηη

r

xx

r

xx

,r,r

1

g

1

g

1

11

()()

()()()()

2

1

2

3

2

1

2

1

2

1

3

1

2

1

2

1

2

1

3

1

2

1

xxxx2x

2

1

xxxx2x

2

1

xxxx2x

2

1

⎥

⎦

⎤

⎢

⎣

⎡

−++−+

⎥

⎦

⎤

⎢

⎣

⎡

−++−

⎥

⎦

⎤

⎢

⎣

⎡

−++−

=

ηη

η

(3.91)

() ()

(

)

()

(

)

()

=

−−

=

ηη

r

xx

r

xx

,r,r

1

2

g

2

1

2

g

2

22

()()

()()()()

2

1

2

3

2

1

2

1

2

1

3

1

2

1

2

1

2

3

2

1

2

xxxx2x

2

1

xxxx2x

2

1

xxxx2x

2

1

⎥

⎦

⎤

⎢

⎣

⎡

−++−+

⎥

⎦

⎤

⎢

⎣

⎡

−++−

⎥

⎦

⎤

⎢

⎣

⎡

−++−

=

ηη

η

(3.92)

() ()

(

)

()

(

)

()

=

−−

=

ηη

r

xx

r

xx

,r,r

1

2

g

2

1

g

1

21

90

()

(

)

(

)( )

()()()()

2

1

2

3

2

1

2

1

2

1

3

1

2

1

2

3

2

1

2

1

2

1

3

1

2

1

xxxx2x

2

1

xxxx2x

2

1

xxxx2x

2

1

xxxx2x

2

1

⎥

⎦

⎤

⎢

⎣

⎡

−++−+

⎥

⎦

⎤

⎢

⎣

⎡

−++−

⎥

⎦

⎤

⎢

⎣

⎡

−++−

⎥

⎦

⎤

⎢

⎣

⎡

−++−

=

ηη

(3.93)

Caso 2:

()

(

)

(

)

2

g

2

1

g

1

xxxxr −+−=

η

()

(

)

(

)

[

]

()

()()

[]

⎪

⎩

⎪

⎨

⎧

−++−=−

1

2

3

2

3

2

1

2

g

2

1

3

1

3

1

2

1

2

1

g

1

xxxx2x

2

1

xx

xxxx2x

2

1

xx

ηη

()

(

)

(

)

[

]

()()

[]

2

1

2

3

2

3

2

1

2

1

3

1

3

1

2

1

xxxx2x

2

1

r

−++−+

+−++−

=

η

ηη

(3.94)

() ()

(

)

()

(

)

()

=

−−

=

ηη

r

xx

r

xx

,r,r

2

1

g

1

2

1

g

1

11

(

)

(

)

[]

()()

[]

()()

[]

2

1

2

3

2

3

2

1

2

1

3

1

3

1

2

1

2

1

3

1

3

1

2

1

xxxx2xxxxx2x

xxxx2x

−++−+−++−

−++−

=

ηη

η

(3.95)

() ()

(

)

()

(

)

()

=

−−

=

ηη

r

xx

r

xx

,r,r

2

g

2

g

2

22

(

)

(

)

[

]

()()

[]

()()

[]

2

1

2

3

2

3

2

1

2

1

3

1

3

1

2

1

2

1

2

3

2

3

2

1

2

xxxx2xxxxx2x

xxxx2x

−++−+−++−

−++−

=

ηη

η

(3.96)

() ()

(

)

()

(

)

()

=

−−

=

ηη

r

xx

r

xx

,r,r

2

g

2

1

g

1

21

()()

[]

()()

[]

()()

[]

()()

[]

2

1

2

3

2

3

2

1

2

1

3

1

3

1

2

1

2

3

2

3

2

1

2

1

3

1

3

1

2

1

xxxx2xxxxx2x

−++−+−++−

−++−−++−

=

ηη

(3.97)

Caso 3:

()

(

)

(

)

2

3

2

g

2

3

1

g

1

xxxxr −+−=

η

()

()()

()

()()

⎪

⎩

⎪

⎨

⎧

⎥

⎦

⎤

⎢

⎣

⎡

−++−−=−

⎥

⎦

⎤

⎢

⎣

⎡

−++−−=−

2

3

2

3

2

1

2

3

2

g

2

1

3

1

3

1

2

1

3

1

g

1

xxxx2x

2

1

1xx

xxxx2x

2

1

1xx

ηη

91

()

()()

2

3

2

3

2

1

2

1

3

1

3

1

2

1

xxxx2x

2

1

xxxx2x

2

1

1r

⎥

⎦

⎤

⎢

⎣

⎡

−++−+

+

⎥

⎦

⎤

⎢

⎣

⎡

−++−

−=

η

ηη

(3.98)

() ()

(

)

()

(

)

()

=

−−

=

ηη

r

xx

r

xx

,r,r

3

1

g

1

3

1

g

1

11

()()

()()()()

2

3

2

3

2

1

2

1

3

1

3

1

2

1

2

1

3

1

3

1

2

1

xxxx2x

2

1

xxxx2x

2

1

xxxx2x

2

1

⎥

⎦

⎤

⎢

⎣

⎡

−++−+

⎥

⎦

⎤

⎢

⎣

⎡

−++−

⎥

⎦

⎤

⎢

⎣

⎡

−++−

=

ηη

η

(3.99)

() ()

(

)

()

(

)

()

=

−−

=

ηη

r

xx

r

xx

,r,r

3

2

g

2

3

2

g

2

22

()()

()()()()

2

3

2

3

2

1

2

1

3

1

3

1

2

1

2

3

2

3

2

1

2

xxxx2x

2

1

xxxx2x

2

1

xxxx2x

2

1

⎥

⎦

⎤

⎢

⎣

⎡

−++−+

⎥

⎦

⎤

⎢

⎣

⎡

−++−

⎥

⎦

⎤

⎢

⎣

⎡

−++−

=

ηη

η

(3.100)

() ()

(

)

()

(

)

()

=

−−

=

ηη

r

xx

r

xx

,r,r

3

2

g

2

3

1

g

1

21

()

(

)

(

)( )

()()()()

2

3

2

3

2

1

2

1

3

1

3

1

2

1

2

3

2

3

2

1

2

1

3

1

3

1

2

1

xxxx2x

2

1

xxxx2x

2

1

xxxx2x

2

1

xxxx2x

2

1

⎥

⎦

⎤

⎢

⎣

⎡

−++−+

⎥

⎦

⎤

⎢

⎣

⎡

−++−

⎥

⎦

⎤

⎢

⎣

⎡

−++−

⎥

⎦

⎤

⎢

⎣

⎡

−++−

=

ηη

(3.101)

Substituindo esses resultados na equação (3.89) resulta

Caso 1:

(

)

(

)

ηηη

f1r +=

() ()

(

)

()

++−=

∫

+

−

1

nij

cn

ij

dJf1lnN1CTEQ

ηηδηηη

α

() () () ()

∫

+

−

+

1

nji

dJ,r,rN2CTE

ηηηηη

α

()

() () ()() ()

+

⎥

⎦

⎤

⎢

⎣

⎡

++−=

∫∫

+

−

+

−

1

nij

1

nij

cn

ij

dJflnNdJ1lnN1CTEQ

ηηδηηηηδηη

αα

α

() () () ()

∫

+

−

+

1

nji

dJ,r,rN2CTE

ηηηηη

α

()

() () ()() ()

∫∫

+

−

+

−

−+−=

1

nij

1

nij

cn

ij

dJflnN1CTEdJ1lnN1CTEQ

ηηδηηηηδηη

αα

α

92

() () () ()

∫

+

−

+

1

nji

dJ,r,rN2CTE

ηηηηη

α

(3.102)

Caso 2:

() ()

ηηη

f

2

1

r =

() () ()

+

⎟

⎠

⎞

⎜

⎝

⎛

−=

∫

+

−

1

nij

cn

ij

dJf

2

1

lnN1CTEQ

ηηδηηη

α

() () () ()

∫

+

−

+

1

nji

dJ,r,rN2CTE

ηηηηη

α

()

() () () ()

+

⎥

⎦

⎤

⎢

⎣

⎡

⎟

⎠

⎞

⎜

⎝

⎛

+−=

∫∫

+

−

+

−

1

nij

1

nij

cn

ij

dJf

2

1

lnNdJlnN1CTEQ

ηηδηηηηδηη

αα

α

() () () ()

∫

+

−

+

1

nji

dJ,r,rN2CTE

ηηηηη

α

()

() () () ()

∫∫

+

−

+

−

⎟

⎠

⎞

⎜

⎝

⎛

−−=

1

nij

1

nij

cn

ij

dJf

2

1

lnN1CTEdJlnN1CTEQ

ηηδηηηηδηη

αα

α

() () () ()

∫

+

−

+

1

nji

dJ,r,rN2CTE

ηηηηη

α

()

(

)

() ()

()

∫∫

+

−

−−=

1

0

nij

0

1

nij

cn

ij

dJlnN1CTEdJlnN1CTEQ

ηηδηηηηδηη

αα

α

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

∫∫

+

−

+

−

+

⎟

⎠

⎞

⎜

⎝

⎛

−

1

nji

1

nij

dJ,r,rN2CTEdJf

2

1

lnN1CTE

ηηηηηηηδηη

αα

(3.103)

Caso 3:

(

)

(

)

ηηη

f1r −=

() ()

(

)

()

+−−=

∫

+

−

1

nij

cn

ij

dJf1lnN1CTEQ

ηηδηηη

α

() () () ()

∫

+

−

+

1

nji

dJ,r,rN2CTE

ηηηηη

α

()

() () ()() ()

+

⎥

⎦

⎤

⎢

⎣

⎡

+−−=

∫∫

+

−

+

−

1

nij

1

nij

cn

ij

dJflnNdJ1lnN1CTEQ

ηηδηηηηδηη

αα

α

() () () ()

∫

+

−

+

1

nji

dJ,r,rN2CTE

ηηηηη

α

()

() () ()() ()

∫∫

+

−

+

−

−−−=

1

nij

1

nij

cn

ij

dJflnN1CTEdJ1lnN1CTEQ

ηηδηηηηδηη

αα

α

() () () ()

∫

+

−

+

1

nji

dJ,r,rN2CTE

ηηηηη

α

(3.104)

93

As singularidades na primeira integral em (3.102), na primeira e na segunda

integrais em (3.103) e na primeira integral em (3.104) podem ser contornadas através

das transformações de coordenadas apresentadas na Figura 3.16, ou seja,

Caso 1:

2

1

1c

η

ξ

+

=⇒= (3.105)

Caso 2:

η

ξ

η

ξ

−

=

⇒= 'e2c (3.106)

Caso 3:

2

1

3c

η

ξ

−

=⇒= (3.107)

Consequentemente, novas funções de interpolação serão obtidas para cada

caso da posição do ponte fonte, isto é,

Caso 1:

()

(

)( )

112N

1

−−=

ξ

(

)

(

)

ξ

−

=

14N

2

(

)( )

12N

3

−=

ξ

(3.108)

Caso 2:

() ()

1

2

1

N

1

−=

ξξξ

(

)

2

1N

ξξ

−=

() ()

1

2

1

N

3

+=

ξξξ

(3.109)

() ()

1''

2

1

'N

1

+=

ξξξ

(

)

(

)

2

'1'N

ξξ

−=

() ()

1''

2

1

'N

3

−=

ξξξ

(3.110)

Caso 3:

() ( )

12N

1

−=

ξ

(

)

(

)

ξ

−

=

14N

2

(

)

(

)( )

112N

3

−−

=

ξ

(3.111)

Também são encontradas novas relações (3.83) e (3.84) para cada caso da

posição do ponto fonte.

Caso 1:

()

(

)

()( )

⎪

⎩

⎪

⎨

⎧

+−+−++−=

1

2

3

2

1

2

3

2

1

2

g

2

1

3

1

2

1

3

1

2

1

2

g

1

xxx4x3xx2x2x

ξξ

(3.112)

()

(

)

(

)

[

]

()( )

[]

2

3

2

1

2

3

2

1

2

3

1

2

1

3

1

2

1

2

1

n

xx4x3xx2x4

d

dx

d

dx

J

−+−++−+

+−+−++−

=

⎟

⎠

⎞

⎜

⎝

⎛

+

⎟

⎠

⎞

⎜

⎝

⎛

=

ξ

ξξ

ξ

(3.113)

94

Caso 2:

(

)

(

)

()()

⎪

⎩

⎪

⎨

⎧

+−++−=

2

1

2

3

2

3

2

1

2

g

2

1

3

1

3

1

2

1

2

g

1

xxx

2

1

xx2x

2

1

x

xxx

2

1

xx2x

2

1

x

ξξ

(3.114)

()

()()

()

2

1

2

3

2

3

2

1

2

1

3

1

3

1

2

1

2

1

n

xx

2

1

xx2x

xx

2

1

xx2x

d

dx

d

dx

J

⎥

⎦

⎤

⎢

⎣

⎡

−++−+

+

⎥

⎦

⎤

⎢

⎣

⎡

−++−

=

⎟

⎠

⎞

⎜

⎝

⎛

+

⎟

⎠

⎞

⎜

⎝

⎛

=

ξ

ξξ

ξ

(3.115)

()

(

)

(

)

()

()()

⎪

⎩

⎪

⎨

⎧

+−−+−=

2

1

2

3

2

3

2

1

2

g

2

1

3

1

3

1

2

1

2

g

1

xxx'

2

1

xx2x'

2

1

x

xxx'

2

1

xx2x'

2

1

x

ξξ

(3.116)

()

()()

()

2

1

2

3

2

3

2

1

2

1

3

1

3

1

2

1

2

1

n

xx

2

1

xx2x'

xx

2

1

xx2x'

'd

dx

'd

dx

'J

⎥

⎦

⎤

⎢

⎣

⎡

−−+−+

+

⎥

⎦

⎤

⎢

⎣

⎡

−−+−

=

⎟

⎠

⎞

⎜

⎝

⎛

+

⎟

⎠

⎞

⎜

⎝

⎛

=

ξ

ξξ

ξ

(3.117)

Caso 3:

()

(

)

()( )

⎪

⎩

⎪

⎨

⎧

+−+−++−=

3

2

3

2

1

2

3

2

1

2

g

2

3

1

3

1

2

1

3

1

2

1

2

g

1

xx3x4xxx2x2x

ξξ

(3.118)

()

(

)

(

)

[

]

()( )

[]

2

3

2

1

2

3

2

1

2

3

1

2

1

3

1

2

1

2

1

n

x3x4xxx2x4

d

dx

d

dx

J

−+−++−+

+−+−++−

=

⎟

⎠

⎞

⎜

⎝

⎛

+

⎟

⎠

⎞

⎜

⎝

⎛

=

ξ

ξξ

ξ

(3.119)

Substituindo (3.105), (3.106) e (3.107) em (3.102), (3.103) e (3.104),

respectivamente, obtém-se

Caso 1:

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

∫∫

+

−

+

−−=

1

nij

1

0

nij

cn

ij

dJflnN1CTEdJ2lnN1CTEQ

ηηδηηξξδξξ

αα

α

() () () ()

∫

+

−

+

1

nji

dJ,r,rN2CTE

ηηηηη

α

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

∫∫

+

−

+

−−=

1

nij

1

0

nij

cn

ij

dJf2lnN1CTEdJlnN1CTEQ

ηηδηηξξδξξ

αα

α

() () () ()

∫

+

−

+

1

nji

dJ,r,rN2CTE

ηηηηη

α

95

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

∫∫

+

−

+

−

⎟

⎠

⎞

⎜

⎝

⎛

=

1

nij

1

0

nij

cn

ij

dJf2lnN1CTEdJ

1

lnN1CTEQ

ηηδηηξξδ

ξ

αα

α

() () () ()

∫

+

−

+

1

nji

dJ,r,rN2CTE

ηηηηη

α

(3.120)

Caso 2:

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

∫∫

+

−−=

1

0

nij

0

1

nij

cn

ij

dJlnN1CTE'd'J'ln'N1CTEQ

ξξδξξξξδξξ

αα

α

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

∫∫

+

−

+

−

+

⎟

⎠

⎞

⎜

⎝

⎛

−

1

nji

1

nij

dJ,r,rN2CTEdJf

2

1

lnN1CTE

ηηηηηηηδηη

αα

() ( ) () () ()

∫∫

+

−

⎟

⎠

⎞

⎜

⎝

⎛

=

1

0

nij

0

1

nij

cn

ij

dJlnN1CTE'd'J

'

1

ln'N1CTEQ

ξξδξξξξδ

ξ

αα

α

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

∫∫

+

−

+

−

+

⎟

⎠

⎞

⎜

⎝

⎛

−

1

nji

1

nij

dJ,r,rN2CTEdJf

2

1

lnN1CTE

ηηηηηηηδηη

αα

() () () ()

∫∫

+

⎟

⎠

⎞

⎜

⎝

⎛

+

⎟

⎠

⎞

⎜

⎝

⎛

=

1

0

nij

0

1

nij

cn

ij

dJ

1

lnN1CTE'd'J

'

1

ln'N1CTEQ

ξξδ

ξ

ξξξδ

ξ

αα

α

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

∫∫

+

−

+

−

+

⎟

⎠

⎞

⎜

⎝

⎛

−

1

nji

1

nij

dJ,r,rN2CTEdJf

2

1

lnN1CTE

ηηηηηηηδηη

αα

(3.121)

Caso 3:

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

∫∫

+

−+

−−=

1

nij

0

1

nij

cn

ij

dJflnN1CTEdJ2lnN1CTEQ

ηηδηηξξδξξ

αα

α

() () () ()

∫

+

−

+

1

nji

dJ,r,rN2CTE

ηηηηη

α

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

∫∫

+

−+

−−=

1

nij

0

1

nij

cn

ij

dJf2lnN1CTEdJlnN1CTEQ

ηηδηηξξδξξ

αα

α

() () () ()

∫

+

−

+

1

nji

dJ,r,rN2CTE

ηηηηη

α

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

∫∫

+

−+

−

⎟

⎠

⎞

⎜

⎝

⎛

=

1

nij

0

1

nij

cn

ij

dJf2lnN1CTEdJ

1

lnN1CTEQ

ηηδηηξξδ

ξ

αα

α

() () () ()

∫

+

−

+

1

nji

dJ,r,rN2CTE

ηηηηη

α

(3.122)

A primeira integral em (3.120), a primeira e a segunda integrais em (3.121) e

a primeira integral em

(3.122) são avaliadas por meio de uma quadratura de Gauss

modificada, denominada quadratura de Gauss logarítmica (para valores calculados e

tabelados de pontos e pesos de Gauss para a quadratura Gaussiana logarítmica ver,

96

por exemplo, BREBBIA e DOMINGUEZ

15

). Como mencionado para o caso de

elementos lineares, essas integrais poderiam, alternativamente, ser calculadas

analiticamente.

3.8.2.4 Elementos Quadráticos – Coeficientes

α

cn

ij

P

A expressão (3.39) pode ser reescrita como

() ()

[]

()

∫

+

−

=

1

n

cn

ij

cn

ij

dJrTNP

ηηηη

α

(3.123)

sendo que

()

[]

()()

(

)

()

[

]

()

⎭

⎬

⎫

⎩

⎨

⎧

−−++−

∂

−

=

ijjijiij

cn

ij

,rn,rn21,r,r221

n

r

r14

1

rT

νδν

η

ηνπ

η

ou

()

[]

()

[]

η

rg

r

1

rT

cn

ij

=

onde g [r (η)] é, nitidamente, não singular.

Verificando o comportamento da expressão (3.123), para cada caso da

posição do ponto fonte, observa-se:

Caso 1:

(

)

(

)

ηηη

f1r +=

()

[]

() ()

()

[]

()

=

+

⎟

⎠

⎞

⎜

⎝

⎛

−==

∫∫

+

−

+

−

1

n

1

n1

1cn

ij

dJrg

f1

1

2

1

dJrg

r

1

NP

ηηη

ηη

ηηηηη

η

(

)

()

[]

()

∫

+

−

+

−

=

1

n

dJrg

f1

1

2

1

ηηη

ηη

η

()

[]

() ( )( )()

()

[]

()

=

+

−+==

∫∫

+

−

+

−

1

n

1

n2

2cn

ij

dJrg

f1

1

11dJrg

r

1

NP

ηηη

ηη

ηηηηη

η

(

)

()

[]

()

∫

+

−

=

1

n

dJrg

f

1

ηηη

η

15

BREBBIA, C. A., DOMINGUEZ, J. Boundary Elements, An Introductory Course. Boston:

Computational Mechanics Publications, 1989, p. 271.

97

()

[]

() ()

()

[]

()

=

+

⎟

⎠

⎞

⎜

⎝

⎛

+==

∫∫

+

−

+

−

1

n

1

n3

3cn

ij

dJrg

f1

1

2

1

dJrg

r

1

NP

ηηη

ηη

ηηηηη

η

()

[]

()

∫

+

−

=

1

n

dJrg

f2

1

ηηη

η

(Conclusão: apenas a primeira expressão apresenta singularidade).

Caso 2:

() ()

ηηη

f

2

1

r =

()

[]

() ()

()

[]

()

=

⎟

⎠

⎞

⎜

⎝

⎛

⎟

⎠

⎞

⎜

⎝

⎛

−==

∫∫

+

−

+

−

1

n

1

n1

1cn

ij

dJrg

f

2

1

2

1

dJrg

r

1

NP

ηηη

ηη

ηηηηη

η

(

)

()

[]

()

∫

+

−

=

1

n

dJrg

f

1

ηηη

η

()

[]

() ( )( )()

()

[]

()

=

⎟

⎠

⎞

⎜

⎝

⎛

−+==

∫∫

+

−

+

−

1

n

1

n2

2cn

ij

dJrg

f

2

1

11dJrg

r

1

NP

ηηη

ηη

ηηηηη

η

(

)

(

)

()

[]

()

∫

+

−

+

=

1

n

dJrg

f

11

2

ηηη

ηη

η

()

[]

() ()

()

[]

()

=

⎟

⎠

⎞

⎜

⎝

⎛

⎟

⎠

⎞

⎜

⎝

⎛

+==

∫∫

+

−

+

−

1

n

1

n3

3cn

ij

dJrg

f

2

1

2

1

dJrg

r

1

NP

ηηη

ηη

ηηηηη

η

(

)

()

[]

()

∫

+

−

+

=

1

n

dJrg

f

1

ηηη

η

(Conclusão: apenas a segunda expressão apresenta singularidade).

Caso 3:

(

)

(

)

ηηη

f1r −=

()

[]

() ()

()

[]

()

=

−

⎟

⎠

⎞

⎜

⎝

⎛

−==

∫∫

+

−

+

−

1

n

1

n1

1cn

ij

dJrg

f1

1

2

1

dJrg

r

1

NP

ηηη

ηη

ηηηηη

η

()

[]

()

∫

+

−

=

1

n

dJrg

f2

1

ηηη

η

()

[]

() ( )( )()

()

[]

()

=

−

−+==

∫∫

+

−

+

−

1

n

1

n2

2cn

ij

dJrg

f1

1

11dJrg

r

1

NP

ηηη

ηη

ηηηηη

η

(

)

()

[]

()

∫

+

−

+

=

1

n

dJrg

f

1

ηηη

η

98

()

[]

() ()

()

[]

()

=

−

⎟

⎠

⎞

⎜

⎝

⎛

+==

∫∫

+

−

+

−

1

n

1

n3

3cn

ij

dJrg

f1

1

2

1

dJrg

r

1

NP

ηηη

ηη

ηηηηη

η

(

)

()

[]

()

∫

+

−

+

=

1

n

dJrg

f1

1

2

1

ηηη

ηη

η

(Conclusão: apenas a terceira expressão apresenta singularidade).

Os coeficientes

α

cn

ij

P singulares são exatamente os coeficientes que

compõem as submatrizes diagonais de H, as quais são calculadas indiretamente

através de consideração de movimento de corpo rígido, como discutido

anteriormente. A conclusão é que, para efeitos de implementação computacional, a

avaliação dos coeficientes

α

cn

ij

P , quando c é um dos nós do elemento n, pode ser

feita através da quadratura de Gauss convencional, pois os coeficientes singulares,

calculados imprecisamente através desse processo, são desprezados e

recalculados em função dos demais coeficientes, resultantes de integrações não

singulares e, portanto, precisamente avaliados.

3.8.3 Deslocamentos e Tensões nos Pontos Internos

É possível obter os valores dos deslocamentos e tensões nos pontos

internos a partir dos valores dos deslocamentos e forças de superfície do contorno

utilizando-se as identidades de Somigliana para deslocamentos, (3.6), e tensões,

(3.8), que podem ser discretizadas de maneira similar à expressão (3.43).

Para um ponto qualquer

X

d

∈

Ω, a versão discretizada de (3.6) é

∑∑∑

===

−=

Nn

1

j

d

ij

Ne

1n

m

1

n

j

dn

ij

d

i

uHtGu

γ

γγ

α

αα

(3.124)

e a versão discretizada de (3.8) é

∑∑∑

===

−=

Nn

1

k

d

kij

Ne

1n

m

1

n

k

dn

kij

d

ij

uStD

γ

γγ

α

αα

σ

(3.125)

99

4 CÓDIGOS COMPUTACIONAIS

4.1 Elementos Constantes

As considerações matemáticas relativas ao desenvolvimento do Método de

Elementos de Contorno utilizando elementos de contorno constantes, apresentadas

no capítulo precedente, foram implementadas no código denominado

MEC_ELASTIC_CONST, escrito em linguagem de programação Fortran 95, que consta

do CD anexo a este trabalho. Esse código é formado pelos seguintes arquivos fonte:

a) MEC_ELASTIC_CONST.f90;

b) MEC_ELASTIC_CONST_SUBROTINAS.f90.

4.1.1 Organização e Funcionamento do Código

A organização e o funcionamento desse código podem ser resumidos nas

seguintes etapas:

a) declaração do nome do programa e dos módulos e bibliotecas utilizadas;

100

b) declaração das variáveis globais e definições;

c) atribuição do número dos arquivos de entrada e saída de dados;

d) leitura do nome dos arquivos de entrada e saída de dados e abertura dos

mesmos (em caso de erro o programa encerra sua execução alertando o

usuário);

e) chamada da sub-rotina INPUT1 para leitura dos dados contidos no arquivo

de entrada de dados: leitura da linha de título do problema e linha de

parâmetros globais do problema;

f) alocação de memória para as variáveis do problema (em caso de erro o

programa encerra sua execução alertando o usuário);

g) chamada da sub-rotina INPUT2 para leitura dos demais dados contidos no

arquivo de entrada de dados;

h) chamada da sub-rotina MATRIZES_G_H para cálculo dos coeficientes das

matrizes G e H e formação do sistema AY = F;

 chamada da sub-rotina INTEGRACAO_LOCAL para cálculo dos coeficientes

das submatrizes diagonais de G;

 chamada da sub-rotina INTEGRACAO_EXTERNA para cálculo dos

coeficientes das submatrizes posicionadas fora das diagonais de G e H;

i) chamada da sub-rotina SOLVER para solução do sistema de equações

através do método de Gauss (em caso de erro o programa encerra sua

execução alertando o usuário);

j) chamada da sub-rotina

SOLUCAO_CONTORNO para reordenação dos valores

calculados no contorno do problema;

k) chamada da sub-rotina

PONTOS_INTERNOS para cálculo dos deslocamentos

e tensões nos pontos internos;

 chamada da sub-rotina INTEGRACAO_EXTERNA para cálculo dos

coeficientes necessários para a obtenção dos deslocamentos nos

pontos internos;

 chamada da sub-rotina TENSOES_INTERNAS para cálculo dos coeficientes

necessários para a obtenção das tensões nos pontos internos;

l) chamada da sub-rotina OUTPUT para transcrição dos resultados no arquivo

de saída de dados;

m) fechamento dos arquivos de entrada e saída de dados.

101

4.1.2 Modelo de Arquivo de Entrada de Dados

Um exemplo de arquivo de entrada de dados para o código

MEC_ELASTIC_CONST é apresentado na Figura 4.1 (e ilustrado na Figura 4.2),

constituindo-se dos seguintes campos:

*** CAVIDADE CIRCULAR SOB PRESSAO CONSTANTE (12 ELEMENTOS CONSTANTES) ***

0 0 12 7 10 20500.0 0.3 1.0E-10 1.0E-3 1.0E-3

1 -0.7765 -2.8978

2 -2.1213 -2.1213

3 -2.8978 -0.7765

4 -2.8978 0.7765

5 -2.1213 2.1213

6 -0.7765 2.8978

7 0.7765 2.8978

8 2.1213 2.1213

9 2.8978 0.7765

10 2.8978 -0.7765

11 2.1213 -2.1213

12 0.7765 -2.8978

1 1 -50.0000 1 -86.6025

2 1 -86.6025 1 -50.0000

3 1 -100.0000 0 0.0000

4 1 -86.6025 1 50.0000

5 1 -50.0000 1 86.6025

6 0 0.0000 1 100.0000

7 1 50.0000 1 86.6025

8 1 86.6025 1 50.0000

9 1 100.0000 0 0.0000

10 1 86.6025 1 -50.0000

11 1 50.0000 1 -86.6025

12 0 0.0000 1 -100.0000

1 4.0000 0.0000

2 5.0000 0.0000

3 6.0000 0.0000

4 7.0000 0.0000

5 8.0000 0.0000

6 9.0000 0.0000

7 10.0000 0.0000

FIGURA 4.1 – Modelo de arquivo de entrada de dados do código MEC_ELASTIC_CONST

a) linha de título do problema: contém uma breve descrição do problema

analisado (máximo 100 caracteres);

b) linha de parâmetros globais do problema: índice do tipo de domínio do

problema (

0 para domínio infinito e 1 para domínio finito); índice do estado

plano adotado (

0 para Estado Plano de Deformações e 1 para Estado

Plano de Tensões); número de elementos de contorno; número de pontos

internos; número de pontos de Gauss para integração numérica (valores

possíveis: 2, 4, 8 ou 10); módulo de elasticidade longitudinal ou módulo de

Young; coeficiente de Poisson; valor de tolerância adotado na solução do

sistema de equações para determinar coeficientes nulos; valor de

tolerância adotado para o arredondamento dos valores de deslocamentos

102

calculados; valor de tolerância adotado para o arredondamento dos

valores de tensões e forças de superfície calculados;

c) bloco de coordenadas dos pontos extremos dos elementos de contorno:

número do ponto; coordenada X; coordenada Y;

d) bloco de condições de contorno prescritas para os nós de contorno:

número do nó; índice do tipo de condição de contorno prescrita para a

direção X (0 para deslocamentos prescritos e 1 para forças de superfície

prescritas); valor da condição de contorno prescrita para a direção X;

índice do tipo de condição de contorno prescrita para a direção Y (0 para

deslocamentos prescritos e 1 para forças de superfície prescritas); valor

da condição de contorno prescrita para a direção Y;

e) bloco de coordenadas dos pontos internos: número do ponto; coordenada

X; coordenada Y.

FIGURA 4.2 – Ilustração da malha do modelo de arquivo de entrada de dados

do código MEC_ELASTIC_CONST

4.2 Elementos Lineares

As considerações matemáticas relativas ao desenvolvimento do Método de

Elementos de Contorno utilizando elementos de contorno lineares, apresentadas no

103

capítulo precedente, foram implementadas no código denominado MEC_ELASTIC_LIN,

escrito em linguagem de programação Fortran 95, que consta do CD anexo a este

trabalho. Esse código é formado pelos seguintes arquivos fonte:

a) MEC_ELASTIC_LIN.f90;

b) MEC_ELASTIC_LIN_SUBROTINAS.f90.

4.2.1 Organização e Funcionamento do Código

A organização e o funcionamento desse código podem ser resumidos nas

seguintes etapas:

a) declaração do nome do programa e dos módulos e bibliotecas utilizadas;

b) declaração das variáveis globais e definições;

c) atribuição do número dos arquivos de entrada e saída de dados;

d) leitura do nome dos arquivos de entrada e saída de dados e abertura dos

mesmos (em caso de erro o programa encerra sua execução alertando o

usuário);

e) chamada da sub-rotina INPUT1 para leitura dos dados contidos no arquivo

de entrada de dados: leitura da linha de título do problema e linha de

parâmetros globais do problema;

f) alocação de memória para as variáveis do problema (em caso de erro o

programa encerra sua execução alertando o usuário);

g) chamada da sub-rotina INPUT2 para leitura dos demais dados contidos no

arquivo de entrada de dados;

h) chamada da sub-rotina MATRIZES_G_H para cálculo dos coeficientes das

matrizes G e H e formação do sistema AY = F;

 chamada da sub-rotina INTEGRACAO_LOCAL_P para cálculo dos

coeficientes das submatrizes P quando o ponto fonte é um dos nós do

elemento que está sendo integrado;

 chamada da sub-rotina

INTEGRACAO_LOCAL_Q para cálculo dos

coeficientes das submatrizes

Q quando o ponto fonte é um dos nós do

elemento que está sendo integrado;

104

 chamada da sub-rotina INTEGRACAO_EXTERNA para cálculo dos

coeficientes das submatrizes P e Q quando o ponto fonte não é um dos

nós do elemento que está sendo integrado;

i) chamada da sub-rotina SOLVER para solução do sistema de equações

através do método de Gauss (em caso de erro o programa encerra sua

execução alertando o usuário);

j) chamada da sub-rotina SOLUCAO_CONTORNO para reordenação dos valores

calculados no contorno do problema;

k) chamada da sub-rotina PONTOS_INTERNOS para cálculo dos deslocamentos

e tensões nos pontos internos;

 chamada da sub-rotina INTEGRACAO_EXTERNA para cálculo dos

coeficientes necessários para a obtenção dos deslocamentos nos

pontos internos;

 chamada da sub-rotina TENSOES_INTERNAS para cálculo dos coeficientes

necessários para a obtenção das tensões nos pontos internos;

l) chamada da sub-rotina OUTPUT para transcrição dos resultados no arquivo

de saída de dados;

m) fechamento dos arquivos de entrada e saída de dados.

4.2.2 Modelo de Arquivo de Entrada de Dados

Um exemplo de arquivo de entrada de dados para o código MEC_ELASTIC_LIN

é apresentado na Figura 4.3 (e ilustrado na Figura 4.4), constituindo-se dos

seguintes campos:

a) linha de título do problema: contém uma breve descrição do problema

analisado (máximo 100 caracteres);

b) linha de parâmetros globais do problema: índice do tipo de domínio do

problema (0 para domínio infinito e 1 para domínio finito); índice do estado

plano adotado (0 para Estado Plano de Deformações e 1 para Estado

Plano de Tensões); número de elementos de contorno; número de pontos

internos; número de pontos de Gauss para integração numérica (valores

105

possíveis: 2, 4, 8 ou 10); módulo de elasticidade longitudinal ou módulo de

Young; coeficiente de Poisson; valor de tolerância adotado na solução do

sistema de equações para determinar coeficientes nulos; valor de

tolerância adotado para o arredondamento dos valores de deslocamentos

calculados; valor de tolerância adotado para o arredondamento dos

valores de tensões e forças de superfície calculados;

*** CAVIDADE CIRCULAR SOB PRESSAO CONSTANTE (12 ELEMENTOS LINEARES) ***

0 0 12 7 10 20500.0 0.3 1.0E-10 1.0E-3 1.0E-3

1 -1.5000 -2.5981

2 -2.5981 -1.5000

3 -3.0000 0.0000

5 -2.5981 1.5000

5 -1.5000 2.5981

6 0.0000 3.0000

7 1.5000 2.5981

8 2.5981 1.5000

9 3.0000 0.0000

10 2.5981 -1.5000

11 1.5000 -2.5981

12 0.0000 -3.0000

1 1 -50.0000 1 -86.6025 1 -86.6025 1 -50.0000

2 1 -86.6025 1 -50.0000 1 -100.0000 0 0.0000

3 1 -100.0000 0 0.0000 1 -86.6025 1 50.0000

4 1 -86.6025 1 50.0000 1 -50.0000 1 86.6025

5 1 -50.0000 1 86.6025 0 0.0000 1 100.0000

6 0 0.0000 1 100.0000 1 50.0000 1 86.6025

7 1 50.0000 1 86.6025 1 86.6025 1 50.0000

8 1 86.6025 1 50.0000 1 100.0000 0 0.0000

9 1 100.0000 0 0.0000 1 86.6025 1 -50.0000

10 1 86.6025 1 -50.0000 1 50.0000 1 -86.6025

11 1 50.0000 1 -86.6025 0 0.0000 1 -100.0000

12 0 0.0000 1 -100.0000 1 -50.0000 1 -86.6025

1 4.0000 0.0000

2 5.0000 0.0000

3 6.0000 0.0000

4 7.0000 0.0000

5 8.0000 0.0000

6 9.0000 0.0000

7 10.0000 0.0000

FIGURA 4.3 – Modelo de arquivo de entrada de dados do código MEC_ELASTIC_LIN

c) bloco de coordenadas dos nós de contorno: número do nó; coordenada X;

coordenada Y;

d) bloco de condições de contorno prescritas para os elementos de

contorno: número do elemento; índice do tipo de condição de contorno

prescrita para a direção X do primeiro nó do elemento (0 para

deslocamentos prescritos e 1 para forças de superfície prescritas); valor

da condição de contorno prescrita para a direção X do primeiro nó do

elemento; índice do tipo de condição de contorno prescrita para a direção

Y do primeiro nó do elemento (0 para deslocamentos prescritos e 1 para

forças de superfície prescritas); valor da condição de contorno prescrita

para a direção

Y do primeiro nó do elemento; índice do tipo de condição

de contorno prescrita para a direção

X do segundo nó do elemento (0 para

106

deslocamentos prescritos e 1 para forças de superfície prescritas); valor

da condição de contorno prescrita para a direção X do segundo nó do

elemento; índice do tipo de condição de contorno prescrita para a direção

Y do segundo nó do elemento (0 para deslocamentos prescritos e 1 para

forças de superfície prescritas); valor da condição de contorno prescrita

para a direção Y do segundo nó do elemento;

e) bloco de coordenadas dos pontos internos: número do ponto; coordenada

X; coordenada Y.

FIGURA 4.4 – Ilustração da malha do modelo de arquivo de entrada de dados

do código MEC_ELASTIC_LIN

4.3 Elementos Quadráticos

As considerações matemáticas relativas ao desenvolvimento do Método de

Elementos de Contorno utilizando elementos de contorno quadráticos, apresentadas

no capítulo precedente, foram implementadas no código denominado

MEC_ELASTIC_QUAD, escrito em linguagem de programação Fortran 95, que consta do

CD anexo a este trabalho. Esse código é formado pelos seguintes arquivos fonte:

a)

MEC_ELASTIC_QUAD.f90;

b) MEC_ELASTIC_QUAD_SUBROTINAS.f90.

107

4.3.1 Organização e Funcionamento do Código

A organização e o funcionamento desse código podem ser resumidos nas

seguintes etapas:

a) declaração do nome do programa e dos módulos e bibliotecas utilizadas;

b) declaração das variáveis globais e definições;

c) atribuição do número dos arquivos de entrada e saída de dados;

d) leitura do nome dos arquivos de entrada e saída de dados e abertura dos

mesmos (em caso de erro o programa encerra sua execução alertando o

usuário);

e) chamada da sub-rotina INPUT1 para leitura dos dados contidos no arquivo

de entrada de dados: leitura da linha de título do problema e linha de

parâmetros globais do problema;

f) alocação de memória para as variáveis do problema (em caso de erro o

programa encerra sua execução alertando o usuário);

g) chamada da sub-rotina INPUT2 para leitura dos demais dados contidos no

arquivo de entrada de dados;

h) chamada da sub-rotina MATRIZES_G_H para cálculo dos coeficientes das

matrizes G e H e formação do sistema AY = F;

 chamada da sub-rotina INTEGRACAO_LOCAL_P para cálculo dos

coeficientes das submatrizes

P quando o ponto fonte é um dos nós do

elemento que está sendo integrado;

 chamada da sub-rotina

INTEGRACAO_LOCAL_Q para cálculo dos

coeficientes das submatrizes Q quando o ponto fonte é um dos nós do

elemento que está sendo integrado;

 chamada da sub-rotina INTEGRACAO_EXTERNA para cálculo dos

coeficientes das submatrizes P e Q quando o ponto fonte não é um dos

nós do elemento que está sendo integrado;

i) chamada da sub-rotina SOLVER para solução do sistema de equações

através do método de Gauss (em caso de erro o programa encerra sua

execução alertando o usuário);

108

j) chamada da sub-rotina SOLUCAO_CONTORNO para reordenação dos valores

calculados no contorno do problema;

k) chamada da sub-rotina PONTOS_INTERNOS para cálculo dos deslocamentos

e tensões nos pontos internos;

 chamada da sub-rotina INTEGRACAO_EXTERNA para cálculo dos

coeficientes necessários para a obtenção dos deslocamentos nos

pontos internos;

 chamada da sub-rotina TENSOES_INTERNAS para cálculo dos coeficientes

necessários para a obtenção das tensões nos pontos internos;

l) chamada da sub-rotina OUTPUT para transcrição dos resultados no arquivo

de saída de dados;

m) fechamento dos arquivos de entrada e saída de dados.

4.3.2 Modelo de Arquivo de Entrada de Dados

Um exemplo de arquivo de entrada de dados para o código

MEC_ELASTIC_QUAD é apresentado na Figura 4.5 (e ilustrado na Figura 4.6),