Análise dinâmica de uma coluna de perfuração …...de Moraes Lobo, Daniel An alise din^amica de uma coluna de perfura˘c~ao durante a mudan˘ca de caracter sticas da rocha com

ANALISE DINAMICA DE UMA COLUNA DE PERFURACAO DURANTE A

MUDANCA DE CARACTERISTICAS DA ROCHA COM QUANTIFICACAO

DE INCERTEZAS

Daniel de Moraes Lobo

Projeto de Graduacao apresentado ao Curso

de Engenharia Mecanica da Escola Politecnica,

Universidade Federal do Rio de Janeiro, como

parte dos requisitos necessarios a obtencao do

tıtulo de Engenheiro.

Orientadores: Thiago Gamboa Ritto

Daniel Alves Castello

Rio de Janeiro

Julho de 2016

UNIVERSIDADE FEDERAL DO RIO DE JANEIRO

Departamento de Engenharia Mecnica

DEM/POLI/UFRJ



DE INCERTEZAS


PROJETO FINAL SUBMETIDO AO CORPO DOCENTE DO DEPARTAMENTO

DE ENGENHARIA MECANICA DA ESCOLA POLITECNICA DA

UNIVERSIDADE FEDERAL DO RIO DE JANEIRO COMO PARTE

DOS REQUISITOS NECESSARIOS PARA A OBTENCAO DO GRAU DE

ENGENHEIRO MECANICO.

Aprovada por:

Prof. Thiago Gamboa Ritto, D.Sc.

Prof. Daniel Alves Castello, D.Sc.

Prof. Lavinia Maria Sanabio Alves Borges, D.Sc.

Prof. Fernando Augusto de Noronha Castro Pinto, Dr.Ing.

RIO DE JANEIRO, RJ – BRASIL

JULHO DE 2016

de Moraes Lobo, Daniel

Analise dinamica de uma coluna de perfuracao durante

a mudanca de caracterısticas da rocha com quantificacao

de incertezas/ Daniel de Moraes Lobo. – Rio de Janeiro:

UFRJ/Escola Politecnica, 2016.

XIII, 72 p.: il.; 29, 7cm.



Projeto de Graduacao – UFRJ/ Escola Politecnica/

Curso de Engenharia Mecanica, 2016.

Referencias Bibliograficas: p. 41 – 42.

1. Coluna de perfuracao. 2. Vibracao. 3. Stick-Slip.

4. Mudanca de rocha. I. Gamboa Ritto, Thiago et al..

II. Universidade Federal do Rio de Janeiro, UFRJ, Curso

de Engenharia Mecanica. III. Analise dinamica de uma

coluna de perfuracao durante a mudanca de caracterısticas

da rocha com quantificacao de incertezas.

iii

”Voce precisa fazer aquilo que

pensa que nao e capaz de fazer.”

Eleanor Roosevelt

iv

Agradecimentos

O meu sucesso hoje atingido nao seria possıvel sem o apoio da minha famılia.

A minha noiva por todo o amor, carinho e apoio incondicional. Obrigado por

ter acreditado em mim sempre, por ter me consolado nos momentos difıceis e come-

morado minhas conquistas. Eu nao conseguiria metade do que consegui sem o seu

apoio.

Aos meu pais pelo amor, incentivo e suporte. Gostaria de agradecer tambem por

todo o investimento feito na minha educacao desde pequeno e pela valorizacao da

minha formacao. O apoio deles foi fundamental para mim.

Aos meus orientadores, Thiago Ritto e Daniel Castello, pela confianca, orientacao

e incentivo durante a construcao desse trabalho. Gostaria de agradecer principal-

mente pela orientacao profissional que me permitiu identificar um futuro com mais

clareza.

As professoras Lavınia e Anna Carla pelas orientacoes dadas durante esse meu

percurso e por sempre me receberem bem.

Ao professo Silvio Carlos, do qual fui monitor por 2 anos. Obrigado pelas ori-

entacoes academicas, pessoais e profissionais. Agradeco ainda o suporte, confianca

e incentivo que me foi dado durante todo esse perıodo trabalhando juntos.

Aos amigos que fiz durante meu estagio na Vallourec. Seria inviavel listar todos

os nomes importantes para mim aqui, entao agradeco a todos pelas orientacoes,

trabalhos e oportunidades que me foram dadas durante o estagio. Ainda assim,

gostaria de agradecer especialmente ao Pedro Filgueiras que teve papel fundamental

para meu desenvolvimento dentro da empresa, pela confianca depositada em mim

nesse tempo e por todas as oportunidade que me foram abertas durante esse perıodo.

Agradecp especialmente tambem aos demais estagiarios da equipe a qual pertencia.

Mais do que colegas de trabalho, sei que fiz amigos; amigos de verdade. Obrigado!

v

Resumo do Projeto de Graduacao apresentado a Escola Politecnica/UFRJ como

parte dos requisitos necessarios para a obtencao do grau de Engenheiro Mecanico



DE INCERTEZAS


Julho/2016



Programa: Engenharia Mecanica

O processo de perfuracao de pocos de petroleo e um dos mais caros e complexos

presentes na industria. O entendimento dos mecanismos complexos de vibracao da

coluna de perfuracao e crucial para melhorar o controle da operacao e performance.

Neste trabalho, foi proposto um modelo simples com um grau de liberdade para

avaliar a vibracao torcional da coluna. Devido a diversas fontes de incerteza no pro-

cesso, o metodo de Monte Carlo foi usado de forma a considerar essas incertezas nas

analises. A resposta dinamica do sistema foi simulada e um estudo de sensibilidade

das variaveis peso sobre a broca (WOB) e velocidade da mesa rotativa foi realizado,

introduzindo uma proposta de analise atraves de graficos de bifurcacao. Incertezas

foram incluıdas no parametro de peso sobre a broca e analisou-se a probabilidade de

Stick-Slip para cada velocidade de rotacao da mesa, mostrando que velocidades mais

altas resultam em menores probabilidade de Stick-Slip. Em seguida, foram aplica-

das funcoes aos torques de atrito broca-rocha com o objetivo de simular o primeiro

contato broca-rocha e a alteracao de caracterısticas da rocha durante o processo de

perfuracao, observando-se comportamentos caracterısticos para cada funcao. Por

fim, incertezas foram adicionadas a analise para avaliar o impacto para cada funcao

e mostrou-se que a maior influencia ocorre no caso da funcao degrau.

vi

Abstract of Undergraduate Project presented to POLI/UFRJ as a partial fulfillment

of the requirements for the degree of Mechanical Engineer

DYNAMIC ANALYSIS OF A DRILLING STRING DURING THE CHANGE OF

ROCK PROPERTIES WITH UNCERTAINTIES QUANTIFICATION


July/2016

Advisors: Thiago Gamboa Ritto


Department: Mechanical Engineering

The oil well drilling process is one of the most expensive and complex in the in-

dustry. The understanding of the complex mechanisms of the drill string vibration

is crucial to improve the control of the operation and performance. In this study,

a simple model was proposed with one degree of freedom to analyze the torsional

vibration of the string. Due to several sources of uncertainty in the process, the

Monte Carlo method was used in order to consider these uncertainties in the anal-

ysis. The dynamic response of the system was simulated and the sensitivity of the

variables weight on bit (WOB) and rotary table speed was studied by introducing

a proposal for the use of fork graphics. Uncertainties were included in the weight

on bit parameter and the probability of Stick-Slip for each rotary table speed was

analyzed, showing that higher speeds result in smaller probabilities of Stick-Slip.

Then functions were applied to drill-rock friction torques in order to simulate the

first drill-rock contact and the rock changing during the drilling process, observing

specific behaviors for each function. Finally, they were added to the uncertainty

analysis to assess the impact for each function, showing that the most influence

occurs in the case of the step function.

vii

Sumario

Lista de Figuras x

Lista de Tabelas xiii

1 Introducao 1

1.1 Motivacao . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1

1.2 Objetivo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4

1.3 Organizacao . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4

2 Modelagem matematica do processo de perfuracao 6

2.1 Processo de perfuracao . . . . . . . . . . . . . . . . . . . . . . . . . . 6

2.2 Modelo computacional do sistema . . . . . . . . . . . . . . . . . . . . 8

2.3 Consideracao sobre incertezas . . . . . . . . . . . . . . . . . . . . . . 11

3 Resultados e Discussoes 14

3.1 Resposta determinıstica do sistema: domınio do tempo e bifurcacao . 14

3.2 Bifurcacao Estocastica . . . . . . . . . . . . . . . . . . . . . . . . . . 23

3.3 Alteracao do contato broca-rocha no tempo . . . . . . . . . . . . . . 25

3.4 Resposta determinıstica sem contato inicial . . . . . . . . . . . . . . . 27

3.5 Resposta determinıstica com mudanca de caracterısticas da rocha . . 30

3.6 Resposta estocastica com mudanca de caracterısticas da rocha . . . . 33

4 Conclusao e Trabalhos futuros 38

4.1 Conclusao . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38

4.2 Trabalhos futuros . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40

Referencias Bibliograficas 41

viii

A Codigos em MATLAB 43

A.1 Programa principal . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43

A.2 Representacao dos modelos utilizados . . . . . . . . . . . . . . . . . . 68

A.2.1 Modelo simples - sem mudanca de rocha . . . . . . . . . . . . 68

A.2.2 Modelo com multiprocessamento . . . . . . . . . . . . . . . . 69

A.2.3 Modelo com mudanca de rocha tipo degrau . . . . . . . . . . . 70

A.2.4 Modelo com mudanca de rocha tipo linear . . . . . . . . . . . 71

A.2.5 Modelo com mudanca de rocha tipo tangente hiperbolica . . . 72

ix

Lista de Figuras

1.1 Evolucao do Petroleo. Adaptado de [1]. . . . . . . . . . . . . . . . . . 2

1.2 Tipos de viracao em uma coluna de perfuracao. Adaptado de [2]. . . 3

2.1 Esquema geral de uma sonda de perfuracao. Adaptado de [3]. . . . . 7

2.2 Modelo mecanico proposto para a coluna de perfuracao. . . . . . . . . 9

2.3 Esquema do modelo de atrito seco na broca. Adaptado de [4]. . . . . 10

2.4 Esquema do metodo de monte carlo. Baseado em [5]. . . . . . . . . . 12

3.1 Resposta determinıstica do sistema para os regimes com e sem Stick-

Slip em (a) velocidade angular da broca e (b) angulo da broca. Linha

cheia: Ω = 19 rad/s; Linha tracejada: Ω = 9,5 rad/s. . . . . . . . . . 15

3.2 Modelo de atrito seco. . . . . . . . . . . . . . . . . . . . . . . . . . . 16

3.3 Resposta do torque de atrito seco na broca (a) durante o tempo total

de simulacao e (b) detalhado entre 21 e 28 segundos. Linha cheia: Ω

= 19 rad/s; Linha tracejada: Ω = 9,5 rad/s. . . . . . . . . . . . . . . 17

3.4 Sensibilidade da resposta sem stick-slip do sistema para uma variacao

de 5% na velocidade de rotacao da mesa. Resultados em (a) veloci-

dade de rotacao da broca e (b) angulo da broca. . . . . . . . . . . . . 18

3.5 Sensibilidade da resposta com stick-slip do sistema para uma variacao

de 5% na velocidade de rotacao da mesa. Resultados em (a) veloci-

dade de rotacao da broca e (b) angulo da broca. . . . . . . . . . . . . 19

3.6 Analise por bifurcacao da resposta do sistema mediante variacao da

velocidade de rotacao da mesa. . . . . . . . . . . . . . . . . . . . . . 20

3.7 Sensibilidade da resposta sem stick-slip do sistema para uma variacao

de 5% no peso sobre a broca. Resultados em (a) velocidade de rotacao

da broca e (b) angulo da broca. . . . . . . . . . . . . . . . . . . . . . 21

x

3.8 Sensibilidade da resposta com stick-slip do sistema para uma variacao

de 5% no peso sobre a broca. Resultados em (a) velocidade de rotacao

da broca e (b) angulo da broca. . . . . . . . . . . . . . . . . . . . . . 22

3.9 Analise por bifurcacao da resposta do sistema mediante variacao do

peso sobre a broca. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22

3.10 Analise de estabilidade de acordo com velocidade de rotacao da mesa

e peso sobre a broca. Azul: sem Stick-Slip; Vermelho: com Stick-Slip. 23

3.11 Analise por bifurcacao estocastica para determinar a probabilidade

de ocorrer o fenomeno de Stick-Slip. . . . . . . . . . . . . . . . . . . . 25

3.12 Alteracao, durante o primeiro contato broca-rocha, do (a) torque de

atrito estatico e (b) torque de atrito dinamico. . . . . . . . . . . . . . 26

3.13 Alteracao, durante o processo de perfuracao, do (a) torque de atrito

estatico e (b) torque de atrito dinamico. . . . . . . . . . . . . . . . . 27

3.14 Velocidade da broca com rotacao inicial livre de atrito e aplicacoes

do tipo degrau, linear e tanh para o peso sobre a broca. . . . . . . . . 28

3.15 Detalhe da velocidade da broca com rotacao inicial livre de atrito e

aplicacoes do tipo degrau, linear e tanh para o peso sobre a broca. . . 29

3.16 Torque de atrito com rotacao inicial livre de atrito e aplicacoes do

tipo degrau, linear e tanh para o peso sobre a broca. . . . . . . . . . 29

3.17 Detalhe do torque de atrito com rotacao inicial livre de atrito e

aplicacoes do tipo degrau, linear e tanh para o peso sobre a broca. . . 30

3.18 Velocidade da broca com mudanca de rocha seguindo uma funcao

degrau, linear e tanh. . . . . . . . . . . . . . . . . . . . . . . . . . . . 31

3.19 Detalhe da velocidade da broca com mudanca de rocha seguindo uma

funcao degrau, linear e tanh. . . . . . . . . . . . . . . . . . . . . . . . 32

3.20 Torque de atrito com mudanca de rocha seguindo uma funcao degrau,

linear e tanh. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32

3.21 Detalhe do torque de atrito com mudanca de rocha seguindo uma

funcao degrau, linear e tanh. . . . . . . . . . . . . . . . . . . . . . . . 33

3.22 Envelope de confianca de 99% para a velocidade da broca com mu-

danca de rocha seguindo uma funcao linear. . . . . . . . . . . . . . . 35

xi


danca de rocha seguindo uma funcao tangente hiperbolica. . . . . . . 35

3.24 Envelope de confianca de 99% para o torque de atrito com mudanca

de rocha seguindo uma funcao (a) linear e (b) tangente hiperbolica. . 36


danca de rocha seguindo uma funcao degrau. . . . . . . . . . . . . . . 37

3.26 Envelope de confianca de 99% para o torque de atrito com mudanca

de rocha seguindo uma funcao degrau. . . . . . . . . . . . . . . . . . 37

xii

Lista de Tabelas

2.1 Propriedades da coluna de perfuracao . . . . . . . . . . . . . . . . . . 11

xiii

Capıtulo 1

Introducao

Inicialmente, o petroleo foi utilizado pelas civilizacoes antigas, como os meso-

potamicos, gregos, persas e egıpcios para aquecimento e iluminacao de casas e es-

tradas, embalsamador de corpos, impermeabilizador de construcoes, entre outros.

O petroleo, desde a sua primeira perfuracao, no dia 27 de agosto de 1859, no estado

americano da Pensilvania, por Edwin Laurentine Drake (conhecido como coronel

Drake), foi tendo sua importancia ampliada com o passar do tempo e o desenvol-

vimento das tecnologias. Na era pos-Drake, observou-se o seu grande potencial e,

assim, passou-se a utiliza-lo com maior amplitude e a pesquisar novas formas e

aplicabilidades. [6]

Hoje em dia, o Petroleo e considerado um dos recursos naturais mais importantes

para o desenvolvimento humano. Derivados do petroleo como a gasolina e o diesel

sao utilizados para mover importantes maquinas como caminhoes, carros, avioes,

navios, entre outros. Alem disso o petroleo ainda e utilizado como uma importante

fonte de energia, inclusive para a producao de energia eletrica. O Petroleo ainda e

utilizado na producao de bens consumo presentes no dia-a-dia da populacao mundial,

como plasticos, ceras, fertilizantes, entre muitos outros. [7]

1.1 Motivacao

Devido a alta e crescente importancia do petroleo para a sociedade, esse recurso

continua a ser explorado hoje em dia e as buscas por novas reservas continuam.

Conforme o passar dos anos, desde a primeira descoberta do petroleo, as reservas

1

descobertas foram se tornando cada vez mais difıceis de serem exploradas, neces-

sitando um investimento alto no desenvolvimento de novas tecnologias. Conforme

pode ser visto na figura 1.1, as profundidades dos novos pocos vem aumentando com

o passar dos anos e, a partir de 2007, quando da descoberta dos primeiros pocos

no pre-sal localizados na bacia de santos, essa profundidade aumentou consideravel-

mente, chegando a incrıveis 4 mil metros de rocha perfurada em um lamina d’agua

de mais de 2 mil metros.

Figura 1.1: Evolucao do Petroleo. Adaptado de [1].

O processo de perfuracao e certamente um dos processos mais importantes na

construcao de um poco de petroleo e, de acordo com KHULIEF et al. [8], a vibracao

da coluna de perfuracao e uma das principais causas para a reducao da performance

do processo. Experiencias de campo ainda revelaram que e crucial o entendimento

dos mecanismos complexos de vibracao apresentados pelo sistema de perduracao

com o intuito de controlar melhor a operacao e melhorar a performance do processo.

[8]

Na coluna de perfuracao, podem ser percebidos tres tipos de vibracao: axial,

lateral e torcional (figura 1.2). Uma das formas mais recorrentes e destrutivas e

a torcional que quando atinge seu estado de maior severidade e conhecida como

Stick-Slip.[8]

O termo Stick-Slip e utilizado pela industria para descrever o fenomeno que

acontece quando a broca apresenta ciclos de velocidade repetidos com valores de

mınimo proximo ou igual a zero e maximos elevados, enquanto que a mesa rotativa

2

continua a girar normalmente. Quando a broca esta travada, a energia torcional

aumente ate atingir um patamar que a broca nao consegue mais resistir, voltando

a girar de forma brusca com altos valores de velocidade angular ate que a energia

torcional reduza e a broca ”trave”novamente, gerando um ciclo. De acordo com

BRETT [9], seus experimentos mostraram que o fenomeno de Stick-Slip apenas

acontece quando a broca esta em contato com a rocha a ser perfurada, ou seja, o

fenomeno de Stick-Slip nao se manifestou antes do contato broca-rocha. Portanto, e

evidente que um dos principais mecanismos de Stick-Slip e a interacao broca-rocha.

KHULIEF et al. [8] tambem diz que as principais razoes para as vibracoes da coluna

de perfuracao sao o contato da broca com a formacao rochosa e o contato da coluna

de perfuracao com a parede do poco.

Figura 1.2: Tipos de viracao em uma coluna de perfuracao. Adaptado de [2].

Sabe-se que o petroleo e gerado a partir do acumulo de materia organica presente

em rochas sedimentares que, quando submetidas a certas condicoes de pressao e

temperatura, formam o petroleo que por sua vez e acumulado em uma rocha com

certos nıveis de porosidade e permeabilidade que e chamada de rocha reservatorio.

[6]

3

1.2 Objetivo

Os objetivos desse trabalho sao: i. propor um modelo simplificado para descrever

o movimento torcional da coluna de perfuracao, baseando-se no modelo proposto

por NAVARRO e LOPEZ [4]; ii. Analisar a resposta determinıstica da dinamica

do sistema modelado e propor um novo tipo de analise por bifurcacao; iii. Analisar

a resposta estocastica do sistema para variacoes nos torques de atrito atraves do

grafico de bifurcacao estocastica; iv. Avaliar a resposta determinıstica quando ha

alteracao no contato broca-rocha e; v. Analisar a resposta estocastica do sistema no

caso de alteracao do contato broca-rocha durante a perfuracao.

1.3 Organizacao

O primeiro capıtulo tem como objetivo introduzir o assunto ao leitor e apresentar a

motivacao do trabalho com uma contextualizacao do problema abordado, o objetivo

do trabalho e a organizacao deste documento.

O segundo capıtulo comeca com uma explicacao breve do processo de perfuracao

com foco na atividade de perfuracao em si, sem aprofundar nos processos laterais

da operacao. O entendimento do processo de perfuracao permite o desenvolvimento

de um modelos simplificado para descrever a coluna de perfuracao e suas condicoes

de contorno. Em seguida e apresentada a modelagem proposta para a mudanca de

caracterısticas da rocha e por ultimo a consideracao das incertezas e abordado.

O terceiro capıtulo comeca apresentando a resposta dinamica do sistema sim-

plificado desenvolvido no capıtulo anterior. Primeiramente a resposta dinamica do

sistema e analisada, seguida de uma analise de sensibilidade dos parametros do mo-

delo. A seguir, uma analise de estabilidade e apresentada com a quantificacao de

incertezas nos parametros de interacao broca-rocha. Apos isso, a resposta do sis-

tema e analisada quando a rotacao inicial nao possui contato com a rocha de forma

a avaliar o comportamento do sistema no momento do contato com a rocha. O

capıtulo entao segue com a analise da resposta do sistema mediante a alteracao da

rocha durante a perfuracao e, em seguida, apresenta a resposta estocastica para esse

caso.

O capıtulo cinco conclui o trabalho e contem algumas sugestoes para trabalhos

4

futuros.

Em anexo pode ser encontrado o programa, codificado em MATLAB, utilizado

para o desenvolvimento desse trabalho.

5

Capıtulo 2

Modelagem matematica do

processo de perfuracao

Com o intuito de possibilitar a realizacao de analises e estudos de um sistema fısico, e

importante que o modelo represente de forma simplificada os fenomenos pertinentes

a serem estudados. Outro ponto importante e a necessidade desse sistema ser o mais

simples possıvel, mas que ainda assim forneca um resultado satisfatorio. Apesar do

modelo usado nesse estudo ser altamente simplificado, ele engloba os fenomenos

mais importantes da vibracao de colunas de perfuracao.

2.1 Processo de perfuracao

Um dos principais sistemas fısicos em um processo de perfuracao e a coluna de per-

furacao que e uma estrutura tubular que liga a sonda de perfuracao ate a rocha. Essa

coluna ainda permite a circulacao da lama de perfuracao, transmissao da rotacao,

sensoriamento, entre outras funcoes. O esquema geral de uma sonda de perfuracao

pode ser visto na figura 2.1.

6

Figura 2.1: Esquema geral de uma sonda de perfuracao. Adaptado de [3].

A coluna de perfuracao consiste do Bottom Hole Assembly (BHA) e dos drillpipes

que sao conectados entre si formando um longo tubo que pode chegar a quilometros

de comprimento. O BHA e composto pela broca, estabilizadores que previnem o

desbalanceamento da coluna, e uma seria de tubos relativamente pesados conhecidos

como drill collars.

Outro elemento muito importante no processo de perfuracao e a lama de per-

furacao que possui, entre outras, as funcoes de limpeza, resfriamento e lubrificacao

da broca. A lama de perfuracao e injetada no topo da coluna de perfuracao, percor-

rendo por dentro dos drill pipes ate chegar na rocha sendo perfurada, onde cumpre

sua funcao e retorna a sonda de perfuracao pelo espaco anular entre os drill pipes e

a parede do poco.

Na extremidade inferior da coluna de perfuracao, encontra-se a broca responsavel

pela perfuracao da rocha. A rotacao e originada na sonda de perfuracao por um

motor eletrico e transmitida atraves dos drill pipes ate a broca. O mecanismo de

rotacao pode ser de dois tipos: mesa rotativa ou top drive.

Durante a perfuracao, o acompanhamento e a analise do sistema mecanico da

7

coluna de perfuracao se da pela interpretacao dos parametros de superfıcie impostos

e medidos. Em poucos cenarios de exploracao sao utilizados sensores proximos a

broca com envio de dados em tempo real para a superfıcie [3]. Normalmente o

operador controla tres parametros: rotacao da mesa, peso sobre a broca (WOB,

abreviacao de Weight on Bit) e vazao do fluido de perfuracao. O peso sobre a broca

e um dos principais fatores na interacao broca-rocha.

2.2 Modelo computacional do sistema

O modelo matematico exposto a seguir foi baseado em NAVARRO e SUAREZ [4]. A

coluna de perfuracao foi modelada como um pendulo torcional simples cuja rotacao

e imposta por um motor eletrico, cuja velocidade angular foi modelada como cons-

tante, e a interacao broca-rocha foi modelada atraves de um modelo de atrito seco.

Os drill pipes sao representados atraves de uma mola torcional linear de rigidez

k e atraves de um amortecimento c que conecta duas inercias, Jc e Jb, representando

a inercia da mesa rotativa e a inercia dos drill pipes junto com o BHA. Apenas a

inercia Jb e representada no esquema da figura 2.2 pois a inercia Jc nao e utilizada

na representacao matematica do modelo proposto devido a hipotese de rotacao cons-

tante da mesa rotativa. A inercia Jb e comumente calculada como sendo a inercia

do BHA mais 1/3 da inercia dos drill pipes [9].

Para a modelagem da interacao broca-rocha considerou-se um modelo de atrito

seco que causa um torque reativo Tfb aliado a um amortecimento cb que representa

o atrito viscoso devido a lama de perfuracao, principalmente.

8

Figura 2.2: Modelo mecanico proposto para a coluna de perfuracao.

Com o intuito de simplificacao do modelo, algumas hipotese foram consideradas:

(i) O poco e a coluna de perfuracao estao perfeitamente na vertical; (ii) Nao existe

nenhum movimento lateral; (iii) O atrito nas conexoes e entre a coluna e a parede do

poco foram desconsiderados; (iv) A lama de perfuracao e simplificada por um atrito

viscoso na broca; (v) O escoamento da lama de perfuracao e considerado laminar.

Finalmente, levando em conta todas as consideracoes feitas ate esse ponto, e

possıvel escrever a equacao de movimento que descreve o modelo mecanico proposto:

Jbθb + c(θb − Ω) + k(θb − Ωt) = −cbθb − Tfb(θb) (2.1)

A equacao acima tambem pode ser representada atraves do espaco de estados,

especialmente util na resolucao de equacoes diferenciais nao-lineares de ordem su-

perior. Considerando y1 = θb e y2 = θb, a representacao no espaco de estados fica

da seguinte maneira:

y =

y1y2

=

y2

− k

Jby1 −

c

Jby2 +

k

JbΩt+

c

JbΩ− Tfb

Jb(y1)−

cbJby2

(2.2)

O modelo utilizado para representar a interacao broca-rocha foi sugerido por NA-

VARRO e SUAREZ [4] e e uma variacao do efeito Stribeck junto com o modelo de

atrito seco. O parametro Tfb e portanto dependente da velocidade da broca.

9

Figura 2.3: Esquema do modelo de atrito seco na broca. Adaptado de [4].

Na figura 2.3 esta representado o modelo utilizado para o torque resistivo da

interacao broca-rocha. Atraves dessa figura e possıvel identificar o torque Tsb que

caracteriza o atrito estatico e o torque Tcb que caracteriza o atrito dinamico. Os

torques de atrito estatico e dinamico, Tsb e Tcb, dependem do raio da broca (Rb),

do peso sobre a broca (Wob) e dos coeficientes de atrito estatico, no caso do Tsb, e

dinamico, no caso do Tcb. O decaimento presente entre os atritos estatico e dinamico

e modelado por uma funcao exponencial com uma constante positiva γb.

Tfb =

Teb , θb < Dv e |Teb| < Tsb

Tsbsign(Teb) , θb < Dv e |Teb| > Tsb

Tcb + (Tsb − Tcb)e−γbθbsign(θb) , θb > Dv

(2.3)

Teb = c(Ω− θb) + k(Ωt− θb)− cbθb (2.4)

Tsb = RbWobµsb (2.5)

Tcb = RbWobµcb (2.6)

O primeiro caso da equacao 2.3 e chamado de Stick e representa a situacao onde

a broca esta presa solidaria a rocha, ou seja, quando a velocidade angular da broca

10

e nula. O terceiro caso e conhecido como Slip e representa a situacao onde a broca

esta perfurando a rocha, ou seja, quando a velocidade angular da broca e diferente

de zero. O segundo caso descreve a transicao entre os regimes Stick e Slip.

A equacao 2.4 representa o torque reativo no regime Stick que e igual ao torque

aplicado pois o sistema esta em equilıbrio. Como a mesa rotativa continua a girar,

o sistema acumula energia e o torque aplicado aumenta ate que o torque reativo em

equilıbrio atinge um valor maximo definido pelo torque estatico e a broca comeca a

girar, entrando no regime do atrito dinamico (Slip).

O parametro Dv especifica uma vizinhanca ao redor de θb = 0 suficientemente

pequena para compensar erros numericos. Os valores dos parametros utilizados ate

aqui foram retirados de NAVARRO e SUAREZ [4] e estao descritos na tabela 2.1.

Tabela 2.1: Propriedades da coluna de perfuracao

Variavel Valor

Jb 0.0318 kg.m2

c 0.0001 N.m.s/rad

cb 0.03 N.m.s/rad

k 0.073 N.m/rad

Tcb 5 N.m

Tsb 8 N.m

Dv 10−6s

γb 0.9

Ω 19 rad/s

2.3 Consideracao sobre incertezas

A industria do petroleo enfrenta desafios mais complexos a cada dia e e necessario

que as analises realizadas se tornem casa vez mais precisas. Com isso metodos

estocasticos se tornam cada vez mais comuns nas solucoes de problemas.

Os metodos estocasticos permitem a consideracao de incertezas nos modelos

matematicos desenvolvidos, quantificando o impacto destas no comportamento do

sistema e assim, aumentando a confiabilidade nas analises dos resultados. Essas

11

incertezas podem ser provenientes da falta de informacoes ou a ocorrencias nao

previstas no projeto de um poco. Os parametros modelados como variaveis aleatorias

nesse trabalho foram escolhidos devido ao alto grau de dificuldade em identifica-los

de maneira precisa. Esses metodos sao comumente utilizados quando e difıcil ou

impossıvel obter uma formulacao fechada para a resposta do sistema como e o caso

da analise dinamica nao linear de colunas de perfuracao.

As analises estocasticas presentes neste trabalho foram realizadas utilizando o

Metodo de Monte Carlo (MC) que consiste na simulacao de variaveis aleatorias

independentes para os parametros de entrada do sistema afim de realizar analises

probabilısticas de uma resposta desconhecida. As variaveis aleatorias simuladas

devem seguir a distribuicao de probabilidade das variaveis de entrada que represen-

tam. Neste trabalho foram consideradas as seguintes distribuicoes de probabilidade

contınuas: uniforme e normal.

Figura 2.4: Esquema do metodo de monte carlo. Baseado em [5].

A distribuicao contınua uniforme e definida na equacao 2.7, onde a e b repre-

sentam o limite inferior e superior, respectivamente. Essa distribuicao foi utilizada

para descrever as variaveis aleatorias Tcb e Tsb na analise de bifurcacao estocastica.

f(x|a, b) =

1b−a , a ≤ x ≤ b

0 , caso contrario

(2.7)

A distribuicao contınua normal e apresentada na equacao 2.8, onde µ e σ representam

a media e o desvio padrao respectivamente. Essa distribuicao foi utilizada para

descrever os parametros envolvidos na alteracao dos torques de atrito estatico e

12

dinamico na analise estocastica da resposta do sistema.

f(x|µ, σ) =1√

2πσ2e(−

(x−µ2)2σ2

),−∞ < x <∞, σ > 0 (2.8)

No caso da utilizacao da distribuicao contınua normal, a dispersao das variaveis

aleatorias e comumente avaliada atraves de um parametro chamado coeficiente de

variacao (CV), definido na equacao 2.9. Esse parametro permite uma analise mais

compreensıvel do que o desvio padrao ja que e expressado em porcentagem.

CV =

√E([X − E(X)]2)

E(X)=σ

µ(2.9)

13

Capıtulo 3

Resultados e Discussoes

Esse capıtulo possui o objetivo de discutir a simulacao do modelo para diversas

situacoes. Com isso, e possıvel identificar as principais caracterısticas do modelo

utilizado, as diferencas entre os regimes com e sem Stick-slip, alem de fazer uma

analise de sensibilidade das varaveis de interacao broca-rocha e da velocidade de

rotacao da mesa rotativa. Alem disso, e apresentada uma analise estocastica in-

cluindo incertezas nos parametros de interacao broca-rocha com o intuito de avaliar

a probabilidade do regime com Stick-slip se manifestar. Por fim sao comentados os

resultados obtidos quando ha alteracao do contato broca-rocha no tempo.

3.1 Resposta determinıstica do sistema: domınio

do tempo e bifurcacao

A resposta dinamica do sistema, na ausencia de Stick-slip, deve apresentar um regime

inicial com velocidade nula devido ao atrito estatico presente na interacao broca-

rocha, caracterizado pelo fenomeno de Stick. Apos a energia torcional acumulada

na coluna vencer o atrito estatico, o movimento se inicia dando lugar ao atrito

dinamico apresentando o fenomeno de Slip. Com o inıcio do movimento, um regime

transitorio se manifesta devido a elasticidade e amortecimento do sistema seguido

de um regime permanente cujo modulo deve ser igual a rotacao da mesa rotativa ja

que esta e responsavel por impor a rotacao da coluna de perfuracao.

No caso com Stick-slip, o regime inicial com velocidade nula tambem esta pre-

sente porem o inicio do movimento se da de forma brusca devido ao alto acumulo de

14

energia torcional anterior ao movimento, resultando em altas velocidades de rotacao

na broca. Essas altas rotacoes fazem com que o angulo da broca se aproxime rapida-

mente do angulo da mesa rotativa, reduzindo consideravelmente a energia torcional

da coluna, causando assim uma nova parada da broca. Esse ciclo entao se repete

ate que haja uma intervencao.

A figura 3.1 exemplifica os dois casos explicitados anteriormente. Os dados

utilizados foram apresentados na tabela 2.1, com excecao da velocidade de rotacao

da mesa no caso com Stick-slip que foi reduzida em 50% com o objetivo de induzir

o Stick-slip.

E possıvel notar que mesmo com uma reducao de 50% na velocidade de rotacao

da mesa, a broca ainda atinge velocidades bem proximas do caso sem reducao devido

ao Stick-slip. Essa reducao tambem causou um delay no comeco da rotacao pois o

acumulo de energia na coluna demora mais tempo para atingir o nıvel crıtico e vencer

o atrito estatico.

0 10 20 30 400

10

20

30

40

50

60

70

tempo (s)

θb(rad/s)

sem Stick−Slipcom Stick−Slip

(a)

0 10 20 30 400

200

400

600

800

1000

tempo (s)

θb(rad)


(b)

Figura 3.1: Resposta determinıstica do sistema para os regimes com e sem Stick-Slip

em (a) velocidade angular da broca e (b) angulo da broca. Linha cheia: Ω = 19

rad/s; Linha tracejada: Ω = 9,5 rad/s.

A figura 3.2 representa o modelo de atrito com os parametros descritos na tabela

2.1. Conforme descrito na secao 2.2, quando a velocidade de rotacao da broca e

proxima de zero, o torque de atrito e igual ao torque reativo da coluna, ate um valor

maximo definido por Tsb. Nesse caso o torque de atrito apresenta um comportamento

linear no tempo pois a velocidade de rotacao da mesa e constante. Para velocidades

nao-nulas na broca, o torque de atrito tende exponencialmente para Tcb.

15

−10 −5 0 5 10−10

−5

0

5

10

θb (rad/s)

Tfb(N

.m)

← Tsb

Tcb

Figura 3.2: Modelo de atrito seco.

O comportamento supracitado pode ser claramente observado na figura 3.3a no

caso sem Stick-slip, pois o torque de atrito sobe linearmente com o tempo ate o valor

maximo de 8N.m e decresce exponencialmente para o valor de 5N.m. O decrescimo

exponencial acontece muito rapidamente devido ao aumento rapido da velocidade.

No caso com Stick-slip, percebe o mesmo comportamento linear seguido de um

decrescimo exponencial. Conforme dito anteriormente, no regime Stick-slip existe

um ciclo de repeticoes de travamento e destravamento da broca. Esse ciclo de

repeticoes pode ser visto na figura 3.3a.

A figura 3.3b apresenta um ciclo de repeticao em detalhe. Percebe-se que ha um

aumento exponencial do torque de atrito devido a parada da broca, que representa

o caminho inverso no grafico da figura 3.2 e e caracterıstico do modelo escolhido. A

obtencao desse valor maximo na parada da broca acontece devido a parte Slip do

modelo de atrito, mais especificamente no inıcio da parte exponencial, e nao possui

relacao com a parte Stick do modelo. Isso se justifica pois a diferenca angular da

mesa rotativa para a broca nao e capaz de produzir um torque reativo Teb proximo

ao torque maximo Tsb logo no inıcio da parada da broca, que tambem explica o salto

do torque de atrito para um valor bem abaixo.

Apos a broca parar completamente, o ciclo continua com a subida linear do

torque de atrito ate o valor maximo Tsb, quando a broca reinicia o movimento e o

torque de atrito decresce exponencialmente para Tcb, reiniciando o ciclo.

16

0 10 20 30 400

2

4

6

8

10

tempo (s)

Tfb(N

.m)


(a)

21 22 23 24 25 26 27 280

2

4

6

8

10

tempo (s)

Tfb(N

.m)

com Stick−Slip

(b)

Figura 3.3: Resposta do torque de atrito seco na broca (a) durante o tempo total

de simulacao e (b) detalhado entre 21 e 28 segundos. Linha cheia: Ω = 19 rad/s;

Linha tracejada: Ω = 9,5 rad/s.

Na secao anterior foi apresentada a resposta dinamica do sistema com e sem

Stick-Slip. No caso com Stick-Slip, foi necessario reduzir a velocidade de rotacao da

mesa em 50% para produzir tal efeito. Isso ilustra a mudanca de comportamento do

sistema quando ha variacoes nas variaveis de entrada, o que demonstra a importancia

de estudar a sensibilidade do sistema para as variaveis que possuem maior incerteza

ou possibilidade de modificacao durante o processo de perfuracao.

Os graficos presentes nas figuras 3.4 e 3.5 apresentam a resposta do sistema para

uma variacao de 5% na velocidade de rotacao da mesa rotativa nos casos sem e com

Stick-Slip, respectivamente. Ambos os casos foram analisados de forma que nao

houvesse mudanca de regime dentro da variacao de 5%. A variacao de 5% no caso

com Stick-Slip foi obtida em cima da reducao de 50% realizada para induzir esse

regime, ou seja, uma variacao de -5% nesse regime representa 47.5% da velocidade

de rotacao da tabela 2.1.

O grafico da figura 3.4 mostra que a variacao da velocidade de rotacao da mesa

induz uma alteracao igual na resposta permanente da velocidade de rotacao da

broca. Isso acontece pois a velocidade da broca tende para a velocidade de rotacao

da mesa rotativa no regime sem Stick-Slip. Esse mesmo comportamento pode ser

observado no grafico do angulo da broca que muda a inclinacao da reta quando a

velocidade se torna constante.

A resposta transitoria apresenta um pequeno deslocamento no tempo devido a

17

alteracao do tempo necessario para a coluna acumular energia suficiente para vencer

o atrito estatico. O pico de velocidade apresentou uma variacao de 1,4% no modulo

em relacao a velocidade de referencia da mesa.

0 20 40 60 800

10

20

30

40

50

60

tempo (s)

θb(rad/s)

ΩΩ + 5%Ω − 5%

(a)

0 20 40 60 800

500

1000

1500

tempo (s)

θb(rad)

ΩΩ + 5%Ω − 5%

(b)

Figura 3.4: Sensibilidade da resposta sem stick-slip do sistema para uma variacao

de 5% na velocidade de rotacao da mesa. Resultados em (a) velocidade de rotacao

da broca e (b) angulo da broca.

A figura 3.5, referente ao regime com Stick-Slip, mostra uma diferenca temporal

da resposta que vai aumentando com o tempo. Isso pode ser explicado pois ve-

locidades inferiores da mesa rotativa demandam um ciclo de Stick-Slip com maior

duracao. Essa diferenca temporal vai acumulando com o tempo e ficando cada vez

maior. O pico apresentou uma variacao de 0,7% no modulo em relacao a rotacao de

referencia utilizada nesse regime.

18

0 20 40 60 800

10

20

30

40

50

60

tempo (s)

θb(rad/s)

ΩΩ + 5%Ω − 5%

(a)

0 20 40 60 800

100

200

300

400

500

600

tempo (s)

θb(rad)

ΩΩ + 5%Ω − 5%

(b)

Figura 3.5: Sensibilidade da resposta com stick-slip do sistema para uma variacao

de 5% na velocidade de rotacao da mesa. Resultados em (a) velocidade de rotacao

da broca e (b) angulo da broca.

O grafico de bifurcacao apresentado na figura 3.6 tem, por sua vez, o objetivo de

avaliar a mudanca de regime mediante a variacao da velocidade de rotacao da mesa

assim como avaliar a sensibilidade do sistema a esse parametro. O eixo das abscissas

apresenta a razao da velocidade real da mesa em relacao a velocidade de referencia,

presente na tabela 2.1. O eixo das ordenadas representa os valores mınimo e maximo

para a velocidade da broca no regime permanente. Essa analise foi inspirada nos

graficos contido no trabalho de CAYRES [10].

O grafico no inıcio apresenta duas curvas que se encontram quando a velocidade

da mesa atinge cerca de 90% da velocidade de referencia. A primeira parte do

grafico, onde ha duas curvas, representa o regime com Stick-Slip pois a velocidade

mınima da broca e zero, a broca esta travada (Stick). A segunda parte, onde ha

apenas uma curva, representa o regime sem Stick-Slip pois o maximo e mınimo sao

equivalentes, representando um regime permanente estavel.

A curva de maximos sobe com o aumento da velocidade da mesa rotativa. No

regime com Stick-Slip, a broca atinge velocidades bem superiores a da mesa rotativa.

Por outro lado, no caso sem Stick-Slip, a broca atinge a mesma velocidade da mesa

rotativa e conforme a velocidade da mesa aumenta, a velocidade da broca aumentar

igualmente.

19

0.5 0.75 1 1.25 1.5

0

0.5

1

1.5

2

2.5

3

3.5

Ω

Ωref

min,m

ax(

θb

Ωref)(rad/s)

Figura 3.6: Analise por bifurcacao da resposta do sistema mediante variacao da

velocidade de rotacao da mesa.

Outro parametro muito importante no processo de perfuracao que e controlado

pelo operador, alem de possuir diversas fontes de incerteza, e o peso sobre a broca

(Wob). O peso sobre a broca esta presente no modelo atraves dos torques de atrito

Tsb e Tcb nas equacoes 2.5 e 2.6, respectivamente. Devido a dependencia linear de

tais variaveis, uma variacao no peso sobre a broca pode ser modelada diretamente

nos torques de atrito.

O grafico da figura 3.7 mostra a sensibilidade do sistema no regime sem Stick-Slip.

E possıvel notar que, ao contrario da velocidade da mesa, a variacao no peso sobre

a broca nao provoca uma variacao da velocidade da broca em regime permanente.

Isso pode ser confirmado no grafico 3.7b pois nao ha variacao da inclinacao da reta

no regime permanente.

Existe ainda uma leve variacao de 3,5% no modulo do primeiro pico do grafico

3.7a. Alem disso a resposta e levemente deslocada no tempo pois a variacao no

torque de atrito faz com que seja demandada uma energia acumulada maior para

vencer o atrito broca-rocha.

20

0 20 40 60 800

10

20

30

40

50

60

tempo (s)

θb(rad/s)

WobWob + 5%Wob − 5%

(a)

0 20 40 60 800

500

1000

1500

tempo (s)

θb(rad)


(b)

Figura 3.7: Sensibilidade da resposta sem stick-slip do sistema para uma variacao

de 5% no peso sobre a broca. Resultados em (a) velocidade de rotacao da broca e

(b) angulo da broca.

No caso com Stick-Slip (figura 3.8), os graficos foram gerados considerando tor-

ques de atrito iguais ao dobro dos da referencia. Isso foi feito com o intuito de

induzir o fenomeno Stick-Slip.

Percebe-se um aumento significativo (72%) no modulo do pico de velocidade

em relacao a resposta sem Stick-Slip. Isso mostra a severidade desse fenomeno.

Alem disso, existe uma variacao temporal entre as respostas. Essa variacao possui

a mesma causa explicada no caso sem Stick-Slip, sendo mais evidente aqui devido

aos ciclos de repeticao do Stick-Slip.

21

0 20 40 60 800

20

40

60

80

100

tempo (s)

θb(rad/s)


(a)

0 20 40 60 800

500

1000

1500

tempo (s)

θb(rad)


(b)

Figura 3.8: Sensibilidade da resposta com stick-slip do sistema para uma variacao

de 5% no peso sobre a broca. Resultados em (a) velocidade de rotacao da broca e

(b) angulo da broca.

A figura 3.9 mostra o grafico de bifurcacao para variacoes do peso sobre a broca.

A primeira parte do grafico, referente ao regime sem Stick-Slip, mostra que a veloci-

dade em regime permanente permanece constante nao obstante a variacao do peso

sobre a broca. Isso confirma o fato constatado na figura 3.7a. A segunda parte do

grafico, referente ao regime com Stick-Slip, mostra o agravamento do fenomeno de

Stick-Slip quanto maior e o peso sobre a broca.

0.5 1 1.5 2 2.5

0

1

2

3

4

5

6

7

WobWobref

min,m

ax(

θb

Ωref)(rad/s)

Figura 3.9: Analise por bifurcacao da resposta do sistema mediante variacao do peso

sobre a broca.

O grafico da figura 3.10 mostra quando ocorre o fenomeno de Stick-Slip para

22

diferentes conbinacoes de rotacao da mesa e peso sobre a broca. A parte vermelha

do grafico representa situacoes com Stick-Slip e a parte azul representa situacoes

sem Stick-Slip. Nota-se que para valores baixos de velocidade da mesa e valores

altos de peso sobre a broca, o fenomeno de Stick-Slip se manifesta.

Esse grafico e obtido atraves do calculo de um parametro chamado ”Stick-Slip

Severity”que por sua vez e calculado da seguinte forma:

SSS =(θbmax − θbmin)

2 · Ω(3.1)

onde θbmax e θbmax sao as velocidades maxima e mınima da broca no regime perma-

nente, respectivamente. Caso o valor de SSS supere 0.8, foi considerado que ocorreu

Stick-Slip; caso contrario, nao houve.

0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.50.5

0.6

0.7

0.8

0.9

1

1.1

1.2

1.3

1.4

1.5

Ω

Ωref

Wob

Wob r

ef

Figura 3.10: Analise de estabilidade de acordo com velocidade de rotacao da mesa

e peso sobre a broca. Azul: sem Stick-Slip; Vermelho: com Stick-Slip.

3.2 Bifurcacao Estocastica

Na secao 2.3 ficou clara a importancia de modelar as incertezas das variaveis im-

portantes ao problema. Na secao 3.1 foi feito um estudo de sensibilidade dessas

variaveis. Esse estudo foi realizado de forma determinıstica para variacoes de 5%

nos parametros.

O estudo apresentado a seguir foi realizado considerando o peso sobre a broca

como variavel aleatoria, ou seja, essa variavel passa a ter uma funcao de probabili-

dade. A funcao de probabilidade escolhida para modelar tal variavel foi a Uniforme

23

com limites 10% acima e abaixo do valor de referencia mostrado na tabela 2.1. Todas

as demais variaveis foram considerada determinısticas. Esse estudo possui o objetivo

de avaliar a probabilidade de haver instabilidade, ou seja, ocorrer Stick-Slip, para

diferentes velocidades da mesa rotativa.

Tsb ∝ Tsbref · Unif(0.9, 1.1) (3.2)

Tcb ∝ Tcbref · Unif(0.9, 1.1) (3.3)

Utilizou-se o metodo de Monte Carlo para simular a resposta do sistema atraves

de 1000 simulacoes para velocidades de rotacao da mesa com passo de 0,5% da

velocidade de referencia. Em cada simulacao verificou o SSS a fim de determinar

o numero de simulacoes cujo SSS ficou maior que 0,8, caracterizando o numero de

”falhas”. A probabilidade de instabilidade para essa determinada velocidade de

rotacao da mesa foi entao calculada da seguinte maneira:

Pinstabilidade =nº de simulacoes com SSS > 0, 8

nº total de simulacoes [103](3.4)

O grafico da figura 3.11 mostra a probabilidade de instabilidade em funcao da

velocidade da mesa rotativa. Como esperado, a probabilidade de instabilidade reduz

com o aumento da velocidade da mesa. Esse comportamento e compatıvel com a

analise determinıstica que evidenciou o fato do Stick-Slip nao acontecer para velo-

cidade altas da mesa rotativa. A probabilidade de instabilidade atinge os valores

extremos de 0% e 100% devido a escolha da funcao uniforme de probabilidade que,

por sua vez, possui valores limites e agregam essa propriedade ao sistema. O fato

do valor de referencia de 19 rad/s para a velocidade de rotacao da mesa n]ao repre-

sentar 50% de probabilidade e devido ao fato da referencia nao estar no ponto de

instabilidade.

24

12 14 16 18 20 22 0%

20%

40%

60%

80%

100%

Ω (rad/s)

Pro

babi

lidad

e de

Inst

abili

dade

Figura 3.11: Analise por bifurcacao estocastica para determinar a probabilidade de

ocorrer o fenomeno de Stick-Slip.

3.3 Alteracao do contato broca-rocha no tempo

Levando em consideracao a modelagem escolhida e as simulacoes mostradas anteri-

ormente, fica claro que os parametros principais na interacao broca-rocha sao Tsb e

Tcb que representam os atritos estatico e dinamico, respectivamente.

Anteriormente a analise da alteracao de caracterısticas da rocha ao longo do

processo de perfuracao, achou-se oportuna a realizacao de simulacoes do primeiro

contato da broca com a rocha, apos o inıcio de sua rotacao com as superfıcies

separadas. Isso foi realizado devido a ser o procedimento padrao na operacao de

perfuracao. Esse primeiro contato depende principalmente da forma de aplicacao

do peso sobre a broca. Essa aplicacao foi simulada atraves dos parametros Tsb e Tcb,

utilizando tres funcoes diferentes: degrau, linear e tangente hiperbolica.

Antes do inıcio do contato, e importante que o sistema estivesse em regime

permanente para permitir uma analise mais criteriosa. Com esse objetivo, utilizou

um tempo de 20 segundos de simulacao para o inıcio do contato. O contato entao e

aplicado seguindo cada funcao escolhida ate um tempo de 60 segundos. Os torques de

atrito possuem valor nulo (sem contato) ate o tempo inicial de 20 segundos e possuem

valor igual a referencia acima de 60 segundo de simulacao. Esse comportamento esta

representado na figura 3.12.

25

0 20 40 60 800

1

2

3

4

5

6

7

8

tempo (s)

Tsb(N

.m)

DegrauLinearTanh

(a)

0 20 40 60 800

1

2

3

4

5

tempo (s)

Tcb(N

.m)

DegrauLinearTanh

(b)

Figura 3.12: Alteracao, durante o primeiro contato broca-rocha, do (a) torque de

atrito estatico e (b) torque de atrito dinamico.

Apos essa analise inicial, sera avaliado o impacto da alteracao de caracterısticas

da rocha durante o processo de perfuracao. Para isso, foram consideradas as mesmas

tres funcoes para descrever esse impacto atraves dos parametros Tsb e Tcb. Esses

parametro foram escolhidos pois a alteracao da rocha sendo perfurada implica em

um novo coeficiente de atrito que possui influencia de diversas propriedades da rocha

como geometria, dureza, rugosidade, entre outros.

Foi utilizado um tempo de 30 segundos de simulacao para o inıcio da alteracao

dos parametros com o intuito de aguardar o regime permanente se consolidar. Essa

alteracao acontece seguindo cada funcao especificada ate o tempo de 70 segundos.

Foi escolhido um aumento de 2N.m para cada torque de atrito, mantendo o valor

final igual a referencia. A figura 3.13 mostra as funcoes utilizadas.

26

0 20 40 60 80

5

5.5

6

6.5

7

7.5

8

8.5

tempo (s)

Tsb(N

.m)

DegrauLinearTanh

(a)

0 20 40 60 802.5

3

3.5

4

4.5

5

5.5

tempo (s)

Tcb(N

.m)

DegrauLinearTanh

(b)

Figura 3.13: Alteracao, durante o processo de perfuracao, do (a) torque de atrito

estatico e (b) torque de atrito dinamico.

3.4 Resposta determinıstica sem contato inicial

A simulacao da resposta do sistema com rotacao inicial livre e mostrada na figura

3.14. O tempo entre as linhas azuis e referente ao contato progressivo da broca com

a rocha.

A resposta do sistema no caso da funcao linear apresentou uma nova velocidade

de regime permanente durante a aplicacao do contato. Esse novo regime perma-

nente foi atingido apos uma pequena oscilacao. Esse comportamento e devido ao

modelo escolhido para o sistema ja que, no caso linear, a solucao particular para o

angulo da broca e linear (considerando que nao existe Stick-Slip), sendo constante

na velocidade da broca e, quando somada a solucao homogenea, provoca uma nova

velocidade de regime permanente. E importante ressaltar que a mesa rotativa per-

manece girando na mesma velocidade enquanto que a broca passa a girar em uma

velocidade menor, ou seja, a energia torcional comecar a ser acumulada progres-

sivamente devido a torcao da coluna. Essa energia torcional sendo acumulada e

suficiente para vencer o atrito e manter a broca girando e perfurando.

Quando o contato e aplicado seguindo uma funcao tangente hiperbolica, a ve-

locidade nao atinge um patamar novo como no caso linear. A funcao hiperbolica

possui a caracterıstica de ser bem suave e apresentar um aumento progressivo da

sua derivada seguido de uma diminuicao tambem progressiva. Essa caracterıstica

27

influencia diretamente a resposta do sistema. E possıvel perceber isso na figura

3.14 pois a velocidade da broca reduz gradativamente para compensar o atrito que

aumenta devagar no inıcio ate chegar um valor mınimo que representa o ponto de

inflexao da funcao tangente hiperbolica e maxima e volta a subir gradativamente

quando o atrito se aproxima do valor final.

0 20 40 60 800

10

20

30

40

50

60

tempo (s)

θb(rad/s)

DegrauLinearTanh

Figura 3.14: Velocidade da broca com rotacao inicial livre de atrito e aplicacoes do

tipo degrau, linear e tanh para o peso sobre a broca.

No caso da funcao degrau, o contato e feito instantaneamente. Nesse caso,

percebe-se um travamento inicial da broca pois, antes do contato, a energia torcional

acumulada na coluna e muito pequena ja que havia uma resistencia bem menor.

Quando o contato e realizado de forma brusca, a energia torcional acumulada nao e

suficiente para vencer o atrito e continuar girando, causando o travamento da broca.

Em outras palavras, a desaceleracao da broca acontece mais rapidamente do que a

coluna acumula energia torcional para vencer esse atrito.

Apesar da funcao de aplicacao do atrito ser degrau, a velocidade da broca nao cai

de forma degrau. Segundo pode ser observado na figura 3.15, a broca e desacelerada

pelo atrito entre os tempos de 20 e 20,13 segundos. Apesar disso, o tempo de

caimento e bem pequeno e pode ser considerado instantaneo.

28

19.9 20 20.1 20.2 20.30

10

20

30

40

50

60

tempo (s)

θb(rad/s)

DegrauLinearTanh

Figura 3.15: Detalhe da velocidade da broca com rotacao inicial livre de atrito e

aplicacoes do tipo degrau, linear e tanh para o peso sobre a broca.

O grafico da figura 3.16 mostra o torque de atrito ao longo do tempo. Para

as funcoes linear e tangente hiperbolica, percebe-se que o grafico possui a mesa

forma da funcao e isso acontece pois nao ha alteracao de regime Stick-Slip durante

o contato. Portanto, nos casos linear e tangente hiperbolica, o grafico representa o

torque de atrito dinamico (Tcb).

0 20 40 60 800

1

2

3

4

5

6

7

8

tempo (s)

Tfb(N

.m)

DegrauLinearTanh

Figura 3.16: Torque de atrito com rotacao inicial livre de atrito e aplicacoes do tipo

degrau, linear e tanh para o peso sobre a broca.

29

No caso da funcao degrau, a figura 3.17 mostra que o torque de atrito sobe

instantaneamente para o torque de atrito dinamico final. Porem, conforme mostrado

anteriormente, a velocidade da broca cai ate o valor nulo e o torque de atrito portanto

torna-se estatico e entao decai conforme o modelo escolhido. O torque de atrito nao

cai para zero pois ja existe uma energia torcional na coluna, apesar de pequena. A

figura 3.16 mostra entao que a coluna comeca a acumular energia ate vencer o atrito

estatico e sair da inercia, dando lugar para o atrito dinamico.

19.9 20 20.1 20.2 20.30

1

2

3

4

5

6

7

8

tempo (s)

Tfb(N

.m)

DegrauLinearTanh

Figura 3.17: Detalhe do torque de atrito com rotacao inicial livre de atrito e

aplicacoes do tipo degrau, linear e tanh para o peso sobre a broca.

3.5 Resposta determinıstica com mudanca de ca-

racterısticas da rocha

Esse item abordara a mudanca de caracterısticas da rocha durante o processo de

perfuracao, ou seja, a alteracao da interacao broca-rocha quando a broca ja esta em

contato com a rocha. A metodologia de analise dos graficos e semelhante ao item

anterior.

O grafico da figura 3.18 mostra que o efeito da alteracao do atrito nas funcoes

linear e tangente hiperbolica e significantemente menor que o efeito mostrado no

item anterior. Isso e explicado pois a alteracao dos torques de atrito possuem menor

amplitude e mesmo tempo de duracao, permitindo uma melhor compensacao da

30

energia torcional.

Contudo, o efeito dessa alteracao possui a mesma forma e justificativa que foram

explicadas no item anterior. Levando em conta que nao ha manifestacao de Stick-

Slip, a minimizacao do efeito pode ser atribuıda ao fato da alteracao dos parametros

de atrito possuir menor amplitude.

0 20 40 60 800

10

20

30

40

50

60

tempo (s)

θb(rad/s)

DegrauLinearTanh

Figura 3.18: Velocidade da broca com mudanca de rocha seguindo uma funcao

degrau, linear e tanh.

No caso da alteracao seguindo uma funcao degrau, o comportamento do sistema e

bem similar ao observado no estudo com rotacao inicial livre. Isso pode ser explicado

pois a energia torcional necessaria para vencer o atrito estatico nos dois caso e a

mesma, pois o atrito estatico final nos dois casos tambem e igual. Portanto, a

liberacao dessa energia acumulada gera rotacoes de mesma amplitude.

A broca desacelera com o torque de atrito e possui um comportamento bem

similar ao explicado no item anterior, conforme pode ser constatado na figura 3.19.

Nesse caso, o travamento da broca leva mais tempo pois a alteracao do torque de

atrito e menor, permitindo uma maior compensacao do atrito por meio da energia

torcional acumulada.

31

29.9 30 30.1 30.2 30.3 30.4 30.50

10

20

30

40

50

60

tempo (s)

θb(rad/s)

DegrauLinearTanh

Figura 3.19: Detalhe da velocidade da broca com mudanca de rocha seguindo uma

funcao degrau, linear e tanh.

O grafico da figura 3.20 mostra o torque de atrito para todas as funcoes. Nova-

mente, para os casos linear e tangente hiperbolica, o grafico demonstra comporta-

mento identico a mudanca modelada ja que o regime de atrito permanece dinamico

todo o tempo.

0 20 40 60 800

1

2

3

4

5

6

7

8

tempo (s)

Tfb(N

.m)

DegrauLinearTanh

Figura 3.20: Torque de atrito com mudanca de rocha seguindo uma funcao degrau,

linear e tanh.

No caso da funcao degrau, o comportamento e bastante semelhante ao do estudo

32

com rotacao inicial livre. A figura 3.21 mostra que o torque aumenta instantanea-

mente para o torque de atrito dinamico final e permanece assim ate o travamento

da broca, onde o atrito cai rapidamente. Apos isso, a energia torcional e acumulada

ate superar o atrito estatico e voltar para o atrito dinamico final.

29.9 30 30.1 30.2 30.3 30.4 30.50

1

2

3

4

5

6

7

8

tempo (s)

Tfb(N

.m)

DegrauLinearTanh

Figura 3.21: Detalhe do torque de atrito com mudanca de rocha seguindo uma

funcao degrau, linear e tanh.

3.6 Resposta estocastica com mudanca de carac-

terısticas da rocha

Os estudos feitos sobre mudanca de rocha consideraram parametros determinısticos

para descrever essa mudanca. O estudo presente nessa secao visa avaliar a con-

sequencia da incerteza desses parametros na resposta do sistema.

Com o intuito de simular a resposta do sistema levando em conta variaveis

aleatorias, foi utilizado o metodo de Monte Carlo descrito no item 2.3. As variaveis

aleatorias foram considerada normais com coeficiente de variacao igual a 2% e mo-

33

deladas da seguinte forma:

Tcbi ∝ 3 ·Normal(1; 0, 02) (3.5)

Tcbf ∝ 5 ·Normal(1; 0, 02) (3.6)

Tsbi ∝ 6 ·Normal(1; 0, 02) (3.7)

Tsbf ∝ 8 ·Normal(1; 0, 02) (3.8)

ti ∝ 30 ·Normal(1; 0, 02) (3.9)

tf ∝ 70 ·Normal(1; 0, 02) (3.10)

onde Tcbi e Tcbf sao os torques de atrito dinamico inicial e final, respectivamente;

Tsbi e Tsbf sao os torques de atrito estatico inicial e final, respectivamente, e; ti e tf

sao os tempos inicial e final da mudanca de rocha, respectivamente.

Os graficos apresentados a seguir sao conhecidos por envelopes de confianca.

A linha cheia no grafico representa a media dos valores obtidos para determinada

reposta em cada instante de tempo. Por outro lado, as linhas tracejadas delimitam

uma regiao onde 90% dos valores obtidos estao contidos.

O valor de 90% de confianca foi escolhido pois nao apresenta casos com Stick-Slip

em seu domınio. Como os parametros sao modelados como funcoes normais e essas

por sua vez possuem um domınio de (−∞,∞), podem surgir valores que induzem

o Stick-Slip. Observar a ocorrencia desse fenomeno foge ao objetivo desse item, ja

sendo abordado em itens anteriores e, devido a isso, considerou-se os 99% para o

envelope de confianca.

As figuras 3.22 e 3.23 mostram o envelope de confianca para a velocidade da

broca quando a mudanca de rocha segue uma funcao linear e tangente hiperbolica,

respectivamente. Percebe-se nesses graficos que as incertezas impostas nao oferecem

alteracoes expressivas na resposta do sistema. A explicacao desse comportamento

remete a analise feita na figura 3.9 que mostrou que alteracoes nos fatores de atrito

apenas impactam a resposta da velocidade da broca se houver Stick-Slip. Como nos

casos linear e tangente hiperbolica nao ha Stick-Slip, os modulos das velocidades

permanecem muito proximos. Existe apenas um deslocamento no tempo devido as

incertezas nos tempos de inıcio e fim da mudanca de rocha, mas pode ser consideravel

desprezıvel e praticamente imperceptıvel.

34

0 20 40 60 800

10

20

30

40

50

60

tempo (s)

θb(rad/s)

Figura 3.22: Envelope de confianca de 99% para a velocidade da broca com mudanca

de rocha seguindo uma funcao linear.

0 20 40 60 800

10

20

30

40

50

60

tempo (s)

θb(rad/s)


de rocha seguindo uma funcao tangente hiperbolica.

No caso das funcoes linear e tangente hiperbolica, a figura 3.24 mostra os graficos

de envelope de confianca para o torque de atrito. As incertezas embutidas podem

ser claramente identificadas. No caso dessas funcoes, apos vencer o atrito inicial, a

35

broca permanece no regime de atrito dinamico e a incerteza nesse parametro, antes

e depois da mudanca de rocha, pode ser percebida nos graficos.

A incerteza nos tempos inicial e final da mudanca de rocha tambem e percebida

atraves de um leve deslocamento no tempo, quando comparadas as 3 linhas de

cada grafico. As incertezas modeladas nao sao grandes o suficiente para perceber a

dispersao dos valores de torque estatico desses graficos.

0 20 40 60 800

2

4

6

8

tempo (s)

Tfb(N

.m)

(a)

0 20 40 60 800

2

4

6

8

tempo (s)

Tfb(N

.m)

(b)

Figura 3.24: Envelope de confianca de 99% para o torque de atrito com mudanca

de rocha seguindo uma funcao (a) linear e (b) tangente hiperbolica.

No caso da funcao degrau, percebe-se uma dispersao muito maior no resultado

pois o envelope de confianca fica significativamente mais espesso. Isso acontece

devido a manifestacao do Stick-Slip no inıcio da mudanca de rocha. Retornando

novamente a figura 3.9, pode-se atribuir uma maior dispersao no caso degrau pois a

amplitude da resposta e afetada quando existe Stick-Slip. Alem disso tambem esta

presente um leve deslocamento no tempo devido as incertezas nos tempos de inıcio

e fim da mudanca de rocha.

36

0 20 40 60 800

10

20

30

40

50

60

tempo (s)

θb(rad/s)


de rocha seguindo uma funcao degrau.

A figura 3.26 mostra o grafico de envelope de confianca para o torque no caso

da funcao degrau. O envelope de confianca tambem demonstra maior dispersao

no inıcio da alteracao da rocha. Essa dispersao e atribuıda ao Stick-Slip pois esse

fenomeno engloba todos os regimes de atrito e acaba sofrendo maior influencia da

incerteza desses parametros.

0 20 40 60 800

2

4

6

8

tempo (s)

Tfb(N

.m)

Figura 3.26: Envelope de confianca de 99% para o torque de atrito com mudanca

de rocha seguindo uma funcao degrau.

37

Capıtulo 4

Conclusao e Trabalhos futuros

4.1 Conclusao

O objetivo do presente trabalho foi de introduzir novas analises sobre o problema da

dinamica de colunas de perfuracao. Primeiramente foi proposto um modelo simples

para descrever a coluna de perfuracao de forma a permite as analises posteriores.

Esse modelo simplificado foi baseado em NAVARRO e SUAREZ [4] com a diferenca

de imposicao de velocidade constante na mesa rotativa.

Apos a definicao do modelo, foram realizadas algumas simulacoes convencionais

de forma a analisar as respostas classicas para esse tipo de estrutura e assim validou-

se o modelo utilizado.

A seguir foi realizado um estudo que avaliou a sensibilidade dos parametros do

modelo na resposta do sistema. Esse estudo introduziu uma nova forma de analise

grafica dessa sensibilidade atraves do que foi intitulado graficos de bifurcacao.

Essa analise mostrou que a velocidade da mesa rotativa influencia a resposta

do sistema independente da presenca de Stick-Slip. Por outro lado, o peso sobre a

broca, descrito atraves dos torques de atrito, apenas influencia a resposta do sistema

no caso com Stick-Slip. Para ambos os parametros (velocidade da mesa rotativa e

peso sobre a broca), a analise mostrou que a broca atinge velocidades bem superiores

na presenca de Stick-Slip.

Para finalizar essa analise, foi introduzida uma nova analise do sistema conside-

rando incertezas nos torques de atrito. Essa analise visou avaliar a probabilidade

de ocorrencia de Stick-Slip para diversas velocidades da mesa rotativa. Percebeu-se

38

que velocidades baixas da mesa rotativa geram uma probabilidade alta de Stick-

Slip, enquanto que altas velocidade da mesa resultam em baixas probabilidades de

Stick-Slip.

Por fim, conforme objetivo do presente trabalho, foi analisado o impacto da

mudanca de rocha no processo de perfuracao. Esse impacto foi avaliado atraves dos

torques de atrito que sao responsaveis pela modelagem da interacao broca-rocha.

A mudanca dos torques de atrito no tempo foi modelado atraves de tres funcoes:

linear, tangente hiperbolica e degrau.

Primeiramente foi analisado o caso onde a broca entra em contato com a rocha

somente apos atingir o regime permanente em rotacao. Esse caso foi avaliado pois

representa as operacoes reais de perfuracao de pocos. A seguir, foi avaliada o caso

onde ha alteracao da rocha durante o processo de perfuracao.

Em ambos os casos, os resultados mostraram que a variacao da velocidade da

broca foi mais suave com a funcao tangente hiperbolica. A funcao linear, por sua

vez, apresentou uma suavidade menor apesar de sua variacao possuir menor ampli-

tude. Por fim, a funcao degrau impos uma resposta tao brusca quanto a excitacao.

Inclusive, a funcao degrau induziu um princıpio de Stick-Slip que nao se manteve

ate o final.

No caso da alteracao da rocha durante o processo de perfuracao, a resposta

do sistema foi mais suave para as funcoes linear e tangente hiperbolica. Isso foi

explicado pois a variacao dos parametros de atrito foram menores. O mesmo nao

aconteceu para a funcao degrau que demonstrou um impacto equivalente em ambos

os casos. Esse acontecimento foi explicado pelo fato dos parametros de atrito finais

serem os mesmos e pela presenca de Stick-Slip.

Como analise final, foram avaliadas as incertezas dos parametros que definem

a alteracao da rocha durante o processo de perfuracao. A principal conclusao a

ser tirada dessa ultima analise e que essas incertezas nao afetaram a velocidade da

broca nas mudancas linear e tangente hiperbolica pois essas nao apresentaram Stick-

Slip. No caso degrau existe Stick-Slip e as incertezas afetam a resposta de maneira

significativa.

39

4.2 Trabalhos futuros

O presente trabalho possuiu o objetivo de introduzir novas metodologias de analise

para esse problema comum na industria do petroleo. Os proximos passos para

tal trabalho seria refinar o modelo e calibra-lo para que seja capaz de descrever

satisfatoriamente o caso real.

Apos isso, seria interessante refazer as analises aqui desenvolvidas e avaliar qual

funcao melhor descreve a mudanca dos parametros de atrito na alteracao da rocha

e como esses parametros de atrito sao alterados dependendo do tipo de rochas e

condicoes que o poco esta submetido. Dessa forma espera-se ser capaz de avaliar os

riscos de forma mais precisa a fim de mitigar gastos excessivos e otimizar o processo

de perfuracao de pocos.

Outra sugestao seria incluir oscilacoes dos parametros de atrito durante o tempo

de simulacao. Dentro de uma mesma rocha, os parametro de atrito podem variar e

seria interessante avaliar a resposta do sistema a essas oscilacoes.

Por fim, sugere-se a utilizacao das analises aqui desenvolvidas para otimizar a

performance e o controle desse tipo de operacao.

40

Referencias Bibliograficas

[1] “Petrobras 60 anos”, http://infograficos.oglobo.globo.com/economia/

petrobras-60-anos-1.html, Online. Acessado em: 2016-06-13.

[2] FRANCA, L., Perfuracao Percurssiva-Rotativa Auto-Excitada em Rochas Du-

ras , Ph.D. Thesis, Pontifıcia Universidade Catolica, Rio de Janeiro, Brasil,

2004.

[3] PERCY, J. G., ANALISE DE INCERTEZAS EM VIBRACOES LATERAIS E

DE TORCAO ACOPLADAS EM COLUNAS DE PERFURACAO , Mas-

ter’s Thesis, UFRJ/COPPE, Rio de Janeiro, Brasil, 2014.

[4] NAVARRO-LOPEZ, E. M., SUAREZ, R., “Practical approach to modelling and

controlling stick-slip oscillations in oilwell drillstrings”. In: Control Ap-

plications, 2004. Proceedings of the 2004 IEEE International Conference

on, v. 2, pp. 1454–1460 Vol.2, Sept 2004.

[5] SEUNG-KYUM CHOI, RAMANA V. GRANDHI, R. A. C., Reliability-based

Structural Design. Springer-Verlag London.

[6] “Geologia - Texto produzido pelos alunos da Primeira Turma do Curso de

Engenharia de Petroleo da UFC com base nos livros ”O Universo

da Industria Petrolıfera- Da Pesquisa a Refinacao”, de Jose Salgado

Gomes. Fundacao Calouste-Gulbenkian, Portugal, 2ª Edicao, 2011 e

”Fundamentos de Engenharia de Petroleo”, de Jose Eduardo Thomas.

Editora Interciencia, 2001.” http://www.petroleo.ufc.br/index.php?

option=com_content&task=view&id=394&Itemid=56, Online. Acessado

em: 2016-06-13.

41

http://infograficos.oglobo.globo.com/economia/petrobras-60-anos-1.html

http://infograficos.oglobo.globo.com/economia/petrobras-60-anos-1.html

http://www.petroleo.ufc.br/index.php?option=com_content&task=view&id=394&Itemid=56

http://www.petroleo.ufc.br/index.php?option=com_content&task=view&id=394&Itemid=56

[7] “Petroleum”, http://www.scienceclarified.com/Oi-Ph/Petroleum.html,

Online. Acessado em: 2016-06-13.

[8] KHULIEF, Y., AL-SULAIMAN, F., BASHMAL, S., “Vibration analysis of

drillstrings with self-excited stick–slip oscillations”, Journal of Sound and

Vibration, v. 299, n. 3, pp. 540 – 558, 2007.

[9] BRETT, J. F., “The Genesis of Bit-Induced Torsional Drillstring Vibrations”,

Society of Petroleum Engineers .

[10] CAYRES, B. C., Numerical and Experimental Analysis of the Nonlinear Torsi-

onal Dynamics of a Drilling System, Master’s Thesis, Pontifıcia Univer-

sidade Catolica do Rio de Janeiro, Rio de Janeiro, Brasil, 2013.

[11] L.C. CUNHA LIMA, R.R. AGUIAR, T. R., HBAIEB, S., “Analysis of the

torsional stability of a simplified drillstring”. In: DINAME 2015 - Pro-

ceedings of the XVII International Symposium on Dynamic Problems of

Mechanics , Feb 2015.

[12] RITTO, T., SOIZE, C., SAMPAIO, R., “Non-linear dynamics of a drill-string

with uncertain model of the bit–rock interaction”, International Journal

of Non-Linear Mechanics , v. 44, n. 8, pp. 865 – 876, 2009.

42

http://www.scienceclarified.com/Oi-Ph/Petroleum.html

Apendice A

Codigos em MATLAB

A.1 Programa principal

1 %% Declarando variaveis globais

2 global k Jb c Vrot Gb Dv cb VrotR TsbR TcbR Tsb Tcb Tsb1 Tcb1 Tsb2 ...

Tcb2 ti tf

3

4 % Variaveis caracteristicas do sistema

5 Jb=0.0318; % Inercia da coluna

6 c=0.0001; % constante de amortecimento da coluna

7 cb=0.03; % constate do atrito viscoso da broca

8 k=0.073; % rigidez da coluna

9 Dv=10ˆ−6; % Correcao de intervalo do modelo de atrito seco ...

(problemas de convergencia)

10 Gb=0.9; % Constante gama do modelo de atrito seco

11

12 % Variaveis controlaveis

13 VrotR = 19; % Velocidade de rotacao da mesa

14 TsbR = 8; % Pico de torque de referencia

15 TcbR = 5; % Torque de atrito referencia

16

17 % Vetor temporal

18 ta = 0:0.001:120;

19

20 %% Analise Dinamica Deterministica do Sistema

43

21

22 % Propriedades de interacao com a rocha iguais a referencia

23 Tsb = TsbR;

24 Tcb = TcbR;

25

26 % Resposta Estavel

27 % Velocidade de rotacao da mesa igual a referencia

28 Vrot = VrotR;

29 % Solucionando o sistema

30 [ts,ys] = ode23tb(@modelo simples,ta,[0 0]);

31

32 % Resposta Instavel

33 % Velocidade de rotacao da mesa igual a metade da referencia

34 Vrot = 0.5*VrotR;


36 [te,ye] = ode23tb(@modelo simples,ta,[0 0]);

37

38 % Grafico Angulo da broca

39 figure

40 plot(ts,ys(:,1),'k',te,ye(:,1),'k−−','LineWidth',2)

41 legend('sem Stick−Slip','com Stick−Slip','Location','northwest')

42 xlabel ('tempo (s)','Interpreter','latex','FontSize',26)

43 ylabel ('$\theta b$(rad)','Interpreter','latex','FontSize',26)

44 ylim([0 2500])

45 xlim([0 40])

46 set(gca,'fontsize',20)

47 print −depsc Resposta angulo

48

49 % Grafico Velocidade da broca

50 figure

51 plot(ts,ys(:,2),'k',te,ye(:,2),'k−−','LineWidth',2)



54 ylabel ('$\dot\theta b$(rad/s)','Interpreter','latex','FontSize',26)

55 ylim([0 70])

56 xlim([0 40])


58 print −depsc Resposta velocidade

59

44

60 % Calculo do torque de atrito na broca para a resposta estavel

61 Vrot = VrotR;

62 Tfbs = zeros(length(ys),1);

63 for i = 1:length(ys)

64 Teb = c*(Vrot−ys(i,1))+k*(Vrot*ts(i)−ys(i,1))−cb*ys(i,2); % ...

Toque Stick

65

66 if and(ys(i,2)<Dv,abs(Teb)<Tsb) % Stick

67 Tfbs(i)=Teb;

68 elseif and(ys(i,2)<Dv,abs(Teb)>Tsb) % Stick−Slip ...

Transition

69 Tfbs(i)=Tsb*sign(Teb);

70 else % Slip

71 Tfbs(i)=Tcb+(Tsb−Tcb)*exp(−Gb*abs(ys(i,2)))*sign(ys(i,2));

72 end

73 end

74

75 % Calculo do torque de atrito na broca para a resposta instavel

76 Vrot = VrotR*0.5;

77 Tfbe = zeros(length(ye),1);

78 for i = 1:length(ye)

79 Teb = c*(Vrot−ye(i,1))+k*(Vrot*te(i)−ye(i,1))−cb*ye(i,2); % ...

Toque Stick

80

81 if and(ye(i,2)<Dv,abs(Teb)<Tsb) % Stick

82 Tfbe(i)=Teb;

83 elseif and(ye(i,2)<Dv,abs(Teb)>Tsb) % Stick−Slip ...

Transition

84 Tfbe(i)=Tsb*sign(Teb);

85 else % Slip

86 Tfbe(i)=Tcb+(Tsb−Tcb)*exp(−Gb*abs(ye(i,2)))*sign(ye(i,2));

87 end

88 end

89

90 % Grafico Torque de atrito na broca

91 figure

92 plot(ts,Tfbs,'k',te,Tfbe,'k−−','LineWidth',2)



45

95 ylabel ('$T fb$(N.m)','Interpreter','latex','FontSize',26)

96 ylim([0 10])

97 xlim([0 40])


99 print −depsc Resposta torque

100

101 % Grafico detalhado do Torque de atrito na broca

102 figure

103 plot(te,Tfbe,'k−−','LineWidth',2)

104 legend('com Stick−Slip','Location','northwest')



107 ylim([0 10])

108 xlim([21 28])


110 print −depsc Resposta torque detalhado

111

112 %% Grafico do modelo de atrito seco

113 y2 = −10:0.01:10;

114 Tfb=TcbR+(TsbR−TcbR)*exp(−Gb*abs(y2));

115 Tfb=Tfb.*sign(y2);

116 figure

117 plot(y2,Tfb,'k','LineWidth',2)

118 xlabel ('$\dot\theta b$ (rad/s)','Interpreter','latex','FontSize',26)

119 ylabel ('$T fb$ (N.m)','Interpreter','latex','FontSize',26)

120 ylim([−10 10])

121 text(0,8,'$\leftarrow T sb$','Interpreter','latex','FontSize',26)

122 text(8,4,'$T cb$','Interpreter','latex','FontSize',26)


124 print −depsc Modelo atrito

125

126 %% Analise de sensibilidade para velocidade de rotacao da mesa

127

128 % Torques de atrito iguais a referencia

129 Tcb = TcbR; % Torque de atrito dinamico

130 Tsb = TsbR; % Pico de torque de atrito

131

132 % Resposta Estavel

133 Vrot = VrotR; % Velocidade de rotacao da mesa

46

134 [t1,y1] = ode23tb(@modelo simples,ta,[0 0]);

135

136 Vrot = VrotR*1.05; % Velocidade de rotacao da mesa + 5%


138

139 Vrot = VrotR*0.95; % Velocidade de rotacao da mesa − 5%


141

142 % Grafico angulo da broca

143 figure

144 plot(t1,y1(:,1),'k−',t2,y2(:,1),'k−−',t3,y3(:,1),'k:','LineWidth',2)

145 legend('\Omega','\Omega + 5%','\Omega − 5%','Location','northwest')

146 xlabel('tempo (s)','Interpreter','latex','FontSize',26)

147 ylabel ('$\theta b$ (rad)','Interpreter','latex','FontSize',26)

148 ylim([0 1500])


150 print −depsc Sensibilidade multiVrot est angulo

151

152 % Grafico velocidade da broca

153 figure


155 legend('\Omega','\Omega + 5%','\Omega − 5%','Location','northeast')


157 ylabel ('$\dot\theta b$ (rad/s)','Interpreter','latex','FontSize',26)

158 ylim([0 60])


160 print −depsc Sensibilidade multiVrot est velocidade

161

162 % Resposta Instavel

163 Vrot = VrotR*0.5; % Velocidade de rotacao da mesa


165

166 Vrot = VrotR*0.5*1.05; % Velocidade de rotacao da mesa + 5%


168

169 Vrot = VrotR*0.5*0.95; % Velocidade de rotacao da mesa − 5%


171


47

173 figure





178 ylim([0 600])


180 print −depsc Sensibilidade multiVrot inst angulo

181


183 figure





188 ylim([0 60])


190 print −depsc Sensibilidade multiVrot inst velocidade

191

192 %% Analise de sensibilidade para peso sobre a broca

193

194 % Velocidade da mesa igual a referencia

195 Vrot = VrotR;

196

197 % Resposta estavel

198 Tcb = TcbR; % Torque de atrito dinamico

199 Tsb = TsbR; % Pico de torque de atrito


201

202 Tcb = TcbR*1.05; % Torque de atrito dinamico + 5%

203 Tsb = TsbR*1.05; % Pico de torque de atrito − 5%


205

206 Tcb = TcbR*0.95; % Torque de atrito dinamico + 5%

207 Tsb = TsbR*0.95; % Pico de torque de atrito − 5%


209


211 figure

48


213 legend('Wob','Wob + 5%','Wob − 5%','Location','northwest')



216 ylim([0 1500])


218 print −depsc Sensibilidade multiWob est angulo

219


221 figure


223 legend('Wob','Wob + 5%','Wob − 5%','Location','northeast')



226 ylim([0 60])


228 print −depsc Sensibilidade multiWob est velocidade

229

230 % Resposta instavel

231 Tcb = TcbR*2; % Torque de atrito dinamico

232 Tsb = TsbR*2; % Pico de torque de atrito


234

235 Tcb = TcbR*2*1.05; % Torque de atrito dinamico

236 Tsb = TsbR*2*1.05; % Pico de torque de atrito


238

239 Tcb = TcbR*2*0.95; % Torque de atrito dinamico

240 Tsb = TsbR*2*0.95; % Pico de torque de atrito


242


244 figure





249 ylim([0 1500])


49

251 print −depsc Sensibilidade multiWob inst angulo

252


254 figure





259 ylim([0 110])


261 print −depsc Sensibilidade multiWob inst velocidade

262

263 %% Analise de sensibilidade por bifurcacao para peso sobre a broca

264

265 %Velocidade da mesa igual a referencia

266 Vrot = VrotR;

267

268 % Criacao de variaveis de suporte

269 x = 0.5:0.005:2.5;

270 V topo = zeros(length(x),1);

271 V topo2 = zeros(length(x),1);

272 i=1;

273

274 % Calculo dos maximos e minimos

275 for f = x

276 Tcb = TcbR*f; % Torque de atrito dinamico

277 Tsb = TsbR*f; % Pico de torque de atrito


279 V topo(i) = max(ys(60000:end,2))/19;

280 V topo2(i) = min(ys(60000:end,2))/19;

281 i = i + 1;

282 end

283

284 % Geracao do grafico

285 figure

286 plot(x,V topo,'k',x,V topo2,'k','LineWidth',2)

287 xlabel('$\fracWobWob ref$','Interpreter','latex','FontSize',26)

288 ylabel('$min,max(\frac\dot\theta b\Omega ref)$ ...

(rad/s)','Interpreter','latex','FontSize',26)

50

289 ylim([−0.50 7])

290 xlim([min(x) max(x)])

291 set(gca,'fontsize',20,'XTick',[0.5:0.5:2.5],'YTick',[0:1:7])

292 print −depsc Sensibilidade maxWob

293

294

295 %% Analise de sensibilidade por bifurcacao para rotacao da mesa

296

297 % Torques iguais a referencia

298 Tcb = TcbR;

299 Tsb = TsbR;

300

301 % Criacao de varieveis de suporte

302 x = 0.5:0.005:1.5;

303 V topo = zeros(length(x),1);

304 V topo2 = zeros(length(x),1);

305 i=1;

306

307 % Calculo dos maximos e minimos

308 for f = x

309 Vrot=VrotR*f; % Velocidade de rotacao da mesa


311 V topo(i) = max(ys(60000:end,2))/19;

312 V topo2(i) = min(ys(60000:end,2))/19;

313 i = i + 1;

314 end

315


317 figure

318 plot(x,V topo,'k',x,V topo2,'k','LineWidth',2)

319 xlabel('$\frac\Omega\Omega ref$','Interpreter','latex','FontSize',26)

320 ylabel('$min,max(\frac\dot\theta b\Omega ref)$ ...

(rad/s)','Interpreter','latex','FontSize',26)

321 ylim([−0.50 3.5])

322 xlim([min(x) max(x)])

323 set(gca,'fontsize',20,'XTick',[0.5:0.25:1.5],'YTick',[0:0.5:3.5])

324 print −depsc Sensibilidade maxVrot

325

326

51

327 %% Analise de Instabilidade por mapeamento do SSS

328


330 SSS = zeros(441,1);

331 x chart = zeros(441,1);

332 y chart = zeros(441,1);

333 filled chart = zeros(21,2);

334 i=0;

335 j=1;

336

337 % Calculo do SSS

338 for p2 = 0.5:0.05:1.5

339 filled chart(j,1) = 1.5;

340 flag = 0;

341 for p1 = 0.5:0.05:1.5

342 if flag == 1

343 continue

344 end

345 i = i + 1;

346 Tcb = TcbR*p1;

347 Tsb = TsbR*p1;

348 Vrot = VrotR*p2;


350 SSS(i) = (max(ys(60000:end,2))−min(ys(60000:end,2)))/(2*Vrot);

351 if SSS(i) > 0.8

352 flag=1;

353 filled chart(j,1) = p1;

354 filled chart(j,2) = 1.5 − filled chart(j,1);

355 end

356 end

357 j = j + 1;

358 end

359


361 figure

362 ax=0.5:0.05:1.5;

363 area(ax,filled chart,0.5)

364 xlabel ...

('$\frac\Omega\Omega ref$','Interpreter','latex','FontSize',26)

52

365 ylabel ('$\fracWobWob ref$','Interpreter','latex','FontSize',26)

366 ylim([0.5 1.5])

367 set(gca,'fontsize',20,'XTick',0.5:0.1:1.5,'YTick',0.5:0.1:1.5)

368 print −depsc SSS

369

370 %% Bifurcacao Estocastica

371


373 Var Tcb = 0.1;

374 Var Tsb = 0.1;

375 n sim = 1000;

376 Vrot array = 0.6:0.005:1.2;

377

378 % Igualar variaveis as referencias

379 VrotRl = VrotR;

380 TcbRl = TcbR;

381 TsbRl = TsbR;

382

383 % Simulacoes utilizando processamento paralelo

384 n = zeros(length(Vrot array),1);

385 matlabpool(4)

386 parfor j=1:length(Vrot array)

387 Tcb sim = (1−Var Tcb)*TcbRl + 2*Var Tcb*TcbRl.*rand(n sim,1);

388 Tsb sim = (1−Var Tsb)*TsbRl + 2*Var Tsb*TsbRl.*rand(n sim,1);

389 Vrotl = VrotRl*Vrot array(j);

390 SSS = zeros(n sim,1);

391 for i=1:n sim

392 Tcbl = Tcb sim(i);

393 Tsbl = Tsb sim(i);

394 disp([i Vrotl Tsbl Tcbl])

395 [ts,ys] = ode23tb(@(t,y) modelo multiproc(t,y,[Vrotl Tsbl ...

Tcbl]),ta,[0 0]);

396 SSS(i) = (max(ys(60000:end,2))−min(ys(60000:end,2)))/(2*Vrotl);

397 end

398 n(j) = sum(SSS>0.8)/n sim;

399 end

400 matlabpool close

401

402 % Criacao do grafico

53

403 figure

404 plot(Vrot array.*VrotR,n,'k','LineWidth',2)

405 xlabel('$\Omega$ (rad/s)','Interpreter','latex','FontSize',26)

406 ylabel('Probabilidade de Instabilidade','FontSize',26)

407 xlim([min(Vrot array.*VrotR) max(Vrot array.*VrotR)])


409

410 % Converte o eixo y em percentagem

411 a=[cellstr(num2str(get(gca,'ytick')'*100))];

412 pct = char(ones(size(a,1),1)*'%');

413 new yticks = [char(a),pct];

414 set(gca,'yticklabel',new yticks)

415

416 print −depsc bifurcacao estocastica

417

418 %% Analise Dinamica do Sistema com mudanca de parametros de atrito seco

419

420 Analise = 'alt Tcb Tsb'; % Tipo de analise a ser realizada (output ...

file)

421 Tcb1 = 3; % Torque de atrito dinamico para t < t1

422 TcbR = 5; % Torque de atrito dinamico referencia

423 Tsb1 = 6; % Pico de torque de atrito para t < t1

424 TsbR = 8; % Pico de torque de atrito referencia

425 ti = 30; % Tempo inicial para transicao das ...

propriedades de atrito

426 tf = 70; % Tempo final para transicao das ...


427 Vrot = VrotR; % Velocidade de rotacao da mesa igual a ...

referencia

428 Tcb2 = TcbR; % Torque de atrito dinamico

429 Tsb2 = TsbR; % Torque de atrito maximo

430 ta = 0:0.001:90; % Vetor temporal

431


433 [t degrau,y degrau] = ode23tb(@modelo degrau,ta,[0 0]);

434 [t linear,y linear] = ode23tb(@modelo linear,ta,[0 0]);

435 [t tanh,y tanh] = ode23tb(@modelo tanh,ta,[0 0]);

436

437 % Calculo do torque de atrito

54

438 Tfb = zeros(length(ta),3);

439 y = [y degrau y linear y tanh];

440

441 for j = 1:3

442

443 for i = 1:length(ta)

444 Teb = ...

c*(Vrot−y(i,2*j−1))+k*(Vrot*ta(i)−y(i,2*j−1))−cb*y(i,2*j); ...

% Toque Stick

445

446 if j == 3

447 Tcb = Tcb1 + ...

(Tcb2−Tcb1)*(tanh((((ta(i)−ti)/(tf−ti))*5)−2.5)+1)/2;

448 Tsb = Tsb1 + ...

(Tsb2−Tsb1)*(tanh((((ta(i)−ti)/(tf−ti))*5)−2.5)+1)/2;

449 elseif (ti>ta(i))

450 Tsb = Tsb1;

451 Tcb = Tcb1;

452 elseif (ti<ta(i)) && (ta(i)≤tf)

453 if j == 1

454 Tsb = Tsb2;

455 Tcb = Tcb2;

456 elseif j == 2

457 Tcb = Tcb1 + ((ta(i)−ti)/(tf−ti))*(Tcb2−Tcb1);

458 Tsb = Tsb1 + ((ta(i)−ti)/(tf−ti))*(Tsb2−Tsb1);

459 end

460 elseif ta(i) > tf

461 Tsb = Tsb2;

462 Tcb = Tcb2;

463 end

464

465

466 if and(y(i,2*j)<Dv,abs(Teb)<Tsb2) % Stick

467 Tfb(i,j)=Teb;

468 elseif and(y(i,2*j)<Dv,abs(Teb)>Tsb2) % ...

Stick−Slip Transition

469 Tfb(i,j)=Tsb*sign(Teb);

470 else % Slip

471 Tfb(i,j)=Tcb+(Tsb−Tcb)*exp(−Gb*abs(y(i,2*j)))*sign(y(i,2*j));

55

472 end

473 end

474 end

475

476 % Grafico Angulo da broca

477 figure

478 plot(t degrau,y degrau(:,1),'k',t linear,y linear(:,1),'k−−',t tanh,y tanh(:,1),'k:','LineWidth',2)

479 legend('Degrau','Linear','Tanh','Location','northwest')



482 ylim([0 ...

1.1*max([max(y degrau(:,1)),max(y linear(:,1)),max(y tanh(:,1))])])

483 xlim([0 max(ta)])

484 y1=get(gca,'ylim');

485 hold on

486 plot([ti ti],y1)

487 plot([tf tf],y1)

488 hold off


490 print(strcat('Resposta angulo ',Analise),'−depsc')

491

492 % Grafico Velocidade da broca

493 figure





498 ylim([0 ...




501 hold on



504 hold off


506 print(strcat('Resposta velocidade ',Analise),'−depsc')

507

508 % Grafico Torque de atrito na broca

56

509 figure

510 plot(ta,Tfb(:,1),'k',ta,Tfb(:,2),'k−−',ta,Tfb(:,3),'k:','LineWidth',2)




514 ylim([0 1.1*max(max(Tfb))])



517 hold on



520 hold off


522 print(strcat('Resposta torque ',Analise),'−depsc')

523

524 % Grafico Velocidade da broca (detalhe)

525 figure





530 ylim([0 ...


531 xlim([29.9 30.5])


533 hold on



536 hold off


538 print(strcat('Resposta velocidade detalhe ',Analise),'−depsc')

539

540 % Grafico Torque de atrito na broca (detalhe)

541 figure

542 plot(ta,Tfb(:,1),'k',ta,Tfb(:,2),'k−−',ta,Tfb(:,3),'k:','LineWidth',2)




546 ylim([0 1.1*max(max(Tfb))])

57

547 xlim([29.9 30.5])


549 hold on



552 hold off


554 print(strcat('Resposta torque detalhe ',Analise),'−depsc')

555

556 %% Analise de sensibilidade da reposta do sistema quando existe ...

variacao de propriedades da interacao broca−rocha

557


referencia

559 ta = 0:0.001:90; % Vetor temporal

560


562 n sim = 1000;

563 a = 0.95;

564 b = 1.05;

565 t degrau = zeros(length(ta),n sim);

566 y degrau = zeros(length(ta),n sim*2);

567 t linear = zeros(length(ta),n sim);

568 y linear = zeros(length(ta),n sim*2);

569 t tanh = zeros(length(ta),n sim);

570 y tanh = zeros(length(ta),n sim*2);

571 Tcb1 l = zeros(length(ta),n sim);

572 Tsb1 l = zeros(length(ta),n sim);

573 TcbR l = zeros(length(ta),n sim);

574 TsbR l = zeros(length(ta),n sim);

575 ti l = zeros(length(ta),n sim);

576 tf l = zeros(length(ta),n sim);

577

578 % Simulacoes

579 h = waitbar(0,'calculating . . .');


581 for j = 1:n sim

582 Tcb1 l(j) = 3*normrnd(1,0.02); % Torque de atrito dinamico ...

para t < t1

58

583 Tcb1 = Tcb1 l(j);

584 TcbR l(j) = 5*normrnd(1,0.02); % Torque de atrito dinamico ...

referencia

585 TcbR = TcbR l(j);

586 Tsb1 l(j) = 6*normrnd(1,0.02); % Pico de torque de atrito ...

para t < t1

587 Tsb1 = Tsb1 l(j);

588 TsbR l(j) = 8*normrnd(1,0.02); % Pico de torque de atrito ...

referencia

589 TsbR = TsbR l(j);

590 ti l(j) = 30*normrnd(1,0.02); % Tempo inicial para ...

transicao das propriedades de atrito

591 ti = ti l(j);

592 tf l(j) = 70*normrnd(1,0.02); % Tempo final para transicao ...

das propriedades de atrito

593 tf = tf l(j);



596

597 index = j;

598

599 [t degrau(:,j),y degrau(:,j*2−1:j*2)] = ...

ode23tb(@modelo degrau,ta,[0 0]);

600 [t linear(:,j),y linear(:,j*2−1:j*2)] = ...

ode23tb(@modelo linear,ta,[0 0]);

601 [t tanh(:,j),y tanh(:,j*2−1:j*2)] = ode23tb(@modelo tanh,ta,[0 0]);

602

603 waitbar(j/n sim,h,sprintf('calculating . . . (%d/%d)',j,n sim));

604 end

605 close(h)

606

607 % Calculo do torque de atrito

608 h = waitbar(0,'calculating . . .');

609 for j=1:3

610 if j==1

611 y = y degrau;

612 elseif j == 2

613 y = y linear;

614 else

59

615 y = y tanh;

616 end

617

618 for n = 1:n sim

619 for i = 1:length(ta)

620 Teb = ...

c*(Vrot−y(i,n*2−1))+k*(Vrot*ta(i)−y(i,n*2−1))−cb*y(i,n*2); ...

% Toque Stick

621

622 if j == 3

623 Tcb = Tcb1 l(n) + ...

(TcbR l(n)−Tcb1 l(n))*(tanh((((ta(i)−ti l(n))/(tf l(n)−ti l(n)))*5)−2.5)+1)/2;

624 Tsb = Tsb1 l(n) + ...

(TsbR l(n)−Tsb1 l(n))*(tanh((((ta(i)−ti l(n))/(tf l(n)−ti l(n)))*5)−2.5)+1)/2;

625 elseif (ti>ta(i))

626 Tsb = Tsb1 l(n);

627 Tcb = Tcb1 l(n);

628 elseif (ti l(n)<ta(i)) && (ta(i)≤tf l(n))

629 if j == 1

630 Tsb = TsbR l(n);

631 Tcb = TcbR l(n);

632 elseif j == 2

633 Tcb = Tcb1 l(n) + ...

((ta(i)−ti l(n))/(tf l(n)−ti l(n)))*(TcbR l(n)−Tcb1 l(n));

634 Tsb = Tsb1 l(n) + ...

((ta(i)−ti l(n))/(tf l(n)−ti l(n)))*(TsbR l(n)−Tsb1 l(n));

635 end

636 elseif ta(i) > tf l(n)

637 Tsb = TsbR l(n);

638 Tcb = TcbR l(n);

639 end

640

641

642 if and(y(i,2*j)<Dv,abs(Teb)<TsbR l(n)) ...

% Stick

643 Tfb(i,(j−1)*n sim+n)=Teb;

644 elseif and(y(i,2*j)<Dv,abs(Teb)>TsbR l(n)) ...

% Stick−Slip Transition

645 Tfb(i,(j−1)*n sim+n)=Tsb*sign(Teb);

60

646 else % Slip

647 Tfb(i,(j−1)*n sim+n)=Tcb+(Tsb−Tcb)*exp(−Gb*abs(y(i,n*2)))*sign(y(i,n*2));

648 end

649

650 end

651 waitbar((n)/(n sim),h,sprintf('calculating . . . ...

(%d/%d)',((j−1)*n sim+n),(n sim)));

652 end

653 end

654 close(h)

655

656 % Grafico Angulo da broca (Degrau)

657 figure

658 plot(ta,mean(transp(y degrau(:,1:2:end))),'k',ta,prctile(transp(y degrau(:,1:2:end)),5),'k−−',ta,prctile(transp(y degrau(:,1:2:end)),95),'k−−','LineWidth',2)



661 ylim([0 ...




664 hold on



667 hold off


669 print('Resposta angulo alt Tcb Tsb estoc degrau','−depsc')

670

671 % Grafico Angulo da broca (Linear)

672 figure

673 plot(ta,mean(transp(y linear(:,1:2:end))),'k',ta,prctile(transp(y linear(:,1:2:end)),5),'k−−',ta,prctile(transp(y linear(:,1:2:end)),95),'k−−','LineWidth',2)



676 ylim([0 ...




679 hold on



61

682 hold off


684 print('Resposta angulo alt Tcb Tsb estoc linear','−depsc')

685

686 % Grafico Angulo da broca (Tanh)

687 figure

688 plot(ta,mean(transp(y tanh(:,1:2:end))),'k',ta,prctile(transp(y tanh(:,1:2:end)),5),'k−−',ta,prctile(transp(y tanh(:,1:2:end)),95),'k−−','LineWidth',2)



691 ylim([0 ...




694 hold on



697 hold off


699 print('Resposta angulo alt Tcb Tsb estoc tanh','−depsc')

700

701 % Grafico Velocidade da broca (Degrau)

702 figure

703 plot(ta,mean(transp(y degrau(:,2:2:end))),'k',ta,prctile(transp(y degrau(:,2:2:end)),5),'k−−',ta,prctile(transp(y degrau(:,2:2:end)),95),'k−−','LineWidth',2)



706 ylim([0 ...




709 hold on



712 hold off


714 print('Resposta velocidade alt Tcb Tsb estoc degrau','−depsc')

715

716 % Grafico Velocidade da broca (Linear)

717 figure

718 plot(ta,mean(transp(y linear(:,2:2:end))),'k',ta,prctile(transp(y linear(:,2:2:end)),5),'k−−',ta,prctile(transp(y linear(:,2:2:end)),95),'k−−','LineWidth',2)

62



721 ylim([0 ...




724 hold on



727 hold off


729 print('Resposta velocidade alt Tcb Tsb estoc linear','−depsc')

730

731 % Grafico Velocidade da broca (Tanh)

732 figure

733 plot(ta,mean(transp(y tanh(:,2:2:end))),'k',ta,prctile(transp(y tanh(:,2:2:end)),5),'k−−',ta,prctile(transp(y tanh(:,2:2:end)),95),'k−−','LineWidth',2)



736 ylim([0 ...




739 hold on



742 hold off


744 print('Resposta velocidade alt Tcb Tsb estoc tanh','−depsc')

745

746 % Grafico Torque de atrito na broca (Degrau)

747 figure

748 plot(ta,mean(transp(Tfb(:,1:1000))),'k',ta,prctile(transp(Tfb(:,1:1000)),5),'k−−',ta,prctile(transp(Tfb(:,1:1000)),95),'k−−','LineWidth',2)



751 ylim([0 1.1*max(max(Tfb(:,1:1000)))])



754 hold on


63


757 hold off


759 print('Resposta torque alt Tcb Tsb estoc degrau','−depsc')

760

761 % Grafico Torque de atrito na broca (Linear)

762 figure




766 ylim([0 1.1*max(max(Tfb(:,1:1000)))])



769 hold on



772 hold off


774 print('Resposta torque alt Tcb Tsb estoc linear','−depsc')

775

776 % Grafico Torque de atrito na broca (Tanh)

777 figure




781 ylim([0 1.1*max(max(Tfb(:,1:1000)))])



784 hold on



787 hold off


789 print('Resposta torque alt Tcb Tsb estoc tanh','−depsc')

790

791 %% Graficos para alteracao dos torques de atrito no tempo

792

793 % Mudanca de caracteristicas da rocha


64









referencia



803 tb = 0:0.1:90;

804 tb1 = tb;

805 tb2 = tb;

806 tb3 = tb;

807 j=1;

808 Tsb degrau = zeros(length(tb),1);

809 Tcb degrau = zeros(length(tb),1);

810 Tsb linear = zeros(length(tb),1);

811 Tcb linear = zeros(length(tb),1);

812 Tsb tanh = zeros(length(tb),1);

813 Tcb tanh = zeros(length(tb),1);

814

815 for t = tb

816 if t < ti

817 Tsb degrau(j) = Tsb1;

818 Tcb degrau(j) = Tcb1;

819 elseif t ≥ ti



822 end

823

824 if t < ti

825 Tsb linear(j) = Tsb1;

826 Tcb linear(j) = Tcb1;

827 elseif (ti≤t) && (t≤tf)

828 Tcb linear(j) = Tcb1 + ((t−ti)/(tf−ti))*(Tcb2−Tcb1);

829 Tsb linear(j) = Tsb1 + ((t−ti)/(tf−ti))*(Tsb2−Tsb1);

830 elseif t > tf

65



833 end

834

835 Tcb tanh(j) = Tcb1 + ...

(Tcb2−Tcb1)*(tanh((((t−ti)/(tf−ti))*5)−2.5)+1)/2;

836 Tsb tanh(j) = Tsb1 + ...

(Tsb2−Tsb1)*(tanh((((t−ti)/(tf−ti))*5)−2.5)+1)/2;

837

838 j = j + 1;

839 end

840

841 %Grafico Tsb

842 figure

843 plot(tb1,Tsb degrau,'k',tb2,Tsb linear,'k−−',tb3,Tsb tanh,'k:','LineWidth',2)


845 ylabel ('$T sb$(N.m)','Interpreter','latex','FontSize',26)

846 ylim([0.8*Tsb1 1.1*Tsb2])

847 xlim([0 max(tb)])



850 print('Tsb tempo alt Tsb Tcb','−depsc')

851

852 % Grafico Tcb

853 figure

854 plot(tb1,Tcb degrau,'k',tb2,Tcb linear,'k−−',tb3,Tcb tanh,'k:','LineWidth',2)



857 ylabel ('$T cb$(N.m)','Interpreter','latex','FontSize',26)

858 ylim([0.8*Tcb1 1.1*Tcb2])



861 print('Tcb tempo alt Tsb Tcb','−depsc')

862

863 % Rotacao Inicial Livre





66






referencia



873 tb = 0:0.1:90;

874 j=1;

875 Tsb degrau = zeros(length(tb),1);

876 Tcb degrau = zeros(length(tb),1);

877 Tsb linear = zeros(length(tb),1);

878 Tcb linear = zeros(length(tb),1);

879 Tsb tanh = zeros(length(tb),1);

880 Tcb tanh = zeros(length(tb),1);

881

882 for t = tb

883 if t < ti



886 elseif t ≥ ti



889 end

890

891 if t < ti




895 Tcb linear(j) = Tcb1 + ((t−ti)/(tf−ti))*(Tcb2−Tcb1);

896 Tsb linear(j) = Tsb1 + ((t−ti)/(tf−ti))*(Tsb2−Tsb1);

897 elseif t > tf



900 end

901

902 Tcb tanh(j) = Tcb1 + ...

(Tcb2−Tcb1)*(tanh((((t−ti)/(tf−ti))*5)−2.5)+1)/2;

67

903 Tsb tanh(j) = Tsb1 + ...

(Tsb2−Tsb1)*(tanh((((t−ti)/(tf−ti))*5)−2.5)+1)/2;

904

905 j = j + 1;

906 end

907

908 % Grafico Tsb

909 figure

910 plot(tb,Tsb degrau,'k',tb,Tsb linear,'k−−',tb,Tsb tanh,'k:','LineWidth',2)



913 ylabel ('$T sb$(N.m)','Interpreter','latex','FontSize',26)

914 ylim([0.8*Tsb1 1.1*Tsb2])



917 print('Tsb tempo rot livre','−depsc')

918

919 % Grafico Tcb

920 figure

921 plot(tb,Tcb degrau,'k',tb,Tcb linear,'k−−',tb,Tcb tanh,'k:','LineWidth',2)



924 ylabel ('$T cb$(N.m)','Interpreter','latex','FontSize',26)

925 ylim([0.8*Tcb1 1.1*Tcb2])



928 print('Tcb tempo rot livre','−depsc')

A.2 Representacao dos modelos utilizados

A.2.1 Modelo simples - sem mudanca de rocha

1 function dy=modelo simples(t,y)

2

3 global k Jb c Vrot Gb Dv cb Tsb Tcb

4

5 Teb = c*(Vrot−y(2))+k*(Vrot*t−y(1))−cb*y(2); % Torque Stick

68

6

7 if and(y(2)<Dv,abs(Teb)<Tsb) % Stick

8 Tfb=Teb;

9 elseif and(y(2)<Dv,abs(Teb)>Tsb) % Stick−Slip ...

Transition

10 Tfb=Tsb*sign(Teb);

11 else % Slip

12 Tfb=Tcb+(Tsb−Tcb)*exp(−Gb*abs(y(2)));

13 Tfb=Tfb.*sign(y(2));

14 end

15

16 % Discretizacao em espaco de estados

17 dy = zeros(2,1);

18 dy(1)=y(2);

19 dy(2)=(−k*(y(1)−Vrot*t)−c*(y(2)−Vrot)−cb*y(2)−Tfb)/Jb;

20 end

A.2.2 Modelo com multiprocessamento

1 function dy=modelo multiproc(t,y,p)

2

3 Jb=0.0318; % Inercia da coluna

4 c=0.0001; % constante de amortecimento da coluna

5 cb=0.03; % constate do atrito viscoso da broca

6 k=0.073; % rigidez da coluna

7 Dv=10ˆ−6; % Correcao de intervalo do modelo de atrito ...

seco (problemas de convergencia)

8 Gb=0.9; % Constante gama do modelo de atrito seco

9 Vrot = p(1); % Velocidade de rotacao da mesa

10 Tsb = p(2); % Pico de torque de atrito

11 Tcb = p(3); % Torque de atrito dinamico

12

13

14 Teb = c*(Vrot−y(2))+k*(Vrot*t−y(1))−cb*y(2); % Toque Stick

15


17 Tfb=Teb;


69

Transition


20 else % Slip



23 end

24


26 dy = zeros(2,1);

27 dy(1)=y(2);


29 end

A.2.3 Modelo com mudanca de rocha tipo degrau

1 function dy=modelo degrau(t,y)

2

3 global k Jb c Vrot Gb Dv cb Tsb1 Tsb2 Tcb1 Tcb2 ti

4


6

7 % Mudanca degrau

8 if t < ti

9 Tsb = Tsb1;

10 Tcb = Tcb1;

11 elseif t ≥ ti

12 Tsb = Tsb2;

13 Tcb = Tcb2;

14 end

15


17 Tfb=Teb;


Transition


20 else % Slip



70

23 end

24


26 dy = zeros(2,1);

27 dy(1)=y(2);


29 end

A.2.4 Modelo com mudanca de rocha tipo linear

1 function dy=modelo linear(t,y)

2

3 global k Jb c Vrot Gb Dv cb Tsb1 Tsb2 Tcb1 Tcb2 ti tf

4


6

7 %Mundanca Linear

8 if t < ti

9 Tsb = Tsb1;

10 Tcb = Tcb1;


12 Tcb = Tcb1 + ((t−ti)/(tf−ti))*(Tcb2−Tcb1);

13 Tsb = Tsb1 + ((t−ti)/(tf−ti))*(Tsb2−Tsb1);

14 elseif t > tf

15 Tsb = Tsb2;

16 Tcb = Tcb2;

17 end

18


20 Tfb=Teb;


Transition


23 else % Slip



26 end

27

71


29 dy = zeros(2,1);

30 dy(1)=y(2);


32 end

A.2.5 Modelo com mudanca de rocha tipo tangente hi-

perbolica

1 function dy=modelo tanh(t,y)

2

3 global k Jb c Vrot Gb Dv cb Tsb1 Tsb2 Tcb1 Tcb2 ti tf

4


6

7 % Mudanca Tangente Hiperbolica

8 Tcb = Tcb1 + (Tcb2−Tcb1)*(tanh((((t−ti)/(tf−ti))*5)−2.5)+1)/2;

9 Tsb = Tsb1 + (Tsb2−Tsb1)*(tanh((((t−ti)/(tf−ti))*5)−2.5)+1)/2;

10


12 Tfb=Teb;


Transition


15 else % Slip



18 end

19


21 dy = zeros(2,1);

22 dy(1)=y(2);


24 end

72

Documents

Análise dinâmica de uma coluna de perfuração …...de Moraes Lobo, Daniel An alise din^amica de uma coluna de perfura˘c~ao durante a mudan˘ca de caracter sticas da rocha com